23 #ifndef OPM_STANDARDWELL_EVAL_HEADER_INCLUDED
24 #define OPM_STANDARDWELL_EVAL_HEADER_INCLUDED
26 #include <opm/simulators/wells/StandardWellGeneric.hpp>
28 #include <opm/material/densead/DynamicEvaluation.hpp>
36 class ConvergenceReport;
41 class WellContributions;
42 template<
class Flu
idSystem,
class Indices,
class Scalar>
class WellInterfaceIndices;
45 template<
class Flu
idSystem,
class Indices,
class Scalar>
50 static constexpr
int numWellConservationEq = Indices::numPhases + Indices::numSolvents;
52 static constexpr
int numWellControlEq = 1;
55 static constexpr
int numStaticWellEq = numWellConservationEq + numWellControlEq;
59 static constexpr
int Bhp = numStaticWellEq - numWellControlEq;
79 static const int WQTotal = 0;
84 static const bool waterEnabled = Indices::waterEnabled;
85 static const bool gasEnabled = Indices::gasEnabled;
86 static const bool oilEnabled = Indices::oilEnabled;
88 static constexpr
bool has_wfrac_variable = Indices::waterEnabled && Indices::oilEnabled;
89 static constexpr
bool has_gfrac_variable = Indices::gasEnabled && Indices::numPhases > 1;
90 static constexpr
int WFrac = has_wfrac_variable ? 1 : -1000;
91 static constexpr
int GFrac = has_gfrac_variable ? has_wfrac_variable + 1 : -1000;
92 static constexpr
int SFrac = !Indices::enableSolvent ? -1000 : 3;
95 using EvalWell = DenseAd::DynamicEvaluation<Scalar, numStaticWellEq + Indices::numEq + 1>;
96 using Eval = DenseAd::Evaluation<Scalar, Indices::numEq>;
97 using BVectorWell =
typename StandardWellGeneric<Scalar>::BVectorWell;
99 #if HAVE_CUDA || HAVE_OPENCL
109 void initPrimaryVariablesEvaluation()
const;
111 const EvalWell& getBhp()
const
113 return primary_variables_evaluation_[Bhp];
116 const EvalWell& getWQTotal()
const
118 return primary_variables_evaluation_[WQTotal];
121 EvalWell extendEval(
const Eval& in)
const;
122 EvalWell getQs(
const int compIdx)
const;
123 EvalWell wellSurfaceVolumeFraction(
const int compIdx)
const;
124 EvalWell wellVolumeFraction(
const unsigned compIdx)
const;
125 EvalWell wellVolumeFractionScaled(
const int phase)
const;
129 static double relaxationFactorFractionsProducer(
const std::vector<double>& primary_variables,
130 const BVectorWell& dwells);
132 void assembleControlEq(
const WellState& well_state,
134 const Schedule& schedule,
135 const SummaryState& summaryState,
139 void computeAccumWell();
143 void computeConnectionDensities(
const std::vector<double>& perfComponentRates,
144 const std::vector<double>& b_perf,
145 const std::vector<double>& rsmax_perf,
146 const std::vector<double>& rvmax_perf,
147 const std::vector<double>& surf_dens_perf,
151 const std::vector<double>& B_avg,
152 const double maxResidualAllowed,
153 const double tol_wells,
154 const double relaxed_tolerance_flow,
155 const bool relax_tolerance,
156 std::vector<double>& res,
159 void init(std::vector<double>& perf_depth,
160 const std::vector<double>& depth_arg,
162 const bool has_polymermw);
165 void processFractions()
const;
167 void updatePrimaryVariables(
const WellState& well_state,
170 void updatePrimaryVariablesPolyMW(
const BVectorWell& dwells)
const;
172 void updateWellStateFromPrimaryVariables(
WellState& well_state,
175 void updatePrimaryVariablesNewton(
const BVectorWell& dwells,
176 const double dFLimit,
177 const double dBHPLimit)
const;
179 void updateWellStateFromPrimaryVariablesPolyMW(
WellState& well_state)
const;
186 int numWellEq_ = numStaticWellEq;
190 mutable std::vector<double> primary_variables_;
193 mutable std::vector<EvalWell> primary_variables_evaluation_;
196 std::vector<double> F0_;
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: StandardWellEval.hpp:47
Definition: StandardWellGeneric.hpp:49
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