|
||||||||||
PREV CLASS NEXT CLASS | FRAMES NO FRAMES All Classes | |||||||||
SUMMARY: NESTED | FIELD | CONSTR | METHOD | DETAIL: FIELD | CONSTR | METHOD |
java.lang.Object | +--optimization.Minpack_f77
This class contains Java translations of the MINPACK nonlinear least squares routines. As of November 2000, it does not yet contain the MINPACK solvers of systems of nonlinear equations. They should be added in the Spring of 2001.
The original FORTRAN MINPACK package was produced by Burton S. Garbow, Kenneth E. Hillstrom, and Jorge J. More as part of the Argonne National Laboratory MINPACK project, March 1980.
IMPORTANT: The "_f77" suffixes indicate that these routines use FORTRAN style indexing. For example, you will see
for (i = 1; i <= n; i++)rather than
for (i = 0; i < n; i++)To use the "_f77" routines you will have to declare your vectors and matrices to be one element larger (e.g., v[101] rather than v[100], and a[101][101] rather than a[100][100]), and you will have to fill elements 1 through n rather than elements 0 through n - 1. Versions of these programs that use C/Java style indexing might eventually be available. They would end with the suffix "_j".
This class was translated by a statistician from FORTRAN versions of lmder and lmdif. It is NOT an official translation. It wastes memory by failing to use the first elements of vectors. When public domain Java optimization routines become available from the people who produced MINPACK, then THE CODE PRODUCED BY THE NUMERICAL ANALYSTS SHOULD BE USED.
Meanwhile, if you have suggestions for improving this code, please contact Steve Verrill at steve@ws13.fpl.fs.fed.us.
Constructor Summary | |
Minpack_f77()
|
Method Summary | |
static double |
enorm_f77(int n,
double[] x)
The enorm_f77 method calculates the Euclidean norm of a vector. |
static void |
fdjac2_f77(Lmdif_fcn nlls,
int m,
int n,
double[] x,
double[] fvec,
double[][] fjac,
int[] iflag,
double epsfcn,
double[] wa)
The fdjac2 method computes a forward-difference approximation to the m by n Jacobian matrix associated with a specified problem of m functions in n variables. |
static void |
lmder_f77(Lmder_fcn nlls,
int m,
int n,
double[] x,
double[] fvec,
double[][] fjac,
double ftol,
double xtol,
double gtol,
int maxfev,
double[] diag,
int mode,
double factor,
int nprint,
int[] info,
int[] nfev,
int[] njev,
int[] ipvt,
double[] qtf)
The lmder_f77 method minimizes the sum of the squares of m nonlinear functions in n variables by a modification of the Levenberg-Marquardt algorithm. |
static void |
lmder1_f77(Lmder_fcn nlls,
int m,
int n,
double[] x,
double[] fvec,
double[][] fjac,
double tol,
int[] info,
int[] ipvt)
The lmder1_f77 method minimizes the sum of the squares of m nonlinear functions in n variables by a modification of the Levenberg-Marquardt algorithm. |
static void |
lmdif_f77(Lmdif_fcn nlls,
int m,
int n,
double[] x,
double[] fvec,
double ftol,
double xtol,
double gtol,
int maxfev,
double epsfcn,
double[] diag,
int mode,
double factor,
int nprint,
int[] info,
int[] nfev,
double[][] fjac,
int[] ipvt,
double[] qtf)
The lmdif_f77 method minimizes the sum of the squares of m nonlinear functions in n variables by a modification of the Levenberg-Marquardt algorithm. |
static void |
lmdif1_f77(Lmdif_fcn nlls,
int m,
int n,
double[] x,
double[] fvec,
double tol,
int[] info)
The lmdif1_f77 method minimizes the sum of the squares of m nonlinear functions in n variables by a modification of the Levenberg-Marquardt algorithm. |
static void |
lmpar_f77(int n,
double[][] r,
int[] ipvt,
double[] diag,
double[] qtb,
double delta,
double[] par,
double[] x,
double[] sdiag,
double[] wa1,
double[] wa2)
Given an m by n matrix A, an n by n nonsingular diagonal matrix D, an m-vector b, and a positive number delta, the problem is to determine a value for the parameter par such that if x solves the system |
static void |
qrfac_f77(int m,
int n,
double[][] a,
boolean pivot,
int[] ipvt,
double[] rdiag,
double[] acnorm,
double[] wa)
The qrfac_f77 method uses Householder transformations with column pivoting (optional) to compute a QR factorization of the m by n matrix A. |
static void |
qrsolv_f77(int n,
double[][] r,
int[] ipvt,
double[] diag,
double[] qtb,
double[] x,
double[] sdiag,
double[] wa)
Given an m by n matrix A, an n by n diagonal matrix D, and an m-vector b, the problem is to determine an x which solves the system |
Methods inherited from class java.lang.Object |
clone, equals, finalize, getClass, hashCode, notify, notifyAll, toString, wait, wait, wait |
Constructor Detail |
public Minpack_f77()
Method Detail |
public static void lmder1_f77(Lmder_fcn nlls, int m, int n, double[] x, double[] fvec, double[][] fjac, double tol, int[] info, int[] ipvt)
The lmder1_f77 method minimizes the sum of the squares of m nonlinear functions in n variables by a modification of the Levenberg-Marquardt algorithm. This is done by using the more general least-squares solver lmder_f77. The user must provide a method which calculates the functions and the Jacobian.
Translated by Steve Verrill on November 17, 2000 from the FORTRAN MINPACK source produced by Garbow, Hillstrom, and More.
nlls
- A class that implements the Lmder_fcn interface
(see the definition in Lmder_fcn.java). See
LmderTest_f77.java for an example of such a class.
The class must define a method, fcn, that must
have the form
public static void fcn(int m, int n, double x[],
double fvec[], double fjac[][], int iflag[])
If iflag[1] equals 1, fcn calculates the values of
the m functions [residuals] at x and returns this
vector in fvec. If iflag[1] equals 2, fcn calculates
the Jacobian at x and returns this matrix in fjac
(and does not alter fvec).
The value of iflag[1] should not be changed by fcn unless
the user wants to terminate execution of lmder_f77.
In this case set iflag[1] to a negative integer.m
- A positive integer set to the number of functions
[number of observations]n
- A positive integer set to the number of variables
[number of parameters]. n must not exceed m.x
- On input, it contains the initial estimate of
the solution vector [the least squares parameters].
On output it contains the final estimate of the
solution vector.fvec
- An output vector that contains the m functions
[residuals] evaluated at x.fjac
- An output m by n array. The upper n by n submatrix
of fjac contains an upper triangular matrix R with
diagonal elements of nonincreasing magnitude such that
t t t P (jac jac)P = R R,where P is a permutation matrix and jac is the final calculated Jacobian. Column j of P is column ipvt[j] of the identity matrix. The lower trapezoidal part of fjac contains information generated during the computation of R.
tol
- tol is a nonnegative input variable. Termination occurs
when the algorithm estimates either that the relative
error in the sum of squares is at most tol or that
the relative error between x and the solution is at
most tol.info
- An integer output variable. If the user has
terminated execution, info is set to the (negative)
value of iflag[1]. See description of fcn. Otherwise,
info is set as follows.
info = 0 improper input parameters.
info = 1 algorithm estimates that the relative error
in the sum of squares is at most tol.
info = 2 algorithm estimates that the relative error
between x and the solution is at most tol.
info = 3 conditions for info = 1 and info = 2 both hold.
info = 4 fvec is orthogonal to the columns of the
Jacobian to machine precision.
info = 5 number of calls to fcn with iflag[1] = 1 has
reached 100*(n+1).
info = 6 tol is too small. No further reduction in
the sum of squares is possible.
info = 7 tol is too small. No further improvement in
the approximate solution x is possible.ipvt
- An integer output array of length n. ipvt
defines a permutation matrix P such that jac*P = QR,
where jac is the final calculated Jacobian, Q is
orthogonal (not stored), and R is upper triangular
with diagonal elements of nonincreasing magnitude.
Column j of P is column ipvt[j] of the identity matrix.public static void lmder_f77(Lmder_fcn nlls, int m, int n, double[] x, double[] fvec, double[][] fjac, double ftol, double xtol, double gtol, int maxfev, double[] diag, int mode, double factor, int nprint, int[] info, int[] nfev, int[] njev, int[] ipvt, double[] qtf)
The lmder_f77 method minimizes the sum of the squares of m nonlinear functions in n variables by a modification of the Levenberg-Marquardt algorithm. The user must provide a method which calculates the functions and the Jacobian.
Translated by Steve Verrill on November 3, 2000 from the FORTRAN MINPACK source produced by Garbow, Hillstrom, and More.
nlls
- A class that implements the Lmder_fcn interface
(see the definition in Lmder_fcn.java). See
LmderTest_f77.java for an example of such a class.
The class must define a method, fcn, that must
have the form
public static void fcn(int m, int n, double x[],
double fvec[], double fjac[][], int iflag[])
If iflag[1] equals 1, fcn calculates the values of
the m functions [residuals] at x and returns this
vector in fvec. If iflag[1] equals 2, fcn calculates
the Jacobian at x and returns this matrix in fjac
(and does not alter fvec).
The value of iflag[1] should not be changed by fcn unless
the user wants to terminate execution of lmder_f77.
In this case set iflag[1] to a negative integer.m
- A positive integer set to the number of functions
[number of observations]n
- A positive integer set to the number of variables
[number of parameters]. n must not exceed m.x
- On input, it contains the initial estimate of
the solution vector [the least squares parameters].
On output it contains the final estimate of the
solution vector.fvec
- An output vector that contains the m functions
[residuals] evaluated at x.fjac
- An output m by n array. The upper n by n submatrix
of fjac contains an upper triangular matrix R with
diagonal elements of nonincreasing magnitude such that
t t t P (jac jac)P = R R,where P is a permutation matrix and jac is the final calculated Jacobian. Column j of P is column ipvt[j] of the identity matrix. The lower trapezoidal part of fjac contains information generated during the computation of R.
ftol
- A nonnegative input variable. Termination
occurs when both the actual and predicted relative
reductions in the sum of squares are at most ftol.
Therefore, ftol measures the relative error desired
in the sum of squares.xtol
- A nonnegative input variable. Termination
occurs when the relative error between two consecutive
iterates is at most xtol. Therefore, xtol measures the
relative error desired in the approximate solution.gtol
- A nonnegative input variable. Termination
occurs when the cosine of the angle between fvec and
any column of the Jacobian is at most gtol in absolute
value. Therefore, gtol measures the orthogonality
desired between the function vector and the columns
of the Jacobian.maxfev
- A positive integer input variable. Termination
occurs when the number of calls to fcn with iflag[1] = 1
has reached maxfev.diag
- An vector of length n. If mode = 1 (see
below), diag is internally set. If mode = 2, diag
must contain positive entries that serve as
multiplicative scale factors for the variables.mode
- If mode = 1, the
variables will be scaled internally. If mode = 2,
the scaling is specified by the input diag. Other
values of mode are equivalent to mode = 1.factor
- A positive input variable used in determining the
initial step bound. This bound is set to the product of
factor and the euclidean norm of diag*x if nonzero, or else
to factor itself. In most cases factor should lie in the
interval (.1,100). 100 is a generally recommended value.nprint
- An integer input variable that enables controlled
printing of iterates if it is positive. In this case,
fcn is called with iflag[1] = 0 at the beginning of the first
iteration and every nprint iterations thereafter and
immediately prior to return, with x, fvec, and fjac
available for printing. fvec and fjac should not be
altered. If nprint is not positive, no special calls
of fcn with iflag[1] = 0 are made.info
- An integer output variable. If the user has
terminated execution, info is set to the (negative)
value of iflag[1]. See description of fcn. Otherwise,
info is set as follows.
info = 0 improper input parameters.
info = 1 both actual and predicted relative reductions
in the sum of squares are at most ftol.
info = 2 relative error between two consecutive iterates
is at most xtol.
info = 3 conditions for info = 1 and info = 2 both hold.
info = 4 the cosine of the angle between fvec and any
column of the Jacobian is at most gtol in
absolute value.
info = 5 number of calls to fcn with iflag[1] = 1 has
reached maxfev.
info = 6 ftol is too small. no further reduction in
the sum of squares is possible.
info = 7 xtol is too small. no further improvement in
the approximate solution x is possible.
info = 8 gtol is too small. fvec is orthogonal to the
columns of the Jacobian to machine precision.nfev
- An integer output variable set to the number of
calls to fcn with iflag[1] = 1.njev
- An integer output variable set to the number of
calls to fcn with iflag[1] = 2.ipvt
- An integer output array of length n. ipvt
defines a permutation matrix P such that jac*P = QR,
where jac is the final calculated Jacobian, Q is
orthogonal (not stored), and R is upper triangular
with diagonal elements of nonincreasing magnitude.
column j of P is column ipvt[j] of the identity matrix.qtf
- An output array of length n which contains
the first n elements of the vector (Q transpose)fvec.public static double enorm_f77(int n, double[] x)
The enorm_f77 method calculates the Euclidean norm of a vector.
Translated by Steve Verrill on November 14, 2000 from the FORTRAN MINPACK source produced by Garbow, Hillstrom, and More.
n
- The length of the vector, x.x
- The vector whose Euclidean norm is to be calculated.public static void qrfac_f77(int m, int n, double[][] a, boolean pivot, int[] ipvt, double[] rdiag, double[] acnorm, double[] wa)
The qrfac_f77 method uses Householder transformations with column pivoting (optional) to compute a QR factorization of the m by n matrix A. That is, qrfac_f77 determines an orthogonal matrix Q, a permutation matrix P, and an upper trapezoidal matrix R with diagonal elements of nonincreasing magnitude, such that AP = QR.
Translated by Steve Verrill on November 17, 2000 from the FORTRAN MINPACK source produced by Garbow, Hillstrom, and More.
m
- The number of rows of A.n
- The number of columns of A.a
- A is an m by n array. On input A contains the matrix for
which the QR factorization is to be computed. On output
the strict upper trapezoidal part of A contains the strict
upper trapezoidal part of R, and the lower trapezoidal
part of A contains a factored form of Q.pivot
- pivot is a logical input variable. If pivot is set true,
then column pivoting is enforced. If pivot is set false,
then no column pivoting is done.ipvt
- ipvt is an integer output array. ipvt
defines the permutation matrix P such that A*P = Q*R.
Column j of P is column ipvt[j] of the identity matrix.
If pivot is false, ipvt is not referenced.rdiag
- rdiag is an output array of length n which contains the
diagonal elements of R.acnorm
- acnorm is an output array of length n which contains the
norms of the corresponding columns of the input matrix A.wa
- wa is a work array of length n.public static void qrsolv_f77(int n, double[][] r, int[] ipvt, double[] diag, double[] qtb, double[] x, double[] sdiag, double[] wa)
Given an m by n matrix A, an n by n diagonal matrix D, and an m-vector b, the problem is to determine an x which solves the system
Ax = b , Dx = 0 ,in the least squares sense.
This method completes the solution of the problem if it is provided with the necessary information from the QR factorization, with column pivoting, of A. That is, if AP = QR, where P is a permutation matrix, Q has orthogonal columns, and R is an upper triangular matrix with diagonal elements of nonincreasing magnitude, then qrsolv_f77 expects the full upper triangle of R, the permutation matrix P, and the first n components of (Q transpose)b. The system
Ax = b, Dx = 0, is then equivalent to t t Rz = Q b, P DPz = 0 ,where x = Pz. If this system does not have full rank, then a least squares solution is obtained. On output qrsolv_f77 also provides an upper triangular matrix S such that
t t t P (A A + DD)P = S S .S is computed within qrsolv_f77 and may be of separate interest.
Translated by Steve Verrill on November 17, 2000 from the FORTRAN MINPACK source produced by Garbow, Hillstrom, and More.
n
- The order of r.r
- r is an n by n array. On input the full upper triangle
must contain the full upper triangle of the matrix R.
On output the full upper triangle is unaltered, and the
strict lower triangle contains the strict upper triangle
(transposed) of the upper triangular matrix S.ipvt
- ipvt is an integer input array of length n which defines the
permutation matrix P such that AP = QR. Column j of P
is column ipvt[j] of the identity matrix.diag
- diag is an input array of length n which must contain the
diagonal elements of the matrix D.qtb
- qtb is an input array of length n which must contain the first
n elements of the vector (Q transpose)b.x
- x is an output array of length n which contains the least
squares solution of the system Ax = b, Dx = 0.sdiag
- sdiag is an output array of length n which contains the
diagonal elements of the upper triangular matrix S.wa
- wa is a work array of length n.public static void lmpar_f77(int n, double[][] r, int[] ipvt, double[] diag, double[] qtb, double delta, double[] par, double[] x, double[] sdiag, double[] wa1, double[] wa2)
Given an m by n matrix A, an n by n nonsingular diagonal matrix D, an m-vector b, and a positive number delta, the problem is to determine a value for the parameter par such that if x solves the system
A*x = b , sqrt(par)*D*x = 0in the least squares sense, and dxnorm is the Euclidean norm of D*x, then either par is zero and
(dxnorm-delta) <= 0.1*delta ,or par is positive and
abs(dxnorm-delta) <= 0.1*delta .This method (lmpar_f77) completes the solution of the problem if it is provided with the necessary information from the QR factorization, with column pivoting, of A. That is, if AP = QR, where P is a permutation matrix, Q has orthogonal columns, and R is an upper triangular matrix with diagonal elements of nonincreasing magnitude, then lmpar_f77 expects the full upper triangle of R, the permutation matrix P, and the first n components of (Q transpose)b. On output lmpar_f77 also provides an upper triangular matrix S such that
t t t P (A A + par*DD)P = S S .S is employed within lmpar_f77 and may be of separate interest.
Only a few iterations are generally needed for convergence of the algorithm. If, however, the limit of 10 iterations is reached, then the output par will contain the best value obtained so far.
Translated by Steve Verrill on November 17, 2000 from the FORTRAN MINPACK source produced by Garbow, Hillstrom, and More.
n
- The order of r.r
- r is an n by n array. On input the full upper triangle
must contain the full upper triangle of the matrix R.
On output the full upper triangle is unaltered, and the
strict lower triangle contains the strict upper triangle
(transposed) of the upper triangular matrix S.ipvt
- ipvt is an integer input array of length n which defines the
permutation matrix P such that AP = QR. Column j of P
is column ipvt[j] of the identity matrix.diag
- diag is an input array of length n which must contain the
diagonal elements of the matrix D.qtb
- qtb is an input array of length n which must contain the first
n elements of the vector (Q transpose)b.delta
- delta is a positive input variable which specifies an upper
bound on the Euclidean norm of Dx.par
- par is a nonnegative variable. On input par contains an
initial estimate of the Levenberg-Marquardt parameter.
On output par contains the final estimate.x
- x is an output array of length n which contains the least
squares solution of the system Ax = b, sqrt(par)*Dx = 0,
for the output par.sdiag
- sdiag is an output array of length n which contains the
diagonal elements of the upper triangular matrix S.wa1
- wa1 is a work array of length n.wa2
- wa2 is a work array of length n.public static void lmdif1_f77(Lmdif_fcn nlls, int m, int n, double[] x, double[] fvec, double tol, int[] info)
The lmdif1_f77 method minimizes the sum of the squares of m nonlinear functions in n variables by a modification of the Levenberg-Marquardt algorithm. This is done by using the more general least-squares solver lmdif. The user must provide a method that calculates the functions. The Jacobian is then calculated by a forward-difference approximation.
Translated by Steve Verrill on November 24, 2000 from the FORTRAN MINPACK source produced by Garbow, Hillstrom, and More.
nlls
- A class that implements the Lmdif_fcn interface
(see the definition in Lmdif_fcn.java). See
LmdifTest_f77.java for an example of such a class.
The class must define a method, fcn, that must
have the form
public static void fcn(int m, int n, double x[],
double fvec[], int iflag[])
The value of iflag[1] should not be changed by fcn unless
the user wants to terminate execution of lmdif_f77.
In this case set iflag[1] to a negative integer.m
- A positive integer set to the number of functions
[number of observations]n
- A positive integer set to the number of variables
[number of parameters]. n must not exceed m.x
- On input, it contains the initial estimate of
the solution vector [the least squares parameters].
On output it contains the final estimate of the
solution vector.fvec
- An output vector that contains the m functions
[residuals] evaluated at x.tol
- tol is a nonnegative input variable. Termination occurs
when the algorithm estimates either that the relative
error in the sum of squares is at most tol or that
the relative error between x and the solution is at
most tol.info
- An integer output variable. If the user has
terminated execution, info is set to the (negative)
value of iflag[1]. See description of fcn. Otherwise,
info is set as follows.
info = 0 improper input parameters.
info = 1 algorithm estimates that the relative error
in the sum of squares is at most tol.
info = 2 algorithm estimates that the relative error
between x and the solution is at most tol.
info = 3 conditions for info = 1 and info = 2 both hold.
info = 4 fvec is orthogonal to the columns of the
Jacobian to machine precision.
info = 5 number of calls to fcn has
reached or exceeded 200*(n+1).
info = 6 tol is too small. No further reduction in
the sum of squares is possible.
info = 7 tol is too small. No further improvement in
the approximate solution x is possible.public static void lmdif_f77(Lmdif_fcn nlls, int m, int n, double[] x, double[] fvec, double ftol, double xtol, double gtol, int maxfev, double epsfcn, double[] diag, int mode, double factor, int nprint, int[] info, int[] nfev, double[][] fjac, int[] ipvt, double[] qtf)
The lmdif_f77 method minimizes the sum of the squares of m nonlinear functions in n variables by a modification of the Levenberg-Marquardt algorithm. The user must provide a method that calculates the functions. The Jacobian is then calculated by a forward-difference approximation.
Translated by Steve Verrill on November 20, 2000 from the FORTRAN MINPACK source produced by Garbow, Hillstrom, and More.
nlls
- A class that implements the Lmdif_fcn interface
(see the definition in Lmdif_fcn.java). See
LmdifTest_f77.java for an example of such a class.
The class must define a method, fcn, that must
have the form
public static void fcn(int m, int n, double x[],
double fvec[], int iflag[])
The value of iflag[1] should not be changed by fcn unless
the user wants to terminate execution of lmdif_f77.
In this case set iflag[1] to a negative integer.m
- A positive integer set to the number of functions
[number of observations]n
- A positive integer set to the number of variables
[number of parameters]. n must not exceed m.x
- On input, it contains the initial estimate of
the solution vector [the least squares parameters].
On output it contains the final estimate of the
solution vector.fvec
- An output vector that contains the m functions
[residuals] evaluated at x.ftol
- A nonnegative input variable. Termination
occurs when both the actual and predicted relative
reductions in the sum of squares are at most ftol.
Therefore, ftol measures the relative error desired
in the sum of squares.xtol
- A nonnegative input variable. Termination
occurs when the relative error between two consecutive
iterates is at most xtol. Therefore, xtol measures the
relative error desired in the approximate solution.gtol
- A nonnegative input variable. Termination
occurs when the cosine of the angle between fvec and
any column of the Jacobian is at most gtol in absolute
value. Therefore, gtol measures the orthogonality
desired between the function vector and the columns
of the Jacobian.maxfev
- A positive integer input variable. Termination
occurs when the number of calls to fcn is at least
maxfev by the end of an iteration.epsfcn
- An input variable used in determining a suitable
step length for the forward-difference approximation. This
approximation assumes that the relative errors in the
functions are of the order of epsfcn. If epsfcn is less
than the machine precision, it is assumed that the relative
errors in the functions are of the order of the machine
precision.diag
- An vector of length n. If mode = 1 (see
below), diag is internally set. If mode = 2, diag
must contain positive entries that serve as
multiplicative scale factors for the variables.mode
- If mode = 1, the
variables will be scaled internally. If mode = 2,
the scaling is specified by the input diag. Other
values of mode are equivalent to mode = 1.factor
- A positive input variable used in determining the
initial step bound. This bound is set to the product of
factor and the euclidean norm of diag*x if nonzero, or else
to factor itself. In most cases factor should lie in the
interval (.1,100). 100 is a generally recommended value.nprint
- An integer input variable that enables controlled
printing of iterates if it is positive. In this case,
fcn is called with iflag[1] = 0 at the beginning of the first
iteration and every nprint iterations thereafter and
immediately prior to return, with x and fvec
available for printing. If nprint is not positive,
no special calls of fcn with iflag[1] = 0 are made.info
- An integer output variable. If the user has
terminated execution, info is set to the (negative)
value of iflag[1]. See description of fcn. Otherwise,
info is set as follows.
info = 0 improper input parameters.
info = 1 both actual and predicted relative reductions
in the sum of squares are at most ftol.
info = 2 relative error between two consecutive iterates
is at most xtol.
info = 3 conditions for info = 1 and info = 2 both hold.
info = 4 the cosine of the angle between fvec and any
column of the Jacobian is at most gtol in
absolute value.
info = 5 number of calls to fcn with iflag[1] = 1 has
reached maxfev.
info = 6 ftol is too small. no further reduction in
the sum of squares is possible.
info = 7 xtol is too small. no further improvement in
the approximate solution x is possible.
info = 8 gtol is too small. fvec is orthogonal to the
columns of the Jacobian to machine precision.nfev
- An integer output variable set to the number of
calls to fcn.fjac
- An output m by n array. The upper n by n submatrix
of fjac contains an upper triangular matrix R with
diagonal elements of nonincreasing magnitude such that
t t t
P (jac *jac)P = R R,
where P is a permutation matrix and jac is the final
calculated Jacobian. Column j of P is column ipvt[j]
(see below) of the identity matrix. The lower trapezoidal
part of fjac contains information generated during
the computation of R.ipvt
- An integer output array of length n. ipvt
defines a permutation matrix P such that jac*P = QR,
where jac is the final calculated Jacobian, Q is
orthogonal (not stored), and R is upper triangular
with diagonal elements of nonincreasing magnitude.
column j of P is column ipvt[j] of the identity matrix.qtf
- An output array of length n which contains
the first n elements of the vector (Q transpose)fvec.public static void fdjac2_f77(Lmdif_fcn nlls, int m, int n, double[] x, double[] fvec, double[][] fjac, int[] iflag, double epsfcn, double[] wa)
The fdjac2 method computes a forward-difference approximation to the m by n Jacobian matrix associated with a specified problem of m functions in n variables.
Translated by Steve Verrill on November 24, 2000 from the FORTRAN MINPACK source produced by Garbow, Hillstrom, and More.
nlls
- A class that implements the Lmdif_fcn interface
(see the definition in Lmdif_fcn.java). See
LmdifTest_f77.java for an example of such a class.
The class must define a method, fcn, that must
have the form
public static void fcn(int m, int n, double x[],
double fvec[], int iflag[])
The value of iflag[1] should not be changed by fcn unless
the user wants to terminate execution of fdjac2_f77.
In this case iflag[1] should be set to a negative integer.m
- A positive integer set to the number of functions
[number of observations]n
- A positive integer set to the number of variables
[number of parameters]. n must not exceed m.x
- An input array.fvec
- An input array that contains the functions
evaluated at x.fjac
- An output m by n array that contains the
approximation to the Jacobian matrix evaluated at x.iflag
- An integer variable that can be used to terminate
the execution of fdjac2. See the description of nlls.epsfcn
- An input variable used in determining a suitable
step length for the forward-difference approximation. This
approximation assumes that the relative errors in the
functions are of the order of epsfcn. If epsfcn is less
than the machine precision, it is assumed that the relative
errors in the functions are of the order of the machine
precision.wa
- A work array.
|
||||||||||
PREV CLASS NEXT CLASS | FRAMES NO FRAMES All Classes | |||||||||
SUMMARY: NESTED | FIELD | CONSTR | METHOD | DETAIL: FIELD | CONSTR | METHOD |