www.mooseframework.org
BDF2.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 "BDF2.h"
11 #include "NonlinearSystem.h"
12 
13 registerMooseObject("MooseApp", BDF2);
14 
17 {
19  params.addClassDescription(
20  "Second order backward differentiation formula time integration scheme.");
21  return params;
22 }
23 
24 BDF2::BDF2(const InputParameters & parameters)
25  : TimeIntegrator(parameters),
26  _weight(declareRestartableData<std::vector<Real>>("weight")),
27  _solution_older(_sys.solutionState(2))
28 {
29  _weight.resize(3);
30 }
31 
32 void
34 {
35  if (_t_step > 1)
36  {
37  Real sum = _dt + _dt_old;
38  _weight[0] = 1. + _dt / sum;
39  _weight[1] = -sum / _dt_old;
40  _weight[2] = _dt * _dt / _dt_old / sum;
41  }
42 }
43 
44 void
46 {
47  if (!_sys.solutionUDot())
48  mooseError("BDF2: Time derivative of solution (`u_dot`) is not stored. Please set "
49  "uDotRequested() to true in FEProblemBase befor requesting `u_dot`.");
50 
52  if (_t_step == 1)
53  {
54  u_dot = *_solution;
55  _du_dot_du = 1. / _dt;
56  }
57  else
58  {
59  u_dot.zero();
60  _du_dot_du = _weight[0] / _dt;
61  }
63  u_dot.close();
64 }
65 
66 void
68  const dof_id_type & dof,
69  DualReal & /*ad_u_dotdot*/) const
70 {
71  auto ad_sln = ad_u_dot;
72  if (_t_step != 1)
73  ad_u_dot = 0;
74  computeTimeDerivativeHelper(ad_u_dot, ad_sln, _solution_old(dof), _solution_older(dof));
75 }
76 
77 void
79 {
80  residual += _Re_time;
81  residual += _Re_non_time;
82  residual.close();
83 }
virtual NumericVector< Number > * solutionUDot()=0
virtual void computeTimeDerivatives() override
Computes the time derivative and the Jacobian of the time derivative.
Definition: BDF2.C:45
DualNumber< Real, DNDerivativeType, true > DualReal
Definition: DualReal.h:49
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...
virtual void zero()=0
BDF2 time integrator.
Definition: BDF2.h:18
const NumericVector< Number > & _solution_older
The older solution.
Definition: BDF2.h:44
std::vector< Real > & _weight
Definition: BDF2.h:41
virtual void preStep() override
Definition: BDF2.C:33
void computeTimeDerivativeHelper(T &u_dot, const T2 &u, const T3 &u_old, const T4 &u_older) const
Helper function that actually does the math for computing the time derivative.
Definition: BDF2.h:53
registerMooseObject("MooseApp", BDF2)
virtual void close()=0
BDF2(const InputParameters &parameters)
Definition: BDF2.C:24
static InputParameters validParams()
Definition: BDF2.C:16
Real & _du_dot_du
Derivative of time derivative with respect to current solution: .
const NumericVector< Number > *const & _solution
solution vectors
virtual void postResidual(NumericVector< Number > &residual) override
Callback to the TimeIntegrator called immediately after the residuals are computed in NonlinearSystem...
Definition: BDF2.C:78
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Base class for time integrators.
void computeADTimeDerivatives(DualReal &ad_u_dot, const dof_id_type &dof, DualReal &ad_u_dotdot) const override
method for computing local automatic differentiation time derivatives
Definition: BDF2.C:67
void mooseError(Args &&... args) const
Emits an error prefixed with object name and type.
NumericVector< Number > & _Re_time
residual vector for time contributions
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...
const NumericVector< Number > & _solution_old
static InputParameters validParams()
uint8_t dof_id_type