My Project
Loading...
Searching...
No Matches
WellState.hpp
1/*
2 Copyright 2014 SINTEF ICT, Applied Mathematics.
3 Copyright 2017 IRIS AS
4
5 This file is part of the Open Porous Media project (OPM).
6
7 OPM is free software: you can redistribute it and/or modify
8 it under the terms of the GNU General Public License as published by
9 the Free Software Foundation, either version 3 of the License, or
10 (at your option) any later version.
11
12 OPM is distributed in the hope that it will be useful,
13 but WITHOUT ANY WARRANTY; without even the implied warranty of
14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 GNU General Public License for more details.
16
17 You should have received a copy of the GNU General Public License
18 along with OPM. If not, see <http://www.gnu.org/licenses/>.
19*/
20
21#ifndef OPM_WELLSTATEFULLYIMPLICITBLACKOIL_HEADER_INCLUDED
22#define OPM_WELLSTATEFULLYIMPLICITBLACKOIL_HEADER_INCLUDED
23
24#include <dune/common/version.hh>
25#include <dune/common/parallel/mpihelper.hh>
26
27#include <opm/common/ErrorMacros.hpp>
28
29#include <opm/input/eclipse/Schedule/Events.hpp>
30
31#include <opm/output/data/Wells.hpp>
32
33#include <opm/simulators/wells/ALQState.hpp>
34#include <opm/simulators/wells/GlobalWellInfo.hpp>
35#include <opm/simulators/wells/PerfData.hpp>
36#include <opm/simulators/wells/SegmentState.hpp>
37#include <opm/simulators/wells/SingleWellState.hpp>
38#include <opm/simulators/wells/WellContainer.hpp>
39
40#include <opm/simulators/utils/BlackoilPhases.hpp>
41#include <opm/simulators/utils/ParallelCommunication.hpp>
42
43#include <algorithm>
44#include <cstddef>
45#include <cstdint>
46#include <functional>
47#include <map>
48#include <optional>
49#include <string>
50#include <utility>
51#include <vector>
52
53namespace Opm
54{
55
56template<class Scalar> class ParallelWellInfo;
57template<class Scalar> struct PerforationData;
58template<class Scalar> class ConnFracStatistics;
59class Schedule;
60enum class WellStatus;
61
64template<class Scalar>
66{
67public:
68 static const std::uint64_t event_mask = ScheduleEvents::WELL_STATUS_CHANGE
69 | ScheduleEvents::PRODUCTION_UPDATE
70 | ScheduleEvents::INJECTION_UPDATE;
71
72 // TODO: same definition with WellInterface, eventually they should go to a common header file.
73 static const int Water = BlackoilPhases::Aqua;
74 static const int Oil = BlackoilPhases::Liquid;
75 static const int Gas = BlackoilPhases::Vapour;
76
77 // Only usable for testing purposes
78 explicit WellState(const ParallelWellInfo<Scalar>& pinfo);
79
80 explicit WellState(const PhaseUsage& pu)
81 : phase_usage_(pu)
82 {}
83
84 static WellState serializationTestObject(const ParallelWellInfo<Scalar>& pinfo);
85
86 std::size_t size() const
87 {
88 return this->wells_.size();
89 }
90
91 std::vector<std::string> wells() const
92 {
93 return this->wells_.wells();
94 }
95
96 int numWells() const
97 {
98 return this->size();
99 }
100
101 const ParallelWellInfo<Scalar>& parallelWellInfo(std::size_t well_index) const;
102
106 void init(const std::vector<Scalar>& cellPressures,
107 const std::vector<Scalar>& cellTemperatures,
108 const Schedule& schedule,
109 const std::vector<Well>& wells_ecl,
110 const std::vector<std::reference_wrapper<ParallelWellInfo<Scalar>>>& parallel_well_info,
111 const int report_step,
112 const WellState* prevState,
113 const std::vector<std::vector<PerforationData<Scalar>>>& well_perf_data,
115 const bool enableDistributedWells);
116
117 void resize(const std::vector<Well>& wells_ecl,
118 const std::vector<std::reference_wrapper<ParallelWellInfo<Scalar>>>& parallel_well_info,
119 const Schedule& schedule,
120 const bool handle_ms_well,
121 const std::size_t numCells,
122 const std::vector<std::vector<PerforationData<Scalar>>>& well_perf_data,
124
125 void setCurrentWellRates(const std::string& wellName,
126 const std::vector<Scalar>& new_rates)
127 {
128 auto& [owner, rates] = this->well_rates.at(wellName);
129 if (owner)
130 rates = new_rates;
131 }
132
133 const std::vector<Scalar>& currentWellRates(const std::string& wellName) const;
134
135 bool hasWellRates(const std::string& wellName) const
136 {
137 return this->well_rates.find(wellName) != this->well_rates.end();
138 }
139
140 void clearWellRates()
141 {
142 this->well_rates.clear();
143 }
144
145 void gatherVectorsOnRoot(const std::vector<data::Connection>& from_connections,
146 std::vector<data::Connection>& to_connections,
147 const Parallel::Communication& comm) const;
148
149 data::Wells
150 report(const int* globalCellIdxMap,
151 const std::function<bool(const int)>& wasDynamicallyClosed) const;
152
153 void reportConnections(std::vector<data::Connection>& connections,
154 const PhaseUsage &pu,
155 std::size_t well_index,
156 const int* globalCellIdxMap) const;
157
159 void initWellStateMSWell(const std::vector<Well>& wells_ecl,
161
162 static void calculateSegmentRates(const ParallelWellInfo<Scalar>& pw_info,
163 const std::vector<std::vector<int>>& segment_inlets,
164 const std::vector<std::vector<int>>& segment_perforations,
165 const std::vector<Scalar>& perforation_rates,
166 const int np,
167 const int segment,
168 std::vector<Scalar>& segment_rates);
169
170
171 void communicateGroupRates(const Parallel::Communication& comm);
172
173 void updateGlobalIsGrup(const Parallel::Communication& comm);
174 void updateEfficiencyScalingFactor(const std::string& wellName,
175 const Scalar value);
176
177 bool isInjectionGrup(const std::string& name) const
178 {
179 return this->global_well_info.value().in_injecting_group(name);
180 }
181
182 bool isProductionGrup(const std::string& name) const
183 {
184 return this->global_well_info.value().in_producing_group(name);
185 }
186
187 bool isOpen(const std::string& name) const
188 {
189 return this->global_well_info.value().is_open(name);
190 }
191
192 Scalar getGlobalEfficiencyScalingFactor(const std::string& name) const
193 {
194 return this->global_well_info.value().efficiency_scaling_factor(name);
195 }
196
197 Scalar getALQ(const std::string& name) const
198 {
199 return this->alq_state.get(name);
200 }
201
202 void setALQ(const std::string& name, Scalar value)
203 {
204 this->alq_state.set(name, value);
205 }
206
207 int gliftGetDebugCounter()
208 {
209 return this->alq_state.get_debug_counter();
210 }
211
212 void gliftSetDebugCounter(int value)
213 {
214 return this->alq_state.set_debug_counter(value);
215 }
216
217 int gliftUpdateDebugCounter()
218 {
219 return this->alq_state.update_debug_counter();
220 }
221
222 bool gliftCheckAlqOscillation(const std::string &name) const
223 {
224 return this->alq_state.oscillation(name);
225 }
226
227 int gliftGetAlqDecreaseCount(const std::string &name)
228 {
229 return this->alq_state.get_decrement_count(name);
230 }
231
232 int gliftGetAlqIncreaseCount(const std::string &name)
233 {
234 return this->alq_state.get_increment_count(name);
235 }
236
237 void gliftUpdateAlqIncreaseCount(const std::string &name, bool increase)
238 {
239 this->alq_state.update_count(name, increase);
240 }
241
242 void gliftTimeStepInit()
243 {
244 this->alq_state.reset_count();
245 }
246
247 // If the ALQ has changed since the previous time step,
248 // reset current_alq and update default_alq. ALQ is used for
249 // constant lift gas injection and for gas lift optimization
250 // (THP controlled wells).
251 void updateWellsDefaultALQ(const Schedule& schedule,
252 const int report_step,
254
255 int wellNameToGlobalIdx(const std::string& name)
256 {
257 return this->global_well_info.value().well_index(name);
258 }
259
260 std::string globalIdxToWellName(const int index)
261 {
262 return this->global_well_info.value().well_name(index);
263 }
264
265 bool wellIsOwned(std::size_t well_index,
266 const std::string& wellName) const;
267
268 bool wellIsOwned(const std::string& wellName) const;
269
270 bool isRank0() const {
271 return this->global_well_info.value().isRank0();
272 }
273
274 void updateStatus(int well_index, WellStatus status);
275
276 void openWell(int well_index);
277 void shutWell(int well_index);
278 void stopWell(int well_index);
279
280 void switchToProducer(const std::string& name) {
281 alq_state.insert(name);
282 }
283
285 int numPhases() const
286 {
287 return this->phase_usage_.num_phases;
288 }
289
290 const PhaseUsage& phaseUsage() const
291 {
292 return this->phase_usage_;
293 }
294
296 std::vector<Scalar>& wellRates(std::size_t well_index)
297 { return this->wells_[well_index].surface_rates; }
298 const std::vector<Scalar>& wellRates(std::size_t well_index) const
299 { return this->wells_[well_index].surface_rates; }
300
301 const std::string& name(std::size_t well_index) const
302 {
303 return this->wells_.well_name(well_index);
304 }
305
306 std::optional<std::size_t> index(const std::string& well_name) const
307 {
308 return this->wells_.well_index(well_name);
309 }
310
311 const SingleWellState<Scalar>& operator[](std::size_t well_index) const
312 {
313 return this->wells_[well_index];
314 }
315
316 const SingleWellState<Scalar>& operator[](const std::string& well_name) const
317 {
318 return this->wells_[well_name];
319 }
320
321 SingleWellState<Scalar>& operator[](std::size_t well_index)
322 {
323 return this->wells_[well_index];
324 }
325
326 SingleWellState<Scalar>& operator[](const std::string& well_name)
327 {
328 return this->wells_[well_name];
329 }
330
331 const SingleWellState<Scalar>& well(std::size_t well_index) const
332 {
333 return this->operator[](well_index);
334 }
335
336 const SingleWellState<Scalar>& well(const std::string& well_name) const
337 {
338 return this->operator[](well_name);
339 }
340
341 SingleWellState<Scalar>& well(std::size_t well_index)
342 {
343 return this->operator[](well_index);
344 }
345
346 SingleWellState<Scalar>& well(const std::string& well_name)
347 {
348 return this->operator[](well_name);
349 }
350
351 bool has(const std::string& well_name) const
352 {
353 return this->wells_.has(well_name);
354 }
355
356 bool operator==(const WellState&) const;
357
358 template<class Serializer>
359 void serializeOp(Serializer& serializer)
360 {
361 serializer(alq_state);
362 serializer(well_rates);
363 if (serializer.isSerializing()) {
364 serializer(wells_.size());
365 } else {
366 std::size_t size = 0;
367 serializer(size);
368 if (size != wells_.size()) {
369 OPM_THROW(std::runtime_error, "Error deserializing WellState: size mismatch");
370 }
371 }
372 for (auto& w : wells_) {
373 serializer(w);
374 }
375 serializer(permanently_inactive_well_names_);
376 }
377
378 bool is_permanently_inactive_well(const std::string& wname) const {
379 return std::find(this->permanently_inactive_well_names_.begin(), this->permanently_inactive_well_names_.end(), wname) != this->permanently_inactive_well_names_.end();
380 }
381
382private:
383 bool enableDistributedWells_ = false;
384
385 PhaseUsage phase_usage_;
386
387 // The wells_ variable is essentially a map of all the wells on the current
388 // process. Observe that since a well can be split over several processes a
389 // well might appear in the WellContainer on different processes.
391
392 // The members alq_state, global_well_info and well_rates are map like
393 // structures which will have entries for *all* the wells in the system.
394
395 // Use of std::optional<> here is a technical crutch, the
396 // WellStateFullyImplicitBlackoil class should be default constructible,
397 // whereas the GlobalWellInfo is not.
398 std::optional<GlobalWellInfo<Scalar>> global_well_info;
399 ALQState<Scalar> alq_state;
400
401 // The well_rates variable is defined for all wells on all processors. The
402 // bool in the value pair is whether the current process owns the well or
403 // not.
404 std::map<std::string, std::pair<bool, std::vector<Scalar>>> well_rates;
405
406 // Keep track of permanently inactive well names
407 std::vector<std::string> permanently_inactive_well_names_;
408
409 data::Segment
410 reportSegmentResults(const int well_id,
411 const int seg_ix,
412 const int seg_no) const;
413
414
420 void base_init(const std::vector<Scalar>& cellPressures,
421 const std::vector<Scalar>& cellTemperatures,
422 const std::vector<Well>& wells_ecl,
423 const std::vector<std::reference_wrapper<ParallelWellInfo<Scalar>>>& parallel_well_info,
424 const std::vector<std::vector<PerforationData<Scalar>>>& well_perf_data,
426
427 void initSingleWell(const std::vector<Scalar>& cellPressures,
428 const std::vector<Scalar>& cellTemperatures,
429 const Well& well,
430 const std::vector<PerforationData<Scalar>>& well_perf_data,
433
434 void initSingleProducer(const Well& well,
436 Scalar pressure_first_connection,
437 const std::vector<PerforationData<Scalar>>& well_perf_data,
439
440 void initSingleInjector(const Well& well,
442 Scalar pressure_first_connection,
444 const std::vector<PerforationData<Scalar>>& well_perf_data,
446
447 static void calculateSegmentRatesBeforeSum(const ParallelWellInfo<Scalar>& pw_info,
448 const std::vector<std::vector<int>>& segment_inlets,
449 const std::vector<std::vector<int>>& segment_perforations,
450 const std::vector<Scalar>& perforation_rates,
451 const int np,
452 const int segment,
453 std::vector<Scalar>& segment_rates);
454
455 void reportConnectionFactors(const std::size_t well_index,
456 std::vector<data::Connection>& connections) const;
457
458 void reportConnectionPressuresAndRates(const std::size_t well_index,
459 const PhaseUsage& pu,
460 std::vector<data::Connection>& connections) const;
461
462 void reportConnectionFilterCake(const std::size_t well_index,
463 std::vector<data::Connection>& connections) const;
464
465 void reportFractureStatistics(const std::vector<ConnFracStatistics<Scalar>>& stats,
466 std::vector<data::Connection>& connections) const;
467};
468
469} // namespace Opm
470
471#endif // OPM_WELLSTATEFULLYIMPLICITBLACKOIL_HEADER_INCLUDED
Class encapsulating some information about parallel wells.
Definition ParallelWellInfo.hpp:195
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition WellState.hpp:66
void initWellStateMSWell(const std::vector< Well > &wells_ecl, const WellState *prev_well_state)
init the MS well related.
Definition WellState.cpp:660
int numPhases() const
The number of phases present.
Definition WellState.hpp:285
void init(const std::vector< Scalar > &cellPressures, const std::vector< Scalar > &cellTemperatures, const Schedule &schedule, const std::vector< Well > &wells_ecl, const std::vector< std::reference_wrapper< ParallelWellInfo< Scalar > > > &parallel_well_info, const int report_step, const WellState *prevState, const std::vector< std::vector< PerforationData< Scalar > > > &well_perf_data, const SummaryState &summary_state, const bool enableDistributedWells)
Allocate and initialize if wells is non-null.
Definition WellState.cpp:270
std::vector< Scalar > & wellRates(std::size_t well_index)
One rate per well and phase.
Definition WellState.hpp:296
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
Static data associated with a well perforation.
Definition PerforationData.hpp:30
Definition BlackoilPhases.hpp:46