www.mooseframework.org
ExplicitRK2.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 "ExplicitRK2.h"
11 #include "NonlinearSystem.h"
12 #include "FEProblem.h"
13 #include "PetscSupport.h"
14 
15 template <>
18 {
20 
21  return params;
22 }
23 
25  : TimeIntegrator(parameters),
26  _stage(1),
27  _residual_old(_nl.addVector("residual_old", false, GHOSTED))
28 {
29  mooseInfo("ExplicitRK2-derived TimeIntegrators (ExplicitMidpoint, Heun, Ralston) and other "
30  "multistage TimeIntegrators are known not to work with "
31  "Materials/AuxKernels that accumulate 'state' and should be used with caution.");
32 }
33 
34 void
36 {
37  if (_dt == _dt_old)
39  else
41 }
42 
43 void
45 {
46  // Since advanceState() is called in between stages 2 and 3, this
47  // changes the meaning of "_solution_old". In the second stage,
48  // "_solution_older" is actually the original _solution_old.
49  if (!_sys.solutionUDot())
50  mooseError("ExplicitRK2: Time derivative of solution (`u_dot`) is not stored. Please set "
51  "uDotRequested() to true in FEProblemBase befor requesting `u_dot`.");
52 
53  NumericVector<Number> & u_dot = *_sys.solutionUDot();
54  u_dot = *_solution;
56 
57  _du_dot_du = 1. / _dt;
58  u_dot.close();
59 }
60 
61 void
62 ExplicitRK2::computeADTimeDerivatives(DualReal & ad_u_dot, const dof_id_type & dof) const
63 {
65 }
66 
67 void
69 {
70  Real time_new = _fe_problem.time();
71  Real time_old = _fe_problem.timeOld();
72  Real time_stage2 = time_old + a() * _dt;
73 
74  // Reset iteration counts
77 
78  // There is no work to do for the first stage (Y_1 = y_n). The
79  // first solve therefore happens in the second stage. Note that the
80  // non-time Kernels (which should be marked implicit=false) are
81  // evaluated at the old solution during this stage.
83  _console << "1st solve\n";
84  _stage = 2;
85  _fe_problem.timeOld() = time_old;
86  _fe_problem.time() = time_stage2;
90 
91  // Abort time step immediately on stage failure - see TimeIntegrator doc page
92  if (!_fe_problem.converged())
93  return;
94 
95  // Advance solutions old->older, current->old. Also moves Material
96  // properties and other associated state forward in time.
98 
99  // The "update" stage (which we call stage 3) requires an additional
100  // solve with the mass matrix.
102  _console << "2nd solve\n";
103  _stage = 3;
104  _fe_problem.timeOld() = time_stage2;
105  _fe_problem.time() = time_new;
109 
110  // Reset time at beginning of step to its original value
111  _fe_problem.timeOld() = time_old;
112 }
113 
114 void
115 ExplicitRK2::postResidual(NumericVector<Number> & residual)
116 {
117  if (_stage == 1)
118  {
119  // If postResidual() is called before solve(), _stage==1 and we don't
120  // need to do anything.
121  }
122  else if (_stage == 2)
123  {
124  // In the standard RK notation, the stage 2 residual is given by:
125  //
126  // R := M*(Y_2 - y_n)/dt - a*f(t_n, Y_1) = 0
127  //
128  // where:
129  // .) M is the mass matrix.
130  // .) f(t_n, Y_1) is the residual we are currently computing,
131  // since this method is intended to be used with "implicit=false"
132  // kernels.
133  // .) M*(Y_2 - y_n)/dt corresponds to the residual of the time kernels.
134  // .) The minus signs are "baked in" to the non-time residuals, so
135  // they do not appear here.
136  // .) The current non-time residual is saved for the next stage.
138  _residual_old.close();
139 
140  residual.add(1., _Re_time);
141  residual.add(a(), _residual_old);
142  residual.close();
143  }
144  else if (_stage == 3)
145  {
146  // In the standard RK notation, the update step residual is given by:
147  //
148  // R := M*(y_{n+1} - y_n)/dt - f(t_n+dt/2, Y_2) = 0
149  //
150  // where:
151  // .) M is the mass matrix.
152  // .) f(t_n+dt/2, Y_2) is the residual from stage 2.
153  // .) The minus sign is already baked in to the non-time
154  // residuals, so it does not appear here.
155  // .) Although this is an update step, we have to do a "solve"
156  // using the mass matrix.
157  residual.add(1., _Re_time);
158  residual.add(b1(), _residual_old);
159  residual.add(b2(), _Re_non_time);
160  residual.close();
161  }
162  else
163  mooseError("ExplicitRK2::postResidual(): _stage = ", _stage, ", only _stage = 1-3 is allowed.");
164 }
virtual NumericVector< Number > * solutionUDot()=0
void initPetscOutput()
Reinitialize petsc output for proper linear/nonlinear iteration display.
InputParameters validParams< ExplicitRK2 >()
Definition: ExplicitRK2.C:17
virtual Real & time() const
NonlinearSystemBase & getNonlinearSystemBase()
FEProblemBase & _fe_problem
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...
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 advanceState()
Advance all of the state holding vectors / datastructures so that we can move to the next timestep...
virtual bool converged() override
virtual void solve() override
Solves the time step and sets the number of nonlinear and linear iterations.
Definition: ExplicitRK2.C:68
virtual Real a() const =0
The method coefficients.
InputParameters validParams< TimeIntegrator >()
unsigned int _n_linear_iterations
Total number of linear iterations over all stages of the time step.
unsigned int _stage
Definition: ExplicitRK2.h:79
const NumericVector< Number > & _solution_older
Real & _du_dot_du
solution vector for
NumericVector< Number > & _residual_old
Buffer to store non-time residual from the first stage.
Definition: ExplicitRK2.h:82
const NumericVector< Number > *const & _solution
solution vectors
Base class for time integrators.
virtual Real b1() const =0
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)
virtual System & system() override
Get the reference to the libMesh system.
virtual void preSolve() override
Definition: ExplicitRK2.C:35
ExplicitRK2(const InputParameters &parameters)
Definition: ExplicitRK2.C:24
virtual void postResidual(NumericVector< Number > &residual) override
Callback to the TimeIntegrator called immediately after the residuals are computed in NonlinearSystem...
Definition: ExplicitRK2.C:115
NumericVector< Number > & _Re_time
residual vector for time contributions
virtual void computeTimeDerivatives() override
Definition: ExplicitRK2.C:44
virtual Real & timeOld() const
const NumericVector< Number > & _solution_old
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.
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.
Definition: ExplicitRK2.h:94
virtual Real b2() const =0
void computeADTimeDerivatives(DualReal &ad_u_dot, const dof_id_type &dof) const override
method for computing local automatic differentiation time derivatives
Definition: ExplicitRK2.C:62