My Project
Loading...
Searching...
No Matches
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#include <map>
35
36namespace Opm {
37
38 class Deck;
39 namespace EclIO { class EclFile; }
40 struct NNCdata;
41 class UnitSystem;
42 class ZcornMapper;
43
55 class EclipseGrid : public GridDims {
56 public:
58 explicit EclipseGrid(const std::string& filename);
59
60 /*
61 These constructors will make a copy of the src grid, with
62 zcorn and or actnum have been adjustments.
63 */
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 double top = 0.0);
70 explicit EclipseGrid(const GridDims& gd);
71
72 EclipseGrid(const std::array<int, 3>& dims ,
73 const std::vector<double>& coord ,
74 const std::vector<double>& zcorn ,
75 const int * actnum = nullptr);
76
77
80 explicit EclipseGrid(const Deck& deck, const int * actnum = nullptr);
81
82 static bool hasGDFILE(const Deck& deck);
83 static bool hasRadialKeywords(const Deck& deck);
84 static bool hasSpiderKeywords(const Deck& deck);
85 static bool hasCylindricalKeywords(const Deck& deck);
86 static bool hasCornerPointKeywords(const Deck&);
87 static bool hasCartesianKeywords(const Deck&);
88 size_t getNumActive( ) const;
89 bool allActive() const;
90
91 size_t activeIndex(size_t i, size_t j, size_t k) const;
92 size_t activeIndex(size_t globalIndex) const;
93
94 size_t getActiveIndex(size_t i, size_t j, size_t k) const {
95 return activeIndex(i, j, k);
96 }
97
98 void save(const std::string& filename, bool formatted, const std::vector<Opm::NNCdata>& nnc, const Opm::UnitSystem& units) const;
99 /*
100 Observe that the there is a getGlobalIndex(i,j,k)
101 implementation in the base class. This method - translating
102 from an active index to a global index must be implemented
103 in the current class.
104 */
105 size_t getGlobalIndex(size_t active_index) const;
106 size_t getGlobalIndex(size_t i, size_t j, size_t k) const;
107
108 /*
109 For RADIAL grids you can *optionally* use the keyword
110 'CIRCLE' to denote that period boundary conditions should be
111 applied in the 'THETA' direction; this will only apply if
112 the theta keywords entered sum up to exactly 360 degrees!
113 */
114
115 bool circle() const;
116 bool isPinchActive() const;
117 double getPinchThresholdThickness() const;
118 PinchMode getPinchOption() const;
119 PinchMode getMultzOption() const;
120 PinchMode getPinchGapMode() const;
121 double getPinchMaxEmptyGap() const;
122
123 MinpvMode getMinpvMode() const;
124 const std::vector<double>& getMinpvVector( ) const;
125
126 /*
127 Will return a vector of nactive elements. The method will
128 behave differently depending on the lenght of the
129 input_vector:
130
131 nx*ny*nz: only the values corresponding to active cells
132 are copied out.
133
134 nactive: The input vector is copied straight out again.
135
136 ??? : Exception.
137 */
138
139 template<typename T>
140 std::vector<T> compressedVector(const std::vector<T>& input_vector) const {
141 if( input_vector.size() == this->getNumActive() ) {
142 return input_vector;
143 }
144
145 if (input_vector.size() != getCartesianSize())
146 throw std::invalid_argument("Input vector must have full size");
147
148 {
149 std::vector<T> compressed_vector( this->getNumActive() );
150 const auto& active_map = this->getActiveMap( );
151
152 for (size_t i = 0; i < this->getNumActive(); ++i)
153 compressed_vector[i] = input_vector[ active_map[i] ];
154
155 return compressed_vector;
156 }
157 }
158
159
162 const std::vector<int>& getActiveMap() const;
163 std::array<double, 3> getCellCenter(size_t i,size_t j, size_t k) const;
164 std::array<double, 3> getCellCenter(size_t globalIndex) const;
165 std::array<double, 3> getCornerPos(size_t i,size_t j, size_t k, size_t corner_index) const;
166 const std::vector<double>& activeVolume() const;
167 double getCellVolume(size_t globalIndex) const;
168 double getCellVolume(size_t i , size_t j , size_t k) const;
169 double getCellThickness(size_t globalIndex) const;
170 double getCellThickness(size_t i , size_t j , size_t k) const;
171 std::array<double, 3> getCellDims(size_t i,size_t j, size_t k) const;
172 std::array<double, 3> getCellDims(size_t globalIndex) const;
173 bool cellActive( size_t globalIndex ) const;
174 bool cellActive( size_t i , size_t j, size_t k ) const;
175
176 std::array<double, 3> getCellDimensions(size_t i, size_t j, size_t k) const {
177 return getCellDims(i, j, k);
178 }
179
180 bool isCellActive(size_t i, size_t j, size_t k) const {
181 return cellActive(i, j, k);
182 }
183
189 bool isValidCellGeomtry(const std::size_t globalIndex,
190 const UnitSystem& usys) const;
191
192 double getCellDepth(size_t i,size_t j, size_t k) const;
193 double getCellDepth(size_t globalIndex) const;
194 ZcornMapper zcornMapper() const;
195
196 const std::vector<double>& getCOORD() const;
197 const std::vector<double>& getZCORN() const;
198 const std::vector<int>& getACTNUM( ) const;
199
200 /*
201 The fixupZCORN method is run as part of constructiong the grid. This will adjust the
202 z-coordinates to ensure that cells do not overlap. The return value is the number of
203 points which have been adjusted. The number of zcorn nodes that has ben fixed is
204 stored in private member zcorn_fixed.
205 */
206
207 size_t fixupZCORN();
208 size_t getZcornFixed() { return zcorn_fixed; };
209
210 // resetACTNUM with no arguments will make all cells in the grid active.
211
212 void resetACTNUM();
213 void resetACTNUM( const std::vector<int>& actnum);
214
216 void setMINPVV(const std::vector<double>& minpvv);
217
218 bool equal(const EclipseGrid& other) const;
219 static bool hasDVDEPTHZKeywords(const Deck&);
220
221 /*
222 For ALugrid we can *only* use the keyword <DXV, DXYV, DZV, DEPTHZ> so to
223 initialize a Regular Cartesian Grid; further we need equidistant mesh
224 spacing in each direction to initialize ALuGrid (mandatory for
225 mesh refinement!).
226 */
227
228 static bool hasEqualDVDEPTHZ(const Deck&);
229 static bool allEqual(const std::vector<double> &v);
230
231 private:
232 std::vector<double> m_minpvVector;
233 MinpvMode m_minpvMode;
234 std::optional<double> m_pinch;
235 PinchMode m_pinchoutMode;
236 PinchMode m_multzMode;
237 PinchMode m_pinchGapMode;
238 double m_pinchMaxEmptyGap;
239
240 mutable std::optional<std::vector<double>> active_volume;
241
242 bool m_circle = false;
243
244 size_t zcorn_fixed = 0;
245 bool m_useActnumFromGdfile = false;
246
247 // Input grid data.
248 mutable std::optional<std::vector<double>> m_input_zcorn;
249 mutable std::optional<std::vector<double>> m_input_coord;
250
251 std::vector<double> m_zcorn;
252 std::vector<double> m_coord;
253
254
255 std::vector<int> m_actnum;
256 std::optional<MapAxes> m_mapaxes;
257
258 // Mapping to/from active cells.
259 int m_nactive {};
260 std::vector<int> m_active_to_global;
261 std::vector<int> m_global_to_active;
262 // Numerical aquifer cells, needs to be active
263 std::unordered_set<size_t> m_aquifer_cells;
264 // Keep track of aquifer cell depths
265 std::map<size_t, double> m_aquifer_cell_depths;
266
267 // Radial grids need this for volume calculations.
268 std::optional<std::vector<double>> m_thetav;
269 std::optional<std::vector<double>> m_rv;
270
271 void updateNumericalAquiferCells(const Deck&);
272 double computeCellGeometricDepth(size_t globalIndex) const;
273
274 void initGridFromEGridFile(Opm::EclIO::EclFile& egridfile, std::string fileName);
275 void resetACTNUM( const int* actnum);
276
277 void initBinaryGrid(const Deck& deck);
278
279 void initCornerPointGrid(const std::vector<double>& coord ,
280 const std::vector<double>& zcorn ,
281 const int * actnum);
282
283 bool keywInputBeforeGdfile(const Deck& deck, const std::string& keyword) const;
284
285 void initCylindricalGrid(const Deck&);
286 void initSpiderwebGrid(const Deck&);
287 void initSpiderwebOrCylindricalGrid(const Deck&, const bool);
288 void initCartesianGrid(const Deck&);
289 void initDTOPSGrid(const Deck&);
290 void initDVDEPTHZGrid(const Deck&);
291 void initGrid(const Deck&, const int* actnum);
292 void initCornerPointGrid(const Deck&);
293 void assertCornerPointKeywords(const Deck&);
294
295 static bool hasDTOPSKeywords(const Deck&);
296 static void assertVectorSize(const std::vector<double>& vector, size_t expectedSize, const std::string& msg);
297
298 static std::vector<double> createTOPSVector(const std::array<int, 3>& dims, const std::vector<double>& DZ, const Deck&);
299 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&);
300 static void scatterDim(const std::array<int, 3>& dims , size_t dim , const std::vector<double>& DV , std::vector<double>& D);
301
302
303 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;
304 std::vector<double> makeZcornDzTops(const std::vector<double>& dz, const std::vector<double>& tops) const ;
305 std::vector<double> makeZcornDzvDepthz(const std::vector<double>& dzv, const std::vector<double>& depthz) const;
306 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;
307
308 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;
309 void getCellCorners(const std::size_t globalIndex,
310 std::array<double,8>& X,
311 std::array<double,8>& Y,
312 std::array<double,8>& Z) const;
313
314 };
315
317 public:
318 CoordMapper(size_t nx, size_t ny);
319 size_t size() const;
320
321
322 /*
323 dim = 0,1,2 for x, y and z coordinate respectively.
324 layer = 0,1 for k=0 and k=nz layers respectively.
325 */
326
327 size_t index(size_t i, size_t j, size_t dim, size_t layer) const;
328 private:
329 size_t nx;
330 size_t ny;
331 };
332
333
334
336 public:
337 ZcornMapper(size_t nx, size_t ny, size_t nz);
338 size_t index(size_t i, size_t j, size_t k, int c) const;
339 size_t index(size_t g, int c) const;
340 size_t size() const;
341
342 /*
343 The fixupZCORN method will take a zcorn vector as input and
344 run through it. If the following situation is detected:
345
346 /| /|
347 / | / |
348 / | / |
349 / | / |
350 / | ==> / |
351 / | / |
352 ---/------x /---------x
353 | /
354 |/
355
356 */
357 size_t fixupZCORN( std::vector<double>& zcorn);
358 bool validZCORN( const std::vector<double>& zcorn) const;
359 private:
360 std::array<size_t,3> dims;
361 std::array<size_t,3> stride;
362 std::array<size_t,8> cell_shift;
363 };
364}
365
366#endif // OPM_PARSER_ECLIPSE_GRID_HPP
Definition EclipseGrid.hpp:316
Definition Deck.hpp:49
Definition EclFile.hpp:36
About cell information and dimension: The actual grid information is held in a pointer to an ERT ecl_...
Definition EclipseGrid.hpp:55
void setMINPVV(const std::vector< double > &minpvv)
Sets MINPVV if MINPV and MINPORV are not used.
Definition EclipseGrid.cpp:1936
size_t getGlobalIndex(size_t active_index) const
Observe: the input argument is assumed to be in the space [0,num_active).
Definition EclipseGrid.cpp:569
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 EclipseGrid.cpp:1891
bool isValidCellGeomtry(const std::size_t globalIndex, const UnitSystem &usys) const
Whether or not given cell has a valid cell geometry.
Definition EclipseGrid.cpp:1704
Definition GridDims.hpp:31
Definition UnitSystem.hpp:34
Definition EclipseGrid.hpp:335
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30