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
58
59#include <dune/common/fmatrix.hh>
60#include <dune/istl/bcrsmatrix.hh>
61#include <dune/istl/matrixmatrix.hh>
62
63#include <opm/material/densead/Evaluation.hpp>
64
65#include <vector>
66
67namespace Opm
68{
69
70class WellInjectionProperties;
71class WellProductionProperties;
72
73template<typename TypeTag>
74class WellInterface : public WellInterfaceIndices<GetPropType<TypeTag, Properties::FluidSystem>,
75 GetPropType<TypeTag, Properties::Indices>>
76{
79public:
84 using IndexTraits = typename FluidSystem::IndexTraitsType;
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
104
108
110
111 static constexpr bool has_solvent = getPropValue<TypeTag, Properties::EnableSolvent>();
112 static constexpr bool has_zFraction = getPropValue<TypeTag, Properties::EnableExtbo>();
113 static constexpr bool has_polymer = getPropValue<TypeTag, Properties::EnablePolymer>();
114 static constexpr bool has_energy = getPropValue<TypeTag, Properties::EnableEnergy>();
115 static const bool has_temperature = getPropValue<TypeTag, Properties::EnableTemperature>();
116 // flag for polymer molecular weight related
117 static constexpr bool has_polymermw = getPropValue<TypeTag, Properties::EnablePolymerMW>();
118 static constexpr bool has_foam = getPropValue<TypeTag, Properties::EnableFoam>();
119 static constexpr bool has_brine = getPropValue<TypeTag, Properties::EnableBrine>();
120 static constexpr bool has_watVapor = getPropValue<TypeTag, Properties::EnableVapwat>();
121 static constexpr bool has_disgas_in_water = getPropValue<TypeTag, Properties::EnableDisgasInWater>();
122 static constexpr bool has_saltPrecip = getPropValue<TypeTag, Properties::EnableSaltPrecipitation>();
123 static constexpr bool has_micp = getPropValue<TypeTag, Properties::EnableMICP>();
124
125 // For the conversion between the surface volume rate and reservoir voidage rate
126 using FluidState = BlackOilFluidState<Eval,
130 Indices::compositionSwitchIdx >= 0,
132 has_brine,
135 Indices::numPhases >;
137 WellInterface(const Well& well,
138 const ParallelWellInfo<Scalar>& pw_info,
139 const int time_step,
140 const ModelParameters& param,
141 const RateConverterType& rate_converter,
142 const int pvtRegionIdx,
143 const int num_conservation_quantities,
144 const int num_phases,
145 const int index_of_well,
146 const std::vector<PerforationData<Scalar>>& perf_data);
147
149 virtual ~WellInterface() = default;
150
151 virtual void init(const std::vector<Scalar>& depth_arg,
152 const Scalar gravity_arg,
153 const std::vector<Scalar>& B_avg,
154 const bool changed_to_open_this_step);
155
157 const WellStateType& well_state,
158 const std::vector<Scalar>& B_avg,
159 DeferredLogger& deferred_logger,
160 const bool relax_tolerance) const = 0;
161
162 virtual void solveEqAndUpdateWellState(const Simulator& simulator,
163 WellStateType& well_state,
164 DeferredLogger& deferred_logger) = 0;
165
166 void assembleWellEq(const Simulator& simulator,
167 const double dt,
168 WellStateType& well_state,
169 const GroupState<Scalar>& group_state,
170 DeferredLogger& deferred_logger);
171
172 void assembleWellEqWithoutIteration(const Simulator& simulator,
173 const double dt,
174 WellStateType& well_state,
175 const GroupState<Scalar>& group_state,
176 DeferredLogger& deferred_logger);
177
178 // TODO: better name or further refactoring the function to make it more clear
179 void prepareWellBeforeAssembling(const Simulator& simulator,
180 const double dt,
181 WellStateType& well_state,
182 const GroupState<Scalar>& group_state,
183 DeferredLogger& deferred_logger);
184
185
186 virtual void computeWellRatesWithBhp(const Simulator& ebosSimulator,
187 const Scalar& bhp,
188 std::vector<Scalar>& well_flux,
189 DeferredLogger& deferred_logger) const = 0;
190
191 virtual std::optional<Scalar>
193 const SummaryState& summary_state,
194 const Scalar alq_value,
195 DeferredLogger& deferred_logger,
196 bool iterate_if_no_solution) const = 0;
197
201 const BVector& x,
202 WellStateType& well_state,
203 DeferredLogger& deferred_logger) = 0;
204
206 virtual void apply(const BVector& x, BVector& Ax) const = 0;
207
209 virtual void apply(BVector& r) const = 0;
210
211 // TODO: before we decide to put more information under mutable, this function is not const
212 virtual void computeWellPotentials(const Simulator& simulator,
213 const WellStateType& well_state,
214 std::vector<Scalar>& well_potentials,
215 DeferredLogger& deferred_logger) = 0;
216
217 virtual void updateWellStateWithTarget(const Simulator& simulator,
218 const GroupState<Scalar>& group_state,
219 WellStateType& well_state,
220 DeferredLogger& deferred_logger) const;
221
222 virtual void computeWellRatesWithBhpIterations(const Simulator& simulator,
223 const Scalar& bhp,
224 std::vector<Scalar>& well_flux,
225 DeferredLogger& deferred_logger) const = 0;
226
227 bool wellUnderZeroRateTarget(const Simulator& simulator,
228 const WellStateType& well_state,
229 DeferredLogger& deferred_logger) const;
230
231 bool wellUnderZeroGroupRateTarget(const Simulator& simulator,
232 const WellStateType& well_state,
233 DeferredLogger& deferred_logger,
234 std::optional<bool> group_control = std::nullopt) const;
235
236 bool stoppedOrZeroRateTarget(const Simulator& simulator,
237 const WellStateType& well_state,
238 DeferredLogger& deferred_logger) const;
239
240 bool updateWellStateWithTHPTargetProd(const Simulator& simulator,
241 WellStateType& well_state,
242 DeferredLogger& deferred_logger) const;
243
244 enum class IndividualOrGroup { Individual, Group, Both };
245 bool updateWellControl(const Simulator& simulator,
246 const IndividualOrGroup iog,
247 WellStateType& well_state,
248 const GroupState<Scalar>& group_state,
249 DeferredLogger& deferred_logger) /* const */;
250
252 WellStateType& well_state,
253 const GroupState<Scalar>& group_state,
254 const Well::InjectionControls& inj_controls,
255 const Well::ProductionControls& prod_controls,
256 const Scalar WQTotal,
257 DeferredLogger& deferred_logger,
258 const bool fixed_control = false,
259 const bool fixed_status = false);
260
261 virtual void updatePrimaryVariables(const Simulator& simulator,
262 const WellStateType& well_state,
263 DeferredLogger& deferred_logger) = 0;
264
265 virtual void calculateExplicitQuantities(const Simulator& simulator,
266 const WellStateType& well_state,
267 DeferredLogger& deferred_logger) = 0; // should be const?
268
269 virtual void updateProductivityIndex(const Simulator& simulator,
270 const WellProdIndexCalculator<Scalar>& wellPICalc,
271 WellStateType& well_state,
272 DeferredLogger& deferred_logger) const = 0;
273
274 // Add well contributions to matrix
276
278 const BVector& x,
279 const int pressureVarIndex,
280 const bool use_well_weights,
281 const WellStateType& well_state) const = 0;
282
283 void addCellRates(RateVector& rates, int cellIdx) const;
284
285 Scalar volumetricSurfaceRateForConnection(int cellIdx, int phaseIdx) const;
286
287 // TODO: theoretically, it should be a const function
288 // Simulator is not const is because that assembleWellEq is non-const Simulator
289 void wellTesting(const Simulator& simulator,
290 const double simulation_time,
291 /* const */ WellStateType& well_state,
292 const GroupState<Scalar>& group_state,
293 WellTestState& welltest_state,
294 GLiftEclWells& ecl_well_map,
295 std::map<std::string, double>& open_times,
296 DeferredLogger& deferred_logger);
297
298 void checkWellOperability(const Simulator& simulator,
299 const WellStateType& well_state,
300 DeferredLogger& deferred_logger);
301
303 WellStateType& well_state,
304 const GroupState<Scalar>& group_state,
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 WellStateType& well_state,
312 DeferredLogger& deferred_logger);
313
314 bool updateWellOperabilityFromWellEq(const Simulator& simulator,
315 const WellStateType& 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 WellStateType& well_state) const = 0;
321
324 virtual std::vector<Scalar>
326 DeferredLogger& deferred_logger) const = 0;
327
331 void initializeProducerWellState(const Simulator& simulator,
332 WellStateType& well_state,
333 DeferredLogger& deferred_logger) const;
334
335 void solveWellEquation(const Simulator& simulator,
336 WellStateType& 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 SingleWellStateType& ws) const;
349
350 void updateConnectionDFactor(const Simulator& simulator,
351 SingleWellStateType& ws) const;
352
354 SingleWellStateType& 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 WellStateType& 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 WellStateType& well_state) const;
384
385 // check whether the well is operable under BHP limit with current reservoir condition
386 virtual void checkOperabilityUnderBHPLimit(const WellStateType& 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 WellStateType& 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 WellStateType& 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 WellStateType& well_state,
412 const GroupState<Scalar>& group_state,
413 DeferredLogger& deferred_logger) = 0;
414
415 virtual void updateIPRImplicit(const Simulator& simulator,
416 WellStateType& well_state,
417 DeferredLogger& deferred_logger) = 0;
418
419 bool iterateWellEquations(const Simulator& simulator,
420 const double dt,
421 WellStateType& 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 WellStateType& 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 WellStateType& 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 WellStateType& well_state,
444 DeferredLogger& deferred_logger);
445
446 bool solveWellWithZeroRate(const Simulator& simulator,
447 const double dt,
448 WellStateType& well_state,
449 DeferredLogger& deferred_logger);
450
451 bool solveWellForTesting(const Simulator& simulator,
452 WellStateType& 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 WellStateType& well_state,
460 const GroupState<Scalar>& group_state,
461 GLiftEclWells& ecl_well_map,
462 DeferredLogger& deferred_logger);
463
464 Eval getPerfCellPressure(const FluidState& fs) const;
465
466 // get the mobility for specific perforation
467 template<class Value, class Callback>
468 void getMobility(const Simulator& simulator,
469 const int local_perf_index,
470 std::vector<Value>& mob,
471 Callback& extendEval,
472 [[maybe_unused]] DeferredLogger& deferred_logger) const;
473
474 void computeConnLevelProdInd(const FluidState& fs,
475 const std::function<Scalar(const Scalar)>& connPICalc,
476 const std::vector<Scalar>& mobility,
477 Scalar* connPI) const;
478
479 void computeConnLevelInjInd(const FluidState& fs,
480 const Phase preferred_phase,
481 const std::function<Scalar(const Scalar)>& connIICalc,
482 const std::vector<Scalar>& mobility,
483 Scalar* connII,
484 DeferredLogger& deferred_logger) const;
485
486 Scalar computeConnectionDFactor(const int perf,
487 const IntensiveQuantities& intQuants,
488 const SingleWellStateType& ws) const;
489};
490
491} // namespace Opm
492
493#include "WellInterface_impl.hpp"
494
495#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:64
Definition: GasLiftSingleWell.hpp:39
Definition: GroupState.hpp:41
Class encapsulating some information about parallel wells.
Definition: ParallelWellInfo.hpp:198
Definition: RateConverter.hpp:70
Manages the initializing and running of time dependent problems.
Definition: simulator.hh:84
Definition: SingleWellState.hpp:43
Definition: WellInterfaceFluidSystem.hpp:49
static constexpr int Oil
Definition: WellInterfaceFluidSystem.hpp:62
static constexpr int Water
Definition: WellInterfaceFluidSystem.hpp:61
static constexpr int Gas
Definition: WellInterfaceFluidSystem.hpp:63
int pvtRegionIdx() const
Definition: WellInterfaceGeneric.hpp:127
const std::vector< FluidSystem::Scalar > & wellIndex() const
Definition: WellInterfaceGeneric.hpp:149
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:76
bool solveWellWithBhp(const Simulator &simulator, const double dt, const Scalar bhp, WellStateType &well_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:696
virtual std::vector< Scalar > computeCurrentWellRates(const Simulator &simulator, DeferredLogger &deferred_logger) const =0
void checkWellOperability(const Simulator &simulator, const WellStateType &well_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:986
Scalar woxygen() const
Definition: WellInterface_impl.hpp:165
static constexpr bool has_brine
Definition: WellInterface.hpp:119
IndividualOrGroup
Definition: WellInterface.hpp:244
GetPropType< TypeTag, Properties::Simulator > Simulator
Definition: WellInterface.hpp:81
Scalar computeConnectionDFactor(const int perf, const IntensiveQuantities &intQuants, const SingleWellStateType &ws) const
Definition: WellInterface_impl.hpp:1813
void solveWellEquation(const Simulator &simulator, WellStateType &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:776
virtual void computeWellRatesWithBhp(const Simulator &ebosSimulator, const Scalar &bhp, std::vector< Scalar > &well_flux, DeferredLogger &deferred_logger) const =0
bool updateWellControlAndStatusLocalIteration(const Simulator &simulator, WellStateType &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:272
virtual void addWellContributions(SparseMatrixAdapter &) const =0
typename WellInterfaceFluidSystem< FluidSystem >::RateConverterType RateConverterType
Definition: WellInterface.hpp:103
bool iterateWellEquations(const Simulator &simulator, const double dt, WellStateType &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:500
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:135
Scalar wfoam() const
Definition: WellInterface_impl.hpp:127
virtual void updateIPR(const Simulator &simulator, DeferredLogger &deferred_logger) const =0
std::optional< Scalar > estimateOperableBhp(const Simulator &ebos_simulator, const double dt, WellStateType &well_state, const SummaryState &summary_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:663
std::vector< RateVector > connectionRates_
Definition: WellInterface.hpp:367
void computeConnLevelProdInd(const FluidState &fs, const std::function< Scalar(const Scalar)> &connPICalc, const std::vector< Scalar > &mobility, Scalar *connPI) const
Definition: WellInterface_impl.hpp:2008
void addCellRates(RateVector &rates, int cellIdx) const
Definition: WellInterface_impl.hpp:950
Dune::BCRSMatrix< Opm::MatrixBlock< Scalar, 1, 1 > > PressureMatrix
Definition: WellInterface.hpp:97
void gliftBeginTimeStepWellTestUpdateALQ(const Simulator &simulator, WellStateType &well_state, const GroupState< Scalar > &group_state, GLiftEclWells &ecl_well_map, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:1013
void assembleWellEq(const Simulator &simulator, const double dt, WellStateType &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:838
Scalar volumetricSurfaceRateForConnection(int cellIdx, int phaseIdx) const
Definition: WellInterface_impl.hpp:966
virtual void init(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 getMobility(const Simulator &simulator, const int local_perf_index, std::vector< Value > &mob, Callback &extendEval, DeferredLogger &deferred_logger) const
Definition: WellInterface_impl.hpp:1899
void assembleWellEqWithoutIteration(const Simulator &simulator, const double dt, WellStateType &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:854
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:115
GetPropType< TypeTag, Properties::IntensiveQuantities > IntensiveQuantities
Definition: WellInterface.hpp:86
GetPropType< TypeTag, Properties::Scalar > Scalar
Definition: WellInterface.hpp:82
virtual void checkOperabilityUnderTHPLimit(const Simulator &simulator, const WellStateType &well_state, DeferredLogger &deferred_logger)=0
std::vector< Scalar > initialWellRateFractions(const Simulator &ebosSimulator, const WellStateType &well_state) const
Definition: WellInterface_impl.hpp:1575
bool updateWellStateWithTHPTargetProd(const Simulator &simulator, WellStateType &well_state, DeferredLogger &deferred_logger) const
Definition: WellInterface_impl.hpp:1977
virtual bool iterateWellEqWithControl(const Simulator &simulator, const double dt, const WellInjectionControls &inj_controls, const WellProductionControls &prod_controls, WellStateType &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)=0
static constexpr bool has_watVapor
Definition: WellInterface.hpp:120
void updateConnectionDFactor(const Simulator &simulator, SingleWellStateType &ws) const
Definition: WellInterface_impl.hpp:1793
virtual void updateWaterThroughput(const double dt, WellStateType &well_state) const =0
GetPropType< TypeTag, Properties::MaterialLaw > MaterialLaw
Definition: WellInterface.hpp:87
virtual Scalar getRefDensity() const =0
typename FluidSystem::IndexTraitsType IndexTraits
Definition: WellInterface.hpp:84
virtual void computeWellPotentials(const Simulator &simulator, const WellStateType &well_state, std::vector< Scalar > &well_potentials, DeferredLogger &deferred_logger)=0
Eval getPerfCellPressure(const FluidState &fs) const
Definition: WellInterface_impl.hpp:1884
Dune::BlockVector< VectorBlockType > BVector
Definition: WellInterface.hpp:96
void initializeProducerWellState(const Simulator &simulator, WellStateType &well_state, DeferredLogger &deferred_logger) const
Definition: WellInterface_impl.hpp:1625
std::vector< Scalar > B_avg_
Definition: WellInterface.hpp:368
bool changed_to_stopped_this_step_
Definition: WellInterface.hpp:369
static constexpr bool has_polymermw
Definition: WellInterface.hpp:117
virtual void calculateExplicitQuantities(const Simulator &simulator, const WellStateType &well_state, DeferredLogger &deferred_logger)=0
static constexpr bool has_polymer
Definition: WellInterface.hpp:113
virtual void apply(BVector &r) const =0
r = r - C D^-1 Rw
GetPropType< TypeTag, Properties::Grid > Grid
Definition: WellInterface.hpp:80
bool wellUnderZeroGroupRateTarget(const Simulator &simulator, const WellStateType &well_state, DeferredLogger &deferred_logger, std::optional< bool > group_control=std::nullopt) const
Definition: WellInterface_impl.hpp:1543
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
typename Base::ModelParameters ModelParameters
Definition: WellInterface.hpp:109
virtual void updatePrimaryVariables(const Simulator &simulator, const WellStateType &well_state, DeferredLogger &deferred_logger)=0
bool updateWellControl(const Simulator &simulator, const IndividualOrGroup iog, WellStateType &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:189
virtual void computeWellRatesWithBhpIterations(const Simulator &simulator, const Scalar &bhp, std::vector< Scalar > &well_flux, DeferredLogger &deferred_logger) const =0
GetPropType< TypeTag, Properties::FluidSystem > FluidSystem
Definition: WellInterface.hpp:83
static constexpr bool has_solvent
Definition: WellInterface.hpp:111
virtual void assembleWellEqWithoutIteration(const Simulator &simulator, const double dt, const WellInjectionControls &inj_controls, const WellProductionControls &prod_controls, WellStateType &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)=0
static constexpr bool has_disgas_in_water
Definition: WellInterface.hpp:121
static constexpr bool has_saltPrecip
Definition: WellInterface.hpp:122
virtual void checkOperabilityUnderBHPLimit(const WellStateType &well_state, const Simulator &simulator, DeferredLogger &deferred_logger)=0
static constexpr bool has_foam
Definition: WellInterface.hpp:118
GetPropType< TypeTag, Properties::RateVector > RateVector
Definition: WellInterface.hpp:89
virtual void recoverWellSolutionAndUpdateWellState(const Simulator &simulator, const BVector &x, WellStateType &well_state, DeferredLogger &deferred_logger)=0
void updateConnectionTransmissibilityFactor(const Simulator &simulator, SingleWellStateType &ws) const
Definition: WellInterface_impl.hpp:1855
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:2042
typename GasLiftGroupInfo< Scalar, IndexTraits >::GLiftEclWells GLiftEclWells
Definition: WellInterface.hpp:91
std::unique_ptr< GasLiftSingleWell > initializeGliftWellTest_(const Simulator &simulator, WellStateType &well_state, const GroupState< Scalar > &group_state, GLiftEclWells &ecl_well_map, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:2076
virtual void addWellPressureEquations(PressureMatrix &mat, const BVector &x, const int pressureVarIndex, const bool use_well_weights, const WellStateType &well_state) const =0
virtual void updateIPRImplicit(const Simulator &simulator, WellStateType &well_state, DeferredLogger &deferred_logger)=0
virtual void solveEqAndUpdateWellState(const Simulator &simulator, WellStateType &well_state, DeferredLogger &deferred_logger)=0
void wellTesting(const Simulator &simulator, const double simulation_time, WellStateType &well_state, const GroupState< Scalar > &group_state, WellTestState &welltest_state, GLiftEclWells &ecl_well_map, std::map< std::string, double > &open_times, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:380
Scalar wsalt() const
Definition: WellInterface_impl.hpp:141
static constexpr bool has_micp
Definition: WellInterface.hpp:123
Dune::FieldMatrix< Scalar, Indices::numEq, Indices::numEq > MatrixBlockType
Definition: WellInterface.hpp:94
void prepareWellBeforeAssembling(const Simulator &simulator, const double dt, WellStateType &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:874
typename Base::Eval Eval
Definition: WellInterface.hpp:95
bool wellUnderZeroRateTarget(const Simulator &simulator, const WellStateType &well_state, DeferredLogger &deferred_logger) const
Definition: WellInterface_impl.hpp:1525
virtual ConvergenceReport getWellConvergence(const Simulator &simulator, const WellStateType &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_conservation_quantities, const int num_phases, const int index_of_well, const std::vector< PerforationData< Scalar > > &perf_data)
Constructor.
Definition: WellInterface_impl.hpp:58
WellState< Scalar, IndexTraits > WellStateType
Definition: WellInterface.hpp:99
bool solveWellForTesting(const Simulator &simulator, WellStateType &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:754
void updateWellOperability(const Simulator &simulator, const WellStateType &well_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:1075
static constexpr bool has_energy
Definition: WellInterface.hpp:114
bool solveWellWithOperabilityCheck(const Simulator &simulator, const double dt, const Well::InjectionControls &inj_controls, const Well::ProductionControls &prod_controls, WellStateType &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:563
Scalar wpolymer() const
Definition: WellInterface_impl.hpp:111
SingleWellState< Scalar, IndexTraits > SingleWellStateType
Definition: WellInterface.hpp:100
virtual void updateProductivityIndex(const Simulator &simulator, const WellProdIndexCalculator< Scalar > &wellPICalc, WellStateType &well_state, DeferredLogger &deferred_logger) const =0
bool updateWellOperabilityFromWellEq(const Simulator &simulator, const WellStateType &well_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:1112
GetPropType< TypeTag, Properties::Indices > Indices
Definition: WellInterface.hpp:85
Scalar wurea() const
Definition: WellInterface_impl.hpp:177
virtual bool iterateWellEqWithSwitching(const Simulator &simulator, const double dt, const WellInjectionControls &inj_controls, const WellProductionControls &prod_controls, WellStateType &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger, const bool fixed_control=false, const bool fixed_status=false)=0
bool solveWellWithZeroRate(const Simulator &simulator, const double dt, WellStateType &well_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:733
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:153
GetPropType< TypeTag, Properties::SparseMatrixAdapter > SparseMatrixAdapter
Definition: WellInterface.hpp:88
virtual void updateWellStateWithTarget(const Simulator &simulator, const GroupState< Scalar > &group_state, WellStateType &well_state, DeferredLogger &deferred_logger) const
Definition: WellInterface_impl.hpp:1131
static constexpr bool has_zFraction
Definition: WellInterface.hpp:112
bool stoppedOrZeroRateTarget(const Simulator &simulator, const WellStateType &well_state, DeferredLogger &deferred_logger) const
Definition: WellInterface_impl.hpp:1562
Definition: WellProdIndexCalculator.hpp:37
Definition: WellState.hpp:66
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