(9)

To solve a system of ordinary differential equations. More...

Public Member Functions

 ODESolver ()
 Default constructor.
 
 ODESolver (size_t nb_eq)
 Constructor providing the number of equations.
 
 ODESolver (TimeScheme s, real_t time_step=theTimeStep, real_t final_time=theFinalTime, size_t nb_eq=1)
 Constructor using time discretization data. More...
 
 ODESolver (MyODE &my_ode, TimeScheme s, real_t time_step=theTimeStep, real_t final_time=theFinalTime, size_t nb_eq=1)
 Constructor using time discretization data and a user defined class to define the ODE. More...
 
 ~ODESolver ()
 Destructor.
 
void set (TimeScheme s, real_t time_step=theTimeStep, real_t final_time=theFinalTime)
 Define data of the differential equation or system. More...
 
void setNbEq (size_t n)
 Set the number of equations [Default: 1]. More...
 
void setCoef (real_t a0, real_t a1, real_t a2, real_t f)
 Define coefficients in the case of a scalar differential equation. More...
 
void setCoef (string a0, string a1, string a2, string f)
 Define coefficients in the case of a scalar differential equation. More...
 
void setLinear ()
 Claim that ODE is linear. More...
 
void setF (string f)
 Set time derivative, given as an algebraic expression, for a nonlinear ODE. More...
 
void setF (string f, int i)
 Set time derivative, given as an algebraic expression, for a nonlinear ODE. More...
 
void setDF (string df, int i, int j)
 Set time derivative of the function defining the ODE. More...
 
void setdFdt (string df, int i)
 Set time derivative of the function defining the ODE. More...
 
void setRK4RHS (real_t f)
 Set intermediate right-hand side vector for the Runge-Kutta method. More...
 
void setRK4RHS (Vect< real_t > &f)
 Set intermediate right-hand side vector for the Runge-Kutta method. More...
 
void setInitial (Vect< real_t > &u)
 Set initial condition for a first-oder system of differential equations. More...
 
void setInitial (real_t u, int i)
 Set initial condition for a first-oder system of differential equations. More...
 
void setInitial (Vect< real_t > &u, Vect< real_t > &v)
 Set initial condition for a second-order system of differential equations. More...
 
void setInitialRHS (Vect< real_t > &f)
 Set initial RHS for a system of differential equations. More...
 
void setInitial (real_t u, real_t v)
 Set initial condition for a second-order ordinary differential equation. More...
 
void setInitial (real_t u)
 Set initial condition for a first-order ordinary differential equation. More...
 
void setInitialRHS (real_t f)
 Set initial right-hand side for a single differential equation. More...
 
void setMatrices (DMatrix< real_t > &A0, DMatrix< real_t > &A1)
 Define matrices for a system of first-order ODEs. More...
 
void setMatrices (DMatrix< real_t > &A0, DMatrix< real_t > &A1, DMatrix< real_t > &A2)
 Define matrices for a system of second-order ODEs. More...
 
void seODEVectors (Vect< real_t > &a0, Vect< real_t > &a1)
 Define matrices for an implicit nonlinear system of first-order ODEs. More...
 
void seODEVectors (Vect< real_t > &a0, Vect< real_t > &a1, Vect< real_t > &a2)
 Define matrices for an implicit nonlinear system of second-order ODEs. More...
 
void setRHS (Vect< real_t > &b)
 Set right-hand side vector for a system of ODE. More...
 
void setRHS (real_t f)
 Set right-hand side for a linear ODE. More...
 
void setRHS (string f)
 Set right-hand side value for a linear ODE.
 
void setNewmarkParameters (real_t beta, real_t gamma)
 Define parameters for the Newmarxk scheme. More...
 
void setConstantMatrix ()
 Say that matrix problem is constant. More...
 
void setNonConstantMatrix ()
 Say that matrix problem is variable. More...
 
void setLinearSolver (Iteration s=DIRECT_SOLVER, Preconditioner p=DIAG_PREC)
 Set linear solver data. More...
 
