diff --git a/RELEASE-NOTES.md b/RELEASE-NOTES.md index b02804ac2e..e63620b1df 100644 --- a/RELEASE-NOTES.md +++ b/RELEASE-NOTES.md @@ -60,6 +60,9 @@ The Axom project release numbers follow [Semantic Versioning](http://semver.org/ - Updates to [MFEM version 4.8.0][https://github.com/mfem/mfem/releases/tag/v4.8] - Readers in Quest were moved from a `quest/readers` directory to `quest/io`. - Sina: Renames a Fortran module to `sina_hdf5_config` (from `hdf5_config`) +- Spin: Uses `axom::FlatMap` in `SparseOctreeLevel` implementation. We have observed a performance regression + of about 20% during InOutOctree construction and queries over STL surface meshes relative to the previous sparsehash + implementation. Please reach out to Axom developers if this affects you while we work on fixes for these. ### Fixed - Core: prevent incorrect instantiations of `axom::Array` from a host-only compile, when Axom is compiled diff --git a/src/axom/core/Macros.hpp b/src/axom/core/Macros.hpp index fbf81566a6..c211846d5e 100644 --- a/src/axom/core/Macros.hpp +++ b/src/axom/core/Macros.hpp @@ -301,35 +301,29 @@ * \note Can be used in test files after including `gtest/gtest.h` */ // Specifically, we expose the `TestBody()` method as public instead of private -#define AXOM_TYPED_TEST(CaseName, TestName) \ - static_assert(sizeof(GTEST_STRINGIFY_(TestName)) > 1, \ - "test-name must not be empty"); \ - template \ - class GTEST_TEST_CLASS_NAME_(CaseName, TestName) \ - : public CaseName { \ - public: \ - typedef CaseName TestFixture; \ - typedef gtest_TypeParam_ TypeParam; \ - void TestBody() override; \ - }; \ - GTEST_INTERNAL_ATTRIBUTE_MAYBE_UNUSED static bool \ - gtest_##CaseName##_##TestName##_registered_ = \ - ::testing::internal::TypeParameterizedTest< \ - CaseName, \ - ::testing::internal::TemplateSel, \ - GTEST_TYPE_PARAMS_( \ - CaseName)>::Register("", \ - ::testing::internal::CodeLocation( \ - __FILE__, __LINE__), \ - GTEST_STRINGIFY_(CaseName), \ - GTEST_STRINGIFY_(TestName), 0, \ - ::testing::internal::GenerateNames< \ - GTEST_NAME_GENERATOR_(CaseName), \ - GTEST_TYPE_PARAMS_(CaseName)>()); \ - template \ - void GTEST_TEST_CLASS_NAME_(CaseName, \ - TestName)::TestBody() - +#define AXOM_TYPED_TEST(CaseName, TestName) \ + static_assert(sizeof(GTEST_STRINGIFY_(TestName)) > 1, "test-name must not be empty"); \ + template \ + class GTEST_TEST_CLASS_NAME_(CaseName, TestName) : public CaseName \ + { \ + public: \ + typedef CaseName TestFixture; \ + typedef gtest_TypeParam_ TypeParam; \ + void TestBody() override; \ + }; \ + GTEST_INTERNAL_ATTRIBUTE_MAYBE_UNUSED static bool gtest_##CaseName##_##TestName##_registered_ = \ + ::testing::internal::TypeParameterizedTest< \ + CaseName, \ + ::testing::internal::TemplateSel, \ + GTEST_TYPE_PARAMS_( \ + CaseName)>::Register("", \ + ::testing::internal::CodeLocation(__FILE__, __LINE__), \ + GTEST_STRINGIFY_(CaseName), \ + GTEST_STRINGIFY_(TestName), \ + 0, \ + ::testing::internal::GenerateNames()); \ + template \ + void GTEST_TEST_CLASS_NAME_(CaseName, TestName)::TestBody() #endif // AXOM_MACROS_HPP_ diff --git a/src/axom/quest/InOutOctree.hpp b/src/axom/quest/InOutOctree.hpp index 901de99eab..06f33d9125 100644 --- a/src/axom/quest/InOutOctree.hpp +++ b/src/axom/quest/InOutOctree.hpp @@ -30,7 +30,6 @@ #include #include #include -#include #ifndef DUMP_VTK_MESH // #define DUMP_VTK_MESH @@ -177,9 +176,7 @@ class InOutOctree : public spin::SpatialOctree setVertexWeldThreshold(DEFAULT_VERTEX_WELD_THRESHOLD); } - /** - * \brief Generate the spatial index over the surface mesh - */ + /// \brief Generate the spatial index over the surface mesh void generateIndex(); /** @@ -229,9 +226,7 @@ class InOutOctree : public spin::SpatialOctree */ void insertVertex(VertexIndex idx, int startingLevel = 0); - /** - * \brief Insert all mesh cells into the octree, generating a PM octree - */ + /// \brief Insert all mesh cells into the octree, generating a PM octree void insertMeshCells(); /** @@ -448,6 +443,7 @@ void InOutOctree::generateIndex() // STEP 1 -- Add mesh vertices to octree { AXOM_ANNOTATE_SCOPE("insert vertices"); + timer.reset(); timer.start(); int numMeshVerts = m_meshWrapper.numMeshVertices(); for(int idx = 0; idx < numMeshVerts; ++idx) @@ -464,6 +460,7 @@ void InOutOctree::generateIndex() // STEP 1(b) -- Update the mesh vertices and cells with after vertex welding from octree { AXOM_ANNOTATE_SCOPE("update surface mesh vertices"); + timer.reset(); timer.start(); updateSurfaceMeshVertices(); timer.stop(); @@ -496,6 +493,7 @@ void InOutOctree::generateIndex() // STEP 2 -- Add mesh cells (segments in 2D; triangles in 3D) to octree { AXOM_ANNOTATE_SCOPE("insert mesh cells"); + timer.reset(); timer.start(); insertMeshCells(); timer.stop(); @@ -509,6 +507,7 @@ void InOutOctree::generateIndex() // -- Black (in), White(out), Gray(Intersects surface) { AXOM_ANNOTATE_SCOPE("color octree leaves"); + timer.reset(); timer.start(); colorOctreeLeaves(); @@ -519,8 +518,8 @@ void InOutOctree::generateIndex() "\t--Coloring octree leaves took {:.3Lf} seconds.", timer.elapsed())); -// -- Print some stats about the octree #ifdef DUMP_OCTREE_INFO + // -- Print some stats about the octree SLIC_INFO_ROOT("** Octree stats after inserting cells"); { AXOM_ANNOTATE_SCOPE("dump stats after inserting cells"); @@ -533,9 +532,11 @@ void InOutOctree::generateIndex() AXOM_ANNOTATE_SCOPE("validate after inserting cells"); checkValid(); } + // CLEANUP -- Finally, fix up the surface mesh after octree operations { AXOM_ANNOTATE_SCOPE("regenerate surface mesh"); + timer.reset(); timer.start(); m_meshWrapper.regenerateSurfaceMesh(); timer.stop(); diff --git a/src/axom/quest/detail/inout/BlockData.hpp b/src/axom/quest/detail/inout/BlockData.hpp index 5b4b560d1b..64326bc082 100644 --- a/src/axom/quest/detail/inout/BlockData.hpp +++ b/src/axom/quest/detail/inout/BlockData.hpp @@ -60,13 +60,13 @@ class InOutBlockData */ InOutBlockData() : m_idx(LEAF_BLOCK_UNCOLORED) { } - /** \brief Constructor from a given index */ + /// \brief Constructor from a given index explicit InOutBlockData(int dataIdx) : m_idx(dataIdx) { } - /** \brief Copy constructor for an InOutBlockData instance */ + /// \brief Copy constructor for an InOutBlockData instance InOutBlockData(const InOutBlockData& other) : m_idx(other.m_idx) { } - /** \brief Assignment operator for an InOutBlockData instance */ + /// \brief Assignment operator for an InOutBlockData instance InOutBlockData& operator=(const InOutBlockData& other) { this->m_idx = other.m_idx; @@ -81,17 +81,16 @@ class InOutBlockData */ bool isLeaf() const { return m_idx > INTERNAL_BLOCK; } - /** \brief Marks the associated block as internal */ + /// \brief Marks the associated block as internal void setInternal() { m_idx = INTERNAL_BLOCK; } - /** \brief Marks the associated block as a non-block (i.e. not in the tree) */ + /// \brief Marks the associated block as a non-block (i.e. not in the tree) void setNonBlock() { m_idx = NON_BLOCK; } /** * \brief Predicate to determine if the associated block is in the tree * - * \return True, if the block is in the tree (internal or leaf), False - * otherwise + * \return True, if the block is in the tree (internal or leaf), False otherwise */ bool isBlock() const { return m_idx != NON_BLOCK; } @@ -111,9 +110,7 @@ class InOutBlockData * */ bool hasData() const { return m_idx >= 0; } - /** - * Returns the index of the data associated with the block - */ + /// Returns the index of the data associated with the block const int& dataIndex() const { //SLIC_ASSERT(hasData()); @@ -134,24 +131,24 @@ class InOutBlockData m_idx = idx; } - /** Marks the block as Black (the entire domain is inside the surface) */ + /// Marks the block as Black (the entire domain is inside the surface) void setBlack() { SLIC_ASSERT(isLeaf()); m_idx = LEAF_BLOCK_BLACK; } - /** Marks the block as Black (the entire domain is outside the surface) */ + /// Marks the block as Black (the entire domain is outside the surface) void setWhite() { SLIC_ASSERT(isLeaf()); m_idx = LEAF_BLOCK_WHITE; } - /** Sets the data associated with the block to the given index idx */ + /// Sets the data associated with the block to the given index idx void setData(int idx) { m_idx = idx; } - /** Marks the block as uncolored */ + /// Marks the block as uncolored void setUncoloredLeaf() { SLIC_ASSERT(isLeaf()); @@ -187,13 +184,14 @@ class InOutBlockData return Undetermined; } - /** Predicate to determine if the associated block has a color + /** + * \brief Predicate to determine if the associated block has a color * \return True if the block has a color, false otherwise * \sa color() */ bool isColored() const { return color() != Undetermined; } - /** Friend function to compare equality of two InOutBlockData instances */ + /// Friend function to compare equality of two InOutBlockData instances friend bool operator==(const InOutBlockData& lhs, const InOutBlockData& rhs) { return lhs.m_idx == rhs.m_idx; @@ -204,7 +202,7 @@ class InOutBlockData }; /** - * Free function to print an InOutBlockData to an output stream + * \brief Free function to print an InOutBlockData to an output stream * \param os The output stream to write to * \param iob The InOUtBlockData instance that we are writing */ @@ -274,9 +272,7 @@ class DynamicGrayBlockData using CellList = std::vector; public: - /** - * \brief Default constructor for an InOutLeafData - */ + /// \brief Default constructor for an InOutLeafData DynamicGrayBlockData() : m_vertIndex(NO_VERTEX), m_isLeaf(true) { } /** @@ -287,18 +283,14 @@ class DynamicGrayBlockData */ DynamicGrayBlockData(VertexIndex vInd, bool isLeaf) : m_vertIndex(vInd), m_isLeaf(isLeaf) { } - /** - * \brief Copy constructor for an DynamicGrayBlockData instance - */ + /// \brief Copy constructor for an DynamicGrayBlockData instance DynamicGrayBlockData(const DynamicGrayBlockData& other) : m_vertIndex(other.m_vertIndex) , m_cells(other.m_cells) , m_isLeaf(other.m_isLeaf) { } - /** - * \brief Assignment operator for an InOutLeafData instance - */ + /// \brief Assignment operator for an InOutLeafData instance DynamicGrayBlockData& operator=(const DynamicGrayBlockData& other) { this->m_vertIndex = other.m_vertIndex; @@ -334,52 +326,48 @@ class DynamicGrayBlockData } public: // Functions related to whether this is a leaf - /** Predicate to determine if the associated block is a leaf in the octree */ + /// Predicate to determine if the associated block is a leaf in the octree bool isLeaf() const { return m_isLeaf; } - /** Sets a flag to determine whether the associated block is a leaf or internal */ + /// Sets a flag to determine whether the associated block is a leaf or internal void setLeafFlag(bool isLeaf) { m_isLeaf = isLeaf; } public: // Functions related to the associated vertex - /** - * \brief Checks whether there is a vertex associated with this leaf - */ + /// \brief Checks whether there is a vertex associated with this leaf bool hasVertex() const { return m_vertIndex >= 0; } - /** Sets the vertex associated with this leaf */ + /// Sets the vertex associated with this leaf void setVertex(VertexIndex vInd) { m_vertIndex = vInd; } - /** Clears the associated vertex index */ + /// Clears the associated vertex index void clearVertex() { m_vertIndex = NO_VERTEX; } - /** Accessor for the index of the vertex associated with this leaf */ + /// Accessor for the index of the vertex associated with this leaf VertexIndex& vertexIndex() { return m_vertIndex; } - /** Constant accessor for the index of the vertex associated with this leaf */ + /// Constant accessor for the index of the vertex associated with this leaf const VertexIndex& vertexIndex() const { return m_vertIndex; } public: // Functions related to the associated triangles - /** Check whether this Leaf has any associated triangles */ + /// Check whether this Leaf has any associated triangles bool hasCells() const { return !m_cells.empty(); } /** - * Reserves space for a given number of triangles + * \brief Reserves space for a given number of triangles * \param count The number of triangles for which to reserve space */ void reserveCells(int count) { m_cells.reserve(count); } - /** Find the number of triangles associated with this leaf */ + /// Find the number of triangles associated with this leaf int numCells() const { return static_cast(m_cells.size()); } - /** Associates the surface triangle with the given index with this block */ + /// Associates the surface triangle with the given index with this block void addCell(CellIndex tInd) { m_cells.push_back(tInd); } - /** Returns a const reference to the list of triangle indexes associated with - the block */ + /// Returns a const reference to the list of triangle indexes associated with the block const CellList& cells() const { return m_cells; } - /** Returns a reference to the list of triangle indexes associated with the - block */ + /// Returns a reference to the list of triangle indexes associated with the block CellList& cells() { return m_cells; } private: @@ -388,9 +376,7 @@ class DynamicGrayBlockData bool m_isLeaf; }; -/** - * Free function to print a DynamicGrayBlockData instance to an output stream - */ +/// Free function to print a DynamicGrayBlockData instance to an output stream inline std::ostream& operator<<(std::ostream& os, const DynamicGrayBlockData& bData) { os << "DynamicGrayBlockData{"; diff --git a/src/axom/quest/detail/inout/InOutOctreeMeshDumper.hpp b/src/axom/quest/detail/inout/InOutOctreeMeshDumper.hpp index 5f7b7b00f9..6e65f4a7e4 100644 --- a/src/axom/quest/detail/inout/InOutOctreeMeshDumper.hpp +++ b/src/axom/quest/detail/inout/InOutOctreeMeshDumper.hpp @@ -23,6 +23,8 @@ #include "axom/fmt.hpp" +#include + namespace axom { namespace quest diff --git a/src/axom/quest/examples/containment_driver.cpp b/src/axom/quest/examples/containment_driver.cpp index e6c583f8a7..1282fdc5d4 100644 --- a/src/axom/quest/examples/containment_driver.cpp +++ b/src/axom/quest/examples/containment_driver.cpp @@ -5,7 +5,7 @@ /*! * \file containment_driver.cpp - * \brief Basic demo of point containment acceleration structure over surfaces. + * \brief Demo of InOutOctree point containment acceleration structure over surfaces. */ // Axom includes @@ -401,10 +401,7 @@ class ContainmentDriver nVerts, nCells)); - if(dimension() == 3) - { - SLIC_INFO("Edge length range: " << meshEdgeLenRange); - } + SLIC_INFO_IF(dimension() == 3, "Edge length range: " << meshEdgeLenRange); SLIC_INFO("Cell area range is: " << meshCellAreaRange); if(dimension() == 3) diff --git a/src/axom/quest/examples/delaunay_triangulation.cpp b/src/axom/quest/examples/delaunay_triangulation.cpp index d23d18900d..f027259f68 100644 --- a/src/axom/quest/examples/delaunay_triangulation.cpp +++ b/src/axom/quest/examples/delaunay_triangulation.cpp @@ -172,6 +172,7 @@ void run_delaunay(const Input& params) // Check that the Delaunay complex is valid SLIC_INFO("Checking validity of Delaunay complex and underlying mesh..."); { + timer.reset(); timer.start(); dt.getMeshData()->isValid(true); dt.isValid(true); diff --git a/src/axom/quest/examples/quest_inout_interface.cpp b/src/axom/quest/examples/quest_inout_interface.cpp index ad3d98b763..740e4a78f6 100644 --- a/src/axom/quest/examples/quest_inout_interface.cpp +++ b/src/axom/quest/examples/quest_inout_interface.cpp @@ -288,6 +288,7 @@ int main(int argc, char** argv) // -- Run the queries (the z-coordinate is ignored for 2D queries) int numInside = 0; SLIC_INFO(axom::fmt::format("Querying mesh with {} query points...", nQueryPoints)); + timer.reset(); timer.start(); for(auto& pt : queryPoints) { diff --git a/src/axom/quest/examples/scattered_interpolation.cpp b/src/axom/quest/examples/scattered_interpolation.cpp index 825e21daf4..83d238415b 100644 --- a/src/axom/quest/examples/scattered_interpolation.cpp +++ b/src/axom/quest/examples/scattered_interpolation.cpp @@ -904,6 +904,7 @@ int main(int argc, char** argv) } // Find the simplices containing each of the query points + timer.reset(); timer.start(); switch(params.dimension) { @@ -923,6 +924,7 @@ int main(int argc, char** argv) numQueryPts / timer.elapsedTimeInSec())); // Perform the interpolation on each of the fields + timer.reset(); timer.start(); switch(params.dimension) { @@ -948,6 +950,7 @@ int main(int argc, char** argv) numQueryPts * numFields / timer.elapsedTimeInSec())); // Check interpolation using local interpolation from weights + timer.reset(); timer.start(); switch(params.dimension) { diff --git a/src/axom/spin/Brood.hpp b/src/axom/spin/Brood.hpp index e348f73e17..705b124d20 100644 --- a/src/axom/spin/Brood.hpp +++ b/src/axom/spin/Brood.hpp @@ -15,16 +15,13 @@ namespace spin { /** * \class - * \brief Helper class to handle subindexing of block data within octree - * siblings + * \brief Helper class to handle subindexing of block data within octree siblings * * \note A brood is a collection of siblings that are generated simultaneously. - * \note This class converts a grid point at the given level into a brood index - * of the point. + * \note This class converts a grid point at the given level into a brood index of the point. * - * The base brood is the MortonIndex of the grid point's octree parent - * and its offset index is obtained by interleaving the least significant - * bit of its coordinates. + * The base brood is the MortonIndex of the grid point's octree parent and its offset index + * is obtained by interleaving the least significant bit of its coordinates. */ template struct Brood @@ -59,8 +56,7 @@ struct Brood /** \brief Accessor for the offset of the point within the brood */ const int& offset() const { return m_offset; } - /** \brief Reconstruct a grid point from a brood's Morton index and an offset - */ + /// \brief Reconstruct a grid point from a brood's Morton index and an offset static GridPt reconstructGridPt(MortonIndexType morton, int offset) { return static_cast( @@ -68,10 +64,8 @@ struct Brood } private: - MortonIndexType m_broodIdx; /** MortonIndex of the base point of all - blocks within the brood */ - int m_offset; /** Index of the block within the brood. - Value is in [0, 2^DIM) */ + MortonIndexType m_broodIdx; /// MortonIndex of the base point of all blocks within the brood + int m_offset; /// Index of the block within the brood. Value is in [0, 2^DIM) }; /** @@ -98,18 +92,17 @@ struct Brood { for(int i = 0; i < GridPt::DIMENSION; ++i) { - m_offset |= (pt[i] & 1) << i; // interleave the least significant - // bits + m_offset |= (pt[i] & 1) << i; // interleave the least significant bits } } - /** \brief Accessor for the base point of the entire brood */ + /// \brief Accessor for the base point of the entire brood const GridPt& base() const { return m_broodPt; } - /** \brief Accessor for the index of the point within the brood */ + /// \brief Accessor for the index of the point within the brood const int& offset() const { return m_offset; } - /** \brief Reconstruct a grid point from a brood's base point and an offset */ + /// \brief Reconstruct a grid point from a brood's base point and an offset static GridPt reconstructGridPt(const GridPt& pt, int offset) { // shift and add offset to each coordinate @@ -123,9 +116,8 @@ struct Brood } private: - GridPt m_broodPt; /** Base point of all blocks within the brood */ - int m_offset; /** Index of the block within the brood. Value is in - [0, 2^DIM) */ + GridPt m_broodPt; /// Base point of all blocks within the brood + int m_offset; /// Index of the block within the brood. Value is in [0, 2^DIM) }; } // end namespace spin diff --git a/src/axom/spin/DenseOctreeLevel.hpp b/src/axom/spin/DenseOctreeLevel.hpp index d0040a2658..5a3c0fc311 100644 --- a/src/axom/spin/DenseOctreeLevel.hpp +++ b/src/axom/spin/DenseOctreeLevel.hpp @@ -54,8 +54,7 @@ class DenseOctreeLevel : public OctreeLevel * \brief Concrete instance of the BlockIteratorHelper class defined in the * OctreeLevel base class. * - * \note ParenType must be either BlockIteratorHelper or - * ConstBlockIteratorHelper, + * \note ParenType must be either BlockIteratorHelper or ConstBlockIteratorHelper, * both are defined in the OctreeLevel base class */ template @@ -80,7 +79,7 @@ class DenseOctreeLevel : public OctreeLevel } } - /** Increment to next block of the level */ + /// Increment to next block of the level void increment() { // Note, must skip blocks that are not in the tree @@ -95,40 +94,36 @@ class DenseOctreeLevel : public OctreeLevel } while(m_currentIdx < m_endIdx && !data()->isBlock()); } - /** Access to point associated with the block pointed to by the iterator */ + /// Access to point associated with the block pointed to by the iterator GridPt pt() const { // Reconstruct the grid point from its brood representation return BroodType::reconstructGridPt(m_currentIdx, m_offset); } - /** Accessor for data associated with the iterator's block */ + /// Accessor for data associated with the iterator's block BlockDataType* data() { return &m_octreeLevel->m_data[m_currentIdx][m_offset]; } - /** Const accessor for data associated with the iterator's block */ + /// Const accessor for data associated with the iterator's block const BlockDataType* data() const { return &m_octreeLevel->m_data[m_currentIdx][m_offset]; } - /** \brief Predicate to determine if two block iterators are the same */ + /// \brief Predicate to determine if two block iterators are the same bool equal(const BaseBlockItType* other) { const self* pother = dynamic_cast(other); - return (pother != nullptr) && - (m_currentIdx == pother->m_currentIdx) // iterators - // are the same - && (m_offset == pother->m_offset); // brood - // indices are - // the same + return (pother != nullptr) && (m_currentIdx == pother->m_currentIdx) // iterators are the same + && (m_offset == pother->m_offset); // brood indices are the same } private: - OctreeLevelType* m_octreeLevel; + OctreeLevelType* m_octreeLevel {nullptr}; MortonIndexType m_currentIdx, m_endIdx; int m_offset; bool m_isLevelZero; }; public: - /** \brief Default constructor for an octree level */ + /// \brief Default constructor for an octree level DenseOctreeLevel(int level = -1) : Base(level), m_blockCount(0) { if(level < 0) @@ -182,8 +177,8 @@ class DenseOctreeLevel : public OctreeLevel /** * \brief Factory function to return a ConstGridBlockIterHelper for this level * - * \param begin A boolean to determine if this is to be - * a begin (true) or end (false) iterator + * \param begin A boolean to determine if this is to be + * a begin (true) or end (false) iterator */ ConstBaseBlockIteratorHelper* getIteratorHelper(bool begin) const { @@ -232,14 +227,14 @@ class DenseOctreeLevel : public OctreeLevel } } - /** \brief Accessor for the data associated with pt */ + /// \brief Accessor for the data associated with pt BlockDataType& operator[](const GridPt& pt) { const BroodType brood(pt); return m_data[brood.base()][brood.offset()]; } - /** \brief Const accessor for the data associated with pt */ + /// \brief Const accessor for the data associated with pt const BlockDataType& operator[](const GridPt& pt) const { SLIC_ASSERT_MSG(hasBlock(pt), @@ -249,28 +244,28 @@ class DenseOctreeLevel : public OctreeLevel return m_data[brood.base()][brood.offset()]; } - /** \brief Access the data associated with the entire brood */ + /// \brief Access the data associated with the entire brood BroodData& getBroodData(const GridPt& pt) { return m_data[BroodType::MortonizerType::mortonize(pt)]; } - /** \brief Const access to data associated with the entire brood */ + /// \brief Const access to data associated with the entire brood const BroodData& getBroodData(const GridPt& pt) const { return m_data[BroodType::MortonizerType::mortonize(pt)]; } - /** \brief Predicate to check if there are any blocks in this octree level */ + /// \brief Predicate to check if there are any blocks in this octree level bool empty() const { return m_blockCount == 0; } - /** \brief Returns the number of blocks (internal and leaf) in the level */ + /// \brief Returns the number of blocks (internal and leaf) in the level int numBlocks() const { return m_blockCount; } - /** \brief Returns the number of internal blocks in the level */ + /// \brief Returns the number of internal blocks in the level int numInternalBlocks() const { return numBlocks() - numLeafBlocks(); } - /** \brief Returns the number of leaf blocks in the level */ + /// \brief Returns the number of leaf blocks in the level int numLeafBlocks() const { if(empty()) @@ -294,12 +289,10 @@ class DenseOctreeLevel : public OctreeLevel } /** - * \brief Helper function to determine the status of - * an octree block within this octree level + * \brief Helper function to determine the status of an octree block within this octree level * * \param pt The grid point of the block index that we are testing - * \return The status of the grid point pt (e.g. LeafBlock, InternalBlock, - *...) + * \return The status of the grid point pt (e.g. LeafBlock, InternalBlock, ...) */ TreeBlockStatus blockStatus(const GridPt& pt) const { @@ -319,7 +312,7 @@ class DenseOctreeLevel : public OctreeLevel DISABLE_MOVE_AND_ASSIGNMENT(DenseOctreeLevel); private: - BroodData* m_data; + BroodData* m_data {nullptr}; int m_broodCapacity; int m_blockCount; }; diff --git a/src/axom/spin/MortonIndex.hpp b/src/axom/spin/MortonIndex.hpp index e17c2226ce..16cfab452b 100644 --- a/src/axom/spin/MortonIndex.hpp +++ b/src/axom/spin/MortonIndex.hpp @@ -9,9 +9,8 @@ * \brief Classes and functions to convert between points on an integer grid and * their unidimensional MortonIndex. * - * Also has some utility functions for 'mortonizing' and 'demortonizing' points - * and a PointHash functor class that can be used as a std::hash for - * unordered_maps + * Also has some utility functions for 'mortonizing' and 'demortonizing' points and a PointHash + * functor class that can be used as a std::hash for unordered_maps and axom::FlatMap */ #ifndef AXOM_SPIN_MORTON_INDEX_HPP_ @@ -19,7 +18,7 @@ #include "axom/config.hpp" #include "axom/core/Types.hpp" -#include "axom/core/Macros.hpp" // defines AXOM_STATIC_ASSERT +#include "axom/core/Macros.hpp" #include "axom/core/NumericLimits.hpp" #include "axom/primal/geometry/Point.hpp" @@ -36,78 +35,51 @@ namespace template struct NumReps { - enum - { - value = 5 - }; + static constexpr int value = 5; }; template <> struct NumReps { - enum - { - value = 5 - }; + static constexpr int value = 5; }; template <> struct NumReps { - enum - { - value = 5 - }; + static constexpr int value = 5; }; template <> struct NumReps { - enum - { - value = 4 - }; + static constexpr int value = 4; }; template <> struct NumReps { - enum - { - value = 4 - }; + static constexpr int value = 4; }; template <> struct NumReps { - enum - { - value = 3 - }; + static constexpr int value = 3; }; template <> struct NumReps { - enum - { - value = 3 - }; + static constexpr int value = 3; }; template <> struct NumReps { - enum - { - value = 2 - }; + static constexpr int value = 2; }; template <> struct NumReps { - enum - { - value = 2 - }; + static constexpr int value = 2; }; } // namespace @@ -232,39 +204,19 @@ struct Mortonizer using self = Mortonizer; using Base = MortonBase; - // Magic numbers in 2D - AXOM_HOST_DEVICE static MortonIndexType GetB(int i) - { - constexpr MortonIndexType B[] = {static_cast(0x5555555555555555), // 0101'0101 - static_cast(0x3333333333333333), // 0011'0011 - static_cast(0x0F0F0F0F0F0F0F0F), // 0000'1111 - static_cast(0x00FF00FF00FF00FF), // 0x8 1x8 - static_cast(0x0000FFFF0000FFFF), // 0x16 1x16 - static_cast(0x00000000FFFFFFFF)}; // 0x32 1x32; - return B[i]; - } - - AXOM_HOST_DEVICE static int GetS(int i) - { - constexpr int S[] = {1, 2, 4, 8, 16, 32}; - return S[i]; - } - - enum - { - /*! The dimension of the Mortonizer */ - NDIM = 2, + /*! The dimension of the Mortonizer */ + static constexpr int NDIM = 2; - /*! The number of bits in a CoordType */ - COORD_BITS = axom::numeric_limits::digits, + /*! The number of bits in a CoordType */ + static constexpr int COORD_BITS = axom::numeric_limits::digits; - /*! The number of bits in a MortonIndex */ - MORTON_BITS = axom::numeric_limits::digits, + /*! The number of bits in a MortonIndex */ + static constexpr int MORTON_BITS = axom::numeric_limits::digits; - /*! The number of representable Morton bits per dimension */ - MB_PER_DIM = MORTON_BITS / NDIM, + /*! The number of representable Morton bits per dimension */ + static constexpr int MB_PER_DIM = MORTON_BITS / NDIM; - /*! + /*! * The maximum number of unique bits from each coordinate of type CoordType * that can be represented in a MortonIndex. * @@ -274,20 +226,37 @@ struct Mortonizer * coordinate. * */ - MAX_UNIQUE_BITS = (MB_PER_DIM < COORD_BITS) ? MB_PER_DIM : COORD_BITS, + static constexpr int MAX_UNIQUE_BITS = (MB_PER_DIM < COORD_BITS) ? MB_PER_DIM : COORD_BITS; - /*! + /*! * The number of iterations required for converting from MortonIndexes to * CoordType using the bit interleaving algorithm in MortonBase. */ - CONTRACT_MAX_ITER = NumReps::value, + static constexpr int CONTRACT_MAX_ITER = NumReps::value; - /*! + /*! * The number of iterations required for converting between CoordTypes * and MortonIndexes using the bit interleaving algorithm in MortonBase. */ - EXPAND_MAX_ITER = NumReps::value - }; + static constexpr int EXPAND_MAX_ITER = NumReps::value; + + // Magic numbers in 2D + AXOM_HOST_DEVICE static MortonIndexType GetB(int i) + { + constexpr MortonIndexType B[] = {static_cast(0x5555555555555555), // 0101'0101 + static_cast(0x3333333333333333), // 0011'0011 + static_cast(0x0F0F0F0F0F0F0F0F), // 0000'1111 + static_cast(0x00FF00FF00FF00FF), // 0x8 1x8 + static_cast(0x0000FFFF0000FFFF), // 0x16 1x16 + static_cast(0x00000000FFFFFFFF)}; // 0x32 1x32; + return B[i]; + } + + AXOM_HOST_DEVICE static int GetS(int i) + { + constexpr int S[] = {1, 2, 4, 8, 16, 32}; + return S[i]; + } /*! * \brief A function to convert a 2D point to a Morton index @@ -372,8 +341,41 @@ struct Mortonizer using self = Mortonizer; using Base = MortonBase; - // Magic numbers in 3D from C. Ericson's Real Time Collision Detection book + /*! The dimension of the Mortonizer */ + constexpr static int NDIM = 3; + + /*! The number of bits in a CoordType */ + constexpr static int COORD_BITS = axom::numeric_limits::digits; + + /*! The number of bits in a MortonIndex */ + constexpr static int MORTON_BITS = axom::numeric_limits::digits; + + /*! The number of representable morton bits per dimension */ + constexpr static int MB_PER_DIM = MORTON_BITS / NDIM; + + /*! + * The maximum number of unique bits from each coordinate of type CoordType + * that can be represented in a MortonIndex. + * \note If we are use Mortonizer as a (one-way) hash function, + * it is ok to use more bits. But, if we would like to be + * able to reverse the MortonIndex, then we cannot safely use + * more than MAX_UNIQUE_BITS per coordinate. + */ + constexpr static int MAX_UNIQUE_BITS = (MB_PER_DIM < COORD_BITS) ? MB_PER_DIM : COORD_BITS; + + /*! + * The number of iterations required for converting from MortonIndexes + * to CoordType using the bit interleaving algorithm in MortonBase. + */ + constexpr static int CONTRACT_MAX_ITER = NumReps::value; + /*! + * The number of iterations required for converting between CoordTypes + * and MortonIndexes using the bit interleaving algorithm in MortonBase. + */ + constexpr static int EXPAND_MAX_ITER = NumReps::value - 1; + + // Magic numbers in 3D from C. Ericson's Real Time Collision Detection book AXOM_HOST_DEVICE static MortonIndexType GetB(int i) { constexpr MortonIndexType B[] = { @@ -392,43 +394,6 @@ struct Mortonizer return S[i]; } - enum - { - /*! The dimension of the Mortonizer */ - NDIM = 3, - - /*! The number of bits in a CoordType */ - COORD_BITS = axom::numeric_limits::digits, - - /*! The number of bits in a MortonIndex */ - MORTON_BITS = axom::numeric_limits::digits, - - /*! The number of representable morton bits per dimension */ - MB_PER_DIM = MORTON_BITS / NDIM, - - /*! - * The maximum number of unique bits from each coordinate of type CoordType - * that can be represented in a MortonIndex. - * \note If we are use Mortonizer as a (one-way) hash function, - * it is ok to use more bits. But, if we would like to be - * able to reverse the MortonIndex, then we cannot safely use - * more than MAX_UNIQUE_BITS per coordinate. - */ - MAX_UNIQUE_BITS = (MB_PER_DIM < COORD_BITS) ? MB_PER_DIM : COORD_BITS, - - /*! - * The number of iterations required for converting from MortonIndexes - * to CoordType using the bit interleaving algorithm in MortonBase. - */ - CONTRACT_MAX_ITER = NumReps::value, - - /*! - * The number of iterations required for converting between CoordTypes - * and MortonIndexes using the bit interleaving algorithm in MortonBase. - */ - EXPAND_MAX_ITER = NumReps::value - 1 - }; - /*! * \brief A function to convert a 3D point to a Morton index * @@ -548,6 +513,7 @@ template struct PointHash { using MortonIndex = std::size_t; + using result_type = std::size_t; /*! * \brief Mortonizes a coordinate (viewed as a 1D point) diff --git a/src/axom/spin/OctreeBase.hpp b/src/axom/spin/OctreeBase.hpp index 2077183fc3..bb17c5eab8 100644 --- a/src/axom/spin/OctreeBase.hpp +++ b/src/axom/spin/OctreeBase.hpp @@ -57,25 +57,16 @@ class BlockData return *this; } - /** - * \brief Predicate to determine if this BlockData represents a leaf block in - * the tree - */ + /// \brief Predicate to determine if this BlockData represents a leaf block in the tree bool isLeaf() const { return m_id >= 0; } - /** - * \brief Sets the data associated with this block - */ + /// \brief Sets the data associated with this block void setData(int blockID) { m_id = blockID; } - /** - * Marks the block as not in the octree - */ + /// Marks the block as not in the octree void setNonBlock() { m_id = NON_BLOCK; } - /** - * Predicate to check if the associated block is in the octree - */ + /// Predicate to check if the associated block is in the octree bool isBlock() const { return m_id != NON_BLOCK; } /** @@ -86,9 +77,7 @@ class BlockData const int& dataIndex() const { return m_id; } - /** - * \brief Sets the block type to internal - */ + /// \brief Sets the block type to internal void setInternal() { if(isLeaf()) @@ -97,9 +86,7 @@ class BlockData } } - /** - * Equality operator for comparing two BlockData instances - */ + /// Equality operator for comparing two BlockData instances friend bool operator==(const BlockData& lhs, const BlockData& rhs) { return lhs.m_id == rhs.m_id; @@ -151,15 +138,13 @@ class OctreeBase * \brief Inner class encapsulating the index of an octree block. * * Each block index is represented as a point on an integer grid (the minimum - * point of the block's extent) - * at a given level of resolution. + * point of the block's extent) at a given level of resolution. * * Each level of resolution is a regular grid with \f$ 2^{level} \f$ * grid points along each dimension. The root block (at level 0) * covers the entire domain. * An octree block at level \f$ \ell \f$ has \f$ 2^{DIM} \f$ children - * at level \f$ \ell + 1 \f$ - * covering its domain. + * at level \f$ \ell + 1 \f$ covering its domain. */ class BlockIndex { @@ -181,14 +166,10 @@ class OctreeBase using FaceNeighborIndexSet = slam::OrderedSet; public: - /** - * \brief Default constructor - */ + /// \brief Default constructor BlockIndex() : m_pt(GridPt()), m_lev(0) { } - /** - * \brief Constructor from a point and a level - */ + /// \brief Constructor from a point and a level BlockIndex(const GridPt& pt, int level) : m_pt(pt), m_lev(level) { } /** @@ -219,26 +200,19 @@ class OctreeBase */ int& level() { return m_lev; } - /** - * \brief The level of the block index's parent - */ + /// \brief The level of the block index's parent int parentLevel() const { return m_lev - 1; } - /** - * \brief The level of the block index's child - */ + /// \brief The level of the block index's child int childLevel() const { return m_lev + 1; } - /** - * \brief Returns the grid point of the block's parent - */ + /// \brief Returns the grid point of the block's parent GridPt parentPt() const { return GridPt(m_pt.array() / 2); } /** * \brief Returns the grid point of the block's child at index childIndex * - * \param [in] childIndex The index of the child whose grid point we are - * finding + * \param [in] childIndex The index of the child whose grid point we are finding * \pre \f$ 0 \le childIndex < \f$ Octree::NUM_CHILDREN */ GridPt childPt(int childIndex) const @@ -258,10 +232,7 @@ class OctreeBase return cPoint; } - /** - * \brief Returns a grid point at the specified offset - * from the current block index's point - */ + /// \brief Returns a grid point at the specified offset from the current block index's point GridPt neighborPt(const GridPt& offset) const { GridPt nPoint(m_pt); @@ -284,8 +255,7 @@ class OctreeBase /** * \brief Returns the child BlockIndex of this block * - * \param [in] childIndex The index of the child whose grid point we are - * finding + * \param [in] childIndex The index of the child whose grid point we are finding * \pre \f$ 0 \le childIndex < \f$ Octree::NUM_CHILDREN */ BlockIndex child(int childIndex) const { return BlockIndex(childPt(childIndex), childLevel()); } @@ -295,10 +265,8 @@ class OctreeBase * * \pre 0 <= neighborIndex < 2 * DIM * \note The face neighbors indices cycle through the dimensions, two per - * dimension, - * e.g. Neighbor 0 is at offset (-1, 0,0,...,0), neighbor 1 is at offset - *(1,0,0,..,0) - * and neighbor 2 is at offset (0,-1, 0, 0, ..., 0) etc... + * dimension, e.g. Neighbor 0 is at offset (-1, 0,0,...,0), + * neighbor 1 is at offset (1,0,0,..,0) and neighbor 2 is at offset (0,-1, 0, 0, ..., 0) etc... */ BlockIndex faceNeighbor(int neighborIndex) const { @@ -346,8 +314,8 @@ class OctreeBase /** * \brief Checks the validity of the index. * - * A block index is valid when its level is \f$ \ge 0 \f$ - * and it is inBounds + * A block index is valid when its level is \f$ \ge 0 \f$ and it is inBounds + * * \returns true if the block index is valid, else false * \see inBounds() */ @@ -375,16 +343,12 @@ class OctreeBase } /** - * \brief Predicate to determine if the block instance is a descendant of - * ancestor block + * \brief Predicate to determine if the block instance is a descendant of ancestor block * * \param ancestor The potential ancestor of the block - * \note A block is an ancestor of another block if neither block is an - * invalid_index() - * and the block's are equivalent after 0 or more calls to - * BlockIndex::parent() - * \return True, if the block instance is a descendant of the ancestor - * block + * \note A block is an ancestor of another block if neither block is an invalid_index() + * and the block's are equivalent after 0 or more calls to BlockIndex::parent() + * \return True, if the block instance is a descendant of the ancestor block */ bool isDescendantOf(const BlockIndex& ancestor) const { @@ -426,15 +390,10 @@ class OctreeBase */ static BlockIndex invalid_index() { return BlockIndex(GridPt(), -1); } - /** - * \brief The number of children that an octree block can have - */ + /// \brief The number of children that an octree block can have static int numChildren() { return ChildIndexSet().size(); } - /** - * \brief The number of face neighbors that an octree block - * can have (ignoring boundaries) - */ + /// \brief The number of face neighbors that an octree block can have (ignoring boundaries) static int numFaceNeighbors() { return FaceNeighborIndexSet().size(); } private: @@ -463,11 +422,7 @@ class OctreeBase using Sparse64OctLevPtr = Sparse64OctLevType*; using SparsePtOctLevPtr = SparsePtOctLevType*; - /** - * \brief Simple utility to check if a pointer of type BasePtrType - * - * can be cast to a pointer of type DerivedPtrType - */ + /// \brief Simple utility to check if a pointer of type BasePtrType can be cast to a pointer of type DerivedPtrType template bool checkCast(BasePtrType base) const { @@ -475,9 +430,7 @@ class OctreeBase } public: - /** - * Sets up an octree containing only the root block - */ + /// Sets up an octree containing only the root block OctreeBase() : m_leavesLevelMap(&m_levels) { for(int i = 0; i < maxLeafLevel(); ++i) @@ -514,9 +467,7 @@ class OctreeBase (*m_leavesLevelMap[rootBlock.level()]).addAllChildren(rootBlock.pt()); } - /** - * \brief OctreeBase destructor - */ + /// \brief OctreeBase destructor ~OctreeBase() { for(int i = 0; i < maxLeafLevel(); ++i) @@ -526,22 +477,17 @@ class OctreeBase } } - /** - * \brief The max level for leaf blocks of the octree - */ + /// \brief The max level for leaf blocks of the octree int maxLeafLevel() const { return m_levels.size(); } - /** - * \brief The max level for internal blocks of the octree - */ + /// \brief The max level for internal blocks of the octree int maxInternalLevel() const { return m_levels.size() - 1; } public: //@{ /** - * \brief Utility function to find the number of (possible) - * grid cells at a given level or resolution + * \brief Utility function to find the number of (possible) grid cells at a given level or resolution * * \param [in] level The level or resolution. * \pre \f$ 0 \le lev @@ -558,8 +504,7 @@ class OctreeBase /** * Auxiliary function to return the root of the octree - * \note The root block has no parent. - * Its parent is an invalid BlockIndex. + * \note The root block has no parent. Its parent is an invalid BlockIndex. * I.e. octree.parent( octree.root()).isValid() = false. */ static BlockIndex root() { return BlockIndex(); } @@ -619,15 +564,10 @@ class OctreeBase // @} - /** - * \brief Accessor for a reference to the octree level instance at level lev - */ + /// \brief Accessor for a reference to the octree level instance at level lev OctreeLevelType& getOctreeLevel(int lev) { return *m_leavesLevelMap[lev]; } - /** - * \brief Const accessor for a reference to the octree level instance at level - * lev - */ + /// \brief Const accessor for a reference to the octree level instance at level lev const OctreeLevelType& getOctreeLevel(int lev) const { return *m_leavesLevelMap[lev]; } public: @@ -644,8 +584,7 @@ class OctreeBase * * \param [in] pt The grid point to check * \param [in] lev The level of the grid point - * \returns true if the associated block is a leaf in the octree, false - * otherwise + * \returns true if the associated block is a leaf in the octree, false otherwise */ bool isLeaf(const GridPt& pt, int lev) const { @@ -653,12 +592,10 @@ class OctreeBase } /** - * \brief Determine whether the octree contains a leaf block associated with - * this BlockIndex + * \brief Determine whether the octree contains a leaf block associated with this BlockIndex * * \param [in] block The BlockIndex of the tree to check - * \returns true if the associated block is a leaf in the octree, false - * otherwise + * \returns true if the associated block is a leaf in the octree, false otherwise */ bool isLeaf(const BlockIndex& block) const { @@ -671,8 +608,7 @@ class OctreeBase * * \param [in] pt The grid point to check * \param [in] lev The level of the grid point - * \returns true if the associated block is an internal block of the octree, - * false otherwise + * \returns true if the associated block is an internal block of the octree, false otherwise */ bool isInternal(const GridPt& pt, int lev) const { @@ -680,8 +616,7 @@ class OctreeBase } /** - * \brief Determine whether the octree contains an internal block associated - * with this BlockIndex + * \brief Determine whether the octree contains an internal block associated with this BlockIndex * * \param [in] block The BlockIndex of the tree to check * \returns true if the associated block is an internal block of the octree, @@ -740,8 +675,7 @@ class OctreeBase * associated with this BlockIndex * * \param [in] block The BlockIndex of the tree to check - * \returns true if the associated block is a block of the octree, false - * otherwise + * \returns true if the associated block is a block of the octree, false otherwise */ bool hasBlock(const BlockIndex& block) const { return this->hasBlock(block.pt(), block.level()); } @@ -750,8 +684,7 @@ class OctreeBase * level lev is a possible block in this octree * * \note A block index is out of bounds if its level is not in the tree, or - * its grid point is out of the - * range of possible grid points for its level + * its grid point is out of the range of possible grid points for its level */ bool inBounds(const GridPt& pt, int lev) const { @@ -763,8 +696,7 @@ class OctreeBase * possible block in this octree * * \note A block index is out of bounds if its level is not in the tree, or - * its grid point is out of the - * range of possible grid points for its level + * its grid point is out of the range of possible grid points for its level */ bool inBounds(const BlockIndex& block) const { @@ -877,11 +809,9 @@ class OctreeBase * \param checkInBounds A flag to determine if we should check that * the block lies within the octree bounds (default=true) * \post The returned block, if valid, is blk or one of its ancestor blocks. - * \return The blockIndex of the finest octree leaf covering blk, if it - * exists, + * \return The blockIndex of the finest octree leaf covering blk, if it exists, * BlockIndex::invalid_index otherwise (e.g. blk is an internal block of - * the tree - * or is out of bounds) + * the tree or is out of bounds) */ BlockIndex coveringLeafBlock(const BlockIndex& blk, bool checkInBounds = true) const { @@ -893,8 +823,7 @@ class OctreeBase switch(blockStatus(blk)) { - case BlockNotInTree: // Find its nearest ancestor in the tree (it will be - // a leaf) + case BlockNotInTree: // Find its nearest ancestor in the tree (it will be a leaf) { BlockIndex ancBlk = blk.parent(); while(!this->hasBlock(ancBlk)) @@ -912,20 +841,17 @@ class OctreeBase } SLIC_ASSERT_MSG(false, - "OctreeBase::coveringLeafBlock -- Should never get past " - "the switch statement. " - << " Perhaps a new case was added to the TreeBlock enum"); + "OctreeBase::coveringLeafBlock -- Should never get past the switch statement. " + " Perhaps a new case was added to the TreeBlock enum"); return BlockIndex::invalid_index(); } protected: /** - * \brief Helper function to determine the status of a BlockIndex within an - * octree instance + * \brief Helper function to determine the status of a BlockIndex within an octree instance * * \note This function is meant to help with implementing basic octree - * functionality - * and is not meant to be exposed in the public API + * functionality and is not meant to be exposed in the public API * \param pt The grid point of the block index that we are testing * \param lev The level of the block index that we are testing */ @@ -965,11 +891,9 @@ class OctreeBase } /** - * \brief Helper function to determine the status of a BlockIndex within an - * octree instance + * \brief Helper function to determine the status of a BlockIndex within an octree instance * - * \note This function is meant to help with implementing basic octree - * functionality + * \note This function is meant to help with implementing basic octree functionality * and is not meant to be exposed in the public API * \param blk The block index we are testing */ diff --git a/src/axom/spin/OctreeLevel.hpp b/src/axom/spin/OctreeLevel.hpp index 5408ffc5cc..845eaa15e7 100644 --- a/src/axom/spin/OctreeLevel.hpp +++ b/src/axom/spin/OctreeLevel.hpp @@ -32,10 +32,7 @@ namespace axom { namespace spin { -/** - * \brief Helper enumeration for status of a BlockIndex within an OctreeLevel - * instance - */ +/// \brief Helper enumeration for status of a BlockIndex within an OctreeLevel instance enum TreeBlockStatus { BlockNotInTree, /** Status of blocks that are not in the tree */ @@ -45,17 +42,14 @@ enum TreeBlockStatus /** * \class - * \brief An abstract base class to represent a sparse level of blocks within an - * octree. + * \brief An abstract base class to represent a sparse level of blocks within an octree. * * Each block is associated with an integer grid point whose coordinates * have values between 0 and 2^L (where L = this->level() is the encoded level). - * The OctreeLevel associates data of (templated) type BlockDataType with each - * such block. - * \note For efficiency, the data is stored within a brood (collection of - * siblings that are created simultaneously). - * \note BlockDataType must define a predicate function with the signature: bool - * isLeaf() const; + * The OctreeLevel associates data of (templated) type BlockDataType with each such block. + * \note For efficiency, the data is stored within a brood + * (collection of siblings that are created simultaneously). + * \note BlockDataType must define a predicate function with the signature: bool isLeaf() const; */ template class OctreeLevel @@ -75,63 +69,54 @@ class OctreeLevel BROOD_SIZE = 1 << DIM }; - /** A brood is a collection of sibling blocks that are generated - simultaneously */ + /// A brood is a collection of sibling blocks that are generated simultaneously using BroodData = NumericArray; - /** Predeclare the BlockIterator type */ + /// Predeclare the BlockIterator type template class BlockIterator; protected: - /** - * \brief A virtual base class to help with iteration of an OctreeLevel's - * blocks - */ + /// \brief A virtual base class to help with iteration of an OctreeLevel's blocks class BlockIteratorHelper { public: - /** Virtual destructor */ + /// Virtual destructor virtual ~BlockIteratorHelper() { } - /** \brief A function to increment to the next Block in the level */ + /// \brief A function to increment to the next Block in the level virtual void increment() = 0; - /** \brief Predicate to determine if two BlockIteratorHelpers are equivalent - */ + /// \brief Predicate to determine if two BlockIteratorHelpers are equivalent virtual bool equal(const BlockIteratorHelper* other) = 0; - /** Accessor for the point associated with the current octree block */ + /// Accessor for the point associated with the current octree block virtual GridPt pt() const = 0; - /** Accessor for the data associated with the current octree block */ + /// Accessor for the data associated with the current octree block virtual BlockDataType* data() = 0; - /** Const accessor for the data associated with the current octree block */ + /// Const accessor for the data associated with the current octree block virtual const BlockDataType* data() const = 0; }; - /** - * \brief A virtual base class to help with constant iteration of an - * OctreeLevel's blocks - */ + /// \brief A virtual base class to help with constant iteration of an OctreeLevel's blocks class ConstBlockIteratorHelper { public: - /** Virtual destructor */ + /// Virtual destructor virtual ~ConstBlockIteratorHelper() { } - /** \brief A function to increment to the next Block in the level */ + /// \brief A function to increment to the next Block in the level virtual void increment() = 0; - /** \brief Predicate to determine if two BlockIteratorHelpers are equivalent - */ + /// \brief Predicate to determine if two BlockIteratorHelpers are equivalent virtual bool equal(const ConstBlockIteratorHelper* other) = 0; - /** Accessor for the point associated with the current octree block */ + /// Accessor for the point associated with the current octree block virtual GridPt pt() const = 0; - /** Const accessor for the data associated with the current octree block */ + /// Const accessor for the data associated with the current octree block virtual const BlockDataType* data() const = 0; }; @@ -141,10 +126,10 @@ class OctreeLevel BlockIterator; public: - /** \brief Constructor of an OctreeLevel at level lev */ + /// \brief Constructor of an OctreeLevel at level \a lev OctreeLevel(int level = -1) : m_level(level) { } - /** \brief Virtual destructor of an OctreeLevel */ + /// \brief Virtual destructor of an OctreeLevel virtual ~OctreeLevel() {}; /** @@ -161,7 +146,7 @@ class OctreeLevel */ GridPt maxGridCell() const { return GridPt(maxCoord()); } - /** Accessor for the instance's level */ + /// Accessor for the instance's level int level() const { return m_level; } /** @@ -204,8 +189,7 @@ class OctreeLevel // In C++17, inheriting from std::iterator was deprecated. // We provide these typedefs for class BlockIterator to avoid inheriting - // from std::iterator and causing warnings for those compiling to C++17 - // or newer. + // from std::iterator and causing warnings for those compiling to C++17 or newer. using value_type = DataType; using difference_type = std::ptrdiff_t; using pointer = DataType*; @@ -215,8 +199,7 @@ class OctreeLevel BlockIterator(OctreeLevel* octLevel, bool begin = false) : m_octLevel(octLevel) { SLIC_ASSERT(octLevel != nullptr); - m_iterHelper = octLevel->getIteratorHelper(begin); // factory - // function + m_iterHelper = octLevel->getIteratorHelper(begin); // factory function } ~BlockIterator() @@ -228,34 +211,33 @@ class OctreeLevel } } - /** \brief A const dereference function to access the data */ + /// \brief A const dereference function to access the data DataType& operator*() const { return *m_iterHelper->data(); } - /** \brief A const pointer dereference function to access the data */ + /// \brief A const pointer dereference function to access the data DataType* operator->() const { return m_iterHelper->data(); } - /** \brief Const accessor for the iterator's current grid point */ + /// \brief Const accessor for the iterator's current grid point GridPt pt() const { return m_iterHelper->pt(); } - /** \brief Equality test against another iterator */ + /// \brief Equality test against another iterator bool operator==(const iter& other) const { - return (m_octLevel == other.m_octLevel) // point to same - // object + return (m_octLevel == other.m_octLevel) // point to same object && m_iterHelper->equal(other.m_iterHelper); } - /** \brief Inequality test against another iterator */ + /// \brief Inequality test against another iterator bool operator!=(const iter& other) const { return !operator==(other); } - /** \brief Increment the iterator to the next point */ + /// \brief Increment the iterator to the next point iter& operator++() { m_iterHelper->increment(); return *this; } - /** \brief Increment the iterator to the next point */ + /// \brief Increment the iterator to the next point iter operator++(int) { iter res = *this; @@ -264,16 +246,12 @@ class OctreeLevel } private: - OctreeLevel* m_octLevel; /** Pointer to the iterator's - container class */ - IterHelper* m_iterHelper; /// Instance of iterator helper class + OctreeLevel* m_octLevel {nullptr}; /// Pointer to the iterator's container class + IterHelper* m_iterHelper {nullptr}; /// Instance of iterator helper class }; public: - /** - * \brief Predicate to check whether the block associated with the - * given GridPt pt is a leaf block - */ + /// \brief Predicate to check whether the block associated with the given GridPt pt is a leaf block bool isLeaf(const GridPt& pt) const { return this->blockStatus(pt) == LeafBlock; } /** @@ -282,49 +260,53 @@ class OctreeLevel */ bool isInternal(const GridPt& pt) const { return this->blockStatus(pt) == InternalBlock; } - /** \brief Begin iterator to points and data in tree level */ + /// \brief Begin iterator to points and data in tree level BlockIter begin() { return BlockIter(this, true); } - /** \brief Const begin iterator to points and data in tree level */ + /// \brief Const begin iterator to points and data in tree level ConstBlockIter begin() const { return ConstBlockIter(this, true); } - /** \brief End iterator to points and data in tree level */ + /// \brief End iterator to points and data in tree level BlockIter end() { return BlockIter(this, false); } - /** \brief Const end iterator to points and data in tree level */ + /// \brief Const end iterator to points and data in tree level ConstBlockIter end() const { return ConstBlockIter(this, false); } - /** \brief Virtual function to check the status of a block + /** + * \brief Virtual function to check the status of a block * (e.g. Leaf, Internal, NotInTree) */ virtual TreeBlockStatus blockStatus(const GridPt& pt) const = 0; - /** \brief Virtual predicate to determine if the OctreeLevel is empty */ + /// \brief Virtual predicate to determine if the OctreeLevel is empty virtual bool empty() const = 0; - /** \brief Virtual predicate to determine if the OctreeLevel has + /** + * \brief Virtual predicate to determine if the OctreeLevel has * a block with the given grid point pt */ virtual bool hasBlock(const GridPt& pt) const = 0; - /** \brief Virtual function to add all children of the given + /** + * \brief Virtual function to add all children of the given * grid point pt to the OctreeLevel */ virtual void addAllChildren(const GridPt& pt) = 0; - /** \brief Virtual const accessor for the data associated with grid point pt - */ + /// \brief Virtual const accessor for the data associated with grid point pt virtual const BlockDataType& operator[](const GridPt& pt) const = 0; - /** \brief Virtual accessor for the data associated with grid point pt */ + /// \brief Virtual accessor for the data associated with grid point pt virtual BlockDataType& operator[](const GridPt& pt) = 0; - /** \brief Virtual accessor for the data associated with + /** + * \brief Virtual accessor for the data associated with * all children of the given grid point (i.e. the brood) */ virtual BroodData& getBroodData(const GridPt& pt) = 0; - /** \brief Virtual const accessor for the data associated + /** + * \brief Virtual const accessor for the data associated * with all children of the given grid point (i.e. the brood) */ virtual const BroodData& getBroodData(const GridPt& pt) const = 0; @@ -344,9 +326,7 @@ class OctreeLevel */ virtual ConstBlockIteratorHelper* getIteratorHelper(bool) const = 0; - /** \brief Virtual function to compute the number of - * blocks (internal and leaf) in the level - */ + /// \brief Virtual function to compute the number of blocks (internal and leaf) in the level virtual int numBlocks() const = 0; /** \brief Virtual function to compute the number of @@ -354,8 +334,7 @@ class OctreeLevel */ virtual int numInternalBlocks() const = 0; - /** \brief Virtual function to compute the number of leaf blocks in the level - */ + /// \brief Virtual function to compute the number of leaf blocks in the level virtual int numLeafBlocks() const = 0; protected: diff --git a/src/axom/spin/SparseOctreeLevel.hpp b/src/axom/spin/SparseOctreeLevel.hpp index 6230cd36d8..043c85d6f9 100644 --- a/src/axom/spin/SparseOctreeLevel.hpp +++ b/src/axom/spin/SparseOctreeLevel.hpp @@ -18,22 +18,14 @@ #include -#if defined(AXOM_USE_SPARSEHASH) - #include "axom/sparsehash/dense_hash_map" -#else - #include -#endif - namespace axom { namespace spin { /** - * \brief Traits class to manage types for different point representations in a - * SparseOctreeLevel + * \brief Traits class to manage types for different point representations in a SparseOctreeLevel * - * The general case is meant for Representations types that are unsigned - * integers + * The general case is meant for Representations types that are unsigned integers * and uses a Morton-based index as the hashmap key. */ template @@ -49,15 +41,11 @@ struct BroodRepresentationTraits "RepresentationType must be unsigned"); // Requires a uint for RepresentationType with 8-,16-,32-, or 64- bits -#if defined(AXOM_USE_SPARSEHASH) - using MapType = axom::google::dense_hash_map; -#else - using MapType = std::unordered_map; -#endif + using MapType = axom::FlatMap; using BroodType = Brood; - /** Simple function to convert a point to its representation type */ + /// Simple function to convert a point to its representation type static PointRepresentationType convertPoint(const GridPt& pt) { return BroodType::MortonizerType::mortonize(pt); @@ -68,25 +56,14 @@ struct BroodRepresentationTraits * * \note sparsehash's maps require setting some default keys */ - static void initializeMap(MapType& map) - { -#if defined(AXOM_USE_SPARSEHASH) - const PointRepresentationType maxVal = axom::numeric_limits::max(); - map.set_empty_key(maxVal); - map.set_deleted_key(maxVal - 1); -#else - AXOM_UNUSED_VAR(map); -#endif - } + static void initializeMap(MapType& map) { AXOM_UNUSED_VAR(map); } }; /** - * \brief Traits class to manage types for different point representations in a - * SparseOctreeLevel + * \brief Traits class to manage types for different point representations in a SparseOctreeLevel * - * This is a specialization meant for point representation - * that use an integer grid point. The underlying hashmap uses a Morton-based - * hash function. + * This is a specialization meant for point representation that use an integer grid point. + * The underlying hashmap uses a Morton-based hash function. */ template struct BroodRepresentationTraits> @@ -97,17 +74,13 @@ struct BroodRepresentationTraits::value, "CoordType must be integral"); -#if defined(AXOM_USE_SPARSEHASH) - using MapType = axom::google::dense_hash_map; -#else - using MapType = std::unordered_map; -#endif + using MapType = axom::FlatMap; using BroodType = Brood; - /** Simple function to convert a point to its representation type - * \note This is a pass through function - * since the representation and grid point types are the same + /** + * \brief Simple function to convert a point to its representation type + * \note This is a pass through function since the representation and grid point types are the same */ static const PointRepresentationType& convertPoint(const GridPt& pt) { @@ -119,19 +92,7 @@ struct BroodRepresentationTraits::max(); - GridPt maxPt(maxCoord); - map.set_empty_key(maxPt); - - maxPt[DIM - 1]--; - map.set_deleted_key(maxPt); -#else - AXOM_UNUSED_VAR(map); -#endif - } + static void initializeMap(MapType& map) { AXOM_UNUSED_VAR(map); } }; /** @@ -188,13 +149,12 @@ class SparseOctreeLevel : public OctreeLevel using BaseBlockItType = ParentType; IteratorHelper(OctreeLevelType* octLevel, bool begin) - : m_offset(0) + : m_currentIter(begin ? octLevel->m_map.begin() : octLevel->m_map.end()) + , m_offset(0) , m_isLevelZero(octLevel->level() == 0) - { - m_currentIter = begin ? octLevel->m_map.begin() : octLevel->m_map.end(); - } + { } - /** Increment to next block in the level */ + /// Increment to next block in the level void increment() { ++m_offset; @@ -206,15 +166,15 @@ class SparseOctreeLevel : public OctreeLevel } } - /** Accessor for point associated with iterator's block */ + /// Accessor for point associated with iterator's block GridPt pt() const { return BroodType::reconstructGridPt(m_currentIter->first, m_offset); } - /** Accessor for data associated with the iterator's block */ + /// Accessor for data associated with the iterator's block BlockDataType* data() { return &m_currentIter->second[m_offset]; } - /** Const accessor for data associated with the iterator's block */ + /// Const accessor for data associated with the iterator's block const BlockDataType* data() const { return &m_currentIter->second[m_offset]; } - /** \brief Predicate to determine if two block iterators are the same */ + /// \brief Predicate to determine if two block iterators are the same bool equal(const BaseBlockItType* other) { const self* pother = dynamic_cast(other); @@ -231,7 +191,7 @@ class SparseOctreeLevel : public OctreeLevel }; public: - /** \brief Default constructor for an octree level */ + /// \brief Default constructor for an octree level SparseOctreeLevel(int level = -1) : Base(level) { BroodTraits::initializeMap(m_map); } /** @@ -243,8 +203,7 @@ class SparseOctreeLevel : public OctreeLevel BaseBlockIteratorHelper* getIteratorHelper(bool begin) { return new IterHelper(this, begin); } /** - * \brief Factory function to return a ConstSparseBlockIterHelper for this - * level + * \brief Factory function to return a ConstSparseBlockIterHelper for this level * * \param begin A boolean to determine if this is to be * a begin (true) or end (false) iterator @@ -280,8 +239,7 @@ class SparseOctreeLevel : public OctreeLevel << pt << " into octree level " << this->m_level << ". Point was out of bounds -- " << "each coordinate must be between 0 and " << this->maxCoord() << "."); - BroodData& bd = getBroodData(pt); // Adds entire brood at once - // (default constructed) + BroodData& bd = getBroodData(pt); // Adds entire brood at once (default constructed) if(this->level() == 0) { for(int j = 1; j < Base::BROOD_SIZE; ++j) @@ -291,14 +249,14 @@ class SparseOctreeLevel : public OctreeLevel } } - /** \brief Accessor for the data associated with pt */ + /// \brief Accessor for the data associated with pt BlockDataType& operator[](const GridPt& pt) { const BroodType brood(pt); return m_map[brood.base()][brood.offset()]; } - /** \brief Const accessor for the data associated with pt */ + /// \brief Const accessor for the data associated with pt const BlockDataType& operator[](const GridPt& pt) const { SLIC_ASSERT_MSG(hasBlock(pt), @@ -310,10 +268,10 @@ class SparseOctreeLevel : public OctreeLevel return blockIt->second[brood.offset()]; } - /** \brief Access the data associated with the entire brood */ + /// \brief Access the data associated with the entire brood BroodData& getBroodData(const GridPt& pt) { return m_map[BroodTraits::convertPoint(pt)]; } - /** \brief Const access to data associated with the entire brood */ + /// \brief Const access to data associated with the entire brood const BroodData& getBroodData(const GridPt& pt) const { SLIC_ASSERT_MSG(hasBlock(pt), @@ -324,10 +282,10 @@ class SparseOctreeLevel : public OctreeLevel return blockIt->second; } - /** \brief Predicate to check if there are any blocks in this octree level */ + /// \brief Predicate to check if there are any blocks in this octree level bool empty() const { return m_map.empty(); } - /** \brief Returns the number of blocks (internal and leaf) in the level */ + /// \brief Returns the number of blocks (internal and leaf) in the level int numBlocks() const { if(empty()) @@ -338,10 +296,10 @@ class SparseOctreeLevel : public OctreeLevel return ((this->m_level == 0) ? 1 : (static_cast(m_map.size()) * Base::BROOD_SIZE)); } - /** \brief Returns the number of internal blocks in the level */ + /// \brief Returns the number of internal blocks in the level int numInternalBlocks() const { return numBlocks() - numLeafBlocks(); } - /** \brief Returns the number of leaf blocks in the level */ + /// \brief Returns the number of leaf blocks in the level int numLeafBlocks() const { if(empty()) @@ -369,8 +327,7 @@ class SparseOctreeLevel : public OctreeLevel * octree block within this octree level * * \param pt The grid point of the block index that we are testing - * \return The status of the grid point pt (e.g. LeafBlock, InternalBlock, - *...) + * \return The status of the grid point pt (e.g. LeafBlock, InternalBlock, ...) */ TreeBlockStatus blockStatus(const GridPt& pt) const { diff --git a/src/axom/spin/SpatialOctree.hpp b/src/axom/spin/SpatialOctree.hpp index 4a916315a2..7b3dab0a34 100644 --- a/src/axom/spin/SpatialOctree.hpp +++ b/src/axom/spin/SpatialOctree.hpp @@ -108,15 +108,12 @@ class SpatialOctree : public OctreeBase * * \param [in] pt The query point in space * \param [in] startingLevel (Optional) starting level for the query - * \pre pt must be in the bounding box of the octree (i.e. - * boundingBox.contains(pt) == true ) + * \pre pt must be in the bounding box of the octree (i.e. boundingBox.contains(pt) == true ) * \note The collection of leaves covers the bounding box, and the interiors - * of the leaves do not - * intersect, so every point in the bounding box should be located in a unique - * leaf block. + * of the leaves do not intersect, so every point in the bounding box + * should be located in a unique leaf block. * \note We are assuming a half-open interval on the bounding boxes. - * \return The block index (i.e. grid point and level) of the leaf block - * containing the query point + * \return The block index (i.e. grid point and level) of the leaf block containing the query point */ BlockIndex findLeafBlock(const SpacePt& pt, int startingLevel = -1) const { @@ -124,8 +121,7 @@ class SpatialOctree : public OctreeBase m_boundingBox.contains(pt), "SpatialOctree::findLeafNode -- Did not find " << pt << " in bounding box " << m_boundingBox); - // Perform binary search on levels to find the leaf block containing the - // point + // Perform binary search on levels to find the leaf block containing the point CoordType minLev = 0; CoordType maxLev = this->maxLeafLevel(); CoordType lev = (startingLevel == -1) ? maxLev >> 1 : startingLevel; @@ -156,14 +152,12 @@ class SpatialOctree : public OctreeBase } /** - * \brief Utility function to find the quantized grid cell at level lev for - * query point pt + * \brief Utility function to find the quantized grid cell at level lev for query point pt * * \param [in] pt The point at which we are querying. * \param [in] lev The level or resolution. * \pre \f$ 0 \le lev < octree.maxLeafLevel() \f$ - * \post Each coordinate of the returned gridPt is in range - * \f$ [0, 2^{lev}) \f$ + * \post Each coordinate of the returned gridPt is in range \f$ [0, 2^{lev}) \f$ * \return The grid point of the block covering this point at this level * \internal * \todo KW: Should this function be protected? Is it generally useful? @@ -194,10 +188,8 @@ class SpatialOctree : public OctreeBase DISABLE_MOVE_AND_ASSIGNMENT(SpatialOctree); protected: - SpaceVectorLevelMap m_deltaLevelMap; // The width of a cell at each - // level or resolution - SpaceVectorLevelMap m_invDeltaLevelMap; // Its inverse is useful for - // quantizing + SpaceVectorLevelMap m_deltaLevelMap; // The width of a cell at each level or resolution + SpaceVectorLevelMap m_invDeltaLevelMap; // Its inverse is useful for quantizing GeometricBoundingBox m_boundingBox; }; diff --git a/src/axom/spin/docs/sphinx/helper.rst b/src/axom/spin/docs/sphinx/helper.rst index 4b043a25ed..83a04188b8 100644 --- a/src/axom/spin/docs/sphinx/helper.rst +++ b/src/axom/spin/docs/sphinx/helper.rst @@ -38,9 +38,9 @@ The ``Mortonizer`` (along with its associated class ``MortonBase``) implements the Morton index, an operation that associates each point in N-D space with a point on a space-filling curve [#f1]_. The ``PointHash`` class adapts the ``Mortonizer`` to provide a hashing functionality for use with -``std::unordered_map`` or similar container classes. +``std::unordered_map``, ``axom::FlatMap`` or similar container classes. -The math of the Morton index works with integers. Thus the ``Mortonizer`` and +The math of the Morton index works with integers. Thus, the ``Mortonizer`` and its dependent class ``PointHash`` will not work with floating point coordinates. The following code example shows how the cells of a ``RectangularLattice``, which have integer coordinates, can be used with a hash table. diff --git a/src/axom/spin/examples/spin_introduction.cpp b/src/axom/spin/examples/spin_introduction.cpp index 50ad8e5f43..804d5ba218 100644 --- a/src/axom/spin/examples/spin_introduction.cpp +++ b/src/axom/spin/examples/spin_introduction.cpp @@ -159,8 +159,7 @@ void demoMorton() // Make the map from grid point to DataContainer MapType map; - // For several query points, create a DataContainer if necessary and register - // the point. + // For several query points, create a DataContainer if necessary and register the point. std::vector pts = generatePoints(); for(RLSpacePt p : pts) { @@ -317,8 +316,7 @@ void findNeighborCandidates(TriangleType& t1, int i, UniformGridType* ugrid, std const std::vector bToCheck = ugrid->getBinsForBbox(bbox); size_t checkcount = bToCheck.size(); - // Load all the triangles in these bins whose indices are - // greater than i into a vector. + // Load all the triangles in these bins whose indices are greater than i into a vector. for(size_t curb = 0; curb < checkcount; ++curb) { axom::ArrayView ntlist = ugrid->getBinContents(bToCheck[curb]); @@ -389,8 +387,7 @@ void showImplicitGrid() const int numElts = static_cast(tris.size()); IGridT grid(bbox, &res, numElts); - // load the bounding box of each triangle, along with its index, - // into the ImplicitGrid. + // load the bounding box of each triangle, along with its index, into the ImplicitGrid. for(int i = 0; i < numElts; ++i) { grid.insert(findBbox(tris[i]), i); @@ -599,7 +596,7 @@ void driveBVHTree() { std::vector tris; makeTreeTriangles(tris); - Point2DType ppoint = Point2DType::make_point(0.45, 0.25); + Point2DType ppoint {0.45, 0.25}; std::vector intersections, candidates; BVH2DType* tree = buildBVHTree(tris); diff --git a/src/axom/spin/internal/linear_bvh/build_radix_tree.hpp b/src/axom/spin/internal/linear_bvh/build_radix_tree.hpp index 407c43c177..756e1ef3f9 100644 --- a/src/axom/spin/internal/linear_bvh/build_radix_tree.hpp +++ b/src/axom/spin/internal/linear_bvh/build_radix_tree.hpp @@ -42,44 +42,42 @@ namespace internal namespace linear_bvh { //------------------------------------------------------------------------------ -//Returns 30 bit morton code for coordinates for -// x, y, and z are expecting to be between [0,1] +//Returns 30 bit morton code for coordinates point is expected to be between [0,1] template -static inline AXOM_HOST_DEVICE std::int32_t morton32_encode(const primal::Vector& point) +static inline AXOM_HOST_DEVICE std::uint32_t morton32_encode(const primal::Vector& point) { using PointType = primal::Point; + using axom::utilities::clampVal; //for a float, take the first 10 bits. Note, 2^10 = 1024 constexpr int NUM_BITS_PER_DIM = 32 / Dims; constexpr FloatType FLOAT_TO_INT = 1 << NUM_BITS_PER_DIM; constexpr FloatType FLOAT_CEILING = FLOAT_TO_INT - 1; + constexpr FloatType FLOAT_ZERO {0.}; std::int32_t int_coords[Dims]; for(int i = 0; i < Dims; i++) { - int_coords[i] = static_cast( - fmin(fmax(point[i] * FLOAT_TO_INT, static_cast(0)), FLOAT_CEILING)); + int_coords[i] = + static_cast(clampVal(point[i] * FLOAT_TO_INT, FLOAT_ZERO, FLOAT_CEILING)); } - return static_cast(convertPointToMorton(PointType(int_coords))); + return convertPointToMorton(PointType(int_coords)); } //------------------------------------------------------------------------------ -//Returns 30 bit morton code for coordinates for -//coordinates in the unit cude -static inline AXOM_HOST_DEVICE std::int64_t morton64_encode(axom::float32 x, - axom::float32 y, - axom::float32 z = 0.0) +//Returns 63 bit morton code for a 3D point with coordinates in the unit cube +static inline AXOM_HOST_DEVICE std::uint64_t morton64_encode(float x, float y, float z = 0.0) { using PointType = primal::Point; + using axom::utilities::clampVal; //take the first 21 bits. Note, 2^21= 2097152.0f constexpr float F_21 = 2097152.0f; - const PointType integer_pt = {static_cast(fmin(fmax(x * F_21, 0.0f), F_21)), - static_cast(fmin(fmax(y * F_21, 0.0f), F_21)), - static_cast(fmin(fmax(z * F_21, 0.0f), F_21))}; - - return static_cast(convertPointToMorton(integer_pt)); + return convertPointToMorton( + PointType {static_cast(clampVal(x * F_21, 0.f, F_21)), + static_cast(clampVal(y * F_21, 0.f, F_21)), + static_cast(clampVal(z * F_21, 0.f, F_21))}); } template @@ -293,7 +291,7 @@ void build_tree(RadixTree& data) // http://research.nvidia.com/sites/default/files/publications/karras2012hpg_paper.pdf // Pointers and vars are redeclared because I have a faint memory - // of a huge amount of pain and suffering due so cuda + // of a huge amount of pain and suffering due to cuda // lambda captures of pointers inside a struct. Bad memories // of random segfaults ........ be warned const std::int32_t inner_size = data.m_inner_size;