23#ifndef OPM_STANDARDWELL_HEADER_INCLUDED
24#define OPM_STANDARDWELL_HEADER_INCLUDED
26#include <opm/simulators/timestepping/ConvergenceReport.hpp>
28#include <opm/simulators/wells/VFPInjProperties.hpp>
29#include <opm/simulators/wells/VFPProdProperties.hpp>
30#include <opm/simulators/wells/WellInterface.hpp>
31#include <opm/simulators/wells/WellProdIndexCalculator.hpp>
32#include <opm/simulators/wells/ParallelWellInfo.hpp>
34#include <opm/models/blackoil/blackoilpolymermodules.hh>
35#include <opm/models/blackoil/blackoilsolventmodules.hh>
36#include <opm/models/blackoil/blackoilextbomodules.hh>
37#include <opm/models/blackoil/blackoilfoammodules.hh>
38#include <opm/models/blackoil/blackoilbrinemodules.hh>
39#include <opm/models/blackoil/blackoilmicpmodules.hh>
41#include <opm/material/densead/Evaluation.hpp>
42#include <opm/input/eclipse/Schedule/ScheduleTypes.hpp>
44#include <opm/simulators/wells/StandardWellEval.hpp>
46#include <dune/common/dynvector.hh>
47#include <dune/common/dynmatrix.hh>
55 template<
typename TypeTag>
58 GetPropType<TypeTag, Properties::Indices>,
59 GetPropType<TypeTag, Properties::Scalar>>
65 GetPropType<TypeTag, Properties::Indices>,
66 GetPropType<TypeTag, Properties::Scalar>>;
71 using typename Base::Simulator;
72 using typename Base::IntensiveQuantities;
73 using typename Base::FluidSystem;
74 using typename Base::MaterialLaw;
76 using typename Base::Indices;
77 using typename Base::RateConverterType;
78 using typename Base::SparseMatrixAdapter;
79 using typename Base::FluidState;
80 using typename Base::RateVector;
82 using Base::has_solvent;
83 using Base::has_zFraction;
84 using Base::has_polymer;
85 using Base::has_polymermw;
87 using Base::has_brine;
88 using Base::has_energy;
91 using PolymerModule = BlackOilPolymerModule<TypeTag>;
92 using FoamModule = BlackOilFoamModule<TypeTag>;
93 using BrineModule = BlackOilBrineModule<TypeTag>;
94 using typename Base::PressureMatrix;
97 static constexpr int numWellConservationEq = Indices::numPhases + Indices::numSolvents;
99 static constexpr int numWellControlEq = 1;
102 static constexpr int numStaticWellEq = numWellConservationEq + numWellControlEq;
107 static constexpr int Bhp = numStaticWellEq - numWellControlEq;
109 using StdWellEval::WQTotal;
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,
141 const bool changed_to_open_this_step)
override;
144 void initPrimaryVariablesEvaluation()
override;
148 const std::vector<double>& B_avg,
150 const bool relax_tolerance =
false)
const override;
153 virtual void apply(
const BVector& x, BVector& Ax)
const override;
155 virtual void apply(BVector& r)
const override;
167 std::vector<double>& well_potentials,
170 void updatePrimaryVariables(
const SummaryState& summary_state,
174 void solveEqAndUpdateWellState(
const SummaryState& summary_state,
178 virtual void calculateExplicitQuantities(
const Simulator& ebosSimulator,
182 virtual void updateProductivityIndex(
const Simulator& ebosSimulator,
187 virtual void addWellContributions(SparseMatrixAdapter& mat)
const override;
189 virtual void addWellPressureEquations(PressureMatrix& mat,
191 const int pressureVarIndex,
192 const bool use_well_weights,
193 const WellState& well_state)
const override;
196 bool iterateWellEqWithControl(
const Simulator& ebosSimulator,
198 const Well::InjectionControls& inj_controls,
199 const Well::ProductionControls& prod_controls,
211 double computeWellRatesAndBhpWithThpAlqProd(
const Simulator &ebos_simulator,
212 const SummaryState &summary_state,
214 std::vector<double> &potentials,
217 void computeWellRatesWithThpAlqProd(
218 const Simulator &ebos_simulator,
219 const SummaryState &summary_state,
221 std::vector<double> &potentials,
224 std::optional<double> computeBhpAtThpLimitProdWithAlq(
225 const Simulator& ebos_simulator,
226 const SummaryState& summary_state,
227 const double alq_value,
230 virtual void computeWellRatesWithBhp(
231 const Simulator& ebosSimulator,
233 std::vector<double>& well_flux,
237 using Base::phaseUsage;
238 using Base::vfp_properties_;
243 void computeConnLevelProdInd(
const FluidState& fs,
244 const std::function<
double(
const double)>& connPICalc,
245 const std::vector<EvalWell>& mobility,
246 double* connPI)
const;
248 void computeConnLevelInjInd(
const typename StandardWell<TypeTag>::FluidState& fs,
249 const Phase preferred_phase,
250 const std::function<
double(
const double)>& connIICalc,
251 const std::vector<EvalWell>& mobility,
260 void updateWellState(
const SummaryState& summary_state,
261 const BVectorWell& dwells,
267 void computePropertiesForWellConnectionPressures(
const Simulator& ebosSimulator,
269 std::vector<double>& b_perf,
270 std::vector<double>& rsmax_perf,
271 std::vector<double>& rvmax_perf,
272 std::vector<double>& rvwmax_perf,
273 std::vector<double>& rswmax_perf,
274 std::vector<double>& surf_dens_perf)
const;
276 void computeWellConnectionDensitesPressures(
const Simulator& ebosSimulator,
278 const std::vector<double>& b_perf,
279 const std::vector<double>& rsmax_perf,
280 const std::vector<double>& rvmax_perf,
281 const std::vector<double>& rvwmax_perf,
282 const std::vector<double>& rswmax_perf,
283 const std::vector<double>& surf_dens_perf,
286 void computeWellConnectionPressures(
const Simulator& ebosSimulator,
290 void computePerfRateEval(
const IntensiveQuantities& intQuants,
291 const std::vector<EvalWell>& mob,
296 std::vector<EvalWell>& cq_s,
297 double& perf_dis_gas_rate,
298 double& perf_dis_gas_rate_in_water,
299 double& perf_vap_oil_rate,
300 double& perf_vap_wat_rate,
303 void computePerfRateScalar(
const IntensiveQuantities& intQuants,
304 const std::vector<Scalar>& mob,
309 std::vector<Scalar>& cq_s,
312 template<
class Value>
313 void computePerfRate(
const std::vector<Value>& mob,
314 const Value& pressure,
320 std::vector<Value>& b_perfcells_dense,
324 const Value& skin_pressure,
325 const std::vector<Value>& cmix_s,
326 std::vector<Value>& cq_s,
327 double& perf_dis_gas_rate,
328 double& perf_dis_gas_rate_in_water,
329 double& perf_vap_oil_rate,
330 double& perf_vap_wat_rate,
333 void computeWellRatesWithBhpIterations(
const Simulator& ebosSimulator,
335 std::vector<double>& well_flux,
338 std::vector<double> computeWellPotentialWithTHP(
339 const Simulator& ebosSimulator,
343 bool updateWellStateWithTHPTargetProd(
const Simulator& ebos_simulator,
347 virtual double getRefDensity()
const override;
350 void getMobilityEval(
const Simulator& ebosSimulator,
352 std::vector<EvalWell>& mob,
356 void getMobilityScalar(
const Simulator& ebosSimulator,
358 std::vector<Scalar>& mob,
362 void updateWaterMobilityWithPolymer(
const Simulator& ebos_simulator,
364 std::vector<EvalWell>& mob_water,
367 void updatePrimaryVariablesNewton(
const BVectorWell& dwells,
368 const bool stop_or_zero_rate_target,
371 void updateWellStateFromPrimaryVariables(
const bool stop_or_zero_rate_target,
375 virtual void assembleWellEqWithoutIteration(
const Simulator& ebosSimulator,
377 const Well::InjectionControls& inj_controls,
378 const Well::ProductionControls& prod_controls,
383 void assembleWellEqWithoutIterationImpl(
const Simulator& ebosSimulator,
385 const Well::InjectionControls& inj_controls,
386 const Well::ProductionControls& prod_controls,
391 void calculateSinglePerf(
const Simulator& ebosSimulator,
394 std::vector<RateVector>& connectionRates,
395 std::vector<EvalWell>& cq_s,
396 EvalWell& water_flux_s,
397 EvalWell& cq_s_zfrac_effective,
401 virtual void checkOperabilityUnderBHPLimit(
const WellState& well_state,
const Simulator& ebos_simulator,
DeferredLogger& deferred_logger)
override;
404 virtual void checkOperabilityUnderTHPLimit(
const Simulator& ebos_simulator,
const WellState& well_state,
DeferredLogger& deferred_logger)
override;
407 virtual void updateIPR(
const Simulator& ebos_simulator,
DeferredLogger& deferred_logger)
const override;
411 bool allDrawDownWrongDirection(
const Simulator& ebos_simulator)
const;
414 bool canProduceInjectWithCurrentBhp(
const Simulator& ebos_simulator,
423 bool openCrossFlowAvoidSingularity(
const Simulator& ebos_simulator)
const;
429 EvalWell pskin(
const double throuhgput,
430 const EvalWell& water_velocity,
431 const EvalWell& poly_inj_conc,
435 EvalWell pskinwater(
const double throughput,
436 const EvalWell& water_velocity,
440 EvalWell wpolymermw(
const double throughput,
441 const EvalWell& water_velocity,
445 void handleInjectivityRate(
const Simulator& ebosSimulator,
447 std::vector<EvalWell>& cq_s)
const;
450 void handleInjectivityEquations(
const Simulator& ebosSimulator,
453 const EvalWell& water_flux_s,
456 virtual void updateWaterThroughput(
const double dt,
WellState& well_state)
const override;
459 void checkConvergenceExtraEqs(
const std::vector<double>& res,
463 void updateConnectionRatePolyMW(
const EvalWell& cq_s_poly,
464 const IntensiveQuantities& int_quants,
467 std::vector<RateVector>& connectionRates,
471 std::optional<double> computeBhpAtThpLimitProd(
const WellState& well_state,
472 const Simulator& ebos_simulator,
473 const SummaryState& summary_state,
476 std::optional<double> computeBhpAtThpLimitInj(
const Simulator& ebos_simulator,
477 const SummaryState& summary_state,
484#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:38
Definition: DeferredLogger.hpp:57
Definition: GroupState.hpp:34
Class encapsulating some information about parallel wells.
Definition: ParallelWellInfo.hpp:184
Definition: StandardWellEval.hpp:49
Definition: StandardWell.hpp:60
virtual bool jacobianContainsWellContributions() const override
Wether the Jacobian will also have well contributions in it.
Definition: StandardWell.hpp:205
void recoverWellSolutionAndUpdateWellState(const SummaryState &summary_state, const BVector &x, WellState &well_state, DeferredLogger &deferred_logger) override
using the solution x to recover the solution xw for wells and applying xw to update Well State
Definition: StandardWell_impl.hpp:1659
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:2544
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:1407
virtual void apply(const BVector &x, BVector &Ax) const override
Ax = Ax - C D^-1 B x.
Definition: StandardWell_impl.hpp:1627
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:1885
const std::string & name() const
Well name.
Definition: WellInterfaceGeneric.cpp:146
Definition: WellInterface.hpp:74
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:60
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