diff --git a/opm/simulators/utils/PartiallySupportedFlowKeywords.cpp b/opm/simulators/utils/PartiallySupportedFlowKeywords.cpp index dd49ca63c77..39638eaa5aa 100644 --- a/opm/simulators/utils/PartiallySupportedFlowKeywords.cpp +++ b/opm/simulators/utils/PartiallySupportedFlowKeywords.cpp @@ -107,7 +107,7 @@ partiallySupported() { "GECON", { - {7,{true, allow_values {"NONE"}, "GECON(WORKOVER): Workover procedures not implemented"}}, + {7,{true, allow_values {"NONE", "WELL"}, "GECON(WORKOVER): Workover procedures only implemented for WELL and NONE"}}, {8,{true, allow_values {"NO"}, "GECON(ENDRUN): End run not implemented"}}, }, }, diff --git a/opm/simulators/wells/BlackoilWellModelGeneric.cpp b/opm/simulators/wells/BlackoilWellModelGeneric.cpp index 2cab4703dd4..784c56cfb30 100644 --- a/opm/simulators/wells/BlackoilWellModelGeneric.cpp +++ b/opm/simulators/wells/BlackoilWellModelGeneric.cpp @@ -458,8 +458,11 @@ checkGEconLimits( if (checker.minOilRate() || checker.minGasRate()) { checker.closeWells(); } - else if (checker.waterCut() || checker.GOR() || checker.WGR()) { - checker.doWorkOver(); + else { + const auto ratio_details = checker.ratioViolation(); + if (ratio_details.has_value()) { + checker.doWorkOver(ratio_details.value()); + } } if (checker.endRun() && (checker.numProducersOpenInitially() >= 1) && (checker.numProducersOpen() == 0)) diff --git a/opm/simulators/wells/GroupEconomicLimitsChecker.cpp b/opm/simulators/wells/GroupEconomicLimitsChecker.cpp index 5a2bc8bcb81..681938999a1 100644 --- a/opm/simulators/wells/GroupEconomicLimitsChecker.cpp +++ b/opm/simulators/wells/GroupEconomicLimitsChecker.cpp @@ -24,7 +24,9 @@ #include #include +#include #include +#include #include @@ -32,13 +34,17 @@ #include #include +#include +#include #include -#include #include -#include +#include #include +#include +#include +#include namespace Opm { @@ -106,10 +112,19 @@ closeWells() template void GroupEconomicLimitsChecker:: -doWorkOver() +doWorkOver(const RatioDetails& ratio_details) { - if (this->gecon_props_.workover() != GroupEconProductionLimits::EconWorkover::NONE) { - throwNotImplementedError("workover procedure"); + switch (this->gecon_props_.workover()) { + case GroupEconProductionLimits::EconWorkover::NONE: + break; + case GroupEconProductionLimits::EconWorkover::WELL: + this->closeWorstOffendingRatioWell(ratio_details); + break; + default: + const auto workover_name = GroupEconProductionLimits::econWorkoverToString(this->gecon_props_.workover()); + throwNotImplementedError(fmt::format( + "workover procedure {}", + workover_name)); } } @@ -124,37 +139,25 @@ endRun() } template -bool GroupEconomicLimitsChecker:: +std::optional::RatioDetails> +GroupEconomicLimitsChecker:: GOR() { - auto oil_phase_idx = this->phase_idx_reverse_map_[IndexTraits::oilPhaseIdx]; - auto gas_phase_idx = this->phase_idx_reverse_map_[IndexTraits::gasPhaseIdx]; - auto oil_rate = this->production_rates_[oil_phase_idx]; - auto gas_rate = this->production_rates_[gas_phase_idx]; - Scalar gor; - if (gas_rate <= 0.0) { - gor = 0.0; - } - else if (oil_rate <= 0.0) { - gor = 1e100; - } - else { - gor = gas_rate / oil_rate; - } - if (auto max_gor = this->gecon_props_.maxGasOilRatio(); max_gor) { - if (gor > *max_gor) { + const auto details = this->groupRatioDetails(RatioViolation::GOR); + if (details.has_value()) { + if (details->ratio > details->limit) { if (this->debug_) { const std::string msg = fmt::format( "GOR={} is greater than maximum: {}", - gor, *max_gor); + details->ratio, details->limit); displayDebugMessage(msg); } - addPrintMessage(" Gas/oil ratio = {:.2f} {} which is greater than the minimum economic value = {:.2f} {}", - gor, *max_gor, UnitSystem::measure::gas_oil_ratio); - return true; + addPrintMessage(" Gas/oil ratio = {:.2f} {} which is greater than the maximum economic value = {:.2f} {}", + details->ratio, details->limit, details->measure); + return details; } } - return false; + return std::nullopt; } template @@ -225,78 +228,69 @@ numProducersOpenInitially() return 1; } +template +std::optional::RatioDetails> +GroupEconomicLimitsChecker:: +ratioViolation() +{ + if (auto water_cut_details = this->waterCut(); water_cut_details.has_value()) + { + return water_cut_details; + } + else if (auto gor_details = this->GOR(); gor_details.has_value()) + { + return gor_details; + } + else if (auto wgr_details = this->WGR(); wgr_details.has_value()) + { + return wgr_details; + } + + return std::nullopt; +} + template -bool GroupEconomicLimitsChecker:: +std::optional::RatioDetails> +GroupEconomicLimitsChecker:: waterCut() { - auto oil_phase_idx = this->phase_idx_reverse_map_[IndexTraits::oilPhaseIdx]; - auto water_phase_idx = this->phase_idx_reverse_map_[IndexTraits::waterPhaseIdx]; - auto oil_rate = this->production_rates_[oil_phase_idx]; - auto water_rate = this->production_rates_[water_phase_idx]; - auto liquid_rate = oil_rate + water_rate; - Scalar water_cut; - if (liquid_rate == 0.0) { - water_cut = 0.0; - } - else { - if (water_rate < 0.0) { - water_cut = 0.0; - } - else if (oil_rate < 0.0) { - water_cut = 1.0; - } - else { - water_cut = water_rate / liquid_rate; - } - } - if (auto max_water_cut = this->gecon_props_.maxWaterCut(); max_water_cut) { - if (water_cut > *max_water_cut) { + const auto details = this->groupRatioDetails(RatioViolation::WATER_CUT); + if (details.has_value()) { + if (details->ratio > details->limit) { if (this->debug_) { const std::string msg = fmt::format( "water_cut={} is greater than maximum: {}", - water_cut, *max_water_cut); + details->ratio, details->limit); displayDebugMessage(msg); } addPrintMessage(" Water cut = {:.2f} {} which is greater than the maximum economic value = {:.2f} {}", - water_cut, *max_water_cut, UnitSystem::measure::water_cut); - return true; + details->ratio, details->limit, details->measure); + return details; } } - return false; + return std::nullopt; } template -bool GroupEconomicLimitsChecker:: +std::optional::RatioDetails> +GroupEconomicLimitsChecker:: WGR() { - auto water_phase_idx = this->phase_idx_reverse_map_[IndexTraits::waterPhaseIdx]; - auto gas_phase_idx = this->phase_idx_reverse_map_[IndexTraits::gasPhaseIdx]; - auto water_rate = this->production_rates_[water_phase_idx]; - auto gas_rate = this->production_rates_[gas_phase_idx]; - Scalar wgr; - if (water_rate <= 0.0) { - wgr = 0.0; - } - else if (gas_rate <= 0.0) { - wgr = 1e100; - } - else { - wgr = water_rate / gas_rate; - } - if (auto max_wgr = this->gecon_props_.maxWaterGasRatio(); max_wgr) { - if (wgr > *max_wgr) { + const auto details = this->groupRatioDetails(RatioViolation::WGR); + if (details.has_value()) { + if (details->ratio > details->limit) { if (this->debug_) { const std::string msg = fmt::format( "WGR={} is greater than maximum: {}", - wgr, *max_wgr); + details->ratio, details->limit); displayDebugMessage(msg); } addPrintMessage(" Water/gas ratio = {:.2f} {} which is greater than the maximum economic value = {:.2f} {}", - wgr, *max_wgr, UnitSystem::measure::gas_oil_ratio); // Same units - return true; + details->ratio, details->limit, details->measure); // Same units + return details; } } - return false; + return std::nullopt; } /**************************************** @@ -407,6 +401,258 @@ throwNotImplementedError(const std::string& error) const OPM_DEFLOG_THROW(std::runtime_error, msg, this->deferred_logger_); } +template +void GroupEconomicLimitsChecker:: +collectProducerWells(const Group& group, std::vector& well_names) const +{ + for (const std::string& group_name : group.groups()) { + const auto& sub_group = this->schedule_.getGroup(group_name, this->report_step_idx_); + this->collectProducerWells(sub_group, well_names); + } + for (const std::string& well_name : group.wells()) { + const auto& well_ecl = this->schedule_.getWell(well_name, this->report_step_idx_); + if (well_ecl.isProducer()) { + well_names.push_back(well_name); + } + } +} + +template +Scalar GroupEconomicLimitsChecker:: +computeWellRatio(const std::string& well_name, + const RatioViolation ratio_violation) const +{ + constexpr Scalar invalid_ratio = std::numeric_limits::lowest(); + + // The checks below avoid accessing WellState data for wells that are already closed, + // not present on this rank, or not owned by this rank (important for parallel runs). + if (this->well_test_state_.well_is_closed(well_name)) { + return invalid_ratio; + } + + const auto well_index = this->well_state_.index(well_name); + if (!well_index.has_value()) { + return invalid_ratio; + } + if (!this->well_state_.wellIsOwned(well_index.value(), well_name)) { + return invalid_ratio; + } + + const auto& ws = this->well_state_.well(well_index.value()); + if (ws.status == Well::Status::SHUT) { + return invalid_ratio; + } + + const auto& pu = this->well_model_.phaseUsage(); + constexpr Scalar big_value = std::numeric_limits::max(); + + switch (ratio_violation) { + case RatioViolation::WATER_CUT: { + if (!pu.phaseIsActive(IndexTraits::oilPhaseIdx) || + !pu.phaseIsActive(IndexTraits::waterPhaseIdx)) + { + return invalid_ratio; + } + const auto oil_pos = pu.canonicalToActivePhaseIdx(IndexTraits::oilPhaseIdx); + const auto water_pos = pu.canonicalToActivePhaseIdx(IndexTraits::waterPhaseIdx); + // Surface rates are negative for producers in WellState; flip sign here. + const Scalar oil_rate = -ws.surface_rates[oil_pos]; + const Scalar water_rate = -ws.surface_rates[water_pos]; + const Scalar liquid = oil_rate + water_rate; + if (liquid <= 0.0) return Scalar{0}; + if (water_rate <= 0.0) return Scalar{0}; + if (oil_rate <= 0.0) return Scalar{1}; + return water_rate / liquid; + } + case RatioViolation::GOR: { + if (!pu.phaseIsActive(IndexTraits::oilPhaseIdx) || + !pu.phaseIsActive(IndexTraits::gasPhaseIdx)) + { + return invalid_ratio; + } + const auto oil_pos = pu.canonicalToActivePhaseIdx(IndexTraits::oilPhaseIdx); + const auto gas_pos = pu.canonicalToActivePhaseIdx(IndexTraits::gasPhaseIdx); + // Surface rates are negative for producers in WellState; flip sign here. + const Scalar oil_rate = -ws.surface_rates[oil_pos]; + const Scalar gas_rate = -ws.surface_rates[gas_pos]; + if (gas_rate <= 0.0) return Scalar{0}; + if (oil_rate <= 0.0) return big_value; + return gas_rate / oil_rate; + } + case RatioViolation::WGR: { + if (!pu.phaseIsActive(IndexTraits::gasPhaseIdx) || + !pu.phaseIsActive(IndexTraits::waterPhaseIdx)) + { + return invalid_ratio; + } + const auto gas_pos = pu.canonicalToActivePhaseIdx(IndexTraits::gasPhaseIdx); + const auto water_pos = pu.canonicalToActivePhaseIdx(IndexTraits::waterPhaseIdx); + // Surface rates are negative for producers in WellState; flip sign here. + const Scalar gas_rate = -ws.surface_rates[gas_pos]; + const Scalar water_rate = -ws.surface_rates[water_pos]; + if (water_rate <= 0.0) return Scalar{0}; + if (gas_rate <= 0.0) return big_value; + return water_rate / gas_rate; + } + case RatioViolation::NONE: + default: + return invalid_ratio; + } +} + +template +std::optional::RatioDetails> +GroupEconomicLimitsChecker:: +groupRatioDetails(const RatioViolation ratio_violation) const +{ + constexpr Scalar big_value = std::numeric_limits::max(); + switch (ratio_violation) { + case RatioViolation::WGR: { + const auto max_wgr = this->gecon_props_.maxWaterGasRatio(); + if (!max_wgr.has_value()) { + return std::nullopt; + } + const auto water_phase_idx = this->phase_idx_reverse_map_.at(IndexTraits::waterPhaseIdx); + const auto gas_phase_idx = this->phase_idx_reverse_map_.at(IndexTraits::gasPhaseIdx); + const Scalar water_rate = this->production_rates_[water_phase_idx]; + const Scalar gas_rate = this->production_rates_[gas_phase_idx]; + const Scalar ratio = (water_rate <= 0.0) ? Scalar{0} + : (gas_rate <= 0.0) ? Scalar{big_value} + : water_rate / gas_rate; + return RatioDetails{RatioViolation::WGR, "water-gas ratio", ratio, static_cast(max_wgr.value()), UnitSystem::measure::gas_oil_ratio}; + } + case RatioViolation::GOR: { + const auto max_gor = this->gecon_props_.maxGasOilRatio(); + if (!max_gor.has_value()) { + return std::nullopt; + } + const auto oil_phase_idx = this->phase_idx_reverse_map_.at(IndexTraits::oilPhaseIdx); + const auto gas_phase_idx = this->phase_idx_reverse_map_.at(IndexTraits::gasPhaseIdx); + const Scalar oil_rate = this->production_rates_[oil_phase_idx]; + const Scalar gas_rate = this->production_rates_[gas_phase_idx]; + const Scalar ratio = (gas_rate <= 0.0) ? Scalar{0} + : (oil_rate <= 0.0) ? Scalar{big_value} + : gas_rate / oil_rate; + return RatioDetails{RatioViolation::GOR, "gas-oil ratio", ratio, static_cast(max_gor.value()), UnitSystem::measure::gas_oil_ratio}; + } + case RatioViolation::WATER_CUT: { + const auto max_water_cut = this->gecon_props_.maxWaterCut(); + if (!max_water_cut.has_value()) { + return std::nullopt; + } + const auto oil_phase_idx = this->phase_idx_reverse_map_.at(IndexTraits::oilPhaseIdx); + const auto water_phase_idx = this->phase_idx_reverse_map_.at(IndexTraits::waterPhaseIdx); + const Scalar oil_rate = this->production_rates_[oil_phase_idx]; + const Scalar water_rate = this->production_rates_[water_phase_idx]; + const Scalar liquid_rate = oil_rate + water_rate; + const Scalar ratio = (liquid_rate == 0.0) ? Scalar{0} + : (water_rate < 0.0) ? Scalar{0} + : (oil_rate < 0.0) ? Scalar{1} + : water_rate / liquid_rate; + return RatioDetails{RatioViolation::WATER_CUT, "water cut", ratio, static_cast(max_water_cut.value()), UnitSystem::measure::water_cut}; + } + case RatioViolation::NONE: + default: + return std::nullopt; + } +} + +template +void GroupEconomicLimitsChecker:: +closeWorstOffendingRatioWell(const RatioDetails& ratio_details) +{ + if (ratio_details.violation == RatioViolation::NONE) { + // Should not happen: doWorkOver() is only called after a ratio + // violation was reported, but guard against misuse. + return; + } + + std::vector producer_wells; + this->collectProducerWells(this->group_, producer_wells); + + // Determine the worst-offending well consistently across all ranks. + // Each rank reports its locally-owned ratio (or "lowest" for non-owned). + // Reduce all ratios in one bulk MPI max-reduction, then select the + // globally worst ratio locally from the reduced array. + std::string worst_well; + Scalar worst_ratio = std::numeric_limits::lowest(); + std::vector global_ratios; + global_ratios.reserve(producer_wells.size()); + for (const std::string& well_name : producer_wells) { + global_ratios.push_back(this->computeWellRatio(well_name, ratio_details.violation)); + } + + if (!global_ratios.empty()) { + this->well_model_.comm().max(global_ratios.data(), global_ratios.size()); + for (std::size_t i = 0; i < global_ratios.size(); ++i) { + if (global_ratios[i] > worst_ratio) { + worst_ratio = global_ratios[i]; + worst_well = producer_wells[i]; + } + } + } + + if (worst_well.empty() || + worst_ratio == std::numeric_limits::lowest()) + { + // No candidate producer found (e.g. all already closed). + if (this->debug_) { + displayDebugMessage( + "GECON WELL workover: no open producer well found to close."); + } + return; + } + + + // Perform the actual close on the owning rank(s). The well's automatic + // shut-in instruction (WELSPECS item 9) selects between SHUT and STOP for + // the log message; the underlying WellTestState::close_well call is the + // same in either case (the well status is updated subsequently). + const auto& well_ecl = + this->schedule_.getWell(worst_well, this->report_step_idx_); + const std::string action = well_ecl.getAutomaticShutIn() ? "shut" : "stopped"; + + bool well_was_closed_locally = false; + if (this->well_model_.hasLocalWell(worst_well)) { + if (!this->well_test_state_.well_is_closed(worst_well)) { + this->well_test_state_.close_well( + worst_well, WellTestConfig::Reason::GROUP, this->simulation_time_); + this->well_model_.updateClosedWellsThisStep(worst_well); + well_was_closed_locally = true; + } + } + + if (this->well_model_.comm().max(static_cast(well_was_closed_locally)) == 0) { + // Nothing to do (well was already closed everywhere). + return; + } + + if (this->well_model_.comm().rank() == 0) { + const Scalar well_ratio_display = this->unit_system_.from_si(ratio_details.measure, worst_ratio); + const Scalar group_ratio_display = this->unit_system_.from_si(ratio_details.measure, ratio_details.ratio); + const Scalar limit_display = this->unit_system_.from_si(ratio_details.measure, ratio_details.limit); + const std::string unit_name = this->unit_system_.name(ratio_details.measure); + + const std::string msg = fmt::format( + "{}\nwell {} ({} {:.4e} {}) is {} at time {:.2f} {} (date = {}),\n" + "because the group {} {} {:.4e} {}" + " exceeds the limit {:.4e} {}.\n{}", + this->message_separator(), + worst_well, + ratio_details.ratio_type, well_ratio_display, unit_name, + action, + this->unit_system_.from_si(UnitSystem::measure::time, this->simulation_time_), + this->unit_system_.name(UnitSystem::measure::time), + this->date_string_, + this->group_.name(), + ratio_details.ratio_type, + group_ratio_display, unit_name, + limit_display, unit_name, + this->message_separator()); + this->deferred_logger_.info(msg); + } +} + template class GroupEconomicLimitsChecker; #if FLOW_INSTANTIATE_FLOAT diff --git a/opm/simulators/wells/GroupEconomicLimitsChecker.hpp b/opm/simulators/wells/GroupEconomicLimitsChecker.hpp index 3c855803c57..4c0a9d4351e 100644 --- a/opm/simulators/wells/GroupEconomicLimitsChecker.hpp +++ b/opm/simulators/wells/GroupEconomicLimitsChecker.hpp @@ -25,7 +25,9 @@ #include #include +#include #include +#include namespace Opm { @@ -40,6 +42,22 @@ template class GroupEconomicLimitsChecker { public: + enum class RatioViolation { + NONE, + WATER_CUT, + GOR, + WGR + }; + + struct RatioDetails + { + RatioViolation violation; + std::string ratio_type; + Scalar ratio; + Scalar limit; + UnitSystem::measure measure; + }; + GroupEconomicLimitsChecker(const BlackoilWellModelGeneric& well_model, WellTestState& well_test_state, const Group& group, @@ -49,10 +67,11 @@ class GroupEconomicLimitsChecker void closeWells(); bool minGasRate(); bool minOilRate(); - bool waterCut(); - bool GOR(); - bool WGR(); - void doWorkOver(); + std::optional waterCut(); + std::optional GOR(); + std::optional WGR(); + std::optional ratioViolation(); + void doWorkOver(const RatioDetails& ratio_details); bool endRun(); int numProducersOpenInitially(); int numProducersOpen(); @@ -72,6 +91,24 @@ class GroupEconomicLimitsChecker bool closeWellsRecursive(const Group& group, int level = 0); void throwNotImplementedError(const std::string& error) const; + /// Collect names of all producer wells belonging to \p group (recursively). + void collectProducerWells(const Group& group, + std::vector& well_names) const; + + /// Compute the relevant ratio (per \p ratio_violation) for a single well. + /// Returns \c std::numeric_limits::lowest() if the well is not owned + /// by the current rank or is already shut, so that an MPI max-reduction picks + /// up the value from the owning rank. + Scalar computeWellRatio(const std::string& well_name, + const RatioViolation ratio_violation) const; + + /// Implements the WELL workover procedure: identify the worst-offending + /// producer in the group hierarchy and close it (shut/stop based on + /// WELSPECS automatic shut-in policy). + void closeWorstOffendingRatioWell(const RatioDetails& ratio_details); + + std::optional groupRatioDetails(const RatioViolation ratio_violation) const; + const BlackoilWellModelGeneric& well_model_; const Group& group_; const double simulation_time_; diff --git a/regressionTests.cmake b/regressionTests.cmake index dcc997edcce..aa930ea7977 100644 --- a/regressionTests.cmake +++ b/regressionTests.cmake @@ -2667,6 +2667,23 @@ add_test_compareECLFiles( satellite ) +add_test_compareECLFiles( + CASENAME + gecon_ratio_well + FILENAME + GECON-02 + SIMULATOR + flow + DEV_SIMULATOR + flow_blackoil + ABS_TOL + ${abs_tol} + REL_TOL + ${rel_tol} + DIR + gecon +) + if(BUILD_FLOW_POLY_GRID) add_test_compareECLFiles( CASENAME