Is it possible to code a function like movsum to perform as fast as the diff function when using gpu arrays?

3 次查看（过去 30 天）
Nathan Zechar
Nathan Zechar 2021-12-22
评论： Nathan Zechar 2021-12-23
I'm attempting to find the fastest way to average several points on 3D data. In other parts of a program I'm writting, I've used the diff function which is much faster than coding the difference between adjacent points along a dimension.
Since movsum appears to be closely related to diff, I thought it would be just as fast. But it's not. At least it isn't when using a gpu array. Perhaps then I need to understand how diff is written so I can write my own mx file to speed up the averaging process? Or would I need to understand both how cuda is implemented and mx files? Here is how I timed the performance of the movsum and diff functions with gpu arrays using the "Run and Time" button.
%% Compare movsum to diff
clear all
[Nx,Ny,Nz] = deal(150);
A = rand(Nx,Ny,Nz);
A = gpuArray(A);
dev = gpuDevice();
for i = 1:100
A0 = movsum(A,2,3);
A0 = A0(:,:,2:Nz);
wait(dev);
A1 = A(:,:,1:Nz-1)+A(:,:,2:Nz);
wait(dev);
end
for i = 1:100
A2 = diff(A,1,3);
wait(dev);
A3 = A(:,:,2:Nz)-A(:,:,1:Nz-1);
wait(dev);
end
As for what I'd ultimately like to do. Here is the code I used to average points upon 3D data, and how the same process can be done using movsum. Using the "Run and Time" button there is a large difference in performance. Because the diff function is more computationally efficient than hardcoding differences, surely there must be a way to gain that efficicency when calculating adjacent sums, yes?
[Nx,Ny,Nz] = deal(100);
A = rand(Nx,Ny,Nz);
A = gpuArray(A);
dev = gpuDevice();
for i = 1:100
%% Code I'm using to average points
A1 = (A(1:Nx-1,:,1:Nz-1) +...
A(1:Nx-1,:,2:Nz) +...
A(2:Nx ,:,1:Nz-1) +...
A(2:Nx ,:,2:Nz))*0.25;
wait(dev);
%% Idea I would like to use to speed up code
A0 = movsum(movsum(A,2,3),2,1)*0.25;
A0 = A0(2:Nx,:,2:Nz);
wait(dev);
end
If movsum performed similiar to diff when using gpuArray, then I believe I could speed up my code significantly. Is this a situation where I would need to learn
  2 个评论
Matt J
Matt J 2021-12-22
编辑：Matt J 2021-12-22
I am uncertain about your method of timing. By my measurements (on the GTX 1080Ti), your current method is only about a factor of 2 away from the speed offered by diff().
% function test
T=60;
[Nx,Ny,Nz] = deal(150);
A = gpuArray.rand(Nx,Ny,Nz);
dev = gpuDevice();
tic;
for i = 1:T
B0 = movsum(movsum(A,2,3,'End','d'),2,1,'End','d')*0.25;
end
wait(dev);
toc;%Elapsed time is 0.270494 seconds.
tic;
e=reshape(ones(2),2,1,2)/4;
for i = 1:T
B1 = convn(A,e,'valid');
end
wait(dev);
toc;%Elapsed time is 0.152638 seconds.
tic;
for i = 1:T
B2 = (A(1:Nx-1,:,1:Nz-1) +...
A(1:Nx-1,:,2:Nz) +...
A(2:Nx ,:,1:Nz-1) +...
A(2:Nx ,:,2:Nz))*0.25;
end
wait(dev);
toc;%Elapsed time is 0.042915 seconds.
tic;
for i = 1:T
D = diff(diff(A,1,1),1,3);
end
wait(dev);
toc;%Elapsed time is 0.019765 seconds.
Nathan Zechar
Nathan Zechar 2021-12-23
Interesting, I'll have to dive in to what you did with reshape. I believe you've answered other questions using the reshape approach. You're certainly a master at it! I've tried implementing it in other portions of my code. Sometimes it works faster, sometimes it's just a fun thought exercise.
And yes, as you can see, diff is just blazingly fast. It might offer better than 2x performance as 3D matricies become larger, I would have to test it though.
As for the uncertaintly in my method of timing... I was looking up other posts which focused on timing GPU computations. I hope I'm doing it correctly by adding "wait(dev)" in the code as I thought that's what was needed for the most accurate results. If there is a better way to time this, or a standard I should be following, I'm ready to learn.

请先登录，再进行评论。

请先登录，再回答此问题。

回答（1 个）

Matt J
Matt J 2021-12-22
编辑：Matt J 2021-12-22
Perhaps then I need to understand how diff is written so I can write my own mx file to speed up the averaging process?
I'm not a Mathworks developer so this is all my speculation, but it stands to reason that accelerating diff() on the GPU is going to be easier than for other routines like movsum. This is because each thread always requires exactly 2 pieces of data. The CUDA kernel threads can read all elements of A(1:Nx-1,:,:) from GPU global memory in one synchronized step, and then a A(2:Nx,:,:) in a second synchronized step. No memory sharing is required between the threads, all data access is embarassingly parallel and naturally coalesced, and the memory consumption per thread as so small that the fastest local memory registers can be used.
Conversely, for movsum and conv, the code needs to handle moving windows of arbitrary sizes, and so the parallelization strategy needs to be more complicated (and hence is slower).
You could almost certainly write your own optimized movsum CUDA kernel using the same principles as with diff(), as long as, like diff, the moving window size is always 2. Or at least, as long as the window size is small, fixed, and a priori known.
  1 个评论
Nathan Zechar
Nathan Zechar 2021-12-23
Yeah, I envision it should be as simple as changing a plus sign to a minus sign in order to mimic exactly what diff() is doing. If I knew how diff() was coded I could make so much more of my code operate so much faster knowing whatever principels were used.

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Logical 的更多信息

标签

产品


版本

R2020b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by