void setMaxIter (int max_it)
 Set maximal number of iterations. More...
 
void setTolerance (real_t toler)
 Set tolerance value for convergence. More...
 
real_t runOneTimeStep ()
 Run one time step. More...
 
void run (bool opt=false)
 Run the time stepping procedure. More...
 
size_t getNbEq () const
 Return number of equations.
 
LinearSolvergetLSolver ()
 Return LinearSolver instance.
 
real_t getTimeDerivative (int i=1) const
 Get time derivative of solution. More...
 
void getTimeDerivative (Vect< real_t > &y) const
 Get time derivative of solution (for a system) More...
 
real_t get () const
 Return solution in the case of a scalar equation.
 
void getPhase (Vect< real_t > &x, Vect< real_t > &v, size_t i=1)
 Get phase portrait vectors. More...
 

Detailed Description

To solve a system of ordinary differential equations.

The class ODESolver enables solving by a numerical scheme a system or ordinary differential equations taking one of the forms:

  • A linear system of differential equations of the first-order:
    A1(t)u'(t) + A0(t)u(t) = f(t)
  • A linear system of differential equations of the second-order:
    A2(t)u''(t) + A1(t)u'(t) + A0(t)u(t) = f(t)
  • A system of ordinary differential equations of the form:
    u'(t) = f(t,u(t))

The following time integration schemes can be used:

  • Forward Euler scheme (value: FORWARD_EULER) for first-order systems
  • Backward Euler scheme (value: BACKWARD_EULER) for first-order linear systems
  • Crank-Nicolson (value: CRANK_NICOLSON) for first-order linear systems
  • Heun (value: HEUN) for first-order systems
  • 2nd Order Adams-Bashforth (value: AB2) for first-order systems
  • 4-th order Runge-Kutta (value: RK4) for first-order systems
  • 2nd order Backward Differentiation Formula (value: BDF2) for linear first-order systems
  • Newmark (value: NEWMARK) for linear second-order systems with constant matrices
Author
Rachid Touzani

Constructor & Destructor Documentation

◆ ODESolver() [1/2]

ODESolver ( TimeScheme  s,
real_t  time_step = theTimeStep,
real_t  final_time = theFinalTime,
size_t  nb_eq = 1 
)

Constructor using time discretization data.

Parameters
[in]sChoice of the scheme: To be chosen in the enumerated variable Scheme (see the presentation of the class)
[in]time_stepValue of the time step. This value will be modified if an adaptive method is used. The default value for this parameter if the value given by the global variable theTimeStep
[in]final_timeValue of the final time (time starts at 0). The default value for this parameter is the value given by the global variable theFinalTime
[in]nb_eqNumber of differential equations (size of the system) [Default: 1]

◆ ODESolver() [2/2]

ODESolver ( MyODE my_ode,
TimeScheme  s,
real_t  time_step = theTimeStep,
real_t  final_time = theFinalTime,
size_t  nb_eq = 1 
)

Constructor using time discretization data and a user defined class to define the ODE.

Parameters
[in]my_odeReference to instance of user defined class. This class inherits from abstract class MyODE. It must contain the member function double Function(double t, const double& y) which returns the value of the nonlinear function, as a vector, for a given solution vector x. The user defined class must contain, if the iterative scheme requires it the member function Vect<double> Gradient(const Vect<real_t>& x) which returns the gradient as a n*n vector, each index (i,j) containing the j-th partial derivative of the i-th function.
[in]sChoice of the scheme: To be chosen in the enumerated variable Scheme (see the presentation of the class)
[in]time_stepValue of the time step. This value will be modified if an adaptive method is used. The default value for this parameter if the value given by the global variable theTimeStep
[in]final_timeValue of the final time (time starts at 0). The default value for this parameter is the value given by the global variable theFinalTime
[in]nb_eqNumber of differential equations (size of the system) [Default: 1]

Member Function Documentation

