Line data    Source code

       1             : //* This file is part of the MOOSE framework
       2             : //* https://mooseframework.inl.gov
       3             : //*
       4             : //* All rights reserved, see COPYRIGHT for full restrictions
       5             : //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
       6             : //*
       7             : //* Licensed under LGPL 2.1, please see LICENSE for details
       8             : //* https://www.gnu.org/licenses/lgpl-2.1.html
       9             : 
      10             : #pragma once
      11             : 
      12             : #include "TimeIntegrator.h"
      13             : 
      14             : /**
      15             :  * Base class for three different explicit second-order Runge-Kutta
      16             :  * time integration methods:
      17             :  * - Explicit midpoint method (ExplicitMidpoint.C)
      18             :  * - Heun's method (Heun.C)
      19             :  * - Ralston's method (Ralston.C)
      20             :  *
      21             :  * Each of these methods is characterized by the following generic
      22             :  * Butcher tableau:
      23             :  * 0   | 0
      24             :  * a   | a   0
      25             :  * ---------------------
      26             :  *     | b1  b2
      27             :  *
      28             :  * where a, b1, and b2 are constants that define the different
      29             :  * methods.
      30             :  *
      31             :  * The stability function for each of these methods is:
      32             :  * R(z) = z*(z + 2)/2 + 1
      33             :  *
      34             :  * The methods are all explicit, so they are neither A-stable nor L-stable,
      35             :  * lim R(z), z->oo = oo
      36             :  *
      37             :  * See also:
      38             :  *   https://en.wikipedia.org/wiki/Runge%E2%80%93Kutta_methods
      39             :  *   https://en.wikipedia.org/wiki/Midpoint_method
      40             :  *   https://en.wikipedia.org/wiki/Heun%27s_method
      41             :  *
      42             :  * All three methods require two mass matrix (linear) system solves
      43             :  * per timestep. Although strictly speaking these are "two stage"
      44             :  * methods, we treat the "update" step as a third stage, since in
      45             :  * finite element analysis the update step requires a mass matrix
      46             :  * solve.
      47             :  *
      48             :  * IMPORTANT: To use the explicit TimeIntegrators derived from this
      49             :  * method, you must generally add "implicit=false" to the Kernels,
      50             :  * Materials, etc. used in your simulation, so that MOOSE evaluates
      51             :  * them correctly!  An important exception are TimeDerivative kernels,
      52             :  * which should never be marked "implicit=false".
      53             :  */
      54             : class ExplicitRK2 : public TimeIntegrator
      55             : {
      56             : public:
      57             :   static InputParameters validParams();
      58             : 
      59             :   ExplicitRK2(const InputParameters & parameters);
      60             : 
      61             :   virtual void preSolve() override;
      62           0 :   virtual int order() override { return 2; }
      63             : 
      64             :   virtual void computeTimeDerivatives() override;
      65             :   void computeADTimeDerivatives(ADReal & ad_u_dot,
      66             :                                 const dof_id_type & dof,
      67             :                                 ADReal & ad_u_dotdot) const override;
      68             :   virtual void solve() override;
      69             :   virtual void postResidual(NumericVector<Number> & residual) override;
      70        4720 :   virtual bool overridesSolve() const override { return true; }
      71             : 
      72             : protected:
      73             :   /**
      74             :    * Helper function that actually does the math for computing the time derivative
      75             :    */
      76             :   template <typename T, typename T2, typename T3>
      77             :   void computeTimeDerivativeHelper(T & u_dot, const T2 & u_old, const T3 & u_older) const;
      78             : 
      79             :   unsigned int _stage;
      80             : 
      81             :   /// Buffer to store non-time residual from the first stage.
      82             :   NumericVector<Number> * _residual_old;
      83             : 
      84             :   /// The older solution
      85             :   const NumericVector<Number> & _solution_older;
      86             : 
      87             :   /// The method coefficients.  See the table above for description.
      88             :   /// These are pure virtual in the base class, and must be overridden
      89             :   /// in subclasses to implement different schemes.
      90             :   virtual Real a() const = 0;
      91             :   virtual Real b1() const = 0;
      92             :   virtual Real b2() const = 0;
      93             : };
      94             : 
      95             : template <typename T, typename T2, typename T3>
      96             : void
      97        7944 : ExplicitRK2::computeTimeDerivativeHelper(T & u_dot, const T2 & u_old, const T3 & u_older) const
      98             : {
      99        7944 :   if (_stage < 3)
     100        2832 :     u_dot -= u_old;
     101             :   else
     102        5112 :     u_dot -= u_older;
     103             : 
     104        7944 :   u_dot *= 1. / _dt;
     105        7944 : }