From 899ab9b13dcf015e32f73f953793ace75673f256 Mon Sep 17 00:00:00 2001 From: Kai Bao Date: Tue, 31 Jan 2023 10:54:33 +0100 Subject: [PATCH] fix the running when solution does not change between time steps The PID time step control will not pose any constraint on the time step size by returning a largest possible time step size. --- .../timestepping/TimeStepControl.cpp | 19 ++++++++++++++----- 1 file changed, 14 insertions(+), 5 deletions(-) diff --git a/opm/simulators/timestepping/TimeStepControl.cpp b/opm/simulators/timestepping/TimeStepControl.cpp index 1e7ad40bc..09a494b4d 100644 --- a/opm/simulators/timestepping/TimeStepControl.cpp +++ b/opm/simulators/timestepping/TimeStepControl.cpp @@ -27,11 +27,15 @@ #include #include #include +#include #include +#include #include #include +#include + namespace Opm { //////////////////////////////////////////////////////// @@ -141,17 +145,22 @@ namespace Opm errors_[ 2 ] = error; for( int i=0; i<2; ++i ) { assert(std::isfinite(errors_[i])); - assert(errors_[i]>0); } - if( error > tol_ ) + if( errors_[2] > tol_ ) { // adjust dt by given tolerance const double newDt = dt * tol_ / error; - if( verbose_ ) - std::cout << "Computed step size (tol): " << unit::convert::to( newDt, unit::day ) << " (days)" << std::endl; + if ( verbose_ ) + OpmLog::info(fmt::format("Computed step size (tol): {} days", unit::convert::to( newDt, unit::day ))); return newDt; } + else if (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 PID time step control "); + return std::numeric_limits::max(); + } else { // values taking from turek time stepping paper @@ -162,7 +171,7 @@ namespace Opm std::pow( tol_ / errors_[ 2 ], kI ) * std::pow( errors_[0]*errors_[0]/errors_[ 1 ]/errors_[ 2 ], kD )); if( verbose_ ) - std::cout << "Computed step size (pow): " << unit::convert::to( newDt, unit::day ) << " (days)" << std::endl; + OpmLog::info(fmt::format("Computed step size (pow): {} days", unit::convert::to( newDt, unit::day ))); return newDt; } }