MFEM v4.8.0
Finite element discretization library
Loading...
Searching...
No Matches
mfem::ODESolver Class Referenceabstract

Abstract class for solving systems of ODEs: dx/dt = f(x,t) More...

#include <ode.hpp>

Inheritance diagram for mfem::ODESolver:
[legend]
Collaboration diagram for mfem::ODESolver:
[legend]

Public Member Functions

 ODESolver ()
 
virtual void Init (TimeDependentOperator &f_)
 Associate a TimeDependentOperator with the ODE solver.
 
virtual void Step (Vector &x, real_t &t, real_t &dt)=0
 Perform a time step from time t [in] to time t [out] based on the requested step size dt [in].
 
virtual void Run (Vector &x, real_t &t, real_t &dt, real_t tf)
 Perform time integration from time t [in] to time tf [in].
 
virtual int GetStateSize ()
 Returns how many State vectors the ODE requires.
 
virtual ~ODESolver ()
 

Static Public Member Functions

static MFEM_EXPORT std::unique_ptr< ODESolverSelect (const int ode_solver_type)
 
static MFEM_EXPORT std::unique_ptr< ODESolverSelectExplicit (const int ode_solver_type)
 
static MFEM_EXPORT std::unique_ptr< ODESolverSelectImplicit (const int ode_solver_type)
 

Static Public Attributes

static MFEM_EXPORT std::string ExplicitTypes
 
static MFEM_EXPORT std::string ImplicitTypes
 
static MFEM_EXPORT std::string Types
 

Protected Attributes

TimeDependentOperatorf
 Pointer to the associated TimeDependentOperator.
 
MemoryType mem_type
 

Detailed Description

Abstract class for solving systems of ODEs: dx/dt = f(x,t)

Definition at line 109 of file ode.hpp.

Constructor & Destructor Documentation

◆ ODESolver()

mfem::ODESolver::ODESolver ( )
inline

Definition at line 117 of file ode.hpp.

◆ ~ODESolver()

virtual mfem::ODESolver::~ODESolver ( )
inlinevirtual

Definition at line 206 of file ode.hpp.

Member Function Documentation

◆ GetStateSize()

virtual int mfem::ODESolver::GetStateSize ( )
inlinevirtual

Returns how many State vectors the ODE requires.

Reimplemented in mfem::ODESolverWithStates.

Definition at line 182 of file ode.hpp.

◆ Init()

void mfem::ODESolver::Init ( TimeDependentOperator & f_)
virtual

◆ Run()

virtual void mfem::ODESolver::Run ( Vector & x,
real_t & t,
real_t & dt,
real_t tf )
inlinevirtual

Perform time integration from time t [in] to time tf [in].

Parameters
[in,out]xApproximate solution.
[in,out]tTime associated with the approximate solution x.
[in,out]dtTime step size.
[in]tfRequested final time.

The default implementation makes consecutive calls to Step() until reaching tf. The following rules describe the common behavior of the method:

  • The input x [in] is the approximate solution for the input time t [in].
  • The input dt [in] is the initial time step size.
  • The output dt [out] is the last time step taken by the method which may be smaller or larger than the input dt [in] value, e.g. because of time step control.
  • The output value of t [out] is not smaller than tf [in].

Reimplemented in mfem::PetscODESolver.

Definition at line 176 of file ode.hpp.

◆ Select()

std::unique_ptr< ODESolver > mfem::ODESolver::Select ( const int ode_solver_type)
static

Function for selecting the desired ODESolver (Explicit and Implicit) Returns an ODESolver pointer based on an type Caller gets ownership of the object and is responsible for its deletion

Definition at line 34 of file ode.cpp.

◆ SelectExplicit()

std::unique_ptr< ODESolver > mfem::ODESolver::SelectExplicit ( const int ode_solver_type)
static

Function for selecting the desired Explicit ODESolver Returns an ODESolver pointer based on an type Caller gets ownership of the object and is responsible for its deletion

Definition at line 46 of file ode.cpp.

◆ SelectImplicit()

