MFEM
v4.0
Finite element discretization library
|
Wrapper for SUNDIALS' CVODE library – Multi-step time integration. More...
#include <sundials.hpp>
Public Member Functions | |
CVODESolver (int lmm, int iter) | |
Construct a serial CVODESolver, a wrapper for SUNDIALS' CVODE solver. More... | |
CVODESolver (MPI_Comm comm, int lmm, int iter) | |
Construct a parallel CVODESolver, a wrapper for SUNDIALS' CVODE solver. More... | |
void | SetSStolerances (double reltol, double abstol) |
Set the scalar relative and scalar absolute tolerances. More... | |
void | SetLinearSolver (SundialsODELinearSolver &ls_spec) |
void | SetStepMode (int itask) |
CVode supports two modes, specified by itask: CV_NORMAL (default) and CV_ONE_STEP. More... | |
void | SetMaxOrder (int max_order) |
Set the maximum order of the linear multistep method. More... | |
void | SetMaxStep (double dt_max) |
Set the maximum time step of the linear multistep 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 CVODE to integrate over [t, t + dt], with the specified step mode. More... | |
void | PrintInfo () const |
Print CVODE statistics. More... | |
virtual | ~CVODESolver () |
Destroy the associated CVODE 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... | |
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... | |
Protected Attributes inherited from mfem::ODESolver | |
TimeDependentOperator * | f |
Pointer to the associated TimeDependentOperator. More... | |
MemoryType | mem_type |
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... | |
Static Protected Attributes inherited from mfem::SundialsSolver | |
static const double | default_rel_tol = 1e-4 |
static const double | default_abs_tol = 1e-9 |
Wrapper for SUNDIALS' CVODE library – Multi-step time integration.
Definition at line 133 of file sundials.hpp.
mfem::CVODESolver::CVODESolver | ( | int | lmm, |
int | iter | ||
) |
Construct a serial CVODESolver, a wrapper for SUNDIALS' CVODE solver.
[in] | lmm | Specifies the linear multistep method, the options are CV_ADAMS (explicit methods) or CV_BDF (implicit methods). |
[in] | iter | Specifies type of nonlinear solver iteration, the options are CV_FUNCTIONAL (usually with CV_ADAMS) or CV_NEWTON (usually with CV_BDF). For parameter desciption, see the CVodeCreate documentation (cvode.h). |
Definition at line 175 of file sundials.cpp.
mfem::CVODESolver::CVODESolver | ( | MPI_Comm | comm, |
int | lmm, | ||
int | iter | ||
) |
Construct a parallel CVODESolver, a wrapper for SUNDIALS' CVODE solver.
[in] | comm | The MPI communicator used to partition the ODE system. |
[in] | lmm | Specifies the linear multistep method, the options are CV_ADAMS (explicit methods) or CV_BDF (implicit methods). |
[in] | iter | Specifies type of nonlinear solver iteration, the options are CV_FUNCTIONAL (usually with CV_ADAMS) or CV_NEWTON (usually with CV_BDF). For parameter desciption, see the CVodeCreate documentation (cvode.h). |
Definition at line 195 of file sundials.cpp.
|
virtual |
Destroy the associated CVODE memory.
Definition at line 417 of file sundials.cpp.
|
virtual |
Set the ODE right-hand-side operator.
The start time of CVODE is initialized from the current time of f_.
Reimplemented from mfem::ODESolver.
Definition at line 292 of file sundials.cpp.
void mfem::CVODESolver::PrintInfo | ( | ) | const |
Print CVODE statistics.
Definition at line 400 of file sundials.cpp.
void mfem::CVODESolver::SetLinearSolver | ( | SundialsODELinearSolver & | ls_spec | ) |
Set a custom Jacobian system solver for the CV_NEWTON option usually used with implicit CV_BDF.
Definition at line 233 of file sundials.cpp.
void mfem::CVODESolver::SetMaxOrder | ( | int | max_order | ) |
Set the maximum order of the linear multistep method.
The default is 12 (CV_ADAMS) or 5 (CV_BDF). CVODE uses adaptive-order integration, based on the local truncation error. Use this if you know a priori that your system is such that higher order integration formulas are unstable.
Definition at line 259 of file sundials.cpp.
|
inline |
Set the maximum time step of the linear multistep method.
Definition at line 183 of file sundials.hpp.
void mfem::CVODESolver::SetSStolerances | ( | double | reltol, |
double | abstol | ||
) |
Set the scalar relative and scalar absolute tolerances.
Definition at line 224 of file sundials.cpp.
void mfem::CVODESolver::SetStepMode | ( | int | itask | ) |
CVode supports two modes, specified by itask: CV_NORMAL (default) and CV_ONE_STEP.
In the CV_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 CV_ONE_STEP mode, it takes one internal step and returns.
Definition at line 254 of file sundials.cpp.
|
virtual |
Use CVODE to integrate over [t, t + dt], with the specified step mode.
Calls CVode(), which is the main driver of the CVODE package.
[in,out] | x | Solution vector to advance. On input/output x=x(t) for t corresponding to the input/output value of t, respectively. |
[in,out] | t | Input: the starting time value. Output: the time value of the solution output, as returned by CVode(). |
[in,out] | dt | Input: desired time step. Output: the last incremental time step used. |
Implements mfem::ODESolver.
Definition at line 356 of file sundials.cpp.