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/input/eclipse/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>
40 template<
class Matrix>
class UMFPack;
46 class ConvergenceReport;
49 class WellContributions;
50 template<
class Flu
idSystem,
class Indices,
class Scalar>
class WellInterfaceIndices;
53 template<
typename Flu
idSystem,
typename Indices,
typename Scalar>
74 static constexpr
bool has_water = (Indices::waterSaturationIdx >= 0);
75 static constexpr
bool has_gas = (Indices::compositionSwitchIdx >= 0);
76 static constexpr
bool has_oil = (Indices::numPhases - has_gas - has_water) > 0;
81 static constexpr
bool has_wfrac_variable = has_water && Indices::numPhases > 1;
82 static constexpr
bool has_gfrac_variable = has_gas && has_oil;
84 static constexpr
int WQTotal = 0;
85 static constexpr
int WFrac = has_wfrac_variable ? 1 : -1000;
86 static constexpr
int GFrac = has_gfrac_variable ? has_wfrac_variable + 1 : -1000;
87 static constexpr
int SPres = has_wfrac_variable + has_gfrac_variable + 1;
90 static constexpr
int numWellEq = Indices::numPhases + 1;
97 using VectorBlockWellType = Dune::FieldVector<Scalar, numWellEq>;
98 using BVectorWell = Dune::BlockVector<VectorBlockWellType>;
100 using VectorBlockType = Dune::FieldVector<Scalar, Indices::numEq>;
101 using BVector = Dune::BlockVector<VectorBlockType>;
104 using DiagMatrixBlockWellType = Dune::FieldMatrix<Scalar, numWellEq, numWellEq>;
105 using DiagMatWell = Dune::BCRSMatrix<DiagMatrixBlockWellType>;
108 using OffDiagMatrixBlockWellType = Dune::FieldMatrix<Scalar, numWellEq, Indices::numEq>;
109 using OffDiagMatWell = Dune::BCRSMatrix<OffDiagMatrixBlockWellType>;
114 using EvalWell = DenseAd::Evaluation<double, Indices::numEq + numWellEq>;
115 using Eval = DenseAd::Evaluation<Scalar, Indices::numEq>;
119 void initMatrixAndVectors(
const int num_cells)
const;
120 void initPrimaryVariablesEvaluation()
const;
122 void assembleControlEq(
const WellState& well_state,
124 const Schedule& schedule,
125 const SummaryState& summaryState,
126 const Well::InjectionControls& inj_controls,
127 const Well::ProductionControls& prod_controls,
131 void assembleDefaultPressureEq(
const int seg,
136 void assembleICDPressureEq(
const int seg,
137 const UnitSystem& unit_system,
142 void assemblePressureEq(
const int seg,
143 const UnitSystem& unit_system,
147 void checkConvergenceControlEq(
const WellState& well_state,
149 const double tolerance_pressure_ms_wells,
150 const double tolerance_wells,
151 const double max_residual_allowed,
156 const std::vector<double>& B_avg,
158 const double max_residual_allowed,
159 const double tolerance_wells,
160 const double relaxed_inner_tolerance_flow_ms_well,
161 const double tolerance_pressure_ms_wells,
162 const double relaxed_inner_tolerance_pressure_ms_well,
163 const bool relax_tolerance)
const;
166 void processFractions(
const int seg)
const;
169 void recoverSolutionWell(
const BVector& x,
170 BVectorWell& xw)
const;
172 void updatePrimaryVariables(
const WellState& well_state)
const;
174 void updateUpwindingSegments();
177 void updatePrimaryVariablesNewton(
const BVectorWell& dwells,
178 const double relaxation_factor,
179 const double DFLimit,
180 const double max_pressure_change)
const;
182 void computeSegmentFluidProperties(
const EvalWell& temperature,
183 const EvalWell& saltConcentration,
184 int pvt_region_index,
187 EvalWell getBhp()
const;
188 EvalWell getFrictionPressureLoss(
const int seg)
const;
189 EvalWell getHydroPressureLoss(
const int seg)
const;
190 EvalWell getQs(
const int comp_idx)
const;
191 EvalWell getSegmentWQTotal(
const int seg)
const;
192 EvalWell getSegmentPressure(
const int seg)
const;
193 EvalWell getSegmentRate(
const int seg,
const int comp_idx)
const;
194 EvalWell getSegmentRateUpwinding(
const int seg,
195 const size_t comp_idx)
const;
196 EvalWell getSegmentSurfaceVolume(
const EvalWell& temperature,
197 const EvalWell& saltConcentration,
198 const int pvt_region_index,
199 const int seg_idx)
const;
200 EvalWell getWQTotal()
const;
203 std::pair<bool, std::vector<Scalar> >
204 getFiniteWellResiduals(
const std::vector<Scalar>& B_avg,
207 double getControlTolerance(
const WellState& well_state,
208 const double tolerance_wells,
209 const double tolerance_pressure_ms_wells,
212 double getResidualMeasureValue(
const WellState& well_state,
213 const std::vector<double>& residuals,
214 const double tolerance_wells,
215 const double tolerance_pressure_ms_wells,
218 void handleAccelerationPressureLoss(
const int seg,
222 EvalWell pressureDropAutoICD(
const int seg,
223 const UnitSystem& unit_system)
const;
226 EvalWell pressureDropSpiralICD(
const int seg)
const;
229 EvalWell pressureDropValve(
const int seg)
const;
235 void updateWellStateFromPrimaryVariables(
WellState& well_state,
241 EvalWell volumeFraction(
const int seg,
242 const unsigned compIdx)
const;
245 EvalWell volumeFractionScaled(
const int seg,
246 const int comp_idx)
const;
249 EvalWell surfaceVolumeFraction(
const int seg,
250 const int comp_idx)
const;
253 EvalWell extendEval(
const Eval& in)
const;
259 mutable OffDiagMatWell duneB_;
260 mutable OffDiagMatWell duneC_;
262 mutable DiagMatWell duneD_;
270 mutable BVectorWell resWell_;
274 mutable std::vector<std::array<double, numWellEq> > primary_variables_;
277 mutable std::vector<std::array<EvalWell, numWellEq> > primary_variables_evaluation_;
280 std::vector<int> upwinding_segments_;
284 std::vector<EvalWell> segment_densities_;
287 std::vector<EvalWell> segment_mass_rates_;
290 std::vector<EvalWell> segment_viscosities_;
292 std::vector<std::vector<EvalWell>> segment_phase_densities_;
293 std::vector<std::vector<EvalWell>> segment_phase_fractions_;
294 std::vector<std::vector<EvalWell>> segment_phase_viscosities_;
297 std::vector<double> cell_perforation_depth_diffs_;
300 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:55
void addWellContribution(WellContributions &wellContribs) const
add the contribution (C, D, B matrices) of this Well to the WellContributions object
Definition: MultisegmentWellEval.cpp:1827
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:231
std::shared_ptr< Dune::UMFPack< DiagMatWell > > duneDSolver_
solver for diagonal matrix
Definition: MultisegmentWellEval.hpp:267
Definition: MultisegmentWellGeneric.hpp:42
This class serves to eliminate the need to include the WellContributions into the matrix (with –matri...
Definition: WellContributions.hpp:53
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:27