Next: , Previous: Iteration of the Minimization Algorithm, Up: Nonlinear Least-Squares Fitting


38.5 Search Stopping Parameters

A minimization procedure should stop when one of the following conditions is true:

The handling of these conditions is under user control. The functions below allow the user to test the current estimate of the best-fit parameters in several standard ways.

— Function: int gsl_multifit_test_delta (const gsl_vector * dx, const gsl_vector * x, double epsabs, double epsrel)

This function tests for the convergence of the sequence by comparing the last step dx with the absolute error epsabs and relative error epsrel to the current position x. The test returns GSL_SUCCESS if the following condition is achieved,

          |dx_i| < epsabs + epsrel |x_i|

for each component of x and returns GSL_CONTINUE otherwise.

— Function: int gsl_multifit_test_gradient (const gsl_vector * g, double epsabs)

This function tests the residual gradient g against the absolute error bound epsabs. Mathematically, the gradient should be exactly zero at the minimum. The test returns GSL_SUCCESS if the following condition is achieved,

          \sum_i |g_i| < epsabs

and returns GSL_CONTINUE otherwise. This criterion is suitable for situations where the precise location of the minimum, x, is unimportant provided a value can be found where the gradient is small enough.

— Function: int gsl_multifit_gradient (const gsl_matrix * J, const gsl_vector * f, gsl_vector * g)

This function computes the gradient g of \Phi(x) = (1/2) ||F(x)||^2 from the Jacobian matrix J and the function values f, using the formula g = J^T f.