◆ getPhase()

void getPhase ( Vect< real_t > &  x,
Vect< real_t > &  v,
size_t  i = 1 
)

Get phase portrait vectors.

This function gets vectors containing solution and its time derivative to determine phase portraite of the ode. The function is valid for a single ode only.

Parameters
xVector containing solution of the ode at each computed time step
vVector containing discrete time derivative of the ode at each computed time step
iComponent for which the phase is extracted [Dafault: 1]

◆ getTimeDerivative() [1/2]

real_t getTimeDerivative ( int  i = 1) const

Get time derivative of solution.

Return approximate time derivative of solution in the case of a single equation

Parameters
[in]iIndex of component whose time derivative is sought
Returns
Time derivative of the i-th component of the solution
Remarks
If we are solving one equation, this parameter is not used.

◆ getTimeDerivative() [2/2]

void getTimeDerivative ( Vect< real_t > &  y) const

Get time derivative of solution (for a system)

Get approximate time derivative of solution in the case of an ODE system

Parameters
[out]yVector containing time derivative of solution

◆ run()

void run ( bool  opt = false)

Run the time stepping procedure.

Parameters
[in]optFlag to say if problem matrix is constant while time stepping (true) or not (Default value is false)
Note
This argument is not used if the time stepping scheme is explicit

◆ runOneTimeStep()

real_t runOneTimeStep ( )

Run one time step.

Returns
Value of new time step if this one is updated

◆ seODEVectors() [1/2]

void seODEVectors ( Vect< real_t > &  a0,
Vect< real_t > &  a1 
)

Define matrices for an implicit nonlinear system of first-order ODEs.

The system has the nonlinear implicit form a1(u)' + a0(u) = 0 Vectors a0, a1 are given as references to class Vect.

Parameters
[in]a0Reference to vector in front of the 0-th order term (no time derivative)
[in]a1Reference to vector in front of the 1-st order term (first time derivative)
Remarks
This function has to be called at each time step

◆ seODEVectors() [2/2]

void seODEVectors ( Vect< real_t > &  a0,
Vect< real_t > &  a1,
Vect< real_t > &  a2 
)

Define matrices for an implicit nonlinear system of second-order ODEs.

The system has the nonlinear implicit form a2(u)'' + a1(u)' + a0(u) = 0 Vectors a0, a1, a2 are given as references to class Vect.

Parameters
[in]a0Reference to vector in front of the 0-th order term (no time derivative)
[in]a1Reference to vector in front of the 1-st order term (first time derivative)
[in]a2Reference to vector in front of the 2-nd order term (second time derivative)
Remarks
This function has to be called at each time step

◆ set()

void set ( TimeScheme  s,
real_t  time_step = theTimeStep,
real_t  final_time = theFinalTime 
)

Define data of the differential equation or system.

Parameters
[in]sChoice of the scheme: To be chosen in the enumerated variable Scheme (see the presentation of the class)
[in]time_stepValue of the time step. This value will be modified if an adaptive method is used. The default value for this parameter if the value given by the global variable theTimeStep
[in]final_timeValue of the final time (time starts at 0). The default value for this parameter is the value given by the global variable theFinalTime

◆ setCoef() [1/2]

void setCoef ( real_t  a0,
real_t  a1,
real_t  a2,
real_t  f 
)

Define coefficients in the case of a scalar differential equation.

This function enables giving coefficients of the differential equation as an algebraic expression of time t (see the function fparse)

Parameters
[in]a0Coefficient of the 0-th order term
[in]a1Coefficient of the 1-st order term
[in]a2Coefficient of the 2-nd order term
[in]fValue of the right-hand side
Note
Naturally, the equation is of the first order if a2=0

◆ setCoef() [2/2]

void setCoef ( string  a0,
string  a1,
string  a2,
string  f 
)

Define coefficients in the case of a scalar differential equation.

