www.mooseframework.org
AB2PredictorCorrector.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 "AB2PredictorCorrector.h"
16 #include "AdamsPredictor.h"
17 #include "Problem.h"
18 #include "FEProblem.h"
19 #include "MooseApp.h"
20 #include "NonlinearSystem.h"
21 #include "AuxiliarySystem.h"
22 #include "TimeIntegrator.h"
23 #include "Conversion.h"
24 
25 #include "libmesh/nonlinear_solver.h"
26 #include "libmesh/numeric_vector.h"
27 
28 // C++ Includes
29 #include <iomanip>
30 #include <iostream>
31 #include <fstream>
32 
33 template <>
36 {
38  params.addRequiredParam<Real>("e_tol", "Target error tolerance.");
39  params.addRequiredParam<Real>("e_max", "Maximum acceptable error.");
40  params.addRequiredParam<Real>("dt", "Initial time step size");
41  params.addParam<Real>("max_increase", 1.0e9, "Maximum ratio that the time step can increase.");
42  params.addParam<int>(
43  "steps_between_increase", 1, "the number of time steps before recalculating dt");
44  params.addParam<int>("start_adapting", 2, "when to start taking adaptive time steps");
45  params.addParam<Real>("scaling_parameter", .8, "scaling parameter for dt selection");
46  return params;
47 }
48 
50  : TimeStepper(parameters),
51  _u1(_fe_problem.getNonlinearSystemBase().addVector("u1", true, GHOSTED)),
52  _aux1(_fe_problem.getAuxiliarySystem().addVector("aux1", true, GHOSTED)),
53  _pred1(_fe_problem.getNonlinearSystemBase().addVector("pred1", true, GHOSTED)),
54  _dt_full(declareRestartableData<Real>("dt_full", 0)),
55  _error(declareRestartableData<Real>("error", 0)),
56  _e_tol(getParam<Real>("e_tol")),
57  _e_max(getParam<Real>("e_max")),
58  _max_increase(getParam<Real>("max_increase")),
59  _steps_between_increase(getParam<int>("steps_between_increase")),
60  _dt_steps_taken(declareRestartableData<int>("dt_steps_taken", 0)),
61  _start_adapting(getParam<int>("start_adapting")),
62  _my_dt_old(declareRestartableData<Real>("my_dt_old", 0)),
63  _infnorm(declareRestartableData<Real>("infnorm", 0)),
64  _scaling_parameter(getParam<Real>("scaling_parameter"))
65 {
66  Real predscale = 1.;
67  InputParameters params = _app.getFactory().getValidParams("AdamsPredictor");
68  params.set<Real>("scale") = predscale;
69  _fe_problem.addPredictor("AdamsPredictor", "adamspredictor", params);
70 }
71 
72 void
74 {
76 }
77 
78 void
80 {
81  // save dt
82  _dt_full = _dt;
83 }
84 
85 void
87 {
90 
93  if (_converged)
94  {
95  _u1 = *nl.currentSolution();
96  _u1.close();
97 
98  _aux1 = *aux.currentSolution();
99  _aux1.close();
100  if (_t_step >= _start_adapting)
101  {
102  // Calculate error if past the first solve
104 
105  _infnorm = _u1.linfty_norm();
106  _e_max = 1.1 * _e_tol * _infnorm;
107  _console << "Time Error Estimate: " << _error << std::endl;
108  }
109  else
110  {
111  // First time step is problematic, sure we converged but what does that mean? We don't know.
112  // Nor can we calculate the error on the first time step.
113  }
114  }
115 }
116 
117 bool
119 {
120  if (!_converged)
121  {
122  _dt_steps_taken = 0;
123  return false;
124  }
125  if (_error < _e_max)
126  {
127  return true;
128  }
129  else
130  {
131  _console << "Marking last solve not converged " << _error << " " << _e_max << std::endl;
132  _dt_steps_taken = 0;
133  return false;
134  }
135 }
136 
137 Real
139 {
140  if (_t_step <= _start_adapting)
141  return _dt;
142 
143  _my_dt_old = _dt;
144 
145  _dt_steps_taken += 1;
147  {
148 
149  Real new_dt = _dt_full * _scaling_parameter * std::pow(_infnorm * _e_tol / _error, 1.0 / 3.0);
150 
151  if (new_dt / _dt_full > _max_increase)
152  new_dt = _dt_full * _max_increase;
153  _dt_steps_taken = 0;
154  return new_dt;
155  }
156 
157  return _dt;
158 }
159 
160 Real
162 {
163  return getParam<Real>("dt");
164 }
165 
166 Real
167 AB2PredictorCorrector::estimateTimeError(NumericVector<Number> & solution)
168 {
171  auto scheme = Moose::stringToEnum<Moose::TimeIntegratorType>(ti->name());
172  Real dt_old = _my_dt_old;
173  if (dt_old == 0)
174  dt_old = _dt;
175 
176  switch (scheme)
177  {
179  {
180  _pred1 *= -1;
181  _pred1 += solution;
182  Real calc = _dt * _dt * .5;
183  _pred1 *= calc;
184  return _pred1.l2_norm();
185  }
187  {
188  _pred1 -= solution;
189  _pred1 *= (_dt) / (3.0 * (_dt + dt_old));
190  return _pred1.l2_norm();
191  }
192  case Moose::TI_BDF2:
193  {
194  _pred1 *= -1.0;
195  _pred1 += solution;
196  Real topcalc = 2.0 * (_dt + dt_old) * (_dt + dt_old);
197  Real bottomcalc = 6.0 * _dt * _dt + 12.0 * _dt * dt_old + 5.0 * dt_old * dt_old;
198  _pred1 *= topcalc / bottomcalc;
199 
200  return _pred1.l2_norm();
201  }
202  default:
203  break;
204  }
205  return -1;
206 }
const std::string & name() const
Get the name of the object.
Definition: MooseObject.h:47
virtual const NumericVector< Number > *& currentSolution() override
The solution vector that is currently being operated on.
AB2PredictorCorrector(const InputParameters &parameters)
Real & _error
global relative time discretization error estimate
virtual void addPredictor(const std::string &type, const std::string &name, InputParameters parameters)
virtual Real computeInitialDT() override
Called to compute _current_dt for the first timestep.
int _steps_between_increase
steps to take before increasing dt
TimeIntegrator * getTimeIntegrator()
NonlinearSystemBase & getNonlinearSystemBase()
InputParameters getValidParams(const std::string &name)
Get valid parameters for the object.
Definition: Factory.C:26
Real _e_tol
error tolerance
NumericVector< Number > & _aux1
T & set(const std::string &name, bool quiet_mode=false)
Returns a writable reference to the named parameters.
Base class for time stepping.
Definition: TimeStepper.h:31
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
NumericVector< Number > & _pred1
virtual void preSolve() override
virtual const NumericVector< Number > *& currentSolution() override
The solution vector that is currently being operated on.
Real & _dt_full
dt of the big step
void addRequiredParam(const std::string &name, const std::string &doc_string)
This method adds a parameter and documentation string to the InputParameters object that will be extr...
Factory & getFactory()
Retrieve the Factory associated with this App.
Definition: MooseApp.h:253
virtual void step() override
Take a time step.
NonlinearSystemBase * nl
Nonlinear system to be solved.
FEProblemBase & _fe_problem
Definition: TimeStepper.h:124
virtual bool converged() override
virtual Real estimateTimeError(NumericVector< Number > &sol)
Real _e_max
maximal error
Real _max_increase
maximum increase ratio
Real _scaling_parameter
scaling_parameter for time step selection, default is 0.8
AuxiliarySystem & getAuxiliarySystem()
bool _converged
Whether or not the previous solve converged.
Definition: TimeStepper.h:144
NumericVector< Number > & _u1
InputParameters validParams< TimeStepper >()
Definition: TimeStepper.C:22
Base class for time integrators.
virtual NumericVector< Number > & solutionPredictor()
Definition: Predictor.h:50
MooseApp & _app
The MooseApp this object is associated with.
Definition: MooseObject.h:108
virtual void preExecute()
Definition: TimeStepper.C:65
Real & _infnorm
infinity norm of the solution vector
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...
int & _t_step
Definition: TimeStepper.h:131
virtual bool converged() override
If the time step converged.
virtual Real computeDT() override
Called to compute _current_dt for a normal step.
const ConsoleStream _console
An instance of helper class to write streams to the Console objects.
virtual void preExecute() override
virtual void solve() override
InputParameters validParams< AB2PredictorCorrector >()
int & _dt_steps_taken
steps taken at current dt
Real & _dt
Definition: TimeStepper.h:132
A system that holds auxiliary variables.