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 StandardWell(const Well& well,
123 const ParallelWellInfo<Scalar>& pw_info,
124 const int time_step,
125 const ModelParameters& param,
126 const RateConverterType& rate_converter,
127 const int pvtRegionIdx,
128 const int num_components,
129 const int num_phases,
130 const int index_of_well,
131 const std::vector<PerforationData<Scalar>>& perf_data);
132
133 virtual void init(const PhaseUsage* phase_usage_arg,
134 const std::vector<Scalar>& depth_arg,
135 const Scalar gravity_arg,
136 const std::vector<Scalar>& B_avg,
137 const bool changed_to_open_this_step) override;
138
140 virtual ConvergenceReport getWellConvergence(const Simulator& simulator,
141 const WellState<Scalar>& well_state,
142 const std::vector<Scalar>& B_avg,
143 DeferredLogger& deferred_logger,
144 const bool relax_tolerance) const override;
145
147 virtual void apply(const BVector& x, BVector& Ax) const override;
149 virtual void apply(BVector& r) const override;
150
154 const BVector& x,
155 WellState<Scalar>& well_state,
156 DeferredLogger& deferred_logger) override;
157
159 void computeWellPotentials(const Simulator& simulator,
160 const WellState<Scalar>& well_state,
161 std::vector<Scalar>& well_potentials,
162 DeferredLogger& deferred_logger) /* const */ override;
163
164 void updatePrimaryVariables(const Simulator& simulator,
165 const WellState<Scalar>& well_state,
166 DeferredLogger& deferred_logger) override;
167
168 void solveEqAndUpdateWellState(const Simulator& simulator,
169 WellState<Scalar>& well_state,
170 DeferredLogger& deferred_logger) override;
171
172 void calculateExplicitQuantities(const Simulator& simulator,
173 const WellState<Scalar>& well_state,
174 DeferredLogger& deferred_logger) override; // should be const?
175
176 void updateProductivityIndex(const Simulator& simulator,
177 const WellProdIndexCalculator<Scalar>& wellPICalc,
178 WellState<Scalar>& well_state,
179 DeferredLogger& deferred_logger) const override;
180
181 Scalar connectionDensity(const int globalConnIdx,
182 const int openConnIdx) const override;
183
184 void addWellContributions(SparseMatrixAdapter& mat) const override;
185
187 const BVector& x,
188 const int pressureVarIndex,
189 const bool use_well_weights,
190 const WellState<Scalar>& well_state) const override;
191
192 // iterate well equations with the specified control until converged
193 bool iterateWellEqWithControl(const Simulator& simulator,
194 const double dt,
195 const Well::InjectionControls& inj_controls,
196 const Well::ProductionControls& prod_controls,
197 WellState<Scalar>& well_state,
198 const GroupState<Scalar>& group_state,
199 DeferredLogger& deferred_logger) override;
200
201 // iterate well equations including control switching
202 bool iterateWellEqWithSwitching(const Simulator& simulator,
203 const double dt,
204 const Well::InjectionControls& inj_controls,
205 const Well::ProductionControls& prod_controls,
206 WellState<Scalar>& well_state,
207 const GroupState<Scalar>& group_state,
208 DeferredLogger& deferred_logger,
209 const bool fixed_control = false,
210 const bool fixed_status = false) override;
211
212 /* returns BHP */
214 const SummaryState &summary_state,
215 DeferredLogger& deferred_logger,
216 std::vector<Scalar>& potentials,
217 Scalar alq) const;
218
219 void computeWellRatesWithThpAlqProd(const Simulator& ebos_simulator,
220 const SummaryState& summary_state,
221 DeferredLogger& deferred_logger,
222 std::vector<Scalar>& potentials,
223 Scalar alq) const;
224
225 std::optional<Scalar>
226 computeBhpAtThpLimitProdWithAlq(const Simulator& ebos_simulator,
227 const SummaryState& summary_state,
228 const Scalar alq_value,
229 DeferredLogger& deferred_logger,
230 bool iterate_if_no_solution) const override;
231
232 void updateIPRImplicit(const Simulator& simulator,
233 WellState<Scalar>& well_state,
234 DeferredLogger& deferred_logger) override;
235
236 void computeWellRatesWithBhp(const Simulator& ebosSimulator,
237 const Scalar& bhp,
238 std::vector<Scalar>& well_flux,
239 DeferredLogger& deferred_logger) const override;
240
241 // NOTE: These cannot be protected since they are used by GasLiftRuntime
242 using Base::phaseUsage;
244
245 std::vector<Scalar>
246 computeCurrentWellRates(const Simulator& ebosSimulator,
247 DeferredLogger& deferred_logger) const override;
248
249 std::vector<Scalar> getPrimaryVars() const override;
250
251 int setPrimaryVars(typename std::vector<Scalar>::const_iterator it) override;
252
253 protected:
255
256 // updating the well_state based on well solution dwells
257 void updateWellState(const Simulator& simulator,
258 const BVectorWell& dwells,
259 WellState<Scalar>& well_state,
260 DeferredLogger& deferred_logger);
261
262 using WellConnectionProps = typename StdWellEval::StdWellConnections::Properties;
263
264 // Compute connection level PVT properties needed to calulate the
265 // pressure difference between well connections.
268 const WellState<Scalar>& well_state) const;
269
271 const WellState<Scalar>& well_state,
272 const WellConnectionProps& props,
273 DeferredLogger& deferred_logger);
274
275 void computeWellConnectionPressures(const Simulator& simulator,
276 const WellState<Scalar>& well_state,
277 DeferredLogger& deferred_logger);
278
279 template<class Value>
280 void computePerfRate(const IntensiveQuantities& intQuants,
281 const std::vector<Value>& mob,
282 const Value& bhp,
283 const std::vector<Scalar>& Tw,
284 const int perf,
285 const bool allow_cf,
286 std::vector<Value>& cq_s,
287 PerforationRates<Scalar>& perf_rates,
288 DeferredLogger& deferred_logger) const;
289
290 template<class Value>
291 void computePerfRate(const std::vector<Value>& mob,
292 const Value& pressure,
293 const Value& bhp,
294 const Value& rs,
295 const Value& rv,
296 const Value& rvw,
297 const Value& rsw,
298 std::vector<Value>& b_perfcells_dense,
299 const std::vector<Scalar>& Tw,
300 const int perf,
301 const bool allow_cf,
302 const Value& skin_pressure,
303 const std::vector<Value>& cmix_s,
304 std::vector<Value>& cq_s,
305 PerforationRates<Scalar>& perf_rates,
306 DeferredLogger& deferred_logger) const;
307
308 void computeWellRatesWithBhpIterations(const Simulator& ebosSimulator,
309 const Scalar& bhp,
310 std::vector<Scalar>& well_flux,
311 DeferredLogger& deferred_logger) const override;
312
313 std::vector<Scalar>
314 computeWellPotentialWithTHP(const Simulator& ebosSimulator,
315 DeferredLogger& deferred_logger,
316 const WellState<Scalar>& well_state) const;
317
318 bool computeWellPotentialsImplicit(const Simulator& ebos_simulator,
319 const WellState<Scalar>& well_state,
320 std::vector<Scalar>& well_potentials,
321 DeferredLogger& deferred_logger) const;
322
323 Scalar getRefDensity() const override;
324
325 // get the mobility for specific perforation
326 template<class Value>
327 void getMobility(const Simulator& simulator,
328 const int perf,
329 std::vector<Value>& mob,
330 DeferredLogger& deferred_logger) const;
331
332 void updateWaterMobilityWithPolymer(const Simulator& simulator,
333 const int perf,
334 std::vector<EvalWell>& mob_water,
335 DeferredLogger& deferred_logger) const;
336
337 void updatePrimaryVariablesNewton(const BVectorWell& dwells,
338 const bool stop_or_zero_rate_target,
339 DeferredLogger& deferred_logger);
340
342 const SummaryState& summary_state,
343 DeferredLogger& deferred_logger) const;
344
345 void assembleWellEqWithoutIteration(const Simulator& simulator,
346 const double dt,
347 const Well::InjectionControls& inj_controls,
348 const Well::ProductionControls& prod_controls,
349 WellState<Scalar>& well_state,
350 const GroupState<Scalar>& group_state,
351 DeferredLogger& deferred_logger) override;
352
353 void assembleWellEqWithoutIterationImpl(const Simulator& simulator,
354 const double dt,
355 const Well::InjectionControls& inj_controls,
356 const Well::ProductionControls& prod_controls,
357 WellState<Scalar>& well_state,
358 const GroupState<Scalar>& group_state,
359 DeferredLogger& deferred_logger);
360
361 void calculateSinglePerf(const Simulator& simulator,
362 const int perf,
363 WellState<Scalar>& well_state,
364 std::vector<RateVector>& connectionRates,
365 std::vector<EvalWell>& cq_s,
366 EvalWell& water_flux_s,
367 EvalWell& cq_s_zfrac_effective,
368 DeferredLogger& deferred_logger) const;
369
370 // check whether the well is operable under BHP limit with current reservoir condition
372 const Simulator& simulator,
373 DeferredLogger& deferred_logger) override;
374
375 // check whether the well is operable under THP limit with current reservoir condition
376 void checkOperabilityUnderTHPLimit(const Simulator& simulator,
377 const WellState<Scalar>& well_state,
378 DeferredLogger& deferred_logger) override;
379
380 // updating the inflow based on the current reservoir condition
381 void updateIPR(const Simulator& simulator,
382 DeferredLogger& deferred_logger) const override;
383
384 // for a well, when all drawdown are in the wrong direction, then this well will not
385 // be able to produce/inject .
386 bool allDrawDownWrongDirection(const Simulator& simulator) const;
387
388 // turn on crossflow to avoid singular well equations
389 // when the well is banned from cross-flow and the BHP is not properly initialized,
390 // we turn on crossflow to avoid singular well equations. It can result in wrong-signed
391 // well rates, it can cause problem for THP calculation
392 // TODO: looking for better alternative to avoid wrong-signed well rates
393 bool openCrossFlowAvoidSingularity(const Simulator& simulator) const;
394
395 // calculate the skin pressure based on water velocity, throughput and polymer concentration.
396 // throughput is used to describe the formation damage during water/polymer injection.
397 // calculated skin pressure will be applied to the drawdown during perforation rate calculation
398 // to handle the effect from formation damage.
399 EvalWell pskin(const Scalar throughput,
400 const EvalWell& water_velocity,
401 const EvalWell& poly_inj_conc,
402 DeferredLogger& deferred_logger) const;
403
404 // calculate the skin pressure based on water velocity, throughput during water injection.
405 EvalWell pskinwater(const Scalar throughput,
406 const EvalWell& water_velocity,
407 DeferredLogger& deferred_logger) const;
408
409 // calculate the injecting polymer molecular weight based on the througput and water velocity
410 EvalWell wpolymermw(const Scalar throughput,
411 const EvalWell& water_velocity,
412 DeferredLogger& deferred_logger) const;
413
414 // modify the water rate for polymer injectivity study
415 void handleInjectivityRate(const Simulator& simulator,
416 const int perf,
417 std::vector<EvalWell>& cq_s) const;
418
419 // handle the extra equations for polymer injectivity study
420 void handleInjectivityEquations(const Simulator& simulator,
421 const WellState<Scalar>& well_state,
422 const int perf,
423 const EvalWell& water_flux_s,
424 DeferredLogger& deferred_logger);
425
426 void updateWaterThroughput(const double dt,
427 WellState<Scalar>& well_state) const override;
428
429 // checking convergence of extra equations, if there are any
430 void checkConvergenceExtraEqs(const std::vector<Scalar>& res,
431 ConvergenceReport& report) const;
432
433 // updating the connectionRates_ related polymer molecular weight
434 void updateConnectionRatePolyMW(const EvalWell& cq_s_poly,
435 const IntensiveQuantities& int_quants,
436 const WellState<Scalar>& well_state,
437 const int perf,
438 std::vector<RateVector>& connectionRates,
439 DeferredLogger& deferred_logger) const;
440
441 std::optional<Scalar>
443 const Simulator& simulator,
444 const SummaryState& summary_state,
445 DeferredLogger& deferred_logger) const;
446
447 std::optional<Scalar>
448 computeBhpAtThpLimitInj(const Simulator& simulator,
449 const SummaryState& summary_state,
450 DeferredLogger& deferred_logger) const;
451
452 private:
453 Eval connectionRateEnergy(const Scalar maxOilSaturation,
454 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:43
Class encapsulating some information about parallel wells.
Definition: ParallelWellInfo.hpp:195
Definition: StandardWellEval.hpp:46
typename StandardWellEquations< Scalar, Indices::numEq >::BVectorWell BVectorWell
Definition: StandardWellEval.hpp:64
DenseAd::Evaluation< Scalar, Indices::numEq > Eval
Definition: StandardWellEval.hpp:63
Definition: StandardWell.hpp:60
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) override
Definition: StandardWell_impl.hpp:76
void computeWellConnectionPressures(const Simulator &simulator, const WellState< Scalar > &well_state, DeferredLogger &deferred_logger)
Definition: StandardWell_impl.hpp:1343
EvalWell wpolymermw(const Scalar throughput, const EvalWell &water_velocity, DeferredLogger &deferred_logger) const
Definition: StandardWell_impl.hpp:1997
void computeWellPotentials(const Simulator &simulator, const WellState< Scalar > &well_state, std::vector< Scalar > &well_potentials, DeferredLogger &deferred_logger) override
computing the well potentials for group control
Definition: StandardWell_impl.hpp:1719
typename StdWellEval::EvalWell EvalWell
Definition: StandardWell.hpp:119
static constexpr int numStaticWellEq
Definition: StandardWell.hpp:100
bool regularize_
Definition: StandardWell.hpp:254
typename StdWellEval::BVectorWell BVectorWell
Definition: StandardWell.hpp:120
void checkOperabilityUnderBHPLimit(const WellState< Scalar > &well_state, const Simulator &simulator, DeferredLogger &deferred_logger) override
Definition: StandardWell_impl.hpp:951
bool iterateWellEqWithControl(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) override
Definition: StandardWell_impl.hpp:2319
void computeWellRatesWithThpAlqProd(const Simulator &ebos_simulator, const SummaryState &summary_state, DeferredLogger &deferred_logger, std::vector< Scalar > &potentials, Scalar alq) const
Definition: StandardWell_impl.hpp:1702
void addWellContributions(SparseMatrixAdapter &mat) const override
Definition: StandardWell_impl.hpp:1892
std::vector< Scalar > getPrimaryVars() const override
Definition: StandardWell_impl.hpp:2518
void updateIPRImplicit(const Simulator &simulator, WellState< Scalar > &well_state, DeferredLogger &deferred_logger) override
Definition: StandardWell_impl.hpp:875
void updateProductivityIndex(const Simulator &simulator, const WellProdIndexCalculator< Scalar > &wellPICalc, WellState< Scalar > &well_state, DeferredLogger &deferred_logger) const override
Definition: StandardWell_impl.hpp:1217
void updateWaterMobilityWithPolymer(const Simulator &simulator, const int perf, std::vector< EvalWell > &mob_water, DeferredLogger &deferred_logger) const
Definition: StandardWell_impl.hpp:1828
bool iterateWellEqWithSwitching(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, const bool fixed_control=false, const bool fixed_status=false) override
Definition: StandardWell_impl.hpp:2366
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_components, const int num_phases, const int index_of_well, const std::vector< PerforationData< Scalar > > &perf_data)
Definition: StandardWell_impl.hpp:52
void updateWaterThroughput(const double dt, WellState< Scalar > &well_state) const override
Definition: StandardWell_impl.hpp:2026
std::vector< Scalar > computeCurrentWellRates(const Simulator &ebosSimulator, DeferredLogger &deferred_logger) const override
Definition: StandardWell_impl.hpp:2482
void updatePrimaryVariablesNewton(const BVectorWell &dwells, const bool stop_or_zero_rate_target, DeferredLogger &deferred_logger)
Definition: StandardWell_impl.hpp:740
void calculateExplicitQuantities(const Simulator &simulator, const WellState< Scalar > &well_state, DeferredLogger &deferred_logger) override
Definition: StandardWell_impl.hpp:1383
void solveEqAndUpdateWellState(const Simulator &simulator, WellState< Scalar > &well_state, DeferredLogger &deferred_logger) override
Definition: StandardWell_impl.hpp:1361
void recoverWellSolutionAndUpdateWellState(const Simulator &simulator, const BVector &x, WellState< Scalar > &well_state, DeferredLogger &deferred_logger) override
Definition: StandardWell_impl.hpp:1429
void updateWellStateFromPrimaryVariables(WellState< Scalar > &well_state, const SummaryState &summary_state, DeferredLogger &deferred_logger) const
Definition: StandardWell_impl.hpp:763
typename StdWellEval::StdWellConnections::Properties WellConnectionProps
Definition: StandardWell.hpp:262
void computeWellRatesWithBhp(const Simulator &ebosSimulator, const Scalar &bhp, std::vector< Scalar > &well_flux, DeferredLogger &deferred_logger) const override
Definition: StandardWell_impl.hpp:1449
WellConnectionProps computePropertiesForWellConnectionPressures(const Simulator &simulator, const WellState< Scalar > &well_state) const
Definition: StandardWell_impl.hpp:1122
void getMobility(const Simulator &simulator, const int perf, std::vector< Value > &mob, DeferredLogger &deferred_logger) const
Definition: StandardWell_impl.hpp:662
void computeWellRatesWithBhpIterations(const Simulator &ebosSimulator, const Scalar &bhp, std::vector< Scalar > &well_flux, DeferredLogger &deferred_logger) const override
Definition: StandardWell_impl.hpp:1495
void updateIPR(const Simulator &simulator, DeferredLogger &deferred_logger) const override
Definition: StandardWell_impl.hpp:782
std::vector< Scalar > computeWellPotentialWithTHP(const Simulator &ebosSimulator, DeferredLogger &deferred_logger, const WellState< Scalar > &well_state) const
Definition: StandardWell_impl.hpp:1559
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:2205
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:94
void updateWellState(const Simulator &simulator, const BVectorWell &dwells, WellState< Scalar > &well_state, DeferredLogger &deferred_logger)
Definition: StandardWell_impl.hpp:718
void assembleWellEqWithoutIteration(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) override
Definition: StandardWell_impl.hpp:342
void computeWellConnectionDensitesPressures(const Simulator &simulator, const WellState< Scalar > &well_state, const WellConnectionProps &props, DeferredLogger &deferred_logger)
Definition: StandardWell_impl.hpp:1292
virtual void apply(const BVector &x, BVector &Ax) const override
Ax = Ax - C D^-1 B x.
Definition: StandardWell_impl.hpp:1397
std::optional< Scalar > computeBhpAtThpLimitProd(const WellState< Scalar > &well_state, const Simulator &simulator, const SummaryState &summary_state, DeferredLogger &deferred_logger) const
Definition: StandardWell_impl.hpp:2190
void checkConvergenceExtraEqs(const std::vector< Scalar > &res, ConvergenceReport &report) const
Definition: StandardWell_impl.hpp:2128
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:1673
void updateConnectionRatePolyMW(const EvalWell &cq_s_poly, const IntensiveQuantities &int_quants, const WellState< Scalar > &well_state, const int perf, std::vector< RateVector > &connectionRates, DeferredLogger &deferred_logger) const
Definition: StandardWell_impl.hpp:2147
static constexpr int Bhp
Definition: StandardWell.hpp:105
bool openCrossFlowAvoidSingularity(const Simulator &simulator) const
Definition: StandardWell_impl.hpp:1111
void updatePrimaryVariables(const Simulator &simulator, const WellState< Scalar > &well_state, DeferredLogger &deferred_logger) override
Definition: StandardWell_impl.hpp:1794
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 override
check whether the well equations get converged for this well
Definition: StandardWell_impl.hpp:1173
std::optional< Scalar > computeBhpAtThpLimitInj(const Simulator &simulator, const SummaryState &summary_state, DeferredLogger &deferred_logger) const
Definition: StandardWell_impl.hpp:2287
bool allDrawDownWrongDirection(const Simulator &simulator) const
Definition: StandardWell_impl.hpp:1069
EvalWell pskin(const Scalar throughput, const EvalWell &water_velocity, const EvalWell &poly_inj_conc, DeferredLogger &deferred_logger) const
Definition: StandardWell_impl.hpp:1952
static constexpr int numWellConservationEq
Definition: StandardWell.hpp:95
void handleInjectivityEquations(const Simulator &simulator, const WellState< Scalar > &well_state, const int perf, const EvalWell &water_flux_s, DeferredLogger &deferred_logger)
Definition: StandardWell_impl.hpp:2080
int setPrimaryVars(typename std::vector< Scalar >::const_iterator it) override
Definition: StandardWell_impl.hpp:2535
void assembleWellEqWithoutIterationImpl(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: StandardWell_impl.hpp:368
void checkOperabilityUnderTHPLimit(const Simulator &simulator, const WellState< Scalar > &well_state, DeferredLogger &deferred_logger) override
Definition: StandardWell_impl.hpp:1020
void addWellPressureEquations(PressureMatrix &mat, const BVector &x, const int pressureVarIndex, const bool use_well_weights, const WellState< Scalar > &well_state) const override
Definition: StandardWell_impl.hpp:1900
EvalWell pskinwater(const Scalar throughput, const EvalWell &water_velocity, DeferredLogger &deferred_logger) const
Definition: StandardWell_impl.hpp:1920
void handleInjectivityRate(const Simulator &simulator, const int perf, std::vector< EvalWell > &cq_s) const
Definition: StandardWell_impl.hpp:2057
static constexpr int numWellControlEq
Definition: StandardWell.hpp:97
Scalar getRefDensity() const override
Definition: StandardWell_impl.hpp:1817
void calculateSinglePerf(const Simulator &simulator, const int perf, WellState< Scalar > &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:494
bool computeWellPotentialsImplicit(const Simulator &ebos_simulator, const WellState< Scalar > &well_state, std::vector< Scalar > &well_potentials, DeferredLogger &deferred_logger) const
Definition: StandardWell_impl.hpp:1594
Scalar connectionDensity(const int globalConnIdx, const int openConnIdx) const override
Definition: StandardWell_impl.hpp:1779
static constexpr int Oil
Definition: WellInterfaceFluidSystem.hpp:63
static constexpr int Water
Definition: WellInterfaceFluidSystem.hpp:62
static constexpr int Gas
Definition: WellInterfaceFluidSystem.hpp:64
const std::string & name() const
Well name.
int pvtRegionIdx() const
Definition: WellInterfaceGeneric.hpp:124
const VFPProperties< FluidSystem::Scalar > * vfp_properties_
Definition: WellInterfaceGeneric.hpp:392
Definition: WellInterface.hpp:77
static constexpr bool has_brine
Definition: WellInterface.hpp:116
GetPropType< TypeTag, Properties::Simulator > Simulator
Definition: WellInterface.hpp:82
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
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:83
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:114
static constexpr bool has_polymer
Definition: WellInterface.hpp:110
typename Base::ModelParameters ModelParameters
Definition: WellInterface.hpp:106
GetPropType< TypeTag, Properties::FluidSystem > FluidSystem
Definition: WellInterface.hpp:84
static constexpr bool has_solvent
Definition: WellInterface.hpp:108
static constexpr bool has_foam
Definition: WellInterface.hpp:115
GetPropType< TypeTag, Properties::RateVector > RateVector
Definition: WellInterface.hpp:89
static constexpr bool has_micp
Definition: WellInterface.hpp:120
static constexpr bool has_energy
Definition: WellInterface.hpp:111
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:109
Definition: WellProdIndexCalculator.hpp:37
Definition: WellState.hpp:65
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
Definition: BlackoilPhases.hpp:46