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/input/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 Base::has_solvent;
57 using Base::has_polymer;
62 using typename Base::Scalar;
65 using typename Base::BVector;
66 using typename Base::Eval;
68 using typename MSWEval::EvalWell;
69 using typename MSWEval::BVectorWell;
70 using typename MSWEval::DiagMatWell;
71 using typename MSWEval::OffDiagMatrixBlockWellType;
74 using MSWEval::WQTotal;
76 using MSWEval::numWellEq;
77 using typename Base::PressureMatrix;
83 const RateConverterType& rate_converter,
84 const int pvtRegionIdx,
85 const int num_components,
87 const int index_of_well,
88 const std::vector<PerforationData>& perf_data);
90 virtual void init(
const PhaseUsage* phase_usage_arg,
91 const std::vector<double>& depth_arg,
92 const double gravity_arg,
94 const std::vector< Scalar >& B_avg,
95 const bool changed_to_open_this_step)
override;
97 virtual void initPrimaryVariablesEvaluation()
const override;
107 const std::vector<double>& B_avg,
109 const bool relax_tolerance =
false)
const override;
112 virtual void apply(
const BVector& x, BVector& Ax)
const override;
114 virtual void apply(BVector& r)
const override;
125 std::vector<double>& well_potentials,
128 virtual void updatePrimaryVariables(
const WellState& well_state,
DeferredLogger& deferred_logger)
const override;
132 virtual void calculateExplicitQuantities(
const Simulator& ebosSimulator,
136 virtual void updateProductivityIndex(
const Simulator& ebosSimulator,
141 virtual void addWellContributions(SparseMatrixAdapter& jacobian)
const override;
143 virtual void addWellPressureEquations(PressureMatrix& mat,
145 const int pressureVarIndex,
146 const bool use_well_weights,
147 const WellState& well_state)
const override;
152 void computeConnLevelProdInd(
const FluidState& fs,
153 const std::function<
double(
const double)>& connPICalc,
154 const std::vector<Scalar>& mobility,
155 double* connPI)
const;
157 void computeConnLevelInjInd(
const FluidState& fs,
158 const Phase preferred_phase,
159 const std::function<
double(
const double)>& connIICalc,
160 const std::vector<Scalar>& mobility,
164 virtual std::optional<double> computeBhpAtThpLimitProdWithAlq(
165 const Simulator& ebos_simulator,
166 const SummaryState& summary_state,
168 double alq_value)
const override;
171 int number_segments_;
177 WellSegments::CompPressureDrop compPressureDrop()
const;
179 WellSegments::MultiPhaseModel multiphaseModel()
const;
182 std::vector<std::vector<double> > segment_fluid_initial_;
184 mutable int debug_cost_counter_ = 0;
187 void updateWellState(
const BVectorWell& dwells,
190 const double relaxation_factor=1.0)
const;
194 void computeInitialSegmentFluids(
const Simulator& ebos_simulator);
197 void computePerfCellPressDiffs(
const Simulator& ebosSimulator);
199 void computePerfRateScalar(
const IntensiveQuantities& int_quants,
200 const std::vector<Scalar>& mob_perfcells,
204 const Scalar& segment_pressure,
205 const bool& allow_cf,
206 std::vector<Scalar>& cq_s,
209 void computePerfRateEval(
const IntensiveQuantities& int_quants,
210 const std::vector<EvalWell>& mob_perfcells,
214 const EvalWell& segment_pressure,
215 const bool& allow_cf,
216 std::vector<EvalWell>& cq_s,
217 EvalWell& perf_press,
218 double& perf_dis_gas_rate,
219 double& perf_vap_oil_rate,
222 template<
class Value>
223 void computePerfRate(
const Value& pressure_cell,
226 const std::vector<Value>& b_perfcells,
227 const std::vector<Value>& mob_perfcells,
230 const Value& segment_pressure,
231 const Value& segment_density,
232 const bool& allow_cf,
233 const std::vector<Value>& cmix_s,
234 std::vector<Value>& cq_s,
236 double& perf_dis_gas_rate,
237 double& perf_vap_oil_rate,
242 void computeSegmentFluidProperties(
const Simulator& ebosSimulator,
246 void getMobilityEval(
const Simulator& ebosSimulator,
248 std::vector<EvalWell>& mob)
const;
251 void getMobilityScalar(
const Simulator& ebosSimulator,
253 std::vector<Scalar>& mob)
const;
255 void computeWellRatesAtBhpLimit(
const Simulator& ebosSimulator,
256 std::vector<double>& well_flux,
259 virtual void computeWellRatesWithBhp(
const Simulator& ebosSimulator,
261 std::vector<double>& well_flux,
264 void computeWellRatesWithBhpIterations(
const Simulator& ebosSimulator,
266 std::vector<double>& well_flux,
269 std::vector<double> computeWellPotentialWithTHP(
271 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;
309 std::optional<double> computeBhpAtThpLimitProd(
311 const Simulator& ebos_simulator,
312 const SummaryState& summary_state,
315 std::optional<double> computeBhpAtThpLimitInj(
const Simulator& ebos_simulator,
316 const SummaryState& summary_state,
319 double maxPerfPress(
const Simulator& ebos_simulator)
const;
322 virtual void checkOperabilityUnderBHPLimit(
const WellState& well_state,
const Simulator& ebos_simulator,
DeferredLogger& deferred_logger)
override;
325 virtual void checkOperabilityUnderTHPLimit(
const Simulator& ebos_simulator,
const WellState& well_state,
DeferredLogger& deferred_logger)
override;
328 virtual void updateIPR(
const Simulator& ebos_simulator,
DeferredLogger& deferred_logger)
const override;
333 #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: GroupState.hpp:34
Definition: MultisegmentWellEval.hpp:55
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:153
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:237
virtual void apply(const BVector &x, BVector &Ax) const override
Ax = Ax - C D^-1 B x.
Definition: MultisegmentWell_impl.hpp:195
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:2003
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:255
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:172
Class encapsulating some information about parallel wells.
Definition: ParallelWellInfo.hpp:243
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
Definition: BlackoilPhases.hpp:46