libode
Easy to compile, fast ODE integrators as C++ classes
Loading...
Searching...
No Matches
ode::OdeRKCK Class Reference

Explicit 5/4 pair, also with 3rd, 2nd, and 1st order embedded methods, from Cash & Karp. More...

#include <ode_rkck.h>

Inheritance diagram for ode::OdeRKCK:
ode::OdeEmbedded ode::OdeRK ode::OdeERK ode::OdeAdaptive ode::OdeBase

Public Member Functions

 OdeRKCK (unsigned long neq)
 constructs
 
- Public Member Functions inherited from ode::OdeEmbedded
 OdeEmbedded (unsigned long neq, bool need_jac, int lowerord)
 constructs
 
virtual ~OdeEmbedded ()
 destructs
 
double get_facsafe ()
 gets safety factor applied to time step selection
 
double get_facmin ()
 gets minimum allowable fraction change in time step (a number <1)
 
double get_facmax ()
 gets maximum allowable fraction change in time step (a number >1)
 
void set_facsafe (double facsafe)
 sets safety factor applied to time step selection
 
void set_facmin (double facmin)
 sets minimum allowable fraction change in time step (a number <1)
 
void set_facmax (double facmax)
 sets maximum allowable fraction change in time step (a number >1)
 
- Public Member Functions inherited from ode::OdeAdaptive
 OdeAdaptive (unsigned long neq, bool need_jac)
 constructs
 
virtual ~OdeAdaptive ()
 destructs
 
long unsigned get_nrej ()
 gets the count of rejected steps
 
double get_abstol ()
 gets the absolute error tolerance
 
double get_reltol ()
 gets the relative error tolerance
 
double get_dtmax ()
 gets the maximum allowable time step
 
void set_abstol (double tol)
 sets the absolute error tolerance
 
void set_reltol (double tol)
 sets the relative error tolerance
 
void set_tol (double tol)
 sets the absolute and relative error tolerance to the same value
 
void set_dtmax (double dtmax)
 sets the maximum allowable time step
 
void solve_adaptive (double tint, double dt0, bool extras=true)
 integrates for a specified duration of independent variable without output
 
void solve_adaptive (double tint, double dt0, const char *dirout, int inter)
 lots of output, solves and stores every "inter" point along the way
 
void solve_adaptive (double tint, double dt0, unsigned long nsnap, const char *dirout)
 solves and writes evenly spaced snapshots
 
void solve_adaptive (double dt0, double *tsnap, unsigned long nsnap, const char *dirout)
 solves and writes snapshots at times specified in the tsnap array
 
- Public Member Functions inherited from ode::OdeBase
 OdeBase (unsigned long neq, bool need_jac)
 constructs
 
virtual ~OdeBase ()
 destructs
 
const char * get_name ()
 gets the name of the ODE system
 
const char * get_method ()
 gets the name of the solver/method
 
const char * get_dirout ()
 gets output directory string if one has been set
 
bool get_quiet ()
 gets the boolean determining if updates are printed during solves
 
bool get_silent_snap ()
 gets whether to skip writing the solution vector to file when snapping
 
unsigned long get_neq ()
 gets the size of the ODE system
 
double get_t ()
 gets the current value of the independent variable
 
double get_dt ()
 gets the most recent or current time step size
 
double * get_sol ()
 gets a pointer to the whole solution array
 
double get_sol (unsigned long i)
 gets an element of the solution array
 
long unsigned get_nstep ()
 gets the total number of steps taken
 
long unsigned get_neval ()
 gets the total number of ODE system evaluation
 
long unsigned get_icheck ()
 gets the number of steps after which the solution is checked for integrity
 
long unsigned get_nJac ()
 gets the total number of Jacobian evaluations performed
 
void set_t (double t)
 sets the "time," or independent variable used to track progress
 
void set_sol (unsigned long i, double x)
 sets an element of the solution array
 
void set_sol (double *sol)
 copies an array into the solution array
 
void set_name (std::string name)
 sets the name of the ODE system
 
void set_name (const char *name)
 sets the name of the ODE system
 
void set_quiet (bool quiet)
 sets the boolean determining if updates are printed during solves
 
void set_silent_snap (bool silent_snap)
 sets whether to skip writing the solution vector to file when snapping
 
void set_icheck (unsigned long icheck)
 sets the number of steps after which the solution is checked for integrity
 
void step (double dt, bool extra=true)
 increments the step counter and the time, checks the solution integrity if needed, stores the time step in the object, and executes after_step() if extra is true
 
void solve_fixed (double tint, double dt, bool extras=true)
 integrates for a specified duration of independent variable without output
 
void solve_fixed (double tint, double dt, const char *dirout, int inter=1)
 lots of output, solves and stores every "inter" point along the way
 
void solve_fixed (double tint, double dt, unsigned long nsnap, const char *dirout)
 solves and writes evenly spaced snapshots
 
void solve_fixed (double dt, double *tsnap, unsigned long nsnap, const char *dirout)
 solves and writes snapshots at times specified in the tsnap array
 
void reset (double t, double *sol)
 reset to a specified time and initial condition array
 

Additional Inherited Members

- Protected Member Functions inherited from ode::OdeEmbedded
double error (double abstol, double reltol)
 calculates error estimate with lower and higher order solutions
 
double facopt (double err)
 calculates factor for "optimal" next time step
 
