23 #ifndef OPM_STANDARDWELL_HEADER_INCLUDED
24 #define OPM_STANDARDWELL_HEADER_INCLUDED
26 #include <opm/simulators/timestepping/ConvergenceReport.hpp>
28 #include <opm/simulators/wells/StandardWellGeneric.hpp>
29 #include <opm/simulators/wells/VFPInjProperties.hpp>
30 #include <opm/simulators/wells/VFPProdProperties.hpp>
31 #include <opm/simulators/wells/WellInterface.hpp>
32 #include <opm/simulators/wells/WellProdIndexCalculator.hpp>
33 #include <opm/simulators/wells/ParallelWellInfo.hpp>
34 #include <opm/simulators/wells/GasLiftSingleWell.hpp>
35 #include <opm/simulators/wells/GasLiftGroupInfo.hpp>
37 #include <opm/models/blackoil/blackoilpolymermodules.hh>
38 #include <opm/models/blackoil/blackoilsolventmodules.hh>
39 #include <opm/models/blackoil/blackoilextbomodules.hh>
40 #include <opm/models/blackoil/blackoilfoammodules.hh>
41 #include <opm/models/blackoil/blackoilbrinemodules.hh>
42 #include <opm/models/blackoil/blackoilmicpmodules.hh>
44 #include <opm/material/densead/DynamicEvaluation.hpp>
45 #include <opm/parser/eclipse/EclipseState/Runspec.hpp>
46 #include <opm/parser/eclipse/EclipseState/Schedule/ScheduleTypes.hpp>
48 #include <opm/simulators/wells/StandardWellEval.hpp>
50 #include <dune/common/dynvector.hh>
51 #include <dune/common/dynmatrix.hh>
55 #include <fmt/format.h>
60 template<
typename TypeTag>
63 GetPropType<TypeTag, Properties::Indices>,
64 GetPropType<TypeTag, Properties::Scalar>>
70 GetPropType<TypeTag, Properties::Indices>,
71 GetPropType<TypeTag, Properties::Scalar>>;
76 using typename Base::Simulator;
77 using typename Base::IntensiveQuantities;
78 using typename Base::FluidSystem;
79 using typename Base::MaterialLaw;
81 using typename Base::Indices;
82 using typename Base::RateConverterType;
83 using typename Base::SparseMatrixAdapter;
84 using typename Base::FluidState;
85 using typename Base::RateVector;
87 using typename Base::GLiftOptWells;
88 using typename Base::GLiftProdWells;
89 using typename Base::GLiftWellStateMap;
90 using typename Base::GLiftSyncGroups;
92 using Base::has_solvent;
93 using Base::has_zFraction;
94 using Base::has_polymer;
95 using Base::has_polymermw;
97 using Base::has_brine;
98 using Base::has_energy;
101 using PolymerModule = BlackOilPolymerModule<TypeTag>;
102 using FoamModule = BlackOilFoamModule<TypeTag>;
103 using BrineModule = BlackOilBrineModule<TypeTag>;
106 static constexpr
int numWellConservationEq = Indices::numPhases + Indices::numSolvents;
108 static constexpr
int numWellControlEq = 1;
111 static constexpr
int numStaticWellEq = numWellConservationEq + numWellControlEq;
116 static constexpr
int Bhp = numStaticWellEq - numWellControlEq;
118 using typename Base::Scalar;
126 using typename Base::BVector;
128 using Eval =
typename StdWellEval::Eval;
129 using EvalWell =
typename StdWellEval::EvalWell;
130 using BVectorWell =
typename StdWellEval::BVectorWell;
136 const RateConverterType& rate_converter,
137 const int pvtRegionIdx,
138 const int num_components,
139 const int num_phases,
140 const int index_of_well,
141 const std::vector<PerforationData>& perf_data);
143 virtual void init(
const PhaseUsage* phase_usage_arg,
144 const std::vector<double>& depth_arg,
145 const double gravity_arg,
147 const std::vector< Scalar >& B_avg)
override;
150 virtual void initPrimaryVariablesEvaluation()
const override;
154 const std::vector<double>& B_avg,
156 const bool relax_tolerance =
false)
const override;
159 virtual void apply(
const BVector& x, BVector& Ax)
const override;
161 virtual void apply(BVector& r)
const override;
172 std::vector<double>& well_potentials,
175 virtual void updatePrimaryVariables(
const WellState& well_state,
DeferredLogger& deferred_logger)
const override;
179 virtual void calculateExplicitQuantities(
const Simulator& ebosSimulator,
183 virtual void updateProductivityIndex(
const Simulator& ebosSimulator,
188 virtual void addWellContributions(SparseMatrixAdapter& mat)
const override;
191 bool iterateWellEqWithControl(
const Simulator& ebosSimulator,
193 const Well::InjectionControls& inj_controls,
194 const Well::ProductionControls& prod_controls,
205 virtual void gasLiftOptimizationStage1 (
208 const Simulator& ebosSimulator,
210 GLiftProdWells &prod_wells,
211 GLiftOptWells &glift_wells,
212 GLiftWellStateMap &state_map,
214 GLiftSyncGroups &sync_groups
218 double computeWellRatesAndBhpWithThpAlqProd(
const Simulator &ebos_simulator,
219 const SummaryState &summary_state,
221 std::vector<double> &potentials,
224 void computeWellRatesWithThpAlqProd(
225 const Simulator &ebos_simulator,
226 const SummaryState &summary_state,
228 std::vector<double> &potentials,
232 std::optional<double> computeBhpAtThpLimitProdWithAlq(
233 const Simulator& ebos_simulator,
234 const SummaryState& summary_state,
236 double alq_value)
const;
239 void computeWellRatesWithBhp(
240 const Simulator& ebosSimulator,
242 std::vector<double>& well_flux,
246 using Base::phaseUsage;
247 using Base::vfp_properties_;
252 void computeConnLevelProdInd(
const FluidState& fs,
253 const std::function<
double(
const double)>& connPICalc,
254 const std::vector<EvalWell>& mobility,
255 double* connPI)
const;
257 void computeConnLevelInjInd(
const typename StandardWell<TypeTag>::FluidState& fs,
258 const Phase preferred_phase,
259 const std::function<
double(
const double)>& connIICalc,
260 const std::vector<EvalWell>& mobility,
266 Eval getPerfCellPressure(
const FluidState& fs)
const;
269 void recoverSolutionWell(
const BVector& x, BVectorWell& xw)
const;
272 void updateWellState(
const BVectorWell& dwells,
278 void computePropertiesForWellConnectionPressures(
const Simulator& ebosSimulator,
280 std::vector<double>& b_perf,
281 std::vector<double>& rsmax_perf,
282 std::vector<double>& rvmax_perf,
283 std::vector<double>& surf_dens_perf)
const;
285 void computeWellConnectionDensitesPressures(
const Simulator& ebosSimulator,
287 const std::vector<double>& b_perf,
288 const std::vector<double>& rsmax_perf,
289 const std::vector<double>& rvmax_perf,
290 const std::vector<double>& surf_dens_perf);
292 void computeWellConnectionPressures(
const Simulator& ebosSimulator,
295 void computePerfRateEval(
const IntensiveQuantities& intQuants,
296 const std::vector<EvalWell>& mob,
301 std::vector<EvalWell>& cq_s,
302 double& perf_dis_gas_rate,
303 double& perf_vap_oil_rate,
306 void computePerfRateScalar(
const IntensiveQuantities& intQuants,
307 const std::vector<Scalar>& mob,
312 std::vector<Scalar>& cq_s,
315 template<
class Value>
316 void computePerfRate(
const std::vector<Value>& mob,
317 const Value& pressure,
321 std::vector<Value>& b_perfcells_dense,
325 const Value& skin_pressure,
326 const std::vector<Value>& cmix_s,
327 std::vector<Value>& cq_s,
328 double& perf_dis_gas_rate,
329 double& perf_vap_oil_rate,
332 void computeWellRatesWithBhpIterations(
const Simulator& ebosSimulator,
334 std::vector<double>& well_flux,
337 std::vector<double> computeWellPotentialWithTHP(
338 const Simulator& ebosSimulator,
343 virtual double getRefDensity()
const override;
346 void getMobilityEval(
const Simulator& ebosSimulator,
348 std::vector<EvalWell>& mob,
352 void getMobilityScalar(
const Simulator& ebosSimulator,
354 std::vector<Scalar>& mob,
358 void updateWaterMobilityWithPolymer(
const Simulator& ebos_simulator,
360 std::vector<EvalWell>& mob_water,
363 void updatePrimaryVariablesNewton(
const BVectorWell& dwells,
368 void updateExtraPrimaryVariables(
const BVectorWell& dwells)
const;
373 virtual void assembleWellEqWithoutIteration(
const Simulator& ebosSimulator,
375 const Well::InjectionControls& inj_controls,
376 const Well::ProductionControls& prod_controls,
381 void assembleWellEqWithoutIterationImpl(
const Simulator& ebosSimulator,
387 void calculateSinglePerf(
const Simulator& ebosSimulator,
390 std::vector<RateVector>& connectionRates,
391 std::vector<EvalWell>& cq_s,
392 EvalWell& water_flux_s,
393 EvalWell& cq_s_zfrac_effective,
397 virtual void checkOperabilityUnderBHPLimitProducer(
const WellState& well_state,
const Simulator& ebos_simulator,
DeferredLogger& deferred_logger)
override;
400 virtual void checkOperabilityUnderTHPLimitProducer(
const Simulator& ebos_simulator,
const WellState& well_state,
DeferredLogger& deferred_logger)
override;
403 virtual void updateIPR(
const Simulator& ebos_simulator,
DeferredLogger& deferred_logger)
const override;
407 bool allDrawDownWrongDirection(
const Simulator& ebos_simulator)
const;
410 bool canProduceInjectWithCurrentBhp(
const Simulator& ebos_simulator,
419 bool openCrossFlowAvoidSingularity(
const Simulator& ebos_simulator)
const;
425 EvalWell pskin(
const double throuhgput,
426 const EvalWell& water_velocity,
427 const EvalWell& poly_inj_conc,
431 EvalWell pskinwater(
const double throughput,
432 const EvalWell& water_velocity,
436 EvalWell wpolymermw(
const double throughput,
437 const EvalWell& water_velocity,
441 void handleInjectivityRate(
const Simulator& ebosSimulator,
443 std::vector<EvalWell>& cq_s)
const;
446 void handleInjectivityEquations(
const Simulator& ebosSimulator,
449 const EvalWell& water_flux_s,
452 virtual void updateWaterThroughput(
const double dt,
WellState& well_state)
const override;
455 void checkConvergenceExtraEqs(
const std::vector<double>& res,
459 void updateConnectionRatePolyMW(
const EvalWell& cq_s_poly,
460 const IntensiveQuantities& int_quants,
463 std::vector<RateVector>& connectionRates,
467 std::optional<double> computeBhpAtThpLimitProd(
const WellState& well_state,
468 const Simulator& ebos_simulator,
469 const SummaryState& summary_state,
472 std::optional<double> computeBhpAtThpLimitInj(
const Simulator& ebos_simulator,
473 const SummaryState& summary_state,
480 #include "StandardWell_impl.hpp"
Facility for converting component rates at surface conditions to phase (voidage) rates at reservoir c...
Represents the convergence status of the whole simulator, to make it possible to query and store the ...
Definition: ConvergenceReport.hpp:36
Definition: DeferredLogger.hpp:57
Definition: GasLiftGroupInfo.hpp:46
Definition: GasLiftSingleWell.hpp:46
Definition: GroupState.hpp:34
Class encapsulating some information about parallel wells.
Definition: ParallelWellInfo.hpp:252
Definition: StandardWellEval.hpp:47
Definition: StandardWell.hpp:65
virtual bool jacobianContainsWellContributions() const override
Wether the Jacobian will also have well contributions in it.
Definition: StandardWell.hpp:200
virtual std::vector< double > computeCurrentWellRates(const Simulator &ebosSimulator, DeferredLogger &deferred_logger) const override
Compute well rates based on current reservoir conditions and well variables.
Definition: StandardWell_impl.hpp:2426
virtual ConvergenceReport getWellConvergence(const WellState &well_state, const std::vector< double > &B_avg, DeferredLogger &deferred_logger, const bool relax_tolerance=false) const override
check whether the well equations get converged for this well
Definition: StandardWell_impl.hpp:1344
virtual void recoverWellSolutionAndUpdateWellState(const BVector &x, WellState &well_state, DeferredLogger &deferred_logger) const override
using the solution x to recover the solution xw for wells and applying xw to update Well State
Definition: StandardWell_impl.hpp:1642
virtual void apply(const BVector &x, BVector &Ax) const override
Ax = Ax - C D^-1 B x.
Definition: StandardWell_impl.hpp:1578
virtual void computeWellPotentials(const Simulator &ebosSimulator, const WellState &well_state, std::vector< double > &well_potentials, DeferredLogger &deferred_logger) override
computing the well potentials for group control
Definition: StandardWell_impl.hpp:1858
const std::string & name() const
Well name.
Definition: WellInterfaceGeneric.cpp:108
Definition: WellInterface.hpp:71
Collect per-connection static information to enable calculating connection-level or well-level produc...
Definition: WellProdIndexCalculator.hpp:36
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition: WellState.hpp:56
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:26
Solver parameters for the BlackoilModel.
Definition: BlackoilModelParametersEbos.hpp:327
bool matrix_add_well_contributions_
Whether to add influences of wells between cells to the matrix and preconditioner matrix.
Definition: BlackoilModelParametersEbos.hpp:414
Definition: BlackoilPhases.hpp:45