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