Public Member Functions | Protected Member Functions | Protected Attributes

Abstract base class for all kinds of Runge-Kutta schemes for integrating ODEs. More...

#include <integrator_runge_kutta.hpp>

Inheritance diagram for IntegratorRK:
Inheritance graph
[legend]

List of all members.

Public Member Functions

virtual Integratorclone () const =0
virtual returnValue freezeAll ()
virtual returnValue freezeMesh ()
virtual int getNumberOfRejectedSteps () const
virtual int getNumberOfSteps () const
virtual double getStepSize () const
virtual returnValue init (const DifferentialEquation &rhs_)
returnValue init (const DifferentialEquation &rhs_, const Transition &trs_)
 IntegratorRK ()
 IntegratorRK (int dim_, double power_)
 IntegratorRK (const DifferentialEquation &rhs_, int dim_, double power_)
 IntegratorRK (const IntegratorRK &arg)
virtual IntegratorRKoperator= (const IntegratorRK &arg)
virtual returnValue setDxInitialization (double *dx0)
virtual returnValue step (int number)
virtual returnValue stop ()
virtual returnValue unfreeze ()
virtual ~IntegratorRK ()

Protected Member Functions

void allocateMemory ()
void constructAll (const IntegratorRK &arg)
void deleteAll ()
double determineEta45 ()
double determineEta45 (int number)
void determineEtaGForward (int number)
void determineEtaGForward2 (int number)
void determineEtaHBackward (int number)
void determineEtaHBackward2 (int number)
virtual returnValue evaluate (const DVector &x0, const DVector &xa, const DVector &p, const DVector &u, const DVector &w, const Grid &t_)
virtual returnValue evaluateSensitivities ()
virtual int getDim () const
virtual returnValue getProtectedBackwardSensitivities (DVector &Dx_x0, DVector &Dx_p, DVector &Dx_u, DVector &Dx_w, int order) const
virtual returnValue getProtectedForwardSensitivities (DMatrix *Dx, int order) const
virtual returnValue getProtectedX (DVector *xEnd) const
virtual void initializeButcherTableau ()=0
void initializeVariables ()
void interpolate (int jj, double *e1, double *d1, double *e2, VariablesGrid &poly)
void logCurrentIntegratorStep (const DVector &currentX=emptyConstVector)
void printIntermediateResults ()
virtual returnValue setBackwardSeed2 (const DVector &seed)
returnValue setForwardSeed2 (const DVector &xSeed, const DVector &pSeed, const DVector &uSeed, const DVector &wSeed)
virtual returnValue setProtectedBackwardSeed (const DVector &seed, const int &order)
virtual returnValue setProtectedForwardSeed (const DVector &xSeed, const DVector &pSeed, const DVector &uSeed, const DVector &wSeed, const int &order)

Protected Attributes

double ** A
double * b4
double * b5
DVector bseed
DVector bseed2
double * c
int dim
double err_power
double * eta4
double * eta4_
double * eta5
double * eta5_
double * etaG
double * etaG2
double * etaG3
double * etaH
double * etaH2
double * etaH3
DVector fseed
DVector fseed2
double * G
double * G2
double * G3
double * H
double * H2
double * H3
double ** k
double ** k2
double ** l
double ** l2
int maxAlloc
double t
double * x

Detailed Description

Abstract base class for all kinds of Runge-Kutta schemes for integrating ODEs.

The class IntegratorRK serves as an abstract base class for all kinds of Runge-Kutta schemes for integrating ordinary differential equations (ODEs).

Author:
Boris Houska, Hans Joachim Ferreau

Definition at line 54 of file integrator_runge_kutta.hpp.


Constructor & Destructor Documentation

Default constructor.

Definition at line 51 of file integrator_runge_kutta.cpp.

IntegratorRK::IntegratorRK ( int  dim_,
double  power_ 
)

Default constructor.

Definition at line 58 of file integrator_runge_kutta.cpp.

IntegratorRK::IntegratorRK ( const DifferentialEquation rhs_,
int  dim_,
double  power_ 
)

Default constructor.

Definition at line 80 of file integrator_runge_kutta.cpp.

Copy constructor (deep copy).

Definition at line 103 of file integrator_runge_kutta.cpp.

