5

由于 matlab 在执行 for 循环时很慢，我通常避免对所有代码进行 for 循环，并将它们转换为矩阵计算，这样会很快。但这里有一个问题我找不到聪明的方法：

我有anxn矩阵

A=[a1,a2,a3,...,an],

这里 a1,a2,a3....an 是矩阵的列。

另一个 nxn 矩阵

B=[b1,b2,b3,...,bn],

同样 b1,b2,b3... 也是 B 的列。

还有anxn矩阵M。

我想计算nxn矩阵

C=[c1,c2,c3,...,cn], 

thus (M+diag(ai))*ci = bi.

即

ci = (M+diag(ai))\bi.

我知道没有 for 循环的一种方法是：

C(:)=( blkdiag(M)+diag(A(:)) )\B(:).

但这会做太多的计算而不是需要。

有什么聪明的解决方案吗？您可以假设计算中没有奇点问题。

4

1 回答 1

2

自从 Matlab...euhm，R2008a 以来，“for 循环在 Matlab 中很慢”的说法通常不再正确？（有人请填写这个：）

无论如何，试试这个：

clc
clear all

M = rand(50);
a = rand(50);
b = rand(50);

% simple loop approach
tic
c = zeros(size(b));
for ii = 1:size(b,2)
    c(:,ii) = ( M+diag(a(:,ii)) ) \ b(:,ii);
end
toc

% not-so-simple vectorized approach
tic
MM = repmat({M}, 50,1);
c2 = (blkdiag(MM{:})+diag(a(:)))\b(:);
toc

norm(c(:)-c2(:))

结果：

Elapsed time is 0.011226 seconds.  % loop
Elapsed time is 1.049130 seconds.  % no-loop
ans =
     5.091221148787843e-10         % results are indeed "equal"

可能有更好的方法来矢量化操作，但我怀疑它会比 JIT 循环版本快得多。

有些问题不适合矢量化。我认为这是一个。

于 2012-12-03T11:36:07.957 回答