www.mooseframework.org
ExplicitTVDRK2.h
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
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 
46 {
47 public:
49 
51 
52  virtual void preSolve() override;
53  virtual int order() override { return 2; }
54 
55  virtual void computeTimeDerivatives() override;
56  void computeADTimeDerivatives(DualReal & ad_u_dot,
57  const dof_id_type & dof,
58  DualReal & ad_u_dotdot) const override;
59  virtual void solve() override;
60  virtual void postResidual(NumericVector<Number> & residual) override;
61 
62 protected:
66  template <typename T, typename T2, typename T3>
67  void computeTimeDerivativeHelper(T & u_dot, const T2 & u_old, const T3 & u_older) const;
68 
69  unsigned int _stage;
70 
73 
76 };
77 
78 template <typename T, typename T2, typename T3>
79 void
80 ExplicitTVDRK2::computeTimeDerivativeHelper(T & u_dot, const T2 & u_old, const T3 & u_older) const
81 {
82  if (_stage < 3)
83  {
84  u_dot -= u_old;
85  u_dot *= 1. / _dt;
86  }
87  else
88  {
89  u_dot *= 2.;
90  u_dot -= u_old;
91  u_dot -= u_older;
92  u_dot *= 0.5 / _dt;
93  }
94 }
virtual void computeTimeDerivatives() override
Computes the time derivative and the Jacobian of the time derivative.
void computeTimeDerivativeHelper(T &u_dot, const T2 &u_old, const T3 &u_older) const
Helper function that actually does the math for computing the time derivative.
unsigned int _stage
DualNumber< Real, DNDerivativeType, true > DualReal
Definition: DualReal.h:49
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
virtual int order() override
virtual void preSolve() override
const NumericVector< Number > & _solution_older
The older solution.
NumericVector< Number > & _residual_old
Buffer to store non-time residual from the first stage.
Base class for time integrators.
ExplicitTVDRK2(const InputParameters &parameters)
virtual void postResidual(NumericVector< Number > &residual) override
Callback to the TimeIntegrator called immediately after the residuals are computed in NonlinearSystem...
const InputParameters & parameters() const
Get the parameters of the object.
Explicit TVD (total-variation-diminishing) second-order Runge-Kutta time integration methods: ...
virtual void solve() override
Solves the time step and sets the number of nonlinear and linear iterations.
static InputParameters validParams()
void computeADTimeDerivatives(DualReal &ad_u_dot, const dof_id_type &dof, DualReal &ad_u_dotdot) const override
method for computing local automatic differentiation time derivatives
uint8_t dof_id_type