Destructor.

Definition at line 110 of file integrator_runge_kutta.cpp.


Member Function Documentation

void IntegratorRK::allocateMemory ( ) [protected]

This routine is protected and sets up all
variables (i.e. allocates memory etc.).
Note that this routine assumes that the
dimensions are already set correctly and is
thus for internal use only.

Definition at line 151 of file integrator_runge_kutta.cpp.

virtual Integrator* IntegratorRK::clone ( ) const [pure virtual]

The (virtual) copy constructor

Implements Integrator.

Implemented in IntegratorDiscretizedODE, IntegratorRK23, IntegratorRK45, IntegratorRK78, and IntegratorRK12.

void IntegratorRK::constructAll ( const IntegratorRK arg) [protected]

Implementation of the copy constructor.

Definition at line 407 of file integrator_runge_kutta.cpp.

void IntegratorRK::deleteAll ( ) [protected]

Implementation of the delete operator.

Definition at line 291 of file integrator_runge_kutta.cpp.

double IntegratorRK::determineEta45 ( ) [protected]

computes eta4 and eta5 (only for internal use)

Returns:
The error estimate.

Definition at line 1625 of file integrator_runge_kutta.cpp.

double IntegratorRK::determineEta45 ( int  number) [protected]

computes eta4 and eta5 (only for internal use)

Returns:
The error estimate.

Definition at line 1689 of file integrator_runge_kutta.cpp.

void IntegratorRK::determineEtaGForward ( int  number) [protected]

computes etaG in forward direction (only for internal use)

Definition at line 1748 of file integrator_runge_kutta.cpp.

void IntegratorRK::determineEtaGForward2 ( int  number) [protected]

computes etaG and etaG2 in forward direction
(only for internal use)

Definition at line 1779 of file integrator_runge_kutta.cpp.

void IntegratorRK::determineEtaHBackward ( int  number) [protected]

computes etaH in backward direction (only for internal use)

Definition at line 1817 of file integrator_runge_kutta.cpp.

void IntegratorRK::determineEtaHBackward2 ( int  number) [protected]

computes etaH2 in backward direction (only for internal use)

Definition at line 1851 of file integrator_runge_kutta.cpp.

returnValue IntegratorRK::evaluate ( const DVector x0,
const DVector xa,
const DVector p,
const DVector u,
const DVector w,
const Grid t_ 
) [protected, virtual]

Starts integration: cf. integrate(...) for
more details.

Parameters:
x0the initial state
xathe initial algebraic state
pthe parameters
uthe controls
wthe disturbance
t_the time interval

Implements Integrator.

Definition at line 645 of file integrator_runge_kutta.cpp.

< Integrates forward and/or backward depending on the specified seeds.

Returns:
SUCCESSFUL_RETURN
RET_NO_SEED_ALLOCATED

Implements Integrator.

Definition at line 1075 of file integrator_runge_kutta.cpp.

Freezes the mesh as well as all intermediate values. This function
is necessary for the case that automatic differentiation in backward
mode should is used. (Note: This function might for large right hand
sides lead to memory problems as all intemediate values will be
stored!)

Returns:
SUCCESSFUL_RETURN
RET_ALREADY_FROZEN

Implements Integrator.

Definition at line 621 of file integrator_runge_kutta.cpp.

Freezes the mesh: Storage of the step sizes. If the integrator is
freezed the mesh will be stored when calling the function integrate
for the first time. If the function integrate is called more than
once the same mesh will be reused (i.e. the step size control will
be turned off). Note that the mesh should be frozen if any kind of
sensitivity generation is used.

Returns:
SUCCESSFUL_RETURN
RET_ALREADY_FROZEN

Implements Integrator.

Definition at line 605 of file integrator_runge_kutta.cpp.

int IntegratorRK::getDim ( ) const [protected, virtual]

Returns the dimension of the Differential Equation

Implements Integrator.

Definition at line 2059 of file integrator_runge_kutta.cpp.

int IntegratorRK::getNumberOfRejectedSteps ( ) const [virtual]

Returns the number of rejected Steps.

Returns:
The requested number of rejected steps.

Implements Integrator.

Definition at line 1603 of file integrator_runge_kutta.cpp.

int IntegratorRK::getNumberOfSteps ( ) const [virtual]

