1 #include "cminpack.h"
2 #include "cminpackP.h"
3
4 __cminpack_attr__
__cminpack_func__(lmder1)5 int __cminpack_func__(lmder1)(__cminpack_decl_fcnder_mn__ void *p, int m, int n, real *x,
6 real *fvec, real *fjac, int ldfjac, real tol,
7 int *ipvt, real *wa, int lwa)
8 {
9 /* Initialized data */
10
11 const real factor = 100.;
12
13 /* Local variables */
14 int mode, nfev, njev;
15 real ftol, gtol, xtol;
16 int maxfev, nprint;
17 int info;
18
19 /* ********** */
20
21 /* subroutine lmder1 */
22
23 /* the purpose of lmder1 is to minimize the sum of the squares of */
24 /* m nonlinear functions in n variables by a modification of the */
25 /* levenberg-marquardt algorithm. this is done by using the more */
26 /* general least-squares solver lmder. the user must provide a */
27 /* subroutine which calculates the functions and the jacobian. */
28
29 /* the subroutine statement is */
30
31 /* subroutine lmder1(fcn,m,n,x,fvec,fjac,ldfjac,tol,info, */
32 /* ipvt,wa,lwa) */
33
34 /* where */
35
36 /* fcn is the name of the user-supplied subroutine which */
37 /* calculates the functions and the jacobian. fcn must */
38 /* be declared in an external statement in the user */
39 /* calling program, and should be written as follows. */
40
41 /* subroutine fcn(m,n,x,fvec,fjac,ldfjac,iflag) */
42 /* integer m,n,ldfjac,iflag */
43 /* double precision x(n),fvec(m),fjac(ldfjac,n) */
44 /* ---------- */
45 /* if iflag = 1 calculate the functions at x and */
46 /* return this vector in fvec. do not alter fjac. */
47 /* if iflag = 2 calculate the jacobian at x and */
48 /* return this matrix in fjac. do not alter fvec. */
49 /* ---------- */
50 /* return */
51 /* end */
52
53 /* the value of iflag should not be changed by fcn unless */
54 /* the user wants to terminate execution of lmder1. */
55 /* in this case set iflag to a negative integer. */
56
57 /* m is a positive integer input variable set to the number */
58 /* of functions. */
59
60 /* n is a positive integer input variable set to the number */
61 /* of variables. n must not exceed m. */
62
63 /* x is an array of length n. on input x must contain */
64 /* an initial estimate of the solution vector. on output x */
65 /* contains the final estimate of the solution vector. */
66
67 /* fvec is an output array of length m which contains */
68 /* the functions evaluated at the output x. */
69
70 /* fjac is an output m by n array. the upper n by n submatrix */
71 /* of fjac contains an upper triangular matrix r with */
72 /* diagonal elements of nonincreasing magnitude such that */
73
74 /* t t t */
75 /* p *(jac *jac)*p = r *r, */
76
77 /* where p is a permutation matrix and jac is the final */
78 /* calculated jacobian. column j of p is column ipvt(j) */
79 /* (see below) of the identity matrix. the lower trapezoidal */
80 /* part of fjac contains information generated during */
81 /* the computation of r. */
82
83 /* ldfjac is a positive integer input variable not less than m */
84 /* which specifies the leading dimension of the array fjac. */
85
86 /* tol is a nonnegative input variable. termination occurs */
87 /* when the algorithm estimates either that the relative */
88 /* error in the sum of squares is at most tol or that */
89 /* the relative error between x and the solution is at */
90 /* most tol. */
91
92 /* info is an integer output variable. if the user has */
93 /* terminated execution, info is set to the (negative) */
94 /* value of iflag. see description of fcn. otherwise, */
95 /* info is set as follows. */
96
97 /* info = 0 improper input parameters. */
98
99 /* info = 1 algorithm estimates that the relative error */
100 /* in the sum of squares is at most tol. */
101
102 /* info = 2 algorithm estimates that the relative error */
103 /* between x and the solution is at most tol. */
104
105 /* info = 3 conditions for info = 1 and info = 2 both hold. */
106
107 /* info = 4 fvec is orthogonal to the columns of the */
108 /* jacobian to machine precision. */
109
110 /* info = 5 number of calls to fcn with iflag = 1 has */
111 /* reached 100*(n+1). */
112
113 /* info = 6 tol is too small. no further reduction in */
114 /* the sum of squares is possible. */
115
116 /* info = 7 tol is too small. no further improvement in */
117 /* the approximate solution x is possible. */
118
119 /* ipvt is an integer output array of length n. ipvt */
120 /* defines a permutation matrix p such that jac*p = q*r, */
121 /* where jac is the final calculated jacobian, q is */
122 /* orthogonal (not stored), and r is upper triangular */
123 /* with diagonal elements of nonincreasing magnitude. */
124 /* column j of p is column ipvt(j) of the identity matrix. */
125
126 /* wa is a work array of length lwa. */
127
128 /* lwa is a positive integer input variable not less than 5*n+m. */
129
130 /* subprograms called */
131
132 /* user-supplied ...... fcn */
133
134 /* minpack-supplied ... lmder */
135
136 /* argonne national laboratory. minpack project. march 1980. */
137 /* burton s. garbow, kenneth e. hillstrom, jorge j. more */
138
139 /* ********** */
140
141 /* check the input parameters for errors. */
142
143 if (n <= 0 || m < n || ldfjac < m || tol < 0. || lwa < n * 5 + m) {
144 return 0;
145 }
146
147 /* call lmder. */
148
149 maxfev = (n + 1) * 100;
150 ftol = tol;
151 xtol = tol;
152 gtol = 0.;
153 mode = 1;
154 nprint = 0;
155 info = __cminpack_func__(lmder)(__cminpack_param_fcnder_mn__ p, m, n, x, fvec, fjac, ldfjac,
156 ftol, xtol, gtol, maxfev, wa, mode, factor, nprint,
157 &nfev, &njev, ipvt, &wa[n], &wa[(n << 1)], &
158 wa[n * 3], &wa[(n << 2)], &wa[n * 5]);
159 if (info == 8) {
160 info = 4;
161 }
162 return info;
163
164 /* last card of subroutine lmder1. */
165
166 } /* lmder1_ */
167
168