非対称実スパースソルバー - サンプルソースコード Fortran

ホーム > 製品 > nAG数値計算ライブラリ > nAG Fortranライブラリ > サンプルソースコード集 > 非対称実スパースソルバー

Keyword: スパースソルバー, sparse solver, RGMRES, 疎行列, SSOR

概要

本サンプルは以下に示す非対称実スパース連立方程式を解くFortranによるプログラムです。連立方程式は係数行列と右辺ベクトルとして与えます。　

係数行列（スパース）　右辺ベクトル

以下のプログラム例ではRGMRES法 (restarted generalized minimal residual method) を用いていますが、その他にCGS法 (conjugate gradient squared method)、及びBICGSTAB法 (bi-conjugate gradient stabilised method) に対応しています。（入力データの'RGMRES'の部分を'CGS'、'BiCGSTAB'にそれぞれ置き換える事で手法を選ぶことが可能です）

またプリコンディショナーにはSSOR (symmetric successive-over-relaxation) を用いていますが、その他にJacobi の指定及びプリコンディショナーを使用しない設定が可能です。（入力データの'S'の部分をそれぞれ'J'、'N'に置き換えることで指定可能です）

※本サンプルはnAG Fortranライブラリに含まれるルーチン f11def()のExampleコードです。本サンプル及びルーチンの詳細情報はf11def のマニュアルページをご参照ください。
ご相談やお問い合わせはこちらまで

入力データ

（本ルーチンの詳細はf11def のマニュアルページを参照）

このデータをダウンロード

F11DEF Example Program Data
  5   1               N, M
 16                   NNZ
  'RGMRES' 'S'        METHOD, PRECON
  1.05                OMEGA
  1.D-10 1000         TOL, MAXITN
  2.   1    1
  1.   1    2
 -1.   1    4
 -3.   2    2
 -2.   2    3
  1.   2    5
  1.   3    1
  5.   3    3
  3.   3    4
  1.   3    5
 -2.   4    1
 -3.   4    4
 -1.   4    5
  4.   5    2
 -2.   5    3
 -6.   5    5         A(I), IROW(I), ICOL(I), I=1,...,NNZ
  0.  -7.  33.
-19. -28.             B(I), I=1,...,N
  0.   0.   0.
  0.   0.             X(I), I=1,...,N

１行目はタイトル行で読み飛ばされます。
２行目は係数行列の大きさ（n=5,5x5なので5)、RGMRES法の再開サブスペース空間の次元 m (=1) を指定しています。数字より後ろは行末まで読み飛ばされます。
３行目では係数行列（スパース）に含まれる非ゼロ要素の数(nnz=16)を与えます。数字より後ろは行末まで読み飛ばされます。
４行目は用いる計算手法 (method='RGMRES') とプリコンディショナー (precon='S'（SSOR）) を指定しています。プリコンディショナーの指定より後ろは行末まで読み飛ばされます。
※ここで他の手法（'CGS' もしくは 'BiCGSTAB'）やプリコンディショナー（'J'（ヤコビ）もしくは 'N'（プリコンディショナー無し））の指定も可能です。
５行目ではSSORプリコンディショナーのリラックスパラメータであるω (=1.05) を与えています。SSOR以外のプリコンディショナーを用いる場合にはこの数値は無視されます。数字より後ろは行末まで読み飛ばされます。
６行目では計算精度 tol (=1.D-10) 、最大反復数 maxitn (=1000) をそれぞれ与えています。２つの数値より後ろは行末まで読み飛ばされます。
７行目～２２行目にはスパース行列の非ゼロ係数の値を指定しています。それぞれの行は値(a)、行番号(irow)、列番号(icol)で全部で16個（３行目で与えた値）の非ゼロ係数が指定されています。行番号及び列番号は1から最大5（２行目で指定した値）までの値を取ります。　例えば行列の一番左の一番上の値が3.5であった場合には
3.5　　1　　1
と与えます。
２３行目と２４行目は右辺ベクトルの値(b)を指定しています。全部で５個（２行目で指定した値）の数値が指定されています。
２５行目と２６行目は解の初期値(x)を与えています。全部で５個（２行目で指定した値）の数値が指定されています。ここでは初期値としてすべて 0 を与えています。

