22 #ifndef OPM_MULTISEGMENTWELL_HEADER_INCLUDED 23 #define OPM_MULTISEGMENTWELL_HEADER_INCLUDED 26 #include <opm/autodiff/WellInterface.hpp> 31 template<
typename TypeTag>
38 using typename Base::Simulator;
39 using typename Base::IntensiveQuantities;
40 using typename Base::FluidSystem;
42 using typename Base::MaterialLaw;
43 using typename Base::BlackoilIndices;
48 using Base::has_solvent;
49 using Base::has_polymer;
59 static const bool gasoil =
numEq == 2 && (BlackoilIndices::compositionSwitchIdx >= 0);
60 static const int GTotal = 0;
61 static const int WFrac = gasoil? -1000: 1;
62 static const int GFrac = gasoil? 1 : 2;
63 static const int SPres = gasoil? 2 : 3;
68 using typename Base::Scalar;
72 using typename Base::Mat;
73 using typename Base::BVector;
74 using typename Base::Eval;
81 typedef Dune::FieldVector<Scalar, numWellEq> VectorBlockWellType;
82 typedef Dune::BlockVector<VectorBlockWellType> BVectorWell;
85 typedef Dune::FieldMatrix<Scalar, numWellEq, numWellEq > DiagMatrixBlockWellType;
86 typedef Dune::BCRSMatrix <DiagMatrixBlockWellType> DiagMatWell;
89 typedef Dune::FieldMatrix<Scalar, numWellEq, numEq> OffDiagMatrixBlockWellType;
90 typedef Dune::BCRSMatrix<OffDiagMatrixBlockWellType> OffDiagMatWell;
99 virtual void init(
const PhaseUsage* phase_usage_arg,
100 const std::vector<bool>* active_arg,
101 const std::vector<double>& depth_arg,
102 const double gravity_arg,
103 const int num_cells);
106 virtual void initPrimaryVariablesEvaluation()
const;
108 virtual void assembleWellEq(Simulator& ebosSimulator,
122 virtual void apply(
const BVector& x, BVector& Ax)
const;
124 virtual void apply(BVector& r)
const;
134 std::vector<double>& well_potentials);
136 virtual void updatePrimaryVariables(
const WellState& well_state)
const;
138 virtual void solveEqAndUpdateWellState(
WellState& well_state);
140 virtual void calculateExplicitQuantities(
const Simulator& ebosSimulator,
147 int numberOfPerforations()
const;
150 int number_segments_;
153 WellSegment::CompPressureDropEnum compPressureDrop()
const;
155 WellSegment::MultiPhaseModelEnum multiphaseModel()
const;
158 const SegmentSet& segmentSet()
const;
161 using Base::well_ecl_;
162 using Base::number_of_perforations_;
163 using Base::current_step_;
164 using Base::index_of_well_;
165 using Base::number_of_phases_;
169 using Base::well_cells_;
171 using Base::well_index_;
172 using Base::well_type_;
173 using Base::first_perf_;
174 using Base::saturation_table_number_;
175 using Base::well_efficiency_factor_;
176 using Base::gravity_;
177 using Base::well_controls_;
178 using Base::perf_depth_;
182 using Base::phaseUsage;
184 using Base::numComponents;
185 using Base::flowPhaseToEbosPhaseIdx;
186 using Base::flowPhaseToEbosCompIdx;
187 using Base::getAllowCrossFlow;
198 std::vector<std::vector<int> > segment_perforations_;
201 std::vector<std::vector<int> > segment_inlets_;
205 int segmentNumberToIndex(
const int segment_number)
const;
209 mutable OffDiagMatWell duneB_;
210 mutable OffDiagMatWell duneC_;
212 mutable DiagMatWell duneD_;
215 mutable BVectorWell resWell_;
219 mutable std::vector<std::array<double, numWellEq> > primary_variables_;
222 mutable std::vector<std::array<EvalWell, numWellEq> > primary_variables_evaluation_;
225 std::vector<double> cell_perforation_depth_diffs_;
228 std::vector<double> cell_perforation_pressure_diffs_;
233 std::vector<double> perforation_segment_depth_diffs_;
236 std::vector<std::vector<double> > segment_comp_initial_;
240 std::vector<EvalWell> segment_densities_;
243 std::vector<EvalWell> segment_viscosities_;
246 std::vector<EvalWell> segment_mass_rates_;
248 std::vector<double> segment_depth_diffs_;
250 void initMatrixAndVectors(
const int num_cells)
const;
258 void recoverSolutionWell(
const BVector& x, BVectorWell& xw)
const;
261 void updateWellState(
const BVectorWell& dwells,
262 const bool inner_iteration,
268 void initSegmentRatesWithWellRates(
WellState& well_state)
const;
271 void computeInitialComposition();
274 void computePerfCellPressDiffs(
const Simulator& ebosSimulator);
278 EvalWell volumeFraction(
const int seg,
const int comp_idx)
const;
281 EvalWell volumeFractionScaled(
const int seg,
const int comp_idx)
const;
284 EvalWell surfaceVolumeFraction(
const int seg,
const int comp_idx)
const;
286 void computePerfRate(
const IntensiveQuantities& int_quants,
287 const std::vector<EvalWell>& mob_perfcells,
290 const EvalWell& segment_pressure,
291 const bool& allow_cf,
292 std::vector<EvalWell>& cq_s)
const;
295 EvalWell extendEval(
const Eval& in)
const;
299 void computeSegmentFluidProperties(
const Simulator& ebosSimulator);
301 EvalWell getSegmentPressure(
const int seg)
const;
303 EvalWell getSegmentRate(
const int seg,
const int comp_idx)
const;
305 EvalWell getSegmentGTotal(
const int seg)
const;
308 void getMobility(
const Simulator& ebosSimulator,
310 std::vector<EvalWell>& mob)
const;
312 void assembleControlEq()
const;
314 void assemblePressureEq(
const int seg)
const;
317 EvalWell getHydroPressureLoss(
const int seg)
const;
320 EvalWell getFrictionPressureLoss(
const int seg)
const;
322 void handleAccelerationPressureLoss(
const int seg)
const;
325 void processFractions(
const int seg)
const;
327 void updateWellStateFromPrimaryVariables(
WellState& well_state)
const;
329 double scalingFactor(
const int comp_idx)
const;
331 bool frictionalPressureLossConsidered()
const;
333 bool accelerationalPressureLossConsidered()
const;
336 void iterateWellEquations(Simulator& ebosSimulator,
340 void assembleWellEqWithoutIteration(Simulator& ebosSimulator,
348 #include "MultisegmentWell_impl.hpp" 350 #endif // OPM_MULTISEGMENTWELL_HEADER_INCLUDED virtual void recoverWellSolutionAndUpdateWellState(const BVector &x, WellState &well_state) const
using the solution x to recover the solution xw for wells and applying xw to update Well State ...
Definition: MultisegmentWell_impl.hpp:533
static const int numEq
the number of reservior equations
Definition: WellInterface.hpp:74
virtual ConvergenceReport getWellConvergence(const std::vector< double > &B_avg) const
check whether the well equations get converged for this well
Definition: MultisegmentWell_impl.hpp:409
a struct to collect information about the convergence checking
Definition: WellInterface.hpp:117
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: AdditionalObjectDeleter.hpp:22
int numberOfSegments() const
number of segments for this well int number_of_segments_;
Definition: MultisegmentWell_impl.hpp:824
Solver parameters for the BlackoilModel.
Definition: BlackoilModelParameters.hpp:31
Definition: WellInterface.hpp:60
const std::string & name() const
Well name.
Definition: WellInterface_impl.hpp:138
static const int numWellEq
the number of well equations // TODO: it should have a more general strategy for it ...
Definition: MultisegmentWell.hpp:66
virtual void apply(const BVector &x, BVector &Ax) const
Ax = Ax - C D^-1 B x.
Definition: MultisegmentWell_impl.hpp:490
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition: WellStateFullyImplicitBlackoil.hpp:44
virtual void computeWellPotentials(const Simulator &ebosSimulator, const WellState &well_state, std::vector< double > &well_potentials)
computing the well potentials for group control
Definition: MultisegmentWell_impl.hpp:548
virtual void updateWellStateWithTarget(const int current, WellState &well_state) const
updating the well state based the control mode specified with current
Definition: MultisegmentWell_impl.hpp:251
Definition: MultisegmentWell.hpp:32