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