www.mooseframework.org
NewmarkBeta.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 "NewmarkBeta.h"
11 #include "NonlinearSystem.h"
12 #include "FEProblemBase.h"
13 
14 registerMooseObject("MooseApp", NewmarkBeta);
15 
16 template <>
19 {
21  params.addClassDescription(
22  "Computes the first and second time derivative of variable using Newmark-Beta method.");
23  params.addRangeCheckedParam<Real>("beta", 0.25, "beta > 0.0", "beta value");
24  params.addRangeCheckedParam<Real>("gamma", 0.5, "gamma >= 0.5", "gamma value");
25  return params;
26 }
27 
29  : TimeIntegrator(parameters),
30  _beta(getParam<Real>("beta")),
31  _gamma(getParam<Real>("gamma")),
32  _du_dotdot_du(_sys.duDotDotDu())
33 {
37 
38  if (_gamma > 2.0 * _beta)
39  mooseError("NewmarkBeta: For Newmark method to be unconditionally stable, gamma should lie "
40  "between 0.5 and 2.0*beta.");
41 
42  if (_gamma != 0.5)
43  mooseWarning("NewmarkBeta: For gamma > 0.5, Newmark method is only first order accurate. "
44  "Please use either HHT time integration method or set gamma = 0.5 for second "
45  "order solution accuracy with time.");
46 }
47 
49 
50 void
52 {
53  // compute second derivative
54  // according to Newmark-Beta method
55  // u_dotdot = first_term - second_term - third_term
56  // first_term = (u - u_old) / beta / dt ^ 2
57  // second_term = u_dot_old / beta / dt
58  // third_term = u_dotdot_old * (1 / 2 / beta - 1)
59  if (!_sys.solutionUDot())
60  mooseError("NewmarkBeta: Time derivative of solution (`u_dot`) is not stored. Please set "
61  "uDotRequested() to true in FEProblemBase befor requesting `u_dot`.");
62 
63  if (!_sys.solutionUDotDot())
64  mooseError("NewmarkBeta: Second time derivative of solution (`u_dotdot`) is not stored. Please "
65  "set uDotDotRequested() to true in FEProblemBase befor requesting `u_dotdot`.");
66 
67  if (!_sys.solutionUDotOld())
68  mooseError("NewmarkBeta: Old time derivative of solution (`u_dot_old`) is not stored. Please "
69  "set uDotOldRequested() to true in FEProblemBase befor requesting `u_dot_old`.");
70 
71  if (!_sys.solutionUDotDotOld())
72  mooseError("NewmarkBeta: Old second time derivative of solution (`u_dotdot_old`) is not "
73  "stored. Please set uDotDotOldRequested() to true in FEProblemBase befor requesting "
74  "`u_dotdot_old`.");
75 
76  NumericVector<Number> & u_dot = *_sys.solutionUDot();
77  NumericVector<Number> & u_dotdot = *_sys.solutionUDotDot();
78  NumericVector<Number> & u_dot_old = *_sys.solutionUDotOld();
79  NumericVector<Number> & u_dotdot_old = *_sys.solutionUDotDotOld();
80 
81  u_dotdot = *_solution;
82 
83  computeTimeDerivativeHelper(u_dot, _solution_old, u_dot_old, u_dotdot, u_dotdot_old);
84 
85  // make sure _u_dotdot and _u_dot are in good state
86  u_dotdot.close();
87  u_dot.close();
88 
89  // used for Jacobian calculations
90  _du_dotdot_du = 1.0 / _beta / _dt / _dt;
91  _du_dot_du = _gamma / _beta / _dt;
92 }
93 
94 void
95 NewmarkBeta::computeADTimeDerivatives(DualReal & ad_u_dot, const dof_id_type & dof) const
96 {
97  const auto & u_old = _solution_old(dof);
98  const auto & u_dot_old = (*_sys.solutionUDotOld())(dof);
99  const auto & u_dotdot_old = (*_sys.solutionUDotDotOld())(dof);
100 
101  auto u_dotdot = ad_u_dot;
102 
103  computeTimeDerivativeHelper(ad_u_dot, u_old, u_dot_old, u_dotdot, u_dotdot_old);
104 }
105 
106 void
107 NewmarkBeta::postResidual(NumericVector<Number> & residual)
108 {
109  residual += _Re_time;
110  residual += _Re_non_time;
111  residual.close();
112 }
virtual void setUDotDotOldRequested(const bool u_dotdot_old_requested)
Set boolean flag to true to store old solution second time derivative.
virtual NumericVector< Number > * solutionUDot()=0
Real & _du_dotdot_du
solution vector for
Definition: NewmarkBeta.h:53
virtual ~NewmarkBeta()
Definition: NewmarkBeta.C:48
void mooseWarning(Args &&... args) const
Definition: MooseObject.h:155
virtual void setUDotDotRequested(const bool u_dotdot_requested)
Set boolean flag to true to store solution second time derivative.
FEProblemBase & _fe_problem
Newmark-Beta time integration method.
Definition: NewmarkBeta.h:23
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 postResidual(NumericVector< Number > &residual) override
Callback to the TimeIntegrator called immediately after the residuals are computed in NonlinearSystem...
Definition: NewmarkBeta.C:107
DualNumber< Real, NumberArray< AD_MAX_DOFS_PER_ELEM, Real > > DualReal
Definition: DualReal.h:29
virtual NumericVector< Number > * solutionUDotDotOld()=0
void mooseError(Args &&... args) const
Definition: MooseObject.h:147
Real _gamma
Newmark time integration parameter-gamma.
Definition: NewmarkBeta.h:50
virtual void setUDotOldRequested(const bool u_dot_old_requested)
Set boolean flag to true to store old solution time derivative.
Real _beta
Newmark time integration parameter-beta.
Definition: NewmarkBeta.h:47
InputParameters validParams< TimeIntegrator >()
registerMooseObject("MooseApp", NewmarkBeta)
NewmarkBeta(const InputParameters &parameters)
Definition: NewmarkBeta.C:28
InputParameters validParams< NewmarkBeta >()
Definition: NewmarkBeta.C:18
virtual void computeADTimeDerivatives(DualReal &ad_u_dot, const dof_id_type &dof) const override
method for computing local automatic differentiation time derivatives
Definition: NewmarkBeta.C:95
Real & _du_dot_du
solution vector for
const NumericVector< Number > *const & _solution
solution vectors
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:58
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...
virtual NumericVector< Number > * solutionUDotOld()=0
const NumericVector< Number > & _solution_old
virtual NumericVector< Number > * solutionUDotDot()=0
void addRangeCheckedParam(const std::string &name, const T &value, const std::string &parsed_function, const std::string &doc_string)
virtual void computeTimeDerivatives() override
Definition: NewmarkBeta.C:51