std::unique_ptr< ODESolver > mfem::ODESolver::SelectImplicit ( const int ode_solver_type)
static

Function for selecting the desired Implicit ODESolver Returns an ODESolver pointer based on an type Caller gets ownership of the object and is responsible for its deletion

Definition at line 70 of file ode.cpp.

◆ Step()

virtual void mfem::ODESolver::Step ( Vector & x,
real_t & t,
real_t & dt )
pure virtual

Perform a time step from time t [in] to time t [out] based on the requested step size dt [in].

Parameters
[in,out]xApproximate solution.
[in,out]tTime associated with the approximate solution x.
[in,out]dtTime step size.

The following rules describe the common behavior of the method:

  • The input x [in] is the approximate solution for the input time t [in].
  • The input dt [in] is the desired time step size, defining the desired target time: t [target] = t [in] + dt [in].
  • The output x [out] is the approximate solution for the output time t [out].
  • The output dt [out] is the last time step taken by the method which may be smaller or larger than the input dt [in] value, e.g. because of time step control.
  • The method may perform more than one time step internally; in this case dt [out] is the last internal time step size.
  • The output value of t [out] may be smaller or larger than t [target], however, it is not smaller than t [in] + dt [out], if at least one internal time step was performed.
  • The value x [out] may be obtained by interpolation using internally stored data.
  • In some cases, the contents of x [in] may not be used, e.g. when x [out] from a previous Step() call was obtained by interpolation.
  • In consecutive calls to this method, the output t [out] of one Step() call has to be the same as the input t [in] to the next Step() call.
  • If the previous rule has to be broken, e.g. to restart a time stepping sequence, then the ODE solver must be re-initialized by calling Init() between the two Step() calls.

Implemented in mfem::AdamsBashforthSolver, mfem::AdamsMoultonSolver, mfem::ARKStepSolver, mfem::BackwardEulerSolver, mfem::ESDIRK32Solver, mfem::ESDIRK33Solver, mfem::ExplicitRKSolver, mfem::ForwardEulerSolver, mfem::GeneralizedAlphaSolver, mfem::ImplicitMidpointSolver, mfem::PetscODESolver, mfem::RK2Solver, mfem::RK3SSPSolver, mfem::RK4Solver, mfem::SDIRK23Solver, mfem::SDIRK33Solver, mfem::SDIRK34Solver, and mfem::TrapezoidalRuleSolver.

Member Data Documentation

◆ ExplicitTypes

std::string mfem::ODESolver::ExplicitTypes
static
Initial value:
=
"\n\tExplicit solver: \n\t"
" RK : 1 - Forward Euler, 2 - RK2(0.5), 3 - RK3 SSP, 4 - RK4, 6 - RK6,\n\t"
" AB : 11 - AB1, 12 - AB2, 13 - AB3, 14 - AB4, 15 - AB5\n"

Definition at line 185 of file ode.hpp.

◆ f

TimeDependentOperator* mfem::ODESolver::f
protected

Pointer to the associated TimeDependentOperator.

Definition at line 113 of file ode.hpp.

◆ ImplicitTypes

std::string mfem::ODESolver::ImplicitTypes
static
Initial value:
=
"\n\tImplicit solver: \n\t"
" (L-Stab): 21 - Backward Euler, 22 - SDIRK23(2), 23 - SDIRK33,\n\t"
" (A-Stab): 32 - Implicit Midpoint, 33 - SDIRK23, 34 - SDIRK34,\n\t"
" GA : 40 -- 50 - Generalized-alpha,\n\t"
" AM : 51 - AM1, 52 - AM2, 53 - AM3, 54 - AM4\n"

Definition at line 186 of file ode.hpp.

◆ mem_type

MemoryType mfem::ODESolver::mem_type
protected

Definition at line 114 of file ode.hpp.

◆ Types

std::string mfem::ODESolver::Types
static
Initial value:
static MFEM_EXPORT std::string ImplicitTypes
Definition ode.hpp:186
static MFEM_EXPORT std::string ExplicitTypes
Definition ode.hpp:185

Definition at line 187 of file ode.hpp.


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