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