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/Well/WellTestState.hpp>
29#include <opm/input/eclipse/Schedule/Group/GuideRate.hpp>
30
31#include <opm/simulators/utils/DeferredLoggingErrorHelpers.hpp>
32
33#include <opm/simulators/wells/PerforationData.hpp>
34#include <opm/simulators/wells/WellProdIndexCalculator.hpp>
35#include <opm/simulators/wells/WGState.hpp>
36
37#include <functional>
38#include <map>
39#include <memory>
40#include <string>
41#include <unordered_map>
42#include <unordered_set>
43#include <vector>
44
45namespace Opm {
46 class DeferredLogger;
47 class EclipseState;
48 class GasLiftSingleWellGeneric;
49 class GasLiftWellState;
50 class GasLiftGroupInfo;
51 class Group;
52 class GuideRateConfig;
53 class ParallelWellInfo;
54 class RestartValue;
55 class Schedule;
56 class SummaryConfig;
57 class VFPProperties;
58 class WellInterfaceGeneric;
59 class WellState;
60} // namespace Opm
61
62namespace Opm { namespace data {
63 struct GroupData;
64 struct GroupGuideRates;
65 class GroupAndNetworkValues;
66 struct NodeData;
67}} // namespace Opm::data
68
69namespace Opm {
70
73{
74public:
75 // --------- Types ---------
76 using GLiftOptWells = std::map<std::string, std::unique_ptr<GasLiftSingleWellGeneric>>;
77 using GLiftProdWells = std::map<std::string, const WellInterfaceGeneric*>;
78 using GLiftWellStateMap = std::map<std::string, std::unique_ptr<GasLiftWellState>>;
79
81 const SummaryState& summaryState,
82 const EclipseState& eclState,
84 const Parallel::Communication& comm);
85
86 virtual ~BlackoilWellModelGeneric() = default;
87
88 int numLocalWells() const;
89 int numPhases() const;
90
92 bool wellsActive() const;
93 bool hasWell(const std::string& wname) const;
94
95 // whether there exists any multisegment well open on this process
96 bool anyMSWellOpenLocal() const;
97
98 const Well& getWellEcl(const std::string& well_name) const;
99 std::vector<Well> getLocalWells(const int timeStepIdx) const;
100 const Schedule& schedule() const { return schedule_; }
101 const PhaseUsage& phaseUsage() const { return phase_usage_; }
102 const GroupState& groupState() const { return this->active_wgstate_.group_state; }
103 std::vector<const WellInterfaceGeneric*> genericWells() const
104 { return {well_container_generic_.begin(), well_container_generic_.end()}; }
105
106 /*
107 Immutable version of the currently active wellstate.
108 */
109 const WellState& wellState() const
110 {
111 return this->active_wgstate_.well_state;
112 }
113
114 /*
115 Mutable version of the currently active wellstate.
116 */
117 WellState& wellState()
118 {
119 return this->active_wgstate_.well_state;
120 }
121
122 GroupState& groupState() { return this->active_wgstate_.group_state; }
123
124 WellTestState& wellTestState() { return this->active_wgstate_.well_test_state; }
125
126 const WellTestState& wellTestState() const { return this->active_wgstate_.well_test_state; }
127
128
129 double wellPI(const int well_index) const;
130 double wellPI(const std::string& well_name) const;
131
132 void updateEclWells(const int timeStepIdx,
133 const std::unordered_set<std::string>& wells,
134 const SummaryState& st);
135
136 void initFromRestartFile(const RestartValue& restartValues,
138 const size_t numCells,
139 bool handle_ms_well);
140
141 void prepareDeserialize(int report_step,
142 const size_t numCells,
143 bool handle_ms_well);
144
145 /*
146 Will assign the internal member last_valid_well_state_ to the
147 current value of the this->active_well_state_. The state stored
148 with storeWellState() can then subsequently be recovered with the
149 resetWellState() method.
150 */
151 void commitWGState()
152 {
153 this->last_valid_wgstate_ = this->active_wgstate_;
154 }
155
156 data::GroupAndNetworkValues groupAndNetworkData(const int reportStepIdx) const;
157
159 bool hasTHPConstraints() const;
160
163 bool forceShutWellByName(const std::string& wellname,
164 const double simulation_time);
165
166 const std::vector<PerforationData>& perfData(const int well_idx) const
167 { return well_perf_data_[well_idx]; }
168
169 const Parallel::Communication& comm() const { return comm_; }
170
171 const SummaryState& summaryState() const { return summaryState_; }
172
173 const GuideRate& guideRate() const { return guideRate_; }
174
175 bool reportStepStarts() const { return report_step_starts_; }
176
177 bool shouldBalanceNetwork(const int reportStepIndex,
178 const int iterationIdx) const;
179
180 bool shouldIterateNetwork(const int reportStepIndex,
181 const std::size_t recursion_level,
182 const double network_imbalance) const;
183
184 template<class Serializer>
185 void serializeOp(Serializer& serializer)
186 {
187 serializer(initial_step_);
188 serializer(report_step_starts_);
189 serializer(last_run_wellpi_);
190 serializer(local_shut_wells_);
191 serializer(closed_this_step_);
192 serializer(guideRate_);
193 serializer(node_pressures_);
194 serializer(active_wgstate_);
195 serializer(last_valid_wgstate_);
196 serializer(nupcol_wgstate_);
197 serializer(last_glift_opt_time_);
198 serializer(switched_prod_groups_);
199 serializer(switched_inj_groups_);
200 }
201
202 bool operator==(const BlackoilWellModelGeneric& rhs) const
203 {
204 return this->initial_step_ == rhs.initial_step_ &&
205 this->report_step_starts_ == rhs.report_step_starts_ &&
206 this->last_run_wellpi_ == rhs.last_run_wellpi_ &&
207 this->local_shut_wells_ == rhs.local_shut_wells_ &&
208 this->closed_this_step_ == rhs.closed_this_step_ &&
209 this->node_pressures_ == rhs.node_pressures_ &&
210 this->active_wgstate_ == rhs.active_wgstate_ &&
211 this->last_valid_wgstate_ == rhs.last_valid_wgstate_ &&
212 this->nupcol_wgstate_ == rhs.nupcol_wgstate_ &&
213 this->last_glift_opt_time_ == rhs.last_glift_opt_time_ &&
214 this->switched_prod_groups_ == rhs.switched_prod_groups_ &&
215 this->switched_inj_groups_ == rhs.switched_inj_groups_;
216 }
217
218protected:
219
220 /*
221 The dynamic state of the well model is maintained with an instance
222 of the WellState class. Currently we have
223 three different wellstate instances:
224
225 1. The currently active wellstate is in the active_well_state_
226 member. That is the state which is mutated by the simulator.
227
228 2. In the case timestep fails to converge and we must go back and
229 try again with a smaller timestep we need to recover the last
230 valid wellstate. This is maintained with the
231 last_valid_well_state_ member and the functions
232 commitWellState() and resetWellState().
233
234 3. For the NUPCOL functionality we should either use the
235 currently active wellstate or a wellstate frozen at max
236 nupcol iterations. This is handled with the member
237 nupcol_well_state_ and the initNupcolWellState() function.
238 */
239
240
241 /*
242 Will return the last good wellstate. This is typcially used when
243 initializing a new report step where the Schedule object might
244 have introduced new wells. The wellstate returned by
245 prevWellState() must have been stored with the commitWellState()
246 function first.
247 */
248 const WellState& prevWellState() const
249 {
250 return this->last_valid_wgstate_.well_state;
251 }
252
253 const WGState& prevWGState() const
254 {
255 return this->last_valid_wgstate_;
256 }
257 /*
258 Will return the currently active nupcolWellState; must initialize
259 the internal nupcol wellstate with initNupcolWellState() first.
260 */
261 const WellState& nupcolWellState() const
262 {
263 return this->nupcol_wgstate_.well_state;
264 }
265
266 /*
267 Will store a copy of the input argument well_state in the
268 last_valid_well_state_ member, that state can then be recovered
269 with a subsequent call to resetWellState().
270 */
271 void commitWGState(WGState wgstate)
272 {
273 this->last_valid_wgstate_ = std::move(wgstate);
274 }
275
276 /*
277 Will update the internal variable active_well_state_ to whatever
278 was stored in the last_valid_well_state_ member. This function
279 works in pair with commitWellState() which should be called first.
280 */
281 void resetWGState()
282 {
283 this->active_wgstate_ = this->last_valid_wgstate_;
284 }
285
286 /*
287 Will store the current active wellstate in the nupcol_well_state_
288 member. This can then be subsequently retrieved with accessor
289 nupcolWellState().
290 */
291 void updateNupcolWGState()
292 {
293 this->nupcol_wgstate_ = this->active_wgstate_;
294 }
295
298 std::vector<std::reference_wrapper<ParallelWellInfo>> createLocalParallelWellInfo(const std::vector<Well>& wells);
299
300 void initializeWellProdIndCalculators();
301 void initializeWellPerfData();
302
303 bool wasDynamicallyShutThisTimeStep(const int well_index) const;
304
305 std::pair<bool, double> updateNetworkPressures(const int reportStepIdx);
306
307 void updateWsolvent(const Group& group,
308 const int reportStepIdx,
309 const WellState& wellState);
310 void setWsolvent(const Group& group,
311 const int reportStepIdx,
312 double wsolvent);
313 virtual void calcRates(const int fipnum,
314 const int pvtreg,
315 const std::vector<double>& production_rates,
316 std::vector<double>& resv_coeff) = 0;
317 virtual void calcInjRates(const int fipnum,
318 const int pvtreg,
319 std::vector<double>& resv_coeff) = 0;
320
321 void assignShutConnections(data::Wells& wsrpt,
322 const int reportStepIndex) const;
323 void assignGroupControl(const Group& group,
324 data::GroupData& gdata) const;
325 void assignGroupValues(const int reportStepIdx,
326 std::map<std::string, data::GroupData>& gvalues) const;
327 void assignNodeValues(std::map<std::string, data::NodeData>& nodevalues) const;
328
329 void calculateEfficiencyFactors(const int reportStepIdx);
330
331 void checkGconsaleLimits(const Group& group,
332 WellState& well_state,
333 const int reportStepIdx,
335
336 bool checkGroupHigherConstraints(const Group& group,
338 const int reportStepIdx);
339
340 void updateAndCommunicateGroupData(const int reportStepIdx,
341 const int iterationIdx);
342
343 void inferLocalShutWells();
344
345 void setRepRadiusPerfLength();
346
347 void gliftDebug(const std::string& msg,
349
350 void gliftDebugShowALQ(DeferredLogger& deferred_logger);
351
352 void gasLiftOptimizationStage2(DeferredLogger& deferred_logger,
353 GLiftProdWells& prod_wells,
354 GLiftOptWells& glift_wells,
356 GLiftWellStateMap& map,
357 const int episodeIndex);
358
359 virtual void computePotentials(const std::size_t widx,
361 std::string& exc_msg,
362 ExceptionType::ExcEnum& exc_type,
364
365 // Calculating well potentials for each well
366 void updateWellPotentials(const int reportStepIdx,
367 const bool onlyAfterEvent,
370
371 // create the well container
372 virtual void createWellContainer(const int time_step) = 0;
373 virtual void initWellContainer(const int reportStepIdx) = 0;
374
375 virtual void calculateProductivityIndexValuesShutWells(const int reportStepIdx,
377 virtual void calculateProductivityIndexValues(DeferredLogger& deferred_logger) = 0;
378
379 void runWellPIScaling(const int timeStepIdx,
381
384
385 std::vector<int> getCellsForConnections(const Well& well) const;
386
387 std::vector<std::string> getWellsForTesting(const int timeStepIdx,
388 const double simulationTime);
389
390 Schedule& schedule_;
391 const SummaryState& summaryState_;
392 const EclipseState& eclState_;
393 const Parallel::Communication& comm_;
394
395 PhaseUsage phase_usage_;
396 bool terminal_output_{false};
397 bool wells_active_{false};
398 bool initial_step_{};
399 bool report_step_starts_{};
400
401 std::optional<int> last_run_wellpi_{};
402
403 std::vector<Well> wells_ecl_;
404 std::vector<std::vector<PerforationData>> well_perf_data_;
405 std::function<bool(const Well&)> not_on_process_{};
406
407 // a vector of all the wells.
408 std::vector<WellInterfaceGeneric*> well_container_generic_{};
409
410 std::vector<int> local_shut_wells_{};
411
412 std::vector<ParallelWellInfo> parallel_well_info_;
413 std::vector<std::reference_wrapper<ParallelWellInfo>> local_parallel_well_info_;
414
415 std::vector<WellProdIndexCalculator> prod_index_calc_;
416
417 std::vector<int> pvt_region_idx_;
418
419 mutable std::unordered_set<std::string> closed_this_step_;
420
421 GuideRate guideRate_;
422 std::unique_ptr<VFPProperties> vfp_properties_{};
423 std::map<std::string, double> node_pressures_; // Storing network pressures for output.
424
425 /*
426 The various wellState members should be accessed and modified
427 through the accessor functions wellState(), prevWellState(),
428 commitWellState(), resetWellState(), nupcolWellState() and
429 updateNupcolWellState().
430 */
431 WGState active_wgstate_;
432 WGState last_valid_wgstate_;
433 WGState nupcol_wgstate_;
434
435 bool glift_debug = false;
436
437 double last_glift_opt_time_ = -1.0;
438
439 std::map<std::string, std::string> switched_prod_groups_;
440 std::map<std::pair<std::string, Opm::Phase>, std::string> switched_inj_groups_;
441
442
443private:
444 WellInterfaceGeneric* getGenWell(const std::string& well_name);
445};
446
447
448} // namespace Opm
449
450#endif
Definition AquiferInterface.hpp:35
Class for handling the blackoil well model.
Definition BlackoilWellModelGeneric.hpp:73
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:255
bool wellsActive() const
return true if wells are available in the reservoir
Definition BlackoilWellModelGeneric.cpp:131
bool hasTHPConstraints() const
Return true if any well has a THP constraint.
Definition BlackoilWellModelGeneric.cpp:903
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:910
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