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

Question

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

这是一个Nx1 列向量。正如我之前所说，我一次只存储一个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