https://mooseframework.inl.gov
ActuallyExplicitEuler.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 // MOOSE includes
11 #include "ActuallyExplicitEuler.h"
12 #include "NonlinearSystem.h"
13 #include "FEProblem.h"
14 
15 // libMesh includes
16 #include "libmesh/nonlinear_solver.h"
17 
18 using namespace libMesh;
19 
21 
24 {
26 
27  params.addClassDescription(
28  "Implementation of Explicit/Forward Euler without invoking any of the nonlinear solver");
29 
30  params.addParam<bool>("use_constant_mass",
31  false,
32  "If set to true, will only compute the mass matrix in the first time step, "
33  "and keep using it throughout the simulation.");
34 
35  return params;
36 }
37 
39  : ExplicitTimeIntegrator(parameters), _constant_mass(getParam<bool>("use_constant_mass"))
40 {
41 }
42 
43 void
45 {
46  if (!_sys.solutionUDot())
47  mooseError("ActuallyExplicitEuler: Time derivative of solution (`u_dot`) is not stored. Please "
48  "set uDotRequested() to true in FEProblemBase before requesting `u_dot`.");
49 
51  u_dot = *_solution;
53  u_dot.close();
55 }
56 
57 void
59  const dof_id_type & dof,
60  ADReal & /*ad_u_dotdot*/) const
61 {
63 }
64 
65 void
67 {
68  // Reset iteration counts
71 
73 
74  // Set time to the time at which to evaluate the residual
77 
78  // Compute the residual
79  _explicit_residual->zero();
82 
83  // Move the residual to the RHS
84  *_explicit_residual *= -1.0;
85 
86  // Compute the mass matrix
87  auto & mass_matrix = _nonlinear_implicit_system->get_system_matrix();
88  if (!_constant_mass || (_constant_mass && _t_step == 1))
91 
92  // Perform the linear solve
93  bool converged = performExplicitSolve(mass_matrix);
94 
95  // Update the solution
98 
99  // Constraints may be solved in an uncoupled way. For example, momentum-balance equations may be
100  // solved node-wise and then the solution (e.g. velocities or positions)can be applied to those
101  // nodes without solving for such constraints on a system level. This strategy is being used for
102  // node-face contact in explicit dynamics.
104 
105  // Enforce contraints on the solution
110 
112 
113  _nonlinear_implicit_system->nonlinear_solver->converged = converged;
114 }
115 
116 void
118 {
119  residual += *_Re_time;
120  residual += *_Re_non_time;
121  residual.close();
122 
123  // Reset time to the time at which to evaluate nodal BCs, which comes next
125 }
virtual void computeJacobianTag(const NumericVector< libMesh::Number > &soln, libMesh::SparseMatrix< libMesh::Number > &jacobian, TagID tag)
Form a Jacobian matrix for a given tag.
std::unique_ptr< NonlinearSolver< Number > > nonlinear_solver
void overwriteNodeFace(NumericVector< Number > &soln)
Called from explicit time stepping to overwrite boundary positions (explicit dynamics).
virtual Real & time() const
NumericVector< Real > * _solution_update
Solution vector for the linear solve.
Real _current_time
Save off current time to reset it back and forth.
virtual void postResidual(NumericVector< Number > &residual) override
Callback to the NonLinearTimeIntegratorInterface called immediately after the residuals are computed ...
FEProblemBase & _fe_problem
Reference to the problem.
NonlinearSystemBase * _nl
Pointer to the nonlinear system, can happen that we dont have any.
TagID _Ke_time_tag
For computing the mass matrix.
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...
static InputParameters validParams()
void setSolution(const NumericVector< Number > &soln)
Set the solution to a given vector.
Definition: SolverSystem.C:67
The following methods are specializations for using the libMesh::Parallel::packed_range_* routines fo...
DualNumber< Real, DNDerivativeType, true > ADReal
Definition: ADRealForward.h:46
NumericVector< Number > * _Re_time
residual vector for time contributions
libMesh::NonlinearImplicitSystem * _nonlinear_implicit_system
libMesh nonlinear implicit system, if applicable; otherwise, nullptr
static InputParameters validParams()
NumericVector< Real > * _explicit_residual
Residual used for the RHS.
virtual NumericVector< Number > * solutionUDot()
Definition: SystemBase.h:252
std::unique_ptr< NumericVector< Number > > solution
unsigned int _n_linear_iterations
Total number of linear iterations over all stages of the time step.
registerMooseObject("MooseApp", ActuallyExplicitEuler)
virtual bool performExplicitSolve(SparseMatrix< Number > &mass_matrix)
Solves a linear system using the chosen solve type.
unsigned int number() const
Gets the number of this system.
Definition: SystemBase.C:1159
virtual void close()=0
void computeResidual(libMesh::NonlinearImplicitSystem &sys, const NumericVector< libMesh::Number > &soln, NumericVector< libMesh::Number > &residual)
This function is called by Libmesh to form a residual.
const NumericVector< Number > *const & _solution
virtual void solve() override
Solves the time step and sets the number of nonlinear and linear iterations.
ActuallyExplicitEuler(const InputParameters &parameters)
void computeTimeDerivativeHelper(T &u_dot, const T2 &u_old) const
Helper function that actually does the math for computing the time derivative.
void mooseError(Args &&... args) const
Emits an error prefixed with object name and type.
std::unique_ptr< NumericVector< Number > > current_local_solution
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...
Implements a truly explicit (no nonlinear solve) first-order, forward Euler time integration scheme...
virtual Real & timeOld() const
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
unsigned int _n_nonlinear_iterations
Total number of nonlinear iterations over all stages of the time step.
NumericVector< Number > * _Re_non_time
residual vector for non-time contributions
NumericVector< Number > & solutionOld()
Definition: SystemBase.h:196
int & _t_step
The current time step number.
virtual void computeTimeDerivatives() override
Computes the time derivative and the Jacobian of the time derivative.
const DofMap & get_dof_map() const
const SparseMatrix< Number > & get_system_matrix() const
Base class for explicit time integrators that are implemented without using a nonlinear solver...
uint8_t dof_id_type
void computeDuDotDu()
Compute _du_dot_du.
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
void enforce_constraints_exactly(const System &system, NumericVector< Number > *v=nullptr, bool homogeneous=false) const