www.mooseframework.org
CrankNicolson.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 "CrankNicolson.h"
11 #include "NonlinearSystem.h"
12 #include "FEProblem.h"
13 
15 
16 template <>
19 {
21 
22  return params;
23 }
24 
26  : TimeIntegrator(parameters), _residual_old(_nl.addVector("residual_old", false, GHOSTED))
27 {
28 }
29 
30 void
32 {
33  if (!_sys.solutionUDot())
34  mooseError("CrankNicolson: Time derivative of solution (`u_dot`) is not stored. Please set "
35  "uDotRequested() to true in FEProblemBase befor requesting `u_dot`.");
36 
37  NumericVector<Number> & u_dot = *_sys.solutionUDot();
38  u_dot = *_solution;
40  u_dot.close();
41 
42  _du_dot_du = 2. / _dt;
43 }
44 
45 void
46 CrankNicolson::computeADTimeDerivatives(DualReal & ad_u_dot, const dof_id_type & dof) const
47 {
49 }
50 
51 void
53 {
54  if (!_sys.solutionUDot())
55  mooseError("CrankNicolson: Time derivative of solution (`u_dot`) is not stored. Please set "
56  "uDotRequested() to true in FEProblemBase befor requesting `u_dot`.");
57 
58  // time derivative is assumed to be zero on initial
59  NumericVector<Number> & u_dot = *_sys.solutionUDot();
60  u_dot.zero();
61  _du_dot_du = 0;
62 
63  // compute residual for the initial time step
64  // Note: we can not directly pass _residual_old in computeResidualType because
65  // the function will call postResidual, which will cause _residual_old
66  // to be added on top of itself prohibited by PETSc.
67  // Objects executed on initial have been executed by FEProblem,
68  // so we can and should directly call NonlinearSystem residual evaluation.
70  _residual_old = _nl.RHS();
71 }
72 
73 void
74 CrankNicolson::postResidual(NumericVector<Number> & residual)
75 {
76  residual += _Re_time;
77  residual += _Re_non_time;
78  residual += _residual_old;
79 }
80 
81 void
83 {
84  // shift the residual in time
86 }
NonlinearSystemBase & _nl
virtual NumericVector< Number > * solutionUDot()=0
virtual void postResidual(NumericVector< Number > &residual) override
Callback to the TimeIntegrator called immediately after the residuals are computed in NonlinearSystem...
Definition: CrankNicolson.C:74
void computeTimeDerivativeHelper(T &u_dot, const T2 &u_old) const
Helper function that actually does the math for computing the time derivative.
Definition: CrankNicolson.h:51
virtual void init() override
Called only before the very first timestep (t_step = 0) Never called again (not even during recover/r...
Definition: CrankNicolson.C:52
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...
DualNumber< Real, NumberArray< AD_MAX_DOFS_PER_ELEM, Real > > DualReal
Definition: DualReal.h:29
virtual void postStep() override
Callback to the TimeIntegrator called at the very end of time step.
Definition: CrankNicolson.C:82
void mooseError(Args &&... args) const
Definition: MooseObject.h:147
NumericVector< Number > & _residual_old
Definition: CrankNicolson.h:46
virtual TagID nonTimeVectorTag() override
InputParameters validParams< TimeIntegrator >()
void computeADTimeDerivatives(DualReal &ad_u_dot, const dof_id_type &dof) const override
method for computing local automatic differentiation time derivatives
Definition: CrankNicolson.C:46
registerMooseObject("MooseApp", CrankNicolson)
CrankNicolson(const InputParameters &parameters)
Definition: CrankNicolson.C:25
Real & _du_dot_du
solution vector for
virtual NumericVector< Number > & RHS()=0
const NumericVector< Number > *const & _solution
solution vectors
Base class for time integrators.
NumericVector< Number > & _Re_time
residual vector for time contributions
const NumericVector< Number > & _solution_old
void computeResidualTag(NumericVector< Number > &residual, TagID tag_id)
Computes residual for a given tag.
virtual void computeTimeDerivatives() override
Definition: CrankNicolson.C:31
Crank-Nicolson time integrator.
Definition: CrankNicolson.h:27
InputParameters validParams< CrankNicolson >()
Definition: CrankNicolson.C:18