WellInterface.hpp
Go to the documentation of this file.
1/*
2 Copyright 2017 SINTEF Digital, Mathematics and Cybernetics.
3 Copyright 2017 Statoil ASA.
4 Copyright 2017 IRIS
5 Copyright 2019 Norce
6
7 This file is part of the Open Porous Media project (OPM).
8
9 OPM is free software: you can redistribute it and/or modify
10 it under the terms of the GNU General Public License as published by
11 the Free Software Foundation, either version 3 of the License, or
12 (at your option) any later version.
13
14 OPM is distributed in the hope that it will be useful,
15 but WITHOUT ANY WARRANTY; without even the implied warranty of
16 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17 GNU General Public License for more details.
18
19 You should have received a copy of the GNU General Public License
20 along with OPM. If not, see <http://www.gnu.org/licenses/>.
21*/
22
23#ifndef OPM_WELLINTERFACE_HEADER_INCLUDED
24#define OPM_WELLINTERFACE_HEADER_INCLUDED
25
26// NOTE: GasLiftSingleWell.hpp includes StandardWell.hpp which includes ourself
27// (WellInterface.hpp), so we need to forward declare GasLiftSingleWell
28// for it to be defined in this file. Similar for BlackoilWellModel
29namespace Opm {
30 template<typename TypeTag> class GasLiftSingleWell;
31 template<typename TypeTag> class BlackoilWellModel;
32}
33
34#include <opm/common/OpmLog/OpmLog.hpp>
35#include <opm/common/ErrorMacros.hpp>
36#include <opm/common/Exceptions.hpp>
37
38#include <opm/input/eclipse/Schedule/Well/WellTestState.hpp>
39
40#include <opm/material/fluidstates/BlackOilFluidState.hpp>
41
43
45
54
56
59
60#include <dune/common/fmatrix.hh>
61#include <dune/istl/bcrsmatrix.hh>
62#include <dune/istl/matrixmatrix.hh>
63
64#include <opm/material/densead/Evaluation.hpp>
65
66#include <vector>
67
68namespace Opm
69{
70
71class WellInjectionProperties;
72class WellProductionProperties;
73
74template<typename TypeTag>
75class WellInterface : public WellInterfaceIndices<GetPropType<TypeTag, Properties::FluidSystem>,
76 GetPropType<TypeTag, Properties::Indices>>
77{
80public:
92
93 using VectorBlockType = Dune::FieldVector<Scalar, Indices::numEq>;
94 using MatrixBlockType = Dune::FieldMatrix<Scalar, Indices::numEq, Indices::numEq>;
95 using Eval = typename Base::Eval;
96 using BVector = Dune::BlockVector<VectorBlockType>;
97 using PressureMatrix = Dune::BCRSMatrix<Opm::MatrixBlock<Scalar, 1, 1>>;
98
101
105
107
108 static constexpr bool has_solvent = getPropValue<TypeTag, Properties::EnableSolvent>();
109 static constexpr bool has_zFraction = getPropValue<TypeTag, Properties::EnableExtbo>();
110 static constexpr bool has_polymer = getPropValue<TypeTag, Properties::EnablePolymer>();
111 static constexpr bool has_energy = getPropValue<TypeTag, Properties::EnableEnergy>();
112 static const bool has_temperature = getPropValue<TypeTag, Properties::EnableTemperature>();
113 // flag for polymer molecular weight related
114 static constexpr bool has_polymermw = getPropValue<TypeTag, Properties::EnablePolymerMW>();
115 static constexpr bool has_foam = getPropValue<TypeTag, Properties::EnableFoam>();
116 static constexpr bool has_brine = getPropValue<TypeTag, Properties::EnableBrine>();
117 static constexpr bool has_watVapor = getPropValue<TypeTag, Properties::EnableVapwat>();
118 static constexpr bool has_disgas_in_water = getPropValue<TypeTag, Properties::EnableDisgasInWater>();
119 static constexpr bool has_saltPrecip = getPropValue<TypeTag, Properties::EnableSaltPrecipitation>();
120 static constexpr bool has_micp = getPropValue<TypeTag, Properties::EnableMICP>();
121
122 // For the conversion between the surface volume rate and reservoir voidage rate
123 using FluidState = BlackOilFluidState<Eval,
127 Indices::compositionSwitchIdx >= 0,
129 has_brine,
132 Indices::numPhases >;
134 WellInterface(const Well& well,
135 const ParallelWellInfo<Scalar>& pw_info,
136 const int time_step,
137 const ModelParameters& param,
138 const RateConverterType& rate_converter,
139 const int pvtRegionIdx,
140 const int num_components,
141 const int num_phases,
142 const int index_of_well,
143 const std::vector<PerforationData<Scalar>>& perf_data);
144
146 virtual ~WellInterface() = default;
147
148 virtual void init(const PhaseUsage* phase_usage_arg,
149 const std::vector<Scalar>& depth_arg,
150 const Scalar gravity_arg,
151 const std::vector<Scalar>& B_avg,
152 const bool changed_to_open_this_step);
153
155 const WellState<Scalar>& well_state,
156 const std::vector<Scalar>& B_avg,
157 DeferredLogger& deferred_logger,
158 const bool relax_tolerance) const = 0;
159
160 virtual void solveEqAndUpdateWellState(const Simulator& simulator,
161 WellState<Scalar>& well_state,
162 DeferredLogger& deferred_logger) = 0;
163
164 void assembleWellEq(const Simulator& simulator,
165 const double dt,
166 WellState<Scalar>& well_state,
167 const GroupState<Scalar>& group_state,
168 DeferredLogger& deferred_logger);
169
170 void assembleWellEqWithoutIteration(const Simulator& simulator,
171 const double dt,
172 WellState<Scalar>& well_state,
173 const GroupState<Scalar>& group_state,
174 DeferredLogger& deferred_logger);
175
176 // TODO: better name or further refactoring the function to make it more clear
177 void prepareWellBeforeAssembling(const Simulator& simulator,
178 const double dt,
179 WellState<Scalar>& well_state,
180 const GroupState<Scalar>& group_state,
181 DeferredLogger& deferred_logger);
182
183
184 virtual void computeWellRatesWithBhp(const Simulator& ebosSimulator,
185 const Scalar& bhp,
186 std::vector<Scalar>& well_flux,
187 DeferredLogger& deferred_logger) const = 0;
188
189 virtual std::optional<Scalar>
191 const SummaryState& summary_state,
192 const Scalar alq_value,
193 DeferredLogger& deferred_logger,
194 bool iterate_if_no_solution) const = 0;
195
199 const BVector& x,
200 WellState<Scalar>& well_state,
201 DeferredLogger& deferred_logger) = 0;
202
204 virtual void apply(const BVector& x, BVector& Ax) const = 0;
205
207 virtual void apply(BVector& r) const = 0;
208
209 // TODO: before we decide to put more information under mutable, this function is not const
210 virtual void computeWellPotentials(const Simulator& simulator,
211 const WellState<Scalar>& well_state,
212 std::vector<Scalar>& well_potentials,
213 DeferredLogger& deferred_logger) = 0;
214
215 virtual void updateWellStateWithTarget(const Simulator& simulator,
216 const GroupState<Scalar>& group_state,
217 WellState<Scalar>& well_state,
218 DeferredLogger& deferred_logger) const;
219
220 virtual void computeWellRatesWithBhpIterations(const Simulator& simulator,
221 const Scalar& bhp,
222 std::vector<Scalar>& well_flux,
223 DeferredLogger& deferred_logger) const = 0;
224
225 bool wellUnderZeroRateTarget(const Simulator& simulator,
226 const WellState<Scalar>& well_state,
227 DeferredLogger& deferred_logger) const;
228
229 bool wellUnderZeroGroupRateTarget(const Simulator& simulator,
230 const WellState<Scalar>& well_state,
231 DeferredLogger& deferred_logger,
232 std::optional<bool> group_control = std::nullopt) const;
233
234 bool stoppedOrZeroRateTarget(const Simulator& simulator,
235 const WellState<Scalar>& well_state,
236 DeferredLogger& deferred_logger) const;
237
238 bool updateWellStateWithTHPTargetProd(const Simulator& simulator,
239 WellState<Scalar>& well_state,
240 DeferredLogger& deferred_logger) const;
241
242 enum class IndividualOrGroup { Individual, Group, Both };
243 bool updateWellControl(const Simulator& simulator,
244 const IndividualOrGroup iog,
245 WellState<Scalar>& well_state,
246 const GroupState<Scalar>& group_state,
247 DeferredLogger& deferred_logger) /* const */;
248
250 WellState<Scalar>& well_state,
251 const GroupState<Scalar>& group_state,
252 const Well::InjectionControls& inj_controls,
253 const Well::ProductionControls& prod_controls,
254 const Scalar WQTotal,
255 DeferredLogger& deferred_logger,
256 const bool fixed_control = false,
257 const bool fixed_status = false);
258
259 virtual void updatePrimaryVariables(const Simulator& simulator,
260 const WellState<Scalar>& well_state,
261 DeferredLogger& deferred_logger) = 0;
262
263 virtual void calculateExplicitQuantities(const Simulator& simulator,
264 const WellState<Scalar>& well_state,
265 DeferredLogger& deferred_logger) = 0; // should be const?
266
267 virtual void updateProductivityIndex(const Simulator& simulator,
268 const WellProdIndexCalculator<Scalar>& wellPICalc,
269 WellState<Scalar>& well_state,
270 DeferredLogger& deferred_logger) const = 0;
271
272 // Add well contributions to matrix
274
276 const BVector& x,
277 const int pressureVarIndex,
278 const bool use_well_weights,
279 const WellState<Scalar>& well_state) const = 0;
280
281 void addCellRates(RateVector& rates, int cellIdx) const;
282
283 Scalar volumetricSurfaceRateForConnection(int cellIdx, int phaseIdx) const;
284
285 // TODO: theoretically, it should be a const function
286 // Simulator is not const is because that assembleWellEq is non-const Simulator
287 void wellTesting(const Simulator& simulator,
288 const double simulation_time,
289 /* const */ WellState<Scalar>& well_state,
290 const GroupState<Scalar>& group_state,
291 WellTestState& welltest_state,
292 const PhaseUsage& phase_usage,
293 GLiftEclWells& ecl_well_map,
294 std::map<std::string, double>& open_times,
295 DeferredLogger& deferred_logger);
296
297 void checkWellOperability(const Simulator& simulator,
298 const WellState<Scalar>& well_state,
299 DeferredLogger& deferred_logger);
300
302 WellState<Scalar>& well_state,
303 const GroupState<Scalar>& group_state,
304 const PhaseUsage& phase_usage,
305 GLiftEclWells& ecl_well_map,
306 DeferredLogger& deferred_logger);
307
308 // check whether the well is operable under the current reservoir condition
309 // mostly related to BHP limit and THP limit
310 void updateWellOperability(const Simulator& simulator,
311 const WellState<Scalar>& well_state,
312 DeferredLogger& deferred_logger);
313
314 bool updateWellOperabilityFromWellEq(const Simulator& simulator,
315 const WellState<Scalar>& well_state,
316 DeferredLogger& deferred_logger);
317
318 // update perforation water throughput based on solved water rate
319 virtual void updateWaterThroughput(const double dt,
320 WellState<Scalar>& well_state) const = 0;
321
324 virtual std::vector<Scalar>
326 DeferredLogger& deferred_logger) const = 0;
327
331 void initializeProducerWellState(const Simulator& simulator,
332 WellState<Scalar>& well_state,
333 DeferredLogger& deferred_logger) const;
334
335 void solveWellEquation(const Simulator& simulator,
336 WellState<Scalar>& well_state,
337 const GroupState<Scalar>& group_state,
338 DeferredLogger& deferred_logger);
339
340 const std::vector<RateVector>& connectionRates() const
341 {
342 return connectionRates_;
343 }
344
345 std::vector<Scalar> wellIndex(const int perf,
346 const IntensiveQuantities& intQuants,
347 const Scalar trans_mult,
348 const SingleWellState<Scalar>& ws) const;
349
350 void updateConnectionDFactor(const Simulator& simulator,
351 SingleWellState<Scalar>& ws) const;
352
354 SingleWellState<Scalar>& ws) const;
355
356 virtual bool iterateWellEqWithSwitching(const Simulator& simulator,
357 const double dt,
358 const WellInjectionControls& inj_controls,
359 const WellProductionControls& prod_controls,
360 WellState<Scalar>& well_state,
361 const GroupState<Scalar>& group_state,
362 DeferredLogger& deferred_logger,
363 const bool fixed_control = false,
364 const bool fixed_status = false) = 0;
365protected:
366 // simulation parameters
367 std::vector<RateVector> connectionRates_;
368 std::vector<Scalar> B_avg_;
371
372 Scalar wpolymer() const;
373 Scalar wfoam() const;
374 Scalar wsalt() const;
375 Scalar wmicrobes() const;
376 Scalar woxygen() const;
377 Scalar wurea() const;
378
379 virtual Scalar getRefDensity() const = 0;
380
381 std::vector<Scalar>
382 initialWellRateFractions(const Simulator& ebosSimulator,
383 const WellState<Scalar>& well_state) const;
384
385 // check whether the well is operable under BHP limit with current reservoir condition
386 virtual void checkOperabilityUnderBHPLimit(const WellState<Scalar>& well_state,
387 const Simulator& simulator,
388 DeferredLogger& deferred_logger) = 0;
389
390 // check whether the well is operable under THP limit with current reservoir condition
391 virtual void checkOperabilityUnderTHPLimit(const Simulator& simulator,
392 const WellState<Scalar>& well_state,
393 DeferredLogger& deferred_logger) = 0;
394
395 virtual void updateIPR(const Simulator& simulator,
396 DeferredLogger& deferred_logger) const = 0;
397
398 virtual void assembleWellEqWithoutIteration(const Simulator& simulator,
399 const double dt,
400 const WellInjectionControls& inj_controls,
401 const WellProductionControls& prod_controls,
402 WellState<Scalar>& well_state,
403 const GroupState<Scalar>& group_state,
404 DeferredLogger& deferred_logger) = 0;
405
406 // iterate well equations with the specified control until converged
407 virtual bool iterateWellEqWithControl(const Simulator& simulator,
408 const double dt,
409 const WellInjectionControls& inj_controls,
410 const WellProductionControls& prod_controls,
411 WellState<Scalar>& well_state,
412 const GroupState<Scalar>& group_state,
413 DeferredLogger& deferred_logger) = 0;
414
415 virtual void updateIPRImplicit(const Simulator& simulator,
416 WellState<Scalar>& well_state,
417 DeferredLogger& deferred_logger) = 0;
418
419 bool iterateWellEquations(const Simulator& simulator,
420 const double dt,
421 WellState<Scalar>& well_state,
422 const GroupState<Scalar>& group_state,
423 DeferredLogger& deferred_logger);
424
425 bool solveWellWithOperabilityCheck(const Simulator& simulator,
426 const double dt,
427 const Well::InjectionControls& inj_controls,
428 const Well::ProductionControls& prod_controls,
429 WellState<Scalar>& well_state,
430 const GroupState<Scalar>& group_state,
431 DeferredLogger& deferred_logger);
432
433 std::optional<Scalar>
434 estimateOperableBhp(const Simulator& ebos_simulator,
435 const double dt,
436 WellState<Scalar>& well_state,
437 const SummaryState& summary_state,
438 DeferredLogger& deferred_logger);
439
440 bool solveWellWithBhp(const Simulator& simulator,
441 const double dt,
442 const Scalar bhp,
443 WellState<Scalar>& well_state,
444 DeferredLogger& deferred_logger);
445
446 bool solveWellWithZeroRate(const Simulator& simulator,
447 const double dt,
448 WellState<Scalar>& well_state,
449 DeferredLogger& deferred_logger);
450
451 bool solveWellForTesting(const Simulator& simulator,
452 WellState<Scalar>& well_state,
453 const GroupState<Scalar>& group_state,
454 DeferredLogger& deferred_logger);
455
456
457 template<class GasLiftSingleWell>
458 std::unique_ptr<GasLiftSingleWell> initializeGliftWellTest_(const Simulator& simulator,
459 WellState<Scalar>& well_state,
460 const GroupState<Scalar>& group_state,
461 const PhaseUsage& phase_usage,
462 GLiftEclWells& ecl_well_map,
463 DeferredLogger& deferred_logger);
464
465 Eval getPerfCellPressure(const FluidState& fs) const;
466
467 // get the mobility for specific perforation
468 template<class Value, class Callback>
469 void getMobility(const Simulator& simulator,
470 const int local_perf_index,
471 std::vector<Value>& mob,
472 Callback& extendEval,
473 [[maybe_unused]] DeferredLogger& deferred_logger) const;
474
475 void computeConnLevelProdInd(const FluidState& fs,
476 const std::function<Scalar(const Scalar)>& connPICalc,
477 const std::vector<Scalar>& mobility,
478 Scalar* connPI) const;
479
480 void computeConnLevelInjInd(const FluidState& fs,
481 const Phase preferred_phase,
482 const std::function<Scalar(const Scalar)>& connIICalc,
483 const std::vector<Scalar>& mobility,
484 Scalar* connII,
485 DeferredLogger& deferred_logger) const;
486
487 Scalar computeConnectionDFactor(const int perf,
488 const IntensiveQuantities& intQuants,
489 const SingleWellState<Scalar>& ws) const;
490};
491
492} // namespace Opm
493
494#include "WellInterface_impl.hpp"
495
496#endif // OPM_WELLINTERFACE_HEADER_INCLUDED
Declares the properties required by the black oil model.
Definition: ConvergenceReport.hpp:38
Definition: DeferredLogger.hpp:57
std::map< std::string, std::pair< const Well *, int > > GLiftEclWells
Definition: GasLiftGroupInfo.hpp:70
Definition: GasLiftSingleWell.hpp:37
Definition: GroupState.hpp:43
Class encapsulating some information about parallel wells.
Definition: ParallelWellInfo.hpp:195
Definition: RateConverter.hpp:71
Manages the initializing and running of time dependent problems.
Definition: simulator.hh:84
Definition: SingleWellState.hpp:42
Definition: WellInterfaceFluidSystem.hpp:51
static constexpr int Oil
Definition: WellInterfaceFluidSystem.hpp:63
static constexpr int Water
Definition: WellInterfaceFluidSystem.hpp:62
static constexpr int Gas
Definition: WellInterfaceFluidSystem.hpp:64
int pvtRegionIdx() const
Definition: WellInterfaceGeneric.hpp:124
const std::vector< FluidSystem::Scalar > & wellIndex() const
Definition: WellInterfaceGeneric.hpp:146
Definition: WellInterfaceIndices.hpp:34
typename WellInterfaceFluidSystem< GetPropType< TypeTag, Properties::FluidSystem > >::ModelParameters ModelParameters
Definition: WellInterfaceIndices.hpp:41
DenseAd::Evaluation< Scalar, Indices::numEq > Eval
Definition: WellInterfaceIndices.hpp:40
Definition: WellInterface.hpp:77
bool updateWellControl(const Simulator &simulator, const IndividualOrGroup iog, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:191
virtual std::vector< Scalar > computeCurrentWellRates(const Simulator &simulator, DeferredLogger &deferred_logger) const =0
Scalar woxygen() const
Definition: WellInterface_impl.hpp:167
static constexpr bool has_brine
Definition: WellInterface.hpp:116
IndividualOrGroup
Definition: WellInterface.hpp:242
GetPropType< TypeTag, Properties::Simulator > Simulator
Definition: WellInterface.hpp:82
virtual ConvergenceReport getWellConvergence(const Simulator &simulator, const WellState< Scalar > &well_state, const std::vector< Scalar > &B_avg, DeferredLogger &deferred_logger, const bool relax_tolerance) const =0
WellInterface(const Well &well, const ParallelWellInfo< Scalar > &pw_info, const int time_step, const ModelParameters &param, const RateConverterType &rate_converter, const int pvtRegionIdx, const int num_components, const int num_phases, const int index_of_well, const std::vector< PerforationData< Scalar > > &perf_data)
Constructor.
Definition: WellInterface_impl.hpp:58
virtual void recoverWellSolutionAndUpdateWellState(const Simulator &simulator, const BVector &x, WellState< Scalar > &well_state, DeferredLogger &deferred_logger)=0
virtual void computeWellPotentials(const Simulator &simulator, const WellState< Scalar > &well_state, std::vector< Scalar > &well_potentials, DeferredLogger &deferred_logger)=0
void assembleWellEq(const Simulator &simulator, const double dt, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:841
virtual void computeWellRatesWithBhp(const Simulator &ebosSimulator, const Scalar &bhp, std::vector< Scalar > &well_flux, DeferredLogger &deferred_logger) const =0
virtual void addWellContributions(SparseMatrixAdapter &) const =0
typename WellInterfaceFluidSystem< FluidSystem >::RateConverterType RateConverterType
Definition: WellInterface.hpp:100
const std::vector< RateVector > & connectionRates() const
Definition: WellInterface.hpp:340
BlackOilFluidState< Eval, FluidSystem, has_temperature, has_energy, Indices::compositionSwitchIdx >=0, has_watVapor, has_brine, has_saltPrecip, has_disgas_in_water, Indices::numPhases > FluidState
Definition: WellInterface.hpp:132
virtual void assembleWellEqWithoutIteration(const Simulator &simulator, const double dt, const WellInjectionControls &inj_controls, const WellProductionControls &prod_controls, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)=0
Scalar wfoam() const
Definition: WellInterface_impl.hpp:129
virtual void updateIPR(const Simulator &simulator, DeferredLogger &deferred_logger) const =0
virtual void checkOperabilityUnderBHPLimit(const WellState< Scalar > &well_state, const Simulator &simulator, DeferredLogger &deferred_logger)=0
std::vector< RateVector > connectionRates_
Definition: WellInterface.hpp:367
bool wellUnderZeroRateTarget(const Simulator &simulator, const WellState< Scalar > &well_state, DeferredLogger &deferred_logger) const
Definition: WellInterface_impl.hpp:1521
void computeConnLevelProdInd(const FluidState &fs, const std::function< Scalar(const Scalar)> &connPICalc, const std::vector< Scalar > &mobility, Scalar *connPI) const
Definition: WellInterface_impl.hpp:2004
bool updateWellStateWithTHPTargetProd(const Simulator &simulator, WellState< Scalar > &well_state, DeferredLogger &deferred_logger) const
Definition: WellInterface_impl.hpp:1973
void addCellRates(RateVector &rates, int cellIdx) const
Definition: WellInterface_impl.hpp:953
virtual void calculateExplicitQuantities(const Simulator &simulator, const WellState< Scalar > &well_state, DeferredLogger &deferred_logger)=0
Dune::BCRSMatrix< Opm::MatrixBlock< Scalar, 1, 1 > > PressureMatrix
Definition: WellInterface.hpp:97
Scalar volumetricSurfaceRateForConnection(int cellIdx, int phaseIdx) const
Definition: WellInterface_impl.hpp:969
virtual void addWellPressureEquations(PressureMatrix &mat, const BVector &x, const int pressureVarIndex, const bool use_well_weights, const WellState< Scalar > &well_state) const =0
void getMobility(const Simulator &simulator, const int local_perf_index, std::vector< Value > &mob, Callback &extendEval, DeferredLogger &deferred_logger) const
Definition: WellInterface_impl.hpp:1895
virtual void apply(const BVector &x, BVector &Ax) const =0
Ax = Ax - C D^-1 B x.
static const bool has_temperature
Definition: WellInterface.hpp:112
virtual void updatePrimaryVariables(const Simulator &simulator, const WellState< Scalar > &well_state, DeferredLogger &deferred_logger)=0
bool iterateWellEquations(const Simulator &simulator, const double dt, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:503
GetPropType< TypeTag, Properties::IntensiveQuantities > IntensiveQuantities
Definition: WellInterface.hpp:86
virtual void updateProductivityIndex(const Simulator &simulator, const WellProdIndexCalculator< Scalar > &wellPICalc, WellState< Scalar > &well_state, DeferredLogger &deferred_logger) const =0
GetPropType< TypeTag, Properties::Scalar > Scalar
Definition: WellInterface.hpp:83
void checkWellOperability(const Simulator &simulator, const WellState< Scalar > &well_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:989
std::optional< Scalar > estimateOperableBhp(const Simulator &ebos_simulator, const double dt, WellState< Scalar > &well_state, const SummaryState &summary_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:666
void updateWellOperability(const Simulator &simulator, const WellState< Scalar > &well_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:1080
static constexpr bool has_watVapor
Definition: WellInterface.hpp:117
void wellTesting(const Simulator &simulator, const double simulation_time, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, WellTestState &welltest_state, const PhaseUsage &phase_usage, GLiftEclWells &ecl_well_map, std::map< std::string, double > &open_times, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:381
GetPropType< TypeTag, Properties::MaterialLaw > MaterialLaw
Definition: WellInterface.hpp:87
bool solveWellWithOperabilityCheck(const Simulator &simulator, const double dt, const Well::InjectionControls &inj_controls, const Well::ProductionControls &prod_controls, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:566
virtual void updateWellStateWithTarget(const Simulator &simulator, const GroupState< Scalar > &group_state, WellState< Scalar > &well_state, DeferredLogger &deferred_logger) const
Definition: WellInterface_impl.hpp:1136
virtual Scalar getRefDensity() const =0
Eval getPerfCellPressure(const FluidState &fs) const
Definition: WellInterface_impl.hpp:1880
Dune::BlockVector< VectorBlockType > BVector
Definition: WellInterface.hpp:96
std::vector< Scalar > B_avg_
Definition: WellInterface.hpp:368
bool changed_to_stopped_this_step_
Definition: WellInterface.hpp:369
virtual void solveEqAndUpdateWellState(const Simulator &simulator, WellState< Scalar > &well_state, DeferredLogger &deferred_logger)=0
static constexpr bool has_polymermw
Definition: WellInterface.hpp:114
static constexpr bool has_polymer
Definition: WellInterface.hpp:110
virtual void apply(BVector &r) const =0
r = r - C D^-1 Rw
GetPropType< TypeTag, Properties::Grid > Grid
Definition: WellInterface.hpp:81
void updateConnectionTransmissibilityFactor(const Simulator &simulator, SingleWellState< Scalar > &ws) const
Definition: WellInterface_impl.hpp:1851
virtual ~WellInterface()=default
Virtual destructor.
virtual std::optional< Scalar > computeBhpAtThpLimitProdWithAlq(const Simulator &ebos_simulator, const SummaryState &summary_state, const Scalar alq_value, DeferredLogger &deferred_logger, bool iterate_if_no_solution) const =0
bool stoppedOrZeroRateTarget(const Simulator &simulator, const WellState< Scalar > &well_state, DeferredLogger &deferred_logger) const
Definition: WellInterface_impl.hpp:1558
void initializeProducerWellState(const Simulator &simulator, WellState< Scalar > &well_state, DeferredLogger &deferred_logger) const
Definition: WellInterface_impl.hpp:1621
typename Base::ModelParameters ModelParameters
Definition: WellInterface.hpp:106
virtual void computeWellRatesWithBhpIterations(const Simulator &simulator, const Scalar &bhp, std::vector< Scalar > &well_flux, DeferredLogger &deferred_logger) const =0
virtual void updateIPRImplicit(const Simulator &simulator, WellState< Scalar > &well_state, DeferredLogger &deferred_logger)=0
GetPropType< TypeTag, Properties::FluidSystem > FluidSystem
Definition: WellInterface.hpp:84
static constexpr bool has_solvent
Definition: WellInterface.hpp:108
static constexpr bool has_disgas_in_water
Definition: WellInterface.hpp:118
void prepareWellBeforeAssembling(const Simulator &simulator, const double dt, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:877
static constexpr bool has_saltPrecip
Definition: WellInterface.hpp:119
static constexpr bool has_foam
Definition: WellInterface.hpp:115
GetPropType< TypeTag, Properties::RateVector > RateVector
Definition: WellInterface.hpp:89
virtual void updateWaterThroughput(const double dt, WellState< Scalar > &well_state) const =0
void computeConnLevelInjInd(const FluidState &fs, const Phase preferred_phase, const std::function< Scalar(const Scalar)> &connIICalc, const std::vector< Scalar > &mobility, Scalar *connII, DeferredLogger &deferred_logger) const
Definition: WellInterface_impl.hpp:2039
void solveWellEquation(const Simulator &simulator, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:779
virtual void init(const PhaseUsage *phase_usage_arg, const std::vector< Scalar > &depth_arg, const Scalar gravity_arg, const std::vector< Scalar > &B_avg, const bool changed_to_open_this_step)
Definition: WellInterface_impl.hpp:95
void updateConnectionDFactor(const Simulator &simulator, SingleWellState< Scalar > &ws) const
Definition: WellInterface_impl.hpp:1789
virtual bool iterateWellEqWithControl(const Simulator &simulator, const double dt, const WellInjectionControls &inj_controls, const WellProductionControls &prod_controls, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)=0
std::vector< Scalar > initialWellRateFractions(const Simulator &ebosSimulator, const WellState< Scalar > &well_state) const
Definition: WellInterface_impl.hpp:1571
Scalar wsalt() const
Definition: WellInterface_impl.hpp:143
static constexpr bool has_micp
Definition: WellInterface.hpp:120
bool updateWellOperabilityFromWellEq(const Simulator &simulator, const WellState< Scalar > &well_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:1117
virtual bool iterateWellEqWithSwitching(const Simulator &simulator, const double dt, const WellInjectionControls &inj_controls, const WellProductionControls &prod_controls, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger, const bool fixed_control=false, const bool fixed_status=false)=0
Dune::FieldMatrix< Scalar, Indices::numEq, Indices::numEq > MatrixBlockType
Definition: WellInterface.hpp:94
void gliftBeginTimeStepWellTestUpdateALQ(const Simulator &simulator, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, const PhaseUsage &phase_usage, GLiftEclWells &ecl_well_map, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:1016
typename GasLiftGroupInfo< Scalar >::GLiftEclWells GLiftEclWells
Definition: WellInterface.hpp:91
typename Base::Eval Eval
Definition: WellInterface.hpp:95
static constexpr bool has_energy
Definition: WellInterface.hpp:111
Scalar wpolymer() const
Definition: WellInterface_impl.hpp:113
Scalar computeConnectionDFactor(const int perf, const IntensiveQuantities &intQuants, const SingleWellState< Scalar > &ws) const
Definition: WellInterface_impl.hpp:1809
bool solveWellWithBhp(const Simulator &simulator, const double dt, const Scalar bhp, WellState< Scalar > &well_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:699
bool solveWellWithZeroRate(const Simulator &simulator, const double dt, WellState< Scalar > &well_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:736
bool solveWellForTesting(const Simulator &simulator, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:757
GetPropType< TypeTag, Properties::Indices > Indices
Definition: WellInterface.hpp:85
Scalar wurea() const
Definition: WellInterface_impl.hpp:179
bool wellUnderZeroGroupRateTarget(const Simulator &simulator, const WellState< Scalar > &well_state, DeferredLogger &deferred_logger, std::optional< bool > group_control=std::nullopt) const
Definition: WellInterface_impl.hpp:1539
bool thp_update_iterations
Definition: WellInterface.hpp:370
Dune::FieldVector< Scalar, Indices::numEq > VectorBlockType
Definition: WellInterface.hpp:93
Scalar wmicrobes() const
Definition: WellInterface_impl.hpp:155
void assembleWellEqWithoutIteration(const Simulator &simulator, const double dt, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:857
GetPropType< TypeTag, Properties::SparseMatrixAdapter > SparseMatrixAdapter
Definition: WellInterface.hpp:88
bool updateWellControlAndStatusLocalIteration(const Simulator &simulator, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, const Well::InjectionControls &inj_controls, const Well::ProductionControls &prod_controls, const Scalar WQTotal, DeferredLogger &deferred_logger, const bool fixed_control=false, const bool fixed_status=false)
Definition: WellInterface_impl.hpp:274
std::unique_ptr< GasLiftSingleWell > initializeGliftWellTest_(const Simulator &simulator, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, const PhaseUsage &phase_usage, GLiftEclWells &ecl_well_map, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:2075
virtual void checkOperabilityUnderTHPLimit(const Simulator &simulator, const WellState< Scalar > &well_state, DeferredLogger &deferred_logger)=0
static constexpr bool has_zFraction
Definition: WellInterface.hpp:109
Definition: WellProdIndexCalculator.hpp:37
Definition: WellState.hpp:65
Declares the properties required by the black oil model.
Definition: blackoilboundaryratevector.hh:39
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition: propertysystem.hh:233
Static data associated with a well perforation.
Definition: PerforationData.hpp:30
Definition: BlackoilPhases.hpp:46