◆ Lmdif1_r()

Sub Lmdif1_r	(	M As	Long,
		N As	Long,
		X() As	Double,
		Fvec() As	Double,
		Tol As	Double,
		Info As	Long,
		XX() As	Double,
		YY() As	Double,
		IRev As	Long,
		Optional Info2 As	Long
	)

非線形最小二乗法 (レーベンバーグ・マルカート法) (ヤコビ行列不要) (シンプルドライバ) (リバースコミュニケーション版)

目的: 本ルーチンはM個のN変数非線形関数の二乗和の最小点をレーベンバーグ・マルカート法により求める.
min Σfi(x1, x2, ..., xN)^2 (ただし, Σは i = 1 〜 M)

IRevに従って関数値をユーザーが与える必要がある. ヤコビ行列はルーチン内で前進差分により計算されるため, ユーザーがヤコビ行列を求める必要はない.

Lmdif1_rは, FTol = Tol, XTol = Tol, GTol = 0, Mode = 1 としてLmdif_rを呼び出すのと同じである.

引数

[in]	M	関数の数. (M > 0)
[in]	N	変数の数. (0 < N <= M)
[in,out]	X()	配列 X(LX - 1) (LX >= N) [in] 初期近似解. [out] 求められた解ベクトル.
[out]	Fvec()	配列 Fvec(LFvec - 1) (LFvec >= M) 求められた解ベクトルX()における関数値.
[in]	Tol	二乗和および解の相対誤差の許容値. (Tol >= 0)
[out]	Info	= 0: 正常終了. (サブコードをInfo2に返す) = -1: パラメータ M の誤り. (M < N) = -2: パラメータ N の誤り. (N <= 0) = -3: パラメータ X() の誤り. (配列X()の大きさが不足) = -4: パラメータ Fvec() の誤り. (配列Fvec()の大きさが不足) = -5: パラメータ Tol の誤り. (Tol < 0) = -7: パラメータ XX() の誤り. (配列XX()の大きさが不足) = -8: パラメータ YY() の誤り. (配列YY()の大きさが不足) = -9: パラメータ IRev の誤り. = 1: 関数呼び出し(IFlag = 1または2)回数がMaxfevに達した. = 2: 残差二乗和が減少しなくなった. (Tolが小さすぎる) = 3: 解が改善されなくなった. (Tolが小さすぎる)
[out]	XX()	配列 XX(LXX - 1) (LXX >= N) IRev = 1 または 2の場合, 関数値を求めるべき点を返す.
[in]	YY()	配列 YY(LYY - 1) (LYY >= M) IRev = 1 または 2の場合, 再呼び出し時に関数値fi(XX()) (i = 1〜M)を与えること.
[in,out]	IRev	リバースコミュニケーションの制御変数. [in] 最初の呼び出し時に 0 に設定しておくこと. 2回目以降の呼び出し時には値を変更してはならない. [out] 0 以外の時には下記処理を行ってから再び本ルーチンを呼び出すこと. = 0: 処理終了. 正常終了かどうかはInfoをチェックすること. = 1 または 2: XX()における関数値を求めYY()に設定すること. YY()以外の変数を変更してはならない.
[out]	Info2	(省略可) Info = 0 で戻ったときのサブコード. = 1: 二乗和の相対減少値およびその予測値の両方がTol以下になった. = 2: 連続する2回の反復の相対誤差がTol以下になった. = 3: 上の2つ共満たした. = 4: Fvecがヤコビ行列の列とマシンイプシロン以内で直交した.

出典: netlib/minpack

使用例: 次のデータをモデル関数 f(x) = c1*(1 - exp(-c2*x)) で近似する. 2つのパラメータc1, c2は非線形最小二乗法により定める.
f(x) x

10.07 77.6

29.61 239.9

50.76 434.8

81.78 760.0

初期値は, c1 = 500, c2 = 0.0001 とする.
Sub FLmdif(M As Long, N As Long, X() As Double, Fvec() As Double, IFlag As Long)

Dim Xdata(3) As Double, Ydata(3) As Double, I As Long

Ydata(0) = 10.07: Xdata(0) = 77.6

Ydata(1) = 29.61: Xdata(1) = 239.9

Ydata(2) = 50.76: Xdata(2) = 434.8

Ydata(3) = 81.78: Xdata(3) = 760

If IFlag = 1 Or IFlag = 2 Then

For I = 0 To M - 1

Fvec(I) = Ydata(I) - X(0) * (1 - Exp(-Xdata(I) * X(1)))

Next

End If

End Sub

Sub Ex_Lmdif1_r()

Const M As Long = 4, N As Long = 2

Dim X(N - 1) As Double, Fvec(M - 1) As Double, Tol As Double, Info As Long

Dim XX(N - 1) As Double, YY(M - 1) As Double, IRev As Long, IFlag As Long

Tol = 0.00000001 '1.0e-8

X(0) = 500: X(1) = 0.0001

IRev = 0

Do

Call Lmdif1_r(M, N, X(), Fvec(), Tol, Info, XX(), YY(), IRev)

If IRev = 1 Or IRev = 2 Then

IFlag = IRev

Call FLmdif(M, N, XX(), YY(), IFlag)

End If

Loop While IRev <> 0

Debug.Print "C1, C2 =", X(0), X(1)

Debug.Print "Info =", Info

End Sub

Lmdif1_r
Sub Lmdif1_r(M As Long, N As Long, X() As Double, Fvec() As Double, Tol As Double, Info As Long, XX() As Double, YY() As Double, IRev As Long, Optional Info2 As Long)
非線形最小二乗法 (レーベンバーグ・マルカート法) (ヤコビ行列不要) (シンプルドライバ) (リバースコミュニケーション版)

実行結果: C1, C2 = 241.084897132837 5.44942231312974E-04

Info = 0