Class ode_bsimp_gsl (o2scl)¶
-
template<class
func_t
= ode_funct, classjac_func_t
= ode_jac_funct, classvec_t
= boost::numeric::ublas::vector<double>, classmat_t
= boost::numeric::ublas::matrix<double>>
classo2scl
::
ode_bsimp_gsl
¶ Bulirsch-Stoer implicit ODE stepper (GSL)
Implicit Bulirsch-Stoer method of Bader and Deuflhard ( [Bader83] ).
There is an example for the usage of this class in
examples/ex_stiff.cpp
documented in the ex_stiff_sect section.- Note
The variable
h_next
was defined in the original GSL version has been removed here, as it was unused by the stepper routine.- Note
At the moment, this class retains the GSL approach to handling non-integer return values in the user-specified derivative function. If the user-specified function returns
exc_efailed
, then the stepper attempts to decrease the stepsize to continue. If the user-specified function returns a non-zero value other thanexc_efailed
, or if the Jacobian evaluation returns any non-zero value, then the stepper aborts and returns the error value without calling the error handler. This behavior may change in the future.
- Idea for Future:
More detailed documentation about the algorithm
- Idea for Future:
Figure out if this should be a child of ode_step or astep. The function step_local() is actually its own ODE stepper and could be reimplemented as an object of type ode_step.
- Idea for Future:
I don’t like setting yerr to GSL_POSINF, there should be a better way to force an adaptive stepper which is calling this stepper to readjust the stepsize.
- Idea for Future:
The functions deuf_kchoice() and compute_weights() can be moved out of this header file.
- Idea for Future:
Rework internal arrays
- Idea for Future:
Rework linear solver to be amenable to using a sparse matrix solver
Workspace for extrapolation step
Workspace for the basic stepper
-
size_t
order
¶ Order of last step.
-
int
verbose
¶ Verbose parameter.
-
size_t
deuf_kchoice
(double eps2, size_t dimension)¶ Calculate a choice for the “order” of the method, using the Deuflhard criteria.
Used in the allocate() function.
-
int
poly_extrap
(ubmatrix &dloc, const double x[], const unsigned int i_step, const double x_i, const vec_t &y_i, vec_t &y_0, vec_t &y_0_err, ubarr &work)¶ Polynomial extrapolation.
Compute the step of index
i_step
using polynomial extrapolation to evaulate functions by fitting a polynomial to estimates(x_i,y_i)
and output the result toy_0
andy_0_err
.The index
i_step
begins with zero.
-
int
step_local
(const double t0, const double h_total, const unsigned int n_step, const ubarr &y, const ubarr &yp_loc, const vec_t &dfdt_loc, const mat_t &dfdy_loc, vec_t &y_out)¶ Basic implicit Bulirsch-Stoer step.
Divide the step
h_total
inton_step
smaller steps and do the Bader-Deuflhard semi-implicit iteration. This function starts att0
with function valuesy
, derivativesyp_loc
, and information from the Jacobian to compute the final valuey_out
.
-
int
allocate
(size_t n)¶ Allocate memory for a system of size
n
.
-
void
free
()¶ Free allocated memory.
-
ode_bsimp_gsl
()¶
-
~ode_bsimp_gsl
()¶
-
int
reset
()¶ Reset stepper.
-
int
step
(double x, double h, size_t n, vec_t &y, vec_t &dydx, vec_t &yout, vec_t &yerr, vec_t &dydx_out, func_t &derivs, jac_func_t &jac)¶ Perform an integration step.
Given initial value of the n-dimensional function in
y
and the derivative indydx
(which must generally be computed beforehand) at the pointx
, take a step of sizeh
giving the result inyout
, the uncertainty inyerr
, and the new derivative indydx_out
(atx+h
) using functionderivs
to calculate derivatives. Implementations which do not calculateyerr
and/ordydx_out
do not reference these variables so that a blankvec_t
can be given. All of the implementations allowyout=y
anddydx_out=dydx
if necessary.
Public Types
-
typedef boost::numeric::ublas::unbounded_array<double>
ubarr
¶
-
typedef boost::numeric::ublas::vector<double>
ubvector
¶
-
typedef boost::numeric::ublas::matrix<double>
ubmatrix
¶
Protected Attributes
-
size_t
dim
¶ Size of allocated vectors.
-
jac_func_t *
jfuncp
¶ Jacobian.
-
double
ex_wk
[sequence_max
]¶ Workspace for extrapolation.
(This state variable was named ‘x’ in GSL.)