My Project
Runspec.hpp
1 /*
2  Copyright 2016 Statoil ASA.
3 
4  This file is part of the Open Porous Media project (OPM).
5 
6  OPM is free software: you can redistribute it and/or modify it under the terms
7  of the GNU General Public License as published by the Free Software
8  Foundation, either version 3 of the License, or (at your option) any later
9  version.
10 
11  OPM is distributed in the hope that it will be useful, but WITHOUT ANY
12  WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
13  A PARTICULAR PURPOSE. See the GNU General Public License for more details.
14 
15  You should have received a copy of the GNU General Public License along with
16  OPM. If not, see <http://www.gnu.org/licenses/>.
17 */
18 
19 #ifndef OPM_RUNSPEC_HPP
20 #define OPM_RUNSPEC_HPP
21 
22 #include <iosfwd>
23 #include <string>
24 #include <optional>
25 
26 #include <opm/common/OpmLog/KeywordLocation.hpp>
27 #include <opm/input/eclipse/EclipseState/Tables/Tabdims.hpp>
28 #include <opm/input/eclipse/EclipseState/Tables/Regdims.hpp>
29 #include <opm/input/eclipse/EclipseState/EndpointScaling.hpp>
30 #include <opm/input/eclipse/Schedule/UDQ/UDQParams.hpp>
31 #include <opm/input/eclipse/Schedule/Action/Actdims.hpp>
32 
33 namespace Opm {
34 class Deck;
35 
36 
37 enum class Phase {
38  OIL = 0,
39  GAS = 1,
40  WATER = 2,
41  SOLVENT = 3,
42  POLYMER = 4,
43  ENERGY = 5,
44  POLYMW = 6,
45  FOAM = 7,
46  BRINE = 8,
47  ZFRACTION = 9
48 
49  // If you add more entries to this enum, remember to update NUM_PHASES_IN_ENUM below.
50 };
51 
52 constexpr int NUM_PHASES_IN_ENUM = static_cast<int>(Phase::ZFRACTION) + 1; // Used to get correct size of the bitset in class Phases.
53 
54 Phase get_phase( const std::string& );
55 std::ostream& operator<<( std::ostream&, const Phase& );
56 
57 class Phases {
58  public:
59  Phases() noexcept = default;
60  Phases( bool oil, bool gas, bool wat, bool solvent = false, bool polymer = false, bool energy = false,
61  bool polymw = false, bool foam = false, bool brine = false, bool zfraction = false ) noexcept;
62 
63  static Phases serializeObject();
64 
65  bool active( Phase ) const noexcept;
66  size_t size() const noexcept;
67 
68  bool operator==(const Phases& data) const;
69 
70  template<class Serializer>
71  void serializeOp(Serializer& serializer)
72  {
73  if (serializer.isSerializing())
74  serializer(bits.to_ulong());
75  else {
76  unsigned long Bits = 0;
77  serializer(Bits);
78  bits = std::bitset<NUM_PHASES_IN_ENUM>(Bits);
79  }
80  }
81 
82  private:
83  std::bitset< NUM_PHASES_IN_ENUM > bits;
84 };
85 
86 
87 class Welldims {
88 public:
89  Welldims() = default;
90  explicit Welldims(const Deck& deck);
91 
92  static Welldims serializeObject();
93 
94  int maxConnPerWell() const
95  {
96  return this->nCWMax;
97  }
98 
99  int maxWellsPerGroup() const
100  {
101  return this->nWGMax;
102  }
103 
104  int maxGroupsInField() const
105  {
106  return this->nGMax;
107  }
108 
109  int maxWellsInField() const
110  {
111  return this->nWMax;
112  }
113 
114  int maxWellListsPrWell() const
115  {
116  return this->nWlistPrWellMax;
117  }
118 
119  int maxDynamicWellLists() const
120  {
121  return this->nDynWlistMax;
122  }
123 
124  const std::optional<KeywordLocation>& location() const {
125  return this->m_location;
126  }
127 
128  static bool rst_cmp(const Welldims& full_dims, const Welldims& rst_dims) {
129  return full_dims.maxConnPerWell() == rst_dims.maxConnPerWell() &&
130  full_dims.maxWellsPerGroup() == rst_dims.maxWellsPerGroup() &&
131  full_dims.maxGroupsInField() == rst_dims.maxGroupsInField() &&
132  full_dims.maxWellsInField() == rst_dims.maxWellsInField() &&
133  full_dims.maxWellListsPrWell() == rst_dims.maxWellListsPrWell() &&
134  full_dims.maxDynamicWellLists() == rst_dims.maxDynamicWellLists();
135  }
136 
137  bool operator==(const Welldims& data) const {
138  return this->location() == data.location() &&
139  rst_cmp(*this, data);
140  }
141 
142 
143  template<class Serializer>
144  void serializeOp(Serializer& serializer)
145  {
146  serializer(nWMax);
147  serializer(nCWMax);
148  serializer(nWGMax);
149  serializer(nGMax);
150  serializer(nWlistPrWellMax);
151  serializer(nDynWlistMax);
152  serializer(m_location);
153  }
154 
155 private:
156  int nWMax { 0 };
157  int nCWMax { 0 };
158  int nWGMax { 0 };
159  int nGMax { 0 };
160  int nWlistPrWellMax { 1 };
161  int nDynWlistMax { 1 };
162  std::optional<KeywordLocation> m_location;
163 };
164 
166 public:
167  WellSegmentDims();
168  explicit WellSegmentDims(const Deck& deck);
169 
170  static WellSegmentDims serializeObject();
171 
172 
173  int maxSegmentedWells() const
174  {
175  return this->nSegWellMax;
176  }
177 
178  int maxSegmentsPerWell() const
179  {
180  return this->nSegmentMax;
181  }
182 
183  int maxLateralBranchesPerWell() const
184  {
185  return this->nLatBranchMax;
186  }
187 
188  bool operator==(const WellSegmentDims& data) const;
189 
190  template<class Serializer>
191  void serializeOp(Serializer& serializer)
192  {
193  serializer(nSegWellMax);
194  serializer(nSegmentMax);
195  serializer(nLatBranchMax);
196  }
197 
198 private:
199  int nSegWellMax;
200  int nSegmentMax;
201  int nLatBranchMax;
202 };
203 
204 class NetworkDims {
205 public:
206  NetworkDims();
207  explicit NetworkDims(const Deck& deck);
208 
209  static NetworkDims serializeObject();
210 
211  int maxNONodes() const
212  {
213  return this->nMaxNoNodes;
214  }
215 
216  int maxNoBranches() const
217  {
218  return this->nMaxNoBranches;
219  }
220 
221  int maxNoBranchesConToNode() const
222  {
223  return this->nMaxNoBranchesConToNode;
224  }
225 
226  bool active() const;
227 
228  bool operator==(const NetworkDims& data) const;
229 
230  template<class Serializer>
231  void serializeOp(Serializer& serializer)
232  {
233  serializer(nMaxNoNodes);
234  serializer(nMaxNoBranches);
235  serializer(nMaxNoBranchesConToNode);
236  }
237 
238 private:
239  int nMaxNoNodes;
240  int nMaxNoBranches;
241  int nMaxNoBranchesConToNode;
242 };
243 
245 public:
247  explicit AquiferDimensions(const Deck& deck);
248 
249  static AquiferDimensions serializeObject();
250 
251  int maxAnalyticAquifers() const
252  {
253  return this->maxNumAnalyticAquifers;
254  }
255 
256  int maxAnalyticAquiferConnections() const
257  {
258  return this->maxNumAnalyticAquiferConn;
259  }
260 
261  template <class Serializer>
262  void serializeOp(Serializer& serializer)
263  {
264  serializer(this->maxNumAnalyticAquifers);
265  serializer(this->maxNumAnalyticAquiferConn);
266  }
267 
268 private:
269  int maxNumAnalyticAquifers;
270  int maxNumAnalyticAquiferConn;
271 };
272 
273 bool operator==(const AquiferDimensions& lhs, const AquiferDimensions& rhs);
274 
276 {
277 public:
278  EclHysterConfig() = default;
279  explicit EclHysterConfig(const Deck& deck);
280 
281  static EclHysterConfig serializeObject();
282 
286  //void setActive(bool yesno);
287 
291  bool active() const;
292 
299  int pcHysteresisModel() const;
300 
307  int krHysteresisModel() const;
308 
309  bool operator==(const EclHysterConfig& data) const;
310 
311  template<class Serializer>
312  void serializeOp(Serializer& serializer)
313  {
314  serializer(activeHyst);
315  serializer(pcHystMod);
316  serializer(krHystMod);
317  }
318 
319 private:
320  // enable hysteresis at all
321  bool activeHyst { false };
322 
323  // the capillary pressure and the relperm hysteresis models to be used
324  int pcHystMod { 0 };
325  int krHystMod { 0 };
326 };
327 
329 public:
330  enum class ThreePhaseOilKrModel {
331  Default,
332  Stone1,
333  Stone2
334  };
335 
336  enum class KeywordFamily {
337  Family_I, // SGOF, SWOF, SLGOF
338  Family_II, // SGFN, SOF{2,3}, SWFN
339  Undefined,
340  };
341 
342  SatFuncControls();
343  explicit SatFuncControls(const Deck& deck);
344  explicit SatFuncControls(const double tolcritArg,
345  const ThreePhaseOilKrModel model,
346  const KeywordFamily family);
347 
348  static SatFuncControls serializeObject();
349 
350  double minimumRelpermMobilityThreshold() const
351  {
352  return this->tolcrit;
353  }
354 
355  ThreePhaseOilKrModel krModel() const
356  {
357  return this->krmodel;
358  }
359 
360  KeywordFamily family() const
361  {
362  return this->satfunc_family;
363  }
364 
365  bool operator==(const SatFuncControls& rhs) const;
366 
367  template<class Serializer>
368  void serializeOp(Serializer& serializer)
369  {
370  serializer(tolcrit);
371  serializer(krmodel);
372  serializer(satfunc_family);
373  }
374 
375 private:
376  double tolcrit;
377  ThreePhaseOilKrModel krmodel = ThreePhaseOilKrModel::Default;
378  KeywordFamily satfunc_family = KeywordFamily::Undefined;
379 };
380 
381 
382 class Nupcol {
383 public:
384  Nupcol();
385  explicit Nupcol(int min_value);
386  void update(int value);
387  int value() const;
388 
389  static Nupcol serializeObject();
390  bool operator==(const Nupcol& data) const;
391 
392  template<class Serializer>
393  void serializeOp(Serializer& serializer) {
394  serializer(this->nupcol_value);
395  serializer(this->min_nupcol);
396  }
397 
398 private:
399  int min_nupcol;
400  int nupcol_value;
401 };
402 
403 
404 class Tracers {
405 public:
406 
407  Tracers() = default;
408 
409  explicit Tracers(const Deck& );
410  int water_tracers() const;
411 
412  template<class Serializer>
413  void serializeOp(Serializer& serializer) {
414  serializer(this->m_oil_tracers);
415  serializer(this->m_water_tracers);
416  serializer(this->m_gas_tracers);
417  serializer(this->m_env_tracers);
418  serializer(this->diffusion_control);
419  serializer(this->max_iter);
420  serializer(this->min_iter);
421  }
422 
423  static Tracers serializeObject();
424  bool operator==(const Tracers& data) const;
425 
426 private:
427  int m_oil_tracers;
428  int m_water_tracers;
429  int m_gas_tracers;
430  int m_env_tracers;
431  bool diffusion_control;
432  int max_iter;
433  int min_iter;
434  // The TRACERS keyword has some additional options which seem quite arcane,
435  // for now not included here.
436 };
437 
438 
439 class Runspec {
440 public:
441  Runspec() = default;
442  explicit Runspec( const Deck& );
443 
444  static Runspec serializeObject();
445 
446  std::time_t start_time() const noexcept;
447  const UDQParams& udqParams() const noexcept;
448  const Phases& phases() const noexcept;
449  const Tabdims& tabdims() const noexcept;
450  const Regdims& regdims() const noexcept;
451  const EndpointScaling& endpointScaling() const noexcept;
452  const Welldims& wellDimensions() const noexcept;
453  const WellSegmentDims& wellSegmentDimensions() const noexcept;
454  const NetworkDims& networkDimensions() const noexcept;
455  const AquiferDimensions& aquiferDimensions() const noexcept;
456  int eclPhaseMask( ) const noexcept;
457  const EclHysterConfig& hysterPar() const noexcept;
458  const Actdims& actdims() const noexcept;
459  const SatFuncControls& saturationFunctionControls() const noexcept;
460  const Nupcol& nupcol() const noexcept;
461  const Tracers& tracers() const;
462  bool co2Storage() const noexcept;
463  bool micp() const noexcept;
464 
465  bool operator==(const Runspec& data) const;
466  static bool rst_cmp(const Runspec& full_state, const Runspec& rst_state);
467 
468  template<class Serializer>
469  void serializeOp(Serializer& serializer)
470  {
471  serializer(this->m_start_time);
472  active_phases.serializeOp(serializer);
473  m_tabdims.serializeOp(serializer);
474  m_regdims.serializeOp(serializer);
475  endscale.serializeOp(serializer);
476  welldims.serializeOp(serializer);
477  wsegdims.serializeOp(serializer);
478  netwrkdims.serializeOp(serializer);
479  aquiferdims.serializeOp(serializer);
480  udq_params.serializeOp(serializer);
481  hystpar.serializeOp(serializer);
482  m_actdims.serializeOp(serializer);
483  m_sfuncctrl.serializeOp(serializer);
484  m_nupcol.serializeOp(serializer);
485  serializer(m_co2storage);
486  serializer(m_micp);
487  }
488 
489 private:
490  std::time_t m_start_time;
491  Phases active_phases;
492  Tabdims m_tabdims;
493  Regdims m_regdims;
494  EndpointScaling endscale;
495  Welldims welldims;
496  WellSegmentDims wsegdims;
497  NetworkDims netwrkdims;
498  AquiferDimensions aquiferdims;
499  UDQParams udq_params;
500  EclHysterConfig hystpar;
501  Actdims m_actdims;
502  SatFuncControls m_sfuncctrl;
503  Nupcol m_nupcol;
504  Tracers m_tracers;
505  bool m_co2storage;
506  bool m_micp;
507 };
508 
509 
510 }
511 
512 #endif // OPM_RUNSPEC_HPP
Definition: Actdims.hpp:30
Definition: Runspec.hpp:244
Definition: Deck.hpp:63
Definition: Runspec.hpp:276
int pcHysteresisModel() const
Return the type of the hysteresis model which is used for capillary pressure.
bool active() const
Specify whether hysteresis is enabled or not.
int krHysteresisModel() const
Return the type of the hysteresis model which is used for relative permeability.
Definition: EndpointScaling.hpp:28
Definition: Runspec.hpp:204
Definition: Runspec.hpp:382
Definition: Runspec.hpp:57
Definition: Regdims.hpp:36
Definition: Runspec.hpp:439
Definition: Runspec.hpp:328
Definition: Serializer.hpp:38
Definition: Tabdims.hpp:36
Definition: Runspec.hpp:404
Definition: UDQParams.hpp:31
Definition: Runspec.hpp:165
Definition: Runspec.hpp:87
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition: Exceptions.hpp:29