From 67a0044cea9a84c9d27d31d4d7d60c49f919b944 Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Mon, 3 Nov 2025 17:50:11 -0800 Subject: [PATCH 01/35] Material view iterator for unibuffer matset. --- src/axom/bump/tests/bump_views.cpp | 40 +++++++++- src/axom/bump/views/MaterialView.hpp | 112 +++++++++++++++++++++++++++ 2 files changed, 151 insertions(+), 1 deletion(-) diff --git a/src/axom/bump/tests/bump_views.cpp b/src/axom/bump/tests/bump_views.cpp index f1a2c1d963..14a449f465 100644 --- a/src/axom/bump/tests/bump_views.cpp +++ b/src/axom/bump/tests/bump_views.cpp @@ -571,6 +571,8 @@ struct test_braid2d_mat // _bump_views_matsetview_end // clang-format on test_matsetview(nzones, matsetView, allocatorID); + // Test iterators. + test_matsetview_iterators(matsetView, allocatorID); } else if(mattype == "multibuffer") { @@ -658,6 +660,41 @@ struct test_braid2d_mat EXPECT_EQ(results[i], resultsHost[i]); } } + + template + static void test_matsetview_iterators(MatsetView matsetView, int allocatorID) + { + // Allocate results array on device. + const int nResults = matsetView.numberOfZones(); + axom::Array resultsArrayDevice(nResults, nResults, allocatorID); + auto resultsView = resultsArrayDevice.view(); + + axom::for_all( + matsetView.numberOfZones(), + AXOM_LAMBDA(axom::IndexType index) { + typename MatsetView::IDList ids {}; + typename MatsetView::VFList vfs {}; + matsetView.zoneMaterials(index, ids, vfs); + + const auto end = matsetView.endZone(index); + int i = 0; + int eq_count = 0; + for(auto it = matsetView.beginZone(index); it != end; it++, i++) + { + eq_count += (vfs[i] == it.volume_fraction() && + ids[i] == it.material_id()) ? 1 : 0; + } + resultsView[index] = (eq_count == ids.size()) ? 1 : 0; + }); + + // Get containsView data to the host and compare results + std::vector resultsHost(nResults); + axom::copy(resultsHost.data(), resultsView.data(), sizeof(int) * nResults); + for(int i = 0; i < nResults; i++) + { + EXPECT_EQ(resultsHost[i], 1); + } + } }; // Unibuffer @@ -665,6 +702,7 @@ TEST(bump_views, matset_unibuffer_seq) { test_braid2d_mat::test("uniform", "unibuffer", "uniform2d_unibuffer"); } +#if 0 #if defined(AXOM_USE_OPENMP) TEST(bump_views, matset_unibuffer_omp) { @@ -826,7 +864,7 @@ TEST(bump_views, matset_multibuffer) EXPECT_EQ(vf, 1.); }); } - +#endif //------------------------------------------------------------------------------ int main(int argc, char *argv[]) { diff --git a/src/axom/bump/views/MaterialView.hpp b/src/axom/bump/views/MaterialView.hpp index d97549ff6f..f14daa1493 100644 --- a/src/axom/bump/views/MaterialView.hpp +++ b/src/axom/bump/views/MaterialView.hpp @@ -163,6 +163,118 @@ class UnibufferMaterialView return false; } + /*! + * \brief An iterator class for iterating over the mat/vf data in a zone. + */ + class iterator + { + // Let the material view call the iterator constructor. + friend class UnibufferMaterialView; + public: + /// Get the current material id for the iterator. + MaterialID AXOM_HOST_DEVICE material_id() const + { + return m_material_ids[m_indices[m_currentIndex]]; + } + /// Get the current volume fraction for the iterator. + FloatType AXOM_HOST_DEVICE volume_fraction() const + { + return m_volume_fractions[m_indices[m_currentIndex]]; + } + axom::IndexType AXOM_HOST_DEVICE size() const + { + return m_indices.size(); + } + void AXOM_HOST_DEVICE operator ++() + { + SLIC_ASSERT(m_currentIndex < m_indices.size()); + m_currentIndex++; + } + void AXOM_HOST_DEVICE operator ++(int) + { + SLIC_ASSERT(m_currentIndex < m_indices.size()); + m_currentIndex++; + } + bool AXOM_HOST_DEVICE operator == (const iterator &rhs) const + { + return m_currentIndex == rhs.m_currentIndex && + m_material_ids == rhs.m_material_ids && + m_volume_fractions == rhs.m_volume_fractions && + m_indices == rhs.m_indices; + } + bool AXOM_HOST_DEVICE operator < (const iterator &rhs) const + { + return m_currentIndex < rhs.m_currentIndex && + (m_material_ids == rhs.m_material_ids && + m_volume_fractions == rhs.m_volume_fractions && + m_indices == rhs.m_indices); + } + bool AXOM_HOST_DEVICE operator > (const iterator &rhs) const + { + return m_currentIndex > rhs.m_currentIndex && + (m_material_ids == rhs.m_material_ids && + m_volume_fractions == rhs.m_volume_fractions && + m_indices == rhs.m_indices); + } + bool AXOM_HOST_DEVICE operator != (const iterator &rhs) const + { + return m_currentIndex != rhs.m_currentIndex || + m_material_ids != rhs.m_material_ids || + m_volume_fractions != rhs.m_volume_fractions || + m_indices != rhs.m_indices; + } + private: + AXOM_HOST_DEVICE iterator(axom::ArrayView material_ids, + axom::ArrayView volume_fractions, + axom::ArrayView indices, + axom::IndexType currentIndex = 0) : + m_material_ids(material_ids), + m_volume_fractions(volume_fractions), + m_indices(indices), + m_currentIndex(currentIndex) + { + } + + axom::ArrayView m_material_ids; + axom::ArrayView m_volume_fractions; + axom::ArrayView m_indices; + axom::IndexType m_currentIndex; + }; + + /*! + * \brief Return the iterator for the beginning of a zone's material data. + * + * \param zi The zone index being queried. + * + * \return The iterator for the beginning of a zone's material data. + */ + iterator AXOM_HOST_DEVICE beginZone(ZoneIndex zi) const + { + SLIC_ASSERT(zi < static_cast(numberOfZones())); + + return iterator(m_material_ids, + m_volume_fractions, + axom::ArrayView(m_indices.data() + m_offsets[zi], m_sizes[zi]), + 0); + } + + /*! + * \brief Return the iterator for the end of a zone's material data. + * + * \param zi The zone index being queried. + * + * \return The iterator for the end of a zone's material data. + */ + iterator AXOM_HOST_DEVICE endZone(ZoneIndex zi) const + { + SLIC_ASSERT(zi < static_cast(numberOfZones())); + + return iterator(m_material_ids, + m_volume_fractions, + axom::ArrayView(m_indices.data() + m_offsets[zi], m_sizes[zi]), + m_sizes[zi]); + } + private: axom::ArrayView m_material_ids; axom::ArrayView m_volume_fractions; From 951df9b23e68fb177727e53dfb9418d594856774 Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Tue, 4 Nov 2025 11:30:08 -0800 Subject: [PATCH 02/35] Added MultiBufferMaterialView::iterator. --- src/axom/bump/tests/bump_views.cpp | 9 +- src/axom/bump/views/MaterialView.hpp | 128 +++++++++++++++++++++++++++ 2 files changed, 134 insertions(+), 3 deletions(-) diff --git a/src/axom/bump/tests/bump_views.cpp b/src/axom/bump/tests/bump_views.cpp index 14a449f465..96e4ee4495 100644 --- a/src/axom/bump/tests/bump_views.cpp +++ b/src/axom/bump/tests/bump_views.cpp @@ -578,7 +578,11 @@ struct test_braid2d_mat { axom::bump::views::dispatch_material_multibuffer( deviceMesh["matsets/mat"], - [&](auto matsetView) { test_matsetview(nzones, matsetView, allocatorID); }); + [&](auto matsetView) { + test_matsetview(nzones, matsetView, allocatorID); + // Test iterators. + test_matsetview_iterators(matsetView, allocatorID); + }); } else if(mattype == "element_dominant") { @@ -702,7 +706,6 @@ TEST(bump_views, matset_unibuffer_seq) { test_braid2d_mat::test("uniform", "unibuffer", "uniform2d_unibuffer"); } -#if 0 #if defined(AXOM_USE_OPENMP) TEST(bump_views, matset_unibuffer_omp) { @@ -745,7 +748,7 @@ TEST(bump_views, matset_multibuffer_hip) test_braid2d_mat::test("uniform", "multibuffer", "uniform2d_multibuffer"); } #endif - +#if 0 // Element-dominant TEST(bump_views, matset_element_dominant_seq) { diff --git a/src/axom/bump/views/MaterialView.hpp b/src/axom/bump/views/MaterialView.hpp index f14daa1493..a08cf8c559 100644 --- a/src/axom/bump/views/MaterialView.hpp +++ b/src/axom/bump/views/MaterialView.hpp @@ -416,6 +416,134 @@ class MultiBufferMaterialView return found; } + /*! + * \brief An iterator class for iterating over the mat/vf data in a zone. + */ + class iterator + { + // Let the material view call the iterator constructor. + friend class MultiBufferMaterialView; + public: + /// Get the current material id for the iterator. + MaterialID AXOM_HOST_DEVICE material_id() const + { + SLIC_ASSERT(m_currentIndex < m_view->m_size); + return m_view->m_matnos[m_currentIndex]; + } + /// Get the current volume fraction for the iterator. + FloatType AXOM_HOST_DEVICE volume_fraction() const + { + SLIC_ASSERT(m_currentIndex < m_view->m_size); + const auto &curIndices = m_view->m_indices[m_currentIndex]; + const auto &curValues = m_view->m_values[m_currentIndex]; + const auto idx = curIndices[m_zoneIndex]; + return curValues[idx]; + } + axom::IndexType AXOM_HOST_DEVICE size() const + { + return m_view->numberOfMaterials(m_zoneIndex); + } + void AXOM_HOST_DEVICE operator ++() + { + m_currentIndex += (m_currentIndex < m_view->m_size) ? 1 : 0; + advance(); + } + void AXOM_HOST_DEVICE operator ++(int) + { + m_currentIndex += (m_currentIndex < m_view->m_size) ? 1 : 0; + advance(); + } + bool AXOM_HOST_DEVICE operator == (const iterator &rhs) const + { + return m_currentIndex == rhs.m_currentIndex && + m_zoneIndex == rhs.m_zoneIndex && + m_view == rhs.m_view; + } + bool AXOM_HOST_DEVICE operator < (const iterator &rhs) const + { + return m_currentIndex < rhs.m_currentIndex && + m_zoneIndex < rhs.m_zoneIndex && + m_view == rhs.m_view; + } + bool AXOM_HOST_DEVICE operator > (const iterator &rhs) const + { + return m_currentIndex > rhs.m_currentIndex && + m_zoneIndex > rhs.m_zoneIndex && + m_view == rhs.m_view; + } + bool AXOM_HOST_DEVICE operator != (const iterator &rhs) const + { + return m_currentIndex != rhs.m_currentIndex || + m_zoneIndex != rhs.m_zoneIndex || + m_view != rhs.m_view; + } + private: + /// Constructor + AXOM_HOST_DEVICE iterator(const MultiBufferMaterialView *view, + ZoneIndex zoneIndex, + axom::IndexType currentIndex = 0) : + m_view(view), + m_zoneIndex(zoneIndex), + m_currentIndex(currentIndex) + { + } + + /// Advance to the next valid material slot for the zone. + void AXOM_HOST_DEVICE advance() + { + while(m_currentIndex < m_view->m_size) + { + const auto &curIndices = m_view->m_indices[m_currentIndex]; + const auto &curValues = m_view->m_values[m_currentIndex]; + + if(m_zoneIndex < static_cast(curIndices.size())) + { + const auto idx = curIndices[m_zoneIndex]; + if(curValues[idx] > 0) + { + break; + } + } + m_currentIndex++; + } + } + + const MultiBufferMaterialView *m_view; + ZoneIndex m_zoneIndex; + axom::IndexType m_currentIndex; + }; + // Let the iterator access members. + friend class iterator; + + /*! + * \brief Return the iterator for the beginning of a zone's material data. + * + * \param zi The zone index being queried. + * + * \return The iterator for the beginning of a zone's material data. + */ + iterator AXOM_HOST_DEVICE beginZone(ZoneIndex zi) const + { + SLIC_ASSERT(zi < static_cast(numberOfZones())); + + auto it = iterator(this, zi, 0); + it.advance(); + return it; + } + + /*! + * \brief Return the iterator for the end of a zone's material data. + * + * \param zi The zone index being queried. + * + * \return The iterator for the end of a zone's material data. + */ + iterator AXOM_HOST_DEVICE endZone(ZoneIndex zi) const + { + SLIC_ASSERT(zi < static_cast(numberOfZones())); + return iterator(this, zi, m_size); + } + private: AXOM_HOST_DEVICE axom::IndexType indexOfMaterialID(MaterialID mat) const From 77d5a003d3a9394e4e9ea1f44cccb4b19ef9403b Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Tue, 4 Nov 2025 13:34:23 -0800 Subject: [PATCH 03/35] Added ElementDominantMaterialView::iterator --- src/axom/bump/tests/bump_views.cpp | 10 ++- src/axom/bump/views/MaterialView.hpp | 118 +++++++++++++++++++++++++++ 2 files changed, 125 insertions(+), 3 deletions(-) diff --git a/src/axom/bump/tests/bump_views.cpp b/src/axom/bump/tests/bump_views.cpp index 96e4ee4495..faa81f71ac 100644 --- a/src/axom/bump/tests/bump_views.cpp +++ b/src/axom/bump/tests/bump_views.cpp @@ -588,7 +588,11 @@ struct test_braid2d_mat { axom::bump::views::dispatch_material_element_dominant( deviceMesh["matsets/mat"], - [&](auto matsetView) { test_matsetview(nzones, matsetView, allocatorID); }); + [&](auto matsetView) { + test_matsetview(nzones, matsetView, allocatorID); + // Test iterators. + test_matsetview_iterators(matsetView, allocatorID); + }); } else if(mattype == "material_dominant") { @@ -748,7 +752,7 @@ TEST(bump_views, matset_multibuffer_hip) test_braid2d_mat::test("uniform", "multibuffer", "uniform2d_multibuffer"); } #endif -#if 0 + // Element-dominant TEST(bump_views, matset_element_dominant_seq) { @@ -772,7 +776,7 @@ TEST(bump_views, matset_element_dominant_hip) test_braid2d_mat::test("uniform", "element_dominant", "uniform2d_element_dominant"); } #endif - +#if 0 // Material-dominant TEST(bump_views, matset_material_dominant_seq) { diff --git a/src/axom/bump/views/MaterialView.hpp b/src/axom/bump/views/MaterialView.hpp index a08cf8c559..4d5add0f98 100644 --- a/src/axom/bump/views/MaterialView.hpp +++ b/src/axom/bump/views/MaterialView.hpp @@ -685,6 +685,124 @@ class ElementDominantMaterialView return found; } + /*! + * \brief An iterator class for iterating over the mat/vf data in a zone. + */ + class iterator + { + // Let the material view call the iterator constructor. + friend class ElementDominantMaterialView; + public: + /// Get the current material id for the iterator. + MaterialID AXOM_HOST_DEVICE material_id() const + { + SLIC_ASSERT(m_currentIndex < m_view->m_volume_fractions.size()); + return m_view->m_matnos[m_currentIndex]; + } + /// Get the current volume fraction for the iterator. + FloatType AXOM_HOST_DEVICE volume_fraction() const + { + SLIC_ASSERT(m_currentIndex < m_view->m_volume_fractions.size()); + return m_view->m_volume_fractions[m_currentIndex][m_zoneIndex]; + } + axom::IndexType AXOM_HOST_DEVICE size() const + { + return m_view->numberOfMaterials(m_zoneIndex); + } + void AXOM_HOST_DEVICE operator ++() + { + m_currentIndex += (m_currentIndex < m_view->m_volume_fractions.size()) ? 1 : 0; + advance(); + } + void AXOM_HOST_DEVICE operator ++(int) + { + m_currentIndex += (m_currentIndex < m_view->m_volume_fractions.size()) ? 1 : 0; + advance(); + } + bool AXOM_HOST_DEVICE operator == (const iterator &rhs) const + { + return m_currentIndex == rhs.m_currentIndex && + m_zoneIndex == rhs.m_zoneIndex && + m_view == rhs.m_view; + } + bool AXOM_HOST_DEVICE operator < (const iterator &rhs) const + { + return m_currentIndex < rhs.m_currentIndex && + m_zoneIndex < rhs.m_zoneIndex && + m_view == rhs.m_view; + } + bool AXOM_HOST_DEVICE operator > (const iterator &rhs) const + { + return m_currentIndex > rhs.m_currentIndex && + m_zoneIndex > rhs.m_zoneIndex && + m_view == rhs.m_view; + } + bool AXOM_HOST_DEVICE operator != (const iterator &rhs) const + { + return m_currentIndex != rhs.m_currentIndex || + m_zoneIndex != rhs.m_zoneIndex || + m_view != rhs.m_view; + } + private: + /// Constructor + AXOM_HOST_DEVICE iterator(const ElementDominantMaterialView *view, + ZoneIndex zoneIndex, + axom::IndexType currentIndex = 0) : + m_view(view), + m_zoneIndex(zoneIndex), + m_currentIndex(currentIndex) + { + } + + /// Advance to the next valid material slot for the zone. + void AXOM_HOST_DEVICE advance() + { + while(m_currentIndex < m_view->m_volume_fractions.size()) + { + if(m_view->m_volume_fractions[m_currentIndex][m_zoneIndex] > 0) + { + break; + } + m_currentIndex++; + } + } + + const ElementDominantMaterialView *m_view; + ZoneIndex m_zoneIndex; + axom::IndexType m_currentIndex; + }; + // Let the iterator access members. + friend class iterator; + + /*! + * \brief Return the iterator for the beginning of a zone's material data. + * + * \param zi The zone index being queried. + * + * \return The iterator for the beginning of a zone's material data. + */ + iterator AXOM_HOST_DEVICE beginZone(ZoneIndex zi) const + { + SLIC_ASSERT(zi < static_cast(numberOfZones())); + + auto it = iterator(this, zi, 0); + it.advance(); + return it; + } + + /*! + * \brief Return the iterator for the end of a zone's material data. + * + * \param zi The zone index being queried. + * + * \return The iterator for the end of a zone's material data. + */ + iterator AXOM_HOST_DEVICE endZone(ZoneIndex zi) const + { + SLIC_ASSERT(zi < static_cast(numberOfZones())); + return iterator(this, zi, m_volume_fractions.size()); + } + private: #if !defined(AXOM_DEVICE_CODE) void checkSizes() const From 53e40cfb14764f5867a699f3da95cdc5014fabfd Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Tue, 4 Nov 2025 13:42:37 -0800 Subject: [PATCH 04/35] Delete default constructors for matset iterator types. --- src/axom/bump/views/MaterialView.hpp | 7 +++++++ 1 file changed, 7 insertions(+) diff --git a/src/axom/bump/views/MaterialView.hpp b/src/axom/bump/views/MaterialView.hpp index 4d5add0f98..7ffcc428e7 100644 --- a/src/axom/bump/views/MaterialView.hpp +++ b/src/axom/bump/views/MaterialView.hpp @@ -224,6 +224,9 @@ class UnibufferMaterialView m_indices != rhs.m_indices; } private: + DISABLE_DEFAULT_CTOR(iterator); + + /// Constructor AXOM_HOST_DEVICE iterator(axom::ArrayView material_ids, axom::ArrayView volume_fractions, axom::ArrayView indices, @@ -478,6 +481,8 @@ class MultiBufferMaterialView m_view != rhs.m_view; } private: + DISABLE_DEFAULT_CTOR(iterator); + /// Constructor AXOM_HOST_DEVICE iterator(const MultiBufferMaterialView *view, ZoneIndex zoneIndex, @@ -744,6 +749,8 @@ class ElementDominantMaterialView m_view != rhs.m_view; } private: + DISABLE_DEFAULT_CTOR(iterator); + /// Constructor AXOM_HOST_DEVICE iterator(const ElementDominantMaterialView *view, ZoneIndex zoneIndex, From b302ed4a1547817c581e18d39b4486709e3ff572 Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Tue, 4 Nov 2025 17:43:37 -0800 Subject: [PATCH 05/35] Material view iterator improvements. --- src/axom/bump/tests/bump_views.cpp | 34 ++-- src/axom/bump/utilities/conduit_memory.hpp | 2 +- src/axom/bump/views/MaterialView.hpp | 196 ++++++++++++++++----- 3 files changed, 176 insertions(+), 56 deletions(-) diff --git a/src/axom/bump/tests/bump_views.cpp b/src/axom/bump/tests/bump_views.cpp index faa81f71ac..ee11e8fb50 100644 --- a/src/axom/bump/tests/bump_views.cpp +++ b/src/axom/bump/tests/bump_views.cpp @@ -571,8 +571,6 @@ struct test_braid2d_mat // _bump_views_matsetview_end // clang-format on test_matsetview(nzones, matsetView, allocatorID); - // Test iterators. - test_matsetview_iterators(matsetView, allocatorID); } else if(mattype == "multibuffer") { @@ -580,8 +578,6 @@ struct test_braid2d_mat deviceMesh["matsets/mat"], [&](auto matsetView) { test_matsetview(nzones, matsetView, allocatorID); - // Test iterators. - test_matsetview_iterators(matsetView, allocatorID); }); } else if(mattype == "element_dominant") @@ -590,8 +586,6 @@ struct test_braid2d_mat deviceMesh["matsets/mat"], [&](auto matsetView) { test_matsetview(nzones, matsetView, allocatorID); - // Test iterators. - test_matsetview_iterators(matsetView, allocatorID); }); } else if(mattype == "material_dominant") @@ -667,11 +661,15 @@ struct test_braid2d_mat { EXPECT_EQ(results[i], resultsHost[i]); } + + // Test iterators. + test_matsetview_iterators(matsetView, allocatorID); } template static void test_matsetview_iterators(MatsetView matsetView, int allocatorID) { + using ZoneIndex = typename MatsetView::ZoneIndex; // Allocate results array on device. const int nResults = matsetView.numberOfZones(); axom::Array resultsArrayDevice(nResults, nResults, allocatorID); @@ -684,15 +682,31 @@ struct test_braid2d_mat typename MatsetView::VFList vfs {}; matsetView.zoneMaterials(index, ids, vfs); + // Get the end iterator for the zone. const auto end = matsetView.endZone(index); - int i = 0; + int eq_count = 0; + int count = 0; + + // Make sure the iterator is for the right zone. + eq_count += (end.zoneIndex() == static_cast(index)) ? 1 : 0; + count++; + + // Make sure incrementing the last iterator has no effect. + auto end2 = end; + end2++; + eq_count += (end == end2) ? 1 : 0; + count++; + + // Make sure the iterator order is the same as for the values we got from zoneMaterials(). + int i = 0; for(auto it = matsetView.beginZone(index); it != end; it++, i++) { eq_count += (vfs[i] == it.volume_fraction() && ids[i] == it.material_id()) ? 1 : 0; + count++; } - resultsView[index] = (eq_count == ids.size()) ? 1 : 0; + resultsView[index] = (eq_count == count) ? 1 : 0; }); // Get containsView data to the host and compare results @@ -776,7 +790,7 @@ TEST(bump_views, matset_element_dominant_hip) test_braid2d_mat::test("uniform", "element_dominant", "uniform2d_element_dominant"); } #endif -#if 0 + // Material-dominant TEST(bump_views, matset_material_dominant_seq) { @@ -871,7 +885,7 @@ TEST(bump_views, matset_multibuffer) EXPECT_EQ(vf, 1.); }); } -#endif + //------------------------------------------------------------------------------ int main(int argc, char *argv[]) { diff --git a/src/axom/bump/utilities/conduit_memory.hpp b/src/axom/bump/utilities/conduit_memory.hpp index 838e00fd95..fd0810fbd3 100644 --- a/src/axom/bump/utilities/conduit_memory.hpp +++ b/src/axom/bump/utilities/conduit_memory.hpp @@ -142,7 +142,7 @@ void copy(conduit::Node &dest, const conduit::Node &src) else { const int allocatorID = axom::getAllocatorIDFromPointer(src.data_ptr()); - bool deviceAllocated = isDeviceAllocator(allocatorID); + bool deviceAllocated = (allocatorID == INVALID_ALLOCATOR_ID) ? false : isDeviceAllocator(allocatorID); if(deviceAllocated || (!src.dtype().is_string() && src.dtype().number_of_elements() > 1)) { // Allocate the node's memory in the right place. diff --git a/src/axom/bump/views/MaterialView.hpp b/src/axom/bump/views/MaterialView.hpp index 7ffcc428e7..44eb186d49 100644 --- a/src/axom/bump/views/MaterialView.hpp +++ b/src/axom/bump/views/MaterialView.hpp @@ -174,51 +174,41 @@ class UnibufferMaterialView /// Get the current material id for the iterator. MaterialID AXOM_HOST_DEVICE material_id() const { + SLIC_ASSERT(m_currentIndex < size()); return m_material_ids[m_indices[m_currentIndex]]; } /// Get the current volume fraction for the iterator. FloatType AXOM_HOST_DEVICE volume_fraction() const { + SLIC_ASSERT(m_currentIndex < size()); return m_volume_fractions[m_indices[m_currentIndex]]; } axom::IndexType AXOM_HOST_DEVICE size() const { return m_indices.size(); } + ZoneIndex AXOM_HOST_DEVICE zoneIndex() const { return m_zoneIndex; } + void AXOM_HOST_DEVICE operator ++() { - SLIC_ASSERT(m_currentIndex < m_indices.size()); - m_currentIndex++; + m_currentIndex += (m_currentIndex < size()) ? 1 : 0; } void AXOM_HOST_DEVICE operator ++(int) { - SLIC_ASSERT(m_currentIndex < m_indices.size()); - m_currentIndex++; + m_currentIndex += (m_currentIndex < size()) ? 1 : 0; } bool AXOM_HOST_DEVICE operator == (const iterator &rhs) const { return m_currentIndex == rhs.m_currentIndex && + m_zoneIndex == rhs.m_zoneIndex && m_material_ids == rhs.m_material_ids && m_volume_fractions == rhs.m_volume_fractions && m_indices == rhs.m_indices; } - bool AXOM_HOST_DEVICE operator < (const iterator &rhs) const - { - return m_currentIndex < rhs.m_currentIndex && - (m_material_ids == rhs.m_material_ids && - m_volume_fractions == rhs.m_volume_fractions && - m_indices == rhs.m_indices); - } - bool AXOM_HOST_DEVICE operator > (const iterator &rhs) const - { - return m_currentIndex > rhs.m_currentIndex && - (m_material_ids == rhs.m_material_ids && - m_volume_fractions == rhs.m_volume_fractions && - m_indices == rhs.m_indices); - } bool AXOM_HOST_DEVICE operator != (const iterator &rhs) const { return m_currentIndex != rhs.m_currentIndex || + m_zoneIndex != rhs.m_zoneIndex || m_material_ids != rhs.m_material_ids || m_volume_fractions != rhs.m_volume_fractions || m_indices != rhs.m_indices; @@ -227,13 +217,15 @@ class UnibufferMaterialView DISABLE_DEFAULT_CTOR(iterator); /// Constructor - AXOM_HOST_DEVICE iterator(axom::ArrayView material_ids, - axom::ArrayView volume_fractions, - axom::ArrayView indices, + AXOM_HOST_DEVICE iterator(const axom::ArrayView &material_ids, + const axom::ArrayView &volume_fractions, + const axom::ArrayView &indices, + ZoneIndex zoneIndex, axom::IndexType currentIndex = 0) : m_material_ids(material_ids), m_volume_fractions(volume_fractions), m_indices(indices), + m_zoneIndex(zoneIndex), m_currentIndex(currentIndex) { } @@ -241,6 +233,7 @@ class UnibufferMaterialView axom::ArrayView m_material_ids; axom::ArrayView m_volume_fractions; axom::ArrayView m_indices; + ZoneIndex m_zoneIndex; axom::IndexType m_currentIndex; }; @@ -258,6 +251,7 @@ class UnibufferMaterialView return iterator(m_material_ids, m_volume_fractions, axom::ArrayView(m_indices.data() + m_offsets[zi], m_sizes[zi]), + zi, 0); } @@ -275,6 +269,7 @@ class UnibufferMaterialView return iterator(m_material_ids, m_volume_fractions, axom::ArrayView(m_indices.data() + m_offsets[zi], m_sizes[zi]), + zi, m_sizes[zi]); } @@ -446,6 +441,7 @@ class MultiBufferMaterialView { return m_view->numberOfMaterials(m_zoneIndex); } + ZoneIndex AXOM_HOST_DEVICE zoneIndex() const { return m_zoneIndex; } void AXOM_HOST_DEVICE operator ++() { m_currentIndex += (m_currentIndex < m_view->m_size) ? 1 : 0; @@ -462,18 +458,6 @@ class MultiBufferMaterialView m_zoneIndex == rhs.m_zoneIndex && m_view == rhs.m_view; } - bool AXOM_HOST_DEVICE operator < (const iterator &rhs) const - { - return m_currentIndex < rhs.m_currentIndex && - m_zoneIndex < rhs.m_zoneIndex && - m_view == rhs.m_view; - } - bool AXOM_HOST_DEVICE operator > (const iterator &rhs) const - { - return m_currentIndex > rhs.m_currentIndex && - m_zoneIndex > rhs.m_zoneIndex && - m_view == rhs.m_view; - } bool AXOM_HOST_DEVICE operator != (const iterator &rhs) const { return m_currentIndex != rhs.m_currentIndex || @@ -714,6 +698,7 @@ class ElementDominantMaterialView { return m_view->numberOfMaterials(m_zoneIndex); } + ZoneIndex AXOM_HOST_DEVICE zoneIndex() const { return m_zoneIndex; } void AXOM_HOST_DEVICE operator ++() { m_currentIndex += (m_currentIndex < m_view->m_volume_fractions.size()) ? 1 : 0; @@ -730,18 +715,6 @@ class ElementDominantMaterialView m_zoneIndex == rhs.m_zoneIndex && m_view == rhs.m_view; } - bool AXOM_HOST_DEVICE operator < (const iterator &rhs) const - { - return m_currentIndex < rhs.m_currentIndex && - m_zoneIndex < rhs.m_zoneIndex && - m_view == rhs.m_view; - } - bool AXOM_HOST_DEVICE operator > (const iterator &rhs) const - { - return m_currentIndex > rhs.m_currentIndex && - m_zoneIndex > rhs.m_zoneIndex && - m_view == rhs.m_view; - } bool AXOM_HOST_DEVICE operator != (const iterator &rhs) const { return m_currentIndex != rhs.m_currentIndex || @@ -998,6 +971,139 @@ class MaterialDominantMaterialView return found; } + /*! + * \brief An iterator class for iterating over the mat/vf data in a zone. + */ + class iterator + { + // Let the material view call the iterator constructor. + friend class MaterialDominantMaterialView; + public: + /// Get the current material id for the iterator. + MaterialID AXOM_HOST_DEVICE material_id() const + { + SLIC_ASSERT(m_miIndex < m_view->m_size); + return m_view->m_matnos[m_miIndex]; + } + /// Get the current volume fraction for the iterator. + FloatType AXOM_HOST_DEVICE volume_fraction() const + { + SLIC_ASSERT(m_miIndex < m_view->m_size && m_index < m_view->m_element_ids[m_miIndex].size()); + return m_view->m_volume_fractions[m_miIndex][m_index]; + } + ZoneIndex AXOM_HOST_DEVICE zoneIndex() const { return m_zoneIndex; } + axom::IndexType AXOM_HOST_DEVICE size() const + { + return m_view->numberOfMaterials(m_zoneIndex); + } + void AXOM_HOST_DEVICE operator ++() + { + advance(true); + } + void AXOM_HOST_DEVICE operator ++(int) + { + advance(true); + } + bool AXOM_HOST_DEVICE operator == (const iterator &rhs) const + { + return m_miIndex == rhs.m_miIndex && + m_index == rhs.m_index && + m_zoneIndex == rhs.m_zoneIndex && + m_view == rhs.m_view; + } + bool AXOM_HOST_DEVICE operator != (const iterator &rhs) const + { + return m_miIndex != rhs.m_miIndex || + m_index != rhs.m_index || + m_zoneIndex != rhs.m_zoneIndex || + m_view != rhs.m_view; + } + private: + DISABLE_DEFAULT_CTOR(iterator); + + /// Constructor + AXOM_HOST_DEVICE iterator(const MaterialDominantMaterialView *view, + ZoneIndex zoneIndex, + axom::IndexType miIndex, + axom::IndexType index) : + m_view(view), + m_zoneIndex(zoneIndex), + m_miIndex(miIndex), + m_index(index) + { + } + + /// Advance to the next valid material slot for the zone. + void AXOM_HOST_DEVICE advance(bool doIncrement) + { + if(doIncrement) + { + if(m_miIndex < m_view->m_size) + { + m_index = 0; + m_miIndex++; + } + if(m_miIndex == m_view->m_size) + { + m_index = m_view->m_element_ids[m_view->m_size - 1].size(); + } + } + + // Look for the next m_miIndex,m_index pair that contains material for the selected zone index. + for(; m_miIndex < m_view->m_size; m_miIndex++) + { + const auto &element_ids = m_view->m_element_ids[m_miIndex]; + const auto sz = element_ids.size(); + for(; m_index < sz; m_index++) + { + if(element_ids[m_index] == m_zoneIndex) + { + return; + } + } + m_index = (m_miIndex + 1 == m_view->m_size) ? m_index : 0; + } + } + + const MaterialDominantMaterialView *m_view; + ZoneIndex m_zoneIndex; + axom::IndexType m_miIndex; + axom::IndexType m_index; + }; + // Let the iterator access members. + friend class iterator; + + /*! + * \brief Return the iterator for the beginning of a zone's material data. + * + * \param zi The zone index being queried. + * + * \return The iterator for the beginning of a zone's material data. + */ + iterator AXOM_HOST_DEVICE beginZone(ZoneIndex zi) const + { + SLIC_ASSERT(zi < static_cast(numberOfZones())); + + auto it = iterator(this, zi, 0, 0); + it.advance(false); + return it; + } + + /*! + * \brief Return the iterator for the end of a zone's material data. + * + * \param zi The zone index being queried. + * + * \return The iterator for the end of a zone's material data. + */ + iterator AXOM_HOST_DEVICE endZone(ZoneIndex zi) const + { + SLIC_ASSERT(zi < static_cast(numberOfZones())); + const axom::IndexType miIndex = m_size; + const axom::IndexType index = (m_size > 0) ? (m_volume_fractions[m_size - 1].size()) : 0; + return iterator(this, zi, miIndex, index); + } + private: AXOM_HOST_DEVICE axom::IndexType indexOfMaterialID(MaterialID mat) const From 0047a14f0e8af17fe259d28a54e518d7b14b553e Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Tue, 4 Nov 2025 17:44:13 -0800 Subject: [PATCH 06/35] Extend mir_concentric_circles with a traversal method that times matset view iterators. --- .../examples/concentric_circles/runMIR.hpp | 47 +++++++++++++++++++ 1 file changed, 47 insertions(+) diff --git a/src/axom/mir/examples/concentric_circles/runMIR.hpp b/src/axom/mir/examples/concentric_circles/runMIR.hpp index bcf033338c..4208676d94 100644 --- a/src/axom/mir/examples/concentric_circles/runMIR.hpp +++ b/src/axom/mir/examples/concentric_circles/runMIR.hpp @@ -10,6 +10,47 @@ #include "axom/bump.hpp" #include "axom/mir.hpp" // for Mir classes & functions +template +void test_matset_traveral(MatsetView matsetView) +{ + AXOM_ANNOTATE_SCOPE("test_matset_traversal"); + double vf1, vf2; + { + AXOM_ANNOTATE_SCOPE("zoneMaterials"); + axom::ReduceSum vfSum(0.); + axom::for_all(matsetView.numberOfZones(), AXOM_LAMBDA(axom::IndexType zoneIndex) + { + typename MatsetView::IDList ids; + typename MatsetView::VFList vfs; + matsetView.zoneMaterials(zoneIndex, ids, vfs); + double sum = 0.; + for(axom::IndexType i = 0; i < vfs.size(); i++) + { + sum += vfs[i]; + } + vfSum += sum; + }); + vf1 = vfSum.get(); + } + { + AXOM_ANNOTATE_SCOPE("iterators"); + axom::ReduceSum vfSum(0.); + axom::for_all(matsetView.numberOfZones(), AXOM_LAMBDA(axom::IndexType zoneIndex) + { + const auto end = matsetView.endZone(zoneIndex); + double sum = 0.; + for(auto it = matsetView.beginZone(zoneIndex); it != end; it++) + { + sum += it.volume_fraction(); + } + vfSum += sum; + }); + vf2 = vfSum.get(); + } + + std::cout << "test_matset_traversal: vf1=" << vf1 << ", vf2=" << vf2 << ", nzones=" << matsetView.numberOfZones() << std::endl; +} + template int runMIR(const conduit::Node &hostMesh, const conduit::Node &options, conduit::Node &hostResult) { @@ -98,6 +139,12 @@ int runMIR(const conduit::Node &hostMesh, const conduit::Node &options, conduit: MIR m(topologyView, coordsetView, matsetView); m.execute(deviceMesh, options, deviceResult); } + else if(method == "traversal") + { + using namespace axom::bump::views; + auto matsetView = make_unibuffer_matset::view(n_matset); + test_matset_traveral(matsetView); + } else { SLIC_ERROR(axom::fmt::format("Unsupported MIR method {}", method)); From 4aac6c99e70456f074358f282e2e47665a624066 Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Tue, 4 Nov 2025 18:12:12 -0800 Subject: [PATCH 07/35] Far faster UnibufferMaterialView::iterator --- src/axom/bump/views/MaterialView.hpp | 55 ++++++++++++---------------- 1 file changed, 24 insertions(+), 31 deletions(-) diff --git a/src/axom/bump/views/MaterialView.hpp b/src/axom/bump/views/MaterialView.hpp index 44eb186d49..74d24f94bb 100644 --- a/src/axom/bump/views/MaterialView.hpp +++ b/src/axom/bump/views/MaterialView.hpp @@ -175,67 +175,66 @@ class UnibufferMaterialView MaterialID AXOM_HOST_DEVICE material_id() const { SLIC_ASSERT(m_currentIndex < size()); - return m_material_ids[m_indices[m_currentIndex]]; + return m_view->m_material_ids[m_index]; } /// Get the current volume fraction for the iterator. FloatType AXOM_HOST_DEVICE volume_fraction() const { SLIC_ASSERT(m_currentIndex < size()); - return m_volume_fractions[m_indices[m_currentIndex]]; + return m_view->m_volume_fractions[m_index]; } axom::IndexType AXOM_HOST_DEVICE size() const { - return m_indices.size(); + return m_view->m_sizes[m_zoneIndex]; } ZoneIndex AXOM_HOST_DEVICE zoneIndex() const { return m_zoneIndex; } void AXOM_HOST_DEVICE operator ++() { - m_currentIndex += (m_currentIndex < size()) ? 1 : 0; + advance(true); } void AXOM_HOST_DEVICE operator ++(int) { - m_currentIndex += (m_currentIndex < size()) ? 1 : 0; + advance(true); } bool AXOM_HOST_DEVICE operator == (const iterator &rhs) const { return m_currentIndex == rhs.m_currentIndex && m_zoneIndex == rhs.m_zoneIndex && - m_material_ids == rhs.m_material_ids && - m_volume_fractions == rhs.m_volume_fractions && - m_indices == rhs.m_indices; + m_view == rhs.m_view; } bool AXOM_HOST_DEVICE operator != (const iterator &rhs) const { return m_currentIndex != rhs.m_currentIndex || m_zoneIndex != rhs.m_zoneIndex || - m_material_ids != rhs.m_material_ids || - m_volume_fractions != rhs.m_volume_fractions || - m_indices != rhs.m_indices; + m_view != rhs.m_view; } private: DISABLE_DEFAULT_CTOR(iterator); /// Constructor - AXOM_HOST_DEVICE iterator(const axom::ArrayView &material_ids, - const axom::ArrayView &volume_fractions, - const axom::ArrayView &indices, + AXOM_HOST_DEVICE iterator(const UnibufferMaterialView *view, ZoneIndex zoneIndex, axom::IndexType currentIndex = 0) : - m_material_ids(material_ids), - m_volume_fractions(volume_fractions), - m_indices(indices), + m_view(view), m_zoneIndex(zoneIndex), - m_currentIndex(currentIndex) + m_currentIndex(currentIndex), + m_index(0) { } - axom::ArrayView m_material_ids; - axom::ArrayView m_volume_fractions; - axom::ArrayView m_indices; + void advance(bool doIncrement) + { + m_currentIndex += (doIncrement && m_currentIndex < size()) ? 1 : 0; + m_index = m_view->m_indices[m_view->m_offsets[m_zoneIndex] + m_currentIndex]; + } + + const UnibufferMaterialView *m_view; ZoneIndex m_zoneIndex; axom::IndexType m_currentIndex; + axom::IndexType m_index; // not considered in ==, != }; + friend class iterator; /*! * \brief Return the iterator for the beginning of a zone's material data. @@ -248,11 +247,9 @@ class UnibufferMaterialView { SLIC_ASSERT(zi < static_cast(numberOfZones())); - return iterator(m_material_ids, - m_volume_fractions, - axom::ArrayView(m_indices.data() + m_offsets[zi], m_sizes[zi]), - zi, - 0); + auto it = iterator(this, zi, 0); + it.advance(false); + return it; } /*! @@ -266,11 +263,7 @@ class UnibufferMaterialView { SLIC_ASSERT(zi < static_cast(numberOfZones())); - return iterator(m_material_ids, - m_volume_fractions, - axom::ArrayView(m_indices.data() + m_offsets[zi], m_sizes[zi]), - zi, - m_sizes[zi]); + return iterator(this, zi, m_sizes[zi]); } private: From b81fbdeb2eae21d2124d88c31a05a8c4830507d7 Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Tue, 4 Nov 2025 18:12:34 -0800 Subject: [PATCH 08/35] Updated command line option doc. --- src/axom/mir/examples/concentric_circles/MIRApplication.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/axom/mir/examples/concentric_circles/MIRApplication.cpp b/src/axom/mir/examples/concentric_circles/MIRApplication.cpp index 92e86f4dcc..531d47a837 100644 --- a/src/axom/mir/examples/concentric_circles/MIRApplication.cpp +++ b/src/axom/mir/examples/concentric_circles/MIRApplication.cpp @@ -42,7 +42,7 @@ int MIRApplication::initialize(int argc, char **argv) app.add_option("--gridsize", gridSize) ->check(axom::CLI::PositiveNumber) ->description("The number of zones along an axis."); - app.add_option("--method", method)->description("The MIR method name (equiz, elvira)"); + app.add_option("--method", method)->description("The MIR method (or operation) name (equiz, elvira, traversal)"); app.add_option("--numcircles", numCircles) ->check(axom::CLI::PositiveNumber) ->description("The number of circles to use for material creation."); From 47fae970d367d46ebdbd2c4e240d38bb7dbb063e Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Tue, 4 Nov 2025 18:13:31 -0800 Subject: [PATCH 09/35] Added comment --- src/axom/bump/views/MaterialView.hpp | 1 + 1 file changed, 1 insertion(+) diff --git a/src/axom/bump/views/MaterialView.hpp b/src/axom/bump/views/MaterialView.hpp index 74d24f94bb..e7e1462574 100644 --- a/src/axom/bump/views/MaterialView.hpp +++ b/src/axom/bump/views/MaterialView.hpp @@ -234,6 +234,7 @@ class UnibufferMaterialView axom::IndexType m_currentIndex; axom::IndexType m_index; // not considered in ==, != }; + // Let the iterator access members. friend class iterator; /*! From 6aa28973bced0a335246eddb350c05a106e26514 Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Tue, 4 Nov 2025 18:17:16 -0800 Subject: [PATCH 10/35] make style --- src/axom/bump/tests/bump_views.cpp | 11 +- src/axom/bump/utilities/conduit_memory.hpp | 3 +- src/axom/bump/views/MaterialView.hpp | 164 +++++++----------- src/axom/core/memory_management.cpp | 3 +- .../concentric_circles/MIRApplication.cpp | 3 +- .../examples/concentric_circles/runMIR.hpp | 49 +++--- 6 files changed, 102 insertions(+), 131 deletions(-) diff --git a/src/axom/bump/tests/bump_views.cpp b/src/axom/bump/tests/bump_views.cpp index ee11e8fb50..9b4f7e10e8 100644 --- a/src/axom/bump/tests/bump_views.cpp +++ b/src/axom/bump/tests/bump_views.cpp @@ -576,17 +576,13 @@ struct test_braid2d_mat { axom::bump::views::dispatch_material_multibuffer( deviceMesh["matsets/mat"], - [&](auto matsetView) { - test_matsetview(nzones, matsetView, allocatorID); - }); + [&](auto matsetView) { test_matsetview(nzones, matsetView, allocatorID); }); } else if(mattype == "element_dominant") { axom::bump::views::dispatch_material_element_dominant( deviceMesh["matsets/mat"], - [&](auto matsetView) { - test_matsetview(nzones, matsetView, allocatorID); - }); + [&](auto matsetView) { test_matsetview(nzones, matsetView, allocatorID); }); } else if(mattype == "material_dominant") { @@ -702,8 +698,7 @@ struct test_braid2d_mat int i = 0; for(auto it = matsetView.beginZone(index); it != end; it++, i++) { - eq_count += (vfs[i] == it.volume_fraction() && - ids[i] == it.material_id()) ? 1 : 0; + eq_count += (vfs[i] == it.volume_fraction() && ids[i] == it.material_id()) ? 1 : 0; count++; } resultsView[index] = (eq_count == count) ? 1 : 0; diff --git a/src/axom/bump/utilities/conduit_memory.hpp b/src/axom/bump/utilities/conduit_memory.hpp index fd0810fbd3..52db3f139e 100644 --- a/src/axom/bump/utilities/conduit_memory.hpp +++ b/src/axom/bump/utilities/conduit_memory.hpp @@ -142,7 +142,8 @@ void copy(conduit::Node &dest, const conduit::Node &src) else { const int allocatorID = axom::getAllocatorIDFromPointer(src.data_ptr()); - bool deviceAllocated = (allocatorID == INVALID_ALLOCATOR_ID) ? false : isDeviceAllocator(allocatorID); + bool deviceAllocated = + (allocatorID == INVALID_ALLOCATOR_ID) ? false : isDeviceAllocator(allocatorID); if(deviceAllocated || (!src.dtype().is_string() && src.dtype().number_of_elements() > 1)) { // Allocate the node's memory in the right place. diff --git a/src/axom/bump/views/MaterialView.hpp b/src/axom/bump/views/MaterialView.hpp index e7e1462574..5a9607c9f2 100644 --- a/src/axom/bump/views/MaterialView.hpp +++ b/src/axom/bump/views/MaterialView.hpp @@ -170,6 +170,7 @@ class UnibufferMaterialView { // Let the material view call the iterator constructor. friend class UnibufferMaterialView; + public: /// Get the current material id for the iterator. MaterialID AXOM_HOST_DEVICE material_id() const @@ -183,45 +184,34 @@ class UnibufferMaterialView SLIC_ASSERT(m_currentIndex < size()); return m_view->m_volume_fractions[m_index]; } - axom::IndexType AXOM_HOST_DEVICE size() const - { - return m_view->m_sizes[m_zoneIndex]; - } + axom::IndexType AXOM_HOST_DEVICE size() const { return m_view->m_sizes[m_zoneIndex]; } ZoneIndex AXOM_HOST_DEVICE zoneIndex() const { return m_zoneIndex; } - void AXOM_HOST_DEVICE operator ++() - { - advance(true); - } - void AXOM_HOST_DEVICE operator ++(int) - { - advance(true); - } - bool AXOM_HOST_DEVICE operator == (const iterator &rhs) const + void AXOM_HOST_DEVICE operator++() { advance(true); } + void AXOM_HOST_DEVICE operator++(int) { advance(true); } + bool AXOM_HOST_DEVICE operator==(const iterator &rhs) const { - return m_currentIndex == rhs.m_currentIndex && - m_zoneIndex == rhs.m_zoneIndex && - m_view == rhs.m_view; + return m_currentIndex == rhs.m_currentIndex && m_zoneIndex == rhs.m_zoneIndex && + m_view == rhs.m_view; } - bool AXOM_HOST_DEVICE operator != (const iterator &rhs) const + bool AXOM_HOST_DEVICE operator!=(const iterator &rhs) const { - return m_currentIndex != rhs.m_currentIndex || - m_zoneIndex != rhs.m_zoneIndex || - m_view != rhs.m_view; + return m_currentIndex != rhs.m_currentIndex || m_zoneIndex != rhs.m_zoneIndex || + m_view != rhs.m_view; } + private: DISABLE_DEFAULT_CTOR(iterator); /// Constructor AXOM_HOST_DEVICE iterator(const UnibufferMaterialView *view, ZoneIndex zoneIndex, - axom::IndexType currentIndex = 0) : - m_view(view), - m_zoneIndex(zoneIndex), - m_currentIndex(currentIndex), - m_index(0) - { - } + axom::IndexType currentIndex = 0) + : m_view(view) + , m_zoneIndex(zoneIndex) + , m_currentIndex(currentIndex) + , m_index(0) + { } void advance(bool doIncrement) { @@ -232,7 +222,7 @@ class UnibufferMaterialView const UnibufferMaterialView *m_view; ZoneIndex m_zoneIndex; axom::IndexType m_currentIndex; - axom::IndexType m_index; // not considered in ==, != + axom::IndexType m_index; // not considered in ==, != }; // Let the iterator access members. friend class iterator; @@ -415,6 +405,7 @@ class MultiBufferMaterialView { // Let the material view call the iterator constructor. friend class MultiBufferMaterialView; + public: /// Get the current material id for the iterator. MaterialID AXOM_HOST_DEVICE material_id() const @@ -431,45 +422,40 @@ class MultiBufferMaterialView const auto idx = curIndices[m_zoneIndex]; return curValues[idx]; } - axom::IndexType AXOM_HOST_DEVICE size() const - { - return m_view->numberOfMaterials(m_zoneIndex); - } + axom::IndexType AXOM_HOST_DEVICE size() const { return m_view->numberOfMaterials(m_zoneIndex); } ZoneIndex AXOM_HOST_DEVICE zoneIndex() const { return m_zoneIndex; } - void AXOM_HOST_DEVICE operator ++() + void AXOM_HOST_DEVICE operator++() { m_currentIndex += (m_currentIndex < m_view->m_size) ? 1 : 0; advance(); } - void AXOM_HOST_DEVICE operator ++(int) + void AXOM_HOST_DEVICE operator++(int) { m_currentIndex += (m_currentIndex < m_view->m_size) ? 1 : 0; advance(); } - bool AXOM_HOST_DEVICE operator == (const iterator &rhs) const + bool AXOM_HOST_DEVICE operator==(const iterator &rhs) const { - return m_currentIndex == rhs.m_currentIndex && - m_zoneIndex == rhs.m_zoneIndex && - m_view == rhs.m_view; + return m_currentIndex == rhs.m_currentIndex && m_zoneIndex == rhs.m_zoneIndex && + m_view == rhs.m_view; } - bool AXOM_HOST_DEVICE operator != (const iterator &rhs) const + bool AXOM_HOST_DEVICE operator!=(const iterator &rhs) const { - return m_currentIndex != rhs.m_currentIndex || - m_zoneIndex != rhs.m_zoneIndex || - m_view != rhs.m_view; + return m_currentIndex != rhs.m_currentIndex || m_zoneIndex != rhs.m_zoneIndex || + m_view != rhs.m_view; } + private: DISABLE_DEFAULT_CTOR(iterator); /// Constructor AXOM_HOST_DEVICE iterator(const MultiBufferMaterialView *view, ZoneIndex zoneIndex, - axom::IndexType currentIndex = 0) : - m_view(view), - m_zoneIndex(zoneIndex), - m_currentIndex(currentIndex) - { - } + axom::IndexType currentIndex = 0) + : m_view(view) + , m_zoneIndex(zoneIndex) + , m_currentIndex(currentIndex) + { } /// Advance to the next valid material slot for the zone. void AXOM_HOST_DEVICE advance() @@ -675,6 +661,7 @@ class ElementDominantMaterialView { // Let the material view call the iterator constructor. friend class ElementDominantMaterialView; + public: /// Get the current material id for the iterator. MaterialID AXOM_HOST_DEVICE material_id() const @@ -688,45 +675,40 @@ class ElementDominantMaterialView SLIC_ASSERT(m_currentIndex < m_view->m_volume_fractions.size()); return m_view->m_volume_fractions[m_currentIndex][m_zoneIndex]; } - axom::IndexType AXOM_HOST_DEVICE size() const - { - return m_view->numberOfMaterials(m_zoneIndex); - } + axom::IndexType AXOM_HOST_DEVICE size() const { return m_view->numberOfMaterials(m_zoneIndex); } ZoneIndex AXOM_HOST_DEVICE zoneIndex() const { return m_zoneIndex; } - void AXOM_HOST_DEVICE operator ++() + void AXOM_HOST_DEVICE operator++() { m_currentIndex += (m_currentIndex < m_view->m_volume_fractions.size()) ? 1 : 0; advance(); } - void AXOM_HOST_DEVICE operator ++(int) + void AXOM_HOST_DEVICE operator++(int) { m_currentIndex += (m_currentIndex < m_view->m_volume_fractions.size()) ? 1 : 0; advance(); } - bool AXOM_HOST_DEVICE operator == (const iterator &rhs) const + bool AXOM_HOST_DEVICE operator==(const iterator &rhs) const { - return m_currentIndex == rhs.m_currentIndex && - m_zoneIndex == rhs.m_zoneIndex && - m_view == rhs.m_view; + return m_currentIndex == rhs.m_currentIndex && m_zoneIndex == rhs.m_zoneIndex && + m_view == rhs.m_view; } - bool AXOM_HOST_DEVICE operator != (const iterator &rhs) const + bool AXOM_HOST_DEVICE operator!=(const iterator &rhs) const { - return m_currentIndex != rhs.m_currentIndex || - m_zoneIndex != rhs.m_zoneIndex || - m_view != rhs.m_view; + return m_currentIndex != rhs.m_currentIndex || m_zoneIndex != rhs.m_zoneIndex || + m_view != rhs.m_view; } + private: DISABLE_DEFAULT_CTOR(iterator); /// Constructor AXOM_HOST_DEVICE iterator(const ElementDominantMaterialView *view, ZoneIndex zoneIndex, - axom::IndexType currentIndex = 0) : - m_view(view), - m_zoneIndex(zoneIndex), - m_currentIndex(currentIndex) - { - } + axom::IndexType currentIndex = 0) + : m_view(view) + , m_zoneIndex(zoneIndex) + , m_currentIndex(currentIndex) + { } /// Advance to the next valid material slot for the zone. void AXOM_HOST_DEVICE advance() @@ -972,6 +954,7 @@ class MaterialDominantMaterialView { // Let the material view call the iterator constructor. friend class MaterialDominantMaterialView; + public: /// Get the current material id for the iterator. MaterialID AXOM_HOST_DEVICE material_id() const @@ -986,32 +969,20 @@ class MaterialDominantMaterialView return m_view->m_volume_fractions[m_miIndex][m_index]; } ZoneIndex AXOM_HOST_DEVICE zoneIndex() const { return m_zoneIndex; } - axom::IndexType AXOM_HOST_DEVICE size() const - { - return m_view->numberOfMaterials(m_zoneIndex); - } - void AXOM_HOST_DEVICE operator ++() - { - advance(true); - } - void AXOM_HOST_DEVICE operator ++(int) - { - advance(true); - } - bool AXOM_HOST_DEVICE operator == (const iterator &rhs) const + axom::IndexType AXOM_HOST_DEVICE size() const { return m_view->numberOfMaterials(m_zoneIndex); } + void AXOM_HOST_DEVICE operator++() { advance(true); } + void AXOM_HOST_DEVICE operator++(int) { advance(true); } + bool AXOM_HOST_DEVICE operator==(const iterator &rhs) const { - return m_miIndex == rhs.m_miIndex && - m_index == rhs.m_index && - m_zoneIndex == rhs.m_zoneIndex && - m_view == rhs.m_view; + return m_miIndex == rhs.m_miIndex && m_index == rhs.m_index && + m_zoneIndex == rhs.m_zoneIndex && m_view == rhs.m_view; } - bool AXOM_HOST_DEVICE operator != (const iterator &rhs) const + bool AXOM_HOST_DEVICE operator!=(const iterator &rhs) const { - return m_miIndex != rhs.m_miIndex || - m_index != rhs.m_index || - m_zoneIndex != rhs.m_zoneIndex || - m_view != rhs.m_view; + return m_miIndex != rhs.m_miIndex || m_index != rhs.m_index || + m_zoneIndex != rhs.m_zoneIndex || m_view != rhs.m_view; } + private: DISABLE_DEFAULT_CTOR(iterator); @@ -1019,13 +990,12 @@ class MaterialDominantMaterialView AXOM_HOST_DEVICE iterator(const MaterialDominantMaterialView *view, ZoneIndex zoneIndex, axom::IndexType miIndex, - axom::IndexType index) : - m_view(view), - m_zoneIndex(zoneIndex), - m_miIndex(miIndex), - m_index(index) - { - } + axom::IndexType index) + : m_view(view) + , m_zoneIndex(zoneIndex) + , m_miIndex(miIndex) + , m_index(index) + { } /// Advance to the next valid material slot for the zone. void AXOM_HOST_DEVICE advance(bool doIncrement) diff --git a/src/axom/core/memory_management.cpp b/src/axom/core/memory_management.cpp index 5b8c592bae..c9ef6b0841 100644 --- a/src/axom/core/memory_management.cpp +++ b/src/axom/core/memory_management.cpp @@ -45,7 +45,8 @@ int getSharedMemoryAllocatorID() const auto name = axom::fmt::format("SHARED::{}", UMPIRE_DEFAULT_SHARED_MEMORY_RESOURCE); auto traits {umpire::get_default_resource_traits(name)}; traits.scope = umpire::MemoryResourceTraits::shared_scope::node; - auto axom_node_allocator {rm.makeResource(axom::fmt::format("{}::axom_node_allocator", name), traits)}; + auto axom_node_allocator { + rm.makeResource(axom::fmt::format("{}::axom_node_allocator", name), traits)}; auto axom_shared_allocator { rm.makeAllocator(allocatorName, axom_node_allocator)}; allocator_id = axom_shared_allocator.getId(); diff --git a/src/axom/mir/examples/concentric_circles/MIRApplication.cpp b/src/axom/mir/examples/concentric_circles/MIRApplication.cpp index 531d47a837..4d48f6c319 100644 --- a/src/axom/mir/examples/concentric_circles/MIRApplication.cpp +++ b/src/axom/mir/examples/concentric_circles/MIRApplication.cpp @@ -42,7 +42,8 @@ int MIRApplication::initialize(int argc, char **argv) app.add_option("--gridsize", gridSize) ->check(axom::CLI::PositiveNumber) ->description("The number of zones along an axis."); - app.add_option("--method", method)->description("The MIR method (or operation) name (equiz, elvira, traversal)"); + app.add_option("--method", method) + ->description("The MIR method (or operation) name (equiz, elvira, traversal)"); app.add_option("--numcircles", numCircles) ->check(axom::CLI::PositiveNumber) ->description("The number of circles to use for material creation."); diff --git a/src/axom/mir/examples/concentric_circles/runMIR.hpp b/src/axom/mir/examples/concentric_circles/runMIR.hpp index 4208676d94..ed81e5a06c 100644 --- a/src/axom/mir/examples/concentric_circles/runMIR.hpp +++ b/src/axom/mir/examples/concentric_circles/runMIR.hpp @@ -18,37 +18,40 @@ void test_matset_traveral(MatsetView matsetView) { AXOM_ANNOTATE_SCOPE("zoneMaterials"); axom::ReduceSum vfSum(0.); - axom::for_all(matsetView.numberOfZones(), AXOM_LAMBDA(axom::IndexType zoneIndex) - { - typename MatsetView::IDList ids; - typename MatsetView::VFList vfs; - matsetView.zoneMaterials(zoneIndex, ids, vfs); - double sum = 0.; - for(axom::IndexType i = 0; i < vfs.size(); i++) - { - sum += vfs[i]; - } - vfSum += sum; - }); + axom::for_all( + matsetView.numberOfZones(), + AXOM_LAMBDA(axom::IndexType zoneIndex) { + typename MatsetView::IDList ids; + typename MatsetView::VFList vfs; + matsetView.zoneMaterials(zoneIndex, ids, vfs); + double sum = 0.; + for(axom::IndexType i = 0; i < vfs.size(); i++) + { + sum += vfs[i]; + } + vfSum += sum; + }); vf1 = vfSum.get(); } { AXOM_ANNOTATE_SCOPE("iterators"); axom::ReduceSum vfSum(0.); - axom::for_all(matsetView.numberOfZones(), AXOM_LAMBDA(axom::IndexType zoneIndex) - { - const auto end = matsetView.endZone(zoneIndex); - double sum = 0.; - for(auto it = matsetView.beginZone(zoneIndex); it != end; it++) - { - sum += it.volume_fraction(); - } - vfSum += sum; - }); + axom::for_all( + matsetView.numberOfZones(), + AXOM_LAMBDA(axom::IndexType zoneIndex) { + const auto end = matsetView.endZone(zoneIndex); + double sum = 0.; + for(auto it = matsetView.beginZone(zoneIndex); it != end; it++) + { + sum += it.volume_fraction(); + } + vfSum += sum; + }); vf2 = vfSum.get(); } - std::cout << "test_matset_traversal: vf1=" << vf1 << ", vf2=" << vf2 << ", nzones=" << matsetView.numberOfZones() << std::endl; + std::cout << "test_matset_traversal: vf1=" << vf1 << ", vf2=" << vf2 + << ", nzones=" << matsetView.numberOfZones() << std::endl; } template From 1da959f4a737b150b61c98cc7d2acb6ba06beaed Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Wed, 5 Nov 2025 17:05:41 -0800 Subject: [PATCH 11/35] Added missing AXOM_HOST_DEVICE macro. --- src/axom/bump/views/MaterialView.hpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/axom/bump/views/MaterialView.hpp b/src/axom/bump/views/MaterialView.hpp index 5a9607c9f2..6c58fe3c70 100644 --- a/src/axom/bump/views/MaterialView.hpp +++ b/src/axom/bump/views/MaterialView.hpp @@ -213,7 +213,7 @@ class UnibufferMaterialView , m_index(0) { } - void advance(bool doIncrement) + void AXOM_HOST_DEVICE advance(bool doIncrement) { m_currentIndex += (doIncrement && m_currentIndex < size()) ? 1 : 0; m_index = m_view->m_indices[m_view->m_offsets[m_zoneIndex] + m_currentIndex]; From fe1a41f512693ac3390029d635f9fd9905ce38dd Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Wed, 5 Nov 2025 18:09:20 -0800 Subject: [PATCH 12/35] Fix an iterator. --- src/axom/bump/views/MaterialView.hpp | 6 +++++- 1 file changed, 5 insertions(+), 1 deletion(-) diff --git a/src/axom/bump/views/MaterialView.hpp b/src/axom/bump/views/MaterialView.hpp index 6c58fe3c70..b317aac74e 100644 --- a/src/axom/bump/views/MaterialView.hpp +++ b/src/axom/bump/views/MaterialView.hpp @@ -216,7 +216,11 @@ class UnibufferMaterialView void AXOM_HOST_DEVICE advance(bool doIncrement) { m_currentIndex += (doIncrement && m_currentIndex < size()) ? 1 : 0; - m_index = m_view->m_indices[m_view->m_offsets[m_zoneIndex] + m_currentIndex]; + const auto idx = m_view->m_offsets[m_zoneIndex] + m_currentIndex; + if(idx < m_view->m_indices.size()) + { + m_index = m_view->m_indices[idx]; + } } const UnibufferMaterialView *m_view; From 500639562fcde7c7ca7cb203196504fdd9910f7c Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Wed, 5 Nov 2025 18:09:41 -0800 Subject: [PATCH 13/35] Use iterator in MatsetSlicer. --- src/axom/bump/MatsetSlicer.hpp | 11 ++++------- 1 file changed, 4 insertions(+), 7 deletions(-) diff --git a/src/axom/bump/MatsetSlicer.hpp b/src/axom/bump/MatsetSlicer.hpp index e44145b544..6916333237 100644 --- a/src/axom/bump/MatsetSlicer.hpp +++ b/src/axom/bump/MatsetSlicer.hpp @@ -112,15 +112,12 @@ class MatsetSlicer const auto size = static_cast(sizesView[index]); const auto offset = offsetsView[index]; - typename MatsetView::IDList ids; - typename MatsetView::VFList vfs; - deviceMatsetView.zoneMaterials(deviceSelectedZonesView[index], ids, vfs); - - for(int i = 0; i < size; i++) + auto zoneMat = deviceMatsetView.beginZone(deviceSelectedZonesView[index]); + for(int i = 0; i < size; i++, zoneMat++) { const auto destIndex = offset + i; - materialIdsView[destIndex] = ids[i]; - volumeFractionsView[destIndex] = vfs[i]; + materialIdsView[destIndex] = zoneMat.material_id(); + volumeFractionsView[destIndex] = zoneMat.volume_fraction(); indicesView[destIndex] = destIndex; } }); From 674b300c51feca1719c174a4f84b2344eb4f61a5 Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Wed, 5 Nov 2025 18:13:24 -0800 Subject: [PATCH 14/35] Use matset iterator in MergeMeshes. --- src/axom/bump/MergeMeshes.hpp | 13 +++++-------- 1 file changed, 5 insertions(+), 8 deletions(-) diff --git a/src/axom/bump/MergeMeshes.hpp b/src/axom/bump/MergeMeshes.hpp index eb0ca81e3c..9d5ce9f08d 100644 --- a/src/axom/bump/MergeMeshes.hpp +++ b/src/axom/bump/MergeMeshes.hpp @@ -1811,8 +1811,6 @@ class MergeMeshesAndMatsets : public MergeMeshes axom::IndexType nzones, axom::IndexType zOffset) const { - using IDList = typename decltype(matsetView)::IDList; - using VFList = typename decltype(matsetView)::VFList; using MatID = typename decltype(matsetView)::IndexType; // Make some maps for renumbering material numbers. @@ -1844,19 +1842,18 @@ class MergeMeshesAndMatsets : public MergeMeshes nzones, AXOM_LAMBDA(axom::IndexType zoneIndex) { // Get this zone's materials. - IDList ids; - VFList vfs; - matsetView.zoneMaterials(zoneIndex, ids, vfs); + auto zoneMat = matsetView.beginZone(zoneIndex); + const auto nmats = zoneMat.size(); // Store the materials in the new material. const auto zoneStart = offsetsView[zOffset + zoneIndex]; - for(axom::IndexType mi = 0; mi < ids.size(); mi++) + for(axom::IndexType mi = 0; mi < nmats; mi++, zoneMat++) { const auto destIndex = zoneStart + mi; - volumeFractionsView[destIndex] = vfs[mi]; + volumeFractionsView[destIndex] = zoneMat.volume_fraction(); // Get the index of the material number in the local map. - const auto mapIndex = axom::utilities::binary_search(localView, ids[mi]); + const auto mapIndex = axom::utilities::binary_search(localView, zoneMat.material_id()); SLIC_ASSERT(mapIndex != -1); // We'll store the all materials number. const auto allMatno = allView[mapIndex]; From be4f7c9c5cdaf0caa5514ed6cca48cfa41024679 Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Wed, 5 Nov 2025 18:16:00 -0800 Subject: [PATCH 15/35] Use matset iterator in TopologyMapper. --- src/axom/bump/TopologyMapper.hpp | 8 +++----- 1 file changed, 3 insertions(+), 5 deletions(-) diff --git a/src/axom/bump/TopologyMapper.hpp b/src/axom/bump/TopologyMapper.hpp index 01c53374f8..84543bb3d4 100644 --- a/src/axom/bump/TopologyMapper.hpp +++ b/src/axom/bump/TopologyMapper.hpp @@ -653,11 +653,9 @@ class TopologyMapper // Get the src material - there should just be one because we assume // that a clean matset is being mapped. - typename SrcMatsetView::IDList zoneMatIds; - typename SrcMatsetView::VFList zoneMatVFs; - srcMatsetView.zoneMaterials(srcZone, zoneMatIds, zoneMatVFs); - SLIC_ASSERT(zoneMatIds.size() == 1); - const auto mat = zoneMatIds[0]; + auto zoneMat = srcMatsetView.beginZone(srcZone); + SLIC_ASSERT(zoneMat.size() == 1); + const auto mat = zoneMat.material_id(); #if defined(AXOM_DEBUG_TOPOLOGY_MAPPER) && !defined(AXOM_DEVICE_CODE) std::cout << "\tintersection:" << std::endl From c05797234e401f94f28dcb98c9f95b745577a87e Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Wed, 5 Nov 2025 19:25:25 -0800 Subject: [PATCH 16/35] Use iterators in test. --- src/axom/mir/tests/mir_elvira3d.cpp | 12 ++++-------- 1 file changed, 4 insertions(+), 8 deletions(-) diff --git a/src/axom/mir/tests/mir_elvira3d.cpp b/src/axom/mir/tests/mir_elvira3d.cpp index 52707f9824..9e27336301 100644 --- a/src/axom/mir/tests/mir_elvira3d.cpp +++ b/src/axom/mir/tests/mir_elvira3d.cpp @@ -300,20 +300,16 @@ struct test_Elvira3D axom::for_all( matsetView.numberOfZones(), AXOM_LAMBDA(axom::IndexType zi) { - // Get the materials for this zone. - typename MatsetView::IDList ids; - typename MatsetView::VFList vfs; - matsetView.zoneMaterials(zi, ids, vfs); - // Add the material volumes to the total volumes. - for(axom::IndexType i = 0; i < ids.size(); i++) + const auto end = matsetView.endZone(zi); + for(auto zoneMat = matsetView.beginZone(zi); zoneMat != end; zoneMat++) { - auto index = axom::utilities::binary_search(sortedIdsView, ids[i]); + auto index = axom::utilities::binary_search(sortedIdsView, zoneMat.material_id()); // RelWithDebInfo workaround - "sortedIdsView.size()" substitutes lambda capture device failure for "nmats" SLIC_ASSERT(index >= 0 && index < sortedIdsView.size()); // Use an atomic to sum the value. - axom::atomicAdd(totalVolumeView.data() + index, zoneVolumes[zi] * vfs[i]); + axom::atomicAdd(totalVolumeView.data() + index, zoneVolumes[zi] * zoneMat.volume_fraction()); } }); From de79be5aef059b0bfe3e872d7d656e5dac9f901a Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Wed, 5 Nov 2025 19:26:44 -0800 Subject: [PATCH 17/35] Add a new ArrayView-based zoneMaterials method in matset views. Use it inside ElviraAlgorithm. --- src/axom/bump/views/MaterialView.hpp | 82 ++++++++++++++++++++++++++++ src/axom/mir/ElviraAlgorithm.hpp | 24 ++++---- 2 files changed, 93 insertions(+), 13 deletions(-) diff --git a/src/axom/bump/views/MaterialView.hpp b/src/axom/bump/views/MaterialView.hpp index b317aac74e..cc07b7c440 100644 --- a/src/axom/bump/views/MaterialView.hpp +++ b/src/axom/bump/views/MaterialView.hpp @@ -136,6 +136,24 @@ class UnibufferMaterialView } } + AXOM_HOST_DEVICE + axom::IndexType zoneMaterials(ZoneIndex zi, axom::ArrayView &ids, axom::ArrayView &vfs) const + { + SLIC_ASSERT(zi < static_cast(numberOfZones())); + + const auto sz = numberOfMaterials(zi); + SLIC_ASSERT(sz <= ids.size()); + const auto offset = m_offsets[zi]; + for(axom::IndexType i = 0; i < sz; i++) + { + const auto idx = m_indices[offset + i]; + + ids[i] = m_material_ids[idx]; + vfs[i] = m_volume_fractions[idx]; + } + return sz; + } + AXOM_HOST_DEVICE bool zoneContainsMaterial(ZoneIndex zi, MaterialID mat) const { @@ -375,6 +393,29 @@ class MultiBufferMaterialView } } + AXOM_HOST_DEVICE + axom::IndexType zoneMaterials(ZoneIndex zi, axom::ArrayView &ids, axom::ArrayView &vfs) const + { + axom::IndexType n = 0; + for(axom::IndexType i = 0; i < m_size; i++) + { + const auto &curIndices = m_indices[i]; + const auto &curValues = m_values[i]; + + if(zi < static_cast(curIndices.size())) + { + const auto idx = curIndices[zi]; + if(curValues[idx] > 0) + { + ids[n] = m_matnos[i]; + vfs[n] = curValues[idx]; + n++; + } + } + } + return n; + } + AXOM_HOST_DEVICE bool zoneContainsMaterial(ZoneIndex zi, MaterialID mat) const { @@ -635,6 +676,24 @@ class ElementDominantMaterialView } } + AXOM_HOST_DEVICE + axom::IndexType zoneMaterials(ZoneIndex zi, axom::ArrayView &ids, axom::ArrayView &vfs) const + { + axom::IndexType n = 0; + for(axom::IndexType i = 0; i < m_volume_fractions.size(); i++) + { + const auto ¤tVF = m_volume_fractions[i]; + SLIC_ASSERT(zi < currentVF.size()); + if(currentVF[zi] > 0) + { + ids[n] = m_matnos[i]; + vfs[n] = currentVF[zi]; + n++; + } + } + return n; + } + AXOM_HOST_DEVICE bool zoneContainsMaterial(ZoneIndex zi, MaterialID mat) const { @@ -920,6 +979,29 @@ class MaterialDominantMaterialView } } + AXOM_HOST_DEVICE + axom::IndexType zoneMaterials(ZoneIndex zi, axom::ArrayView &ids, axom::ArrayView &vfs) const + { + axom::IndexType n = 0; + for(axom::IndexType mi = 0; mi < m_size; mi++) + { + const auto &element_ids = m_element_ids[mi]; + const auto &volume_fractions = m_volume_fractions[mi]; + const auto sz = element_ids.size(); + for(axom::IndexType i = 0; i < sz; i++) + { + if(element_ids[i] == zi) + { + ids[n] = m_matnos[mi]; + vfs[n] = volume_fractions[i]; + n++; + break; + } + } + } + return n; + } + AXOM_HOST_DEVICE bool zoneContainsMaterial(ZoneIndex zi, MaterialID mat) const { diff --git a/src/axom/mir/ElviraAlgorithm.hpp b/src/axom/mir/ElviraAlgorithm.hpp index ea8c9a6735..c6496ea084 100644 --- a/src/axom/mir/ElviraAlgorithm.hpp +++ b/src/axom/mir/ElviraAlgorithm.hpp @@ -614,11 +614,15 @@ class ElviraAlgorithm : public axom::mir::MIRAlgorithm axom::Array fragmentVFStencil(numFragmentsStencil, numFragmentsStencil, allocatorID); auto fragmentVFStencilView = fragmentVFStencil.view(); - // Sorted material ids for each zone. + // Sorted material ids / vfs for each zone. axom::Array sortedMaterialIds(numFragments, numFragments, allocatorID); + axom::Array sortedMaterialVfs(numFragments, + numFragments, + allocatorID); auto sortedMaterialIdsView = sortedMaterialIds.view(); + auto sortedMaterialVfsView = sortedMaterialVfs.view(); // Coordinate stencil data for each zone. const auto nzonesStencil = nzones * StencilSize; @@ -641,22 +645,16 @@ class ElviraAlgorithm : public axom::mir::MIRAlgorithm // Where to begin writing this zone's fragment data. const auto offset = matOffsetView[szIndex]; - // Get materials for this zone from the matset. - typename MatsetView::IDList ids; - typename MatsetView::VFList vfs; - deviceMatsetView.zoneMaterials(matZoneIndex, ids, vfs); + // Determine the views for this zone's material data. + axom::ArrayView ids(sortedMaterialIdsView.data() + offset, matCount); + axom::ArrayView vfs(sortedMaterialVfsView.data() + offset, matCount); + // Get materials for this zone from the matset, directly into the "sorted" views. + [[maybe_unused]] auto ids_size = deviceMatsetView.zoneMaterials(matZoneIndex, ids, vfs); // Reverse sort the materials by the volume fraction so the larger VFs are first. - SLIC_ASSERT(ids.size() == matCount); + SLIC_ASSERT(ids_size == matCount); axom::utilities::reverse_sort_multiple(vfs.data(), ids.data(), matCount); - // Save sorted ids in sortedMaterialIdsView. - for(axom::IndexType m = 0; m < matCount; m++) - { - const auto fragmentIndex = offset + m; - sortedMaterialIdsView[fragmentIndex] = ids[m]; - } - // Retrieve the stencil data from neighbor zones. auto logical = deviceTopologyView.indexing().IndexToLogicalIndex(zoneIndex); for(int si = 0; si < StencilSize; si++) From 6d709435f2fbf60e8cd8d19f8ccb8260118cf990 Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Thu, 6 Nov 2025 11:23:23 -0800 Subject: [PATCH 18/35] RELEASE-NOTES.md --- RELEASE-NOTES.md | 6 ++++++ 1 file changed, 6 insertions(+) diff --git a/RELEASE-NOTES.md b/RELEASE-NOTES.md index b8959062d6..c60f670edd 100644 --- a/RELEASE-NOTES.md +++ b/RELEASE-NOTES.md @@ -33,6 +33,12 @@ The Axom project release numbers follow [Semantic Versioning](http://semver.org/ - The maximum number of vertices allowed in polygon primitives can now be passed as a template argument to `axom::bump::TopologyMapper`, `axom::bump::PrimalAdaptor`, and `axom::mir::ElviraAlgorithm`. +- Material views in `axom::bump::views` were enhanced with `iterator` classes that + enable traversal of material data for zones so kernels do not need to use large fixed size + buffers inside kernels. +- Material views in `axom::bump::views` were enhanced with an overloaded `zoneMaterials()` + method that allows data to be gathered into `axom::ArrayView` objects. + ### Fixed ### Deprecated From 032fd01f1442e0aad7c523ab3d827c663b86fcb3 Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Thu, 6 Nov 2025 11:43:12 -0800 Subject: [PATCH 19/35] Added another zoneMaterials test --- src/axom/bump/tests/bump_views.cpp | 20 +++++++++++++++++++- 1 file changed, 19 insertions(+), 1 deletion(-) diff --git a/src/axom/bump/tests/bump_views.cpp b/src/axom/bump/tests/bump_views.cpp index 712dcbcc8a..d5aba7d9e3 100644 --- a/src/axom/bump/tests/bump_views.cpp +++ b/src/axom/bump/tests/bump_views.cpp @@ -702,7 +702,25 @@ struct test_braid2d_mat eq_count += (vfs[i] == it.volume_fraction() && ids[i] == it.material_id()) ? 1 : 0; count++; } - resultsView[index] = (eq_count == count) ? 1 : 0; + + // Test ArrayView version of zoneMaterials(). + using IndexType = typename MatsetView::IndexType; + using FloatType = typename MatsetView::FloatType; + constexpr int ARRAY_SIZE = 10; + IndexType idStorage[ARRAY_SIZE]; + FloatType vfStorage[ARRAY_SIZE]; + axom::ArrayView idView(idStorage, ARRAY_SIZE); + axom::ArrayView vfView(vfStorage, ARRAY_SIZE); + const auto nmats = matsetView.zoneMaterials(index, idView, vfView); + eq_count += (nmats == ids.size()) ? 1 : 0; + count++; + for(axom::IndexType j = 0; j < nmats; j++) + { + eq_count += (vfs[j] == vfView[j] && ids[j] == idView[j]) ? 1 : 0; + count++; + } + + resultsView[index] = (eq_count == count) ? 1 : 0; }); // Get containsView data to the host and compare results From 833ea28ac47077c9194c81b1ccccb2344ddb7622 Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Thu, 6 Nov 2025 18:25:25 -0800 Subject: [PATCH 20/35] Added heavily_mixed MIR example. --- src/axom/mir/examples/CMakeLists.txt | 1 + .../mir/examples/heavily_mixed/CMakeLists.txt | 61 +++ .../examples/heavily_mixed/HMApplication.cpp | 446 ++++++++++++++++++ .../examples/heavily_mixed/HMApplication.hpp | 73 +++ .../heavily_mixed/mir_heavily_mixed.cpp | 18 + .../mir/examples/heavily_mixed/runMIR.hpp | 121 +++++ .../examples/heavily_mixed/runMIR_cuda.cpp | 34 ++ .../mir/examples/heavily_mixed/runMIR_hip.cpp | 34 ++ .../mir/examples/heavily_mixed/runMIR_omp.cpp | 32 ++ .../mir/examples/heavily_mixed/runMIR_seq.cpp | 22 + 10 files changed, 842 insertions(+) create mode 100644 src/axom/mir/examples/heavily_mixed/CMakeLists.txt create mode 100644 src/axom/mir/examples/heavily_mixed/HMApplication.cpp create mode 100644 src/axom/mir/examples/heavily_mixed/HMApplication.hpp create mode 100644 src/axom/mir/examples/heavily_mixed/mir_heavily_mixed.cpp create mode 100644 src/axom/mir/examples/heavily_mixed/runMIR.hpp create mode 100644 src/axom/mir/examples/heavily_mixed/runMIR_cuda.cpp create mode 100644 src/axom/mir/examples/heavily_mixed/runMIR_hip.cpp create mode 100644 src/axom/mir/examples/heavily_mixed/runMIR_omp.cpp create mode 100644 src/axom/mir/examples/heavily_mixed/runMIR_seq.cpp diff --git a/src/axom/mir/examples/CMakeLists.txt b/src/axom/mir/examples/CMakeLists.txt index 538202e59d..fc7d11078d 100644 --- a/src/axom/mir/examples/CMakeLists.txt +++ b/src/axom/mir/examples/CMakeLists.txt @@ -4,4 +4,5 @@ # SPDX-License-Identifier: (BSD-3-Clause) add_subdirectory(concentric_circles) +add_subdirectory(heavily_mixed) add_subdirectory(tutorial_simple) diff --git a/src/axom/mir/examples/heavily_mixed/CMakeLists.txt b/src/axom/mir/examples/heavily_mixed/CMakeLists.txt new file mode 100644 index 0000000000..bd910aefab --- /dev/null +++ b/src/axom/mir/examples/heavily_mixed/CMakeLists.txt @@ -0,0 +1,61 @@ +# Copyright (c) 2017-2025, Lawrence Livermore National Security, LLC and +# other Axom Project Developers. See the top-level COPYRIGHT file for details. +# +# SPDX-License-Identifier: (BSD-3-Clause) + +set( mir_example_dependencies + core + slic + mir + bump + ) + +# Speed up compilation by separating out the execspaces into different runMIR_xxx.cpp files. +axom_add_executable( + NAME mir_heavily_mixed + SOURCES mir_heavily_mixed.cpp HMApplication.cpp runMIR_seq.cpp runMIR_omp.cpp runMIR_cuda.cpp runMIR_hip.cpp + OUTPUT_DIR ${EXAMPLE_OUTPUT_DIRECTORY} + DEPENDS_ON ${mir_example_dependencies} + FOLDER axom/mir/examples + ) + +if(AXOM_ENABLE_TESTS) + set (_policies "seq") + if(RAJA_FOUND AND UMPIRE_FOUND) + blt_list_append(TO _policies ELEMENTS "omp" IF AXOM_ENABLE_OPENMP) + blt_list_append(TO _policies ELEMENTS "cuda" IF AXOM_ENABLE_CUDA) + blt_list_append(TO _policies ELEMENTS "hip" IF AXOM_ENABLE_HIP) + endif() + + foreach(_policy ${_policies}) + set(_testname "mir_heavily_mixed_2D_${_policy}") + axom_add_test( + NAME ${_testname} + COMMAND mir_heavily_mixed + --dimension 2 + --dims 50,50 + --refinement 10 + --nmats 10 + --policy ${_policy} + --disable-write + ) + + set_tests_properties(${_testname} PROPERTIES + PASS_REGULAR_EXPRESSION "Material interface reconstruction time:") + + set(_testname "mir_heavily_mixed_3D_${_policy}") + axom_add_test( + NAME ${_testname} + COMMAND mir_heavily_mixed + --dimension 3 + --dims 50,50,10 + --refinement 10 + --nmats 10 + --policy ${_policy} + --disable-write + ) + + set_tests_properties(${_testname} PROPERTIES + PASS_REGULAR_EXPRESSION "Material interface reconstruction time:") + endforeach() +endif() diff --git a/src/axom/mir/examples/heavily_mixed/HMApplication.cpp b/src/axom/mir/examples/heavily_mixed/HMApplication.cpp new file mode 100644 index 0000000000..404b08167a --- /dev/null +++ b/src/axom/mir/examples/heavily_mixed/HMApplication.cpp @@ -0,0 +1,446 @@ +// Copyright (c) 2017-2025, Lawrence Livermore National Security, LLC and +// other Axom Project Developers. See the top-level COPYRIGHT file for details. +// +// SPDX-License-Identifier: (BSD-3-Clause) + +#include "axom/config.hpp" +#include "axom/core.hpp" // for axom macros +#include "axom/slic.hpp" +#include "axom/bump.hpp" +#include "axom/mir.hpp" // for Mir classes & functions +#include "runMIR.hpp" +#include "HMApplication.hpp" + +#include +#include + +#include + +using RuntimePolicy = axom::runtime_policy::Policy; + +namespace detail +{ +/*! + * \brief Turn a field on a fine mesh into a matset on the coarse mesh. + * + * \param topoName The name of the topology. + * \param dims The number of cells in each logical dimension. + * \param refinement The amount of refinement between coarse/fine meshes. + * \param n_coarse The coarse mesh (it gets the matset). + * \param n_field The field used for matset creation. + * \param nmats The number of materials to make. + */ +void heavily_mixed_matset(const std::string &topoName, int dims[3], int refinement, conduit::Node &n_coarse, const conduit::Node &n_field, int nmats) +{ + const auto fine = n_field.as_float64_accessor(); + int nzones = dims[0] * dims[1] * dims[2]; + int nslots = nzones * nmats; + std::vector vfs(nslots, 0.); + + // break the data range into nmats parts. + const double matSize = 1000. / nmats; //fine.max() / nmats; + + const int rdims[] = {dims[0] * refinement, dims[1] * refinement, dims[2] * refinement}; + + for(int k = 0; k < dims[2]; k++) + { + for(int j = 0; j < dims[1]; j++) + { + for(int i = 0; i < dims[0]; i++) + { + const int zoneIndex = k * dims[1] * dims[0] + j * dims[0] + i; + const int kr = k * refinement; + const int jr = j * refinement; + const int ir = i * refinement; + for(int jj = 0; jj < refinement; jj++) + for(int ii = 0; ii < refinement; ii++) + { + const int fine_index = (kr * rdims[0] * rdims[1]) + ((jr + jj) * rdims[0]) + (ir + ii); + const int matid = axom::utilities::clampVal(static_cast(fine[fine_index] / matSize), 0, nmats-1); + const int matslot = zoneIndex * nmats + matid; + vfs[matslot] += 1. / (refinement * refinement); + } + } + } + } + + std::vector material_ids; + std::vector volume_fractions; + std::vector indices; + std::vector sizes; + std::vector offsets; + for(int k = 0; k < dims[2]; k++) + { + for(int j = 0; j < dims[1]; j++) + { + for(int i = 0; i < dims[0]; i++) + { + const int zoneIndex = k * dims[0] * dims[1] + j * dims[0] + i; + + int size = 0; + offsets.push_back(indices.size()); + for(int m = 0; m < nmats; m++) + { + int matslot = zoneIndex * nmats + m; + if(vfs[matslot] > 0) + { + indices.push_back(material_ids.size()); + material_ids.push_back(m + 1); + volume_fractions.push_back(vfs[matslot]); + size++; + } + } + sizes.push_back(size); + } + } + } + conduit::Node &n_matset = n_coarse["matsets/mat"]; + n_matset["topology"] = topoName; + conduit::Node &n_material_map = n_matset["material_map"]; + for(int i = 0; i < nmats; i++) + { + int matno = i + 1; + std::stringstream ss; + ss << "mat" << matno; + n_material_map[ss.str()] = matno; + } + n_matset["material_ids"].set(material_ids); + n_matset["volume_fractions"].set(volume_fractions); + n_matset["indices"].set(indices); + n_matset["sizes"].set(sizes); + n_matset["offsets"].set(offsets); + + n_coarse["fields/nmats/association"] = "element"; + n_coarse["fields/nmats/topology"] = topoName; + n_coarse["fields/nmats/values"].set(sizes); +} + +template +void heavily_mixed(conduit::Node &n_mesh, int dims[3], int refinement, int nmats) +{ + const int rdims[] = {refinement * dims[0], refinement * dims[1], refinement * dims[2]}; + + // Default window + const conduit::float64 x_min = -0.6; + const conduit::float64 x_max = 0.6; + const conduit::float64 y_min = -0.5; + const conduit::float64 y_max = 0.5; + const conduit::float64 c_re = -0.5125; + const conduit::float64 c_im = 0.5213; + + conduit::blueprint::mesh::examples::julia(dims[0], + dims[1], + x_min, + x_max, + y_min, + y_max, + c_re, + c_im, + n_mesh); + if(dims[2] > 1) + { + // Add another dimension to the coordset. + const conduit::float64 z_min = 0.; + const conduit::float64 z_max = x_max - x_min; + std::vector z; + for(int i = 0; i <= dims[2]; i++) + { + const auto t = static_cast(i) / dims[2]; + const auto zc = axom::utilities::lerp(z_min, z_max, t); + z.push_back(zc); + } + n_mesh["coordsets/coords/values/z"].set(z); + + // Destination window + const conduit::float64 s = 0.9; + const conduit::float64 x1_min = x_min * s; + const conduit::float64 x1_max = x_max * s; + const conduit::float64 y1_min = y_min * s; + const conduit::float64 y1_max = y_max * s; + + conduit::Node n_field; + n_field.set(conduit::DataType::int32(rdims[0] * rdims[1] * rdims[2])); + conduit::int32 *destPtr = n_field.as_int32_ptr(); + axom::for_all(rdims[2], AXOM_LAMBDA(int k) + { + const auto t = static_cast(k) / (dims[2] - 1); + // Interpolate the window + const conduit::float64 x0 = axom::utilities::lerp(x_min, x1_min, t); + const conduit::float64 x1 = axom::utilities::lerp(x_max, x1_max, t); + const conduit::float64 y0 = axom::utilities::lerp(y_min, y1_min, t); + const conduit::float64 y1 = axom::utilities::lerp(y_max, y1_max, t); + conduit::Node n_rmesh; + conduit::blueprint::mesh::examples::julia(rdims[0], + rdims[1], + x0, + x1, + y0, + y1, + c_re, + c_im, + n_rmesh); + const conduit::Node &n_src_field = n_rmesh["fields/iters/values"]; + const conduit::int32 *srcPtr = n_src_field.as_int32_ptr(); + conduit::int32 *currentDestPtr = destPtr + k * rdims[0] * rdims[1]; + axom::copy(currentDestPtr, srcPtr, rdims[0] * rdims[1] * sizeof(conduit::int32)); + SLIC_INFO(axom::fmt::format("Made slice {}/{}", k+1, rdims[2])); + }); + + // Make a matset based on the higher resolution julia field. + heavily_mixed_matset("topo", dims, refinement, n_mesh, n_field, nmats); + } + else + { + // Generate the same julia set at higher resolution to use as materials. + conduit::Node n_rmesh; + conduit::blueprint::mesh::examples::julia(rdims[0], + rdims[1], + x_min, + x_max, + y_min, + y_max, + c_re, + c_im, + n_rmesh); + + // Make a matset based on the higher resolution julia field. + const conduit::Node &n_field = n_rmesh["fields/iters/values"]; + heavily_mixed_matset("topo", dims, refinement, n_mesh, n_field, nmats); + } +} + +} // end namespace detail + +//-------------------------------------------------------------------------------- +HMApplication::HMApplication() + : m_handler(true) + , m_dims{100, 100, 1} + , m_numMaterials(40) + , m_refinement(40) + , m_numTrials(1) + , m_writeFiles(true) + , m_outputFilePath("output") + , m_method("elvira") + , m_policy(RuntimePolicy::seq) + , m_annotationMode("report") +{ } + +//-------------------------------------------------------------------------------- +int HMApplication::initialize(int argc, char **argv) +{ + axom::CLI::App app; + app.add_flag("--handler", m_handler) + ->description("Install a custom error handler that loops forever.") + ->capture_default_str(); + + std::vector dims; + app.add_option("--dims", dims, "Dimensions in x,y,z (z optional)") + ->expected(2, 3) + ->check(axom::CLI::PositiveNumber); + + app.add_option("--method", m_method) + ->description("The MIR method (or operation) name (equiz, elvira, traversal)"); + app.add_option("--materials", m_numMaterials) + ->check(axom::CLI::PositiveNumber) + ->description("The number of materials to create."); + app.add_option("--refinement", m_refinement) + ->check(axom::CLI::PositiveNumber) + ->description("The refinement within a zone when creating materials."); + app.add_option("--output", m_outputFilePath) + ->description("The file path for HDF5/YAML output files"); + bool disable_write = !m_writeFiles; + app.add_flag("--disable-write", disable_write)->description("Disable writing data files"); + app.add_option("--trials", m_numTrials) + ->check(axom::CLI::PositiveNumber) + ->description("The number of MIR trials to run on the mesh."); + +#if defined(AXOM_USE_CALIPER) + app.add_option("--caliper", m_annotationMode) + ->description( + "caliper annotation mode. Valid options include 'none' and 'report'. " + "Use 'help' to see full list.") + ->capture_default_str() + ->check(axom::utilities::ValidCaliperMode); +#endif + + std::stringstream pol_sstr; + pol_sstr << "Set MIR runtime policy method."; +#if defined(AXOM_USE_RAJA) && defined(AXOM_USE_UMPIRE) + pol_sstr << "\nSet to 'seq' or 0 to use the RAJA sequential policy."; + #ifdef AXOM_USE_OPENMP + pol_sstr << "\nSet to 'omp' or 1 to use the RAJA OpenMP policy."; + #endif + #ifdef AXOM_USE_CUDA + pol_sstr << "\nSet to 'cuda' or 2 to use the RAJA CUDA policy."; + #endif + #ifdef AXOM_USE_HIP + pol_sstr << "\nSet to 'hip' or 3 to use the RAJA HIP policy."; + #endif +#endif + app.add_option("-p, --policy", m_policy, pol_sstr.str()) + ->capture_default_str() + ->transform(axom::CLI::CheckedTransformer(axom::runtime_policy::s_nameToPolicy)); + + // Parse command line options. + int retval = 0; + try + { + app.parse(argc, argv); + m_writeFiles = !disable_write; + } + catch(axom::CLI::CallForHelp &e) + { + std::cout << app.help() << std::endl; + retval = -1; + } + catch(axom::CLI::ParseError &e) + { + // Handle other parsing errors + std::cerr << e.what() << std::endl; + retval = -2; + } + + for(size_t i = 0; i < axom::utilities::min(dims.size(), static_cast(3)); i++) + { + m_dims[i] = dims[i]; + } + + return retval; +} + +//-------------------------------------------------------------------------------- +int HMApplication::execute() +{ + axom::slic::SimpleLogger logger(axom::slic::message::Info); + + if(m_handler) + { + conduit::utils::set_error_handler(conduit_debug_err_handler); + } +#if defined(AXOM_USE_CALIPER) + axom::utilities::raii::AnnotationsWrapper annotations_raii_wrapper(m_annotationMode); +#endif + int retval = 0; + try + { + retval = runMIR(); + } + catch(std::invalid_argument const &e) + { + SLIC_WARNING("Bad input. " << e.what()); + retval = -2; + } + catch(std::out_of_range const &e) + { + SLIC_WARNING("Integer overflow. " << e.what()); + retval = -3; + } + return retval; +} + +//-------------------------------------------------------------------------------- +int HMApplication::runMIR() +{ + // Initialize a mesh for testing MIR + auto timer = axom::utilities::Timer(true); + conduit::Node mesh; + { + AXOM_ANNOTATE_SCOPE("generate"); + int dims[3]; + dims[0] = m_dims[0]; + dims[1] = m_dims[1]; + dims[2] = m_dims[2]; + SLIC_INFO(axom::fmt::format("dims: {},{},{}", dims[0], dims[1], dims[2])); + SLIC_INFO(axom::fmt::format("refinement: {}", m_refinement)); + SLIC_INFO(axom::fmt::format("numMaterials: {}", m_numMaterials)); +#if defined(AXOM_USE_RAJA) && defined(AXOM_USE_UMPIRE) && defined(AXOM_USE_OPENMP) + using CPUExecSpace = axom::OMP_EXEC; +#else + using CPUExecSpace = axom::SEQ_EXEC; +#endif + detail::heavily_mixed(mesh, dims, m_refinement, m_numMaterials); + } + timer.stop(); + SLIC_INFO("Mesh init time: " << timer.elapsedTimeInMilliSec() << " ms."); + + // Output initial mesh. + if(m_writeFiles) + { + saveMesh(mesh, "heavily_mixed"); + } + + // Begin material interface reconstruction + timer.reset(); + timer.start(); + conduit::Node options, resultMesh; + options["matset"] = "mat"; + options["method"] = m_method; // pass method via options. + options["trials"] = m_numTrials; + + const int dimension = (m_dims[2] > 1) ? 3 : 2; + int retval = 0; + if(m_policy == RuntimePolicy::seq) + { + retval = runMIR_seq(dimension, mesh, options, resultMesh); + } +#if defined(AXOM_USE_RAJA) && defined(AXOM_USE_UMPIRE) + #if defined(AXOM_USE_OPENMP) + else if(m_policy == RuntimePolicy::omp) + { + retval = runMIR_omp(dimension, mesh, options, resultMesh); + } + #endif + #if defined(AXOM_USE_CUDA) + else if(m_policy == RuntimePolicy::cuda) + { + constexpr int CUDA_BLOCK_SIZE = 256; + using cuda_exec = axom::CUDA_EXEC; + retval = runMIR_cuda(dimension, mesh, options, resultMesh); + } + #endif + #if defined(AXOM_USE_HIP) + else if(m_policy == RuntimePolicy::hip) + { + retval = runMIR_hip(dimension, mesh, options, resultMesh); + } + #endif +#endif + else + { + retval = -1; + SLIC_ERROR("Unhandled policy."); + } + timer.stop(); + SLIC_INFO("Material interface reconstruction time: " << timer.elapsedTimeInMilliSec() << " ms."); + + // Output results + if(m_writeFiles) + { + AXOM_ANNOTATE_SCOPE("save_output"); + saveMesh(resultMesh, m_outputFilePath); + } + + return retval; +} + +//-------------------------------------------------------------------------------- +void HMApplication::adjustMesh(conduit::Node &) { } + +//-------------------------------------------------------------------------------- +void HMApplication::saveMesh(const conduit::Node &n_mesh, const std::string &path) +{ +#if defined(CONDUIT_RELAY_IO_HDF5_ENABLED) + std::string protocol("hdf5"); +#else + std::string protocol("yaml"); +#endif + conduit::relay::io::blueprint::save_mesh(n_mesh, path, protocol); +} + +//-------------------------------------------------------------------------------- +void HMApplication::conduit_debug_err_handler(const std::string &s1, const std::string &s2, int i1) +{ + SLIC_ERROR(axom::fmt::format("Error from Conduit: s1={}, s2={}, i1={}", s1, s2, i1)); + // This is on purpose. + while(1); +} diff --git a/src/axom/mir/examples/heavily_mixed/HMApplication.hpp b/src/axom/mir/examples/heavily_mixed/HMApplication.hpp new file mode 100644 index 0000000000..2f8fc7787c --- /dev/null +++ b/src/axom/mir/examples/heavily_mixed/HMApplication.hpp @@ -0,0 +1,73 @@ +// Copyright (c) 2017-2025, Lawrence Livermore National Security, LLC and +// other Axom Project Developers. See the top-level COPYRIGHT file for details. +// +// SPDX-License-Identifier: (BSD-3-Clause) +#ifndef AXOM_MIR_EXAMPLES_HEAVILY_MIXED_APPLICATION_HPP +#define AXOM_MIR_EXAMPLES_HEAVILY_MIXED_APPLICATION_HPP +#include "axom/config.hpp" +#include "axom/core.hpp" + +#include + +#include + +/*! + * \brief The heavily mixed materials application. + */ +class HMApplication +{ +public: + /// Constructor + HMApplication(); + + /*! + * \brief Initialize the application from command line args. + * \return 0 on success; less than zero otherwise. + */ + int initialize(int argc, char **argv); + + /*! + * \brief Execute the main application logic. + * \return 0 on success; less than zero otherwise. + */ + int execute(); + +protected: + /*! + * \brief Invoke the MIR appropriate for the selected runtime policy. + * \return 0 on success; less than zero otherwise. + */ + int runMIR(); + + /*! + * \brief Make any adjustments to the mesh. + */ + virtual void adjustMesh(conduit::Node &); + + /*! + * \brief Save the mesh to a file. + * + * \param path The filepath where the file will be saved. + * \param n_mesh The mesh to be saved. + */ + virtual void saveMesh(const conduit::Node &n_mesh, const std::string &path); + + /*! + * \brief A static error handler for Conduit. + */ + static void conduit_debug_err_handler(const std::string &s1, const std::string &s2, int i1); + + bool m_handler; + axom::StackArray m_dims; + int m_numMaterials; + int m_refinement; + int m_numTrials; + bool m_writeFiles; + std::string m_outputFilePath; + std::string m_method; + axom::runtime_policy::Policy m_policy; + std::string m_annotationMode; + std::string m_protocol; +}; + +#endif diff --git a/src/axom/mir/examples/heavily_mixed/mir_heavily_mixed.cpp b/src/axom/mir/examples/heavily_mixed/mir_heavily_mixed.cpp new file mode 100644 index 0000000000..3df9133de0 --- /dev/null +++ b/src/axom/mir/examples/heavily_mixed/mir_heavily_mixed.cpp @@ -0,0 +1,18 @@ +// Copyright (c) 2017-2025, Lawrence Livermore National Security, LLC and +// other Axom Project Developers. See the top-level COPYRIGHT file for details. +// +// SPDX-License-Identifier: (BSD-3-Clause) + +#include "HMApplication.hpp" + +int main(int argc, char **argv) +{ + HMApplication app; + int retval = app.initialize(argc, argv); + if(retval == 0) + { + retval = app.execute(); + } + + return retval; +} diff --git a/src/axom/mir/examples/heavily_mixed/runMIR.hpp b/src/axom/mir/examples/heavily_mixed/runMIR.hpp new file mode 100644 index 0000000000..adff6b363f --- /dev/null +++ b/src/axom/mir/examples/heavily_mixed/runMIR.hpp @@ -0,0 +1,121 @@ +// Copyright (c) 2017-2025, Lawrence Livermore National Security, LLC and +// other Axom Project Developers. See the top-level COPYRIGHT file for details. +// +// SPDX-License-Identifier: (BSD-3-Clause) +#ifndef AXOM_MIR_EXAMPLES_HEAVILY_MIXED_RUNMIR_HPP +#define AXOM_MIR_EXAMPLES_HEAVILY_MIXED_RUNMIR_HPP +#include "axom/config.hpp" +#include "axom/core.hpp" +#include "axom/slic.hpp" +#include "axom/bump.hpp" +#include "axom/mir.hpp" + +template +int runMIR(const conduit::Node &hostMesh, const conduit::Node &options, conduit::Node &hostResult) +{ + AXOM_ANNOTATE_SCOPE("runMIR"); + + namespace utils = axom::bump::utilities; + + // Pick the method out of the options. + std::string method("equiz"); + if(options.has_child("method")) + { + method = options["method"].as_string(); + } + SLIC_INFO(axom::fmt::format("Using policy {} for {} {}D", + axom::execution_space::name(), + method, + NDIMS)); + + // Get the number of times we want to run MIR. + int trials = 1; + if(options.has_child("trials")) + { + trials = std::max(1, options["trials"].to_int()); + } + + // TODO: set MAXMATERIALS to 1, remove warning. + // Check materials. + constexpr int MAXMATERIALS = 100; + auto materialInfo = axom::bump::views::materials(hostMesh["matsets/mat"]); + if(materialInfo.size() >= MAXMATERIALS) + { + SLIC_WARNING( + axom::fmt::format("To use more than {} materials, recompile with " + "larger MAXMATERIALS value.", + MAXMATERIALS)); + return -4; + } + + conduit::Node deviceMesh; + { + AXOM_ANNOTATE_SCOPE("host->device"); + utils::copy(deviceMesh, hostMesh); + } + + const conduit::Node &n_coordset = deviceMesh["coordsets/coords"]; + const conduit::Node &n_topology = deviceMesh["topologies/topo"]; + const conduit::Node &n_matset = deviceMesh["matsets/mat"]; + conduit::Node deviceResult; + for(int trial = 0; trial < trials; trial++) + { + deviceResult.reset(); + + // Make views + using namespace axom::bump::views; + auto coordsetView = make_rectilinear_coordset::view(n_coordset); + using CoordsetView = decltype(coordsetView); + + auto topologyView = make_rectilinear_topology::view(n_topology); + using TopologyView = decltype(topologyView); + + auto matsetView = make_unibuffer_matset::view(n_matset); + using MatsetView = decltype(matsetView); + + if(method == "equiz") + { + using MIR = axom::mir::EquiZAlgorithm; + MIR m(topologyView, coordsetView, matsetView); + m.execute(deviceMesh, options, deviceResult); + } + else if(method == "elvira") + { + using IndexingPolicy = typename TopologyView::IndexingPolicy; + using MIR = axom::mir::ElviraAlgorithm; + MIR m(topologyView, coordsetView, matsetView); + m.execute(deviceMesh, options, deviceResult); + } + else + { + SLIC_ERROR(axom::fmt::format("Unsupported MIR method {}", method)); + } + } + + { + AXOM_ANNOTATE_SCOPE("device->host"); + utils::copy(hostResult, deviceResult); + } + + return 0; +} + +// Prototypes. +int runMIR_seq(int dimension, + const conduit::Node &mesh, + const conduit::Node &options, + conduit::Node &result); +int runMIR_omp(int dimension, + const conduit::Node &mesh, + const conduit::Node &options, + conduit::Node &result); +int runMIR_cuda(int dimension, + const conduit::Node &mesh, + const conduit::Node &options, + conduit::Node &result); +int runMIR_hip(int dimension, + const conduit::Node &mesh, + const conduit::Node &options, + conduit::Node &result); + +#endif diff --git a/src/axom/mir/examples/heavily_mixed/runMIR_cuda.cpp b/src/axom/mir/examples/heavily_mixed/runMIR_cuda.cpp new file mode 100644 index 0000000000..e7ef3f0188 --- /dev/null +++ b/src/axom/mir/examples/heavily_mixed/runMIR_cuda.cpp @@ -0,0 +1,34 @@ +// Copyright (c) 2017-2025, Lawrence Livermore National Security, LLC and +// other Axom Project Developers. See the top-level COPYRIGHT file for details. +// +// SPDX-License-Identifier: (BSD-3-Clause) +#include "runMIR.hpp" + +#if defined(AXOM_USE_RAJA) && defined(AXOM_USE_UMPIRE) && defined(AXOM_USE_CUDA) +int runMIR_cuda(int dimension, + const conduit::Node &mesh, + const conduit::Node &options, + conduit::Node &result) +{ + constexpr int CUDA_BLOCK_SIZE = 256; + using cuda_exec = axom::CUDA_EXEC; + int retval = 0; + if(dimension == 3) + { + retval = runMIR(mesh, options, result); + } + else + { + retval = runMIR(mesh, options, result); + } + return retval; +} +#else +int runMIR_cuda(int AXOM_UNUSED_PARAM(dimension), + const conduit::Node &AXOM_UNUSED_PARAM(mesh), + const conduit::Node &AXOM_UNUSED_PARAM(options), + conduit::Node &AXOM_UNUSED_PARAM(result)) +{ + return 0; +} +#endif diff --git a/src/axom/mir/examples/heavily_mixed/runMIR_hip.cpp b/src/axom/mir/examples/heavily_mixed/runMIR_hip.cpp new file mode 100644 index 0000000000..b03bff1b8f --- /dev/null +++ b/src/axom/mir/examples/heavily_mixed/runMIR_hip.cpp @@ -0,0 +1,34 @@ +// Copyright (c) 2017-2025, Lawrence Livermore National Security, LLC and +// other Axom Project Developers. See the top-level COPYRIGHT file for details. +// +// SPDX-License-Identifier: (BSD-3-Clause) +#include "runMIR.hpp" + +#if defined(AXOM_USE_RAJA) && defined(AXOM_USE_UMPIRE) && defined(AXOM_USE_HIP) +int runMIR_hip(int dimension, + const conduit::Node &mesh, + const conduit::Node &options, + conduit::Node &result) +{ + constexpr int HIP_BLOCK_SIZE = 64; + using hip_exec = axom::HIP_EXEC; + int retval = 0; + if(dimension == 3) + { + retval = runMIR(mesh, options, result); + } + else + { + retval = runMIR(mesh, options, result); + } + return retval; +} +#else +int runMIR_hip(int AXOM_UNUSED_PARAM(dimension), + const conduit::Node &AXOM_UNUSED_PARAM(mesh), + const conduit::Node &AXOM_UNUSED_PARAM(options), + conduit::Node &AXOM_UNUSED_PARAM(result)) +{ + return 0; +} +#endif diff --git a/src/axom/mir/examples/heavily_mixed/runMIR_omp.cpp b/src/axom/mir/examples/heavily_mixed/runMIR_omp.cpp new file mode 100644 index 0000000000..6b57981575 --- /dev/null +++ b/src/axom/mir/examples/heavily_mixed/runMIR_omp.cpp @@ -0,0 +1,32 @@ +// Copyright (c) 2017-2025, Lawrence Livermore National Security, LLC and +// other Axom Project Developers. See the top-level COPYRIGHT file for details. +// +// SPDX-License-Identifier: (BSD-3-Clause) +#include "runMIR.hpp" + +#if defined(AXOM_USE_RAJA) && defined(AXOM_USE_UMPIRE) && defined(AXOM_USE_OPENMP) +int runMIR_omp(int dimension, + const conduit::Node &mesh, + const conduit::Node &options, + conduit::Node &result) +{ + int retval = 0; + if(dimension == 3) + { + retval = runMIR(mesh, options, result); + } + else + { + retval = runMIR(mesh, options, result); + } + return retval; +} +#else +int runMIR_omp(int AXOM_UNUSED_PARAM(dimension), + const conduit::Node &AXOM_UNUSED_PARAM(mesh), + const conduit::Node &AXOM_UNUSED_PARAM(options), + conduit::Node &AXOM_UNUSED_PARAM(result)) +{ + return 0; +} +#endif diff --git a/src/axom/mir/examples/heavily_mixed/runMIR_seq.cpp b/src/axom/mir/examples/heavily_mixed/runMIR_seq.cpp new file mode 100644 index 0000000000..af38756a9c --- /dev/null +++ b/src/axom/mir/examples/heavily_mixed/runMIR_seq.cpp @@ -0,0 +1,22 @@ +// Copyright (c) 2017-2025, Lawrence Livermore National Security, LLC and +// other Axom Project Developers. See the top-level COPYRIGHT file for details. +// +// SPDX-License-Identifier: (BSD-3-Clause) +#include "runMIR.hpp" + +int runMIR_seq(int dimension, + const conduit::Node &mesh, + const conduit::Node &options, + conduit::Node &result) +{ + int retval = 0; + if(dimension == 3) + { + retval = runMIR(mesh, options, result); + } + else + { + retval = runMIR(mesh, options, result); + } + return retval; +} From f4e87d119decc266947ef4a118c8415acb496164 Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Thu, 6 Nov 2025 18:29:44 -0800 Subject: [PATCH 21/35] CMake testing fixes. --- .../examples/concentric_circles/MIRApplication.cpp | 1 + src/axom/mir/examples/heavily_mixed/CMakeLists.txt | 12 +++++------- 2 files changed, 6 insertions(+), 7 deletions(-) diff --git a/src/axom/mir/examples/concentric_circles/MIRApplication.cpp b/src/axom/mir/examples/concentric_circles/MIRApplication.cpp index 4d48f6c319..9d33a486b7 100644 --- a/src/axom/mir/examples/concentric_circles/MIRApplication.cpp +++ b/src/axom/mir/examples/concentric_circles/MIRApplication.cpp @@ -173,6 +173,7 @@ int MIRApplication::runMIR() } // Begin material interface reconstruction + timer.reset(); timer.start(); conduit::Node options, resultMesh; options["matset"] = "mat"; diff --git a/src/axom/mir/examples/heavily_mixed/CMakeLists.txt b/src/axom/mir/examples/heavily_mixed/CMakeLists.txt index bd910aefab..4729454873 100644 --- a/src/axom/mir/examples/heavily_mixed/CMakeLists.txt +++ b/src/axom/mir/examples/heavily_mixed/CMakeLists.txt @@ -32,10 +32,9 @@ if(AXOM_ENABLE_TESTS) axom_add_test( NAME ${_testname} COMMAND mir_heavily_mixed - --dimension 2 - --dims 50,50 - --refinement 10 - --nmats 10 + --dims 100 100 + --refinement 20 + --materials 10 --policy ${_policy} --disable-write ) @@ -47,10 +46,9 @@ if(AXOM_ENABLE_TESTS) axom_add_test( NAME ${_testname} COMMAND mir_heavily_mixed - --dimension 3 - --dims 50,50,10 + --dims 50 50 10 --refinement 10 - --nmats 10 + --materials 10 --policy ${_policy} --disable-write ) From 86088154bd33166919ed0d60acd04795b6823000 Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Fri, 7 Nov 2025 17:01:42 -0800 Subject: [PATCH 22/35] Update release notes --- RELEASE-NOTES.md | 2 ++ 1 file changed, 2 insertions(+) diff --git a/RELEASE-NOTES.md b/RELEASE-NOTES.md index c60f670edd..28cfceea95 100644 --- a/RELEASE-NOTES.md +++ b/RELEASE-NOTES.md @@ -38,6 +38,8 @@ The Axom project release numbers follow [Semantic Versioning](http://semver.org/ buffers inside kernels. - Material views in `axom::bump::views` were enhanced with an overloaded `zoneMaterials()` method that allows data to be gathered into `axom::ArrayView` objects. +- A new `heavily_mixed` example program was added in `axom::mir` to demonstrate running MIR on + meshes with heavily mixed zones. ### Fixed From 5ef3e68d1f6fd770e4507bd667dbe42d9dc1a7f7 Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Fri, 7 Nov 2025 17:02:54 -0800 Subject: [PATCH 23/35] Debugging code for bad clip shapes. --- src/axom/mir/ElviraAlgorithm.hpp | 42 +++++++++++++++++++++++--------- 1 file changed, 30 insertions(+), 12 deletions(-) diff --git a/src/axom/mir/ElviraAlgorithm.hpp b/src/axom/mir/ElviraAlgorithm.hpp index ee8a277586..81e9e03cca 100644 --- a/src/axom/mir/ElviraAlgorithm.hpp +++ b/src/axom/mir/ElviraAlgorithm.hpp @@ -35,13 +35,13 @@ #include // Uncomment to save inputs and outputs. -//#define AXOM_ELVIRA_DEBUG +// #define AXOM_ELVIRA_DEBUG // Uncomment to debug make fragments. -//#define AXOM_ELVIRA_DEBUG_MAKE_FRAGMENTS +// #define AXOM_ELVIRA_DEBUG_MAKE_FRAGMENTS // Uncomment to gather ELVIRA data and save to YAML file. -//#define AXOM_ELVIRA_GATHER_INFO +// #define AXOM_ELVIRA_GATHER_INFO #if defined(AXOM_ELVIRA_DEBUG) #include @@ -918,6 +918,11 @@ class ElviraAlgorithm : public axom::mir::MIRAlgorithm // Get the starting shape. const auto inputShape = deviceShapeView.getShape(zoneIndex); +#if defined(AXOM_ELVIRA_DEBUG_MAKE_FRAGMENTS) && !defined(AXOM_DEVICE_CODE) + // Get the shape's bounding box and enlarge it a little. + auto inputShapeBBox = axom::primal::compute_bounding_box(inputShape); + inputShapeBBox.scale(1.05); +#endif // Get the zone's actual volume. const double zoneVol = utils::ComputeShapeAmount::execute(inputShape); @@ -1000,27 +1005,40 @@ class ElviraAlgorithm : public axom::mir::MIRAlgorithm // Emit clippedShape as material matId buildView.addShape(zoneIndex, fragmentIndex, clippedShape, matId, pt, planeOffset, normalPtr); +#if defined(AXOM_ELVIRA_DEBUG_MAKE_FRAGMENTS) && !defined(AXOM_DEVICE_CODE) + // Examine clippedShape's bounding box. It should NEVER be larger than the + // original inputShape's bounding box. If so, there was probably an error + // in clipping. + const auto clippedShapeBBox = axom::primal::compute_bounding_box(clippedShape); + if(!inputShapeBBox.contains(clippedShapeBBox)) + { + SLIC_ERROR("\tclip: BAD CLIPPED SHAPE IN ZONE " << zoneIndex + << "\n\t\tinputShape=" << inputShape + << "\n\t\tinputShapeBBox=" << inputShapeBBox + << "\n\t\tclippedShape=" << clippedShape + << "\n\t\tclippedShapeBBox=" << clippedShapeBBox + ); + } +#endif + // Clip in the other direction to get the remaining fragment for the next material. if(m == 0) { #if defined(AXOM_ELVIRA_DEBUG_MAKE_FRAGMENTS) && !defined(AXOM_DEVICE_CODE) - SLIC_DEBUG("\tclip: P=" << P << ", before=" << inputShape); -#endif - remaining = axom::primal::clip(inputShape, P); -#if defined(AXOM_ELVIRA_DEBUG_MAKE_FRAGMENTS) && !defined(AXOM_DEVICE_CODE) - SLIC_DEBUG("\tclip: after=" << clippedShape); + SLIC_DEBUG("\tclip: before=" << inputShape << ", P=" << P); #endif + remaining = axom::primal::clip(inputShape, P, detail::clip_precision::eps); } else { #if defined(AXOM_ELVIRA_DEBUG_MAKE_FRAGMENTS) && !defined(AXOM_DEVICE_CODE) - SLIC_DEBUG("\tclip: P=" << P << ", before=" << remaining); + SLIC_DEBUG("\tclip: before=" << remaining << ", P=" << P); #endif - remaining = axom::primal::clip(remaining, P); + remaining = axom::primal::clip(remaining, P, detail::clip_precision::eps); + } #if defined(AXOM_ELVIRA_DEBUG_MAKE_FRAGMENTS) && !defined(AXOM_DEVICE_CODE) - SLIC_DEBUG("\tclip: after=" << remaining); + SLIC_DEBUG("\tclip: after=" << remaining); #endif - } } // Emit the last leftover fragment. From b91d136e45787020a9bbeecf3ef9090e7d190b66 Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Fri, 7 Nov 2025 17:03:40 -0800 Subject: [PATCH 24/35] Logging changes in example applications --- src/axom/mir/examples/concentric_circles/MIRApplication.cpp | 1 + src/axom/mir/examples/heavily_mixed/HMApplication.cpp | 6 +++--- 2 files changed, 4 insertions(+), 3 deletions(-) diff --git a/src/axom/mir/examples/concentric_circles/MIRApplication.cpp b/src/axom/mir/examples/concentric_circles/MIRApplication.cpp index 9d33a486b7..94ff2dafd8 100644 --- a/src/axom/mir/examples/concentric_circles/MIRApplication.cpp +++ b/src/axom/mir/examples/concentric_circles/MIRApplication.cpp @@ -109,6 +109,7 @@ int MIRApplication::initialize(int argc, char **argv) int MIRApplication::execute() { axom::slic::SimpleLogger logger(axom::slic::message::Info); + axom::slic::setLoggingMsgLevel(axom::slic::message::Debug); if(handler) { diff --git a/src/axom/mir/examples/heavily_mixed/HMApplication.cpp b/src/axom/mir/examples/heavily_mixed/HMApplication.cpp index 404b08167a..2ac03460ea 100644 --- a/src/axom/mir/examples/heavily_mixed/HMApplication.cpp +++ b/src/axom/mir/examples/heavily_mixed/HMApplication.cpp @@ -100,9 +100,8 @@ void heavily_mixed_matset(const std::string &topoName, int dims[3], int refineme for(int i = 0; i < nmats; i++) { int matno = i + 1; - std::stringstream ss; - ss << "mat" << matno; - n_material_map[ss.str()] = matno; + const std::string name = axom::fmt::format("mat{:02d}", matno); + n_material_map[name] = matno; } n_matset["material_ids"].set(material_ids); n_matset["volume_fractions"].set(volume_fractions); @@ -312,6 +311,7 @@ int HMApplication::initialize(int argc, char **argv) int HMApplication::execute() { axom::slic::SimpleLogger logger(axom::slic::message::Info); + axom::slic::setLoggingMsgLevel(axom::slic::message::Debug); if(m_handler) { From 68d0ab5518d950530fc06c9caae681e98b4da635 Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Fri, 7 Nov 2025 17:04:03 -0800 Subject: [PATCH 25/35] Fix view type --- src/axom/mir/examples/heavily_mixed/runMIR.hpp | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/src/axom/mir/examples/heavily_mixed/runMIR.hpp b/src/axom/mir/examples/heavily_mixed/runMIR.hpp index adff6b363f..b87527ba5c 100644 --- a/src/axom/mir/examples/heavily_mixed/runMIR.hpp +++ b/src/axom/mir/examples/heavily_mixed/runMIR.hpp @@ -35,7 +35,6 @@ int runMIR(const conduit::Node &hostMesh, const conduit::Node &options, conduit: trials = std::max(1, options["trials"].to_int()); } - // TODO: set MAXMATERIALS to 1, remove warning. // Check materials. constexpr int MAXMATERIALS = 100; auto materialInfo = axom::bump::views::materials(hostMesh["matsets/mat"]); @@ -70,7 +69,7 @@ int runMIR(const conduit::Node &hostMesh, const conduit::Node &options, conduit: auto topologyView = make_rectilinear_topology::view(n_topology); using TopologyView = decltype(topologyView); - auto matsetView = make_unibuffer_matset::view(n_matset); + auto matsetView = make_unibuffer_matset::view(n_matset); using MatsetView = decltype(matsetView); if(method == "equiz") From ea443fd183fcdc90817f00f3f45cb3e08b5a7c7e Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Fri, 7 Nov 2025 17:04:58 -0800 Subject: [PATCH 26/35] Adjust clip precision based on type. --- src/axom/mir/detail/elvira_detail.hpp | 17 ++++++++++++++++- 1 file changed, 16 insertions(+), 1 deletion(-) diff --git a/src/axom/mir/detail/elvira_detail.hpp b/src/axom/mir/detail/elvira_detail.hpp index 41c5cf0307..2b3cba7fb6 100644 --- a/src/axom/mir/detail/elvira_detail.hpp +++ b/src/axom/mir/detail/elvira_detail.hpp @@ -72,6 +72,21 @@ class ELVIRAOptions : public axom::bump::Options namespace detail { +//------------------------------------------------------------------------------ +template +struct clip_precision +{ + // The default precision used in axom::primal::clip + static constexpr double eps = 1.e-10; +}; + +template <> +struct clip_precision +{ + // Make the axom::primal::clip more precise for double. + static constexpr double eps = 1.e-15; +}; + //------------------------------------------------------------------------------ /*! * \brief Compute a range given by 2 points given a shape and a normal. The @@ -162,7 +177,7 @@ AXOM_HOST_DEVICE inline ClipResultType clipToVolume(const ShapeType &shape, const auto P = axom::primal::Plane(clipNormal, pt, false); // Clip the shape at the current plane. - clippedShape = axom::primal::clip(shape, P); + clippedShape = axom::primal::clip(shape, P, clip_precision::eps); // Find the volume of the clipped shape. const double fragmentVolume = utils::ComputeShapeAmount::execute(clippedShape); From 9905b8adea7d9e830d6b1c2e8a049052bb72b468 Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Fri, 7 Nov 2025 17:05:53 -0800 Subject: [PATCH 27/35] make style --- src/axom/bump/tests/bump_views.cpp | 2 +- src/axom/bump/views/MaterialView.hpp | 16 +++- src/axom/mir/ElviraAlgorithm.hpp | 16 ++-- .../examples/heavily_mixed/HMApplication.cpp | 79 ++++++++----------- src/axom/mir/tests/mir_elvira3d.cpp | 3 +- 5 files changed, 58 insertions(+), 58 deletions(-) diff --git a/src/axom/bump/tests/bump_views.cpp b/src/axom/bump/tests/bump_views.cpp index d5aba7d9e3..07471d4587 100644 --- a/src/axom/bump/tests/bump_views.cpp +++ b/src/axom/bump/tests/bump_views.cpp @@ -720,7 +720,7 @@ struct test_braid2d_mat count++; } - resultsView[index] = (eq_count == count) ? 1 : 0; + resultsView[index] = (eq_count == count) ? 1 : 0; }); // Get containsView data to the host and compare results diff --git a/src/axom/bump/views/MaterialView.hpp b/src/axom/bump/views/MaterialView.hpp index cc07b7c440..a32d50ecaa 100644 --- a/src/axom/bump/views/MaterialView.hpp +++ b/src/axom/bump/views/MaterialView.hpp @@ -137,7 +137,9 @@ class UnibufferMaterialView } AXOM_HOST_DEVICE - axom::IndexType zoneMaterials(ZoneIndex zi, axom::ArrayView &ids, axom::ArrayView &vfs) const + axom::IndexType zoneMaterials(ZoneIndex zi, + axom::ArrayView &ids, + axom::ArrayView &vfs) const { SLIC_ASSERT(zi < static_cast(numberOfZones())); @@ -394,7 +396,9 @@ class MultiBufferMaterialView } AXOM_HOST_DEVICE - axom::IndexType zoneMaterials(ZoneIndex zi, axom::ArrayView &ids, axom::ArrayView &vfs) const + axom::IndexType zoneMaterials(ZoneIndex zi, + axom::ArrayView &ids, + axom::ArrayView &vfs) const { axom::IndexType n = 0; for(axom::IndexType i = 0; i < m_size; i++) @@ -677,7 +681,9 @@ class ElementDominantMaterialView } AXOM_HOST_DEVICE - axom::IndexType zoneMaterials(ZoneIndex zi, axom::ArrayView &ids, axom::ArrayView &vfs) const + axom::IndexType zoneMaterials(ZoneIndex zi, + axom::ArrayView &ids, + axom::ArrayView &vfs) const { axom::IndexType n = 0; for(axom::IndexType i = 0; i < m_volume_fractions.size(); i++) @@ -980,7 +986,9 @@ class MaterialDominantMaterialView } AXOM_HOST_DEVICE - axom::IndexType zoneMaterials(ZoneIndex zi, axom::ArrayView &ids, axom::ArrayView &vfs) const + axom::IndexType zoneMaterials(ZoneIndex zi, + axom::ArrayView &ids, + axom::ArrayView &vfs) const { axom::IndexType n = 0; for(axom::IndexType mi = 0; mi < m_size; mi++) diff --git a/src/axom/mir/ElviraAlgorithm.hpp b/src/axom/mir/ElviraAlgorithm.hpp index 81e9e03cca..a3898d5cd0 100644 --- a/src/axom/mir/ElviraAlgorithm.hpp +++ b/src/axom/mir/ElviraAlgorithm.hpp @@ -682,8 +682,10 @@ class ElviraAlgorithm : public axom::mir::MIRAlgorithm const auto offset = matOffsetView[szIndex]; // Determine the views for this zone's material data. - axom::ArrayView ids(sortedMaterialIdsView.data() + offset, matCount); - axom::ArrayView vfs(sortedMaterialVfsView.data() + offset, matCount); + axom::ArrayView ids(sortedMaterialIdsView.data() + offset, + matCount); + axom::ArrayView vfs(sortedMaterialVfsView.data() + offset, + matCount); // Get materials for this zone from the matset, directly into the "sorted" views. [[maybe_unused]] auto ids_size = deviceMatsetView.zoneMaterials(matZoneIndex, ids, vfs); @@ -1012,12 +1014,10 @@ class ElviraAlgorithm : public axom::mir::MIRAlgorithm const auto clippedShapeBBox = axom::primal::compute_bounding_box(clippedShape); if(!inputShapeBBox.contains(clippedShapeBBox)) { - SLIC_ERROR("\tclip: BAD CLIPPED SHAPE IN ZONE " << zoneIndex - << "\n\t\tinputShape=" << inputShape - << "\n\t\tinputShapeBBox=" << inputShapeBBox - << "\n\t\tclippedShape=" << clippedShape - << "\n\t\tclippedShapeBBox=" << clippedShapeBBox - ); + SLIC_ERROR("\tclip: BAD CLIPPED SHAPE IN ZONE " + << zoneIndex << "\n\t\tinputShape=" << inputShape << "\n\t\tinputShapeBBox=" + << inputShapeBBox << "\n\t\tclippedShape=" << clippedShape + << "\n\t\tclippedShapeBBox=" << clippedShapeBBox); } #endif diff --git a/src/axom/mir/examples/heavily_mixed/HMApplication.cpp b/src/axom/mir/examples/heavily_mixed/HMApplication.cpp index 2ac03460ea..3e6369f554 100644 --- a/src/axom/mir/examples/heavily_mixed/HMApplication.cpp +++ b/src/axom/mir/examples/heavily_mixed/HMApplication.cpp @@ -30,7 +30,12 @@ namespace detail * \param n_field The field used for matset creation. * \param nmats The number of materials to make. */ -void heavily_mixed_matset(const std::string &topoName, int dims[3], int refinement, conduit::Node &n_coarse, const conduit::Node &n_field, int nmats) +void heavily_mixed_matset(const std::string &topoName, + int dims[3], + int refinement, + conduit::Node &n_coarse, + const conduit::Node &n_field, + int nmats) { const auto fine = n_field.as_float64_accessor(); int nzones = dims[0] * dims[1] * dims[2]; @@ -38,7 +43,7 @@ void heavily_mixed_matset(const std::string &topoName, int dims[3], int refineme std::vector vfs(nslots, 0.); // break the data range into nmats parts. - const double matSize = 1000. / nmats; //fine.max() / nmats; + const double matSize = 1000. / nmats; //fine.max() / nmats; const int rdims[] = {dims[0] * refinement, dims[1] * refinement, dims[2] * refinement}; @@ -53,13 +58,14 @@ void heavily_mixed_matset(const std::string &topoName, int dims[3], int refineme const int jr = j * refinement; const int ir = i * refinement; for(int jj = 0; jj < refinement; jj++) - for(int ii = 0; ii < refinement; ii++) - { - const int fine_index = (kr * rdims[0] * rdims[1]) + ((jr + jj) * rdims[0]) + (ir + ii); - const int matid = axom::utilities::clampVal(static_cast(fine[fine_index] / matSize), 0, nmats-1); - const int matslot = zoneIndex * nmats + matid; - vfs[matslot] += 1. / (refinement * refinement); - } + for(int ii = 0; ii < refinement; ii++) + { + const int fine_index = (kr * rdims[0] * rdims[1]) + ((jr + jj) * rdims[0]) + (ir + ii); + const int matid = + axom::utilities::clampVal(static_cast(fine[fine_index] / matSize), 0, nmats - 1); + const int matslot = zoneIndex * nmats + matid; + vfs[matslot] += 1. / (refinement * refinement); + } } } } @@ -127,15 +133,7 @@ void heavily_mixed(conduit::Node &n_mesh, int dims[3], int refinement, int nmats const conduit::float64 c_re = -0.5125; const conduit::float64 c_im = 0.5213; - conduit::blueprint::mesh::examples::julia(dims[0], - dims[1], - x_min, - x_max, - y_min, - y_max, - c_re, - c_im, - n_mesh); + conduit::blueprint::mesh::examples::julia(dims[0], dims[1], x_min, x_max, y_min, y_max, c_re, c_im, n_mesh); if(dims[2] > 1) { // Add another dimension to the coordset. @@ -160,30 +158,23 @@ void heavily_mixed(conduit::Node &n_mesh, int dims[3], int refinement, int nmats conduit::Node n_field; n_field.set(conduit::DataType::int32(rdims[0] * rdims[1] * rdims[2])); conduit::int32 *destPtr = n_field.as_int32_ptr(); - axom::for_all(rdims[2], AXOM_LAMBDA(int k) - { - const auto t = static_cast(k) / (dims[2] - 1); - // Interpolate the window - const conduit::float64 x0 = axom::utilities::lerp(x_min, x1_min, t); - const conduit::float64 x1 = axom::utilities::lerp(x_max, x1_max, t); - const conduit::float64 y0 = axom::utilities::lerp(y_min, y1_min, t); - const conduit::float64 y1 = axom::utilities::lerp(y_max, y1_max, t); - conduit::Node n_rmesh; - conduit::blueprint::mesh::examples::julia(rdims[0], - rdims[1], - x0, - x1, - y0, - y1, - c_re, - c_im, - n_rmesh); - const conduit::Node &n_src_field = n_rmesh["fields/iters/values"]; - const conduit::int32 *srcPtr = n_src_field.as_int32_ptr(); - conduit::int32 *currentDestPtr = destPtr + k * rdims[0] * rdims[1]; - axom::copy(currentDestPtr, srcPtr, rdims[0] * rdims[1] * sizeof(conduit::int32)); - SLIC_INFO(axom::fmt::format("Made slice {}/{}", k+1, rdims[2])); - }); + axom::for_all( + rdims[2], + AXOM_LAMBDA(int k) { + const auto t = static_cast(k) / (dims[2] - 1); + // Interpolate the window + const conduit::float64 x0 = axom::utilities::lerp(x_min, x1_min, t); + const conduit::float64 x1 = axom::utilities::lerp(x_max, x1_max, t); + const conduit::float64 y0 = axom::utilities::lerp(y_min, y1_min, t); + const conduit::float64 y1 = axom::utilities::lerp(y_max, y1_max, t); + conduit::Node n_rmesh; + conduit::blueprint::mesh::examples::julia(rdims[0], rdims[1], x0, x1, y0, y1, c_re, c_im, n_rmesh); + const conduit::Node &n_src_field = n_rmesh["fields/iters/values"]; + const conduit::int32 *srcPtr = n_src_field.as_int32_ptr(); + conduit::int32 *currentDestPtr = destPtr + k * rdims[0] * rdims[1]; + axom::copy(currentDestPtr, srcPtr, rdims[0] * rdims[1] * sizeof(conduit::int32)); + SLIC_INFO(axom::fmt::format("Made slice {}/{}", k + 1, rdims[2])); + }); // Make a matset based on the higher resolution julia field. heavily_mixed_matset("topo", dims, refinement, n_mesh, n_field, nmats); @@ -208,12 +199,12 @@ void heavily_mixed(conduit::Node &n_mesh, int dims[3], int refinement, int nmats } } -} // end namespace detail +} // end namespace detail //-------------------------------------------------------------------------------- HMApplication::HMApplication() : m_handler(true) - , m_dims{100, 100, 1} + , m_dims {100, 100, 1} , m_numMaterials(40) , m_refinement(40) , m_numTrials(1) diff --git a/src/axom/mir/tests/mir_elvira3d.cpp b/src/axom/mir/tests/mir_elvira3d.cpp index 9e27336301..ec46a03774 100644 --- a/src/axom/mir/tests/mir_elvira3d.cpp +++ b/src/axom/mir/tests/mir_elvira3d.cpp @@ -309,7 +309,8 @@ struct test_Elvira3D SLIC_ASSERT(index >= 0 && index < sortedIdsView.size()); // Use an atomic to sum the value. - axom::atomicAdd(totalVolumeView.data() + index, zoneVolumes[zi] * zoneMat.volume_fraction()); + axom::atomicAdd(totalVolumeView.data() + index, + zoneVolumes[zi] * zoneMat.volume_fraction()); } }); From 4703bf184b286ee97b6ec12113e0560ea46d475f Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Fri, 7 Nov 2025 19:15:19 -0800 Subject: [PATCH 28/35] Added a check for the last fragment's bounding box --- src/axom/mir/ElviraAlgorithm.hpp | 18 ++++++++++++++++-- 1 file changed, 16 insertions(+), 2 deletions(-) diff --git a/src/axom/mir/ElviraAlgorithm.hpp b/src/axom/mir/ElviraAlgorithm.hpp index a3898d5cd0..ea08ce7ed1 100644 --- a/src/axom/mir/ElviraAlgorithm.hpp +++ b/src/axom/mir/ElviraAlgorithm.hpp @@ -1025,14 +1025,14 @@ class ElviraAlgorithm : public axom::mir::MIRAlgorithm if(m == 0) { #if defined(AXOM_ELVIRA_DEBUG_MAKE_FRAGMENTS) && !defined(AXOM_DEVICE_CODE) - SLIC_DEBUG("\tclip: before=" << inputShape << ", P=" << P); + SLIC_DEBUG("\tclip: before=" << inputShape << ", P=" << P << ", pt=" << pt); #endif remaining = axom::primal::clip(inputShape, P, detail::clip_precision::eps); } else { #if defined(AXOM_ELVIRA_DEBUG_MAKE_FRAGMENTS) && !defined(AXOM_DEVICE_CODE) - SLIC_DEBUG("\tclip: before=" << remaining << ", P=" << P); + SLIC_DEBUG("\tclip: before=" << remaining << ", P=" << P << ", pt=" << pt); #endif remaining = axom::primal::clip(remaining, P, detail::clip_precision::eps); } @@ -1053,6 +1053,20 @@ class ElviraAlgorithm : public axom::mir::MIRAlgorithm lastNormal[d] = -normal[d]; } buildView.addShape(zoneIndex, fragmentIndex, remaining, matId, pt, -planeOffset, lastNormal); + +#if defined(AXOM_ELVIRA_DEBUG_MAKE_FRAGMENTS) && !defined(AXOM_DEVICE_CODE) + // Examine remaining's bounding box. It should NEVER be larger than the + // original inputShape's bounding box. If so, there was probably an error + // in clipping. + const auto remainingBBox = axom::primal::compute_bounding_box(remaining); + if(!inputShapeBBox.contains(remainingBBox)) + { + SLIC_ERROR("\tclip: BAD CLIPPED SHAPE IN ZONE " + << zoneIndex << "\n\t\tinputShape=" << inputShape + << "\n\t\tinputShapeBBox=" << inputShapeBBox << "\n\t\tremaining=" << remaining + << "\n\t\tremainingBBox=" << remainingBBox); + } +#endif }); reportErrors(__LINE__); } From 0654350ed4aa431063c2e6ca8961eb181a161460 Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Wed, 12 Nov 2025 11:03:35 -0800 Subject: [PATCH 29/35] Renamed iterator to const_iterator and improved its docs. Also changed a small comparison in mir_concentric_circles. --- src/axom/bump/views/MaterialView.hpp | 130 ++++++++++-------- src/axom/mir/detail/elvira_detail.hpp | 16 ++- .../examples/concentric_circles/runMIR.hpp | 12 +- 3 files changed, 87 insertions(+), 71 deletions(-) diff --git a/src/axom/bump/views/MaterialView.hpp b/src/axom/bump/views/MaterialView.hpp index a32d50ecaa..5dbeb6b345 100644 --- a/src/axom/bump/views/MaterialView.hpp +++ b/src/axom/bump/views/MaterialView.hpp @@ -184,11 +184,13 @@ class UnibufferMaterialView } /*! - * \brief An iterator class for iterating over the mat/vf data in a zone. + * \brief An iterator class for iterating over read-only data in a zone. + * The iterator can access material ids and volume fractions for one + * material at a time in the associated zone. */ - class iterator + class const_iterator { - // Let the material view call the iterator constructor. + // Let the material view call the const_iterator constructor. friend class UnibufferMaterialView; public: @@ -209,24 +211,24 @@ class UnibufferMaterialView void AXOM_HOST_DEVICE operator++() { advance(true); } void AXOM_HOST_DEVICE operator++(int) { advance(true); } - bool AXOM_HOST_DEVICE operator==(const iterator &rhs) const + bool AXOM_HOST_DEVICE operator==(const const_iterator &rhs) const { return m_currentIndex == rhs.m_currentIndex && m_zoneIndex == rhs.m_zoneIndex && m_view == rhs.m_view; } - bool AXOM_HOST_DEVICE operator!=(const iterator &rhs) const + bool AXOM_HOST_DEVICE operator!=(const const_iterator &rhs) const { return m_currentIndex != rhs.m_currentIndex || m_zoneIndex != rhs.m_zoneIndex || m_view != rhs.m_view; } private: - DISABLE_DEFAULT_CTOR(iterator); + DISABLE_DEFAULT_CTOR(const_iterator); /// Constructor - AXOM_HOST_DEVICE iterator(const UnibufferMaterialView *view, - ZoneIndex zoneIndex, - axom::IndexType currentIndex = 0) + AXOM_HOST_DEVICE const_iterator(const UnibufferMaterialView *view, + ZoneIndex zoneIndex, + axom::IndexType currentIndex = 0) : m_view(view) , m_zoneIndex(zoneIndex) , m_currentIndex(currentIndex) @@ -248,8 +250,8 @@ class UnibufferMaterialView axom::IndexType m_currentIndex; axom::IndexType m_index; // not considered in ==, != }; - // Let the iterator access members. - friend class iterator; + // Let the const_iterator access members. + friend class const_iterator; /*! * \brief Return the iterator for the beginning of a zone's material data. @@ -258,11 +260,11 @@ class UnibufferMaterialView * * \return The iterator for the beginning of a zone's material data. */ - iterator AXOM_HOST_DEVICE beginZone(ZoneIndex zi) const + const_iterator AXOM_HOST_DEVICE beginZone(ZoneIndex zi) const { SLIC_ASSERT(zi < static_cast(numberOfZones())); - auto it = iterator(this, zi, 0); + auto it = const_iterator(this, zi, 0); it.advance(false); return it; } @@ -274,11 +276,11 @@ class UnibufferMaterialView * * \return The iterator for the end of a zone's material data. */ - iterator AXOM_HOST_DEVICE endZone(ZoneIndex zi) const + const_iterator AXOM_HOST_DEVICE endZone(ZoneIndex zi) const { SLIC_ASSERT(zi < static_cast(numberOfZones())); - return iterator(this, zi, m_sizes[zi]); + return const_iterator(this, zi, m_sizes[zi]); } private: @@ -448,11 +450,13 @@ class MultiBufferMaterialView } /*! - * \brief An iterator class for iterating over the mat/vf data in a zone. + * \brief An iterator class for iterating over read-only data in a zone. + * The iterator can access material ids and volume fractions for one + * material at a time in the associated zone. */ - class iterator + class const_iterator { - // Let the material view call the iterator constructor. + // Let the material view call the const_iterator constructor. friend class MultiBufferMaterialView; public: @@ -483,24 +487,24 @@ class MultiBufferMaterialView m_currentIndex += (m_currentIndex < m_view->m_size) ? 1 : 0; advance(); } - bool AXOM_HOST_DEVICE operator==(const iterator &rhs) const + bool AXOM_HOST_DEVICE operator==(const const_iterator &rhs) const { return m_currentIndex == rhs.m_currentIndex && m_zoneIndex == rhs.m_zoneIndex && m_view == rhs.m_view; } - bool AXOM_HOST_DEVICE operator!=(const iterator &rhs) const + bool AXOM_HOST_DEVICE operator!=(const const_iterator &rhs) const { return m_currentIndex != rhs.m_currentIndex || m_zoneIndex != rhs.m_zoneIndex || m_view != rhs.m_view; } private: - DISABLE_DEFAULT_CTOR(iterator); + DISABLE_DEFAULT_CTOR(const_iterator); /// Constructor - AXOM_HOST_DEVICE iterator(const MultiBufferMaterialView *view, - ZoneIndex zoneIndex, - axom::IndexType currentIndex = 0) + AXOM_HOST_DEVICE const_iterator(const MultiBufferMaterialView *view, + ZoneIndex zoneIndex, + axom::IndexType currentIndex = 0) : m_view(view) , m_zoneIndex(zoneIndex) , m_currentIndex(currentIndex) @@ -530,8 +534,8 @@ class MultiBufferMaterialView ZoneIndex m_zoneIndex; axom::IndexType m_currentIndex; }; - // Let the iterator access members. - friend class iterator; + // Let the const_iterator access members. + friend class const_iterator; /*! * \brief Return the iterator for the beginning of a zone's material data. @@ -540,11 +544,11 @@ class MultiBufferMaterialView * * \return The iterator for the beginning of a zone's material data. */ - iterator AXOM_HOST_DEVICE beginZone(ZoneIndex zi) const + const_iterator AXOM_HOST_DEVICE beginZone(ZoneIndex zi) const { SLIC_ASSERT(zi < static_cast(numberOfZones())); - auto it = iterator(this, zi, 0); + auto it = const_iterator(this, zi, 0); it.advance(); return it; } @@ -556,10 +560,10 @@ class MultiBufferMaterialView * * \return The iterator for the end of a zone's material data. */ - iterator AXOM_HOST_DEVICE endZone(ZoneIndex zi) const + const_iterator AXOM_HOST_DEVICE endZone(ZoneIndex zi) const { SLIC_ASSERT(zi < static_cast(numberOfZones())); - return iterator(this, zi, m_size); + return const_iterator(this, zi, m_size); } private: @@ -724,11 +728,13 @@ class ElementDominantMaterialView } /*! - * \brief An iterator class for iterating over the mat/vf data in a zone. + * \brief An iterator class for iterating over read-only data in a zone. + * The iterator can access material ids and volume fractions for one + * material at a time in the associated zone. */ - class iterator + class const_iterator { - // Let the material view call the iterator constructor. + // Let the material view call the const_iterator constructor. friend class ElementDominantMaterialView; public: @@ -756,24 +762,24 @@ class ElementDominantMaterialView m_currentIndex += (m_currentIndex < m_view->m_volume_fractions.size()) ? 1 : 0; advance(); } - bool AXOM_HOST_DEVICE operator==(const iterator &rhs) const + bool AXOM_HOST_DEVICE operator==(const const_iterator &rhs) const { return m_currentIndex == rhs.m_currentIndex && m_zoneIndex == rhs.m_zoneIndex && m_view == rhs.m_view; } - bool AXOM_HOST_DEVICE operator!=(const iterator &rhs) const + bool AXOM_HOST_DEVICE operator!=(const const_iterator &rhs) const { return m_currentIndex != rhs.m_currentIndex || m_zoneIndex != rhs.m_zoneIndex || m_view != rhs.m_view; } private: - DISABLE_DEFAULT_CTOR(iterator); + DISABLE_DEFAULT_CTOR(const_iterator); /// Constructor - AXOM_HOST_DEVICE iterator(const ElementDominantMaterialView *view, - ZoneIndex zoneIndex, - axom::IndexType currentIndex = 0) + AXOM_HOST_DEVICE const_iterator(const ElementDominantMaterialView *view, + ZoneIndex zoneIndex, + axom::IndexType currentIndex = 0) : m_view(view) , m_zoneIndex(zoneIndex) , m_currentIndex(currentIndex) @@ -796,8 +802,8 @@ class ElementDominantMaterialView ZoneIndex m_zoneIndex; axom::IndexType m_currentIndex; }; - // Let the iterator access members. - friend class iterator; + // Let the const_iterator access members. + friend class const_iterator; /*! * \brief Return the iterator for the beginning of a zone's material data. @@ -806,11 +812,11 @@ class ElementDominantMaterialView * * \return The iterator for the beginning of a zone's material data. */ - iterator AXOM_HOST_DEVICE beginZone(ZoneIndex zi) const + const_iterator AXOM_HOST_DEVICE beginZone(ZoneIndex zi) const { SLIC_ASSERT(zi < static_cast(numberOfZones())); - auto it = iterator(this, zi, 0); + auto it = const_iterator(this, zi, 0); it.advance(); return it; } @@ -822,10 +828,10 @@ class ElementDominantMaterialView * * \return The iterator for the end of a zone's material data. */ - iterator AXOM_HOST_DEVICE endZone(ZoneIndex zi) const + const_iterator AXOM_HOST_DEVICE endZone(ZoneIndex zi) const { SLIC_ASSERT(zi < static_cast(numberOfZones())); - return iterator(this, zi, m_volume_fractions.size()); + return const_iterator(this, zi, m_volume_fractions.size()); } private: @@ -1042,11 +1048,13 @@ class MaterialDominantMaterialView } /*! - * \brief An iterator class for iterating over the mat/vf data in a zone. + * \brief An iterator class for iterating over read-only data in a zone. + * The iterator can access material ids and volume fractions for one + * material at a time in the associated zone. */ - class iterator + class const_iterator { - // Let the material view call the iterator constructor. + // Let the material view call the const_iterator constructor. friend class MaterialDominantMaterialView; public: @@ -1066,25 +1074,25 @@ class MaterialDominantMaterialView axom::IndexType AXOM_HOST_DEVICE size() const { return m_view->numberOfMaterials(m_zoneIndex); } void AXOM_HOST_DEVICE operator++() { advance(true); } void AXOM_HOST_DEVICE operator++(int) { advance(true); } - bool AXOM_HOST_DEVICE operator==(const iterator &rhs) const + bool AXOM_HOST_DEVICE operator==(const const_iterator &rhs) const { return m_miIndex == rhs.m_miIndex && m_index == rhs.m_index && m_zoneIndex == rhs.m_zoneIndex && m_view == rhs.m_view; } - bool AXOM_HOST_DEVICE operator!=(const iterator &rhs) const + bool AXOM_HOST_DEVICE operator!=(const const_iterator &rhs) const { return m_miIndex != rhs.m_miIndex || m_index != rhs.m_index || m_zoneIndex != rhs.m_zoneIndex || m_view != rhs.m_view; } private: - DISABLE_DEFAULT_CTOR(iterator); + DISABLE_DEFAULT_CTOR(const_iterator); /// Constructor - AXOM_HOST_DEVICE iterator(const MaterialDominantMaterialView *view, - ZoneIndex zoneIndex, - axom::IndexType miIndex, - axom::IndexType index) + AXOM_HOST_DEVICE const_iterator(const MaterialDominantMaterialView *view, + ZoneIndex zoneIndex, + axom::IndexType miIndex, + axom::IndexType index) : m_view(view) , m_zoneIndex(zoneIndex) , m_miIndex(miIndex) @@ -1128,8 +1136,8 @@ class MaterialDominantMaterialView axom::IndexType m_miIndex; axom::IndexType m_index; }; - // Let the iterator access members. - friend class iterator; + // Let the const_iterator access members. + friend class const_iterator; /*! * \brief Return the iterator for the beginning of a zone's material data. @@ -1138,11 +1146,11 @@ class MaterialDominantMaterialView * * \return The iterator for the beginning of a zone's material data. */ - iterator AXOM_HOST_DEVICE beginZone(ZoneIndex zi) const + const_iterator AXOM_HOST_DEVICE beginZone(ZoneIndex zi) const { SLIC_ASSERT(zi < static_cast(numberOfZones())); - auto it = iterator(this, zi, 0, 0); + auto it = const_iterator(this, zi, 0, 0); it.advance(false); return it; } @@ -1154,12 +1162,12 @@ class MaterialDominantMaterialView * * \return The iterator for the end of a zone's material data. */ - iterator AXOM_HOST_DEVICE endZone(ZoneIndex zi) const + const_iterator AXOM_HOST_DEVICE endZone(ZoneIndex zi) const { SLIC_ASSERT(zi < static_cast(numberOfZones())); const axom::IndexType miIndex = m_size; const axom::IndexType index = (m_size > 0) ? (m_volume_fractions[m_size - 1].size()) : 0; - return iterator(this, zi, miIndex, index); + return const_iterator(this, zi, miIndex, index); } private: diff --git a/src/axom/mir/detail/elvira_detail.hpp b/src/axom/mir/detail/elvira_detail.hpp index 2b3cba7fb6..205a3b92c6 100644 --- a/src/axom/mir/detail/elvira_detail.hpp +++ b/src/axom/mir/detail/elvira_detail.hpp @@ -104,14 +104,18 @@ AXOM_HOST_DEVICE inline void computeRange(const ShapeType &shape, const axom::primal::Vector &normal, axom::primal::Point range[2]) { +#if 1 + // Use the vertex mean as the center of the shape + const auto center = shape.vertexMean(); +#else // Compute the shape bounding box. const auto bbox = axom::primal::compute_bounding_box(shape); - - const axom::primal::Plane P(normal, bbox.getCentroid(), false); + const auto center = bbox.getCentroid(); +#endif + const axom::primal::Plane P(normal, center, false); + range[0] = range[1] = center; // Compute distances from all points in shape to plane. - const auto centroid = bbox.getCentroid(); - range[0] = range[1] = centroid; double dist[2] = {0., 0.}; for(axom::IndexType ip = 0; ip < shape.numVertices(); ip++) { @@ -119,12 +123,12 @@ AXOM_HOST_DEVICE inline void computeRange(const ShapeType &shape, if(d < dist[0]) { dist[0] = d; - range[0] = centroid + (normal * d); + range[0] = center + (normal * d); } if(d > dist[1]) { dist[1] = d; - range[1] = centroid + (normal * d); + range[1] = center + (normal * d); } } } diff --git a/src/axom/mir/examples/concentric_circles/runMIR.hpp b/src/axom/mir/examples/concentric_circles/runMIR.hpp index ed81e5a06c..59406a32a9 100644 --- a/src/axom/mir/examples/concentric_circles/runMIR.hpp +++ b/src/axom/mir/examples/concentric_circles/runMIR.hpp @@ -11,7 +11,7 @@ #include "axom/mir.hpp" // for Mir classes & functions template -void test_matset_traveral(MatsetView matsetView) +void test_matset_traversal(MatsetView matsetView) { AXOM_ANNOTATE_SCOPE("test_matset_traversal"); double vf1, vf2; @@ -50,8 +50,12 @@ void test_matset_traveral(MatsetView matsetView) vf2 = vfSum.get(); } - std::cout << "test_matset_traversal: vf1=" << vf1 << ", vf2=" << vf2 - << ", nzones=" << matsetView.numberOfZones() << std::endl; + const double eps = 1.e-10; + SLIC_INFO(axom::fmt::format("test_matset_traversal: vf1={}, vf2={}, nzones={}, result={}", + vf1, + vf2, + matsetView.numberOfZones(), + (axom::utilities::abs(vf1 - vf2) < eps) ? "pass" : "fail")); } template @@ -146,7 +150,7 @@ int runMIR(const conduit::Node &hostMesh, const conduit::Node &options, conduit: { using namespace axom::bump::views; auto matsetView = make_unibuffer_matset::view(n_matset); - test_matset_traveral(matsetView); + test_matset_traversal(matsetView); } else { From b49ffc170bc67100286d96ec9eef85f97b02c73f Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Wed, 12 Nov 2025 11:59:54 -0800 Subject: [PATCH 30/35] Tweak to the RELEASE-NOTES. --- RELEASE-NOTES.md | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/RELEASE-NOTES.md b/RELEASE-NOTES.md index 28cfceea95..45605ba8b5 100644 --- a/RELEASE-NOTES.md +++ b/RELEASE-NOTES.md @@ -33,9 +33,9 @@ The Axom project release numbers follow [Semantic Versioning](http://semver.org/ - The maximum number of vertices allowed in polygon primitives can now be passed as a template argument to `axom::bump::TopologyMapper`, `axom::bump::PrimalAdaptor`, and `axom::mir::ElviraAlgorithm`. -- Material views in `axom::bump::views` were enhanced with `iterator` classes that +- Material views in `axom::bump::views` were enhanced with `const_iterator` classes that enable traversal of material data for zones so kernels do not need to use large fixed size - buffers inside kernels. + buffers to gather that data inside kernels. - Material views in `axom::bump::views` were enhanced with an overloaded `zoneMaterials()` method that allows data to be gathered into `axom::ArrayView` objects. - A new `heavily_mixed` example program was added in `axom::mir` to demonstrate running MIR on From 69daa62fcfc215a29090ad4c3d3ef67239b8b1de Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Wed, 12 Nov 2025 16:03:50 -0800 Subject: [PATCH 31/35] Revert a change I did not mean to commit. --- src/axom/mir/detail/elvira_detail.hpp | 16 ++++++---------- 1 file changed, 6 insertions(+), 10 deletions(-) diff --git a/src/axom/mir/detail/elvira_detail.hpp b/src/axom/mir/detail/elvira_detail.hpp index 205a3b92c6..2b3cba7fb6 100644 --- a/src/axom/mir/detail/elvira_detail.hpp +++ b/src/axom/mir/detail/elvira_detail.hpp @@ -104,18 +104,14 @@ AXOM_HOST_DEVICE inline void computeRange(const ShapeType &shape, const axom::primal::Vector &normal, axom::primal::Point range[2]) { -#if 1 - // Use the vertex mean as the center of the shape - const auto center = shape.vertexMean(); -#else // Compute the shape bounding box. const auto bbox = axom::primal::compute_bounding_box(shape); - const auto center = bbox.getCentroid(); -#endif - const axom::primal::Plane P(normal, center, false); - range[0] = range[1] = center; + + const axom::primal::Plane P(normal, bbox.getCentroid(), false); // Compute distances from all points in shape to plane. + const auto centroid = bbox.getCentroid(); + range[0] = range[1] = centroid; double dist[2] = {0., 0.}; for(axom::IndexType ip = 0; ip < shape.numVertices(); ip++) { @@ -123,12 +119,12 @@ AXOM_HOST_DEVICE inline void computeRange(const ShapeType &shape, if(d < dist[0]) { dist[0] = d; - range[0] = center + (normal * d); + range[0] = centroid + (normal * d); } if(d > dist[1]) { dist[1] = d; - range[1] = center + (normal * d); + range[1] = centroid + (normal * d); } } } From d149b90b9fc4fce3cb405a38e64a5293095d10ea Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Mon, 17 Nov 2025 14:07:08 -0800 Subject: [PATCH 32/35] Added default MAXMATERIALS value in matset views. --- src/axom/bump/views/MaterialView.hpp | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/src/axom/bump/views/MaterialView.hpp b/src/axom/bump/views/MaterialView.hpp index 5dbeb6b345..5a58658f11 100644 --- a/src/axom/bump/views/MaterialView.hpp +++ b/src/axom/bump/views/MaterialView.hpp @@ -75,7 +75,7 @@ MaterialInformation materials(const conduit::Node &matset); \endverbatim */ -template +template class UnibufferMaterialView { public: @@ -316,7 +316,7 @@ class UnibufferMaterialView \endverbatim */ -template +template class MultiBufferMaterialView { public: @@ -611,7 +611,7 @@ class MultiBufferMaterialView \endverbatim */ -template +template class ElementDominantMaterialView { public: @@ -897,7 +897,7 @@ class ElementDominantMaterialView \note This matset type does not seem so GPU friendly since there is some work to do for some of the queries. */ -template +template class MaterialDominantMaterialView { public: From 43d28b5c74be22976b132c05c8ce07960713fbc9 Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Tue, 18 Nov 2025 12:49:46 -0800 Subject: [PATCH 33/35] Fix OMP race condition with Umpire. --- src/axom/mir/examples/heavily_mixed/HMApplication.cpp | 8 +++++--- 1 file changed, 5 insertions(+), 3 deletions(-) diff --git a/src/axom/mir/examples/heavily_mixed/HMApplication.cpp b/src/axom/mir/examples/heavily_mixed/HMApplication.cpp index 3e6369f554..b0a14eecfd 100644 --- a/src/axom/mir/examples/heavily_mixed/HMApplication.cpp +++ b/src/axom/mir/examples/heavily_mixed/HMApplication.cpp @@ -338,9 +338,11 @@ int HMApplication::runMIR() { AXOM_ANNOTATE_SCOPE("generate"); int dims[3]; - dims[0] = m_dims[0]; - dims[1] = m_dims[1]; - dims[2] = m_dims[2]; + // NOTE: Use axom::copy to copy m_dims into dims. This way, the Umpire resource + // manager gets created now so by the time we have multiple threads using + // axom::copy, there is no race condition. + axom::copy(dims, m_dims, sizeof(int) * 3); + SLIC_INFO(axom::fmt::format("dims: {},{},{}", dims[0], dims[1], dims[2])); SLIC_INFO(axom::fmt::format("refinement: {}", m_refinement)); SLIC_INFO(axom::fmt::format("numMaterials: {}", m_numMaterials)); From 7e584bf9eec97e18549c9187aead685aa7fea728 Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Tue, 9 Dec 2025 16:30:59 -0800 Subject: [PATCH 34/35] Small changes and asserts. --- src/axom/bump/TopologyMapper.hpp | 7 +++++++ src/axom/bump/tests/bump_topology_mapper.cpp | 20 ++++++++++++------- src/axom/bump/views/Shapes.hpp | 14 +++++++++++++ .../UnstructuredTopologyMixedShapeView.hpp | 2 +- .../views/dispatch_unstructured_topology.hpp | 1 + 5 files changed, 36 insertions(+), 8 deletions(-) diff --git a/src/axom/bump/TopologyMapper.hpp b/src/axom/bump/TopologyMapper.hpp index b1342766be..df2f331b6d 100644 --- a/src/axom/bump/TopologyMapper.hpp +++ b/src/axom/bump/TopologyMapper.hpp @@ -623,8 +623,12 @@ class TopologyMapper axom::for_all( targetSelectionView.size(), AXOM_LAMBDA(axom::IndexType index) { + SLIC_ASSERT(index >= 0 && index < targetSelectionView.size()); + // Get the target zone as a primal shape. const axom::IndexType zi = targetSelectionView[index]; + SLIC_ASSERT(zi >= 0 && zi < targetView.numberOfZones()); + const auto targetBBox = targetView.getBoundingBox(zi); const auto targetShape = targetView.getShape(zi); #if defined(AXOM_DEBUG_TOPOLOGY_MAPPER) && !defined(AXOM_DEVICE_CODE) @@ -638,7 +642,10 @@ class TopologyMapper // Handle intersection in-depth of the bounding boxes intersected. auto handleIntersection = [&](std::int32_t currentNode, const std::int32_t *leafNodes) { const auto srcBboxIndex = leafNodes[currentNode]; + SLIC_ASSERT(srcBboxIndex >= 0 && srcBboxIndex < srcSelectionView.size()); + const auto srcZone = srcSelectionView[srcBboxIndex]; + SLIC_ASSERT(srcZone >= 0 && srcZone < srcView.numberOfZones()); #if defined(AXOM_DEBUG_TOPOLOGY_MAPPER) && !defined(AXOM_DEVICE_CODE) std::cout << "handleIntersection: targetZone=" << zi << ", srcZone=" << srcZone << std::endl; diff --git a/src/axom/bump/tests/bump_topology_mapper.cpp b/src/axom/bump/tests/bump_topology_mapper.cpp index 3fe205eade..c50068d6a4 100644 --- a/src/axom/bump/tests/bump_topology_mapper.cpp +++ b/src/axom/bump/tests/bump_topology_mapper.cpp @@ -461,18 +461,24 @@ class test_TopologyMapper auto srcCoordset = views::make_explicit_coordset::view(n_dev["coordsets/epm_coords"]); using SrcCoordsetView = decltype(srcCoordset); - using SrcTopologyView = views::UnstructuredTopologyMixedShapeView; axom::Array shapeValues, shapeIds; const conduit::Node &n_srcTopo = n_dev["topologies/epm"]; + EXPECT_EQ(n_srcTopo["type"].as_string(), "unstructured"); + EXPECT_EQ(n_srcTopo["elements/shape"].as_string(), "mixed"); + const int allocatorID = axom::execution_space::allocatorID(); auto shapeMap = views::buildShapeMap(n_srcTopo, shapeValues, shapeIds, allocatorID); - SrcTopologyView srcTopo( - utils::make_array_view(n_srcTopo["elements/connectivity"]), - utils::make_array_view(n_srcTopo["elements/shapes"]), - utils::make_array_view(n_srcTopo["elements/sizes"]), - utils::make_array_view(n_srcTopo["elements/offsets"]), - shapeMap); + const auto srcConnView = utils::make_array_view(n_srcTopo["elements/connectivity"]); + const auto srcShapesView = utils::make_array_view(n_srcTopo["elements/shapes"]); + const auto srcSizesView = utils::make_array_view(n_srcTopo["elements/sizes"]); + const auto srcOffsetsView = utils::make_array_view(n_srcTopo["elements/offsets"]); + // Check sizes with the current version of this test. + EXPECT_TRUE(srcSizesView.size() == 54); + EXPECT_TRUE(srcSizesView.size() == srcOffsetsView.size()); + EXPECT_TRUE(srcSizesView.size() == srcShapesView.size()); + + SrcTopologyView srcTopo(srcConnView, srcShapesView, srcSizesView, srcOffsetsView, shapeMap); const conduit::Node &n_srcMatset = n_dev["matsets/epm_matset"]; auto srcMatset = views::make_unibuffer_matset::view(n_srcMatset); diff --git a/src/axom/bump/views/Shapes.hpp b/src/axom/bump/views/Shapes.hpp index fc4c6010e5..1d650aabb5 100644 --- a/src/axom/bump/views/Shapes.hpp +++ b/src/axom/bump/views/Shapes.hpp @@ -1146,6 +1146,20 @@ inline int shapeNameToID(const std::string &name) return id; } +/*! + * \brief Determine whether a value is a valid ShapeID. + * + * \param shapeID The value to test for validity. + * + * \return True if the value is a valid ShapeID; False otherwise. + */ +template +AXOM_HOST_DEVICE +constexpr bool isValidShapeID(T shapeID) +{ + return shapeID >= static_cast(Point_ShapeID) && shapeID <= static_cast(Mixed_ShapeID); +} + } // end namespace views } // end namespace bump } // end namespace axom diff --git a/src/axom/bump/views/UnstructuredTopologyMixedShapeView.hpp b/src/axom/bump/views/UnstructuredTopologyMixedShapeView.hpp index c36df98bbb..42ae367dcd 100644 --- a/src/axom/bump/views/UnstructuredTopologyMixedShapeView.hpp +++ b/src/axom/bump/views/UnstructuredTopologyMixedShapeView.hpp @@ -181,7 +181,7 @@ class UnstructuredTopologyMixedShapeView const ConnectivityView shapeData(m_connectivity.data() + m_offsets[zoneIndex], m_sizes[zoneIndex]); const auto shapeID = m_shapeMap[m_shapes[zoneIndex]]; - SLIC_ASSERT(shapeID >= Point_ShapeID && shapeID <= Mixed_ShapeID); + SLIC_ASSERT(isValidShapeID(shapeID)); return ShapeType(shapeID, shapeData); } diff --git a/src/axom/bump/views/dispatch_unstructured_topology.hpp b/src/axom/bump/views/dispatch_unstructured_topology.hpp index 5e343426a9..84c37478a2 100644 --- a/src/axom/bump/views/dispatch_unstructured_topology.hpp +++ b/src/axom/bump/views/dispatch_unstructured_topology.hpp @@ -63,6 +63,7 @@ struct make_unstructured_single_shape_topology { sizesView = utils::make_array_view(n_topo["elements/sizes"]); offsetsView = utils::make_array_view(n_topo["elements/offsets"]); + SLIC_ASSERT(sizesView.size() == offsetsView.size()); } else { From 18f864ca07e419c741e19b189339795d998df295 Mon Sep 17 00:00:00 2001 From: Brad Whitlock Date: Tue, 9 Dec 2025 16:58:28 -0800 Subject: [PATCH 35/35] make style --- src/axom/bump/tests/bump_topology_mapper.cpp | 9 ++++++--- src/axom/bump/views/Shapes.hpp | 3 +-- 2 files changed, 7 insertions(+), 5 deletions(-) diff --git a/src/axom/bump/tests/bump_topology_mapper.cpp b/src/axom/bump/tests/bump_topology_mapper.cpp index c50068d6a4..c49c80a1c3 100644 --- a/src/axom/bump/tests/bump_topology_mapper.cpp +++ b/src/axom/bump/tests/bump_topology_mapper.cpp @@ -469,10 +469,13 @@ class test_TopologyMapper const int allocatorID = axom::execution_space::allocatorID(); auto shapeMap = views::buildShapeMap(n_srcTopo, shapeValues, shapeIds, allocatorID); - const auto srcConnView = utils::make_array_view(n_srcTopo["elements/connectivity"]); - const auto srcShapesView = utils::make_array_view(n_srcTopo["elements/shapes"]); + const auto srcConnView = + utils::make_array_view(n_srcTopo["elements/connectivity"]); + const auto srcShapesView = + utils::make_array_view(n_srcTopo["elements/shapes"]); const auto srcSizesView = utils::make_array_view(n_srcTopo["elements/sizes"]); - const auto srcOffsetsView = utils::make_array_view(n_srcTopo["elements/offsets"]); + const auto srcOffsetsView = + utils::make_array_view(n_srcTopo["elements/offsets"]); // Check sizes with the current version of this test. EXPECT_TRUE(srcSizesView.size() == 54); EXPECT_TRUE(srcSizesView.size() == srcOffsetsView.size()); diff --git a/src/axom/bump/views/Shapes.hpp b/src/axom/bump/views/Shapes.hpp index 1d650aabb5..ee8b172680 100644 --- a/src/axom/bump/views/Shapes.hpp +++ b/src/axom/bump/views/Shapes.hpp @@ -1154,8 +1154,7 @@ inline int shapeNameToID(const std::string &name) * \return True if the value is a valid ShapeID; False otherwise. */ template -AXOM_HOST_DEVICE -constexpr bool isValidShapeID(T shapeID) +AXOM_HOST_DEVICE constexpr bool isValidShapeID(T shapeID) { return shapeID >= static_cast(Point_ShapeID) && shapeID <= static_cast(Mixed_ShapeID); }