23#ifndef OPM_BLACKOILWELLMODEL_GENERIC_HEADER_INCLUDED
24#define OPM_BLACKOILWELLMODEL_GENERIC_HEADER_INCLUDED
26#include <opm/output/data/GuideRateValue.hpp>
28#include <opm/input/eclipse/Schedule/Group/GuideRate.hpp>
29#include <opm/input/eclipse/Schedule/Well/PAvg.hpp>
30#include <opm/input/eclipse/Schedule/Well/PAvgCalculator.hpp>
31#include <opm/input/eclipse/Schedule/Well/PAvgCalculatorCollection.hpp>
32#include <opm/input/eclipse/Schedule/Well/WellTestState.hpp>
34#include <opm/simulators/utils/DeferredLoggingErrorHelpers.hpp>
36#include <opm/simulators/wells/ParallelPAvgDynamicSourceData.hpp>
37#include <opm/simulators/wells/ParallelWBPCalculation.hpp>
38#include <opm/simulators/wells/PerforationData.hpp>
39#include <opm/simulators/wells/WellFilterCake.hpp>
40#include <opm/simulators/wells/WellProdIndexCalculator.hpp>
41#include <opm/simulators/wells/WGState.hpp>
49#include <unordered_map>
50#include <unordered_set>
56 template<
class Scalar>
class GasLiftGroupInfo;
57 template<
class Scalar>
class GasLiftSingleWellGeneric;
58 template<
class Scalar>
class GasLiftWellState;
60 class GuideRateConfig;
61 template<
class Scalar>
class ParallelWellInfo;
64 struct SimulatorUpdate;
66 template<
class Scalar>
class VFPProperties;
67 template<
class Scalar>
class WellInterfaceGeneric;
68 template<
class Scalar>
class WellState;
71namespace Opm {
namespace data {
73 struct GroupGuideRates;
74 class GroupAndNetworkValues;
82class BlackoilWellModelGeneric
86 using GLiftOptWells = std::map<std::string, std::unique_ptr<GasLiftSingleWellGeneric<Scalar>>>;
87 using GLiftProdWells = std::map<std::string, const WellInterfaceGeneric<Scalar>*>;
88 using GLiftWellStateMap = std::map<std::string, std::unique_ptr<GasLiftWellState<Scalar>>>;
90 BlackoilWellModelGeneric(Schedule& schedule,
91 const SummaryState& summaryState,
92 const EclipseState& eclState,
93 const PhaseUsage& phase_usage,
94 const Parallel::Communication& comm);
96 virtual ~BlackoilWellModelGeneric() =
default;
98 int numLocalWells()
const;
99 int numLocalWellsEnd()
const;
100 int numLocalNonshutWells()
const;
101 int numPhases()
const;
105 bool hasWell(
const std::string& wname)
const;
111 bool anyMSWellOpenLocal()
const;
113 const Well& getWellEcl(
const std::string& well_name)
const;
114 std::vector<Well> getLocalWells(
const int timeStepIdx)
const;
115 const Schedule& schedule()
const {
return schedule_; }
116 const PhaseUsage& phaseUsage()
const {
return phase_usage_; }
117 const GroupState<Scalar>& groupState()
const {
return this->active_wgstate_.group_state; }
118 std::vector<const WellInterfaceGeneric<Scalar>*> genericWells()
const
119 {
return {well_container_generic_.begin(), well_container_generic_.end()}; }
124 const WellState<Scalar>& wellState()
const
126 return this->active_wgstate_.well_state;
132 WellState<Scalar>& wellState()
134 return this->active_wgstate_.well_state;
141 const WellState<Scalar>& nupcolWellState()
const
143 return this->nupcol_wgstate_.well_state;
145 GroupState<Scalar>& groupState() {
return this->active_wgstate_.group_state; }
147 WellTestState& wellTestState() {
return this->active_wgstate_.well_test_state; }
149 const WellTestState& wellTestState()
const {
return this->active_wgstate_.well_test_state; }
151 Scalar wellPI(
const int well_index)
const;
152 Scalar wellPI(
const std::string& well_name)
const;
154 void updateEclWells(
const int timeStepIdx,
155 const SimulatorUpdate& sim_update,
156 const SummaryState& st);
158 void initFromRestartFile(
const RestartValue& restartValues,
159 std::unique_ptr<WellTestState> wtestState,
160 const std::size_t numCells,
161 bool handle_ms_well);
163 void prepareDeserialize(
int report_step,
164 const std::size_t numCells,
165 bool handle_ms_well);
175 this->last_valid_wgstate_ = this->active_wgstate_;
178 data::GroupAndNetworkValues groupAndNetworkData(
const int reportStepIdx)
const;
194 const double simulation_time);
196 const std::vector<PerforationData<Scalar>>& perfData(
const int well_idx)
const
197 {
return well_perf_data_[well_idx]; }
199 const Parallel::Communication& comm()
const {
return comm_; }
201 const EclipseState& eclipseState()
const {
return eclState_; }
203 const SummaryState& summaryState()
const {
return summaryState_; }
205 const GuideRate& guideRate()
const {
return guideRate_; }
207 bool reportStepStarts()
const {
return report_step_starts_; }
209 bool shouldBalanceNetwork(
const int reportStepIndex,
210 const int iterationIdx)
const;
212 void updateClosedWellsThisStep(
const std::string& well_name)
const
214 this->closed_this_step_.insert(well_name);
216 bool wasDynamicallyShutThisTimeStep(
const std::string& well_name)
const;
218 template<
class Serializer>
219 void serializeOp(Serializer& serializer)
221 serializer(initial_step_);
222 serializer(report_step_starts_);
223 serializer(last_run_wellpi_);
224 serializer(local_shut_wells_);
225 serializer(closed_this_step_);
226 serializer(guideRate_);
227 serializer(node_pressures_);
228 serializer(prev_inj_multipliers_);
229 serializer(active_wgstate_);
230 serializer(last_valid_wgstate_);
231 serializer(nupcol_wgstate_);
232 serializer(last_glift_opt_time_);
233 serializer(switched_prod_groups_);
234 serializer(switched_inj_groups_);
235 serializer(closed_offending_wells_);
238 bool operator==(
const BlackoilWellModelGeneric& rhs)
const
240 return this->initial_step_ == rhs.initial_step_ &&
241 this->report_step_starts_ == rhs.report_step_starts_ &&
242 this->last_run_wellpi_ == rhs.last_run_wellpi_ &&
243 this->local_shut_wells_ == rhs.local_shut_wells_ &&
244 this->closed_this_step_ == rhs.closed_this_step_ &&
245 this->node_pressures_ == rhs.node_pressures_ &&
246 this->prev_inj_multipliers_ == rhs.prev_inj_multipliers_ &&
247 this->active_wgstate_ == rhs.active_wgstate_ &&
248 this->last_valid_wgstate_ == rhs.last_valid_wgstate_ &&
249 this->nupcol_wgstate_ == rhs.nupcol_wgstate_ &&
250 this->last_glift_opt_time_ == rhs.last_glift_opt_time_ &&
251 this->switched_prod_groups_ == rhs.switched_prod_groups_ &&
252 this->switched_inj_groups_ == rhs.switched_inj_groups_ &&
253 this->closed_offending_wells_ == rhs.closed_offending_wells_;
284 const WellState<Scalar>& prevWellState()
const
286 return this->last_valid_wgstate_.well_state;
289 const WGState<Scalar>& prevWGState()
const
291 return this->last_valid_wgstate_;
299 void commitWGState(WGState<Scalar> wgstate)
301 this->last_valid_wgstate_ = std::move(wgstate);
311 this->active_wgstate_ = this->last_valid_wgstate_;
319 void updateNupcolWGState()
321 this->nupcol_wgstate_ = this->active_wgstate_;
326 std::vector<std::reference_wrapper<ParallelWellInfo<Scalar>>>
329 void initializeWellProdIndCalculators();
330 void initializeWellPerfData();
332 bool wasDynamicallyShutThisTimeStep(
const int well_index)
const;
334 Scalar updateNetworkPressures(
const int reportStepIdx);
336 void updateWsolvent(
const Group& group,
337 const int reportStepIdx,
338 const WellState<Scalar>& wellState);
339 void setWsolvent(
const Group& group,
340 const int reportStepIdx,
342 virtual void calcResvCoeff(
const int fipnum,
344 const std::vector<Scalar>& production_rates,
345 std::vector<Scalar>& resv_coeff) = 0;
346 virtual void calcInjResvCoeff(
const int fipnum,
348 std::vector<Scalar>& resv_coeff) = 0;
350 void assignShutConnections(data::Wells& wsrpt,
351 const int reportStepIndex)
const;
352 void assignWellTargets(data::Wells& wsrpt)
const;
353 void assignProductionWellTargets(
const Well& well, data::WellControlLimits& limits)
const;
354 void assignInjectionWellTargets(
const Well& well, data::WellControlLimits& limits)
const;
355 void assignGroupControl(
const Group& group,
356 data::GroupData& gdata)
const;
357 void assignGroupValues(
const int reportStepIdx,
358 std::map<std::string, data::GroupData>& gvalues)
const;
359 void assignNodeValues(std::map<std::string, data::NodeData>& nodevalues,
360 const int reportStepIdx)
const;
362 void calculateEfficiencyFactors(
const int reportStepIdx);
364 void checkGconsaleLimits(
const Group& group,
365 WellState<Scalar>& well_state,
366 const int reportStepIdx,
367 DeferredLogger& deferred_logger);
369 void checkGEconLimits(
const Group& group,
370 const double simulation_time,
371 const int report_step_idx,
372 DeferredLogger& deferred_logger);
374 bool checkGroupHigherConstraints(
const Group& group,
375 DeferredLogger& deferred_logger,
376 const int reportStepIdx);
378 void updateAndCommunicateGroupData(
const int reportStepIdx,
379 const int iterationIdx);
381 void inferLocalShutWells();
383 void setRepRadiusPerfLength();
385 void gliftDebug(
const std::string& msg,
386 DeferredLogger& deferred_logger)
const;
388 void gliftDebugShowALQ(DeferredLogger& deferred_logger);
390 void gasLiftOptimizationStage2(DeferredLogger& deferred_logger,
391 GLiftProdWells& prod_wells,
392 GLiftOptWells& glift_wells,
393 GasLiftGroupInfo<Scalar>& group_info,
394 GLiftWellStateMap& map,
395 const int episodeIndex);
397 virtual void computePotentials(
const std::size_t widx,
398 const WellState<Scalar>& well_state_copy,
399 std::string& exc_msg,
400 ExceptionType::ExcEnum& exc_type,
401 DeferredLogger& deferred_logger) = 0;
404 void updateWellPotentials(
const int reportStepIdx,
405 const bool onlyAfterEvent,
406 const SummaryConfig& summaryConfig,
407 DeferredLogger& deferred_logger);
411 void updateInjMult(DeferredLogger& deferred_logger);
412 void updateInjFCMult(DeferredLogger& deferred_logger);
414 void updateFiltrationModelsPostStep(
const double dt,
415 const std::size_t water_index,
416 DeferredLogger& deferred_logger);
418 void updateFiltrationModelsPreStep(DeferredLogger& deferred_logger);
421 virtual void createWellContainer(
const int time_step) = 0;
422 virtual void initWellContainer(
const int reportStepIdx) = 0;
424 virtual void calculateProductivityIndexValuesShutWells(
const int reportStepIdx,
425 DeferredLogger& deferred_logger) = 0;
426 virtual void calculateProductivityIndexValues(DeferredLogger& deferred_logger) = 0;
428 void runWellPIScaling(
const int reportStepIdx,
429 DeferredLogger& local_deferredLogger);
434 std::vector<int> getCellsForConnections(
const Well& well)
const;
435 std::vector<std::vector<int>> getMaxWellConnections()
const;
437 std::vector<std::string> getWellsForTesting(
const int timeStepIdx,
438 const double simulationTime);
440 using WellTracerRates = std::map<std::pair<std::string, std::string>, Scalar>;
441 void assignWellTracerRates(data::Wells& wsrpt,
442 const WellTracerRates& wellTracerRates)
const;
443 using MswTracerRates = std::map<std::tuple<std::string, std::string, std::size_t>, Scalar>;
444 void assignMswTracerRates(data::Wells& wsrpt,
445 const MswTracerRates& mswTracerRates)
const;
446 void assignMassGasRate(data::Wells& wsrpt,
447 const Scalar& gasDensity)
const;
450 const SummaryState& summaryState_;
451 const EclipseState& eclState_;
452 const Parallel::Communication& comm_;
455 bool terminal_output_{
false};
456 bool wells_active_{
false};
457 bool network_active_{
false};
458 bool initial_step_{};
459 bool report_step_starts_{};
461 std::optional<int> last_run_wellpi_{};
463 std::vector<Well> wells_ecl_;
464 std::vector<std::vector<PerforationData<Scalar>>> well_perf_data_;
475 : local_(numConns, -1)
477 this->global_.reserve(numConns);
478 this->open_.reserve(numConns);
489 const bool connIsOpen)
491 this->local_[connIdx] =
492 static_cast<int>(this->global_.size());
494 this->global_.push_back(connIdx);
496 const auto open_conn_idx = connIsOpen
497 ? this->num_open_conns_++
500 this->open_.push_back(open_conn_idx);
508 const std::vector<int>&
local()
const
520 return this->global_[connIdx];
530 int open(
const int connIdx)
const
532 return this->open_[connIdx];
539 std::vector<int> local_{};
543 std::vector<int> global_{};
546 std::vector<int> open_{};
549 int num_open_conns_{0};
552 std::vector<ConnectionIndexMap> conn_idx_map_{};
553 std::function<bool(
const Well&)> not_on_process_{};
556 std::vector<WellInterfaceGeneric<Scalar>*> well_container_generic_{};
558 std::vector<int> local_shut_wells_{};
560 std::vector<ParallelWellInfo<Scalar>> parallel_well_info_;
561 std::vector<std::reference_wrapper<ParallelWellInfo<Scalar>>> local_parallel_well_info_;
563 std::vector<WellProdIndexCalculator<Scalar>> prod_index_calc_;
564 mutable ParallelWBPCalculation<Scalar> wbpCalculationService_;
566 std::vector<int> pvt_region_idx_;
568 mutable std::unordered_set<std::string> closed_this_step_;
570 GuideRate guideRate_;
571 std::unique_ptr<VFPProperties<Scalar>> vfp_properties_{};
572 std::map<std::string, Scalar> node_pressures_;
575 std::unordered_map<std::string, std::vector<Scalar>> prev_inj_multipliers_;
578 std::unordered_map<std::string, WellFilterCake<Scalar>> filter_cake_;
586 WGState<Scalar> active_wgstate_;
587 WGState<Scalar> last_valid_wgstate_;
588 WGState<Scalar> nupcol_wgstate_;
590 bool glift_debug =
false;
592 double last_glift_opt_time_ = -1.0;
594 bool wellStructureChangedDynamically_{
false};
597 std::map<std::string, std::string> switched_prod_groups_;
598 std::map<std::pair<std::string, Phase>, std::string> switched_inj_groups_;
600 std::map<std::string, std::pair<std::string, std::string>> closed_offending_wells_;
603 WellInterfaceGeneric<Scalar>* getGenWell(
const std::string& well_name);
605 template <
typename Iter,
typename Body>
606 void wellUpdateLoop(Iter first, Iter last,
const int timeStepIdx, Body&& body);
608 void updateEclWellsConstraints(
const int timeStepIdx,
609 const SimulatorUpdate& sim_update,
610 const SummaryState& st);
612 void updateEclWellsCTFFromAction(
const int timeStepIdx,
613 const SimulatorUpdate& sim_update);
Connection index mappings.
Definition BlackoilWellModelGeneric.hpp:468
int open(const int connIdx) const
Get open connection ID of local (on-rank) connection.
Definition BlackoilWellModelGeneric.hpp:530
const std::vector< int > & local() const
Get local connection IDs/indices of every existing well connection.
Definition BlackoilWellModelGeneric.hpp:508
int global(const int connIdx) const
Get global connection ID of local (on-rank) connection.
Definition BlackoilWellModelGeneric.hpp:518
ConnectionIndexMap(const std::size_t numConns)
Constructor.
Definition BlackoilWellModelGeneric.hpp:474
void addActiveConnection(const int connIdx, const bool connIsOpen)
Enumerate/map new active connection.
Definition BlackoilWellModelGeneric.hpp:488
void updateNetworkActiveState(const int report_step)
Checks if network is active (at least one network well on prediction).
Definition BlackoilWellModelGeneric.cpp:1280
bool wellsActive() const
return true if wells are available in the reservoir
Definition BlackoilWellModelGeneric.cpp:158
bool networkActive() const
return true if network is active (at least one network well in prediction mode)
Definition BlackoilWellModelGeneric.cpp:165
bool forceShutWellByName(const std::string &wellname, const double simulation_time)
Shut down any single well Returns true if the well was actually found and shut.
Definition BlackoilWellModelGeneric.cpp:1320
std::vector< std::reference_wrapper< ParallelWellInfo< Scalar > > > createLocalParallelWellInfo(const std::vector< Well > &wells)
Create the parallel well information.
Definition BlackoilWellModelGeneric.cpp:291
bool needPreStepNetworkRebalance(const int report_step) const
Checks if there are reasons to perform a pre-step network re-balance.
Definition BlackoilWellModelGeneric.cpp:1301
bool hasTHPConstraints() const
Return true if any well has a THP constraint.
Definition BlackoilWellModelGeneric.cpp:1273
virtual int compressedIndexForInterior(int cartesian_cell_idx) const =0
get compressed index for interior cells (-1, otherwise
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
Definition BlackoilPhases.hpp:46