StandardWells.hpp
1 /*
2  Copyright 2016 SINTEF ICT, Applied Mathematics.
3  Copyright 2016 Statoil ASA.
4 
5  This file is part of the Open Porous Media project (OPM).
6 
7  OPM is free software: you can redistribute it and/or modify
8  it under the terms of the GNU General Public License as published by
9  the Free Software Foundation, either version 3 of the License, or
10  (at your option) any later version.
11 
12  OPM is distributed in the hope that it will be useful,
13  but WITHOUT ANY WARRANTY; without even the implied warranty of
14  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15  GNU General Public License for more details.
16 
17  You should have received a copy of the GNU General Public License
18  along with OPM. If not, see <http://www.gnu.org/licenses/>.
19 */
20 
21 
22 #ifndef OPM_STANDARDWELLS_HEADER_INCLUDED
23 #define OPM_STANDARDWELLS_HEADER_INCLUDED
24 
25 #include <dune/common/parallel/mpihelper.hh>
26 
27 #include <opm/common/OpmLog/OpmLog.hpp>
28 
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>
33 
34 #include <cassert>
35 #include <tuple>
36 
37 #include <opm/parser/eclipse/EclipseState/Schedule/Schedule.hpp>
38 
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>
46 
47 namespace Opm {
48 
49 
51  class StandardWells {
52  public:
53  struct WellOps {
54  explicit WellOps(const Wells* wells);
55  Eigen::SparseMatrix<double> w2p; // well -> perf (scatter)
56  Eigen::SparseMatrix<double> p2w; // perf -> well (gather)
57  std::vector<int> well_cells; // the set of perforated cells
58  };
59 
60  // --------- Types ---------
61  using ADB = AutoDiffBlock<double>;
62  using Vector = ADB::V;
63  using Communication =
64  Dune::CollectiveCommunication<typename Dune::MPIHelper
65  ::MPICommunicator>;
66 
67  // copied from BlackoilModelBase
68  // should put to somewhere better
69  using DataBlock = Eigen::Array<double,
70  Eigen::Dynamic,
71  Eigen::Dynamic,
72  Eigen::RowMajor>;
73  // --------- Public methods ---------
74  StandardWells(const Wells* wells_arg, WellCollection* well_collection);
75 
76  void init(const BlackoilPropsAdFromDeck* fluid_arg,
77  const std::vector<bool>* active_arg,
78  const std::vector<PhasePresence>* pc_arg,
79  const VFPProperties* vfp_properties_arg,
80  const double gravity_arg,
81  const Vector& depth_arg);
82 
83  const WellOps& wellOps() const;
84 
85  int numPhases() const { return wells().number_of_phases; };
86 
87  const Wells& wells() const;
88 
89  const Wells* wellsPointer() const;
90 
92  bool wellsActive() const;
93  void setWellsActive(const bool wells_active);
95  bool localWellsActive() const;
96 
97  int numWellVars() const;
98 
100  Vector& wellPerforationDensities(); // mutable version kept for BlackoilMultisegmentModel
101  const Vector& wellPerforationDensities() const;
102 
104  Vector& wellPerforationPressureDiffs(); // mutable version kept for BlackoilMultisegmentModel
105  const Vector& wellPerforationPressureDiffs() const;
106 
107  template <class ReservoirResidualQuant, 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;
112 
113  template <class SolutionState>
114  void computeWellFlux(const SolutionState& state,
115  const std::vector<ADB>& mob_perfcells,
116  const std::vector<ADB>& b_perfcells,
117  Vector& aliveWells,
118  std::vector<ADB>& cq_s) const;
119 
120  template <class SolutionState, class WellState>
121  void updatePerfPhaseRatesAndPressures(const std::vector<ADB>& cq_s,
122  const SolutionState& state,
123  WellState& xw) const;
124 
125  template <class WellState>
126  void updateWellState(const Vector& dwells,
127  const double dpmaxrel,
128  WellState& well_state);
129 
130  template <class WellState>
131  void updateWellControls(WellState& xw) const;
132 
133  // TODO: should LinearisedBlackoilResidual also be a template class?
134  template <class SolutionState>
135  void addWellFluxEq(const std::vector<ADB>& cq_s,
136  const SolutionState& state,
137  LinearisedBlackoilResidual& residual);
138 
139  // TODO: some parameters, like gravity, maybe it is better to put in the member list
140  template <class SolutionState, class WellState>
141  void addWellControlEq(const SolutionState& state,
142  const WellState& xw,
143  const Vector& aliveWells,
144  LinearisedBlackoilResidual& residual);
145 
146  template <class SolutionState, class WellState>
147  void computeWellConnectionPressures(const SolutionState& state,
148  const WellState& xw);
149 
150  // state0 is non-constant, while it will not be used outside of the function
151  template <class SolutionState, class WellState>
152  void
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;
158 
159  template <class SolutionState>
160  void
161  variableStateExtractWellsVars(const std::vector<int>& indices,
162  std::vector<ADB>& vars,
163  SolutionState& state) const;
164 
165  void
166  variableStateWellIndices(std::vector<int>& indices,
167  int& next) const;
168 
169  std::vector<int>
170  variableWellStateIndices() const;
171 
172  template <class WellState>
173  void
174  variableWellStateInitials(const WellState& xw,
175  std::vector<Vector>& vars0) const;
176 
179  void setStoreWellPerforationFluxesFlag(const bool store_fluxes);
180 
184  const Vector& getStoredWellPerforationFluxes() const;
185 
187  template<class WellState>
188  void
189  updateListEconLimited(const Schedule& schedule,
190  const int current_step,
191  const Wells* wells,
192  const WellState& well_state,
193  DynamicListEconLimited& list_econ_limited) const;
194 
195 
196  WellCollection* wellCollection() const;
197 
198  void calculateEfficiencyFactors();
199 
200  const Vector& wellPerfEfficiencyFactors() const;
201 
202  protected:
203  bool wells_active_;
204  const Wells* wells_;
205  const WellOps wops_;
206  // It will probably need to be updated during running time.
207  WellCollection* well_collection_;
208 
209  // The efficiency factor for each connection
210  // It is specified based on wells and groups
211  // By default, they should all be one.
212  Vector well_perforation_efficiency_factors_;
213 
214  const BlackoilPropsAdFromDeck* fluid_;
215  const std::vector<bool>* active_;
216  const std::vector<PhasePresence>* phase_condition_;
217  const VFPProperties* vfp_properties_;
218  double gravity_;
219  // the depth of the all the cell centers
220  // for standard Wells, it the same with the perforation depth
221  Vector perf_cell_depth_;
222 
223  Vector well_perforation_densities_;
224  Vector well_perforation_pressure_diffs_;
225 
226  bool store_well_perforation_fluxes_;
227  Vector well_perforation_fluxes_;
228 
229  // protected methods
230  template <class SolutionState, class WellState>
231  void computePropertiesForWellConnectionPressures(const SolutionState& state,
232  const WellState& xw,
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);
237 
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,
245  const double grav);
246 
247 
248  template <class WellState>
249  bool checkRateEconLimits(const WellEconProductionLimits& econ_production_limits,
250  const WellState& well_state,
251  const int well_number) const;
252 
253  using WellMapType = typename WellState::WellMapType;
254  using WellMapEntryType = typename WellState::mapentry_t;
255 
256  // a tuple type for ratio limit check.
257  // first value indicates whether ratio limit is violated, when the ratio limit is not violated, the following three
258  // values should not be used.
259  // second value indicates whehter there is only one connection left.
260  // third value indicates the indx of the worst-offending connection.
261  // the last value indicates the extent of the violation for the worst-offending connection, which is defined by
262  // the ratio of the actual value to the value of the violated limit.
263  using RatioCheckTuple = std::tuple<bool, bool, int, double>;
264 
265  enum ConnectionIndex {
266  INVALIDCONNECTION = -10000
267  };
268 
269 
270  template <class WellState>
271  RatioCheckTuple checkRatioEconLimits(const WellEconProductionLimits& econ_production_limits,
272  const WellState& well_state,
273  const WellMapEntryType& map_entry) const;
274 
275  template <class WellState>
276  RatioCheckTuple checkMaxWaterCutLimit(const WellEconProductionLimits& econ_production_limits,
277  const WellState& well_state,
278  const WellMapEntryType& map_entry) const;
279 
280  template <class WellState>
281  void updateWellStateWithTarget(const WellControls* wc,
282  const int current,
283  const int well_index,
284  WellState& xw) const;
285 
286  };
287 
288 
289 } // namespace Opm
290 
291 #include "StandardWells_impl.hpp"
292 
293 #endif
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
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
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: AdditionalObjectDeleter.hpp:22
Class for handling the standard well model.
Definition: StandardWells.hpp:51
Definition: StandardWells.hpp:53
const Vector & getStoredWellPerforationFluxes() const
Retrieves the stored fluxes.
Definition: StandardWells_impl.hpp:1195
bool localWellsActive() const
return true if wells are available on this process
Definition: StandardWells_impl.hpp:148
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
bool wellsActive() const
return true if wells are available in the reservoir
Definition: StandardWells_impl.hpp:130
Eigen::Array< double, Eigen::Dynamic, 1 > V
Underlying type for values.
Definition: AutoDiffBlock.hpp:99
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