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

Multidimensional root-finding algorithm using Powell's Hybrid method (GSL) More...

#include <mroot_hybrids.h>

Inheritance diagram for o2scl::mroot_hybrids< func_t, vec_t, mat_t, jfunc_t >:
o2scl::mroot< func_t, vec_t, jfunc_t > o2scl::mroot_hybrids_base o2scl::mroot_hybrids_arma_qr_econ< func_t, vec_t, mat_t, jfunc_t > o2scl::mroot_hybrids_eigen< func_t, vec_t, mat_t, jfunc_t >

Public Member Functions

virtual int set_jacobian (jacobian< func_t, vec_t, mat_t > &j)
 Set the automatic Jacobian object.
 
int iterate ()
 Perform an iteration. More...
 
void allocate (size_t n)
 Allocate the memory.
 
virtual const char * type ()
 Return the type,"mroot_hybrids".
 
virtual int msolve_de (size_t nn, vec_t &xx, func_t &ufunc, jfunc_t &dfunc)
 Solve func with derivatives dfunc using x as an initial guess, returning x. More...
 
virtual int msolve (size_t nn, vec_t &xx, func_t &ufunc)
 Solve ufunc using xx as an initial guess, returning xx.
 
int set (size_t nn, vec_t &ax, func_t &ufunc)
 Set the function, the parameters, and the initial guess.
 
int set_de (size_t nn, vec_t &ax, func_t &ufunc, jfunc_t &dfunc)
 Set the function, the Jacobian, the parameters, and the initial guess.
 
- Public Member Functions inherited from o2scl::mroot< func_t, vec_t, jfunc_t >
template<class vec2_t , class vec3_t >
int print_iter (size_t n, const vec2_t &x, const vec3_t &y, int iter, double value=0.0, double limit=0.0, std::string comment="")
 Print out iteration information. More...
 

Public Attributes

bool shrink_step
 If true, iterate() will shrink the step-size automatically if the function returns a non-zero value (default true) More...
 
bool extra_finite_check
 If true, double check that the input function values are finite (default true)
 
bool int_scaling
 If true, use the internal scaling method (default true)
 
jacobian_gsl< func_t, vec_t, mat_t > def_jac
 Default automatic Jacobian object.
 
vec_t f
 The value of the function at the present iteration. More...
 
vec_t x
 The present solution.
 
- Public Attributes inherited from o2scl::mroot< func_t, vec_t, jfunc_t >
double tol_rel
 The maximum value of the functions for success (default 1.0e-8)
 
double tol_abs
 The minimum allowable stepsize (default 1.0e-12)
 
int verbose
 Output control (default 0)
 
int ntrial
 Maximum number of iterations (default 100)
 
int last_ntrial
 The number of iterations for in the most recent minimization.
 
bool err_nonconv
 If true, call the error handler if msolve() or msolve_de() does not converge (default true)
 

Protected Member Functions

virtual int solve_set (size_t nn, vec_t &xx, func_t &ufunc)
 Finish the solution after set() or set_de() has been called.
 

Protected Attributes

int iter
 Number of iterations.
 
size_t ncfail
 Compute the number of failures.
 
size_t ncsuc
 Compute the number of successes.
 
size_t nslow1
 The number of times the actual reduction is less than 0.001.
 
size_t nslow2
 The number of times the actual reduction is less than 0.1.
 
double fnorm
 The norm of the current function value.
 
double delta
 The limit of the Nuclidean norm.
 
mat_t J
 Jacobian.
 
mat_t q
 Q matrix from QR decomposition.
 
mat_t r
 R matrix from QR decomposition.
 
ubvector diag
 The diagonal elements.
 
ubvector qtf
 The value of $ Q^T f $.
 
ubvector newton
 The Newton direction.
 
ubvector gradient
 The gradient direction.
 
ubvector df
 The change in the function value.
 
ubvector qtdf
 The value of $ Q^T \cdot \mathrm{df} $.
 
ubvector rdx
 The value of $ R \cdot \mathrm{dx} $.
 
ubvector w
 The value of $ w=(Q^T df-R dx)/|dx| $.
 
ubvector v
 The value of $ v=D^2 dx/|dx| $.
 
jfunc_t * jac
 The user-specified Jacobian.
 
jacobian< func_t, vec_t, mat_t > * ajac
 The automatic Jacobian.
 
vec_t dx
 The value of the derivative.
 
vec_t x_trial
 Trial root.
 
vec_t f_trial
 Trial function value.
 
size_t dim
 The number of equations and unknowns.
 
