当前位置：文江博客话题详情

将 3D 矩阵与 2D 矩阵相乘

发布于 2024-08-11 04:53:55 字数 185 浏览 4 评论 0原文

假设我有一个 AxBxC 矩阵 X 和一个 BxD 矩阵 Y。

是否有一种非循环方法可以将每个 C AxB 矩阵与 Y 相乘？

原文

分享到QQ

分享到微博

如果你对这篇内容有疑问，欢迎到本站社区发帖提问参与讨论，获取更多帮助，或者扫码二维码加入 Web 技术交流群。

发布评论

需要登录才能够评论，你可以免费注册一个本站的账号。

苦妄 2024-08-18 04:53:55

作为个人喜好，我希望我的代码尽可能简洁和可读。

这就是我会做的，尽管它不满足您的“无循环”要求：

for m = 1:C

    Z(:,:,m) = X(:,:,m)*Y;

end

这会产生一个 A x D x C 矩阵 Z。

当然，您始终可以使用 Z = Zeros(A,D,C); 预先分配 Z 以加快速度。

As a personal preference, I like my code to be as succinct and readable as possible.

Here's what I would have done, though it doesn't meet your 'no-loops' requirement:

for m = 1:C

    Z(:,:,m) = X(:,:,m)*Y;

end

This results in an A x D x C matrix Z.

And of course, you can always pre-allocate Z to speed things up by using Z = zeros(A,D,C);.

回复收藏 0 原文

可遇━不可求 2024-08-18 04:53:55

您可以使用函数 NUM2CELL 在一行中完成此操作将矩阵 X 分解为元胞数组和 CELLFUN 跨单元格操作：

Z = cellfun(@(x) x*Y,num2cell(X,[1 2]),'UniformOutput',false);

结果 Z 是一个 1-by-C 单元格数组，其中每个单元格包含一个 A-by -D 矩阵。如果您希望 Z 成为 A-by-D-by-C 矩阵，您可以使用 CAT 函数：

Z = cat(3,Z{:});

注意： 我使用的旧解决方案MAT2CELL 而不是 NUM2CELL，这并不简洁：

[A,B,C] = size(X);
Z = cellfun(@(x) x*Y,mat2cell(X,A,B,ones(1,C)),'UniformOutput',false);

You can do this in one line using the functions NUM2CELL to break the matrix X into a cell array and CELLFUN to operate across the cells:

Z = cellfun(@(x) x*Y,num2cell(X,[1 2]),'UniformOutput',false);

The result Z is a 1-by-C cell array where each cell contains an A-by-D matrix. If you want Z to be an A-by-D-by-C matrix, you can use the CAT function:

Z = cat(3,Z{:});

NOTE: My old solution used MAT2CELL instead of NUM2CELL, which wasn't as succinct:

[A,B,C] = size(X);
Z = cellfun(@(x) x*Y,mat2cell(X,A,B,ones(1,C)),'UniformOutput',false);

回复收藏 0 原文

撧情箌佬 2024-08-18 04:53:55

这是一个单行解决方案（如果要拆分为第三维，则需要两行）：

A = 2;
B = 3;
C = 4;
D = 5;

X = rand(A,B,C);
Y = rand(B,D);

%# calculate result in one big matrix
Z = reshape(reshape(permute(X, [2 1 3]), [A B*C]), [B A*C])' * Y;

