1#ifndef DUNE_FEM_SOLVER_RUNGEKUTTA_TIMESTEPCONTROL_HH
2#define DUNE_FEM_SOLVER_RUNGEKUTTA_TIMESTEPCONTROL_HH
9#include <dune/common/exceptions.hh>
75 static const std::string verboseTypeTable[] = {
"none",
"noconv",
"cfl",
"full" };
104 virtual bool cflFactor(
const double imOpTimeStepEstimate,
105 const double exOpTimeStepEstimate,
106 const int numberOfLinearIterations,
108 double &factor)
const
110 const int iter = numberOfLinearIterations;
112 bool changed =
false;
117 if( imOpTimeStepEstimate <= exOpTimeStepEstimate )
137 virtual void initTimeStepEstimate (
const double dtEstExpl,
const double dtEstImpl,
double &dtEst,
double &cfl )
const
144 if( (dtEstImpl > 0) && (dtEstExpl > dtEstImpl) )
145 cfl = dtEstExpl / dtEstImpl;
154 const std::string names [] = {
"ImplicitEuler",
"CrankNicolson",
"DIRK23",
"DIRK34",
"SDIRK22" };
196 double estimate = std::numeric_limits< double >::max();
202 template<
class Monitor >
203 void reduceTimeStep (
double helmholtzEstimate,
double sourceTermEstimate,
const Monitor &monitor )
206 DUNE_THROW( Dune::InvalidStateException,
"ImplicitRungeKuttaSolver must be initialized before first solve." );
209 parameters().
cflFactor( helmholtzEstimate, sourceTermEstimate, monitor.linearSolverIterations_,
false, factor );
211 if( !((factor >= std::numeric_limits< double >::min()) && (factor < 1.0)) )
212 DUNE_THROW( Dune::InvalidStateException,
"invalid cfl factor: " << factor );
218 Dune::derr <<
"Implicit Runge-Kutta step failed to converge." << std::endl;
219 Dune::derr <<
"New cfl number: " <<
cfl_
220 <<
", iterations per time step: ILS = " << monitor.linearSolverIterations_
221 <<
", INLS = " << monitor.newtonIterations_
229 template<
class Monitor >
230 void timeStepEstimate (
double helmholtzEstimate,
double sourceTermEstimate,
const Monitor &monitor )
233 DUNE_THROW( Dune::InvalidStateException,
"ImplicitRungeKuttaSolver must be initialized before first solve." );
238 parameters().
cflFactor( helmholtzEstimate, sourceTermEstimate, monitor.linearSolverIterations_,
true, factor );
239 if( !((factor >= std::numeric_limits< double >::min()) && (factor <= std::numeric_limits< double >::max())) )
240 DUNE_THROW( Dune::InvalidStateException,
"invalid cfl factor: " << factor );
242 const double oldCfl =
cfl_;
249 Dune::derr <<
"New cfl number: " <<
cfl_
250 <<
", iterations per time step: ILS = " << monitor.linearSolverIterations_
251 <<
", INLS = " << monitor.newtonIterations_
313 :
BaseType( timeProvider, parameter ),
317 if( parameter.getValue(
"fem.ode.pidcontrol",
bool(
false) ) )
319 tol_ = parameter.getValue(
"fem.ode.pidtolerance",
tol_ );
326 template<
class Monitor >
327 void timeStepEstimate (
double helmholtzEstimate,
double sourceTermEstimate,
const Monitor &monitor )
330 DUNE_THROW( Dune::InvalidStateException,
"ImplicitRungeKuttaSolver must be initialized before first solve." );
344 std::cout <<
"Set dt = " << dtEst << std::endl;
349 Dune::derr <<
"New dt: " << dtEst
350 <<
", iterations per time step: ILS = " << monitor.linearSolverIterations_
351 <<
", INLS = " << monitor.newtonIterations_
361 template <
class Monitor >
365 const double error = monitor.error_;
366 std::cout << error <<
" error " << std::endl;
367 if(
std::abs( error ) < 1e-12 )
return 10. * dt;
370 for(
int i=0; i<2; ++i )
381 const double newDt = dt *
tol_ / error;
384 else if( error > 1e-12 )
387 const double kP = 0.075 ;
388 const double kI = 0.175 ;
389 const double kD = 0.01 ;
393 std::cout <<
"newDt = " << newDt << std::endl;
Dune::Fem::Double abs(const Dune::Fem::Double &a)
Definition: double.hh:942
double min(const Dune::Fem::Double &v, const double p)
Definition: double.hh:953
Dune::Fem::Double pow(const Dune::Fem::Double &a, const Dune::Fem::Double &b)
Definition: double.hh:947
Definition: multistep.hh:17
static ParameterContainer & container()
Definition: io/parameter.hh:193
Definition: io/parameter.hh:551
int getEnum(const std::string &key, const std::string(&values)[n]) const
Definition: reader.hh:225
T getValue(const std::string &key) const
get mandatory parameter
Definition: reader.hh:159
static int rank()
Definition: mpimanager.hh:401
Definition: timestepcontrol.hh:23
virtual double cflMax() const
Definition: timestepcontrol.hh:84
double initialDeltaT(double dt) const
Definition: timestepcontrol.hh:89
virtual int maxLinearIterations() const
Definition: timestepcontrol.hh:149
virtual double tolerance() const
tolerance for the non-linear solver (should be larger than the tolerance for the linear solver
Definition: timestepcontrol.hh:62
ImplicitRungeKuttaSolverParameters(const Dune::Fem::ParameterReader ¶meter=Dune::Fem::Parameter::container())
Definition: timestepcontrol.hh:51
virtual int verbose() const
verbosity level ( none, noconv, cfl, full )
Definition: timestepcontrol.hh:73
virtual ~ImplicitRungeKuttaSolverParameters()
Definition: timestepcontrol.hh:56
const int maxIter_
Definition: timestepcontrol.hh:36
virtual double cflStart() const
Definition: timestepcontrol.hh:79
ImplicitRungeKuttaSolverParameters(const std::string keyPrefix, const Dune::Fem::ParameterReader ¶meter=Dune::Fem::Parameter::container())
Definition: timestepcontrol.hh:43
virtual void initTimeStepEstimate(const double dtEstExpl, const double dtEstImpl, double &dtEst, double &cfl) const
Definition: timestepcontrol.hh:137
const int minIter_
Definition: timestepcontrol.hh:33
@ cflVerbosity
Definition: timestepcontrol.hh:25
@ noConvergenceVerbosity
Definition: timestepcontrol.hh:24
@ fullVerbosity
Definition: timestepcontrol.hh:25
@ noVerbosity
Definition: timestepcontrol.hh:24
const std::string keyPrefix_
Definition: timestepcontrol.hh:29
const double sigma_
Definition: timestepcontrol.hh:38
virtual bool cflFactor(const double imOpTimeStepEstimate, const double exOpTimeStepEstimate, const int numberOfLinearIterations, bool converged, double &factor) const
return multiplication factor for the current cfl number
Definition: timestepcontrol.hh:104
virtual int iterations() const
Definition: timestepcontrol.hh:67
Dune::Fem::ParameterReader parameter_
Definition: timestepcontrol.hh:40
const Dune::Fem::ParameterReader & parameter() const
Definition: timestepcontrol.hh:58
virtual int selectedSolver(const int order) const
return number of selected solver (default = order of solver)
Definition: timestepcontrol.hh:152
Definition: timestepcontrol.hh:166
TimeProviderType & timeProvider_
Definition: timestepcontrol.hh:265
double cfl_
Definition: timestepcontrol.hh:267
ImplicitRungeKuttaSolverParameters ParameterType
Definition: timestepcontrol.hh:171
void initialTimeStepSize(double helmholtzEstimate, double sourceTermEstimate)
Definition: timestepcontrol.hh:194
ImplicitRungeKuttaTimeStepControl(TimeProviderType &timeProvider, const Dune::Fem::ParameterReader ¶meter=Dune::Fem::Parameter::container())
Definition: timestepcontrol.hh:182
void reduceTimeStep(double helmholtzEstimate, double sourceTermEstimate, const Monitor &monitor)
Definition: timestepcontrol.hh:203
Dune::Fem::TimeProviderBase TimeProviderType
Definition: timestepcontrol.hh:170
bool computeError() const
Definition: timestepcontrol.hh:256
const ParameterType & parameters() const
Definition: timestepcontrol.hh:259
void timeStepEstimate(double helmholtzEstimate, double sourceTermEstimate, const Monitor &monitor)
Definition: timestepcontrol.hh:230
std::shared_ptr< const ParameterType > parameters_
Definition: timestepcontrol.hh:266
double timeStepSize() const
Definition: timestepcontrol.hh:192
int verbose_
Definition: timestepcontrol.hh:268
double time() const
Definition: timestepcontrol.hh:191
double cflMax_
Definition: timestepcontrol.hh:267
ImplicitRungeKuttaTimeStepControl(TimeProviderType &timeProvider, const ParameterType ¶meters)
Definition: timestepcontrol.hh:173
bool initialized_
Definition: timestepcontrol.hh:269
PID time step control.
Definition: timestepcontrol.hh:288
double tol_
Definition: timestepcontrol.hh:402
PIDTimeStepControl(TimeProviderType &timeProvider, const Dune::Fem::ParameterReader ¶meter=Dune::Fem::Parameter::container())
Definition: timestepcontrol.hh:312
PIDTimeStepControl(TimeProviderType &timeProvider, const ParameterType ¶meters)
Definition: timestepcontrol.hh:300
BaseType::ParameterType ParameterType
Definition: timestepcontrol.hh:298
Dune::Fem::TimeProviderBase TimeProviderType
Definition: timestepcontrol.hh:297
double pidTimeStepControl(const double dt, const Monitor &monitor)
Definition: timestepcontrol.hh:362
std::vector< double > errors_
Definition: timestepcontrol.hh:401
void timeStepEstimate(double helmholtzEstimate, double sourceTermEstimate, const Monitor &monitor)
Definition: timestepcontrol.hh:327
bool computeError() const
Definition: timestepcontrol.hh:324
general base for time providers
Definition: timeprovider.hh:36
double time() const
obtain the current time
Definition: timeprovider.hh:94
void provideTimeStepEstimate(const double dtEstimate)
set time step estimate to minimum of given value and internal time step estiamte
Definition: timeprovider.hh:142
double deltaT() const
obtain the size of the current time step
Definition: timeprovider.hh:113
void invalidateTimeStep()
count current time step a not valid
Definition: timeprovider.hh:158