User Tools

Site Tools


solving nonlinear equations with lmfit

Levenberg-Marquardt least-squares minimization can also be used to solve systems of nonlinear equations.

There are n equations in n unkowns p_j. We write each equation in the form f_i(p)=0. The resulting system of equations, f(p)=0, can be solved by minimizing ||f(p)||.

In the following example, we search for an intersection of the unit circle x^2+y^2=1 with the standard parabola y=x^2. The dimension of this problem is n=2, and we have p=(x, y). Note that x and y have completely different meaning than in other typical applications of lmfit.

/* demo/nonlin1.c */
 
#include "lmmin.h"
#include <stdio.h>
#include <stdlib.h>
 
/* componentwise computation of f(p) */
void evaluate_nonlin1(
    const double *p, int n, const void *data, double *f, int *info )
{
    f[0] = p[0]*p[0] + p[1]*p[1] - 1; /* unit circle    x^2+y^2=1 */
    f[1] = p[1] - p[0]*p[0];          /* standard parabola  y=x^2 */
}
 
 
int main( int argc, char **argv )
{
    int n = 2;   /* dimension of the problem */
    double p[2]; /* parameter vector p=(x,y) */
 
    /* auxiliary parameters */
    lm_control_struct control = lm_control_double;
    lm_princon_struct princon = lm_princon_std;
    lm_status_struct  status;
 
    princon.form  =  1;
    princon.flags = 15;
 
    /* get start values from command line */
    if( argc!=3 ){
        fprintf( stderr, "usage: nonlin1 x_start y_start\n" );
        exit(-1);
    }
    p[0] = atof( argv[1] );
    p[1] = atof( argv[2] );
 
    /* the minimization */
    printf( "Minimization:\n" );
    lmmin( n, p, n, NULL, evaluate_nonlin1,
           lm_printout_std, &control, &princon, &status );
 
    /* print results */
    printf( "\n" );
    printf( "lmmin status after %d function evaluations:\n  %s\n",
            status.nfev, lm_infmsg[status.info] );
 
    printf( "\n" );
    printf("Solution:\n");
    printf("  x = %19.11f\n", p[0]);
    printf("  y = %19.11f\n", p[1]);
    printf("  d = %19.11f => ", status.fnorm);
 
    /* convergence of lmfit is not enough to ensure validity of the solution */
    if( status.fnorm >= control.ftol )
        printf( "not a valid solution, try other starting values\n" );
    else
        printf( "valid, though not the only solution: "
                "try other starting values\n" );
 
    return 0;
}

The program takes the starting values x_0, y_0 as command-line arguments. As obvious from the symmetry of problem, the solution depends on starting values.