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 <type_traits>
50#include <unordered_map>
51#include <unordered_set>
52#include <vector>
53
54namespace Opm {
55 class DeferredLogger;
56 class EclipseState;
57 class GasLiftSingleWellGeneric;
58 class GasLiftWellState;
59 class GasLiftGroupInfo;
60 class Group;
61 class GuideRateConfig;
62 class ParallelWellInfo;
63 class RestartValue;
64 class Schedule;
65 class SummaryConfig;
66 class VFPProperties;
67 class WellInterfaceGeneric;
68 class WellState;
69} // namespace Opm
70
71namespace Opm { namespace data {
72 struct GroupData;
73 struct GroupGuideRates;
74 class GroupAndNetworkValues;
75 struct NodeData;
76}} // namespace Opm::data
77
78namespace Opm {
79
82{
83public:
84 // --------- Types ---------
85 using GLiftOptWells = std::map<std::string, std::unique_ptr<GasLiftSingleWellGeneric>>;
86 using GLiftProdWells = std::map<std::string, const WellInterfaceGeneric*>;
87 using GLiftWellStateMap = std::map<std::string, std::unique_ptr<GasLiftWellState>>;
88
90 const SummaryState& summaryState,
91 const EclipseState& eclState,
93 const Parallel::Communication& comm);
94
95 virtual ~BlackoilWellModelGeneric() = default;
96
97 int numLocalWells() const;
98 int numLocalWellsEnd() const;
99 int numLocalNonshutWells() const;
100 int numPhases() const;
101
103 bool wellsActive() const;
104 bool hasWell(const std::string& wname) const;
105
106 // whether there exists any multisegment well open on this process
107 bool anyMSWellOpenLocal() const;
108
109 const Well& getWellEcl(const std::string& well_name) const;
110 std::vector<Well> getLocalWells(const int timeStepIdx) const;
111 const Schedule& schedule() const { return schedule_; }
112 const PhaseUsage& phaseUsage() const { return phase_usage_; }
113 const GroupState& groupState() const { return this->active_wgstate_.group_state; }
114 std::vector<const WellInterfaceGeneric*> genericWells() const
115 { return {well_container_generic_.begin(), well_container_generic_.end()}; }
116
117 /*
118 Immutable version of the currently active wellstate.
119 */
120 const WellState& wellState() const
121 {
122 return this->active_wgstate_.well_state;
123 }
124
125 /*
126 Mutable version of the currently active wellstate.
127 */
128 WellState& wellState()
129 {
130 return this->active_wgstate_.well_state;
131 }
132
133 GroupState& groupState() { return this->active_wgstate_.group_state; }
134
135 WellTestState& wellTestState() { return this->active_wgstate_.well_test_state; }
136
137 const WellTestState& wellTestState() const { return this->active_wgstate_.well_test_state; }
138
139
140 double wellPI(const int well_index) const;
141 double wellPI(const std::string& well_name) const;
142
143 void updateEclWells(const int timeStepIdx,
144 const std::unordered_set<std::string>& wells,
145 const SummaryState& st);
146
147 void initFromRestartFile(const RestartValue& restartValues,
149 const std::size_t numCells,
150 bool handle_ms_well);
151
152 void prepareDeserialize(int report_step,
153 const std::size_t numCells,
154 bool handle_ms_well);
155
156 /*
157 Will assign the internal member last_valid_well_state_ to the
158 current value of the this->active_well_state_. The state stored
159 with storeWellState() can then subsequently be recovered with the
160 resetWellState() method.
161 */
162 void commitWGState()
163 {
164 this->last_valid_wgstate_ = this->active_wgstate_;
165 }
166
167 data::GroupAndNetworkValues groupAndNetworkData(const int reportStepIdx) const;
168
170 bool hasTHPConstraints() const;
171
173 bool needRebalanceNetwork(const int report_step) const;
174
177 bool forceShutWellByName(const std::string& wellname,
178 const double simulation_time);
179
180 const std::vector<PerforationData>& perfData(const int well_idx) const
181 { return well_perf_data_[well_idx]; }
182
183 const Parallel::Communication& comm() const { return comm_; }
184
185 const EclipseState& eclipseState() const { return eclState_; }
186
187 const SummaryState& summaryState() const { return summaryState_; }
188
189 const GuideRate& guideRate() const { return guideRate_; }
190
191 bool reportStepStarts() const { return report_step_starts_; }
192
193 bool shouldBalanceNetwork(const int reportStepIndex,
194 const int iterationIdx) const;
195
196 void updateClosedWellsThisStep(const std::string& well_name) const {
197 this->closed_this_step_.insert(well_name);
198 }
199
200 template<class Serializer>
201 void serializeOp(Serializer& serializer)
202 {
203 serializer(initial_step_);
204 serializer(report_step_starts_);
205 serializer(last_run_wellpi_);
206 serializer(local_shut_wells_);
207 serializer(closed_this_step_);
208 serializer(guideRate_);
209 serializer(node_pressures_);
210 serializer(prev_inj_multipliers_);
211 serializer(active_wgstate_);
212 serializer(last_valid_wgstate_);
213 serializer(nupcol_wgstate_);
214 serializer(last_glift_opt_time_);
215 serializer(switched_prod_groups_);
216 serializer(switched_inj_groups_);
217 }
218
219 bool operator==(const BlackoilWellModelGeneric& rhs) const
220 {
221 return this->initial_step_ == rhs.initial_step_ &&
222 this->report_step_starts_ == rhs.report_step_starts_ &&
223 this->last_run_wellpi_ == rhs.last_run_wellpi_ &&
224 this->local_shut_wells_ == rhs.local_shut_wells_ &&
225 this->closed_this_step_ == rhs.closed_this_step_ &&
226 this->node_pressures_ == rhs.node_pressures_ &&
227 this->prev_inj_multipliers_ == rhs.prev_inj_multipliers_ &&
228 this->active_wgstate_ == rhs.active_wgstate_ &&
229 this->last_valid_wgstate_ == rhs.last_valid_wgstate_ &&
230 this->nupcol_wgstate_ == rhs.nupcol_wgstate_ &&
231 this->last_glift_opt_time_ == rhs.last_glift_opt_time_ &&
232 this->switched_prod_groups_ == rhs.switched_prod_groups_ &&
233 this->switched_inj_groups_ == rhs.switched_inj_groups_;
234 }
235
236protected:
237
238 /*
239 The dynamic state of the well model is maintained with an instance
240 of the WellState class. Currently we have
241 three different wellstate instances:
242
243 1. The currently active wellstate is in the active_well_state_
244 member. That is the state which is mutated by the simulator.
245
246 2. In the case timestep fails to converge and we must go back and
247 try again with a smaller timestep we need to recover the last
248 valid wellstate. This is maintained with the
249 last_valid_well_state_ member and the functions
250 commitWellState() and resetWellState().
251
252 3. For the NUPCOL functionality we should either use the
253 currently active wellstate or a wellstate frozen at max
254 nupcol iterations. This is handled with the member
255 nupcol_well_state_ and the initNupcolWellState() function.
256 */
257
258
259 /*
260 Will return the last good wellstate. This is typcially used when
261 initializing a new report step where the Schedule object might
262 have introduced new wells. The wellstate returned by
263 prevWellState() must have been stored with the commitWellState()
264 function first.
265 */
266 const WellState& prevWellState() const
267 {
268 return this->last_valid_wgstate_.well_state;
269 }
270
271 const WGState& prevWGState() const
272 {
273 return this->last_valid_wgstate_;
274 }
275 /*
276 Will return the currently active nupcolWellState; must initialize
277 the internal nupcol wellstate with initNupcolWellState() first.
278 */
279 const WellState& nupcolWellState() const
280 {
281 return this->nupcol_wgstate_.well_state;
282 }
283
284 /*
285 Will store a copy of the input argument well_state in the
286 last_valid_well_state_ member, that state can then be recovered
287 with a subsequent call to resetWellState().
288 */
289 void commitWGState(WGState wgstate)
290 {
291 this->last_valid_wgstate_ = std::move(wgstate);
292 }
293
294 /*
295 Will update the internal variable active_well_state_ to whatever
296 was stored in the last_valid_well_state_ member. This function
297 works in pair with commitWellState() which should be called first.
298 */
299 void resetWGState()
300 {
301 this->active_wgstate_ = this->last_valid_wgstate_;
302 }
303
304 /*
305 Will store the current active wellstate in the nupcol_well_state_
306 member. This can then be subsequently retrieved with accessor
307 nupcolWellState().
308 */
309 void updateNupcolWGState()
310 {
311 this->nupcol_wgstate_ = this->active_wgstate_;
312 }
313
316 std::vector<std::reference_wrapper<ParallelWellInfo>> createLocalParallelWellInfo(const std::vector<Well>& wells);
317
318 void initializeWellProdIndCalculators();
319 void initializeWellPerfData();
320
321 bool wasDynamicallyShutThisTimeStep(const int well_index) const;
322
323 double updateNetworkPressures(const int reportStepIdx);
324
325 void updateWsolvent(const Group& group,
326 const int reportStepIdx,
327 const WellState& wellState);
328 void setWsolvent(const Group& group,
329 const int reportStepIdx,
330 double wsolvent);
331 virtual void calcRates(const int fipnum,
332 const int pvtreg,
333 const std::vector<double>& production_rates,
334 std::vector<double>& resv_coeff) = 0;
335 virtual void calcInjRates(const int fipnum,
336 const int pvtreg,
337 std::vector<double>& resv_coeff) = 0;
338
339 void assignShutConnections(data::Wells& wsrpt,
340 const int reportStepIndex) const;
341 void assignGroupControl(const Group& group,
342 data::GroupData& gdata) const;
343 void assignGroupValues(const int reportStepIdx,
344 std::map<std::string, data::GroupData>& gvalues) const;
345 void assignNodeValues(std::map<std::string, data::NodeData>& nodevalues) const;
346
347 void calculateEfficiencyFactors(const int reportStepIdx);
348
349 void checkGconsaleLimits(const Group& group,
350 WellState& well_state,
351 const int reportStepIdx,
353
354 void checkGEconLimits(const Group& group,
355 const double simulation_time,
356 const int report_step_idx,
358
359 bool checkGroupHigherConstraints(const Group& group,
361 const int reportStepIdx);
362
363 void updateAndCommunicateGroupData(const int reportStepIdx,
364 const int iterationIdx);
365
366 void inferLocalShutWells();
367
368 void setRepRadiusPerfLength();
369
370 void gliftDebug(const std::string& msg,
372
373 void gliftDebugShowALQ(DeferredLogger& deferred_logger);
374
375 void gasLiftOptimizationStage2(DeferredLogger& deferred_logger,
376 GLiftProdWells& prod_wells,
377 GLiftOptWells& glift_wells,
379 GLiftWellStateMap& map,
380 const int episodeIndex);
381
382 virtual void computePotentials(const std::size_t widx,
384 std::string& exc_msg,
385 ExceptionType::ExcEnum& exc_type,
387
388 // Calculating well potentials for each well
389 void updateWellPotentials(const int reportStepIdx,
390 const bool onlyAfterEvent,
393
394 void initInjMult();
395
396
397 void updateInjMult(DeferredLogger& deferred_logger);
398 void updateInjFCMult(DeferredLogger& deferred_logger);
399
400 void updateFiltrationParticleVolume(const double dt, const std::size_t water_index);
401
402 // create the well container
403 virtual void createWellContainer(const int time_step) = 0;
404 virtual void initWellContainer(const int reportStepIdx) = 0;
405
406 virtual void calculateProductivityIndexValuesShutWells(const int reportStepIdx,
408 virtual void calculateProductivityIndexValues(DeferredLogger& deferred_logger) = 0;
409
410 void runWellPIScaling(const int timeStepIdx,
412
415
416 std::vector<int> getCellsForConnections(const Well& well) const;
417 std::vector<std::vector<int>> getMaxWellConnections() const;
418
419 std::vector<std::string> getWellsForTesting(const int timeStepIdx,
420 const double simulationTime);
421
422 using WellTracerRates = std::map<std::pair<std::string, std::string>, double>;
423 void assignWellTracerRates(data::Wells& wsrpt,
424 const WellTracerRates& wellTracerRates) const;
425
426 Schedule& schedule_;
427 const SummaryState& summaryState_;
428 const EclipseState& eclState_;
429 const Parallel::Communication& comm_;
430
431 PhaseUsage phase_usage_;
432 bool terminal_output_{false};
433 bool wells_active_{false};
434 bool initial_step_{};
435 bool report_step_starts_{};
436
437 std::optional<int> last_run_wellpi_{};
438
439 std::vector<Well> wells_ecl_;
440 std::vector<std::vector<PerforationData>> well_perf_data_;
441
444 {
445 public:
450 explicit ConnectionIndexMap(const std::size_t numConns)
451 : local_(numConns, -1)
452 {
453 this->global_.reserve(numConns);
454 this->open_.reserve(numConns);
455 }
456
465 const bool connIsOpen)
466 {
467 this->local_[connIdx] =
468 static_cast<int>(this->global_.size());
469
470 this->global_.push_back(connIdx);
471
472 const auto open_conn_idx = connIsOpen
473 ? this->num_open_conns_++
474 : -1;
475
476 this->open_.push_back(open_conn_idx);
477 }
478
484 const std::vector<int>& local() const
485 {
486 return this->local_;
487 }
488
494 int global(const int connIdx) const
495 {
496 return this->global_[connIdx];
497 }
498
506 int open(const int connIdx) const
507 {
508 return this->open_[connIdx];
509 }
510
511 private:
515 std::vector<int> local_{};
516
519 std::vector<int> global_{};
520
522 std::vector<int> open_{};
523
525 int num_open_conns_{0};
526 };
527
528 std::vector<ConnectionIndexMap> conn_idx_map_{};
529 std::function<bool(const Well&)> not_on_process_{};
530
531 // a vector of all the wells.
532 std::vector<WellInterfaceGeneric*> well_container_generic_{};
533
534 std::vector<int> local_shut_wells_{};
535
536 std::vector<ParallelWellInfo> parallel_well_info_;
537 std::vector<std::reference_wrapper<ParallelWellInfo>> local_parallel_well_info_;
538
539 std::vector<WellProdIndexCalculator> prod_index_calc_;
540 mutable ParallelWBPCalculation wbpCalculationService_;
541
542 std::vector<int> pvt_region_idx_;
543
544 mutable std::unordered_set<std::string> closed_this_step_;
545
546 GuideRate guideRate_;
547 std::unique_ptr<VFPProperties> vfp_properties_{};
548 std::map<std::string, double> node_pressures_; // Storing network pressures for output.
549
550 // previous injection multiplier, it is used in the injection multiplier calculation for WINJMULT keyword
551 std::unordered_map<std::string, std::vector<double>> prev_inj_multipliers_;
552
553 // Handling for filter cake injection multipliers
554 std::unordered_map<std::string, WellFilterCake> filter_cake_;
555
556 /*
557 The various wellState members should be accessed and modified
558 through the accessor functions wellState(), prevWellState(),
559 commitWellState(), resetWellState(), nupcolWellState() and
560 updateNupcolWellState().
561 */
562 WGState active_wgstate_;
563 WGState last_valid_wgstate_;
564 WGState nupcol_wgstate_;
565
566 bool glift_debug = false;
567
568 double last_glift_opt_time_ = -1.0;
569
570 std::map<std::string, std::string> switched_prod_groups_;
571 std::map<std::pair<std::string, Opm::Phase>, std::string> switched_inj_groups_;
572
573private:
574 WellInterfaceGeneric* getGenWell(const std::string& well_name);
575};
576
577
578} // namespace Opm
579
580#endif
Definition AquiferInterface.hpp:35
Connection index mappings.
Definition BlackoilWellModelGeneric.hpp:444
int open(const int connIdx) const
Get open connection ID of local (on-rank) connection.
Definition BlackoilWellModelGeneric.hpp:506
const std::vector< int > & local() const
Get local connection IDs/indices of every existing well connection.
Definition BlackoilWellModelGeneric.hpp:484
int global(const int connIdx) const
Get global connection ID of local (on-rank) connection.
Definition BlackoilWellModelGeneric.hpp:494
void addActiveConnection(const int connIdx, const bool connIsOpen)
Enumerate/map new active connection.
Definition BlackoilWellModelGeneric.hpp:464
ConnectionIndexMap(const std::size_t numConns)
Constructor.
Definition BlackoilWellModelGeneric.hpp:450
Class for handling the blackoil well model.
Definition BlackoilWellModelGeneric.hpp:82
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:261
bool needRebalanceNetwork(const int report_step) const
Whether it is necessary to re-balance network.
Definition BlackoilWellModelGeneric.cpp:983
bool wellsActive() const
return true if wells are available in the reservoir
Definition BlackoilWellModelGeneric.cpp:137
bool hasTHPConstraints() const
Return true if any well has a THP constraint.
Definition BlackoilWellModelGeneric.cpp:976
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:1007
Definition DeferredLogger.hpp:57
Definition GasLiftGroupInfo.hpp:50
Definition GroupState.hpp:34
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition WellState.hpp:60
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition BlackoilPhases.hpp:27
Definition BlackoilPhases.hpp:46
Definition WGState.hpp:37