xref: /aosp_15_r20/external/lmfit/man/lmmin.pod (revision 5ddc57e5d924f146ab5fd87df586563e2270da38)
1*5ddc57e5SXin Li=pod
2*5ddc57e5SXin Li
3*5ddc57e5SXin Li=begin html
4*5ddc57e5SXin Li
5*5ddc57e5SXin Li<link rel="stylesheet" href="podstyle.css" type="text/css" />
6*5ddc57e5SXin Li
7*5ddc57e5SXin Li=end html
8*5ddc57e5SXin Li
9*5ddc57e5SXin Li=head1 NAME
10*5ddc57e5SXin Li
11*5ddc57e5SXin Lilmmin - Levenberg-Marquardt least-squares minimization
12*5ddc57e5SXin Li
13*5ddc57e5SXin Li
14*5ddc57e5SXin Li=head1 SYNOPSIS
15*5ddc57e5SXin Li
16*5ddc57e5SXin LiB<#include <lmmin.h>>
17*5ddc57e5SXin Li
18*5ddc57e5SXin LiB<void lmmin( const int> I<n_par>B<, double *>I<par>B<, const int> I<m_dat>B<,
19*5ddc57e5SXin Li            constS< >void *>I<data>B<,
20*5ddc57e5SXin Li            void *>I<evaluate>B<(
21*5ddc57e5SXin Li                 constS< >double *>I<par>B<, const int >I<m_dat>B<,
22*5ddc57e5SXin Li                 constS< >void *>I<data>B<, double *>I<fvec>B<, int *>I<userbreak>B<),
23*5ddc57e5SXin Li            constS< >lm_control_struct *>I<control>B<,
24*5ddc57e5SXin Li            lm_status_struct *>I<status>B< );>
25*5ddc57e5SXin Li
26*5ddc57e5SXin LiB<extern const lm_control_struct lm_control_double;>
27*5ddc57e5SXin Li
28*5ddc57e5SXin LiB<extern const lm_control_struct lm_control_float;>
29*5ddc57e5SXin Li
30*5ddc57e5SXin LiB<extern const char *lm_infmsg[];>
31*5ddc57e5SXin Li
32*5ddc57e5SXin LiB<extern const char *lm_shortmsg[];>
33*5ddc57e5SXin Li
34*5ddc57e5SXin Li=head1 DESCRIPTION
35*5ddc57e5SXin Li
36*5ddc57e5SXin LiB<lmmin()> determines a vector I<par> that minimizes the sum of squared elements of a vector I<fvec> that is computed by a user-supplied function I<evaluate>().
37*5ddc57e5SXin LiOn success, I<par> represents a local minimum, not necessarily a global one; it may depend on its starting value.
38*5ddc57e5SXin Li
39*5ddc57e5SXin LiFor applications in curve fitting, the wrapper function B<lmcurve(3)> offers a simplified API.
40*5ddc57e5SXin Li
41*5ddc57e5SXin LiThe Levenberg-Marquardt minimization starts with a steepest-descent exploration of the parameter space, and achieves rapid convergence by crossing over into the Newton-Gauss method.
42*5ddc57e5SXin Li
43*5ddc57e5SXin LiFunction arguments:
44*5ddc57e5SXin Li
45*5ddc57e5SXin Li=over
46*5ddc57e5SXin Li
47*5ddc57e5SXin Li=item I<n_par>
48*5ddc57e5SXin Li
49*5ddc57e5SXin LiNumber of free variables.
50*5ddc57e5SXin LiLength of parameter vector I<par>.
51*5ddc57e5SXin Li
52*5ddc57e5SXin Li=item I<par>
53*5ddc57e5SXin Li
54*5ddc57e5SXin LiParameter vector.
55*5ddc57e5SXin LiOn input, it must contain a reasonable guess.
56*5ddc57e5SXin LiOn output, it contains the solution found to minimize ||I<fvec>||.
57*5ddc57e5SXin Li
58*5ddc57e5SXin Li=item I<m_dat>
59*5ddc57e5SXin Li
60*5ddc57e5SXin LiLength of vector I<fvec>.
61*5ddc57e5SXin LiMust statisfy I<n_par> <= I<m_dat>.
62*5ddc57e5SXin Li
63*5ddc57e5SXin Li=item I<data>
64*5ddc57e5SXin Li
65*5ddc57e5SXin LiThis pointer is ignored by the fit algorithm,
66*5ddc57e5SXin Liexcept for appearing as an argument in all calls to the user-supplied
67*5ddc57e5SXin Liroutine I<evaluate>.
68*5ddc57e5SXin Li
69*5ddc57e5SXin Li=item I<evaluate>
70*5ddc57e5SXin Li
71*5ddc57e5SXin LiPointer to a user-supplied function that computes I<m_dat> elements of vector I<fvec> for a given parameter vector I<par>.
72*5ddc57e5SXin LiIf I<evaluate> return with *I<userbreak> set to a negative value, B<lmmin()> will interrupt the fitting and terminate.
73*5ddc57e5SXin Li
74*5ddc57e5SXin Li=item I<control>
75*5ddc57e5SXin Li
76*5ddc57e5SXin LiParameter collection for tuning the fit procedure.
77*5ddc57e5SXin LiIn most cases, the default &I<lm_control_double> is adequate.
78*5ddc57e5SXin LiIf I<f> is only computed with single-precision accuracy,
79*5ddc57e5SXin LiI<&lm_control_float> should be used.
80*5ddc57e5SXin LiSee also below, NOTES on initializing parameter records.
81*5ddc57e5SXin Li
82*5ddc57e5SXin LiI<control> has the following members (for more details, see the source file I<lmstruct.h>):
83*5ddc57e5SXin Li
84*5ddc57e5SXin Li=over
85*5ddc57e5SXin Li
86*5ddc57e5SXin Li=item B<double> I<control.ftol>
87*5ddc57e5SXin Li
88*5ddc57e5SXin LiRelative error desired in the sum of squares.
89*5ddc57e5SXin LiRecommended setting: somewhat above machine precision; less if I<fvec> is computed with reduced accuracy.
90*5ddc57e5SXin Li
91*5ddc57e5SXin Li=item B<double> I<control.xtol>
92*5ddc57e5SXin Li
93*5ddc57e5SXin LiRelative error between last two approximations.
94*5ddc57e5SXin LiRecommended setting: as I<ftol>.
95*5ddc57e5SXin Li
96*5ddc57e5SXin Li=item B<double> I<control.gtol>
97*5ddc57e5SXin Li
98*5ddc57e5SXin LiA measure for degeneracy.
99*5ddc57e5SXin LiRecommended setting: as I<ftol>.
100*5ddc57e5SXin Li
101*5ddc57e5SXin Li=item B<double> I<control.epsilon>
102*5ddc57e5SXin Li
103*5ddc57e5SXin LiStep used to calculate the Jacobian.
104*5ddc57e5SXin LiRecommended setting: as I<ftol>, but definitely less than the accuracy of I<fvec>.
105*5ddc57e5SXin Li
106*5ddc57e5SXin Li=item B<double> I<control.stepbound>
107*5ddc57e5SXin Li
108*5ddc57e5SXin LiInitial bound to steps in the outer loop, generally between 0.01 and 100; recommended value is 100.
109*5ddc57e5SXin Li
110*5ddc57e5SXin Li=item B<int> I<control.patience>
111*5ddc57e5SXin Li
112*5ddc57e5SXin LiUsed to set the maximum number of function evaluations to patience*n_par.
113*5ddc57e5SXin Li
114*5ddc57e5SXin Li=item B<int> I<control.scale_diag>
115*5ddc57e5SXin Li
116*5ddc57e5SXin LiLogical switch (0 or 1).
117*5ddc57e5SXin LiIf 1, then scale parameters to their initial value.
118*5ddc57e5SXin LiThis is the recommended setting.
119*5ddc57e5SXin Li
120*5ddc57e5SXin Li=item B<FILE*> I<control.msgfile>
121*5ddc57e5SXin Li
122*5ddc57e5SXin LiProgress messages will be written to this file.
123*5ddc57e5SXin LiTypically I<stdout> or I<stderr>.
124*5ddc57e5SXin LiThe value I<NULL> will be interpreted as I<stdout>.
125*5ddc57e5SXin Li
126*5ddc57e5SXin Li=item B<int> I<control.verbosity>
127*5ddc57e5SXin Li
128*5ddc57e5SXin LiIf nonzero, some progress information from within the LM algorithm
129*5ddc57e5SXin Liis written to control.stream.
130*5ddc57e5SXin Li
131*5ddc57e5SXin Li=item B<int> I<control.n_maxpri>
132*5ddc57e5SXin Li
133*5ddc57e5SXin Li-1, or maximum number of parameters to print.
134*5ddc57e5SXin Li
135*5ddc57e5SXin Li=item B<int> I<control.m_maxpri>
136*5ddc57e5SXin Li
137*5ddc57e5SXin Li-1, or maximum number of residuals to print.
138*5ddc57e5SXin Li
139*5ddc57e5SXin Li=back
140*5ddc57e5SXin Li
141*5ddc57e5SXin Li=item I<status>
142*5ddc57e5SXin Li
143*5ddc57e5SXin LiA record used to return information about the minimization process:
144*5ddc57e5SXin Li
145*5ddc57e5SXin Li=over
146*5ddc57e5SXin Li
147*5ddc57e5SXin Li=item B<double> I<status.fnorm>
148*5ddc57e5SXin Li
149*5ddc57e5SXin LiNorm of the vector I<fvec>;
150*5ddc57e5SXin Li
151*5ddc57e5SXin Li=item B<int> I<status.nfev>
152*5ddc57e5SXin Li
153*5ddc57e5SXin LiActual number of iterations;
154*5ddc57e5SXin Li
155*5ddc57e5SXin Li=item B<int> I<status.outcome>
156*5ddc57e5SXin Li
157*5ddc57e5SXin LiStatus of minimization;
158*5ddc57e5SXin Lifor the corresponding text message, print I<lm_infmsg>B<[>I<status.outcome>B<]>;
159*5ddc57e5SXin Lifor a short code, print I<lm_shortmsg>B<[>I<status.outcome>B<]>.
160*5ddc57e5SXin Li
161*5ddc57e5SXin Li=item B<int> I<status.userbreak>
162*5ddc57e5SXin Li
163*5ddc57e5SXin LiSet when termination has been forced by the user-supplied routine I<evaluate>.
164*5ddc57e5SXin Li
165*5ddc57e5SXin Li=back
166*5ddc57e5SXin Li
167*5ddc57e5SXin Li=back
168*5ddc57e5SXin Li
169*5ddc57e5SXin Li
170*5ddc57e5SXin Li=head1 NOTES
171*5ddc57e5SXin Li
172*5ddc57e5SXin Li=head2 Initializing parameter records.
173*5ddc57e5SXin Li
174*5ddc57e5SXin LiThe parameter record I<control> should always be initialized
175*5ddc57e5SXin Lifrom supplied default records:
176*5ddc57e5SXin Li
177*5ddc57e5SXin Li    lm_control_struct control = lm_control_double; /* or _float */
178*5ddc57e5SXin Li
179*5ddc57e5SXin LiAfter this, parameters may be overwritten:
180*5ddc57e5SXin Li
181*5ddc57e5SXin Li    control.patience = 500; /* allow more iterations */
182*5ddc57e5SXin Li    control.verbosity = 15; /* for verbose monitoring */
183*5ddc57e5SXin Li
184*5ddc57e5SXin LiAn application written this way is guaranteed to work even if new parameters
185*5ddc57e5SXin Liare added to I<lm_control_struct>.
186*5ddc57e5SXin Li
187*5ddc57e5SXin LiConversely, addition of parameters is not considered an API change; it may happen without increment of the major version number.
188*5ddc57e5SXin Li
189*5ddc57e5SXin Li=head1 EXAMPLES
190*5ddc57e5SXin Li
191*5ddc57e5SXin Li=head2 Fitting a surface
192*5ddc57e5SXin Li
193*5ddc57e5SXin LiFit a data set y(t) by a function f(t;p) where t is a two-dimensional vector:
194*5ddc57e5SXin Li
195*5ddc57e5SXin Li    #include "lmmin.h"
196*5ddc57e5SXin Li    #include <stdio.h>
197*5ddc57e5SXin Li
198*5ddc57e5SXin Li    /* fit model: a plane p0 + p1*tx + p2*tz */
199*5ddc57e5SXin Li    double f( double tx, double tz, const double *p )
200*5ddc57e5SXin Li    {
201*5ddc57e5SXin Li        return p[0] + p[1]*tx + p[2]*tz;
202*5ddc57e5SXin Li    }
203*5ddc57e5SXin Li
204*5ddc57e5SXin Li    /* data structure to transmit data arays and fit model */
205*5ddc57e5SXin Li    typedef struct {
206*5ddc57e5SXin Li        double *tx, *tz;
207*5ddc57e5SXin Li        double *y;
208*5ddc57e5SXin Li        double (*f)( double tx, double tz, const double *p );
209*5ddc57e5SXin Li    } data_struct;
210*5ddc57e5SXin Li
211*5ddc57e5SXin Li    /* function evaluation, determination of residues */
212*5ddc57e5SXin Li    void evaluate_surface( const double *par, int m_dat,
213*5ddc57e5SXin Li        const void *data, double *fvec, int *userbreak )
214*5ddc57e5SXin Li    {
215*5ddc57e5SXin Li        /* for readability, explicit type conversion */
216*5ddc57e5SXin Li        data_struct *D;
217*5ddc57e5SXin Li        D = (data_struct*)data;
218*5ddc57e5SXin Li
219*5ddc57e5SXin Li        int i;
220*5ddc57e5SXin Li        for ( i = 0; i < m_dat; i++ )
221*5ddc57e5SXin Li    	fvec[i] = D->y[i] - D->f( D->tx[i], D->tz[i], par );
222*5ddc57e5SXin Li    }
223*5ddc57e5SXin Li
224*5ddc57e5SXin Li    int main()
225*5ddc57e5SXin Li    {
226*5ddc57e5SXin Li        /* parameter vector */
227*5ddc57e5SXin Li        int n_par = 3; /* number of parameters in model function f */
228*5ddc57e5SXin Li        double par[3] = { -1, 0, 1 }; /* arbitrary starting value */
229*5ddc57e5SXin Li
230*5ddc57e5SXin Li        /* data points */
231*5ddc57e5SXin Li        int m_dat = 4;
232*5ddc57e5SXin Li        double tx[4] = { -1, -1,  1,  1 };
233*5ddc57e5SXin Li        double tz[4] = { -1,  1, -1,  1 };
234*5ddc57e5SXin Li        double y[4]  = {  0,  1,  1,  2 };
235*5ddc57e5SXin Li
236*5ddc57e5SXin Li        data_struct data = { tx, tz, y, f };
237*5ddc57e5SXin Li
238*5ddc57e5SXin Li        /* auxiliary parameters */
239*5ddc57e5SXin Li        lm_status_struct status;
240*5ddc57e5SXin Li        lm_control_struct control = lm_control_double;
241*5ddc57e5SXin Li        control.verbosity = 3;
242*5ddc57e5SXin Li
243*5ddc57e5SXin Li        /* perform the fit */
244*5ddc57e5SXin Li        printf( "Fitting:\n" );
245*5ddc57e5SXin Li        lmmin( n_par, par, m_dat, (const void*) &data, evaluate_surface,
246*5ddc57e5SXin Li               &control, &status );
247*5ddc57e5SXin Li
248*5ddc57e5SXin Li        /* print results */
249*5ddc57e5SXin Li        printf( "\nResults:\n" );
250*5ddc57e5SXin Li        printf( "status after %d function evaluations:\n  %s\n",
251*5ddc57e5SXin Li                status.nfev, lm_infmsg[status.outcome] );
252*5ddc57e5SXin Li
253*5ddc57e5SXin Li        printf("obtained parameters:\n");
254*5ddc57e5SXin Li        int i;
255*5ddc57e5SXin Li        for ( i=0; i<n_par; ++i )
256*5ddc57e5SXin Li    	printf("  par[%i] = %12g\n", i, par[i]);
257*5ddc57e5SXin Li        printf("obtained norm:\n  %12g\n", status.fnorm );
258*5ddc57e5SXin Li
259*5ddc57e5SXin Li        printf("fitting data as follows:\n");
260*5ddc57e5SXin Li        double ff;
261*5ddc57e5SXin Li        for ( i=0; i<m_dat; ++i ){
262*5ddc57e5SXin Li            ff = f(tx[i], tz[i], par);
263*5ddc57e5SXin Li            printf( "  t[%2d]=%12g,%12g y=%12g fit=%12g residue=%12g\n",
264*5ddc57e5SXin Li                    i, tx[i], tz[i], y[i], ff, y[i] - ff );
265*5ddc57e5SXin Li        }
266*5ddc57e5SXin Li
267*5ddc57e5SXin Li        return 0;
268*5ddc57e5SXin Li    }
269*5ddc57e5SXin Li
270*5ddc57e5SXin Li=head2 More examples
271*5ddc57e5SXin Li
272*5ddc57e5SXin LiFor more examples, see the homepage and directories demo/ and test/ in the source distribution.
273*5ddc57e5SXin Li
274*5ddc57e5SXin Li=head1 COPYING
275*5ddc57e5SXin Li
276*5ddc57e5SXin LiCopyright (C):
277*5ddc57e5SXin Li   1980-1999 University of Chicago
278*5ddc57e5SXin Li   2004-2015 Joachim Wuttke, Forschungszentrum Juelich GmbH
279*5ddc57e5SXin Li
280*5ddc57e5SXin LiSoftware: FreeBSD License
281*5ddc57e5SXin Li
282*5ddc57e5SXin LiDocumentation: Creative Commons Attribution Share Alike
283*5ddc57e5SXin Li
284*5ddc57e5SXin Li
285*5ddc57e5SXin Li=head1 SEE ALSO
286*5ddc57e5SXin Li
287*5ddc57e5SXin Li=begin html
288*5ddc57e5SXin Li
289*5ddc57e5SXin Li<a href="http://apps.jcns.fz-juelich.de/man/lmcurve.html"><b>lmcurve</b>(3)</a>
290*5ddc57e5SXin Li
291*5ddc57e5SXin Li=end html
292*5ddc57e5SXin Li
293*5ddc57e5SXin Li=begin man
294*5ddc57e5SXin Li
295*5ddc57e5SXin Li\fBlmcurve\fR(3)
296*5ddc57e5SXin Li.PP
297*5ddc57e5SXin Li
298*5ddc57e5SXin Li=end man
299*5ddc57e5SXin Li
300*5ddc57e5SXin LiHomepage: http://apps.jcns.fz-juelich.de/lmfit
301*5ddc57e5SXin Li
302*5ddc57e5SXin Li=head1 BUGS
303*5ddc57e5SXin Li
304*5ddc57e5SXin LiPlease send bug reports and suggestions to the author <[email protected]>.
305