www.mooseframework.org
ExplicitRK2.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 "ExplicitRK2.h"
16 #include "NonlinearSystem.h"
17 #include "FEProblem.h"
18 #include "PetscSupport.h"
19 
20 template <>
23 {
25 
26  return params;
27 }
28 
30  : TimeIntegrator(parameters),
31  _stage(1),
32  _residual_old(_nl.addVector("residual_old", false, GHOSTED))
33 {
34 }
35 
37 
38 void
40 {
41  if (_dt == _dt_old)
43  else
45 }
46 
47 void
49 {
50  // Since advanceState() is called in between stages 2 and 3, this
51  // changes the meaning of "_solution_old". In the second stage,
52  // "_solution_older" is actually the original _solution_old.
53  _u_dot = *_solution;
54  if (_stage < 3)
56  else
58 
59  _u_dot *= 1. / _dt;
60  _du_dot_du = 1. / _dt;
61  _u_dot.close();
62 }
63 
64 void
66 {
67  Real time_new = _fe_problem.time();
68  Real time_old = _fe_problem.timeOld();
69  Real time_stage2 = time_old + a() * _dt;
70 
71  // There is no work to do for the first stage (Y_1 = y_n). The
72  // first solve therefore happens in the second stage. Note that the
73  // non-time Kernels (which should be marked implicit=false) are
74  // evaluated at the old solution during this stage.
76  _console << "1st solve\n";
77  _stage = 2;
78  _fe_problem.timeOld() = time_old;
79  _fe_problem.time() = time_stage2;
81 
82  // Advance solutions old->older, current->old. Also moves Material
83  // properties and other associated state forward in time.
85 
86  // The "update" stage (which we call stage 3) requires an additional
87  // solve with the mass matrix.
89  _console << "2nd solve\n";
90  _stage = 3;
91  _fe_problem.timeOld() = time_stage2;
92  _fe_problem.time() = time_new;
94 
95  // Reset time at beginning of step to its original value
96  _fe_problem.timeOld() = time_old;
97 }
98 
99 void
100 ExplicitRK2::postStep(NumericVector<Number> & residual)
101 {
102  if (_stage == 1)
103  {
104  // If postStep() is called before solve(), _stage==1 and we don't
105  // need to do anything.
106  }
107  else if (_stage == 2)
108  {
109  // In the standard RK notation, the stage 2 residual is given by:
110  //
111  // R := M*(Y_2 - y_n)/dt - a*f(t_n, Y_1) = 0
112  //
113  // where:
114  // .) M is the mass matrix.
115  // .) f(t_n, Y_1) is the residual we are currently computing,
116  // since this method is intended to be used with "implicit=false"
117  // kernels.
118  // .) M*(Y_2 - y_n)/dt corresponds to the residual of the time kernels.
119  // .) The minus signs are "baked in" to the non-time residuals, so
120  // they do not appear here.
121  // .) The current non-time residual is saved for the next stage.
123  _residual_old.close();
124 
125  residual.add(1., _Re_time);
126  residual.add(a(), _residual_old);
127  residual.close();
128  }
129  else if (_stage == 3)
130  {
131  // In the standard RK notation, the update step residual is given by:
132  //
133  // R := M*(y_{n+1} - y_n)/dt - f(t_n+dt/2, Y_2) = 0
134  //
135  // where:
136  // .) M is the mass matrix.
137  // .) f(t_n+dt/2, Y_2) is the residual from stage 2.
138  // .) The minus sign is already baked in to the non-time
139  // residuals, so it does not appear here.
140  // .) Although this is an update step, we have to do a "solve"
141  // using the mass matrix.
142  residual.add(1., _Re_time);
143  residual.add(b1(), _residual_old);
144  residual.add(b2(), _Re_non_time);
145  residual.close();
146  }
147  else
148  mooseError("ExplicitRK2::postStep(): _stage = ", _stage, ", only _stage = 1-3 is allowed.");
149 }
virtual Real b2() const =0
virtual Real & timeOld() const
void initPetscOutput()
Reinitialize petsc output for proper linear/nonlinear iteration display.
InputParameters validParams< ExplicitRK2 >()
Definition: ExplicitRK2.C:22
NumericVector< Number > & _u_dot
solution vector for u^dot
NonlinearSystemBase & getNonlinearSystemBase()
FEProblemBase & _fe_problem
virtual void solve()
Definition: ExplicitRK2.C:65
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...
virtual Real b1() const =0
virtual void computeTimeDerivatives()
Definition: ExplicitRK2.C:48
virtual void advanceState()
Advance all of the state holding vectors / datastructures so that we can move to the next timestep...
virtual ~ExplicitRK2()
Definition: ExplicitRK2.C:36
InputParameters validParams< TimeIntegrator >()
virtual void postStep(NumericVector< Number > &residual)
Callback to the TimeIntegrator called immediately after the residuals are computed in NonlinearSystem...
Definition: ExplicitRK2.C:100
unsigned int _stage
Definition: ExplicitRK2.h:79
virtual Real a() const =0
The method coefficients.
const NumericVector< Number > & _solution_older
Real & _du_dot_du
solution vector for
NumericVector< Number > & _residual_old
Buffer to store non-time residual from the first stage.
Definition: ExplicitRK2.h:82
Base class for time integrators.
void setConstJacobian(bool state)
Set flag that Jacobian is constant (for optimization purposes)
virtual System & system() override
Get the reference to the libMesh system.
ExplicitRK2(const InputParameters &parameters)
Definition: ExplicitRK2.C:29
NumericVector< Number > & _Re_time
residual vector for time contributions
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
virtual void preSolve()
Definition: ExplicitRK2.C:39