22 #ifndef OPM_STANDARDWELLS_HEADER_INCLUDED
23 #define OPM_STANDARDWELLS_HEADER_INCLUDED
25 #include <dune/common/parallel/mpihelper.hh>
27 #include <opm/common/OpmLog/OpmLog.hpp>
29 #include <opm/common/utility/platform_dependent/disable_warnings.h>
30 #include <Eigen/Eigen>
31 #include <Eigen/Sparse>
32 #include <opm/common/utility/platform_dependent/reenable_warnings.h>
37 #include <opm/parser/eclipse/EclipseState/Schedule/Schedule.hpp>
39 #include <opm/core/wells.h>
40 #include <opm/core/wells/DynamicListEconLimited.hpp>
41 #include <opm/core/wells/WellCollection.hpp>
42 #include <opm/autodiff/AutoDiffBlock.hpp>
43 #include <opm/autodiff/AutoDiffHelpers.hpp>
44 #include <opm/autodiff/BlackoilPropsAdFromDeck.hpp>
45 #include <opm/simulators/WellSwitchingLogger.hpp>
54 explicit WellOps(
const Wells* wells);
55 Eigen::SparseMatrix<double> w2p;
56 Eigen::SparseMatrix<double> p2w;
57 std::vector<int> well_cells;
64 Dune::CollectiveCommunication<
typename Dune::MPIHelper
69 using DataBlock = Eigen::Array<double,
74 StandardWells(
const Wells* wells_arg, WellCollection* well_collection);
77 const std::vector<bool>* active_arg,
78 const std::vector<PhasePresence>* pc_arg,
80 const double gravity_arg,
81 const Vector& depth_arg);
85 int numPhases()
const {
return wells().number_of_phases; };
87 const Wells& wells()
const;
89 const Wells* wellsPointer()
const;
93 void setWellsActive(
const bool wells_active);
97 int numWellVars()
const;
107 template <
class ReservoirRes
idualQuant,
class SolutionState>
108 void extractWellPerfProperties(
const SolutionState& state,
109 const std::vector<ReservoirResidualQuant>& rq,
110 std::vector<ADB>& mob_perfcells,
111 std::vector<ADB>& b_perfcells)
const;
113 template <
class SolutionState>
114 void computeWellFlux(
const SolutionState& state,
115 const std::vector<ADB>& mob_perfcells,
116 const std::vector<ADB>& b_perfcells,
118 std::vector<ADB>& cq_s)
const;
120 template <
class SolutionState,
class WellState>
121 void updatePerfPhaseRatesAndPressures(
const std::vector<ADB>& cq_s,
122 const SolutionState& state,
123 WellState& xw)
const;
125 template <
class WellState>
126 void updateWellState(
const Vector& dwells,
127 const double dpmaxrel,
128 WellState& well_state);
130 template <
class WellState>
131 void updateWellControls(WellState& xw)
const;
134 template <
class SolutionState>
135 void addWellFluxEq(
const std::vector<ADB>& cq_s,
136 const SolutionState& state,
140 template <
class SolutionState,
class WellState>
141 void addWellControlEq(
const SolutionState& state,
143 const Vector& aliveWells,
146 template <
class SolutionState,
class WellState>
147 void computeWellConnectionPressures(
const SolutionState& state,
148 const WellState& xw);
151 template <
class SolutionState,
class WellState>
153 computeWellPotentials(
const std::vector<ADB>& mob_perfcells,
154 const std::vector<ADB>& b_perfcells,
155 const WellState& well_state,
156 SolutionState& state0,
157 std::vector<double>& well_potentials)
const;
159 template <
class SolutionState>
161 variableStateExtractWellsVars(
const std::vector<int>& indices,
162 std::vector<ADB>& vars,
163 SolutionState& state)
const;
166 variableStateWellIndices(std::vector<int>& indices,
170 variableWellStateIndices()
const;
172 template <
class WellState>
174 variableWellStateInitials(
const WellState& xw,
175 std::vector<Vector>& vars0)
const;
187 template<
class WellState>
190 const int current_step,
192 const WellState& well_state,
193 DynamicListEconLimited& list_econ_limited)
const;
196 WellCollection* wellCollection()
const;
198 void calculateEfficiencyFactors();
200 const Vector& wellPerfEfficiencyFactors()
const;
207 WellCollection* well_collection_;
212 Vector well_perforation_efficiency_factors_;
215 const std::vector<bool>* active_;
216 const std::vector<PhasePresence>* phase_condition_;
221 Vector perf_cell_depth_;
223 Vector well_perforation_densities_;
224 Vector well_perforation_pressure_diffs_;
226 bool store_well_perforation_fluxes_;
227 Vector well_perforation_fluxes_;
230 template <
class SolutionState,
class WellState>
231 void computePropertiesForWellConnectionPressures(
const SolutionState& state,
233 std::vector<double>& b_perf,
234 std::vector<double>& rsmax_perf,
235 std::vector<double>& rvmax_perf,
236 std::vector<double>& surf_dens_perf);
238 template <
class WellState>
239 void computeWellConnectionDensitesPressures(
const WellState& xw,
240 const std::vector<double>& b_perf,
241 const std::vector<double>& rsmax_perf,
242 const std::vector<double>& rvmax_perf,
243 const std::vector<double>& surf_dens_perf,
244 const std::vector<double>& depth_perf,
248 template <
class WellState>
249 bool checkRateEconLimits(
const WellEconProductionLimits& econ_production_limits,
250 const WellState& well_state,
251 const int well_number)
const;
253 using WellMapType =
typename WellState::WellMapType;
254 using WellMapEntryType =
typename WellState::mapentry_t;
263 using RatioCheckTuple = std::tuple<bool, bool, int, double>;
265 enum ConnectionIndex {
266 INVALIDCONNECTION = -10000
270 template <
class WellState>
271 RatioCheckTuple checkRatioEconLimits(
const WellEconProductionLimits& econ_production_limits,
272 const WellState& well_state,
273 const WellMapEntryType& map_entry)
const;
275 template <
class WellState>
276 RatioCheckTuple checkMaxWaterCutLimit(
const WellEconProductionLimits& econ_production_limits,
277 const WellState& well_state,
278 const WellMapEntryType& map_entry)
const;
280 template <
class WellState>
281 void updateWellStateWithTarget(
const WellControls* wc,
283 const int well_index,
284 WellState& xw)
const;
291 #include "StandardWells_impl.hpp"
A thin wrapper class that holds one VFPProdProperties and one VFPInjProperties object.
Definition: VFPProperties.hpp:37
Vector & wellPerforationDensities()
Density of each well perforation.
Definition: StandardWells_impl.hpp:184
bool localWellsActive() const
return true if wells are available on this process
Definition: StandardWells_impl.hpp:148
bool wellsActive() const
return true if wells are available in the reservoir
Definition: StandardWells_impl.hpp:130
Vector & wellPerforationPressureDiffs()
Diff to bhp for each well perforation.
Definition: StandardWells_impl.hpp:204
void setStoreWellPerforationFluxesFlag(const bool store_fluxes)
If set, computeWellFlux() will additionally store the total reservoir volume perforation fluxes...
Definition: StandardWells_impl.hpp:1185
Class for handling the standard well model.
Definition: StandardWells.hpp:51
Definition: StandardWells.hpp:53
Residual structure of the fully implicit solver.
Definition: LinearisedBlackoilResidual.hpp:47
This class implements the AD-adapted fluid interface for three-phase black-oil.
Definition: BlackoilPropsAdFromDeck.hpp:61
void updateListEconLimited(const Schedule &schedule, const int current_step, const Wells *wells, const WellState &well_state, DynamicListEconLimited &list_econ_limited) const
upate the dynamic lists related to economic limits
Definition: StandardWells_impl.hpp:1209
Eigen::Array< double, Eigen::Dynamic, 1 > V
Underlying type for values.
Definition: AutoDiffBlock.hpp:99
const Vector & getStoredWellPerforationFluxes() const
Retrieves the stored fluxes.
Definition: StandardWells_impl.hpp:1195