为什么 gpuArray 上的 bsxfun 这么慢？

Question

我有 MATLAB 2016a 在 Win 10 机器上运行，配备库存 i5 2500K 和 2 GTX 970。我是 GPU 计算的新手，我正在探索如何使用我的 GPU 加快计算速度。

所以我运行以下简单代码：

clear;

A = randn(1000,1);
B = randn(100,1);
n = 10000;

gA = gpuArray(A);
gB = gpuArray(B);

myfunc = @(a,b)(a.*b);

tic;
for i = 1:n
    C = bsxfun(myfunc,A,B');
end
disp(toc);

tic;
for i = 1:n
    C = gather(bsxfun(myfunc,gA,gB'));
end
disp(toc);

我分别得到 8.2（秒）和 321.3864（秒）。

clear;

A = randn(1000,1);
B = randn(100,1);
n = 10000;

gA = gpuArray(A);
gB = gpuArray(B);

myfunc = @(a,b)(a.*b);

tic;
parfor i = 1:n
    C = bsxfun(myfunc,A,B');
end
disp(toc);

tic;
parfor i = 1:n
    C = gather(bsxfun(myfunc,gA,gB'));
end
disp(toc);

（区别：for --> parfor）。我得到了 2.7（秒）和 6.3（秒）。

为什么 GPU 方法在这两种情况下都比较慢？ 在我的工作中，myfunc 要复杂得多。我已经定义了它，以便它可以很好地与非 GPU bsxfun 一起工作，但是当我像上面所做的那样进行 GPU 化时，我遇到了错误 Use of functional workspace is not supported.（在我的工作中，myfunc 是在内部和parfor 循环的开始。）你能否解释一下这个错误是什么意思？

Answer 1

首先我要说的是，GPU 并不是可以以某种方式提高任何计算速率的神奇对象。它们是适用于某些工作的工具，但存在需要考虑的局限性。 GPU 的经验法则是数学运算比内存访问“便宜”，因此例如，如果您每次需要重新计算某个数组而不是将其保存到临时变量并访问它，那么为 GPU 编写的代码可能会运行得更好.底线 - GPU 编码需要一些不同的想法，这些事情超出了当前答案的范围。

---

这里列出了可以改进的地方：

1。随机数生成：

在 GPU 上生成随机数的效率要高得多，更不用说它为您节省了昂贵的通信开销。 MATLAB 为我们提供了several convenience functions to establish arrays on a GPU。换句话说，

A = randn(1000,1);
gA = gpuArray(A);

可以替换为：

gA = gpuArray.randn(1000,1);

2。为bsxfun重新定义现有函数：

没有必要这样做。看看the list of builtin functions supported by bsxfun：.*或times已经是其中之一！因此，您可以替换：

myfunc = @(a,b)(a.*b);
...
bsxfun(myfunc,A,B');

与：

bsxfun(@times,A,B.');

（或在 MATLAB 版本中 >= R2016b：A.*B.'）。

此外，最好在脚本文件中将自定义函数定义为nested function 并使用@myFunc 调用它，即：

function main
...
bsxfun(@myFunc,A,B')

% later in the same file, or in a completely different one:
function out = myFunc(a,b)
out = ...

3。使用ctranspose 而不是transpose：

这在这里解释得很好。长话短说：您应该养成使用.' 进行转置和使用' 进行复共轭转置的习惯。

4。定时函数执行：

长话短说：tic 和 toc 通常不是一个好的指示，请改用 timeit。

5。隐式创建并行池：

这是一个相当小的注释：在第二个^nd 代码 sn-p 中，您使用 parfor 而不首先调用 parpool。这意味着如果该阶段没有创建池，则创建时间（几秒）将添加到tic/toc 报告的时间中。为避免这种情况，请遵循"Explicit is better than implicit" 的编程原则，只需事先调用parpool。

6。比较苹果和苹果：

这两行代码做的工作量不同：

C = bsxfun(myfunc,A,B');
C = gather(bsxfun(myfunc,gA,gB'));

这是为什么呢？因为第 2^nd 版本还必须将 bsxfun 的结果从 GPU 内存传输到 RAM - 这不是 免费 （就运行时间而言）。在本示例中，这意味着您要在每次迭代中添加约 800KB 的数据传输。我假设你的实际问题有更大的矩阵，所以你明白这种开销很快就会变得很严重。

7。保留不需要的变量：

另一个小评论：而不是做：

parfor i = 1:n % or "for"
    C = bsxfun(myfunc,A,B');
end

你可以这样做：

parfor i = 1:n % or "for"
    [~] = bsxfun(myfunc,A,B');
end

---

至于错误，我无法在我的 R2016b 上重现它，但这听起来像是与捕获不兼容有关的一些问题（即在匿名函数中使用的变量快照的机制是在创建）与parfor 所需的切片。我不知道您到底做错了什么，但在我看来，您不应该在 parfor 迭代中定义函数。也许这些帖子会有所帮助：1、2。

Answer 2

我使用 Tesla K20 和 Quadro K620 在台式机上运行您的代码，在笔记本电脑上使用 GTX 9XXM 芯片。在这两种情况下，我都没有得到任何像你的第一组数字那样糟糕的东西，而且两张 GeForce 卡（嗯，你的卡和我的芯片）的统计数据也没有显示出如此大的差异。也许您的第一组时间安排因一些开销而出现偏差？当然，您的第二组时间安排似乎证明了这一点，因为时间改善了太多。

不过，我可以回答一般性问题。你为什么打电话给gather？大部分成本都在数据传输上，参见：

>> gputimeit(@()gA.*gB') 答案= 2.6918e-04 >> gputimeit(@()gather(gA.*gB')) 答案= 0.0011

因此，10 万个元素 (0.8MB) 的数据传输时间为 11 毫秒，而实际计算的成本为 0.3 毫秒。 （请注意，我没有使用bsxfun，因为 MATLAB R2016b 会自动对元素操作进行维度扩展，从而消除了对它的需要。）

除非您需要显示它或将其写入磁盘（或在其上运行不支持 GPU 计算的操作），否则您永远不需要收集 gpuArray，所以不要，将其留在设备上. CPU 并没有被这个数据传输问题所困扰，所以相比之下看起来当然不错。

另一点是 GeForce 卡的双精度性能很差 - 例如，您的 GTX 970 报告的单精度性能为 3494 GFlops，而双精度性能则为 109 GFlops。它们针对您看到的显示进行了优化。但是，在这种情况下切换到单精度不会有太大的不同，因为操作首先是数据传输限制，然后是内存带宽限制。计算时间并没有真正涉及。

就并行化而言，我质疑您的数字，因为对于仅两个 GPU 而言，这种改进太好了。但是，您有两个 GPU，因此它们可以并行进行数据传输（和计算），因此您可以获得改进。不过，似乎不足以抵消您的情况下的开销。