Public Member Functions | List of all members
RKDP45 Class Reference

4/5th Order Dormand–Prince Runge-Kutta ODE solver. More...

Inheritance diagram for RKDP45:
Inheritance graph
[legend]
Collaboration diagram for RKDP45:
Collaboration graph
[legend]

Public Member Functions

 TypeName ("RKDP45")
 Runtime type information. More...
 
 RKDP45 (const ODESystem &ode, const dictionary &dict)
 Construct from ODE. More...
 
scalar solve (const scalar x0, const scalarField &y0, const scalarField &dydx0, const scalar dx, scalarField &y) const
 Solve a single step dx and return the error. More...
 
void solve (scalar &x, scalarField &y, scalar &dxTry) const
 Solve the ODE system and the update the state. More...
 
- Public Member Functions inherited from ODESolver
 TypeName ("ODESolver")
 Runtime type information. More...
 
 declareRunTimeSelectionTable (autoPtr, ODESolver, dictionary,(const ODESystem &ode, const dictionary &dict),(ode, dict))
 
 ODESolver (const ODESystem &ode, const dictionary &dict)
 Construct for given ODESystem. More...
 
 ODESolver (const ODESystem &ode, const scalarField &absTol, const scalarField &relTol)
 Construct for given ODESystem specifying tolerances. More...
 
virtual ~ODESolver ()
 Destructor. More...
 
scalarFieldabsTol ()
 
scalarFieldrelTol ()
 
virtual void solve (scalar &x, scalarField &y, stepState &step) const
 Solve the ODE system as far as possible upto dxTry. More...
 
virtual void solve (const scalar xStart, const scalar xEnd, scalarField &y, scalar &dxEst) const
 Solve the ODE system from xStart to xEnd, update the state. More...
 
- Public Member Functions inherited from adaptiveSolver
 adaptiveSolver (const ODESystem &ode, const dictionary &dict)
 Construct from ODESystem. More...
 
virtual ~adaptiveSolver ()
 Destructor. More...
 
void solve (const ODESystem &ode, scalar &x, scalarField &y, scalar &dxTry) const
 Solve the ODE system and the update the state. More...
 

Additional Inherited Members

- Static Public Member Functions inherited from ODESolver
static autoPtr< ODESolverNew (const ODESystem &ode, const dictionary &dict)
 Select null constructed. More...
 
- Protected Member Functions inherited from ODESolver
scalar normalizeError (const scalarField &y0, const scalarField &y, const scalarField &err) const
 Return the nomalized scalar error. More...
 
 ODESolver (const ODESolver &)
 Disallow default bitwise copy construct. More...
 
void operator= (const ODESolver &)
 Disallow default bitwise assignment. More...
 
- Protected Attributes inherited from ODESolver
const ODESystemodes_
 Reference to ODESystem. More...
 
label n_
 Size of the ODESystem. More...
 
scalarField absTol_
 Absolute convergence tolerance per step. More...
 
scalarField relTol_
 Relative convergence tolerance per step. More...
 
label maxSteps_
 The maximum number of sub-steps allowed for the integration step. More...
 

Detailed Description

4/5th Order Dormand–Prince Runge-Kutta ODE solver.

References:

    "A family of embedded Runge-Kutta formulae"
    Dormand, J. R.,
    Prince, P. J.,
    Journal of Computational and Applied Mathematics,
    6 (1), 1980: pp. 19-26.

    "Solving Ordinary Differential Equations I: Nonstiff Problems,
     second edition",
    Hairer, E.,
    Nørsett, S.,
    Wanner, G.,
    Springer-Verlag, Berlin. 1993, ISBN 3-540-56670-8.
See also
Foam::RKF45 Foam::RKCK45
Source files

Definition at line 69 of file RKDP45.H.

Constructor & Destructor Documentation

RKDP45 ( const ODESystem ode,
const dictionary dict 
)

Construct from ODE.

Definition at line 79 of file RKDP45.C.

References RKDP45::solve().

Here is the call graph for this function:

Member Function Documentation

TypeName ( "RKDP45"  )

Runtime type information.

Foam::scalar solve ( const scalar  x0,
const scalarField y0,
const scalarField dydx0,
const scalar  dx,
scalarField y 
) const
virtual

Solve a single step dx and return the error.

Implements adaptiveSolver.

Definition at line 96 of file RKDP45.C.

References ODESystem::derivatives(), forAll, ODESolver::normalizeError(), and ODESolver::odes_.

Referenced by RKDP45::RKDP45().

Here is the call graph for this function:

Here is the caller graph for this function:

void solve ( scalar &  x,
scalarField y,
scalar &  dxTry 
) const
virtual

Solve the ODE system and the update the state.

Reimplemented from ODESolver.

Definition at line 164 of file RKDP45.C.

References ODESolver::odes_, and adaptiveSolver::solve().

Here is the call graph for this function:


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