| =pod |
| |
| =begin html |
| |
| <link rel="stylesheet" href="podstyle.css" type="text/css" /> |
| |
| =end html |
| |
| =head1 NAME |
| |
| lmmin - Levenberg-Marquardt least-squares minimization |
| |
| |
| =head1 SYNOPSIS |
| |
| B<#include <lmmin.h>> |
| |
| B<void lmmin( const int> I<n_par>B<, double *>I<par>B<, const int> I<m_dat>B<, |
| constS< >void *>I<data>B<, |
| void *>I<evaluate>B<( |
| constS< >double *>I<par>B<, const int >I<m_dat>B<, |
| constS< >void *>I<data>B<, double *>I<fvec>B<, int *>I<userbreak>B<), |
| constS< >lm_control_struct *>I<control>B<, |
| lm_status_struct *>I<status>B< );> |
| |
| B<extern const lm_control_struct lm_control_double;> |
| |
| B<extern const lm_control_struct lm_control_float;> |
| |
| B<extern const char *lm_infmsg[];> |
| |
| B<extern const char *lm_shortmsg[];> |
| |
| =head1 DESCRIPTION |
| |
| B<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>(). |
| On success, I<par> represents a local minimum, not necessarily a global one; it may depend on its starting value. |
| |
| For applications in curve fitting, the wrapper function B<lmcurve(3)> offers a simplified API. |
| |
| The 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. |
| |
| Function arguments: |
| |
| =over |
| |
| =item I<n_par> |
| |
| Number of free variables. |
| Length of parameter vector I<par>. |
| |
| =item I<par> |
| |
| Parameter vector. |
| On input, it must contain a reasonable guess. |
| On output, it contains the solution found to minimize ||I<fvec>||. |
| |
| =item I<m_dat> |
| |
| Length of vector I<fvec>. |
| Must statisfy I<n_par> <= I<m_dat>. |
| |
| =item I<data> |
| |
| This pointer is ignored by the fit algorithm, |
| except for appearing as an argument in all calls to the user-supplied |
| routine I<evaluate>. |
| |
| =item I<evaluate> |
| |
| Pointer to a user-supplied function that computes I<m_dat> elements of vector I<fvec> for a given parameter vector I<par>. |
| If I<evaluate> return with *I<userbreak> set to a negative value, B<lmmin()> will interrupt the fitting and terminate. |
| |
| =item I<control> |
| |
| Parameter collection for tuning the fit procedure. |
| In most cases, the default &I<lm_control_double> is adequate. |
| If I<f> is only computed with single-precision accuracy, |
| I<&lm_control_float> should be used. |
| See also below, NOTES on initializing parameter records. |
| |
| I<control> has the following members (for more details, see the source file I<lmstruct.h>): |
| |
| =over |
| |
| =item B<double> I<control.ftol> |
| |
| Relative error desired in the sum of squares. |
| Recommended setting: somewhat above machine precision; less if I<fvec> is computed with reduced accuracy. |
| |
| =item B<double> I<control.xtol> |
| |
| Relative error between last two approximations. |
| Recommended setting: as I<ftol>. |
| |
| =item B<double> I<control.gtol> |
| |
| A measure for degeneracy. |
| Recommended setting: as I<ftol>. |
| |
| =item B<double> I<control.epsilon> |
| |
| Step used to calculate the Jacobian. |
| Recommended setting: as I<ftol>, but definitely less than the accuracy of I<fvec>. |
| |
| =item B<double> I<control.stepbound> |
| |
| Initial bound to steps in the outer loop, generally between 0.01 and 100; recommended value is 100. |
| |
| =item B<int> I<control.patience> |
| |
| Used to set the maximum number of function evaluations to patience*n_par. |
| |
| =item B<int> I<control.scale_diag> |
| |
| Logical switch (0 or 1). |
| If 1, then scale parameters to their initial value. |
| This is the recommended setting. |
| |
| =item B<FILE*> I<control.msgfile> |
| |
| Progress messages will be written to this file. |
| Typically I<stdout> or I<stderr>. |
| The value I<NULL> will be interpreted as I<stdout>. |
| |
| =item B<int> I<control.verbosity> |
| |
| If nonzero, some progress information from within the LM algorithm |
| is written to control.stream. |
| |
| =item B<int> I<control.n_maxpri> |
| |
| -1, or maximum number of parameters to print. |
| |
| =item B<int> I<control.m_maxpri> |
| |
| -1, or maximum number of residuals to print. |
| |
| =back |
| |
| =item I<status> |
| |
| A record used to return information about the minimization process: |
| |
| =over |
| |
| =item B<double> I<status.fnorm> |
| |
| Norm of the vector I<fvec>; |
| |
| =item B<int> I<status.nfev> |
| |
| Actual number of iterations; |
| |
| =item B<int> I<status.outcome> |
| |
| Status of minimization; |
| for the corresponding text message, print I<lm_infmsg>B<[>I<status.outcome>B<]>; |
| for a short code, print I<lm_shortmsg>B<[>I<status.outcome>B<]>. |
| |
| =item B<int> I<status.userbreak> |
| |
| Set when termination has been forced by the user-supplied routine I<evaluate>. |
| |
| =back |
| |
| =back |
| |
| |
| =head1 NOTES |
| |
| =head2 Initializing parameter records. |
| |
| The parameter record I<control> should always be initialized |
| from supplied default records: |
| |
| lm_control_struct control = lm_control_double; /* or _float */ |
| |
| After this, parameters may be overwritten: |
| |
| control.patience = 500; /* allow more iterations */ |
| control.verbosity = 15; /* for verbose monitoring */ |
| |
| An application written this way is guaranteed to work even if new parameters |
| are added to I<lm_control_struct>. |
| |
| Conversely, addition of parameters is not considered an API change; it may happen without increment of the major version number. |
| |
| =head1 EXAMPLES |
| |
| =head2 Fitting a surface |
| |
| Fit a data set y(t) by a function f(t;p) where t is a two-dimensional vector: |
| |
| #include "lmmin.h" |
| #include <stdio.h> |
| |
| /* fit model: a plane p0 + p1*tx + p2*tz */ |
| double f( double tx, double tz, const double *p ) |
| { |
| return p[0] + p[1]*tx + p[2]*tz; |
| } |
| |
| /* data structure to transmit data arays and fit model */ |
| typedef struct { |
| double *tx, *tz; |
| double *y; |
| double (*f)( double tx, double tz, const double *p ); |
| } data_struct; |
| |
| /* function evaluation, determination of residues */ |
| void evaluate_surface( const double *par, int m_dat, |
| const void *data, double *fvec, int *userbreak ) |
| { |
| /* for readability, explicit type conversion */ |
| data_struct *D; |
| D = (data_struct*)data; |
| |
| int i; |
| for ( i = 0; i < m_dat; i++ ) |
| fvec[i] = D->y[i] - D->f( D->tx[i], D->tz[i], par ); |
| } |
| |
| int main() |
| { |
| /* parameter vector */ |
| int n_par = 3; /* number of parameters in model function f */ |
| double par[3] = { -1, 0, 1 }; /* arbitrary starting value */ |
| |
| /* data points */ |
| int m_dat = 4; |
| double tx[4] = { -1, -1, 1, 1 }; |
| double tz[4] = { -1, 1, -1, 1 }; |
| double y[4] = { 0, 1, 1, 2 }; |
| |
| data_struct data = { tx, tz, y, f }; |
| |
| /* auxiliary parameters */ |
| lm_status_struct status; |
| lm_control_struct control = lm_control_double; |
| control.verbosity = 3; |
| |
| /* perform the fit */ |
| printf( "Fitting:\n" ); |
| lmmin( n_par, par, m_dat, (const void*) &data, evaluate_surface, |
| &control, &status ); |
| |
| /* print results */ |
| printf( "\nResults:\n" ); |
| printf( "status after %d function evaluations:\n %s\n", |
| status.nfev, lm_infmsg[status.outcome] ); |
| |
| printf("obtained parameters:\n"); |
| int i; |
| for ( i=0; i<n_par; ++i ) |
| printf(" par[%i] = %12g\n", i, par[i]); |
| printf("obtained norm:\n %12g\n", status.fnorm ); |
| |
| printf("fitting data as follows:\n"); |
| double ff; |
| for ( i=0; i<m_dat; ++i ){ |
| ff = f(tx[i], tz[i], par); |
| printf( " t[%2d]=%12g,%12g y=%12g fit=%12g residue=%12g\n", |
| i, tx[i], tz[i], y[i], ff, y[i] - ff ); |
| } |
| |
| return 0; |
| } |
| |
| =head2 More examples |
| |
| For more examples, see the homepage and directories demo/ and test/ in the source distribution. |
| |
| =head1 COPYING |
| |
| Copyright (C): |
| 1980-1999 University of Chicago |
| 2004-2015 Joachim Wuttke, Forschungszentrum Juelich GmbH |
| |
| Software: FreeBSD License |
| |
| Documentation: Creative Commons Attribution Share Alike |
| |
| |
| =head1 SEE ALSO |
| |
| =begin html |
| |
| <a href="http://apps.jcns.fz-juelich.de/man/lmcurve.html"><b>lmcurve</b>(3)</a> |
| |
| =end html |
| |
| =begin man |
| |
| \fBlmcurve\fR(3) |
| .PP |
| |
| =end man |
| |
| Homepage: http://apps.jcns.fz-juelich.de/lmfit |
| |
| =head1 BUGS |
| |
| Please send bug reports and suggestions to the author <j.wuttke@fz-juelich.de>. |