MFEM  v3.3
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
Public Types | Public Member Functions | Protected Attributes | List of all members
mfem::ARKODESolver Class Reference

Wrapper for SUNDIALS' ARKODE library – Runge-Kutta time integration. More...

#include <sundials.hpp>

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

Public Types

enum  Type { EXPLICIT, IMPLICIT }
 Types of ARKODE solvers. More...
 

Public Member Functions

 ARKODESolver (Type type=EXPLICIT)
 Construct a serial ARKODESolver, a wrapper for SUNDIALS' ARKODE solver. More...
 
 ARKODESolver (MPI_Comm comm, Type type=EXPLICIT)
 Construct a parallel ARKODESolver, a wrapper for SUNDIALS' ARKODE solver. More...
 
void SetSStolerances (double reltol, double abstol)
 Specify the scalar relative and scalar absolute tolerances. More...
 
void SetLinearSolver (SundialsODELinearSolver &ls_spec)
 Set a custom Jacobian system solver for implicit methods. More...
 
void SetStepMode (int itask)
 ARKode supports two modes, specified by itask: ARK_NORMAL (default) and ARK_ONE_STEP. More...
 
void SetOrder (int order)
 Chooses integration order for all explicit / implicit / IMEX methods. More...
 
void SetIRKTableNum (int table_num)
 Choose a specific Butcher table for implicit RK method. More...
 
void SetERKTableNum (int table_num)
 Choose a specific Butcher table for explicit RK method. More...
 
void SetFixedStep (double dt)
 Use a fixed time step size, instead of performing any form of temporal adaptivity. More...
 
void SetMaxStep (double dt_max)
 Set the maximum time step of the Runge-Kutta method. More...
 
virtual void Init (TimeDependentOperator &f_)
 Set the ODE right-hand-side operator. More...
 
virtual void Step (Vector &x, double &t, double &dt)
 Use ARKODE to integrate over [t, t + dt], with the specified step mode. More...
 
void PrintInfo () const
 Print ARKODE statistics. More...
 
virtual ~ARKODESolver ()
 Destroy the associated ARKODE memory. More...
 
- Public Member Functions inherited from mfem::ODESolver
 ODESolver ()
 
virtual void Run (Vector &x, double &t, double &dt, double tf)
 Perform time integration from time t [in] to time tf [in]. More...
 
virtual ~ODESolver ()
 
- Public Member Functions inherited from mfem::SundialsSolver
void * SundialsMem () const
 Access the underlying SUNDIALS object. More...
 
int GetFlag () const
 Return the flag returned by the last call to a SUNDIALS function. More...
 

Protected Attributes

bool use_implicit
 
int irk_table
 
int erk_table
 
- Protected Attributes inherited from mfem::ODESolver
TimeDependentOperatorf
 Pointer to the associated TimeDependentOperator. More...
 
- Protected Attributes inherited from mfem::SundialsSolver
void * sundials_mem
 Pointer to the SUNDIALS mem object. More...
 
int flag
 Flag returned by the last call to SUNDIALS. More...
 
N_Vector y
 Auxiliary N_Vector. More...
 

Additional Inherited Members

- Protected Member Functions inherited from mfem::SundialsSolver
bool Parallel () const
 
bool Parallel () const
 
 SundialsSolver ()
 
 SundialsSolver (void *mem)
 
- Static Protected Member Functions inherited from mfem::SundialsSolver
static int ODEMult (realtype t, const N_Vector y, N_Vector ydot, void *td_oper)
 Callback function used in CVODESolver and ARKODESolver. More...
 
- Static Protected Attributes inherited from mfem::SundialsSolver
static const double default_rel_tol = 1e-4
 
static const double default_abs_tol = 1e-9
 

Detailed Description

Wrapper for SUNDIALS' ARKODE library – Runge-Kutta time integration.

Definition at line 220 of file sundials.hpp.

Member Enumeration Documentation

Types of ARKODE solvers.

Enumerator
EXPLICIT 
IMPLICIT 

Definition at line 228 of file sundials.hpp.

Constructor & Destructor Documentation

mfem::ARKODESolver::ARKODESolver ( Type  type = EXPLICIT)

Construct a serial ARKODESolver, a wrapper for SUNDIALS' ARKODE solver.

Parameters
[in]typeSpecifies the Type of ARKODE solver to construct.