virtual void adapt (double abstol, double reltol)
 does the calculations to determine isrej and dtopt
 
virtual bool is_rejected ()
 simply returns isrej
 
virtual double dt_adapt ()
 simply returns dtopt
 
- Protected Member Functions inherited from ode::OdeAdaptive
void solve_adaptive_ (double tint, double dt0, bool extra=true)
 integrates without output or any counters, trackers, extra functions...
 
bool solve_done_adaptive (double tend)
 determines whether an adaptive solve is finished
 
bool step_adaptive_ (double dt, bool extra=true)
 executes a single time and calls all necessary adapting functions
 
double dt_adapt_ (double tend)
 wrapper around dt_adapt() to perform additional checks
 
- Protected Member Functions inherited from ode::OdeBase
void solve_fixed_ (double tint, double dt, bool extra=true)
 integrates without output or any counters, trackers, extra functions...
 
virtual void ode_fun (double *solin, double *fout)=0
 evaluates the system of ODEs in autonomous form and must be defined by a derived class
 
virtual void ode_jac (double *solin, double **Jout)
 evaluates the system's Jacobian matrix, also in autonomous form, and can either be defined in a derived class or left to numerical approximation
 
void ode_fun_ (double *solin, double *fout)
 wrapper, calls ode_fun() and increases the neval counter by one
 
void ode_jac_ (double *solin, double **Jout)
 wrapper, calls ode_jac() and increments nJac;
 
virtual void before_solve ()
 does any extra stuff before starting a solve
 
virtual void after_step (double t)
 does any extra stuff after each step
 
virtual void after_capture (double t)
 does any extra stuff only when a step is captured
 
virtual void after_snap (std::string dirout, long isnap, double t)
 does any extra stuff after each snap
 
virtual void after_solve ()
 does any extra stuff after completing a solve
 
void snap (std::string dirout, long isnap, double tsnap)
 writes the current value of the solution to a binary file
 
bool solve_done (double dt, double tend)
 checks if the solution is within a single time step of the end point
 
void check_sol_integrity ()
 checks solution for nans and infs, exiting the program if they're found
 
void check_pre_solve (double tint, double dt)
 checks that a solve can be performed with given tend and dt values
 
void check_pre_snaps (double dt, double *tsnap, unsigned long nsnap)
 checks that snap times are monotonically increasing and > current time
 
- Protected Attributes inherited from ode::OdeEmbedded
double facsafe_
 safety factor applied to time step selection
 
double facmin_
 minimum allowable fraction change in time step
 
double facmax_
 maximum allowable fraction change in time step
 
double * solemb_
 embedded solution array
 
- Protected Attributes inherited from ode::OdeAdaptive
long unsigned nrej_
 counter for rejected steps
 
double abstol_
 absolute error tolerance
 
double reltol_
 absolute error tolerance
 
double dtmax_
 maximum allowable time step
 
- Protected Attributes inherited from ode::OdeBase
std::string name_
 the "name" of the system, which is used for output
 
std::string method_
 the "name" of the solver/method, as in "Euler" or "RK4"
 
std::string dirout_
 output directory if one is being used by a solver
 
bool quiet_
 whether stuff should be printed during a solve
 
bool silent_snap_
 whether to skip writing the solution vector to file when snapping but still execute after_snap()
 
unsigned long neq_
 number of equations in the system of ODEs
 
double t_
 time, initialized to zero
 
double dt_
 time step is stored and updated during solves
 
double * sol_
 array for the solution, changing over time
 
long unsigned nstep_
 number of time steps
 
long unsigned neval_
 function evaluation counter, must be incremented in step() when defined
 
long unsigned icheck_
 interval of steps after which to check for nans and infs (zero to ignore)
 
double ** Jac_
 storage for the ODE system's Jacobian matrix, only allocated for the methods that need it
 
long unsigned nJac_
 counter for jacobian evaluations
 
double absjacdel_
 absolute adjustment fraction for numerical Jacobian, if needed
 
double reljacdel_
 relative adjustment fraction for numerical Jacobian, if needed
 

Detailed Description

Explicit 5/4 pair, also with 3rd, 2nd, and 1st order embedded methods, from Cash & Karp.

This class implements a 5th order method developed by Cash and Karp (reference below) which includes embeddes solutions of order 1, 2, 3, and 4. It's a complete family of solutions up to order 5. The lower order solutions can be used to vary the solution's order and control the step size in sophisticated ways, canceling a step early if the error estimate is large. The fact that function evaluations are made fairly evenly through the interval of a time step is also advantageous for surveying potential trouble and controlling the step.

So far, the sophisticated time stepping is not implemented, so this solver simply uses the 4th and 5th order solutions for adaptation.

+ J. R. Cash, A. H. Karp. "A variable order Runge-Kutta method for initial value problems with rapidly varying right-hand sides", ACM Transactions on Mathematical Software 16: 201-222, 1990. doi:10.1145/79505.79507
+ https://en.wikipedia.org/wiki/Cash%E2%80%93Karp_method

Definition at line 21 of file ode_rkck.h.

Constructor & Destructor Documentation

◆ OdeRKCK()

ode::OdeRKCK::OdeRKCK ( unsigned long neq)

constructs

Parameters
[in]neqsize of ODE sytem

Definition at line 7 of file ode_rkck.cc.


The documentation for this class was generated from the following files: