19 #ifndef OPM_IO_EGRID_HPP
20 #define OPM_IO_EGRID_HPP
22 #include <opm/io/eclipse/EclFile.hpp>
23 #include <opm/common/utility/FileSystem.hpp>
33 namespace Opm {
namespace EclIO {
38 explicit EGrid(
const std::string& filename, std::string grid_name =
"global");
40 int global_index(
int i,
int j,
int k)
const;
41 int active_index(
int i,
int j,
int k)
const;
43 const std::array<int, 3>& dimension()
const {
return nijk; }
45 std::array<int, 3> ijk_from_active_index(
int actInd)
const;
46 std::array<int, 3> ijk_from_global_index(
int globInd)
const;
48 void getCellCorners(
int globindex, std::array<double,8>& X, std::array<double,8>& Y, std::array<double,8>& Z);
49 void getCellCorners(
const std::array<int, 3>& ijk, std::array<double,8>& X, std::array<double,8>& Y, std::array<double,8>& Z);
51 std::vector<std::array<float, 3>> getXYZ_layer(
int layer,
bool bottom=
false);
52 std::vector<std::array<float, 3>> getXYZ_layer(
int layer,
const std::array<int, 4>& box,
bool bottom=
false);
54 int activeCells()
const {
return nactive; }
55 int totalNumberOfCells()
const {
return nijk[0] * nijk[1] * nijk[2]; }
57 void load_grid_data();
59 bool is_radial()
const {
return m_radial; }
61 const std::vector<int>& hostCellsGlobalIndex()
const {
return host_cells; }
62 std::vector<std::array<int, 3>> hostCellsIJK();
65 using NNCentry = std::tuple<int, int, int, int, int, int, float>;
66 std::vector<NNCentry> get_nnc_ijk();
68 const std::vector<std::string>& list_of_lgrs()
const {
return lgr_names; }
70 const std::vector<float>& get_mapaxes()
const {
return m_mapaxes; }
71 const std::string& get_mapunits()
const {
return m_mapunits; }
74 Opm::filesystem::path inputFileName, initFileName;
75 std::string m_grid_name;
78 std::vector<float> m_mapaxes;
79 std::string m_mapunits;
81 std::array<int, 3> nijk;
82 std::array<int, 3> host_nijk;
85 mutable bool m_nncs_loaded;
87 std::vector<int> act_index;
88 std::vector<int> glob_index;
90 std::vector<float> coord_array;
91 std::vector<float> zcorn_array;
93 std::vector<int> nnc1_array;
94 std::vector<int> nnc2_array;
95 std::vector<float> transnnc_array;
96 std::vector<int> host_cells;
98 std::vector<std::string> lgr_names;
100 int zcorn_array_index;
101 int coord_array_index;
102 int actnum_array_index;
103 int nnc1_array_index;
104 int nnc2_array_index;
106 std::vector<float> get_zcorn_from_disk(
int layer,
bool bottom);
108 void getCellCorners(
const std::array<int, 3>& ijk,
const std::vector<float>& zcorn_layer,
109 std::array<double,4>& X, std::array<double,4>& Y, std::array<double,4>& Z);
Definition: EclFile.hpp:35
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition: Exceptions.hpp:29