My Project
TableManager.hpp
1 /*
2  Copyright 2015 Statoil ASA.
3  Copyright 2018 IRIS
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 #ifndef OPM_TABLE_MANAGER_HPP
22 #define OPM_TABLE_MANAGER_HPP
23 
24 #include <cassert>
25 #include <optional>
26 #include <set>
27 
28 #include <opm/input/eclipse/EclipseState/Tables/DenT.hpp>
29 #include <opm/input/eclipse/EclipseState/Tables/JouleThomson.hpp>
30 #include <opm/input/eclipse/EclipseState/Tables/PvtgTable.hpp>
31 #include <opm/input/eclipse/EclipseState/Tables/PvtgwTable.hpp>
32 #include <opm/input/eclipse/EclipseState/Tables/PvtgwoTable.hpp>
33 #include <opm/input/eclipse/EclipseState/Tables/PvtoTable.hpp>
34 #include <opm/input/eclipse/EclipseState/Tables/PvtsolTable.hpp>
35 #include <opm/input/eclipse/EclipseState/Tables/RocktabTable.hpp>
36 #include <opm/input/eclipse/EclipseState/Tables/Rock2dTable.hpp>
37 #include <opm/input/eclipse/EclipseState/Tables/Rock2dtrTable.hpp>
38 #include <opm/input/eclipse/EclipseState/Tables/PlyshlogTable.hpp>
39 #include <opm/input/eclipse/EclipseState/Tables/PvtwsaltTable.hpp>
40 #include <opm/input/eclipse/EclipseState/Tables/RwgsaltTable.hpp>
41 #include <opm/input/eclipse/EclipseState/Tables/BrineDensityTable.hpp>
42 #include <opm/input/eclipse/EclipseState/Tables/SolventDensityTable.hpp>
43 #include <opm/input/eclipse/EclipseState/Tables/StandardCond.hpp>
44 #include <opm/input/eclipse/EclipseState/Tables/FlatTable.hpp>
45 #include <opm/input/eclipse/EclipseState/Tables/SorwmisTable.hpp>
46 #include <opm/input/eclipse/EclipseState/Tables/SgcwmisTable.hpp>
47 #include <opm/input/eclipse/EclipseState/Tables/MiscTable.hpp>
48 #include <opm/input/eclipse/EclipseState/Tables/PmiscTable.hpp>
49 #include <opm/input/eclipse/EclipseState/Tables/MsfnTable.hpp>
50 #include <opm/input/eclipse/EclipseState/Tables/JFunc.hpp>
51 #include <opm/input/eclipse/EclipseState/Tables/Tabdims.hpp>
52 #include <opm/input/eclipse/EclipseState/Tables/TableContainer.hpp>
53 #include <opm/input/eclipse/EclipseState/Tables/Aqudims.hpp>
54 #include <opm/input/eclipse/EclipseState/Tables/PlymwinjTable.hpp>
55 #include <opm/input/eclipse/EclipseState/Tables/SkprwatTable.hpp>
56 #include <opm/input/eclipse/EclipseState/Tables/SkprpolyTable.hpp>
57 #include <opm/input/eclipse/EclipseState/Tables/Eqldims.hpp>
58 #include <opm/input/eclipse/EclipseState/Tables/Regdims.hpp>
59 #include <opm/input/eclipse/EclipseState/Tables/TLMixpar.hpp>
60 
61 namespace Opm {
62 
63  class Deck;
64 
65  class TableManager {
66  public:
67  explicit TableManager( const Deck& deck );
68  TableManager() = default;
69 
70  static TableManager serializationTestObject();
71 
72  const TableContainer& getTables( const std::string& tableName ) const;
73  const TableContainer& operator[](const std::string& tableName) const;
74  bool hasTables( const std::string& tableName ) const;
75 
76  const Tabdims& getTabdims() const;
77  const Eqldims& getEqldims() const;
78  const Aqudims& getAqudims() const;
79  const Regdims& getRegdims() const;
80  const TLMixpar& getTLMixpar() const;
81  /*
82  WIll return max{ Tabdims::NTFIP , Regdims::NTFIP }.
83  */
84  size_t numFIPRegions() const;
85 
86  const TableContainer& getSwofTables() const;
87  const TableContainer& getSgwfnTables() const;
88  const TableContainer& getSof2Tables() const;
89  const TableContainer& getSof3Tables() const;
90  const TableContainer& getSgofTables() const;
91  const TableContainer& getSlgofTables() const;
92  const TableContainer& getSwfnTables() const;
93  const TableContainer& getSgfnTables() const;
94  const TableContainer& getSsfnTables() const;
95  const TableContainer& getRsvdTables() const;
96  const TableContainer& getRvvdTables() const;
97  const TableContainer& getRvwvdTables() const;
98  const TableContainer& getPbvdTables() const;
99  const TableContainer& getPdvdTables() const;
100  const TableContainer& getSaltvdTables() const;
101  const TableContainer& getSaltpvdTables() const;
102  const TableContainer& getSaltsolTables() const;
103  const TableContainer& getPermfactTables() const;
104  const TableContainer& getEnkrvdTables() const;
105  const TableContainer& getEnptvdTables() const;
106  const TableContainer& getImkrvdTables() const;
107  const TableContainer& getImptvdTables() const;
108  const TableContainer& getPvdgTables() const;
109  const TableContainer& getPvdoTables() const;
110  const TableContainer& getPvdsTables() const;
111  const TableContainer& getSpecheatTables() const;
112  const TableContainer& getSpecrockTables() const;
113  const TableContainer& getWatvisctTables() const;
114  const TableContainer& getOilvisctTables() const;
115  const TableContainer& getGasvisctTables() const;
116  const TableContainer& getRtempvdTables() const;
117  const TableContainer& getRocktabTables() const;
118  const TableContainer& getPlyadsTables() const;
119  const TableContainer& getPlyviscTables() const;
120  const TableContainer& getPlydhflfTables() const;
121  const TableContainer& getPlymaxTables() const;
122  const TableContainer& getPlyrockTables() const;
123  const TableContainer& getPlyshlogTables() const;
124  const TableContainer& getAqutabTables() const;
125  const TableContainer& getFoamadsTables() const;
126  const TableContainer& getFoammobTables() const;
127 
128  const TableContainer& getSorwmisTables() const;
129  const TableContainer& getSgcwmisTables() const;
130  const TableContainer& getMiscTables() const;
131  const TableContainer& getPmiscTables() const;
132  const TableContainer& getMsfnTables() const;
133  const TableContainer& getTlpmixpaTables() const;
134 
135  const JFunc& getJFunc() const;
136 
137  const std::vector<PvtgTable>& getPvtgTables() const;
138  const std::vector<PvtgwTable>& getPvtgwTables() const;
139  const std::vector<PvtgwoTable>& getPvtgwoTables() const;
140  const std::vector<PvtoTable>& getPvtoTables() const;
141  const std::vector<PvtsolTable>& getPvtsolTables() const;
142  const std::vector<Rock2dTable>& getRock2dTables() const;
143  const std::vector<Rock2dtrTable>& getRock2dtrTables() const;
144  const TableContainer& getRockwnodTables() const;
145  const TableContainer& getOverburdTables() const;
146 
147  const DenT& WatDenT() const;
148  const DenT& GasDenT() const;
149  const DenT& OilDenT() const;
150  const JouleThomson& WatJT() const;
151  const JouleThomson& GasJT() const;
152  const JouleThomson& OilJT() const;
153  const StandardCond& stCond() const;
154  std::size_t gas_comp_index() const;
155  const PvtwTable& getPvtwTable() const;
156  const std::vector<PvtwsaltTable>& getPvtwSaltTables() const;
157  const std::vector<RwgsaltTable>& getRwgSaltTables() const;
158  const std::vector<BrineDensityTable>& getBrineDensityTables() const;
159  const std::vector<SolventDensityTable>& getSolventDensityTables() const;
160 
161  const PvcdoTable& getPvcdoTable() const;
162  const DensityTable& getDensityTable() const;
163  const DiffCoeffTable& getDiffusionCoefficientTable() const;
164  const PlyvmhTable& getPlyvmhTable() const;
165  const RockTable& getRockTable() const;
166  const ViscrefTable& getViscrefTable() const;
167  const PlmixparTable& getPlmixparTable() const;
168  const ShrateTable& getShrateTable() const;
169  const Stone1exTable& getStone1exTable() const;
170  const WatdentTable& getWatdentTable() const;
171  const SgofletTable& getSgofletTable() const;
172  const SwofletTable& getSwofletTable() const;
173  const std::map<int, PlymwinjTable>& getPlymwinjTables() const;
174  const std::map<int, SkprwatTable>& getSkprwatTables() const;
175  const std::map<int, SkprpolyTable>& getSkprpolyTables() const;
176  const std::map<std::string, TableContainer>& getSimpleTables() const;
177 
179  bool useImptvd() const;
180 
182  bool useEnptvd() const;
183 
185  bool useEqlnum() const;
186 
188  bool useShrate() const;
189 
191  bool useJFunc() const;
192 
193  double rtemp() const;
194 
195  double salinity() const;
196 
197  bool operator==(const TableManager& data) const;
198 
199  template<class Serializer>
200  void serializeOp(Serializer& serializer)
201  {
202  auto simpleTables = m_simpleTables;
203  auto split = splitSimpleTable(simpleTables);
204  serializer(simpleTables);
205  serializer(split.plyshMax);
206  serializer(split.plyshMap);
207  serializer(split.rockMax);
208  serializer(split.rockMap);
209  serializer(m_pvtgTables);
210  serializer(m_pvtgwTables);
211  serializer(m_pvtgwoTables);
212  serializer(m_pvtoTables);
213  serializer(m_pvtsolTables);
214  serializer(m_rock2dTables);
215  serializer(m_rock2dtrTables);
216  serializer(m_pvtwTable);
217  serializer(m_pvcdoTable);
218  serializer(m_densityTable);
219  serializer(m_diffCoeffTable);
220  serializer(m_plyvmhTable);
221  serializer(m_rockTable);
222  serializer(m_plmixparTable);
223  serializer(m_shrateTable);
224  serializer(m_stone1exTable);
225  serializer(m_viscrefTable);
226  serializer(m_watdentTable);
227  serializer(m_sgofletTable);
228  serializer(m_swofletTable);
229  serializer(m_pvtwsaltTables);
230  serializer(m_rwgsaltTables);
231  serializer(m_bdensityTables);
232  serializer(m_sdensityTables);
233  serializer(m_plymwinjTables);
234  serializer(m_skprwatTables);
235  serializer(m_skprpolyTables);
236  serializer(m_tabdims);
237  serializer(m_regdims);
238  serializer(m_eqldims);
239  serializer(m_aqudims);
240  serializer(hasImptvd);
241  serializer(hasEnptvd);
242  serializer(hasEqlnum);
243  serializer(hasShrate);
244  serializer(jfunc);
245  serializer(oilDenT);
246  serializer(gasDenT);
247  serializer(watDenT);
248  serializer(oilJT);
249  serializer(gasJT);
250  serializer(watJT);
251  serializer(stcond);
252  serializer(m_gas_comp_index);
253  serializer(m_rtemp);
254  serializer(m_salinity);
255  serializer(m_tlmixpar);
256  if (!serializer.isSerializing()) {
257  m_simpleTables = simpleTables;
258  if (split.plyshMax > 0) {
259  TableContainer container(split.plyshMax);
260  for (const auto& it : split.plyshMap) {
261  container.addTable(it.first, it.second);
262  }
263  m_simpleTables.insert(std::make_pair("PLYSHLOG", container));
264  }
265  if (split.rockMax > 0) {
266  TableContainer container(split.rockMax);
267  for (const auto& it : split.rockMap) {
268  container.addTable(it.first, it.second);
269  }
270  m_simpleTables.insert(std::make_pair("ROCKTAB", container));
271  }
272  }
273  }
274 
275  private:
276  TableContainer& forceGetTables( const std::string& tableName , size_t numTables);
277 
278  void complainAboutAmbiguousKeyword(const Deck& deck, const std::string& keywordName);
279 
280  void addTables( const std::string& tableName , size_t numTables);
281  void initSimpleTables(const Deck& deck);
282  void initRTempTables(const Deck& deck);
283  void initDims(const Deck& deck);
284  void initRocktabTables(const Deck& deck);
285  void initGasvisctTables(const Deck& deck);
286 
287  void initPlymaxTables(const Deck& deck);
288  void initPlyrockTables(const Deck& deck);
289  void initPlyshlogTables(const Deck& deck);
290 
291  void initPlymwinjTables(const Deck& deck);
292  void initSkprwatTables(const Deck& deck);
293  void initSkprpolyTables(const Deck& deck);
294 
295  //void initRockTables(const Deck& deck, const std::string& keywordName);
296 
297  template <class TableType>
298  void initRockTables(const Deck& deck, const std::string& keywordName, std::vector<TableType>& rocktable );
299 
300  template <class TableType>
301  void initPvtwsaltTables(const Deck& deck, std::vector<TableType>& pvtwtables );
302 
303  template <class TableType>
304  void initRwgsaltTables(const Deck& deck, std::vector<TableType>& rwgtables );
305 
306  template <class TableType>
307  void initBrineTables(const Deck& deck, std::vector<TableType>& brinetables );
308 
309  void initSolventTables(const Deck& deck, std::vector<SolventDensityTable>& solventtables);
310 
314  template <class TableType>
315  void initSimpleTableContainerWithJFunc(const Deck& deck,
316  const std::string& keywordName,
317  const std::string& tableName,
318  size_t numTables);
319 
320  template <class TableType>
321  void initSimpleTableContainer(const Deck& deck,
322  const std::string& keywordName,
323  const std::string& tableName,
324  size_t numTables);
325 
326  template <class TableType>
327  void initSimpleTableContainer(const Deck& deck,
328  const std::string& keywordName,
329  size_t numTables);
330 
331  template <class TableType>
332  void initSimpleTableContainerWithJFunc(const Deck& deck,
333  const std::string& keywordName,
334  size_t numTables);
335 
336  template <class TableType>
337  void initSimpleTable(const Deck& deck,
338  const std::string& keywordName,
339  std::vector<TableType>& tableVector);
340 
341  template <class TableType>
342  void initFullTables(const Deck& deck,
343  const std::string& keywordName,
344  std::vector<TableType>& tableVector);
345 
346  void checkPVTOMonotonicity(const Deck& deck) const;
347 
348  void logPVTOMonotonicityFailure(const Deck& deck,
349  const std::size_t tableID,
350  const std::vector<PvtoTable::FlippedFVF>& flipped_Bo) const;
351 
352  std::map<std::string , TableContainer> m_simpleTables;
353  std::vector<PvtgTable> m_pvtgTables;
354  std::vector<PvtgwTable> m_pvtgwTables;
355  std::vector<PvtgwoTable> m_pvtgwoTables;
356  std::vector<PvtoTable> m_pvtoTables;
357  std::vector<PvtsolTable> m_pvtsolTables;
358  std::vector<Rock2dTable> m_rock2dTables;
359  std::vector<Rock2dtrTable> m_rock2dtrTables;
360  PvtwTable m_pvtwTable;
361  PvcdoTable m_pvcdoTable;
362  DensityTable m_densityTable;
363  DiffCoeffTable m_diffCoeffTable;
364  PlyvmhTable m_plyvmhTable;
365  RockTable m_rockTable;
366  PlmixparTable m_plmixparTable;
367  ShrateTable m_shrateTable;
368  Stone1exTable m_stone1exTable;
369  ViscrefTable m_viscrefTable;
370  WatdentTable m_watdentTable;
371  SgofletTable m_sgofletTable;
372  SwofletTable m_swofletTable;
373  std::vector<PvtwsaltTable> m_pvtwsaltTables;
374  std::vector<RwgsaltTable> m_rwgsaltTables;
375  std::vector<BrineDensityTable> m_bdensityTables;
376  std::vector<SolventDensityTable> m_sdensityTables;
377  std::map<int, PlymwinjTable> m_plymwinjTables;
378  std::map<int, SkprwatTable> m_skprwatTables;
379  std::map<int, SkprpolyTable> m_skprpolyTables;
380 
381  Tabdims m_tabdims;
382  Regdims m_regdims;
383  Eqldims m_eqldims;
384  Aqudims m_aqudims;
385  TLMixpar m_tlmixpar;
386 
387  bool hasImptvd = false;// if deck has keyword IMPTVD
388  bool hasEnptvd = false;// if deck has keyword ENPTVD
389  bool hasEqlnum = false;// if deck has keyword EQLNUM
390  bool hasShrate = false;// if deck has keyword SHRATE
391  std::optional<JFunc> jfunc;
392 
393  DenT oilDenT;
394  DenT gasDenT;
395  DenT watDenT;
396  JouleThomson oilJT;
397  JouleThomson gasJT;
398  JouleThomson watJT;
399  StandardCond stcond;
400  std::size_t m_gas_comp_index = 77;
401  double m_rtemp {288.7056}; // 60 Fahrenheit in Kelvin
402  double m_salinity {0.0};
403 
404  struct SplitSimpleTables {
405  size_t plyshMax = 0;
406  size_t rockMax = 0;
407  std::map<size_t, std::shared_ptr<PlyshlogTable>> plyshMap;
408  std::map<size_t, std::shared_ptr<RocktabTable>> rockMap;
409  };
410 
411  SplitSimpleTables splitSimpleTable(std::map<std::string,TableContainer>& simpleTables);
412  };
413 }
414 
415 
416 #endif
417 
Definition: Aqudims.hpp:34
Definition: Deck.hpp:63
Definition: DenT.hpp:30
Definition: Eqldims.hpp:32
Definition: JFunc.hpp:28
Definition: JouleThomson.hpp:30
Definition: Regdims.hpp:36
Class for (de-)serializing.
Definition: Serializer.hpp:75
bool isSerializing() const
Returns true if we are currently doing a serialization operation.
Definition: Serializer.hpp:174
Definition: TLMixpar.hpp:55
Definition: Tabdims.hpp:36
Definition: TableContainer.hpp:31
Definition: TableManager.hpp:65
bool useEqlnum() const
deck has keyword "EQLNUM" — Equilibriation region numbers
bool useJFunc() const
deck has keyword "JFUNC" — Use Leverett's J Function for capillary pressure
bool useEnptvd() const
deck has keyword "ENPTVD" — Saturation end-point versus depth tables
bool useImptvd() const
deck has keyword "IMPTVD" — Imbition end-point versus depth tables
bool useShrate() const
deck has keyword "SHRATE"
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition: Exceptions.hpp:29
Definition: FlatTable.hpp:130
Definition: FlatTable.hpp:185
Definition: FlatTable.hpp:321
Definition: FlatTable.hpp:355
Definition: FlatTable.hpp:296
Definition: FlatTable.hpp:223
Definition: FlatTable.hpp:259
Definition: FlatTable.hpp:575
Definition: FlatTable.hpp:380
Definition: StandardCond.hpp:24
Definition: FlatTable.hpp:405
Definition: FlatTable.hpp:565
Definition: FlatTable.hpp:461
Definition: FlatTable.hpp:492