Loading [MathJax]/extensions/TeX/AMSmath.js
3.5-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
Public Member Functions | List of all members
Dumux::Experimental::MultiStageMethod< Scalar > Class Template Referenceabstract

Abstract interface for one-step multi-stage method parameters in Shu/Osher form. More...

#include <dumux/timestepping/multistagemethods.hh>

Inheritance diagram for Dumux::Experimental::MultiStageMethod< Scalar >:

Description

template<class Scalar>
class Dumux::Experimental::MultiStageMethod< Scalar >

Abstract interface for one-step multi-stage method parameters in Shu/Osher form.

forward declaration

This implementation is based on the Shu/Osher form from: Chi W. Shu and Stanley Osher. Efficient implementation of essentially non- oscillatory shock-capturing schemes. J. Comput. Phys., 77:439–471, 1988. https://doi.org/10.1016/0021-9991(88)90177-5. To this end Eq. (2.12) is extended for implicit schemes.

We consider the general PDE form

\begin{equation} \frac{\partial M(x)}{\partial t} - R(x, t) = 0, \end{equation}

where M(x) is the temporal operator/residual in terms of the solution x , and R(x) is the discrete spatial operator/residual. M(x) usually corresponds to the conserved quantity (e.g. mass), and R(x) contains the rest of the residual. We can then construct m -stage time discretization methods. Integrating from time t^n to t^{n+1} with time step size \Delta t^n, we solve:

\begin{aligned} x^{(0)} &= u^n\\ \sum_{k=0}^i \left[ \alpha_{ik} M\left(x^{(k)}, t^n + d_k\Delta t^n\right) + \beta_{ik}\Delta t^n R \left(x^{(k)}, t^n+d_k\Delta t^n \right)\right] &= 0 & \forall i \in \{1,\ldots,m\} \\ x^{n+1} &= x^{(m)} \end{aligned}

where x^{(k)} denotes the intermediate solution at stage k . Dependent on the number of stages m , and the coefficients \alpha, \beta, d, schemes with different properties and order of accuracy can be constructed.

Public Member Functions

virtual bool implicit () const =0
 
virtual std::size_t numStages () const =0
 
virtual Scalar temporalWeight (std::size_t i, std::size_t k) const =0
 weights of the temporal operator residual ( \alpha_{ik} ) More...
 
virtual Scalar spatialWeight (std::size_t i, std::size_t k) const =0
 weights of the spatial operator residual ( \beta_{ik} ) More...
 
virtual Scalar timeStepWeight (std::size_t k) const =0
 time step weights for each stage ( d_k ) More...
 
virtual std::string name () const =0
 
virtual ~MultiStageMethod ()=default
 

Constructor & Destructor Documentation

◆ ~MultiStageMethod()

template<class Scalar >
virtual Dumux::Experimental::MultiStageMethod< Scalar >::~MultiStageMethod ( )
virtualdefault

Member Function Documentation

◆ implicit()

template<class Scalar >
virtual bool Dumux::Experimental::MultiStageMethod< Scalar >::implicit ( ) const
pure virtual

◆ name()

template<class Scalar >
virtual std::string Dumux::Experimental::MultiStageMethod< Scalar >::name ( ) const
pure virtual

◆ numStages()

template<class Scalar >
virtual std::size_t Dumux::Experimental::MultiStageMethod< Scalar >::numStages ( ) const
pure virtual

◆ spatialWeight()

template<class Scalar >
virtual Scalar Dumux::Experimental::MultiStageMethod< Scalar >::spatialWeight ( std::size_t  i,
std::size_t  k 
) const
pure virtual

weights of the spatial operator residual ( \beta_{ik} )

Implemented in Dumux::Experimental::MultiStage::RungeKuttaExplicitFourthOrder< Scalar >, and Dumux::Experimental::MultiStage::Theta< Scalar >.

◆ temporalWeight()

template<class Scalar >
virtual Scalar Dumux::Experimental::MultiStageMethod< Scalar >::temporalWeight ( std::size_t  i,
std::size_t  k 
) const
pure virtual

weights of the temporal operator residual ( \alpha_{ik} )

Implemented in Dumux::Experimental::MultiStage::RungeKuttaExplicitFourthOrder< Scalar >, and Dumux::Experimental::MultiStage::Theta< Scalar >.

◆ timeStepWeight()

template<class Scalar >
virtual Scalar Dumux::Experimental::MultiStageMethod< Scalar >::timeStepWeight ( std::size_t  k) const
pure virtual

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