Scilab Website | Contribute with GitLab | Mailing list archives | ATOMS toolboxes
Scilab Online Help
2023.0.0 - 日本語


qmr

プリコンディショナ付きのQuasi Minimal Residual法

呼び出し手順

[x, flag, err, iter, res] = qmr(A, b, x0, M1, M2, maxi, tol)
[x, flag, err, iter, res] = qmr(A,Ap,b,x0,M1,M1p,M2,M2p,maxi,tol) // deprecated

Parameters

A
Square dense or sparse matrix of size n-by-n, or function:

  • If A is a function which returns A*x, it must have the following header :

    function y = A ( x )
  • If A is a function which returns A*x or A'*x depending on a option t, it must have the following header :

    function y = A(x, t)

    • If t = "notransp": the function returns A*x.
    • If t = "transp": the function returns A'*x.

Ap
function returning A'*x. It must have the following header :

function y = Ap(x)

b
右辺ベクトル.

x0
初期推定ベクトル (デフォルト: zeros(n,1)).

M1
左プリコンディショナ: 行列またはM1*xを返す関数 (前者のデフォルト値: eye(n,n)). If M1 is a function, it returns either,
  • only M1*x, or
  • M1*x or M1'*x, depending on t.

M1p
M1が関数の場合のみ指定する 必要があります. この場合, M1p は M1'*xを返す関数です.

M2
右プリコンディショナ: 行列またはM2*xを 返す関数 (前者のデフォルト値: eye(n,n)) If M2 is a function, it returns either
  • only M2*x, or
  • M2*x or M2'*x depending on an option t.

M2p
M2が関数の場合のみ指定する 必要があります. この場合, M2pはM2'*xを返す関数です.

maxi
最大反復回数 (デフォルト: n).

tol
許容誤差 (デフォルト: 1000*%eps).

x
解ベクトル.

flag
  • flag=0: gmr は maxi回の反復の間に 許容誤差内に収束しました
  • flag=1: 指定したmaxi回の反復の間に 収束しませんでした
  • -7 < flag < 0: A breakdown occurred because one of the scalar quantities calculated was equal to zero.

res
残差ベクトル.

err
最終残差ノルム.

iter
実行した反復回数.

説明

プリコンディショナ付きのQuasi Minimal Residual法により, 線形システムAx=bを解きます.

例

If A is a matrix:

A = [ 94   0   0   0    0   28  0   0   32  0
       0   59  13  5    0   0   0   10  0   0
       0   13  72  34   2   0   0   0   0   65
       0   5   34  114  0   0   0   0   0   55
       0   0   2   0    70  0   28  32  12  0
       28  0   0   0    0   87  20  0   33  0
       0   0   0   0    28  20  71  39  0   0
       0   10  0   0    32  0   39  46  8   0
       32  0   0   0    12  33  0   8   82  11
       0   0   65  55   0   0   0   0   11  100];
b = ones(10,1);
[x,flag,err,iter,res] = qmr(A, b)

[x,flag,err,iter,res] = qmr(A, b, zeros(10,1), eye(10,10), eye(10,10), 10, 1d-12)

If A is a function:

function y=Atimesx(x, t)
    A = [ 94  0   0   0    0   28  0   0   32  0
          0   59  13  5    0   0   0   10  0   0
          0   13  72  34   2   0   0   0   0   65
          0   5   34  114  0   0   0   0   0   55
          0   0   2   0    70  0   28  32  12  0
          28  0   0   0    0   87  20  0   33  0
          0   0   0   0    28  20  71  39  0   0
          0   10  0   0    32  0   39  46  8   0
          32  0   0   0    12  33  0   8   82  11
          0   0   65  55   0   0   0   0   11  100];
     if (t == 'notransp') then
        y = A*x;
    elseif (t ==  'transp') then
        y = A'*x;
    end
endfunction
b = ones(10,1);

[x,flag,err,iter,res] = qmr(Atimesx, b)

[x,flag,err,iter,res] = qmr(Atimesx, b, zeros(10,1), eye(10,10), eye(10,10), 10, 1d-12)

OR

function y=funA(x)
    A = [ 94  0   0   0    0   28  0   0   32  0
          0   59  13  5    0   0   0   10  0   0
          0   13  72  34   2   0   0   0   0   65
          0   5   34  114  0   0   0   0   0   55
          0   0   2   0    70  0   28  32  12  0
          28  0   0   0    0   87  20  0   33  0
          0   0   0   0    28  20  71  39  0   0
          0   10  0   0    32  0   39  46  8   0
          32  0   0   0    12  33  0   8   82  11
          0   0   65  55   0   0   0   0   11  100];
     y = A*x