Parameters
[in]a0Coefficient of the 0-th order term
[in]a1Coefficient of the 1-st order term
[in]a2Coefficient of the 2-nd order term
[in]fValue of the right-hand side
Note
Naturally, the equation if of the first order if a2=0

◆ setConstantMatrix()

void setConstantMatrix ( )

Say that matrix problem is constant.

This is useful if the linear system is solved by a factorization method but has no effect otherwise

◆ setDF()

void setDF ( string  df,
int  i,
int  j 
)

Set time derivative of the function defining the ODE.

This function enables prescribing the value of the 1-st derivative for a 1st order ODE or the 2nd one for a 2nd-order ODE. It is to be used for nonlinear ODEs of the form y'(t) = f(t,y(t)) or y''(t) = f(t,y(t),y'(t))
In the case of a system of ODEs, this function can be called once for each equation, given in the order of the unknowns

◆ setdFdt()

void setdFdt ( string  df,
int  i 
)

Set time derivative of the function defining the ODE.

This function enables prescribing the value of the 1-st derivative for a 1st order ODE or the 2nd one for a 2nd-order ODE. It is to be used for nonlinear ODEs of the form y'(t) = f(t,y(t)) or y''(t) = f(t,y(t),y'(t))
In the case of a system of ODEs, this function can be called once for each equation, given in the order of the unknowns

◆ setF() [1/2]

void setF ( string  f)

Set time derivative, given as an algebraic expression, for a nonlinear ODE.

This function enables prescribing the value of the 1-st derivative for a 1st order ODE or the 2nd one for a 2nd-order ODE. It is to be used for nonlinear ODEs of the form y'(t) = f(t,y(t)) or y''(t) = f(t,y(t),y'(t))
In the case of a system of ODEs, this function can be called once for each equation, given in the order of the unknowns

Parameters
[in]fExpression of the function

◆ setF() [2/2]

void setF ( string  f,
int  i 
)

Set time derivative, given as an algebraic expression, for a nonlinear ODE.

This function enables prescribing the value of the 1-st derivative for a 1st order ODE or the 2nd one for a 2nd-order ODE. It is to be used for nonlinear ODEs of the form y'(t) = f(t,y(t)) or y''(t) = f(t,y(t),y'(t))
This function is to be used for the i-th equation of a system of ODEs

Parameters
[in]fExpression of the function
[in]iIndex of equation. Must be not larger than the number of equations

◆ setInitial() [1/5]

void setInitial ( Vect< real_t > &  u)

Set initial condition for a first-oder system of differential equations.

Parameters
[in]uVector containing initial condition for the unknown

◆ setInitial() [2/5]

void setInitial ( real_t  u,
int  i 
)

Set initial condition for a first-oder system of differential equations.

Parameters
[in]uInitial condition for an unknown
[in]iIndex of the unknown

◆ setInitial() [3/5]

void setInitial ( Vect< real_t > &  u,
Vect< real_t > &  v 
)

Set initial condition for a second-order system of differential equations.

Giving the right-hand side at initial time is somtimes required for high order methods like Runge-Kutta

Parameters
[in]uVector containing initial condition for the unknown
[in]vVector containing initial condition for the time derivative of the unknown

◆ setInitial() [4/5]

void setInitial ( real_t  u,
real_t  v 
)

Set initial condition for a second-order ordinary differential equation.

Parameters
[in]uInitial condition (unknown) value
[in]vInitial condition (time derivative of the unknown) value

◆ setInitial() [5/5]

void setInitial ( real_t  u)

Set initial condition for a first-order ordinary differential equation.

Parameters
[in]uInitial condition (unknown) value

◆ setInitialRHS() [1/2]

void setInitialRHS ( Vect< real_t > &  f)

Set initial RHS for a system of differential equations.

Giving the right-hand side at initial time is somtimes required for high order methods like Runge-Kutta

Parameters
[in]fVector containing right-hand side at initial time. This vector is helpful for high order methods

◆ setInitialRHS() [2/2]

void setInitialRHS ( real_t  f)

