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/Well/WellTestState.hpp>
29 #include <opm/input/eclipse/Schedule/Group/GuideRate.hpp>
31 #include <opm/simulators/utils/DeferredLoggingErrorHelpers.hpp>
33 #include <opm/simulators/wells/PerforationData.hpp>
34 #include <opm/simulators/wells/WellProdIndexCalculator.hpp>
35 #include <opm/simulators/wells/WGState.hpp>
41 #include <unordered_map>
42 #include <unordered_set>
48 class GasLiftSingleWellGeneric;
49 class GasLiftWellState;
51 class GuideRateConfig;
52 class ParallelWellInfo;
57 class WellInterfaceGeneric;
61 namespace Opm {
namespace data {
63 struct GroupGuideRates;
64 class GroupAndNetworkValues;
75 using GLiftOptWells = std::map<std::string, std::unique_ptr<GasLiftSingleWellGeneric>>;
76 using GLiftProdWells = std::map<std::string, const WellInterfaceGeneric*>;
77 using GLiftWellStateMap = std::map<std::string, std::unique_ptr<GasLiftWellState>>;
80 const SummaryState& summaryState,
81 const EclipseState& eclState,
83 const Parallel::Communication& comm);
87 int numLocalWells()
const;
88 int numPhases()
const;
92 bool hasWell(
const std::string& wname);
95 bool anyMSWellOpenLocal()
const;
97 const Well& getWellEcl(
const std::string& well_name)
const;
98 std::vector<Well> getLocalWells(
const int timeStepIdx)
const;
99 const Schedule& schedule()
const {
return schedule_; }
100 const PhaseUsage& phaseUsage()
const {
return phase_usage_; }
101 const GroupState& groupState()
const {
return this->active_wgstate_.group_state; }
108 return this->active_wgstate_.well_state;
116 return this->active_wgstate_.well_state;
119 GroupState& groupState() {
return this->active_wgstate_.group_state; }
121 WellTestState& wellTestState() {
return this->active_wgstate_.well_test_state; }
123 const WellTestState& wellTestState()
const {
return this->active_wgstate_.well_test_state; }
126 double wellPI(
const int well_index)
const;
127 double wellPI(
const std::string& well_name)
const;
129 void updateEclWells(
const int timeStepIdx,
130 const std::unordered_set<std::string>& wells,
131 const SummaryState& st);
134 void loadRestartData(
const data::Wells& rst_wells,
135 const data::GroupAndNetworkValues& grpNwrkValues,
137 const bool handle_ms_well,
140 void initFromRestartFile(
const RestartValue& restartValues,
141 WellTestState wtestState,
142 const size_t numCells,
143 bool handle_ms_well);
153 this->last_valid_wgstate_ = this->active_wgstate_;
156 data::GroupAndNetworkValues groupAndNetworkData(
const int reportStepIdx)
const;
164 const double simulation_time);
198 return this->last_valid_wgstate_.well_state;
201 const WGState& prevWGState()
const
203 return this->last_valid_wgstate_;
211 return this->nupcol_wgstate_.well_state;
219 void commitWGState(
WGState wgstate)
221 this->last_valid_wgstate_ = std::move(wgstate);
231 this->active_wgstate_ = this->last_valid_wgstate_;
239 void updateNupcolWGState()
241 this->nupcol_wgstate_ = this->active_wgstate_;
248 void initializeWellProdIndCalculators();
249 void initializeWellPerfData();
251 bool wasDynamicallyShutThisTimeStep(
const int well_index)
const;
253 std::pair<bool, double> updateNetworkPressures(
const int reportStepIdx);
255 void updateWsolvent(
const Group& group,
256 const int reportStepIdx,
258 void setWsolvent(
const Group& group,
259 const int reportStepIdx,
261 virtual void calcRates(
const int fipnum,
263 std::vector<double>& resv_coeff) = 0;
264 virtual void calcInjRates(
const int fipnum,
266 std::vector<double>& resv_coeff) = 0;
268 data::GuideRateValue getGuideRateValues(
const Group& group)
const;
269 data::GuideRateValue getGuideRateValues(
const Well& well)
const;
270 data::GuideRateValue getGuideRateInjectionGroupValues(
const Group& group)
const;
271 void getGuideRateValues(
const GuideRate::RateVector& qs,
273 const std::string& wgname,
274 data::GuideRateValue& grval)
const;
276 void assignWellGuideRates(data::Wells& wsrpt,
277 const int reportStepIdx)
const;
278 void assignShutConnections(data::Wells& wsrpt,
279 const int reportStepIndex)
const;
280 void assignGroupControl(
const Group& group,
281 data::GroupData& gdata)
const;
282 void assignGroupGuideRates(
const Group& group,
283 const std::unordered_map<std::string, data::GroupGuideRates>& groupGuideRates,
284 data::GroupData& gdata)
const;
285 void assignGroupValues(
const int reportStepIdx,
286 std::map<std::string, data::GroupData>& gvalues)
const;
287 void assignNodeValues(std::map<std::string, data::NodeData>& nodevalues)
const;
289 void loadRestartConnectionData(
const std::vector<data::Rates::opt>& phs,
290 const data::Well& rst_well,
291 const std::vector<PerforationData>& old_perf_data,
294 void loadRestartSegmentData(
const std::string& well_name,
295 const std::vector<data::Rates::opt>& phs,
296 const data::Well& rst_well,
299 void loadRestartWellData(
const std::string& well_name,
300 const bool handle_ms_well,
301 const std::vector<data::Rates::opt>& phs,
302 const data::Well& rst_well,
303 const std::vector<PerforationData>& old_perf_data,
306 void loadRestartGroupData(
const std::string& group,
307 const data::GroupData& value);
309 void loadRestartGuideRates(
const int report_step,
310 const GuideRateModel::Target target,
311 const data::Wells& rst_wells);
313 void loadRestartGuideRates(
const int report_step,
314 const GuideRateConfig& config,
315 const std::map<std::string, data::GroupData>& rst_groups);
317 std::unordered_map<std::string, data::GroupGuideRates>
318 calculateAllGroupGuiderates(
const int reportStepIdx)
const;
320 void calculateEfficiencyFactors(
const int reportStepIdx);
322 bool checkGroupConstraints(
const Group& group,
323 const int reportStepIdx,
326 std::pair<Group::InjectionCMode, double> checkGroupInjectionConstraints(
const Group& group,
327 const int reportStepIdx,
328 const Phase& phase)
const;
329 std::pair<Group::ProductionCMode, double> checkGroupProductionConstraints(
const Group& group,
330 const int reportStepIdx,
333 void checkGconsaleLimits(
const Group& group,
335 const int reportStepIdx,
338 bool checkGroupHigherConstraints(
const Group& group,
340 const int reportStepIdx);
342 bool updateGroupIndividualControl(
const Group& group,
344 const int reportStepIdx);
346 void actionOnBrokenConstraints(
const Group& group,
347 const Group::ExceedAction& exceed_action,
348 const Group::ProductionCMode& newControl,
350 void actionOnBrokenConstraints(
const Group& group,
351 const Group::InjectionCMode& newControl,
352 const Phase& controlPhase,
355 void updateAndCommunicateGroupData(
const int reportStepIdx,
356 const int iterationIdx);
358 void inferLocalShutWells();
360 void setRepRadiusPerfLength();
362 void gliftDebug(
const std::string& msg,
368 GLiftProdWells& prod_wells,
369 GLiftOptWells& glift_wells,
370 GLiftWellStateMap& map,
371 const int episodeIndex);
373 virtual void computePotentials(
const std::size_t widx,
375 std::string& exc_msg,
376 ExceptionType::ExcEnum& exc_type,
380 void updateWellPotentials(
const int reportStepIdx,
381 const bool onlyAfterEvent,
382 const SummaryConfig& summaryConfig,
385 bool guideRateUpdateIsNeeded(
const int reportStepIdx)
const;
388 virtual void createWellContainer(
const int time_step) = 0;
389 virtual void initWellContainer(
const int reportStepIdx) = 0;
391 virtual void calculateProductivityIndexValuesShutWells(
const int reportStepIdx,
393 virtual void calculateProductivityIndexValues(
DeferredLogger& deferred_logger) = 0;
395 void runWellPIScaling(
const int timeStepIdx,
403 const SummaryState& summaryState_;
404 const EclipseState& eclState_;
405 const Parallel::Communication& comm_;
408 bool terminal_output_{
false};
409 bool wells_active_{
false};
410 bool initial_step_{};
411 bool report_step_starts_{};
413 std::optional<int> last_run_wellpi_{};
415 std::vector<Well> wells_ecl_;
416 std::vector<std::vector<PerforationData>> well_perf_data_;
417 std::function<bool(
const Well&)> not_on_process_{};
420 std::vector<WellInterfaceGeneric*> well_container_generic_{};
422 std::vector<int> local_shut_wells_{};
424 std::vector<ParallelWellInfo> parallel_well_info_;
425 std::vector<std::reference_wrapper<ParallelWellInfo>> local_parallel_well_info_;
427 std::vector<WellProdIndexCalculator> prod_index_calc_;
429 std::vector<int> pvt_region_idx_;
431 mutable std::unordered_set<std::string> closed_this_step_;
433 GuideRate guideRate_;
434 std::unique_ptr<VFPProperties> vfp_properties_{};
435 std::map<std::string, double> node_pressures_;
443 WGState active_wgstate_;
444 WGState last_valid_wgstate_;
445 WGState nupcol_wgstate_;
447 bool glift_debug =
false;
449 double last_glift_opt_time_ = -1.0;
451 std::map<std::string, std::string> switched_prod_groups_;
452 std::map<std::pair<std::string, Opm::Phase>, std::string> switched_inj_groups_;
456 WellInterfaceGeneric* getGenWell(
const std::string& well_name);
Class for handling the blackoil well model.
Definition: BlackoilWellModelGeneric.hpp:72
virtual int compressedIndexForInterior(int cartesian_cell_idx) const =0
get compressed index for interior cells (-1, otherwise
std::vector< std::reference_wrapper< ParallelWellInfo > > createLocalParallelWellInfo(const std::vector< Well > &wells)
Create the parallel well information.
Definition: BlackoilWellModelGeneric.cpp:677
bool wellsActive() const
return true if wells are available in the reservoir
Definition: BlackoilWellModelGeneric.cpp:384
bool hasTHPConstraints() const
Return true if any well has a THP constraint.
Definition: BlackoilWellModelGeneric.cpp:2003
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:2016
Definition: DeferredLogger.hpp:57
Definition: GroupState.hpp:34
Definition: SingleWellState.hpp:38
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition: WellState.hpp:56
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:27
Definition: BlackoilPhases.hpp:46
Definition: WGState.hpp:35