#include <gsl_mroot_hybrids.h>
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. 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 mroot::tolf 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 Multidimensional solver example .
Definition at line 352 of file gsl_mroot_hybrids.h.
Public Member Functions | |
virtual int | set_jacobian (jacobian< param_t, 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, param_t &pa, 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, param_t &pa, func_t &ufunc) |
Solve ufunc using xx as an initial guess, returning xx . | |
int | set (size_t nn, vec_t &ax, func_t &ufunc, param_t &pa) |
Set the function, the parameters, and the initial guess. | |
int | set_de (size_t nn, vec_t &ax, func_t &ufunc, jfunc_t &dfunc, param_t &pa) |
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< param_t, 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 | |
void | compute_Rg (size_t N, const gsl_matrix *r, const gsl_vector *gradient, vec_t &Rg) |
Desc. | |
void | compute_wv (size_t n, const gsl_vector *qtdf, const gsl_vector *rdx, const vec_t &dxx, const gsl_vector *diag, double pnorm, gsl_vector *w, gsl_vector *v) |
Desc. | |
void | compute_rdx (size_t N, const gsl_matrix *r, const vec_t &dxx, gsl_vector *rdx) |
Desc. | |
double | scaled_enorm_tvec (size_t n, const gsl_vector *d, const vec_t &ff) |
Desc. | |
double | compute_delta (size_t n, gsl_vector *diag, vec_t &xx) |
Desc. | |
double | enorm_tvec (const vec_t &ff) |
Desc. | |
int | compute_trial_step_tvec (size_t N, vec_t &xl, vec_t &dxl, vec_t &xx_trial) |
Desc. | |
int | compute_df_tvec (size_t n, const vec_t &ff_trial, const vec_t &fl, gsl_vector *dfl) |
Desc. | |
void | compute_diag_tvec (size_t n, const mat_t &J, gsl_vector *diag) |
Desc. | |
void | compute_qtf_tvec (size_t N, const gsl_matrix *q, const vec_t &ff, gsl_vector *qtf) |
Desc. | |
void | update_diag_tvec (size_t n, const mat_t &J, gsl_vector *diag) |
Desc. | |
void | scaled_addition_tvec (size_t N, double alpha, gsl_vector *newton, double beta, gsl_vector *gradient, vec_t &pp) |
Desc. | |
int | dogleg (size_t n, const gsl_matrix *r, const gsl_vector *qtf, const gsl_vector *diag, double delta, gsl_vector *newton, gsl_vector *gradient, vec_t &p) |
Take a dogleg step. | |
int | solve_set (size_t nn, vec_t &xx, param_t &pa, func_t &ufunc) |
Finish the solution after set() or set_de() has been called. | |
Protected Attributes | |
jfunc_t * | jac |
The user-specified Jacobian. | |
jacobian< param_t, func_t, vec_t, mat_t > * | ajac |
The automatic Jacobian. | |
alloc_t | ao |
Memory allocator for objects of type alloc_vec_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. | |
o2scl_hybrid_state_t< vec_t, alloc_vec_t, alloc_t, mat_t, alloc_mat_t, mat_alloc_t > | state |
The solver state. | |
param_t * | params |
The function parameters. | |
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. |
int iterate | ( | ) | [inline] |
Perform an iteration.
At the end of the iteration, the current value of the solution is stored in x.
Definition at line 746 of file gsl_mroot_hybrids.h.
virtual int msolve_de | ( | size_t | nn, | |
vec_t & | xx, | |||
param_t & | pa, | |||
func_t & | ufunc, | |||
jfunc_t & | dfunc | |||
) | [inline, virtual] |
Solve func
with derivatives dfunc
using x
as an initial guess, returning x
.
Reimplemented from mroot.
Definition at line 1014 of file gsl_mroot_hybrids.h.
bool shrink_step |
If true, iterate() will shrink the step-size automatically if the function returns a non-zero value (default true).
The original GSL behavior can be obtained by setting this to false
.
Definition at line 711 of file gsl_mroot_hybrids.h.
Documentation generated with Doxygen and provided under the GNU Free Documentation License. See License Information for details.
Project hosting provided by
,
O2scl Sourceforge Project Page