bool jac_given
 True if the jacobian has been given.
 
func_t * fnewp
 The user-specified function.
 
bool set_called
 True if "set" has been called.
 

Private Member Functions

 mroot_hybrids (const mroot_hybrids< func_t, vec_t, mat_t, jfunc_t > &)
 
mroot_hybrids< func_t, vec_t, mat_t, jfunc_t > & operator= (const mroot_hybrids< func_t, vec_t, mat_t, jfunc_t > &)
 
- Private Member Functions inherited from o2scl::mroot_hybrids_base
double compute_actual_reduction (double fnorm0, double fnorm1)
 Compute the actual reduction.
 
double compute_predicted_reduction (double fnorm0, double fnorm1)
 Compute the predicted reduction phi1p=|Q^T f + R dx|.
 
template<class vec2_t , class mat_t >
void compute_Rg (size_t N, const mat_t &r2, const ubvector &gradient2, vec2_t &Rg)
 Compute $ R \cdot g $ where g is the gradient.
 
template<class vec2_t >
void compute_wv (size_t n, const ubvector &qtdf2, const ubvector &rdx2, const vec2_t &dx2, const ubvector &diag2, double pnorm, ubvector &w2, ubvector &v2)
 Compute w and v.
 
template<class vec2_t , class mat_t >
void compute_rdx (size_t N, const mat_t &r2, const vec2_t &dx2, ubvector &rdx2)
 Compute $ R \cdot \mathrm{dx} $.
 
template<class vec2_t , class vec3_t >
double scaled_enorm (size_t n, const vec2_t &d, const vec3_t &ff)
 Compute the norm of the vector $ \vec{v} $ defined by $ v_i = d_i ff_i $.
 
template<class vec2_t >
double compute_delta (size_t n, ubvector &diag2, vec2_t &x2)
 Compute delta.
 
template<class vec2_t >
double enorm (size_t N, const vec2_t &ff)
 Compute the Euclidean norm of ff. More...
 
double enorm_sum (size_t n, const ubvector &a, const ubvector &b)
 Compute the Euclidean norm of the sum of a and b.
 
template<class vec2_t >
int compute_trial_step (size_t N, vec2_t &xl, vec2_t &dxl, vec2_t &xx_trial)
 Compute a trial step and put the result in xx_trial. More...
 
template<class vec2_t >
int compute_df (size_t n, const vec2_t &ff_trial, const vec2_t &fl, ubvector &dfl)
 Compute the change in the function value. More...
 
template<class mat2_t >
void compute_diag (size_t n, const mat2_t &J2, ubvector &diag2)
 Compute diag, the norm of the columns of the Jacobian.
 
template<class vec2_t , class vec3_t , class vec4_t >
void compute_qtf (size_t N, const vec2_t &q2, const vec3_t &f2, vec4_t &qtf2)
 Compute $ Q^{T} f $. More...
 
template<class mat2_t >
void update_diag (size_t n, const mat2_t &J2, ubvector &diag2)
 Update diag.
 
template<class vec2_t >
void scaled_addition (size_t N, double alpha, ubvector &newton2, double beta, ubvector &gradient2, vec2_t &pp)
 Form appropriate convex combination of the Gauss-Newton direction and the scaled gradient direction. More...
 
template<class mat_t >
int newton_direction (const size_t N, const mat_t &r2, const ubvector &qtf2, ubvector &p)
 Compute the Gauss-Newton direction.
 
template<class mat_t >
void gradient_direction (const size_t M, const size_t N, const mat_t &r2, const ubvector &qtf2, const ubvector &diag2, ubvector &g)
 Compute the gradient direction.
 
void minimum_step (const size_t N, double gnorm, const ubvector &diag2, ubvector &g)
 Compute the point at which the gradient is minimized.
 
template<class vec2_t , class mat_t >
int dogleg (size_t n, const mat_t &r2, const ubvector &qtf2, const ubvector &diag2, double delta2, ubvector &newton2, ubvector &gradient2, vec2_t &p)
 Take a dogleg step. More...
 

Additional Inherited Members

- Private Types inherited from o2scl::mroot_hybrids_base
typedef boost::numeric::ublas::vector< double > ubvector
 
typedef boost::numeric::ublas::matrix< double > ubmatrix
 

Detailed Description

template<class func_t = mm_funct11, class vec_t = boost::numeric::ublas::vector<double>, class mat_t = boost::numeric::ublas::matrix<double>, class jfunc_t = jac_funct11>
class o2scl::mroot_hybrids< func_t, vec_t, mat_t, jfunc_t >

