22 #ifndef OPM_MULTISEGMENTWELL_EVAL_HEADER_INCLUDED
23 #define OPM_MULTISEGMENTWELL_EVAL_HEADER_INCLUDED
25 #include <opm/simulators/wells/MultisegmentWellGeneric.hpp>
27 #include <opm/material/densead/Evaluation.hpp>
29 #include <opm/parser/eclipse/EclipseState/Schedule/Well/Well.hpp>
31 #include <dune/common/fmatrix.hh>
32 #include <dune/common/fvector.hh>
33 #include <dune/istl/bcrsmatrix.hh>
34 #include <dune/istl/bvector.hh>
35 #include <dune/istl/umfpack.hh>
43 class ConvergenceReport;
46 class WellContributions;
47 template<
class Flu
idSystem,
class Indices,
class Scalar>
class WellInterfaceIndices;
50 template<
typename Flu
idSystem,
typename Indices,
typename Scalar>
54 #if HAVE_CUDA || HAVE_OPENCL
73 static constexpr
bool has_water = (Indices::waterSaturationIdx >= 0);
74 static constexpr
bool has_gas = (Indices::compositionSwitchIdx >= 0);
75 static constexpr
bool has_oil = (Indices::numPhases - has_gas - has_water) > 0;
80 static constexpr
bool has_wfrac_variable = has_water && Indices::numPhases > 1;
81 static constexpr
bool has_gfrac_variable = has_gas && has_oil;
83 static constexpr
int GTotal = 0;
84 static constexpr
int WFrac = has_wfrac_variable ? 1 : -1000;
85 static constexpr
int GFrac = has_gfrac_variable ? has_wfrac_variable + 1 : -1000;
86 static constexpr
int SPres = has_wfrac_variable + has_gfrac_variable + 1;
89 static constexpr
int numWellEq = Indices::numPhases + 1;
96 using VectorBlockWellType = Dune::FieldVector<Scalar, numWellEq>;
97 using BVectorWell = Dune::BlockVector<VectorBlockWellType>;
99 using VectorBlockType = Dune::FieldVector<Scalar, Indices::numEq>;
100 using BVector = Dune::BlockVector<VectorBlockType>;
103 using DiagMatrixBlockWellType = Dune::FieldMatrix<Scalar, numWellEq, numWellEq>;
104 using DiagMatWell = Dune::BCRSMatrix<DiagMatrixBlockWellType>;
107 using OffDiagMatrixBlockWellType = Dune::FieldMatrix<Scalar, numWellEq, Indices::numEq>;
108 using OffDiagMatWell = Dune::BCRSMatrix<OffDiagMatrixBlockWellType>;
113 using EvalWell = DenseAd::Evaluation<double, Indices::numEq + numWellEq>;
114 using Eval = DenseAd::Evaluation<Scalar, Indices::numEq>;
118 void initMatrixAndVectors(
const int num_cells)
const;
119 void initPrimaryVariablesEvaluation()
const;
121 void assembleControlEq(
const WellState& well_state,
123 const Schedule& schedule,
124 const SummaryState& summaryState,
125 const Well::InjectionControls& inj_controls,
126 const Well::ProductionControls& prod_controls,
130 void assembleDefaultPressureEq(
const int seg,
135 void assembleICDPressureEq(
const int seg,
136 const UnitSystem& unit_system,
141 void assemblePressureEq(
const int seg,
142 const UnitSystem& unit_system,
146 void checkConvergenceControlEq(
const WellState& well_state,
148 const double tolerance_pressure_ms_wells,
149 const double tolerance_wells,
150 const double max_residual_allowed,
155 const std::vector<double>& B_avg,
157 const double max_residual_allowed,
158 const double tolerance_wells,
159 const double relaxed_inner_tolerance_flow_ms_well,
160 const double tolerance_pressure_ms_wells,
161 const double relaxed_inner_tolerance_pressure_ms_well,
162 const bool relax_tolerance)
const;
165 void processFractions(
const int seg)
const;
168 void recoverSolutionWell(
const BVector& x,
169 BVectorWell& xw)
const;
171 void updatePrimaryVariables(
const WellState& well_state)
const;
173 void updateUpwindingSegments();
176 void updateWellState(
const BVectorWell& dwells,
177 const double relaxation_factor,
178 const double DFLimit,
179 const double max_pressure_change)
const;
181 void computeSegmentFluidProperties(
const EvalWell& temperature,
182 const EvalWell& saltConcentration,
183 int pvt_region_index);
185 EvalWell getBhp()
const;
186 EvalWell getFrictionPressureLoss(
const int seg)
const;
187 EvalWell getHydroPressureLoss(
const int seg)
const;
188 EvalWell getQs(
const int comp_idx)
const;
189 EvalWell getSegmentGTotal(
const int seg)
const;
190 EvalWell getSegmentPressure(
const int seg)
const;
191 EvalWell getSegmentRate(
const int seg,
const int comp_idx)
const;
192 EvalWell getSegmentRateUpwinding(
const int seg,
193 const size_t comp_idx)
const;
194 EvalWell getSegmentSurfaceVolume(
const EvalWell& temperature,
195 const EvalWell& saltConcentration,
196 const int pvt_region_index,
197 const int seg_idx)
const;
198 EvalWell getWQTotal()
const;
201 std::vector<Scalar> getWellResiduals(
const std::vector<Scalar>& B_avg,
204 double getControlTolerance(
const WellState& well_state,
205 const double tolerance_wells,
206 const double tolerance_pressure_ms_wells,
209 double getResidualMeasureValue(
const WellState& well_state,
210 const std::vector<double>& residuals,
211 const double tolerance_wells,
212 const double tolerance_pressure_ms_wells,
215 void handleAccelerationPressureLoss(
const int seg,
219 EvalWell pressureDropAutoICD(
const int seg,
220 const UnitSystem& unit_system)
const;
223 EvalWell pressureDropSpiralICD(
const int seg)
const;
226 EvalWell pressureDropValve(
const int seg)
const;
232 void updateWellStateFromPrimaryVariables(
WellState& well_state,
238 EvalWell volumeFraction(
const int seg,
239 const unsigned compIdx)
const;
242 EvalWell volumeFractionScaled(
const int seg,
243 const int comp_idx)
const;
246 EvalWell surfaceVolumeFraction(
const int seg,
247 const int comp_idx)
const;
250 EvalWell extendEval(
const Eval& in)
const;
256 mutable OffDiagMatWell duneB_;
257 mutable OffDiagMatWell duneC_;
259 mutable DiagMatWell duneD_;
267 mutable BVectorWell resWell_;
271 mutable std::vector<std::array<double, numWellEq> > primary_variables_;
274 mutable std::vector<std::array<EvalWell, numWellEq> > primary_variables_evaluation_;
277 std::vector<int> upwinding_segments_;
281 std::vector<EvalWell> segment_densities_;
284 std::vector<EvalWell> segment_mass_rates_;
287 std::vector<EvalWell> segment_viscosities_;
289 std::vector<std::vector<EvalWell>> segment_phase_densities_;
290 std::vector<std::vector<EvalWell>> segment_phase_fractions_;
291 std::vector<std::vector<EvalWell>> segment_phase_viscosities_;
294 std::vector<double> cell_perforation_depth_diffs_;
297 std::vector<double> cell_perforation_pressure_diffs_;
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
ConvergenceReport getWellConvergence(const WellState &well_state, const std::vector< double > &B_avg, DeferredLogger &deferred_logger, const double max_residual_allowed, const double tolerance_wells, const double relaxed_inner_tolerance_flow_ms_well, const double tolerance_pressure_ms_wells, const double relaxed_inner_tolerance_pressure_ms_well, const bool relax_tolerance) const
check whether the well equations get converged for this well
Definition: MultisegmentWellEval.cpp:229
std::shared_ptr< Dune::UMFPack< DiagMatWell > > duneDSolver_
solver for diagonal matrix
Definition: MultisegmentWellEval.hpp:264
Definition: MultisegmentWellGeneric.hpp:41
This class serves to eliminate the need to include the WellContributions into the matrix (with –matri...
Definition: WellContributions.hpp:61
Definition: WellInterfaceIndices.hpp:35
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