2

我正在使用 Matlab 来查找 Jacobi 迭代矩阵的光谱半径，其中A=[4 2 1;1 3 1;1 1 4].

在 5 次迭代后，我似乎无法输入正确的命令来获取错误的大小。有人能帮我吗？

以下是迄今为止我在 Matlab 中输入的命令列表：

A=[4 2 1;1 3 1;1 1 4]

A =

 4     2     1
 1     3     1
 1     1     4

 D=diagonal(diagonal(A));L=(A,-1);U=(A,1);

 b=([3 -1 4])

 x0j=zeros([0 0 0]);

 x=D\(-(U+L)*x0j+b);r=b-A*x    %Jacobi iteration.
------------------------------------------------------------------------------
 Error using  * 
 Inputs must be 2-D, o enter code here r at least one input must be scalar.
 To compute element wise TIMES, use TIMES (.*) instead.
4

1 回答 1

5

矩阵的谱半径是其特征值模的最大值。它可以简单地使用计算max(abs(eig(·)))。

但是，正如其他人所注意到的那样，您的整个代码似乎很混乱，实际上不是有效的 Matlab 代码，所以您的问题并不是真正计算光谱半径，是吗？该算法非常简单且易于实现：

% diagonal part of A and rest
D = diag(diag(A));
R = A - D;

% iteration matrix and offset
T = - inv(D) * R;
C = inv(D) * b;

% spectral radius condition
rho = max(abs(eig(T)));
if rho >= 1
    error('no convergence')
end

% initial guess
x = randn(size(b));

% iteration
while norm(A * x - b) > 1e-15
    x = T * x + C;
end

请注意，我inv(D)以前直接遵循 Wikipedia 中的描述，但是可以使用diag(1 ./ diag(D)).

我真的不明白为什么需要R分成上下两部分。我想这与数值效率有关，但是，Matlab 已经是一种非常有效的矩阵计算高级语言。所以实际上，当你可以简单地编写时，实际上没有必要在其中显式地实现 Jacobi 算法A \ b——我猜是为了教育目的。

于 2013-11-15T13:53:09.390 回答