https://mooseframework.inl.gov
ExplicitTVDRK2.C
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 #include "ExplicitTVDRK2.h"
11 #include "NonlinearSystemBase.h"
12 #include "FEProblem.h"
13 #include "PetscSupport.h"
14 
16 
19 {
21  params.addClassDescription("Explicit TVD (total-variation-diminishing) second-order Runge-Kutta "
22  "time integration method.");
23  return params;
24 }
25 
27  : TimeIntegrator(parameters),
28  _stage(1),
29  _residual_old(addVector("residual_old", false, libMesh::GHOSTED)),
30  _solution_older(_sys.solutionState(2))
31 {
32  mooseInfo("ExplicitTVDRK2 and other multistage TimeIntegrators are known not to work with "
33  "Materials/AuxKernels that accumulate 'state' and should be used with caution.");
34 }
35 
36 void
38 {
39  if (_dt == _dt_old)
41  else
43 }
44 
45 void
47 {
48  // Since advanceState() is called in between stages 2 and 3, this
49  // changes the meaning of "_solution_old". In the second stage,
50  // "_solution_older" is actually the original _solution_old.
51  if (!_sys.solutionUDot())
52  mooseError("ExplicitTVDRK2: Time derivative of solution (`u_dot`) is not stored. Please set "
53  "uDotRequested() to true in FEProblemBase befor requesting `u_dot`.");
54 
56  u_dot = *_solution;
58 
60  u_dot.close();
61 }
62 
63 void
65  const dof_id_type & dof,
66  ADReal & /*ad_u_dotdot*/) const
67 {
69 }
70 
71 void
73 {
74  Real time_new = _fe_problem.time();
75  Real time_old = _fe_problem.timeOld();
76  Real time_stage2 = time_old + _dt;
77 
78  // Reset numbers of iterations
81 
82  // There is no work to do for the first stage (Y_1 = y_n). The
83  // first solve therefore happens in the second stage. Note that the
84  // non-time Kernels (which should be marked implicit=false) are
85  // evaluated at the old solution during this stage.
87  _console << "1st solve" << std::endl;
88  _stage = 2;
89  _fe_problem.timeOld() = time_old;
90  _fe_problem.time() = time_stage2;
91  _nl->system().solve();
94 
95  // Abort time step immediately on stage failure - see TimeIntegrator doc page
97  return;
98 
99  // Advance solutions old->older, current->old. Also moves Material
100  // properties and other associated state forward in time.
102 
103  // The "update" stage (which we call stage 3) requires an additional
104  // solve with the mass matrix.
106  _console << "2nd solve" << std::endl;
107  _stage = 3;
108  _fe_problem.timeOld() = time_stage2;
109  _fe_problem.time() = time_new;
110  _nl->system().solve();
113 
114  // Error if solve didn't work, since we can't undo the advanceState
115  if (!_fe_problem.converged(_nl->number()))
116  mooseError("Aborting as ",
117  type(),
118  " has undefined behavior if it attempts to re-do a timestep after the first stage.");
119 
120  // Reset time at beginning of step to its original value
121  _fe_problem.timeOld() = time_old;
122 }
123 
124 void
126 {
127  if (_stage == 1)
128  {
129  // If postResidual() is called before solve(), _stage==1 and we don't
130  // need to do anything.
131  }
132  else if (_stage == 2)
133  {
134  // In the standard RK notation, the stage 2 residual is given by:
135  //
136  // R := M*(Y_2 - y_n)/dt - f(t_n, Y_1) = 0
137  //
138  // where:
139  // .) M is the mass matrix.
140  // .) f(t_n, Y_1) is the residual we are currently computing,
141  // since this method is intended to be used with "implicit=false"
142  // kernels.
143  // .) M*(Y_2 - y_n)/dt corresponds to the residual of the time kernels.
144  // .) The minus signs are "baked in" to the non-time residuals, so
145  // they do not appear here.
146  // .) The current non-time residual is saved for the next stage.
148  _residual_old->close();
149 
150  residual.add(1.0, *_Re_time);
151  residual.add(1.0, *_residual_old);
152  residual.close();
153  }
154  else if (_stage == 3)
155  {
156  // In the standard RK notation, the update step residual is given by:
157  //
158  // R := M*(2*y_{n+1} - Y_2 - y_n)/(2*dt) - (1/2)*f(t_n+dt/2, Y_2) = 0
159  //
160  // where:
161  // .) M is the mass matrix.
162  // .) f(t_n+dt/2, Y_2) is the residual from stage 2.
163  // .) The minus sign is already baked in to the non-time
164  // residuals, so it does not appear here.
165  // .) Although this is an update step, we have to do a "solve"
166  // using the mass matrix.
167  residual.add(1.0, *_Re_time);
168  residual.add(0.5, *_Re_non_time);
169  residual.close();
170  }
171  else
172  mooseError(
173  "ExplicitTVDRK2::postResidual(): _stage = ", _stage, ", only _stage = 1-3 is allowed.");
174 }
virtual void computeTimeDerivatives() override
Computes the time derivative and the Jacobian of the time derivative.
virtual void initPetscOutputAndSomeSolverSettings()
Reinitialize PETSc output for proper linear/nonlinear iteration display.
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.
virtual Real & time() const
FEProblemBase & _fe_problem
Reference to the problem.
void mooseInfo(Args &&... args) const
unsigned int _stage
NonlinearSystemBase * _nl
Pointer to the nonlinear system, can happen that we dont have any.
SystemBase & _sys
Reference to the system this time integrator operates on.
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
virtual void preSolve() override
The following methods are specializations for using the libMesh::Parallel::packed_range_* routines fo...
DualNumber< Real, DNDerivativeType, true > ADReal
Definition: ADRealForward.h:46
Real & _dt
The current time step size.
Real & _dt_old
The previous time step size.
NumericVector< Number > * _Re_time
residual vector for time contributions
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
virtual void advanceState()
Advance all of the state holding vectors / datastructures so that we can move to the next timestep...
const NumericVector< Number > & _solution_older
The older solution.
GHOSTED
virtual NumericVector< Number > * solutionUDot()
Definition: SystemBase.h:252
registerMooseObject("MooseApp", ExplicitTVDRK2)
virtual bool converged(const unsigned int sys_num)
Eventually we want to convert this virtual over to taking a solver system number argument.
Definition: SubProblem.h:113
unsigned int _n_linear_iterations
Total number of linear iterations over all stages of the time step.
const std::string & type() const
Get the type of this class.
Definition: MooseBase.h:51
unsigned int number() const
Gets the number of this system.
Definition: SystemBase.C:1159
virtual void solve()
virtual void close()=0
const NumericVector< Number > *const & _solution
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Base class for time integrators.
unsigned int getNumLinearIterationsLastSolve() const
Gets the number of linear iterations in the most recent solve.
void setConstJacobian(bool state)
Set flag that Jacobian is constant (for optimization purposes)
ExplicitTVDRK2(const InputParameters &parameters)
void mooseError(Args &&... args) const
Emits an error prefixed with object name and type.
virtual void postResidual(NumericVector< Number > &residual) override
Callback to the NonLinearTimeIntegratorInterface called immediately after the residuals are computed ...
void addClassDescription(const std::string &doc_string)
This method adds a description of the class that will be displayed in the input file syntax dump...
virtual Real & timeOld() const
NumericVector< Number > * _residual_old
Buffer to store non-time residual from the first stage.
const NumericVector< Number > & _solution_old
unsigned int _n_nonlinear_iterations
Total number of nonlinear iterations over all stages of the time step.
const ConsoleStream _console
An instance of helper class to write streams to the Console objects.
unsigned int getNumNonlinearIterationsLastSolve() const
Gets the number of nonlinear iterations in the most recent solve.
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.
virtual void add(const numeric_index_type i, const Number value)=0
NumericVector< Number > * _Re_non_time
residual vector for non-time contributions
static InputParameters validParams()
static InputParameters validParams()
uint8_t dof_id_type
void computeDuDotDu()
Compute _du_dot_du.
virtual libMesh::System & system() override
Get the reference to the libMesh system.