Next: , Previous: , Up: Nonlinear Least-Squares Fitting   [Index]


39.4 Initializing the Solver

Function: gsl_multifit_fsolver * gsl_multifit_fsolver_alloc (const gsl_multifit_fsolver_type * T, size_t n, size_t p)

This function returns a pointer to a newly allocated instance of a solver of type T for n observations and p parameters. The number of observations n must be greater than or equal to parameters p.

If there is insufficient memory to create the solver then the function returns a null pointer and the error handler is invoked with an error code of GSL_ENOMEM.

Function: gsl_multifit_fdfsolver * gsl_multifit_fdfsolver_alloc (const gsl_multifit_fdfsolver_type * T, size_t n, size_t p)

This function returns a pointer to a newly allocated instance of a derivative solver of type T for n observations and p parameters. For example, the following code creates an instance of a Levenberg-Marquardt solver for 100 data points and 3 parameters,

const gsl_multifit_fdfsolver_type * T 
    = gsl_multifit_fdfsolver_lmder;
gsl_multifit_fdfsolver * s 
    = gsl_multifit_fdfsolver_alloc (T, 100, 3);

The number of observations n must be greater than or equal to parameters p.

If there is insufficient memory to create the solver then the function returns a null pointer and the error handler is invoked with an error code of GSL_ENOMEM.

Function: gsl_multifit_fdfridge * gsl_multifit_fdfridge_alloc (const gsl_multifit_fdfsolver_type * T, size_t n, size_t p)

This function returns a pointer to a newly allocated instance of a derivative solver of type T for n observations and p parameters. The solver will automatically form the augmented system \tilde{f}(x) and \tilde{J} for ridge (Tikhonov) regression. If there is insufficient memory to create the solver then the function returns a null pointer and the error handler is invoked with an error code of GSL_ENOMEM.

Function: int gsl_multifit_fsolver_set (gsl_multifit_fsolver * s, gsl_multifit_function * f, const gsl_vector * x)

This function initializes, or reinitializes, an existing solver s to use the function f and the initial guess x.

Function: int gsl_multifit_fdfsolver_set (gsl_multifit_fdfsolver * s, gsl_multifit_function_fdf * fdf, const gsl_vector * x)
Function: int gsl_multifit_fdfsolver_wset (gsl_multifit_fdfsolver * s, gsl_multifit_function_fdf * fdf, const gsl_vector * x, const gsl_vector * wts)

These functions initialize, or reinitialize, an existing solver s to use the function and derivative fdf and the initial guess x.

Optionally, a weight vector wts can be given to perform a weighted nonlinear regression. Here, the weighting matrix is W = diag(w_1,w_2,...,w_n). The wts vector is referenced throughout the iteration so it should not be freed by the caller until the iteration terminates.

Function: int gsl_multifit_fdfridge_set (gsl_multifit_fdfridge * s, gsl_multifit_function_fdf * fdf, const gsl_vector * x, const double lambda)
Function: int gsl_multifit_fdfridge_wset (gsl_multifit_fdfridge * s, gsl_multifit_function_fdf * fdf, const gsl_vector * x, const double lambda, const gsl_vector * wts)

This function initializes, or reinitializes, an existing ridge solver s to use the function and derivative fdf and the initial guess x. Here, the regularization matrix is set to L = \lambda I, with \lambda specified in lambda.

Optionally, a weight vector wts can be given to perform a weighted nonlinear regression. Here, the weighting matrix is W = diag(w_1,w_2,...,w_n). The wts vector is referenced throughout the iteration so it should not be freed by the caller until the iteration terminates.

Function: int gsl_multifit_fdfridge_set2 (gsl_multifit_fdfridge * s, gsl_multifit_function_fdf * fdf, const gsl_vector * x, const gsl_vector * lambda)
Function: int gsl_multifit_fdfridge_wset2 (gsl_multifit_fdfridge * s, gsl_multifit_function_fdf * fdf, const gsl_vector * x, const gsl_vector * lambda, const gsl_vector * wts)

This function initializes, or reinitializes, an existing ridge solver s to use the function and derivative fdf and the initial guess x. Here, the regularization matrix is set to L = diag(\lambda_1,\lambda_2,...,\lambda_p), where the \lambda_i are given in lambda.

Optionally, a weight vector wts can be given to perform a weighted nonlinear regression. Here, the weighting matrix is W = diag(w_1,w_2,...,w_n). The wts vector is referenced throughout the iteration so it should not be freed by the caller until the iteration terminates.

Function: int gsl_multifit_fdfridge_set3 (gsl_multifit_fdfridge * s, gsl_multifit_function_fdf * fdf, const gsl_vector * x, const gsl_matrix * L)
Function: int gsl_multifit_fdfridge_wset3 (gsl_multifit_fdfridge * s, gsl_multifit_function_fdf * fdf, const gsl_vector * x, const gsl_matrix * L, const gsl_vector * wts)

This function initializes, or reinitializes, an existing ridge solver s to use the function and derivative fdf and the initial guess x. Here, the regularization matrix is set to L, which must have p columns but may have any number of rows.

Optionally, a weight vector wts can be given to perform a weighted nonlinear regression. Here, the weighting matrix is W = diag(w_1,w_2,...,w_n). The wts vector is referenced throughout the iteration so it should not be freed by the caller until the iteration terminates.

Function: void gsl_multifit_fsolver_free (gsl_multifit_fsolver * s)
Function: void gsl_multifit_fdfsolver_free (gsl_multifit_fdfsolver * s)
Function: void gsl_multifit_fdfridge_free (gsl_multifit_fdfridge * s)

These functions free all the memory associated with the solver s.

Function: const char * gsl_multifit_fsolver_name (const gsl_multifit_fsolver * s)
Function: const char * gsl_multifit_fdfsolver_name (const gsl_multifit_fdfsolver * s)
Function: const char * gsl_multifit_fdfridge_name (const gsl_multifit_fdfridge * s)

These functions return a pointer to the name of the solver. For example,

printf ("s is a '%s' solver\n", 
        gsl_multifit_fdfsolver_name (s));

would print something like s is a 'lmder' solver.


Next: , Previous: , Up: Nonlinear Least-Squares Fitting   [Index]