NAME

lmmin - Levenberg-Marquardt least-squares minimization (simple/legacy API without error estimates)

SYNOPSIS

#include <lmmin.h>

void lmmin( const int n_par, double *par, const int m_dat,
const void *y, const void *data,
void *evaluate(
const double *par, const int m_dat,
const void *data, double *fvec, int *userbreak),
const lm_control_struct *control,
lm_status_struct *status );

extern const lm_control_struct lm_control_double;

extern const lm_control_struct lm_control_float;

extern const char *lm_infmsg[];

extern const char *lm_shortmsg[];

DESCRIPTION

lmmin() determines a vector par that minimizes the sum of squared elements of fvec-y. The vector fvec is computed by a user-supplied function evaluate(); the vector y contains user-provided values. On success, par represents a local minimum, not necessarily a global one; it may depend on its starting value.

This is a simple wrapper of the function lmmin2(3), which also returns error estimates. Conversely, the function lmcurve(3) provides an even simpler wrapper, for use in curve fitting.

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:

n_par: Number of free variables. Length of parameter vector par.
par: Parameter vector. On input, it must contain a reasonable guess. On output, it contains the solution found to minimize ||fvec||.
m_dat: Length of vector fvec. Must statisfy n_par <= m_dat.
y: Input vector of length m_dat. May also be the null pointer; in this case, lmmin() minimizes the squared sum of fvec instead of fvec-y.
data: This pointer is ignored by the fit algorithm, except for appearing as an argument in all calls to the user-supplied routine evaluate.
evaluate: Pointer to a user-supplied function that computes m_dat elements of vector fvec for a given parameter vector par. If evaluate return with *userbreak set to a negative value, lmmin() will interrupt the fitting and terminate.
control: Parameter collection for tuning the fit procedure. In most cases, the default &lm_control_double is adequate. If f is only computed with single-precision accuracy, &lm_control_float should be used. See also below, NOTES on initializing parameter records.
control has the following members (for more details, see the source file lmstruct.h):

double control.ftol: Relative error desired in the sum of squares. Recommended setting: somewhat above machine precision; less if fvec is computed with reduced accuracy.
double control.xtol: Relative error between last two approximations. Recommended setting: as ftol.
double control.gtol: A measure for degeneracy. Recommended setting: as ftol.
double control.epsilon: Step used to calculate the Jacobian. Recommended setting: as ftol, but definitely less than the accuracy of fvec.
double control.stepbound: Initial bound to steps in the outer loop, generally between 0.01 and 100; recommended value is 100.
int control.patience: Used to set the maximum number of function evaluations to patience*n_par.
int control.scale_diag: Logical switch (0 or 1). If 1, then scale parameters to their initial value. This is the recommended setting.
FILE* control.msgfile: Progress messages will be written to this file. Typically stdout or stderr. The value NULL will be interpreted as stdout.
int control.verbosity: If nonzero, some progress information from within the LM algorithm is written to control.stream.
int control.n_maxpri: -1, or maximum number of parameters to print.
int control.m_maxpri: -1, or maximum number of residuals to print.

status: A record used to return information about the minimization process:

double status.fnorm: Norm of the vector fvec;
int status.nfev: Actual number of iterations;
int status.outcome: Status of minimization; for the corresponding text message, print lm_infmsg[status.outcome]; for a short code, print lm_shortmsg[status.outcome].
int status.userbreak: Set when termination has been forced by the user-supplied routine evaluate.

NOTES

Initializing parameter records.

The parameter record 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 lm_control_struct.

Conversely, addition of parameters is not considered an API change; it may happen without increment of the major version number.

EXAMPLES

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;
    }

More examples

For more examples, see the homepage and directories demo/ and test/ in the source distribution.

COPYING

Software: FreeBSD License

Documentation: Creative Commons Attribution Share Alike

BUGS

Please send bug reports and suggestions to the author <j.wuttke@fz-juelich.de>.