www.mooseframework.org
LStableDirk2.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* DO NOT MODIFY THIS HEADER */
3 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
4 /* */
5 /* (c) 2010 Battelle Energy Alliance, LLC */
6 /* ALL RIGHTS RESERVED */
7 /* */
8 /* Prepared by Battelle Energy Alliance, LLC */
9 /* Under Contract No. DE-AC07-05ID14517 */
10 /* With the U. S. Department of Energy */
11 /* */
12 /* See COPYRIGHT for full restrictions */
13 /****************************************************************/
14 
15 #include "LStableDirk2.h"
16 #include "NonlinearSystem.h"
17 #include "FEProblem.h"
18 #include "PetscSupport.h"
19 
20 template <>
23 {
25  return params;
26 }
27 
29  : TimeIntegrator(parameters),
30  _stage(1),
31  _residual_stage1(_nl.addVector("residual_stage1", false, GHOSTED)),
32  _residual_stage2(_nl.addVector("residual_stage2", false, GHOSTED)),
33  _alpha(1. - 0.5 * std::sqrt(2))
34 {
35 }
36 
38 
39 void
41 {
42  // We are multiplying by the method coefficients in postStep(), so
43  // the time derivatives are of the same form at every stage although
44  // the current solution varies depending on the stage.
45  _u_dot = *_solution;
47  _u_dot *= 1. / _dt;
48  _u_dot.close();
49  _du_dot_du = 1. / _dt;
50 }
51 
52 void
54 {
55  // Time at end of step
56  Real time_new = _fe_problem.time();
57 
58  // Time at beginning of step
59  Real time_old = _fe_problem.timeOld();
60 
61  // Time at stage 1
62  Real time_stage1 = time_old + _alpha * _dt;
63 
64  // Compute first stage
66  _console << "1st stage\n";
67  _stage = 1;
68  _fe_problem.time() = time_stage1;
70 
71  // Compute second stage
73  _console << "2nd stage\n";
74  _stage = 2;
75  _fe_problem.timeOld() = time_stage1;
76  _fe_problem.time() = time_new;
78 
79  // Reset time at beginning of step to its original value
80  _fe_problem.timeOld() = time_old;
81 }
82 
83 void
84 LStableDirk2::postStep(NumericVector<Number> & residual)
85 {
86  if (_stage == 1)
87  {
88  // In the standard RK notation, the stage 1 residual is given by:
89  //
90  // R := (Y_1 - y_n)/dt - alpha*f(t_n + alpha*dt, Y_1) = 0
91  //
92  // where:
93  // .) f(t_n + alpha*dt, Y_1) corresponds to the residuals of the
94  // non-time kernels. We'll save this as "_residual_stage1" to use
95  // later.
96  // .) (Y_1 - y_n)/dt corresponds to the residual of the time kernels.
97  // .) The minus sign in front of alpha is already "baked in" to
98  // the non-time residuals, so it does not appear here.
100  _residual_stage1.close();
101 
102  residual.add(1., _Re_time);
103  residual.add(_alpha, _residual_stage1);
104  residual.close();
105  }
106  else if (_stage == 2)
107  {
108  // In the standard RK notation, the stage 2 residual is given by:
109  //
110  // R := (Y_2 - y_n)/dt - (1-alpha)*f(t_n + alpha*dt, Y_1) - alpha*f(t_n + dt, Y_2) = 0
111  //
112  // where:
113  // .) f(t_n + alpha*dt, Y_1) has already been saved as _residual_stage1.
114  // .) f(t_n + dt, Y_2) will now be saved as "_residual_stage2".
115  // .) (Y_2 - y_n)/dt corresponds to the residual of the time kernels.
116  // .) The minus signs are once again "baked in" to the non-time
117  // residuals, so they do not appear here.
118  //
119  // The solution at the end of stage 2, i.e. Y_2, is also the final solution.
121  _residual_stage2.close();
122 
123  residual.add(1., _Re_time);
124  residual.add(1. - _alpha, _residual_stage1);
125  residual.add(_alpha, _residual_stage2);
126  residual.close();
127  }
128  else
129  mooseError("LStableDirk2::postStep(): Member variable _stage can only have values 1 or 2.");
130 }
virtual Real & timeOld() const
void initPetscOutput()
Reinitialize petsc output for proper linear/nonlinear iteration display.
NumericVector< Number > & _u_dot
solution vector for u^dot
NonlinearSystemBase & getNonlinearSystemBase()
InputParameters validParams< LStableDirk2 >()
Definition: LStableDirk2.C:22
FEProblemBase & _fe_problem
LStableDirk2(const InputParameters &parameters)
Definition: LStableDirk2.C:28
virtual ~LStableDirk2()
Definition: LStableDirk2.C:37
NumericVector< Number > & _Re_non_time
residual vector for non-time contributions
const Real _alpha
Definition: LStableDirk2.h:71
virtual Real & time() const
virtual void solve()
Definition: LStableDirk2.C:53
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
virtual void postStep(NumericVector< Number > &residual)
Callback to the TimeIntegrator called immediately after the residuals are computed in NonlinearSystem...
Definition: LStableDirk2.C:84
InputParameters validParams< TimeIntegrator >()
Real & _du_dot_du
solution vector for
Base class for time integrators.
virtual System & system() override
Get the reference to the libMesh system.
unsigned int _stage
Indicates the current stage (1 or 2).
Definition: LStableDirk2.h:62
NumericVector< Number > & _Re_time
residual vector for time contributions
void mooseError(Args &&...args) const
Definition: MooseObject.h:80
const NumericVector< Number > & _solution_old
virtual void computeTimeDerivatives()
Definition: LStableDirk2.C:40
NumericVector< Number > & _residual_stage2
Buffer to store non-time residual from second stage solve.
Definition: LStableDirk2.h:68
const ConsoleStream _console
An instance of helper class to write streams to the Console objects.
const NumericVector< Number > *& _solution
solution vectors
NumericVector< Number > & _residual_stage1
Buffer to store non-time residual from first stage solve.
Definition: LStableDirk2.h:65