◆ ZSor()

Sub ZSor	(	N As	Long,
		Val() As	Complex,
		Ptr() As	Long,
		Ind() As	Long,
		B() As	Complex,
		X() As	Complex,
		Optional Info As	Long,
		Optional Iter As	Long,
		Optional Res As	Double,
		Optional Format As	Long = `0`,
		Optional Omega As	Double = `1.5`,
		Optional MaxIter As	Long = `500`,
		Optional Tol As	Double = `1.0E-10`,
		Optional Base As	Long = `-1`
	)

逐次的過剰緩和(SOR)法による連立一次方程式 Ax = b の解 (複素行列) (ドライバ)

目的: 反復法(逐次的過剰緩和(SOR)法)により連立一次方程式 Ax = b の解を求める. 行列 A は CSR または CSC 形式で表される.

引数

[in]	N	行列 A の次数. (N >= 0) (N = 0 の場合, 処理を行わずに戻る)
[in]	Val()	配列 Val(LVal - 1) (LVal >= Nnz) (Nnz は行列 A の非ゼロ要素数) 行列 A の非ゼロ要素の値.
[in]	Ptr()	配列 Ptr(LPtr - 1) (LPtr >= N + 1) 行列 A の列ポインタ(CSC の場合)または行ポインタ(CSR の場合).
[in]	Ind()	配列 Ind(LInd - 1) (LInd >= Nnz) 行列 A の行インデクス(CSC の場合)または列インデクス(CSR の場合).
[in]	B()	配列 B(LB - 1) (LB >= N) 右辺ベクトル b.
[in,out]	X()	配列 X(LX - 1) (LX >= N) [in] 解の初期推定値. [out] 求められた近似解.
[out]	Info	(省略可) = 0: 正常終了. = i < 0: (-i)番目の入力パラメータの誤り. = 11: 最大反復回数を超えた. = 12: 行列が特異である(対角要素が0). = 15: 内部処理のエラー.
[out]	Iter	(省略可) 収束時の反復回数.
[out]	Res	(省略可) 最終的な残差ノルム norm(b - A*x) の値.
[in]	Format	(省略可) 行列の格納形式. (省略時 = 0) = 0: CSR 形式. = 1: CSC 形式.
[in]	Omega	(省略可) 緩和パラメータω. (0 < ω < 2) (省略時 = 1.5)
[in]	MaxIter	(省略可) 最大反復回数. (MaxIter > 0) (省略時 = 500) MaxIter <= 0 であれば省略時の既定値とみなす.
[in]	Tol	(省略可) 収束判定基準値. (省略時 = 1.0e-10) norm(b - Ax) <= Tolnorm(b) であれば収束とみなす. Tol < eps (マシンイプシロン) であれば Tol = eps とみなす.
[in]	Base	(省略可) Ptr() および Ind() のインデクス形式. = 0: 0-ベース(C形式): 開始インデクス値が 0. = 1: 1-ベース(Fortran形式): 開始インデクス値が 1. (省略時: Ptr(0) = 1 であれば 1, そうでなければ 0 とみなす)

使用例: 連立一次方程式 Ax = B を解く. ただし,
( 4 0 1 0.7 0 )

( 2i 4 0 1 0.7 )

A = ( 0 2i 4 0 1 )

( 0 0 2i 4 0 )

( 0 0 0 2i 4 )

( 5.7 )

( 5.7 + 2i )

B = ( 5 + 2i )

( 4 + 2i )

( 4 + 2i )

とする.
Sub Ex_ZSor()

Const N = 5, Nnz = 14, Omega = 1.05

Dim A(Nnz - 1) As Complex, Ia(N) As Long, Ja(Nnz - 1) As Long

Dim B(N - 1) As Complex, X(N - 1) As Complex

Dim Iter As Long, Res As Double, Info As Long

A(0) = Cmplx(4): A(1) = Cmplx(1): A(2) = Cmplx(0.7): A(3) = Cmplx(0, 2): A(4) = Cmplx(4): A(5) = Cmplx(1): A(6) = Cmplx(0.7): A(7) = Cmplx(0, 2): A(8) = Cmplx(4): A(9) = Cmplx(1): A(10) = Cmplx(0, 2): A(11) = Cmplx(4): A(12) = Cmplx(0, 2): A(13) = Cmplx(4)

Ia(0) = 0: Ia(1) = 3: Ia(2) = 7: Ia(3) = 10: Ia(4) = 12: Ia(5) = 14

Ja(0) = 0: Ja(1) = 2: Ja(2) = 3: Ja(3) = 0: Ja(4) = 1: Ja(5) = 3: Ja(6) = 4: Ja(7) = 1: Ja(8) = 2: Ja(9) = 4: Ja(10) = 2: Ja(11) = 3: Ja(12) = 3: Ja(13) = 4

B(0) = Cmplx(5.7): B(1) = Cmplx(5.7, 2): B(2) = Cmplx(5, 2): B(3) = Cmplx(4, 2): B(4) = Cmplx(4, 2)

Call ZSor(N, A(), Ia(), Ja(), B(), X(), Info, Iter, Res, , Omega)

Debug.Print "X ="

Debug.Print "(" + CStr(Creal(X(0))) + "," + CStr(Cimag(X(0))) + ")"

Debug.Print "(" + CStr(Creal(X(1))) + "," + CStr(Cimag(X(1))) + ")"

Debug.Print "(" + CStr(Creal(X(2))) + "," + CStr(Cimag(X(2))) + ")"

Debug.Print "(" + CStr(Creal(X(3))) + "," + CStr(Cimag(X(3))) + ")"

Debug.Print "(" + CStr(Creal(X(4))) + "," + CStr(Cimag(X(4))) + ")"

Debug.Print "Iter =" + Str(Iter) + ", Res =" + Str(Res) + ", Info =" + Str(Info)

End Sub

Cmplx
Function Cmplx(R As Double, Optional I As Double=0) As Complex
複素数の作成

Cimag
Function Cimag(A As Complex) As Double
複素数の虚数部

Creal
Function Creal(A As Complex) As Double
複素数の実数部

ZSor
Sub ZSor(N As Long, Val() As Complex, Ptr() As Long, Ind() As Long, B() As Complex, X() As Complex, Optional Info As Long, Optional Iter As Long, Optional Res As Double, Optional Format As Long=0, Optional Omega As Double=1.5, Optional MaxIter As Long=500, Optional Tol As Double=1.0E-10, Optional Base As Long=-1)
逐次的過剰緩和(SOR)法による連立一次方程式 Ax = b の解 (複素行列) (ドライバ)

実行結果: X =

( .999999999924489,-7.08452827945102E-11)

( .999999999965658, 7.93322103076301E-11)

( 1.00000000005495, 8.34202476797265E-12)

( .999999999994847,-2.33985620845386E-11)

( .999999999991225, 5.71383743476044E-12)

Iter = 11, Res = 4.1830004025201E-10, Info = 0