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>
35 #include <opm/models/blackoil/blackoilpolymermodules.hh>
36 #include <opm/models/blackoil/blackoilsolventmodules.hh>
37 #include <opm/models/blackoil/blackoilextbomodules.hh>
38 #include <opm/models/blackoil/blackoilfoammodules.hh>
39 #include <opm/models/blackoil/blackoilbrinemodules.hh>
40 #include <opm/models/blackoil/blackoilmicpmodules.hh>
42 #include <opm/material/densead/DynamicEvaluation.hpp>
43 #include <opm/input/eclipse/EclipseState/Runspec.hpp>
44 #include <opm/input/eclipse/Schedule/ScheduleTypes.hpp>
46 #include <opm/simulators/wells/StandardWellEval.hpp>
48 #include <dune/common/dynvector.hh>
49 #include <dune/common/dynmatrix.hh>
53 #include <fmt/format.h>
58 template<
typename TypeTag>
61 GetPropType<TypeTag, Properties::Indices>,
62 GetPropType<TypeTag, Properties::Scalar>>
68 GetPropType<TypeTag, Properties::Indices>,
69 GetPropType<TypeTag, Properties::Scalar>>;
74 using typename Base::Simulator;
75 using typename Base::IntensiveQuantities;
76 using typename Base::FluidSystem;
77 using typename Base::MaterialLaw;
79 using typename Base::Indices;
80 using typename Base::RateConverterType;
81 using typename Base::SparseMatrixAdapter;
82 using typename Base::FluidState;
83 using typename Base::RateVector;
85 using Base::has_solvent;
86 using Base::has_zFraction;
87 using Base::has_polymer;
88 using Base::has_polymermw;
90 using Base::has_brine;
91 using Base::has_energy;
94 using PolymerModule = BlackOilPolymerModule<TypeTag>;
95 using FoamModule = BlackOilFoamModule<TypeTag>;
96 using BrineModule = BlackOilBrineModule<TypeTag>;
97 using typename Base::PressureMatrix;
100 static constexpr
int numWellConservationEq = Indices::numPhases + Indices::numSolvents;
102 static constexpr
int numWellControlEq = 1;
105 static constexpr
int numStaticWellEq = numWellConservationEq + numWellControlEq;
110 static constexpr
int Bhp = numStaticWellEq - numWellControlEq;
112 using typename Base::Scalar;
120 using typename Base::BVector;
122 using Eval =
typename StdWellEval::Eval;
123 using EvalWell =
typename StdWellEval::EvalWell;
124 using BVectorWell =
typename StdWellEval::BVectorWell;
125 using DiagMatrixBlockWellType =
typename StdWellEval::DiagMatrixBlockWellType;
131 const RateConverterType& rate_converter,
132 const int pvtRegionIdx,
133 const int num_components,
134 const int num_phases,
135 const int index_of_well,
136 const std::vector<PerforationData>& perf_data);
138 virtual void init(
const PhaseUsage* phase_usage_arg,
139 const std::vector<double>& depth_arg,
140 const double gravity_arg,
142 const std::vector< Scalar >& B_avg,
143 const bool changed_to_open_this_step)
override;
146 virtual void initPrimaryVariablesEvaluation()
const override;
150 const std::vector<double>& B_avg,
152 const bool relax_tolerance =
false)
const override;
155 virtual void apply(
const BVector& x, BVector& Ax)
const override;
157 virtual void apply(BVector& r)
const override;
168 std::vector<double>& well_potentials,
171 virtual void updatePrimaryVariables(
const WellState& well_state,
DeferredLogger& deferred_logger)
const override;
175 virtual void calculateExplicitQuantities(
const Simulator& ebosSimulator,
179 virtual void updateProductivityIndex(
const Simulator& ebosSimulator,
184 virtual void addWellContributions(SparseMatrixAdapter& mat)
const override;
186 virtual void addWellPressureEquations(PressureMatrix& mat,
188 const int pressureVarIndex,
189 const bool use_well_weights,
190 const WellState& well_state)
const override;
193 bool iterateWellEqWithControl(
const Simulator& ebosSimulator,
195 const Well::InjectionControls& inj_controls,
196 const Well::ProductionControls& prod_controls,
208 double computeWellRatesAndBhpWithThpAlqProd(
const Simulator &ebos_simulator,
209 const SummaryState &summary_state,
211 std::vector<double> &potentials,
214 void computeWellRatesWithThpAlqProd(
215 const Simulator &ebos_simulator,
216 const SummaryState &summary_state,
218 std::vector<double> &potentials,
221 virtual std::optional<double> computeBhpAtThpLimitProdWithAlq(
222 const Simulator& ebos_simulator,
223 const SummaryState& summary_state,
225 double alq_value)
const override;
227 virtual void computeWellRatesWithBhp(
228 const Simulator& ebosSimulator,
230 std::vector<double>& well_flux,
234 using Base::phaseUsage;
235 using Base::vfp_properties_;
240 void computeConnLevelProdInd(
const FluidState& fs,
241 const std::function<
double(
const double)>& connPICalc,
242 const std::vector<EvalWell>& mobility,
243 double* connPI)
const;
245 void computeConnLevelInjInd(
const typename StandardWell<TypeTag>::FluidState& fs,
246 const Phase preferred_phase,
247 const std::function<
double(
const double)>& connIICalc,
248 const std::vector<EvalWell>& mobility,
257 void recoverSolutionWell(
const BVector& x, BVectorWell& xw)
const;
260 void updateWellState(
const BVectorWell& dwells,
266 void computePropertiesForWellConnectionPressures(
const Simulator& ebosSimulator,
268 std::vector<double>& b_perf,
269 std::vector<double>& rsmax_perf,
270 std::vector<double>& rvmax_perf,
271 std::vector<double>& rvwmax_perf,
272 std::vector<double>& surf_dens_perf)
const;
274 void computeWellConnectionDensitesPressures(
const Simulator& ebosSimulator,
276 const std::vector<double>& b_perf,
277 const std::vector<double>& rsmax_perf,
278 const std::vector<double>& rvmax_perf,
279 const std::vector<double>& rvwmax_perf,
280 const std::vector<double>& surf_dens_perf,
283 void computeWellConnectionPressures(
const Simulator& ebosSimulator,
287 void computePerfRateEval(
const IntensiveQuantities& intQuants,
288 const std::vector<EvalWell>& mob,
293 std::vector<EvalWell>& cq_s,
294 double& perf_dis_gas_rate,
295 double& perf_vap_oil_rate,
296 double& perf_vap_wat_rate,
299 void computePerfRateScalar(
const IntensiveQuantities& intQuants,
300 const std::vector<Scalar>& mob,
305 std::vector<Scalar>& cq_s,
308 template<
class Value>
309 void computePerfRate(
const std::vector<Value>& mob,
310 const Value& pressure,
315 std::vector<Value>& b_perfcells_dense,
319 const Value& skin_pressure,
320 const std::vector<Value>& cmix_s,
321 std::vector<Value>& cq_s,
322 double& perf_dis_gas_rate,
323 double& perf_vap_oil_rate,
324 double& perf_vap_wat_rate,
327 void computeWellRatesWithBhpIterations(
const Simulator& ebosSimulator,
329 std::vector<double>& well_flux,
332 std::vector<double> computeWellPotentialWithTHP(
333 const Simulator& ebosSimulator,
338 virtual double getRefDensity()
const override;
341 void getMobilityEval(
const Simulator& ebosSimulator,
343 std::vector<EvalWell>& mob,
347 void getMobilityScalar(
const Simulator& ebosSimulator,
349 std::vector<Scalar>& mob,
353 void updateWaterMobilityWithPolymer(
const Simulator& ebos_simulator,
355 std::vector<EvalWell>& mob_water,
358 void updatePrimaryVariablesNewton(
const BVectorWell& dwells,
363 void updateExtraPrimaryVariables(
const BVectorWell& dwells)
const;
368 virtual void assembleWellEqWithoutIteration(
const Simulator& ebosSimulator,
370 const Well::InjectionControls& inj_controls,
371 const Well::ProductionControls& prod_controls,
376 void assembleWellEqWithoutIterationImpl(
const Simulator& ebosSimulator,
382 void calculateSinglePerf(
const Simulator& ebosSimulator,
385 std::vector<RateVector>& connectionRates,
386 std::vector<EvalWell>& cq_s,
387 EvalWell& water_flux_s,
388 EvalWell& cq_s_zfrac_effective,
392 virtual void checkOperabilityUnderBHPLimit(
const WellState& well_state,
const Simulator& ebos_simulator,
DeferredLogger& deferred_logger)
override;
395 virtual void checkOperabilityUnderTHPLimit(
const Simulator& ebos_simulator,
const WellState& well_state,
DeferredLogger& deferred_logger)
override;
398 virtual void updateIPR(
const Simulator& ebos_simulator,
DeferredLogger& deferred_logger)
const override;
402 bool allDrawDownWrongDirection(
const Simulator& ebos_simulator)
const;
405 bool canProduceInjectWithCurrentBhp(
const Simulator& ebos_simulator,
414 bool openCrossFlowAvoidSingularity(
const Simulator& ebos_simulator)
const;
420 EvalWell pskin(
const double throuhgput,
421 const EvalWell& water_velocity,
422 const EvalWell& poly_inj_conc,
426 EvalWell pskinwater(
const double throughput,
427 const EvalWell& water_velocity,
431 EvalWell wpolymermw(
const double throughput,
432 const EvalWell& water_velocity,
436 void handleInjectivityRate(
const Simulator& ebosSimulator,
438 std::vector<EvalWell>& cq_s)
const;
441 void handleInjectivityEquations(
const Simulator& ebosSimulator,
444 const EvalWell& water_flux_s,
447 virtual void updateWaterThroughput(
const double dt,
WellState& well_state)
const override;
450 void checkConvergenceExtraEqs(
const std::vector<double>& res,
454 void updateConnectionRatePolyMW(
const EvalWell& cq_s_poly,
455 const IntensiveQuantities& int_quants,
458 std::vector<RateVector>& connectionRates,
462 std::optional<double> computeBhpAtThpLimitProd(
const WellState& well_state,
463 const Simulator& ebos_simulator,
464 const SummaryState& summary_state,
467 std::optional<double> computeBhpAtThpLimitInj(
const Simulator& ebos_simulator,
468 const SummaryState& summary_state,
475 #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: GroupState.hpp:34
Class encapsulating some information about parallel wells.
Definition: ParallelWellInfo.hpp:243
Definition: StandardWellEval.hpp:47
Definition: StandardWell.hpp:63
virtual bool jacobianContainsWellContributions() const override
Wether the Jacobian will also have well contributions in it.
Definition: StandardWell.hpp:202
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:2709
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:1465
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:1777
virtual void apply(const BVector &x, BVector &Ax) const override
Ax = Ax - C D^-1 B x.
Definition: StandardWell_impl.hpp:1713
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:1963
const std::string & name() const
Well name.
Definition: WellInterfaceGeneric.cpp:134
Definition: WellInterface.hpp:72
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:27
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:46