StandardWell.hpp
Go to the documentation of this file.
1/*
2 Copyright 2017 SINTEF Digital, Mathematics and Cybernetics.
3 Copyright 2017 Statoil ASA.
4 Copyright 2016 - 2017 IRIS AS.
5
6 This file is part of the Open Porous Media project (OPM).
7
8 OPM is free software: you can redistribute it and/or modify
9 it under the terms of the GNU General Public License as published by
10 the Free Software Foundation, either version 3 of the License, or
11 (at your option) any later version.
12
13 OPM is distributed in the hope that it will be useful,
14 but WITHOUT ANY WARRANTY; without even the implied warranty of
15 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16 GNU General Public License for more details.
17
18 You should have received a copy of the GNU General Public License
19 along with OPM. If not, see <http://www.gnu.org/licenses/>.
20*/
21
22
23#ifndef OPM_STANDARDWELL_HEADER_INCLUDED
24#define OPM_STANDARDWELL_HEADER_INCLUDED
25
34
41
42#include <opm/material/densead/Evaluation.hpp>
43#include <opm/input/eclipse/Schedule/ScheduleTypes.hpp>
44
46
47#include <dune/common/dynvector.hh>
48#include <dune/common/dynmatrix.hh>
49
50#include <memory>
51#include <optional>
52
53namespace Opm
54{
55
56 template<typename TypeTag>
57 class StandardWell : public WellInterface<TypeTag>
58 , public StandardWellEval<GetPropType<TypeTag, Properties::FluidSystem>,
59 GetPropType<TypeTag, Properties::Indices>>
60 {
61
62 public:
66
67 // TODO: some functions working with AD variables handles only with values (double) without
68 // dealing with derivatives. It can be beneficial to make functions can work with either AD or scalar value.
69 // And also, it can also be beneficial to make these functions hanle different types of AD variables.
70 using typename Base::Simulator;
71 using typename Base::IntensiveQuantities;
72 using typename Base::FluidSystem;
73 using typename Base::MaterialLaw;
74 using typename Base::ModelParameters;
75 using typename Base::Indices;
76 using typename Base::RateConverterType;
77 using typename Base::SparseMatrixAdapter;
78 using typename Base::FluidState;
79 using typename Base::RateVector;
80
85 using Base::has_foam;
86 using Base::has_brine;
87 using Base::has_energy;
88 using Base::has_micp;
89
92 using typename Base::PressureMatrix;
93
94 // number of the conservation equations
95 static constexpr int numWellConservationEq = Indices::numPhases + Indices::numSolvents;
96 // number of the well control equations
97 static constexpr int numWellControlEq = 1;
98 // number of the well equations that will always be used
99 // based on the solution strategy, there might be other well equations be introduced
101
102 // the index for Bhp in primary variables and also the index of well control equation
103 // they both will be the last one in their respective system.
104 // TODO: we should have indices for the well equations and well primary variables separately
105 static constexpr int Bhp = numStaticWellEq - numWellControlEq;
106
108
109 using typename Base::Scalar;
110
111 using Base::name;
112 using Base::Water;
113 using Base::Oil;
114 using Base::Gas;
115
116 using typename Base::BVector;
117
118 using Eval = typename StdWellEval::Eval;
121
122 using IndexTraits = typename FluidSystem::IndexTraitsType;
124
125 StandardWell(const Well& well,
126 const ParallelWellInfo<Scalar>& pw_info,
127 const int time_step,
128 const ModelParameters& param,
129 const RateConverterType& rate_converter,
130 const int pvtRegionIdx,
131 const int num_conservation_quantities,
132 const int num_phases,
133 const int index_of_well,
134 const std::vector<PerforationData<Scalar>>& perf_data);
135
136 virtual void init(const std::vector<Scalar>& depth_arg,
137 const Scalar gravity_arg,
138 const std::vector<Scalar>& B_avg,
139 const bool changed_to_open_this_step) override;
140
142 virtual ConvergenceReport getWellConvergence(const Simulator& simulator,
143 const WellStateType& well_state,
144 const std::vector<Scalar>& B_avg,
145 DeferredLogger& deferred_logger,
146 const bool relax_tolerance) const override;
147
149 virtual void apply(const BVector& x, BVector& Ax) const override;
151 virtual void apply(BVector& r) const override;
152
156 const BVector& x,
157 WellStateType& well_state,
158 DeferredLogger& deferred_logger) override;
159
161 void computeWellPotentials(const Simulator& simulator,
162 const WellStateType& well_state,
163 std::vector<Scalar>& well_potentials,
164 DeferredLogger& deferred_logger) /* const */ override;
165
166 void updatePrimaryVariables(const Simulator& simulator,
167 const WellStateType& well_state,
168 DeferredLogger& deferred_logger) override;
169
170 void solveEqAndUpdateWellState(const Simulator& simulator,
171 WellStateType& well_state,
172 DeferredLogger& deferred_logger) override;
173
174 void calculateExplicitQuantities(const Simulator& simulator,
175 const WellStateType& well_state,
176 DeferredLogger& deferred_logger) override; // should be const?
177
178 void updateProductivityIndex(const Simulator& simulator,
179 const WellProdIndexCalculator<Scalar>& wellPICalc,
180 WellStateType& well_state,
181 DeferredLogger& deferred_logger) const override;
182
183 Scalar connectionDensity(const int globalConnIdx,
184 const int openConnIdx) const override;
185
186 void addWellContributions(SparseMatrixAdapter& mat) const override;
187
189 const BVector& x,
190 const int pressureVarIndex,
191 const bool use_well_weights,
192 const WellStateType& well_state) const override;
193
194 // iterate well equations with the specified control until converged
195 bool iterateWellEqWithControl(const Simulator& simulator,
196 const double dt,
197 const Well::InjectionControls& inj_controls,
198 const Well::ProductionControls& prod_controls,
199 WellStateType& well_state,
200 const GroupState<Scalar>& group_state,
201 DeferredLogger& deferred_logger) override;
202
203 // iterate well equations including control switching
204 bool iterateWellEqWithSwitching(const Simulator& simulator,
205 const double dt,
206 const Well::InjectionControls& inj_controls,
207 const Well::ProductionControls& prod_controls,
208 WellStateType& well_state,
209 const GroupState<Scalar>& group_state,
210 DeferredLogger& deferred_logger,
211 const bool fixed_control = false,
212 const bool fixed_status = false) override;
213
214 /* returns BHP */
216 const SummaryState &summary_state,
217 DeferredLogger& deferred_logger,
218 std::vector<Scalar>& potentials,
219 Scalar alq) const;
220
221 void computeWellRatesWithThpAlqProd(const Simulator& ebos_simulator,
222 const SummaryState& summary_state,
223 DeferredLogger& deferred_logger,
224 std::vector<Scalar>& potentials,
225 Scalar alq) const;
226
227 std::optional<Scalar>
228 computeBhpAtThpLimitProdWithAlq(const Simulator& ebos_simulator,
229 const SummaryState& summary_state,
230 const Scalar alq_value,
231 DeferredLogger& deferred_logger,
232 bool iterate_if_no_solution) const override;
233
234 void updateIPRImplicit(const Simulator& simulator,
235 WellStateType& well_state,
236 DeferredLogger& deferred_logger) override;
237
238 void computeWellRatesWithBhp(const Simulator& ebosSimulator,
239 const Scalar& bhp,
240 std::vector<Scalar>& well_flux,
241 DeferredLogger& deferred_logger) const override;
242
243 // NOTE: These cannot be protected since they are used by GasLiftRuntime
245
246 std::vector<Scalar>
247 computeCurrentWellRates(const Simulator& ebosSimulator,
248 DeferredLogger& deferred_logger) const override;
249
250 std::vector<Scalar> getPrimaryVars() const override;
251
252 int setPrimaryVars(typename std::vector<Scalar>::const_iterator it) override;
253
254 protected:
256
257 // updating the well_state based on well solution dwells
258 void updateWellState(const Simulator& simulator,
259 const BVectorWell& dwells,
260 WellStateType& well_state,
261 DeferredLogger& deferred_logger);
262
263 using WellConnectionProps = typename StdWellEval::StdWellConnections::Properties;
264
265 // Compute connection level PVT properties needed to calulate the
266 // pressure difference between well connections.
269 const WellStateType& well_state) const;
270
272 const WellStateType& well_state,
273 const WellConnectionProps& props,
274 DeferredLogger& deferred_logger);
275
276 void computeWellConnectionPressures(const Simulator& simulator,
277 const WellStateType& well_state,
278 DeferredLogger& deferred_logger);
279
280 template<class Value>
281 void computePerfRate(const IntensiveQuantities& intQuants,
282 const std::vector<Value>& mob,
283 const Value& bhp,
284 const std::vector<Scalar>& Tw,
285 const int perf,
286 const bool allow_cf,
287 std::vector<Value>& cq_s,
288 PerforationRates<Scalar>& perf_rates,
289 DeferredLogger& deferred_logger) const;
290
291 template<class Value>
292 void computePerfRate(const std::vector<Value>& mob,
293 const Value& pressure,
294 const Value& bhp,
295 const Value& rs,
296 const Value& rv,
297 const Value& rvw,
298 const Value& rsw,
299 std::vector<Value>& b_perfcells_dense,
300 const std::vector<Scalar>& Tw,
301 const int perf,
302 const bool allow_cf,
303 const Value& skin_pressure,
304 const std::vector<Value>& cmix_s,
305 std::vector<Value>& cq_s,
306 PerforationRates<Scalar>& perf_rates,
307 DeferredLogger& deferred_logger) const;
308
309 void computeWellRatesWithBhpIterations(const Simulator& ebosSimulator,
310 const Scalar& bhp,
311 std::vector<Scalar>& well_flux,
312 DeferredLogger& deferred_logger) const override;
313
314 std::vector<Scalar>
315 computeWellPotentialWithTHP(const Simulator& ebosSimulator,
316 DeferredLogger& deferred_logger,
317 const WellStateType& well_state) const;
318
319 bool computeWellPotentialsImplicit(const Simulator& ebos_simulator,
320 const WellStateType& well_state,
321 std::vector<Scalar>& well_potentials,
322 DeferredLogger& deferred_logger) const;
323
324 Scalar getRefDensity() const override;
325
326 // get the mobility for specific perforation
327 template<class Value>
328 void getMobility(const Simulator& simulator,
329 const int perf,
330 std::vector<Value>& mob,
331 DeferredLogger& deferred_logger) const;
332
333 void updateWaterMobilityWithPolymer(const Simulator& simulator,
334 const int perf,
335 std::vector<EvalWell>& mob_water,
336 DeferredLogger& deferred_logger) const;
337
338 void updatePrimaryVariablesNewton(const BVectorWell& dwells,
339 const bool stop_or_zero_rate_target,
340 DeferredLogger& deferred_logger);
341
343 const SummaryState& summary_state,
344 DeferredLogger& deferred_logger) const;
345
346 void assembleWellEqWithoutIteration(const Simulator& simulator,
347 const double dt,
348 const Well::InjectionControls& inj_controls,
349 const Well::ProductionControls& prod_controls,
350 WellStateType& well_state,
351 const GroupState<Scalar>& group_state,
352 DeferredLogger& deferred_logger) override;
353
354 void assembleWellEqWithoutIterationImpl(const Simulator& simulator,
355 const double dt,
356 const Well::InjectionControls& inj_controls,
357 const Well::ProductionControls& prod_controls,
358 WellStateType& well_state,
359 const GroupState<Scalar>& group_state,
360 DeferredLogger& deferred_logger);
361
362 void calculateSinglePerf(const Simulator& simulator,
363 const int perf,
364 WellStateType& well_state,
365 std::vector<RateVector>& connectionRates,
366 std::vector<EvalWell>& cq_s,
367 EvalWell& water_flux_s,
368 EvalWell& cq_s_zfrac_effective,
369 DeferredLogger& deferred_logger) const;
370
371 // check whether the well is operable under BHP limit with current reservoir condition
372 void checkOperabilityUnderBHPLimit(const WellStateType& well_state,
373 const Simulator& simulator,
374 DeferredLogger& deferred_logger) override;
375
376 // check whether the well is operable under THP limit with current reservoir condition
377 void checkOperabilityUnderTHPLimit(const Simulator& simulator,
378 const WellStateType& well_state,
379 DeferredLogger& deferred_logger) override;
380
381 // updating the inflow based on the current reservoir condition
382 void updateIPR(const Simulator& simulator,
383 DeferredLogger& deferred_logger) const override;
384
385 // for a well, when all drawdown are in the wrong direction, then this well will not
386 // be able to produce/inject .
387 bool allDrawDownWrongDirection(const Simulator& simulator) const;
388
389 // turn on crossflow to avoid singular well equations
390 // when the well is banned from cross-flow and the BHP is not properly initialized,
391 // we turn on crossflow to avoid singular well equations. It can result in wrong-signed
392 // well rates, it can cause problem for THP calculation
393 // TODO: looking for better alternative to avoid wrong-signed well rates
394 bool openCrossFlowAvoidSingularity(const Simulator& simulator) const;
395
396 // calculate the skin pressure based on water velocity, throughput and polymer concentration.
397 // throughput is used to describe the formation damage during water/polymer injection.
398 // calculated skin pressure will be applied to the drawdown during perforation rate calculation
399 // to handle the effect from formation damage.
400 EvalWell pskin(const Scalar throughput,
401 const EvalWell& water_velocity,
402 const EvalWell& poly_inj_conc,
403 DeferredLogger& deferred_logger) const;
404
405 // calculate the skin pressure based on water velocity, throughput during water injection.
406 EvalWell pskinwater(const Scalar throughput,
407 const EvalWell& water_velocity,
408 DeferredLogger& deferred_logger) const;
409
410 // calculate the injecting polymer molecular weight based on the througput and water velocity
411 EvalWell wpolymermw(const Scalar throughput,
412 const EvalWell& water_velocity,
413 DeferredLogger& deferred_logger) const;
414
415 // modify the water rate for polymer injectivity study
416 void handleInjectivityRate(const Simulator& simulator,
417 const int perf,
418 std::vector<EvalWell>& cq_s) const;
419
420 // handle the extra equations for polymer injectivity study
421 void handleInjectivityEquations(const Simulator& simulator,
422 const WellStateType& well_state,
423 const int perf,
424 const EvalWell& water_flux_s,
425 DeferredLogger& deferred_logger);
426
427 void updateWaterThroughput(const double dt,
428 WellStateType& well_state) const override;
429
430 // checking convergence of extra equations, if there are any
431 void checkConvergenceExtraEqs(const std::vector<Scalar>& res,
432 ConvergenceReport& report) const;
433
434 // updating the connectionRates_ related polymer molecular weight
435 void updateConnectionRatePolyMW(const EvalWell& cq_s_poly,
436 const IntensiveQuantities& int_quants,
437 const WellStateType& well_state,
438 const int perf,
439 std::vector<RateVector>& connectionRates,
440 DeferredLogger& deferred_logger) const;
441
442 std::optional<Scalar>
443 computeBhpAtThpLimitProd(const WellStateType& well_state,
444 const Simulator& simulator,
445 const SummaryState& summary_state,
446 DeferredLogger& deferred_logger) const;
447
448 std::optional<Scalar>
449 computeBhpAtThpLimitInj(const Simulator& simulator,
450 const SummaryState& summary_state,
451 DeferredLogger& deferred_logger) const;
452
453 private:
454 Eval connectionRateEnergy(const std::vector<EvalWell>& cq_s,
455 const IntensiveQuantities& intQuants,
456 DeferredLogger& deferred_logger) const;
457 };
458
459}
460
461#include "StandardWell_impl.hpp"
462
463#endif // OPM_STANDARDWELL_HEADER_INCLUDED
Contains the classes required to extend the black-oil model by brine.
Contains the classes required to extend the black-oil model by solvent component. For details,...
Contains the classes required to extend the black-oil model to include the effects of foam.
Contains the classes required to extend the black-oil model by MICP.
Contains the classes required to extend the black-oil model by polymer.
Contains the classes required to extend the black-oil model by solvents.
Contains the high level supplements required to extend the black oil model to include the effects of ...
Definition: blackoilfoammodules.hh:58
Contains the high level supplements required to extend the black oil model by polymer.
Definition: blackoilpolymermodules.hh:64
Definition: ConvergenceReport.hpp:38
Definition: DeferredLogger.hpp:57
Definition: GroupState.hpp:41
Class encapsulating some information about parallel wells.
Definition: ParallelWellInfo.hpp:198
Definition: StandardWellEval.hpp:46
typename StandardWellEquations< Scalar, IndexTraits, Indices::numEq >::BVectorWell BVectorWell
Definition: StandardWellEval.hpp:65
DenseAd::Evaluation< Scalar, Indices::numEq > Eval
Definition: StandardWellEval.hpp:64
Definition: StandardWell.hpp:60
void calculateExplicitQuantities(const Simulator &simulator, const WellStateType &well_state, DeferredLogger &deferred_logger) override
Definition: StandardWell_impl.hpp:1380
EvalWell wpolymermw(const Scalar throughput, const EvalWell &water_velocity, DeferredLogger &deferred_logger) const
Definition: StandardWell_impl.hpp:1994
typename StdWellEval::EvalWell EvalWell
Definition: StandardWell.hpp:119
void updateWellStateFromPrimaryVariables(WellStateType &well_state, const SummaryState &summary_state, DeferredLogger &deferred_logger) const
Definition: StandardWell_impl.hpp:760
static constexpr int numStaticWellEq
Definition: StandardWell.hpp:100
bool regularize_
Definition: StandardWell.hpp:255
WellConnectionProps computePropertiesForWellConnectionPressures(const Simulator &simulator, const WellStateType &well_state) const
Definition: StandardWell_impl.hpp:1119
typename StdWellEval::BVectorWell BVectorWell
Definition: StandardWell.hpp:120
std::optional< Scalar > computeBhpAtThpLimitProd(const WellStateType &well_state, const Simulator &simulator, const SummaryState &summary_state, DeferredLogger &deferred_logger) const
Definition: StandardWell_impl.hpp:2186
void computeWellRatesWithThpAlqProd(const Simulator &ebos_simulator, const SummaryState &summary_state, DeferredLogger &deferred_logger, std::vector< Scalar > &potentials, Scalar alq) const
Definition: StandardWell_impl.hpp:1699
void addWellContributions(SparseMatrixAdapter &mat) const override
Definition: StandardWell_impl.hpp:1889
std::vector< Scalar > getPrimaryVars() const override
Definition: StandardWell_impl.hpp:2518
void addWellPressureEquations(PressureMatrix &mat, const BVector &x, const int pressureVarIndex, const bool use_well_weights, const WellStateType &well_state) const override
Definition: StandardWell_impl.hpp:1897
void updateWaterMobilityWithPolymer(const Simulator &simulator, const int perf, std::vector< EvalWell > &mob_water, DeferredLogger &deferred_logger) const
Definition: StandardWell_impl.hpp:1825
void assembleWellEqWithoutIteration(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) override
Definition: StandardWell_impl.hpp:341
std::vector< Scalar > computeCurrentWellRates(const Simulator &ebosSimulator, DeferredLogger &deferred_logger) const override
Definition: StandardWell_impl.hpp:2482
void calculateSinglePerf(const Simulator &simulator, const int perf, WellStateType &well_state, std::vector< RateVector > &connectionRates, std::vector< EvalWell > &cq_s, EvalWell &water_flux_s, EvalWell &cq_s_zfrac_effective, DeferredLogger &deferred_logger) const
Definition: StandardWell_impl.hpp:493
void updatePrimaryVariablesNewton(const BVectorWell &dwells, const bool stop_or_zero_rate_target, DeferredLogger &deferred_logger)
Definition: StandardWell_impl.hpp:737
void computeWellConnectionPressures(const Simulator &simulator, const WellStateType &well_state, DeferredLogger &deferred_logger)
Definition: StandardWell_impl.hpp:1340
StandardWell(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)
Definition: StandardWell_impl.hpp:52
typename StdWellEval::StdWellConnections::Properties WellConnectionProps
Definition: StandardWell.hpp:263
void updateConnectionRatePolyMW(const EvalWell &cq_s_poly, const IntensiveQuantities &int_quants, const WellStateType &well_state, const int perf, std::vector< RateVector > &connectionRates, DeferredLogger &deferred_logger) const
Definition: StandardWell_impl.hpp:2144
bool iterateWellEqWithSwitching(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, const bool fixed_control=false, const bool fixed_status=false) override
Definition: StandardWell_impl.hpp:2362
void computeWellRatesWithBhp(const Simulator &ebosSimulator, const Scalar &bhp, std::vector< Scalar > &well_flux, DeferredLogger &deferred_logger) const override
Definition: StandardWell_impl.hpp:1446
void getMobility(const Simulator &simulator, const int perf, std::vector< Value > &mob, DeferredLogger &deferred_logger) const
Definition: StandardWell_impl.hpp:659
void computeWellRatesWithBhpIterations(const Simulator &ebosSimulator, const Scalar &bhp, std::vector< Scalar > &well_flux, DeferredLogger &deferred_logger) const override
Definition: StandardWell_impl.hpp:1492
void updateIPR(const Simulator &simulator, DeferredLogger &deferred_logger) const override
Definition: StandardWell_impl.hpp:779
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 override
Definition: StandardWell_impl.hpp:2201
void updateIPRImplicit(const Simulator &simulator, WellStateType &well_state, DeferredLogger &deferred_logger) override
Definition: StandardWell_impl.hpp:872
void computeWellConnectionDensitesPressures(const Simulator &simulator, const WellStateType &well_state, const WellConnectionProps &props, DeferredLogger &deferred_logger)
Definition: StandardWell_impl.hpp:1289
void computePerfRate(const IntensiveQuantities &intQuants, const std::vector< Value > &mob, const Value &bhp, const std::vector< Scalar > &Tw, const int perf, const bool allow_cf, std::vector< Value > &cq_s, PerforationRates< Scalar > &perf_rates, DeferredLogger &deferred_logger) const
Definition: StandardWell_impl.hpp:93
void updateWellState(const Simulator &simulator, const BVectorWell &dwells, WellStateType &well_state, DeferredLogger &deferred_logger)
Definition: StandardWell_impl.hpp:715
void handleInjectivityEquations(const Simulator &simulator, const WellStateType &well_state, const int perf, const EvalWell &water_flux_s, DeferredLogger &deferred_logger)
Definition: StandardWell_impl.hpp:2077
virtual void apply(const BVector &x, BVector &Ax) const override
Ax = Ax - C D^-1 B x.
Definition: StandardWell_impl.hpp:1394
virtual ConvergenceReport getWellConvergence(const Simulator &simulator, const WellStateType &well_state, const std::vector< Scalar > &B_avg, DeferredLogger &deferred_logger, const bool relax_tolerance) const override
check whether the well equations get converged for this well
Definition: StandardWell_impl.hpp:1170
void checkConvergenceExtraEqs(const std::vector< Scalar > &res, ConvergenceReport &report) const
Definition: StandardWell_impl.hpp:2125
typename StdWellEval::Eval Eval
Definition: StandardWell.hpp:118
Scalar computeWellRatesAndBhpWithThpAlqProd(const Simulator &ebos_simulator, const SummaryState &summary_state, DeferredLogger &deferred_logger, std::vector< Scalar > &potentials, Scalar alq) const
Definition: StandardWell_impl.hpp:1670
static constexpr int Bhp
Definition: StandardWell.hpp:105
bool openCrossFlowAvoidSingularity(const Simulator &simulator) const
Definition: StandardWell_impl.hpp:1108
std::optional< Scalar > computeBhpAtThpLimitInj(const Simulator &simulator, const SummaryState &summary_state, DeferredLogger &deferred_logger) const
Definition: StandardWell_impl.hpp:2283
typename FluidSystem::IndexTraitsType IndexTraits
Definition: StandardWell.hpp:122
bool allDrawDownWrongDirection(const Simulator &simulator) const
Definition: StandardWell_impl.hpp:1066
EvalWell pskin(const Scalar throughput, const EvalWell &water_velocity, const EvalWell &poly_inj_conc, DeferredLogger &deferred_logger) const
Definition: StandardWell_impl.hpp:1949
static constexpr int numWellConservationEq
Definition: StandardWell.hpp:95
int setPrimaryVars(typename std::vector< Scalar >::const_iterator it) override
Definition: StandardWell_impl.hpp:2535
void checkOperabilityUnderTHPLimit(const Simulator &simulator, const WellStateType &well_state, DeferredLogger &deferred_logger) override
Definition: StandardWell_impl.hpp:1017
void computeWellPotentials(const Simulator &simulator, const WellStateType &well_state, std::vector< Scalar > &well_potentials, DeferredLogger &deferred_logger) override
computing the well potentials for group control
Definition: StandardWell_impl.hpp:1716
bool computeWellPotentialsImplicit(const Simulator &ebos_simulator, const WellStateType &well_state, std::vector< Scalar > &well_potentials, DeferredLogger &deferred_logger) const
Definition: StandardWell_impl.hpp:1591
void updateWaterThroughput(const double dt, WellStateType &well_state) const override
Definition: StandardWell_impl.hpp:2023
void checkOperabilityUnderBHPLimit(const WellStateType &well_state, const Simulator &simulator, DeferredLogger &deferred_logger) override
Definition: StandardWell_impl.hpp:948
void recoverWellSolutionAndUpdateWellState(const Simulator &simulator, const BVector &x, WellStateType &well_state, DeferredLogger &deferred_logger) override
Definition: StandardWell_impl.hpp:1426
void assembleWellEqWithoutIterationImpl(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: StandardWell_impl.hpp:367
bool iterateWellEqWithControl(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) override
Definition: StandardWell_impl.hpp:2315
EvalWell pskinwater(const Scalar throughput, const EvalWell &water_velocity, DeferredLogger &deferred_logger) const
Definition: StandardWell_impl.hpp:1917
void solveEqAndUpdateWellState(const Simulator &simulator, WellStateType &well_state, DeferredLogger &deferred_logger) override
Definition: StandardWell_impl.hpp:1358
void handleInjectivityRate(const Simulator &simulator, const int perf, std::vector< EvalWell > &cq_s) const
Definition: StandardWell_impl.hpp:2054
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) override
Definition: StandardWell_impl.hpp:76
std::vector< Scalar > computeWellPotentialWithTHP(const Simulator &ebosSimulator, DeferredLogger &deferred_logger, const WellStateType &well_state) const
Definition: StandardWell_impl.hpp:1556
void updateProductivityIndex(const Simulator &simulator, const WellProdIndexCalculator< Scalar > &wellPICalc, WellStateType &well_state, DeferredLogger &deferred_logger) const override
Definition: StandardWell_impl.hpp:1214
void updatePrimaryVariables(const Simulator &simulator, const WellStateType &well_state, DeferredLogger &deferred_logger) override
Definition: StandardWell_impl.hpp:1791
static constexpr int numWellControlEq
Definition: StandardWell.hpp:97
Scalar getRefDensity() const override
Definition: StandardWell_impl.hpp:1814
Scalar connectionDensity(const int globalConnIdx, const int openConnIdx) const override
Definition: StandardWell_impl.hpp:1776
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 VFPProperties< FluidSystem::Scalar, FluidSystem::IndexTraitsType > * vfp_properties_
Definition: WellInterfaceGeneric.hpp:395
Definition: WellInterface.hpp:76
static constexpr bool has_brine
Definition: WellInterface.hpp:119
GetPropType< TypeTag, Properties::Simulator > Simulator
Definition: WellInterface.hpp:81
typename WellInterfaceFluidSystem< FluidSystem >::RateConverterType RateConverterType
Definition: WellInterface.hpp:103
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
Dune::BCRSMatrix< Opm::MatrixBlock< Scalar, 1, 1 > > PressureMatrix
Definition: WellInterface.hpp:97
GetPropType< TypeTag, Properties::IntensiveQuantities > IntensiveQuantities
Definition: WellInterface.hpp:86
GetPropType< TypeTag, Properties::Scalar > Scalar
Definition: WellInterface.hpp:82
GetPropType< TypeTag, Properties::MaterialLaw > MaterialLaw
Definition: WellInterface.hpp:87
Dune::BlockVector< VectorBlockType > BVector
Definition: WellInterface.hpp:96
static constexpr bool has_polymermw
Definition: WellInterface.hpp:117
static constexpr bool has_polymer
Definition: WellInterface.hpp:113
typename Base::ModelParameters ModelParameters
Definition: WellInterface.hpp:109
GetPropType< TypeTag, Properties::FluidSystem > FluidSystem
Definition: WellInterface.hpp:83
static constexpr bool has_solvent
Definition: WellInterface.hpp:111
static constexpr bool has_foam
Definition: WellInterface.hpp:118
GetPropType< TypeTag, Properties::RateVector > RateVector
Definition: WellInterface.hpp:89
static constexpr bool has_micp
Definition: WellInterface.hpp:123
static constexpr bool has_energy
Definition: WellInterface.hpp:114
GetPropType< TypeTag, Properties::Indices > Indices
Definition: WellInterface.hpp:85
GetPropType< TypeTag, Properties::SparseMatrixAdapter > SparseMatrixAdapter
Definition: WellInterface.hpp:88
static constexpr bool has_zFraction
Definition: WellInterface.hpp:112
Definition: WellProdIndexCalculator.hpp:37
Definition: WellState.hpp:66
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: PerforationData.hpp:41