Public Member Functions | Public Attributes | Protected Member Functions | Protected Attributes | List of all members
o2scl::fit_nonlin< func_t, vec_t, mat_t > Class Template Reference

Non-linear least-squares fitting class (GSL) More...

#include <fit_nonlin.h>

Inheritance diagram for o2scl::fit_nonlin< func_t, vec_t, mat_t >:
o2scl::fit_nonlin_b< vec_t, mat_t > o2scl::fit_base< func_t, vec_t, mat_t >

Public Member Functions

virtual int print_iter (size_t nv, vec_t &x, vec_t &dx, int iter2, double l_epsabs, double l_epsrel)
 Print the progress in the current iteration.
 
void resize (size_t n, size_t p)
 Allocate memory with n data points and p parameters.
 
int set (size_t npar, vec_t &parms, func_t &fitfun)
 Set the initial values of the parameters and the fitting function to use for the next call to iterate()
 
int iterate ()
 Perform an iteration.
 
virtual int fit (size_t npar, vec_t &parms, mat_t &covar, double &chi2, func_t &fitfun)
 Fit the data specified in (xdat,ydat) to the function fitfun with the parameters in par. More...
 
virtual const char * type ()
 Return string denoting type ("fit_nonlin")
 
- Public Member Functions inherited from o2scl::fit_nonlin_b< vec_t, mat_t >
int test_delta_f (size_t nparm, vec_t &dx, vec_t &x, double l_epsabs, double l_epsrel)
 Test if converged.
 
int test_gradient_f (size_t nparm, vec_t &g, double l_epsabs)
 Test if converged.
 
- Public Member Functions inherited from o2scl::fit_base< func_t, vec_t, mat_t >
virtual int print_iter (size_t nv, vec_t &x, double y, int iter, double value=0.0, double limit=0.0)
 Print out iteration information. More...
 

Public Attributes

bool err_nonconv
 If true, call the error handler if fit() does not converge (default true)
 
vec_t dx_
 The last step taken in parameter space.
 
vec_t f_
 Desc.
 
bool test_gradient
 If true, test the gradient also (default false)
 
bool use_scaled
 Use the scaled routine if true (default true)
 
- Public Attributes inherited from o2scl::fit_nonlin_b< vec_t, mat_t >
double tol_rel_covar
 The relative tolerance for the computation of the covariance matrix (default 0)
 
- Public Attributes inherited from o2scl::fit_base< func_t, vec_t, mat_t >
size_t ntrial
 Maximum number of iterations (default 500)
 
double tol_abs
 Absolute tolerance (default 1.0e-4)
 
double tol_rel
 (default 1.0e-4)
 
int verbose
 An integer describing the verbosity of the output.
 
size_t n_dat
 The number of data points.
 
size_t n_par
 The number of parameters.
 

Protected Member Functions

void free ()
 Free allocated memory.
 
- Protected Member Functions inherited from o2scl::fit_nonlin_b< vec_t, mat_t >
double compute_actual_reduction (double fnorm0, double fnorm1)
 Desc.
 
size_t count_nsing (const size_t ncols, const mat_t &r2)
 Desc.
 
void compute_newton_direction (size_t n, const mat_t &r2, const permutation &perm2, const vec_t &qtf2, vec_t &x)
 Desc.
 
void compute_newton_bound (size_t nd, size_t np, const mat_t &r2, const vec_t &x, double dxnorm, const permutation &perm, const vec_t &diag, vec_t &w)
 Desc.
 
void compute_gradient_direction (size_t n, const mat_t &r, const permutation &p, const vec_t &qtf2, const vec_t &diag, vec_t &g)
 Desc.
 
void update_diag (size_t n, const mat_t &J, vec_t &diag2)
 Desc.
 
double scaled_enorm (const vec_t &d, size_t n, const vec_t &f)
 Euclidean norm of vector f of length n, scaled by vector d.
 
double compute_delta (vec_t &diag2, size_t n, const vec_t &x)
 Desc.
 
void compute_rptdx (const mat_t &r2, const permutation &p, size_t N, vec_t &dx, vec_t &rptdx2)
 Desc.
 
