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