1
subroutine lmstr1(fcn,m,n,x,fvec,fjac,ldfjac,tol,info,ipvt,wa,
3
integer m,n,ldfjac,info,lwa
6
double precision x(n),fvec(m),fjac(ldfjac,n),wa(lwa)
12
c the purpose of lmstr1 is to minimize the sum of the squares of
13
c m nonlinear functions in n variables by a modification of
14
c the levenberg-marquardt algorithm which uses minimal storage.
15
c this is done by using the more general least-squares solver
16
c lmstr. the user must provide a subroutine which calculates
17
c the functions and the rows of the jacobian.
19
c the subroutine statement is
21
c subroutine lmstr1(fcn,m,n,x,fvec,fjac,ldfjac,tol,info,
26
c fcn is the name of the user-supplied subroutine which
27
c calculates the functions and the rows of the jacobian.
28
c fcn must be declared in an external statement in the
29
c user calling program, and should be written as follows.
31
c subroutine fcn(m,n,x,fvec,fjrow,iflag)
33
c double precision x(n),fvec(m),fjrow(n)
35
c if iflag = 1 calculate the functions at x and
36
c return this vector in fvec.
37
c if iflag = i calculate the (i-1)-st row of the
38
c jacobian at x and return this vector in fjrow.
43
c the value of iflag should not be changed by fcn unless
44
c the user wants to terminate execution of lmstr1.
45
c in this case set iflag to a negative integer.
47
c m is a positive integer input variable set to the number
50
c n is a positive integer input variable set to the number
51
c of variables. n must not exceed m.
53
c x is an array of length n. on input x must contain
54
c an initial estimate of the solution vector. on output x
55
c contains the final estimate of the solution vector.
57
c fvec is an output array of length m which contains
58
c the functions evaluated at the output x.
60
c fjac is an output n by n array. the upper triangle of fjac
61
c contains an upper triangular matrix r such that
64
c p *(jac *jac)*p = r *r,
66
c where p is a permutation matrix and jac is the final
67
c calculated jacobian. column j of p is column ipvt(j)
68
c (see below) of the identity matrix. the lower triangular
69
c part of fjac contains information generated during
70
c the computation of r.
72
c ldfjac is a positive integer input variable not less than n
73
c which specifies the leading dimension of the array fjac.
75
c tol is a nonnegative input variable. termination occurs
76
c when the algorithm estimates either that the relative
77
c error in the sum of squares is at most tol or that
78
c the relative error between x and the solution is at
81
c info is an integer output variable. if the user has
82
c terminated execution, info is set to the (negative)
83
c value of iflag. see description of fcn. otherwise,
84
c info is set as follows.
86
c info = 0 improper input parameters.
88
c info = 1 algorithm estimates that the relative error
89
c in the sum of squares is at most tol.
91
c info = 2 algorithm estimates that the relative error
92
c between x and the solution is at most tol.
94
c info = 3 conditions for info = 1 and info = 2 both hold.
96
c info = 4 fvec is orthogonal to the columns of the
97
c jacobian to machine precision.
99
c info = 5 number of calls to fcn with iflag = 1 has
102
c info = 6 tol is too small. no further reduction in
103
c the sum of squares is possible.
105
c info = 7 tol is too small. no further improvement in
106
c the approximate solution x is possible.
108
c ipvt is an integer output array of length n. ipvt
109
c defines a permutation matrix p such that jac*p = q*r,
110
c where jac is the final calculated jacobian, q is
111
c orthogonal (not stored), and r is upper triangular.
112
c column j of p is column ipvt(j) of the identity matrix.
114
c wa is a work array of length lwa.
116
c lwa is a positive integer input variable not less than 5*n+m.
120
c user-supplied ...... fcn
122
c minpack-supplied ... lmstr
124
c argonne national laboratory. minpack project. march 1980.
125
c burton s. garbow, dudley v. goetschel, kenneth e. hillstrom,
129
integer maxfev,mode,nfev,njev,nprint
130
double precision factor,ftol,gtol,xtol,zero
131
data factor,zero /1.0d2,0.0d0/
134
c check the input parameters for errors.
136
if (n .le. 0 .or. m .lt. n .or. ldfjac .lt. n .or. tol .lt. zero
137
* .or. lwa .lt. 5*n + m) go to 10
147
call lmstr(fcn,m,n,x,fvec,fjac,ldfjac,ftol,xtol,gtol,maxfev,
148
* wa(1),mode,factor,nprint,info,nfev,njev,ipvt,wa(n+1),
149
* wa(2*n+1),wa(3*n+1),wa(4*n+1),wa(5*n+1))
150
if (info .eq. 8) info = 4
154
c last card of subroutine lmstr1.