symmlq (MATLAB Function 参考文献)

対称 LQ 法

表示

x = symmlq(A,b)
symmlq(A,b,tol)
symmlq(A,b,tol,maxit)
symmlq(A,b,tol,maxit,M)
symmlq(A,b,tol,maxit,M1,M2)
symmlq(A,b,tol,maxit,M1,M2,x0)
symmlq(afun,b,tol,maxit,m1fun,m2fun,x0,p1,p2,...)
[x,flag] = symmlq(A,b,...)
[x,flag,relres] = symmlq(A,b,...)
[x,flag,relres,iter] = symmlq(A,b,...)
[x,flag,relres,iter,resvec] = symmlq(A,b,...)
[x,flag,relres,iter,resvec,resveccg] = symmlq(A,b,...)

詳細

x = symmlq(A,b) は、線形方程式 A*x=b を x に関して解きます。n 行 n 列の係数行列 A は、対称である必要がありますが、必ずしも正定である必要はありません。列ベクトル b は、長さ n をもっています。A は、afun(x) が、A*x を戻すような関数 afun です。

symmlq が収束する場合、その効果に関するメッセージが表示されます。symmlq が、設定された最大繰り返し回数に達しても収束しなかったり、または、何らかの理由で、計算を停止すると、ワーニングメッセージが、相対残差 norm(b-A*x)/norm(b) と停止したときの繰り返し回数を表示します。

symmlq(A,b,tol) は、使用する方法のトレランスを設定します。tol が [] の場合、symmlq は、デフォルトの 1e-6 を使います。

symmlq(A,b,tol,maxit) は、繰り返しの最大回数を設定します。maxit が [] の場合、symmlq は、デフォルトの min(n,20) を使います。

symmlq(A,b,tol,maxit,M) と symmlq(A,b,tol,maxit,M1,M2) は、対称正定前置条件子 M、または、M = M1*M2 を使って、システム inv(sqrt(M))*A*inv(sqrt(M))*y = inv(sqrt(M))*b を y に対して、効率的に解き、x = inv(sqrt(M))*yを戻します。M が [] の場合、symmlq は、前置条件子を適用しません。M は、M\x を戻す関数です。

symmlq(A,b,tol,maxit,M1,M2,x0) は、初期推定を指定するものです。x0 が [] の場合、symmlq は、デフォルトのすべての要素がゼロのベクトルを使用します。

symmlq(afun,b,tol,maxit,m1fun,m2fun,x0,p1,p2,...) は、パラメータ p1,p2,... を関数 afun(x,p1,p2,...), m1fun(x,p1,p2,...), m2fun(x,p1,p2,...) に渡します。

[x,flag] = symmlq(A,b,tol,maxit,M1,M2,x0,p1,p2,...) は、収束に関するフラグも戻します。

フラグ
収束

0
symmlq は、maxit 以内の繰り返し回数で、希望するトレランス tol 内に収束します。

1
symmlq は、maxit の繰り返しを実行しましたが、収束しませんでした。

2
前置条件子 M は、条件数が良くありません。

3
symmlq の出力結果に変化が見られません(連続する２回の繰り返しで、変化が認められない)。

4

symmlq を実行中にスカラ量の一つが、計算を続けるには、あまりにも大き過ぎるか、小さ過ぎます。

5

前置条件子 M は、対称で正定な行列でありません。

フラグ	収束
`0`	`symmlq` は、`maxit` 以内の繰り返し回数で、希望するトレランス `tol` 内に収束します。
`1`	`symmlq` は、`maxit` の繰り返しを実行しましたが、収束しませんでした。
`2`	前置条件子 `M` は、条件数が良くありません。
`3`	`symmlq` の出力結果に変化が見られません(連続する２回の繰り返しで、変化が認められない)。
`4`	`symmlq` を実行中にスカラ量の一つが、計算を続けるには、あまりにも大き過ぎるか、小さ過ぎます。
`5`	前置条件子 `M` は、対称で正定な行列でありません。

flag が 0 でない場合は常に、戻される解 x は、すべての繰り返しでの計算の中の最小ノルム残差を戻します。flag 出力が指定されている場合は、メッセージが表示されません。

[x,flag,relres] = symmlq(A,b,tol,maxit,M1,M2,x0,p1,p2,...) は、相対残差 norm(b-A*x)/norm(b) も戻します。flag が0 の場合、relres <= tol になります。

[x,flag,relres,iter] = symmlq(A,b,tol,maxit,M1,M2,x0,p1,p2,...) は、x が計算されたときの繰り返し回数も出力します。ここで、0 <= iter <= maxit です。

[x,flag,relres,iter,resvec] = symmlq(A,b,tol,maxit,M1,M2,x0, p1,p2,...) は、norm(b-A*x0) を含んで、各繰り返しで、symmlq の残差ノルムの推定から構成されるベクトルも出力します。

[x,flag,relres,iter,resvec,resveccg] = symmlq(A,b,tol,maxit,M1,M2, x0,p1,p2,...) は、各繰り返しで、共役勾配残差ノルムの推定からなるベクトルも出力します。

例題

例題 1.

n = 100; 
on = ones(n,1); 
A = spdiags([-2*on 4*on -2*on],-1:1,n,n);
b = sum(A,2); 
tol = 1e-10; 
maxit = 50; M1 = spdiags(4*on,0,n,n);

x = symmlq(A,b,tol,maxit,M1,[],[]);
symmlq converged at iteration 49 to a solution with relative 
residual 4.3e-015

また、つぎの行列-ベクトルの積の関数を使います。

function y = afun(x,n)
   y = 4 * x;
   y(2:n) = y(2:n) - 2 * x(1:n-1);
   y(1:n-1) = y(1:n-1) - 2 * x(2:n);

symmlq に入力

x1 = symmlq(@afun,b,tol,maxit,M1,[],[],n);

例題 2.

pcg でうまく機能しなかった対称な不定行列を使います。

A = diag([20:-1:1,-1:-1:-20]);
b = sum(A,2);      % 真の解は、すべてが１のベクトルです。
x = pcg(A,b);      % 最初の計算でエラーを出力します。
pcg stopped at iteration 1 without converging to the desired
tolerance 1e-006 because a scalar quantity became too small or 
too large to continue computing. 
The iterate returned (number 0) has relative residual 1

しかし、symmlq は、対称な不定行列 A を取り扱うこともできます。

x = symmlq(A,b,1e-6,40);
symmlq converged at iteration 39 to a solution with relative 
residual 1.3e-007

参考

bicg, bicgstab, cgs, lsqr, gmres, minres, pcg, qmr

@ (関数ハンドル), / (スラッシュ)

参考文献

Barrett, R., M. Berry, T. F. Chan, et al., Templates for the Solution of Linear Systems: Building Blocks for Iterative Methods, SIAM, Philadelphia, 1994.

Paige, C. C. and M. A., "Solution of Sparse Indefinite Systems of Linear Equations." SIAM J. Numer. Anal., Vol.12, 1975, pp. 617-629.

symbfact symmmd