Returns the number of accepted Steps.

Returns:
The requested number of accepted steps.

Implements Integrator.

Definition at line 1598 of file integrator_runge_kutta.cpp.

returnValue IntegratorRK::getProtectedBackwardSensitivities ( DVector Dx_x0,
DVector Dx_p,
DVector Dx_u,
DVector Dx_w,
int  order 
) const [protected, virtual]

Returns the result for the backward sensitivities at the time tend.

Parameters:
Dx_x0backward sensitivities w.r.t. the initial states
Dx_pbackward sensitivities w.r.t. the parameters
Dx_ubackward sensitivities w.r.t. the controls
Dx_wbackward sensitivities w.r.t. the disturbance
orderthe order of the derivative

Returns:
SUCCESSFUL_RETURN
RET_INPUT_OUT_OF_RANGE

Implements Integrator.

Definition at line 1509 of file integrator_runge_kutta.cpp.

returnValue IntegratorRK::getProtectedForwardSensitivities ( DMatrix Dx,
int  order 
) const [protected, virtual]

Returns the result for the forward sensitivities at the time tend.

Returns:
SUCCESSFUL_RETURN
RET_INPUT_OUT_OF_RANGE
Parameters:
Dxthe result for the forward sensitivi- ties
orderthe order

Implements Integrator.

Definition at line 1475 of file integrator_runge_kutta.cpp.

returnValue IntegratorRK::getProtectedX ( DVector xEnd) const [protected, virtual]

Returns the result for the state at the time tend.

Returns:
SUCCESSFUL_RETURN
Parameters:
xEndthe result for the states at the time tend.

Implements Integrator.

Definition at line 1461 of file integrator_runge_kutta.cpp.

double IntegratorRK::getStepSize ( ) const [virtual]

Returns the current step size

Implements Integrator.

Definition at line 1609 of file integrator_runge_kutta.cpp.

The initialization routine which takes the right-hand side of
the differential equation to be integrated.

Parameters:
rhsthe right-hand side of the ODE/DAE.

Returns:
SUCCESSFUL_RETURN if all dimension checks succeed.
otherwise: integrator dependent error message.

Implements Integrator.

Reimplemented in IntegratorDiscretizedODE, and IntegratorRK12.

Definition at line 116 of file integrator_runge_kutta.cpp.

returnValue IntegratorRK::init ( const DifferentialEquation rhs_,
const Transition trs_ 
) [inline]

The initialization routine which takes the right-hand side of
the differential equation to be integrated. In addition a
transition function can be set which is evaluated at the end
of the integration interval.

Parameters:
rhsthe right-hand side of the ODE/DAE.
trsthe transition to be evaluated at the end.

Returns:
SUCCESSFUL_RETURN if all dimension checks succeed.
otherwise: integrator dependent error message.

Reimplemented from Integrator.

virtual void IntegratorRK::initializeButcherTableau ( ) [protected, pure virtual]

This routine initializes the coefficients of the Butcher Tableau.

Implemented in IntegratorRK12, IntegratorRK23, IntegratorRK45, and IntegratorRK78.

void IntegratorRK::initializeVariables ( ) [protected]

This routine is protected and is basically used
to set all pointer-valued member to the NULL pointer.
In addition some dimensions are initialized with 0 as
a default value.

Definition at line 134 of file integrator_runge_kutta.cpp.

void IntegratorRK::interpolate ( int  jj,
double *  e1,
double *  d1,
double *  e2,
VariablesGrid poly 
) [protected]

Definition at line 2042 of file integrator_runge_kutta.cpp.

void IntegratorRK::logCurrentIntegratorStep ( const DVector currentX = emptyConstVector) [protected]
IntegratorRK & IntegratorRK::operator= ( const IntegratorRK arg) [virtual]

Assignment operator (deep copy).

Definition at line 395 of file integrator_runge_kutta.cpp.

prints intermediate results for the case that the PrintLevel is
HIGH.

Definition at line 1891 of file integrator_runge_kutta.cpp.

returnValue IntegratorRK::setBackwardSeed2 ( const DVector seed) [protected, virtual]

Initializes a second backward seed. (only for internal use)

Parameters:
seedthe seed matrix

