23 #ifndef OPM_STANDARDWELL_HEADER_INCLUDED 24 #define OPM_STANDARDWELL_HEADER_INCLUDED 27 #include <opm/autodiff/WellInterface.hpp> 32 template<
typename TypeTag>
41 using typename Base::Simulator;
43 using typename Base::IntensiveQuantities;
44 using typename Base::FluidSystem;
45 using typename Base::MaterialLaw;
47 using typename Base::BlackoilIndices;
48 using typename Base::PolymerModule;
56 static const bool gasoil = numEq == 2 && (BlackoilIndices::compositionSwitchIdx >= 0);
57 static const int XvarWell = 0;
58 static const int WFrac = gasoil? -1000: 1;
59 static const int GFrac = gasoil? 1: 2;
60 static const int SFrac = 3;
62 using typename Base::Scalar;
66 using Base::has_solvent;
67 using Base::has_polymer;
72 static const int numWellEq = GET_PROP_VALUE(TypeTag, EnablePolymer)? numEq-1 : numEq;
73 using typename Base::Mat;
74 using typename Base::BVector;
75 using typename Base::Eval;
82 typedef Dune::FieldVector<Scalar, numWellEq> VectorBlockWellType;
83 typedef Dune::BlockVector<VectorBlockWellType> BVectorWell;
86 typedef Dune::FieldMatrix<Scalar, numWellEq, numWellEq > DiagMatrixBlockWellType;
87 typedef Dune::BCRSMatrix <DiagMatrixBlockWellType> DiagMatWell;
90 typedef Dune::FieldMatrix<Scalar, numWellEq, numEq> OffDiagMatrixBlockWellType;
91 typedef Dune::BCRSMatrix<OffDiagMatrixBlockWellType> OffDiagMatWell;
93 typedef DenseAd::Evaluation<double, numEq + numWellEq> EvalWell;
96 static const int contiSolventEqIdx = BlackoilIndices::contiSolventEqIdx;
97 static const int contiPolymerEqIdx = BlackoilIndices::contiPolymerEqIdx;
98 static const int solventSaturationIdx = BlackoilIndices::solventSaturationIdx;
99 static const int polymerConcentrationIdx = BlackoilIndices::polymerConcentrationIdx;
104 virtual void init(
const PhaseUsage* phase_usage_arg,
105 const std::vector<bool>* active_arg,
106 const std::vector<double>& depth_arg,
107 const double gravity_arg,
108 const int num_cells);
111 virtual void initPrimaryVariablesEvaluation()
const;
113 virtual void assembleWellEq(Simulator& ebosSimulator,
127 virtual void apply(
const BVector& x, BVector& Ax)
const;
129 virtual void apply(BVector& r)
const;
139 std::vector<double>& well_potentials) ;
141 virtual void updatePrimaryVariables(
const WellState& well_state)
const;
143 virtual void solveEqAndUpdateWellState(
WellState& well_state);
145 virtual void calculateExplicitQuantities(
const Simulator& ebosSimulator,
150 using Base::getAllowCrossFlow;
151 using Base::phaseUsage;
153 using Base::flowPhaseToEbosPhaseIdx;
154 using Base::flowPhaseToEbosCompIdx;
155 using Base::numComponents;
156 using Base::wsolvent;
157 using Base::wpolymer;
158 using Base::wellHasTHPConstraints;
159 using Base::mostStrictBhpFromBhpLimits;
162 using Base::vfp_properties_;
163 using Base::gravity_;
165 using Base::well_efficiency_factor_;
166 using Base::first_perf_;
167 using Base::ref_depth_;
168 using Base::perf_depth_;
169 using Base::well_cells_;
170 using Base::number_of_perforations_;
171 using Base::number_of_phases_;
172 using Base::saturation_table_number_;
173 using Base::comp_frac_;
174 using Base::well_index_;
175 using Base::index_of_well_;
176 using Base::well_controls_;
177 using Base::well_type_;
179 using Base::perf_rep_radius_;
180 using Base::perf_length_;
181 using Base::bore_diameters_;
184 std::vector<double> perf_densities_;
186 std::vector<double> perf_pressure_diffs_;
189 BVectorWell resWell_;
192 OffDiagMatWell duneB_;
193 OffDiagMatWell duneC_;
195 DiagMatWell invDuneD_;
198 mutable BVectorWell Bx_;
199 mutable BVectorWell invDrw_;
203 mutable std::vector<double> primary_variables_;
206 mutable std::vector<EvalWell> primary_variables_evaluation_;
209 std::vector<double> F0_;
214 EvalWell getBhp()
const;
217 EvalWell getQs(
const int comp_idx)
const;
219 EvalWell wellVolumeFractionScaled(
const int phase)
const;
221 EvalWell wellVolumeFraction(
const int phase)
const;
223 EvalWell wellSurfaceVolumeFraction(
const int phase)
const;
225 EvalWell extendEval(
const Eval& in)
const;
227 bool crossFlowAllowed(
const Simulator& ebosSimulator)
const;
230 void recoverSolutionWell(
const BVector& x, BVectorWell& xw)
const;
233 void updateWellState(
const BVectorWell& dwells,
238 void computePropertiesForWellConnectionPressures(
const Simulator& ebosSimulator,
240 std::vector<double>& b_perf,
241 std::vector<double>& rsmax_perf,
242 std::vector<double>& rvmax_perf,
243 std::vector<double>& surf_dens_perf)
const;
247 void computeConnectionDensities(
const std::vector<double>& perfComponentRates,
248 const std::vector<double>& b_perf,
249 const std::vector<double>& rsmax_perf,
250 const std::vector<double>& rvmax_perf,
251 const std::vector<double>& surf_dens_perf);
253 void computeConnectionPressureDelta();
255 void computeWellConnectionDensitesPressures(
const WellState& xw,
256 const std::vector<double>& b_perf,
257 const std::vector<double>& rsmax_perf,
258 const std::vector<double>& rvmax_perf,
259 const std::vector<double>& surf_dens_perf);
262 void computeAccumWell();
264 void computeWellConnectionPressures(
const Simulator& ebosSimulator,
268 void computePerfRate(
const IntensiveQuantities& intQuants,
269 const std::vector<EvalWell>& mob_perfcells_dense,
270 const double Tw,
const EvalWell& bhp,
const double& cdp,
271 const bool& allow_cf, std::vector<EvalWell>& cq_s)
const;
275 void computeWellRatesWithBhp(
const Simulator& ebosSimulator,
277 std::vector<double>& well_flux)
const;
279 std::vector<double> computeWellPotentialWithTHP(
const Simulator& ebosSimulator,
280 const double initial_bhp,
281 const std::vector<double>& initial_potential)
const;
283 template <
class ValueType>
284 ValueType calculateBhpFromThp(
const std::vector<ValueType>& rates,
const int control_index)
const;
286 double calculateThpFromBhp(
const std::vector<double>& rates,
const int control_index,
const double bhp)
const;
289 void getMobility(
const Simulator& ebosSimulator,
291 std::vector<EvalWell>& mob)
const;
293 double scalingFactor(
const int comp_idx)
const;
298 #include "StandardWell_impl.hpp" 300 #endif // OPM_STANDARDWELL_HEADER_INCLUDED virtual void updateWellStateWithTarget(const int current, WellState &xw) const
updating the well state based the control mode specified with current
Definition: StandardWell_impl.hpp:1013
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: StandardWell_impl.hpp:1594
Definition: StandardWell.hpp:33
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
Solver parameters for the BlackoilModel.
Definition: BlackoilModelParameters.hpp:31
virtual ConvergenceReport getWellConvergence(const std::vector< double > &B_avg) const
check whether the well equations get converged for this well
Definition: StandardWell_impl.hpp:1366
Definition: WellInterface.hpp:60
const std::string & name() const
Well name.
Definition: WellInterface_impl.hpp:138
virtual void computeWellPotentials(const Simulator &ebosSimulator, const WellState &well_state, std::vector< double > &well_potentials)
computing the well potentials for group control
Definition: StandardWell_impl.hpp:1744
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition: WellStateFullyImplicitBlackoil.hpp:44
virtual void apply(const BVector &x, BVector &Ax) const
Ax = Ax - C D^-1 B x.
Definition: StandardWell_impl.hpp:1538