Object-oriented Scientific Computing Library: Version 0.910
Public Member Functions | Data Fields | Protected Member Functions | Protected Attributes
gsl_mroot_hybrids< func_t, vec_t, alloc_vec_t, alloc_t, mat_t, alloc_mat_t, mat_alloc_t, jfunc_t > Class Template Reference

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

#include <gsl_mroot_hybrids.h>

Inheritance diagram for gsl_mroot_hybrids< func_t, vec_t, alloc_vec_t, alloc_t, mat_t, alloc_mat_t, mat_alloc_t, jfunc_t >:
mroot< func_t, vec_t, jfunc_t > hybrids_base

Detailed Description

template<class func_t = mm_funct<>, class vec_t = ovector_base, class alloc_vec_t = ovector, class alloc_t = ovector_alloc, class mat_t = omatrix_base, class alloc_mat_t = omatrix, class mat_alloc_t = omatrix_alloc, class jfunc_t = jac_funct<vec_t,mat_t>>
class gsl_mroot_hybrids< func_t, vec_t, alloc_vec_t, alloc_t, mat_t, alloc_mat_t, mat_alloc_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 free() and 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.

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. To return to the default GSL behavior, set shrink_step to false.

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

There is an example for the usage of the multidimensional solver classes given in examples/ex_mroot.cpp, see the Multi-dimensional solver .

Todo:
Are all the set_all() statements really necessary? Do they need to be executed even if memory hasn't been recently allocated?
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:
It might be good to make sure that we're directly testing to make sure that GSL and O2SCL are more or less identical.
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 485 of file gsl_mroot_hybrids.h.

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.
int allocate (size_t n)
 Allocate the memory.
int free ()
 Free the allocated memory.
virtual const char * type ()
 Return the type,"gsl_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.
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.

Data Fields

bool shrink_step
 If true, iterate() will shrink the step-size automatically if the function returns a non-zero value (default true)
bool int_scaling
 If true, use the internal scaling method (default true)
simple_jacobian< func_t, vec_t,
mat_t, alloc_vec_t, alloc_t > 
def_jac
 Default automatic Jacobian object.
alloc_vec_t f
 The value of the function at the present iteration.
alloc_vec_t x
 The present solution.

Protected Member Functions

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.
alloc_mat_t J
 Jacobian.
omatrix q
 Q matrix from QR decomposition.
omatrix r
 R matrix from QR decomposition.
ovector tau
 The tau vector from QR decomposition.
ovector diag
 The diagonal elements.
ovector qtf
 The value of $ Q^T f $.
ovector newton
 The Newton direction.
ovector gradient
 The gradient direction.
ovector df
 The change in the function value.
ovector qtdf
 The value of $ Q^T \cdot \mathrm{df} $.
ovector rdx
 The value of $ R \cdot \mathrm{dx} $.
ovector w
 The value of $ w=(Q^T df-R dx)/|dx| $.
ovector 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.
alloc_t ao
 Memory allocator for objects of type alloc_vec_t.
mat_alloc_t am
 Memory allocator for objects of type alloc_mat_t.
alloc_vec_t dx
 The value of the derivative.
alloc_vec_t x_trial
 Trial root.
alloc_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.

Member Function Documentation

template<class func_t = mm_funct<>, class vec_t = ovector_base, class alloc_vec_t = ovector, class alloc_t = ovector_alloc, class mat_t = omatrix_base, class alloc_mat_t = omatrix, class mat_alloc_t = omatrix_alloc, class jfunc_t = jac_funct<vec_t,mat_t>>
int gsl_mroot_hybrids< func_t, vec_t, alloc_vec_t, alloc_t, mat_t, alloc_mat_t, mat_alloc_t, jfunc_t >::iterate ( ) [inline]

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

Definition at line 675 of file gsl_mroot_hybrids.h.

template<class func_t = mm_funct<>, class vec_t = ovector_base, class alloc_vec_t = ovector, class alloc_t = ovector_alloc, class mat_t = omatrix_base, class alloc_mat_t = omatrix, class mat_alloc_t = omatrix_alloc, class jfunc_t = jac_funct<vec_t,mat_t>>
virtual int gsl_mroot_hybrids< func_t, vec_t, alloc_vec_t, alloc_t, mat_t, alloc_mat_t, mat_alloc_t, jfunc_t >::msolve_de ( size_t  nn,
vec_t &  xx,
func_t &  ufunc,
jfunc_t &  dfunc 
) [inline, virtual]

Reimplemented from mroot< func_t, vec_t, jfunc_t >.

Definition at line 972 of file gsl_mroot_hybrids.h.


Field Documentation

template<class func_t = mm_funct<>, class vec_t = ovector_base, class alloc_vec_t = ovector, class alloc_t = ovector_alloc, class mat_t = omatrix_base, class alloc_mat_t = omatrix, class mat_alloc_t = omatrix_alloc, class jfunc_t = jac_funct<vec_t,mat_t>>
bool gsl_mroot_hybrids< func_t, vec_t, alloc_vec_t, alloc_t, mat_t, alloc_mat_t, mat_alloc_t, jfunc_t >::shrink_step

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

Definition at line 644 of file gsl_mroot_hybrids.h.

template<class func_t = mm_funct<>, class vec_t = ovector_base, class alloc_vec_t = ovector, class alloc_t = ovector_alloc, class mat_t = omatrix_base, class alloc_mat_t = omatrix, class mat_alloc_t = omatrix_alloc, class jfunc_t = jac_funct<vec_t,mat_t>>
alloc_vec_t gsl_mroot_hybrids< func_t, vec_t, alloc_vec_t, alloc_t, mat_t, alloc_mat_t, mat_alloc_t, jfunc_t >::f

Definition at line 665 of file gsl_mroot_hybrids.h.


The documentation for this class was generated from the following file:
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines

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

Get Object-oriented Scientific Computing
Lib at SourceForge.net. Fast, secure and Free Open Source software
downloads.