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::GTotal;
76 using MSWEval::numWellEq;
82 const RateConverterType& rate_converter,
83 const int pvtRegionIdx,
84 const int num_components,
86 const int index_of_well,
87 const std::vector<PerforationData>& perf_data);
89 virtual void init(
const PhaseUsage* phase_usage_arg,
90 const std::vector<double>& depth_arg,
91 const double gravity_arg,
93 const std::vector< Scalar >& B_avg)
override;
95 virtual void initPrimaryVariablesEvaluation()
const override;
105 const std::vector<double>& B_avg,
107 const bool relax_tolerance =
false)
const override;
110 virtual void apply(
const BVector& x, BVector& Ax)
const override;
112 virtual void apply(BVector& r)
const override;
123 std::vector<double>& well_potentials,
126 virtual void updatePrimaryVariables(
const WellState& well_state,
DeferredLogger& deferred_logger)
const override;
130 virtual void calculateExplicitQuantities(
const Simulator& ebosSimulator,
134 virtual void updateProductivityIndex(
const Simulator& ebosSimulator,
139 virtual void addWellContributions(SparseMatrixAdapter& jacobian)
const override;
144 void computeConnLevelProdInd(
const FluidState& fs,
145 const std::function<
double(
const double)>& connPICalc,
146 const std::vector<Scalar>& mobility,
147 double* connPI)
const;
149 void computeConnLevelInjInd(
const FluidState& fs,
150 const Phase preferred_phase,
151 const std::function<
double(
const double)>& connIICalc,
152 const std::vector<Scalar>& mobility,
156 virtual std::optional<double> computeBhpAtThpLimitProdWithAlq(
157 const Simulator& ebos_simulator,
158 const SummaryState& summary_state,
160 double alq_value)
const override;
163 int number_segments_;
166 WellSegments::CompPressureDrop compPressureDrop()
const;
168 WellSegments::MultiPhaseModel multiphaseModel()
const;
171 std::vector<std::vector<double> > segment_fluid_initial_;
173 mutable int debug_cost_counter_ = 0;
176 void updateWellState(
const BVectorWell& dwells,
179 const double relaxation_factor=1.0)
const;
183 void computeInitialSegmentFluids(
const Simulator& ebos_simulator);
186 void computePerfCellPressDiffs(
const Simulator& ebosSimulator);
188 void computePerfRateScalar(
const IntensiveQuantities& int_quants,
189 const std::vector<Scalar>& mob_perfcells,
193 const Scalar& segment_pressure,
194 const bool& allow_cf,
195 std::vector<Scalar>& cq_s,
198 void computePerfRateEval(
const IntensiveQuantities& int_quants,
199 const std::vector<EvalWell>& mob_perfcells,
203 const EvalWell& segment_pressure,
204 const bool& allow_cf,
205 std::vector<EvalWell>& cq_s,
206 EvalWell& perf_press,
207 double& perf_dis_gas_rate,
208 double& perf_vap_oil_rate,
211 template<
class Value>
212 void computePerfRate(
const Value& pressure_cell,
215 const std::vector<Value>& b_perfcells,
216 const std::vector<Value>& mob_perfcells,
219 const Value& segment_pressure,
220 const Value& segment_density,
221 const bool& allow_cf,
222 const std::vector<Value>& cmix_s,
223 std::vector<Value>& cq_s,
225 double& perf_dis_gas_rate,
226 double& perf_vap_oil_rate,
231 void computeSegmentFluidProperties(
const Simulator& ebosSimulator,
235 void getMobilityEval(
const Simulator& ebosSimulator,
237 std::vector<EvalWell>& mob)
const;
240 void getMobilityScalar(
const Simulator& ebosSimulator,
242 std::vector<Scalar>& mob)
const;
244 void computeWellRatesAtBhpLimit(
const Simulator& ebosSimulator,
245 std::vector<double>& well_flux,
248 virtual void computeWellRatesWithBhp(
const Simulator& ebosSimulator,
250 std::vector<double>& well_flux,
253 void computeWellRatesWithBhpIterations(
const Simulator& ebosSimulator,
255 std::vector<double>& well_flux,
258 std::vector<double> computeWellPotentialWithTHP(
260 const Simulator& ebos_simulator,
263 virtual double getRefDensity()
const override;
265 virtual bool iterateWellEqWithControl(
const Simulator& ebosSimulator,
267 const Well::InjectionControls& inj_controls,
268 const Well::ProductionControls& prod_controls,
273 virtual void assembleWellEqWithoutIteration(
const Simulator& ebosSimulator,
275 const Well::InjectionControls& inj_controls,
276 const Well::ProductionControls& prod_controls,
281 virtual void updateWaterThroughput(
const double dt,
WellState& well_state)
const override;
283 EvalWell getSegmentSurfaceVolume(
const Simulator& ebos_simulator,
const int seg_idx)
const;
290 bool openCrossFlowAvoidSingularity(
const Simulator& ebos_simulator)
const;
294 bool allDrawDownWrongDirection(
const Simulator& ebos_simulator)
const;
298 std::optional<double> computeBhpAtThpLimitProd(
300 const Simulator& ebos_simulator,
301 const SummaryState& summary_state,
304 std::optional<double> computeBhpAtThpLimitInj(
const Simulator& ebos_simulator,
305 const SummaryState& summary_state,
308 double maxPerfPress(
const Simulator& ebos_simulator)
const;
311 virtual void checkOperabilityUnderBHPLimit(
const WellState& well_state,
const Simulator& ebos_simulator,
DeferredLogger& deferred_logger)
override;
314 virtual void checkOperabilityUnderTHPLimit(
const Simulator& ebos_simulator,
const WellState& well_state,
DeferredLogger& deferred_logger)
override;
317 virtual void updateIPR(
const Simulator& ebos_simulator,
DeferredLogger& deferred_logger)
const override;
322 #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: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:1909
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: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