Definition at line 1029 of file integrator_runge_kutta.cpp.

returnValue IntegratorRK::setDxInitialization ( double *  dx0) [virtual]

Sets an initial guess for the differential state derivatives
(consistency condition)

Returns:
SUCCESSFUL_RETURN
Parameters:
dx0initial guess for the differential state derivatives

Implements Integrator.

Definition at line 1615 of file integrator_runge_kutta.cpp.

returnValue IntegratorRK::setForwardSeed2 ( const DVector xSeed,
const DVector pSeed,
const DVector uSeed,
const DVector wSeed 
) [protected]

Initializes a second forward seed. (only for internal use)

Parameters:
xSeedthe seed w.r.t the initial states
pSeedthe seed w.r.t the parameters
uSeedthe seed w.r.t the controls
wSeedthe seed w.r.t the disturbances

Definition at line 909 of file integrator_runge_kutta.cpp.

returnValue IntegratorRK::setProtectedBackwardSeed ( const DVector seed,
const int &  order 
) [protected, virtual]

Define a backward seed

Returns:
SUCCESFUL_RETURN
RET_INPUT_OUT_OF_RANGE
Parameters:
seedthe seed matrix
orderthe order of the seed.

Implements Integrator.

Definition at line 981 of file integrator_runge_kutta.cpp.

returnValue IntegratorRK::setProtectedForwardSeed ( const DVector xSeed,
const DVector pSeed,
const DVector uSeed,
const DVector wSeed,
const int &  order 
) [protected, virtual]

Define a forward seed.

Returns:
SUCCESFUL RETURN
RET_INPUT_OUT_OF_RANGE
Parameters:
xSeedthe seed w.r.t the initial states
pSeedthe seed w.r.t the parameters
uSeedthe seed w.r.t the controls
wSeedthe seed w.r.t the disturbances
orderthe order of the seed.

Implements Integrator.

Definition at line 835 of file integrator_runge_kutta.cpp.

returnValue IntegratorRK::step ( int  number) [virtual]

Executes the next single step. This function can be used to
call the integrator step wise. Note that this function is e.g.
useful in real-time simulations where after each step a time
out limit has to be checked. This function will usually return

Returns:
RET_FINAL_STEP_NOT_PERFORMED_YET
for the case that the final step has not been
performed yet, i.e. the integration routine is not yet
at the time tend.
Otherwise it will either return
SUCCESSFUL_RETURN
or any other error message that might occur during
an integration step.

In most real situations you can define the maximum number of
step sizes to be 1 before calling the function integrate
Then the function integrate should return after one step with
the message
RET_MAXIMUM_NUMBER_OF_STEPS_EXCEEDED. After that you can call
step() until the final time is reached.
(You can use the PrintLevel NONE to avoid that the message
RET_MAXIMUM_NUMBER_OF_STEPS_EXCEEDED is printed.)

Parameters:
numberthe step number

Reimplemented in IntegratorDiscretizedODE, and IntegratorRK12.

Definition at line 1194 of file integrator_runge_kutta.cpp.

Stops the integration even if the final time has not been
reached yet. This function will also give all memory free.
In particular, the function unfreeze() will be called.
(This function is designed for the usage in real-time
contexts in order to deal with error messages without
deleting and re-initializing the integrator.)

Returns:
SUCCESSFUL_RETURN

Definition at line 1455 of file integrator_runge_kutta.cpp.

Unfreezes the mesh: Gives the memory free that has previously
been allocated by "freeze". If you use the function
integrate after unfreezing the usual step size control will be
switched on.

Returns:
SUCCESSFUL_RETURN
RET_MESH_ALREADY_UNFROZED

Implements Integrator.

Definition at line 635 of file integrator_runge_kutta.cpp.


Member Data Documentation

double** IntegratorRK::A [protected]

the coefficient A of the Butcher Tableau.

Definition at line 436 of file integrator_runge_kutta.hpp.

double* IntegratorRK::b4 [protected]

the 4th order coefficients of the Butcher Tableau.

Definition at line 437 of file integrator_runge_kutta.hpp.

double* IntegratorRK::b5 [protected]

the 5th order coefficients of the Butcher Tableau.

Definition at line 438 of file integrator_runge_kutta.hpp.