endfunction

function y=funAp(x)
    A = [ 94  0   0   0    0   28  0   0   32  0
          0   59  13  5    0   0   0   10  0   0
          0   13  72  34   2   0   0   0   0   65
          0   5   34  114  0   0   0   0   0   55
          0   0   2   0    70  0   28  32  12  0
          28  0   0   0    0   87  20  0   33  0
          0   0   0   0    28  20  71  39  0   0
          0   10  0   0    32  0   39  46  8   0
          32  0   0   0    12  33  0   8   82  11
          0   0   65  55   0   0   0   0   11  100];
     y = A'*x
endfunction

b = ones(10,1);

[x,flag,err,iter,res] = qmr(funA, funAp, b)

[x,flag,err,iter,res] = qmr(funA, funAp, b, zeros(10,1), eye(10,10), eye(10,10), 10, 1d-12)

If A is a matrix, M1 and M2 are functions:

A = [ 94   0   0   0    0   28  0   0   32  0
       0   59  13  5    0   0   0   10  0   0
       0   13  72  34   2   0   0   0   0   65
       0   5   34  114  0   0   0   0   0   55
       0   0   2   0    70  0   28  32  12  0
       28  0   0   0    0   87  20  0   33  0
       0   0   0   0    28  20  71  39  0   0
       0   10  0   0    32  0   39  46  8   0
       32  0   0   0    12  33  0   8   82  11
       0   0   65  55   0   0   0   0   11  100];

b = ones(10,1);

function y=M1timesx(x, t)
    M1 = eye(10,10);
    if(t=="notransp") then
        y = M1*x;
    elseif (t=="transp") then
        y = M1'*x;
    end
endfunction

function y=M2timesx(x, t)
    M2 = eye(10,10);
    if(t=="notransp") then
        y = M2*x;
    elseif (t=="transp") then
        y = M2'*x;
    end
endfunction

[x, flag, err, iter, res] = qmr(A, b, zeros(10,1), M1timesx, M2timesx, 10, 1d-12)

OR

A = [ 94   0   0   0    0   28  0   0   32  0
       0   59  13  5    0   0   0   10  0   0
       0   13  72  34   2   0   0   0   0   65
       0   5   34  114  0   0   0   0   0   55
       0   0   2   0    70  0   28  32  12  0
       28  0   0   0    0   87  20  0   33  0
       0   0   0   0    28  20  71  39  0   0
       0   10  0   0    32  0   39  46  8   0
       32  0   0   0    12  33  0   8   82  11
       0   0   65  55   0   0   0   0   11  100];

b = ones(10,1);

function y=funM1(x)
    M1 = eye(10,10);
    y = M1*x;
endfunction

function y=funM1p(x)
    M1 = eye(10,10);
    y = M1'*x;
endfunction

function y=funM2(x)
    M2 = eye(10,10);
    y = M2*x;
endfunction

function y=funM2p(x)
    M2 = eye(10,10);
    y = M2'*x;
endfunction

[x,flag,err,iter,res] = qmr(A, b, zeros(10,1), funM1, funM1p, funM2, funM2p, 10, 1d-12)

If A, M1, M2 are functions:

// See functions defined above in previous examples. Then,

[x,flag,err,iter,res] = qmr(funA, funAp, b, zeros(10,1), funM1, funM1p, funM2, funM2p, 10, 1d-12)
// or
[x,flag,err,iter,res] = qmr(Atimesx, b, zeros(10,1), M1timesx, M2timesx, 10, 1d-12)

参照

  • gmres — Generalized Minimum RESidual 法
  • conjgrad — 共役勾配ソルバ

履歴

バージョン記述
5.4.0 qmr(A, Ap) のコールは廃止されました. qmr(A) を代わりに使用してください.
Report an issue
<< gmres Linear Equations (Iterative Solvers) Sparse Matrix Manipulation >>
Copyright (c) 2022-2024 (Dassault Systèmes)
Copyright (c) 2017-2022 (ESI Group)
Copyright (c) 2011-2017 (Scilab Enterprises)
Copyright (c) 1989-2012 (INRIA)
Copyright (c) 1989-2007 (ENPC)
with contributors 		Last updated:
Tue Mar 07 09:28:47 CET 2023