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 GasLiftGroupInfo;
58 class GasLiftSingleWellGeneric;
59 class GasLiftWellState;
60 class Group;
61 class GuideRateConfig;
62 class ParallelWellInfo;
63 class RestartValue;
64 class Schedule;
65 struct SimulatorUpdate;
66 class SummaryConfig;
67 class VFPProperties;
68 class WellInterfaceGeneric;
69 class WellState;
70} // namespace Opm
71
72namespace Opm { namespace data {
73 struct GroupData;
74 struct GroupGuideRates;
75 class GroupAndNetworkValues;
76 struct NodeData;
77}} // namespace Opm::data
78
79namespace Opm {
80
83{
84public:
85 // --------- Types ---------
86 using GLiftOptWells = std::map<std::string, std::unique_ptr<GasLiftSingleWellGeneric>>;
87 using GLiftProdWells = std::map<std::string, const WellInterfaceGeneric*>;
88 using GLiftWellStateMap = std::map<std::string, std::unique_ptr<GasLiftWellState>>;
89
90 BlackoilWellModelGeneric(Schedule& schedule,
91 const SummaryState& summaryState,
92 const EclipseState& eclState,
93 const PhaseUsage& phase_usage,
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& groupState() const { return this->active_wgstate_.group_state; }
118 std::vector<const WellInterfaceGeneric*> 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& wellState() const
125 {
126 return this->active_wgstate_.well_state;
127 }
128
129 /*
130 Mutable version of the currently active wellstate.
131 */
132 WellState& 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& nupcolWellState() const
142 {
143 return this->nupcol_wgstate_.well_state;
144 }
145 GroupState& 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
152 double wellPI(const int well_index) const;
153 double wellPI(const std::string& well_name) const;
154
155 void updateEclWells(const int timeStepIdx,
156 const SimulatorUpdate& sim_update,
157 const SummaryState& st);
158
159 void initFromRestartFile(const RestartValue& restartValues,
160 WellTestState wtestState,
161 const std::size_t numCells,
162 bool handle_ms_well);
163
164 void prepareDeserialize(int report_step,
165 const std::size_t numCells,
166 bool handle_ms_well);
167
168 /*
169 Will assign the internal member last_valid_well_state_ to the
170 current value of the this->active_well_state_. The state stored
171 with storeWellState() can then subsequently be recovered with the
172 resetWellState() method.
173 */
174 void commitWGState()
175 {
176 this->last_valid_wgstate_ = this->active_wgstate_;
177 }
178
179 data::GroupAndNetworkValues groupAndNetworkData(const int reportStepIdx) const;
180
182 bool hasTHPConstraints() const;
183
185 void updateNetworkActiveState(const int report_step);
186
190 bool needPreStepNetworkRebalance(const int report_step) const;
191
194 bool forceShutWellByName(const std::string& wellname,
195 const double simulation_time);
196
197 const std::vector<PerforationData>& perfData(const int well_idx) const
198 { return well_perf_data_[well_idx]; }
199
200 const Parallel::Communication& comm() const { return comm_; }
201
202 const EclipseState& eclipseState() const { return eclState_; }
203
204 const SummaryState& summaryState() const { return summaryState_; }
205
206 const GuideRate& guideRate() const { return guideRate_; }
207
208 bool reportStepStarts() const { return report_step_starts_; }
209
210 bool shouldBalanceNetwork(const int reportStepIndex,
211 const int iterationIdx) const;
212
213 void updateClosedWellsThisStep(const std::string& well_name) const {
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 /*
259 The dynamic state of the well model is maintained with an instance
260 of the WellState class. Currently we have
261 three different wellstate instances:
262
263 1. The currently active wellstate is in the active_well_state_
264 member. That is the state which is mutated by the simulator.
265
266 2. In the case timestep fails to converge and we must go back and
267 try again with a smaller timestep we need to recover the last
268 valid wellstate. This is maintained with the
269 last_valid_well_state_ member and the functions
270 commitWellState() and resetWellState().
271
272 3. For the NUPCOL functionality we should either use the
273 currently active wellstate or a wellstate frozen at max
274 nupcol iterations. This is handled with the member
275 nupcol_well_state_ and the initNupcolWellState() function.
276 */
277
278
279 /*
280 Will return the last good wellstate. This is typcially used when
281 initializing a new report step where the Schedule object might
282 have introduced new wells. The wellstate returned by
283 prevWellState() must have been stored with the commitWellState()
284 function first.
285 */
286 const WellState& prevWellState() const
287 {
288 return this->last_valid_wgstate_.well_state;
289 }
290
291
292 const WGState& prevWGState() const
293 {
294 return this->last_valid_wgstate_;
295 }
296
297
298
299 /*
300 Will store a copy of the input argument well_state in the
301 last_valid_well_state_ member, that state can then be recovered
302 with a subsequent call to resetWellState().
303 */
304 void commitWGState(WGState wgstate)
305 {
306 this->last_valid_wgstate_ = std::move(wgstate);
307 }
308
309 /*
310 Will update the internal variable active_well_state_ to whatever
311 was stored in the last_valid_well_state_ member. This function
312 works in pair with commitWellState() which should be called first.
313 */
314 void resetWGState()
315 {
316 this->active_wgstate_ = this->last_valid_wgstate_;
317 }
318
319 /*
320 Will store the current active wellstate in the nupcol_well_state_
321 member. This can then be subsequently retrieved with accessor
322 nupcolWellState().
323 */
324 void updateNupcolWGState()
325 {
326 this->nupcol_wgstate_ = this->active_wgstate_;
327 }
328
331 std::vector<std::reference_wrapper<ParallelWellInfo>> createLocalParallelWellInfo(const std::vector<Well>& wells);
332
333 void initializeWellProdIndCalculators();
334 void initializeWellPerfData();
335
336 bool wasDynamicallyShutThisTimeStep(const int well_index) const;
337
338 double updateNetworkPressures(const int reportStepIdx);
339
340 void updateWsolvent(const Group& group,
341 const int reportStepIdx,
342 const WellState& wellState);
343 void setWsolvent(const Group& group,
344 const int reportStepIdx,
345 double wsolvent);
346 virtual void calcRates(const int fipnum,
347 const int pvtreg,
348 const std::vector<double>& production_rates,
349 std::vector<double>& resv_coeff) = 0;
350 virtual void calcInjRates(const int fipnum,
351 const int pvtreg,
352 std::vector<double>& resv_coeff) = 0;
353
354 void assignShutConnections(data::Wells& wsrpt,
355 const int reportStepIndex) const;
356 void assignGroupControl(const Group& group,
357 data::GroupData& gdata) const;
358 void assignGroupValues(const int reportStepIdx,
359 std::map<std::string, data::GroupData>& gvalues) const;
360 void assignNodeValues(std::map<std::string, data::NodeData>& nodevalues, const int reportStepIdx) const;
361
362 void calculateEfficiencyFactors(const int reportStepIdx);
363
364 void checkGconsaleLimits(const Group& group,
365 WellState& well_state,
366 const int reportStepIdx,
367 DeferredLogger& deferred_logger);
368
369 void checkGEconLimits(const Group& group,
370 const double simulation_time,
371 const int report_step_idx,
372 DeferredLogger& deferred_logger);
373
374 bool checkGroupHigherConstraints(const Group& group,
375 DeferredLogger& deferred_logger,
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,
386 DeferredLogger& deferred_logger) const;
387
388 void gliftDebugShowALQ(DeferredLogger& deferred_logger);
389
390 void gasLiftOptimizationStage2(DeferredLogger& deferred_logger,
391 GLiftProdWells& prod_wells,
392 GLiftOptWells& glift_wells,
393 GasLiftGroupInfo& group_info,
394 GLiftWellStateMap& map,
395 const int episodeIndex);
396
397 virtual void computePotentials(const std::size_t widx,
398 const WellState& well_state_copy,
399 std::string& exc_msg,
400 ExceptionType::ExcEnum& exc_type,
401 DeferredLogger& deferred_logger) = 0;
402
403 // Calculating well potentials for each well
404 void updateWellPotentials(const int reportStepIdx,
405 const bool onlyAfterEvent,
406 const SummaryConfig& summaryConfig,
407 DeferredLogger& deferred_logger);
408
409 void initInjMult();
410
411
412 void updateInjMult(DeferredLogger& deferred_logger);
413 void updateInjFCMult(DeferredLogger& deferred_logger);
414
415 void updateFiltrationParticleVolume(const double dt, const std::size_t water_index);
416
417 // create the well container
418 virtual void createWellContainer(const int time_step) = 0;
419 virtual void initWellContainer(const int reportStepIdx) = 0;
420
421 virtual void calculateProductivityIndexValuesShutWells(const int reportStepIdx,
422 DeferredLogger& deferred_logger) = 0;
423 virtual void calculateProductivityIndexValues(DeferredLogger& deferred_logger) = 0;
424
425 void runWellPIScaling(const int reportStepIdx,
426 DeferredLogger& local_deferredLogger);
427
429 virtual int compressedIndexForInterior(int cartesian_cell_idx) const = 0;
430
431 std::vector<int> getCellsForConnections(const Well& well) const;
432 std::vector<std::vector<int>> getMaxWellConnections() const;
433
434 std::vector<std::string> getWellsForTesting(const int timeStepIdx,
435 const double simulationTime);
436
437 using WellTracerRates = std::map<std::pair<std::string, std::string>, double>;
438 void assignWellTracerRates(data::Wells& wsrpt,
439 const WellTracerRates& wellTracerRates) const;
440
441 Schedule& schedule_;
442 const SummaryState& summaryState_;
443 const EclipseState& eclState_;
444 const Parallel::Communication& comm_;
445
446 PhaseUsage phase_usage_;
447 bool terminal_output_{false};
448 bool wells_active_{false};
449 bool network_active_{false};
450 bool initial_step_{};
451 bool report_step_starts_{};
452
453 std::optional<int> last_run_wellpi_{};
454
455 std::vector<Well> wells_ecl_;
456 std::vector<std::vector<PerforationData>> well_perf_data_;
457
460 {
461 public:
466 explicit ConnectionIndexMap(const std::size_t numConns)
467 : local_(numConns, -1)
468 {
469 this->global_.reserve(numConns);
470 this->open_.reserve(numConns);
471 }
472
480 void addActiveConnection(const int connIdx,
481 const bool connIsOpen)
482 {
483 this->local_[connIdx] =
484 static_cast<int>(this->global_.size());
485
486 this->global_.push_back(connIdx);
487
488 const auto open_conn_idx = connIsOpen
489 ? this->num_open_conns_++
490 : -1;
491
492 this->open_.push_back(open_conn_idx);
493 }
494
500 const std::vector<int>& local() const
501 {
502 return this->local_;
503 }
504
510 int global(const int connIdx) const
511 {
512 return this->global_[connIdx];
513 }
514
522 int open(const int connIdx) const
523 {
524 return this->open_[connIdx];
525 }
526
527 private:
531 std::vector<int> local_{};
532
535 std::vector<int> global_{};
536
538 std::vector<int> open_{};
539
541 int num_open_conns_{0};
542 };
543
544 std::vector<ConnectionIndexMap> conn_idx_map_{};
545 std::function<bool(const Well&)> not_on_process_{};
546
547 // a vector of all the wells.
548 std::vector<WellInterfaceGeneric*> well_container_generic_{};
549
550 std::vector<int> local_shut_wells_{};
551
552 std::vector<ParallelWellInfo> parallel_well_info_;
553 std::vector<std::reference_wrapper<ParallelWellInfo>> local_parallel_well_info_;
554
555 std::vector<WellProdIndexCalculator> prod_index_calc_;
556 mutable ParallelWBPCalculation wbpCalculationService_;
557
558 std::vector<int> pvt_region_idx_;
559
560 mutable std::unordered_set<std::string> closed_this_step_;
561
562 GuideRate guideRate_;
563 std::unique_ptr<VFPProperties> vfp_properties_{};
564 std::map<std::string, double> node_pressures_; // Storing network pressures for output.
565
566 // previous injection multiplier, it is used in the injection multiplier calculation for WINJMULT keyword
567 std::unordered_map<std::string, std::vector<double>> prev_inj_multipliers_;
568
569 // Handling for filter cake injection multipliers
570 std::unordered_map<std::string, WellFilterCake> filter_cake_;
571
572 /*
573 The various wellState members should be accessed and modified
574 through the accessor functions wellState(), prevWellState(),
575 commitWellState(), resetWellState(), nupcolWellState() and
576 updateNupcolWellState().
577 */
578 WGState active_wgstate_;
579 WGState last_valid_wgstate_;
580 WGState nupcol_wgstate_;
581
582 bool glift_debug = false;
583
584 double last_glift_opt_time_ = -1.0;
585
586 bool wellStructureChangedDynamically_{false};
587
588 // Store maps of group name and new group controls for output
589 std::map<std::string, std::string> switched_prod_groups_;
590 std::map<std::pair<std::string, Opm::Phase>, std::string> switched_inj_groups_;
591 // Store map of group name and close offending well for output
592 std::map<std::string, std::pair<std::string, std::string>> closed_offending_wells_;
593
594
595private:
596 WellInterfaceGeneric* getGenWell(const std::string& well_name);
597};
598
599
600} // namespace Opm
601
602#endif
Connection index mappings.
Definition BlackoilWellModelGeneric.hpp:460
int open(const int connIdx) const
Get open connection ID of local (on-rank) connection.
Definition BlackoilWellModelGeneric.hpp:522
const std::vector< int > & local() const
Get local connection IDs/indices of every existing well connection.
Definition BlackoilWellModelGeneric.hpp:500
int global(const int connIdx) const
Get global connection ID of local (on-rank) connection.
Definition BlackoilWellModelGeneric.hpp:510
void addActiveConnection(const int connIdx, const bool connIsOpen)
Enumerate/map new active connection.
Definition BlackoilWellModelGeneric.hpp:480
ConnectionIndexMap(const std::size_t numConns)
Constructor.
Definition BlackoilWellModelGeneric.hpp:466
Class for handling the blackoil well model.
Definition BlackoilWellModelGeneric.hpp:83
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:272
bool networkActive() const
return true if network is active (at least one network well in prediction mode)
Definition BlackoilWellModelGeneric.cpp:152
void updateNetworkActiveState(const int report_step)
Checks if network is active (at least one network well on prediction).
Definition BlackoilWellModelGeneric.cpp:1051
bool needPreStepNetworkRebalance(const int report_step) const
Checks if there are reasons to perform a pre-step network re-balance.
Definition BlackoilWellModelGeneric.cpp:1072
bool wellsActive() const
return true if wells are available in the reservoir
Definition BlackoilWellModelGeneric.cpp:145
bool hasTHPConstraints() const
Return true if any well has a THP constraint.
Definition BlackoilWellModelGeneric.cpp:1044
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:1091
Definition DeferredLogger.hpp:57
Definition GasLiftGroupInfo.hpp:45
Definition GroupState.hpp:34
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition WellState.hpp:61
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition BlackoilPhases.hpp:27
Definition BlackoilPhases.hpp:46
Definition WGState.hpp:37