My Project
BlackoilWellModelGeneric.hpp
1 /*
2  Copyright 2016 SINTEF ICT, Applied Mathematics.
3  Copyright 2016 - 2017 Statoil ASA.
4  Copyright 2017 Dr. Blatt - HPC-Simulation-Software & Services
5  Copyright 2016 - 2018 IRIS AS
6 
7  This file is part of the Open Porous Media project (OPM).
8 
9  OPM is free software: you can redistribute it and/or modify
10  it under the terms of the GNU General Public License as published by
11  the Free Software Foundation, either version 3 of the License, or
12  (at your option) any later version.
13 
14  OPM is distributed in the hope that it will be useful,
15  but WITHOUT ANY WARRANTY; without even the implied warranty of
16  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17  GNU General Public License for more details.
18 
19  You should have received a copy of the GNU General Public License
20  along with OPM. If not, see <http://www.gnu.org/licenses/>.
21 */
22 
23 
24 #ifndef OPM_BLACKOILWELLMODEL_GENERIC_HEADER_INCLUDED
25 #define OPM_BLACKOILWELLMODEL_GENERIC_HEADER_INCLUDED
26 
27 #include <functional>
28 #include <map>
29 #include <memory>
30 #include <string>
31 #include <unordered_map>
32 #include <unordered_set>
33 #include <vector>
34 
35 #include <dune/common/version.hh>
36 #include <opm/output/data/GuideRateValue.hpp>
37 #include <opm/parser/eclipse/EclipseState/Schedule/Well/WellTestState.hpp>
38 #include <opm/parser/eclipse/EclipseState/Schedule/Group/GuideRate.hpp>
39 
40 #include <opm/simulators/utils/DeferredLoggingErrorHelpers.hpp>
41 #include <opm/simulators/wells/ParallelWellInfo.hpp>
42 #include <opm/simulators/wells/PerforationData.hpp>
43 #include <opm/simulators/wells/WellInterfaceGeneric.hpp>
44 #include <opm/simulators/wells/WellProdIndexCalculator.hpp>
45 #include <opm/simulators/wells/WGState.hpp>
46 
47 namespace Opm {
48 
49 namespace data {
50 struct GroupData;
51 struct GroupGuideRates;
52 class GroupAndNetworkValues;
53 struct NodeData;
54 }
55 
56 class DeferredLogger;
57 class EclipseState;
58 class GasLiftSingleWellGeneric;
59 class GasLiftWellState;
60 class Group;
61 class RestartValue;
62 class Schedule;
63 class SummaryConfig;
64 class VFPProperties;
65 class WellState;
66 
69 {
70 public:
71  // --------- Types ---------
72  using GLiftOptWells = std::map<std::string,std::unique_ptr<GasLiftSingleWellGeneric>>;
73  using GLiftProdWells = std::map<std::string,const WellInterfaceGeneric*>;
74  using GLiftWellStateMap = std::map<std::string,std::unique_ptr<GasLiftWellState>>;
75 
76  BlackoilWellModelGeneric(Schedule& schedule,
77  const SummaryState& summaryState,
78  const EclipseState& eclState,
79  const PhaseUsage& phase_usage,
80  const Parallel::Communication& comm);
81 
82  virtual ~BlackoilWellModelGeneric() = default;
83 
84  int numLocalWells() const;
85  int numPhases() const;
86 
88  bool wellsActive() const;
89  bool hasWell(const std::string& wname);
91  bool localWellsActive() const;
92  // whether there exists any multisegment well open on this process
93  bool anyMSWellOpenLocal() const;
94 
95  const Well& getWellEcl(const std::string& well_name) const;
96  std::vector<Well> getLocalWells(const int timeStepIdx) const;
97  const Schedule& schedule() const { return schedule_; }
98  const PhaseUsage& phaseUsage() const { return phase_usage_; }
99  const GroupState& groupState() const { return this->active_wgstate_.group_state; }
100 
101  /*
102  Immutable version of the currently active wellstate.
103  */
104  const WellState& wellState() const
105  {
106  return this->active_wgstate_.well_state;
107  }
108 
109  /*
110  Mutable version of the currently active wellstate.
111  */
112  WellState& wellState()
113  {
114  return this->active_wgstate_.well_state;
115  }
116 
117  GroupState& groupState() { return this->active_wgstate_.group_state; }
118 
119  WellTestState& wellTestState() { return this->active_wgstate_.well_test_state; }
120 
121  const WellTestState& wellTestState() const { return this->active_wgstate_.well_test_state; }
122 
123 
124  double wellPI(const int well_index) const;
125  double wellPI(const std::string& well_name) const;
126 
127  void updateEclWells(const int timeStepIdx,
128  const std::unordered_set<std::string>& wells);
129 
130 
131  void loadRestartData(const data::Wells& rst_wells,
132  const data::GroupAndNetworkValues& grpNwrkValues,
133  const PhaseUsage& phases,
134  const bool handle_ms_well,
135  WellState& well_state);
136 
137  void initFromRestartFile(const RestartValue& restartValues,
138  WellTestState wtestState,
139  const size_t numCells,
140  bool handle_ms_well);
141 
142  void setWellsActive(const bool wells_active);
143 
144  /*
145  Will assign the internal member last_valid_well_state_ to the
146  current value of the this->active_well_state_. The state stored
147  with storeWellState() can then subsequently be recovered with the
148  resetWellState() method.
149  */
150  void commitWGState()
151  {
152  this->last_valid_wgstate_ = this->active_wgstate_;
153  }
154 
155  data::GroupAndNetworkValues groupAndNetworkData(const int reportStepIdx) const;
156 
158  bool hasTHPConstraints() const;
159 
162  bool forceShutWellByNameIfPredictionMode(const std::string& wellname,
163  const double simulation_time);
164 
165 protected:
166 
167  /*
168  The dynamic state of the well model is maintained with an instance
169  of the WellState class. Currently we have
170  three different wellstate instances:
171 
172  1. The currently active wellstate is in the active_well_state_
173  member. That is the state which is mutated by the simulator.
174 
175  2. In the case timestep fails to converge and we must go back and
176  try again with a smaller timestep we need to recover the last
177  valid wellstate. This is maintained with the
178  last_valid_well_state_ member and the functions
179  commitWellState() and resetWellState().
180 
181  3. For the NUPCOL functionality we should either use the
182  currently active wellstate or a wellstate frozen at max
183  nupcol iterations. This is handled with the member
184  nupcol_well_state_ and the initNupcolWellState() function.
185  */
186 
187 
188  /*
189  Will return the last good wellstate. This is typcially used when
190  initializing a new report step where the Schedule object might
191  have introduced new wells. The wellstate returned by
192  prevWellState() must have been stored with the commitWellState()
193  function first.
194  */
195  const WellState& prevWellState() const
196  {
197  return this->last_valid_wgstate_.well_state;
198  }
199 
200  const WGState& prevWGState() const
201  {
202  return this->last_valid_wgstate_;
203  }
204  /*
205  Will return the currently active nupcolWellState; must initialize
206  the internal nupcol wellstate with initNupcolWellState() first.
207  */
208  const WellState& nupcolWellState() const
209  {
210  return this->nupcol_wgstate_.well_state;
211  }
212 
213  /*
214  Will store a copy of the input argument well_state in the
215  last_valid_well_state_ member, that state can then be recovered
216  with a subsequent call to resetWellState().
217  */
218  void commitWGState(WGState wgstate)
219  {
220  this->last_valid_wgstate_ = std::move(wgstate);
221  }
222 
223  /*
224  Will update the internal variable active_well_state_ to whatever
225  was stored in the last_valid_well_state_ member. This function
226  works in pair with commitWellState() which should be called first.
227  */
228  void resetWGState()
229  {
230  this->active_wgstate_ = this->last_valid_wgstate_;
231  }
232 
233  /*
234  Will store the current active wellstate in the nupcol_well_state_
235  member. This can then be subsequently retrieved with accessor
236  nupcolWellState().
237  */
238  void updateNupcolWGState()
239  {
240  this->nupcol_wgstate_ = this->active_wgstate_;
241  }
242 
245  std::vector<std::reference_wrapper<ParallelWellInfo>> createLocalParallelWellInfo(const std::vector<Well>& wells);
246 
247  void initializeWellProdIndCalculators();
248  void initializeWellPerfData();
249 
250  bool wasDynamicallyShutThisTimeStep(const int well_index) const;
251 
252  void updateNetworkPressures(const int reportStepIdx);
253 
254  void updateWsolvent(const Group& group,
255  const int reportStepIdx,
256  const WellState& wellState);
257  void setWsolvent(const Group& group,
258  const int reportStepIdx,
259  double wsolvent);
260  virtual void calcRates(const int fipnum,
261  const int pvtreg,
262  std::vector<double>& resv_coeff) = 0;
263  virtual void calcInjRates(const int fipnum,
264  const int pvtreg,
265  std::vector<double>& resv_coeff) = 0;
266 
267  data::GuideRateValue getGuideRateValues(const Group& group) const;
268  data::GuideRateValue getGuideRateValues(const Well& well) const;
269  data::GuideRateValue getGuideRateInjectionGroupValues(const Group& group) const;
270  void getGuideRateValues(const GuideRate::RateVector& qs,
271  const bool is_inj,
272  const std::string& wgname,
273  data::GuideRateValue& grval) const;
274 
275  void assignWellGuideRates(data::Wells& wsrpt) const;
276  void assignShutConnections(data::Wells& wsrpt,
277  const int reportStepIndex) const;
278  void assignGroupControl(const Group& group,
279  data::GroupData& gdata) const;
280  void assignGroupGuideRates(const Group& group,
281  const std::unordered_map<std::string, data::GroupGuideRates>& groupGuideRates,
282  data::GroupData& gdata) const;
283  void assignGroupValues(const int reportStepIdx,
284  std::map<std::string, data::GroupData>& gvalues) const;
285  void assignNodeValues(std::map<std::string, data::NodeData>& nodevalues) const;
286 
287  std::unordered_map<std::string, data::GroupGuideRates>
288  calculateAllGroupGuiderates(const int reportStepIdx) const;
289 
290  void calculateEfficiencyFactors(const int reportStepIdx);
291 
292  bool checkGroupConstraints(const Group& group,
293  const int reportStepIdx,
294  DeferredLogger& deferred_logger) const;
295 
296  std::pair<Group::InjectionCMode, double> checkGroupInjectionConstraints(const Group& group,
297  const int reportStepIdx,
298  const Phase& phase) const;
299  std::pair<Group::ProductionCMode, double> checkGroupProductionConstraints(const Group& group,
300  const int reportStepIdx,
301  DeferredLogger& deferred_logger) const;
302 
303  void checkGconsaleLimits(const Group& group,
304  WellState& well_state,
305  const int reportStepIdx,
306  DeferredLogger& deferred_logger);
307 
308  bool checkGroupHigherConstraints(const Group& group,
309  DeferredLogger& deferred_logger,
310  const int reportStepIdx,
311  std::set<std::string>& switched_groups);
312 
313  bool updateGroupIndividualControl(const Group& group,
314  DeferredLogger& deferred_logger,
315  const int reportStepIdx,
316  std::set<std::string>& switched_groups);
317 
318  bool updateGroupIndividualControls(DeferredLogger& deferred_logger,
319  std::set<std::string>& switched_groups,
320  const int reportStepIdx,
321  const int iterationIdx);
322 
323  bool updateGroupHigherControls(DeferredLogger& deferred_logger,
324  const int reportStepIdx,
325  std::set<std::string>& switched_groups);
326 
327  void actionOnBrokenConstraints(const Group& group,
328  const Group::ExceedAction& exceed_action,
329  const Group::ProductionCMode& newControl,
330  DeferredLogger& deferred_logger);
331  void actionOnBrokenConstraints(const Group& group,
332  const Group::InjectionCMode& newControl,
333  const Phase& controlPhase,
334  DeferredLogger& deferred_logger);
335 
336  void updateAndCommunicateGroupData(const int reportStepIdx,
337  const int iterationIdx);
338 
339  void inferLocalShutWells();
340 
341  void setRepRadiusPerfLength();
342 
343  void gliftDebug(const std::string& msg,
344  DeferredLogger& deferred_logger) const;
345 
346  void gliftDebugShowALQ(DeferredLogger& deferred_logger);
347 
348  void gasLiftOptimizationStage2(DeferredLogger& deferred_logger,
349  GLiftProdWells& prod_wells,
350  GLiftOptWells& glift_wells,
351  GLiftWellStateMap& map,
352  const int episodeIndex);
353 
354  virtual void computePotentials(const std::size_t widx,
355  const WellState& well_state_copy,
356  std::string& exc_msg,
357  ExceptionType::ExcEnum& exc_type,
358  DeferredLogger& deferred_logger) = 0;
359 
360  // Calculating well potentials for each well
361  void updateWellPotentials(const int reportStepIdx,
362  const bool onlyAfterEvent,
363  const SummaryConfig& summaryConfig,
364  DeferredLogger& deferred_logger);
365 
366  // create the well container
367  virtual void createWellContainer(const int time_step) = 0;
368  virtual void initWellContainer() = 0;
369 
370  virtual void calculateProductivityIndexValuesShutWells(const int reportStepIdx,
371  DeferredLogger& deferred_logger) = 0;
372  virtual void calculateProductivityIndexValues(DeferredLogger& deferred_logger) = 0;
373 
374  void runWellPIScaling(const int timeStepIdx,
375  DeferredLogger& local_deferredLogger);
376 
377  Schedule& schedule_;
378  const SummaryState& summaryState_;
379  const EclipseState& eclState_;
380  const Parallel::Communication& comm_;
381 
382  PhaseUsage phase_usage_;
383  bool terminal_output_{false};
384  bool wells_active_{false};
385  bool initial_step_{};
386  bool report_step_starts_{};
387 
388  std::optional<int> last_run_wellpi_{};
389 
390  std::vector<Well> wells_ecl_;
391  std::vector<std::vector<PerforationData>> well_perf_data_;
392  std::function<bool(const Well&)> not_on_process_{};
393 
394  // a vector of all the wells.
395  std::vector<WellInterfaceGeneric*> well_container_generic_{};
396 
397  std::vector<int> local_shut_wells_{};
398 
399  std::vector<ParallelWellInfo> parallel_well_info_;
400  std::vector<std::reference_wrapper<ParallelWellInfo>> local_parallel_well_info_;
401 
402  std::vector<WellProdIndexCalculator> prod_index_calc_;
403 
404  // Map from logically cartesian cell indices to compressed ones.
405  // Cells not in the interior are not mapped. This deactivates
406  // these for distributed wells and makes the distribution non-overlapping.
407  std::vector<int> cartesian_to_compressed_;
408 
409  std::vector<int> pvt_region_idx_;
410 
411  mutable std::unordered_set<std::string> closed_this_step_;
412 
413  GuideRate guideRate_;
414  std::unique_ptr<VFPProperties> vfp_properties_{};
415  std::map<std::string, double> node_pressures_; // Storing network pressures for output.
416 
417  /*
418  The various wellState members should be accessed and modified
419  through the accessor functions wellState(), prevWellState(),
420  commitWellState(), resetWellState(), nupcolWellState() and
421  updateNupcolWellState().
422  */
423  WGState active_wgstate_;
424  WGState last_valid_wgstate_;
425  WGState nupcol_wgstate_;
426 
427  bool glift_debug = false;
428 
429  private:
430  WellInterfaceGeneric* getGenWell(const std::string& well_name);
431 };
432 
433 
434 } // namespace Opm
435 
436 #endif
Class for handling the blackoil well model.
Definition: BlackoilWellModelGeneric.hpp:69
std::vector< std::reference_wrapper< ParallelWellInfo > > createLocalParallelWellInfo(const std::vector< Well > &wells)
Create the parallel well information.
Definition: BlackoilWellModelGeneric.cpp:305
bool forceShutWellByNameIfPredictionMode(const std::string &wellname, const double simulation_time)
Shut down any single well, but only if it is in prediction mode.
Definition: BlackoilWellModelGeneric.cpp:1589
bool wellsActive() const
return true if wells are available in the reservoir
Definition: BlackoilWellModelGeneric.cpp:105
bool hasTHPConstraints() const
Return true if any well has a THP constraint.
Definition: BlackoilWellModelGeneric.cpp:1576
bool localWellsActive() const
return true if wells are available on this process
Definition: BlackoilWellModelGeneric.cpp:112
Definition: DeferredLogger.hpp:57
Definition: GroupState.hpp:34
Definition: WellInterfaceGeneric.hpp:51
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition: WellState.hpp:56
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:26
Definition: BlackoilPhases.hpp:45
Definition: WGState.hpp:35