int qrsolv (size_t n, mat_t &r2, const permutation &p, const double lambda, const vec_t &diag2, const vec_t &qtb, vec_t &x, vec_t &sdiag2, vec_t &wa)
 Compute the solution to a least squares system. More...
 
void compute_newton_correction (size_t n, const mat_t &r2, const vec_t &sdiag2, const permutation &p, vec_t &x, double dxnorm, const vec_t &diag2, vec_t &w2)
 Desc.
 
void lmpar (mat_t &r2, const permutation &perm2, const vec_t &qtf2, const vec_t &diag2, double delta2, double *par_inout, vec_t &newton2, vec_t &gradient2, vec_t &sdiag2, vec_t &x, vec_t &w2, size_t nparm, size_t ndata)
 Determine Levenburg-Marquardt parameter.
 
void compute_trial_step (size_t N, vec_t &x, vec_t &dx, vec_t &trial)
 Compute trial step, $ \mathrm{trial}=\mathrm{x}+\mathrm{dx} $.
 
int compute_diag (size_t nparm, size_t ndata, const mat_t &J, vec_t &diag_vec)
 Compute the root of the sum of the squares of the columns of J. More...
 
int covariance (size_t m, size_t n, const mat_t &J, mat_t &covar, vec_t &norm, mat_t &r, vec_t &tau, permutation &perm, double epsrel)
 Compute the covarance matrix covar given the Jacobian J. More...
 

Protected Attributes

func_t * cff
 Function to fit.
 
vec_t x_trial
 Trial step.
 
vec_t f_trial
 Trial function value.
 
size_t iter
 Number of iterations.
 
double xnorm
 Desc.
 
double fnorm
 Desc.
 
double delta
 Desc.
 
double par
 Desc.
 
mat_t r
 Desc.
 
vec_t tau
 Desc.
 
vec_t diag
 Desc.
 
vec_t qtf
 Desc.
 
vec_t df
 Desc.
 
vec_t rptdx
 Desc.
 
vec_t newton
 Desc.
 
vec_t gradient
 Desc.
 
vec_t sdiag
 Desc.
 
vec_t w
 Desc.
 
vec_t work1
 Desc.
 
permutation perm
 Desc.
 
size_t ndata
 Number of data points.
 
size_t nparm
 Number of parameters.
 
vec_t g_
 Desc.
 
mat_t J_
 Desc.
 
vec_t * x_
 Desc.
 

Detailed Description

template<class func_t = gen_fit_funct<>, class vec_t = boost::numeric::ublas::vector<double>, class mat_t = boost::numeric::ublas::matrix<double>>
class o2scl::fit_nonlin< func_t, vec_t, mat_t >

The GSL-based fitting class using a Levenberg-Marquardt type algorithm. The algorithm stops when

\[ |dx_i| < \mathrm{tol\_abs}+\mathrm{tol\_rel}\times|x_i| \]

where $dx$ is the last step and $x$ is the current position. If test_gradient is true, then additionally fit() requires that

\[ \sum_i |g_i| < \mathrm{tol\_abs} \]

where $g_i$ is the $i$-th component of the gradient of the function $\Phi(x)$ where

\[ \Phi(x) = || F(x) ||^2 \]

Default template arguments

Todo:

Allow the user to specify the derivatives

Fix so that the user can specify automatic scaling of the fitting parameters, where the initial guess are used for scaling so that the fitting parameters are near unity.

Idea for Future:
Some of these member functions (like update_diag()) don't depend on member data and could be possibly be moved to a parent class?

Definition at line 921 of file fit_nonlin.h.

Member Function Documentation

◆ fit()

template<class func_t = gen_fit_funct<>, class vec_t = boost::numeric::ublas::vector<double>, class mat_t = boost::numeric::ublas::matrix<double>>
virtual int o2scl::fit_nonlin< func_t, vec_t, mat_t >::fit ( size_t  npar,
vec_t &  parms,
mat_t &  covar,
double &  chi2,
func_t &  fitfun 
)
inlinevirtual

The covariance matrix for the parameters is returned in covar and the value of $ \chi^2 $ is returned in chi2.

Implements o2scl::fit_base< func_t, vec_t, mat_t >.

Definition at line 1388 of file fit_nonlin.h.


The documentation for this class was generated from the following file:

Documentation generated with Doxygen. Provided under the GNU Free Documentation License (see License Information).