From 2e81e67e737bb0dca923bb061a8fdf97f059c026 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Atgeirr=20Fl=C3=B8=20Rasmussen?= Date: Fri, 28 Mar 2025 16:56:26 +0100 Subject: [PATCH] Updates to clean up. --- .../wells/BlackoilWellModel_impl.hpp | 3 +- .../MultisegmentWellPrimaryVariables.cpp | 35 ++++++----- opm/simulators/wells/SingleWellState.cpp | 14 +++-- opm/simulators/wells/SingleWellState.hpp | 17 ++--- .../wells/StandardWellPrimaryVariables.cpp | 33 +++++----- opm/simulators/wells/WellInterface.hpp | 4 -- opm/simulators/wells/WellInterface_impl.hpp | 63 ++----------------- 7 files changed, 57 insertions(+), 112 deletions(-) diff --git a/opm/simulators/wells/BlackoilWellModel_impl.hpp b/opm/simulators/wells/BlackoilWellModel_impl.hpp index d83cf6f49a..4a04c7a201 100644 --- a/opm/simulators/wells/BlackoilWellModel_impl.hpp +++ b/opm/simulators/wells/BlackoilWellModel_impl.hpp @@ -436,8 +436,7 @@ namespace Opm { // This is done only for producers, as injectors will only have a single // nonzero phase anyway. for (const auto& well : well_container_) { - //const bool zero_target = well->stoppedOrZeroRateTarget(simulator_, this->wellState(), local_deferredLogger); - if (well->isProducer()){// && !zero_target) { + if (well->isProducer()) { well->initializeWellState(simulator_, this->groupState(), summaryState, this->wellState(), local_deferredLogger); } } diff --git a/opm/simulators/wells/MultisegmentWellPrimaryVariables.cpp b/opm/simulators/wells/MultisegmentWellPrimaryVariables.cpp index 7a95993346..74bdb4458a 100644 --- a/opm/simulators/wells/MultisegmentWellPrimaryVariables.cpp +++ b/opm/simulators/wells/MultisegmentWellPrimaryVariables.cpp @@ -111,29 +111,29 @@ update(const WellState& well_state, if (stop_or_zero_rate_target && seg == 0) { value_[seg][WQTotal] = 0; } - assert(ws.initializedFromReservoir()); - // tot to map old fraction to new perforations for now start from scratch. - bool prim_set = ws.multiseg_primaryvar.size() == value_.size(); - if(prim_set){ - //if (std::abs(total_seg_rate) > 0.) { + + // Production wells should have been initialized using reservoir data. + assert(well_.isProducer() == ws.initialized_from_reservoir); + + // todo to map old fraction to new perforations for now start from scratch. + if (ws.primaryvar.size() == value_.size() * numWellEq) { if (has_wfrac_variable) { - value_[seg][WFrac] = ws.multiseg_primaryvar[seg][WFrac]; + value_[seg][WFrac] = ws.primaryvar[seg * numWellEq + WFrac]; } if (has_gfrac_variable) { - value_[seg][GFrac] = ws.multiseg_primaryvar[seg][GFrac]; + value_[seg][GFrac] = ws.primaryvar[seg * numWellEq + GFrac]; } - } else { - if (std::abs(total_seg_rate) > 0.) { + } else if (std::abs(total_seg_rate) > 0.) { if (has_wfrac_variable) { const int water_pos = pu.phase_pos[Water]; value_[seg][WFrac] = well_.scalingFactor(water_pos) * segment_rates[well_.numPhases() * seg + water_pos] / total_seg_rate; } if (has_gfrac_variable) { const int gas_pos = pu.phase_pos[Gas]; - value_[seg][GFrac] = well_.scalingFactor(gas_pos) * segment_rates[well_.numPhases() * seg + gas_pos] / total_seg_rate; + value_[seg][GFrac] = well_.scalingFactor(gas_pos) * segment_rates[well_.numPhases() * seg + gas_pos] / total_seg_rate; } // what about water and gas injection? - } else { // total_seg_rate == 0 + } else { // total_seg_rate == 0 if (well_.isInjector()) { // only single phase injection handled auto phase = well.getInjectionProperties().injectorType; @@ -164,7 +164,6 @@ update(const WellState& well_state, } } } - } } setEvaluationsFromValues(); } @@ -237,16 +236,22 @@ copyToWellState(const MultisegmentWellGeneric& mswell, static constexpr int Gas = BlackoilPhases::Vapour; static constexpr int Oil = BlackoilPhases::Liquid; static constexpr int Water = BlackoilPhases::Aqua; - + const auto pvtReg = std::max(well_.wellEcl().pvt_table_number() - 1, 0); const PhaseUsage& pu = well_.phaseUsage(); const int oil_pos = pu.phase_pos[Oil]; auto& ws = well_state.well(well_.indexOfWell()); - if constexpr ( numWellEq == 4) { - ws.multiseg_primaryvar = value_; + + // Store primary variables + ws.primaryvar.resize(value_.size() * numWellEq); + for (std::size_t seg = 0; seg < value_.size(); ++seg) { + for (int ii = 0; ii < numWellEq; ++ii) { + ws.primaryvar[seg * numWellEq + ii] = value_[seg][ii]; + } } + auto& segments = ws.segments; auto& segment_rates = segments.rates; auto& disgas = segments.dissolved_gas_rate; diff --git a/opm/simulators/wells/SingleWellState.cpp b/opm/simulators/wells/SingleWellState.cpp index 00521c5799..2b1ecb76c8 100644 --- a/opm/simulators/wells/SingleWellState.cpp +++ b/opm/simulators/wells/SingleWellState.cpp @@ -51,7 +51,7 @@ SingleWellState(const std::string& name_, , prev_surface_rates(pu_.num_phases) , perf_data(perf_input.size(), pressure_first_connection, !is_producer, pu_.num_phases) , trivial_group_target(false) - , initialized_from_reservoir_(false) + , initialized_from_reservoir(false) { for (std::size_t perf = 0; perf < perf_input.size(); perf++) { this->perf_data.cell_index[perf] = perf_input[perf].cell_index; @@ -263,14 +263,20 @@ update_producer_targets(const Well& ecl_well, const SummaryState& st) case Well::ProducerCMode::GRUP: case Well::ProducerCMode::THP: case Well::ProducerCMode::BHP: + // Setting all rates to zero, they will be initialized properly in + // a call to WellInterface::initializeWellState() later, which will + // use the reservoir state to find a better initial value. + // This also applies to the ORAT/WRAT/GRAT above, but then only the + // rate not set in the above will be modified in + // WellInterface::initializeWellState(). if (this->pu.phase_used[BlackoilPhases::Liquid]) { - this->surface_rates[pu.phase_pos[BlackoilPhases::Liquid]] = 0;//-1000.0 * Opm::unit::cubic(Opm::unit::meter) / Opm::unit::day; + this->surface_rates[pu.phase_pos[BlackoilPhases::Liquid]] = 0.0; } if (this->pu.phase_used[BlackoilPhases::Aqua]) { - this->surface_rates[pu.phase_pos[BlackoilPhases::Aqua]] = 0;//-1000.0 * Opm::unit::cubic(Opm::unit::meter) / Opm::unit::day; + this->surface_rates[pu.phase_pos[BlackoilPhases::Aqua]] = 0.0; } if (this->pu.phase_used[BlackoilPhases::Vapour]){ - this->surface_rates[pu.phase_pos[BlackoilPhases::Vapour]] = 0;//-100000.0 * Opm::unit::cubic(Opm::unit::meter) / Opm::unit::day; + this->surface_rates[pu.phase_pos[BlackoilPhases::Vapour]] = 0.0; } break; diff --git a/opm/simulators/wells/SingleWellState.hpp b/opm/simulators/wells/SingleWellState.hpp index 34786b3119..55783147a1 100644 --- a/opm/simulators/wells/SingleWellState.hpp +++ b/opm/simulators/wells/SingleWellState.hpp @@ -65,7 +65,7 @@ public: serializer(well_potentials); serializer(productivity_index); serializer(implicit_ipr_a); - serializer(implicit_ipr_b); + serializer(implicit_ipr_b); serializer(surface_rates); serializer(reservoir_rates); serializer(prev_surface_rates); @@ -76,9 +76,8 @@ public: serializer(production_cmode); serializer(filtrate_conc); serializer(perf_data); - serializer(stw_primaryvar); - serializer(multiseg_primaryvar); - serializer(initialized_from_reservoir_); + serializer(primaryvar); + serializer(initialized_from_reservoir); } bool operator==(const SingleWellState&) const; @@ -118,7 +117,8 @@ public: Events events; WellInjectorCMode injection_cmode{WellInjectorCMode::CMODE_UNDEFINED}; WellProducerCMode production_cmode{WellProducerCMode::CMODE_UNDEFINED}; - + std::vector primaryvar; + bool initialized_from_reservoir{false}; /// Special purpose method to support dynamically rescaling a well's /// CTFs through WELPI. @@ -153,13 +153,8 @@ public: Scalar sum_filtrate_rate() const; Scalar sum_filtrate_total() const; - std::vector stw_primaryvar; - std::vector> multiseg_primaryvar; - bool initializedFromReservoir() const { return initialized_from_reservoir_; } - void setInitializedFromReservoir(bool value) { initialized_from_reservoir_ = value; } -private: - bool initialized_from_reservoir_ = false; +private: Scalar sum_connection_rates(const std::vector& connection_rates) const; }; diff --git a/opm/simulators/wells/StandardWellPrimaryVariables.cpp b/opm/simulators/wells/StandardWellPrimaryVariables.cpp index c960df0117..0540b6cd34 100644 --- a/opm/simulators/wells/StandardWellPrimaryVariables.cpp +++ b/opm/simulators/wells/StandardWellPrimaryVariables.cpp @@ -162,30 +162,28 @@ update(const WellState& well_state, value_[WQTotal] = 0.; } } - assert(ws.initializedFromReservoir()); - //if (std::abs(total_well_rate) > 0.) { - - bool prim_set = ws.stw_primaryvar.size()>0; - if(prim_set){ + + // Production wells should have been initialized using reservoir data. + assert(well_.isProducer() == ws.initialized_from_reservoir); + + if (ws.primaryvar.size() > 0) { if constexpr (has_wfrac_variable) { - value_[WFrac] = ws.stw_primaryvar[WFrac]; + value_[WFrac] = ws.primaryvar[WFrac]; } if constexpr (has_gfrac_variable) { - value_[GFrac] = ws.stw_primaryvar[GFrac]; + value_[GFrac] = ws.primaryvar[GFrac]; } if constexpr (Indices::enableSolvent) { - value_[SFrac] = ws.stw_primaryvar[SFrac]; + value_[SFrac] = ws.primaryvar[SFrac]; } - }else{ - assert(prim_set==false); - if (std::abs(total_well_rate) > 0.) { + } else if (std::abs(total_well_rate) > 0.) { if constexpr (has_wfrac_variable) { value_[WFrac] = well_.scalingFactor(pu.phase_pos[Water]) * ws.surface_rates[pu.phase_pos[Water]] / total_well_rate; } if constexpr (has_gfrac_variable) { value_[GFrac] = well_.scalingFactor(pu.phase_pos[Gas]) * - (ws.surface_rates[pu.phase_pos[Gas]] - - (Indices::enableSolvent ? ws.sum_solvent_rates() : 0.0) ) / total_well_rate ; + (ws.surface_rates[pu.phase_pos[Gas]] - + (Indices::enableSolvent ? ws.sum_solvent_rates() : 0.0) ) / total_well_rate ; } if constexpr (Indices::enableSolvent) { value_[SFrac] = well_.scalingFactor(Indices::contiSolventEqIdx) * ws.sum_solvent_rates() / total_well_rate ; @@ -234,8 +232,8 @@ update(const WellState& well_state, OPM_DEFLOG_THROW(std::logic_error, "Expected PRODUCER or INJECTOR type of well", deferred_logger); } } - } - //BHP + + // BHP value_[Bhp] = ws.bhp; setEvaluationsFromValues(); } @@ -348,6 +346,7 @@ copyToWellState(WellState& well_state, static constexpr int Water = BlackoilPhases::Aqua; static constexpr int Oil = BlackoilPhases::Liquid; static constexpr int Gas = BlackoilPhases::Vapour; + const PhaseUsage& pu = well_.phaseUsage(); std::vector F(well_.numPhases(), 0.0); [[maybe_unused]] Scalar F_solvent = 0.0; @@ -410,9 +409,9 @@ copyToWellState(WellState& well_state, F_solvent /= well_.scalingFactor(Indices::contiSolventEqIdx); F[pu.phase_pos[Gas]] += F_solvent; } - + auto& ws = well_state.well(well_.indexOfWell()); - ws.stw_primaryvar = value_; + ws.primaryvar = value_; ws.bhp = value_[Bhp]; // calculate the phase rates based on the primary variables diff --git a/opm/simulators/wells/WellInterface.hpp b/opm/simulators/wells/WellInterface.hpp index c83f3fb1f0..91174a5ded 100644 --- a/opm/simulators/wells/WellInterface.hpp +++ b/opm/simulators/wells/WellInterface.hpp @@ -334,10 +334,6 @@ public: WellState& well_state, DeferredLogger& deferred_logger) const; - // void updateWellStateRates(const Simulator& simulator, - // WellState& well_state, - // DeferredLogger& deferred_logger) const; - void solveWellEquation(const Simulator& simulator, WellState& well_state, const GroupState& group_state, diff --git a/opm/simulators/wells/WellInterface_impl.hpp b/opm/simulators/wells/WellInterface_impl.hpp index e3ebcfb78e..c568b9edd0 100644 --- a/opm/simulators/wells/WellInterface_impl.hpp +++ b/opm/simulators/wells/WellInterface_impl.hpp @@ -1598,57 +1598,6 @@ namespace Opm return scaling_factor; } -#if 0 - - template - void - WellInterface:: - initializeWellState(const Simulator& simulator, - const GroupState& group_state, - const SummaryState& summary_state, - WellState& well_state, - DeferredLogger& deferred_logger) const - { - OPM_TIMEFUNCTION(); - // Check if the rates of this well only are single-phase, do nothing - // if more than one nonzero rate. - auto& ws = well_state.well(this->index_of_well_); - // int nonzero_rate_index = -1; - // const Scalar floating_point_error_epsilon = 1e-14; - // for (int p = 0; p < this->number_of_phases_; ++p) { - // if (std::abs(ws.surface_rates[p]) > floating_point_error_epsilon) { - // if (nonzero_rate_index == -1) { - // nonzero_rate_index = p; - // } else { - // // More than one nonzero rate. - // return; - // } - // } - // } - // Calculate the rates that follow from the current primary variables. - std::vector well_q_s; - Scalar bhp_tmp = this->isInjector() - ? this->well_ecl_.injectionControls(summary_state).bhp_limit - : this->well_ecl_.productionControls(summary_state).bhp_limit; - computeWellRatesWithBhp(simulator, - bhp_tmp, - well_q_s, - deferred_logger); - // constcomputeCurrentWellRates(simulator, deferred_logger); - - //if (nonzero_rate_index == -1) { - // No nonzero rates. - // Use the computed rate directly - for (int p = 0; p < this->number_of_phases_; ++p) { - ws.surface_rates[p] = well_q_s[this->flowPhaseToModelCompIdx(p)]; - } - // set fractions - - ws.setInitializedFromReservoir(true); - //this->updateWellStateWithTarget(simulator, group_state, well_state, deferred_logger); - - } -#endif template @@ -1659,12 +1608,6 @@ namespace Opm const SummaryState& /* summary_state */, WellState& well_state, DeferredLogger& deferred_logger) const - // template - // void - // WellInterface:: - // updateWellStateRates(const Simulator& simulator, - // WellState& well_state, - // DeferredLogger& deferred_logger) const { OPM_TIMEFUNCTION(); // Check if the rates of this well only are single-phase, do nothing @@ -1684,8 +1627,10 @@ namespace Opm } // Calculate the rates that follow from the current primary variables. + // TODO switch to ccomputeWellRatesWithBhp(). std::vector well_q_s = computeCurrentWellRates(simulator, deferred_logger); - + ws.initialized_from_reservoir = true; + if (nonzero_rate_index == -1) { // No nonzero rates. // Use the computed rate directly @@ -1694,7 +1639,7 @@ namespace Opm } return; } - // set fractions + // Set the currently-zero phase flows to be nonzero in proportion to well_q_s. const Scalar initial_nonzero_rate = ws.surface_rates[nonzero_rate_index]; const int comp_idx_nz = this->flowPhaseToModelCompIdx(nonzero_rate_index);