www.mooseframework.org
ActuallyExplicitEuler.h
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 #pragma once
11 
12 #include "TimeIntegrator.h"
13 #include "MeshChangedInterface.h"
14 
15 #include "libmesh/linear_solver.h"
16 
17 // Forward declarations
20 
21 template <>
23 
29 {
30 public:
32 
33  virtual void initialSetup() override;
34  virtual void init() override;
35  virtual void preSolve() override;
36  virtual int order() override { return 1; }
37  virtual void computeTimeDerivatives() override;
38  void computeADTimeDerivatives(DualReal & ad_u_dot, const dof_id_type & dof) const override;
39  virtual void solve() override;
40  virtual void postResidual(NumericVector<Number> & residual) override;
41 
42  virtual void meshChanged() override;
43 
44 protected:
45  enum SolveType
46  {
50  };
51 
56 
60  template <typename T, typename T2>
61  void computeTimeDerivativeHelper(T & u_dot, const T2 & u_old) const;
62 
64 
66  NumericVector<Real> & _explicit_residual;
67 
69  NumericVector<Real> & _explicit_euler_update;
70 
72  NumericVector<Real> & _mass_matrix_diag;
73 
75  NumericVector<Real> * _ones;
76 
79 
81  std::unique_ptr<LinearSolver<Number>> _linear_solver;
82 
84  std::unique_ptr<LumpedPreconditioner> _preconditioner;
85 
88 };
89 
90 template <typename T, typename T2>
91 void
92 ActuallyExplicitEuler::computeTimeDerivativeHelper(T & u_dot, const T2 & u_old) const
93 {
94  u_dot -= u_old;
95  u_dot *= 1. / _dt;
96 }
97 
virtual void initialSetup() override
Called to setup datastructures.
virtual void postResidual(NumericVector< Number > &residual) override
Callback to the TimeIntegrator called immediately after the residuals are computed in NonlinearSystem...
unsigned int TagID
Definition: MooseTypes.h:162
NumericVector< Real > & _mass_matrix_diag
Diagonal of the lumped mass matrix (and its inversion)
NumericVector< Real > * _ones
Just a vector of 1&#39;s to help with creating the lumped mass matrix.
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
DualNumber< Real, NumberArray< AD_MAX_DOFS_PER_ELEM, Real > > DualReal
Definition: DualReal.h:29
bool checkLinearConvergence()
Check for the linear solver convergence.
void computeADTimeDerivatives(DualReal &ad_u_dot, const dof_id_type &dof) const override
method for computing local automatic differentiation time derivatives
Interface for notifications that the mesh has changed.
const InputParameters & parameters() const
Get the parameters of the object.
Definition: MooseObject.h:65
Real _current_time
Save off current time to reset it back and forth.
virtual void init() override
Called only before the very first timestep (t_step = 0) Never called again (not even during recover/r...
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:31
NumericVector< Real > & _explicit_residual
Residual used for the RHS.
virtual void preSolve() override
std::unique_ptr< LumpedPreconditioner > _preconditioner
For solving with lumped preconditioning.
InputParameters validParams< ActuallyExplicitEuler >()
virtual int order() override
virtual void solve() override
Solves the time step and sets the number of nonlinear and linear iterations.
Base class for time integrators.
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.
Implements a truly explicit (no nonlinear solve) first-order, forward Euler time integration scheme...
std::unique_ptr< LinearSolver< Number > > _linear_solver
For solving with the consistent matrix.
virtual void computeTimeDerivatives() override
virtual void meshChanged() override
Called on this object when the mesh changes.
NumericVector< Real > & _explicit_euler_update
Solution vector for the linear solve.
TagID _Ke_time_tag
For computing the mass matrix.
Helper class to apply preconditioner.