22 #ifndef OPM_MULTISEGMENTWELL_HEADER_INCLUDED
23 #define OPM_MULTISEGMENTWELL_HEADER_INCLUDED
25 #include <opm/simulators/wells/WellInterface.hpp>
26 #include <opm/simulators/wells/MultisegmentWellEval.hpp>
28 #include <opm/parser/eclipse/EclipseState/Runspec.hpp>
34 template<
typename TypeTag>
37 GetPropType<TypeTag, Properties::Indices>,
38 GetPropType<TypeTag, Properties::Scalar>>
43 GetPropType<TypeTag, Properties::Indices>,
44 GetPropType<TypeTag, Properties::Scalar>>;
46 using typename Base::Simulator;
47 using typename Base::IntensiveQuantities;
48 using typename Base::FluidSystem;
50 using typename Base::MaterialLaw;
51 using typename Base::Indices;
52 using typename Base::RateConverterType;
53 using typename Base::SparseMatrixAdapter;
54 using typename Base::FluidState;
56 using typename Base::GLiftProdWells;
57 using typename Base::GLiftOptWells;
58 using typename Base::GLiftWellStateMap;
59 using typename Base::GLiftSyncGroups;
61 using Base::has_solvent;
62 using Base::has_polymer;
67 using typename Base::Scalar;
70 using typename Base::BVector;
71 using typename Base::Eval;
73 using typename MSWEval::EvalWell;
74 using typename MSWEval::BVectorWell;
75 using typename MSWEval::DiagMatWell;
76 using typename MSWEval::OffDiagMatrixBlockWellType;
79 using MSWEval::GTotal;
81 using MSWEval::numWellEq;
87 const RateConverterType& rate_converter,
88 const int pvtRegionIdx,
89 const int num_components,
91 const int index_of_well,
92 const std::vector<PerforationData>& perf_data);
94 virtual void init(
const PhaseUsage* phase_usage_arg,
95 const std::vector<double>& depth_arg,
96 const double gravity_arg,
98 const std::vector< Scalar >& B_avg)
override;
100 virtual void initPrimaryVariablesEvaluation()
const override;
102 virtual void gasLiftOptimizationStage1 (
124 const std::vector<double>& B_avg,
126 const bool relax_tolerance =
false)
const override;
129 virtual void apply(
const BVector& x, BVector& Ax)
const override;
131 virtual void apply(BVector& r)
const override;
142 std::vector<double>& well_potentials,
145 virtual void updatePrimaryVariables(
const WellState& well_state,
DeferredLogger& deferred_logger)
const override;
149 virtual void calculateExplicitQuantities(
const Simulator& ebosSimulator,
153 virtual void updateProductivityIndex(
const Simulator& ebosSimulator,
158 virtual void addWellContributions(SparseMatrixAdapter& jacobian)
const override;
163 void computeConnLevelProdInd(
const FluidState& fs,
164 const std::function<
double(
const double)>& connPICalc,
165 const std::vector<Scalar>& mobility,
166 double* connPI)
const;
168 void computeConnLevelInjInd(
const FluidState& fs,
169 const Phase preferred_phase,
170 const std::function<
double(
const double)>& connIICalc,
171 const std::vector<Scalar>& mobility,
176 int number_segments_;
179 WellSegments::CompPressureDrop compPressureDrop()
const;
181 WellSegments::MultiPhaseModel multiphaseModel()
const;
184 std::vector<std::vector<double> > segment_fluid_initial_;
186 mutable int debug_cost_counter_ = 0;
189 void updateWellState(
const BVectorWell& dwells,
192 const double relaxation_factor=1.0)
const;
196 void computeInitialSegmentFluids(
const Simulator& ebos_simulator);
199 void computePerfCellPressDiffs(
const Simulator& ebosSimulator);
201 void computePerfRateScalar(
const IntensiveQuantities& int_quants,
202 const std::vector<Scalar>& mob_perfcells,
206 const Scalar& segment_pressure,
207 const bool& allow_cf,
208 std::vector<Scalar>& cq_s,
211 void computePerfRateEval(
const IntensiveQuantities& int_quants,
212 const std::vector<EvalWell>& mob_perfcells,
216 const EvalWell& segment_pressure,
217 const bool& allow_cf,
218 std::vector<EvalWell>& cq_s,
219 EvalWell& perf_press,
220 double& perf_dis_gas_rate,
221 double& perf_vap_oil_rate,
224 template<
class Value>
225 void computePerfRate(
const Value& pressure_cell,
228 const std::vector<Value>& b_perfcells,
229 const std::vector<Value>& mob_perfcells,
232 const Value& segment_pressure,
233 const Value& segment_density,
234 const bool& allow_cf,
235 const std::vector<Value>& cmix_s,
236 std::vector<Value>& cq_s,
238 double& perf_dis_gas_rate,
239 double& perf_vap_oil_rate,
244 void computeSegmentFluidProperties(
const Simulator& ebosSimulator);
247 void getMobilityEval(
const Simulator& ebosSimulator,
249 std::vector<EvalWell>& mob)
const;
252 void getMobilityScalar(
const Simulator& ebosSimulator,
254 std::vector<Scalar>& mob)
const;
256 void computeWellRatesAtBhpLimit(
const Simulator& ebosSimulator,
257 std::vector<double>& well_flux,
260 void computeWellRatesWithBhp(
const Simulator& ebosSimulator,
262 std::vector<double>& well_flux,
265 void computeWellRatesWithBhpIterations(
const Simulator& ebosSimulator,
267 std::vector<double>& well_flux,
271 computeWellPotentialWithTHP(
const Simulator& ebos_simulator,
274 virtual double getRefDensity()
const override;
276 virtual bool iterateWellEqWithControl(
const Simulator& ebosSimulator,
278 const Well::InjectionControls& inj_controls,
279 const Well::ProductionControls& prod_controls,
284 virtual void assembleWellEqWithoutIteration(
const Simulator& ebosSimulator,
286 const Well::InjectionControls& inj_controls,
287 const Well::ProductionControls& prod_controls,
292 virtual void updateWaterThroughput(
const double dt,
WellState& well_state)
const override;
294 EvalWell getSegmentSurfaceVolume(
const Simulator& ebos_simulator,
const int seg_idx)
const;
301 bool openCrossFlowAvoidSingularity(
const Simulator& ebos_simulator)
const;
305 bool allDrawDownWrongDirection(
const Simulator& ebos_simulator)
const;
308 std::optional<double> computeBhpAtThpLimitProd(
const Simulator& ebos_simulator,
309 const SummaryState& summary_state,
312 std::optional<double> computeBhpAtThpLimitInj(
const Simulator& ebos_simulator,
313 const SummaryState& summary_state,
316 double maxPerfPress(
const Simulator& ebos_simulator)
const;
319 virtual void checkOperabilityUnderBHPLimitProducer(
const WellState& well_state,
const Simulator& ebos_simulator,
DeferredLogger& deferred_logger)
override;
322 virtual void checkOperabilityUnderTHPLimitProducer(
const Simulator& ebos_simulator,
const WellState& well_state,
DeferredLogger& deferred_logger)
override;
325 virtual void updateIPR(
const Simulator& ebos_simulator,
DeferredLogger& deferred_logger)
const override;
330 #include "MultisegmentWell_impl.hpp"
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
Definition: MultisegmentWellEval.hpp:52
Definition: MultisegmentWell.hpp:39
virtual void updateWellStateWithTarget(const Simulator &ebos_simulator, const GroupState &group_state, WellState &well_state, DeferredLogger &deferred_logger) const override
updating the well state based the current control mode
Definition: MultisegmentWell_impl.hpp:142
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: MultisegmentWell_impl.hpp:226
virtual void apply(const BVector &x, BVector &Ax) const override
Ax = Ax - C D^-1 B x.
Definition: MultisegmentWell_impl.hpp:184
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: MultisegmentWell_impl.hpp:1848
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: MultisegmentWell_impl.hpp:244
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: MultisegmentWell_impl.hpp:161
Class encapsulating some information about parallel wells.
Definition: ParallelWellInfo.hpp:252
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
Definition: BlackoilPhases.hpp:45