From 7161773df4dd44cb68bcc59d0ebd59a1b28f5410 Mon Sep 17 00:00:00 2001 From: Antonella Ritorto Date: Thu, 26 Sep 2024 12:15:59 +0200 Subject: [PATCH 1/7] GlobalCellLevel to Leaf Index defined and tested --- CMakeLists_files.cmake | 1 + opm/grid/CpGrid.hpp | 5 + opm/grid/common/CartesianIndexMapper.hpp | 20 ++- opm/grid/cpgrid/CartesianIndexMapper.hpp | 23 ++- opm/grid/cpgrid/CpGrid.cpp | 10 ++ opm/grid/cpgrid/CpGridData.hpp | 17 +- .../polyhedralgrid/cartesianindexmapper.hh | 21 ++- tests/cpgrid/lgr_cartesian_idx_test.cpp | 166 ++++++++++++++++++ 8 files changed, 243 insertions(+), 20 deletions(-) create mode 100644 tests/cpgrid/lgr_cartesian_idx_test.cpp diff --git a/CMakeLists_files.cmake b/CMakeLists_files.cmake index a8dec70b3..f58c8a78c 100644 --- a/CMakeLists_files.cmake +++ b/CMakeLists_files.cmake @@ -122,6 +122,7 @@ if(Boost_VERSION_STRING VERSION_GREATER 1.53) tests/cpgrid/geometry_test.cpp tests/cpgrid/grid_lgr_test.cpp tests/cpgrid/inactiveCell_lgr_test.cpp + tests/cpgrid/lgr_cartesian_idx_test.cpp tests/cpgrid/lookUpCellCentroid_cpgrid_test.cpp tests/cpgrid/lookupdataCpGrid_test.cpp tests/cpgrid/shifted_cart_test.cpp diff --git a/opm/grid/CpGrid.hpp b/opm/grid/CpGrid.hpp index c602d63fa..2ee2c52d2 100644 --- a/opm/grid/CpGrid.hpp +++ b/opm/grid/CpGrid.hpp @@ -968,6 +968,11 @@ namespace Dune /// For nested refinement, we lookup the oldest ancestor, from level zero. void computeGlobalCellLeafGridViewWithLgrs(std::vector& global_cell_leaf); + /// @brief Compute for each level grid, a map from the global_cell_[ cell index in level grid ] to the leaf index of the equivalent cell + /// on the leaf grid view. + /// Notice that cells that vanished and do not appear on the leaf grid view will not be considered. + void mapGlobalCellLevelToLeafIndexSet(); + /// @brief Get the ijk index of a refined corner, given its corner index of a single-cell-refinement. /// /// Given a single-cell, we refine it in {nx, ny, nz} refined children cells (per direction). Then, this single-cell-refinement diff --git a/opm/grid/common/CartesianIndexMapper.hpp b/opm/grid/common/CartesianIndexMapper.hpp index b34561042..17e4401d5 100644 --- a/opm/grid/common/CartesianIndexMapper.hpp +++ b/opm/grid/common/CartesianIndexMapper.hpp @@ -43,12 +43,6 @@ namespace Dune return 0; } - /** \brief return number of cells in the active level zero grid. Only relevant for CpGrid specialization. */ - int compressedLevelZeroSize() const - { - return 0; - } - /** \brief return index of the cells in the logical Cartesian grid */ int cartesianIndex( const int /* compressedElementIndex */) const { @@ -60,11 +54,25 @@ namespace Dune { } + /// Only relevant for CpGrid specialization. + /** \brief return number of cells in the active level zero grid. Only relevant for CpGrid specialization. */ + int compressedLevelZeroSize() const + { + return 0; + } + /** \brief return Cartesian coordinate, i.e. IJK, for a given cell. Only relevant for CpGrid specialization.*/ void cartesianCoordinateLevel(const int /* compressedElementIndexOnLevel */, std::array& /* coordsOnLevel */, int /*level*/) const { } + + /** \brief return index of the cells in the logical Cartesian grid, for refined level grids. Only relevant for CpGrid specialization. */ + int cartesianIndexLevel( const int /* compressedElementIndex */ , const int level) const + { + return 0; + } + /// Only relevant for CpGrid specialization. END }; } // end namespace Opm diff --git a/opm/grid/cpgrid/CartesianIndexMapper.hpp b/opm/grid/cpgrid/CartesianIndexMapper.hpp index 2c50cadb8..fbc3ed511 100644 --- a/opm/grid/cpgrid/CartesianIndexMapper.hpp +++ b/opm/grid/cpgrid/CartesianIndexMapper.hpp @@ -50,11 +50,6 @@ namespace Dune return grid_.globalCell().size(); } - int compressedLevelZeroSize() const - { - return (*grid_.currentData()[0]).size(0); - } - int cartesianIndex( const int compressedElementIndex ) const { assert( compressedElementIndex >= 0 && compressedElementIndex < compressedSize() ); @@ -66,13 +61,29 @@ namespace Dune grid_.getIJK( compressedElementIndex, coords ); } + /// Additional methods realted to LGRs + int compressedLevelZeroSize() const + { + return (*grid_.currentData()[0]).size(0); + } + + int cartesianIndexLevel( const int compressedElementIndex, const int level) const + { + if ((level < 0) || (level > grid_.maxLevel())) { + throw std::invalid_argument("Invalid level.\n"); + } + assert( compressedElementIndex >= 0 && compressedElementIndex < grid_.currentData()[level]->size(0) ); + return grid_.currentData()[level]->globalCell()[compressedElementIndex]; + } + void cartesianCoordinateLevel(const int compressedElementIndexOnLevel, std::array& coordsOnLevel, int level) const { if ((level < 0) || (level > grid_.maxLevel())) { throw std::invalid_argument("Invalid level.\n"); } - (*grid_.currentData()[level]).getIJK( compressedElementIndexOnLevel, coordsOnLevel); + grid_.currentData()[level]->getIJK( compressedElementIndexOnLevel, coordsOnLevel); } + /// Additional methods realted to LGRs. END }; } // end namespace Opm diff --git a/opm/grid/cpgrid/CpGrid.cpp b/opm/grid/cpgrid/CpGrid.cpp index d48b85105..0ae9ed4fb 100644 --- a/opm/grid/cpgrid/CpGrid.cpp +++ b/opm/grid/cpgrid/CpGrid.cpp @@ -679,6 +679,14 @@ void CpGrid::computeGlobalCellLeafGridViewWithLgrs(std::vector& global_cell } } +void CpGrid::mapGlobalCellLevelToLeafIndexSet() +{ + for (const auto& element : elements(leafGridView())) { + const auto& global_cell_level = currentData()[element.level()]->globalCell()[element.getEquivLevelElem().index()]; + (*currentData()[element.level()]).globalCellLevel_to_leafIdx_[global_cell_level] = element.index(); + } +} + void CpGrid::getIJK(const int c, std::array& ijk) const { current_view_data_->getIJK(c, ijk); @@ -1996,6 +2004,8 @@ bool CpGrid::adapt(const std::vector>& cells_per_dim_vec, computeGlobalCellLeafGridViewWithLgrs(global_cell_leaf); (*data[levels + preAdaptMaxLevel +1]).global_cell_.swap(global_cell_leaf); + mapGlobalCellLevelToLeafIndexSet(); + updateCornerHistoryLevels(cornerInMarkedElemWithEquivRefinedCorner, elemLgrAndElemLgrCorner_to_refinedLevelAndRefinedCorner, adaptedCorner_to_elemLgrAndElemLgrCorner, diff --git a/opm/grid/cpgrid/CpGridData.hpp b/opm/grid/cpgrid/CpGridData.hpp index 136da9a38..161de296b 100644 --- a/opm/grid/cpgrid/CpGridData.hpp +++ b/opm/grid/cpgrid/CpGridData.hpp @@ -306,8 +306,13 @@ class CpGridData ijk = getIJK(global_cell_[c], logical_cartesian_size_); } + const std::vector& globalCell() const + { + return global_cell_; + } + /// @brief Extract Cartesian index triplet (i,j,k) given an index between 0 and NXxNYxNZ -1 - /// where NX, NY, and NZ is the total amoung of cells in each direction x-,y-,and z- respectively. + /// where NX, NY, and NZ is the total amoung of cells in each direction x-,y-,and z- respectively. /// /// @param [in] idx Integer between 0 and cells_per_dim[0]*cells_per_dim[1]*cells_per_dim[2]-1 /// @param [in] cells_per_dim @@ -627,6 +632,16 @@ class CpGridData return logical_cartesian_size_; } + /// Get map from global_cell_ to leaf index set + /// + /// Only relevant for CpGrid with LGRs. The maps of + /// all refined level grids is defined in CpGrid::adapt(/*..*/), + /// via CpGrid::mapGlobalCellLevelToLeafIndexSet(). + const std::unordered_map& getGlobalCellLevelToLeafIndexSet() const + { + return globalCellLevel_to_leafIdx_; + } + /// \brief Redistribute a global grid. /// /// The whole grid must be available on all processors. diff --git a/opm/grid/polyhedralgrid/cartesianindexmapper.hh b/opm/grid/polyhedralgrid/cartesianindexmapper.hh index e5c0d848d..284a104fc 100644 --- a/opm/grid/polyhedralgrid/cartesianindexmapper.hh +++ b/opm/grid/polyhedralgrid/cartesianindexmapper.hh @@ -31,7 +31,7 @@ namespace Dune const std::array& cartesianDimensions() const { - return grid_.logicalCartesianSize(); + return grid_.logicalCartesianSize(); } int cartesianSize() const @@ -44,12 +44,6 @@ namespace Dune return grid_.size( 0 ); } - // Only for unifying calls with CartesianIndexMapper where levels are relevant. - int compressedLevelZeroSize() const - { - return grid_.size( 0 ); - } - int cartesianIndex( const int compressedElementIndex ) const { assert( compressedElementIndex >= 0 && compressedElementIndex < compressedSize() ); @@ -73,6 +67,19 @@ namespace Dune coords[ 0 ] = gc ; } + // Only for unifying calls with CartesianIndexMapper where levels are relevant. + int compressedLevelZeroSize() const + { + return grid_.size( 0 ); + } + // Only for unifying calls with CartesianIndexMapper where levels are relevant. + int cartesianIndexLevel( const int compressedElementIndex, const int level ) const + { + if (level) { + throw std::invalid_argument("Invalid level.\n"); + } + return cartesianIndex(compressedElementIndex); + } // Only for unifying calls with CartesianIndexMapper where levels are relevant. void cartesianCoordinateLevel(const int compressedElementIndexOnLevel, std::array& coordsOnLevel, int level) const { diff --git a/tests/cpgrid/lgr_cartesian_idx_test.cpp b/tests/cpgrid/lgr_cartesian_idx_test.cpp new file mode 100644 index 000000000..0bcd62b84 --- /dev/null +++ b/tests/cpgrid/lgr_cartesian_idx_test.cpp @@ -0,0 +1,166 @@ +//=========================================================================== +// +// File: lgr_cartesian_idx_test.cpp +// +// Created: Sep 26 2024 09:15 +// +// Author(s): Antonella Ritorto +// +// $Date$ +// +// $Revision$ +// +//=========================================================================== +/* + Copyright 2024 Equinor ASA. + + This file is part of the Open Porous Media project (OPM). + + OPM is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OPM is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with OPM. If not, see . +*/ +#include "config.h" + +#define BOOST_TEST_MODULE LGRTests +#include +#include +#if BOOST_VERSION / 100000 == 1 && BOOST_VERSION / 100 % 1000 < 71 +#include +#else +#include +#endif +#include +#include +#include +#include +#include +#include +#include + +#include + +#include +#include +#include +#include +#include +#include + +struct Fixture +{ + Fixture() + { + int m_argc = boost::unit_test::framework::master_test_suite().argc; + char** m_argv = boost::unit_test::framework::master_test_suite().argv; + Dune::MPIHelper::instance(m_argc, m_argv); + Opm::OpmLog::setupSimpleDefaultLogging(); + } + + static int rank() + { + int m_argc = boost::unit_test::framework::master_test_suite().argc; + char** m_argv = boost::unit_test::framework::master_test_suite().argv; + return Dune::MPIHelper::instance(m_argc, m_argv).rank(); + } +}; + +BOOST_GLOBAL_FIXTURE(Fixture); + +void checkGlobalCellLgr(Dune::CpGrid& grid) +{ + const Dune::CartesianIndexMapper mapper{grid}; + + for (const auto& element : elements(grid.leafGridView())) + { + // How to get the Cartesian Index of a cell on the leaf grid view. + // 1. The "default" Cartesian Index of a cell on the leaf (which is a mixed grid with coarse and refined cells) is equal to the Cartesian Index of: + // (a) the equivalent cell in level zero (same center and volume, cell not involved in any refinement). + // (b) the parent cell in level zero when the leaf cell is a refined one (it has a equivalent cell that belongs to a refined level grid). + // This default Cartesian Index for leaf cells coincide with the globalCell values. + const auto& cartesian_idx_from_leaf_elem = mapper.cartesianIndex( element.index() ); + const auto& global_cell_idx_leaf = grid.globalCell()[element.index()]; // parent cell index when the cell is a refined one. + BOOST_CHECK_EQUAL(cartesian_idx_from_leaf_elem, global_cell_idx_leaf); + // global_ijk represents the ijk values of the equivalent cell on the level zero, or parent cell if the leaf cell is a refined one. + // Notice that all the refined cells of a same parent cell will get the same global_cell_ value and the same global_ijk (since they + // inherit the parent cell value). + std::array global_ijk = {0,0,0}; + mapper.cartesianCoordinate( element.index(), global_ijk); + + // How to get the Level Cartesian Index of a cell on the leaf grid view. + // Each LGR can be seen as a Cartesian Grid itself, with its own (local) logical_cartesian_size and its own (local) Cartesian indices. + // Given a leaf cell, via the CartesianIndexMapper and its method cartesianIndexLevel(...), we get the local-Cartesian-index. + const auto& cartesian_idx_from_level_elem = mapper.cartesianIndexLevel( element.getEquivLevelElem().index(), element.level() ); + const auto& global_cell_idx_level = grid.currentData()[element.level()]->globalCell()[element.getEquivLevelElem().index()]; + BOOST_CHECK_EQUAL(cartesian_idx_from_level_elem, global_cell_idx_level); + // local_ijk represents the ijk values of the equivalent cell on the level its was born. + std::array local_ijk = {0,0,0}; + mapper.cartesianCoordinateLevel( element.getEquivLevelElem().index(), local_ijk, element.level() ); + + // For leaf cells that were not involved in any refinement, global_ijk and local_ijk must coincide. + if(element.level()==0) + { + BOOST_CHECK_EQUAL( global_ijk[0], local_ijk[0]); + BOOST_CHECK_EQUAL( global_ijk[1], local_ijk[1]); + BOOST_CHECK_EQUAL( global_ijk[2], local_ijk[2]); + } + } + + for (int level = 0; level < grid.maxLevel(); ++level) + { + const auto& globalCellLevel_to_leafIdx = grid.currentData()[level]->getGlobalCellLevelToLeafIndexSet(); + for (const auto& element : elements(grid.levelGridView(level))) + { + if(element.isLeaf()) { + const auto& global_cell_level = grid.currentData()[element.level()]->globalCell()[element.index()]; + globalCellLevel_to_leafIdx.at(global_cell_level); + } + else { + BOOST_CHECK_THROW( globalCellLevel_to_leafIdx.at(element.index()), std::out_of_range); + } + } + } +} + +BOOST_AUTO_TEST_CASE(refine_one_cell) +{ + // Create a grid + Dune::CpGrid grid; + const std::array cell_sizes = {1.0, 1.0, 1.0}; + const std::array grid_dim = {4,3,3}; + grid.createCartesian(grid_dim, cell_sizes); + + const std::array cells_per_dim = {2,2,2}; + const std::array startIJK = {1,0,1}; + const std::array endIJK = {2,1,2};// Single Cell! (with index 13 in level zero grid), LGR dimensions 2x2x2 + const std::string lgr_name = {"LGR1"}; + grid.addLgrsUpdateLeafView({cells_per_dim}, {startIJK}, {endIJK}, {lgr_name}); + + checkGlobalCellLgr(grid); +} + +BOOST_AUTO_TEST_CASE(three_lgrs) +{ + // Create a grid + Dune::CpGrid grid; + const std::array cell_sizes = {1.0, 1.0, 1.0}; + const std::array grid_dim = {4,3,3}; + grid.createCartesian(grid_dim, cell_sizes); + + const std::vector> cells_per_dim_vec = {{2,2,2}, {3,3,3}, {4,4,4}}; + const std::vector> startIJK_vec = {{0,0,0}, {0,0,2}, {3,2,2}}; + const std::vector> endIJK_vec = {{2,1,1}, {1,1,3}, {4,3,3}}; + const std::vector lgr_name_vec = {"LGR1", "LGR2", "LGR3"}; + grid.addLgrsUpdateLeafView(cells_per_dim_vec, startIJK_vec, endIJK_vec, lgr_name_vec); + + checkGlobalCellLgr(grid); +} From 77f4d545ac540ba21402df97d34e7792d1b94420 Mon Sep 17 00:00:00 2001 From: Antonella Ritorto Date: Mon, 30 Sep 2024 15:04:35 +0200 Subject: [PATCH 2/7] Delete data member globalCellLevel_to_leafIdxSet_ from CpGridData --- opm/grid/CpGrid.hpp | 9 ++++----- opm/grid/cpgrid/CpGrid.cpp | 6 ++++-- opm/grid/cpgrid/CpGridData.hpp | 10 ---------- tests/cpgrid/lgr_cartesian_idx_test.cpp | 6 +++--- 4 files changed, 11 insertions(+), 20 deletions(-) diff --git a/opm/grid/CpGrid.hpp b/opm/grid/CpGrid.hpp index 2ee2c52d2..98e2fb823 100644 --- a/opm/grid/CpGrid.hpp +++ b/opm/grid/CpGrid.hpp @@ -968,11 +968,6 @@ namespace Dune /// For nested refinement, we lookup the oldest ancestor, from level zero. void computeGlobalCellLeafGridViewWithLgrs(std::vector& global_cell_leaf); - /// @brief Compute for each level grid, a map from the global_cell_[ cell index in level grid ] to the leaf index of the equivalent cell - /// on the leaf grid view. - /// Notice that cells that vanished and do not appear on the leaf grid view will not be considered. - void mapGlobalCellLevelToLeafIndexSet(); - /// @brief Get the ijk index of a refined corner, given its corner index of a single-cell-refinement. /// /// Given a single-cell, we refine it in {nx, ny, nz} refined children cells (per direction). Then, this single-cell-refinement @@ -1093,6 +1088,10 @@ namespace Dune public: + /// @brief Compute for each level grid, a map from the global_cell_[ cell index in level grid ] to the leaf index of the equivalent cell + /// on the leaf grid view. + /// Notice that cells that vanished and do not appear on the leaf grid view will not be considered. + std::vector> mapGlobalCellLevelToLeafIndexSet(); /// \brief Size of the overlap on the leaf level diff --git a/opm/grid/cpgrid/CpGrid.cpp b/opm/grid/cpgrid/CpGrid.cpp index 0ae9ed4fb..6a5e38000 100644 --- a/opm/grid/cpgrid/CpGrid.cpp +++ b/opm/grid/cpgrid/CpGrid.cpp @@ -679,12 +679,14 @@ void CpGrid::computeGlobalCellLeafGridViewWithLgrs(std::vector& global_cell } } -void CpGrid::mapGlobalCellLevelToLeafIndexSet() +std::vector> CpGrid::mapGlobalCellLevelToLeafIndexSet() { + std::vector> globalCellLevels_to_leafIdx(maxLevel()+1); // Plus level 0 for (const auto& element : elements(leafGridView())) { const auto& global_cell_level = currentData()[element.level()]->globalCell()[element.getEquivLevelElem().index()]; - (*currentData()[element.level()]).globalCellLevel_to_leafIdx_[global_cell_level] = element.index(); + globalCellLevels_to_leafIdx[element.level()][global_cell_level] = element.index(); } + return globalCellLevels_to_leafIdx; } void CpGrid::getIJK(const int c, std::array& ijk) const diff --git a/opm/grid/cpgrid/CpGridData.hpp b/opm/grid/cpgrid/CpGridData.hpp index 161de296b..424c2df59 100644 --- a/opm/grid/cpgrid/CpGridData.hpp +++ b/opm/grid/cpgrid/CpGridData.hpp @@ -632,16 +632,6 @@ class CpGridData return logical_cartesian_size_; } - /// Get map from global_cell_ to leaf index set - /// - /// Only relevant for CpGrid with LGRs. The maps of - /// all refined level grids is defined in CpGrid::adapt(/*..*/), - /// via CpGrid::mapGlobalCellLevelToLeafIndexSet(). - const std::unordered_map& getGlobalCellLevelToLeafIndexSet() const - { - return globalCellLevel_to_leafIdx_; - } - /// \brief Redistribute a global grid. /// /// The whole grid must be available on all processors. diff --git a/tests/cpgrid/lgr_cartesian_idx_test.cpp b/tests/cpgrid/lgr_cartesian_idx_test.cpp index 0bcd62b84..5a70a579d 100644 --- a/tests/cpgrid/lgr_cartesian_idx_test.cpp +++ b/tests/cpgrid/lgr_cartesian_idx_test.cpp @@ -117,15 +117,15 @@ void checkGlobalCellLgr(Dune::CpGrid& grid) for (int level = 0; level < grid.maxLevel(); ++level) { - const auto& globalCellLevel_to_leafIdx = grid.currentData()[level]->getGlobalCellLevelToLeafIndexSet(); + const auto& globalCellLevel_to_leafIdx = grid.mapGlobalCellLevelToLeafIndexSet(); for (const auto& element : elements(grid.levelGridView(level))) { if(element.isLeaf()) { const auto& global_cell_level = grid.currentData()[element.level()]->globalCell()[element.index()]; - globalCellLevel_to_leafIdx.at(global_cell_level); + globalCellLevel_to_leafIdx[element.level()].at(global_cell_level); } else { - BOOST_CHECK_THROW( globalCellLevel_to_leafIdx.at(element.index()), std::out_of_range); + BOOST_CHECK_THROW( globalCellLevel_to_leafIdx[element.level()].at(element.index()), std::out_of_range); } } } From 175cc872126f5a1c40fccf4983c1e1f67c60700c Mon Sep 17 00:00:00 2001 From: Antonella Ritorto Date: Mon, 30 Sep 2024 15:26:50 +0200 Subject: [PATCH 3/7] Map from LeafIdxSet to LevelGlobalCell --- opm/grid/CpGrid.hpp | 5 +++-- opm/grid/cpgrid/CpGrid.cpp | 13 ++++++++++++- 2 files changed, 15 insertions(+), 3 deletions(-) diff --git a/opm/grid/CpGrid.hpp b/opm/grid/CpGrid.hpp index 98e2fb823..2c2db7e66 100644 --- a/opm/grid/CpGrid.hpp +++ b/opm/grid/CpGrid.hpp @@ -1091,8 +1091,9 @@ namespace Dune /// @brief Compute for each level grid, a map from the global_cell_[ cell index in level grid ] to the leaf index of the equivalent cell /// on the leaf grid view. /// Notice that cells that vanished and do not appear on the leaf grid view will not be considered. - std::vector> mapGlobalCellLevelToLeafIndexSet(); - + std::vector> mapGlobalCellLevelToLeafIndexSet() const; + + std::vector> mapLeafIndexSetToGlobalCellLevel() const; /// \brief Size of the overlap on the leaf level unsigned int overlapSize(int) const; diff --git a/opm/grid/cpgrid/CpGrid.cpp b/opm/grid/cpgrid/CpGrid.cpp index 6a5e38000..bdbf44be4 100644 --- a/opm/grid/cpgrid/CpGrid.cpp +++ b/opm/grid/cpgrid/CpGrid.cpp @@ -679,7 +679,7 @@ void CpGrid::computeGlobalCellLeafGridViewWithLgrs(std::vector& global_cell } } -std::vector> CpGrid::mapGlobalCellLevelToLeafIndexSet() +std::vector> CpGrid::mapGlobalCellLevelToLeafIndexSet() const { std::vector> globalCellLevels_to_leafIdx(maxLevel()+1); // Plus level 0 for (const auto& element : elements(leafGridView())) { @@ -689,6 +689,17 @@ std::vector> CpGrid::mapGlobalCellL return globalCellLevels_to_leafIdx; } +std::vector> CpGrid::mapLeafIndexSetToGlobalCellLevel() const +{ + std::vector> leafIdx_to_globalCellLevel(currentData().back()->size(0)); + for (const auto& element : elements(leafGridView())) { + const auto& global_cell_level = currentData()[element.level()]->globalCell()[element.getEquivLevelElem().index()]; + leafIdx_to_globalCellLevel[element.index()] = {element.level(), global_cell_level}; + } + return leafIdx_to_globalCellLevel; +} + + void CpGrid::getIJK(const int c, std::array& ijk) const { current_view_data_->getIJK(c, ijk); From eda07fb13643b7eb72ae2e596b8a63df4e6acf94 Mon Sep 17 00:00:00 2001 From: Antonella Ritorto Date: Mon, 30 Sep 2024 15:45:33 +0200 Subject: [PATCH 4/7] change name maps globalCellLevels->localCartesianIndexSets --- opm/grid/CpGrid.hpp | 8 +++++--- opm/grid/cpgrid/CpGrid.cpp | 19 +++++++++---------- tests/cpgrid/lgr_cartesian_idx_test.cpp | 6 +++--- 3 files changed, 17 insertions(+), 16 deletions(-) diff --git a/opm/grid/CpGrid.hpp b/opm/grid/CpGrid.hpp index 2c2db7e66..3fda62539 100644 --- a/opm/grid/CpGrid.hpp +++ b/opm/grid/CpGrid.hpp @@ -1091,9 +1091,11 @@ namespace Dune /// @brief Compute for each level grid, a map from the global_cell_[ cell index in level grid ] to the leaf index of the equivalent cell /// on the leaf grid view. /// Notice that cells that vanished and do not appear on the leaf grid view will not be considered. - std::vector> mapGlobalCellLevelToLeafIndexSet() const; - - std::vector> mapLeafIndexSetToGlobalCellLevel() const; + /// global_cell_[ cell index in level grid ] coincide with (local) Cartesian Index. + std::vector> mapLocalCartesianIndexSetsToLeafIndexSet() const; + + /// @brief Reverse map: from leaf index cell to { level, local/level Cartesian index of the cell } + std::vector> mapLeafIndexSetToLocalCartesianIndexSets() const; /// \brief Size of the overlap on the leaf level unsigned int overlapSize(int) const; diff --git a/opm/grid/cpgrid/CpGrid.cpp b/opm/grid/cpgrid/CpGrid.cpp index bdbf44be4..0f8333c0c 100644 --- a/opm/grid/cpgrid/CpGrid.cpp +++ b/opm/grid/cpgrid/CpGrid.cpp @@ -679,27 +679,26 @@ void CpGrid::computeGlobalCellLeafGridViewWithLgrs(std::vector& global_cell } } -std::vector> CpGrid::mapGlobalCellLevelToLeafIndexSet() const +std::vector> CpGrid::mapLocalCartesianIndexSetsToLeafIndexSet() const { - std::vector> globalCellLevels_to_leafIdx(maxLevel()+1); // Plus level 0 + std::vector> localCartesianIdxSets_to_leafIdx(maxLevel()+1); // Plus level 0 for (const auto& element : elements(leafGridView())) { const auto& global_cell_level = currentData()[element.level()]->globalCell()[element.getEquivLevelElem().index()]; - globalCellLevels_to_leafIdx[element.level()][global_cell_level] = element.index(); + localCartesianIdxSets_to_leafIdx[element.level()][global_cell_level] = element.index(); } - return globalCellLevels_to_leafIdx; + return localCartesianIdxSets_to_leafIdx; } -std::vector> CpGrid::mapLeafIndexSetToGlobalCellLevel() const +std::vector> CpGrid::mapLeafIndexSetToLocalCartesianIndexSets() const { - std::vector> leafIdx_to_globalCellLevel(currentData().back()->size(0)); + std::vector> leafIdx_to_localCartesianIdxSets(currentData().back()->size(0)); for (const auto& element : elements(leafGridView())) { const auto& global_cell_level = currentData()[element.level()]->globalCell()[element.getEquivLevelElem().index()]; - leafIdx_to_globalCellLevel[element.index()] = {element.level(), global_cell_level}; + leafIdx_to_localCartesianIdxSets[element.index()] = {element.level(), global_cell_level}; } - return leafIdx_to_globalCellLevel; + return leafIdx_to_localCartesianIdxSets; } - void CpGrid::getIJK(const int c, std::array& ijk) const { current_view_data_->getIJK(c, ijk); @@ -2017,7 +2016,7 @@ bool CpGrid::adapt(const std::vector>& cells_per_dim_vec, computeGlobalCellLeafGridViewWithLgrs(global_cell_leaf); (*data[levels + preAdaptMaxLevel +1]).global_cell_.swap(global_cell_leaf); - mapGlobalCellLevelToLeafIndexSet(); + mapLocalCartesianIndexSetsToLeafIndexSet(); updateCornerHistoryLevels(cornerInMarkedElemWithEquivRefinedCorner, elemLgrAndElemLgrCorner_to_refinedLevelAndRefinedCorner, diff --git a/tests/cpgrid/lgr_cartesian_idx_test.cpp b/tests/cpgrid/lgr_cartesian_idx_test.cpp index 5a70a579d..993928ad4 100644 --- a/tests/cpgrid/lgr_cartesian_idx_test.cpp +++ b/tests/cpgrid/lgr_cartesian_idx_test.cpp @@ -117,15 +117,15 @@ void checkGlobalCellLgr(Dune::CpGrid& grid) for (int level = 0; level < grid.maxLevel(); ++level) { - const auto& globalCellLevel_to_leafIdx = grid.mapGlobalCellLevelToLeafIndexSet(); + const auto& localCartesianIdxSets_to_leafIdx = grid.mapLocalCartesianIndexSetsToLeafIndexSet(); for (const auto& element : elements(grid.levelGridView(level))) { if(element.isLeaf()) { const auto& global_cell_level = grid.currentData()[element.level()]->globalCell()[element.index()]; - globalCellLevel_to_leafIdx[element.level()].at(global_cell_level); + localCartesianIdxSets_to_leafIdx[element.level()].at(global_cell_level); } else { - BOOST_CHECK_THROW( globalCellLevel_to_leafIdx[element.level()].at(element.index()), std::out_of_range); + BOOST_CHECK_THROW( localCartesianIdxSets_to_leafIdx[element.level()].at(element.index()), std::out_of_range); } } } From 15a05f4b5a516979bc50c46c279d805408f5cb11 Mon Sep 17 00:00:00 2001 From: Antonella Ritorto Date: Mon, 30 Sep 2024 15:58:17 +0200 Subject: [PATCH 5/7] Call reverse map (leafIdxToLocalCartesianIdx) --- tests/cpgrid/lgr_cartesian_idx_test.cpp | 8 ++++++-- 1 file changed, 6 insertions(+), 2 deletions(-) diff --git a/tests/cpgrid/lgr_cartesian_idx_test.cpp b/tests/cpgrid/lgr_cartesian_idx_test.cpp index 993928ad4..bd1dcf666 100644 --- a/tests/cpgrid/lgr_cartesian_idx_test.cpp +++ b/tests/cpgrid/lgr_cartesian_idx_test.cpp @@ -114,15 +114,19 @@ void checkGlobalCellLgr(Dune::CpGrid& grid) BOOST_CHECK_EQUAL( global_ijk[2], local_ijk[2]); } } + + const auto& localCartesianIdxSets_to_leafIdx = grid.mapLocalCartesianIndexSetsToLeafIndexSet(); + const auto& leafIdx_to_localCartesianIdxSets = grid.mapLeafIndexSetToLocalCartesianIndexSets(); for (int level = 0; level < grid.maxLevel(); ++level) { - const auto& localCartesianIdxSets_to_leafIdx = grid.mapLocalCartesianIndexSetsToLeafIndexSet(); for (const auto& element : elements(grid.levelGridView(level))) { if(element.isLeaf()) { const auto& global_cell_level = grid.currentData()[element.level()]->globalCell()[element.index()]; - localCartesianIdxSets_to_leafIdx[element.level()].at(global_cell_level); + const auto& leaf_idx = localCartesianIdxSets_to_leafIdx[element.level()].at(global_cell_level); + BOOST_CHECK_EQUAL( leafIdx_to_localCartesianIdxSets[leaf_idx][0], element.level()); + BOOST_CHECK_EQUAL( leafIdx_to_localCartesianIdxSets[leaf_idx][1], global_cell_level); } else { BOOST_CHECK_THROW( localCartesianIdxSets_to_leafIdx[element.level()].at(element.index()), std::out_of_range); From 25d8f54c998af78c854910a51945555dd1be538e Mon Sep 17 00:00:00 2001 From: Antonella Ritorto Date: Tue, 1 Oct 2024 09:22:08 +0200 Subject: [PATCH 6/7] Add test with inactive parent cell --- tests/cpgrid/lgr_cartesian_idx_test.cpp | 63 ++++++++++++++++++++++++- 1 file changed, 61 insertions(+), 2 deletions(-) diff --git a/tests/cpgrid/lgr_cartesian_idx_test.cpp b/tests/cpgrid/lgr_cartesian_idx_test.cpp index bd1dcf666..e5f8edae4 100644 --- a/tests/cpgrid/lgr_cartesian_idx_test.cpp +++ b/tests/cpgrid/lgr_cartesian_idx_test.cpp @@ -47,7 +47,11 @@ #include #include +#include #include +#include +#include +#include #include #include @@ -122,14 +126,14 @@ void checkGlobalCellLgr(Dune::CpGrid& grid) { for (const auto& element : elements(grid.levelGridView(level))) { + const auto& global_cell_level = grid.currentData()[element.level()]->globalCell()[element.index()]; if(element.isLeaf()) { - const auto& global_cell_level = grid.currentData()[element.level()]->globalCell()[element.index()]; const auto& leaf_idx = localCartesianIdxSets_to_leafIdx[element.level()].at(global_cell_level); BOOST_CHECK_EQUAL( leafIdx_to_localCartesianIdxSets[leaf_idx][0], element.level()); BOOST_CHECK_EQUAL( leafIdx_to_localCartesianIdxSets[leaf_idx][1], global_cell_level); } else { - BOOST_CHECK_THROW( localCartesianIdxSets_to_leafIdx[element.level()].at(element.index()), std::out_of_range); + BOOST_CHECK_THROW( localCartesianIdxSets_to_leafIdx[element.level()].at(global_cell_level), std::out_of_range); } } } @@ -168,3 +172,58 @@ BOOST_AUTO_TEST_CASE(three_lgrs) checkGlobalCellLgr(grid); } + +BOOST_AUTO_TEST_CASE(inactiveCells_in_lgrs) +{ + + const std::string deckString = + R"( RUNSPEC + DIMENS + 1 1 5 / + GRID + COORD + 0 0 0 + 0 0 1 + 1 0 0 + 1 0 1 + 0 1 0 + 0 1 1 + 1 1 0 + 1 1 1 + / + ZCORN + 4*0 + 8*1 + 8*2 + 8*3 + 8*4 + 4*5 + / + ACTNUM + 0 + 1 + 1 + 1 + 0 + / + PORO + 5*0.15 + /)"; + + const std::vector> cells_per_dim_vec = {{2,2,2}, {3,3,3}}; + const std::vector> startIJK_vec = {{0,0,0}, {0,0,3}}; + const std::vector> endIJK_vec = {{1,1,2}, {1,1,5}}; + // LGR1 cell indices = {0,1}, LGR2 cell indices = {3,4}. + const std::vector lgr_name_vec = {"LGR1", "LGR2"}; + + Opm::Parser parser; + const auto deck = parser.parseString(deckString); + Opm::EclipseState es(deck); + + Dune::CpGrid grid; + grid.processEclipseFormat(&es.getInputGrid(), &es, false, false, false); + + grid.addLgrsUpdateLeafView(cells_per_dim_vec, startIJK_vec, endIJK_vec, lgr_name_vec); + + checkGlobalCellLgr(grid); +} From 2a685178ed295e73f590f21d968a8dea39fa3cbe Mon Sep 17 00:00:00 2001 From: Antonella Ritorto Date: Tue, 1 Oct 2024 11:30:25 +0200 Subject: [PATCH 7/7] Typo fixed --- opm/grid/cpgrid/CartesianIndexMapper.hpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/opm/grid/cpgrid/CartesianIndexMapper.hpp b/opm/grid/cpgrid/CartesianIndexMapper.hpp index fbc3ed511..66487029e 100644 --- a/opm/grid/cpgrid/CartesianIndexMapper.hpp +++ b/opm/grid/cpgrid/CartesianIndexMapper.hpp @@ -61,7 +61,7 @@ namespace Dune grid_.getIJK( compressedElementIndex, coords ); } - /// Additional methods realted to LGRs + /// Additional methods related to LGRs int compressedLevelZeroSize() const { return (*grid_.currentData()[0]).size(0); @@ -83,7 +83,7 @@ namespace Dune } grid_.currentData()[level]->getIJK( compressedElementIndexOnLevel, coordsOnLevel); } - /// Additional methods realted to LGRs. END + /// Additional methods related to LGRs. END }; } // end namespace Opm