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 #ifndef OPM_BLACKOILWELLMODEL_GENERIC_HEADER_INCLUDED
24 #define OPM_BLACKOILWELLMODEL_GENERIC_HEADER_INCLUDED
25 
26 #include <opm/output/data/GuideRateValue.hpp>
27 
28 #include <opm/input/eclipse/Schedule/Well/WellTestState.hpp>
29 #include <opm/input/eclipse/Schedule/Group/GuideRate.hpp>
30 
31 #include <opm/simulators/utils/DeferredLoggingErrorHelpers.hpp>
32 
33 #include <opm/simulators/wells/PerforationData.hpp>
34 #include <opm/simulators/wells/WellProdIndexCalculator.hpp>
35 #include <opm/simulators/wells/WGState.hpp>
36 
37 #include <functional>
38 #include <map>
39 #include <memory>
40 #include <string>
41 #include <unordered_map>
42 #include <unordered_set>
43 #include <vector>
44 
45 namespace Opm {
46  class DeferredLogger;
47  class EclipseState;
48  class GasLiftSingleWellGeneric;
49  class GasLiftWellState;
50  class Group;
51  class GuideRateConfig;
52  class ParallelWellInfo;
53  class RestartValue;
54  class Schedule;
55  class SummaryConfig;
56  class VFPProperties;
57  class WellInterfaceGeneric;
58  class WellState;
59 } // namespace Opm
60 
61 namespace Opm { namespace data {
62  struct GroupData;
63  struct GroupGuideRates;
64  class GroupAndNetworkValues;
65  struct NodeData;
66 }} // namespace Opm::data
67 
68 namespace Opm {
69 
72 {
73 public:
74  // --------- Types ---------
75  using GLiftOptWells = std::map<std::string, std::unique_ptr<GasLiftSingleWellGeneric>>;
76  using GLiftProdWells = std::map<std::string, const WellInterfaceGeneric*>;
77  using GLiftWellStateMap = std::map<std::string, std::unique_ptr<GasLiftWellState>>;
78 
79  BlackoilWellModelGeneric(Schedule& schedule,
80  const SummaryState& summaryState,
81  const EclipseState& eclState,
82  const PhaseUsage& phase_usage,
83  const Parallel::Communication& comm);
84 
85  virtual ~BlackoilWellModelGeneric() = default;
86 
87  int numLocalWells() const;
88  int numPhases() const;
89 
91  bool wellsActive() const;
92  bool hasWell(const std::string& wname);
94  bool localWellsActive() const;
95  // whether there exists any multisegment well open on this process
96  bool anyMSWellOpenLocal() const;
97 
98  const Well& getWellEcl(const std::string& well_name) const;
99  std::vector<Well> getLocalWells(const int timeStepIdx) const;
100  const Schedule& schedule() const { return schedule_; }
101  const PhaseUsage& phaseUsage() const { return phase_usage_; }
102  const GroupState& groupState() const { return this->active_wgstate_.group_state; }
103 
104  /*
105  Immutable version of the currently active wellstate.
106  */
107  const WellState& wellState() const
108  {
109  return this->active_wgstate_.well_state;
110  }
111 
112  /*
113  Mutable version of the currently active wellstate.
114  */
115  WellState& wellState()
116  {
117  return this->active_wgstate_.well_state;
118  }
119 
120  GroupState& groupState() { return this->active_wgstate_.group_state; }
121 
122  WellTestState& wellTestState() { return this->active_wgstate_.well_test_state; }
123 
124  const WellTestState& wellTestState() const { return this->active_wgstate_.well_test_state; }
125 
126 
127  double wellPI(const int well_index) const;
128  double wellPI(const std::string& well_name) const;
129 
130  void updateEclWells(const int timeStepIdx,
131  const std::unordered_set<std::string>& wells,
132  const SummaryState& st);
133 
134 
135  void loadRestartData(const data::Wells& rst_wells,
136  const data::GroupAndNetworkValues& grpNwrkValues,
137  const PhaseUsage& phases,
138  const bool handle_ms_well,
139  WellState& well_state);
140 
141  void initFromRestartFile(const RestartValue& restartValues,
142  WellTestState wtestState,
143  const size_t numCells,
144  bool handle_ms_well);
145 
146  void setWellsActive(const bool wells_active);
147 
148  /*
149  Will assign the internal member last_valid_well_state_ to the
150  current value of the this->active_well_state_. The state stored
151  with storeWellState() can then subsequently be recovered with the
152  resetWellState() method.
153  */
154  void commitWGState()
155  {
156  this->last_valid_wgstate_ = this->active_wgstate_;
157  }
158 
159  data::GroupAndNetworkValues groupAndNetworkData(const int reportStepIdx) const;
160 
162  bool hasTHPConstraints() const;
163 
166  bool forceShutWellByName(const std::string& wellname,
167  const double simulation_time);
168 
169 protected:
170 
171  /*
172  The dynamic state of the well model is maintained with an instance
173  of the WellState class. Currently we have
174  three different wellstate instances:
175 
176  1. The currently active wellstate is in the active_well_state_
177  member. That is the state which is mutated by the simulator.
178 
179  2. In the case timestep fails to converge and we must go back and
180  try again with a smaller timestep we need to recover the last
181  valid wellstate. This is maintained with the
182  last_valid_well_state_ member and the functions
183  commitWellState() and resetWellState().
184 
185  3. For the NUPCOL functionality we should either use the
186  currently active wellstate or a wellstate frozen at max
187  nupcol iterations. This is handled with the member
188  nupcol_well_state_ and the initNupcolWellState() function.
189  */
190 
191 
192  /*
193  Will return the last good wellstate. This is typcially used when
194  initializing a new report step where the Schedule object might
195  have introduced new wells. The wellstate returned by
196  prevWellState() must have been stored with the commitWellState()
197  function first.
198  */
199  const WellState& prevWellState() const
200  {
201  return this->last_valid_wgstate_.well_state;
202  }
203 
204  const WGState& prevWGState() const
205  {
206  return this->last_valid_wgstate_;
207  }
208  /*
209  Will return the currently active nupcolWellState; must initialize
210  the internal nupcol wellstate with initNupcolWellState() first.
211  */
212  const WellState& nupcolWellState() const
213  {
214  return this->nupcol_wgstate_.well_state;
215  }
216 
217  /*
218  Will store a copy of the input argument well_state in the
219  last_valid_well_state_ member, that state can then be recovered
220  with a subsequent call to resetWellState().
221  */
222  void commitWGState(WGState wgstate)
223  {
224  this->last_valid_wgstate_ = std::move(wgstate);
225  }
226 
227  /*
228  Will update the internal variable active_well_state_ to whatever
229  was stored in the last_valid_well_state_ member. This function
230  works in pair with commitWellState() which should be called first.
231  */
232  void resetWGState()
233  {
234  this->active_wgstate_ = this->last_valid_wgstate_;
235  }
236 
237  /*
238  Will store the current active wellstate in the nupcol_well_state_
239  member. This can then be subsequently retrieved with accessor
240  nupcolWellState().
241  */
242  void updateNupcolWGState()
243  {
244  this->nupcol_wgstate_ = this->active_wgstate_;
245  }
246 
249  std::vector<std::reference_wrapper<ParallelWellInfo>> createLocalParallelWellInfo(const std::vector<Well>& wells);
250 
251  void initializeWellProdIndCalculators();
252  void initializeWellPerfData();
253 
254  bool wasDynamicallyShutThisTimeStep(const int well_index) const;
255 
256  void updateNetworkPressures(const int reportStepIdx);
257 
258  void updateWsolvent(const Group& group,
259  const int reportStepIdx,
260  const WellState& wellState);
261  void setWsolvent(const Group& group,
262  const int reportStepIdx,
263  double wsolvent);
264  virtual void calcRates(const int fipnum,
265  const int pvtreg,
266  std::vector<double>& resv_coeff) = 0;
267  virtual void calcInjRates(const int fipnum,
268  const int pvtreg,
269  std::vector<double>& resv_coeff) = 0;
270 
271  data::GuideRateValue getGuideRateValues(const Group& group) const;
272  data::GuideRateValue getGuideRateValues(const Well& well) const;
273  data::GuideRateValue getGuideRateInjectionGroupValues(const Group& group) const;
274  void getGuideRateValues(const GuideRate::RateVector& qs,
275  const bool is_inj,
276  const std::string& wgname,
277  data::GuideRateValue& grval) const;
278 
279  void assignWellGuideRates(data::Wells& wsrpt,
280  const int reportStepIdx) const;
281  void assignShutConnections(data::Wells& wsrpt,
282  const int reportStepIndex) const;
283  void assignGroupControl(const Group& group,
284  data::GroupData& gdata) const;
285  void assignGroupGuideRates(const Group& group,
286  const std::unordered_map<std::string, data::GroupGuideRates>& groupGuideRates,
287  data::GroupData& gdata) const;
288  void assignGroupValues(const int reportStepIdx,
289  std::map<std::string, data::GroupData>& gvalues) const;
290  void assignNodeValues(std::map<std::string, data::NodeData>& nodevalues) const;
291 
292  void loadRestartConnectionData(const std::vector<data::Rates::opt>& phs,
293  const data::Well& rst_well,
294  const std::vector<PerforationData>& old_perf_data,
295  SingleWellState& ws);
296 
297  void loadRestartSegmentData(const std::string& well_name,
298  const std::vector<data::Rates::opt>& phs,
299  const data::Well& rst_well,
300  SingleWellState& ws);
301 
302  void loadRestartWellData(const std::string& well_name,
303  const bool handle_ms_well,
304  const std::vector<data::Rates::opt>& phs,
305  const data::Well& rst_well,
306  const std::vector<PerforationData>& old_perf_data,
307  SingleWellState& ws);
308 
309  void loadRestartGroupData(const std::string& group,
310  const data::GroupData& value);
311 
312  void loadRestartGuideRates(const int report_step,
313  const GuideRateModel::Target target,
314  const data::Wells& rst_wells);
315 
316  void loadRestartGuideRates(const int report_step,
317  const GuideRateConfig& config,
318  const std::map<std::string, data::GroupData>& rst_groups);
319 
320  std::unordered_map<std::string, data::GroupGuideRates>
321  calculateAllGroupGuiderates(const int reportStepIdx) const;
322 
323  void calculateEfficiencyFactors(const int reportStepIdx);
324 
325  bool checkGroupConstraints(const Group& group,
326  const int reportStepIdx,
327  DeferredLogger& deferred_logger) const;
328 
329  std::pair<Group::InjectionCMode, double> checkGroupInjectionConstraints(const Group& group,
330  const int reportStepIdx,
331  const Phase& phase) const;
332  std::pair<Group::ProductionCMode, double> checkGroupProductionConstraints(const Group& group,
333  const int reportStepIdx,
334  DeferredLogger& deferred_logger) const;
335 
336  void checkGconsaleLimits(const Group& group,
337  WellState& well_state,
338  const int reportStepIdx,
339  DeferredLogger& deferred_logger);
340 
341  bool checkGroupHigherConstraints(const Group& group,
342  DeferredLogger& deferred_logger,
343  const int reportStepIdx,
344  std::set<std::string>& switched_groups);
345 
346  bool updateGroupIndividualControl(const Group& group,
347  DeferredLogger& deferred_logger,
348  const int reportStepIdx,
349  std::set<std::string>& switched_groups);
350 
351  bool updateGroupIndividualControls(DeferredLogger& deferred_logger,
352  std::set<std::string>& switched_groups,
353  const int reportStepIdx,
354  const int iterationIdx);
355 
356  bool updateGroupHigherControls(DeferredLogger& deferred_logger,
357  const int reportStepIdx,
358  std::set<std::string>& switched_groups);
359 
360  void actionOnBrokenConstraints(const Group& group,
361  const Group::ExceedAction& exceed_action,
362  const Group::ProductionCMode& newControl,
363  DeferredLogger& deferred_logger);
364  void actionOnBrokenConstraints(const Group& group,
365  const Group::InjectionCMode& newControl,
366  const Phase& controlPhase,
367  DeferredLogger& deferred_logger);
368 
369  void updateAndCommunicateGroupData(const int reportStepIdx,
370  const int iterationIdx);
371 
372  void inferLocalShutWells();
373 
374  void setRepRadiusPerfLength();
375 
376  void gliftDebug(const std::string& msg,
377  DeferredLogger& deferred_logger) const;
378 
379  void gliftDebugShowALQ(DeferredLogger& deferred_logger);
380 
381  void gasLiftOptimizationStage2(DeferredLogger& deferred_logger,
382  GLiftProdWells& prod_wells,
383  GLiftOptWells& glift_wells,
384  GLiftWellStateMap& map,
385  const int episodeIndex);
386 
387  virtual void computePotentials(const std::size_t widx,
388  const WellState& well_state_copy,
389  std::string& exc_msg,
390  ExceptionType::ExcEnum& exc_type,
391  DeferredLogger& deferred_logger) = 0;
392 
393  // Calculating well potentials for each well
394  void updateWellPotentials(const int reportStepIdx,
395  const bool onlyAfterEvent,
396  const SummaryConfig& summaryConfig,
397  DeferredLogger& deferred_logger);
398 
399  bool guideRateUpdateIsNeeded(const int reportStepIdx) const;
400 
401  // create the well container
402  virtual void createWellContainer(const int time_step) = 0;
403  virtual void initWellContainer() = 0;
404 
405  virtual void calculateProductivityIndexValuesShutWells(const int reportStepIdx,
406  DeferredLogger& deferred_logger) = 0;
407  virtual void calculateProductivityIndexValues(DeferredLogger& deferred_logger) = 0;
408 
409  void runWellPIScaling(const int timeStepIdx,
410  DeferredLogger& local_deferredLogger);
411 
413  virtual int compressedIndexForInterior(int cartesian_cell_idx) const = 0;
414 
415 
416  Schedule& schedule_;
417  const SummaryState& summaryState_;
418  const EclipseState& eclState_;
419  const Parallel::Communication& comm_;
420 
421  PhaseUsage phase_usage_;
422  bool terminal_output_{false};
423  bool wells_active_{false};
424  bool initial_step_{};
425  bool report_step_starts_{};
426 
427  std::optional<int> last_run_wellpi_{};
428 
429  std::vector<Well> wells_ecl_;
430  std::vector<std::vector<PerforationData>> well_perf_data_;
431  std::function<bool(const Well&)> not_on_process_{};
432 
433  // a vector of all the wells.
434  std::vector<WellInterfaceGeneric*> well_container_generic_{};
435 
436  std::vector<int> local_shut_wells_{};
437 
438  std::vector<ParallelWellInfo> parallel_well_info_;
439  std::vector<std::reference_wrapper<ParallelWellInfo>> local_parallel_well_info_;
440 
441  std::vector<WellProdIndexCalculator> prod_index_calc_;
442 
443  std::vector<int> pvt_region_idx_;
444 
445  mutable std::unordered_set<std::string> closed_this_step_;
446 
447  GuideRate guideRate_;
448  std::unique_ptr<VFPProperties> vfp_properties_{};
449  std::map<std::string, double> node_pressures_; // Storing network pressures for output.
450 
451  /*
452  The various wellState members should be accessed and modified
453  through the accessor functions wellState(), prevWellState(),
454  commitWellState(), resetWellState(), nupcolWellState() and
455  updateNupcolWellState().
456  */
457  WGState active_wgstate_;
458  WGState last_valid_wgstate_;
459  WGState nupcol_wgstate_;
460 
461  bool glift_debug = false;
462 
463  double last_glift_opt_time_ = -1.0;
464 
465 private:
466  WellInterfaceGeneric* getGenWell(const std::string& well_name);
467 };
468 
469 
470 } // namespace Opm
471 
472 #endif
Class for handling the blackoil well model.
Definition: BlackoilWellModelGeneric.hpp:72
virtual int compressedIndexForInterior(int cartesian_cell_idx) const =0
get compressed index for interior cells (-1, otherwise
std::vector< std::reference_wrapper< ParallelWellInfo > > createLocalParallelWellInfo(const std::vector< Well > &wells)
Create the parallel well information.
Definition: BlackoilWellModelGeneric.cpp:691
bool wellsActive() const
return true if wells are available in the reservoir
Definition: BlackoilWellModelGeneric.cpp:384
bool hasTHPConstraints() const
Return true if any well has a THP constraint.
Definition: BlackoilWellModelGeneric.cpp:2050
bool forceShutWellByName(const std::string &wellname, const double simulation_time)
Shut down any single well Returns true if the well was actually found and shut.
Definition: BlackoilWellModelGeneric.cpp:2063
bool localWellsActive() const
return true if wells are available on this process
Definition: BlackoilWellModelGeneric.cpp:391
Definition: DeferredLogger.hpp:57
Definition: GroupState.hpp:34
Definition: SingleWellState.hpp:38
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:27
Definition: BlackoilPhases.hpp:46
Definition: WGState.hpp:35