The backward seed (only internal use)

Definition at line 460 of file integrator_runge_kutta.hpp.

The backward seed 2 (only internal use)

Definition at line 463 of file integrator_runge_kutta.hpp.

double* IntegratorRK::c [protected]

the time coefficients of the Butcher Tableau.

Definition at line 439 of file integrator_runge_kutta.hpp.

int IntegratorRK::dim [protected]

the dimension of the Butcher Tableau.

Definition at line 435 of file integrator_runge_kutta.hpp.

double IntegratorRK::err_power [protected]

root order of the step size control

Definition at line 454 of file integrator_runge_kutta.hpp.

double* IntegratorRK::eta4 [protected]

the result of order 4

Definition at line 444 of file integrator_runge_kutta.hpp.

double* IntegratorRK::eta4_ [protected]

the result of order 4

Definition at line 446 of file integrator_runge_kutta.hpp.

double* IntegratorRK::eta5 [protected]

the result of order 5

Definition at line 445 of file integrator_runge_kutta.hpp.

double* IntegratorRK::eta5_ [protected]

the result of order 5

Definition at line 447 of file integrator_runge_kutta.hpp.

double* IntegratorRK::etaG [protected]

Sensitivity matrix (only internal use)

Definition at line 466 of file integrator_runge_kutta.hpp.

double* IntegratorRK::etaG2 [protected]

Sensitivity matrix (only internal use)

Definition at line 470 of file integrator_runge_kutta.hpp.

double* IntegratorRK::etaG3 [protected]

Sensitivity matrix (only internal use)

Definition at line 471 of file integrator_runge_kutta.hpp.

double* IntegratorRK::etaH [protected]

Sensitivity matrix (only internal use)

Definition at line 474 of file integrator_runge_kutta.hpp.

double* IntegratorRK::etaH2 [protected]

Sensitivity matrix (only internal use)

Definition at line 478 of file integrator_runge_kutta.hpp.

double* IntegratorRK::etaH3 [protected]

Sensitivity matrix (only internal use)

Definition at line 479 of file integrator_runge_kutta.hpp.

The forward seed (only internal use)

Definition at line 459 of file integrator_runge_kutta.hpp.

The forward seed 2 (only internal use)

Definition at line 462 of file integrator_runge_kutta.hpp.

double* IntegratorRK::G [protected]

Sensitivity matrix (only internal use)

Definition at line 465 of file integrator_runge_kutta.hpp.

double* IntegratorRK::G2 [protected]

Sensitivity matrix (only internal use)

Definition at line 468 of file integrator_runge_kutta.hpp.

double* IntegratorRK::G3 [protected]

Sensitivity matrix (only internal use)

Definition at line 469 of file integrator_runge_kutta.hpp.

double* IntegratorRK::H [protected]

Sensitivity matrix (only internal use)

Definition at line 473 of file integrator_runge_kutta.hpp.

double* IntegratorRK::H2 [protected]

Sensitivity matrix (only internal use)

Definition at line 476 of file integrator_runge_kutta.hpp.

double* IntegratorRK::H3 [protected]

Sensitivity matrix (only internal use)

Definition at line 477 of file integrator_runge_kutta.hpp.

double** IntegratorRK::k [protected]

the intermediate results

Definition at line 448 of file integrator_runge_kutta.hpp.

double** IntegratorRK::k2 [protected]

the intermediate results

Definition at line 449 of file integrator_runge_kutta.hpp.

double** IntegratorRK::l [protected]

the intermediate results

Definition at line 450 of file integrator_runge_kutta.hpp.

double** IntegratorRK::l2 [protected]

the intermediate results

Definition at line 451 of file integrator_runge_kutta.hpp.

int IntegratorRK::maxAlloc [protected]

size of the memory that is allocated to store
the trajectory and the mesh.

Definition at line 484 of file integrator_runge_kutta.hpp.

double IntegratorRK::t [protected]

the actual time

Definition at line 452 of file integrator_runge_kutta.hpp.

double* IntegratorRK::x [protected]

the actual state (only internal use)

Definition at line 453 of file integrator_runge_kutta.hpp.


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


acado
Author(s): Milan Vukov, Rien Quirynen
autogenerated on Sat Jun 8 2019 19:40:24