|
|
◆ N2g()
| Sub N2g |
( |
M As |
Long, |
|
|
N As |
Long, |
|
|
X() As |
Double, |
|
|
F As |
LongPtr, |
|
|
Fj As |
LongPtr, |
|
|
Cov() As |
Double, |
|
|
Rd() As |
Double, |
|
|
Info As |
Long, |
|
|
Optional Itsum As |
LongPtr = NullPtr, |
|
|
Optional Info2 As |
Long, |
|
|
Optional NFcall As |
Long, |
|
|
Optional NFjcall As |
Long, |
|
|
Optional Niter As |
Long, |
|
|
Optional S As |
Double, |
|
|
Optional NFcov As |
Long, |
|
|
Optional NFjcov As |
Long, |
|
|
Optional Rtol As |
Double = -1, |
|
|
Optional Atol As |
Double = -1, |
|
|
Optional Rdreq As |
Long = -1, |
|
|
Optional Covreq As |
Long = -1, |
|
|
Optional MaxFcall As |
Long = -1, |
|
|
Optional MaxIter As |
Long = -1, |
|
|
Optional Dtype As |
Long = -1, |
|
|
Optional Dfac As |
Double = -1, |
|
|
Optional Dtol As |
Double = -1, |
|
|
Optional D0 As |
Double = -1, |
|
|
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 Delta0 As |
Double = -1, |
|
|
Optional Dltfdc As |
Double = -1 |
|
) |
| |
Nonlinear least squares approximation (adaptive algorithm)
- Purpose
- This routine minimizes the sum of the squares of M nonlinear functions in N variables by the adaptive algorithm which combines and augments a Gauss-Newton, Levenberg-Marquardt and other techniques for better convergence.
minimize the sum of fi(x1, x2, ..., xn)^2 (sum for i = 1 to M)
The user must provide a subroutine which calculates the function and Jacobian values.
- Parameters
-
| [in] | M | Number of data. (M > 0) |
| [in] | N | Number of parameters. (0 < N <= M) |
| [in,out] | X() | Array X(LX - 1) (LX >= N)
[in] X must contain an initial estimate of the solution vector.
[out] X contains the final estimate of the solution vector. |
| [in] | F | The user-supplied subroutine which calculates the function fi(x1, x2, ..., xn) defined as follows. Sub F(M As Long, N As Long, X() As Double, Nf As Long, Fvec() As Double)
Calculate the Fi value from given X() and return in Fvec(i-1) (i = 1 to M).
Nf is invocation counter. If given X() is out of bounds, Nf should be set to 0.
The other variables should not be altered.
For same X() value, same Nf will be passed to F and Fj.
End Sub
|
| [in] | Fj | The user-supplied subroutine which calculates the Jacobian of the function fi(x1, x2, ..., xn) defined as follows. Sub Fj(M As Long, N As Long, X() As Double, Nf As Long, Fjac() As Double)
Calculate the Jacobian (dFi/dXj) from given X() and return in Fjac(i-1,j-1) (i = 1 to M, j = 1 to N).
Nf is invocation counter. If given X() is out of bounds, Nf should be set to 0.
The other variables should not be altered.
For same X() value, same Nf will be passed to F and Fj.
End Sub
|
| [out] | Cov() | Array Cov(LCov - 1) (LCov >= N(N + 1)/2)
Covariance matrix. (Lower triangular part stored in column major order)
(Computed only when Rdreq = 1 or 3 and normally terminated) |
| [out] | Rd() | Array Rd(LRd - 1) (LRd >= M)
Regression diagnostic vector. (Rd(i) is an index of the change that would occur in the residual sum of squares if the i-th data is deleted)
(Computed only when Rdreq = 2 or 3 and normally terminated) |
| [out] | Info | = 0: Successful exit. (Sub-code is set to Info2)
= -1: The argument M had an illegal value. (M < N)
= -2: The argument N had an illegal value. (N < 1)
= -3: The argument X() is invalid.
= -6: The argument Cov() is invalid.
= -7: The argument Rd() is invalid.
= 7: Singular convergence. (Hessian near the current iterate appears to be singular)
= 8: False convergence. (Iterate appears to be converging to a noncritical point. Tolerances may be too small)
= 9: Function evaluation limit reached.
= 10: Iteration limit reached.
= 63: F(X) cannot be computed at the initial X.
= 65: The gradient could not be computed at X. |
| [in] | Itsum | (Optional)
The user supplied subroutine to print the intermediate results defined as follows. (default = NUllPtr)
If the address is supplied (if Itsum <> NullPtr), the subroutine is called after every iteration. Sub Itsum(N As Long, X() As Double, NIter As Long, Nf As Long, Nfj As Long, S As Double)
Output the following information in desired format.
N: Number of variables.
X(): Current approximation of the solution vector.
NIter: Iteration counter.
Nf: Number of function calls of F.
Nfj: Number of function calls of Fj.
S: Residual sum of squares at X().
End Sub
Argument values should not be altered. |
| [out] | Info2 | (Optional)
Sub-code for Info = 0.
= 1: X convergence.
= 2: Relative function convergence.
= 3: Both X and relative function convergence.
= 4: Absolute function convergence. |
| [out] | NFcall | (Optional)
Number of function evaluations of F including those used in computing the covariance matrix. |
| [out] | NFjcall | (Optional)
Number of function evaluations of Fj including those used in computing the covariance matrix. |
| [out] | Niter | (Optional)
Number of iterations. |
| [out] | S | (Optional)
Residual sum of squares at the obtained solution vector X(). |
| [out] | NFcov | (Optional)
Number of function evaluations of F to compute the covariance matrix. |
| [out] | NFjcov | (Optional)
Number of function evaluations of Fj to compute the covariance matrix. |
| [in] | Rtol | (Optional)
Relative function convergence tolerance. (Eps <= Rtol <= 0.1) (default = 1e-10)
(Eps shows the machine epsilon hereafter)
(If Rtol < Eps or Rtol > 0.1, the default value will be used) |
| [in] | Atol | (Optional)
Absolute function convergence tolerance. (default = 1e-20)
(If Atol < 0, the default value will be used) |
| [in] | Rdreq | (Optional)
Whether to compute a covariance matrix or regression diagnostic vector. (default = 3)
= 0: Compute neither.
= 1: Compute covariance matrix only.
= 2: Compute regression diagnostic vector only.
= 3: Compute both covariance matrix and regression diagnostic vector.
(For other values, the default value will be used) |
| [in] | Covreq | (Optional)
Which covariance matrix of the form is to be computed. (default = 1)
= 1, -1: sigma * H^(-1) * (J^T * J) * H^(-1)
= 2, -2: sigma * H^(-1)
= 3, -3: sigma * (J^T * J)
where sigma = S / max(1, M-N) (S is residual sum of squares), H is the Hessian matrix and J is the Jacobian matrix.
The sign shows which user's routines are to be used.
< 0: Compute using only F.
> 0: Compute using both F and Fj. (If Fj is not available, only F is to be used)
(For other values, the default value will be used) |
| [in] | MaxFcall | (Optional)
Maximum number of function evaluations of F. (default = 200)
(If MaxFcall <= 0, the default value will be used) |
| [in] | MaxIter | (Optional)
Maximum number of iterations. (default = 150)
(If MaxIter <= 0, the default value will be used) |
| [in] | Dtype | (Optional)
Choice of adaptive scaling. (Dtype = 0, 1 or 2) (default = 1)
= 0: Disable adaptive scaling. (scale factor = 1)
= 1: Enable adaptive scaling during all iterations.
= 2: Enable adaptive scaling during the first iteration and scale factor is left unchanged thereafter.
(For other values, the default value will be used) |
| [in] | Dfac | (Optional)
Factor for adaptive scaling (0 <= Dfac <= 1) (default = 0.6)
A scale factor D(i) is chosen by adaptive scaling so that D(i)*X(i) has about the same magnitude for all i.
Let D1(i) = max(||Ji||, Dfac*D(i)) where ||Ji|| is the 2-norm of the i-th column of Jacobian matrix, then D(i) is chosen as follows.
if D1(i) >= Dtol: D(i) = D1(i)
if D1(i) < Dtol: D(i) = D0
(If Dfac < 0 or Dfac > 1, the default value will be used) |
| [in] | Dtol | (Optional)
Tolerance for adaptive scaling. (Dtol > 0) (default = 1.0e-6)
(If Dtol <= 0, the default value will be used) |
| [in] | D0 | (Optional)
Initial value for adaptive scaling. (D0 > 0) (default = 1)
(If D0 <= 0, the default value will be used) |
| [in] | Tuner1 | (Optional)
Parameter to check for false convergence. (0 <= Tuner1 <= 0.5) (default = 0.1)
(If Tuner1 < 0 or Tuner1 > 0.5, the default value will be used) |
| [in] | Xctol | (Optional)
X convergence tolerance. (0 <= Xctol <= 1) (default = Eps^(1/2))
(If Xctol < 0 or Xctol > 1, the default value will be used) |
| [in] | Xftol | (Optional)
False convergence tolerance. (0 <= Xftol <= 1) (default = 100*Eps)
(If Xftol < 0 or Xftol > 1, the default value will be used) |
| [in] | Lmax0 | (Optional)
Maximum 2-norm allowed for scaled very first step. (Lmax0 > 0) (default = 1)
(If Lmax0 <= 0, the default value will be used) |
| [in] | Lmaxs | (Optional)
|
| [in] | Sctol | (Optional)
Lmaxs and Sctol are the singular convergence test parameters. (Lmaxs > 0, 0 <= Sctol <= 1) (default: Lmaxs = 1, Sctol = 1e-10)
If the function reduction predicted for a step of length bounded by Lmaxs is less than Sctol*abs(f), returns with Info = 7 (f is the function value at the start of the current iteration).
(If Lmaxs <= 0, the default value will be used)
(If Sctol < 0 or Sctol > 1, the default value will be used) |
| [in] | Delta0 | (Optional)
Step size used to compute a covariance matrix is chosen as follows (in the case that Covreq = 1 or 2). (Eps <= Delta0 <= 1) (default = Sqrt(Eps))
Delta0*max(|X(I)|, 1/D(I))*sign(X(I))
(If Delta0 < Eps or Delta0 > 1, the default value will be used) |
| [in] | Dltfdc | (Optional)
Step size used to compute a covariance matrix is chosen as follows (in the case that Covreq = -1 or -2). (Eps <= Dltfdc <= 1) (default = Eps^(1/3))
Dltfdc*max(|X(I)|, 1/D(I))
(If Dltfdc < Eps or Dltfdc > 1, the default value will be used) |
- Reference
- netlib/port
- Example Program
- Approximate the following data with model function f(x) = c1*(1 - exp(-c2*x)). Determine two parameters c1 and c2 by the nonlinear least squares method, and compute the covariance matrix.
f(x) x
10.07 77.6
29.61 239.9
50.76 434.8
81.78 760.0
The initial estimates of the solution are c1 = 500 and c2 = 0.0001. Sub FN2g(M As Long, N As Long, X() As Double, Nf As Long, Fvec() As Double)
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
For I = 0 To M - 1
Fvec(I) = Ydata(I) - X(0) * (1 - Exp(-Xdata(I) * X(1)))
Next
End Sub
Sub JN2g(M As Long, N As Long, X() As Double, Nf As Long, Fjac() As Double)
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
For I = 0 To M - 1
Fjac(I, 0) = Exp(-Xdata(I) * X(1)) - 1
Fjac(I, 1) = -Xdata(I) * X(0) * Exp(-X(1) * Xdata(I))
Next
End Sub
Sub Ex_N2g()
Const M = 4, N = 2
Dim X(N - 1) As Double, Cov(N * (N + 1) / 2 - 1) As Double, Rd(M - 1) As Double
Dim Info As Long
X(0) = 500: X(1) = 0.0001
Call N2g(M, N, X(), AddressOf FN2g, AddressOf JN2g, Cov(), Rd(), Info)
Debug.Print "C1, C2 =", X(0), X(1)
Debug.Print "Cov ="
Debug.Print Cov(0), Cov(1)
Debug.Print Cov(1), Cov(2)
Debug.Print "Info =", Info
End Sub
Function Rd(X As Double, Y As Double, Z As Double, Optional Info As Long) As Double Carlson form of elliptic integral RD(x, y, z)
Sub N2g(M As Long, N As Long, X() As Double, F As LongPtr, Fj As LongPtr, Cov() As Double, Rd() As Double, Info As Long, Optional Itsum As LongPtr=NullPtr, Optional Info2 As Long, Optional NFcall As Long, Optional NFjcall As Long, Optional Niter As Long, Optional S As Double, Optional NFcov As Long, Optional NFjcov As Long, Optional Rtol As Double=-1, Optional Atol As Double=-1, Optional Rdreq As Long=-1, Optional Covreq As Long=-1, Optional MaxFcall As Long=-1, Optional MaxIter As Long=-1, Optional Dtype As Long=-1, Optional Dfac As Double=-1, Optional Dtol As Double=-1, Optional D0 As Double=-1, 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 Delta0 As Double=-1, Optional Dltfdc As Double=-1) Nonlinear least squares approximation (adaptive algorithm)
- Example Results
C1, C2 = 241.084896112856 5.44942234058364E-04
Cov =
20.9246334498711 -5.6146048950131E-05
-5.6146048950131E-05 1.51119765045832E-10
Info = 0
|