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>;
99 static constexpr
int numWellConservationEq = Indices::numPhases + Indices::numSolvents;
101 static constexpr
int numWellControlEq = 1;
104 static constexpr
int numStaticWellEq = numWellConservationEq + numWellControlEq;
109 static constexpr
int Bhp = numStaticWellEq - numWellControlEq;
111 using typename Base::Scalar;
119 using typename Base::BVector;
121 using Eval =
typename StdWellEval::Eval;
122 using EvalWell =
typename StdWellEval::EvalWell;
123 using BVectorWell =
typename StdWellEval::BVectorWell;
129 const RateConverterType& rate_converter,
130 const int pvtRegionIdx,
131 const int num_components,
132 const int num_phases,
133 const int index_of_well,
134 const std::vector<PerforationData>& perf_data);
136 virtual void init(
const PhaseUsage* phase_usage_arg,
137 const std::vector<double>& depth_arg,
138 const double gravity_arg,
140 const std::vector< Scalar >& B_avg)
override;
143 virtual void initPrimaryVariablesEvaluation()
const override;
147 const std::vector<double>& B_avg,
149 const bool relax_tolerance =
false)
const override;
152 virtual void apply(
const BVector& x, BVector& Ax)
const override;
154 virtual void apply(BVector& r)
const override;
165 std::vector<double>& well_potentials,
168 virtual void updatePrimaryVariables(
const WellState& well_state,
DeferredLogger& deferred_logger)
const override;
172 virtual void calculateExplicitQuantities(
const Simulator& ebosSimulator,
176 virtual void updateProductivityIndex(
const Simulator& ebosSimulator,
181 virtual void addWellContributions(SparseMatrixAdapter& mat)
const override;
184 bool iterateWellEqWithControl(
const Simulator& ebosSimulator,
186 const Well::InjectionControls& inj_controls,
187 const Well::ProductionControls& prod_controls,
199 double computeWellRatesAndBhpWithThpAlqProd(
const Simulator &ebos_simulator,
200 const SummaryState &summary_state,
202 std::vector<double> &potentials,
205 void computeWellRatesWithThpAlqProd(
206 const Simulator &ebos_simulator,
207 const SummaryState &summary_state,
209 std::vector<double> &potentials,
212 virtual std::optional<double> computeBhpAtThpLimitProdWithAlq(
213 const Simulator& ebos_simulator,
214 const SummaryState& summary_state,
216 double alq_value)
const override;
218 virtual void computeWellRatesWithBhp(
219 const Simulator& ebosSimulator,
221 std::vector<double>& well_flux,
225 using Base::phaseUsage;
226 using Base::vfp_properties_;
231 void computeConnLevelProdInd(
const FluidState& fs,
232 const std::function<
double(
const double)>& connPICalc,
233 const std::vector<EvalWell>& mobility,
234 double* connPI)
const;
236 void computeConnLevelInjInd(
const typename StandardWell<TypeTag>::FluidState& fs,
237 const Phase preferred_phase,
238 const std::function<
double(
const double)>& connIICalc,
239 const std::vector<EvalWell>& mobility,
246 void recoverSolutionWell(
const BVector& x, BVectorWell& xw)
const;
249 void updateWellState(
const BVectorWell& dwells,
255 void computePropertiesForWellConnectionPressures(
const Simulator& ebosSimulator,
257 std::vector<double>& b_perf,
258 std::vector<double>& rsmax_perf,
259 std::vector<double>& rvmax_perf,
260 std::vector<double>& surf_dens_perf)
const;
262 void computeWellConnectionDensitesPressures(
const Simulator& ebosSimulator,
264 const std::vector<double>& b_perf,
265 const std::vector<double>& rsmax_perf,
266 const std::vector<double>& rvmax_perf,
267 const std::vector<double>& surf_dens_perf,
270 void computeWellConnectionPressures(
const Simulator& ebosSimulator,
274 void computePerfRateEval(
const IntensiveQuantities& intQuants,
275 const std::vector<EvalWell>& mob,
280 std::vector<EvalWell>& cq_s,
281 double& perf_dis_gas_rate,
282 double& perf_vap_oil_rate,
285 void computePerfRateScalar(
const IntensiveQuantities& intQuants,
286 const std::vector<Scalar>& mob,
291 std::vector<Scalar>& cq_s,
294 template<
class Value>
295 void computePerfRate(
const std::vector<Value>& mob,
296 const Value& pressure,
300 std::vector<Value>& b_perfcells_dense,
304 const Value& skin_pressure,
305 const std::vector<Value>& cmix_s,
306 std::vector<Value>& cq_s,
307 double& perf_dis_gas_rate,
308 double& perf_vap_oil_rate,
311 void computeWellRatesWithBhpIterations(
const Simulator& ebosSimulator,
313 std::vector<double>& well_flux,
316 std::vector<double> computeWellPotentialWithTHP(
317 const Simulator& ebosSimulator,
322 virtual double getRefDensity()
const override;
325 void getMobilityEval(
const Simulator& ebosSimulator,
327 std::vector<EvalWell>& mob,
331 void getMobilityScalar(
const Simulator& ebosSimulator,
333 std::vector<Scalar>& mob,
337 void updateWaterMobilityWithPolymer(
const Simulator& ebos_simulator,
339 std::vector<EvalWell>& mob_water,
342 void updatePrimaryVariablesNewton(
const BVectorWell& dwells,
347 void updateExtraPrimaryVariables(
const BVectorWell& dwells)
const;
352 virtual void assembleWellEqWithoutIteration(
const Simulator& ebosSimulator,
354 const Well::InjectionControls& inj_controls,
355 const Well::ProductionControls& prod_controls,
360 void assembleWellEqWithoutIterationImpl(
const Simulator& ebosSimulator,
366 void calculateSinglePerf(
const Simulator& ebosSimulator,
369 std::vector<RateVector>& connectionRates,
370 std::vector<EvalWell>& cq_s,
371 EvalWell& water_flux_s,
372 EvalWell& cq_s_zfrac_effective,
376 virtual void checkOperabilityUnderBHPLimit(
const WellState& well_state,
const Simulator& ebos_simulator,
DeferredLogger& deferred_logger)
override;
379 virtual void checkOperabilityUnderTHPLimit(
const Simulator& ebos_simulator,
const WellState& well_state,
DeferredLogger& deferred_logger)
override;
382 virtual void updateIPR(
const Simulator& ebos_simulator,
DeferredLogger& deferred_logger)
const override;
386 bool allDrawDownWrongDirection(
const Simulator& ebos_simulator)
const;
389 bool canProduceInjectWithCurrentBhp(
const Simulator& ebos_simulator,
398 bool openCrossFlowAvoidSingularity(
const Simulator& ebos_simulator)
const;
404 EvalWell pskin(
const double throuhgput,
405 const EvalWell& water_velocity,
406 const EvalWell& poly_inj_conc,
410 EvalWell pskinwater(
const double throughput,
411 const EvalWell& water_velocity,
415 EvalWell wpolymermw(
const double throughput,
416 const EvalWell& water_velocity,
420 void handleInjectivityRate(
const Simulator& ebosSimulator,
422 std::vector<EvalWell>& cq_s)
const;
425 void handleInjectivityEquations(
const Simulator& ebosSimulator,
428 const EvalWell& water_flux_s,
431 virtual void updateWaterThroughput(
const double dt,
WellState& well_state)
const override;
434 void checkConvergenceExtraEqs(
const std::vector<double>& res,
438 void updateConnectionRatePolyMW(
const EvalWell& cq_s_poly,
439 const IntensiveQuantities& int_quants,
442 std::vector<RateVector>& connectionRates,
446 std::optional<double> computeBhpAtThpLimitProd(
const WellState& well_state,
447 const Simulator& ebos_simulator,
448 const SummaryState& summary_state,
451 std::optional<double> computeBhpAtThpLimitInj(
const Simulator& ebos_simulator,
452 const SummaryState& summary_state,
459 #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:252
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:193
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:2479
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:1364
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:1664
virtual void apply(const BVector &x, BVector &Ax) const override
Ax = Ax - C D^-1 B x.
Definition: StandardWell_impl.hpp:1600
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:1850
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