To solve time stepping problems, i.e. systems of linear ordinary differential equations of the form [A2]{y"} + [A1]{y'} + [A0]{y} = {b}. More...
Public Member Functions | |
TimeStepping () | |
Default constructor. | |
TimeStepping (TimeScheme s, real_t time_step=theTimeStep, real_t final_time=theFinalTime) | |
Constructor using time discretization data. More... | |
~TimeStepping () | |
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 | setPDE (AbsEqua< real_t > &eq) |
Define partial differential equation to solve. More... | |
void | setRK4RHS (Vect< real_t > &f) |
Set intermediate right-hand side vector for the Runge-Kutta method. More... | |
void | setRK3_TVDRHS (Vect< real_t > &f) |
Set intermediate right-hand side vector for the TVD Runge-Kutta 3 method. More... | |
void | setInitial (Vect< real_t > &u) |
Set initial condition for the system of differential equations. More... | |
void | setInitial (Vect< real_t > &u, Vect< real_t > &v) |
Set initial condition for a system of differential equations. More... | |
void | setInitialRHS (Vect< real_t > &f) |
Set initial RHS for a system of differential equations when the used scheme requires it. More... | |
void | setRHS (Vect< real_t > &b) |
Set right-hand side vector. | |
void | setBC (Vect< real_t > &u) |
Set vector containing boundary condition to enforce. | |
void | setNewmarkParameters (real_t beta, real_t gamma) |
Define parameters for the Newmark 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 | setVerbose (int v=0) |
Set verbosity parameter: More... | |
real_t | runOneTimeStep () |
Run one time step. More... | |
void | run (bool opt=false) |
Run the time stepping procedure. More... | |
void | Assembly (const Element &el, real_t *b, real_t *A0, real_t *A1, real_t *A2=NULL) |
Assemble element arrays into global matrix and right-hand side. More... | |
void | SAssembly (const Side &sd, real_t *b, real_t *A=NULL) |
Assemble side arrays into global matrix and right-hand side. More... | |
LinearSolver< real_t > & | getLSolver () |
Return LinearSolver instance. | |
Detailed Description
To solve time stepping problems, i.e. systems of linear ordinary differential equations of the form [A2]{y"} + [A1]{y'} + [A0]{y} = {b}.
Features:
- The system may be first or second order (first and/or second order time derivatives
-
The following time integration schemes can be used:
-
For first order systems: The following schemes are implemented Forward Euler (value: FORWARD_EULER)
Backward Euler (value: BACKWARD_EULER)
Crank-Nicolson (value: CRANK_NICOLSON)
Heun (value: HEUN)
2nd Order Adams-Bashforth (value: AB2)
4-th order Runge-Kutta (value: RK4)
2nd order Backward Differentiation Formula (value: BDF2) - For second order systems: The following schemes are implemented Newmark (value: NEWMARK)
-
For first order systems: The following schemes are implemented Forward Euler (value: FORWARD_EULER)
Constructor & Destructor Documentation
TimeStepping | ( | TimeScheme | s, |
real_t | time_step = theTimeStep , |
||
real_t | final_time = theFinalTime |
||
) |
Constructor using time discretization data.
- Parameters
-
[in] s Choice of the scheme: To be chosen in the enumerated variable TimeScheme (see the presentation of the class) [in] time_step Value 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_time Value of the final time (time starts at 0). The default value for this parameter is the value given by the global variable theFinalTime
Member Function Documentation
void set | ( | TimeScheme | s, |
real_t | time_step = theTimeStep , |
||
real_t | final_time = theFinalTime |
||
) |
Define data of the differential equation or system.
- Parameters
-
[in] s Choice of the scheme: To be chosen in the enumerated variable TimeScheme (see the presentation of the class) [in] time_step Value 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_time Value of the final time (time starts at 0). The default value for this parameter is the value given by the global variable theFinalTime
Define partial differential equation to solve.
The used equation class must have been constructed using the Mesh instance
- Parameters
-
[in] eq Reference to equation instance
Set intermediate right-hand side vector for the Runge-Kutta method.
- Parameters
-
[in] f Vector containing the RHS
Set intermediate right-hand side vector for the TVD Runge-Kutta 3 method.
- Parameters
-
[in] f Vector containing the RHS
Set initial condition for the system of differential equations.
- Parameters
-
[in] u Vector containing initial condition for the unknown
- Remarks
- If a second-order differential equation is to be solved, use the the same function with two initial vectors (one for the unknown, the second for its time derivative)
Set initial condition for a system of differential equations.
- Parameters
-
[in] u Vector containing initial condition for the unknown [in] v Vector containing initial condition for the time derivative of the unknown
- Note
- This function can be used to provide solution at previous time step if a restarting procedure is used.
- This member function is to be used only in the case of a second order system
Set initial RHS for a system of differential equations when the used scheme requires it.
Giving the right-hand side at initial time is somtimes required for high order methods like Runge-Kutta
- Parameters
-
[in] f Vector containing right-hand side at initial time. This vector is helpful for high order methods
- Note
- This function can be used to provide solution at previous time step if a restarting procedure is used.
Define parameters for the Newmark scheme.
- Parameters
-
[in] beta Parameter beta [Default: 0.25
][in] gamma Parameter gamma [Default: 0.5
]
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
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
void setLinearSolver | ( | Iteration | s = DIRECT_SOLVER , |
Preconditioner | p = DIAG_PREC |
||
) |
Set linear solver data.
- Parameters
-
[in] s Solver 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] p Preconditioner 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
void setVerbose | ( | int | v = 0 | ) |
Set verbosity parameter:
- = 0, No output
- = 1, Print step label and time value
- = 2, Print step label, time value, time step and integration scheme
real_t runOneTimeStep | ( | ) |
Run one time step.
- Returns
- Value of new time step if this one is updated
void run | ( | bool | opt = false | ) |
Run the time stepping procedure.
- Parameters
-
[in] opt Flag 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
Assemble element arrays into global matrix and right-hand side.
This member function is to be called from finite element equation classes
- Parameters
-
[in] el Reference to Element class [in] b Pointer to element right-hand side [in] A0 Pointer to matrix of 0-th order term (involving no time derivative) [in] A1 Pointer to matrix of first order term (involving time first derivative) [in] A2 Pointer to matrix of second order term (involving time second derivative) [Default: NULL
]
Assemble side arrays into global matrix and right-hand side.
This member function is to be called from finite element equation classes
- Parameters
-
[in] sd Reference to Side class [in] b Pointer to side right-hand side [in] A Pointer to matrix [Default: NULL
]