My Project
Loading...
Searching...
No Matches
BlackoilWellModelGeneric.hpp
1/*
2 Copyright 2016 SINTEF ICT, Applied Mathematics.
3 Copyright 2016 - 2017 Statoil ASA.
4 Copyright 2017 Dr. Blatt - HPC-Simulation-Software & Services
5 Copyright 2016 - 2018 IRIS AS
6
7 This file is part of the Open Porous Media project (OPM).
8
9 OPM is free software: you can redistribute it and/or modify
10 it under the terms of the GNU General Public License as published by
11 the Free Software Foundation, either version 3 of the License, or
12 (at your option) any later version.
13
14 OPM is distributed in the hope that it will be useful,
15 but WITHOUT ANY WARRANTY; without even the implied warranty of
16 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17 GNU General Public License for more details.
18
19 You should have received a copy of the GNU General Public License
20 along with OPM. If not, see <http://www.gnu.org/licenses/>.
21*/
22
23#ifndef OPM_BLACKOILWELLMODEL_GENERIC_HEADER_INCLUDED
24#define OPM_BLACKOILWELLMODEL_GENERIC_HEADER_INCLUDED
25
26#include <opm/output/data/GuideRateValue.hpp>
27
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>
33
34#include <opm/simulators/utils/DeferredLoggingErrorHelpers.hpp>
35
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>
42
43#include <cstddef>
44#include <functional>
45#include <map>
46#include <memory>
47#include <optional>
48#include <string>
49#include <unordered_map>
50#include <unordered_set>
51#include <vector>
52
53namespace Opm {
54 class DeferredLogger;
55 class EclipseState;
56 template<class Scalar> class GasLiftGroupInfo;
57 template<class Scalar> class GasLiftSingleWellGeneric;
58 template<class Scalar> class GasLiftWellState;
59 class Group;
60 class GuideRateConfig;
61 template<class Scalar> class ParallelWellInfo;
62 class RestartValue;
63 class Schedule;
64 struct SimulatorUpdate;
65 class SummaryConfig;
66 template<class Scalar> class VFPProperties;
67 template<class Scalar> class WellInterfaceGeneric;
68 template<class Scalar> class WellState;
69} // namespace Opm
70
71namespace Opm { namespace data {
72 struct GroupData;
73 struct GroupGuideRates;
75 struct NodeData;
76}} // namespace Opm::data
77
78namespace Opm {
79
81template<class Scalar>
83{
84public:
85 // --------- Types ---------
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>>>;
89
91 const SummaryState& summaryState,
92 const EclipseState& eclState,
94 const Parallel::Communication& comm);
95
96 virtual ~BlackoilWellModelGeneric() = default;
97
98 int numLocalWells() const;
99 int numLocalWellsEnd() const;
100 int numLocalNonshutWells() const;
101 int numPhases() const;
102
104 bool wellsActive() const;
105 bool hasWell(const std::string& wname) const;
106
108 bool networkActive() const;
109
110 // whether there exists any multisegment well open on this process
111 bool anyMSWellOpenLocal() const;
112
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()}; }
120
121 /*
122 Immutable version of the currently active wellstate.
123 */
124 const WellState<Scalar>& wellState() const
125 {
126 return this->active_wgstate_.well_state;
127 }
128
129 /*
130 Mutable version of the currently active wellstate.
131 */
132 WellState<Scalar>& wellState()
133 {
134 return this->active_wgstate_.well_state;
135 }
136
137 /*
138 Will return the currently active nupcolWellState; must initialize
139 the internal nupcol wellstate with initNupcolWellState() first.
140 */
141 const WellState<Scalar>& nupcolWellState() const
142 {
143 return this->nupcol_wgstate_.well_state;
144 }
145 GroupState<Scalar>& groupState() { return this->active_wgstate_.group_state; }
146
147 WellTestState& wellTestState() { return this->active_wgstate_.well_test_state; }
148
149 const WellTestState& wellTestState() const { return this->active_wgstate_.well_test_state; }
150
151 Scalar wellPI(const int well_index) const;
152 Scalar wellPI(const std::string& well_name) const;
153
154 void updateEclWells(const int timeStepIdx,
156 const SummaryState& st);
157
158 void initFromRestartFile(const RestartValue& restartValues,
159 std::unique_ptr<WellTestState> wtestState,
160 const std::size_t numCells,
161 bool handle_ms_well);
162
163 void prepareDeserialize(int report_step,
164 const std::size_t numCells,
165 bool handle_ms_well);
166
167 /*
168 Will assign the internal member last_valid_well_state_ to the
169 current value of the this->active_well_state_. The state stored
170 with storeWellState() can then subsequently be recovered with the
171 resetWellState() method.
172 */
173 void commitWGState()
174 {
175 this->last_valid_wgstate_ = this->active_wgstate_;
176 }
177
178 data::GroupAndNetworkValues groupAndNetworkData(const int reportStepIdx) const;
179
181 bool hasTHPConstraints() const;
182
184 void updateNetworkActiveState(const int report_step);
185
189 bool needPreStepNetworkRebalance(const int report_step) const;
190
193 bool forceShutWellByName(const std::string& wellname,
194 const double simulation_time);
195
196 const std::vector<PerforationData<Scalar>>& perfData(const int well_idx) const
197 { return well_perf_data_[well_idx]; }
198
199 const Parallel::Communication& comm() const { return comm_; }
200
201 const EclipseState& eclipseState() const { return eclState_; }
202
203 const SummaryState& summaryState() const { return summaryState_; }
204
205 const GuideRate& guideRate() const { return guideRate_; }
206
207 bool reportStepStarts() const { return report_step_starts_; }
208
209 bool shouldBalanceNetwork(const int reportStepIndex,
210 const int iterationIdx) const;
211
212 void updateClosedWellsThisStep(const std::string& well_name) const
213 {
214 this->closed_this_step_.insert(well_name);
215 }
216 bool wasDynamicallyShutThisTimeStep(const std::string& well_name) const;
217
218 template<class Serializer>
219 void serializeOp(Serializer& serializer)
220 {
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_);
236 }
237
238 bool operator==(const BlackoilWellModelGeneric& rhs) const
239 {
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_;
254 }
255
256protected:
257 /*
258 The dynamic state of the well model is maintained with an instance
259 of the WellState class. Currently we have
260 three different wellstate instances:
261
262 1. The currently active wellstate is in the active_well_state_
263 member. That is the state which is mutated by the simulator.
264
265 2. In the case timestep fails to converge and we must go back and
266 try again with a smaller timestep we need to recover the last
267 valid wellstate. This is maintained with the
268 last_valid_well_state_ member and the functions
269 commitWellState() and resetWellState().
270
271 3. For the NUPCOL functionality we should either use the
272 currently active wellstate or a wellstate frozen at max
273 nupcol iterations. This is handled with the member
274 nupcol_well_state_ and the initNupcolWellState() function.
275 */
276
277 /*
278 Will return the last good wellstate. This is typcially used when
279 initializing a new report step where the Schedule object might
280 have introduced new wells. The wellstate returned by
281 prevWellState() must have been stored with the commitWellState()
282 function first.
283 */
284 const WellState<Scalar>& prevWellState() const
285 {
286 return this->last_valid_wgstate_.well_state;
287 }
288
289 const WGState<Scalar>& prevWGState() const
290 {
291 return this->last_valid_wgstate_;
292 }
293
294 /*
295 Will store a copy of the input argument well_state in the
296 last_valid_well_state_ member, that state can then be recovered
297 with a subsequent call to resetWellState().
298 */
299 void commitWGState(WGState<Scalar> wgstate)
300 {
301 this->last_valid_wgstate_ = std::move(wgstate);
302 }
303
304 /*
305 Will update the internal variable active_well_state_ to whatever
306 was stored in the last_valid_well_state_ member. This function
307 works in pair with commitWellState() which should be called first.
308 */
309 void resetWGState()
310 {
311 this->active_wgstate_ = this->last_valid_wgstate_;
312 }
313
314 /*
315 Will store the current active wellstate in the nupcol_well_state_
316 member. This can then be subsequently retrieved with accessor
317 nupcolWellState().
318 */
319 void updateNupcolWGState()
320 {
321 this->nupcol_wgstate_ = this->active_wgstate_;
322 }
323
326 std::vector<std::reference_wrapper<ParallelWellInfo<Scalar>>>
327 createLocalParallelWellInfo(const std::vector<Well>& wells);
328
329 void initializeWellProdIndCalculators();
330 void initializeWellPerfData();
331
332 bool wasDynamicallyShutThisTimeStep(const int well_index) const;
333
334 Scalar updateNetworkPressures(const int reportStepIdx);
335
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,
341 Scalar wsolvent);
342 virtual void calcResvCoeff(const int fipnum,
343 const int pvtreg,
344 const std::vector<Scalar>& production_rates,
345 std::vector<Scalar>& resv_coeff) = 0;
346 virtual void calcInjResvCoeff(const int fipnum,
347 const int pvtreg,
348 std::vector<Scalar>& resv_coeff) = 0;
349
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;
361
362 void calculateEfficiencyFactors(const int reportStepIdx);
363
364 void checkGconsaleLimits(const Group& group,
365 WellState<Scalar>& well_state,
366 const int reportStepIdx,
368
369 void checkGEconLimits(const Group& group,
370 const double simulation_time,
371 const int report_step_idx,
373
374 bool checkGroupHigherConstraints(const Group& group,
376 const int reportStepIdx);
377
378 void updateAndCommunicateGroupData(const int reportStepIdx,
379 const int iterationIdx);
380
381 void inferLocalShutWells();
382
383 void setRepRadiusPerfLength();
384
385 void gliftDebug(const std::string& msg,
387
388 void gliftDebugShowALQ(DeferredLogger& deferred_logger);
389
390 void gasLiftOptimizationStage2(DeferredLogger& deferred_logger,
391 GLiftProdWells& prod_wells,
392 GLiftOptWells& glift_wells,
394 GLiftWellStateMap& map,
395 const int episodeIndex);
396
397 virtual void computePotentials(const std::size_t widx,
399 std::string& exc_msg,
400 ExceptionType::ExcEnum& exc_type,
402
403 // Calculating well potentials for each well
404 void updateWellPotentials(const int reportStepIdx,
405 const bool onlyAfterEvent,
406 const SummaryConfig& summaryConfig,
408
409 void initInjMult();
410
411 void updateInjMult(DeferredLogger& deferred_logger);
412 void updateInjFCMult(DeferredLogger& deferred_logger);
413
414 void updateFiltrationModelsPostStep(const double dt,
415 const std::size_t water_index,
417
418 void updateFiltrationModelsPreStep(DeferredLogger& deferred_logger);
419
420 // create the well container
421 virtual void createWellContainer(const int time_step) = 0;
422 virtual void initWellContainer(const int reportStepIdx) = 0;
423
424 virtual void calculateProductivityIndexValuesShutWells(const int reportStepIdx,
426 virtual void calculateProductivityIndexValues(DeferredLogger& deferred_logger) = 0;
427
428 void runWellPIScaling(const int reportStepIdx,
430
433
434 std::vector<int> getCellsForConnections(const Well& well) const;
435 std::vector<std::vector<int>> getMaxWellConnections() const;
436
437 std::vector<std::string> getWellsForTesting(const int timeStepIdx,
438 const double simulationTime);
439
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;
448
449 Schedule& schedule_;
450 const SummaryState& summaryState_;
451 const EclipseState& eclState_;
452 const Parallel::Communication& comm_;
453
454 PhaseUsage phase_usage_;
455 bool terminal_output_{false};
456 bool wells_active_{false};
457 bool network_active_{false};
458 bool initial_step_{};
459 bool report_step_starts_{};
460
461 std::optional<int> last_run_wellpi_{};
462
463 std::vector<Well> wells_ecl_;
464 std::vector<std::vector<PerforationData<Scalar>>> well_perf_data_;
465
468 {
469 public:
474 explicit ConnectionIndexMap(const std::size_t numConns)
475 : local_(numConns, -1)
476 {
477 this->global_.reserve(numConns);
478 this->open_.reserve(numConns);
479 }
480
489 const bool connIsOpen)
490 {
491 this->local_[connIdx] =
492 static_cast<int>(this->global_.size());
493
494 this->global_.push_back(connIdx);
495
496 const auto open_conn_idx = connIsOpen
497 ? this->num_open_conns_++
498 : -1;
499
500 this->open_.push_back(open_conn_idx);
501 }
502
508 const std::vector<int>& local() const
509 {
510 return this->local_;
511 }
512
518 int global(const int connIdx) const
519 {
520 return this->global_[connIdx];
521 }
522
530 int open(const int connIdx) const
531 {
532 return this->open_[connIdx];
533 }
534
535 private:
539 std::vector<int> local_{};
540
543 std::vector<int> global_{};
544
546 std::vector<int> open_{};
547
549 int num_open_conns_{0};
550 };
551
552 std::vector<ConnectionIndexMap> conn_idx_map_{};
553 std::function<bool(const Well&)> not_on_process_{};
554
555 // a vector of all the wells.
556 std::vector<WellInterfaceGeneric<Scalar>*> well_container_generic_{};
557
558 std::vector<int> local_shut_wells_{};
559
560 std::vector<ParallelWellInfo<Scalar>> parallel_well_info_;
561 std::vector<std::reference_wrapper<ParallelWellInfo<Scalar>>> local_parallel_well_info_;
562
563 std::vector<WellProdIndexCalculator<Scalar>> prod_index_calc_;
564 mutable ParallelWBPCalculation<Scalar> wbpCalculationService_;
565
566 std::vector<int> pvt_region_idx_;
567
568 mutable std::unordered_set<std::string> closed_this_step_;
569
570 GuideRate guideRate_;
571 std::unique_ptr<VFPProperties<Scalar>> vfp_properties_{};
572 std::map<std::string, Scalar> node_pressures_; // Storing network pressures for output.
573
574 // previous injection multiplier, it is used in the injection multiplier calculation for WINJMULT keyword
575 std::unordered_map<std::string, std::vector<Scalar>> prev_inj_multipliers_;
576
577 // Handling for filter cake injection multipliers
578 std::unordered_map<std::string, WellFilterCake<Scalar>> filter_cake_;
579
580 /*
581 The various wellState members should be accessed and modified
582 through the accessor functions wellState(), prevWellState(),
583 commitWellState(), resetWellState(), nupcolWellState() and
584 updateNupcolWellState().
585 */
586 WGState<Scalar> active_wgstate_;
587 WGState<Scalar> last_valid_wgstate_;
588 WGState<Scalar> nupcol_wgstate_;
589
590 bool glift_debug = false;
591
592 double last_glift_opt_time_ = -1.0;
593
594 bool wellStructureChangedDynamically_{false};
595
596 // Store maps of group name and new group controls for output
597 std::map<std::string, std::string> switched_prod_groups_;
598 std::map<std::pair<std::string, Phase>, std::string> switched_inj_groups_;
599 // Store map of group name and close offending well for output
600 std::map<std::string, std::pair<std::string, std::string>> closed_offending_wells_;
601
602private:
603 WellInterfaceGeneric<Scalar>* getGenWell(const std::string& well_name);
604
605 template <typename Iter, typename Body>
606 void wellUpdateLoop(Iter first, Iter last, const int timeStepIdx, Body&& body);
607
608 void updateEclWellsConstraints(const int timeStepIdx,
610 const SummaryState& st);
611
612 void updateEclWellsCTFFromAction(const int timeStepIdx,
614};
615
616
617} // namespace Opm
618
619#endif
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
Class for handling the blackoil well model.
Definition BlackoilWellModelGeneric.hpp:83
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
Definition DeferredLogger.hpp:57
Definition GasLiftGroupInfo.hpp:46
Definition GroupState.hpp:38
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition WellState.hpp:62
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
Definition BlackoilPhases.hpp:46
Definition WGState.hpp:41