Matlab中的块矩阵内部产品 [英] Block matrix inner products in Matlab

问题描述

我一直在使用以下自定义函数来将向量与矩阵相乘，其中向量的每个元素都将(3xN)x(3)矩阵内的3x3块相乘:

I have been using the following custom function to perform the multiplication of a vector by a matrix, in which each element of the vector multiplies a 3x3 block within a (3xN)x(3) matrix:

function [B] = BlockScalar(v,A)

 N=size(v,2);
 B=zeros(3*N,3);

 for i=1:N
     B(3*i-2:3*i,:) = v(i).*A(3*i-2:3*i,:);
 end

end

类似地，当我要将3x3矩阵的集合乘以3x3向量的集合时，我使用以下代码

Similarly, when I want to multiply a collection of 3x3 matrices by a collection of 3x3 vectors, I use the following

function [B] = BlockMatrix(A,u)

 N=size(u,2);
 B=zeros(N,3);

 for i=1:N
     B(i,:) = A(3*i-2:3*i,:)*u(:,i);
 end

end

由于我经常调用它们，因此不幸的是，这些操作大大降低了我的代码的运行速度.我想知道上述操作是否有一个更有效的(也许是矢量化的)版本.

Since I call them very often, these unfortunately, slow down the running of my code significantly. I was wondering if there was a more efficient (perhaps vectorised) version of the above operations.

推荐答案

在这两种情况下，您都可以取消for循环(尽管未经测试，但我无法确认这是否一定会加快您的计算速度).

In both instance, you are able to do away with the for-loops (although without testing, I cannot confirm if this will necessarily speed up your computation).

对于第一个功能，您可以执行以下操作:

For the first function, you can do it as follows:

function [B] = BlockScalar(v,A)
% We create a vector N = [1,1,1,2,2,2,3,3,3,...,N,N,N]
N=ceil((1:size(A,1))/3); 

% Use N to index v, and let matlab do the expansion
B = v(N).*A;

end

对于第二个函数，我们可以制作一个块对角矩阵.

For the second function, we can make a block-diagonal matrix.

function [B] = BlockMatrix(A,u)

 N=size(u,2)*3;
 % We use a little meshgrid+sparse magic to convert A to a block matrix
 [X,Y] = meshgrid(1:N,1:3);

% Use sparse matrices to speed up multiplication and save space
 B = reshape(sparse(Y+(ceil((1:N)/3)-1)*3,X,A) * (u(:)),3,size(u,2))';

end

请注意，如果您能够访问各个3x3矩阵，则可以通过使用本机blkdiag来使此操作更快/更简单:

Note that if you are able to access the individual 3x3 matrices, you can potentially make this faster/simpler by using the native blkdiag:

function [B] = BlockMatrix(a,b,c,d,...,u)
% Where A = [a;b;c;d;...];

% We make one of the input matrices sparse to make the whole block matrix sparse
% This saves memory and potentially speeds up multiplication by a lot
% For small enough values of N, however, using sparse may slow things down.
reshape(blkdiag(sparse(a),b,c,d,...) * (u(:)),3,size(u,2))';

end

这篇关于Matlab中的块矩阵内部产品的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！