%'# split into third dimension
Z = permute(reshape(Z',[D A C]),[2 1 3]);

因此现在： Z(:,:,i) 包含 X(:,:,i) 的结果) * Y

解释：

上面的内容可能看起来很混乱，但想法很简单。
首先，我从 X 的第三个维度开始，并沿着第一个维度进行垂直串联：

XX = cat(1, X(:,:,1), X(:,:,2), ..., X(:,:,C))

...困难在于 C 是一个变量，因此你可以' t 使用 cat 或 vertcat 概括该表达式。接下来，我们将其乘以 Y：

ZZ = XX * Y;

最后，我将其拆分回第三维：

Z(:,:,1) = ZZ(1:2, :);
Z(:,:,2) = ZZ(3:4, :);
Z(:,:,3) = ZZ(5:6, :);
Z(:,:,4) = ZZ(7:8, :);

因此您可以看到它只需要一次矩阵乘法，但您必须在之前重塑矩阵以及之后。

Here's a one-line solution (two if you want to split into 3rd dimension):

A = 2;
B = 3;
C = 4;
D = 5;

X = rand(A,B,C);
Y = rand(B,D);

%# calculate result in one big matrix
Z = reshape(reshape(permute(X, [2 1 3]), [A B*C]), [B A*C])' * Y;

%'# split into third dimension
Z = permute(reshape(Z',[D A C]),[2 1 3]);

Hence now: Z(:,:,i) contains the result of X(:,:,i) * Y

Explanation:

The above may look confusing, but the idea is simple.
First I start by take the third dimension of X and do a vertical concatenation along the first dim:

XX = cat(1, X(:,:,1), X(:,:,2), ..., X(:,:,C))

... the difficulty was that C is a variable, hence you can't generalize that expression using cat or vertcat. Next we multiply this by Y:

ZZ = XX * Y;

Finally I split it back into the third dimension:

Z(:,:,1) = ZZ(1:2, :);
Z(:,:,2) = ZZ(3:4, :);
Z(:,:,3) = ZZ(5:6, :);
Z(:,:,4) = ZZ(7:8, :);

So you can see it only requires one matrix multiplication, but you have to reshape the matrix before and after.

回复收藏 0 原文

寻梦旅人 2024-08-18 04:53:55

我正在解决完全相同的问题，着眼于最有效的方法。我看到的大约有三种方法，缺乏使用外部库的方法（即 mtimesx)：

循环遍历 3D 矩阵的切片
repmat-and-permute 魔法
cellfun 乘法

我最近比较了所有三种方法，看看哪种方法最快。我的直觉是（2）将是获胜者。代码如下：

% generate data
A = 20;
B = 30;
C = 40;
D = 50;

X = rand(A,B,C);
Y = rand(B,D);

% ------ Approach 1: Loop (via @Zaid)
tic
Z1 = zeros(A,D,C);
for m = 1:C
    Z1(:,:,m) = X(:,:,m)*Y;
end
toc

% ------ Approach 2: Reshape+Permute (via @Amro)
tic
Z2 = reshape(reshape(permute(X, [2 1 3]), [A B*C]), [B A*C])' * Y;
Z2 = permute(reshape(Z2',[D A C]),[2 1 3]);
toc


% ------ Approach 3: cellfun (via @gnovice)
tic
Z3 = cellfun(@(x) x*Y,num2cell(X,[1 2]),'UniformOutput',false);
Z3 = cat(3,Z3{:});
toc

所有三种方法都产生相同的输出（唷！），但令人惊讶的是，循环是最快的：

Elapsed time is 0.000418 seconds.
Elapsed time is 0.000887 seconds.
Elapsed time is 0.001841 seconds.

请注意，每次试验的时间可能相差很大，有时（2）是最慢的。随着数据的增加，这些差异变得更加显着。但对于更大的数据，(3) 胜过 (2)。循环方法仍然是最好的。

% pretty big data...
A = 200;
B = 300;
C = 400;
D = 500;
Elapsed time is 0.373831 seconds.
Elapsed time is 0.638041 seconds.
Elapsed time is 0.724581 seconds.

% even bigger....
A = 200;
B = 200;
C = 400;
D = 5000;
Elapsed time is 4.314076 seconds.
Elapsed time is 11.553289 seconds.
Elapsed time is 5.233725 seconds.

但是，如果循环尺寸比其他尺寸大得多，则循环方法可能比（2）慢。

A = 2;
B = 3;
C = 400000;
D = 5;
Elapsed time is 0.780933 seconds.
Elapsed time is 0.073189 seconds.
Elapsed time is 2.590697 seconds.

因此，在这种（可能是极端的）情况下，（2）以很大的优势获胜。可能没有一种方法在所有情况下都是最佳的，但循环仍然相当不错，并且在许多情况下都是最好的。在可读性方面也是最好的。循环走开！

I'm approaching the exact same issue, with an eye for the most efficient method. There are roughly three approaches that i see around, short of using outside libraries (i.e., mtimesx):

Loop through slices of the 3D matrix
repmat-and-permute wizardry
cellfun multiplication

I recently compared all three methods to see which was quickest. My intuition was that (2) would be the winner. Here's the code:

% generate data
A = 20;
B = 30;
C = 40;
D = 50;

X = rand(A,B,C);
Y = rand(B,D);

% ------ Approach 1: Loop (via @Zaid)
tic
Z1 = zeros(A,D,C);
for m = 1:C
    Z1(:,:,m) = X(:,:,m)*Y;
end
toc

% ------ Approach 2: Reshape+Permute (via @Amro)
tic
Z2 = reshape(reshape(permute(X, [2 1 3]), [A B*C]), [B A*C])' * Y;
Z2 = permute(reshape(Z2',[D A C]),[2 1 3]);
toc


% ------ Approach 3: cellfun (via @gnovice)
tic
Z3 = cellfun(@(x) x*Y,num2cell(X,[1 2]),'UniformOutput',false);
Z3 = cat(3,Z3{:});
toc

All three approaches produced the same output (phew!), but, surprisingly, the loop was the fastest:

Elapsed time is 0.000418 seconds.
Elapsed time is 0.000887 seconds.
Elapsed time is 0.001841 seconds.

Note that the times can vary quite a lot from one trial to another, and sometimes (2) comes out the slowest. These differences become more dramatic with larger data. But with much bigger data, (3) beats (2). The loop method is still best.

% pretty big data...
A = 200;
B = 300;
C = 400;
D = 500;
Elapsed time is 0.373831 seconds.
Elapsed time is 0.638041 seconds.
Elapsed time is 0.724581 seconds.

% even bigger....
A = 200;
B = 200;
C = 400;
D = 5000;
Elapsed time is 4.314076 seconds.
Elapsed time is 11.553289 seconds.
Elapsed time is 5.233725 seconds.

But the loop method can be slower than (2), if the looped dimension is much larger than the others.

A = 2;
B = 3;
C = 400000;
D = 5;
Elapsed time is 0.780933 seconds.
Elapsed time is 0.073189 seconds.
Elapsed time is 2.590697 seconds.

So (2) wins by a big factor, in this (maybe extreme) case. There may not be an approach that is optimal in all cases, but the loop is still pretty good, and best in many cases. It is also best in terms of readability. Loop away!

回复收藏 0 原文

寻找一个思念的角度 2024-08-18 04:53:55

没有。方法有很多种，但总是直接或间接地循环出现。

只是为了满足我的好奇心，你为什么想要这个？

回复收藏 0 原文

怪我入戏太深 2024-08-18 04:53:55

要回答这个问题，和为了便于阅读，请参阅：

ndmult ，作者：ajuanpi (Juan Pablo Carbajal)，2013 年，GNU GPL

输入

2 数组
暗淡

示例

 nT = 100;
 t = 2*pi*linspace (0,1,nT)’;

 # 2 experiments measuring 3 signals at nT timestamps
 signals = zeros(nT,3,2);
 signals(:,:,1) = [sin(2*t) cos(2*t) sin(4*t).^2];
 signals(:,:,2) = [sin(2*t+pi/4) cos(2*t+pi/4) sin(4*t+pi/6).^2];

 sT(:,:,1) = signals(:,:,1)’;
 sT(:,:,2) = signals(:,:,2)’;
   G = ndmult (signals,sT,[1 2]);

源

原始源。我添加了内嵌评论。

function M = ndmult (A,B,dim)
  dA = dim(1);
  dB = dim(2);

  # reshape A into 2d
  sA = size (A);
  nA = length (sA);
  perA = [1:(dA-1) (dA+1):(nA-1) nA dA](1:nA);
  Ap = permute (A, perA);
  Ap = reshape (Ap, prod (sA(perA(1:end-1))), sA(perA(end)));

  # reshape B into 2d
  sB = size (B);
  nB = length (sB);
  perB = [dB 1:(dB-1) (dB+1):(nB-1) nB](1:nB);
  Bp = permute (B, perB);
  Bp = reshape (Bp, sB(perB(1)), prod (sB(perB(2:end))));

  # multiply
  M = Ap * Bp;

  # reshape back to original format
  s = [sA(perA(1:end-1)) sB(perB(2:end))];
  M = squeeze (reshape (M, s));
endfunction

To answer the question, and for readability, please see:

ndmult, by ajuanpi (Juan Pablo Carbajal), 2013, GNU GPL

Input

2 arrays
dim

Example

 nT = 100;
 t = 2*pi*linspace (0,1,nT)’;

 # 2 experiments measuring 3 signals at nT timestamps
 signals = zeros(nT,3,2);
 signals(:,:,1) = [sin(2*t) cos(2*t) sin(4*t).^2];
 signals(:,:,2) = [sin(2*t+pi/4) cos(2*t+pi/4) sin(4*t+pi/6).^2];

 sT(:,:,1) = signals(:,:,1)’;
 sT(:,:,2) = signals(:,:,2)’;
   G = ndmult (signals,sT,[1 2]);

Source

Original source. I added inline comments.

function M = ndmult (A,B,dim)
  dA = dim(1);
  dB = dim(2);

  # reshape A into 2d
  sA = size (A);
  nA = length (sA);
  perA = [1:(dA-1) (dA+1):(nA-1) nA dA](1:nA);
  Ap = permute (A, perA);
  Ap = reshape (Ap, prod (sA(perA(1:end-1))), sA(perA(end)));

  # reshape B into 2d
  sB = size (B);
  nB = length (sB);
  perB = [dB 1:(dB-1) (dB+1):(nB-1) nB](1:nB);
  Bp = permute (B, perB);
  Bp = reshape (Bp, sB(perB(1)), prod (sB(perB(2:end))));

  # multiply
  M = Ap * Bp;

  # reshape back to original format
  s = [sA(perA(1:end-1)) sB(perB(2:end))];
  M = squeeze (reshape (M, s));
endfunction

回复收藏 0 原文

憧憬巴黎街头的黎明 2024-08-18 04:53:55

我强烈建议您使用 MMX 工具箱MATLAB 的。它可以尽可能快地乘以n维矩阵。

MMX的优点是：

使用简单。
乘以n维矩阵（实际上它可以乘以二维矩阵数组）
它执行其他矩阵运算（转置、二次乘法、Chol分解等等）
它使用< strong>C 编译器和多线程计算以提高速度。

对于这个问题，你只需要编写这个命令：

C=mmx('mul',X,Y);

这里是所有可能方法的基准。有关更多详细信息，请参阅此问题。

    1.6571 # FOR-loop
    4.3110 # ARRAYFUN
    3.3731 # NUM2CELL/FOR-loop/CELL2MAT
    2.9820 # NUM2CELL/CELLFUN/CELL2MAT
    0.0244 # Loop Unrolling
    0.0221 # MMX toolbox  <===================

I highly recommend you use the MMX toolbox of matlab. It can multiply n-dimensional matrices as fast as possible.

The advantages of MMX are:

It is easy to use.
Multiply n-dimensional matrices (actually it can multiply arrays of 2-D matrices)
It performs other matrix operations (transpose, Quadratic Multiply, Chol decomposition and more)
It uses C compiler and multi-thread computation for speed up.

For this problem, you just need to write this command:

C=mmx('mul',X,Y);

here is a benchmark for all possible methods. For more detail refer to this question.

    1.6571 # FOR-loop
    4.3110 # ARRAYFUN
    3.3731 # NUM2CELL/FOR-loop/CELL2MAT
    2.9820 # NUM2CELL/CELLFUN/CELL2MAT
    0.0244 # Loop Unrolling
    0.0221 # MMX toolbox  <===================

回复收藏 0 原文

乖乖兔^ω^ 2024-08-18 04:53:55

我想分享我对以下问题的答案：

1）制作两个张量（任意价）的张量积；

2）使两个张量沿任意维度收缩。

以下是我的第一个和第二个任务的子例程：

1）张量积：

function [C] = tensor(A,B)
   C = squeeze( reshape( repmat(A(:), 1, numel(B)).*B(:).' , [size(A),size(B)] ) );
end

2）收缩：
这里 A 和 B 分别是沿着维度 i 和 j 收缩的张量。当然，这些尺寸的长度应该相等。没有对此进行检查（这会使代码变得模糊），但除此之外它工作得很好。

   function [C] = tensorcontraction(A,B, i,j)
      sa = size(A);
      La = length(sa);
      ia = 1:La;
      ia(i) = [];
      ia = [ia i];

      sb = size(B);
      Lb = length(sb);
      ib = 1:Lb;
      ib(j) = [];
      ib = [j ib];

      % making the i-th dimension the last in A
      A1 = permute(A, ia);
      % making the j-th dimension the first in B
      B1 = permute(B, ib);

      % making both A and B 2D-matrices to make use of the
      % matrix multiplication along the second dimension of A
      % and the first dimension of B
      A2 = reshape(A1, [],sa(i));
      B2 = reshape(B1, sb(j),[]);

      % here's the implicit implication that sa(i) == sb(j),
      % otherwise - crash
      C2 = A2*B2;

      % back to the original shape with the exception
      % of dimensions along which we've just contracted
      sa(i) = [];
      sb(j) = [];
      C = squeeze( reshape( C2, [sa,sb] ) );
   end

有批评吗？

I would like to share my answer to the problems of:

1) making the tensor product of two tensors (of any valence);

2) making the contraction of two tensors along any dimension.

Here are my subroutines for the first and second tasks:

1) tensor product:

function [C] = tensor(A,B)
   C = squeeze( reshape( repmat(A(:), 1, numel(B)).*B(:).' , [size(A),size(B)] ) );
end

2) contraction:
Here A and B are the tensors to be contracted along the dimesions i and j respectively. The lengths of these dimensions should be equal, of course. There's no check for this (this would obscure the code) but apart from this it works well.

   function [C] = tensorcontraction(A,B, i,j)
      sa = size(A);
      La = length(sa);
      ia = 1:La;
      ia(i) = [];
      ia = [ia i];

      sb = size(B);
      Lb = length(sb);
      ib = 1:Lb;
      ib(j) = [];
      ib = [j ib];

      % making the i-th dimension the last in A
      A1 = permute(A, ia);
      % making the j-th dimension the first in B
      B1 = permute(B, ib);

      % making both A and B 2D-matrices to make use of the
      % matrix multiplication along the second dimension of A
      % and the first dimension of B
      A2 = reshape(A1, [],sa(i));
      B2 = reshape(B1, sb(j),[]);

      % here's the implicit implication that sa(i) == sb(j),
      % otherwise - crash
      C2 = A2*B2;

      % back to the original shape with the exception
      % of dimensions along which we've just contracted
      sa(i) = [];
      sb(j) = [];
      C = squeeze( reshape( C2, [sa,sb] ) );
   end

Any critics?

回复收藏 0 原文