Line data    Source code

       1            0 : // Distributed under the MIT License.
       2              : // See LICENSE.txt for details.
       3              : 
       4              : #pragma once
       5              : 
       6              : #include <vector>
       7              : 
       8              : #include "Time/TimeSteppers/ImexTimeStepper.hpp"
       9              : #include "Time/TimeSteppers/RungeKutta.hpp"
      10              : #include "Utilities/Gsl.hpp"
      11              : 
      12              : /// \cond
      13              : class TimeDelta;
      14              : namespace TimeSteppers {
      15              : template <typename T>
      16              : class ConstUntypedHistory;
      17              : }  // namespace TimeSteppers
      18              : /// \endcond
      19              : 
      20              : namespace TimeSteppers {
      21              : /*!
      22              :  * \ingroup TimeSteppersGroup
      23              :  * Intermediate base class implementing a generic IMEX Runge-Kutta
      24              :  * scheme.
      25              :  *
      26              :  * Implements most of the virtual methods of ImexTimeStepper for a
      27              :  * generic Runge-Kutta method.  Derived classes must implement the
      28              :  * requirements of the `RungeKutta` base class, as well as
      29              :  * `imex_order()`, `implicit_stage_order()`, and
      30              :  * `implicit_butcher_tableau()`.
      31              :  */
      32            1 : class ImexRungeKutta : public virtual RungeKutta,
      33              :                        public virtual ImexTimeStepper {
      34              :  public:
      35              :   /// Implicit part of the Butcher tableau.  Most parts of the tableau
      36              :   /// must be the same as the explicit part, and so are omitted.
      37            1 :   struct ImplicitButcherTableau {
      38              :     /*!
      39              :      * The coefficient matrix of the substeps.  We can only reasonably
      40              :      * support EDIRK methods (including special cases such as ESDIRK
      41              :      * and QESDIRK, see \cite Kennedy2016), so the tableau must be
      42              :      * lower-triangular with an empty first row.  As with the explicit
      43              :      * tableau, the initial blank row should be omitted.  For a
      44              :      * stiffly-accurate method, the final row must be the same as the
      45              :      * result coefficients in the explicit tableau.
      46              :      *
      47              :      * More general DIRK methods can be implemented inefficiently by
      48              :      * adding an unused initial substep to convert them to EDIRK form.
      49              :      */
      50            1 :     std::vector<std::vector<double>> substep_coefficients;
      51              :   };
      52              : 
      53              :   /*!
      54              :    * Smallest order of the intermediate result at any substep.  For
      55              :    * the methods supported by this class, this cannot exceed 2.
      56              :    */
      57            1 :   virtual size_t implicit_stage_order() const = 0;
      58              : 
      59            0 :   virtual const ImplicitButcherTableau& implicit_butcher_tableau() const = 0;
      60              : 
      61              :  private:
      62              :   template <typename T>
      63            0 :   void add_inhomogeneous_implicit_terms_impl(
      64              :       gsl::not_null<T*> u, const ConstUntypedHistory<T>& implicit_history,
      65              :       const TimeDelta& time_step) const;
      66              : 
      67              :   template <typename T>
      68            0 :   double implicit_weight_impl(const ConstUntypedHistory<T>& implicit_history,
      69              :                               const TimeDelta& time_step) const;
      70              : 
      71              :   IMEX_TIME_STEPPER_DECLARE_OVERLOADS
      72              : };
      73              : }  // namespace TimeSteppers