20#ifndef DUMUX_TIMESTEPPING_MULTISTAGE_METHODS_HH 
   21#define DUMUX_TIMESTEPPING_MULTISTAGE_METHODS_HH 
   90    virtual std::string 
name () 
const = 0;
 
 
  104template<
class Scalar>
 
  109    : paramAlpha_{{-1.0, 1.0}}
 
  110    , paramBeta_{{1.0-theta, theta}}
 
  111    , paramD_{{0.0, 1.0}}
 
 
  115    { 
return paramBeta_[1] > 0.0; }
 
 
  121    { 
return paramAlpha_[k]; }
 
 
  124    { 
return paramBeta_[k]; }
 
 
  127    { 
return paramD_[k]; }
 
 
  129    std::string 
name ()
 const override 
  130    { 
return "theta scheme"; }
 
 
  133    std::array<Scalar, 2> paramAlpha_;
 
  134    std::array<Scalar, 2> paramBeta_;
 
  135    std::array<Scalar, 2> paramD_;
 
 
  141template<
class Scalar>
 
  147    std::string 
name () const final
 
  148    { 
return "explicit Euler"; }
 
 
 
  154template<
class Scalar>
 
  160    std::string 
name () const final
 
  161    { 
return "implicit Euler"; }
 
 
 
  167template<
class Scalar>
 
  172    : paramAlpha_{{{-1.0, 1.0, 0.0, 0.0, 0.0},
 
  173                   {-1.0, 0.0, 1.0, 0.0, 0.0},
 
  174                   {-1.0, 0.0, 0.0, 1.0, 0.0},
 
  175                   {-1.0, 0.0, 0.0, 0.0, 1.0}}}
 
  176    , paramBeta_{{{0.5, 0.0, 0.0, 0.0, 0.0},
 
  177                  {0.0, 0.5, 0.0, 0.0, 0.0},
 
  178                  {0.0, 0.0, 1.0, 0.0, 0.0},
 
  179                  {1.0/6.0, 1.0/3.0, 1.0/3.0, 1.0/6.0, 0.0}}}
 
  180    , paramD_{{0.0, 0.5, 0.5, 1.0, 1.0}}
 
 
  190    { 
return paramAlpha_[i-1][k]; }
 
 
  193    { 
return paramBeta_[i-1][k]; }
 
 
  196    { 
return paramD_[k]; }
 
 
  198    std::string 
name () const final
 
  199    { 
return "explicit Runge-Kutta 4th order"; }
 
 
  202    std::array<std::array<Scalar, 5>, 4> paramAlpha_;
 
  203    std::array<std::array<Scalar, 5>, 4> paramBeta_;
 
  204    std::array<Scalar, 5> paramD_;
 
 
  211template<
class Scalar>
 
  217        constexpr Scalar alpha = []{
 
  219            Scalar alpha = 0.4358665215; 
 
  220            for (
int i = 0; i < 10; ++i)
 
  221                alpha = alpha - (alpha*(alpha*alpha-3.0*(alpha-0.5))-1.0/6.0)/(3.0*alpha*(alpha-2.0)+1.5);
 
  225        constexpr Scalar tau2 = (1.0+alpha)*0.5;
 
  226        constexpr Scalar b1 = -(6.0*alpha*alpha -16.0*alpha + 1.0)*0.25;
 
  227        constexpr Scalar b2 = (6*alpha*alpha - 20.0*alpha + 5.0)*0.25;
 
  229        paramD_ = {{0.0, alpha, tau2, 1.0}};
 
  231            {-1.0, 1.0, 0.0, 0.0},
 
  232            {-1.0, 0.0, 1.0, 0.0},
 
  233            {-1.0, 0.0, 0.0, 1.0}
 
  236            {0.0, alpha, 0.0, 0.0},
 
  237            {0.0, tau2-alpha, alpha, 0.0},
 
 
  249    { 
return paramAlpha_[i-1][k]; }
 
 
  252    { 
return paramBeta_[i-1][k]; }
 
 
  255    { 
return paramD_[k]; }
 
 
  257    std::string 
name () const final
 
  258    { 
return "diagonally implicit Runge-Kutta 3rd order (Alexander)"; }
 
 
  261    std::array<std::array<Scalar, 4>, 3> paramAlpha_;
 
  262    std::array<std::array<Scalar, 4>, 3> paramBeta_;
 
  263    std::array<Scalar, 4> paramD_;
 
 
 
bool implicit() const final
Definition multistagemethods.hh:242
DIRKThirdOrderAlexander()
Definition multistagemethods.hh:215
Scalar temporalWeight(std::size_t i, std::size_t k) const final
weights of the temporal operator residual ( )
Definition multistagemethods.hh:248
Scalar spatialWeight(std::size_t i, std::size_t k) const final
weights of the spatial operator residual ( )
Definition multistagemethods.hh:251
Scalar timeStepWeight(std::size_t k) const final
time step weights for each stage ( )
Definition multistagemethods.hh:254
std::string name() const final
Definition multistagemethods.hh:257
std::size_t numStages() const final
Definition multistagemethods.hh:245
ExplicitEuler()
Definition multistagemethods.hh:145
std::string name() const final
Definition multistagemethods.hh:147
ImplicitEuler()
Definition multistagemethods.hh:158
std::string name() const final
Definition multistagemethods.hh:160
RungeKuttaExplicitFourthOrder()
Definition multistagemethods.hh:171
Scalar temporalWeight(std::size_t i, std::size_t k) const final
weights of the temporal operator residual ( )
Definition multistagemethods.hh:189
std::string name() const final
Definition multistagemethods.hh:198
bool implicit() const final
Definition multistagemethods.hh:183
std::size_t numStages() const final
Definition multistagemethods.hh:186
Scalar spatialWeight(std::size_t i, std::size_t k) const final
weights of the spatial operator residual ( )
Definition multistagemethods.hh:192
Scalar timeStepWeight(std::size_t k) const final
time step weights for each stage ( )
Definition multistagemethods.hh:195
std::size_t numStages() const final
Definition multistagemethods.hh:117
bool implicit() const final
Definition multistagemethods.hh:114
std::string name() const override
Definition multistagemethods.hh:129
Theta(const Scalar theta)
Definition multistagemethods.hh:108
Scalar spatialWeight(std::size_t, std::size_t k) const final
weights of the spatial operator residual ( )
Definition multistagemethods.hh:123
Scalar temporalWeight(std::size_t, std::size_t k) const final
weights of the temporal operator residual ( )
Definition multistagemethods.hh:120
Scalar timeStepWeight(std::size_t k) const final
time step weights for each stage ( )
Definition multistagemethods.hh:126
Abstract interface for one-step multi-stage method parameters in Shu/Osher form.
Definition multistagemethods.hh:75
virtual std::string name() const =0
virtual std::size_t numStages() const =0
virtual ~MultiStageMethod()=default
virtual Scalar timeStepWeight(std::size_t k) const =0
time step weights for each stage ( )
virtual Scalar temporalWeight(std::size_t i, std::size_t k) const =0
weights of the temporal operator residual ( )
virtual bool implicit() const =0
virtual Scalar spatialWeight(std::size_t i, std::size_t k) const =0
weights of the spatial operator residual ( )
Multi-stage time stepping scheme implementations.
Definition multistagemethods.hh:96
Definition experimental/assembly/cclocalassembler.hh:36