GasLiftSingleWell_impl.hpp
Go to the documentation of this file.
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_IMPL_HEADER_INCLUDED
21#define OPM_GASLIFT_SINGLE_WELL_IMPL_HEADER_INCLUDED
22
23// Improve IDE experience
24#ifndef OPM_GASLIFT_SINGLE_WELL_HEADER_INCLUDED
25#include <config.h>
27#endif
28
29#include <opm/common/TimingMacros.hpp>
30
31#include <opm/input/eclipse/Schedule/GasLiftOpt.hpp>
32#include <opm/input/eclipse/Schedule/Well/Well.hpp>
33
34#include <string>
35#include <vector>
36
37#include <fmt/format.h>
38
39namespace Opm {
40
41template<typename TypeTag>
44 const Simulator& simulator,
45 const SummaryState& summary_state,
46 DeferredLogger& deferred_logger,
47 WellState<Scalar>& well_state,
48 const GroupState<Scalar>& group_state,
49 GasLiftGroupInfo<Scalar>& group_info,
50 GLiftSyncGroups &sync_groups,
51 const Parallel::Communication& comm,
52 bool glift_debug)
53 // The parent class GasLiftSingleWellGeneric contains all stuff
54 // that is not dependent on TypeTag
55 : GasLiftSingleWellGeneric<Scalar>(deferred_logger,
56 well_state,
57 group_state,
58 well.wellEcl(),
59 summary_state,
60 group_info,
61 well.phaseUsage(),
62 simulator.vanguard().schedule(),
63 simulator.episodeIndex(),
64 sync_groups,
65 comm,
66 glift_debug)
67 , simulator_{simulator}
68 , well_{well}
69{
70 const auto& gl_well = *this->gl_well_;
71 if (this->useFixedAlq_(gl_well)) {
72 this->updateWellStateAlqFixedValue_(gl_well);
73 this->optimize_ = false; // lift gas supply is fixed
74 }
75 else {
76 setAlqMaxRate_(gl_well);
77 this->optimize_ = true;
78 }
79
80 setupPhaseVariables_();
81 // get the alq value used for this well for the previous iteration (a
82 // nonlinear iteration in assemble() in BlackoilWellModel).
83 // If gas lift optimization has not been applied to this well yet, the
84 // default value is used.
85 this->orig_alq_ = this->well_state_.well(this->well_name_).alq_state.get();
86 if (this->optimize_) {
87 this->setAlqMinRate_(gl_well);
88 // NOTE: According to item 4 in WLIFTOPT, this value does not
89 // have to be positive.
90 // TODO: Does it make sense to have a negative value?
91 this->alpha_w_ = gl_well.weight_factor();
92 if (this->alpha_w_ <= 0 ) {
93 this->displayWarning_("Nonpositive value for alpha_w ignored");
94 this->alpha_w_ = 1.0;
95 }
96
97 // NOTE: According to item 6 in WLIFTOPT:
98 // "If this value is greater than zero, the incremental gas rate will influence
99 // the calculation of the incremental gradient and may be used
100 // to discourage the allocation of lift gas to wells which produce more gas."
101 // TODO: Does this mean that we should ignore this value if it
102 // is negative?
103 this->alpha_g_ = gl_well.inc_weight_factor();
104
105 // TODO: adhoc value.. Should we keep max_iterations_ as a safety measure
106 // or does it not make sense to have it?
107 this->max_iterations_ = 1000;
108 }
109}
110
111/****************************************
112 * Private methods in alphabetical order
113 ****************************************/
114
115template<typename TypeTag>
116typename GasLiftSingleWell<TypeTag>::BasicRates
118computeWellRates_(Scalar bhp, bool bhp_is_limited, bool debug_output ) const
119{
120 std::vector<Scalar> potentials(this->NUM_PHASES, 0.0);
121 this->well_.computeWellRatesWithBhp(this->simulator_,
122 bhp,
123 potentials,
124 this->deferred_logger_);
125 if (debug_output) {
126 const std::string msg = fmt::format("computed well potentials given bhp {}, "
127 "oil: {}, gas: {}, water: {}", bhp,
128 -potentials[this->oil_pos_], -potentials[this->gas_pos_],
129 -potentials[this->water_pos_]);
130 this->displayDebugMessage_(msg);
131 }
132
133 std::transform(potentials.begin(), potentials.end(),
134 potentials.begin(),
135 [](const auto& potential)
136 { return std::min(Scalar{0}, potential); });
137 return {-potentials[this->oil_pos_],
138 -potentials[this->gas_pos_],
139 -potentials[this->water_pos_],
140 bhp_is_limited
141 };
142}
143
144template<typename TypeTag>
145std::optional<typename GasLiftSingleWell<TypeTag>::Scalar>
146GasLiftSingleWell<TypeTag>::
147computeBhpAtThpLimit_(Scalar alq, bool debug_output) const
148{
149 OPM_TIMEFUNCTION();
150 auto bhp_at_thp_limit = this->well_.computeBhpAtThpLimitProdWithAlq(
151 this->simulator_,
152 this->summary_state_,
153 alq,
154 this->deferred_logger_,
155 /*iterate_if_no_solution */ false);
156 if (bhp_at_thp_limit) {
157 if (*bhp_at_thp_limit < this->controls_.bhp_limit) {
158 if (debug_output && this->debug) {
159 const std::string msg = fmt::format(
160 "Computed bhp ({}) from thp limit is below bhp limit ({}), (ALQ = {})."
161 " Using bhp limit instead",
162 *bhp_at_thp_limit, this->controls_.bhp_limit, alq
163 );
164 this->displayDebugMessage_(msg);
165 }
166 bhp_at_thp_limit = this->controls_.bhp_limit;
167 }
168 //bhp_at_thp_limit = std::max(*bhp_at_thp_limit, this->controls_.bhp_limit);
169 }
170 else {
171 const std::string msg = fmt::format(
172 "Failed in getting converged bhp potential from thp limit (ALQ = {})", alq);
173 this->displayDebugMessage_(msg);
174 }
175 return bhp_at_thp_limit;
176}
177
178template<typename TypeTag>
179void
180GasLiftSingleWell<TypeTag>::
181setupPhaseVariables_()
182{
183 const auto& pu = this->phase_usage_;
184#ifndef NDEBUG
185 bool num_phases_ok = (pu.num_phases == 3);
186#endif
187 if (pu.num_phases == 2) {
188 // NOTE: We support two-phase oil-water flow, by setting the gas flow rate
189 // to zero. This is done by initializing the potential vector to zero:
190 //
191 // std::vector<double> potentials(NUM_PHASES, 0.0);
192 //
193 // see e.g. runOptimizeLoop_() in GasLiftSingleWellGeneric.cpp
194 // In addition the VFP calculations, e.g. to calculate BHP from THP
195 // has been adapted to the two-phase oil-water case, see the comment
196 // in WellInterfaceGeneric.cpp for the method adaptRatesForVFP() for
197 // more information.
198 if ( pu.phase_used[BlackoilPhases::Aqua] == 1
199 && pu.phase_used[BlackoilPhases::Liquid] == 1
200 && pu.phase_used[BlackoilPhases::Vapour] == 0)
201 {
202#ifndef NDEBUG
203 num_phases_ok = true; // two-phase oil-water is also supported
204#endif
205 }
206 else {
207 throw std::logic_error("Two-phase gas lift optimization only supported"
208 " for oil and water");
209 }
210 }
211 assert(num_phases_ok);
212 this->oil_pos_ = pu.phase_pos[this->Oil];
213 this->gas_pos_ = pu.phase_pos[this->Gas];
214 this->water_pos_ = pu.phase_pos[this->Water];
215}
216
217template<typename TypeTag>
218void
219GasLiftSingleWell<TypeTag>::
220setAlqMaxRate_(const GasLiftWell& well)
221{
222 const auto& max_alq_optional = well.max_rate();
223 if (max_alq_optional) {
224 // NOTE: To prevent extrapolation of the VFP tables, any value
225 // entered here must not exceed the largest ALQ value in the well's VFP table.
226 this->max_alq_ = *max_alq_optional;
227 }
228 else { // i.e. WLIFTOPT, item 3 has been defaulted
229 // According to the manual for WLIFTOPT, item 3:
230 // The default value should be set to the largest ALQ
231 // value in the well's VFP table
232 const auto& table = well_.vfpProperties()->getProd()->getTable(
233 this->controls_.vfp_table_number);
234 const auto& alq_values = table.getALQAxis();
235 // Assume the alq_values are sorted in ascending order, so
236 // the last item should be the largest value:
237 this->max_alq_ = alq_values.back();
238 }
239}
240
241template<typename TypeTag>
242bool
243GasLiftSingleWell<TypeTag>::
244checkThpControl_() const
245{
246 const int well_index = this->well_state_.index(this->well_name_).value();
247 const Well::ProducerCMode& control_mode =
248 this->well_state_.well(well_index).production_cmode;
249 bool thp_control = control_mode == Well::ProducerCMode::THP;
250 const auto& well = getWell();
251 thp_control = thp_control || well.thpLimitViolatedButNotSwitched();
252 if (this->debug) {
253 if (!thp_control) {
254 this->displayDebugMessage_("Well is not under THP control, skipping iteration..");
255 }
256 }
257 return thp_control;
258}
259
260}
261
262#endif
Scalar get() const
Definition: DeferredLogger.hpp:57
WellState< GetPropType< TypeTag, Properties::Scalar > > & well_state_
Definition: GasLiftCommon.hpp:54
Definition: GasLiftGroupInfo.hpp:46
Definition: GasLiftSingleWellGeneric.hpp:50
bool optimize_
Definition: GasLiftSingleWellGeneric.hpp:454
std::string well_name_
Definition: GasLiftSingleWellGeneric.hpp:450
GetPropType< TypeTag, Properties::Scalar > alpha_w_
Definition: GasLiftSingleWellGeneric.hpp:440
GetPropType< TypeTag, Properties::Scalar > alpha_g_
Definition: GasLiftSingleWellGeneric.hpp:441
GetPropType< TypeTag, Properties::Scalar > orig_alq_
Definition: GasLiftSingleWellGeneric.hpp:438
int max_iterations_
Definition: GasLiftSingleWellGeneric.hpp:448
const GasLiftWell * gl_well_
Definition: GasLiftSingleWellGeneric.hpp:452
Definition: GasLiftSingleWell.hpp:37
GasLiftSingleWell(const WellInterface< TypeTag > &well, const Simulator &simulator, const SummaryState &summary_state, DeferredLogger &deferred_logger, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, GasLiftGroupInfo< Scalar > &group_info, GLiftSyncGroups &sync_groups, const Parallel::Communication &comm, bool glift_debug)
Definition: GasLiftSingleWell_impl.hpp:43
Definition: GroupState.hpp:43
ALQState< Scalar > alq_state
Definition: SingleWellState.hpp:123
Definition: WellInterface.hpp:77
Definition: WellState.hpp:65
const SingleWellState< Scalar > & well(std::size_t well_index) const
Definition: WellState.hpp:285
Dune::Communication< MPIComm > Communication
Definition: ParallelCommunication.hpp:30
Definition: blackoilboundaryratevector.hh:39
PhaseUsage phaseUsage(const Phases &phases)
Determine the active phases.