#include <gsl_mroot_hybrids.h>
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 Multi-dimensional 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) |
Compute ![]() g is the gradient . | |
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) |
Compute w and v . | |
void | compute_rdx (size_t N, const gsl_matrix *r, const vec_t &dxx, gsl_vector *rdx) |
Compute ![]() | |
double | scaled_enorm_tvec (size_t n, const gsl_vector *d, const vec_t &ff) |
Compute the scaled Euclidean norm. | |
double | compute_delta (size_t n, gsl_vector *diag, vec_t &xx) |
Compute delta. | |
double | enorm_tvec (const vec_t &ff) |
Compute the Euclidiean norm. | |
int | compute_trial_step_tvec (size_t N, vec_t &xl, vec_t &dxl, vec_t &xx_trial) |
Compute a trial step and put the result in xx_trial . | |
int | compute_df_tvec (size_t n, const vec_t &ff_trial, const vec_t &fl, gsl_vector *dfl) |
Compute the change in the function value. | |
void | compute_diag_tvec (size_t n, const mat_t &J, gsl_vector *diag) |
Compute diag , the norm of the columns of the Jacobian. | |
void | compute_qtf_tvec (size_t N, const gsl_matrix *q, const vec_t &ff, gsl_vector *qtf) |
Compute ![]() | |
void | update_diag_tvec (size_t n, const mat_t &J, gsl_vector *diag) |
Update diag . | |
void | scaled_addition_tvec (size_t N, double alpha, gsl_vector *newton, double beta, gsl_vector *gradient, vec_t &pp) |
Form appropriate convex combination of the Gauss-Newton direction and the scaled gradient direction. | |
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. | |
mroot_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 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 | |||
) | [inline, protected] |
Given the QR decomposition of an n
by n
matrix "A", a diagonal matrix diag
, a vector b
, and a positive number delta
, this function determines the convex combination x
of the Gauss-Newton and scaled gradient directions which minimizes in the least squares sense, subject to the restriction that the Euclidean norm of
is smaller than the value of
delta
.
Definition at line 555 of file gsl_mroot_hybrids.h.
int iterate | ( | ) | [inline] |
At the end of the iteration, the current value of the solution is stored in x.
Definition at line 783 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] |
void scaled_addition_tvec | ( | size_t | N, | |
double | alpha, | |||
gsl_vector * | newton, | |||
double | beta, | |||
gsl_vector * | gradient, | |||
vec_t & | pp | |||
) | [inline, protected] |
Using the Gauss-Newton direction given in newton
(a vector of size N), and the gradient direction in gradient
(also a vector of size N), this computes
Definition at line 533 of file gsl_mroot_hybrids.h.
bool shrink_step |
The original GSL behavior can be obtained by setting this to false
.
Definition at line 749 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