diff --git a/opm/input/eclipse/EclipseState/Grid/EclipseGrid.cpp b/opm/input/eclipse/EclipseState/Grid/EclipseGrid.cpp index 75a20358c08..a73ccc1625c 100644 --- a/opm/input/eclipse/EclipseState/Grid/EclipseGrid.cpp +++ b/opm/input/eclipse/EclipseState/Grid/EclipseGrid.cpp @@ -2823,13 +2823,309 @@ namespace Opm { EclipseGrid::perform_refinement(); } - - std::vector EclipseGridLGR::generate_refined_zcorn([[maybe_unused]] const std::vector& parent_coord, - [[maybe_unused]] const std::vector& parent_zcorn, - [[maybe_unused]] const std::array& parent_nxyz) + std::vector EclipseGridLGR::generate_refined_zcorn(const std::vector& parent_coord, + const std::vector& parent_zcorn, + const std::array& parent_nxyz) { + using CoordinateType = std::array; + using PillarType = std::array, 2>; + using ZCornType = double; + + struct Matrix { + std::size_t rows, cols; + std::vector data; + + Matrix() : rows(0), cols(0), data() {} + + Matrix(std::size_t r, std::size_t c) + : rows(r), cols(c), data(r * c, 0) + {} + + ZCornType & operator()(std::size_t i, std::size_t j) { + return data[i * cols + j]; + } + + const ZCornType& operator()(std::size_t i, std::size_t j) const { + return data[i * cols + j]; + } + }; + + + /** + * Bilinear interpolation on a single layer of a quadrilateral, + * selected from an 8-point array (bottom + top layers). + * + * Point order in each layer (0-3): + * + * 2 ---- 3 (top-left → top-right) + * | | + * 0 ---- 1 (bottom-left → bottom-right) + * + * points[0..3] -> bottom layer + * points[4..7] -> top layer + * + * I,J: indices inside the grid cell + * NI,NJ: maximum indices along each direction (size of the cell) + * topLayer: false = use bottom layer, true = use top layer + */ + + auto bilinearInterpolation = [](const std::array& points, + int I, int J, int NI, int NJ, + bool topLayer) -> CoordinateType + { + int offset = topLayer ? 4 : 0; + + // Normalize indices to [0,1] + double x = static_cast(I) / NI; + double y = static_cast(J) / NJ; - const std::vector refined_zcorn; + double w00 = (1 - x) * (1 - y); + double w10 = x * (1 - y); + double w01 = (1 - x) * y; + double w11 = x * y; + + CoordinateType P; + P[0] = w00*points[offset + 0][0] + w10*points[offset + 1][0] + + w01*points[offset + 2][0] + w11*points[offset + 3][0]; + + P[1] = w00*points[offset + 0][1] + w10*points[offset + 1][1] + + w01*points[offset + 2][1] + w11*points[offset + 3][1]; + + P[2] = w00*points[offset + 0][2] + w10*points[offset + 1][2] + + w01*points[offset + 2][2] + w11*points[offset + 3][2]; + + return P; + }; + + auto generate_dest_ijk_refinement_intervals = [this](const std::vector& i_elem, + const std::vector& j_elem, + const std::vector& k_elem) + { + + std::size_t totalX = getNX(); + std::size_t totalY = getNY(); + std::size_t totalZ = getNZ(); + + std::size_t min_i = *std::min_element(i_elem.begin(), i_elem.end()); + std::size_t max_i = *std::max_element(i_elem.begin(), i_elem.end()); + std::size_t min_j = *std::min_element(j_elem.begin(), j_elem.end()); + std::size_t max_j = *std::max_element(j_elem.begin(), j_elem.end()); + std::size_t min_k = *std::min_element(k_elem.begin(), k_elem.end()); + std::size_t max_k = *std::max_element(k_elem.begin(), k_elem.end()); + + std::size_t num_el_i = max_i - min_i + 1; + std::size_t num_el_j = max_j - min_j + 1; + std::size_t num_el_k = max_k - min_k + 1; + + if (num_el_i == 0 || num_el_j == 0 || num_el_k == 0) + { + throw std::invalid_argument("Invalid number of elements."); + } + if (totalX % num_el_i != 0 || totalY % num_el_j != 0 || totalZ % num_el_k != 0) { + // Either handle remainder properly or fail fast + throw std::runtime_error("Refinement factors must exactly divide parent NX,NY and NZ"); + } + + auto dX = totalX/num_el_i; + auto dY = totalY/num_el_j; + auto dZ = totalZ/num_el_k; + + std::vector refined_i_indices; + std::vector refined_j_indices; + std::vector refined_k_indices; + + refined_i_indices.resize(i_elem.size()); + refined_j_indices.resize(j_elem.size()); + refined_k_indices.resize(k_elem.size()); + + std::transform(i_elem.begin(), i_elem.end(), refined_i_indices.begin(), + [min_i, dX](std::size_t i){ + return (i - min_i) * dX; + }); + std::transform(j_elem.begin(), j_elem.end(), refined_j_indices.begin(), + [min_j, dY](std::size_t j){ + return (j - min_j) * dY; + }); + std::transform(k_elem.begin(), k_elem.end(), refined_k_indices.begin(), + [min_k, dZ](std::size_t k){ + return (k - min_k) * dZ; + }); + + return std::make_tuple(refined_i_indices, refined_j_indices, refined_k_indices, dX, dY, dZ); + }; + + CoordMapper parent_coord_mapper(parent_nxyz[0], parent_nxyz[1]); + ZcornMapper parent_zcorn_mapper(parent_nxyz[0], parent_nxyz[1], parent_nxyz[2]); + + CoordMapper local_coord_mapper(this->getNX(), this->getNY()); + ZcornMapper local_zcorn_mapper(this->getNX(), this->getNY(), this->getNZ()); + + auto get_father_pillar = [&parent_coord_mapper, &parent_coord]( std::size_t i, + std::size_t j) + { + PillarType P; + // Top point + P[0][0] = parent_coord[ parent_coord_mapper.index(i, j, 0, 0) ]; // x_top + P[0][1] = parent_coord[ parent_coord_mapper.index(i, j, 1, 0) ]; // y_top + P[0][2] = parent_coord[ parent_coord_mapper.index(i, j, 2, 0) ]; // z_top + // Bottom point + P[1][0] = parent_coord[ parent_coord_mapper.index(i, j, 0, 1) ]; // x_bot + P[1][1] = parent_coord[ parent_coord_mapper.index(i, j, 1, 1) ]; // y_bot + P[1][2] = parent_coord[ parent_coord_mapper.index(i, j, 2, 1) ]; // z_bot + return P; + }; + + auto get_local_pillar = [this, &local_coord_mapper]( std::size_t i, std::size_t j) + { + PillarType P; + // Top point + P[0][0] = m_coord[ local_coord_mapper.index(i, j, 0, 0) ]; // x_top + P[0][1] = m_coord[ local_coord_mapper.index(i, j, 1, 0) ]; // y_top + P[0][2] = m_coord[ local_coord_mapper.index(i, j, 2, 0) ]; // z_top + // Bottom point + P[1][0] = m_coord[ local_coord_mapper.index(i, j, 0, 1) ]; // x_bot + P[1][1] = m_coord[ local_coord_mapper.index(i, j, 1, 1) ]; // y_bot + P[1][2] = m_coord[ local_coord_mapper.index(i, j, 2, 1) ]; // z_bot + return P; + }; + + auto get_zcorn_index = [&parent_zcorn_mapper]( std::size_t i, + std::size_t j, + std::size_t k) + { + std::array father_zcorn_indices; + for (int c = 0; c < 8; c++) { + father_zcorn_indices[c] = parent_zcorn_mapper.index(i, j, k, c); + } + return father_zcorn_indices; + }; + + auto get_zcorn_value = [&parent_zcorn](std::array indices) + { + std::array zcorn_values; + for (int idx = 0; idx < 8; idx++) { + zcorn_values[idx] = parent_zcorn[ indices[idx] ]; + } + return zcorn_values; + }; + + auto pillar_zcorn_to_coord = [](const PillarType& pillar, double zcorn_value) { + double t = (zcorn_value - pillar[0][2]) / (pillar[1][2] - pillar[0][2]); + double x = pillar[0][0] + t * (pillar[1][0] - pillar[0][0]); + double y = pillar[0][1] + t * (pillar[1][1] - pillar[0][1]); + return std::array{x, y, zcorn_value}; + }; + + auto coord_pillar_to_zcorn = [](const PillarType& pillar, const std::array& coord) { + double total_height = pillar[1][2] - pillar[0][2]; + if (total_height == 0.0) { + return pillar[0][2]; // or pillar[1][2], they are the same + } + double t = (coord[2] - pillar[0][2]) / total_height; + return pillar[0][2] + t * total_height; + }; + + auto find_coord_coorners = [&](const std::array& pillars, + const std::array& father_zcorn_values) + { + auto coord_corners = std::array,8>{}; + for (int c = 0; c < 8; c++) { + coord_corners[c] = pillar_zcorn_to_coord(pillars[c % 4], father_zcorn_values[c]); + } + return coord_corners; + }; + + auto interpolate_inner_zcorn = [](std::vector& pillars_matrix) + { + std::size_t nz = pillars_matrix.size(); + for (std::size_t k = 1; k < nz - 1; k++) { + double t = static_cast(k) / (nz - 1); + for (std::size_t j = 0; j < pillars_matrix[0].cols; j++) { + for (std::size_t i = 0; i < pillars_matrix[0].rows; i++) { + pillars_matrix[k](i,j) = (1 - t) * pillars_matrix[0](i,j) + t * pillars_matrix.back()(i,j); + } + } + } + + }; + + + const auto [NX, NY, NZ] = getNXYZ(); + std::vector refined_zcorn; + refined_zcorn.resize(NX*NY*NZ*8); + + const auto [i_list, j_list, k_list ] = VectorUtil::generate_cartesian_product( + low_fatherIJK[0], up_fatherIJK[0], + low_fatherIJK[1], up_fatherIJK[1], + low_fatherIJK[2], up_fatherIJK[2]); + + const auto [dest_i_ref, dest_j_ref, dest_k_ref, dx , dy, dz] = generate_dest_ijk_refinement_intervals(i_list, j_list, k_list); + + + for (std::size_t n = 0; n < i_list.size(); n++) { + const std::size_t i_father = i_list[n]; + const std::size_t j_father = j_list[n]; + const std::size_t k_father = k_list[n]; + + const std::size_t i_ref = dest_i_ref[n]; + const std::size_t j_ref = dest_j_ref[n]; + const std::size_t k_ref = dest_k_ref[n]; + + std::vector local_zcorn_volumes; + local_zcorn_volumes.resize(dz+1); + for (std::size_t k = 0; k <= dz; k++) { + local_zcorn_volumes[k] = Matrix(dx+1, dy+1); + } + + std::array pillars; + pillars[0] = get_father_pillar(i_father, j_father); + pillars[1] = get_father_pillar(i_father+1, j_father); + pillars[2] = get_father_pillar(i_father, j_father+1); + pillars[3] = get_father_pillar(i_father+1, j_father+1); + + auto father_zcorn_indices = get_zcorn_index(i_father, j_father, k_father); + auto father_zcorn_values = get_zcorn_value(father_zcorn_indices); + auto coord_corners = find_coord_coorners(pillars, father_zcorn_values); + + // fill in the refined zcorn values of top and bottom layers + for (std::size_t jj = 0; jj <= dy; jj++) { + for (std::size_t ii = 0; ii <= dx; ii++) { + auto interpolated_coord_bottom = bilinearInterpolation(coord_corners, ii, jj, dx, dy, false); + auto interpolated_coord_top = bilinearInterpolation(coord_corners, ii, jj, dx, dy, true); + auto local_pillar = get_local_pillar(i_ref + ii, j_ref + jj); + local_zcorn_volumes[0](ii,jj) = coord_pillar_to_zcorn(local_pillar, interpolated_coord_bottom); + local_zcorn_volumes.back()(ii,jj) = coord_pillar_to_zcorn(local_pillar, interpolated_coord_top); + } + } + + // interpolate inner layers + interpolate_inner_zcorn(local_zcorn_volumes); + + // assign to refined_zcorn + for (std::size_t kk = 0; kk < dz; kk++) { + for (std::size_t jj = 0; jj < dy; jj++) { + for (std::size_t ii = 0; ii < dx; ii++) { + std::size_t global_i = i_ref + ii; + std::size_t global_j = j_ref + jj; + std::size_t global_k = k_ref + kk; + std::array refined_zcorn_values = {local_zcorn_volumes[kk](ii,jj), + local_zcorn_volumes[kk](ii+1,jj), + local_zcorn_volumes[kk](ii,jj+1), + local_zcorn_volumes[kk](ii+1,jj+1), + local_zcorn_volumes[kk+1](ii,jj), + local_zcorn_volumes[kk+1](ii+1,jj), + local_zcorn_volumes[kk+1](ii,jj+1), + local_zcorn_volumes[kk+1](ii+1,jj+1)}; + + for (std::size_t idx = 0; idx < 8; idx++) { + const std::size_t zcorn_idx = local_zcorn_mapper.index(global_i, global_j, global_k, idx); + refined_zcorn[zcorn_idx] = refined_zcorn_values[idx]; + } + + } + } + } + } return refined_zcorn; } diff --git a/tests/parser/LgrOutputTests.cpp b/tests/parser/LgrOutputTests.cpp index f81800b01be..39c57bcaf66 100644 --- a/tests/parser/LgrOutputTests.cpp +++ b/tests/parser/LgrOutputTests.cpp @@ -24,6 +24,7 @@ #include +#include #include #include #include @@ -34,9 +35,6 @@ #include #include -#include -#include -#include #include #include #include @@ -45,20 +43,95 @@ using namespace Opm; namespace { + template + void check_vec_close(const T& actual, const T& expected, double tol = 1e-5) + { + BOOST_REQUIRE_EQUAL(actual.size(), expected.size()); + for (std::size_t i = 0; i < actual.size(); ++i) + BOOST_CHECK_CLOSE(actual[i], expected[i], tol); + } - auto write_coord_to_txt = [](const std::vector& coord, const std::string& filename) { - std::ofstream outfile(filename); - if (!outfile.is_open()) { - throw std::runtime_error("Failed to open file: " + filename); - } - for (const auto& value : coord) { - outfile << value << "\n"; - } + auto init_deck(const std::string& deck_string) { + Opm::Parser parser; + return parser.parseString(deck_string); + }; + - outfile.close(); +std::tuple, std::vector> final_test_data() +{ + std::vector coord = { + 0, 0, 2537.5, 0, 0, 2601.5, + 101.6, 0, 2537.5, 101.6, 0, 2601.5, + 203.2, 0, 2537.5, 203.2, 0, 2601.5, + 304.8, 0, 2537.5, 304.8, 0, 2601.5, + 406.4, 0, 2537.5, 406.4, 0, 2601.5, + 508, 0, 2537.5, 508, 0, 2601.5, + 609.6, 0, 2537.5, 609.6, 0, 2601.5, + 0, 50.8, 2537.5, 0, 50.8, 2601.5, + 101.6, 50.8, 2537.5, 101.6, 50.8, 2601.5, + 203.2, 50.8, 2537.5, 203.2, 50.8, 2601.5, + 304.8, 50.8, 2537.5, 304.8, 50.8, 2601.5, + 406.4, 50.8, 2537.5, 406.4, 50.8, 2601.5, + 508, 50.8, 2537.5, 508, 50.8, 2601.5, + 609.6, 50.8, 2537.5, 609.6, 50.8, 2601.5, + 0, 101.6, 2537.5, 0, 101.6, 2601.5, + 101.6, 101.6, 2537.5, 101.6, 101.6, 2601.5, + 203.2, 101.6, 2537.5, 203.2, 101.6, 2601.5, + 304.8, 101.6, 2537.5, 304.8, 101.6, 2601.5, + 406.4, 101.6, 2537.5, 406.4, 101.6, 2601.5, + 508, 101.6, 2537.5, 508, 101.6, 2601.5, + 609.6, 101.6, 2537.5, 609.6, 101.6, 2601.5, + 0, 152.4, 2537.5, 0, 152.4, 2601.5, + 101.6, 152.4, 2537.5, 101.6, 152.4, 2601.5, + 203.2, 152.4, 2537.5, 203.2, 152.4, 2601.5, + 304.8, 152.4, 2537.5, 304.8, 152.4, 2601.5, + 406.4, 152.4, 2537.5, 406.4, 152.4, 2601.5, + 508, 152.4, 2537.5, 508, 152.4, 2601.5, + 609.6, 152.4, 2537.5, 609.6, 152.4, 2601.5, + 0, 203.2, 2537.5, 0, 203.2, 2601.5, + 101.6, 203.2, 2537.5, 101.6, 203.2, 2601.5, + 203.2, 203.2, 2537.5, 203.2, 203.2, 2601.5, + 304.8, 203.2, 2537.5, 304.8, 203.2, 2601.5, + 406.4, 203.2, 2537.5, 406.4, 203.2, 2601.5, + 508, 203.2, 2537.5, 508, 203.2, 2601.5, + 609.6, 203.2, 2537.5, 609.6, 203.2, 2601.5, + 0, 254, 2537.5, 0, 254, 2601.5, + 101.6, 254, 2537.5, 101.6, 254, 2601.5, + 203.2, 254, 2537.5, 203.2, 254, 2601.5, + 304.8, 254, 2537.5, 304.8, 254, 2601.5, + 406.4, 254, 2537.5, 406.4, 254, 2601.5, + 508, 254, 2537.5, 508, 254, 2601.5, + 609.6, 254, 2537.5, 609.6, 254, 2601.5, + 0, 304.8, 2537.5, 0, 304.8, 2601.5, + 101.6, 304.8, 2537.5, 101.6, 304.8, 2601.5, + 203.2, 304.8, 2537.5, 203.2, 304.8, 2601.5, + 304.8, 304.8, 2537.5, 304.8, 304.8, 2601.5, + 406.4, 304.8, 2537.5, 406.4, 304.8, 2601.5, + 508, 304.8, 2537.5, 508, 304.8, 2601.5, + 609.6, 304.8, 2537.5, 609.6, 304.8, 2601.5 +}; +std::vector gnint = {144, 288, 288, 288, 288, 288, + 288, 288, 288, 288, 288, 288, + 288, 288, 288, 288, 288, 288, + 288, 288, 144}; +std::vector gnfloat = { + 2537.5, 2538.2, 2539, 2539.8, 2540.5, 2542, 2543.6, 2545.1, 2546.6, 2550.4, 2554.2, 2558, + 2561.8, 2566.4, 2571, 2575.6, 2580.1, 2585.5, 2590.8, 2596.1, 2601.5, }; +int sum_gnint = std::accumulate(gnint.begin(), gnint.end(), 0); +std::vector zcorn; +zcorn.reserve(sum_gnint); +for (size_t i = 0; i < gnint.size(); ++i) { + for (size_t j = 0; j < static_cast(gnint[i]); ++j) { + zcorn.push_back(gnfloat[i]); + } +} + return std::make_tuple(coord, zcorn); +} + + std::array dimensions_collumn_lgr(){ const double i = 500; const double j = 500; @@ -176,7 +249,6 @@ SCHEDULE // LgrCollection is used to initalize LGR Cells in the Eclipse Grid. eclipse_grid_file.init_lgr_cells(lgr_col); // LGR COORD and ZCORN is parsed to EclipseGridLGR children cell. (Simulates the process of recieving the LGR refinement.) - eclipse_grid_file.set_lgr_refinement("LGR1",coord_l,zcorn_l); // Intialize host_cell numbering. eclipse_grid_file.init_children_host_cells(); // Save EclipseGrid. @@ -277,8 +349,6 @@ SCHEDULE eclipse_grid_file.init_lgr_cells(lgr_col); // LGR COORD and ZCORN is parsed to EclipseGridLGR children cell. (Simulates the process of recieving the LGR refinement.) const auto& lgr1 = eclipse_grid_file.getLGRCell("LGR1"); - const auto& coordlgr1 = lgr1.getCOORD(); - eclipse_grid_file.set_lgr_refinement("LGR1",coordlgr1,zcorn_l); // Intialize host_cell numbering. eclipse_grid_file.init_children_host_cells(); // Save EclipseGrid. @@ -302,7 +372,6 @@ SCHEDULE // LgrCollection is used to initalize LGR Cells in the Eclipse Grid. eclipse_grid_OPM.init_lgr_cells(lgr_col); // LGR COORD and ZCORN is parsed to EclipseGridLGR children cell. (Simulates the process of recieving the LGR refinement.) - eclipse_grid_OPM.getLGRCell(0).set_lgr_refinement(coord_l_opm,zcorn_l_opm); // Intialize host_cell numbering. eclipse_grid_OPM.init_children_host_cells(); BOOST_CHECK_EQUAL( coord_g_opm.size() , coord_g.size()); @@ -395,8 +464,6 @@ SCHEDULE // LgrCollection is used to initalize LGR Cells in the Eclipse Grid. eclipse_grid_file.init_lgr_cells(lgr_col); // LGR COORD and ZCORN is parsed to EclipseGridLGR children cell. (Simulates the process of recieving the LGR refinement.) - eclipse_grid_file.getLGRCell(1).set_lgr_refinement(coord_l1,zcorn_l1); - eclipse_grid_file.getLGRCell(0).set_lgr_refinement(coord_l2,zcorn_l2); // Intialize host_cell numbering. eclipse_grid_file.init_children_host_cells(); // Save EclipseGrid. @@ -410,10 +477,6 @@ SCHEDULE Opm::EclipseGrid eclipse_grid_OPM(global_grid_dim, coord_g_opm, zcorn_g_opm); // LgrCollection is used to initalize LGR Cells in the Eclipse Grid. eclipse_grid_OPM.init_lgr_cells(lgr_col); - // LGR COORD and ZCORN is parsed to EclipseGridLGR children cell. (Simulates the process of recieving the LGR refinement.) - eclipse_grid_OPM.getLGRCell(1).set_lgr_refinement(coord_l1_opm,zcorn_l1_opm); - eclipse_grid_OPM.getLGRCell(0).set_lgr_refinement(coord_l2_opm,zcorn_l2_opm); - // Intialize host_cell numbering. eclipse_grid_OPM.init_children_host_cells(); BOOST_CHECK_EQUAL( coord_g_opm.size() , coord_g.size()); @@ -509,13 +572,8 @@ SCHEDULE // LGR COORD and ZCORN is parsed to EclipseGridLGR children cell. (Simulates the process of recieving the LGR refinement.) const auto& lgr1 = eclipse_grid_file.getLGRCell("LGR1"); - const auto& coordlgr1 = lgr1.getCOORD(); const auto& lgr2 = eclipse_grid_file.getLGRCell("LGR2"); - const auto& coordlgr2 = lgr2.getCOORD(); - - eclipse_grid_file.set_lgr_refinement("LGR1",coordlgr1,zcorn_l1); - eclipse_grid_file.set_lgr_refinement("LGR2",coordlgr2,zcorn_l2); // Intialize host_cell numbering. eclipse_grid_file.init_children_host_cells(); // Save EclipseGrid. @@ -598,3 +656,560 @@ SCHEDULE BOOST_CHECK_EQUAL_COLLECTIONS(host_vec.begin(), host_vec.end(), host_vec_sol.begin(), host_vec_sol.end()); } + + + BOOST_AUTO_TEST_CASE(TestLGRDepth) { + const std::string deck_string = R"( +RUNSPEC +TITLE + SPE1 - CASE 1 +DIMENS + 10 10 3 / + +EQLDIMS +/ +TABDIMS +/ +OIL +GAS +WATER +DISGAS +FIELD +START + 1 'JAN' 2015 / +WELLDIMS + 2 1 1 2 / +UNIFOUT +GRID +CARFIN +'LGR1' 5 6 5 6 1 3 6 6 9 / +ENDFIN +CARFIN +'LGR2' 7 8 7 8 1 3 6 6 9 / +ENDFIN +INIT +DX + 300*1000 / +DY + 300*1000 / +DZ + 100*20 100*30 100*50 / +TOPS + 100*8325 / +PORO + 300*0.3 / +PERMX + 100*500 100*50 100*200 / +PERMY + 100*500 100*50 100*200 / +PERMZ + 100*500 100*50 100*200 / +ECHO +PROPS +PVTW + 4017.55 1.038 3.22E-6 0.318 0.0 / +ROCK + 14.7 3E-6 / +SWOF +0.12 0 1 0 +0.18 4.64876033057851E-008 1 0 +0.24 0.000000186 0.997 0 +0.3 4.18388429752066E-007 0.98 0 +0.36 7.43801652892562E-007 0.7 0 +0.42 1.16219008264463E-006 0.35 0 +0.48 1.67355371900826E-006 0.2 0 +0.54 2.27789256198347E-006 0.09 0 +0.6 2.97520661157025E-006 0.021 0 +0.66 3.7654958677686E-006 0.01 0 +0.72 4.64876033057851E-006 0.001 0 +0.78 0.000005625 0.0001 0 +0.84 6.69421487603306E-006 0 0 +0.91 8.05914256198347E-006 0 0 +1 0.00001 0 0 / +SGOF +0 0 1 0 +0.001 0 1 0 +0.02 0 0.997 0 +0.05 0.005 0.980 0 +0.12 0.025 0.700 0 +0.2 0.075 0.350 0 +0.25 0.125 0.200 0 +0.3 0.190 0.090 0 +0.4 0.410 0.021 0 +0.45 0.60 0.010 0 +0.5 0.72 0.001 0 +0.6 0.87 0.0001 0 +0.7 0.94 0.000 0 +0.85 0.98 0.000 0 +0.88 0.984 0.000 0 / +DENSITY + 53.66 64.49 0.0533 / +PVDG +14.700 166.666 0.008000 +264.70 12.0930 0.009600 +514.70 6.27400 0.011200 +1014.7 3.19700 0.014000 +2014.7 1.61400 0.018900 +2514.7 1.29400 0.020800 +3014.7 1.08000 0.022800 +4014.7 0.81100 0.026800 +5014.7 0.64900 0.030900 +9014.7 0.38600 0.047000 / +PVTO +0.0010 14.7 1.0620 1.0400 / +0.0905 264.7 1.1500 0.9750 / +0.1800 514.7 1.2070 0.9100 / +0.3710 1014.7 1.2950 0.8300 / +0.6360 2014.7 1.4350 0.6950 / +0.7750 2514.7 1.5000 0.6410 / +0.9300 3014.7 1.5650 0.5940 / +1.2700 4014.7 1.6950 0.5100 + 9014.7 1.5790 0.7400 / +1.6180 5014.7 1.8270 0.4490 + 9014.7 1.7370 0.6310 / +/ +SOLUTION +EQUIL + 8400 4800 8450 0 8300 0 1 0 0 / +RSVD +8300 1.270 +8450 1.270 / +SUMMARY +FOPR +WGOR + 'PROD' +/ +FGOR +BPR +1 1 1 / +10 10 3 / +/ +BGSAT +1 1 1 / +1 1 2 / +1 1 3 / +10 1 1 / +10 1 2 / +10 1 3 / +10 10 1 / +10 10 2 / +10 10 3 / +/ +WBHP + 'INJ' + 'PROD' +/ +WGIR + 'INJ' + 'PROD' +/ +WGIT + 'INJ' + 'PROD' +/ +WGPR + 'INJ' + 'PROD' +/ +WGPT + 'INJ' + 'PROD' +/ +WOIR + 'INJ' + 'PROD' +/ +WOIT + 'INJ' + 'PROD' +/ +WOPR + 'INJ' + 'PROD' +/ +WOPT + 'INJ' + 'PROD' +/ +WWIR + 'INJ' + 'PROD' +/ +WWIT + 'INJ' + 'PROD' +/ +WWPR + 'INJ' + 'PROD' +/ +WWPT + 'INJ' + 'PROD' +/ +SCHEDULE +RPTSCHED + 'PRES' 'SGAS' 'RS' 'WELLS' / +RPTRST + 'BASIC=1' / +DRSDT + 0 / +WELSPECS + 'PROD' 'G1' 10 10 8400 'OIL' / + 'INJ' 'G1' 1 1 8335 'GAS' / +/ +COMPDAT + 'PROD' 10 10 3 3 'OPEN' 1* 1* 0.5 / + 'INJ' 1 1 1 1 'OPEN' 1* 1* 0.5 / +/ +WCONPROD + 'PROD' 'OPEN' 'ORAT' 20000 4* 1000 / +/ +WCONINJE + 'INJ' 'GAS' 'OPEN' 'RATE' 100000 1* 9014 / +/ +TSTEP +31 28 31 30 31 30 31 31 30 31 30 31 +/ +END +)"; + Opm::UnitSystem units(1); + std::vector vecNNC; + std::array global_grid_dim = {10,10,3}; + std::vector coord_g, zcorn_g, coord_l1, zcorn_l1, coord_l2, zcorn_l2; + // Intialize LgrCollection from string. + LgrCollection lgr_col = read_lgr(deck_string,global_grid_dim[0],global_grid_dim[1],global_grid_dim[2]); + + // Eclipse Grid is intialzied with COORD and ZCORN. + Opm::EclipseGrid eclipse_grid_file(init_deck(deck_string)); + // LgrCollection is used to initalize LGR Cells in the Eclipse Grid. + eclipse_grid_file.init_lgr_cells(lgr_col); + + eclipse_grid_file.save("SPECASE1_CARFIN_TEST.EGRID",false,vecNNC,units); + auto tol = 1e-6; + + // Coarse grid cells + auto l1_cell = eclipse_grid_file.getCellDims(0,0,0); + auto l2_cell = eclipse_grid_file.getCellDims(0,0,1); + auto l3_cell = eclipse_grid_file.getCellDims(0,0,2); + + // LGR grids + auto lgr1 = eclipse_grid_file.getLGRCell("LGR1"); + auto lgr2 = eclipse_grid_file.getLGRCell("LGR2"); + + // --- Level 1 cells in LGR1 --- + auto lgr1_l1_cell_01 = lgr1.getCellDims(0,0,0); + auto lgr1_l1_cell_02 = lgr1.getCellDims(0,0,1); + auto lgr1_l1_cell_03 = lgr1.getCellDims(0,0,2); + + // Check that Z-dimensions are consistent in LGR1 Level 1 + BOOST_CHECK_CLOSE(lgr1_l1_cell_01[2], lgr1_l1_cell_02[2], tol); + BOOST_CHECK_CLOSE(lgr1_l1_cell_02[2], lgr1_l1_cell_03[2], tol); + BOOST_CHECK_CLOSE(lgr1_l1_cell_03[2], lgr1_l1_cell_01[2], tol); + + // --- Level 2 cells in LGR1 --- + auto lgr1_l2_cell_01 = lgr1.getCellDims(0,0,3); + auto lgr1_l2_cell_02 = lgr1.getCellDims(0,0,4); + auto lgr1_l2_cell_03 = lgr1.getCellDims(0,0,5); + + // Check that Z-dimensions are consistent in LGR1 Level 2 + BOOST_CHECK_CLOSE(lgr1_l2_cell_01[2], lgr1_l2_cell_02[2], tol); + BOOST_CHECK_CLOSE(lgr1_l2_cell_02[2], lgr1_l2_cell_03[2], tol); + BOOST_CHECK_CLOSE(lgr1_l2_cell_03[2], lgr1_l2_cell_01[2], tol); + + // --- Level 3 cells in LGR1 --- + auto lgr1_l3_cell_01 = lgr1.getCellDims(0,0,6); + auto lgr1_l3_cell_02 = lgr1.getCellDims(0,0,7); + auto lgr1_l3_cell_03 = lgr1.getCellDims(0,0,8); + + // Check that Z-dimensions are consistent between LGR1 Level 3 + BOOST_CHECK_CLOSE(lgr1_l3_cell_01[2], lgr1_l3_cell_02[2], tol); + BOOST_CHECK_CLOSE(lgr1_l3_cell_02[2], lgr1_l3_cell_03[2], tol); + BOOST_CHECK_CLOSE(lgr1_l3_cell_03[2], lgr1_l3_cell_01[2], tol); + + // CHECKING IF LGR1 and LGR2 LEVEL 1 CELLS ARE THE SAME + auto lgr2_l1_cell_01 = lgr2.getCellDims(0,0,0); + auto lgr2_l1_cell_02 = lgr2.getCellDims(0,0,1); + auto lgr2_l1_cell_03 = lgr2.getCellDims(0,0,2); + BOOST_CHECK_CLOSE(lgr1_l1_cell_01[2], lgr2_l1_cell_01[2], tol); + BOOST_CHECK_CLOSE(lgr1_l1_cell_02[2], lgr2_l1_cell_02[2], tol); + BOOST_CHECK_CLOSE(lgr1_l1_cell_03[2], lgr2_l1_cell_03[2], tol); + + // --- Level 2 cells in LGR2 --- + auto lgr2_l2_cell_01 = lgr2.getCellDims(0,0,3); + auto lgr2_l2_cell_02 = lgr2.getCellDims(0,0,4); + auto lgr2_l2_cell_03 = lgr2.getCellDims(0,0,5); + BOOST_CHECK_CLOSE(lgr1_l2_cell_01[2], lgr2_l2_cell_01[2], tol); + BOOST_CHECK_CLOSE(lgr1_l2_cell_02[2], lgr2_l2_cell_02[2], tol); + BOOST_CHECK_CLOSE(lgr1_l2_cell_03[2], lgr2_l2_cell_03[2], tol); + + // --- Level 3 cells in LGR2 --- + auto lgr2_l3_cell_01 = lgr2.getCellDims(0,0,6); + auto lgr2_l3_cell_02 = lgr2.getCellDims(0,0,7); + auto lgr2_l3_cell_03 = lgr2.getCellDims(0,0,8); + BOOST_CHECK_CLOSE(lgr1_l3_cell_01[2], lgr2_l3_cell_01[2], tol); + BOOST_CHECK_CLOSE(lgr1_l3_cell_02[2], lgr2_l3_cell_02[2], tol); + BOOST_CHECK_CLOSE(lgr1_l3_cell_03[2], lgr2_l3_cell_03[2], tol); + + // compare levels and coarse cell + BOOST_CHECK_CLOSE(l1_cell[2], lgr1_l1_cell_01[2] + lgr1_l1_cell_02[2] + lgr1_l1_cell_03[2], tol); + BOOST_CHECK_CLOSE(l2_cell[2], lgr1_l2_cell_01[2] + lgr1_l2_cell_02[2] + lgr1_l2_cell_03[2], tol); + BOOST_CHECK_CLOSE(l3_cell[2], lgr1_l3_cell_01[2] + lgr1_l3_cell_02[2] + lgr1_l3_cell_03[2], tol); + + } + + + BOOST_AUTO_TEST_CASE(TestFinalDepth) { + const std::string deck_string = R"( +RUNSPEC +TITLE + SPE1 - CASE 1 +DIMENS + 2 2 5 / +EQLDIMS +/ +TABDIMS +/ +OIL +GAS +WATER +DISGAS +FIELD +START + 1 'JAN' 2015 / +WELLDIMS + 2 1 1 2 / +UNIFOUT +GRID +CARFIN +'LGR1' 1 2 1 1 1 5 6 6 20 / +ENDFIN +INIT +DX + 20*1000 / +DY + 20*1000 / +DZ + 4*10 4*20 4*50 4*60 4*70 / +TOPS + 4*8325 / +PORO + 20*0.3 / +PERMX + 20*500/ +PERMY + 20*500/ +PERMZ + 20*500 / +ECHO +PROPS +PVTW + 4017.55 1.038 3.22E-6 0.318 0.0 / +ROCK + 14.7 3E-6 / +SWOF +0.12 0 1 0 +0.18 4.64876033057851E-008 1 0 +0.24 0.000000186 0.997 0 +0.3 4.18388429752066E-007 0.98 0 +0.36 7.43801652892562E-007 0.7 0 +0.42 1.16219008264463E-006 0.35 0 +0.48 1.67355371900826E-006 0.2 0 +0.54 2.27789256198347E-006 0.09 0 +0.6 2.97520661157025E-006 0.021 0 +0.66 3.7654958677686E-006 0.01 0 +0.72 4.64876033057851E-006 0.001 0 +0.78 0.000005625 0.0001 0 +0.84 6.69421487603306E-006 0 0 +0.91 8.05914256198347E-006 0 0 +1 0.00001 0 0 / +SGOF +0 0 1 0 +0.001 0 1 0 +0.02 0 0.997 0 +0.05 0.005 0.980 0 +0.12 0.025 0.700 0 +0.2 0.075 0.350 0 +0.25 0.125 0.200 0 +0.3 0.190 0.090 0 +0.4 0.410 0.021 0 +0.45 0.60 0.010 0 +0.5 0.72 0.001 0 +0.6 0.87 0.0001 0 +0.7 0.94 0.000 0 +0.85 0.98 0.000 0 +0.88 0.984 0.000 0 / +DENSITY + 53.66 64.49 0.0533 / +PVDG +14.700 166.666 0.008000 +264.70 12.0930 0.009600 +514.70 6.27400 0.011200 +1014.7 3.19700 0.014000 +2014.7 1.61400 0.018900 +2514.7 1.29400 0.020800 +3014.7 1.08000 0.022800 +4014.7 0.81100 0.026800 +5014.7 0.64900 0.030900 +9014.7 0.38600 0.047000 / +PVTO +0.0010 14.7 1.0620 1.0400 / +0.0905 264.7 1.1500 0.9750 / +0.1800 514.7 1.2070 0.9100 / +0.3710 1014.7 1.2950 0.8300 / +0.6360 2014.7 1.4350 0.6950 / +0.7750 2514.7 1.5000 0.6410 / +0.9300 3014.7 1.5650 0.5940 / +1.2700 4014.7 1.6950 0.5100 + 9014.7 1.5790 0.7400 / +1.6180 5014.7 1.8270 0.4490 + 9014.7 1.7370 0.6310 / +/ +SOLUTION +EQUIL + 8400 4800 8450 0 8300 0 1 0 0 / +RSVD +8300 1.270 +8450 1.270 / +SUMMARY +FOPR +WGOR + 'PROD' +/ +FGOR +BPR +1 1 1 / +2 2 3 / +/ +BGSAT +1 1 1 / +1 1 2 / +1 1 3 / +2 1 1 / +2 1 2 / +2 1 3 / +2 2 1 / +2 2 2 / +2 2 3 / +/ +WBHP + 'INJ' + 'PROD' +/ +WGIR + 'INJ' + 'PROD' +/ +WGIT + 'INJ' + 'PROD' +/ +WGPR + 'INJ' + 'PROD' +/ +WGPT + 'INJ' + 'PROD' +/ +WOIR + 'INJ' + 'PROD' +/ +WOIT + 'INJ' + 'PROD' +/ +WOPR + 'INJ' + 'PROD' +/ +WOPT + 'INJ' + 'PROD' +/ +WWIR + 'INJ' + 'PROD' +/ +WWIT + 'INJ' + 'PROD' +/ +WWPR + 'INJ' + 'PROD' +/ +WWPT + 'INJ' + 'PROD' +/ +SCHEDULE +RPTSCHED + 'PRES' 'SGAS' 'RS' 'WELLS' / +RPTRST + 'BASIC=1' / +DRSDT + 0 / +WELSPECS + 'PROD' 'G1' 2 2 8400 'OIL' / + 'INJ' 'G1' 1 2 8335 'GAS' / +/ +COMPDAT + 'PROD' 2 2 1 1 'OPEN' 1* 1* 0.5 / + 'INJ' 1 2 1 1 'OPEN' 1* 1* 0.5 / +/ +WCONPROD + 'PROD' 'OPEN' 'ORAT' 20000 4* 1000 / +/ +WCONINJE + 'INJ' 'GAS' 'OPEN' 'RATE' 100000 1* 9014 / +/ +TSTEP +31 28 31 30 31 30 31 31 30 31 30 31 +/ +END +)"; + Opm::UnitSystem units(1); + std::vector vecNNC; + std::array global_grid_dim = {2,2,5}; + std::vector coord_g, zcorn_g, coord_l1, zcorn_l1, coord_l2, zcorn_l2; + // Intialize LgrCollection from string. + LgrCollection lgr_col = read_lgr(deck_string,global_grid_dim[0],global_grid_dim[1],global_grid_dim[2]); + + // Eclipse Grid is intialzied with COORD and ZCORN. + Opm::EclipseGrid eclipse_grid_file(init_deck(deck_string)); + // LgrCollection is used to initalize LGR Cells in the Eclipse Grid. + eclipse_grid_file.init_lgr_cells(lgr_col); + + eclipse_grid_file.save("SPECASE1_CARFIN_TEST_SMALL.EGRID",false,vecNNC,units); + + // LGR1 grids + const auto& lgr1 = eclipse_grid_file.getLGRCell("LGR1"); + const auto lgr1_coord = lgr1.getCOORD(); + const auto lgr1_zcorn = lgr1.getZCORN(); + + auto [expected_coord, expected_zcorn] = final_test_data(); + check_vec_close(lgr1_coord, expected_coord, 1e-2); + check_vec_close(lgr1_zcorn, expected_zcorn, 1e-2); + + const auto l1_cell = eclipse_grid_file.getCellDims(0,0,0); + const auto lgr1_l1_cell = lgr1.getCellDims(0,0,0); + BOOST_CHECK_CLOSE(l1_cell[2], 4*lgr1_l1_cell[2], 1e-6); + + const auto l2_cell = eclipse_grid_file.getCellDims(0,0,1); + const auto lgr1_l2_cell = lgr1.getCellDims(0,0,5); + BOOST_CHECK_CLOSE(l2_cell[2], 4*lgr1_l2_cell[2], 1e-6); + + const auto l3_cell = eclipse_grid_file.getCellDims(0,0,2); + const auto lgr1_l3_cell = lgr1.getCellDims(0,0,9); + BOOST_CHECK_CLOSE(l3_cell[2], 4*lgr1_l3_cell[2], 1e-6); + + const auto l4_cell = eclipse_grid_file.getCellDims(0,0,3); + const auto lgr1_l4_cell = lgr1.getCellDims(0,0,13); + BOOST_CHECK_CLOSE(l4_cell[2], 4*lgr1_l4_cell[2], 1e-6); + + const auto l5_cell = eclipse_grid_file.getCellDims(0,0,4); + const auto lgr1_l5_cell = lgr1.getCellDims(0,0,17); + BOOST_CHECK_CLOSE(l5_cell[2], 4*lgr1_l5_cell[2], 1e-6); + }