如何在 Matlab 中加速嵌套的 for 循环代码

Question

我试图加速以下代码：

xdim=5560;
ydim=6945;
zdim=7;
Nexi=270250785;
Neyi=270260480;

kne=1;
for i=1:xdim
    tic
    for j=1:ydim
        for k=1:zdim
            if j<ydim
                map_ey(j+(k-1)*(ydim-1)+(i-1)*(ydim-1)*zdim+Nexi)=kne;
                kne=kne+1;
            end
        end
    end
    for j=1:ydim
        for k=1:zdim
            if k<zdim
                map_ez(j+(k-1)*ydim+(i-1)*ydim*(zdim-1)+Nexi+Neyi)=kne;
                kne=kne+1;
            end
        end
    end
    for j=1:ydim
         for k=1:zdim
            if i<xdim
                map_ex(j+(k-1)*ydim+(i-1)*ydim*zdim)=kne;
                kne=kne+1;
            end
        end
    end
    toc
end

从上面的代码可以看出，我使用了tic toc，我发现每个i值的每个循环迭代如下：

Elapsed time is 0.378693 seconds.
Elapsed time is 0.378163 seconds.
Elapsed time is 0.380865 seconds.
Elapsed time is 0.378316 seconds.
Elapsed time is 0.377940 seconds.
Elapsed time is 0.379100 seconds.
Elapsed time is 0.378394 seconds.
Elapsed time is 0.383894 seconds.
Elapsed time is 0.397466 seconds.
Elapsed time is 0.415835 seconds.

可以更快吗？请帮忙。提前致谢。

Answer 1

幸运的是，我认为这是 programming puzzle，尽管您不愿意清楚地解释您希望实现的目标，但我还是对您的代码进行了反卷积。这是一个简短的解决方案，测试规模较小：

xdim = 4;
ydim = 5;
zdim = 3;

% avoid magic numbers
Nexi = (xdim-1)*ydim*zdim;
Neyi = xdim*(ydim-1)*zdim;
Nezi = xdim*ydim*(zdim-1);

% preallocate
map_ex2 = zeros(1,Nexi);
map_ey2 = zeros(1,Nexi+Neyi);
map_ez2 = zeros(1,Nexi+Neyi+Nezi);


% index jump between blocks of the same variable
dN = ydim*zdim + (ydim-1)*zdim + (zdim-1)*ydim;

%construct matrices
map_ey2(Nexi+1:Nexi+Neyi) = ...
      permute(reshape(bsxfun(@plus,(1:xdim*zdim).',(0:ydim-2)*dN),...             %'
                      [zdim,xdim,ydim-1]),...
              [2 1 3]);
map_ez2(Nexi+Neyi+1:Nexi+Neyi+Nezi) = ...
      permute(reshape(bsxfun(@plus,(1:(zdim-1)*ydim).'+xdim*zdim,(0:xdim-1)*dN),...    %'
                      [zdim-1,ydim,xdim]),...
              [2 1 3]);
map_ex2(1:Nexi) = ...
      permute(reshape(bsxfun(@plus,(1:(xdim-1)*ydim).'+xdim*zdim+(zdim-1)*ydim,...    %'
                             (0:zdim-1)*dN),...
                      [xdim-1,ydim,zdim]),...
              [2 1 3]);

运行原始循环版本并检查

isequal(map_ex,map_ex2)
isequal(map_ey,map_ey2)
isequal(map_ez,map_ez2)

告诉我这应该和你正在做的一样，只是没有任何循环。

如您所见，我定义了动态变量Nexi、Neyi、Nezi，以避免使用幻数。我建议你采用这种行为。

上述解决方案无疑是混乱的，但您原来的三重循环也是如此。 @Divakar 有可能会以这种方式出现，并为您提供更优雅/高效的解决方案；）

Answer 2

正如您所发现的，for 循环（您实现它们的方式）在 MATLAB 中相对较慢。有时它们是无法避免的，但我们可以做一些简单的事情来尝试优化它们。

首先是查看mlint 警告以寻求帮助（MATLAB 编辑器中的橙色下划线）。如果您注意到，每次循环时都会扩大map_ex、map_ey 和map_ez 的大小。这会导致 MATLAB 不断重新分配用于存储这些变量的内存，这会导致程序随着这些数组大小的增长而变得越来越慢。为避免这种陷阱，请预先分配数组（即在开头创建一个空数组，该数组与您期望的数组一样大）。

在你的情况下，这似乎是以下

map_ex = zeros(1, (xdim - 1) * ydim * zdim);
map_ey = zeros(1, (xdim * (ydim - 1) * zdim) + Nexi);
map_ez = zeros(1, (xdim * ydim * (zdim - 1)) + Nexi + Neyi);

然后您可以按照您当前编写的相同方式执行您的 for 循环，您应该会立即注意到性能提升。

我在查看您的代码时注意到的另一件事是，您在整个过程中都重复了这个主题

for k=1:zdim
    if k < zdim
        % Do something
    end
end

这里的问题是，您每次都在比较 k 和 zdim

for k = 1:(zdim - 1)
    % Do something
end

这删除了一个不必要的 if 语句，该语句在循环中每次都被评估。

我相信您知道，在 MATLAB 中做事的首选方式是通过矢量化而不是 for 循环。如果不了解您在此处尝试实现的具体内容，很难以清晰简洁的方式对代码进行矢量化。我可以矢量化你到目前为止所拥有的，但我觉得可能有更好的方法来做到这一点。如果您能提供更多关于您要解决的问题的信息，也许我可以提供更多帮助。

此外，如果您确实需要性能并且无法避免使用 for 循环，MATLAB 确实提供了mex library，因此您可以在 C/C++ 中实现这种类型的代码，并且仍然可以在 MATLAB 中调用它。