www.mooseframework.org
ExplicitTVDRK2.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 "ExplicitTVDRK2.h"
16 #include "NonlinearSystemBase.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)
55  {
57  _u_dot *= 1. / _dt;
58  }
59  else
60  {
61  _u_dot.scale(2.);
64  _u_dot *= 0.5 / _dt;
65  }
66 
67  _du_dot_du = 1. / _dt;
68  _u_dot.close();
69 }
70 
71 void
73 {
74  Real time_new = _fe_problem.time();
75  Real time_old = _fe_problem.timeOld();
76  Real time_stage2 = time_old + _dt;
77 
78  // There is no work to do for the first stage (Y_1 = y_n). The
79  // first solve therefore happens in the second stage. Note that the
80  // non-time Kernels (which should be marked implicit=false) are
81  // evaluated at the old solution during this stage.
83  _console << "1st solve\n";
84  _stage = 2;
85  _fe_problem.timeOld() = time_old;
86  _fe_problem.time() = time_stage2;
88 
89  // Advance solutions old->older, current->old. Also moves Material
90  // properties and other associated state forward in time.
92 
93  // The "update" stage (which we call stage 3) requires an additional
94  // solve with the mass matrix.
96  _console << "2nd solve\n";
97  _stage = 3;
98  _fe_problem.timeOld() = time_stage2;
99  _fe_problem.time() = time_new;
101 
102  // Reset time at beginning of step to its original value
103  _fe_problem.timeOld() = time_old;
104 }
105 
106 void
107 ExplicitTVDRK2::postStep(NumericVector<Number> & residual)
108 {
109  if (_stage == 1)
110  {
111  // If postStep() is called before solve(), _stage==1 and we don't
112  // need to do anything.
113  }
114  else if (_stage == 2)
115  {
116  // In the standard RK notation, the stage 2 residual is given by:
117  //
118  // R := M*(Y_2 - y_n)/dt - f(t_n, Y_1) = 0
119  //
120  // where:
121  // .) M is the mass matrix.
122  // .) f(t_n, Y_1) is the residual we are currently computing,
123  // since this method is intended to be used with "implicit=false"
124  // kernels.
125  // .) M*(Y_2 - y_n)/dt corresponds to the residual of the time kernels.
126  // .) The minus signs are "baked in" to the non-time residuals, so
127  // they do not appear here.
128  // .) The current non-time residual is saved for the next stage.
130  _residual_old.close();
131 
132  residual.add(1.0, _Re_time);
133  residual.add(1.0, _residual_old);
134  residual.close();
135  }
136  else if (_stage == 3)
137  {
138  // In the standard RK notation, the update step residual is given by:
139  //
140  // R := M*(2*y_{n+1} - Y_2 - y_n)/(2*dt) - (1/2)*f(t_n+dt/2, Y_2) = 0
141  //
142  // where:
143  // .) M is the mass matrix.
144  // .) f(t_n+dt/2, Y_2) is the residual from stage 2.
145  // .) The minus sign is already baked in to the non-time
146  // residuals, so it does not appear here.
147  // .) Although this is an update step, we have to do a "solve"
148  // using the mass matrix.
149  residual.add(1.0, _Re_time);
150  residual.add(0.5, _Re_non_time);
151  residual.close();
152  }
153  else
154  mooseError("ExplicitTVDRK2::postStep(): _stage = ", _stage, ", only _stage = 1-3 is allowed.");
155 }
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()
FEProblemBase & _fe_problem
unsigned int _stage
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 void advanceState()
Advance all of the state holding vectors / datastructures so that we can move to the next timestep...
NumericVector< Number > & _residual_old
Buffer to store non-time residual from the first stage.
InputParameters validParams< TimeIntegrator >()
virtual void preSolve()
virtual void postStep(NumericVector< Number > &residual)
Callback to the TimeIntegrator called immediately after the residuals are computed in NonlinearSystem...
virtual ~ExplicitTVDRK2()
virtual void computeTimeDerivatives()
const NumericVector< Number > & _solution_older
virtual void solve()
Real & _du_dot_du
solution vector for
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.
ExplicitTVDRK2(const InputParameters &parameters)
NumericVector< Number > & _Re_time
residual vector for time contributions
void mooseError(Args &&...args) const
Definition: MooseObject.h:80
const NumericVector< Number > & _solution_old
InputParameters validParams< ExplicitTVDRK2 >()
const ConsoleStream _console
An instance of helper class to write streams to the Console objects.
const NumericVector< Number > *& _solution
solution vectors