Definition at line 383 of file sundials.cpp.

mfem::ARKODESolver::ARKODESolver ( MPI_Comm  comm,
Type  type = EXPLICIT 
)

Construct a parallel ARKODESolver, a wrapper for SUNDIALS' ARKODE solver.

Parameters
[in]commThe MPI communicator used to partition the ODE system.
[in]typeSpecifies the Type of ARKODE solver to construct.

Definition at line 403 of file sundials.cpp.

mfem::ARKODESolver::~ARKODESolver ( )
virtual

Destroy the associated ARKODE memory.

Definition at line 652 of file sundials.cpp.

Member Function Documentation

void mfem::ARKODESolver::Init ( TimeDependentOperator f_)
virtual

Set the ODE right-hand-side operator.

The start time of ARKODE is initialized from the current time of f_.

Note
This method calls ARKodeInit(). Some ARKODE parameters can be set (using the handle returned by SundialsMem()) only after this call.

Reimplemented from mfem::ODESolver.

Definition at line 513 of file sundials.cpp.

void mfem::ARKODESolver::PrintInfo ( ) const

Print ARKODE statistics.

Definition at line 636 of file sundials.cpp.

void mfem::ARKODESolver::SetERKTableNum ( int  table_num)

Choose a specific Butcher table for explicit RK method.

See the documentation for all possible options, stability regions, etc.

Definition at line 481 of file sundials.cpp.

void mfem::ARKODESolver::SetFixedStep ( double  dt)

Use a fixed time step size, instead of performing any form of temporal adaptivity.

Use of this function is not recommended, since there is no assurance of the validity of the computed solutions. It is primarily provided for code-to-code verification testing purposes.

Definition at line 487 of file sundials.cpp.

void mfem::ARKODESolver::SetIRKTableNum ( int  table_num)

Choose a specific Butcher table for implicit RK method.

See the documentation for all possible options, stability regions, etc. For example, table_num = ARK548L2SA_DIRK_8_4_5 is 8-stage 5th order.

Definition at line 475 of file sundials.cpp.

void mfem::ARKODESolver::SetLinearSolver ( SundialsODELinearSolver ls_spec)

Set a custom Jacobian system solver for implicit methods.

Definition at line 441 of file sundials.cpp.

void mfem::ARKODESolver::SetMaxStep ( double  dt_max)
inline

Set the maximum time step of the Runge-Kutta method.

Definition at line 276 of file sundials.hpp.

void mfem::ARKODESolver::SetOrder ( int  order)

Chooses integration order for all explicit / implicit / IMEX methods.

The default is 4, and the allowed ranges are: [2, 8] for explicit; [2, 5] for implicit; [3, 5] for IMEX.

Definition at line 467 of file sundials.cpp.

void mfem::ARKODESolver::SetSStolerances ( double  reltol,
double  abstol 
)

Specify the scalar relative and scalar absolute tolerances.

Definition at line 432 of file sundials.cpp.

void mfem::ARKODESolver::SetStepMode ( int  itask)

ARKode supports two modes, specified by itask: ARK_NORMAL (default) and ARK_ONE_STEP.

In the ARK_NORMAL mode, the solver steps until it reaches or passes tout = t + dt, where t and dt are specified in Step(), and then interpolates to obtain y(tout). In the ARK_ONE_STEP mode, it takes one internal step and returns.

Definition at line 462 of file sundials.cpp.

void mfem::ARKODESolver::Step ( Vector x,
double &  t,
double &  dt 
)
virtual

Use ARKODE to integrate over [t, t + dt], with the specified step mode.

Calls ARKode(), which is the main driver of the ARKODE package.

Parameters
[in,out]xSolution vector to advance. On input/output x=x(t) for t corresponding to the input/output value of t, respectively.
[in,out]tInput: the starting time value. Output: the time value of the solution output, as returned by CVode().
[in,out]dtInput: desired time step. Output: the last incremental time step used.

Implements mfem::ODESolver.

Definition at line 596 of file sundials.cpp.

Member Data Documentation

int mfem::ARKODESolver::erk_table
protected

Definition at line 224 of file sundials.hpp.

int mfem::ARKODESolver::irk_table
protected

Definition at line 224 of file sundials.hpp.

bool mfem::ARKODESolver::use_implicit
protected

Definition at line 223 of file sundials.hpp.


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