MATLAB:使用小矩阵的元素高效构造大型矩阵

Question

我有一个大小为 (n+1) x (n+1) 的小矩阵A，n = O(10) 。我需要N=(n+1)^m B_i构建更大的矩阵，每个矩阵的大小N x m，其中m = O(10)。 幸运的是，我当时只需要存储一个B_i，所以内存不是问题。为了生成每个B_i，我有一个矩阵I索引，大小为 N x m 。 I 包含 1 和 n+1 之间所有可能的整数m -ple。例如，让 n=4 、 m=4 、

I=[1 1 1 1; 
   2 1 1 1; 
   3 1 1 1; 
   4 1 1 1; 
   5 1 1 1; 
   1 2 1 1; 
       .
       .
   5 5 5 5] 

I 的行 i 中的索引是 B_i 的相应列的元素A的列索引。例如，考虑i=3，即B_3=C。然后矩阵C的列1使用A列3中的元素组装，而2列4使用A列1中的元素。同样，I列j中的索引是B_i列j元素A的行索引。再考虑B_3=C：

C(1,1) = A(1,3)
C(2,1) = A(2,3)
C(3,1) = A(3,3)
C(4,1) = A(4,3)
C(5,1) = A(5,3)
C(6,1) = A(1,3)
       .
       .
C(N,1) = A(5,3)
C(1,2) = A(1,1)
C(2,2) = A(1,1)
C(3,2) = A(1,1)
C(4,2) = A(1,1)
C(5,2) = A(1,1)
C(6,2) = A(2,1)
       .
       .
C(N,2) = A(5,1)

等等。使用双 for 循环构建每个B_i是微不足道的。但是，由于我必须构建N B_i，因此生成的三重for循环需要很长时间。你能告诉我一些聪明的 MATLAB 技巧来超快速构建这些矩阵吗？

顺便说一句，我真的不需要B_i，而只需要他们列的元素乘积，即

PolyMD=prod(C,2)

这是一个N x 1列向量。正如我之前所说，我一次只存储一个B_i，所以先计算B_i，然后再PolyMD不是一个大问题。不过，如果有一种快速简单的方法可以直接获得PolyMD，我很想知道。

编辑：这是一个完整的示例，希望使问题更加清晰。

% input
n=2;
m=3;
A=[1.000000000000000e+00    -1.732050807568877e+00     1.414213562373095e+00; 
   1.000000000000000e+00     1.160285448625519e-16    -7.071067811865475e-01;
   1.000000000000000e+00     1.732050807568877e+00     1.414213562373095e+00]; 
I=[     1     1     1;
        2     1     1;
        3     1     1;
        1     2     1;
        2     2     1;
        3     2     1;
        1     3     1;
        2     3     1;
        3     3     1;
        1     1     2;
        2     1     2;
        3     1     2;
        1     2     2;
        2     2     2;
        3     2     2;
        1     3     2;
        2     3     2;
        3     3     2;
        1     1     3;
        2     1     3;
        3     1     3;
        1     2     3;
        2     2     3;
        3     2     3;
        1     3     3;
        2     3     3;
        3     3     3];
% desired output for B_2 (I skip B_1 since it's trivial, and B_3,....B_27 for brevity!)    
B_2=[ -1.732050807568877e+00     1.000000000000000e+00     1.000000000000000e+00;
     1.160285448625519e-16     1.000000000000000e+00     1.000000000000000e+00;
     1.732050807568877e+00     1.000000000000000e+00     1.000000000000000e+00;
    -1.732050807568877e+00     1.000000000000000e+00     1.000000000000000e+00;
     1.160285448625519e-16     1.000000000000000e+00     1.000000000000000e+00;
     1.732050807568877e+00     1.000000000000000e+00     1.000000000000000e+00;
    -1.732050807568877e+00     1.000000000000000e+00     1.000000000000000e+00;
     1.160285448625519e-16     1.000000000000000e+00     1.000000000000000e+00;
     1.732050807568877e+00     1.000000000000000e+00     1.000000000000000e+00;
    -1.732050807568877e+00     1.000000000000000e+00     1.000000000000000e+00;
     1.160285448625519e-16     1.000000000000000e+00     1.000000000000000e+00;
     1.732050807568877e+00     1.000000000000000e+00     1.000000000000000e+00;
    -1.732050807568877e+00     1.000000000000000e+00     1.000000000000000e+00;
     1.160285448625519e-16     1.000000000000000e+00     1.000000000000000e+00;
     1.732050807568877e+00     1.000000000000000e+00     1.000000000000000e+00;
    -1.732050807568877e+00     1.000000000000000e+00     1.000000000000000e+00;
     1.160285448625519e-16     1.000000000000000e+00     1.000000000000000e+00;
     1.732050807568877e+00     1.000000000000000e+00     1.000000000000000e+00;
    -1.732050807568877e+00     1.000000000000000e+00     1.000000000000000e+00;
     1.160285448625519e-16     1.000000000000000e+00     1.000000000000000e+00;
     1.732050807568877e+00     1.000000000000000e+00     1.000000000000000e+00;
    -1.732050807568877e+00     1.000000000000000e+00     1.000000000000000e+00;
     1.160285448625519e-16     1.000000000000000e+00     1.000000000000000e+00;
     1.732050807568877e+00     1.000000000000000e+00     1.000000000000000e+00;
    -1.732050807568877e+00     1.000000000000000e+00     1.000000000000000e+00;
     1.160285448625519e-16     1.000000000000000e+00     1.000000000000000e+00;
     1.732050807568877e+00     1.000000000000000e+00     1.000000000000000e+00];

Answer 1

也许这就是你想要的：

n = 4;
m = 4;
N = (n+1)^m;
A = rand(n+1,n+1); %// example data
I = fliplr(dec2base(0:(n+1)^m-1,n+1)-'0'+1); %// combinations of column indices
for ii = 1:N %// it's better not to use "i" as a variable name
    B_i = A(bsxfun(@plus, I, (I(ii,:)-1)*(n+1))); %// linear indexing into A
    result_i = prod(B_i,2); %// "Poly" is a reserved word. Better not use it.
end