https://mooseframework.inl.gov
LStableDirk3.h
Go to the documentation of this file.
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 
41 {
42 public:
44 
46 
47  virtual int order() override { return 3; }
48  virtual void computeTimeDerivatives() override;
49  virtual void computeADTimeDerivatives(ADReal & ad_u_dot,
50  const dof_id_type & dof,
51  ADReal & ad_u_dotdot) const override;
52  virtual void solve() override;
53  virtual void postResidual(NumericVector<Number> & residual) override;
54  virtual bool overridesSolve() const override { return true; }
55 
56 protected:
60  template <typename T, typename T2>
61  void computeTimeDerivativeHelper(T & u_dot, const T2 & u_old) const;
62 
63  // Indicates the current stage.
64  unsigned int _stage;
65 
66  // Store pointers to the various stage residuals
68 
69  // The parameter of the method, set at construction time and cannot be changed.
70  const Real _gamma; // 0.4358665215084589
71 
72  // Butcher tableau "C" parameters derived from _gamma
73  // 0.4358665215084589, 0.7179332607542295, 1.0000000000000000
74  Real _c[3];
75 
76  // Butcher tableau "A" values derived from _gamma. We only use the
77  // lower triangle of this.
78  // 0.4358665215084589
79  // 0.2820667392457705, 0.4358665215084589
80  // 1.2084966491760099, -0.6443631706844688, 0.4358665215084589
81  Real _a[3][3];
82 };
83 
84 template <typename T, typename T2>
85 void
86 LStableDirk3::computeTimeDerivativeHelper(T & u_dot, const T2 & u_old) const
87 {
88  u_dot -= u_old;
89  u_dot *= 1. / _dt;
90 }
virtual bool overridesSolve() const override
Definition: LStableDirk3.h:54
virtual void postResidual(NumericVector< Number > &residual) override
Callback to the NonLinearTimeIntegratorInterface called immediately after the residuals are computed ...
Definition: LStableDirk3.C:130
NumericVector< Number > * _stage_residuals[3]
Definition: LStableDirk3.h:67
virtual int order() override
Definition: LStableDirk3.h:47
virtual void solve() override
Solves the time step and sets the number of nonlinear and linear iterations.
Definition: LStableDirk3.C:83
Third order diagonally implicit Runge Kutta method (Dirk) with three stages.
Definition: LStableDirk3.h:40
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
static InputParameters validParams()
Definition: LStableDirk3.C:19
DualNumber< Real, DNDerivativeType, true > ADReal
Definition: ADRealForward.h:47
Real & _dt
The current time step size.
void computeTimeDerivativeHelper(T &u_dot, const T2 &u_old) const
Helper function that actually does the math for computing the time derivative.
Definition: LStableDirk3.h:86
Real _a[3][3]
Definition: LStableDirk3.h:81
unsigned int _stage
Definition: LStableDirk3.h:64
virtual void computeTimeDerivatives() override
Computes the time derivative and the Jacobian of the time derivative.
Definition: LStableDirk3.C:58
virtual void computeADTimeDerivatives(ADReal &ad_u_dot, const dof_id_type &dof, ADReal &ad_u_dotdot) const override
method for computing local automatic differentiation time derivatives
Definition: LStableDirk3.C:75
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Base class for time integrators.
const Real _gamma
Definition: LStableDirk3.h:70
LStableDirk3(const InputParameters &parameters)
Definition: LStableDirk3.C:27
const InputParameters & parameters() const
Get the parameters of the object.
Real _c[3]
Definition: LStableDirk3.h:74
uint8_t dof_id_type