My Project
Loading...
Searching...
No Matches
GasLiftSingleWellGeneric.hpp
1/*
2 Copyright 2020 Equinor ASA.
3
4 This file is part of the Open Porous Media project (OPM).
5
6 OPM is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 3 of the License, or
9 (at your option) any later version.
10
11 OPM is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with OPM. If not, see <http://www.gnu.org/licenses/>.
18*/
19
20#ifndef OPM_GASLIFT_SINGLE_WELL_GENERIC_HEADER_INCLUDED
21#define OPM_GASLIFT_SINGLE_WELL_GENERIC_HEADER_INCLUDED
22
23#include <opm/input/eclipse/Schedule/Well/WellProductionControls.hpp>
24
25#include <opm/simulators/wells/GasLiftGroupInfo.hpp>
26#include <opm/simulators/wells/GasLiftCommon.hpp>
27
28#include <opm/simulators/utils/BlackoilPhases.hpp>
29
30#include <optional>
31#include <set>
32#include <stdexcept>
33#include <string>
34#include <tuple>
35#include <utility>
36
37namespace Opm {
38
39class DeferredLogger;
40class GasLiftWell;
41template<class Scalar> class GasLiftWellState;
42class Schedule;
43class SummaryState;
44template<class Scalar> class WellInterfaceGeneric;
45template<class Scalar> class WellState;
46template<class Scalar> class GroupState;
47
48template<class Scalar>
50{
51protected:
52 static constexpr int Water = BlackoilPhases::Aqua;
53 static constexpr int Oil = BlackoilPhases::Liquid;
54 static constexpr int Gas = BlackoilPhases::Vapour;
55 static constexpr int NUM_PHASES = 3;
56 static constexpr Scalar ALQ_EPSILON = 1e-8;
57
58public:
59 using GLiftSyncGroups = std::set<int>;
60 using Rate = typename GasLiftGroupInfo<Scalar>::Rate;
61 using MessageType = typename GasLiftCommon<Scalar>::MessageType;
62
63 struct GradInfo
64 {
65 GradInfo() = default;
66 GradInfo(Scalar grad_,
67 Scalar new_oil_rate_,
68 bool oil_is_limited_,
69 Scalar new_gas_rate_,
70 bool gas_is_limited_,
71 Scalar new_water_rate_,
72 bool water_is_limited_,
73 Scalar alq_,
74 bool alq_is_limited_)
75 : grad{grad_}
76 , new_oil_rate{new_oil_rate_}
77 , oil_is_limited{oil_is_limited_}
78 , new_gas_rate{new_gas_rate_}
79 , gas_is_limited{gas_is_limited_}
80 , new_water_rate{new_water_rate_}
81 , water_is_limited{water_is_limited_}
82 , alq{alq_}
83 , alq_is_limited{alq_is_limited_}
84 {}
85
86 Scalar grad;
87 Scalar new_oil_rate;
88 bool oil_is_limited;
89 Scalar new_gas_rate;
90 bool gas_is_limited;
91 Scalar new_water_rate;
92 bool water_is_limited;
93 Scalar alq;
94 bool alq_is_limited;
95 };
96
97 virtual ~GasLiftSingleWellGeneric() = default;
98
99 const std::string& name() const { return well_name_; }
100
101 std::optional<GradInfo> calcIncOrDecGradient(Scalar oil_rate,
102 Scalar gas_rate,
103 Scalar water_rate,
104 Scalar alq,
105 const std::string& gr_name_dont_limit,
106 bool increase,
107 bool debug_output = true) const;
108
109 std::unique_ptr<GasLiftWellState<Scalar>> runOptimize(const int iteration_idx);
110
111 virtual const WellInterfaceGeneric<Scalar>& getWell() const = 0;
112
113protected:
115 WellState<Scalar>& well_state,
116 const GroupState<Scalar>& group_state,
117 const Well& ecl_well,
118 const SummaryState& summary_state,
119 GasLiftGroupInfo<Scalar>& group_info,
120 const PhaseUsage& phase_usage,
121 const Schedule& schedule,
122 const int report_step_idx,
123 GLiftSyncGroups& sync_groups,
124 const Parallel::Communication& comm,
125 bool glift_debug);
126
127 struct LimitedRates;
129 {
130 BasicRates(const BasicRates& rates) :
131 oil{rates.oil},
132 gas{rates.gas},
133 water{rates.water},
134 bhp_is_limited{rates.bhp_is_limited}
135 {}
136
137 BasicRates(Scalar oil_,
138 Scalar gas_,
139 Scalar water_,
140 bool bhp_is_limited_)
141 : oil{oil_}
142 , gas{gas_}
143 , water{water_}
144 , bhp_is_limited{bhp_is_limited_}
145 {}
146
147 BasicRates& operator=(const BasicRates& rates)
148 {
149 oil = rates.oil;
150 gas = rates.gas;
151 water = rates.water;
152 bhp_is_limited = rates.bhp_is_limited;
153 return *this;
154 }
155
156 // This copy constructor cannot be defined inline here since LimitedRates
157 // has not been defined yet (it is defined below). Instead it is defined in
158 // in the .cpp file
159 BasicRates(const LimitedRates& rates);
160
161 Scalar operator[](Rate rate_type) const
162 {
163 switch (rate_type) {
164 case Rate::oil:
165 return this->oil;
166 case Rate::gas:
167 return this->gas;
168 case Rate::water:
169 return this->water;
170 case Rate::liquid:
171 return this->oil + this->water;
172 default:
173 throw std::runtime_error("This should not happen");
174 }
175 }
176
177 Scalar oil, gas, water;
178 bool bhp_is_limited;
179 };
180
181 struct LimitedRates : public BasicRates
182 {
183 enum class LimitType {well, group, none};
184 LimitedRates(Scalar oil_,
185 Scalar gas_,
186 Scalar water_,
187 bool oil_is_limited_,
188 bool gas_is_limited_,
189 bool water_is_limited_,
190 bool bhp_is_limited_,
191 std::optional<Rate> oil_limiting_target_,
192 std ::optional<Rate> water_limiting_target_)
193 : BasicRates(oil_, gas_, water_, bhp_is_limited_)
194 , oil_is_limited{oil_is_limited_}
195 , gas_is_limited{gas_is_limited_}
196 , water_is_limited{water_is_limited_}
197 , oil_limiting_target{oil_limiting_target_}
198 , water_limiting_target{water_limiting_target_}
199 {
200 set_initial_limit_type_();
201 }
202
203 LimitedRates(const BasicRates& rates,
204 bool oil_is_limited_,
205 bool gas_is_limited_,
206 bool water_is_limited_)
207 : BasicRates(rates)
208 , oil_is_limited{oil_is_limited_}
209 , gas_is_limited{gas_is_limited_}
210 , water_is_limited{water_is_limited_}
211 {
212 set_initial_limit_type_();
213 }
214
215 bool limited() const
216 {
217 return oil_is_limited || gas_is_limited || water_is_limited;
218 }
219
220 // For a given ALQ value, were the rates limited due to group targets
221 // or due to well targets?
222 LimitType limit_type;
223 bool oil_is_limited;
224 bool gas_is_limited;
225 bool water_is_limited;
226 std::optional<Rate> oil_limiting_target;
227 std::optional<Rate> water_limiting_target;
228
229 private:
230 void set_initial_limit_type_()
231 {
232 limit_type = limited() ? LimitType::well : LimitType::none;
233 }
234 };
235
237 {
238 OptimizeState( GasLiftSingleWellGeneric& parent_, bool increase_ )
239 : parent{parent_}
240 , increase{increase_}
241 , it{0}
242 , stop_iteration{false}
243 , bhp{-1}
244 {}
245
247 bool increase;
248 int it;
249 bool stop_iteration;
250 Scalar bhp;
251
252 std::pair<std::optional<Scalar>,bool> addOrSubtractAlqIncrement(Scalar alq);
253 Scalar calcEcoGradient(Scalar oil_rate,
254 Scalar new_oil_rate,
255 Scalar gas_rate,
256 Scalar new_gas_rate);
257
258 bool checkAlqOutsideLimits(Scalar alq, Scalar oil_rate);
259 bool checkEcoGradient(Scalar gradient);
260 bool checkOilRateExceedsTarget(Scalar oil_rate);
261 bool checkRatesViolated(const LimitedRates& rates) const;
262
263 void debugShowIterationInfo(Scalar alq);
264
265 Scalar getBhpWithLimit();
266
267 void warn_(const std::string& msg) { parent.displayWarning_(msg); }
268 };
269
270 bool checkGroupALQrateExceeded(Scalar delta_alq,
271 const std::string& gr_name_dont_limit = "") const;
272 bool checkGroupTotalRateExceeded(Scalar delta_alq,
273 Scalar delta_gas_rate) const;
274
275 std::pair<std::optional<Scalar>, bool>
276 addOrSubtractAlqIncrement_(Scalar alq, bool increase) const;
277
278 Scalar calcEcoGradient_(Scalar oil_rate, Scalar new_oil_rate,
279 Scalar gas_rate, Scalar new_gas_rate, bool increase) const;
280
281 bool checkALQequal_(Scalar alq1, Scalar alq2) const;
282
283 bool checkGroupTargetsViolated(const BasicRates& rates,
284 const BasicRates& new_rates) const;
285 bool checkInitialALQmodified_(Scalar alq, Scalar initial_alq) const;
286
287 virtual bool checkThpControl_() const = 0;
288 virtual std::optional<Scalar > computeBhpAtThpLimit_(Scalar alq,
289 bool debug_output = true) const = 0;
290
291 std::pair<std::optional<Scalar>,Scalar>
292 computeConvergedBhpAtThpLimitByMaybeIncreasingALQ_() const;
293
294 std::pair<std::optional<BasicRates>,Scalar>
295 computeInitialWellRates_() const;
296
297 std::optional<LimitedRates>
298 computeLimitedWellRatesWithALQ_(Scalar alq) const;
299
300 virtual BasicRates computeWellRates_(Scalar bhp,
301 bool bhp_is_limited,
302 bool debug_output = true) const = 0;
303
304 std::optional<BasicRates> computeWellRatesWithALQ_(Scalar alq) const;
305
306 void debugCheckNegativeGradient_(Scalar grad, Scalar alq, Scalar new_alq,
307 Scalar oil_rate, Scalar new_oil_rate,
308 Scalar gas_rate, Scalar new_gas_rate,
309 bool increase) const;
310
311 void debugPrintWellStateRates() const;
312 void debugShowAlqIncreaseDecreaseCounts_();
313 void debugShowBhpAlqTable_();
314 void debugShowLimitingTargets_(const LimitedRates& rates) const;
315 void debugShowProducerControlMode() const;
316 void debugShowStartIteration_(Scalar alq, bool increase, Scalar oil_rate);
317 void debugShowTargets_();
318 void displayDebugMessage_(const std::string& msg) const override;
319 void displayWarning_(const std::string& warning);
320
321 std::pair<Scalar, bool> getBhpWithLimit_(Scalar bhp) const;
322 std::pair<Scalar, bool> getGasRateWithLimit_(const BasicRates& rates) const;
323 std::pair<Scalar, bool> getGasRateWithGroupLimit_(Scalar new_gas_rate,
324 Scalar gas_rate,
325 const std::string& gr_name_dont_limit) const;
326
327 std::pair<std::optional<LimitedRates>,Scalar >
328 getInitialRatesWithLimit_() const;
329
330 LimitedRates getLimitedRatesFromRates_(const BasicRates& rates) const;
331
332 std::tuple<Scalar,Scalar,bool,bool>
333 getLiquidRateWithGroupLimit_(const Scalar new_oil_rate,
334 const Scalar oil_rate,
335 const Scalar new_water_rate,
336 const Scalar water_rate,
337 const std::string& gr_name_dont_limit) const;
338
339 std::pair<Scalar, bool>
340 getOilRateWithGroupLimit_(Scalar new_oil_rate,
341 Scalar oil_rate,
342 const std::string& gr_name_dont_limit) const;
343
344 std::pair<Scalar, bool> getOilRateWithLimit_(const BasicRates& rates) const;
345
346 std::pair<Scalar, std::optional<Rate>>
347 getOilRateWithLimit2_(const BasicRates& rates) const;
348
349 Scalar getProductionTarget_(Rate rate) const;
350 Scalar getRate_(Rate rate_type, const BasicRates& rates) const;
351
352 std::pair<Scalar, std::optional<Rate>>
353 getRateWithLimit_(Rate rate_type, const BasicRates& rates) const;
354
355 std::tuple<Scalar, const std::string*, Scalar>
356 getRateWithGroupLimit_(Rate rate_type,
357 const Scalar new_rate,
358 const Scalar old_rate,
359 const std::string& gr_name_dont_limit) const;
360
361 std::pair<Scalar, bool>
362 getWaterRateWithGroupLimit_(Scalar new_water_rate,
363 Scalar water_rate,
364 const std::string& gr_name_dont_limit) const;
365
366 std::pair<Scalar, bool> getWaterRateWithLimit_(const BasicRates& rates) const;
367
368 std::pair<Scalar, std::optional<Rate>>
369 getWaterRateWithLimit2_(const BasicRates& rates) const;
370
371 BasicRates getWellStateRates_() const;
372 bool hasProductionControl_(Rate rate) const;
373
374 std::pair<LimitedRates, Scalar>
375 increaseALQtoPositiveOilRate_(Scalar alq,
376 const LimitedRates& orig_rates) const;
377
378 std::pair<LimitedRates, Scalar>
379 increaseALQtoMinALQ_(Scalar alq,
380 const LimitedRates& orig_rates) const;
381
382 void logSuccess_(Scalar alq,
383 const int iteration_idx);
384
385 std::pair<LimitedRates, Scalar>
386 maybeAdjustALQbeforeOptimizeLoop_(const LimitedRates& rates,
387 Scalar alq,
388 bool increase) const;
389
390 std::pair<LimitedRates, Scalar>
391 reduceALQtoGroupAlqLimits_(Scalar alq,
392 const LimitedRates& rates) const;
393
394 std::pair<LimitedRates, Scalar>
395 reduceALQtoGroupTarget(Scalar alq,
396 const LimitedRates& rates) const;
397
398 std::pair<LimitedRates, Scalar>
399 reduceALQtoWellTarget_(Scalar alq,
400 const LimitedRates& rates) const;
401
402 std::unique_ptr<GasLiftWellState<Scalar>> runOptimize1_();
403 std::unique_ptr<GasLiftWellState<Scalar>> runOptimize2_();
404 std::unique_ptr<GasLiftWellState<Scalar>> runOptimizeLoop_(bool increase);
405
406 void setAlqMinRate_(const GasLiftWell& well);
407 std::unique_ptr<GasLiftWellState<Scalar>> tryIncreaseLiftGas_();
408 std::unique_ptr<GasLiftWellState<Scalar>> tryDecreaseLiftGas_();
409
410 void updateGroupRates_(const LimitedRates& rates,
411 const LimitedRates& new_rates,
412 Scalar delta_alq) const;
413
415 updateRatesToGroupLimits_(const BasicRates& rates,
416 const LimitedRates& new_rates,
417 const std::string& gr_name = "") const;
418
419 void updateWellStateAlqFixedValue_(const GasLiftWell& well);
420 bool useFixedAlq_(const GasLiftWell& well);
421
422 void debugInfoGroupRatesExceedTarget(Rate rate_type,
423 const std::string& gr_name,
424 Scalar rate,
425 Scalar target) const;
426 void warnMaxIterationsExceeded_();
427
428 const Well& ecl_well_;
429 const SummaryState& summary_state_;
430 GasLiftGroupInfo<Scalar>& group_info_;
431 const PhaseUsage& phase_usage_;
432 GLiftSyncGroups& sync_groups_;
433 const WellProductionControls controls_;
434
435 Scalar increment_;
436 Scalar max_alq_;
437 Scalar min_alq_;
438 Scalar orig_alq_;
439
440 Scalar alpha_w_;
441 Scalar alpha_g_;
442 Scalar eco_grad_;
443
444 int gas_pos_;
445 int oil_pos_;
446 int water_pos_;
447
448 int max_iterations_;
449
450 std::string well_name_;
451
452 const GasLiftWell* gl_well_;
453
454 bool optimize_;
455 bool debug_limit_increase_decrease_;
456 bool debug_abort_if_decrease_and_oil_is_limited_ = false;
457 bool debug_abort_if_increase_and_gas_is_limited_ = false;
458};
459
460} // namespace Opm
461
462#endif // OPM_GASLIFT_SINGLE_WELL_GENERIC_HEADER_INCLUDED
Definition DeferredLogger.hpp:57
Definition GasLiftCommon.hpp:35
Definition GasLiftGroupInfo.hpp:46
Definition GasLiftSingleWellGeneric.hpp:50
Definition WellInterfaceFluidSystem.hpp:43
Definition WellTest.hpp:37
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition WellState.hpp:62
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
Definition GasLiftSingleWellGeneric.hpp:129
Definition GasLiftSingleWellGeneric.hpp:64
Definition GasLiftSingleWellGeneric.hpp:182
Definition GasLiftSingleWellGeneric.hpp:237
Definition BlackoilPhases.hpp:46