www.mooseframework.org
ActuallyExplicitEuler.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 // 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 
19 
22 {
24 
25  params.addClassDescription(
26  "Implementation of Explicit/Forward Euler without invoking any of the nonlinear solver");
27 
28  params.addParam<bool>("use_constant_mass",
29  false,
30  "If set to true, will only compute the mass matrix in the first time step, "
31  "and keep using it throughout the simulation.");
32 
33  return params;
34 }
35 
37  : ExplicitTimeIntegrator(parameters), _constant_mass(getParam<bool>("use_constant_mass"))
38 {
39 }
40 
41 void
43 {
44  if (!_sys.solutionUDot())
45  mooseError("ActuallyExplicitEuler: Time derivative of solution (`u_dot`) is not stored. Please "
46  "set uDotRequested() to true in FEProblemBase before requesting `u_dot`.");
47 
49  u_dot = *_solution;
51  u_dot.close();
52 
53  _du_dot_du = 1.0 / _dt;
54 }
55 
56 void
58  const dof_id_type & dof,
59  DualReal & /*ad_u_dotdot*/) const
60 {
62 }
63 
64 void
66 {
67  // Reset iteration counts
70 
72 
73  // Set time to the time at which to evaluate the residual
76 
77  // Compute the residual
78  _explicit_residual.zero();
80  *_nonlinear_implicit_system->current_local_solution, _explicit_residual, _nl.number());
81 
82  // Move the residual to the RHS
83  _explicit_residual *= -1.0;
84 
85  // Compute the mass matrix
86  auto & mass_matrix = _nonlinear_implicit_system->get_system_matrix();
87  if (!_constant_mass || (_constant_mass && _t_step == 1))
89  *_nonlinear_implicit_system->current_local_solution, mass_matrix, _Ke_time_tag);
90 
91  // Perform the linear solve
92  bool converged = performExplicitSolve(mass_matrix);
93 
94  // Update the solution
97 
98  // Constraints may be solved in an uncoupled way. For example, momentum-balance equations may be
99  // solved node-wise and then the solution (e.g. velocities or positions)can be applied to those
100  // nodes without solving for such constraints on a system level. This strategy is being used for
101  // node-face contact in explicit dynamics.
103 
104  // Enforce contraints on the solution
105  DofMap & dof_map = _nonlinear_implicit_system->get_dof_map();
106  dof_map.enforce_constraints_exactly(*_nonlinear_implicit_system,
107  _nonlinear_implicit_system->solution.get());
108  _nonlinear_implicit_system->update();
109 
110  _nl.setSolution(*_nonlinear_implicit_system->current_local_solution);
111 
112  _nonlinear_implicit_system->nonlinear_solver->converged = converged;
113 }
114 
115 void
117 {
118  residual += _Re_time;
119  residual += _Re_non_time;
120  residual.close();
121 
122  // Reset time to the time at which to evaluate nodal BCs, which comes next
124 }
NonlinearSystemBase & _nl
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
virtual NumericVector< Number > * solutionUDot()=0
void overwriteNodeFace(NumericVector< Number > &soln)
Called from explicit time stepping to overwrite boundary positions (explicit dynamics).
virtual Real & time() const
Real _current_time
Save off current time to reset it back and forth.
virtual void postResidual(NumericVector< Number > &residual) override
Callback to the TimeIntegrator called immediately after the residuals are computed in NonlinearSystem...
virtual void setSolution(const NumericVector< Number > &soln)
FEProblemBase & _fe_problem
DualNumber< Real, DNDerivativeType, true > DualReal
Definition: DualReal.h:49
NumericVector< Number > & _Re_non_time
residual vector for non-time contributions
TagID _Ke_time_tag
For computing the mass matrix.
SystemBase & _sys
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
static InputParameters validParams()
virtual void computeJacobianTag(const NumericVector< Number > &soln, SparseMatrix< Number > &jacobian, TagID tag)
Form a Jacobian matrix for a given tag.
static InputParameters validParams()
NumericVector< Real > & _solution_update
Solution vector for the linear solve.
unsigned int _n_linear_iterations
Total number of linear iterations over all stages of the time step.
registerMooseObject("MooseApp", ActuallyExplicitEuler)
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:1125
virtual void close()=0
Real & _du_dot_du
Derivative of time derivative with respect to current solution: .
const NumericVector< Number > *const & _solution
solution vectors
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.
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...
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 option 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.
void computeResidual(NonlinearImplicitSystem &sys, const NumericVector< Number > &soln, NumericVector< Number > &residual)
This function is called by Libmesh to form a residual.
NumericVector< Number > & solutionOld()
Definition: SystemBase.h:177
virtual void computeTimeDerivatives() override
Computes the time derivative and the Jacobian of the time derivative.
NumericVector< Real > & _explicit_residual
Residual used for the RHS.
Base class for explicit time integrators that are implemented without using a nonlinear solver...
uint8_t dof_id_type
NonlinearImplicitSystem * _nonlinear_implicit_system
Nonlinear implicit system, if applicable; otherwise, nullptr.