3

私は科学計算、特にmatlabの反復法Gauss-SeidelとSORの宿題を行っています。問題は、行列の場合は予期しない結果が得られ(解は収束しない)、別の行列の場合は収束することです。

これがsorのコードです。ここで:

  • A:システムの行列A * x = b
  • Xini:初期反復の配列
  • b:システムに依存しないアレイA * x = b
  • maxiter:最大反復回数
  • tol:許容値;
  • 特に、SOR法は、緩和パラメーターに対応するwと呼ばれる6番目のパラメーターを受け取ります。

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

ガウス・ザイデル法のコードは次のとおりです。

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)

gaussの出力は次のとおりです。

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

メソッドの動作は、両方の行列の条件付けに依存しますか?2番目のマトリックスが最初のマトリックスよりも条件が良いことに気づいたからです。助言がありますか?

4

1 に答える 1

4

Gauss-Seidelに関するwiki記事から:

収束が保証されるのは、行列が対角的に優勢であるか、対称で正定値である場合のみです。

SORはGauss-Seidelに似ているので、SORにも同じ条件が当てはまると思いますが、それを調べてみてください。

最初の行列は、対角的に優勢または対称ではありません。ただし、2番目の行列は対称で正定値です(all(A==A.')all(eig(A)>0))。

Matlabのデフォルトの方法(A\b)を「実際の」ソリューションとして使用し、各反復と「実際の」ソリューションの差のノルムをプロットすると、次の2つのグラフが得られます。最初の行列が収束しないことは明らかですが、2番目の行列は数回の反復後にすでに許容可能な結果を​​生成しています。

アルゴリズムを実際に適用する前に、常にアルゴリズムの制限を理解してください:)

最初のマトリックス-収束は最悪 2番目のマトリックス-収束は良好です

于 2012-09-02T08:43:21.480 に答える