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