ドキュメンテーション

最新のリリースでは、このページがまだ翻訳されていません。 このページの最新版は英語でご覧になれます。

bicgstabl

双共役傾斜安定化 (l) 法

構文

x = bicgstabl(A,b)
x = bicgstabl(afun,b)
x = bicgstabl(A,b,tol)
x = bicgstabl(A,b,tol,maxit)
x = bicgstabl(A,b,tol,maxit,M)
x = bicgstabl(A,b,tol,maxit,M1,M2)
x = bicgstabl(A,b,tol,maxit,M1,M2,x0)
[x,flag] = bicgstabl(A,b,...)
[x,flag,relres] = bicgstabl(A,b,...)
[x,flag,relres,iter] = bicgstabl(A,b,...)
[x,flag,relres,iter,resvec] = bicgstabl(A,b,...)

説明

x = bicgstabl(A,b) は、線形方程式系 A*x=bx について解きます。nn 列の係数行列 A は正方でなければならず、右辺の列ベクトル b は、長さが n でなければなりません。

x = bicgstabl(afun,b) は、行列 A の代わりに、関数ハンドル afun を受け入れます。afun(x) は、ベクトル入力 x を受け入れ、行列とベクトルの積 A*x を返します。以下の構文のすべてにおいて、Aafun で置き換えることが可能です。

x = bicgstabl(A,b,tol) は、このメソッドの許容誤差を指定します。tol が [] の場合、関数 bicgstabl は 既定の設定の 1e-6 を使用します。

x = bicgstabl(A,b,tol,maxit) は、最大反復回数を指定します。maxit が [] の場合は、関数 bicgstabl は既定の設定の min(N,20) を使用します。

x = bicgstabl(A,b,tol,maxit,M) および x = bicgstabl(A,b,tol,maxit,M1,M2) は、前提条件 M または M=M1*M2 を使用し、x 方程式 A*inv(M)*x = b を効率的に解きます。M が [] の場合、前提条件は適用されません。M は、M\x を返す関数ハンドルにもなります。

x = bicgstabl(A,b,tol,maxit,M1,M2,x0) は、初期推定を指定します。x0 が [] の場合、関数 bicgstabl は既定ですべてが 0 のベクトルを使用します。

[x,flag] = bicgstabl(A,b,...) は、収束 flag も返します。

フラグ

収束

0

関数 bicgstabl は希望する許容誤差 tol に反復回数 maxit 以内で収束しました。

1

関数 bicgstabl は、maxit 回の反復に対して、収束しませんでした。

2

前提条件 M は悪条件です。

3

関数 bicgstabl の出力結果に変化が見られません(連続する 2 回の計算で、結果がまったく同じでした)。

4

関数 bicgstabl で計算されたスカラー量の 1 つが計算を続けるには、大きすぎるかまたは小さすぎます。

[x,flag,relres] = bicgstabl(A,b,...) は、相対残差 norm(b-A*x)/norm(b) も返します。flag0 の場合は、relres <= tol になります。

[x,flag,relres,iter] = bicgstabl(A,b,...) は、x が計算された反復回数も返します。ここで、0 <= iter <= maxit です。iterk/4 になる可能性があり (k はいずれかの整数)、特定の 4 分の 1 の反復で収束したことを示します。

[x,flag,relres,iter,resvec] = bicgstabl(A,b,...) は、4 分の 1 の反復ごとに、norm(b-A*x0) を含む残差ノルムのベクトルも返します。

入力または関数と bicgstabl の使用

bicgstabl に直接入力を渡すことができます。

n = 21; 
A = gallery('wilk',n);  
b = sum(A,2);
tol = 1e-12;  
maxit = 15; 
M = diag([10:-1:1 1 1:10]);
x = bicgstabl(A,b,tol,maxit,M);

また、行列-ベクトル積関数を使用することもできます。

function y = afun(x,n)
y = [0; x(1:n-1)] + [((n-1)/2:-1:0)';
(1:(n-1)/2)'].*x+[x(2:n); 0];

前置条件のバックソルブ関数は以下のとおりです。

function y = mfun(r,n)
y = r ./ [((n-1)/2:-1:1)';
1; 
(1:(n-1)/2)'];

関数 bicgstabl への入力として以下の構文を使用します。

x1 = bicgstabl(@(x)afun(x,n),b,tol,maxit,@(x)mfun(x,n));

前提条件子と bicgstabl の使用

この例は、前提条件の使用方法を示します。

479 行 479 列の非対称実スパース行列 west0479 を読み込みます。

load west0479;
A = west0479;

真の解がすべて 1 のベクトルになるように、b を定義します。

b = full(sum(A,2));

許容誤差と最大反復回数を設定します。

tol = 1e-12;
maxit = 20;

bicgstabl を使用して、要求された許容誤差と反復回数で解を求めます。

[x0,fl0,rr0,it0,rv0] = bicgstabl(A,b,tol,maxit);

fl0 は、bicgstabl が要求した反復回数 20 回以内に要求した許容誤差 1e-12 に収束しなかったため、1 となります。実際、bicgstabl の動作が良好でないため、初期推定 (x0 = zeros(size(A,2),1)) が最適解であり、これは it0 = 0 で示される出力です。MATLAB® は残差履歴を rv0 に格納します。

bicgstabl の動作をプロットします。

semilogy(0:0.25:maxit,rv0/norm(b),'-o');
xlabel('Iteration number');
ylabel('Relative residual');

プロットでは、解が収束しないことが示されます。前提条件子を使用して結果を改善できます。

A は非対称であるため、ilu で前提条件子を作成します。

[L,U] = ilu(A,struct('type','ilutp','droptol',1e-5));
Error using ilu
There is a pivot equal to zero.  Consider decreasing the
drop tolerance or consider using the 'udiag' option.

MATLAB は不完全な LU を構成できません。前提条件子では無意味な特異因子が生成されるためです。

エラー メッセージで示されるように、減少させた棄却許容誤差でもう一度試すことができます。

[L,U] = ilu(A,struct('type','ilutp','droptol',1e-6));
[x1,fl1,rr1,it1,rv1] = bicgstabl(A,b,tol,maxit,L,U);

bicgstabl が相対残差 1.0257e-015 (rr1 の値) を導き出すため、fl1 は 0 です。相対残差は、棄却許容誤差 1e-6 で不完全 LU 分解を前提条件としたときの 6 回目の反復 (it1 の値) で指定の許容誤差 1e-12 より小さくなります。出力 rv1(1)norm(b)、出力 rv1(9)norm(b-A*x2) になります。これは、bicgstabl が 4 分の 1 の反復を使用するためです。

初期推定 (反復回数が 0) から始まる各反復での相対残差をプロットして、関数 bicgstabl の進行状況を確認できます。

semilogy(0:0.25:it1,rv1/norm(b),'-o');

h = gca;
h.XTick = 0:0.25:it1;

xlabel('Iteration number');
ylabel('Relative residual');

拡張機能