www.mooseframework.org
AStableDirk4.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 // MOOSE includes
16 #include "AStableDirk4.h"
17 #include "NonlinearSystem.h"
18 #include "FEProblem.h"
19 #include "PetscSupport.h"
20 #include "LStableDirk4.h"
21 
22 template <>
25 {
27  params.addParam<bool>("safe_start", true, "If true, use LStableDirk4 to bootstrap this method.");
28  return params;
29 }
30 
32  : TimeIntegrator(parameters),
33  _stage(1),
34  _gamma(0.5 + std::sqrt(3) / 3. * std::cos(libMesh::pi / 18.)),
35  _safe_start(getParam<bool>("safe_start"))
36 {
37  // Name the stage residuals "residual_stage1", "residual_stage2", etc.
38  for (unsigned int stage = 0; stage < 3; ++stage)
39  {
40  std::ostringstream oss;
41  oss << "residual_stage" << stage + 1;
42  _stage_residuals[stage] = &(_nl.addVector(oss.str(), false, GHOSTED));
43  }
44 
45  // Initialize parameters
46  _c[0] = _gamma;
47  _c[1] = .5;
48  _c[2] = 1.0 - _gamma;
49 
50  _a[0][0] = _gamma;
51  _a[1][0] = .5 - _gamma;
52  _a[1][1] = _gamma;
53  _a[2][0] = 2. * _gamma;
54  _a[2][1] = 1 - 4. * _gamma;
55  _a[2][2] = _gamma;
56 
57  _b[0] = 1. / (24. * (.5 - _gamma) * (.5 - _gamma));
58  _b[1] = 1. - 1. / (12. * (.5 - _gamma) * (.5 - _gamma));
59  _b[2] = _b[0];
60 
61  // If doing a _safe_start, construct the bootstrapping
62  // TimeIntegrator. Note that this method will also add
63  // residual_stage vectors to the system, but since they have the
64  // same name as the vectors we added, they won't be duplicated.
65  if (_safe_start)
66  {
67  Factory & factory = _app.getFactory();
68  InputParameters params = factory.getValidParams("LStableDirk4");
69 
70  // We need to set some parameters that are normally set in
71  // FEProblemBase::addTimeIntegrator() to ensure that the
72  // getCheckedPointerParam() sanity checking is happy. This is why
73  // constructing MOOSE objects "manually" is generally frowned upon.
74  params.set<FEProblemBase *>("_fe_problem_base") = &_fe_problem;
75  params.set<SystemBase *>("_sys") = &_sys;
76 
77  _bootstrap_method = factory.create<LStableDirk4>("LStableDirk4", name() + "_bootstrap", params);
78  }
79 }
80 
82 
83 void
85 {
86  // We are multiplying by the method coefficients in postStep(), so
87  // the time derivatives are of the same form at every stage although
88  // the current solution varies depending on the stage.
89  _u_dot = *_solution;
91  _u_dot *= 1. / _dt;
92  _u_dot.close();
93  _du_dot_du = 1. / _dt;
94 }
95 
96 void
98 {
99  if (_t_step == 1 && _safe_start)
100  _bootstrap_method->solve();
101 
102  else
103  {
104  // Time at end of step
105  Real time_old = _fe_problem.timeOld();
106 
107  // A for-loop would increment _stage too far, so we use an extra
108  // loop counter. The method has three stages and an update stage,
109  // which we treat as just one more (special) stage in the implementation.
110  for (unsigned int current_stage = 1; current_stage < 5; ++current_stage)
111  {
112  // Set the current stage value
113  _stage = current_stage;
114 
115  // This ensures that all the Output objects in the OutputWarehouse
116  // have had solveSetup() called, and sets the default solver
117  // parameters for PETSc.
119 
120  if (current_stage < 4)
121  {
122  _console << "Stage " << _stage << "\n";
123  _fe_problem.time() = time_old + _c[_stage - 1] * _dt;
124  }
125  else
126  {
127  _console << "Update Stage.\n";
128  _fe_problem.time() = time_old + _dt;
129  }
130 
131  // Do the solve
133  }
134  }
135 }
136 
137 void
138 AStableDirk4::postStep(NumericVector<Number> & residual)
139 {
140  if (_t_step == 1 && _safe_start)
141  _bootstrap_method->postStep(residual);
142 
143  else
144  {
145  // Error if _stage got messed up somehow.
146  if (_stage > 4)
147  mooseError("AStableDirk4::postStep(): Member variable _stage can only have values 1-4.");
148 
149  if (_stage < 4)
150  {
151  // In the standard RK notation, the residual of stage 1 of s is given by:
152  //
153  // R := M*(Y_i - y_n)/dt - \sum_{j=1}^s a_{ij} * f(t_n + c_j*dt, Y_j) = 0
154  //
155  // where:
156  // .) M is the mass matrix
157  // .) Y_i is the stage solution
158  // .) dt is the timestep, and is accounted for in the _Re_time residual.
159  // .) f are the "non-time" residuals evaluated for a given stage solution.
160  // .) The minus signs are already "baked in" to the residuals and so do not appear below.
161 
162  // Store this stage's non-time residual. We are calling operator=
163  // here, and that calls close().
165 
166  // Build up the residual for this stage.
167  residual.add(1., _Re_time);
168  for (unsigned int j = 0; j < _stage; ++j)
169  residual.add(_a[_stage - 1][j], *_stage_residuals[j]);
170  residual.close();
171  }
172  else
173  {
174  // The update step is a final solve:
175  //
176  // R := M*(y_{n+1} - y_n)/dt - \sum_{j=1}^s b_j * f(t_n + c_j*dt, Y_j) = 0
177  //
178  // We could potentially fold _b up into an extra row of _a and
179  // just do one more stage, but I think handling it separately like
180  // this is easier to understand and doesn't create too much code
181  // repitition.
182  residual.add(1., _Re_time);
183  for (unsigned int j = 0; j < 3; ++j)
184  residual.add(_b[j], *_stage_residuals[j]);
185  residual.close();
186  }
187  }
188 }
const std::string & name() const
Get the name of the object.
Definition: MooseObject.h:47
NonlinearSystemBase & _nl
virtual void computeTimeDerivatives()
Definition: AStableDirk4.C:84
Fourth-order diagonally implicit Runge Kutta method (Dirk) with five stages.
Definition: LStableDirk4.h:62
virtual Real & timeOld() const
std::shared_ptr< MooseObject > create(const std::string &obj_name, const std::string &name, InputParameters parameters, THREAD_ID tid=0, bool print_deprecated=true)
Build an object (must be registered) - THIS METHOD IS DEPRECATED (Use create<T>()) ...
Definition: Factory.C:46
void initPetscOutput()
Reinitialize petsc output for proper linear/nonlinear iteration display.
std::shared_ptr< LStableDirk4 > _bootstrap_method
Definition: AStableDirk4.h:104
NumericVector< Number > & _u_dot
solution vector for u^dot
NonlinearSystemBase & getNonlinearSystemBase()
Generic factory class for build all sorts of objects.
Definition: Factory.h:152
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
InputParameters getValidParams(const std::string &name)
Get valid parameters for the object.
Definition: Factory.C:26
FEProblemBase & _fe_problem
T & set(const std::string &name, bool quiet_mode=false)
Returns a writable reference to the named parameters.
NumericVector< Number > & _Re_non_time
residual vector for non-time contributions
virtual Real & time() const
SystemBase & _sys
Real _a[3][3]
Definition: AStableDirk4.h:92
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
AStableDirk4(const InputParameters &parameters)
Definition: AStableDirk4.C:31
virtual void postStep(NumericVector< Number > &residual)
Callback to the TimeIntegrator called immediately after the residuals are computed in NonlinearSystem...
Definition: AStableDirk4.C:138
The following methods are specializations for using the libMesh::Parallel::packed_range_* routines fo...
Base class for a system (of equations)
Definition: SystemBase.h:91
Specialization of SubProblem for solving nonlinear equations plus auxiliary equations.
Factory & getFactory()
Retrieve the Factory associated with this App.
Definition: MooseApp.h:253
Real _b[3]
Definition: AStableDirk4.h:96
const Real _gamma
Definition: AStableDirk4.h:81
InputParameters validParams< AStableDirk4 >()
Definition: AStableDirk4.C:24
unsigned int _stage
Definition: AStableDirk4.h:75
NumericVector< Number > * _stage_residuals[3]
Definition: AStableDirk4.h:78
InputParameters validParams< TimeIntegrator >()
Real & _du_dot_du
solution vector for
Base class for time integrators.
virtual void solve()
Definition: AStableDirk4.C:97
virtual System & system() override
Get the reference to the libMesh system.
MooseApp & _app
The MooseApp this object is associated with.
Definition: MooseObject.h:108
virtual ~AStableDirk4()
Definition: AStableDirk4.C:81
Real _c[3]
Definition: AStableDirk4.h:85
NumericVector< Number > & _Re_time
residual vector for time contributions
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...
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