www.mooseframework.org
ImplicitMidpoint.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 "ImplicitMidpoint.h"
11 #include "NonlinearSystem.h"
12 #include "FEProblem.h"
13 #include "PetscSupport.h"
14 
16 
17 template <>
20 {
22 
23  return params;
24 }
25 
27  : TimeIntegrator(parameters),
28  _stage(1),
29  _residual_stage1(_nl.addVector("residual_stage1", false, GHOSTED))
30 {
31  mooseInfo("ImplicitMidpoint and other multistage TimeIntegrators are known not to work with "
32  "Materials/AuxKernels that accumulate 'state' and should be used with caution.");
33 }
34 
35 void
37 {
38  // We are multiplying by the method coefficients in postResidual(), so
39  // the time derivatives are of the same form at every stage although
40  // the current solution varies depending on the stage.
41  if (!_sys.solutionUDot())
42  mooseError("ImplicitMidpoint: Time derivative of solution (`u_dot`) is not stored. Please set "
43  "uDotRequested() to true in FEProblemBase befor requesting `u_dot`.");
44 
45  NumericVector<Number> & u_dot = *_sys.solutionUDot();
46  u_dot = *_solution;
48  u_dot.close();
49  _du_dot_du = 1. / _dt;
50 }
51 
52 void
53 ImplicitMidpoint::computeADTimeDerivatives(DualReal & ad_u_dot, const dof_id_type & dof) const
54 {
56 }
57 
58 void
60 {
61  Real time_new = _fe_problem.time();
62  Real time_old = _fe_problem.timeOld();
63  Real time_half = (time_new + time_old) / 2.;
64 
65  // Reset iteration counts
68 
69  // Compute first stage
71  _console << "1st stage\n";
72  _stage = 1;
73  _fe_problem.time() = time_half;
77 
78  // Abort time step immediately on stage failure - see TimeIntegrator doc page
79  if (!_fe_problem.converged())
80  return;
81 
82  // Compute second stage
84  _console << "2nd stage\n";
85  _stage = 2;
86  _fe_problem.time() = time_new;
90 }
91 
92 void
93 ImplicitMidpoint::postResidual(NumericVector<Number> & residual)
94 {
95  if (_stage == 1)
96  {
97  // In the standard RK notation, the stage 1 residual is given by:
98  //
99  // R := M*(Y_1 - y_n)/dt - (1/2)*f(t_n + dt/2, Y_1) = 0
100  //
101  // where:
102  // .) M is the mass matrix
103  // .) f(t_n + dt/2, Y_1) is saved in _residual_stage1
104  // .) The minus sign is baked in to the non-time residuals, so it does not appear here.
106  _residual_stage1.close();
107 
108  residual.add(1., _Re_time);
109  residual.add(0.5, _Re_non_time);
110  residual.close();
111  }
112  else if (_stage == 2)
113  {
114  // The update step. In the standard RK notation, the update
115  // residual is given by:
116  //
117  // R := M*(y_{n+1} - y_n)/dt - f(t_n + dt/2, Y_1) = 0
118  //
119  // where:
120  // .) M is the mass matrix.
121  // .) f(t_n + dt/2, Y_1) is the residual from stage 1, it has already
122  // been saved as _residual_stage1.
123  // .) The minus signs are "baked in" to the non-time residuals, so
124  // they do not appear here.
125  residual.add(1., _Re_time);
126  residual.add(1., _residual_stage1);
127  residual.close();
128  }
129  else
130  mooseError(
131  "ImplicitMidpoint::postResidual(): _stage = ", _stage, ", only _stage = 1, 2 is allowed.");
132 }
NumericVector< Number > & _residual_stage1
Buffer to store non-time residual from the first stage.
virtual NumericVector< Number > * solutionUDot()=0
void initPetscOutput()
Reinitialize petsc output for proper linear/nonlinear iteration display.
virtual Real & time() const
Second-order Runge-Kutta (implicit midpoint) time integration.
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
ImplicitMidpoint(const InputParameters &parameters)
void mooseError(Args &&... args) const
Definition: MooseObject.h:147
registerMooseObject("MooseApp", ImplicitMidpoint)
void computeADTimeDerivatives(DualReal &ad_u_dot, const dof_id_type &dof) const override
method for computing local automatic differentiation time derivatives
void computeTimeDerivativeHelper(T &u_dot, const T2 &u_old) const
Helper function that actually does the math for computing the time derivative.
virtual void computeTimeDerivatives() override
virtual bool converged() override
InputParameters validParams< TimeIntegrator >()
InputParameters validParams< ImplicitMidpoint >()
unsigned int _n_linear_iterations
Total number of linear iterations over all stages of the time step.
virtual void solve() override
Solves the time step and sets the number of nonlinear and linear iterations.
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.
virtual void postResidual(NumericVector< Number > &residual) override
Callback to the TimeIntegrator called immediately after the residuals are computed in NonlinearSystem...
NumericVector< Number > & _Re_time
residual vector for time contributions
unsigned int _stage
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.