00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022 #ifndef OPM_STANDARDWELLS_HEADER_INCLUDED
00023 #define OPM_STANDARDWELLS_HEADER_INCLUDED
00024
00025 #include <dune/common/parallel/mpihelper.hh>
00026
00027 #include <opm/common/OpmLog/OpmLog.hpp>
00028
00029 #include <opm/common/utility/platform_dependent/disable_warnings.h>
00030 #include <Eigen/Eigen>
00031 #include <Eigen/Sparse>
00032 #include <opm/common/utility/platform_dependent/reenable_warnings.h>
00033
00034 #include <cassert>
00035 #include <tuple>
00036
00037 #include <opm/parser/eclipse/EclipseState/Schedule/Schedule.hpp>
00038
00039 #include <opm/core/wells.h>
00040 #include <opm/core/wells/DynamicListEconLimited.hpp>
00041 #include <opm/core/wells/WellCollection.hpp>
00042 #include <opm/autodiff/AutoDiffBlock.hpp>
00043 #include <opm/autodiff/AutoDiffHelpers.hpp>
00044 #include <opm/autodiff/BlackoilPropsAdFromDeck.hpp>
00045 #include <opm/simulators/WellSwitchingLogger.hpp>
00046
00047 namespace Opm {
00048
00049
00051 class StandardWells {
00052 public:
00053 struct WellOps {
00054 explicit WellOps(const Wells* wells);
00055 Eigen::SparseMatrix<double> w2p;
00056 Eigen::SparseMatrix<double> p2w;
00057 std::vector<int> well_cells;
00058 };
00059
00060
00061 using ADB = AutoDiffBlock<double>;
00062 using Vector = ADB::V;
00063 using Communication =
00064 Dune::CollectiveCommunication<typename Dune::MPIHelper
00065 ::MPICommunicator>;
00066
00067
00068
00069 using DataBlock = Eigen::Array<double,
00070 Eigen::Dynamic,
00071 Eigen::Dynamic,
00072 Eigen::RowMajor>;
00073
00074 StandardWells(const Wells* wells_arg, WellCollection* well_collection);
00075
00076 void init(const BlackoilPropsAdFromDeck* fluid_arg,
00077 const std::vector<bool>* active_arg,
00078 const std::vector<PhasePresence>* pc_arg,
00079 const VFPProperties* vfp_properties_arg,
00080 const double gravity_arg,
00081 const Vector& depth_arg);
00082
00083 const WellOps& wellOps() const;
00084
00085 int numPhases() const { return wells().number_of_phases; };
00086
00087 const Wells& wells() const;
00088
00089 const Wells* wellsPointer() const;
00090
00092 bool wellsActive() const;
00093 void setWellsActive(const bool wells_active);
00095 bool localWellsActive() const;
00096
00097 int numWellVars() const;
00098
00100 Vector& wellPerforationDensities();
00101 const Vector& wellPerforationDensities() const;
00102
00104 Vector& wellPerforationPressureDiffs();
00105 const Vector& wellPerforationPressureDiffs() const;
00106
00107 template <class ReservoirResidualQuant, class SolutionState>
00108 void extractWellPerfProperties(const SolutionState& state,
00109 const std::vector<ReservoirResidualQuant>& rq,
00110 std::vector<ADB>& mob_perfcells,
00111 std::vector<ADB>& b_perfcells) const;
00112
00113 template <class SolutionState>
00114 void computeWellFlux(const SolutionState& state,
00115 const std::vector<ADB>& mob_perfcells,
00116 const std::vector<ADB>& b_perfcells,
00117 Vector& aliveWells,
00118 std::vector<ADB>& cq_s) const;
00119
00120 template <class SolutionState, class WellState>
00121 void updatePerfPhaseRatesAndPressures(const std::vector<ADB>& cq_s,
00122 const SolutionState& state,
00123 WellState& xw) const;
00124
00125 template <class WellState>
00126 void updateWellState(const Vector& dwells,
00127 const double dpmaxrel,
00128 WellState& well_state);
00129
00130 template <class WellState>
00131 void updateWellControls(WellState& xw) const;
00132
00133
00134 template <class SolutionState>
00135 void addWellFluxEq(const std::vector<ADB>& cq_s,
00136 const SolutionState& state,
00137 LinearisedBlackoilResidual& residual);
00138
00139
00140 template <class SolutionState, class WellState>
00141 void addWellControlEq(const SolutionState& state,
00142 const WellState& xw,
00143 const Vector& aliveWells,
00144 LinearisedBlackoilResidual& residual);
00145
00146 template <class SolutionState, class WellState>
00147 void computeWellConnectionPressures(const SolutionState& state,
00148 const WellState& xw);
00149
00150
00151 template <class SolutionState, class WellState>
00152 void
00153 computeWellPotentials(const std::vector<ADB>& mob_perfcells,
00154 const std::vector<ADB>& b_perfcells,
00155 const WellState& well_state,
00156 SolutionState& state0,
00157 std::vector<double>& well_potentials) const;
00158
00159 template <class SolutionState>
00160 void
00161 variableStateExtractWellsVars(const std::vector<int>& indices,
00162 std::vector<ADB>& vars,
00163 SolutionState& state) const;
00164
00165 void
00166 variableStateWellIndices(std::vector<int>& indices,
00167 int& next) const;
00168
00169 std::vector<int>
00170 variableWellStateIndices() const;
00171
00172 template <class WellState>
00173 void
00174 variableWellStateInitials(const WellState& xw,
00175 std::vector<Vector>& vars0) const;
00176
00179 void setStoreWellPerforationFluxesFlag(const bool store_fluxes);
00180
00184 const Vector& getStoredWellPerforationFluxes() const;
00185
00187 template<class WellState>
00188 void
00189 updateListEconLimited(const Schedule& schedule,
00190 const int current_step,
00191 const Wells* wells,
00192 const WellState& well_state,
00193 DynamicListEconLimited& list_econ_limited) const;
00194
00195
00196 WellCollection* wellCollection() const;
00197
00198 void calculateEfficiencyFactors();
00199
00200 const Vector& wellPerfEfficiencyFactors() const;
00201
00202 protected:
00203 bool wells_active_;
00204 const Wells* wells_;
00205 const WellOps wops_;
00206
00207 WellCollection* well_collection_;
00208
00209
00210
00211
00212 Vector well_perforation_efficiency_factors_;
00213
00214 const BlackoilPropsAdFromDeck* fluid_;
00215 const std::vector<bool>* active_;
00216 const std::vector<PhasePresence>* phase_condition_;
00217 const VFPProperties* vfp_properties_;
00218 double gravity_;
00219
00220
00221 Vector perf_cell_depth_;
00222
00223 Vector well_perforation_densities_;
00224 Vector well_perforation_pressure_diffs_;
00225
00226 bool store_well_perforation_fluxes_;
00227 Vector well_perforation_fluxes_;
00228
00229
00230 template <class SolutionState, class WellState>
00231 void computePropertiesForWellConnectionPressures(const SolutionState& state,
00232 const WellState& xw,
00233 std::vector<double>& b_perf,
00234 std::vector<double>& rsmax_perf,
00235 std::vector<double>& rvmax_perf,
00236 std::vector<double>& surf_dens_perf);
00237
00238 template <class WellState>
00239 void computeWellConnectionDensitesPressures(const WellState& xw,
00240 const std::vector<double>& b_perf,
00241 const std::vector<double>& rsmax_perf,
00242 const std::vector<double>& rvmax_perf,
00243 const std::vector<double>& surf_dens_perf,
00244 const std::vector<double>& depth_perf,
00245 const double grav);
00246
00247
00248 template <class WellState>
00249 bool checkRateEconLimits(const WellEconProductionLimits& econ_production_limits,
00250 const WellState& well_state,
00251 const int well_number) const;
00252
00253 using WellMapType = typename WellState::WellMapType;
00254 using WellMapEntryType = typename WellState::mapentry_t;
00255
00256
00257
00258
00259
00260
00261
00262
00263 using RatioCheckTuple = std::tuple<bool, bool, int, double>;
00264
00265 enum ConnectionIndex {
00266 INVALIDCONNECTION = -10000
00267 };
00268
00269
00270 template <class WellState>
00271 RatioCheckTuple checkRatioEconLimits(const WellEconProductionLimits& econ_production_limits,
00272 const WellState& well_state,
00273 const WellMapEntryType& map_entry) const;
00274
00275 template <class WellState>
00276 RatioCheckTuple checkMaxWaterCutLimit(const WellEconProductionLimits& econ_production_limits,
00277 const WellState& well_state,
00278 const WellMapEntryType& map_entry) const;
00279
00280 template <class WellState>
00281 void updateWellStateWithTarget(const WellControls* wc,
00282 const int current,
00283 const int well_index,
00284 WellState& xw) const;
00285
00286 };
00287
00288
00289 }
00290
00291 #include "StandardWells_impl.hpp"
00292
00293 #endif