www.mooseframework.org
LStableDirk4.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 "LStableDirk4.h"
16 #include "NonlinearSystemBase.h"
17 #include "FEProblem.h"
18 #include "PetscSupport.h"
19 
20 template <>
23 {
25  return params;
26 }
27 
28 // Initialize static data
29 const Real LStableDirk4::_c[LStableDirk4::_n_stages] = {.25, 0., .5, 1., 1.};
30 
32  {.25, 0, 0, 0, 0},
33  {-.25, .25, 0, 0, 0},
34  {.125, .125, .25, 0, 0},
35  {-1.5, .75, 1.5, .25, 0},
36  {0, 1. / 6, 2. / 3, -1. / 12, .25}};
37 
39  : TimeIntegrator(parameters), _stage(1)
40 {
41  // Name the stage residuals "residual_stage1", "residual_stage2", etc.
42  for (unsigned int stage = 0; stage < _n_stages; ++stage)
43  {
44  std::ostringstream oss;
45  oss << "residual_stage" << stage + 1;
46  _stage_residuals[stage] = &(_nl.addVector(oss.str(), false, GHOSTED));
47  }
48 }
49 
51 
52 void
54 {
55  // We are multiplying by the method coefficients in postStep(), so
56  // the time derivatives are of the same form at every stage although
57  // the current solution varies depending on the stage.
58  _u_dot = *_solution;
60  _u_dot *= 1. / _dt;
61  _u_dot.close();
62  _du_dot_du = 1. / _dt;
63 }
64 
65 void
67 {
68  // Time at end of step
69  Real time_old = _fe_problem.timeOld();
70 
71  // A for-loop would increment _stage too far, so we use an extra
72  // loop counter.
73  for (unsigned int current_stage = 1; current_stage <= _n_stages; ++current_stage)
74  {
75  // Set the current stage value
76  _stage = current_stage;
77 
78  // This ensures that all the Output objects in the OutputWarehouse
79  // have had solveSetup() called, and sets the default solver
80  // parameters for PETSc.
82 
83  _console << "Stage " << _stage << "\n";
84 
85  // Set the time for this stage
86  _fe_problem.time() = time_old + _c[_stage - 1] * _dt;
87 
88  // Do the solve
90  }
91 }
92 
93 void
94 LStableDirk4::postStep(NumericVector<Number> & residual)
95 {
96  // Error if _stage got messed up somehow.
97  if (_stage > _n_stages)
98  // the explicit cast prevents strange compiler weirdness with the static
99  // const variable and the variadic mooseError function
100  mooseError("LStableDirk4::postStep(): Member variable _stage can only have values 1-",
101  (unsigned int)_n_stages,
102  ".");
103 
104  // In the standard RK notation, the residual of stage 1 of s is given by:
105  //
106  // R := M*(Y_i - y_n)/dt - \sum_{j=1}^s a_{ij} * f(t_n + c_j*dt, Y_j) = 0
107  //
108  // where:
109  // .) M is the mass matrix
110  // .) Y_i is the stage solution
111  // .) dt is the timestep, and is accounted for in the _Re_time residual.
112  // .) f are the "non-time" residuals evaluated for a given stage solution.
113  // .) The minus signs are already "baked in" to the residuals and so do not appear below.
114 
115  // Store this stage's non-time residual. We are calling operator=
116  // here, and that calls close().
118 
119  // Build up the residual for this stage.
120  residual.add(1., _Re_time);
121  for (unsigned int j = 0; j < _stage; ++j)
122  residual.add(_a[_stage - 1][j], *_stage_residuals[j]);
123  residual.close();
124 }
NonlinearSystemBase & _nl
virtual ~LStableDirk4()
Definition: LStableDirk4.C:50
virtual Real & timeOld() const
virtual void postStep(NumericVector< Number > &residual)
Callback to the TimeIntegrator called immediately after the residuals are computed in NonlinearSystem...
Definition: LStableDirk4.C:94
void initPetscOutput()
Reinitialize petsc output for proper linear/nonlinear iteration display.
static const Real _a[_n_stages][_n_stages]
Definition: LStableDirk4.h:91
NumericVector< Number > & _u_dot
solution vector for u^dot
NonlinearSystemBase & getNonlinearSystemBase()
virtual NumericVector< Number > & addVector(const std::string &vector_name, const bool project, const ParallelType type)
Adds a solution length vector to the system.
Definition: SystemBase.C:490
FEProblemBase & _fe_problem
LStableDirk4(const InputParameters &parameters)
Definition: LStableDirk4.C:38
NumericVector< Number > & _Re_non_time
residual vector for non-time contributions
virtual Real & time() const
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
unsigned int _stage
Definition: LStableDirk4.h:75
virtual void computeTimeDerivatives()
Definition: LStableDirk4.C:53
InputParameters validParams< TimeIntegrator >()
virtual void solve()
Definition: LStableDirk4.C:66
static const Real _c[_n_stages]
Definition: LStableDirk4.h:87
Real & _du_dot_du
solution vector for
Base class for time integrators.
virtual System & system() override
Get the reference to the libMesh system.
NumericVector< Number > * _stage_residuals[_n_stages]
Definition: LStableDirk4.h:84
NumericVector< Number > & _Re_time
residual vector for time contributions
static const unsigned int _n_stages
Definition: LStableDirk4.h:81
void mooseError(Args &&...args) const
Definition: MooseObject.h:80
const NumericVector< Number > & _solution_old
const ConsoleStream _console
An instance of helper class to write streams to the Console objects.
const NumericVector< Number > *& _solution
solution vectors
InputParameters validParams< LStableDirk4 >()
Definition: LStableDirk4.C:22