出力結果

（本ルーチンの詳細はf11def のマニュアルページを参照）

この出力例をダウンロード

 F11DEF Example Program Results

 Converged in        13 iterations
 Final residual norm =       5.087E-09

            X
       1.0000E+00
       2.0000E+00
       3.0000E+00
       4.0000E+00
       5.0000E+00

１行目はタイトルです。
３行目に解を求める際の反復数が示されます。（この例では１３回）
４行目は求まった解の残差ノルムが出力されます。
６行目は解のタイトルである x が出力されます。
７行目～１１行目に求まった解が表示されます。

ソースコード

（本ルーチンの詳細はf11def のマニュアルページを参照）

※本サンプルソースコードは科学技術・統計計算ライブラリである「nAG Fortranライブラリ」のルーチンを呼び出します。
サンプルのコンパイル及び実行方法

このソースコードをダウンロード

    PROGRAM f11defe

!      F11DEF Example Program Text

!      Mark 23 Release. nAG Copyright 2011.

!      .. Use Statements ..
       USE nag_library, ONLY : f11def, nag_wp
!      .. Implicit None Statement ..
       IMPLICIT NONE
!      .. Parameters ..
       INTEGER, PARAMETER              :: nin = 5, nout = 6
!      .. Local Scalars ..
       REAL (KIND=nag_wp)              :: omega, rnorm, tol
       INTEGER                         :: i, ifail, itn, l, lwork, m, maxitn,  &
                                          n, nnz
       CHARACTER (8)                   :: method
       CHARACTER (1)                   :: precon
!      .. Local Arrays ..
       REAL (KIND=nag_wp), ALLOCATABLE :: a(:), b(:), work(:), x(:)
       INTEGER, ALLOCATABLE            :: icol(:), irow(:), iwork(:)
!      .. Intrinsic Functions ..
       INTRINSIC                          max
!      .. Executable Statements ..
       WRITE (nout,*) 'F11DEF Example Program Results'
       WRITE (nout,*)
!      Skip heading in data file
       READ (nin,*)

!      Read algorithmic parameters
       READ (nin,*) n, m
       READ (nin,*) nnz
       READ (nin,*) method, precon
       l = n
       IF (precon=='N' .OR. precon=='n') l = 0
       lwork = max(4*n+m*(m+n+5)+l+101,8*n+l+100,2*n*(m+3)+m*(m+2)+l+100, &
          11*n+l+100)

       ALLOCATE (a(nnz),b(n),work(lwork),x(n),icol(nnz),irow(nnz), &
          iwork(2*n+1))
       READ (nin,*) omega
       READ (nin,*) tol, maxitn

!      Read the matrix A

       DO i = 1, nnz
          READ (nin,*) a(i), irow(i), icol(i)
       END DO

!      Read right-hand side vector b and initial approximate solution x

       READ (nin,*) b(1:n)
       READ (nin,*) x(1:n)

!      Solve Ax = b using F11DEF

!      ifail: behaviour on error exit
!             =0 for hard exit, =1 for quiet-soft, =-1 for noisy-soft
       ifail = 0
       CALL f11def(method,precon,n,nnz,a,irow,icol,omega,b,m,tol,maxitn,x, &
          rnorm,itn,work,lwork,iwork,ifail)

       WRITE (nout,'(A,I10,A)') ' Converged in', itn, ' iterations'
       WRITE (nout,'(A,1P,E16.3)') ' Final residual norm =', rnorm
       WRITE (nout,*)

!      Output x

       WRITE (nout,*) '           X'
       WRITE (nout,'(1X,1P,E16.4)') x(1:n)

    END PROGRAM f11defe

非対称実スパースソルバー (Nonsymmetric Sparse Solver)

Fortranによるサンプルソースコード : 使用ルーチン名：f11def

概要

入力データ

出力結果

ソースコード