3

我正在为科学计算做家庭作业,特别是matlab中的迭代方法Gauss-Seidel和SOR,问题是对于一个矩阵给了我意想不到的结果(解决方案不收敛)并且对于另一个矩阵收敛。

这是sor的代码,其中:

  • A:系统矩阵 A * x = b
  • xini:初始迭代数组
  • b:与系统无关的数组 A * x = b
  • maxiter:最大迭代次数
  • tol:公差;
  • 特别是,SOR 方法将接收称为 w 的第六个参数,它对应于松弛参数。

这是sor方法的代码:

 function [x2,iter] = sor(A,xIni, b, maxIter, tol,w)            
     x1 = xIni;
     x2 = x1;
     iter = 0;
     i = 0;
     j = 0;
     n = size(A, 1);

     for iter = 1:maxIter,
         for i = 1:n
             a = w / A(i,i);
             x = 0;
             for j = 1:i-1
                 x = x + (A(i,j) * x2(j));
             end
             for j = i+1:n
                 x = x + (A(i,j) * x1(j));
             end
             x2(i) = (a * (b(i) - x)) + ((1 - w) * x1(i)); 
         end
         x1 = x2;
         if (norm(b - A * x2) < tol);
             break;
         end
     end

这是 Gauss-seidel 方法的代码:

function [x, iter] = Gauss(A, xIni, b, maxIter, tol)

x = xIni;
xnew = x;
iter = 0;
i = 0;
j = 0;
n = size(A,1);

for iter = 1:maxIter,
    for i = 1:n
        a = 1 / A(i,i);
        x1 = 0;
        x2 = 0;
        for j = 1:i-1
            x1 = x1 + (A(i,j) * xnew(j));
        end
        for j = i+1:n
            x2 = x2 + (A(i,j) * x(j));
        end
        xnew(i) = a * (b(i) - x1 - x2);
    end
    x= xnew;
    if ((norm(A*xnew-b)) <= tol);
        break;
    end
end

对于此输入:

A = [1 2 -2; 1 1 1; 2 2 1];
b = [1; 2; 5];

当调用函数 Gauss-Seidel 或 sor 时:

[x, iter] = gauss(A, [0; 0; 0], b, 1000, eps)
[x, iter] = sor(A, [0; 0; 0], b, 1000, eps, 1.5)

高斯的输出是:

x =

  1.0e+304 *

    1.6024
   -1.6030
    0.0011


 iter =

            1000

对于 sor 是:

x =

   NaN
   NaN
   NaN


iter =

        1000

但是对于以下系统能够找到解决方案:

A = [    4 -1  0 -1  0  0;
        -1  4 -1  0 -1  0;
         0 -1  4  0  0 -1;
        -1  0  0  4 -1  0;
         0 -1  0 -1  4 -1;
         0  0 -1  0 -1  4   ]
b = [1 0 0 0 0 0]'

解决方案:

[x, iter] = sor(A, [0; 0; 0], b, 1000, eps, 1.5)
x =

    0.2948
    0.0932
    0.0282
    0.0861
    0.0497
    0.0195


iter =

    52

方法的行为取决于两个矩阵的条件?因为我注意到第二个矩阵比第一个矩阵的条件更好。有什么建议么?

4

1 回答 1

4

来自关于 Gauss-Seidel 的 wiki 文章

仅当矩阵对角占优或对称且正定时,才能保证收敛

由于 SOR 与 Gauss-Seidel 相似,我预计 SOR 也适用相同的条件,但您可能需要查看该条件。

您的第一个矩阵绝对不是对角占优对称的。但是,您的第二个矩阵是对称且正定的(因为all(A==A.')all(eig(A)>0))。

如果您使用 Matlab 的默认方法A\b(很明显,第一个矩阵永远不会收敛,而第二个矩阵在几次迭代后已经产生了可接受的结果。

在将它们应用到野外之前,请务必了解算法的局限性:)

第一个矩阵——收敛很烂 第二个矩阵——收敛性好

于 2012-09-02T08:43:21.480 回答