My Project
EclipseGrid.hpp
1 /*
2  Copyright 2014 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
7  it under the terms of the GNU General Public License as published by
8  the Free Software Foundation, either version 3 of the License, or
9  (at your option) any later version.
10 
11  OPM is distributed in the hope that it will be useful,
12  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  GNU General Public License for more details.
15 
16  You should have received a copy of the GNU General Public License
17  along with OPM. If not, see <http://www.gnu.org/licenses/>.
18 */
19 
20 #ifndef OPM_PARSER_ECLIPSE_GRID_HPP
21 #define OPM_PARSER_ECLIPSE_GRID_HPP
22 
23 #include <opm/input/eclipse/EclipseState/Grid/GridDims.hpp>
24 #include <opm/input/eclipse/EclipseState/Grid/MapAxes.hpp>
25 #include <opm/input/eclipse/EclipseState/Grid/MinpvMode.hpp>
26 #include <opm/input/eclipse/EclipseState/Grid/PinchMode.hpp>
27 
28 #include <array>
29 #include <memory>
30 #include <optional>
31 #include <stdexcept>
32 #include <unordered_set>
33 #include <vector>
34 
35 namespace Opm {
36 
37  class Deck;
38  namespace EclIO { class EclFile; }
39  struct NNCdata;
40  class UnitSystem;
41  class ZcornMapper;
42 
54  class EclipseGrid : public GridDims {
55  public:
56  EclipseGrid() = default;
57  explicit EclipseGrid(const std::string& filename);
58 
59  /*
60  These constructors will make a copy of the src grid, with
61  zcorn and or actnum have been adjustments.
62  */
63  EclipseGrid(const EclipseGrid& src) = default;
64  EclipseGrid(const EclipseGrid& src, const std::vector<int>& actnum);
65  EclipseGrid(const EclipseGrid& src, const double* zcorn, const std::vector<int>& actnum);
66 
67  EclipseGrid(size_t nx, size_t ny, size_t nz,
68  double dx = 1.0, double dy = 1.0, double dz = 1.0);
69 
70  EclipseGrid(const std::array<int, 3>& dims ,
71  const std::vector<double>& coord ,
72  const std::vector<double>& zcorn ,
73  const int * actnum = nullptr);
74 
75 
78  EclipseGrid(const Deck& deck, const int * actnum = nullptr);
79 
80  static bool hasGDFILE(const Deck& deck);
81  static bool hasCylindricalKeywords(const Deck& deck);
82  static bool hasCornerPointKeywords(const Deck&);
83  static bool hasCartesianKeywords(const Deck&);
84  size_t getNumActive( ) const;
85  bool allActive() const;
86 
87  size_t activeIndex(size_t i, size_t j, size_t k) const;
88  size_t activeIndex(size_t globalIndex) const;
89 
90  size_t getActiveIndex(size_t i, size_t j, size_t k) const {
91  return activeIndex(i, j, k);
92  }
93 
94  void save(const std::string& filename, bool formatted, const std::vector<Opm::NNCdata>& nnc, const Opm::UnitSystem& units) const;
95  /*
96  Observe that the there is a getGlobalIndex(i,j,k)
97  implementation in the base class. This method - translating
98  from an active index to a global index must be implemented
99  in the current class.
100  */
101  size_t getGlobalIndex(size_t active_index) const;
102  size_t getGlobalIndex(size_t i, size_t j, size_t k) const;
103 
104  /*
105  For RADIAL grids you can *optionally* use the keyword
106  'CIRCLE' to denote that period boundary conditions should be
107  applied in the 'THETA' direction; this will only apply if
108  the theta keywords entered sum up to exactly 360 degrees!
109  */
110 
111  bool circle( ) const;
112  bool isPinchActive( ) const;
113  double getPinchThresholdThickness( ) const;
114  PinchMode::ModeEnum getPinchOption( ) const;
115  PinchMode::ModeEnum getMultzOption( ) const;
116  PinchMode::ModeEnum getPinchGapMode( ) const;
117 
118  MinpvMode::ModeEnum getMinpvMode() const;
119  const std::vector<double>& getMinpvVector( ) const;
120 
121  /*
122  Will return a vector of nactive elements. The method will
123  behave differently depending on the lenght of the
124  input_vector:
125 
126  nx*ny*nz: only the values corresponding to active cells
127  are copied out.
128 
129  nactive: The input vector is copied straight out again.
130 
131  ??? : Exception.
132  */
133 
134  template<typename T>
135  std::vector<T> compressedVector(const std::vector<T>& input_vector) const {
136  if( input_vector.size() == this->getNumActive() ) {
137  return input_vector;
138  }
139 
140  if (input_vector.size() != getCartesianSize())
141  throw std::invalid_argument("Input vector must have full size");
142 
143  {
144  std::vector<T> compressed_vector( this->getNumActive() );
145  const auto& active_map = this->getActiveMap( );
146 
147  for (size_t i = 0; i < this->getNumActive(); ++i)
148  compressed_vector[i] = input_vector[ active_map[i] ];
149 
150  return compressed_vector;
151  }
152  }
153 
154 
157  const std::vector<int>& getActiveMap() const;
158  std::array<double, 3> getCellCenter(size_t i,size_t j, size_t k) const;
159  std::array<double, 3> getCellCenter(size_t globalIndex) const;
160  std::array<double, 3> getCornerPos(size_t i,size_t j, size_t k, size_t corner_index) const;
161  const std::vector<double>& activeVolume() const;
162  double getCellVolume(size_t globalIndex) const;
163  double getCellVolume(size_t i , size_t j , size_t k) const;
164  double getCellThickness(size_t globalIndex) const;
165  double getCellThickness(size_t i , size_t j , size_t k) const;
166  std::array<double, 3> getCellDims(size_t i,size_t j, size_t k) const;
167  std::array<double, 3> getCellDims(size_t globalIndex) const;
168  bool cellActive( size_t globalIndex ) const;
169  bool cellActive( size_t i , size_t j, size_t k ) const;
170 
171  std::array<double, 3> getCellDimensions(size_t i, size_t j, size_t k) const {
172  return getCellDims(i, j, k);
173  }
174 
175  bool isCellActive(size_t i, size_t j, size_t k) const {
176  return cellActive(i, j, k);
177  }
178 
184  bool isValidCellGeomtry(const std::size_t globalIndex,
185  const UnitSystem& usys) const;
186 
187  double getCellDepth(size_t i,size_t j, size_t k) const;
188  double getCellDepth(size_t globalIndex) const;
189  ZcornMapper zcornMapper() const;
190 
191  const std::vector<double>& getCOORD() const;
192  const std::vector<double>& getZCORN() const;
193  const std::vector<int>& getACTNUM( ) const;
194 
195  /*
196  The fixupZCORN method is run as part of constructiong the grid. This will adjust the
197  z-coordinates to ensure that cells do not overlap. The return value is the number of
198  points which have been adjusted. The number of zcorn nodes that has ben fixed is
199  stored in private member zcorn_fixed.
200  */
201 
202  size_t fixupZCORN();
203  size_t getZcornFixed() { return zcorn_fixed; };
204 
205  // resetACTNUM with no arguments will make all cells in the grid active.
206 
207  void resetACTNUM();
208  void resetACTNUM( const std::vector<int>& actnum);
209 
210  bool equal(const EclipseGrid& other) const;
211 
212  private:
213  std::vector<double> m_minpvVector;
214  MinpvMode::ModeEnum m_minpvMode;
215  std::optional<double> m_pinch;
216  PinchMode::ModeEnum m_pinchoutMode;
217  PinchMode::ModeEnum m_multzMode;
218  PinchMode::ModeEnum m_pinchGapMode;
219 
220  mutable std::optional<std::vector<double>> active_volume;
221 
222  bool m_circle = false;
223 
224  size_t zcorn_fixed = 0;
225  bool m_useActnumFromGdfile = false;
226 
227  // Input grid data.
228  std::vector<double> m_zcorn;
229  std::vector<double> m_coord;
230  std::vector<int> m_actnum;
231  std::optional<MapAxes> m_mapaxes;
232 
233  // Mapping to/from active cells.
234  int m_nactive;
235  std::vector<int> m_active_to_global;
236  std::vector<int> m_global_to_active;
237  // Numerical aquifer cells, needs to be active
238  std::unordered_set<size_t> m_aquifer_cells;
239 
240  // Radial grids need this for volume calculations.
241  std::optional<std::vector<double>> m_thetav;
242  std::optional<std::vector<double>> m_rv;
243 
244  void updateNumericalAquiferCells(const Deck&);
245 
246  void initGridFromEGridFile(Opm::EclIO::EclFile& egridfile, std::string fileName);
247  void resetACTNUM( const int* actnum);
248 
249  void initBinaryGrid(const Deck& deck);
250 
251  void initCornerPointGrid(const std::vector<double>& coord ,
252  const std::vector<double>& zcorn ,
253  const int * actnum);
254 
255  bool keywInputBeforeGdfile(const Deck& deck, const std::string& keyword) const;
256 
257  void initCylindricalGrid(const Deck&);
258  void initSpiderwebGrid(const Deck&);
259  void initSpiderwebOrCylindricalGrid(const Deck&, const bool);
260  void initCartesianGrid(const Deck&);
261  void initDTOPSGrid(const Deck&);
262  void initDVDEPTHZGrid(const Deck&);
263  void initGrid(const Deck&);
264  void initCornerPointGrid(const Deck&);
265  void assertCornerPointKeywords(const Deck&);
266 
267  static bool hasDVDEPTHZKeywords(const Deck&);
268  static bool hasDTOPSKeywords(const Deck&);
269  static void assertVectorSize(const std::vector<double>& vector, size_t expectedSize, const std::string& msg);
270 
271  static std::vector<double> createTOPSVector(const std::array<int, 3>& dims, const std::vector<double>& DZ, const Deck&);
272  static std::vector<double> createDVector(const std::array<int, 3>& dims, std::size_t dim, const std::string& DKey, const std::string& DVKey, const Deck&);
273  static void scatterDim(const std::array<int, 3>& dims , size_t dim , const std::vector<double>& DV , std::vector<double>& D);
274 
275 
276  std::vector<double> makeCoordDxDyDzTops(const std::vector<double>& dx, const std::vector<double>& dy, const std::vector<double>& dz, const std::vector<double>& tops) const;
277  std::vector<double> makeZcornDzTops(const std::vector<double>& dz, const std::vector<double>& tops) const ;
278  std::vector<double> makeZcornDzvDepthz(const std::vector<double>& dzv, const std::vector<double>& depthz) const;
279  std::vector<double> makeCoordDxvDyvDzvDepthz(const std::vector<double>& dxv, const std::vector<double>& dyv, const std::vector<double>& dzv, const std::vector<double>& depthz) const;
280 
281  void getCellCorners(const std::array<int, 3>& ijk, const std::array<int, 3>& dims, std::array<double,8>& X, std::array<double,8>& Y, std::array<double,8>& Z) const;
282  void getCellCorners(const std::size_t globalIndex,
283  std::array<double,8>& X,
284  std::array<double,8>& Y,
285  std::array<double,8>& Z) const;
286 
287  };
288 
289  class CoordMapper {
290  public:
291  CoordMapper(size_t nx, size_t ny);
292  size_t size() const;
293 
294 
295  /*
296  dim = 0,1,2 for x, y and z coordinate respectively.
297  layer = 0,1 for k=0 and k=nz layers respectively.
298  */
299 
300  size_t index(size_t i, size_t j, size_t dim, size_t layer) const;
301  private:
302  size_t nx;
303  size_t ny;
304  };
305 
306 
307 
308  class ZcornMapper {
309  public:
310  ZcornMapper(size_t nx, size_t ny, size_t nz);
311  size_t index(size_t i, size_t j, size_t k, int c) const;
312  size_t index(size_t g, int c) const;
313  size_t size() const;
314 
315  /*
316  The fixupZCORN method will take a zcorn vector as input and
317  run through it. If the following situation is detected:
318 
319  /| /|
320  / | / |
321  / | / |
322  / | / |
323  / | ==> / |
324  / | / |
325  ---/------x /---------x
326  | /
327  |/
328 
329  */
330  size_t fixupZCORN( std::vector<double>& zcorn);
331  bool validZCORN( const std::vector<double>& zcorn) const;
332  private:
333  std::array<size_t,3> dims;
334  std::array<size_t,3> stride;
335  std::array<size_t,8> cell_shift;
336  };
337 }
338 
339 #endif // OPM_PARSER_ECLIPSE_GRID_HPP
Definition: EclipseGrid.hpp:289
Definition: Deck.hpp:63
Definition: EclFile.hpp:35
About cell information and dimension: The actual grid information is held in a pointer to an ERT ecl_...
Definition: EclipseGrid.hpp:54
EclipseGrid(const Deck &deck, const int *actnum=nullptr)
EclipseGrid ignores ACTNUM in Deck, and therefore needs ACTNUM explicitly.
bool isValidCellGeomtry(const std::size_t globalIndex, const UnitSystem &usys) const
Whether or not given cell has a valid cell geometry.
const std::vector< int > & getActiveMap() const
Will return a vector a length num_active; where the value of each element is the corresponding global...
Definition: GridDims.hpp:32
Definition: UnitSystem.hpp:34
Definition: EclipseGrid.hpp:308
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition: Exceptions.hpp:29