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
21
22#ifndef OPM_PARSER_ECLIPSE_GRID_HPP
23#define OPM_PARSER_ECLIPSE_GRID_HPP
24#include <opm/input/eclipse/EclipseState/Grid/GridDims.hpp>
25#include <opm/input/eclipse/EclipseState/Grid/MapAxes.hpp>
26#include <opm/input/eclipse/EclipseState/Grid/MinpvMode.hpp>
27#include <opm/input/eclipse/EclipseState/Grid/PinchMode.hpp>
28
29#include <array>
30#include <memory>
31#include <optional>
32#include <stdexcept>
33#include <string>
34#include <unordered_set>
35#include <vector>
36#include <map>
37
38
39namespace Opm {
40
41 class Deck;
42 namespace EclIO {
43 class EclFile;
44 class EclOutput;
45 }
46 struct NNCdata;
47 class UnitSystem;
48 class ZcornMapper;
49 class EclipseGridLGR;
50 class LgrCollection;
62 class EclipseGrid : public GridDims {
63 public:
65 explicit EclipseGrid(const std::string& filename);
66
67 /*
68 These constructors will make a copy of the src grid, with
69 zcorn and or actnum have been adjustments.
70 */
71 EclipseGrid(const EclipseGrid& src, const std::vector<int>& actnum);
72 EclipseGrid(const EclipseGrid& src, const double* zcorn, const std::vector<int>& actnum);
73
74 EclipseGrid(size_t nx, size_t ny, size_t nz,
75 double dx = 1.0, double dy = 1.0, double dz = 1.0,
76 double top = 0.0);
77 explicit EclipseGrid(const GridDims& gd);
78
79 EclipseGrid(const std::array<int, 3>& dims ,
80 const std::vector<double>& coord ,
81 const std::vector<double>& zcorn ,
82 const int * actnum = nullptr);
83
84 virtual ~EclipseGrid() = default;
85
88 explicit EclipseGrid(const Deck& deck, const int * actnum = nullptr);
89
90 static bool hasGDFILE(const Deck& deck);
91 static bool hasRadialKeywords(const Deck& deck);
92 static bool hasSpiderKeywords(const Deck& deck);
93 static bool hasCylindricalKeywords(const Deck& deck);
94 static bool hasCornerPointKeywords(const Deck&);
95 static bool hasCartesianKeywords(const Deck&);
96 size_t getNumActive( ) const;
97 bool allActive() const;
98
99 size_t activeIndex(size_t i, size_t j, size_t k) const;
100 size_t activeIndex(size_t globalIndex) const;
101
102 size_t getTotalActiveLGR() const;
103 size_t getActiveIndexLGR(const std::string& label, size_t i, size_t j, size_t k) const;
104 size_t getActiveIndexLGR(const std::string& label, size_t localIndex) const;
105
106 size_t activeIndexLGR(const std::string& label, size_t i, size_t j, size_t k) const;
107 size_t activeIndexLGR(const std::string& label, size_t localIndex) const;
108
109 size_t getActiveIndex(size_t i, size_t j, size_t k) const {
110 return activeIndex(i, j, k);
111 }
112
113 size_t getActiveIndex(size_t globalIndex) const {
114 return activeIndex(globalIndex);
115 }
116
117 std::vector<std::string> get_all_lgr_labels() const
118 {
119 return {this->all_lgr_labels.begin() + 1, this->all_lgr_labels.end()};
120 }
121
122 const std::vector<std::string>& get_all_labels() const
123 {
124 return this->all_lgr_labels;
125 }
126
127 std::string get_lgr_tag() const {
128 return this->lgr_label;
129 }
130
131 std::vector<GridDims> get_lgr_children_gridim() const;
132
133
134 void assertIndexLGR(size_t localIndex) const;
135
136 void assertLabelLGR(const std::string& label) const;
137
138 void save(const std::string& filename, bool formatted, const std::vector<Opm::NNCdata>& nnc, const Opm::UnitSystem& units) const;
139 /*
140 Observe that the there is a getGlobalIndex(i,j,k)
141 implementation in the base class. This method - translating
142 from an active index to a global index must be implemented
143 in the current class.
144 */
145 void init_children_host_cells(bool logical = true);
146 void init_children_host_cells_logical(void);
147 void init_children_host_cells_geometrical(void);
148
149 using GridDims::getGlobalIndex;
150 size_t getGlobalIndex(size_t active_index) const;
151
152 /*
153 For RADIAL grids you can *optionally* use the keyword
154 'CIRCLE' to denote that period boundary conditions should be
155 applied in the 'THETA' direction; this will only apply if
156 the theta keywords entered sum up to exactly 360 degrees!
157 */
158 bool circle() const;
159 bool isPinchActive() const;
160 double getPinchThresholdThickness() const;
161 PinchMode getPinchOption() const;
162 PinchMode getMultzOption() const;
163 PinchMode getPinchGapMode() const;
164 double getPinchMaxEmptyGap() const;
165 MinpvMode getMinpvMode() const;
166 const std::vector<double>& getMinpvVector( ) const;
167
168 /*
169 Will return a vector of nactive elements. The method will
170 behave differently depending on the lenght of the
171 input_vector:
172
173 nx*ny*nz: only the values corresponding to active cells
174 are copied out.
175
176 nactive: The input vector is copied straight out again.
177
178 ??? : Exception.
179 */
180
181 template<typename T>
182 std::vector<T> compressedVector(const std::vector<T>& input_vector) const {
183 if( input_vector.size() == this->getNumActive() ) {
184 return input_vector;
185 }
186
187 if (input_vector.size() != getCartesianSize())
188 throw std::invalid_argument("Input vector must have full size");
189
190 {
191 std::vector<T> compressed_vector( this->getNumActive() );
192 const auto& active_map = this->getActiveMap( );
193
194 for (size_t i = 0; i < this->getNumActive(); ++i)
195 compressed_vector[i] = input_vector[ active_map[i] ];
196
197 return compressed_vector;
198 }
199 }
200
201
204 const std::vector<int>& getActiveMap() const;
205
206 void init_lgr_cells(const LgrCollection& lgr_input);
207 void create_lgr_cells_tree(const LgrCollection& );
209 std::tuple<std::array<double, 3>,std::array<double, 3>,std::array<double, 3>>
210 getCellAndBottomCenterNormal(size_t globalIndex) const;
211 std::array<double, 3> getCellCenter(size_t i,size_t j, size_t k) const;
212 std::array<double, 3> getCellCenter(size_t globalIndex) const;
213 std::array<double, 3> getCornerPos(size_t i,size_t j, size_t k, size_t corner_index) const;
214 const std::vector<double>& activeVolume() const;
215 double getCellVolume(size_t globalIndex) const;
216 double getCellVolume(size_t i , size_t j , size_t k) const;
217 double getCellThickness(size_t globalIndex) const;
218 double getCellThickness(size_t i , size_t j , size_t k) const;
219 std::array<double, 3> getCellDims(size_t i,size_t j, size_t k) const;
220 std::array<double, 3> getCellDims(size_t globalIndex) const;
221 bool cellActive( size_t globalIndex ) const;
222 bool cellActive( size_t i , size_t j, size_t k ) const;
223 bool cellActiveAfterMINPV( size_t i , size_t j , size_t k, double cell_porv ) const;
224 bool is_lgr() const {return lgr_grid;};
225 std::array<double, 3> getCellDimensions(size_t i, size_t j, size_t k) const {
226 return getCellDims(i, j, k);
227 }
228 std::array<double,3> getCellDimensionsLGR(const std::size_t i,
229 const std::size_t j,
230 const std::size_t k,
231 const std::string& lgr_tag) const;
232 double getCellDepthLGR(size_t i, size_t j, size_t k, const std::string& lgr_tag) const;
233
234
235 bool isCellActive(size_t i, size_t j, size_t k) const {
236 return cellActive(i, j, k);
237 }
238
244 bool isValidCellGeomtry(const std::size_t globalIndex,
245 const UnitSystem& usys) const;
246
247 double getCellDepth(size_t i,size_t j, size_t k) const;
248 double getCellDepth(size_t globalIndex) const;
249 ZcornMapper zcornMapper() const;
250
251 const std::vector<double>& getCOORD() const;
252 const std::vector<double>& getZCORN() const;
253 const std::vector<int>& getACTNUM( ) const;
254
255 const std::optional<MapAxes>& getMapAxes() const;
256
257 const std::map<size_t, std::array<int,2>>& getAquiferCellTabnums() const;
258
259 /*
260 The fixupZCORN method is run as part of constructiong the grid. This will adjust the
261 z-coordinates to ensure that cells do not overlap. The return value is the number of
262 points which have been adjusted. The number of zcorn nodes that has ben fixed is
263 stored in private member zcorn_fixed.
264 */
265
266 size_t fixupZCORN();
267 size_t getZcornFixed() { return zcorn_fixed; };
268
269 // resetACTNUM with no arguments will make all cells in the grid active.
270
271 void resetACTNUM();
272 void resetACTNUM( const std::vector<int>& actnum);
273
275 void setMINPVV(const std::vector<double>& minpvv);
276
277 bool equal(const EclipseGrid& other) const;
278 static bool hasDVDEPTHZKeywords(const Deck&);
279
280 /*
281 For ALugrid we can *only* use the keyword <DXV, DXYV, DZV, DEPTHZ> so to
282 initialize a Regular Cartesian Grid; further we need equidistant mesh
283 spacing in each direction to initialize ALuGrid (mandatory for
284 mesh refinement!).
285 */
286
287 static bool hasEqualDVDEPTHZ(const Deck&);
288 static bool allEqual(const std::vector<double> &v);
289 EclipseGridLGR& getLGRCell(std::size_t index);
290 const EclipseGridLGR& getLGRCell(const std::string& lgr_tag) const;
291 int getLGR_global_father(std::size_t global_index, const std::string& lgr_tag) const;
292
293 std::vector<EclipseGridLGR> lgr_children_cells;
301 virtual void set_lgr_refinement(const std::string& lgr_tag, const std::vector<double>& coords, const std::vector<double>& zcorn);
302
303 protected:
304 std::size_t lgr_global_counter = 0;
305 std::string lgr_label = "GLOBAL";
306 int lgr_level = 0;
307 int lgr_level_father = 0;
308 std::vector<std::string> lgr_children_labels;
309 std::vector<std::size_t> lgr_active_index;
310 std::vector<std::size_t> lgr_level_active_map;
311 std::vector<std::string> all_lgr_labels;
312 std::map<std::vector<std::size_t>, std::size_t> num_lgr_children_cells;
313 std::vector<double> m_zcorn;
314 std::vector<double> m_coord;
315 std::vector<int> m_actnum;
316 std::vector<std::size_t> m_print_order_lgr_cells;
317 // Input grid data.
318 mutable std::optional<std::vector<double>> m_input_zcorn;
319 mutable std::optional<std::vector<double>> m_input_coord;
320
321 private:
322 std::vector<double> m_minpvVector;
323 MinpvMode m_minpvMode;
324 std::optional<double> m_pinch;
325
326 // Option 4 of PINCH (TOPBOT/ALL), how to calculate TRANS
327 PinchMode m_pinchoutMode;
328 // Option 5 of PINCH (TOP/ALL), how to apply MULTZ
329 PinchMode m_multzMode;
330 // Option 2 of PINCH (GAP/NOGAP)
331 PinchMode m_pinchGapMode;
332 double m_pinchMaxEmptyGap;
333 bool lgr_grid = false;
334 mutable std::optional<std::vector<double>> active_volume;
335
336 bool m_circle = false;
337 size_t zcorn_fixed = 0;
338 bool m_useActnumFromGdfile = false;
339
340 std::optional<MapAxes> m_mapaxes;
341
342 // Mapping to/from active cells.
343 int m_nactive {};
344 std::vector<int> m_active_to_global;
345 std::vector<int> m_global_to_active;
346 // Numerical aquifer cells, needs to be active
347 std::unordered_set<size_t> m_aquifer_cells;
348 // Keep track of aquifer cell depths and (pvtnum,satnum)
349 std::map<size_t, double> m_aquifer_cell_depths;
350 std::map<size_t, std::array<int,2>> m_aquifer_cell_tabnums;
351
352 // Radial grids need this for volume calculations.
353 std::optional<std::vector<double>> m_thetav;
354 std::optional<std::vector<double>> m_rv;
355 void parseGlobalReferenceToChildren(void);
356 int initializeLGRObjectIndices(int);
357 void initializeLGRTreeIndices(void);
358 void propagateParentIndicesToLGRChildren(int);
359 void updateNumericalAquiferCells(const Deck&);
360 double computeCellGeometricDepth(size_t globalIndex) const;
361
362 void initGridFromEGridFile(Opm::EclIO::EclFile& egridfile,
363 const std::string& fileName);
364 void resetACTNUM( const int* actnum);
365
366 void initBinaryGrid(const Deck& deck);
367
368 void initCornerPointGrid(const std::vector<double>& coord ,
369 const std::vector<double>& zcorn ,
370 const int * actnum);
371
372 bool keywInputBeforeGdfile(const Deck& deck, const std::string& keyword) const;
373
374 void initCylindricalGrid(const Deck&);
375 void initSpiderwebGrid(const Deck&);
376 void initSpiderwebOrCylindricalGrid(const Deck&, const bool);
377 void initCartesianGrid(const Deck&);
378 void initDTOPSGrid(const Deck&);
379 void initDVDEPTHZGrid(const Deck&);
380 void initGrid(const Deck&, const int* actnum);
381 void initCornerPointGrid(const Deck&);
382 void assertCornerPointKeywords(const Deck&);
383 void save_all_lgr_labels(const LgrCollection& );
384 static bool hasDTOPSKeywords(const Deck&);
385 static void assertVectorSize(const std::vector<double>& vector, size_t expectedSize, const std::string& msg);
386
387 static std::vector<double> createTOPSVector(const std::array<int, 3>& dims, const std::vector<double>& DZ, const Deck&);
388 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&);
389 static void scatterDim(const std::array<int, 3>& dims , size_t dim , const std::vector<double>& DV , std::vector<double>& D);
390
391
392 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;
393 std::vector<double> makeZcornDzTops(const std::vector<double>& dz, const std::vector<double>& tops) const ;
394 std::vector<double> makeZcornDzvDepthz(const std::vector<double>& dzv, const std::vector<double>& depthz) const;
395 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;
396
397 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;
398 void getCellCorners(const std::size_t globalIndex,
399 std::array<double,8>& X,
400 std::array<double,8>& Y,
401 std::array<double,8>& Z) const;
402 std::array<int,3> getCellSubdivisionRatioLGR(const std::string& lgr_tag,
403 std::array<int,3> acum = {1,1,1}) const;
404
405
406 };
407
410 {
411 public:
412 using vec_size_t = std::vector<std::size_t>;
413 EclipseGridLGR() = default;
414 EclipseGridLGR(const std::string& self_label,
415 const std::string& father_label_,
416 size_t nx,
417 size_t ny,
418 size_t nz,
419 const vec_size_t& father_lgr_index,
420 const std::array<int, 3>& low_fatherIJK_,
421 const std::array<int, 3>& up_fatherIJK_);
422 const vec_size_t& getFatherGlobalID() const;
423 void save(Opm::EclIO::EclOutput&, const std::vector<Opm::NNCdata>&, const Opm::UnitSystem&) const;
424 void save_nnc(Opm::EclIO::EclOutput&) const;
425 void set_lgr_global_counter(std::size_t counter)
426 {
427 lgr_global_counter = counter;
428 }
429 const vec_size_t& get_father_global() const
430 {
431 return father_global;
432 }
433 std::optional<std::reference_wrapper<EclipseGridLGR>>
434 get_child_LGR_cell(const std::string& lgr_tag) const;
435 std::vector<int> save_hostnum(void) const;
436 int get_hostnum(std::size_t global_index) const {return(m_hostnum[global_index]);};
437 void get_label_child_to_top_father(std::vector<std::reference_wrapper<const std::string>>& list) const;
438 void get_global_index_child_to_top_father(std::vector<std::size_t> & list,
439 std::size_t i, std::size_t j, std::size_t k) const;
440 void get_global_index_child_to_top_father(std::vector<std::size_t> & list, std::size_t global_ind) const;
441
442 void set_hostnum(const std::vector<int>&);
443 const std::array<int,3>& get_low_fatherIJK() const{
444 return low_fatherIJK;
445 }
446 const std::string& get_father_label() const{
447 return father_label;
448 }
449
450 const std::array<int,3>& get_up_fatherIJK() const{
451 return up_fatherIJK;
452 }
453
454
462 void set_lgr_refinement(const std::string& lgr_tag,
463 const std::vector<double>& coord,
464 const std::vector<double>& zcorn) override;
465
466 void set_lgr_refinement(const std::vector<double>&, const std::vector<double>&);
467
468 private:
469 void init_father_global();
470 std::string father_label;
471 // references global on the father label
472 vec_size_t father_global;
473 std::array<int, 3> low_fatherIJK {};
474 std::array<int, 3> up_fatherIJK {};
475 std::vector<int> m_hostnum;
476 };
477
478
480 public:
481 CoordMapper(size_t nx, size_t ny);
482 size_t size() const;
483
484
485 /*
486 dim = 0,1,2 for x, y and z coordinate respectively.
487 layer = 0,1 for k=0 and k=nz layers respectively.
488 */
489
490 size_t index(size_t i, size_t j, size_t dim, size_t layer) const;
491 private:
492 size_t nx;
493 size_t ny;
494 };
495
496
497
499 public:
500 ZcornMapper(size_t nx, size_t ny, size_t nz);
501 size_t index(size_t i, size_t j, size_t k, int c) const;
502 size_t index(size_t g, int c) const;
503 size_t size() const;
504
505 /*
506 The fixupZCORN method will take a zcorn vector as input and
507 run through it. If the following situation is detected:
508
509 /| /|
510 / | / |
511 / | / |
512 / | / |
513 / | ==> / |
514 / | / |
515 ---/------x /---------x
516 | /
517 |/
518
519 */
520 size_t fixupZCORN( std::vector<double>& zcorn);
521 bool validZCORN( const std::vector<double>& zcorn) const;
522 private:
523 std::array<size_t,3> dims;
524 std::array<size_t,3> stride;
525 std::array<size_t,8> cell_shift;
526 };
527}
528#endif // OPM_PARSER_ECLIPSE_GRID_HPP
Definition EclipseGrid.hpp:479
Definition Deck.hpp:46
Definition EclFile.hpp:36
Definition EclOutput.hpp:38
Specialized Class to describe LGR refined cells.
Definition EclipseGrid.hpp:410
void set_lgr_refinement(const std::string &lgr_tag, const std::vector< double > &coord, const std::vector< double > &zcorn) override
Sets Local Grid Refinement for the EclipseGridLGR.
Definition EclipseGrid.cpp:2722
About cell information and dimension: The actual grid information is held in a pointer to an ERT ecl_...
Definition EclipseGrid.hpp:62
void setMINPVV(const std::vector< double > &minpvv)
Sets MINPVV if MINPV and MINPORV are not used.
Definition EclipseGrid.cpp:2465
std::array< double, 3 > getCellDimensionsLGR(const std::size_t i, const std::size_t j, const std::size_t k, const std::string &lgr_tag) const
Computes the dimensions of a local grid refinement (LGR) cell.
Definition EclipseGrid.cpp:2037
std::tuple< std::array< double, 3 >, std::array< double, 3 >, std::array< double, 3 > > getCellAndBottomCenterNormal(size_t globalIndex) const
get cell center, and center and normal of bottom face
Definition EclipseGrid.cpp:1674
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:2204
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:589
virtual void set_lgr_refinement(const std::string &lgr_tag, const std::vector< double > &coords, const std::vector< double > &zcorn)
Sets Local Grid Refinement for the EclipseGrid.
Definition EclipseGrid.cpp:2127
bool isValidCellGeomtry(const std::size_t globalIndex, const UnitSystem &usys) const
Whether or not given cell has a valid cell geometry.
Definition EclipseGrid.cpp:1783
Definition GridDims.hpp:31
Definition LgrCollection.hpp:33
Definition UnitSystem.hpp:34
Definition EclipseGrid.hpp:498
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30