1 #ifndef DUNE_FEM_SOLVER_RUNGEKUTTA_TIMESTEPCONTROL_HH 2 #define DUNE_FEM_SOLVER_RUNGEKUTTA_TIMESTEPCONTROL_HH 9 #include <dune/common/exceptions.hh> 44 : keyPrefix_( keyPrefix ),
45 minIter_(
parameter.getValue< int >( keyPrefix_ +
"miniterations", 14 ) ),
46 maxIter_(
parameter.getValue< int >( keyPrefix_ +
"maxiterations" , 16 ) ),
47 sigma_(
parameter.getValue< double >( keyPrefix_ +
"cflincrease" , 1.1 ) ),
75 static const std::string verboseTypeTable[] = {
"none",
"noconv",
"cfl",
"full" };
91 return std::min(
parameter().getValue< double >( keyPrefix_ +
"initialdt", 987654321 ), dt );
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 )
123 else if (iter > maxIter_)
125 factor = (double)maxIter_/(sigma_*(
double)iter);
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" };
156 return parameter().
getEnum( keyPrefix_ +
"solvername", names, order-1 ) + 1;
174 : timeProvider_( timeProvider ),
175 parameters_( parameters.
clone() ),
177 cflMax_( parameters_->
cflMax() ),
178 verbose_( parameters_->
verbose() ),
179 initialized_( false )
183 : timeProvider_( timeProvider ),
184 parameters_(
std::make_shared<ParametersType>(
parameter ) ),
186 cflMax_( parameters_->
cflMax() ),
187 verbose_( parameters_->
verbose() ),
188 initialized_( false )
191 double time ()
const {
return timeProvider_.time(); }
197 parameters().initTimeStepEstimate( sourceTermEstimate, helmholtzEstimate, estimate, cfl_ );
198 timeProvider_.provideTimeStepEstimate( estimate );
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 );
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_
225 timeProvider_.provideTimeStepEstimate( factor * timeStepSize() );
226 timeProvider_.invalidateTimeStep();
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 );
240 DUNE_THROW( Dune::InvalidStateException,
"invalid cfl factor: " << factor );
242 const double oldCfl = cfl_;
243 cfl_ =
std::min( cflMax_, factor * cfl_ );
245 timeProvider_.provideTimeStepEstimate(
std::min( sourceTermEstimate, cfl_ * helmholtzEstimate ) );
249 Dune::derr <<
"New cfl number: " << cfl_
250 <<
", iterations per time step: ILS = " << monitor.linearSolverIterations_
251 <<
", INLS = " << monitor.newtonIterations_
261 assert( parameters_ );
293 using BaseType :: initialized_;
294 using BaseType :: cfl_;
295 using BaseType :: parameters ;
301 : BaseType( timeProvider, parameters ),
308 errors_.resize( 3, tol_ );
320 errors_.resize( 3, 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." );
335 double dtEst = pidTimeStepControl(
std::min( sourceTermEstimate, helmholtzEstimate ), monitor );
344 std::cout <<
"Set dt = " << dtEst << std::endl;
345 timeProvider_.provideTimeStepEstimate( dtEst );
349 Dune::derr <<
"New dt: " << dtEst
350 <<
", iterations per time step: ILS = " << monitor.linearSolverIterations_
351 <<
", INLS = " << monitor.newtonIterations_
357 BaseType::timeStepEstimate( helmholtzEstimate, sourceTermEstimate, monitor );
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 )
372 errors_[ i ] = errors_[i+1];
376 errors_[ 2 ] = error ;
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 ;
390 const double newDt = (dt * std::pow( errors_[ 1 ] / errors_[ 2 ], kP ) *
391 std::pow( tol_ / errors_[ 2 ], kI ) *
392 std::pow( errors_[0]*errors_[0]/errors_[ 1 ]/errors_[ 2 ], kD ));
393 std::cout <<
"newDt = " << newDt << std::endl;
407 #endif // #ifndef DUNE_FEM_SOLVER_RUNGEKUTTA_TIMESTEPCONTROL_HH void initialTimeStepSize(double helmholtzEstimate, double sourceTermEstimate)
Definition: timestepcontrol.hh:194
std::vector< double > errors_
Definition: timestepcontrol.hh:401
ImplicitRungeKuttaSolverParameters ParametersType
Definition: timestepcontrol.hh:171
const int minIter_
Definition: timestepcontrol.hh:33
static int rank()
Definition: mpimanager.hh:116
BaseType::ParametersType ParametersType
Definition: timestepcontrol.hh:298
static constexpr T min(T a)
Definition: utility.hh:81
double pidTimeStepControl(const double dt, const Monitor &monitor)
Definition: timestepcontrol.hh:362
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
bool computeError() const
Definition: timestepcontrol.hh:256
static constexpr T max(T a)
Definition: utility.hh:65
int getEnum(const std::string &key, const std::string(&values)[n]) const
Definition: reader.hh:215
std::shared_ptr< const ParametersType > parameters_
Definition: timestepcontrol.hh:266
Definition: multistep.hh:16
ImplicitRungeKuttaSolverParameters(const std::string keyPrefix, const Dune::Fem::ParameterReader ¶meter=Dune::Fem::Parameter::container())
Definition: timestepcontrol.hh:43
virtual int iterations() const
Definition: timestepcontrol.hh:67
void timeStepEstimate(double helmholtzEstimate, double sourceTermEstimate, const Monitor &monitor)
Definition: timestepcontrol.hh:327
Definition: timestepcontrol.hh:24
double timeStepSize() const
Definition: timestepcontrol.hh:192
ImplicitRungeKuttaSolverParameters(const Dune::Fem::ParameterReader ¶meter=Dune::Fem::Parameter::container())
Definition: timestepcontrol.hh:51
virtual void initTimeStepEstimate(const double dtEstExpl, const double dtEstImpl, double &dtEst, double &cfl) const
Definition: timestepcontrol.hh:137
Dune::Fem::TimeProviderBase TimeProviderType
Definition: timestepcontrol.hh:297
general base for time providers
Definition: timeprovider.hh:35
Double abs(const Double &a)
Definition: double.hh:860
double time() const
Definition: timestepcontrol.hh:191
Dune::Fem::ParameterReader parameter_
Definition: timestepcontrol.hh:40
const double sigma_
Definition: timestepcontrol.hh:38
double tol_
Definition: timestepcontrol.hh:402
const int maxIter_
Definition: timestepcontrol.hh:36
bool initialized_
Definition: timestepcontrol.hh:269
const Dune::Fem::ParameterReader & parameter() const
Definition: timestepcontrol.hh:58
Definition: io/parameter.hh:549
const ParametersType & parameters() const
Definition: timestepcontrol.hh:259
Definition: timestepcontrol.hh:25
Definition: timestepcontrol.hh:165
Definition: timestepcontrol.hh:24
virtual ~ImplicitRungeKuttaSolverParameters()
Definition: timestepcontrol.hh:56
void timeStepEstimate(double helmholtzEstimate, double sourceTermEstimate, const Monitor &monitor)
Definition: timestepcontrol.hh:230
virtual int selectedSolver(const int order) const
return number of selected solver (default = order of solver)
Definition: timestepcontrol.hh:152
static ParameterContainer & container()
Definition: io/parameter.hh:190
PIDTimeStepControl(TimeProviderType &timeProvider, const Dune::Fem::ParameterReader ¶meter=Dune::Fem::Parameter::container())
Definition: timestepcontrol.hh:312
virtual ImplicitRungeKuttaSolverParameters * clone() const
Definition: io/parameter.hh:555
ImplicitRungeKuttaTimeStepControl(TimeProviderType &timeProvider, const Dune::Fem::ParameterReader ¶meter=Dune::Fem::Parameter::container())
Definition: timestepcontrol.hh:182
PID time step control.
Definition: timestepcontrol.hh:287
int verbose_
Definition: timestepcontrol.hh:268
double initialDeltaT(double dt) const
Definition: timestepcontrol.hh:89
void reduceTimeStep(double helmholtzEstimate, double sourceTermEstimate, const Monitor &monitor)
Definition: timestepcontrol.hh:203
ImplicitRungeKuttaTimeStepControl(TimeProviderType &timeProvider, const ParametersType ¶meters)
Definition: timestepcontrol.hh:173
virtual double cflStart() const
Definition: timestepcontrol.hh:79
double cflMax_
Definition: timestepcontrol.hh:267
PIDTimeStepControl(TimeProviderType &timeProvider, const ParametersType ¶meters)
Definition: timestepcontrol.hh:300
const std::string keyPrefix_
Definition: timestepcontrol.hh:29
virtual double tolerance() const
tolerance for the non-linear solver (should be larger than the tolerance for the linear solver ...
Definition: timestepcontrol.hh:62
T getValue(const std::string &key) const
get mandatory parameter
Definition: reader.hh:149
bool computeError() const
Definition: timestepcontrol.hh:324
virtual int maxLinearIterations() const
Definition: timestepcontrol.hh:149
virtual double cflMax() const
Definition: timestepcontrol.hh:84
virtual int verbose() const
verbosity level ( none, noconv, cfl, full )
Definition: timestepcontrol.hh:73
Definition: timestepcontrol.hh:21
TimeProviderType & timeProvider_
Definition: timestepcontrol.hh:265
Dune::Fem::TimeProviderBase TimeProviderType
Definition: timestepcontrol.hh:170
Definition: timestepcontrol.hh:25