Mnfb Mnfb_r Mngb Mngb_r Mnhb Mnhb_r

Mnhb_r() Sub Mnhb_r ( N As  Long, X() As  Double, B() As  Double, Info As  Long, YY As  Double, YYp() As  Double, YYpd() As  Double, IRev As  Long, Optional Iout As  Long = 0, Optional Info2 As  Long, Optional NFcall As  Long, Optional NGcall As  Long, Optional Niter As  Long, Optional Fval As  Double, Optional NFGcal As  Long, Optional Rtol As  Double = -1, Optional Atol As  Double = -1, Optional MaxFcall As  Long = 0, Optional MaxIter As  Long = 0, Optional Dtype As  Long = 0, Optional Dfac As  Double = -1, Optional Dtol As  Double = 0, Optional D0 As  Double = 0, Optional Tuner1 As  Double = -1, Optional Xctol As  Double = -1, Optional Xftol As  Double = -1, Optional Lmax0 As  Double = -1, Optional Lmaxs As  Double = -1, Optional Sctol As  Double = -1, Optional Bias As  Double = -1  ) 多変数非線形関数の最小点 (信頼領域法) (単純制約付き) (二階微分要) (リバースコミュニケーション版) 目的 本ルーチンは非線形関数 f(x1, x2, ..., xn) (2階連続的微分可能な実数関数) の局所的最小点 (xs1, xs2, ..., xsn) を次の制約条件のもとで求める(b0iとb1iはそれぞれ上限値と下限値). b0i <= xsi <= b1i (i = 1 〜 n) ユーザーが解析的に計算した目的関数の微分係数およびヘッセ行列が使われる. 探索アルゴリズムはダブルドッグレッグ信頼領域法を使用する. Mnhb_rはMnhbのリバースコミュニケーション版である. 引数 [in] N 問題の次数(パラメータ数). (N > 0) [in,out] X() 配列 X(LX - 1) (LX >= N) [in] 解ベクトルの初期推定値. [out] IRev = 0: 求められた解ベクトル.   IRev = 1, 2: 関数値または微分係数およびヘッセ行列を求める点.   IRev = 3: 最新の近似解. [in] B() 配列 B(LB1 - 1, LB2 - 1) (LB1 = 2, LB2 >= N) 解ベクトル X の境界制約条件.   B(0, I) <= X(I) <= B(1, I) (I = 0 〜 N-1) [out] Info = 0: 正常終了. (サブコードをInfo2に返す) = -1: パラメータ N の誤り. (N < 1) = -2: パラメータ X() の誤り. = -3: パラメータ B() の誤り. = -6: パラメータ YYp() の誤り. = -7: パラメータ YYpd() の誤り. = 7: 特異収束. (近傍のヘッセ行列が特異になった) = 8: 誤収束. (誤った点での収束と思われる. 目標精度が小さすぎる可能性がある) = 9: 関数評価回数の最大値を超えた. = 10: 反復回数の最大値を超えた. = 63: X()の初期点においてF(X)を求めることができない. = 65: X()において微分係数を求めることができない. [in,out] YY [in] IRev = 1の場合, 再呼び出し時にX()における関数値を与えること. [out] IRev = 3の場合, 中間結果の出力のため現在のX()における関数値を返す. [in,out] YYp() 配列 YYp(LYYp - 1) (LYYp >= N) [in] IRev = 2の場合, 再呼び出し時にX()における微分係数を与えること. [out] IRev = 3の場合, 中間結果の出力のため現在のX()における微分係数値を返す. [in,out] YYpd() 配列 YYpd(LYYpd - 1) (LYYpd >= N(N+1)/2) [in] IRev = 2の場合, 再呼び出し時にX()におけるヘッセ行列の下三角部分を圧縮形式で与えること [out] IRev = 3の場合, 中間結果の出力のため現在のX()におけるヘッセ行列の下三角部分を返す. [in,out] IRev リバースコミュニケーションの制御変数. [in] 最初の呼び出し時に 0 に設定しておくこと. 2回目以降の呼び出し時には値を変更してはならない. [out] 0 以外の時には下記処理を行ってから再び本ルーチンを呼び出すこと. = 0: 処理終了. 正常終了かどうかはInfoをチェックすること. = 1: X()における関数値を求め YYに設定すること. YY以外の変数を変更してはならない. = 2: X()における微分係数を求め YYp()に, また, ヘッセ行列の下三角部分をYYpd()に設定すること. YYp()およびYYpd()以外の変数を変更してはならない. = 3: Iout = 1であれば反復ごとにIRev = 3で戻る. 途中結果(X(), NFcall, NGcall, Niter, Fvalなど)を出力する. どの変数も変更してはならない. [in] Iout (省略可) 計算の途中結果を出力するかどうか指定する. (省略時 = 0) = 0: 出力しない (IRev = 3では戻らない). = 1: 反復ごとにIRev = 3で戻り, 途中結果を出力する. (上記以外の値であればIout = 0とみなす) [out] Info2 (省略可) Info = 0 のときのサブコード. = 1: X値の収束条件を満たした. = 2: 関数値相対収束条件を満たした. = 3: X値および関数値相対収束条件の両方を満たした. = 4: 関数値絶対収束条件を満たした. [out] NFcall (省略可) 関数評価回数(IRev = 1で戻った回数). [out] NGcall (省略可) 微分係数評価回数(IRev = 2で戻った回数). [out] Niter (省略可) 反復回数. [out] Fval (省略可) 求められた解ベクトルX()における関数値. [out] NFGcal (省略可) MnhにおいてサブルーチンFとGに渡されるNf(呼び出しカウンタ). [in] Rtol (省略可) 関数値の目標相対精度. (Eps <= Rtol <= 0.1) (省略時 = 1e-10) (Epsはマシンイプシロン) (Rtol < Eps または Rtol > 0.1 であれば省略時の既定値とみなす) [in] Atol (省略可) 関数値の目標絶対精度. (省略時 = 1e-20) (Atol < 0 であれば省略時の既定値とみなす) [in] MaxFcall (省略可) 関数Fの呼び出し回数の最大値. (省略時 = 200) (MaxFcall <= 0 であれば省略時の既定値とみなす) [in] MaxIter (省略可) 反復回数の最大値. (省略時 = 150) (MaxIter <= 0 であれば省略時の既定値とみなす) [in] Dtype (省略可) 自動スケーリングの設定. (Dtype = 0, 1 または 2) (省略時 = 0) = 0: 自動スケーリングを行わない. (スケール係数 = 1) = 1: 反復ごとに自動スケーリングを行う. = 2: 1回目の反復のみ自動スケーリングを行い, その後スケール係数を変更しない. (上記以外の値であれば省略時の既定値とみなす) [in] Dfac (省略可) 自動スケーリングの係数. (0 <= Dfac <= 1) (省略時 = 0.6) 自動スケーリングではD(i)をスケール係数として, すべてのiについてD(i)*X(i)が同程度の大きさになるように反復ごとにD(i)を調整する. まず, D1(i) = max(Sqr(Abs(H(i,i))), Dfac*D(i))とする(H(i,i)はヘッセ行列の対角要素). 次に, 以下のようにD(i)を調整する.   D1(i) >= Dtolの場合: D(i) = D1(i)   D1(i) < Dtolの場合: D(i) = max(D0, Dtol) (Dfac < 0 または Dfac > 1 であれば省略時の既定値とみなす) [in] Dtol (省略可) 自動スケーリングのしきい値. (Dtol > 0) (省略時 = 1.0e-6) (Dtol <= 0 であれば省略時の既定値とみなす) [in] D0 (省略可) 自動スケーリングの初期値. (D0 > 0) (省略時 = 1) (D0 <= 0 であれば省略時の既定値とみなす) [in] Tuner1 (省略可) 誤収束の判定パラメータ. (0 <= Tuner1 <= 0.5) (省略時 = 0.1) (Tuner1 < 0 または Tuner1 > 0.5 であれば省略時の既定値とみなす) [in] Xctol (省略可) X値の収束判定しきい値. (0 <= Xctol <= 1) (省略時 = Eps^(1/2)) (Xctol < 0 または Xctol > 1 であれば省略時の既定値とみなす) [in] Xftol (省略可) 誤収束の判定しきい値. (0 <= Xftol <= 1) (省略時 = 100*Eps) (Xftol < 0 または Xftol > 1 であれば省略時の既定値とみなす) [in] Lmax0 (省略可) スケーリングされた一番始めのステップ長の最大幅. (Lmax0 > 0) (省略時 = 1) (Lmax0 <= 0 であれば省略時の既定値とみなす) [in] Lmaxs (省略可) [in] Sctol (省略可) LmaxsおよびSctolは特異収束の判定パラメータである. (Lmaxs > 0) (0 <= Sctol <= 1) (省略時: Lmaxs = 1, Sctol = 1e-10) ステップ長 Lmaxs において Sctol*abs(f)よりも関数値が減少しなくなった (fは現在の反復の出発点における関数値). (Lmaxs <= 0 であれば省略時の既定値とみなす) (Sctol < 0 であれば省略時の既定値とみなす) [in] Bias (省略可) ドッグレッグ信頼領域法で使用するパラメータ. (0 <= Bias <= 1) (省略時 = 0.8) (Bias < 0 または Bias > 1 であれば省略時の既定値とみなす) 出典 netlib/port 使用例 次の関数の最小点を求める(ローゼンブロック関数). f(x1, x2) = 100(x2 - x1^2)^2 + (1 - x1)^2 初期値は, (x1, x2) = (-1.2, 1) とする. また, 解ベクトル X の境界制約条件を -2 <= x1 <= 2, 0 <= x2 <= 2 とする. Sub Ex_Mnhb_r() Const N As Long = 2 Dim X(N - 1) As Double, B(1, N - 1) As Double, Info As Long Dim YY As Double, YYp(N - 1) As Double, YYpd(N * (N + 1) / 2) As Double Dim IRev As Long X(0) = -1.2: X(1) = 1 B(0, 0) = -2: B(1, 0) = 2 B(0, 1) = 0: B(1, 1) = 2 IRev = 0 Do Call Mnhb_r(N, X(), B(), Info, YY, YYp(), YYpd(), IRev) If IRev = 1 Then YY = 100 * (X(1) - X(0) ^ 2) ^ 2 + (1 - X(0)) ^ 2 ElseIf IRev = 2 Then YYp(0) = -400 * X(0) * (X(1) - X(0) ^ 2) + 2 * X(0) - 2 YYp(1) = 200 * X(1) - 200 * X(0) ^ 2 YYpd(0) = 1200 * X(0) ^ 2 - 400 * X(1) + 2 YYpd(1) = -400 * X(0) YYpd(2) = 200 End If Loop While IRev <> 0 Debug.Print "X1, X2 =", X(0), X(1) Debug.Print "Info =", Info End Sub 実行結果 X1, X2 = 0.999999999999999 0.999999999999998 Info = 0