XLPack 6.1
Excel VBA 数値計算ライブラリ(基本モジュール) リファレンスマニュアル
読み取り中…
検索中…
一致する文字列を見つけられません

◆ Dgels()

Sub Dgels ( Trans As  String,
M As  Long,
N As  Long,
A() As  Double,
B() As  Double,
Info As  Long,
Optional Nrhs As  Long = 1 
)

優決定または劣決定系連立一次方程式 Ax = b の解 (フルランク)

目的
本ルーチンはm×n実行列Aあるいはその転置からなる優決定または劣決定系連立1次方程式をAのQRあるいはLQ分解を用いて解く. Aはフルランクであること.

下記オプションが提供される:

  1. Trans = "N" かつ m >= n の場合: 優決定系の最小二乗解を求める. すなわち, 次の最小二乗問題を解く.
    || B - A*X || を最小化する.
  2. Trans = "N" かつ m < n の場合: 劣決定系 A * X = B の最小ノルム解を求める.
  3. Trans = "T" かつ m >= n の場合: 劣決定系 A^T * X = B の最小ノルム解を求める.
  4. Trans = "T" かつ m < n の場合: 優決定系の最小二乗解を求める. すなわち, 次の最小二乗問題を解く.
    || B - A^T*X || を最小化する.

いくつかの右辺ベクトル b および解ベクトル x を1回の呼び出しで扱うことができる. これらのベクトルは, m×nrhs右辺行列Bおよびn×nrhs解行列Xの列として格納される.
引数
[in]Trans= "N": 連立一次方程式はAからなる.
= "T": 連立一次方程式はA^Tからなる.
[in]M行列 A の行数. (M >= 0) (M = 0 の場合, 処理を行わずに戻る)
[in]N行列 A の列数. (N >= 0) (N = 0 の場合, 処理を行わずに戻る)
[in,out]A()配列 A(LA1 - 1, LA2 - 1) (LA1 >= M, LA2 >= N)
[in] M×N行列 A.
[out] M >= N: A()はDgeqrfが返すQR分解結果により上書きされる.
  M < N: A()はDgelqfが返すLQ分解結果により上書きされる.
[in,out]B()配列 B(LB1 - 1, LB2 - 1) (LB1 >= max(M, N), LB2 >= Nrhs) (2次元配列) または B(LB - 1) (LB >= max(M, N), Nrhs = 1) (1次元配列)
[in] 右辺ベクトルが列ごとに格納された行列 B. Bは, M×Nrhs行列(Trans = "N"の場合), または, N×Nrhs行列(Trans = "T"の場合)である.
[out] Info = 0の場合、B()は列ごとに格納された解ベクトルにより上書きされる:
  Trans = "N"かつM >= N: B()の行0〜N-1に最小二乗解ベクトルが入る. 各列の解の残差二乗和はその列の行N〜M-1の要素の二乗和で与えられる.
  Trans = "N"かつM < N: B()の行0〜N-1に最小ノルム解ベクトルが入る.
  Trans = "T"かつM >= N: B()の行0〜M-1に最小ノルム解ベクトルが入る.
  Trans = "T"かつM < N: B()の行0〜M-1に最小二乗解ベクトルが入る. 各列の解の残差二乗和はその列の行M〜N-1の要素の二乗和で与えられる.
[out]Info= 0: 正常終了.
= -1: パラメータ Trans の誤り. (Trans <> "T"および"N")
= -2: パラメータ M の誤り. (M < 0)
= -3: パラメータ N の誤り. (N < 0)
= -4: パラメータ A() の誤り.
= -5: パラメータ B() の誤り.
= -7: パラメータ Nrhs の誤り. (Nrhs < 0)
= i > 0: Aの分解結果の三角行列のi番目の対角要素が0になった. 従って, Aはフルランクではない. 最小二乗解は求めることができなかった.
[in]Nrhs(省略可)
右辺の数, すなわち, 行列BおよびXの列数. (Nrhs >= 0) (Nrhs = 0 の場合, 処理を行わずに戻る) (省略時 = 1)
出典
LAPACK
使用例
優決定系連立1次方程式 Ax = B の最小二乗解を求める. また, 分散を求める. ただし,
( -1.06 0.48 -0.04 )
A = ( -1.19 0.73 -0.24 )
( 1.97 -0.89 0.56 )
( 0.68 -0.53 0.08 )
( 0.3884 )
B = ( 0.1120 )
( -0.3644 )
( -0.0002 )
とする.
Sub Ex_Dgels()
Const M = 4, N = 3
Dim A(M - 1, N - 1) As Double, B(M - 1) As Double, Ci(N - 1) As Double
Dim Info As Long
A(0, 0) = -1.06: A(0, 1) = 0.48: A(0, 2) = -0.04
A(1, 0) = -1.19: A(1, 1) = 0.73: A(1, 2) = -0.24
A(2, 0) = 1.97: A(2, 1) = -0.89: A(2, 2) = 0.56
A(3, 0) = 0.68: A(3, 1) = -0.53: A(3, 2) = 0.08
B(0) = 0.3884: B(1) = 0.112: B(2) = -0.3644: B(3) = -0.0002
Call Dgels("N", M, N, A(), B(), Info)
If Info <> 0 Then
Debug.Print "Error in Dgels: Info =", Info
Exit Sub
End If
Debug.Print "X =", B(0), B(1), B(2)
Call Dgecov(0, N, A(), Ci(), Info)
Debug.Print "Var =", Ci(0), Ci(1), Ci(2)
Debug.Print "Info =", Info
End Sub
実行結果
X = -0.820000000000001 -0.940000000000001 0.740000000000001
Var = 6.46959967542666 16.7350408218919 18.7177532773229
Info = 0