https://mooseframework.inl.gov
ImplicitEuler.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 "ImplicitEuler.h"
11 #include "NonlinearSystem.h"
12 
14 
17 {
19  params.addClassDescription("Time integration using the implicit Euler method.");
20  return params;
21 }
22 
23 ImplicitEuler::ImplicitEuler(const InputParameters & parameters) : TimeIntegrator(parameters) {}
24 
26 
27 void
29 {
30  if (!_sys.solutionUDot())
31  mooseError("ImplicitEuler: Time derivative of solution (`u_dot`) is not stored. Please set "
32  "uDotRequested() to true in FEProblemBase befor requesting `u_dot`.");
33 
35  if (!_var_restriction)
36  {
37  u_dot = *_solution;
39  }
40  else
41  {
42  auto u_dot_sub = u_dot.get_subvector(_local_indices);
45  *u_dot_sub = *_solution_sub;
47  u_dot.restore_subvector(std::move(u_dot_sub), _local_indices);
48  // Scatter info needed for ghosts
49  u_dot.close();
50  }
51 
53 }
54 
55 void
57  const dof_id_type & dof,
58  ADReal & /*ad_u_dotdot*/) const
59 {
61 }
62 
63 void
65 {
66  if (!_var_restriction)
67  {
68  residual += *_Re_time;
69  residual += *_Re_non_time;
70  residual.close();
71  }
72  else
73  {
74  auto residual_sub = residual.get_subvector(_local_indices);
75  auto re_time_sub = _Re_time->get_subvector(_local_indices);
76  auto re_non_time_sub = _Re_non_time->get_subvector(_local_indices);
77  *residual_sub += *re_time_sub;
78  *residual_sub += *re_non_time_sub;
79  residual.restore_subvector(std::move(residual_sub), _local_indices);
80  _Re_time->restore_subvector(std::move(re_time_sub), _local_indices);
81  _Re_non_time->restore_subvector(std::move(re_non_time_sub), _local_indices);
82  }
83 }
84 
85 Real
87  const std::vector<Real> & factors) const
88 {
89  mooseAssert(factors.size() == numStatesRequired(),
90  "Either too many or too few states are given!");
91  return factors[0] * _solution_old(dof_id) / _dt;
92 }
93 
94 Real
96 {
97  return factor / _dt;
98 }
bool & _var_restriction
Whether the user has requested that the time integrator be applied to a subset of variables...
virtual void computeTimeDerivatives() override
Computes the time derivative and the Jacobian of the time derivative.
Definition: ImplicitEuler.C:28
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...
registerMooseObject("MooseApp", ImplicitEuler)
DualNumber< Real, DNDerivativeType, true > ADReal
Definition: ADRealForward.h:46
virtual ~ImplicitEuler()
Definition: ImplicitEuler.C:25
Real & _dt
The current time step size.
NumericVector< Number > * _Re_time
residual vector for time contributions
void computeTimeDerivativeHelper(T &u_dot, const T2 &u_old) const
Helper function that actually does the math for computing the time derivative.
Definition: ImplicitEuler.h:47
virtual void create_subvector(NumericVector< Number > &, const std::vector< numeric_index_type > &, bool=true) const
static InputParameters validParams()
Definition: ImplicitEuler.C:16
virtual NumericVector< Number > * solutionUDot()
Definition: SystemBase.h:252
std::vector< dof_id_type > & _local_indices
The local degree of freedom indices this time integrator is being applied to.
virtual std::unique_ptr< NumericVector< Number > > get_subvector(const std::vector< numeric_index_type > &)
std::unique_ptr< NumericVector< Number > > & _solution_old_sub
virtual void close()=0
const NumericVector< Number > *const & _solution
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Base class for time integrators.
ImplicitEuler(const InputParameters &parameters)
Definition: ImplicitEuler.C:23
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: ImplicitEuler.C:56
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...
virtual Real timeDerivativeMatrixContribution(const Real factor) const override
The time derivative&#39;s contribution to the right hand side of a linear system.
Definition: ImplicitEuler.C:95
const NumericVector< Number > & _solution_old
virtual Real timeDerivativeRHSContribution(const dof_id_type dof_id, const std::vector< Real > &factors) const override
The time derivative&#39;s contribution to the right hand side of a linear system.
Definition: ImplicitEuler.C:86
virtual void postResidual(NumericVector< Number > &residual) override
Callback to the NonLinearTimeIntegratorInterface called immediately after the residuals are computed ...
Definition: ImplicitEuler.C:64
NumericVector< Number > * _Re_non_time
residual vector for non-time contributions
static InputParameters validParams()
Implicit Euler&#39;s method.
Definition: ImplicitEuler.h:17
std::unique_ptr< NumericVector< Number > > & _solution_sub
uint8_t dof_id_type
void computeDuDotDu()
Compute _du_dot_du.
virtual unsigned int numStatesRequired() const
Return the number of states this requires in a linear system setting.
virtual void restore_subvector(std::unique_ptr< NumericVector< Number >>, const std::vector< numeric_index_type > &)