www.mooseframework.org
AB2PredictorCorrector.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
10 #include "AB2PredictorCorrector.h"
11 #include "AdamsPredictor.h"
12 #include "Problem.h"
13 #include "FEProblem.h"
14 #include "MooseApp.h"
15 #include "NonlinearSystem.h"
16 #include "AuxiliarySystem.h"
17 #include "TimeIntegrator.h"
18 #include "Conversion.h"
19 
20 #include "libmesh/nonlinear_solver.h"
21 #include "libmesh/numeric_vector.h"
22 
23 // C++ Includes
24 #include <iomanip>
25 #include <iostream>
26 #include <fstream>
27 
29 
30 template <>
33 {
35  params.addRequiredParam<Real>("e_tol", "Target error tolerance.");
36  params.addRequiredParam<Real>("e_max", "Maximum acceptable error.");
37  params.addRequiredParam<Real>("dt", "Initial time step size");
38  params.addParam<Real>("max_increase", 1.0e9, "Maximum ratio that the time step can increase.");
39  params.addParam<int>(
40  "steps_between_increase", 1, "the number of time steps before recalculating dt");
41  params.addParam<int>("start_adapting", 2, "when to start taking adaptive time steps");
42  params.addParam<Real>("scaling_parameter", .8, "scaling parameter for dt selection");
43  return params;
44 }
45 
47  : TimeStepper(parameters),
48  _u1(_fe_problem.getNonlinearSystemBase().addVector("u1", true, GHOSTED)),
49  _aux1(_fe_problem.getAuxiliarySystem().addVector("aux1", true, GHOSTED)),
50  _pred1(_fe_problem.getNonlinearSystemBase().addVector("pred1", true, GHOSTED)),
51  _dt_full(declareRestartableData<Real>("dt_full", 0)),
52  _error(declareRestartableData<Real>("error", 0)),
53  _e_tol(getParam<Real>("e_tol")),
54  _e_max(getParam<Real>("e_max")),
55  _max_increase(getParam<Real>("max_increase")),
56  _steps_between_increase(getParam<int>("steps_between_increase")),
57  _dt_steps_taken(declareRestartableData<int>("dt_steps_taken", 0)),
58  _start_adapting(getParam<int>("start_adapting")),
59  _my_dt_old(declareRestartableData<Real>("my_dt_old", 0)),
60  _infnorm(declareRestartableData<Real>("infnorm", 0)),
61  _scaling_parameter(getParam<Real>("scaling_parameter"))
62 {
63  Real predscale = 1.;
64  InputParameters params = _app.getFactory().getValidParams("AdamsPredictor");
65  params.set<Real>("scale") = predscale;
66  _fe_problem.addPredictor("AdamsPredictor", "adamspredictor", params);
67 }
68 
69 void
71 {
73 }
74 
75 void
77 {
78  // save dt
79  _dt_full = _dt;
80 }
81 
82 void
84 {
87 
90  if (_converged)
91  {
92  _u1 = *nl.currentSolution();
93  _u1.close();
94 
95  _aux1 = *aux.currentSolution();
96  _aux1.close();
97  if (_t_step >= _start_adapting)
98  {
99  // Calculate error if past the first solve
101 
102  _infnorm = _u1.linfty_norm();
103  _e_max = 1.1 * _e_tol * _infnorm;
104  _console << "Time Error Estimate: " << _error << std::endl;
105  }
106  else
107  {
108  // First time step is problematic, sure we converged but what does that mean? We don't know.
109  // Nor can we calculate the error on the first time step.
110  }
111  }
112 }
113 
114 bool
116 {
117  if (!_converged)
118  return false;
119  if (_error < _e_max)
120  return true;
121  else
122  return false;
123 }
124 
125 void
127 {
128  if (!converged())
129  _dt_steps_taken = 0;
130 
131  if (_error >= _e_max)
132  _console << "Marking last solve not converged " << _error << " " << _e_max << std::endl;
133 }
134 
135 Real
137 {
138  if (_t_step <= _start_adapting)
139  return _dt;
140 
141  _my_dt_old = _dt;
142 
143  _dt_steps_taken += 1;
145  {
146 
147  Real new_dt = _dt_full * _scaling_parameter * std::pow(_infnorm * _e_tol / _error, 1.0 / 3.0);
148 
149  if (new_dt / _dt_full > _max_increase)
150  new_dt = _dt_full * _max_increase;
151  _dt_steps_taken = 0;
152  return new_dt;
153  }
154 
155  return _dt;
156 }
157 
158 Real
160 {
161  return getParam<Real>("dt");
162 }
163 
164 Real
165 AB2PredictorCorrector::estimateTimeError(NumericVector<Number> & solution)
166 {
169  auto scheme = Moose::stringToEnum<Moose::TimeIntegratorType>(ti->name());
170  Real dt_old = _my_dt_old;
171  if (dt_old == 0)
172  dt_old = _dt;
173 
174  switch (scheme)
175  {
177  {
178  _pred1 *= -1;
179  _pred1 += solution;
180  Real calc = _dt * _dt * .5;
181  _pred1 *= calc;
182  return _pred1.l2_norm();
183  }
185  {
186  _pred1 -= solution;
187  _pred1 *= (_dt) / (3.0 * (_dt + dt_old));
188  return _pred1.l2_norm();
189  }
190  case Moose::TI_BDF2:
191  {
192  _pred1 *= -1.0;
193  _pred1 += solution;
194  Real topcalc = 2.0 * (_dt + dt_old) * (_dt + dt_old);
195  Real bottomcalc = 6.0 * _dt * _dt + 12.0 * _dt * dt_old + 5.0 * dt_old * dt_old;
196  _pred1 *= topcalc / bottomcalc;
197 
198  return _pred1.l2_norm();
199  }
200  default:
201  break;
202  }
203  return -1;
204 }
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
NonlinearSystemBase & getNonlinearSystemBase()
InputParameters getValidParams(const std::string &name)
Get valid parameters for the object.
Definition: Factory.C:67
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.
TimeIntegrator * getTimeIntegrator()
Definition: SystemBase.h:729
Base class for time stepping.
Definition: TimeStepper.h:26
const NumericVector< Number > *const & currentSolution() const override
The solution vector that is currently being operated on.
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
NumericVector< Number > & _pred1
virtual void preSolve() override
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:286
virtual void step() override
Take a time step.
NonlinearSystemBase * nl
Nonlinear system to be solved.
FEProblemBase & _fe_problem
Definition: TimeStepper.h:119
virtual void postSolve() override
virtual bool converged() override
Real pow(Real x, int e)
Definition: MathUtils.C:211
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()
registerMooseObject("MooseApp", AB2PredictorCorrector)
bool _converged
Whether or not the previous solve converged.
Definition: TimeStepper.h:139
const NumericVector< Number > *const & currentSolution() const override
The solution vector that is currently being operated on.
NumericVector< Number > & _u1
InputParameters validParams< TimeStepper >()
Definition: TimeStepper.C:17
Base class for time integrators.
virtual bool converged() const override
If the time step converged.
virtual NumericVector< Number > & solutionPredictor()
Definition: Predictor.h:44
MooseApp & _app
The MooseApp this object is associated with.
Definition: MooseObject.h:177
virtual void preExecute()
Definition: TimeStepper.C:61
Real & _infnorm
infinity norm of the solution vector
const std::string & name() const
Get the name of the object.
Definition: MooseObject.h:59
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:126
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
A TimeStepper based on the AB2 method.
Real & _dt
Definition: TimeStepper.h:127
A system that holds auxiliary variables.