Set initial right-hand side for a single differential equation.

Parameters
[in]fValue of right-hand side at initial time. This value is helpful for high order methods

◆ setLinear()

void setLinear ( )

Claim that ODE is linear.

Claim that the defined ODE (or system of ODEs) is linear

◆ setLinearSolver()

void setLinearSolver ( Iteration  s = DIRECT_SOLVER,
Preconditioner  p = DIAG_PREC 
)

Set linear solver data.

Parameters
[in]sSolver identification parameter. To be chosen in the enumeration variable Iteration:
DIRECT_SOLVER, CG_SOLVER, CGS_SOLVER, BICG_SOLVER, BICG_STAB_SOLVER, GMRES_SOLVER, QMR_SOLVER [Default: DIRECT_SOLVER]
[in]pPreconditioner identification parameter. To be chosen in the enumeration variable Preconditioner:
IDENT_PREC, DIAG_PREC, ILU_PREC [Default: DIAG_PREC]
Note
The argument p has no effect if the solver is DIRECT_SOLVER

◆ setMatrices() [1/2]

void setMatrices ( DMatrix< real_t > &  A0,
DMatrix< real_t > &  A1 
)

Define matrices for a system of first-order ODEs.

Matrices are given as references to class DMatrix.

Parameters
[in]A0Reference to matrix in front of the 0-th order term (no time derivative)
[in]A1Reference to matrix in front of the 1-st order term (first time derivative)
Remarks
This function has to be called at each time step

◆ setMatrices() [2/2]

void setMatrices ( DMatrix< real_t > &  A0,
DMatrix< real_t > &  A1,
DMatrix< real_t > &  A2 
)

Define matrices for a system of second-order ODEs.

Matrices are given as references to class DMatrix.

Parameters
[in]A0Reference to matrix in front of the 0-th order term (no time derivative)
[in]A1Reference to matrix in front of the 1-st order term (first time derivative)
[in]A2Reference to matrix in front of the 2-nd order term (second time derivative)
Remarks
This function has to be called at each time step

◆ setMaxIter()

void setMaxIter ( int  max_it)

Set maximal number of iterations.

This function is useful for a non linear ODE (or system of ODEs) if an implicit scheme is used

Parameters
[in]max_itMaximal number of iterations [Default: 100]

◆ setNbEq()

void setNbEq ( size_t  n)

Set the number of equations [Default: 1].

This function is to be used if the default constructor was used

◆ setNewmarkParameters()

void setNewmarkParameters ( real_t  beta,
real_t  gamma 
)

Define parameters for the Newmarxk scheme.

Parameters
[in]betaParameter beta [Default: 0.25]
[in]gammaParameter gamma [Default: 0.5]

◆ setNonConstantMatrix()

void setNonConstantMatrix ( )

Say that matrix problem is variable.

This is useful if the linear system is solved by a factorization method but has no effect otherwise

◆ setRHS() [1/2]

void setRHS ( Vect< real_t > &  b)

Set right-hand side vector for a system of ODE.

Parameters
[in]bVect instance containing right-hand side for a linear system of ordinary differential equations

◆ setRHS() [2/2]

void setRHS ( real_t  f)

Set right-hand side for a linear ODE.

Parameters
[in]fValue of the right-hand side for a linear ordinary differential equation

◆ setRK4RHS() [1/2]

void setRK4RHS ( real_t  f)

Set intermediate right-hand side vector for the Runge-Kutta method.

Parameters
[in]fValue of right-hand side

◆ setRK4RHS() [2/2]

void setRK4RHS ( Vect< real_t > &  f)

Set intermediate right-hand side vector for the Runge-Kutta method.

Parameters
[in]fright-hand side vector

◆ setTolerance()

void setTolerance ( real_t  toler)

Set tolerance value for convergence.

This function is useful for a non linear ODE (or system of ODEs) if an implicit scheme is used

Parameters
[in]tolerTolerance value [Default: 1.e-8]