This is a recasted version of the GSL routines which use a modified version of Powell's Hybrid method as implemented in the HYBRJ algorithm in MINPACK (Garbow80). Both the scaled and unscaled options are available by setting int_scaling (the scaled version is the default). If derivatives are not provided, they will be computed automatically. This class provides the GSL-like interface using allocate(), set() (or set_de() in case where derivatives are available), iterate(), and the higher-level interfaces, msolve() and msolve_de(), which perform the solution and the memory allocation automatically. Some additional checking is performed in case the user calls the functions out of order (i.e. set() without allocate()).

The functions msolve() and msolve_de() use the condition $ \sum_i |f_i|<$ mroot::tol_rel to determine if the solver has succeeded.

See the Multi-dimensional solvers section of the User's guide for general information about O2scl solvers. There is an example for the usage of the multidimensional solver classes given in examples/ex_mroot.cpp, see the Multi-dimensional solver example .

Note
The set() and set_de() functions store a pointer to the function object and the user must ensure that the object is still valid for a later call to iterate().

The original GSL algorithm has been modified to shrink the stepsize if a proposed step causes the function to return a non-zero value. This allows the routine to automatically try to avoid regions where the function is not defined. The algorithm is also modified to check that it is not sending non-finite values to the user-specified function. To return to the default GSL behavior, set shrink_step and extra_finite_check to false.

The default method for numerically computing the Jacobian is from jacobian_gsl. This default is identical to the GSL approach, except that the default value of jacobian_gsl::epsmin is non-zero. See jacobian_gsl for more details.

By default convergence failures result in calling the exception handler, but this can be turned off by setting mroot::err_nonconv to false. If mroot::err_nonconv is false, then the functions iterate(), msolve() and msolve_de() will return a non-zero value if convergence fails. Note that the default Jacobian object, def_jac also has a data member jacobian_gsl::err_nonconv which separately handles the case where the one row of the Jacobian is all zero.

Idea for Future:
Is all the setting of vectors and matrices to zero really necessary? Do they need to be executed even if memory hasn't been recently allocated?
Idea for Future:
Convert more ubvectors to vec_t.
Idea for Future:
Some more of the element-wise vector manipulation could be converted to BLAS routines.
Idea for Future:
It's kind of strange that set() sets jac_given to false and set_de() has to reset it to true. Can this be simplified?
Idea for Future:
Many of these minpack functions could be put in their own "minpack_tools" class, or possibly moved to be linear algebra routines instead.
Idea for Future:
There are still some numbers in here which the user could have control over, for example, the nslow2 threshold which indicates failure.

Definition at line 521 of file mroot_hybrids.h.

Member Function Documentation

◆ iterate()

template<class func_t = mm_funct11, class vec_t = boost::numeric::ublas::vector<double>, class mat_t = boost::numeric::ublas::matrix<double>, class jfunc_t = jac_funct11>
int o2scl::mroot_hybrids< func_t, vec_t, mat_t, jfunc_t >::iterate ( )
inline

At the end of the iteration, the current value of the solution is stored in x.

Definition at line 704 of file mroot_hybrids.h.

◆ msolve_de()

template<class func_t = mm_funct11, class vec_t = boost::numeric::ublas::vector<double>, class mat_t = boost::numeric::ublas::matrix<double>, class jfunc_t = jac_funct11>
virtual int o2scl::mroot_hybrids< func_t, vec_t, mat_t, jfunc_t >::msolve_de ( size_t  nn,
vec_t &  xx,
func_t &  ufunc,
jfunc_t &  dfunc 
)
inlinevirtual

Reimplemented from o2scl::mroot< func_t, vec_t, jfunc_t >.

Definition at line 983 of file mroot_hybrids.h.

Member Data Documentation

◆ f

template<class func_t = mm_funct11, class vec_t = boost::numeric::ublas::vector<double>, class mat_t = boost::numeric::ublas::matrix<double>, class jfunc_t = jac_funct11>
vec_t o2scl::mroot_hybrids< func_t, vec_t, mat_t, jfunc_t >::f

Definition at line 694 of file mroot_hybrids.h.

◆ shrink_step

template<class func_t = mm_funct11, class vec_t = boost::numeric::ublas::vector<double>, class mat_t = boost::numeric::ublas::matrix<double>, class jfunc_t = jac_funct11>
bool o2scl::mroot_hybrids< func_t, vec_t, mat_t, jfunc_t >::shrink_step

The original GSL behavior can be obtained by setting this to false.

Definition at line 668 of file mroot_hybrids.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).