3#ifndef OPM_ADAPTIVE_TIME_STEPPING_HPP
4#define OPM_ADAPTIVE_TIME_STEPPING_HPP
6#include <dune/common/version.hh>
7#if DUNE_VERSION_GTE(DUNE_ISTL, 2, 8)
8#include <dune/istl/istlexception.hh>
10#include <dune/istl/ilu.hh>
13#include <opm/common/Exceptions.hpp>
14#include <opm/common/ErrorMacros.hpp>
15#include <opm/common/OpmLog/OpmLog.hpp>
17#include <opm/grid/utility/StopWatch.hpp>
19#include <opm/input/eclipse/Units/Units.hpp>
20#include <opm/input/eclipse/Units/UnitSystem.hpp>
22#include <opm/input/eclipse/Schedule/Tuning.hpp>
28#include <opm/simulators/timestepping/AdaptiveSimulatorTimer.hpp>
29#include <opm/simulators/timestepping/EclTimeSteppingParams.hpp>
30#include <opm/simulators/timestepping/SimulatorReport.hpp>
31#include <opm/simulators/timestepping/SimulatorTimer.hpp>
32#include <opm/simulators/timestepping/TimeStepControl.hpp>
33#include <opm/simulators/timestepping/TimeStepControlInterface.hpp>
35#include <opm/simulators/utils/phaseUsageFromDeck.hpp>
37#include <fmt/format.h>
50namespace Opm::Parameters {
80std::set<std::string> consistentlyFailingWells(
const std::vector<StepReport>&
sr);
82void registerAdaptiveParameters();
88 template<
class TypeTag>
92 template <
class Solver>
95 const Solver& solver_;
97 SolutionTimeErrorSolverWrapper(
const Solver&
solver)
102 double relativeChange()
const
103 {
return solver_.model().relativeChange(); }
111 message =
"Caught Exception: ";
113 OpmLog::debug(message);
125 ,
restartFactor_(Parameters::Get<Parameters::SolverRestartFactor<Scalar>>())
126 ,
growthFactor_(Parameters::Get<Parameters::SolverGrowthFactor<Scalar>>())
127 ,
maxGrowth_(Parameters::Get<Parameters::SolverMaxGrowth<Scalar>>())
128 ,
maxTimeStep_(Parameters::Get<Parameters::SolverMaxTimeStepInDays<Scalar>>() * 24 * 60 * 60)
136 ,
timestepAfterEvent_(Parameters::Get<Parameters::TimeStepAfterEventInDays<Scalar>>() * 24 * 60 * 60)
138 , minTimeStepBeforeShuttingProblematicWells_(Parameters::Get<Parameters::MinTimeStepBeforeShuttingProblematicWellsInDays>() *
unit::
day)
167 , minTimeStepBeforeShuttingProblematicWells_(Parameters::Get<Parameters::MinTimeStepBeforeShuttingProblematicWellsInDays>() *
unit::
day)
172 static void registerParameters()
175 detail::registerAdaptiveParameters();
183 template <
class Solver>
187 const std::vector<int>*
fipnum =
nullptr,
188 const std::function<
bool()>
tuningUpdater = [](){
return false;})
209 auto& simulator =
solver.model().simulator();
210 auto& problem = simulator.problem();
240 OpmLog::debug(
"Overall linear iterations used: " + std::to_string(
substepReport.total_linear_iterations));
245 causeOfFailure =
"Solver convergence failure - Iteration limit reached";
250 catch (
const ConvergenceMonitorFailure&
e) {
263 causeOfFailure =
"Solver convergence failure - Numerical problem encountered";
268 catch (
const std::runtime_error&
e) {
274 catch (
const Dune::ISTLError&
e) {
280 catch (
const Dune::MatrixBlockError&
e) {
297 const auto msg = fmt::format(
298 "Solver failed to converge but timestep "
299 "{} is smaller or equal to {}\n"
300 "which is the minimum threshold given "
301 "by option --solver-min-time-step\n",
304 OpmLog::problem(
msg);
334 if (
solver.model().wellModel().hasTHPConstraints()) {
340 std::ostringstream
ss;
342 OpmLog::info(
ss.str());
356 problem.writeOutput();
358 report.success.output_write_time +=
perfTimer.secsSinceStart();
374 const auto msg = fmt::format(
375 "Solver failed to converge after cutting timestep {} times.",
392 const auto msg = fmt::format(
393 "Solver failed to converge after cutting timestep to {}\n"
394 "which is the minimum threshold given by option --solver-min-time-step\n",
408 const auto msg = fmt::format(
409 "{}\nTimestep chopped to {} days\n",
411 std::to_string(unit::convert::to(
substepTimer.currentStepLength(), unit::day))
413 OpmLog::problem(
msg);
424 std::set<std::string>
failing_wells = detail::consistentlyFailingWells(
solver.model().stepReports());
445 msg =
"\nProblematic well(s) were shut: ";
450 msg +=
"(retrying timestep)\n";
451 OpmLog::problem(
msg);
457 problem.setNextTimeStepSize(
substepTimer.currentStepLength());
463 std::ostringstream
ss;
465 ss <<
"Suggested next step size = " << unit::convert::to(
suggestedNextTimestep_, unit::day) <<
" (days)" << std::endl;
466 OpmLog::debug(
ss.str());
481 void setSuggestedNextStep(
const double x)
502 template<
class Serializer>
507 case TimeStepControlType::HardCodedTimeStep:
510 case TimeStepControlType::PIDAndIterationCount:
513 case TimeStepControlType::SimpleIterationCount:
516 case TimeStepControlType::PID:
533 serializer(minTimeStepBeforeShuttingProblematicWells_);
566 case TimeStepControlType::HardCodedTimeStep:
569 case TimeStepControlType::PIDAndIterationCount:
572 case TimeStepControlType::SimpleIterationCount:
575 case TimeStepControlType::PID:
592 this->minTimeStepBeforeShuttingProblematicWells_ ==
593 rhs.minTimeStepBeforeShuttingProblematicWells_;
597 template<
class Controller>
602 result.restartFactor_ = 1.0;
603 result.growthFactor_ = 2.0;
605 result.maxTimeStep_ = 4.0;
606 result.minTimeStep_ = 5.0;
607 result.ignoreConvergenceFailure_ =
true;
608 result.solverRestartMax_ = 6;
609 result.solverVerbose_ =
true;
610 result.timestepVerbose_ =
true;
611 result.suggestedNextTimestep_ = 7.0;
612 result.fullTimestepInitially_ =
true;
613 result.useNewtonIteration_ =
true;
614 result.minTimeStepBeforeShuttingProblematicWells_ = 9.0;
615 result.timeStepControlType_ = Controller::Type;
616 result.timeStepControl_ = std::make_unique<Controller>(Controller::serializationTestObject());
620 template<
class T,
class Serializer>
621 void allocAndSerialize(Serializer&
serializer)
633 const T* rhs =
static_cast<const T*
>(
Rhs.timeStepControl_.get());
641 std::string
control = Parameters::Get<Parameters::TimeStepControl>();
643 const double tol = Parameters::Get<Parameters::TimeStepControlTolerance>();
648 else if (
control ==
"pid+iteration") {
649 const int iterations = Parameters::Get<Parameters::TimeStepControlTargetIterations>();
650 const double decayDampingFactor = Parameters::Get<Parameters::TimeStepControlDecayDampingFactor>();
651 const double growthDampingFactor = Parameters::Get<Parameters::TimeStepControlGrowthDampingFactor>();
655 else if (
control ==
"pid+newtoniteration") {
656 const int iterations = Parameters::Get<Parameters::TimeStepControlTargetNewtonIterations>();
657 const double decayDampingFactor = Parameters::Get<Parameters::TimeStepControlDecayDampingFactor>();
658 const double growthDampingFactor = Parameters::Get<Parameters::TimeStepControlGrowthDampingFactor>();
667 else if (
control ==
"iterationcount") {
668 const int iterations = Parameters::Get<Parameters::TimeStepControlTargetIterations>();
669 const double decayrate = Parameters::Get<Parameters::TimeStepControlDecayRate>();
670 const double growthrate = Parameters::Get<Parameters::TimeStepControlGrowthRate>();
674 else if (
control ==
"newtoniterationcount") {
675 const int iterations = Parameters::Get<Parameters::TimeStepControlTargetNewtonIterations>();
676 const double decayrate = Parameters::Get<Parameters::TimeStepControlDecayRate>();
677 const double growthrate = Parameters::Get<Parameters::TimeStepControlGrowthRate>();
682 else if (
control ==
"hardcoded") {
683 const std::string
filename = Parameters::Get<Parameters::TimeStepControlFileName>();
689 "Unsupported time step control selected " +
control);
695 using TimeStepController = std::unique_ptr<TimeStepControlInterface>;
712 double minTimeStepBeforeShuttingProblematicWells_;
Defines a type tags and some fundamental properties all models.
Simulation timer for adaptive time stepping.
Definition AdaptiveSimulatorTimer.hpp:41
Definition AdaptiveTimeStepping.hpp:90
double minTimeStep_
minimal allowed time step size before throwing
Definition AdaptiveTimeStepping.hpp:703
double suggestedNextTimestep_
suggested size of next timestep
Definition AdaptiveTimeStepping.hpp:708
bool solverVerbose_
solver verbosity
Definition AdaptiveTimeStepping.hpp:706
double suggestedNextStep() const
Returns the simulator report for the failed substeps of the last report step.
Definition AdaptiveTimeStepping.hpp:478
double growthFactor_
factor to multiply time step when solver recovered from failed convergence
Definition AdaptiveTimeStepping.hpp:700
AdaptiveTimeStepping(const UnitSystem &unitSystem, const double max_next_tstep=-1.0, const bool terminalOutput=true)
contructor taking parameter object
Definition AdaptiveTimeStepping.hpp:121
TimeStepController timeStepControl_
time step control object
Definition AdaptiveTimeStepping.hpp:698
bool fullTimestepInitially_
beginning with the size of the time step from data file
Definition AdaptiveTimeStepping.hpp:709
AdaptiveTimeStepping(double max_next_tstep, const Tuning &tuning, const UnitSystem &unitSystem, const bool terminalOutput=true)
contructor taking parameter object
Definition AdaptiveTimeStepping.hpp:149
double timestepAfterEvent_
suggested size of timestep after an event
Definition AdaptiveTimeStepping.hpp:710
bool useNewtonIteration_
use newton iteration count for adaptive time step control
Definition AdaptiveTimeStepping.hpp:711
bool timestepVerbose_
timestep verbosity
Definition AdaptiveTimeStepping.hpp:707
int solverRestartMax_
how many restart of solver are allowed
Definition AdaptiveTimeStepping.hpp:705
TimeStepControlType timeStepControlType_
type of time step control object
Definition AdaptiveTimeStepping.hpp:697
SimulatorReport step(const SimulatorTimer &simulatorTimer, Solver &solver, const bool isEvent, const std::vector< int > *fipnum=nullptr, const std::function< bool()> tuningUpdater=[](){return false;})
step method that acts like the solver::step method in a sub cycle of time steps
Definition AdaptiveTimeStepping.hpp:184
bool ignoreConvergenceFailure_
continue instead of stop when minimum time step is reached
Definition AdaptiveTimeStepping.hpp:704
double maxTimeStep_
maximal allowed time step size in days
Definition AdaptiveTimeStepping.hpp:702
double restartFactor_
factor to multiply time step with when solver fails to converge
Definition AdaptiveTimeStepping.hpp:699
double maxGrowth_
factor that limits the maximum growth of a time step
Definition AdaptiveTimeStepping.hpp:701
RelativeChangeInterface.
Definition TimeStepControlInterface.hpp:32
Definition SimulatorTimer.hpp:39
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition propertysystem.hh:235
This file provides the infrastructure to retrieve run-time parameters.
The Opm property system, traits with inheritance.
Definition AdaptiveTimeStepping.hpp:57
Definition AdaptiveTimeStepping.hpp:56
Definition AdaptiveTimeStepping.hpp:68
Definition AdaptiveTimeStepping.hpp:67
Definition AdaptiveTimeStepping.hpp:52
Definition AdaptiveTimeStepping.hpp:53
Definition AdaptiveTimeStepping.hpp:54
Definition AdaptiveTimeStepping.hpp:64
Definition AdaptiveTimeStepping.hpp:62
Definition AdaptiveTimeStepping.hpp:66
Definition AdaptiveTimeStepping.hpp:65
Definition AdaptiveTimeStepping.hpp:63
Definition AdaptiveTimeStepping.hpp:60
Definition AdaptiveTimeStepping.hpp:61
Definition AdaptiveTimeStepping.hpp:59
Definition AdaptiveTimeStepping.hpp:58
Definition AdaptiveTimeStepping.hpp:55
Definition SimulatorReport.hpp:100
Definition ConvergenceReport.hpp:447