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