20#ifndef OPM_GASLIFT_SINGLE_WELL_IMPL_HEADER_INCLUDED
21#define OPM_GASLIFT_SINGLE_WELL_IMPL_HEADER_INCLUDED
24#ifndef OPM_GASLIFT_SINGLE_WELL_HEADER_INCLUDED
29#include <opm/common/TimingMacros.hpp>
31#include <opm/input/eclipse/Schedule/GasLiftOpt.hpp>
32#include <opm/input/eclipse/Schedule/Well/Well.hpp>
37#include <fmt/format.h>
41template<
typename TypeTag>
44 const Simulator& simulator,
45 const SummaryState& summary_state,
50 GLiftSyncGroups &sync_groups,
62 simulator.vanguard().schedule(),
63 simulator.episodeIndex(),
67 , simulator_{simulator}
70 const auto& gl_well = *this->
gl_well_;
76 setAlqMaxRate_(gl_well);
80 setupPhaseVariables_();
91 this->
alpha_w_ = gl_well.weight_factor();
103 this->
alpha_g_ = gl_well.inc_weight_factor();
115template<
typename TypeTag>
116typename GasLiftSingleWell<TypeTag>::BasicRates
120 std::vector<Scalar> potentials(this->NUM_PHASES, 0.0);
121 this->well_.computeWellRatesWithBhp(this->simulator_,
124 this->deferred_logger_);
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);
133 std::transform(potentials.begin(), potentials.end(),
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_],
144template<
typename TypeTag>
145std::optional<typename GasLiftSingleWell<TypeTag>::Scalar>
146GasLiftSingleWell<TypeTag>::
147computeBhpAtThpLimit_(Scalar alq,
bool debug_output)
const
150 auto bhp_at_thp_limit = this->well_.computeBhpAtThpLimitProdWithAlq(
152 this->summary_state_,
154 this->deferred_logger_,
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
164 this->displayDebugMessage_(msg);
166 bhp_at_thp_limit = this->controls_.bhp_limit;
171 const std::string msg = fmt::format(
172 "Failed in getting converged bhp potential from thp limit (ALQ = {})", alq);
173 this->displayDebugMessage_(msg);
175 return bhp_at_thp_limit;
178template<
typename TypeTag>
180GasLiftSingleWell<TypeTag>::
181setupPhaseVariables_()
183 const auto& pu = this->phase_usage_;
185 bool num_phases_ok = (pu.num_phases == 3);
187 if (pu.num_phases == 2) {
198 if ( pu.phase_used[BlackoilPhases::Aqua] == 1
199 && pu.phase_used[BlackoilPhases::Liquid] == 1
200 && pu.phase_used[BlackoilPhases::Vapour] == 0)
203 num_phases_ok =
true;
207 throw std::logic_error(
"Two-phase gas lift optimization only supported"
208 " for oil and water");
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];
217template<
typename TypeTag>
219GasLiftSingleWell<TypeTag>::
220setAlqMaxRate_(
const GasLiftWell& well)
222 const auto& max_alq_optional = well.max_rate();
223 if (max_alq_optional) {
226 this->max_alq_ = *max_alq_optional;
232 const auto& table = well_.vfpProperties()->getProd()->getTable(
233 this->controls_.vfp_table_number);
234 const auto& alq_values = table.getALQAxis();
237 this->max_alq_ = alq_values.back();
241template<
typename TypeTag>
243GasLiftSingleWell<TypeTag>::
244checkThpControl_()
const
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();
254 this->displayDebugMessage_(
"Well is not under THP control, skipping iteration..");
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
void displayWarning_(const std::string &warning)
GetPropType< TypeTag, Properties::Scalar > orig_alq_
Definition: GasLiftSingleWellGeneric.hpp:438
bool useFixedAlq_(const GasLiftWell &well)
void updateWellStateAlqFixedValue_(const GasLiftWell &well)
void setAlqMinRate_(const GasLiftWell &well)
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.