1 /* hybrd1.f -- translated by f2c (version 20020621).
2 You must link the resulting object file with the libraries:
3 -lf2c -lm (in that order)
4 */
5
6 #include "minpack.h"
7 #include <math.h>
8 #include "minpackP.h"
9
10 __minpack_attr__
__minpack_func__(hybrd1)11 void __minpack_func__(hybrd1)(__minpack_decl_fcn_nn__ const int *n, real *x, real *
12 fvec, const real *tol, int *info, real *wa, const int *lwa)
13 {
14 /* Initialized data */
15
16 const real factor = 100.;
17
18 /* System generated locals */
19 int i__1;
20
21 /* Local variables */
22 int j, ml, lr, mu, mode, nfev;
23 real xtol;
24 int index;
25 real epsfcn;
26 int maxfev, nprint;
27
28 /* ********** */
29
30 /* subroutine hybrd1 */
31
32 /* the purpose of hybrd1 is to find a zero of a system of */
33 /* n nonlinear functions in n variables by a modification */
34 /* of the powell hybrid method. this is done by using the */
35 /* more general nonlinear equation solver hybrd. the user */
36 /* must provide a subroutine which calculates the functions. */
37 /* the jacobian is then calculated by a forward-difference */
38 /* approximation. */
39
40 /* the subroutine statement is */
41
42 /* subroutine hybrd1(fcn,n,x,fvec,tol,info,wa,lwa) */
43
44 /* where */
45
46 /* fcn is the name of the user-supplied subroutine which */
47 /* calculates the functions. fcn must be declared */
48 /* in an external statement in the user calling */
49 /* program, and should be written as follows. */
50
51 /* subroutine fcn(n,x,fvec,iflag) */
52 /* integer n,iflag */
53 /* double precision x(n),fvec(n) */
54 /* ---------- */
55 /* calculate the functions at x and */
56 /* return this vector in fvec. */
57 /* --------- */
58 /* return */
59 /* end */
60
61 /* the value of iflag should not be changed by fcn unless */
62 /* the user wants to terminate execution of hybrd1. */
63 /* in this case set iflag to a negative integer. */
64
65 /* n is a positive integer input variable set to the number */
66 /* of functions and variables. */
67
68 /* x is an array of length n. on input x must contain */
69 /* an initial estimate of the solution vector. on output x */
70 /* contains the final estimate of the solution vector. */
71
72 /* fvec is an output array of length n which contains */
73 /* the functions evaluated at the output x. */
74
75 /* tol is a nonnegative input variable. termination occurs */
76 /* when the algorithm estimates that the relative error */
77 /* between x and the solution is at most tol. */
78
79 /* info is an integer output variable. if the user has */
80 /* terminated execution, info is set to the (negative) */
81 /* value of iflag. see description of fcn. otherwise, */
82 /* info is set as follows. */
83
84 /* info = 0 improper input parameters. */
85
86 /* info = 1 algorithm estimates that the relative error */
87 /* between x and the solution is at most tol. */
88
89 /* info = 2 number of calls to fcn has reached or exceeded */
90 /* 200*(n+1). */
91
92 /* info = 3 tol is too small. no further improvement in */
93 /* the approximate solution x is possible. */
94
95 /* info = 4 iteration is not making good progress. */
96
97 /* wa is a work array of length lwa. */
98
99 /* lwa is a positive integer input variable not less than */
100 /* (n*(3*n+13))/2. */
101
102 /* subprograms called */
103
104 /* user-supplied ...... fcn */
105
106 /* minpack-supplied ... hybrd */
107
108 /* argonne national laboratory. minpack project. march 1980. */
109 /* burton s. garbow, kenneth e. hillstrom, jorge j. more */
110
111 /* ********** */
112 /* Parameter adjustments */
113 --fvec;
114 --x;
115 --wa;
116
117 /* Function Body */
118 *info = 0;
119
120 /* check the input parameters for errors. */
121
122 if (*n <= 0 || *tol < 0. || *lwa < *n * (*n * 3 + 13) / 2) {
123 /* goto L20; */
124 return;
125 }
126
127 /* call hybrd. */
128
129 maxfev = (*n + 1) * 200;
130 xtol = *tol;
131 ml = *n - 1;
132 mu = *n - 1;
133 epsfcn = 0.;
134 mode = 2;
135 i__1 = *n;
136 for (j = 1; j <= i__1; ++j) {
137 wa[j] = 1.;
138 /* L10: */
139 }
140 nprint = 0;
141 lr = *n * (*n + 1) / 2;
142 index = *n * 6 + lr;
143 __minpack_func__(hybrd)(__minpack_param_fcn_nn__ n, &x[1], &fvec[1], &xtol, &maxfev, &ml, &mu, &epsfcn, &
144 wa[1], &mode, &factor, &nprint, info, &nfev, &wa[index + 1], n, &
145 wa[*n * 6 + 1], &lr, &wa[*n + 1], &wa[(*n << 1) + 1], &wa[*n * 3
146 + 1], &wa[(*n << 2) + 1], &wa[*n * 5 + 1]);
147 if (*info == 5) {
148 *info = 4;
149 }
150 /* L20: */
151 return;
152
153 /* last card of subroutine hybrd1. */
154
155 } /* hybrd1_ */
156
157