https://mooseframework.inl.gov
AB2PredictorCorrector.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://mooseframework.inl.gov
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 
28 using namespace libMesh;
29 
31 
34 {
36  params.addClassDescription(
37  "Implements second order Adams-Bashforth method for timestep calculation.");
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(/*nl_sys=*/0).addVector("u1", true, GHOSTED)),
52  _aux1(_fe_problem.getAuxiliarySystem().addVector("aux1", true, GHOSTED)),
53  _pred1(_fe_problem.getNonlinearSystemBase(/*nl_sys=*/0).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 
91  _fe_problem.solve(/*nl_sys=*/0);
92  _converged = _fe_problem.converged(/*nl_sys=*/0);
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 
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  return false;
122  if (_error < _e_max)
123  return true;
124  else
125  return false;
126 }
127 
128 void
130 {
131  if (!converged())
132  _dt_steps_taken = 0;
133 
134  if (_error >= _e_max)
135  _console << "Marking last solve not converged " << _error << " " << _e_max << std::endl;
136 }
137 
138 Real
140 {
141  if (_t_step <= _start_adapting)
142  return _dt;
143 
144  _my_dt_old = _dt;
145 
146  _dt_steps_taken += 1;
148  {
149 
150  Real new_dt = _dt_full * _scaling_parameter * std::pow(_infnorm * _e_tol / _error, 1.0 / 3.0);
151 
152  if (new_dt / _dt_full > _max_increase)
153  new_dt = _dt_full * _max_increase;
154  _dt_steps_taken = 0;
155  return new_dt;
156  }
157 
158  return _dt;
159 }
160 
161 Real
163 {
164  return getParam<Real>("dt");
165 }
166 
167 Real
169 {
171  const auto & ti =
172  _fe_problem.getNonlinearSystemBase(/*nl_sys=*/0).getTimeIntegrator(/*var_num=*/0);
173 
174  auto scheme = Moose::stringToEnum<Moose::TimeIntegratorType>(ti.type());
175  Real dt_old = _my_dt_old;
176  if (dt_old == 0)
177  dt_old = _dt;
178 
179  switch (scheme)
180  {
182  {
183  _pred1 *= -1;
184  _pred1 += solution;
185  Real calc = _dt * _dt * .5;
186  _pred1 *= calc;
187  return _pred1.l2_norm();
188  }
190  {
191  _pred1 -= solution;
192  _pred1 *= (_dt) / (3.0 * (_dt + dt_old));
193  return _pred1.l2_norm();
194  }
195  case Moose::TI_BDF2:
196  {
197  _pred1 *= -1.0;
198  _pred1 += solution;
199  Real topcalc = 2.0 * (_dt + dt_old) * (_dt + dt_old);
200  Real bottomcalc = 6.0 * _dt * _dt + 12.0 * _dt * dt_old + 5.0 * dt_old * dt_old;
201  _pred1 *= topcalc / bottomcalc;
202 
203  return _pred1.l2_norm();
204  }
205  default:
206  break;
207  }
208  return -1;
209 }
static InputParameters validParams()
Definition: TimeStepper.C:16
AB2PredictorCorrector(const InputParameters &parameters)
Real & _error
global relative time discretization error estimate
virtual Real computeInitialDT() override
Computes time step size for the initial time step.
int _steps_between_increase
steps to take before increasing dt
Real _e_tol
error tolerance
const TimeIntegrator & getTimeIntegrator(const unsigned int var_num) const
Retrieve the time integrator that integrates the given variable&#39;s equation.
Definition: SystemBase.C:1650
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:22
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...
The following methods are specializations for using the libMesh::Parallel::packed_range_* routines fo...
NumericVector< Number > & _pred1
virtual void preSolve() override
Real & _dt_full
dt of the big step
virtual void solve(const unsigned int nl_sys_num)
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 a writable reference to the Factory associated with this App.
Definition: MooseApp.h:434
virtual void step() override
Take a time step.
Nonlinear system to be solved.
FEProblemBase & _fe_problem
Definition: TimeStepper.h:120
virtual void addPredictor(const std::string &type, const std::string &name, InputParameters &parameters)
virtual void postSolve() override
virtual Real l2_norm() const=0
virtual const NumericVector< Number > *const & currentSolution() const override final
The solution vector that is currently being operated on.
Definition: SolverSystem.h:117
virtual Real estimateTimeError(NumericVector< Number > &sol)
virtual bool converged(const unsigned int sys_num)
Eventually we want to convert this virtual over to taking a solver system number argument.
Definition: SubProblem.h:113
Real _e_max
maximal error
NonlinearSystemBase & getNonlinearSystemBase(const unsigned int sys_num)
Real _max_increase
maximum increase ratio
MooseApp & _app
The MOOSE application this is associated with.
Definition: MooseBase.h:84
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:140
virtual void close()=0
NumericVector< Number > & _u1
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual bool converged() const override
If the time step converged.
virtual NumericVector< Number > & solutionPredictor()
Definition: Predictor.h:41
virtual void preExecute()
Definition: TimeStepper.C:70
Real & _infnorm
infinity norm of the solution vector
void addClassDescription(const std::string &doc_string)
This method adds a description of the class that will be displayed in the input file syntax dump...
void addParam(const std::string &name, const S &value, const std::string &doc_string)
These methods add an optional parameter and a documentation string to the InputParameters object...
int & _t_step
Definition: TimeStepper.h:127
virtual Real computeDT() override
Computes time step size after the initial time step.
const ConsoleStream _console
An instance of helper class to write streams to the Console objects.
virtual void preExecute() override
int & _dt_steps_taken
steps taken at current dt
A TimeStepper based on the AB2 method.
Real & _dt
Definition: TimeStepper.h:128
MooseUnits pow(const MooseUnits &, int)
Definition: Units.C:537
A system that holds auxiliary variables.
virtual Real linfty_norm() const=0
void ErrorVector unsigned int
static InputParameters validParams()