From d63e0cbd5093dd0749017032222ef6849328151f Mon Sep 17 00:00:00 2001 From: Erik Hide Saeternes Date: Thu, 15 May 2025 13:56:25 +0200 Subject: [PATCH 1/7] Major updates to general 3rd order controller: different reduction in time step size if tolerance test fails (than if numerical difficulties), possibility to set controller parameters, option to use control error filtering as tolerance test, adjustments/fixes to the controller code --- opm/simulators/flow/NonlinearSolver.hpp | 3 +- .../timestepping/AdaptiveTimeStepping.cpp | 21 ++- .../timestepping/AdaptiveTimeStepping.hpp | 3 + .../AdaptiveTimeStepping_impl.hpp | 7 +- .../timestepping/SimulatorReport.cpp | 3 +- .../timestepping/SimulatorReport.hpp | 2 + .../timestepping/TimeStepControl.cpp | 142 ++++++++++++------ .../timestepping/TimeStepControl.hpp | 35 +++-- .../timestepping/TimeStepControlInterface.hpp | 2 +- 9 files changed, 156 insertions(+), 62 deletions(-) diff --git a/opm/simulators/flow/NonlinearSolver.hpp b/opm/simulators/flow/NonlinearSolver.hpp index 4ac3f3b7c5..0f7e8016f3 100644 --- a/opm/simulators/flow/NonlinearSolver.hpp +++ b/opm/simulators/flow/NonlinearSolver.hpp @@ -176,8 +176,9 @@ struct NonlinearSolverParameters OPM_THROW_NOLOG(TooManyIterations, msg); } auto relativeChange = model_->relativeChange(); - if (timeStepControl && !timeStepControl->timeStepAccepted(relativeChange)) { + if (timeStepControl && !timeStepControl->timeStepAccepted(relativeChange, timer.currentStepLength())) { report.converged = false; + report.time_step_rejected = true; failureReport_ = report; std::string msg = "Relative change in solution for time step was " + std::to_string(relativeChange) + ", which is larger than the tolerance accepted by the timestepping algorithm."; diff --git a/opm/simulators/timestepping/AdaptiveTimeStepping.cpp b/opm/simulators/timestepping/AdaptiveTimeStepping.cpp index da9fdbb881..37c5438d2d 100644 --- a/opm/simulators/timestepping/AdaptiveTimeStepping.cpp +++ b/opm/simulators/timestepping/AdaptiveTimeStepping.cpp @@ -142,6 +142,19 @@ void registerAdaptiveParameters() Parameters::Register ("(Only applicable for the general 3rd order controller.) Include rejection of completed " "time steps if the relative change is larger than the time step control tolerance"); + Parameters::Register + ("(Only applicable for the general 3rd order controller.) Ways to decide if the time step " + "should be rejected. Options: 'standard' and 'control-error-filtering'. The standard " + "version compares relative change to tolerance directly to decide if the time step should " + "be rejected, while the control-error-filtering version compares the relative change " + "between the current and proposed next time steps to the max-reduction-time-step parameter."); + Parameters::Register + ("(Only applicable for the general 3rd order controller, using " + "time-step-control-tolerance-test-version choice 'control-error-filtering') If the relative " + "change in time step size for the next time step is larger than this parameter, the time " + "step will be rejected."); + Parameters::Register + ("TODO!"); } std::tuple, bool> @@ -222,11 +235,17 @@ createController(const UnitSystem& unitSystem) [tol]() { const double safetyFactor = Parameters::Get(); // 0.8 const bool rejectCompletedStep = Parameters::Get(); // false + const std::string toleranceTestVersion = Parameters::Get(); // "standard" + const double maxReductionTimeStep = Parameters::Get(); // 0.1 + const std::string parameters = Parameters::Get(); // "" return RetVal{ TimeStepControlType::General3rdOrder, std::make_unique(tol, safetyFactor, - rejectCompletedStep), + rejectCompletedStep, + toleranceTestVersion, + maxReductionTimeStep, + parameters), false }; }}, diff --git a/opm/simulators/timestepping/AdaptiveTimeStepping.hpp b/opm/simulators/timestepping/AdaptiveTimeStepping.hpp index 4754c4d99c..d65c7c8abb 100644 --- a/opm/simulators/timestepping/AdaptiveTimeStepping.hpp +++ b/opm/simulators/timestepping/AdaptiveTimeStepping.hpp @@ -49,6 +49,9 @@ struct MinTimeStepBeforeShuttingProblematicWellsInDays { static constexpr double struct MinTimeStepBasedOnNewtonIterations { static constexpr double value = 0.0; }; struct TimeStepControlSafetyFactor { static constexpr double value = 0.8; }; struct TimeStepControlRejectCompletedStep { static constexpr bool value = false; }; +struct TimeStepControlToleranceTestVersion { static constexpr auto value = "standard"; }; +struct TimeStepControlMaxReductionTimeStep { static constexpr double value = 0.1; }; +struct TimeStepControlParameters { static constexpr auto value = "0.125;0.25;0.125;0.375;0.125"; }; } // namespace Opm::Parameters diff --git a/opm/simulators/timestepping/AdaptiveTimeStepping_impl.hpp b/opm/simulators/timestepping/AdaptiveTimeStepping_impl.hpp index 39c1eb6e87..40b90c3e47 100644 --- a/opm/simulators/timestepping/AdaptiveTimeStepping_impl.hpp +++ b/opm/simulators/timestepping/AdaptiveTimeStepping_impl.hpp @@ -825,7 +825,12 @@ run() this->substep_timer_.setLastStepFailed(true); checkTimeStepMaxRestartLimit_(restarts); - const double new_time_step = restartFactor_() * dt; + double new_time_step = restartFactor_() * dt; + if (substep_report.time_step_rejected) { + const double tol = Parameters::Get(); + const double safetyFactor = Parameters::Get(); + new_time_step = std::sqrt(safetyFactor * tol / solver_().model().relativeChange()) * dt; + } checkTimeStepMinLimit_(new_time_step); bool wells_shut = false; if (new_time_step > minTimeStepBeforeClosingWells_()) { diff --git a/opm/simulators/timestepping/SimulatorReport.cpp b/opm/simulators/timestepping/SimulatorReport.cpp index d82c7f5188..01cbd99bec 100644 --- a/opm/simulators/timestepping/SimulatorReport.cpp +++ b/opm/simulators/timestepping/SimulatorReport.cpp @@ -35,7 +35,7 @@ namespace Opm return SimulatorReportSingle{1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0, 10.0, 11.0, 12.0, 13, 14, 15, 16, 17, 18, - true, false, 19, 20.0, 21.0, + true, false, false, 19, 20.0, 21.0, 22, 23, 24, 25, 26, 27, 28}; } @@ -60,6 +60,7 @@ namespace Opm this->min_linear_iterations == rhs.min_linear_iterations && this->max_linear_iterations == rhs.max_linear_iterations && this->converged == rhs.converged && + this->time_step_rejected == rhs.time_step_rejected && this->well_group_control_changed == rhs.well_group_control_changed && this->exit_status == rhs.exit_status && this->global_time == rhs.global_time && diff --git a/opm/simulators/timestepping/SimulatorReport.hpp b/opm/simulators/timestepping/SimulatorReport.hpp index c8aba76003..37ecc657ce 100644 --- a/opm/simulators/timestepping/SimulatorReport.hpp +++ b/opm/simulators/timestepping/SimulatorReport.hpp @@ -53,6 +53,7 @@ namespace Opm unsigned int max_linear_iterations = 0; bool converged = false; + bool time_step_rejected = false; bool well_group_control_changed = false; int exit_status = EXIT_SUCCESS; @@ -101,6 +102,7 @@ namespace Opm serializer(min_linear_iterations); serializer(max_linear_iterations); serializer(converged); + serializer(time_step_rejected); serializer(well_group_control_changed); serializer(exit_status); serializer(global_time); diff --git a/opm/simulators/timestepping/TimeStepControl.cpp b/opm/simulators/timestepping/TimeStepControl.cpp index fbe1fc2a65..b947ac088a 100644 --- a/opm/simulators/timestepping/TimeStepControl.cpp +++ b/opm/simulators/timestepping/TimeStepControl.cpp @@ -28,6 +28,7 @@ #include #include #include +#include #include #include @@ -292,57 +293,48 @@ namespace Opm General3rdOrderController::General3rdOrderController( const double tolerance, const double safetyFactor, const bool rejectCompletedStep, + const std::string toleranceTestVersion, + const double maxReductionTimeStep, + const std::string parameters, const bool verbose) : tolerance_( tolerance ) , safetyFactor_( safetyFactor ) , rejectCompletedStep_( rejectCompletedStep ) - , errors_( 3, tolerance_ ) - , timeSteps_ ( 3, 1.0 ) + , toleranceTestVersion_( toleranceTestVersion ) + , maxReductionTimeStep_( maxReductionTimeStep ) , verbose_( verbose ) - {} + { + std::istringstream parameters_stream(parameters); + std::string value; + int counter = 0; + while (std::getline(parameters_stream, value, ';')) { + if (counter < 3) { + beta_[counter] = std::stod(value); + } + else { + alpha_[counter - 3] = std::stod(value); + } + counter++; + } + } General3rdOrderController General3rdOrderController::serializationTestObject() { - General3rdOrderController result(1.0, 0.8, true); + General3rdOrderController result(1.0, 2.0, false, "", 3.0, "", false); result.errors_ = {2.0, 3.0}; return result; } double General3rdOrderController:: - computeTimeStepSize(const double dt, const int /*iterations */, const RelativeChangeInterface& relChange, const AdaptiveSimulatorTimer& substepTimer) const + computeTimeStepSize(const double dt, const int /* iterations */, const RelativeChangeInterface& /* relChange */, const AdaptiveSimulatorTimer& substepTimer) const { - // Shift errors and time steps - for( int i = 0; i < 2; ++i ) + // Use an I controller after report time steps + if (substepTimer.currentStepNum() < 3) { - errors_[i] = errors_[i+1]; - timeSteps_[i] = timeSteps_[i+1]; - } - - // Store new error and time step - const double error = relChange.relativeChange(); - errors_[2] = error; - timeSteps_[2] = dt; - for( int i = 0; i < 2; ++i ) - { - assert(std::isfinite(errors_[i])); - } - - if (errors_[0] == 0 || errors_[1] == 0 || errors_[2] == 0.) - { - if ( verbose_ ) - OpmLog::info("The solution between time steps does not change, there is no time step constraint from the controller."); - return std::numeric_limits::max(); - } - // Use an I controller after report time steps or chopped time steps - else if (substepTimer.currentStepNum() < 3 || substepTimer.lastStepFailed() || counterSinceFailure_ > 0) - { - if (substepTimer.lastStepFailed() || counterSinceFailure_ > 0) - counterSinceFailure_++; - if (counterSinceFailure_ > 1) - counterSinceFailure_ = 0; - const double newDt = dt * std::pow(safetyFactor_ * tolerance_ / errors_[2], 0.35); + controllerVersion_ = "I-controller"; + const double newDt = dt * timeStepFactor(errors_, timeSteps_); if( verbose_ ) OpmLog::info(fmt::format("Computed step size (pow): {} days", unit::convert::to( newDt, unit::day ))); return newDt; @@ -350,25 +342,78 @@ namespace Opm // Use the general third order controller for all other time steps else { - const std::array beta = { 0.125, 0.25, 0.125 }; - const std::array alpha = { 0.375, 0.125 }; - const double newDt = dt * std::pow(safetyFactor_ * tolerance_ / errors_[2], beta[0]) * - std::pow(safetyFactor_ * tolerance_ / errors_[1], beta[1]) * - std::pow(safetyFactor_ * tolerance_ / errors_[0], beta[2]) * - std::pow(timeSteps_[2] / timeSteps_[1], -alpha[0]) * - std::pow(timeSteps_[1] / timeSteps_[0], -alpha[1]); + controllerVersion_ = "3rdOrder"; + const double newDt = dt * timeStepFactor(errors_, timeSteps_); if( verbose_ ) OpmLog::info(fmt::format("Computed step size (pow): {} days", unit::convert::to( newDt, unit::day ))); return newDt; } } - bool General3rdOrderController:: - timeStepAccepted(const double error) const + double General3rdOrderController:: + timeStepFactor(const std::array errors, const std::array timeSteps) const { - if (rejectCompletedStep_ && error > tolerance_) - return false; - return true; + if (controllerVersion_ == "3rdOrder") + { + return std::pow(safetyFactor_ * tolerance_ / errors[2], beta_[0]) * + std::pow(safetyFactor_ * tolerance_ / errors[1], beta_[1]) * + std::pow(safetyFactor_ * tolerance_ / errors[0], beta_[2]) * + std::pow(timeSteps[2] / timeSteps[1], -alpha_[0]) * + std::pow(timeSteps[1] / timeSteps[0], -alpha_[1]); + } + else if (controllerVersion_ == "I-controller") + { + return std::pow(safetyFactor_ * tolerance_ / errors[2], 0.35); + } + } + + bool General3rdOrderController:: + timeStepAccepted(const double error, const double timeStep) const + { + bool acceptTimeStep = true; + + // Reject time step if chosen tolerance test version fails + if (toleranceTestVersion_ == "standard") + { + if (rejectCompletedStep_ && error > tolerance_) { acceptTimeStep = false; } + } + else if (toleranceTestVersion_ == "control-error-filtering") + { + const std::array tempErrors{errors_[0], errors_[1], error}; + const std::array tempTimeSteps{timeSteps_[0], timeSteps_[1], timeStep}; + double stepFactor = timeStepFactor(tempErrors, tempTimeSteps); + if (rejectCompletedStep_ && stepFactor < maxReductionTimeStep_) { acceptTimeStep = false; } + } + else + { + OPM_THROW(std::runtime_error, "Unsupported tolerance test version: " + toleranceTestVersion_); + } + + if (acceptTimeStep) + { + // Shift errors and time steps + for( int i = 0; i < 2; ++i ) + { + errors_[i] = errors_[i+1]; + timeSteps_[i] = timeSteps_[i+1]; + } + + // Store new error and time step + errors_[2] = error; + timeSteps_[2] = timeStep; + + for( int i = 0; i < 2; ++i ) + assert(std::isfinite(errors_[i])); + + if (errors_[0] == 0 || errors_[1] == 0 || errors_[2] == 0.) + { + if ( verbose_ ) + OpmLog::info("The solution between time steps does not change, there is no time step constraint from the controller."); + return std::numeric_limits::max(); + } + } + + return acceptTimeStep; } bool General3rdOrderController::operator==(const General3rdOrderController& ctrl) const @@ -378,6 +423,11 @@ namespace Opm this->rejectCompletedStep_ == ctrl.rejectCompletedStep_ && this->errors_ == ctrl.errors_ && this->timeSteps_ == ctrl.timeSteps_ && + this->beta_ == ctrl.beta_ && + this->alpha_ == ctrl.alpha_ && + this->controllerVersion_ == ctrl.controllerVersion_ && + this->toleranceTestVersion_ == ctrl.toleranceTestVersion_ && + this->maxReductionTimeStep_ == ctrl.maxReductionTimeStep_ && this->verbose_ == ctrl.verbose_; } diff --git a/opm/simulators/timestepping/TimeStepControl.hpp b/opm/simulators/timestepping/TimeStepControl.hpp index 370ee179cb..b21b731817 100644 --- a/opm/simulators/timestepping/TimeStepControl.hpp +++ b/opm/simulators/timestepping/TimeStepControl.hpp @@ -66,7 +66,7 @@ namespace Opm const RelativeChangeInterface& /* relativeChange */, const AdaptiveSimulatorTimer& /* substepTimer */ ) const override; - bool timeStepAccepted(const double /* error */) const override { return true; } + bool timeStepAccepted(const double /* error */, const double /* timeStep */) const override { return true; } template void serializeOp(Serializer& serializer) @@ -119,7 +119,7 @@ namespace Opm const RelativeChangeInterface& relativeChange, const AdaptiveSimulatorTimer& /* substepTimer */ ) const override; - bool timeStepAccepted(const double /* error */) const override { return true; } + bool timeStepAccepted(const double /* error */, const double /* timeStep */) const override { return true; } template void serializeOp(Serializer& serializer) @@ -170,7 +170,7 @@ namespace Opm const RelativeChangeInterface& relativeChange, const AdaptiveSimulatorTimer& /* substepTimer */ ) const override; - bool timeStepAccepted(const double /* error */) const override { return true; } + bool timeStepAccepted(const double /* error */, const double /* timeStep */) const override { return true; } template void serializeOp(Serializer& serializer) @@ -204,16 +204,21 @@ namespace Opm General3rdOrderController( const double tolerance = 1e-3, const double safetyFactor = 0.8, const bool rejectCompletedStep = false, - const bool verbose = false ); + const std::string toleranceTestVersion = "", + const double maxReductionTimeStep = 0.1, + const std::string parameters = "", + const bool verbose = false); static General3rdOrderController serializationTestObject(); double computeTimeStepSize(const double dt, - const int /*iterations */, - const RelativeChangeInterface& relativeChange, + const int /* iterations */, + const RelativeChangeInterface& /* relativeChange */, const AdaptiveSimulatorTimer& substepTimer) const override; - bool timeStepAccepted(const double error) const override; + double timeStepFactor(const std::array errors, const std::array timeSteps) const; + + bool timeStepAccepted(const double error, const double timeStep) const override; template void serializeOp(Serializer& serializer) @@ -223,6 +228,10 @@ namespace Opm serializer(rejectCompletedStep_); serializer(errors_); serializer(timeSteps_); + serializer(beta_); + serializer(alpha_); + serializer(controllerVersion_); + serializer(toleranceTestVersion_); serializer(verbose_); } @@ -233,9 +242,13 @@ namespace Opm const double tolerance_ = 1e-3; const double safetyFactor_ = 0.8; const bool rejectCompletedStep_ = false; - mutable std::vector errors_{}; - mutable std::vector timeSteps_{}; - mutable int counterSinceFailure_ = 0; + mutable std::array errors_{}; + mutable std::array timeSteps_{}; + mutable std::array beta_{}; + mutable std::array alpha_{}; + mutable std::string controllerVersion_ = "I-controller"; + const std::string toleranceTestVersion_ = ""; + const double maxReductionTimeStep_ = 0.1; const bool verbose_ = false; }; @@ -266,7 +279,7 @@ namespace Opm const RelativeChangeInterface& /*relativeChange */, const AdaptiveSimulatorTimer& substepTimer) const override; - bool timeStepAccepted(const double /* error */) const override { return true; } + bool timeStepAccepted(const double /* error */, const double /* timeStep */) const override { return true; } template void serializeOp(Serializer& serializer) diff --git a/opm/simulators/timestepping/TimeStepControlInterface.hpp b/opm/simulators/timestepping/TimeStepControlInterface.hpp index 240eac2baa..6e34443d7d 100644 --- a/opm/simulators/timestepping/TimeStepControlInterface.hpp +++ b/opm/simulators/timestepping/TimeStepControlInterface.hpp @@ -60,7 +60,7 @@ namespace Opm /// \return suggested time step size for the next step virtual double computeTimeStepSize( const double dt, const int iterations, const RelativeChangeInterface& relativeChange , const AdaptiveSimulatorTimer& substepTimer) const = 0; - virtual bool timeStepAccepted(const double error) const = 0; + virtual bool timeStepAccepted(const double error, const double timeStep) const = 0; /// virtual destructor (empty) virtual ~TimeStepControlInterface () {} From f58413225f983268ca08f2b110ae2f573c796103 Mon Sep 17 00:00:00 2001 From: Erik Hide Saeternes Date: Thu, 15 May 2025 14:10:02 +0200 Subject: [PATCH 2/7] Updated explanation for parameters --- opm/simulators/timestepping/AdaptiveTimeStepping.cpp | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) diff --git a/opm/simulators/timestepping/AdaptiveTimeStepping.cpp b/opm/simulators/timestepping/AdaptiveTimeStepping.cpp index 37c5438d2d..462b65ff7c 100644 --- a/opm/simulators/timestepping/AdaptiveTimeStepping.cpp +++ b/opm/simulators/timestepping/AdaptiveTimeStepping.cpp @@ -147,14 +147,14 @@ void registerAdaptiveParameters() "should be rejected. Options: 'standard' and 'control-error-filtering'. The standard " "version compares relative change to tolerance directly to decide if the time step should " "be rejected, while the control-error-filtering version compares the relative change " - "between the current and proposed next time steps to the max-reduction-time-step parameter."); + "in time step size to the max-reduction-time-step parameter."); Parameters::Register - ("(Only applicable for the general 3rd order controller, using " - "time-step-control-tolerance-test-version choice 'control-error-filtering') If the relative " - "change in time step size for the next time step is larger than this parameter, the time " - "step will be rejected."); + ("(Only applicable for the general 3rd order controller, using 'control-error-filtering' " + "as time-step-control-tolerance-test-version) If the (proposed) relative change in time " + "step size is larger than this parameter, the time step will be rejected."); Parameters::Register - ("TODO!"); + ("(Only applicable for the general 3rd order controller.) Parameters for the general " + "3rd order controller. Should be given as \"beta_1;beta_2;beta_3;alpha_2;alpha_3\"."); } std::tuple, bool> From e017c722670d756a89b0d89ea0c0958ccb11911f Mon Sep 17 00:00:00 2001 From: Erik Hide Saeternes Date: Thu, 15 May 2025 14:11:52 +0200 Subject: [PATCH 3/7] Tiny fix to parameter explanation --- opm/simulators/timestepping/AdaptiveTimeStepping.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/opm/simulators/timestepping/AdaptiveTimeStepping.cpp b/opm/simulators/timestepping/AdaptiveTimeStepping.cpp index 462b65ff7c..f919763f57 100644 --- a/opm/simulators/timestepping/AdaptiveTimeStepping.cpp +++ b/opm/simulators/timestepping/AdaptiveTimeStepping.cpp @@ -154,7 +154,7 @@ void registerAdaptiveParameters() "step size is larger than this parameter, the time step will be rejected."); Parameters::Register ("(Only applicable for the general 3rd order controller.) Parameters for the general " - "3rd order controller. Should be given as \"beta_1;beta_2;beta_3;alpha_2;alpha_3\"."); + "3rd order controller. Should be given as 'beta_1;beta_2;beta_3;alpha_2;alpha_3'."); } std::tuple, bool> From dfda88913f73dcd5c18d7a506ecb8ff36ec71e65 Mon Sep 17 00:00:00 2001 From: Erik Hide Saeternes Date: Mon, 19 May 2025 13:13:21 +0200 Subject: [PATCH 4/7] Fixed bugs --- .../timestepping/TimeStepControl.cpp | 35 +++++++++++++------ .../timestepping/TimeStepControl.hpp | 3 +- 2 files changed, 26 insertions(+), 12 deletions(-) diff --git a/opm/simulators/timestepping/TimeStepControl.cpp b/opm/simulators/timestepping/TimeStepControl.cpp index b947ac088a..971244517c 100644 --- a/opm/simulators/timestepping/TimeStepControl.cpp +++ b/opm/simulators/timestepping/TimeStepControl.cpp @@ -333,19 +333,23 @@ namespace Opm // Use an I controller after report time steps if (substepTimer.currentStepNum() < 3) { - controllerVersion_ = "I-controller"; + controllerVersion_ = 0; // "I-controller" const double newDt = dt * timeStepFactor(errors_, timeSteps_); - if( verbose_ ) + if (verbose_) + { OpmLog::info(fmt::format("Computed step size (pow): {} days", unit::convert::to( newDt, unit::day ))); + } return newDt; } // Use the general third order controller for all other time steps else { - controllerVersion_ = "3rdOrder"; + controllerVersion_ = 1; // "3rdOrder"; const double newDt = dt * timeStepFactor(errors_, timeSteps_); - if( verbose_ ) + if (verbose_) + { OpmLog::info(fmt::format("Computed step size (pow): {} days", unit::convert::to( newDt, unit::day ))); + } return newDt; } } @@ -353,7 +357,8 @@ namespace Opm double General3rdOrderController:: timeStepFactor(const std::array errors, const std::array timeSteps) const { - if (controllerVersion_ == "3rdOrder") + // "3rdOrder" + if (controllerVersion_ == 1) { return std::pow(safetyFactor_ * tolerance_ / errors[2], beta_[0]) * std::pow(safetyFactor_ * tolerance_ / errors[1], beta_[1]) * @@ -361,10 +366,8 @@ namespace Opm std::pow(timeSteps[2] / timeSteps[1], -alpha_[0]) * std::pow(timeSteps[1] / timeSteps[0], -alpha_[1]); } - else if (controllerVersion_ == "I-controller") - { - return std::pow(safetyFactor_ * tolerance_ / errors[2], 0.35); - } + // "I-controller" + return std::pow(safetyFactor_ * tolerance_ / errors[2], 0.35); } bool General3rdOrderController:: @@ -375,14 +378,20 @@ namespace Opm // Reject time step if chosen tolerance test version fails if (toleranceTestVersion_ == "standard") { - if (rejectCompletedStep_ && error > tolerance_) { acceptTimeStep = false; } + if (rejectCompletedStep_ && error > tolerance_) + { + acceptTimeStep = false; + } } else if (toleranceTestVersion_ == "control-error-filtering") { const std::array tempErrors{errors_[0], errors_[1], error}; const std::array tempTimeSteps{timeSteps_[0], timeSteps_[1], timeStep}; double stepFactor = timeStepFactor(tempErrors, tempTimeSteps); - if (rejectCompletedStep_ && stepFactor < maxReductionTimeStep_) { acceptTimeStep = false; } + if (rejectCompletedStep_ && stepFactor < maxReductionTimeStep_) + { + acceptTimeStep = false; + } } else { @@ -403,12 +412,16 @@ namespace Opm timeSteps_[2] = timeStep; for( int i = 0; i < 2; ++i ) + { assert(std::isfinite(errors_[i])); + } if (errors_[0] == 0 || errors_[1] == 0 || errors_[2] == 0.) { if ( verbose_ ) + { OpmLog::info("The solution between time steps does not change, there is no time step constraint from the controller."); + } return std::numeric_limits::max(); } } diff --git a/opm/simulators/timestepping/TimeStepControl.hpp b/opm/simulators/timestepping/TimeStepControl.hpp index b21b731817..c4efba75da 100644 --- a/opm/simulators/timestepping/TimeStepControl.hpp +++ b/opm/simulators/timestepping/TimeStepControl.hpp @@ -232,6 +232,7 @@ namespace Opm serializer(alpha_); serializer(controllerVersion_); serializer(toleranceTestVersion_); + serializer(maxReductionTimeStep_); serializer(verbose_); } @@ -246,7 +247,7 @@ namespace Opm mutable std::array timeSteps_{}; mutable std::array beta_{}; mutable std::array alpha_{}; - mutable std::string controllerVersion_ = "I-controller"; + mutable int controllerVersion_ = 0; const std::string toleranceTestVersion_ = ""; const double maxReductionTimeStep_ = 0.1; const bool verbose_ = false; From 108c9649788bea4091fc25b1a24913d3235ad1eb Mon Sep 17 00:00:00 2001 From: Erik Hide Saeternes Date: Mon, 2 Jun 2025 11:08:27 +0200 Subject: [PATCH 5/7] Various minor updates and bug fixes --- .../timestepping/AdaptiveTimeStepping.hpp | 2 +- .../AdaptiveTimeStepping_impl.hpp | 5 +- .../timestepping/TimeStepControl.cpp | 92 ++++++++++--------- .../timestepping/TimeStepControl.hpp | 36 +++++--- .../timestepping/TimeStepControlInterface.hpp | 4 +- 5 files changed, 82 insertions(+), 57 deletions(-) diff --git a/opm/simulators/timestepping/AdaptiveTimeStepping.hpp b/opm/simulators/timestepping/AdaptiveTimeStepping.hpp index d65c7c8abb..726bbf8863 100644 --- a/opm/simulators/timestepping/AdaptiveTimeStepping.hpp +++ b/opm/simulators/timestepping/AdaptiveTimeStepping.hpp @@ -51,7 +51,7 @@ struct TimeStepControlSafetyFactor { static constexpr double value = 0.8; }; struct TimeStepControlRejectCompletedStep { static constexpr bool value = false; }; struct TimeStepControlToleranceTestVersion { static constexpr auto value = "standard"; }; struct TimeStepControlMaxReductionTimeStep { static constexpr double value = 0.1; }; -struct TimeStepControlParameters { static constexpr auto value = "0.125;0.25;0.125;0.375;0.125"; }; +struct TimeStepControlParameters { static constexpr auto value = "0.125;0.25;0.125;0.75;0.25"; }; } // namespace Opm::Parameters diff --git a/opm/simulators/timestepping/AdaptiveTimeStepping_impl.hpp b/opm/simulators/timestepping/AdaptiveTimeStepping_impl.hpp index 40b90c3e47..f06354d054 100644 --- a/opm/simulators/timestepping/AdaptiveTimeStepping_impl.hpp +++ b/opm/simulators/timestepping/AdaptiveTimeStepping_impl.hpp @@ -829,7 +829,10 @@ run() if (substep_report.time_step_rejected) { const double tol = Parameters::Get(); const double safetyFactor = Parameters::Get(); - new_time_step = std::sqrt(safetyFactor * tol / solver_().model().relativeChange()) * dt; + const double temp_time_step = std::sqrt(safetyFactor * tol / solver_().model().relativeChange()) * dt; + if (temp_time_step < dt) { // added in case suggested time step is not a reduction + new_time_step = temp_time_step; + } } checkTimeStepMinLimit_(new_time_step); bool wells_shut = false; diff --git a/opm/simulators/timestepping/TimeStepControl.cpp b/opm/simulators/timestepping/TimeStepControl.cpp index 971244517c..275d6c1283 100644 --- a/opm/simulators/timestepping/TimeStepControl.cpp +++ b/opm/simulators/timestepping/TimeStepControl.cpp @@ -293,28 +293,43 @@ namespace Opm General3rdOrderController::General3rdOrderController( const double tolerance, const double safetyFactor, const bool rejectCompletedStep, - const std::string toleranceTestVersion, + const std::string& toleranceTestVersion, const double maxReductionTimeStep, - const std::string parameters, + const std::string& parameters, const bool verbose) : tolerance_( tolerance ) , safetyFactor_( safetyFactor ) , rejectCompletedStep_( rejectCompletedStep ) - , toleranceTestVersion_( toleranceTestVersion ) , maxReductionTimeStep_( maxReductionTimeStep ) , verbose_( verbose ) { std::istringstream parameters_stream(parameters); std::string value; int counter = 0; - while (std::getline(parameters_stream, value, ';')) { - if (counter < 3) { - beta_[counter] = std::stod(value); + try { + while (std::getline(parameters_stream, value, ';')) + { + if (counter < 3) { + beta_[counter] = std::stod(value); + } + else { + alpha_[counter - 3] = std::stod(value); + } + counter++; } - else { - alpha_[counter - 3] = std::stod(value); - } - counter++; + } + catch (const std::exception& e) { + OPM_THROW(std::runtime_error, "The parameter string for the time step controller is invalid: " + parameters); + } + + if (toleranceTestVersion == "standard") { + toleranceTestVersion_ = ToleranceTestVersions::Standard; + } + else if (toleranceTestVersion == "control-error-filtering") { + toleranceTestVersion_ = ToleranceTestVersions::ControlErrorFiltering; + } + else { + OPM_THROW(std::runtime_error, "Unsupported tolerance test version given: " + toleranceTestVersion); } } @@ -330,10 +345,24 @@ namespace Opm double General3rdOrderController:: computeTimeStepSize(const double dt, const int /* iterations */, const RelativeChangeInterface& /* relChange */, const AdaptiveSimulatorTimer& substepTimer) const { + if (errors_[0] == 0 || errors_[1] == 0 || errors_[2] == 0) + { + if ( verbose_ ) + { + OpmLog::info("The solution between time steps does not change, there is no time step constraint from the controller."); + } + return std::numeric_limits::max(); + } + + for( int i = 0; i < 2; ++i ) + { + assert(std::isfinite(errors_[i])); + } + // Use an I controller after report time steps if (substepTimer.currentStepNum() < 3) { - controllerVersion_ = 0; // "I-controller" + controllerVersion_ = InternalControlVersions::IController; const double newDt = dt * timeStepFactor(errors_, timeSteps_); if (verbose_) { @@ -344,7 +373,7 @@ namespace Opm // Use the general third order controller for all other time steps else { - controllerVersion_ = 1; // "3rdOrder"; + controllerVersion_ = InternalControlVersions::General3rdOrder; const double newDt = dt * timeStepFactor(errors_, timeSteps_); if (verbose_) { @@ -355,10 +384,9 @@ namespace Opm } double General3rdOrderController:: - timeStepFactor(const std::array errors, const std::array timeSteps) const + timeStepFactor(const std::array& errors, const std::array& timeSteps) const { - // "3rdOrder" - if (controllerVersion_ == 1) + if (controllerVersion_ == InternalControlVersions::General3rdOrder) { return std::pow(safetyFactor_ * tolerance_ / errors[2], beta_[0]) * std::pow(safetyFactor_ * tolerance_ / errors[1], beta_[1]) * @@ -366,38 +394,34 @@ namespace Opm std::pow(timeSteps[2] / timeSteps[1], -alpha_[0]) * std::pow(timeSteps[1] / timeSteps[0], -alpha_[1]); } - // "I-controller" + // controllerVersion_ == InternalControlVersions::IController return std::pow(safetyFactor_ * tolerance_ / errors[2], 0.35); } bool General3rdOrderController:: - timeStepAccepted(const double error, const double timeStep) const + timeStepAccepted(const double error, const double timeStepJustCompleted) const { bool acceptTimeStep = true; // Reject time step if chosen tolerance test version fails - if (toleranceTestVersion_ == "standard") + if (toleranceTestVersion_ == ToleranceTestVersions::Standard) { if (rejectCompletedStep_ && error > tolerance_) { acceptTimeStep = false; } } - else if (toleranceTestVersion_ == "control-error-filtering") + else if (toleranceTestVersion_ == ToleranceTestVersions::ControlErrorFiltering) { const std::array tempErrors{errors_[0], errors_[1], error}; - const std::array tempTimeSteps{timeSteps_[0], timeSteps_[1], timeStep}; - double stepFactor = timeStepFactor(tempErrors, tempTimeSteps); + const std::array tempTimeSteps{timeSteps_[0], timeSteps_[1], timeStepJustCompleted}; + const double stepFactor = timeStepFactor(tempErrors, tempTimeSteps); if (rejectCompletedStep_ && stepFactor < maxReductionTimeStep_) { acceptTimeStep = false; } } - else - { - OPM_THROW(std::runtime_error, "Unsupported tolerance test version: " + toleranceTestVersion_); - } - + if (acceptTimeStep) { // Shift errors and time steps @@ -409,21 +433,7 @@ namespace Opm // Store new error and time step errors_[2] = error; - timeSteps_[2] = timeStep; - - for( int i = 0; i < 2; ++i ) - { - assert(std::isfinite(errors_[i])); - } - - if (errors_[0] == 0 || errors_[1] == 0 || errors_[2] == 0.) - { - if ( verbose_ ) - { - OpmLog::info("The solution between time steps does not change, there is no time step constraint from the controller."); - } - return std::numeric_limits::max(); - } + timeSteps_[2] = timeStepJustCompleted; } return acceptTimeStep; diff --git a/opm/simulators/timestepping/TimeStepControl.hpp b/opm/simulators/timestepping/TimeStepControl.hpp index c4efba75da..005e725e77 100644 --- a/opm/simulators/timestepping/TimeStepControl.hpp +++ b/opm/simulators/timestepping/TimeStepControl.hpp @@ -37,6 +37,16 @@ namespace Opm General3rdOrder, }; + enum class ToleranceTestVersions { + Standard, + ControlErrorFiltering, + }; + + enum class InternalControlVersions { + IController, + General3rdOrder, + }; + /////////////////////////////////////////////////////////////////////////////////////////////////////////////// /// /// A simple iteration count based adaptive time step control. @@ -66,7 +76,7 @@ namespace Opm const RelativeChangeInterface& /* relativeChange */, const AdaptiveSimulatorTimer& /* substepTimer */ ) const override; - bool timeStepAccepted(const double /* error */, const double /* timeStep */) const override { return true; } + bool timeStepAccepted(const double /* error */, const double /* timeStepJustCompleted */) const override { return true; } template void serializeOp(Serializer& serializer) @@ -119,7 +129,7 @@ namespace Opm const RelativeChangeInterface& relativeChange, const AdaptiveSimulatorTimer& /* substepTimer */ ) const override; - bool timeStepAccepted(const double /* error */, const double /* timeStep */) const override { return true; } + bool timeStepAccepted(const double /* error */, const double /* timeStepJustCompleted */) const override { return true; } template void serializeOp(Serializer& serializer) @@ -170,7 +180,7 @@ namespace Opm const RelativeChangeInterface& relativeChange, const AdaptiveSimulatorTimer& /* substepTimer */ ) const override; - bool timeStepAccepted(const double /* error */, const double /* timeStep */) const override { return true; } + bool timeStepAccepted(const double /* error */, const double /* timeStepJustCompleted */) const override { return true; } template void serializeOp(Serializer& serializer) @@ -204,9 +214,9 @@ namespace Opm General3rdOrderController( const double tolerance = 1e-3, const double safetyFactor = 0.8, const bool rejectCompletedStep = false, - const std::string toleranceTestVersion = "", + const std::string& toleranceTestVersion = "", const double maxReductionTimeStep = 0.1, - const std::string parameters = "", + const std::string& parameters = "", const bool verbose = false); static General3rdOrderController serializationTestObject(); @@ -216,9 +226,9 @@ namespace Opm const RelativeChangeInterface& /* relativeChange */, const AdaptiveSimulatorTimer& substepTimer) const override; - double timeStepFactor(const std::array errors, const std::array timeSteps) const; + double timeStepFactor(const std::array& errors, const std::array& timeSteps) const; - bool timeStepAccepted(const double error, const double timeStep) const override; + bool timeStepAccepted(const double error, const double timeStepJustCompleted) const override; template void serializeOp(Serializer& serializer) @@ -240,15 +250,15 @@ namespace Opm protected: - const double tolerance_ = 1e-3; + const double tolerance_ = 0.1; const double safetyFactor_ = 0.8; const bool rejectCompletedStep_ = false; mutable std::array errors_{}; mutable std::array timeSteps_{}; - mutable std::array beta_{}; - mutable std::array alpha_{}; - mutable int controllerVersion_ = 0; - const std::string toleranceTestVersion_ = ""; + mutable std::array beta_{0.125, 0.25, 0.125}; + mutable std::array alpha_{0.75, 0.25}; + mutable InternalControlVersions controllerVersion_{InternalControlVersions::IController}; + ToleranceTestVersions toleranceTestVersion_{ToleranceTestVersions::Standard}; const double maxReductionTimeStep_ = 0.1; const bool verbose_ = false; }; @@ -280,7 +290,7 @@ namespace Opm const RelativeChangeInterface& /*relativeChange */, const AdaptiveSimulatorTimer& substepTimer) const override; - bool timeStepAccepted(const double /* error */, const double /* timeStep */) const override { return true; } + bool timeStepAccepted(const double /* error */, const double /* timeStepJustCompleted */) const override { return true; } template void serializeOp(Serializer& serializer) diff --git a/opm/simulators/timestepping/TimeStepControlInterface.hpp b/opm/simulators/timestepping/TimeStepControlInterface.hpp index 6e34443d7d..813b9af469 100644 --- a/opm/simulators/timestepping/TimeStepControlInterface.hpp +++ b/opm/simulators/timestepping/TimeStepControlInterface.hpp @@ -60,7 +60,9 @@ namespace Opm /// \return suggested time step size for the next step virtual double computeTimeStepSize( const double dt, const int iterations, const RelativeChangeInterface& relativeChange , const AdaptiveSimulatorTimer& substepTimer) const = 0; - virtual bool timeStepAccepted(const double error, const double timeStep) const = 0; + /// For the general 3rd order controller, the internal shifting of errors and time steps happens here, and + /// hence this method needs to be called for (after) each time step + virtual bool timeStepAccepted(const double error, const double timeStepJustCompleted) const = 0; /// virtual destructor (empty) virtual ~TimeStepControlInterface () {} From 752fba50f7d114e4e7e2385aff86f7c736fa52de Mon Sep 17 00:00:00 2001 From: Erik Hide Saeternes Date: Tue, 3 Jun 2025 13:55:43 +0200 Subject: [PATCH 6/7] Removed default parameters from constructor --- .../timestepping/AdaptiveTimeStepping.cpp | 14 ++++++++------ opm/simulators/timestepping/TimeStepControl.hpp | 16 +++++++++------- 2 files changed, 17 insertions(+), 13 deletions(-) diff --git a/opm/simulators/timestepping/AdaptiveTimeStepping.cpp b/opm/simulators/timestepping/AdaptiveTimeStepping.cpp index f919763f57..a6150bf1b0 100644 --- a/opm/simulators/timestepping/AdaptiveTimeStepping.cpp +++ b/opm/simulators/timestepping/AdaptiveTimeStepping.cpp @@ -233,11 +233,12 @@ createController(const UnitSystem& unitSystem) }}, {"general3rdorder", [tol]() { - const double safetyFactor = Parameters::Get(); // 0.8 - const bool rejectCompletedStep = Parameters::Get(); // false - const std::string toleranceTestVersion = Parameters::Get(); // "standard" - const double maxReductionTimeStep = Parameters::Get(); // 0.1 - const std::string parameters = Parameters::Get(); // "" + const double safetyFactor = Parameters::Get(); + const bool rejectCompletedStep = Parameters::Get(); + const std::string toleranceTestVersion = Parameters::Get(); + const double maxReductionTimeStep = Parameters::Get(); + const std::string parameters = Parameters::Get(); + const bool verbose = Parameters::Get(); return RetVal{ TimeStepControlType::General3rdOrder, std::make_unique(tol, @@ -245,7 +246,8 @@ createController(const UnitSystem& unitSystem) rejectCompletedStep, toleranceTestVersion, maxReductionTimeStep, - parameters), + parameters, + verbose), false }; }}, diff --git a/opm/simulators/timestepping/TimeStepControl.hpp b/opm/simulators/timestepping/TimeStepControl.hpp index 005e725e77..e7e5da5a55 100644 --- a/opm/simulators/timestepping/TimeStepControl.hpp +++ b/opm/simulators/timestepping/TimeStepControl.hpp @@ -211,13 +211,15 @@ namespace Opm public: static constexpr TimeStepControlType Type = TimeStepControlType::General3rdOrder; - General3rdOrderController( const double tolerance = 1e-3, - const double safetyFactor = 0.8, - const bool rejectCompletedStep = false, - const std::string& toleranceTestVersion = "", - const double maxReductionTimeStep = 0.1, - const std::string& parameters = "", - const bool verbose = false); + General3rdOrderController() = default; + + General3rdOrderController( const double tolerance, + const double safetyFactor, + const bool rejectCompletedStep, + const std::string& toleranceTestVersion, + const double maxReductionTimeStep, + const std::string& parameters, + const bool verbose); static General3rdOrderController serializationTestObject(); From 6cf2d41c99f7ee3e3831f25ab70033851f967118 Mon Sep 17 00:00:00 2001 From: Erik Hide Saeternes Date: Tue, 3 Jun 2025 14:49:39 +0200 Subject: [PATCH 7/7] Fixed tiny bug which caused test failure --- opm/simulators/timestepping/TimeStepControl.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/opm/simulators/timestepping/TimeStepControl.cpp b/opm/simulators/timestepping/TimeStepControl.cpp index 275d6c1283..8093db0268 100644 --- a/opm/simulators/timestepping/TimeStepControl.cpp +++ b/opm/simulators/timestepping/TimeStepControl.cpp @@ -336,7 +336,7 @@ namespace Opm General3rdOrderController General3rdOrderController::serializationTestObject() { - General3rdOrderController result(1.0, 2.0, false, "", 3.0, "", false); + General3rdOrderController result(1.0, 2.0, false, "standard", 3.0, "0.125;0.25;0.125;0.75;0.25", false); result.errors_ = {2.0, 3.0}; return result;