我有两个 3 维数组,前两个维度表示矩阵,最后一个维度通过参数空间计数,举个简单的例子
A = repmat([1,2; 3,4], [1 1 4]);
(但假设 A(:,:,j)
对于每个 j
都是不同的)。如何轻松地执行两个这样的矩阵数组 A
和 B
的每 j
矩阵乘法?
C = A; % pre-allocate, nan(size(A,1), size(B,2)) would be better but slower
for jj = 1:size(A, 3)
C(:,:,jj) = A(:,:,jj) * B(:,:,jj);
end
当然可以完成这项工作,但如果第三维更像是 1e3 元素,这会非常慢,因为它不使用 MATLAB 的矢量化。那么,有没有更快的方法呢?
最佳答案
我现在做了一些计时测试,2x2xN 最快的方法原来是计算矩阵元素:
C = A;
C(1,1,:) = A(1,1,:).*B(1,1,:) + A(1,2,:).*B(2,1,:);
C(1,2,:) = A(1,1,:).*B(1,2,:) + A(1,2,:).*B(2,2,:);
C(2,1,:) = A(2,1,:).*B(1,1,:) + A(2,2,:).*B(2,1,:);
C(2,2,:) = A(2,1,:).*B(1,2,:) + A(2,2,:).*B(2,2,:);
在一般情况下,事实证明 for 循环实际上是最快的(不过不要忘记预先分配 C!)。
不过,如果使用 cellfun 已经得到矩阵元胞数组的结果是最快的选择,它也比遍历单元格元素更快:
C = cellfun(@mtimes, A, B, 'UniformOutput', false);
但是,必须调用 num2cell首先 (Ac = num2cell(A, [1 2])
) 和 cell2mat
对于 3d 数组情况浪费了太多时间。
这是我为一组随机的 2 x 2 x 1e4 做的一些计时:
array-for: 0.057112
arrayfun : 0.14206
num2cell : 0.079468
cell-for : 0.033173
cellfun : 0.025223
cell2mat : 0.010213
explicit : 0.0021338
显式是指直接计算 2 x 2 矩阵元素,见下文。
对于新的随机数组,结果类似,如果之前不需要 num2cell
并且没有 2x2xN 的限制,则 cellfun
是最快的。对于一般的 3d 数组,在三维上循环确实已经是最快的选择了。这是时间代码:
n = 2;
m = 2;
l = 1e4;
A = rand(n,m,l);
B = rand(m,n,l);
% naive for-loop:
tic
%Cf = nan(n,n,l);
Cf = A;
for jl = 1:l
Cf(:,:,jl) = A(:,:,jl) * B(:,:,jl);
end;
disp([' array-for: ' num2str(toc)]);
% using arrayfun:
tic
Ca = arrayfun(@(k) A(:,:,k)*B(:,:,k), 1:size(A,3), 'UniformOutput',false);
Ca = cat(3,Ca{:});
disp([' arrayfun : ' num2str(toc)]);
tic
Ac = num2cell(A, [1 2]);
Bc = num2cell(B, [1 2]);
disp([' num2cell : ' num2str(toc)]);
% cell for-loop:
tic
Cfc = Ac;
for jl = 1:l
Cfc{jl} = Ac{jl} * Bc{jl};
end;
disp([' cell-for : ' num2str(toc)]);
% using cellfun:
tic
Cc = cellfun(@mtimes, Ac, Bc, 'UniformOutput', false);
disp([' cellfun : ' num2str(toc)]);
tic
Cc = cell2mat(Cc);
disp([' cell2mat : ' num2str(toc)]);
tic
Cm = A;
Cm(1,1,:) = A(1,1,:).*B(1,1,:) + A(1,2,:).*B(2,1,:);
Cm(1,2,:) = A(1,1,:).*B(1,2,:) + A(1,2,:).*B(2,2,:);
Cm(2,1,:) = A(2,1,:).*B(1,1,:) + A(2,2,:).*B(2,1,:);
Cm(2,2,:) = A(2,1,:).*B(1,2,:) + A(2,2,:).*B(2,2,:);
disp([' explicit : ' num2str(toc)]);
disp(' ');
关于arrays - MATLAB:如何向量乘两个矩阵数组?,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/6580656/