如何使用Matlab有效地计算矩阵项

Question

我有一个稀疏矩阵B_1的单元格数组myBasis，…， B_n .

我想用Matlab计算矩阵Q(i,j) = trace (B^T_i * B_j)。

因此，我编写了以下代码:

for i=1:n
    for j=1:n
        B=myBasis{i};        
        C=myBasis{j};   
        Q(i,j)=trace(B'*C);
    end
end

当n=1226和B_i有50行和50列时，这已经需要68秒了。

有没有可能加快速度?通常我将for循环从我的matlab代码中排除在c++文件中-但我没有经验如何在c++中处理稀疏单元数组。

Answer 1

1. 正如Inox Q是对称的，因此你只需要显式地计算一半的条目。
2. 计算trace( B.'*C )相当于B(:).'*C(:):
   trace(B.'*C) = sum_i [B.'*C]_ii = sum_i sum_j B_ij * C_ij
   它是元素乘积的和，因此等于B(:).'*C(:)。
   当显式计算trace( B.'*C )时，实际上是预先计算B.'*C的所有 k × k项，以便稍后使用对角线。AFAIK, Matlab并没有优化它的计算，以节省计算所有的条目。

这里有一个方法

for ii = 1:n
    B = myBasis{ii};  
    for jj = ii:n
        C = myBasis{jj};
        t = full( B(:).'*C(:) ); % equivalent to trace(B'*C)!
        Q(ii,jj) = t;
        Q(jj,ii) = t;
    end
end

PS,
在Matlab中，最好不要使用i和j作为变量名。

PPS

,
您应该注意到Matlab中的'算子是而不是矩阵转置，而是厄米共轭，对于实际的转置，您需要使用.'。在大多数情况下，不涉及复数，两个操作符之间没有区别，但是一旦引入了复杂数据，两个操作符之间的混淆会使调试变得相当混乱…

Answer 2

我有几点想法

1)基本的东西:一个"* B = (B * A)和跟踪(A) =跟踪(")。只有这个技巧能让你的计算减少近50%。你的Q(i,j)矩阵是对称的，你只需要计算n(n+1)/2项(而不是n²)

2)为了计算轨迹，你不需要计算B'*C的每一项，只需要计算对角线。尽管如此，我不知道在Matlab中创建一个脚本是否容易，实际上比计算B'*C更快(Matlab在矩阵操作方面非常快)。

但我肯定会实现(1)