Skip to content

Commit

Permalink
Take into account inactive cells, change name local->lgrIjk, add test…
Browse files Browse the repository at this point in the history
… lgr with inac cell
  • Loading branch information
aritorto committed Jan 31, 2025
1 parent 5b01e0a commit bd231ba
Show file tree
Hide file tree
Showing 4 changed files with 190 additions and 31 deletions.
2 changes: 1 addition & 1 deletion CMakeLists_files.cmake
Original file line number Diff line number Diff line change
Expand Up @@ -116,7 +116,7 @@ if(Boost_VERSION_STRING VERSION_GREATER 1.53)
tests/cpgrid/grid_lgr_test.cpp
tests/cpgrid/inactiveCell_lgr_test.cpp
tests/cpgrid/lgr_cartesian_idx_test.cpp
tests/cpgrid/localIJK_test.cpp
tests/cpgrid/lgrIJK_test.cpp
tests/cpgrid/lookUpCellCentroid_cpgrid_test.cpp
tests/cpgrid/lookupdataCpGrid_test.cpp
tests/cpgrid/replace_lgr1_corner_idx_by_lgr2_corner_idx_test.cpp
Expand Down
40 changes: 28 additions & 12 deletions opm/grid/cpgrid/CpGridUtilities.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -26,7 +26,7 @@

