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/core/props/BlackoilPhases.hpp>
28
29#include <opm/common/ErrorMacros.hpp>
30
31#include <opm/input/eclipse/Schedule/Events.hpp>
32
33#include <opm/output/data/Wells.hpp>
34
35#include <opm/simulators/utils/ParallelCommunication.hpp>
36#include <opm/simulators/wells/ALQState.hpp>
37#include <opm/simulators/wells/GlobalWellInfo.hpp>
38#include <opm/simulators/wells/PerfData.hpp>
39#include <opm/simulators/wells/PerforationData.hpp>
40#include <opm/simulators/wells/SegmentState.hpp>
41#include <opm/simulators/wells/SingleWellState.hpp>
42#include <opm/simulators/wells/WellContainer.hpp>
43
44#include <functional>
45#include <map>
46#include <optional>
47#include <string>
48#include <utility>
49#include <vector>
50
51namespace Opm
52{
53
54class ParallelWellInfo;
55class Schedule;
56enum class WellStatus;
57
61{
62public:
63 static const uint64_t event_mask = ScheduleEvents::WELL_STATUS_CHANGE + ScheduleEvents::PRODUCTION_UPDATE + ScheduleEvents::INJECTION_UPDATE;
64 // TODO: same definition with WellInterface, eventually they should go to a common header file.
65 static const int Water = BlackoilPhases::Aqua;
66 static const int Oil = BlackoilPhases::Liquid;
67 static const int Gas = BlackoilPhases::Vapour;
68
69 // Only usable for testing purposes
70 explicit WellState(const ParallelWellInfo& pinfo);
71
72 explicit WellState(const PhaseUsage& pu)
73 : phase_usage_(pu)
74 {}
75
76 static WellState serializationTestObject(const ParallelWellInfo& pinfo);
77
78 std::size_t size() const {
79 return this->wells_.size();
80 }
81
82 std::vector<std::string> wells() const {
83 return this->wells_.wells();
84 }
85
86
87 int numWells() const
88 {
89 return this->size();
90 }
91
92 const ParallelWellInfo& parallelWellInfo(std::size_t well_index) const;
93
94
95
99 void init(const std::vector<double>& cellPressures,
100 const Schedule& schedule,
101 const std::vector<Well>& wells_ecl,
102 const std::vector<std::reference_wrapper<ParallelWellInfo>>& parallel_well_info,
103 const int report_step,
104 const WellState* prevState,
105 const std::vector<std::vector<PerforationData>>& well_perf_data,
106 const SummaryState& summary_state);
107
108 void resize(const std::vector<Well>& wells_ecl,
109 const std::vector<std::reference_wrapper<ParallelWellInfo>>& parallel_well_info,
110 const Schedule& schedule,
111 const bool handle_ms_well,
112 const std::size_t numCells,
113 const std::vector<std::vector<PerforationData>>& well_perf_data,
114 const SummaryState& summary_state);
115
116 void setCurrentWellRates(const std::string& wellName, const std::vector<double>& new_rates ) {
117 auto& [owner, rates] = this->well_rates.at(wellName);
118 if (owner)
119 rates = new_rates;
120 }
121
122 const std::vector<double>& currentWellRates(const std::string& wellName) const;
123
124 bool hasWellRates(const std::string& wellName) const {
125 return this->well_rates.find(wellName) != this->well_rates.end();
126 }
127
128 void clearWellRates()
129 {
130 this->well_rates.clear();
131 }
132
133 void gatherVectorsOnRoot(const std::vector< data::Connection >& from_connections,
134 std::vector< data::Connection >& to_connections,
135 const Parallel::Communication& comm) const;
136
137 data::Wells
138 report(const int* globalCellIdxMap,
139 const std::function<bool(const int)>& wasDynamicallyClosed) const;
140
141 void reportConnections(std::vector<data::Connection>& connections, const PhaseUsage &pu,
142 std::size_t well_index,
143 const int* globalCellIdxMap) const;
144
146 void initWellStateMSWell(const std::vector<Well>& wells_ecl,
147 const WellState* prev_well_state);
148
149 static void calculateSegmentRates(const std::vector<std::vector<int>>& segment_inlets, const std::vector<std::vector<int>>&segment_perforations,
150 const std::vector<double>& perforation_rates, const int np, const int segment, std::vector<double>& segment_rates);
151
152
153 void communicateGroupRates(const Parallel::Communication& comm);
154
155 void updateGlobalIsGrup(const Parallel::Communication& comm);
156
157 bool isInjectionGrup(const std::string& name) const {
158 return this->global_well_info.value().in_injecting_group(name);
159 }
160
161 bool isProductionGrup(const std::string& name) const {
162 return this->global_well_info.value().in_producing_group(name);
163 }
164
165 double getALQ( const std::string& name) const
166 {
167 return this->alq_state.get(name);
168 }
169
170 void setALQ( const std::string& name, double value)
171 {
172 this->alq_state.set(name, value);
173 }
174
175 int gliftGetDebugCounter() {
176 return this->alq_state.get_debug_counter();
177 }
178
179 void gliftSetDebugCounter(int value) {
180 return this->alq_state.set_debug_counter(value);
181 }
182
183 int gliftUpdateDebugCounter() {
184 return this->alq_state.update_debug_counter();
185 }
186
187 bool gliftCheckAlqOscillation(const std::string &name) const {
188 return this->alq_state.oscillation(name);
189 }
190
191 int gliftGetAlqDecreaseCount(const std::string &name) {
192 return this->alq_state.get_decrement_count(name);
193 }
194
195 int gliftGetAlqIncreaseCount(const std::string &name) {
196 return this->alq_state.get_increment_count(name);
197 }
198
199 void gliftUpdateAlqIncreaseCount(const std::string &name, bool increase) {
200 this->alq_state.update_count(name, increase);
201 }
202
203 void gliftTimeStepInit() {
204 this->alq_state.reset_count();
205 }
206
207 // If the ALQ has changed since the previous time step,
208 // reset current_alq and update default_alq. ALQ is used for
209 // constant lift gas injection and for gas lift optimization
210 // (THP controlled wells).
211 void updateWellsDefaultALQ(const std::vector<Well>& wells_ecl, const SummaryState& summary_state);
212
213 int wellNameToGlobalIdx(const std::string &name) {
214 return this->global_well_info.value().well_index(name);
215 }
216
217 std::string globalIdxToWellName(const int index) {
218 return this->global_well_info.value().well_name(index);
219 }
220
221 bool wellIsOwned(std::size_t well_index,
222 const std::string& wellName) const;
223
224 bool wellIsOwned(const std::string& wellName) const;
225
226 void updateStatus(int well_index, WellStatus status);
227
228 void openWell(int well_index);
229 void shutWell(int well_index);
230 void stopWell(int well_index);
231
233 int numPhases() const
234 {
235 return this->phase_usage_.num_phases;
236 }
237
238 const PhaseUsage& phaseUsage() const {
239 return this->phase_usage_;
240 }
241
243 std::vector<double>& wellRates(std::size_t well_index) { return this->wells_[well_index].surface_rates; }
244 const std::vector<double>& wellRates(std::size_t well_index) const { return this->wells_[well_index].surface_rates; }
245
246 const std::string& name(std::size_t well_index) const {
247 return this->wells_.well_name(well_index);
248 }
249
250 std::optional<std::size_t> index(const std::string& well_name) const {
251 return this->wells_.well_index(well_name);
252 }
253
254 const SingleWellState& operator[](std::size_t well_index) const {
255 return this->wells_[well_index];
256 }
257
258 const SingleWellState& operator[](const std::string& well_name) const {
259 return this->wells_[well_name];
260 }
261
262 SingleWellState& operator[](std::size_t well_index) {
263 return this->wells_[well_index];
264 }
265
266 SingleWellState& operator[](const std::string& well_name) {
267 return this->wells_[well_name];
268 }
269
270 const SingleWellState& well(std::size_t well_index) const {
271 return this->operator[](well_index);
272 }
273
274 const SingleWellState& well(const std::string& well_name) const {
275 return this->operator[](well_name);
276 }
277
278 SingleWellState& well(std::size_t well_index) {
279 return this->operator[](well_index);
280 }
281
282 SingleWellState& well(const std::string& well_name) {
283 return this->operator[](well_name);
284 }
285
286 bool has(const std::string& well_name) const {
287 return this->wells_.has(well_name);
288 }
289
290 bool operator==(const WellState&) const;
291
292 template<class Serializer>
293 void serializeOp(Serializer& serializer)
294 {
295 serializer(alq_state);
296 serializer(well_rates);
297 if (serializer.isSerializing()) {
298 serializer(wells_.size());
299 } else {
300 std::size_t size = 0;
301 serializer(size);
302 if (size != wells_.size()) {
303 OPM_THROW(std::runtime_error, "Error deserializing WellState: size mismatch");
304 }
305 }
306 for (auto& w : wells_) {
307 serializer(w);
308 }
309 }
310
311private:
312 PhaseUsage phase_usage_;
313
314 // The wells_ variable is essentially a map of all the wells on the current
315 // process. Observe that since a well can be split over several processes a
316 // well might appear in the WellContainer on different processes.
317 WellContainer<SingleWellState> wells_;
318
319 // The members alq_state, global_well_info and well_rates are map like
320 // structures which will have entries for *all* the wells in the system.
321
322 // Use of std::optional<> here is a technical crutch, the
323 // WellStateFullyImplicitBlackoil class should be default constructible,
324 // whereas the GlobalWellInfo is not.
325 std::optional<GlobalWellInfo> global_well_info;
326 ALQState alq_state;
327
328 // The well_rates variable is defined for all wells on all processors. The
329 // bool in the value pair is whether the current process owns the well or
330 // not.
331 std::map<std::string, std::pair<bool, std::vector<double>>> well_rates;
332
333 data::Segment
334 reportSegmentResults(const int well_id,
335 const int seg_ix,
336 const int seg_no) const;
337
338
344 void base_init(const std::vector<double>& cellPressures,
345 const std::vector<Well>& wells_ecl,
346 const std::vector<std::reference_wrapper<ParallelWellInfo>>& parallel_well_info,
347 const std::vector<std::vector<PerforationData>>& well_perf_data,
348 const SummaryState& summary_state);
349
350 void initSingleWell(const std::vector<double>& cellPressures,
351 const Well& well,
352 const std::vector<PerforationData>& well_perf_data,
353 const ParallelWellInfo& well_info,
354 const SummaryState& summary_state);
355
356 void initSingleProducer(const Well& well,
357 const ParallelWellInfo& well_info,
358 double pressure_first_connection,
359 const std::vector<PerforationData>& well_perf_data,
360 const SummaryState& summary_state);
361
362 void initSingleInjector(const Well& well,
363 const ParallelWellInfo& well_info,
364 double pressure_first_connection,
365 const std::vector<PerforationData>& well_perf_data,
366 const SummaryState& summary_state);
367
368};
369
370} // namespace Opm
371
372
373#endif // OPM_WELLSTATEFULLYIMPLICITBLACKOIL_HEADER_INCLUDED
Class encapsulating some information about parallel wells.
Definition ParallelWellInfo.hpp:184
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition WellState.hpp:61
void init(const std::vector< double > &cellPressures, const Schedule &schedule, const std::vector< Well > &wells_ecl, const std::vector< std::reference_wrapper< ParallelWellInfo > > &parallel_well_info, const int report_step, const WellState *prevState, const std::vector< std::vector< PerforationData > > &well_perf_data, const SummaryState &summary_state)
Allocate and initialize if wells is non-null.
Definition WellState.cpp:241
void initWellStateMSWell(const std::vector< Well > &wells_ecl, const WellState *prev_well_state)
init the MS well related.
Definition WellState.cpp:658
int numPhases() const
The number of phases present.
Definition WellState.hpp:233
std::vector< double > & wellRates(std::size_t well_index)
One rate per well and phase.
Definition WellState.hpp:243
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition BlackoilPhases.hpp:27
Definition BlackoilPhases.hpp:46