https://mooseframework.inl.gov
NewmarkBeta.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://mooseframework.inl.gov
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 "NewmarkBeta.h"
11 #include "NonlinearSystem.h"
12 #include "FEProblemBase.h"
13 
14 registerMooseObject("MooseApp", NewmarkBeta);
15 
18 {
20  params.addClassDescription(
21  "Computes the first and second time derivative of variable using Newmark-Beta method.");
22  params.addRangeCheckedParam<Real>("beta", 0.25, "beta > 0.0", "beta value");
23  params.addRangeCheckedParam<Real>("gamma", 0.5, "gamma >= 0.5", "gamma value");
24  params.addParam<int>("inactive_tsteps",
25  0,
26  "The time derivatives will set to be zero for this number of time steps.");
27  return params;
28 }
29 
31  : TimeIntegrator(parameters),
32  _beta(getParam<Real>("beta")),
33  _gamma(getParam<Real>("gamma")),
34  _inactive_tsteps(getParam<int>("inactive_tsteps")),
35  _du_dotdot_du(_sys.duDotDotDu())
36 {
40 
41  if (_gamma > 2.0 * _beta)
42  mooseError("NewmarkBeta: For Newmark method to be unconditionally stable, gamma should lie "
43  "between 0.5 and 2.0*beta.");
44 
45  if (_gamma != 0.5)
46  mooseWarning("NewmarkBeta: For gamma > 0.5, Newmark method is only first order accurate. "
47  "Please use either HHT time integration method or set gamma = 0.5 for second "
48  "order solution accuracy with time.");
49 }
50 
52 
53 void
55 {
56  if (!_sys.solutionUDot())
57  mooseError("NewmarkBeta: Time derivative of solution (`u_dot`) is not stored. Please set "
58  "uDotRequested() to true in FEProblemBase befor requesting `u_dot`.");
59 
60  if (!_sys.solutionUDotDot())
61  mooseError("NewmarkBeta: Second time derivative of solution (`u_dotdot`) is not stored. Please "
62  "set uDotDotRequested() to true in FEProblemBase befor requesting `u_dotdot`.");
63 
64  if (!_sys.solutionUDotOld())
65  mooseError("NewmarkBeta: Old time derivative of solution (`u_dot_old`) is not stored. Please "
66  "set uDotOldRequested() to true in FEProblemBase befor requesting `u_dot_old`.");
67 
68  if (!_sys.solutionUDotDotOld())
69  mooseError("NewmarkBeta: Old second time derivative of solution (`u_dotdot_old`) is not "
70  "stored. Please set uDotDotOldRequested() to true in FEProblemBase befor requesting "
71  "`u_dotdot_old`.");
72 
76  NumericVector<Number> & u_dotdot_old = *_sys.solutionUDotDotOld();
77 
79  {
80  u_dot.zero();
81  u_dotdot.zero();
82  }
83  else
84  {
85  u_dotdot = *_solution;
86  computeTimeDerivativeHelper(u_dot, _solution_old, u_dot_old, u_dotdot, u_dotdot_old);
87  }
88 
89  // make sure _u_dotdot and _u_dot are in good state
90  u_dotdot.close();
91  u_dot.close();
92 
93  // used for Jacobian calculations
94  _du_dotdot_du = 1.0 / _beta / _dt / _dt;
96 }
97 
98 void
100  const dof_id_type & dof,
101  ADReal & ad_u_dotdot) const
102 {
103  const auto & u_old = _solution_old(dof);
104  const auto & u_dot_old = (*_sys.solutionUDotOld())(dof);
105  const auto & u_dotdot_old = (*_sys.solutionUDotDotOld())(dof);
106 
107  // Seeds ad_u_dotdot with _ad_dof_values and associated derivatives provided via ad_u_dot from
108  // MooseVariableData
109  ad_u_dotdot = ad_u_dot;
110 
111  computeTimeDerivativeHelper(ad_u_dot, u_old, u_dot_old, ad_u_dotdot, u_dotdot_old);
112 }
113 
114 void
116 {
117  residual += *_Re_time;
118  residual += *_Re_non_time;
119  residual.close();
120 }
121 
122 Real
124 {
125  return _gamma / _beta;
126 }
virtual void setUDotDotOldRequested(const bool u_dotdot_old_requested)
Set boolean flag to true to store old solution second time derivative.
Real & _du_dotdot_du
solution vector for
Definition: NewmarkBeta.h:57
virtual ~NewmarkBeta()
Definition: NewmarkBeta.C:51
virtual void setUDotDotRequested(const bool u_dotdot_requested)
Set boolean flag to true to store solution second time derivative.
virtual Real duDotDuCoeff() const override
Definition: NewmarkBeta.C:123
FEProblemBase & _fe_problem
Reference to the problem.
static InputParameters validParams()
Definition: NewmarkBeta.C:17
Newmark-Beta time integration method.
Definition: NewmarkBeta.h:18
SystemBase & _sys
Reference to the system this time integrator operates on.
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
virtual void postResidual(NumericVector< Number > &residual) override
Callback to the NonLinearTimeIntegratorInterface called immediately after the residuals are computed ...
Definition: NewmarkBeta.C:115
virtual NumericVector< Number > * solutionUDotDotOld()
Definition: SystemBase.h:255
DualNumber< Real, DNDerivativeType, true > ADReal
Definition: ADRealForward.h:46
Real & _dt
The current time step size.
NumericVector< Number > * _Re_time
residual vector for time contributions
Real _gamma
Newmark time integration parameter-gamma.
Definition: NewmarkBeta.h:51
void mooseWarning(Args &&... args) const
Emits a warning prefixed with object name and type.
virtual void zero()=0
virtual void setUDotOldRequested(const bool u_dot_old_requested)
Set boolean flag to true to store old solution time derivative.
int _inactive_tsteps
Inactive time steps.
Definition: NewmarkBeta.h:54
Real _beta
Newmark time integration parameter-beta.
Definition: NewmarkBeta.h:48
virtual NumericVector< Number > * solutionUDot()
Definition: SystemBase.h:252
virtual NumericVector< Number > * solutionUDotOld()
Definition: SystemBase.h:254
registerMooseObject("MooseApp", NewmarkBeta)
NewmarkBeta(const InputParameters &parameters)
Definition: NewmarkBeta.C:30
virtual void close()=0
virtual int & timeStep() const
const NumericVector< Number > *const & _solution
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Base class for time integrators.
void computeTimeDerivativeHelper(T &u_dot, const T2 &u_old, const T3 &u_dot_old, T4 &u_dotdot, const T5 &u_dotdot_old) const
Helper function that actually does the math for computing the time derivative.
Definition: NewmarkBeta.h:62
void mooseError(Args &&... args) const
Emits an error prefixed with object name and type.
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...
void addParam(const std::string &name, const S &value, const std::string &doc_string)
These methods add an optional parameter and a documentation string to the InputParameters object...
const NumericVector< Number > & _solution_old
void addRangeCheckedParam(const std::string &name, const T &value, const std::string &parsed_function, const std::string &doc_string)
virtual void computeADTimeDerivatives(ADReal &ad_u_dot, const dof_id_type &dof, ADReal &ad_u_dotdot) const override
method for computing local automatic differentiation time derivatives
Definition: NewmarkBeta.C:99
virtual NumericVector< Number > * solutionUDotDot()
Definition: SystemBase.h:253
NumericVector< Number > * _Re_non_time
residual vector for non-time contributions
virtual void computeTimeDerivatives() override
Computes the time derivative and the Jacobian of the time derivative.
Definition: NewmarkBeta.C:54
static InputParameters validParams()
void ErrorVector unsigned int
uint8_t dof_id_type
void computeDuDotDu()
Compute _du_dot_du.