20#ifndef OPM_ADAPTIVE_TIME_STEPPING_IMPL_HPP
21#define OPM_ADAPTIVE_TIME_STEPPING_IMPL_HPP
24#ifndef OPM_ADAPTIVE_TIME_STEPPING_HPP
30#include <dune/common/timer.hh>
31#include <dune/istl/istlexception.hh>
33#include <opm/common/Exceptions.hpp>
34#include <opm/common/ErrorMacros.hpp>
35#include <opm/common/OpmLog/OpmLog.hpp>
36#include <opm/common/TimingMacros.hpp>
38#include <opm/grid/utility/StopWatch.hpp>
40#include <opm/input/eclipse/Schedule/Tuning.hpp>
42#include <opm/input/eclipse/Units/Units.hpp>
43#include <opm/input/eclipse/Units/UnitSystem.hpp>
55#include <boost/date_time/posix_time/posix_time.hpp>
56#include <fmt/format.h>
57#include <fmt/ranges.h>
65template<
class TypeTag>
69 const double max_next_tstep,
70 const bool terminal_output
72 : time_step_control_{}
73 , restart_factor_{Parameters::
Get<Parameters::SolverRestartFactor<Scalar>>()}
74 , growth_factor_{Parameters::
Get<Parameters::SolverGrowthFactor<Scalar>>()}
75 , max_growth_{Parameters::
Get<Parameters::SolverMaxGrowth<Scalar>>()}
77 Parameters::
Get<Parameters::SolverMaxTimeStepInDays<Scalar>>() * 24 * 60 * 60}
79 unit_system.to_si(UnitSystem::measure::time,
80 Parameters::
Get<Parameters::SolverMinTimeStep<Scalar>>())}
81 , ignore_convergence_failure_{
82 Parameters::
Get<Parameters::SolverContinueOnConvergenceFailure>()}
83 , solver_restart_max_{Parameters::
Get<Parameters::SolverMaxRestarts>()}
84 , solver_verbose_{Parameters::
Get<Parameters::SolverVerbosity>() > 0 && terminal_output}
85 , timestep_verbose_{Parameters::
Get<Parameters::TimeStepVerbosity>() > 0 && terminal_output}
86 , suggested_next_timestep_{
87 (max_next_tstep <= 0 ? Parameters::
Get<Parameters::InitialTimeStepInDays>()
88 : max_next_tstep) * 24 * 60 * 60}
89 , full_timestep_initially_{Parameters::
Get<Parameters::FullTimeStepInitially>()}
90 , timestep_after_event_{
91 Parameters::
Get<Parameters::TimeStepAfterEventInDays<Scalar>>() * 24 * 60 * 60}
92 , use_newton_iteration_{false}
93 , min_time_step_before_shutting_problematic_wells_{
94 Parameters::
Get<Parameters::MinTimeStepBeforeShuttingProblematicWellsInDays>() * unit::day}
106template<
class TypeTag>
109 const Tuning& tuning,
110 const UnitSystem& unit_system,
112 const bool terminal_output
114 : time_step_control_{}
115 , restart_factor_{tuning.TSFCNV}
116 , growth_factor_{tuning.TFDIFF}
117 , max_growth_{tuning.TSFMAX}
118 , max_time_step_{tuning.TSMAXZ}
119 , min_time_step_{tuning.TSMINZ}
120 , ignore_convergence_failure_{true}
121 , solver_restart_max_{Parameters::
Get<Parameters::SolverMaxRestarts>()}
122 , solver_verbose_{Parameters::
Get<Parameters::SolverVerbosity>() > 0 && terminal_output}
123 , timestep_verbose_{Parameters::
Get<Parameters::TimeStepVerbosity>() > 0 && terminal_output}
124 , suggested_next_timestep_{
125 max_next_tstep <= 0 ? Parameters::
Get<Parameters::InitialTimeStepInDays>() * 24 * 60 * 60
127 , full_timestep_initially_{Parameters::
Get<Parameters::FullTimeStepInitially>()}
128 , timestep_after_event_{tuning.TMAXWC}
129 , use_newton_iteration_{false}
130 , min_time_step_before_shutting_problematic_wells_{
131 Parameters::
Get<Parameters::MinTimeStepBeforeShuttingProblematicWellsInDays>() * unit::day}
137template<
class TypeTag>
149 switch (this->time_step_control_type_) {
151 result = castAndComp<HardcodedTimeStepControl>(rhs);
154 result = castAndComp<PIDAndIterationCountTimeStepControl>(rhs);
157 result = castAndComp<SimpleIterationCountTimeStepControl>(rhs);
160 result = castAndComp<PIDTimeStepControl>(rhs);
163 result = castAndComp<General3rdOrderController>(rhs);
179 this->min_time_step_before_shutting_problematic_wells_ ==
183template<
class TypeTag>
188 registerEclTimeSteppingParameters<Scalar>();
193template<
class TypeTag>
194template <
class Solver>
202 SubStepper<Solver> sub_stepper{
203 *
this, simulator_timer, solver, is_event, tuning_updater,
205 return sub_stepper.run();
208template<
class TypeTag>
209template<
class Serializer>
214 serializer(this->time_step_control_type_);
215 switch (this->time_step_control_type_) {
217 allocAndSerialize<HardcodedTimeStepControl>(serializer);
220 allocAndSerialize<PIDAndIterationCountTimeStepControl>(serializer);
223 allocAndSerialize<SimpleIterationCountTimeStepControl>(serializer);
226 allocAndSerialize<PIDTimeStepControl>(serializer);
229 allocAndSerialize<General3rdOrderController>(serializer);
232 serializer(this->restart_factor_);
233 serializer(this->growth_factor_);
234 serializer(this->max_growth_);
235 serializer(this->max_time_step_);
236 serializer(this->min_time_step_);
237 serializer(this->ignore_convergence_failure_);
238 serializer(this->solver_restart_max_);
239 serializer(this->solver_verbose_);
240 serializer(this->timestep_verbose_);
241 serializer(this->suggested_next_timestep_);
242 serializer(this->full_timestep_initially_);
243 serializer(this->timestep_after_event_);
244 serializer(this->use_newton_iteration_);
245 serializer(this->min_time_step_before_shutting_problematic_wells_);
248template<
class TypeTag>
256template<
class TypeTag>
261 return serializationTestObject_<HardcodedTimeStepControl>();
264template<
class TypeTag>
269 return serializationTestObject_<PIDTimeStepControl>();
272template<
class TypeTag>
277 return serializationTestObject_<PIDAndIterationCountTimeStepControl>();
280template<
class TypeTag>
285 return serializationTestObject_<SimpleIterationCountTimeStepControl>();
288template<
class TypeTag>
293 return serializationTestObject_<General3rdOrderController>();
297template<
class TypeTag>
302 this->suggested_next_timestep_ = x;
305template<
class TypeTag>
310 return this->suggested_next_timestep_;
313template<
class TypeTag>
318 return *this->time_step_control_;
322template<
class TypeTag>
329 if (max_next_tstep > 0) {
330 this->suggested_next_timestep_ = max_next_tstep;
334template<
class TypeTag>
337updateTUNING(
double max_next_tstep,
const Tuning& tuning)
339 this->restart_factor_ = tuning.TSFCNV;
340 this->growth_factor_ = tuning.TFDIFF;
341 this->max_growth_ = tuning.TSFMAX;
342 this->max_time_step_ = tuning.TSMAXZ;
343 updateNEXTSTEP(max_next_tstep);
344 this->timestep_after_event_ = tuning.TMAXWC;
351template<
class TypeTag>
352template<
class T,
class Serializer>
357 if (!serializer.isSerializing()) {
358 this->time_step_control_ = std::make_unique<T>();
360 serializer(*
static_cast<T*
>(this->time_step_control_.get()));
363template<
class TypeTag>
366AdaptiveTimeStepping<TypeTag>::
367castAndComp(
const AdaptiveTimeStepping<TypeTag>& Rhs)
const
369 const T* lhs =
static_cast<const T*
>(this->time_step_control_.get());
370 const T* rhs =
static_cast<const T*
>(Rhs.time_step_control_.get());
374template<
class TypeTag>
376AdaptiveTimeStepping<TypeTag>::
377maybeModifySuggestedTimeStepAtBeginningOfReportStep_(
const double original_time_step,
381 if (this->suggested_next_timestep_ < 0) {
382 this->suggested_next_timestep_ = this->restart_factor_ * original_time_step;
385 if (this->full_timestep_initially_) {
386 this->suggested_next_timestep_ = original_time_step;
390 if (is_event && this->timestep_after_event_ > 0) {
391 this->suggested_next_timestep_ = this->timestep_after_event_;
395template<
class TypeTag>
396template<
class Controller>
397AdaptiveTimeStepping<TypeTag>
398AdaptiveTimeStepping<TypeTag>::
399serializationTestObject_()
401 AdaptiveTimeStepping<TypeTag> result;
403 result.restart_factor_ = 1.0;
404 result.growth_factor_ = 2.0;
405 result.max_growth_ = 3.0;
406 result.max_time_step_ = 4.0;
407 result.min_time_step_ = 5.0;
408 result.ignore_convergence_failure_ =
true;
409 result.solver_restart_max_ = 6;
410 result.solver_verbose_ =
true;
411 result.timestep_verbose_ =
true;
412 result.suggested_next_timestep_ = 7.0;
413 result.full_timestep_initially_ =
true;
414 result.use_newton_iteration_ =
true;
415 result.min_time_step_before_shutting_problematic_wells_ = 9.0;
416 result.time_step_control_type_ = Controller::Type;
417 result.time_step_control_ =
418 std::make_unique<Controller>(Controller::serializationTestObject());
427template<
class TypeTag>
429init_(
const UnitSystem& unitSystem)
431 std::tie(time_step_control_type_,
435 if (this->growth_factor_ < 1.0) {
436 OPM_THROW(std::runtime_error,
437 "Growth factor cannot be less than 1.");
447template<
class TypeTag>
448template<
class Solver>
454 const TuningUpdateCallback& tuning_updater)
455 : adaptive_time_stepping_{adaptive_time_stepping}
456 , simulator_timer_{simulator_timer}
458 , is_event_{is_event}
459 , tuning_updater_{tuning_updater}
463template<
class TypeTag>
464template<
class Solver>
465AdaptiveTimeStepping<TypeTag>&
466AdaptiveTimeStepping<TypeTag>::SubStepper<Solver>::
467getAdaptiveTimerStepper()
469 return adaptive_time_stepping_;
472template<
class TypeTag>
473template<
class Solver>
475AdaptiveTimeStepping<TypeTag>::SubStepper<Solver>::
478#ifdef RESERVOIR_COUPLING_ENABLED
479 if (isReservoirCouplingSlave_() && reservoirCouplingSlave_().activated()) {
480 return runStepReservoirCouplingSlave_();
482 else if (isReservoirCouplingMaster_() && reservoirCouplingMaster_().activated()) {
483 return runStepReservoirCouplingMaster_();
486 return runStepOriginal_();
489 return runStepOriginal_();
497#ifdef RESERVOIR_COUPLING_ENABLED
498template<
class TypeTag>
499template<
class Solver>
501AdaptiveTimeStepping<TypeTag>::SubStepper<Solver>::
502isReservoirCouplingMaster_()
const
504 return this->solver_.model().simulator().reservoirCouplingMaster() !=
nullptr;
507template<
class TypeTag>
508template<
class Solver>
510AdaptiveTimeStepping<TypeTag>::SubStepper<Solver>::
511isReservoirCouplingSlave_()
const
513 return this->solver_.model().simulator().reservoirCouplingSlave() !=
nullptr;
517template<
class TypeTag>
518template<
class Solver>
520AdaptiveTimeStepping<TypeTag>::SubStepper<Solver>::
521maybeModifySuggestedTimeStepAtBeginningOfReportStep_(
const double original_time_step)
523 this->adaptive_time_stepping_.maybeModifySuggestedTimeStepAtBeginningOfReportStep_(
524 original_time_step, this->is_event_
531 struct ZeroRelativeChange final : RelativeChangeInterface {
532 double relativeChange()
const override {
return 0.0; }
533 } zero_relative_change;
535 const auto* hardcoded_control =
static_cast<const HardcodedTimeStepControl*
>(
536 this->adaptive_time_stepping_.time_step_control_.get());
537 AdaptiveSimulatorTimer report_step_timer{
538 this->simulator_timer_.startDateTime(),
540 this->simulator_timer_.simulationTimeElapsed(),
542 this->simulator_timer_.reportStepNum(),
546 const double hardcoded_initial_step = hardcoded_control->computeTimeStepSize(
549 zero_relative_change,
551 if (std::isfinite(hardcoded_initial_step) && hardcoded_initial_step > 0.0) {
552 this->adaptive_time_stepping_.setSuggestedNextStep(hardcoded_initial_step);
559template<
class TypeTag>
560template<
class Solver>
562AdaptiveTimeStepping<TypeTag>::SubStepper<Solver>::
563maybeUpdateTuning_(
double elapsed,
double substep_length,
int sub_step_number)
const
565 return this->tuning_updater_(elapsed, substep_length, sub_step_number);
568template<
class TypeTag>
569template<
class Solver>
571AdaptiveTimeStepping<TypeTag>::SubStepper<Solver>::
574 return this->adaptive_time_stepping_.max_time_step_;
577template <
class TypeTag>
578template <
class Solver>
580AdaptiveTimeStepping<TypeTag>::SubStepper<Solver>::
583 const auto elapsed = this->simulator_timer_.simulationTimeElapsed();
584 const auto original_time_step = this->simulator_timer_.currentStepLength();
585 const auto report_step = this->simulator_timer_.reportStepNum();
586 maybeUpdateTuning_(elapsed, suggestedNextTimestep_(), 0);
587 maybeModifySuggestedTimeStepAtBeginningOfReportStep_(original_time_step);
589 AdaptiveSimulatorTimer substep_timer{
590 this->simulator_timer_.startDateTime(),
593 suggestedNextTimestep_(),
597 SubStepIteration<Solver> substepIteration{*
this, substep_timer, original_time_step,
true};
598 return substepIteration.run();
601template <
class TypeTag>
602template <
class Solver>
604AdaptiveTimeStepping<TypeTag>::SubStepper<Solver>::
605suggestedNextTimestep_()
const
607 return this->adaptive_time_stepping_.suggestedNextStep();
611#ifdef RESERVOIR_COUPLING_ENABLED
617template <
class TypeTag>
618template <
class Solver>
620AdaptiveTimeStepping<TypeTag>::SubStepper<Solver>::
621checkIfSlaveIsTerminated_()
623 if (reservoirCouplingSlave_().terminated()) {
624 OPM_THROW(ReservoirCouplingError,
625 "Internal error: attempt to run a coupled substep loop after the slave "
626 "has been terminated by the master process");
645template <
class TypeTag>
646template <
class Solver>
648AdaptiveTimeStepping<TypeTag>::SubStepper<Solver>::
649getRcMasterSyncStepLength_(
double prev_step,
651 double step_end_time)
653 const bool sync_at_report_steps = reservoirCouplingMaster_().syncAtReportSteps();
654 double current_step_length;
655 if (sync_at_report_steps) {
656 current_step_length = prev_step;
658 const double remaining = step_end_time - current_time;
659 current_step_length = std::min({suggestedNextTimestep_(), maxTimeStep_(), remaining});
672 current_step_length = reservoirCouplingMaster_().maybeChopSubStep(current_step_length, current_time);
673 auto num_active = reservoirCouplingMaster_().numActivatedSlaves();
674 OpmLog::info(fmt::format(
675 "\nChoosing next sync time{} between master and {} active slave {}: {:.2f} days",
676 sync_at_report_steps ?
" (RSYNC)" :
"",
677 num_active, (num_active == 1 ?
"process" :
"processes"),
678 current_step_length / unit::day
680 return current_step_length;
683template <
class TypeTag>
684template <
class Solver>
685ReservoirCouplingMaster<typename AdaptiveTimeStepping<TypeTag>::Scalar>&
686AdaptiveTimeStepping<TypeTag>::SubStepper<Solver>::
687reservoirCouplingMaster_()
689 return *(this->solver_.model().simulator().reservoirCouplingMaster());
692template <
class TypeTag>
693template <
class Solver>
694ReservoirCouplingSlave<typename AdaptiveTimeStepping<TypeTag>::Scalar>&
695AdaptiveTimeStepping<TypeTag>::SubStepper<Solver>::
696reservoirCouplingSlave_()
698 return *(this->solver_.model().simulator().reservoirCouplingSlave());
755template <
class TypeTag>
756template <
class Solver>
758AdaptiveTimeStepping<TypeTag>::SubStepper<Solver>::
759runStepReservoirCouplingMaster_()
762 const double original_time_step = this->simulator_timer_.currentStepLength();
763 double current_time{this->simulator_timer_.simulationTimeElapsed()};
764 double step_end_time = current_time + original_time_step;
765 const double report_step_start_time = current_time;
766 int report_step_substep_offset = 0;
773 auto current_step_length = original_time_step;
774 auto report_step_idx = this->simulator_timer_.currentStepNum();
775 if (report_step_idx == 0 && iteration == 0) {
776 reservoirCouplingMaster_().initTimeStepping();
778 SimulatorReport report;
780 reservoirCouplingMaster_().maybeReceiveActivationHandshakeFromSlaves(current_time);
782 reservoirCouplingMaster_().sendDontTerminateSignalToSlaves();
783 reservoirCouplingMaster_().receiveNextReportDateFromSlaves();
784 const bool start_of_report_step = (iteration == 0);
785 if (start_of_report_step) {
786 reservoirCouplingMaster_().initStartOfReportStep(report_step_idx);
787 maybeUpdateTuning_(current_time, suggestedNextTimestep_(), 0);
788 maybeModifySuggestedTimeStepAtBeginningOfReportStep_(original_time_step);
790 current_step_length = getRcMasterSyncStepLength_(
791 current_step_length, current_time, step_end_time);
792 reservoirCouplingMaster_().sendNextTimeStepToSlaves(current_step_length);
793 AdaptiveSimulatorTimer substep_timer{
794 this->simulator_timer_.startDateTime(),
797 suggestedNextTimestep_(),
798 this->simulator_timer_.reportStepNum(),
804 substep_timer.setReportStepStartTime(report_step_start_time);
805 substep_timer.setReportStepTotalTime(step_end_time);
806 substep_timer.setReportStepSubstepOffset(report_step_substep_offset);
808 current_time + current_step_length, step_end_time
813 reservoirCouplingMaster_().setFirstSubstepOfSyncTimestep(
true);
817 reservoirCouplingMaster_().setNeedsSlaveDataReceive(
true);
818 SubStepIteration<Solver> substepIteration{*
this, substep_timer, current_step_length, final_step};
819 const auto sub_steps_report = substepIteration.run();
820 report += sub_steps_report;
821 report_step_substep_offset += substep_timer.currentStepNum();
822 current_time += current_step_length;
831template <
class TypeTag>
832template <
class Solver>
834AdaptiveTimeStepping<TypeTag>::SubStepper<Solver>::
835runStepReservoirCouplingSlave_()
837 checkIfSlaveIsTerminated_();
839 const double original_time_step = this->simulator_timer_.currentStepLength();
840 double current_time{this->simulator_timer_.simulationTimeElapsed()};
841 double step_end_time = current_time + original_time_step;
842 const double report_step_start_time = current_time;
843 int report_step_substep_offset = 0;
844 SimulatorReport report;
845 auto report_step_idx = this->simulator_timer_.currentStepNum();
846 if (report_step_idx == 0 && iteration == 0) {
847 reservoirCouplingSlave_().initTimeStepping();
850 bool start_of_report_step = (iteration == 0);
851 if (reservoirCouplingSlave_().maybeReceiveTerminateSignalFromMaster()) {
855 reservoirCouplingSlave_().sendNextReportDateToMasterProcess();
856 const auto timestep = reservoirCouplingSlave_().receiveNextTimeStepFromMaster();
857 if (start_of_report_step) {
858 maybeUpdateTuning_(current_time, suggestedNextTimestep_(), 0);
859 maybeModifySuggestedTimeStepAtBeginningOfReportStep_(timestep);
861 AdaptiveSimulatorTimer substep_timer{
862 this->simulator_timer_.startDateTime(),
865 suggestedNextTimestep_(),
866 this->simulator_timer_.reportStepNum(),
872 substep_timer.setReportStepStartTime(report_step_start_time);
873 substep_timer.setReportStepTotalTime(step_end_time);
874 substep_timer.setReportStepSubstepOffset(report_step_substep_offset);
876 current_time + timestep, step_end_time
881 reservoirCouplingSlave_().setFirstSubstepOfSyncTimestep(
true);
882 SubStepIteration<Solver> substepIteration{*
this, substep_timer, timestep, final_step};
883 const auto sub_steps_report = substepIteration.run();
884 report += sub_steps_report;
885 report_step_substep_offset += substep_timer.currentStepNum();
886 current_time += timestep;
901template<
class TypeTag>
902template<
class Solver>
903AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
904SubStepIteration(SubStepper<Solver>& substepper,
905 AdaptiveSimulatorTimer& substep_timer,
906 const double original_time_step,
908 : substepper_{substepper}
909 , substep_timer_{substep_timer}
910 , original_time_step_{original_time_step}
911 , final_step_{final_step}
912 , adaptive_time_stepping_{substepper.getAdaptiveTimerStepper()}
916template <
class TypeTag>
917template <
class Solver>
919AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
922 auto& simulator = solver_().model().simulator();
923 auto& problem = simulator.problem();
926 SimulatorReport report;
929 while (!this->substep_timer_.done()) {
933 maybeUpdateTuningAndTimeStep_();
935 const double dt = this->substep_timer_.currentStepLength();
936 if (timeStepVerbose_()) {
940 maybeUpdateLastSubstepOfSyncTimestep_(dt);
941 auto substep_report = runSubStep_();
942 markFirstSubStepAsFinished_();
944 if (substep_report.converged || checkContinueOnUnconvergedSolution_(dt)) {
945 Dune::Timer perfTimer;
948 problem.setSubStepReport(substep_report);
949 auto& full_report = adaptive_time_stepping_.report();
950 full_report += substep_report;
951 problem.setSimulationReport(full_report);
952 problem.endTimeStep();
953 substep_report.pre_post_time += perfTimer.stop();
955 report += substep_report;
957 OPM_TIMEBLOCK(convergenceSucceeded);
958 ++this->substep_timer_;
960 const int iterations = getNumIterations_(substep_report);
961 auto dt_estimate = timeStepControlComputeEstimate_(
962 dt, iterations, this->substep_timer_);
964 assert(dt_estimate > 0);
965 dt_estimate = maybeRestrictTimeStepGrowth_(dt, dt_estimate, restarts);
968 maybeReportSubStep_(substep_report);
969 if (this->final_step_ && this->substep_timer_.done()) {
974 report.success.output_write_time += writeOutput_();
978 checkTimeStepCanAdvance_(dt, dt_estimate);
979 setTimeStep_(dt_estimate);
981 report.success.converged = this->substep_timer_.done();
982 this->substep_timer_.setLastStepFailed(
false);
985 OPM_TIMEBLOCK(convergenceFailed);
986 report += substep_report;
987 this->substep_timer_.setLastStepFailed(
true);
988 checkTimeStepMaxRestartLimit_(restarts);
990 double new_time_step = restartFactor_() * dt;
991 if (substep_report.time_step_rejected) {
992 const double tol = Parameters::Get<Parameters::TimeStepControlTolerance>();
993 const double safetyFactor = Parameters::Get<Parameters::TimeStepControlSafetyFactor>();
994 const double temp_time_step = std::sqrt(safetyFactor * tol / solver_().model().relativeChange()) * dt;
995 if (temp_time_step < dt) {
996 new_time_step = temp_time_step;
999 checkTimeStepMinLimit_(new_time_step);
1000 bool wells_shut =
false;
1001 if (new_time_step > minTimeStepBeforeClosingWells_()) {
1002 chopTimeStep_(new_time_step);
1004 wells_shut = chopTimeStepOrCloseFailingWells_(new_time_step);
1013 problem.setNextTimeStepSize(this->substep_timer_.currentStepLength());
1015 updateSuggestedNextStep_();
1025template<
class TypeTag>
1026template<
class Solver>
1028AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1029checkContinueOnUnconvergedSolution_(
double dt)
const
1031 const bool continue_on_uncoverged_solution = ignoreConvergenceFailure_() && dt <= minTimeStep_();
1032 if (continue_on_uncoverged_solution && solverVerbose_()) {
1034 const auto msg = fmt::format(
1035 "Solver failed to converge but timestep {} is smaller or equal to {}\n"
1036 "which is the minimum threshold given by option --solver-min-time-step\n",
1039 OpmLog::problem(msg);
1041 return continue_on_uncoverged_solution;
1044template<
class TypeTag>
1045template<
class Solver>
1047AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1048checkTimeStepMaxRestartLimit_(
const int restarts)
const
1052 if (restarts >= solverRestartMax_()) {
1053 const auto msg = fmt::format(
1054 fmt::runtime(
"Solver failed to converge after cutting timestep {} times."), restarts
1056 if (solverVerbose_()) {
1060 throw TimeSteppingBreakdown{msg};
1064template<
class TypeTag>
1065template<
class Solver>
1067AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1068checkTimeStepCanAdvance_(
const double current_time_step,
1069 const double new_time_step)
const
1076 if (new_time_step >= minTimeStep_()) {
1084 if (this->substep_timer_.done() || new_time_step >= current_time_step) {
1089 fmt::format(
"Time step control proposed a step of {:.3E} DAYS, below the "
1090 "minimum of {:.3E} DAYS, while every substep converges. The "
1091 "run is not advancing past {:.6E} DAYS.",
1092 new_time_step / 86400.0,
1093 minTimeStep_() / 86400.0,
1094 this->substep_timer_.simulationTimeElapsed() / 86400.0);
1096 if (solverVerbose_()) {
1101 throw TimeSteppingBreakdown{msg};
1104template<
class TypeTag>
1105template<
class Solver>
1107AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1108checkTimeStepMinLimit_(
const double new_time_step)
const
1110 using Meas = UnitSystem::measure;
1113 if (new_time_step < minTimeStep_()) {
1114 std::string msg =
"Solver failed to converge after cutting timestep to ";
1115 if (Parameters::Get<Parameters::EnableTuning>()) {
1116 const UnitSystem& unit_system = solver_().model().simulator().vanguard().eclState().getDeckUnitSystem();
1118 "{:.3E} {}\nwhich is the minimum threshold given by the TUNING keyword\n",
1119 unit_system.from_si(Meas::time, minTimeStep_()),
1120 unit_system.name(Meas::time)
1125 "{:.3E} DAYS\nwhich is the minimum threshold given by option --solver-min-time-step\n",
1126 minTimeStep_() / 86400.0
1129 if (solverVerbose_()) {
1133 throw TimeSteppingBreakdown{msg};
1137template<
class TypeTag>
1138template<
class Solver>
1140AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1141chopTimeStep_(
const double new_time_step)
1143 setTimeStep_(new_time_step);
1144 if (solverVerbose_()) {
1145 const auto msg = fmt::format(fmt::runtime(
"{}\nTimestep chopped to {} days\n"),
1146 this->cause_of_failure_,
1147 unit::convert::to(this->substep_timer_.currentStepLength(), unit::day));
1148 OpmLog::problem(msg);
1152template<
class TypeTag>
1153template<
class Solver>
1155AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1156chopTimeStepOrCloseFailingWells_(
const double new_time_step)
1158 bool wells_shut =
false;
1165 const bool requireRepeatedFailures =
1166 new_time_step > (minTimeStepBeforeClosingWells_() * restartFactor_() * restartFactor_());
1167 const std::set<std::string> failing_wells =
1170 if (failing_wells.empty()) {
1172 chopTimeStep_(new_time_step);
1175 std::vector<std::string> shut_wells;
1176 for (
const auto& well : failing_wells) {
1177 const bool was_shut =
1178 solver_().model().wellModel().forceShutWellByName(well,
1179 this->substep_timer_.simulationTimeElapsed(),
1182 shut_wells.push_back(well);
1186 if (shut_wells.empty()) {
1187 for (
const auto& well : failing_wells) {
1188 const bool was_shut =
1189 solver_().model().wellModel().forceShutWellByName(well,
1190 this->substep_timer_.simulationTimeElapsed(),
1193 shut_wells.push_back(well);
1198 if (shut_wells.empty()) {
1199 chopTimeStep_(new_time_step);
1202 if (solverVerbose_()) {
1203 const std::string msg =
1204 fmt::format(fmt::runtime(
"\nProblematic well(s) were shut: {}"
1205 "(retrying timestep)\n"),
1206 fmt::join(shut_wells,
" "));
1207 OpmLog::problem(msg);
1214template<
class TypeTag>
1215template<
class Solver>
1216boost::posix_time::ptime
1217AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1218currentDateTime_()
const
1220 return simulatorTimer_().currentDateTime();
1223template<
class TypeTag>
1224template<
class Solver>
1226AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1227getNumIterations_(
const SimulatorReportSingle &substep_report)
const
1229 if (useNewtonIteration_()) {
1230 return substep_report.total_newton_iterations;
1233 return substep_report.total_linear_iterations;
1237template<
class TypeTag>
1238template<
class Solver>
1240AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1241growthFactor_()
const
1243 return this->adaptive_time_stepping_.growth_factor_;
1246template<
class TypeTag>
1247template<
class Solver>
1249AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1250ignoreConvergenceFailure_()
const
1252 return adaptive_time_stepping_.ignore_convergence_failure_;
1255template<
class TypeTag>
1256template<
class Solver>
1258AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1259isReservoirCouplingMaster_()
const
1261 return this->substepper_.isReservoirCouplingMaster_();
1264template<
class TypeTag>
1265template<
class Solver>
1267AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1268isReservoirCouplingSlave_()
const
1270 return this->substepper_.isReservoirCouplingSlave_();
1273template<
class TypeTag>
1274template<
class Solver>
1276AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1277markFirstSubStepAsFinished_()
const
1279#ifdef RESERVOIR_COUPLING_ENABLED
1283 if (isReservoirCouplingMaster_()) {
1284 reservoirCouplingMaster_().setFirstSubstepOfSyncTimestep(
false);
1286 else if (isReservoirCouplingSlave_()) {
1287 reservoirCouplingSlave_().setFirstSubstepOfSyncTimestep(
false);
1293template<
class TypeTag>
1294template<
class Solver>
1296AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1299 return this->adaptive_time_stepping_.max_growth_;
1302template<
class TypeTag>
1303template<
class Solver>
1305AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1306maybeReportSubStep_(SimulatorReportSingle substep_report)
const
1308 if (timeStepVerbose_()) {
1309 std::ostringstream ss;
1310 substep_report.reportStep(ss);
1311 OpmLog::info(ss.str());
1315template<
class TypeTag>
1316template<
class Solver>
1318AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1319maybeRestrictTimeStepGrowth_(
const double dt,
double dt_estimate,
const int restarts)
const
1322 dt_estimate = std::min(dt_estimate,
double(maxGrowth_() * dt));
1323 assert(dt_estimate > 0);
1326 dt_estimate = std::min(growthFactor_() * dt, dt_estimate);
1333template<
class TypeTag>
1334template<
class Solver>
1336AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1337maybeUpdateLastSubstepOfSyncTimestep_([[maybe_unused]]
const double dt)
1339#ifdef RESERVOIR_COUPLING_ENABLED
1344 if (isReservoirCouplingSlave_()) {
1346 this->substep_timer_.simulationTimeElapsed() + dt,
1347 this->substep_timer_.totalTime()
1349 reservoirCouplingSlave_().setLastSubstepOfSyncTimestep(is_last);
1357template<
class TypeTag>
1358template<
class Solver>
1360AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1361maybeUpdateTuningAndTimeStep_()
1369 const auto old_value = suggestedNextTimestep_();
1370 if (this->substepper_.maybeUpdateTuning_(this->substep_timer_.simulationTimeElapsed(),
1371 this->substep_timer_.currentStepLength(),
1372 this->substep_timer_.currentStepNum()))
1379 setTimeStep_(suggestedNextTimestep_());
1380 setSuggestedNextStep_(old_value);
1384template<
class TypeTag>
1385template<
class Solver>
1387AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1388minTimeStepBeforeClosingWells_()
const
1390 return this->adaptive_time_stepping_.min_time_step_before_shutting_problematic_wells_;
1393template<
class TypeTag>
1394template<
class Solver>
1396AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1399 return this->adaptive_time_stepping_.min_time_step_;
1402#ifdef RESERVOIR_COUPLING_ENABLED
1403template<
class TypeTag>
1404template<
class Solver>
1405ReservoirCouplingMaster<typename AdaptiveTimeStepping<TypeTag>::Scalar>&
1406AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1407reservoirCouplingMaster_()
const
1409 return this->substepper_.reservoirCouplingMaster_();
1412template<
class TypeTag>
1413template<
class Solver>
1414ReservoirCouplingSlave<typename AdaptiveTimeStepping<TypeTag>::Scalar>&
1415AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1416reservoirCouplingSlave_()
const
1418 return this->substepper_.reservoirCouplingSlave_();
1422template<
class TypeTag>
1423template<
class Solver>
1425AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1426restartFactor_()
const
1428 return this->adaptive_time_stepping_.restart_factor_;
1431template<
class TypeTag>
1432template<
class Solver>
1433SimulatorReportSingle
1434AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1438 SimulatorReportSingle substep_report;
1440 auto handleFailure = [
this, &substep_report]
1441 (
const std::string& failure_reason,
const std::exception& e,
bool log_exception =
true)
1443 substep_report = solver_().failureReport();
1447 this->cause_of_failure_ = failure_reason;
1448 if (
const std::string what = e.what();
1449 !what.empty() && what != failure_reason)
1451 this->cause_of_failure_ +=
" (" + what +
")";
1453 if (log_exception && solverVerbose_()) {
1454 OpmLog::debug(std::string(
"Caught Exception: ") + e.what());
1459 substep_report = solver_().step(this->substep_timer_, &this->adaptive_time_stepping_.timeStepControl());
1460 if (solverVerbose_()) {
1462 OpmLog::debug(
"Overall linear iterations used: "
1466 catch (
const TooManyIterations& e) {
1467 handleFailure(
"Solver convergence failure - Iteration limit reached", e);
1469 catch (
const TimeSteppingBreakdown& e) {
1470 handleFailure(e.what(), e);
1472 catch (
const ConvergenceMonitorFailure& e) {
1473 handleFailure(
"Convergence monitor failure", e,
false);
1475 catch (
const LinearSolverProblem& e) {
1476 handleFailure(
"Linear solver convergence failure", e);
1478 catch (
const NumericalProblem& e) {
1479 handleFailure(
"Solver convergence failure - Numerical problem encountered", e);
1481 catch (
const std::runtime_error& e) {
1482 handleFailure(
"Runtime error encountered", e);
1484 catch (
const Dune::ISTLError& e) {
1485 handleFailure(
"ISTL error - Time step too large", e);
1487 catch (
const Dune::MatrixBlockError& e) {
1488 handleFailure(
"Matrix block error", e);
1491 return substep_report;
1494template<
class TypeTag>
1495template<
class Solver>
1497AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1498setTimeStep_(
double dt_estimate)
1500 this->substep_timer_.provideTimeStepEstimate(dt_estimate);
1503template<
class TypeTag>
1504template<
class Solver>
1506AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1509 return this->substepper_.solver_;
1513template<
class TypeTag>
1514template<
class Solver>
1516AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1517solverRestartMax_()
const
1519 return this->adaptive_time_stepping_.solver_restart_max_;
1522template<
class TypeTag>
1523template<
class Solver>
1525AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1526setSuggestedNextStep_(
double step)
1528 this->adaptive_time_stepping_.setSuggestedNextStep(step);
1531template <
class TypeTag>
1532template <
class Solver>
1533const SimulatorTimer&
1534AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1535simulatorTimer_()
const
1537 return this->substepper_.simulator_timer_;
1540template <
class TypeTag>
1541template <
class Solver>
1543AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1544solverVerbose_()
const
1546 return this->adaptive_time_stepping_.solver_verbose_;
1549template<
class TypeTag>
1550template<
class Solver>
1551boost::posix_time::ptime
1552AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1553startDateTime_()
const
1555 return simulatorTimer_().startDateTime();
1558template <
class TypeTag>
1559template <
class Solver>
1561AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1562suggestedNextTimestep_()
const
1564 return this->adaptive_time_stepping_.suggestedNextStep();
1567template <
class TypeTag>
1568template <
class Solver>
1570AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1571timeStepControlComputeEstimate_(
const double dt,
const int iterations,
1572 const AdaptiveSimulatorTimer& substepTimer)
const
1575 const SolutionTimeErrorSolverWrapper<Solver> relative_change{solver_()};
1576 return this->adaptive_time_stepping_.time_step_control_->computeTimeStepSize(
1577 dt, iterations, relative_change, substepTimer);
1580template <
class TypeTag>
1581template <
class Solver>
1583AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1584timeStepVerbose_()
const
1586 return this->adaptive_time_stepping_.timestep_verbose_;
1596template <
class TypeTag>
1597template <
class Solver>
1599AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1600updateSuggestedNextStep_()
1602 auto suggested_next_step = this->substep_timer_.currentStepLength();
1603 if (! std::isfinite(suggested_next_step)) {
1604 suggested_next_step = this->original_time_step_;
1606 if (timeStepVerbose_()) {
1607 std::ostringstream ss;
1608 this->substep_timer_.report(ss);
1609 ss <<
"Suggested next step size = "
1610 << unit::convert::to(suggested_next_step, unit::day) <<
" (days)" << std::endl;
1611 OpmLog::debug(ss.str());
1613 setSuggestedNextStep_(suggested_next_step);
1616template <
class TypeTag>
1617template <
class Solver>
1619AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1620useNewtonIteration_()
const
1622 return this->adaptive_time_stepping_.use_newton_iteration_;
1625template <
class TypeTag>
1626template <
class Solver>
1628AdaptiveTimeStepping<TypeTag>::SubStepIteration<Solver>::
1631 time::StopWatch perf_timer;
1633 auto& problem = solver_().model().simulator().problem();
1634 problem.writeOutput(
true);
1635 return perf_timer.secsSinceStart();
1642template<
class TypeTag>
1643template<
class Solver>
1644AdaptiveTimeStepping<TypeTag>::
1645SolutionTimeErrorSolverWrapper<Solver>::
1646SolutionTimeErrorSolverWrapper(
const Solver& solver)
1650template<
class TypeTag>
1651template<
class Solver>
1652double AdaptiveTimeStepping<TypeTag>::SolutionTimeErrorSolverWrapper<Solver>::relativeChange()
const
1655 return solver_.model().relativeChange();
Adaptive time-stepping coordinator for the black-oil simulator.
Definition: AdaptiveTimeStepping.hpp:93
double max_growth_
factor that limits the maximum growth of a time step
Definition: AdaptiveTimeStepping.hpp:435
double max_time_step_
maximal allowed time step size in days
Definition: AdaptiveTimeStepping.hpp:436
bool solver_verbose_
solver verbosity
Definition: AdaptiveTimeStepping.hpp:440
int solver_restart_max_
how many restart of solver are allowed
Definition: AdaptiveTimeStepping.hpp:439
double timestep_after_event_
suggested size of timestep after an event
Definition: AdaptiveTimeStepping.hpp:444
void init_(const UnitSystem &unitSystem)
Definition: AdaptiveTimeStepping_impl.hpp:429
void setSuggestedNextStep(const double x)
Set the suggested length for the next substep [s].
Definition: AdaptiveTimeStepping_impl.hpp:300
double suggestedNextStep() const
Definition: AdaptiveTimeStepping_impl.hpp:308
bool operator==(const AdaptiveTimeStepping< TypeTag > &rhs) const
Definition: AdaptiveTimeStepping_impl.hpp:140
static AdaptiveTimeStepping< TypeTag > serializationTestObjectSimple()
Definition: AdaptiveTimeStepping_impl.hpp:283
bool ignore_convergence_failure_
continue instead of stop when minimum time step is reached
Definition: AdaptiveTimeStepping.hpp:438
void serializeOp(Serializer &serializer)
Definition: AdaptiveTimeStepping_impl.hpp:212
void updateTUNING(double max_next_tstep, const Tuning &tuning)
Apply TUNING keyword parameters.
Definition: AdaptiveTimeStepping_impl.hpp:337
TimeStepControlType time_step_control_type_
type of time step control object
Definition: AdaptiveTimeStepping.hpp:431
const TimeStepControlInterface & timeStepControl() const
Definition: AdaptiveTimeStepping_impl.hpp:316
std::function< bool(double elapsed, double substep_length, int sub_step_number)> TuningUpdateCallback
Callback invoked at the start of each substep to apply TUNING, NEXTSTEP (via ACTIONX),...
Definition: AdaptiveTimeStepping.hpp:119
bool full_timestep_initially_
beginning with the size of the time step from data file
Definition: AdaptiveTimeStepping.hpp:443
SimulatorReport step(const SimulatorTimer &simulator_timer, Solver &solver, const bool is_event, const TuningUpdateCallback &tuning_updater)
Run one report step by orchestrating adaptive substepping.
Definition: AdaptiveTimeStepping_impl.hpp:197
double growth_factor_
factor to multiply time step when solver recovered from failed convergence
Definition: AdaptiveTimeStepping.hpp:434
double restart_factor_
factor to multiply time step with when solver fails to converge
Definition: AdaptiveTimeStepping.hpp:433
double min_time_step_
minimal allowed time step size before throwing
Definition: AdaptiveTimeStepping.hpp:437
void updateNEXTSTEP(double max_next_tstep)
Set suggested_next_timestep_ to max_next_tstep iff max_next_tstep > 0.
Definition: AdaptiveTimeStepping_impl.hpp:325
static AdaptiveTimeStepping< TypeTag > serializationTestObjectHardcoded()
Definition: AdaptiveTimeStepping_impl.hpp:259
AdaptiveTimeStepping()=default
TimeStepController time_step_control_
time step control object
Definition: AdaptiveTimeStepping.hpp:432
static AdaptiveTimeStepping< TypeTag > serializationTestObjectPIDIt()
Definition: AdaptiveTimeStepping_impl.hpp:275
double min_time_step_before_shutting_problematic_wells_
< shut problematic wells when time step size in days are less than this
Definition: AdaptiveTimeStepping.hpp:448
static AdaptiveTimeStepping< TypeTag > serializationTestObject3rdOrder()
Definition: AdaptiveTimeStepping_impl.hpp:291
SimulatorReport & report()
Definition: AdaptiveTimeStepping_impl.hpp:251
static AdaptiveTimeStepping< TypeTag > serializationTestObjectPID()
Definition: AdaptiveTimeStepping_impl.hpp:267
static void registerParameters()
Definition: AdaptiveTimeStepping_impl.hpp:186
bool use_newton_iteration_
use newton iteration count for adaptive time step control
Definition: AdaptiveTimeStepping.hpp:445
Definition: SimulatorTimer.hpp:38
Definition: TimeStepControlInterface.hpp:51
auto Get(bool errorIfNotRegistered=true)
Retrieve a runtime parameter.
Definition: parametersystem.hpp:191
void logTimer(const AdaptiveSimulatorTimer &substep_timer)
void registerAdaptiveParameters()
std::set< std::string > consistentlyFailingWells(const std::vector< StepReport > &sr, bool requireRepeatedFailures)
std::tuple< TimeStepControlType, std::unique_ptr< TimeStepControlInterface >, bool > createController(const UnitSystem &unitSystem)
Definition: blackoilbioeffectsmodules.hh:45
std::string to_string(const ConvergenceReport::ReservoirFailure::Type t)
This file provides the infrastructure to retrieve run-time parameters.
static bool compare_gt_or_eq(double a, double b)
Determines if a is greater than b within the specified tolerance.
Definition: SimulatorReport.hpp:125