22 #ifndef OPM_WELLINTERFACE_HEADER_INCLUDED 23 #define OPM_WELLINTERFACE_HEADER_INCLUDED 25 #include <opm/common/OpmLog/OpmLog.hpp> 28 #include <opm/parser/eclipse/EclipseState/Schedule/Well.hpp> 29 #include <opm/core/wells.h> 30 #include <opm/core/well_controls.h> 31 #include <opm/core/props/BlackoilPhases.hpp> 32 #include <opm/core/wells/WellsManager.hpp> 34 #include <opm/autodiff/VFPProperties.hpp> 35 #include <opm/autodiff/VFPInjProperties.hpp> 36 #include <opm/autodiff/VFPProdProperties.hpp> 37 #include <opm/autodiff/WellHelpers.hpp> 38 #include <opm/autodiff/WellStateFullyImplicitBlackoil.hpp> 39 #include <opm/autodiff/BlackoilModelParameters.hpp> 41 #include <opm/simulators/WellSwitchingLogger.hpp> 43 #include<dune/common/fmatrix.hh> 44 #include<dune/istl/bcrsmatrix.hh> 45 #include<dune/istl/matrixmatrix.hh> 47 #include <opm/material/densead/Math.hpp> 48 #include <opm/material/densead/Evaluation.hpp> 59 template<
typename TypeTag>
67 typedef typename GET_PROP_TYPE(TypeTag, Grid) Grid;
68 typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
69 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
70 typedef typename GET_PROP_TYPE(TypeTag, Indices) BlackoilIndices;
71 typedef typename GET_PROP_TYPE(TypeTag, IntensiveQuantities) IntensiveQuantities;
72 typedef typename GET_PROP_TYPE(TypeTag, MaterialLaw) MaterialLaw;
74 static const int numEq = BlackoilIndices::numEq;
75 typedef double Scalar;
77 typedef Dune::FieldVector<Scalar, numEq > VectorBlockType;
78 typedef Dune::FieldMatrix<Scalar, numEq, numEq > MatrixBlockType;
79 typedef Dune::BCRSMatrix <MatrixBlockType> Mat;
80 typedef Dune::BlockVector<VectorBlockType> BVector;
81 typedef DenseAd::Evaluation<double, numEq> Eval;
83 typedef Ewoms::BlackOilPolymerModule<TypeTag> PolymerModule;
85 static const bool has_solvent = GET_PROP_VALUE(TypeTag, EnableSolvent);
86 static const bool has_polymer = GET_PROP_VALUE(TypeTag, EnablePolymer);
95 const std::string&
name()
const;
98 const std::vector<int>&
cells() {
return well_cells_; }
106 void setVFPProperties(
const VFPProperties* vfp_properties_arg);
108 virtual void init(
const PhaseUsage* phase_usage_arg,
109 const std::vector<bool>* active_arg,
110 const std::vector<double>& depth_arg,
111 const double gravity_arg,
112 const int num_cells);
114 virtual void initPrimaryVariablesEvaluation()
const = 0;
119 std::string well_name;
120 std::string phase_name;
122 bool converged =
true;
123 bool nan_residual_found =
false;
124 std::vector<ProblemWell> nan_residual_wells;
126 bool too_large_residual_found =
false;
127 std::vector<ProblemWell> too_large_residual_wells;
130 converged = converged && rhs.converged;
131 nan_residual_found = nan_residual_found || rhs.nan_residual_found;
132 if (rhs.nan_residual_found) {
133 for (
const ProblemWell& well : rhs.nan_residual_wells) {
134 nan_residual_wells.push_back(well);
137 too_large_residual_found = too_large_residual_found || rhs.too_large_residual_found;
138 if (rhs.too_large_residual_found) {
139 for (
const ProblemWell& well : rhs.too_large_residual_wells) {
140 too_large_residual_wells.push_back(well);
147 virtual ConvergenceReport getWellConvergence(
const std::vector<double>& B_avg)
const = 0;
149 virtual void solveEqAndUpdateWellState(WellState& well_state) = 0;
151 virtual void assembleWellEq(Simulator& ebosSimulator,
153 WellState& well_state,
154 bool only_wells) = 0;
156 void updateListEconLimited(
const WellState& well_state,
157 DynamicListEconLimited& list_econ_limited)
const;
159 void setWellEfficiencyFactor(
const double efficiency_factor);
161 void computeRepRadiusPerfLength(
const Grid& grid,
const std::map<int, int>& cartesian_to_compressed);
166 WellState& well_state)
const = 0;
169 virtual void apply(
const BVector& x, BVector& Ax)
const = 0;
172 virtual void apply(BVector& r)
const = 0;
175 virtual void computeWellPotentials(
const Simulator& ebosSimulator,
176 const WellState& well_state,
177 std::vector<double>& well_potentials) = 0;
179 virtual void updateWellStateWithTarget(
const int current,
180 WellState& xw)
const = 0;
182 void updateWellControl(WellState& xw,
183 wellhelpers::WellSwitchingLogger& logger)
const;
185 virtual void updatePrimaryVariables(
const WellState& well_state)
const = 0;
187 virtual void calculateExplicitQuantities(
const Simulator& ebosSimulator,
188 const WellState& xw) = 0;
193 static const int INVALIDCONNECTION = -100000;
195 const Well* well_ecl_;
197 const int current_step_;
203 const ModelParameters& param_;
207 enum WellType well_type_;
210 int number_of_phases_;
214 std::vector<double> comp_frac_;
217 struct WellControls* well_controls_;
220 int number_of_perforations_;
227 std::vector<double> well_index_;
230 std::vector<double> perf_depth_;
235 double well_efficiency_factor_;
238 std::vector<int> well_cells_;
241 std::vector<int> saturation_table_number_;
244 std::vector<double> perf_rep_radius_;
247 std::vector<double> perf_length_;
250 std::vector<double> bore_diameters_;
252 const PhaseUsage* phase_usage_;
254 bool getAllowCrossFlow()
const;
256 const std::vector<bool>* active_;
258 const VFPProperties* vfp_properties_;
262 const std::vector<bool>& active()
const;
264 const PhaseUsage& phaseUsage()
const;
266 int flowPhaseToEbosCompIdx(
const int phaseIdx )
const;
268 int flowPhaseToEbosPhaseIdx(
const int phaseIdx )
const;
271 int numComponents()
const;
273 double wsolvent()
const;
275 double wpolymer()
const;
277 bool checkRateEconLimits(
const WellEconProductionLimits& econ_production_limits,
278 const WellState& well_state)
const;
280 bool wellHasTHPConstraints()
const;
283 const std::vector<double>& compFrac()
const;
285 double mostStrictBhpFromBhpLimits()
const;
294 using RatioCheckTuple = std::tuple<bool, bool, int, double>;
296 RatioCheckTuple checkMaxWaterCutLimit(
const WellEconProductionLimits& econ_production_limits,
297 const WellState& well_state)
const;
299 RatioCheckTuple checkRatioEconLimits(
const WellEconProductionLimits& econ_production_limits,
300 const WellState& well_state)
const;
306 #include "WellInterface_impl.hpp" 308 #endif // OPM_WELLINTERFACE_HEADER_INCLUDED A thin wrapper class that holds one VFPProdProperties and one VFPInjProperties object.
Definition: VFPProperties.hpp:37
virtual void recoverWellSolutionAndUpdateWellState(const BVector &x, WellState &well_state) const =0
using the solution x to recover the solution xw for wells and applying xw to update Well State ...
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
Definition: WellInterface.hpp:60
const std::string & name() const
Well name.
Definition: WellInterface_impl.hpp:138
virtual ~WellInterface()
Virutal destructor.
Definition: WellInterface.hpp:92
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition: WellStateFullyImplicitBlackoil.hpp:44
WellType wellType() const
Well type, INJECTOR or PRODUCER.
Definition: WellInterface_impl.hpp:150
virtual void apply(const BVector &x, BVector &Ax) const =0
Ax = Ax - C D^-1 B x.
WellInterface(const Well *well, const int time_step, const Wells *wells, const ModelParameters ¶m)
Constructor.
Definition: WellInterface_impl.hpp:28
const std::vector< int > & cells()
Well cells.
Definition: WellInterface.hpp:98
Definition: WellInterface.hpp:118
WellControls * wellControls() const
Well controls.
Definition: WellInterface_impl.hpp:162