www.mooseframework.org
LStableDirk4.C
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 #include "LStableDirk4.h"
11 #include "NonlinearSystemBase.h"
12 #include "FEProblem.h"
13 #include "PetscSupport.h"
14 
16 
17 template <>
20 {
22  return params;
23 }
24 
25 // Initialize static data
26 const Real LStableDirk4::_c[LStableDirk4::_n_stages] = {.25, 0., .5, 1., 1.};
27 
29  {.25, 0, 0, 0, 0},
30  {-.25, .25, 0, 0, 0},
31  {.125, .125, .25, 0, 0},
32  {-1.5, .75, 1.5, .25, 0},
33  {0, 1. / 6, 2. / 3, -1. / 12, .25}};
34 
36  : TimeIntegrator(parameters), _stage(1)
37 {
38  mooseInfo("LStableDirk4 and other multistage TimeIntegrators are known not to work with "
39  "Materials/AuxKernels that accumulate 'state' and should be used with caution.");
40 
41  // Name the stage residuals "residual_stage1", "residual_stage2", etc.
42  for (unsigned int stage = 0; stage < _n_stages; ++stage)
43  {
44  std::ostringstream oss;
45  oss << "residual_stage" << stage + 1;
46  _stage_residuals[stage] = &(_nl.addVector(oss.str(), false, GHOSTED));
47  }
48 }
49 
50 void
52 {
53  // We are multiplying by the method coefficients in postResidual(), so
54  // the time derivatives are of the same form at every stage although
55  // the current solution varies depending on the stage.
56  if (!_sys.solutionUDot())
57  mooseError("LStableDirk4: Time derivative of solution (`u_dot`) is not stored. Please set "
58  "uDotRequested() to true in FEProblemBase befor requesting `u_dot`.");
59 
60  NumericVector<Number> & u_dot = *_sys.solutionUDot();
61  u_dot = *_solution;
63  u_dot.close();
64  _du_dot_du = 1. / _dt;
65 }
66 
67 void
68 LStableDirk4::computeADTimeDerivatives(DualReal & ad_u_dot, const dof_id_type & dof) const
69 {
71 }
72 
73 void
75 {
76  // Time at end of step
77  Real time_old = _fe_problem.timeOld();
78 
79  // Reset iteration counts
82 
83  // A for-loop would increment _stage too far, so we use an extra
84  // loop counter.
85  for (unsigned int current_stage = 1; current_stage <= _n_stages; ++current_stage)
86  {
87  // Set the current stage value
88  _stage = current_stage;
89 
90  // This ensures that all the Output objects in the OutputWarehouse
91  // have had solveSetup() called, and sets the default solver
92  // parameters for PETSc.
94 
95  _console << "Stage " << _stage << "\n";
96 
97  // Set the time for this stage
98  _fe_problem.time() = time_old + _c[_stage - 1] * _dt;
99 
100  // Do the solve
102 
103  // Update the iteration counts
106 
107  // Abort time step immediately on stage failure - see TimeIntegrator doc page
108  if (!_fe_problem.converged())
109  return;
110  }
111 }
112 
113 void
114 LStableDirk4::postResidual(NumericVector<Number> & residual)
115 {
116  // Error if _stage got messed up somehow.
117  if (_stage > _n_stages)
118  // the explicit cast prevents strange compiler weirdness with the static
119  // const variable and the variadic mooseError function
120  mooseError("LStableDirk4::postResidual(): Member variable _stage can only have values 1-",
121  (unsigned int)_n_stages,
122  ".");
123 
124  // In the standard RK notation, the residual of stage 1 of s is given by:
125  //
126  // R := M*(Y_i - y_n)/dt - \sum_{j=1}^s a_{ij} * f(t_n + c_j*dt, Y_j) = 0
127  //
128  // where:
129  // .) M is the mass matrix
130  // .) Y_i is the stage solution
131  // .) dt is the timestep, and is accounted for in the _Re_time residual.
132  // .) f are the "non-time" residuals evaluated for a given stage solution.
133  // .) The minus signs are already "baked in" to the residuals and so do not appear below.
134 
135  // Store this stage's non-time residual. We are calling operator=
136  // here, and that calls close().
138 
139  // Build up the residual for this stage.
140  residual.add(1., _Re_time);
141  for (unsigned int j = 0; j < _stage; ++j)
142  residual.add(_a[_stage - 1][j], *_stage_residuals[j]);
143  residual.close();
144 }
NonlinearSystemBase & _nl
Fourth-order diagonally implicit Runge Kutta method (Dirk) with five stages.
Definition: LStableDirk4.h:56
virtual NumericVector< Number > * solutionUDot()=0
void initPetscOutput()
Reinitialize petsc output for proper linear/nonlinear iteration display.
virtual Real & time() const
static const Real _a[_n_stages][_n_stages]
Definition: LStableDirk4.h:92
virtual void computeADTimeDerivatives(DualReal &ad_u_dot, const dof_id_type &dof) const override
method for computing local automatic differentiation time derivatives
Definition: LStableDirk4.C:68
NonlinearSystemBase & getNonlinearSystemBase()
virtual NumericVector< Number > & addVector(const std::string &vector_name, const bool project, const ParallelType type)
Adds a solution length vector to the system.
Definition: SystemBase.C:542
FEProblemBase & _fe_problem
virtual void computeTimeDerivatives() override
Definition: LStableDirk4.C:51
LStableDirk4(const InputParameters &parameters)
Definition: LStableDirk4.C:35
NumericVector< Number > & _Re_non_time
residual vector for non-time contributions
SystemBase & _sys
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
unsigned int _stage
Definition: LStableDirk4.h:76
DualNumber< Real, NumberArray< AD_MAX_DOFS_PER_ELEM, Real > > DualReal
Definition: DualReal.h:29
registerMooseObject("MooseApp", LStableDirk4)
void mooseError(Args &&... args) const
Definition: MooseObject.h:147
virtual void postResidual(NumericVector< Number > &residual) override
Callback to the TimeIntegrator called immediately after the residuals are computed in NonlinearSystem...
Definition: LStableDirk4.C:114
virtual bool converged() override
InputParameters validParams< TimeIntegrator >()
unsigned int _n_linear_iterations
Total number of linear iterations over all stages of the time step.
static const Real _c[_n_stages]
Definition: LStableDirk4.h:88
Real & _du_dot_du
solution vector for
const NumericVector< Number > *const & _solution
solution vectors
Base class for time integrators.
unsigned int getNumLinearIterationsLastSolve() const
Gets the number of linear iterations in the most recent solve.
virtual System & system() override
Get the reference to the libMesh system.
NumericVector< Number > * _stage_residuals[_n_stages]
Definition: LStableDirk4.h:85
NumericVector< Number > & _Re_time
residual vector for time contributions
static const unsigned int _n_stages
Definition: LStableDirk4.h:82
virtual Real & timeOld() const
const NumericVector< Number > & _solution_old
virtual void solve() override
Solves the time step and sets the number of nonlinear and linear iterations.
Definition: LStableDirk4.C:74
void mooseInfo(Args &&... args) const
Definition: MooseObject.h:167
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.
InputParameters validParams< LStableDirk4 >()
Definition: LStableDirk4.C:19
void computeTimeDerivativeHelper(T &u_dot, const T2 &u_old) const
Helper function that actually does the math for computing the time derivative.
Definition: LStableDirk4.h:97