www.mooseframework.org
LStableDirk2.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 "LStableDirk2.h"
11 #include "NonlinearSystem.h"
12 #include "FEProblem.h"
13 #include "PetscSupport.h"
14 
16 
17 template <>
20 {
22  return params;
23 }
24 
26  : TimeIntegrator(parameters),
27  _stage(1),
28  _residual_stage1(_nl.addVector("residual_stage1", false, GHOSTED)),
29  _residual_stage2(_nl.addVector("residual_stage2", false, GHOSTED)),
30  _alpha(1. - 0.5 * std::sqrt(2))
31 {
32  mooseInfo("LStableDirk2 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  // We are multiplying by the method coefficients in postResidual(), so
40  // the time derivatives are of the same form at every stage although
41  // the current solution varies depending on the stage.
42  if (!_sys.solutionUDot())
43  mooseError("LStableDirk2: Time derivative of solution (`u_dot`) is not stored. Please set "
44  "uDotRequested() to true in FEProblemBase befor requesting `u_dot`.");
45 
46  NumericVector<Number> & u_dot = *_sys.solutionUDot();
47  u_dot = *_solution;
49  u_dot.close();
50  _du_dot_du = 1. / _dt;
51 }
52 
53 void
54 LStableDirk2::computeADTimeDerivatives(DualReal & ad_u_dot, const dof_id_type & dof) const
55 {
57 }
58 
59 void
61 {
62  // Time at end of step
63  Real time_new = _fe_problem.time();
64 
65  // Time at beginning of step
66  Real time_old = _fe_problem.timeOld();
67 
68  // Time at stage 1
69  Real time_stage1 = time_old + _alpha * _dt;
70 
71  // Reset iteration counts
74 
75  // Compute first stage
77  _console << "1st stage\n";
78  _stage = 1;
79  _fe_problem.time() = time_stage1;
83 
84  // Abort time step immediately on stage failure - see TimeIntegrator doc page
85  if (!_fe_problem.converged())
86  return;
87 
88  // Compute second stage
90  _console << "2nd stage\n";
91  _stage = 2;
92  _fe_problem.timeOld() = time_stage1;
93  _fe_problem.time() = time_new;
97 
98  // Reset time at beginning of step to its original value
99  _fe_problem.timeOld() = time_old;
100 }
101 
102 void
103 LStableDirk2::postResidual(NumericVector<Number> & residual)
104 {
105  if (_stage == 1)
106  {
107  // In the standard RK notation, the stage 1 residual is given by:
108  //
109  // R := (Y_1 - y_n)/dt - alpha*f(t_n + alpha*dt, Y_1) = 0
110  //
111  // where:
112  // .) f(t_n + alpha*dt, Y_1) corresponds to the residuals of the
113  // non-time kernels. We'll save this as "_residual_stage1" to use
114  // later.
115  // .) (Y_1 - y_n)/dt corresponds to the residual of the time kernels.
116  // .) The minus sign in front of alpha is already "baked in" to
117  // the non-time residuals, so it does not appear here.
119  _residual_stage1.close();
120 
121  residual.add(1., _Re_time);
122  residual.add(_alpha, _residual_stage1);
123  residual.close();
124  }
125  else if (_stage == 2)
126  {
127  // In the standard RK notation, the stage 2 residual is given by:
128  //
129  // R := (Y_2 - y_n)/dt - (1-alpha)*f(t_n + alpha*dt, Y_1) - alpha*f(t_n + dt, Y_2) = 0
130  //
131  // where:
132  // .) f(t_n + alpha*dt, Y_1) has already been saved as _residual_stage1.
133  // .) f(t_n + dt, Y_2) will now be saved as "_residual_stage2".
134  // .) (Y_2 - y_n)/dt corresponds to the residual of the time kernels.
135  // .) The minus signs are once again "baked in" to the non-time
136  // residuals, so they do not appear here.
137  //
138  // The solution at the end of stage 2, i.e. Y_2, is also the final solution.
140  _residual_stage2.close();
141 
142  residual.add(1., _Re_time);
143  residual.add(1. - _alpha, _residual_stage1);
144  residual.add(_alpha, _residual_stage2);
145  residual.close();
146  }
147  else
148  mooseError("LStableDirk2::postResidual(): Member variable _stage can only have values 1 or 2.");
149 }
virtual NumericVector< Number > * solutionUDot()=0
void initPetscOutput()
Reinitialize petsc output for proper linear/nonlinear iteration display.
virtual Real & time() const
virtual void computeTimeDerivatives() override
Definition: LStableDirk2.C:37
NonlinearSystemBase & getNonlinearSystemBase()
InputParameters validParams< LStableDirk2 >()
Definition: LStableDirk2.C:19
FEProblemBase & _fe_problem
LStableDirk2(const InputParameters &parameters)
Definition: LStableDirk2.C:25
NumericVector< Number > & _Re_non_time
residual vector for non-time contributions
const Real _alpha
Definition: LStableDirk2.h:71
SystemBase & _sys
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
DualNumber< Real, NumberArray< AD_MAX_DOFS_PER_ELEM, Real > > DualReal
Definition: DualReal.h:29
void mooseError(Args &&... args) const
Definition: MooseObject.h:147
virtual void solve() override
Solves the time step and sets the number of nonlinear and linear iterations.
Definition: LStableDirk2.C:60
virtual bool converged() override
InputParameters validParams< TimeIntegrator >()
unsigned int _n_linear_iterations
Total number of linear iterations over all stages of the time step.
virtual void postResidual(NumericVector< Number > &residual) override
Callback to the TimeIntegrator called immediately after the residuals are computed in NonlinearSystem...
Definition: LStableDirk2.C:103
void computeTimeDerivativeHelper(T &u_dot, const T2 &u_old) const
Helper function that actually does the math for computing the time derivative.
Definition: LStableDirk2.h:76
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.
unsigned int _stage
Indicates the current stage (1 or 2).
Definition: LStableDirk2.h:62
NumericVector< Number > & _Re_time
residual vector for time contributions
virtual Real & timeOld() const
const NumericVector< Number > & _solution_old
void mooseInfo(Args &&... args) const
Definition: MooseObject.h:167
void computeADTimeDerivatives(DualReal &ad_u_dot, const dof_id_type &dof) const override
method for computing local automatic differentiation time derivatives
Definition: LStableDirk2.C:54
unsigned int _n_nonlinear_iterations
Total number of nonlinear iterations over all stages of the time step.
NumericVector< Number > & _residual_stage2
Buffer to store non-time residual from second stage solve.
Definition: LStableDirk2.h:68
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.
Second order diagonally implicit Runge Kutta method (Dirk) with two stages.
Definition: LStableDirk2.h:43
registerMooseObject("MooseApp", LStableDirk2)
NumericVector< Number > & _residual_stage1
Buffer to store non-time residual from first stage solve.
Definition: LStableDirk2.h:65