namespace Opm {

std::vector<std::array<int,3>> localIJK(const Dune::CpGrid& grid, const std::string& lgr_name)
std::vector<std::array<int,3>> lgrIJK(const Dune::CpGrid& grid, const std::string& lgr_name)
{
// Check lgr_name exists in lgr_names_
const auto& lgr_names = grid.getLgrNameToLevel();
Expand All @@ -35,23 +35,39 @@ std::vector<std::array<int,3>> localIJK(const Dune::CpGrid& grid, const std::str
OPM_THROW(std::runtime_error, "LGR name not found: " + lgr_name);
}

const Opm::LevelCartesianIndexMapper<Dune::CpGrid> levelCartMapp(grid);
const Opm::LevelCartesianIndexMapper<Dune::CpGrid> levelCartMapper(grid);
const auto& level = it->second;

std::vector<std::array<int, 3>> localIJK;
localIJK.resize(grid.levelGridView(level).size(0));
// Determine the logical Cartesian size of the LGR and total cells (including inactive ones)
const auto lgr_dim = grid.currentData()[level]->logicalCartesianSize();
const int lgr_cells = lgr_dim[0] * lgr_dim[1] * lgr_dim[2]; // (including inactive ones)
// Actual size is given by grid.levelGridView(level).size(0)

// To be improved
const auto& leaf_to_lgr_idx = grid.mapLeafIndexSetToLocalCartesianIndexSets();
// Initialize level IJKs with inactive values
static constexpr std::array<int, 3> inactiveIJK = { -1, -1, -1 };
std::vector<std::array<int, 3>> lgrIJK(lgr_cells);
lgrIJK.assign(lgr_cells, inactiveIJK); // Ensures all elements are set to inactive

// Iterate over active elements in the grid. Rewrite ijks at the specified refined level active cells.
//
// Note: to avoid go over all the leaf elements and replace this approach with the
// corresponding grid.levelGridView(level), a map relating level_to_lgr_idx is needed.
for (const auto& element : Dune::elements(grid.leafGridView())) {
if (element.level() == level)
{
std::array<int, 3> local_ijk;
levelCartMapp.cartesianCoordinate( element.getLevelElem().index(), local_ijk, level);
localIJK[ leaf_to_lgr_idx[element.index()][1] ] = local_ijk;
if (element.level() != level) {
continue;
}

std::array<int, 3> ijk;
levelCartMapper.cartesianCoordinate(element.getLevelElem().index(), ijk, level);
// In general, element.getLevelElem().index() != element.getLevelCartesianIdx().
// - element.getLevelElem().index() : integer between 0 and "total active cells - 1" on the level grid.
// - element.getLevelCartesianIdx() : integer between 0 and "lgr_cells - 1", representing an index of
// the underlying Cartesian grid for the level (with/without INACTIVE cells).

const int lgr_cart_index = element.getLevelCartesianIdx();
lgrIJK[lgr_cart_index] = ijk;
}
return localIJK;
return lgrIJK;
}

}
8 changes: 5 additions & 3 deletions opm/grid/cpgrid/CpGridUtilities.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -25,14 +25,16 @@
namespace Opm
{

/// @brief Retrieves the local IJKs of cells within a specified LGR.
/// @brief Retrieves the local grid refinement (LGR) IJKs of cells.
///
/// This method returns a vector of IJK indices corresponding to underlying Cartesian grid of the local
/// grid refinement (LGR) specified by its name. If the specified LGR name is not found, it throws.
/// grid refinement (LGR) specified by its name. If the specified name is not found, an exception is thrown.
/// For inactive cells—i.e., non-existing child cells of inactive coarse cells defined by the CARFIN keyword—
/// the corresponding IJK index is set to {-1, -1, -1}.
///
/// @param lgr_name The name of the LGR whose local IJK coordinates are requested.
/// @return std::vector<std::array<int, 3>> A list of (i, j, k)'s for each cell in the LGR.
std::vector<std::array<int,3>> localIJK(const Dune::CpGrid& grid, const std::string& lgr_name);
std::vector<std::array<int,3>> lgrIJK(const Dune::CpGrid& grid, const std::string& lgr_name);

}

Expand Down
171 changes: 156 additions & 15 deletions tests/cpgrid/localIJK_test.cpp → tests/cpgrid/lgrIJK_test.cpp
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
//===========================================================================
//
// File: localIJK_test.cpp
// File: lgrIJK_test.cpp
//
// Created: Thursday 30.01.2025 08:17:00
//
Expand Down Expand Up @@ -31,7 +31,7 @@
*/
#include "config.h"

#define BOOST_TEST_MODULE LocalIJKTests
#define BOOST_TEST_MODULE LgrIJKTests
#include <boost/test/unit_test.hpp>
#include <boost/version.hpp>
#if BOOST_VERSION / 100000 == 1 && BOOST_VERSION / 100 % 1000 < 71
Expand Down Expand Up @@ -67,7 +67,149 @@ struct Fixture

BOOST_GLOBAL_FIXTURE(Fixture);

BOOST_AUTO_TEST_CASE(localIJK_in_allActiveCellLGRs)
/* Tests lgrIJK() for a grid containing inactive cells within the LGR block.
The test passes if initLgr() is not invoked in opm-common/opm/input/eclipse/EclipseState/EclipseState.cpp.
Currently, inactive cells in LGRs are not supported in "opm-common."
BOOST_AUTO_TEST_CASE(inactive_cell_within_lgr)
{
const std::string deck_string = R"(
RUNSPEC
DIMENS
3 3 1 /
GRID
CARFIN
-- NAME I1-I2 J1-J2 K1-K2 NX NY NZ
'LGR1' 2 3 2 3 1 1 6 6 3/
ENDFIN
DX
9*1000 /
DY
9*1000 /
DZ
9*20 /
TOPS
9*8325 /
ACTNUM
1 1 1
1 1 1
1 1 0
/
PORO
9*0.15 /
PERMX
9*1 /
COPY
PERMX PERMZ /
PERMX PERMY /
/
EDIT
OIL
GAS
TITLE
The title
START
16 JUN 1988 /
PROPS
REGIONS
SOLUTION
SCHEDULE
)";
Opm::Parser parser;
const auto deck = parser.parseString(deck_string);
Opm::EclipseState ecl_state(deck);
Opm::EclipseGrid eclipse_grid = ecl_state.getInputGrid();
Dune::CpGrid grid;
grid.processEclipseFormat(&eclipse_grid, &ecl_state, false, false, false);
// Add LGR and update grid view
//
// LGR1: 3x3x1 child cells in x-,y-, and z-direction per ACTIVE parent cell
const std::vector<std::array<int, 3>> cells_per_dim_vec = { {3, 3, 3} };
// LGR1 starts at (1,1,0) in coarse grid - equivalent to (I1-1, J1-1, K1-1) from its CARFIN block
const std::vector<std::array<int, 3>> startIJK_vec = { {1, 1, 0} };
// LGR1 ends at (3,3,1) in coarse grid - equivalent to (I2, J2, K2) from its CARFIN block
const std::vector<std::array<int, 3>> endIJK_vec = { {3, 3, 1} };
const std::vector<std::string> lgr_name_vec = {"LGR1"};
grid.addLgrsUpdateLeafView(cells_per_dim_vec, startIJK_vec, endIJK_vec, lgr_name_vec);
// Get LGR levels
const auto& lgr_name_to_level = grid.getLgrNameToLevel(); // map std::string to int level: {"GLOBAL", 0}, {"LGR1", 1}, ...
const int lgr1_level = lgr_name_to_level.at("LGR1");
// Get vector of local ijk for each LGR
const std::vector<std::array<int, 3>> ijk_lgr1 = Opm::lgrIJK(grid, "LGR1");
// Note: k = 0, indicating a single-layer grid with dimensions 3x3x1.
// ACTNUM represents the active cell indicator for the parent grid of the LGR (Local Grid Refinement).
// The grid structure is as follows, where '1' denotes an active cell and '0' denotes an inactive cell:
//
// 1 1 1
// 1 1 1
// 1 1 0
//
// The corresponding LGR block parent cell representation appears as:
// 1 1
// 1 0
// Verify the size matches expected elements
BOOST_TEST(grid.levelGridView(lgr1_level).size(0) == 81); // 3 ACTIVE parent cellS into 3x3x3 children each -> 81
BOOST_TEST(ijk_lgr1.size() == 108); // Block size - inlcuding inactive cells (6x6x3)
// LGR1 dimension 6x6x3
// Visual representation of level Cartesian cell indices per k-layer:
//
// k = 2 | 102 103 104 |
// | 96 97 98 | INACTIVE
// | 90 91 92 |
// ----------------------------
// | 84 85 86 | 87 88 89
// | 78 79 80 | 81 82 83
// | 72 73 74 | 75 76 77
//------------------------------------
// k = 1 | 66 67 68 |
// | 60 61 62 | INACTIVE
// | 54 55 56 |
// ----------------------------
// | 48 49 50 | 51 52 53
// | 42 43 44 | 45 46 47
// | 36 37 38 | 39 40 41
//------------------------------------
// k = 0 | 30 31 32 |
// | 24 25 26 | INACTIVE
// | 18 19 20 |
// ----------------------------
// | 12 13 14 | 15 16 17
// | 6 7 8 | 9 10 11
// | 0 1 2 | 3 4 5
// Verify that cell with local (level Cartesian) index 97 in LGR1 has local ijk = {1,4,2}
BOOST_TEST(ijk_lgr1[97][0] == 1);
BOOST_TEST(ijk_lgr1[97][1] == 4);
BOOST_TEST(ijk_lgr1[97][2] == 2);
// Verify that cell with local (level Cartesian) index 88 in LGR1 has local ijk = {4,2,2}
BOOST_TEST(ijk_lgr1[88][0] == 4);
BOOST_TEST(ijk_lgr1[88][1] == 2);
BOOST_TEST(ijk_lgr1[88][2] == 2);
// Accessing inactive (non-existing) cells
for (int i = 0; i < 3; ++i) {
BOOST_TEST(ijk_lgr1[107][i] == -1);
BOOST_TEST(ijk_lgr1[60][i] == -1);
}
// Potential improvement: If lgrIJK() returns a map instead of vector, then
// ijk_lgr1.at( non-exisiting-inactive-cell-cartesian-index ) throws, which
// might be better than accesing an invalid value {-1,-1,-1}.
}
*/

BOOST_AUTO_TEST_CASE(all_active_cells_within_lgrs)
{

const std::string deckString =
Expand Down Expand Up @@ -145,8 +287,8 @@ BOOST_AUTO_TEST_CASE(localIJK_in_allActiveCellLGRs)
const int lgr2_level = lgr_name_to_level.at("LGR2");

// Get vector of local ijk for each LGR
std::vector<std::array<int, 3>> ijk_lgr1 = Opm::localIJK(grid, "LGR1");
std::vector<std::array<int, 3>> ijk_lgr2 = Opm::localIJK(grid, "LGR2");
std::vector<std::array<int, 3>> ijk_lgr1 = Opm::lgrIJK(grid, "LGR1");
std::vector<std::array<int, 3>> ijk_lgr2 = Opm::lgrIJK(grid, "LGR2");

// Ensure it's not empty
BOOST_TEST(!ijk_lgr1.empty()); // 2 active parent cell refined into 2x2x2 refined cells
Expand All @@ -168,7 +310,7 @@ BOOST_AUTO_TEST_CASE(localIJK_in_allActiveCellLGRs)
}

// Invalid LGR should throw an exception
BOOST_CHECK_THROW(Opm::localIJK(grid, "LGR3DOESNOTEXIST"), std::runtime_error);
BOOST_CHECK_THROW(Opm::lgrIJK(grid, "LGR3DOESNOTEXIST"), std::runtime_error);

// LGR1 Dimension: 2x2x4 (Width x Height x Depth)
// Visual representation of level Cartesian cell indices per k-layer:
Expand Down Expand Up @@ -269,8 +411,8 @@ BOOST_AUTO_TEST_CASE(internal_order_differs_from_underlying_cartesian_level_grid
const int lgr2_level = lgr_name_to_level.at("LGR2");

// Get vector of local ijk for each LGR
std::vector<std::array<int, 3>> ijk_lgr1 = Opm::localIJK(grid, "LGR1");
std::vector<std::array<int, 3>> ijk_lgr2 = Opm::localIJK(grid, "LGR2");
std::vector<std::array<int, 3>> ijk_lgr1 = Opm::lgrIJK(grid, "LGR1");
std::vector<std::array<int, 3>> ijk_lgr2 = Opm::lgrIJK(grid, "LGR2");

// Verify the size matches expected elements
BOOST_TEST(ijk_lgr1.size() == grid.levelGridView(lgr1_level).size(0)); // 4 parent cells into 3x3x3 children each -> 108
Expand Down Expand Up @@ -310,12 +452,12 @@ BOOST_AUTO_TEST_CASE(internal_order_differs_from_underlying_cartesian_level_grid
BOOST_TEST(ijk_lgr1[25][1] == 4);
BOOST_TEST(ijk_lgr1[25][2] == 0);

// Verify that cell with local (level Cartesian) index 64 in LGR1 has local ijk = {4,4,1}
// Verify that cell with local (level Cartesian) index 64 in LGR1 has local ijk = {4,4,1}
BOOST_TEST(ijk_lgr1[64][0] == 4);
BOOST_TEST(ijk_lgr1[64][1] == 4);
BOOST_TEST(ijk_lgr1[64][2] == 1);

// Verify that cell with local (level Cartesian) index 98 in LGR1 has local ijk = {2,4,2}
// Verify that cell with local (level Cartesian) index 98 in LGR1 has local ijk = {2,4,2}
BOOST_TEST(ijk_lgr1[98][0] == 2);
BOOST_TEST(ijk_lgr1[98][1] == 4);
BOOST_TEST(ijk_lgr1[98][2] == 2);
Expand Down Expand Up @@ -355,20 +497,19 @@ BOOST_AUTO_TEST_CASE(internal_order_differs_from_underlying_cartesian_level_grid
// | 12 13 14 | 15 16 17
// | 6 7 8 | 9 10 11
// | 0 1 2 | 3 4 5
// Verify that cell with local (level Cartesian) index 139 in LGR2 has local ijk = {1,5,3}

// Verify that cell with local (level Cartesian) index 139 in LGR2 has local ijk = {1,5,3}
BOOST_TEST(ijk_lgr2[139][0] == 1);
BOOST_TEST(ijk_lgr2[139][1] == 5);
BOOST_TEST(ijk_lgr2[139][2] == 3);

// Verify that cell with local (level Cartesian) index 119 in LGR2 has local ijk = {5,1,3}
// Verify that cell with local (level Cartesian) index 119 in LGR2 has local ijk = {5,1,3}
BOOST_TEST(ijk_lgr2[119][0] == 5);
BOOST_TEST(ijk_lgr2[119][1] == 1);
BOOST_TEST(ijk_lgr2[119][2] == 3);

// Verify that cell with local (level Cartesian) index 121 in LGR1 has local ijk = {1,2,3}
// Verify that cell with local (level Cartesian) index 121 in LGR1 has local ijk = {1,2,3}
BOOST_TEST(ijk_lgr2[121][0] == 1);
BOOST_TEST(ijk_lgr2[121][1] == 2);
BOOST_TEST(ijk_lgr2[121][2] == 3);

}

0 comments on commit bd231ba

Please sign in to comment.