From e35263bda007b0b9142e58c340bfb63b2bc524f5 Mon Sep 17 00:00:00 2001 From: Bill Dolinar Date: Sat, 16 Dec 2023 20:06:49 -0700 Subject: [PATCH] Speedup for adjacent face calculation. --- xmsgrid/ugrid/XmUGrid.cpp | 3548 +++++++++++++++++++++---------------- 1 file changed, 1997 insertions(+), 1551 deletions(-) diff --git a/xmsgrid/ugrid/XmUGrid.cpp b/xmsgrid/ugrid/XmUGrid.cpp index 70b8d330..c937ec37 100644 --- a/xmsgrid/ugrid/XmUGrid.cpp +++ b/xmsgrid/ugrid/XmUGrid.cpp @@ -10,7 +10,6 @@ // 1. Precompiled header // 2. My own header -#include "cxxtest/TestSuite.h" #include // 3. Standard library headers @@ -45,6 +44,135 @@ /// XMS Namespace namespace xms { +namespace +{ +//////////////////////////////////////////////////////////////////////////////// +/// \class Span +/// \brief A class that provides a lightweight, read-only reference to a +/// contiguous sequence (an array, part of an array, or a vector), much like +/// std::span in C++20. +/// +/// This class does not own the data it points to, it merely provides a way to +/// access it. +//////////////////////////////////////////////////////////////////////////////// +class Span +{ +public: + //--------------------------------------------------------------------------- + /// \brief Constructs an empty Span. + //--------------------------------------------------------------------------- + Span() : m_data(nullptr), m_size(0) {} + + //--------------------------------------------------------------------------- + /// \brief Constructs a Span from a pointer to data and a size. + /// \param[in] a_data A pointer to the data. + /// \param[in] a_size The number of elements in the data. + //--------------------------------------------------------------------------- + Span(const int* a_data, int a_size) + : m_data(a_data) + , m_size(a_size) + { + } + + //--------------------------------------------------------------------------- + /// \brief Accesses an element of the data. + /// \param[in] a_idx The index of the element to access. + /// \return The element at the specified index. + //--------------------------------------------------------------------------- + int operator[](int a_idx) { return m_data[a_idx]; } + + //--------------------------------------------------------------------------- + /// \brief Accesses an element of the data. + /// \param[in] a_idx The index of the element to access. + /// \return The element at the specified index. + //--------------------------------------------------------------------------- + int operator[](int a_idx) const { return m_data[a_idx]; } + + //--------------------------------------------------------------------------- + /// \brief Returns the size of the span. + /// \return The number of elements in the span. + //--------------------------------------------------------------------------- + [[nodiscard]] int size() const { return m_size; } + + //--------------------------------------------------------------------------- + /// \brief Checks if the span is empty. + /// \return True if the span is empty, false otherwise. + //--------------------------------------------------------------------------- + [[nodiscard]] bool empty() const { return m_size == 0; } + + //--------------------------------------------------------------------------- + /// \brief Returns a pointer to the beginning of the span. + /// \return A pointer to the first element of the span. + //--------------------------------------------------------------------------- + const int* begin() { return m_data; } + + //--------------------------------------------------------------------------- + /// \brief Returns a pointer to the end of the span. + /// \return A pointer to the element following the last element of the span. + //--------------------------------------------------------------------------- + const int* end() { return m_data + m_size; } + + //--------------------------------------------------------------------------- + /// \brief Returns a pointer to the data the span refers to. + /// \return A pointer to the data. + //--------------------------------------------------------------------------- + const int* data() { return m_data; } + +private: + const int* m_data; ///< A pointer to the data. + int m_size; ///< The number of elements in the data. +}; + +//////////////////////////////////////////////////////////////////////////////// +/// \class CellStreamCache +/// \brief A class that provides methods for caching values related to an +/// XmUGrid cell stream. +//////////////////////////////////////////////////////////////////////////////// +class CellStreamCache +{ +public: + void Clear(); + void Update(); + + + // Points + int GetPointCount() const; + const VecPt3d& GetLocations() const; + void SetLocations(const VecPt3d& a_locations); + void GetExtents(Pt3d& a_min, Pt3d& a_max) const; + Pt3d GetPointLocation(int a_pointIdx) const; + bool SetPointLocation(int a_pointIdx, const Pt3d& a_location); + void GetPointAdjacentCells(int a_pointIdx, VecInt& a_adjacentCells) const; + void GetPointAdjacentCells(int a_pointIdx, Span& a_adjacentCells) const; + int GetPointAdjacentCellCount(int a_pointIdx) const; + int GetPointAdjacentCell(int a_pointIdx, int a_cellIdx) const; + + // Cells + int GetCellCount() const; + XmUGridCellType GetCellType(int a_cellIdx) const; + const VecInt& GetCellstream() const; + bool SetCellstream(const VecInt& a_cellstream); + bool GetCellCellstream(int a_cellIdx, VecInt& a_cellstream) const; + Span GetCellCellstream(int a_cellIdx) const; + int GetCellCellstreamIndex(int a_cellIdx) const; + int GetNumberOfItemsForCell(int a_cellIdx) const; + int GetCellDimension(int a_cellIdx) const; + std::vector GetDimensionCounts() const; + +private: + void UpdateCellLinks(); + void UpdatePointLinks(); + + VecPt3d m_locations; ///< UGrid point locations + VecInt m_cellstream; ///< UGrid cell stream. @see SetCellstream, GetCellstream + VecInt m_cellIdxToStreamIdx; ///< Indexes for each cell in the cell stream + VecInt m_pointsToCells; ///< Array of points cells (goes from pointIdx to list of cells) + VecInt m_pointIdxToPointsToCells; ///< Indexes for each point in array of + mutable VecInt m_cellDimensionCounts; ///< Cache for cell dimension counts + ///< points cells +}; +} // namespace + //////////////////////////////////////////////////////////////////////////////// /// Implementation for XmUGrid class XmUGrid::Impl @@ -82,6 +210,7 @@ class XmUGrid::Impl int GetPointAdjacentCellCount(int a_pointIdx) const; VecInt GetPointAdjacentCells(int a_pointIdx) const; // cells associated with point void GetPointAdjacentCells(int a_pointIdx, VecInt& a_adjacentCells) const; + void GetPointAdjacentCells(int a_pointIdx, Span& a_adjacentCells) const; VecInt GetPointsAdjacentCells(const VecInt& a_points) const; void GetPointsAdjacentCells(const int* a_pointIdxs, @@ -139,13 +268,25 @@ class XmUGrid::Impl VecInt GetCell3dFacePoints(int a_cellIdx, int a_faceIdx) const; void GetCell3dFacePoints(int a_cellIdx, int a_faceIdx, VecInt& a_facePtIdxs) const; + void GetCell3dFacePoints(int a_cellIdx, int a_faceIdx, Span& a_facePtIdxs) const; + void CopyCell3dFacePoints(int a_cellIdx, int a_faceIdx, int* a_facePtIdxs) const; VecInt2d GetCell3dFacesPoints(int a_cellIdx) const; int GetCell3dFaceAdjacentCell(int a_cellIdx, int a_faceIdx) const; + int GetCell3dFaceAdjacentCellWithCache(int a_cellIdx, int a_faceIdx) const; + int GetCell3dFaceAdjacentCellNoCache(int a_cellIdx, int a_faceIdx) const; bool GetCell3dFaceAdjacentCell(int a_cellIdx, int a_faceIdx, int& a_neighborCell, int& a_neighborFace) const; + bool GetCell3dFaceAdjacentCellWithCache(int a_cellIdx, + int a_faceIdx, + int& a_neighborCell, + int& a_neighborFace) const; + bool GetCell3dFaceAdjacentCellNoCache(int a_cellIdx, + int a_faceIdx, + int& a_neighborCell, + int& a_neighborFace) const; XmUGridFaceOrientation GetCell3dFaceOrientation(int a_cellIdx, int a_faceIdx) const; XmUGridFaceOrientation FaceOrientation(int a_cellIdx, int a_faceIdx) const; @@ -155,10 +296,6 @@ class XmUGrid::Impl XmUGridFaceOrientation GetOrientationFromArea(int a_cellIdx, int a_faceIdx) const; private: - void UpdateLinks(); // Calls UpdateCellLinks & UpdatePointLinks - void UpdateCellLinks(); - void UpdatePointLinks(); - void SetModified(); bool IsCellValidWithPointChange(int a_cellIdx, @@ -166,43 +303,24 @@ class XmUGrid::Impl const Pt3d& a_newPosition) const; bool IsValidCellIdx(int a_cellIdx) const; - static int DimensionFromCellType(XmUGridCellType a_cellType); - - int GetNumberOfItemsForCell(int a_cellIdx) const; - // Optimized for efficiency - void GetCellCellstream(int a_cellIdx, const int** a_start, int& a_length) const; + Span GetCellCellstream(int a_cellIdx) const; int GetNumberOfPolyhedronEdges(int a_cellIdx) const; - static void GetUniquePointsFromPolyhedronCellstream(const VecInt& a_cellstream, - int a_numCellItems, - int& a_currIdx, - VecInt& a_uniqueGetCellPoints, - VecInt& a_pointLastUsedIdx); - static bool GetUniquePointsFromPolyhedronSingleCellstream(const VecInt& a_cellstream, - VecInt& a_cellPoints); - static void GetUniqueEdgesFromPolyhedronCellstream( - const int* a_start, - int& a_length, - boost::container::flat_set& a_cellEdges, - int& a_currIdx); bool GetPlanViewPolygon2d(int a_cellIdx, VecPt3d& a_polygon) const; bool GetPlanViewPolygon3d(int a_cellIdx, VecPt3d& a_polygon) const; bool IsFaceSide(const VecInt& a_facePts) const; // plan view bool GetCellXySegments(int cellIdx, VecPt3d& a_segments) const; // plan view - void GetEdgesOfFace(const VecInt& a_face, std::vector& a_edges) const; bool DoEdgesCrossWithPointChange(int a_changedPtIdx, const Pt3d& a_newPosition, const std::vector& a_edges) const; - void GetExtentsFromPoints(const VecPt3d& a_locations, Pt3d& a_min, Pt3d& a_max) const; bool GetFaceXySegments(int a_cellIdx, int a_faceIdx, VecPt3d& a_segments) const; // plan view - void CalculateCacheValues() const; + void SetupFaceCacheValues() const; void ClearCacheValues(); int GetCell3dFaceCountNoCache(int a_cellIdx) const; - int GetCell3dFaceAdjacentCellNoCache(int a_cellIdx, int a_faceIdx) const; XmUGridFaceOrientation GetCell3dFaceOrientationNoCache(int a_cellIdx, int a_faceIdx) const; bool GetNextFaceColumn(const VecInt& a_facePoints, int a_starti, @@ -220,22 +338,20 @@ class XmUGrid::Impl NEEDS_CALCULATION = -2 ///< Cached value needs to be calculated }; - VecPt3d m_locations; ///< UGrid point locations - VecInt m_cellstream; ///< UGrid cell stream. @see SetCellstream, GetCellstream - VecInt m_cellIdxToStreamIdx; ///< Indexes for each cell in the cell stream - VecInt m_pointsToCells; ///< Array of points cells (goes from pointIdx to list - ///< of cells) - VecInt m_pointIdxToPointsToCells; ///< Indexes for each point in array of - ///< points cells - bool m_modified = false; ///< Has UGrid been modified since last SetUnmodified call? - bool m_useCache = true; ///< Are we using caching for some calls? - mutable VecInt m_numberOfFaces; ///< Cache for number of cell faces - mutable VecInt m_cellFaceOffset; ///< Cache for offset to m_faceOrientation and m_faceNeighbor - mutable VecInt m_faceOrientation; ///< For vertically prismatic cell is face top, side, bottom - mutable VecInt m_faceNeighbor; ///< Cache for Face neighbor - mutable VecInt m_cellDimensionCounts; ///< Cache for cell dimension counts - XmUGridCellOrdering m_cellOrdering = XMU_CELL_ORDER_INCREASING_DOWN; ///< Cell ordering. When known speeds - ///< up face orientation calculation. + CellStreamCache m_cellStreamCache; ///< Cache for navigating between gridstream points and cells + bool m_modified = false; ///< Has UGrid been modified since last SetUnmodified call? + bool m_useCache = true; ///< Are we using caching for some calls? + mutable VecInt m_numberOfFaces; ///< Cache for number of cell faces (goes from cell index to + ///< number of faces) + mutable VecInt m_cellFaceOffset; ///< Cache for offset to m_faceOrientation, m_faceNeighborCell, + ///< and m_facePointOffset + mutable VecInt m_faceOrientation; ///< For vertically prismatic cell is face top, side, bottom + mutable VecInt m_faceNeighborCell; ///< Cache for Face neighbor cell + mutable VecInt m_facePointOffset; ///< Cache for offset to m_facePoints + mutable VecInt m_facePoints; ///< Cache for Face points + XmUGridCellOrdering m_cellOrdering = + XMU_CELL_ORDER_INCREASING_DOWN; ///< Cell ordering. When known speeds + ///< up face orientation calculation. }; //----- Internal functions ----------------------------------------------------- @@ -349,10 +465,8 @@ const VecEdge& iGetEdgeOffsetTable(int a_cellType) case XMU_PARAMETRIC_TETRA_REGION: case XMU_PARAMETRIC_HEX_REGION: default: - XM_ASSERT(0); - return fg_empty; + break; } - XM_ASSERT(0); return fg_empty; } // iGetEdgeOffsetTable @@ -444,8 +558,7 @@ const VecInt2d& iGetFaceOffsetTable(int a_cellType) case XMU_PARAMETRIC_TETRA_REGION: case XMU_PARAMETRIC_HEX_REGION: default: - XM_ASSERT(0); - return fg_empty; + break; } XM_ASSERT(0); @@ -591,6 +704,129 @@ void iMergeSegmentsToPoly(const VecPt3d& a_segments, VecPt3d& a_polygon) a_polygon.clear(); } } // iMergeSegmentsToPoly +//------------------------------------------------------------------------------ +/// \brief Get the dimension given the cell type (0d, 1d, 2d, or 3d). +/// \param[in] a_cellType the cell type +/// \return the dimension of the cell type +//------------------------------------------------------------------------------ +int iDimensionFromCellType(XmUGridCellType a_cellType) +{ + switch (a_cellType) + { + // invalid + case XMU_INVALID_CELL_TYPE: + return -1; + + // 0D + case XMU_EMPTY_CELL: + case XMU_VERTEX: + case XMU_POLY_VERTEX: + case XMU_CONVEX_POINT_SET: // Special class of cells formed by convex group of + // points + return 0; + + // 1D + case XMU_LINE: + case XMU_POLY_LINE: + + case XMU_QUADRATIC_EDGE: + case XMU_PARAMETRIC_CURVE: + case XMU_HIGHER_ORDER_EDGE: + case XMU_CUBIC_LINE: // Cubic, isoparametric cell + return 1; + // 2D + case XMU_TRIANGLE: + case XMU_TRIANGLE_STRIP: + case XMU_POLYGON: + case XMU_PIXEL: + case XMU_QUAD: + + case XMU_QUADRATIC_TRIANGLE: + case XMU_QUADRATIC_QUAD: + case XMU_BIQUADRATIC_QUAD: + case XMU_QUADRATIC_LINEAR_QUAD: + case XMU_BIQUADRATIC_TRIANGLE: + case XMU_QUADRATIC_POLYGON: + + case XMU_PARAMETRIC_SURFACE: + case XMU_PARAMETRIC_TRI_SURFACE: + case XMU_PARAMETRIC_QUAD_SURFACE: + + case XMU_HIGHER_ORDER_TRIANGLE: + case XMU_HIGHER_ORDER_QUAD: + case XMU_HIGHER_ORDER_POLYGON: + + return 2; + + // assuming the rest are 3D (see below) + default: + break; +#if 0 // Remaining definitions + XMU_TETRA = 10, + XMU_VOXEL = 11, + XMU_HEXAHEDRON = 12, + XMU_WEDGE = 13, + XMU_PYRAMID = 14, + XMU_PENTAGONAL_PRISM = 15, + XMU_HEXAGONAL_PRISM = 16, + + // Quadratic, isoparametric cells + XMU_QUADRATIC_TETRA = 24, + XMU_QUADRATIC_HEXAHEDRON = 25, + XMU_QUADRATIC_WEDGE = 26, + XMU_QUADRATIC_PYRAMID = 27, + XMU_TRIQUADRATIC_HEXAHEDRON = 29, + XMU_QUADRATIC_LINEAR_WEDGE = 31, + XMU_BIQUADRATIC_QUADRATIC_WEDGE = 32, + XMU_BIQUADRATIC_QUADRATIC_HEXAHEDRON = 33, + + // Polyhedron cell (consisting of polygonal faces) + XMU_POLYHEDRON = 42, + + // Higher order cells in parametric form + XMU_PARAMETRIC_TETRA_REGION = 55, + XMU_PARAMETRIC_HEX_REGION = 56, + + // Higher order cells + XMU_HIGHER_ORDER_TETRAHEDRON = 64, + XMU_HIGHER_ORDER_WEDGE = 65, + XMU_HIGHER_ORDER_PYRAMID = 66, + XMU_HIGHER_ORDER_HEXAHEDRON = 67, +#endif + } + return 3; + +} // iDimensionFromCellType +//------------------------------------------------------------------------------ +/// \brief Function to get the extents from a list of points. +/// \param[in] a_locations The point locations to get the extents of. +/// \param[out] a_min Minimum point location. +/// \param[out] a_max Maximum point location. +//------------------------------------------------------------------------------ +void iGetExtentsFromPoints(const VecPt3d& a_locations, Pt3d& a_min, Pt3d& a_max) +{ + a_min.x = a_min.y = a_min.z = xms::XM_DBL_HIGHEST; + a_max.x = a_max.y = a_max.z = xms::XM_DBL_LOWEST; + for (const Pt3d& pt : a_locations) + { + gmAddToExtents(pt, a_min, a_max); + } +} // iGetExtentsFromPoints +//------------------------------------------------------------------------------ +/// \brief Get the edges of a cell given a face +/// \param[in] a_face a vector of point indices of a face +/// \param[out] a_edges a vector of point indices of an edge +//------------------------------------------------------------------------------ +void iGetEdgesOfFace(const VecInt& a_face, std::vector& a_edges) +{ + a_edges.reserve(a_face.size()); + for (int i = 1; i < a_face.size(); ++i) + { + a_edges.emplace_back(a_face[i - 1], a_face[i]); + } + if (a_edges.size() > 1) + a_edges.emplace_back(a_face[a_face.size() - 1], a_face[0]); +} // iGetEdgesOfFace //------------------------------------------------------------------------------ /// \brief Check if a polygon cell in a stream is valid. @@ -706,8 +942,6 @@ int iValidatePolyhedronCell(const VecInt& a_cellStream, int a_points, int a_curr //------------------------------------------------------------------------------ int iValidateCell(const VecInt& a_cellStream, int a_points, int a_currIdx) { - int oldIdx = a_currIdx; - if (a_currIdx + 1 >= a_cellStream.size()) { return BAD_CELL_STREAM; // Cell type is present, but number of points missing @@ -751,412 +985,450 @@ int iValidateCell(const VecInt& a_cellStream, int a_points, int a_currIdx) default: return iValidatePolygonCell(a_cellStream, a_points, a_currIdx, a_cellStream[a_currIdx], false); } - - XM_ASSERT(false); // This should be unreachable. - return BAD_CELL_STREAM; } // iValidateCell +//------------------------------------------------------------------------------ +/// \brief Check a cell stream to make sure it's valid. Compares cell type +/// against expected number of points. Also ensures cells only refer to +/// valid points. +/// \param[in] a_cellstream the cell stream to check +/// \param a_points: The number of points in this grid. Used to ensure the cells +/// only refer to valid points. +/// \return Whether the cell stream is valid. +//------------------------------------------------------------------------------ +bool iIsValidCellstream(const VecInt& a_cellstream, int a_points) +{ + if (a_cellstream.empty()) + { + return true; + } -} // namespace + int currIdx = 0; + while (currIdx < a_cellstream.size()) + { + int step = iValidateCell(a_cellstream, a_points, currIdx); + if (step == BAD_CELL_STREAM) + { + return false; + } + else + { + currIdx += step; + } + } -//////////////////////////////////////////////////////////////////////////////// -/// \class XmUGrid::Impl -/// \brief Implementation for XmUGrid which provides geometry for an -/// unstructured grid. -//////////////////////////////////////////////////////////////////////////////// + return true; +} // iIsValidCellstream //------------------------------------------------------------------------------ -/// \brief Returns the modified flag. Gets set when points or cells get changed. -/// \return the modified flag +/// \brief Get the unique points in a flat set +/// \param[in] a_cellstream the UGrid cell stream +/// \param[in] a_numCellItems the number of cell faces in the polyhedron +/// \param[in] a_currIdx the index of the cell stream; this should reference +/// the number of points in the first face. This variable will be updated +/// to the cell type of the next cell. +/// \param[out] a_uniqueGetCellPoints the unique points of the polyhedron +/// \param[out] a_pointLastUsedIdx a list of integers representing when a point +/// was last found while reading the stream. +/// \note: This function does NOT verify cellstream size!! This function +/// needs to be efficient! //------------------------------------------------------------------------------ -bool XmUGrid::Impl::GetModified() const +void iGetUniquePointsFromPolyhedronCellstream(const VecInt& a_cellstream, + int a_numCellItems, + int& a_currIdx, + VecInt& a_uniqueGetCellPoints, + VecInt& a_pointLastUsedIdx) { - return m_modified; -} // XmUGrid::Impl::GetModified + int stable = a_currIdx; + a_uniqueGetCellPoints.clear(); + int numFaces = a_numCellItems; + for (int faceIdx = 0; faceIdx < numFaces; ++faceIdx) + { + int numFacePoints = a_cellstream[a_currIdx++]; + for (int ptIdx = 0; ptIdx < numFacePoints; ++ptIdx) + { + int pt = a_cellstream[a_currIdx++]; + if (a_pointLastUsedIdx[pt] < stable) + { + a_pointLastUsedIdx[pt] = stable; + a_uniqueGetCellPoints.push_back(pt); + } + } + } +} // iGetUniquePointsFromPolyhedronCellstream //------------------------------------------------------------------------------ -/// \brief Resets the modified flag to false. +/// \brief Get the unique points for cell stream of a single polyhedron cell. +/// \param[in] a_cellstream a single cell stream that is a polyhedron type +/// \param[out] a_cellPoints the points of the cell +/// \return false if invalid stream //------------------------------------------------------------------------------ -void XmUGrid::Impl::SetUnmodified() +bool iGetUniquePointsFromPolyhedronSingleCellstream(const VecInt& a_cellstream, + VecInt& a_cellPoints) { - m_modified = false; -} // XmUGrid::Impl::SetUnmodified -//------------------------------------------------------------------------------ -/// \brief Sets the modified flag to true. + a_cellPoints.clear(); + if (a_cellstream.size() < 2) + { + return false; + } + + auto currItem = a_cellstream.begin(); + int cellType = *currItem++; + if (cellType != XMU_POLYHEDRON) + { + return false; + } + + int numFaces = *currItem++; + for (int faceIdx = 0; faceIdx < numFaces; ++faceIdx) + { + int numFacePoints = *currItem++; + for (int ptIdx = 0; ptIdx < numFacePoints; ++ptIdx) + { + int pt = *currItem++; + if (std::find(a_cellPoints.begin(), a_cellPoints.end(), pt) == a_cellPoints.end()) + { + a_cellPoints.push_back(pt); + } + } + } + + return true; +} // iGetUniquePointsFromPolyhedronSingleCellstream //------------------------------------------------------------------------------ -void XmUGrid::Impl::SetModified() +/// \brief Get the unique edges in a flat set for a given polyhedron +/// \param[in] a_cellstream the UGrid cell stream (integer pointer) +/// \param[in] a_cellEdges Unique cell edges of the polyhedron +/// \param[in] a_currIdx the index of the cell stream; this should reference +/// the number of points in the first face. This variable will be updated +/// to the cell type of the next cell. +/// \note: This function does NOT verify cellstream size!! This function +/// needs to be efficient! +//------------------------------------------------------------------------------ +void iGetUniqueEdgesFromPolyhedronCellstream( + Span a_cellstream, + boost::container::flat_set& a_cellEdges, + int& a_currIdx) { - ClearCacheValues(); - m_modified = true; -} // XmUGrid::Impl::SetModified + int numFaces = a_cellstream[1]; + for (int i(0); i < numFaces; i++) + { + int numPoints = a_cellstream[a_currIdx++]; + for (int pointIdx = 0; pointIdx < numPoints; ++pointIdx) + { + int pt1Idx = pointIdx; + int pt2Idx = (pointIdx + 1) % numPoints; + // The % operator will reset the index back 0 so the final loop will + // provide the last and first point + int pt1 = a_cellstream[a_currIdx + pt1Idx]; + int pt2 = a_cellstream[a_currIdx + pt2Idx]; + // We want unique edges, so we add the lower point index first + if (pt1 < pt2) + a_cellEdges.insert(XmEdge(pt1, pt2)); + else + a_cellEdges.insert(XmEdge(pt2, pt1)); + } + a_currIdx += numPoints; + } +} // iGetUniqueEdgesFromPolyhedronCellstream +//////////////////////////////////////////////////////////////////////////////// +/// \class CellStreamCache +/// \brief Cache for a UGrid's cell stream. +/// \see XmUGrid +//////////////////////////////////////////////////////////////////////////////// //------------------------------------------------------------------------------ -/// \brief Turn on or off use of caching. -/// \param a_useCache Flag to determine if caching will be used. +/// \brief Clear the cache. //------------------------------------------------------------------------------ -void XmUGrid::Impl::SetUseCache(bool a_useCache) +void CellStreamCache::Clear() { - m_useCache = a_useCache; -} // XmUGrid::Impl::SetUseCache -// Points + m_cellIdxToStreamIdx.clear(); + m_pointsToCells.clear(); + m_pointIdxToPointsToCells.clear(); + m_cellDimensionCounts.clear(); +} // CellStreamCache::Clear +//------------------------------------------------------------------------------ +/// \brief Update internal links to navigate between associated points and +/// cells. +//------------------------------------------------------------------------------ +void CellStreamCache::Update() +{ + Clear(); + UpdateCellLinks(); + UpdatePointLinks(); +} // CellStreamCache::Update //------------------------------------------------------------------------------ /// \brief Get the number of points. /// \return the number of points //------------------------------------------------------------------------------ -int XmUGrid::Impl::GetPointCount() const +int CellStreamCache::GetPointCount() const { return (int)m_locations.size(); -} // XmUGrid::Impl::GetPointCount +} // CellStreamCache::GetPointCount //------------------------------------------------------------------------------ /// \brief Get vector of UGrid points. /// \return constant reference to vector of points //------------------------------------------------------------------------------ -const VecPt3d& XmUGrid::Impl::GetLocations() const +const VecPt3d& CellStreamCache::GetLocations() const { return m_locations; -} // XmUGrid::Impl::GetLocations +} // CellStreamCache::GetLocations //------------------------------------------------------------------------------ /// \brief Set UGrid points. /// \param[in] a_locations vector of point locations //------------------------------------------------------------------------------ -void XmUGrid::Impl::SetLocations(const VecPt3d& a_locations) +void CellStreamCache::SetLocations(const VecPt3d& a_locations) { m_locations = a_locations; - SetModified(); -} // XmUGrid::Impl::SetLocations +} // CellStreamCache::SetLocations +//------------------------------------------------------------------------------ +/// \brief Get extents of all points in UGrid +/// \param[out] a_min minimum extent of all points +/// \param[out] a_max maximum extent of all points +//------------------------------------------------------------------------------ +void CellStreamCache::GetExtents(Pt3d& a_min, Pt3d& a_max) const +{ + iGetExtentsFromPoints(m_locations, a_min, a_max); +} // CellStreamCache::GetExtents //------------------------------------------------------------------------------ /// \brief Get the point /// \param[in] a_pointIdx the index of the point /// \return the point or an initialize point if the index is invalid //------------------------------------------------------------------------------ -Pt3d XmUGrid::Impl::GetPointLocation(int a_pointIdx) const +Pt3d CellStreamCache::GetPointLocation(int a_pointIdx) const { if (a_pointIdx >= 0 && a_pointIdx < m_locations.size()) return m_locations[a_pointIdx]; return {}; -} // XmUGrid::Impl::GetPointLocation +} // CellStreamCache::GetPointLocation //------------------------------------------------------------------------------ /// \brief Set the point /// \param[in] a_pointIdx the index of the point /// \param[in] a_location The new location of the specified point /// \return whether the point was successfully set //------------------------------------------------------------------------------ -bool XmUGrid::Impl::SetPointLocation(int a_pointIdx, const Pt3d& a_location) +bool CellStreamCache::SetPointLocation(int a_pointIdx, const Pt3d& a_location) { if (a_pointIdx >= 0 && a_pointIdx < m_locations.size()) { m_locations[a_pointIdx] = a_location; - SetModified(); return true; } return false; -} // XmUGrid::Impl::SetPointLocation - -//------------------------------------------------------------------------------ -/// \brief Get the X, Y location of a point. -/// \param[in] a_pointIdx The index of the point. -/// \return The location of the point with Z set to 0.0. -//------------------------------------------------------------------------------ -Pt3d XmUGrid::Impl::GetPointXy0(int a_pointIdx) const -{ - Pt3d pt = GetPointLocation(a_pointIdx); - pt.z = 0.0; - return pt; -} // XmUGrid::Impl::GetPointXy0 - +} // CellStreamCache::SetPointLocation //------------------------------------------------------------------------------ -/// \brief Convert a vector of point indices into a vector of point 3d -/// \param[in] a_points a vector of point indices -/// \return vector of point 3d +/// \brief Get the cells that are associated with the specified point +/// \param[in] a_pointIdx the index of the point +/// \param[out] a_adjacentCells a vector of the adjacent cell indices +/// \return a vector of the cell indices associated with this point //------------------------------------------------------------------------------ -VecPt3d XmUGrid::Impl::GetPointsLocations(const VecInt& a_points) const +void CellStreamCache::GetPointAdjacentCells(int a_pointIdx, VecInt& a_adjacentCells) const { - VecPt3d point3d; - for (auto point : a_points) + a_adjacentCells.clear(); + int numCells = GetPointAdjacentCellCount(a_pointIdx); + for (int cellIdx = 0; cellIdx < numCells; cellIdx++) { - point3d.push_back(GetPointLocation(point)); + a_adjacentCells.push_back(m_pointsToCells[m_pointIdxToPointsToCells[a_pointIdx] + cellIdx + 1]); } - return point3d; -} // XmUGrid::Impl::GetPointsLocations - +} // CellStreamCache::GetPointAdjacentCells //------------------------------------------------------------------------------ -/// \brief Get extents of all points in UGrid -/// \param[out] a_min minimum extent of all points -/// \param[out] a_max maximum extent of all points +/// \brief Get the cells that are associated with the specified point +/// \param[in] a_pointIdx The index of the point. +/// \return The cell indices associated with this point. //------------------------------------------------------------------------------ -void XmUGrid::Impl::GetExtents(Pt3d& a_min, Pt3d& a_max) const +void CellStreamCache::GetPointAdjacentCells(int a_pointIdx, Span& a_adjacentCells) const { - GetExtentsFromPoints(m_locations, a_min, a_max); -} // XmUGrid::Impl::GetExtents - + const int* cells = &m_pointsToCells[m_pointIdxToPointsToCells[a_pointIdx] + 1]; + int numCells = GetPointAdjacentCellCount(a_pointIdx); + a_adjacentCells = Span(cells, numCells); +} // CellStreamCache::GetPointAdjacentCells //------------------------------------------------------------------------------ /// \brief Get the number of cells that use a point. /// \param a_pointIdx The point to check. /// \return The number of cells that use the point. //------------------------------------------------------------------------------ -int XmUGrid::Impl::GetPointAdjacentCellCount(int a_pointIdx) const +int CellStreamCache::GetPointAdjacentCellCount(int a_pointIdx) const { if (a_pointIdx < 0 || a_pointIdx >= m_locations.size()) return 0; int numCells = m_pointsToCells[m_pointIdxToPointsToCells[a_pointIdx]]; return numCells; -} // XmUGrid::Impl::GetPointAdjacentCellCount - +} // CellStreamCache::GetPointAdjacentCellCount //------------------------------------------------------------------------------ -/// \brief Get the cells that are associated with the specified point -/// \param[in] a_pointIdx the index of the point -/// \return a vector of the cell indices associated with this point +/// \brief Get the cell that is associated with the specified point. +/// \param[in] a_pointIdx the index of the point. +/// \param[in] a_cellIdx the index of the cell. +/// \return the cell index associated with this point. //------------------------------------------------------------------------------ -VecInt XmUGrid::Impl::GetPointAdjacentCells(int a_pointIdx) const +int CellStreamCache::GetPointAdjacentCell(const int a_pointIdx, const int a_cellIdx) const { - VecInt cellsOfPoint; - GetPointAdjacentCells(a_pointIdx, cellsOfPoint); - return cellsOfPoint; -} // XmUGrid::Impl::GetPointAdjacentCells + if (a_pointIdx < 0 || a_pointIdx >= m_locations.size()) + return -1; + int numCells = m_pointsToCells[m_pointIdxToPointsToCells[a_pointIdx]]; + if (a_cellIdx < 0 || a_cellIdx >= numCells) + return -1; + return m_pointsToCells[m_pointIdxToPointsToCells[a_pointIdx] + a_cellIdx + 1]; +} // CellStreamCache::GetPointAdjacentCell //------------------------------------------------------------------------------ -/// \brief Get the cells that are associated with the specified point -/// \param[in] a_pointIdx the index of the point -/// \param[out] a_adjacentCells a vector of the adjacent cell indices -/// \return a vector of the cell indices associated with this point +/// \brief Get the number of cells. +/// \return the number of cells //------------------------------------------------------------------------------ -void XmUGrid::Impl::GetPointAdjacentCells(int a_pointIdx, VecInt& a_adjacentCells) const +int CellStreamCache::GetCellCount() const { - a_adjacentCells.clear(); - int numCells = GetPointAdjacentCellCount(a_pointIdx); - for (int cellIdx = 0; cellIdx < numCells; cellIdx++) - { - a_adjacentCells.push_back(m_pointsToCells[m_pointIdxToPointsToCells[a_pointIdx] + cellIdx + 1]); - } -} // XmUGrid::Impl::GetPointAdjacentCells + if (m_cellIdxToStreamIdx.empty()) + return 0; + else + return (int)m_cellIdxToStreamIdx.size() - 1; +} // CellStreamCache::GetCellCount + //------------------------------------------------------------------------------ -/// \brief Gets the common cells from a vector of points -/// \param[in] a_points a vector of unique points -/// \return a vector of cell indices +/// \brief Get cell stream vector for the entire UGrid. +/// \return constant reference to the cell stream vector //------------------------------------------------------------------------------ -VecInt XmUGrid::Impl::GetPointsAdjacentCells(const VecInt& a_points) const +const VecInt& CellStreamCache::GetCellstream() const { - VecInt commonCells; - if (!a_points.empty()) - GetPointsAdjacentCells(&a_points[0], (int)a_points.size(), commonCells); - return commonCells; -} // XmUGrid::Impl::GetPointsAdjacentCells + return m_cellstream; +} // CellStreamCache::GetCellstream //------------------------------------------------------------------------------ -/// \brief Gets the common cells from a vector of points -/// \param[in] a_pointIdxs an array of point indices -/// \param[in] a_numPointIdxs number of points in array -/// \param[out] a_commonCellIdxs a vector of cell indices +/// \brief Set the ugrid cells for the entire UGrid using a cell stream. +/// \param[in] a_cellstream cells defined as follows: +/// Hexahedrons, polygons, quads, triangles etc: +/// Cell type (ElemTypeEnum), number of points, point numbers. +/// Generally 0-based, CCW, bottom, then top. Not true +/// for pixel or voxel. +/// Polyhedrons: +/// Cell type, number of faces, [num points in face, +/// point numbers (0-based, CCW when looking in)] repeated +/// for each face. +/// \return if successfully set //------------------------------------------------------------------------------ -void XmUGrid::Impl::GetPointsAdjacentCells(const int* a_pointIdxs, - int a_numPointIdxs, - VecInt& a_commonCellIdxs) const +bool CellStreamCache::SetCellstream(const VecInt& a_cellstream) { - a_commonCellIdxs.clear(); - if (a_numPointIdxs == 0) - return; - a_commonCellIdxs = GetPointAdjacentCells(a_pointIdxs[0]); - for (int i = 1; i < a_numPointIdxs; ++i) + if (iIsValidCellstream(a_cellstream, (int)m_locations.size())) { - VecInt tempAssociatedCells = GetPointAdjacentCells(a_pointIdxs[i]); - VecInt remove; - for (int j = 0; j < a_commonCellIdxs.size(); ++j) - { - bool found = false; - for (int tempAssociatedCell : tempAssociatedCells) - { - if (a_commonCellIdxs[j] == tempAssociatedCell) - { - found = true; - break; - } - } - if (!found) - remove.push_back(j); - } - for (int j = (int)remove.size() - 1; j >= 0; --j) - { - a_commonCellIdxs.erase(a_commonCellIdxs.begin() + remove[j]); - } - if (a_commonCellIdxs.empty()) - break; + m_cellstream = a_cellstream; + Update(); + return true; } -} // XmUGrid::Impl::GetPointsAdjacentCells -//------------------------------------------------------------------------------ -/// \brief Gets the cells adjacent to all of a vector of points. -/// \param[in] a_pointIdxs a vector of unique points -/// \param[out] a_adjacentCellIdxs a vector of cell indices -//------------------------------------------------------------------------------ -void XmUGrid::Impl::GetPointsAdjacentCells(const VecInt& a_pointIdxs, - VecInt& a_adjacentCellIdxs) const -{ - if (!a_pointIdxs.empty()) - GetPointsAdjacentCells(&a_pointIdxs[0], (int)a_pointIdxs.size(), a_adjacentCellIdxs); -} // XmUGrid::Impl::GetPointsAdjacentCells + else + { + XM_LOG(xmlog::error, "Invalid cell stream data."); + return false; + } +} // CellStreamCache::SetCellstream //------------------------------------------------------------------------------ -/// \brief Gets the cells adjacent to both of the two points. -/// \param[in] a_pointIdx1 first point index -/// \param[in] a_pointIdx2 second point index -/// \param[out] a_adjacentCellIdxs a vector of cell indices +/// \brief Get the cell type of a specified cell. +/// \param[in] a_cellIdx the index of the cell +/// \return The type of the specified cell or -1 if invalid. //------------------------------------------------------------------------------ -void XmUGrid::Impl::GetPointsAdjacentCells(int a_pointIdx1, - int a_pointIdx2, - VecInt& a_adjacentCellIdxs) const +XmUGridCellType CellStreamCache::GetCellType(int a_cellIdx) const { - int points[] = {a_pointIdx1, a_pointIdx2}; - GetPointsAdjacentCells(points, 2, a_adjacentCellIdxs); -} // XmUGrid::Impl::GetPointsAdjacentCells - -//------------------------------------------------------------------------------ -/// \brief Given a point gets point indices attached to the point by an edge. -/// \param[in] a_pointIdx The point to get adjacent points from. -/// \param[out] a_edgePoints The indices of the adjacent points. -//------------------------------------------------------------------------------ -void XmUGrid::GetPointAdjacentPoints(int a_pointIdx, VecInt& a_edgePoints) const -{ - m_impl->GetPointAdjacentPoints(a_pointIdx, a_edgePoints); -} // XmUGrid::GetPointAdjacentPoints - -//------------------------------------------------------------------------------ -/// \brief Given a point gets point locations attached to the point by an edge. -/// \param[in] a_pointIdx The point to get attached point from. -/// \param[out] a_edgePoints A vector of points attached across edges. -//------------------------------------------------------------------------------ -void XmUGrid::GetPointAdjacentLocations(int a_pointIdx, VecPt3d& a_edgePoints) const -{ - m_impl->GetPointAdjacentLocations(a_pointIdx, a_edgePoints); -} // XmUGrid::GetPointAdjacentLocations -// Cells -//------------------------------------------------------------------------------ -/// \brief Get the number of cells. -/// \return the number of cells -//------------------------------------------------------------------------------ -int XmUGrid::Impl::GetCellCount() const -{ - if (m_cellIdxToStreamIdx.empty()) - return 0; - else - return (int)m_cellIdxToStreamIdx.size() - 1; -} // XmUGrid::Impl::GetCellCount -//------------------------------------------------------------------------------ -/// \brief Get the number of cell points (including polyhedron). -/// \param[in] a_cellIdx the index of the cell -/// \return a vector of point indices -//------------------------------------------------------------------------------ -int XmUGrid::Impl::GetCellPointCount(int a_cellIdx) const -{ - int pointCount = 0; - if (GetCellType(a_cellIdx) == XMU_POLYHEDRON) + if (a_cellIdx < 0 || m_cellIdxToStreamIdx.size() < 2 || + a_cellIdx > m_cellIdxToStreamIdx.size() - 2) { - boost::container::flat_set uniqueIdxs; - const int* cellstream = nullptr; - int streamLength; - GetCellCellstream(a_cellIdx, &cellstream, streamLength); - if (cellstream) - { - int currIdx = 1; - int numFaces = cellstream[currIdx++]; - for (int faceIdx = 0; faceIdx < numFaces; faceIdx++) - { - int numPoints = cellstream[currIdx++]; - for (int pointIdx = 0; pointIdx < numPoints; ++pointIdx) - { - int point = cellstream[currIdx++]; - uniqueIdxs.insert(point); - } - } - } - pointCount = (int)uniqueIdxs.size(); + return XMU_INVALID_CELL_TYPE; } else { - pointCount = GetNumberOfItemsForCell(a_cellIdx); + int cellStart = m_cellIdxToStreamIdx[a_cellIdx]; +#if _DEBUG + if ((m_cellstream[cellStart] > XMU_PYRAMID) && (m_cellstream[cellStart] != XMU_POLYHEDRON) || + (m_cellstream[cellStart] == XMU_TRIANGLE_STRIP)) + { + assert("UNSUPPORTED TYPE!"); + } +#endif + return (XmUGridCellType)m_cellstream[cellStart]; } - - return pointCount; -} // XmUGrid::Impl::GetCellPointCount +} // CellStreamCache::GetCellType //------------------------------------------------------------------------------ -/// \brief Get the points of a cell (including polyhedron) +/// \brief Get cell stream vector for a single cell. /// \param[in] a_cellIdx the index of the cell -/// \return a vector of point indices +/// \param[in] a_cellstream The cellstream of the cell +/// @see SetCellstream for more detail on cell stream definitions. +/// \return whether it was successfull or not //------------------------------------------------------------------------------ -VecInt XmUGrid::Impl::GetCellPoints(int a_cellIdx) const +bool CellStreamCache::GetCellCellstream(int a_cellIdx, VecInt& a_cellstream) const { - VecInt pointsOfCell; - GetCellPoints(a_cellIdx, pointsOfCell); - return pointsOfCell; -} // XmUGrid::Impl::GetCellPoints + a_cellstream.clear(); + if (a_cellIdx < 0 || m_cellIdxToStreamIdx.size() < 2 || + a_cellIdx > m_cellIdxToStreamIdx.size() - 2) + { + return false; + } + int startIndex(m_cellIdxToStreamIdx[a_cellIdx]), endIndex(m_cellIdxToStreamIdx[a_cellIdx + 1]); + a_cellstream.assign(m_cellstream.begin() + startIndex, m_cellstream.begin() + endIndex); + return true; +} // CellStreamCache::GetCellCellstream //------------------------------------------------------------------------------ -/// \brief Get the points of a cell. -/// \param[in] a_cellIdx the index of the cell -/// \param[out] a_cellPoints the points of the cell -/// \return if the cell index is valid +/// \brief Internal function to get start of cell stream for an individual cell. +/// Returns nullptr and zero length for invalid cell index. +/// \param[in] a_cellIdx the index of the cell. +/// \return The cell stream. //------------------------------------------------------------------------------ -bool XmUGrid::Impl::GetCellPoints(int a_cellIdx, VecInt& a_cellPoints) const +Span CellStreamCache::GetCellCellstream(int a_cellIdx) const { - a_cellPoints.clear(); - VecInt cellstream; - if (GetCellCellstream(a_cellIdx, cellstream)) + int cellType = GetCellType(a_cellIdx); + if (cellType == XMU_INVALID_CELL_TYPE) { - XM_ENSURE_TRUE_NO_ASSERT(!cellstream.empty(), false); - int cellType = cellstream[0]; - if (cellType == XMU_POLYHEDRON) - { - GetUniquePointsFromPolyhedronSingleCellstream(cellstream, a_cellPoints); - return true; - } - else if (cellType == XMU_PIXEL) - { - a_cellPoints.push_back(cellstream[2]); - a_cellPoints.push_back(cellstream[3]); - a_cellPoints.push_back(cellstream[5]); - a_cellPoints.push_back(cellstream[4]); - return true; - } - else - { - a_cellPoints.assign(cellstream.begin() + 2, cellstream.end()); - return true; - } + return {}; } - return false; -} // XmUGrid::Impl::GetCellPoints - + else + { + int startIndex = m_cellIdxToStreamIdx[a_cellIdx]; + int endIndex = m_cellIdxToStreamIdx[a_cellIdx + 1]; + return {&m_cellstream[startIndex], endIndex - startIndex}; + } +} // CellStreamCache::GetCellCellstream //------------------------------------------------------------------------------ -/// \brief Get locations of cell points. +/// \brief Get beginning index of cell in cell stream. /// \param[in] a_cellIdx the index of the cell -/// \param[out] a_cellLocations The locations of the cell points +/// \return The index of the cell in the cell stream. //------------------------------------------------------------------------------ -void XmUGrid::Impl::GetCellLocations(int a_cellIdx, VecPt3d& a_cellLocations) const +int CellStreamCache::GetCellCellstreamIndex(int a_cellIdx) const { - VecInt ptIdxs = GetCellPoints(a_cellIdx); - a_cellLocations = GetPointsLocations(ptIdxs); -} // XmUGrid::Impl::GetCellPoints - + if (a_cellIdx < 0 || a_cellIdx >= (int)m_cellIdxToStreamIdx.size()) + { + return -1; + } + return m_cellIdxToStreamIdx[a_cellIdx]; +} // CellStreamCache::GetCellCellstreamIndex //------------------------------------------------------------------------------ -/// \brief Get the cell type of a specified cell. +/// \brief Get number of items given cell. For polyhedron number of items is +/// number of faces. For other cell types it is number of points. /// \param[in] a_cellIdx the index of the cell -/// \return The type of the specified cell or -1 if invalid. +/// \return the number of faces for polyhedron or number of points. //------------------------------------------------------------------------------ -XmUGridCellType XmUGrid::Impl::GetCellType(int a_cellIdx) const +int CellStreamCache::GetNumberOfItemsForCell(int a_cellIdx) const { - if (a_cellIdx < 0 || m_cellIdxToStreamIdx.size() < 2 || - a_cellIdx > m_cellIdxToStreamIdx.size() - 2) + int cellType = GetCellType(a_cellIdx); + if (cellType == XMU_INVALID_CELL_TYPE) { - return XMU_INVALID_CELL_TYPE; + return 0; } - else + else if (cellType != XMU_POLYHEDRON) { - int cellStart = m_cellIdxToStreamIdx[a_cellIdx]; -#if _DEBUG - if ((m_cellstream[cellStart] > XMU_PYRAMID) && (m_cellstream[cellStart] != XMU_POLYHEDRON) || - (m_cellstream[cellStart] == XMU_TRIANGLE_STRIP)) - { - assert("UNSUPPORTED TYPE!"); - } -#endif - return (XmUGridCellType)m_cellstream[cellStart]; + // Number of points is right after cell type + int startIndex = m_cellIdxToStreamIdx[a_cellIdx]; + return m_cellstream[startIndex + 1]; } -} // XmUGrid::Impl::GetCellType + else + { // Polyhedron case, same value is number of faces + int startIndex = m_cellIdxToStreamIdx[a_cellIdx]; + return m_cellstream[startIndex + 1]; + } +} // CellStreamCache::GetNumberOfItemsForCell +//------------------------------------------------------------------------------ +/// \brief Get the dimension of the specified cell. +/// \param[in] a_cellIdx the index of the cell +/// \return the dimension of the cells or -1 if invalid index or invalid +/// dimension +//------------------------------------------------------------------------------ +int CellStreamCache::GetCellDimension(int a_cellIdx) const +{ + return iDimensionFromCellType(GetCellType(a_cellIdx)); +} // CellStreamCache::GetCellDimension //------------------------------------------------------------------------------ /// \brief Count all number of the cells with each dimension (0, 1, 2, 3) /// \return the count of dimensions of all of the cells of the ugrid //------------------------------------------------------------------------------ -std::vector XmUGrid::Impl::GetDimensionCounts() const +std::vector CellStreamCache::GetDimensionCounts() const { if (!m_cellDimensionCounts.empty()) { @@ -1166,1368 +1438,1483 @@ std::vector XmUGrid::Impl::GetDimensionCounts() const { m_cellDimensionCounts.clear(); m_cellDimensionCounts.resize(4, 0); - int itemp = 0; int cellCount = GetCellCount(); - for (int i = 0; i < cellCount; i++) + for (int cellIndex = 0; cellIndex < cellCount; cellIndex++) { - itemp = GetCellDimension(i); - if (itemp >= 0) + auto cellDimension = GetCellDimension(cellIndex); + if (cellDimension >= 0) { - m_cellDimensionCounts[itemp]++; + m_cellDimensionCounts[cellDimension]++; } } return m_cellDimensionCounts; } -} // XmUGrid::Impl::GetDimensionCounts +} // CellStreamCache::GetDimensionCounts //------------------------------------------------------------------------------ -/// \brief Get the dimension of the specified cell. -/// \param[in] a_cellIdx the index of the cell -/// \return the dimension of the cells or -1 if invalid index or invalid -/// dimension +/// \brief Update internal link from cells to cell stream index. //------------------------------------------------------------------------------ -int XmUGrid::Impl::GetCellDimension(int a_cellIdx) const +void CellStreamCache::UpdateCellLinks() { - return DimensionFromCellType(GetCellType(a_cellIdx)); -} // XmUGrid::Impl::GetCellDimension + m_cellIdxToStreamIdx.clear(); -//------------------------------------------------------------------------------ -/// \brief Get the extents of the given cell. -/// \param[in] a_cellIdx The cell index to get the extents of. -/// \param[out] a_min The minimum location. -/// \param[out] a_max The maximum location. -//------------------------------------------------------------------------------ -void XmUGrid::Impl::GetCellExtents(int a_cellIdx, Pt3d& a_min, Pt3d& a_max) const -{ - VecPt3d pts; - GetCellLocations(a_cellIdx, pts); - GetExtentsFromPoints(pts, a_min, a_max); -} // XmUGrid::Impl::GetCellExtents + int currIdx = 0; + if (m_cellstream.empty()) + { + m_cellIdxToStreamIdx.push_back(currIdx); + return; + } + + int numItems = (int)m_cellstream.size(); + while (currIdx < numItems) + { + m_cellIdxToStreamIdx.push_back(currIdx); + + // get cell type + int cellType = m_cellstream[currIdx++]; + if (currIdx >= numItems) + return; + + // get the number of items + int numPoints = m_cellstream[currIdx++]; + if (currIdx >= numItems) + return; + if (cellType == XMU_POLYHEDRON) + { + int numFaces = numPoints; + for (int faceIdx = 0; faceIdx < numFaces; ++faceIdx) + { + int numFacePoints = m_cellstream[currIdx++]; + currIdx += numFacePoints; + } + } + else + { + currIdx += numPoints; + } + } + + m_cellIdxToStreamIdx.push_back(currIdx); +} // CellStreamCache::UpdateCellLinks //------------------------------------------------------------------------------ -/// \brief Get cell stream vector for the entire UGrid. -/// \return constant reference to the cell stream vector +/// \brief Update internal links from points to associated cells. //------------------------------------------------------------------------------ -const VecInt& XmUGrid::Impl::GetCellstream() const +void CellStreamCache::UpdatePointLinks() { - return m_cellstream; -} // XmUGrid::Impl::GetCellstream -//------------------------------------------------------------------------------ -/// \brief Set the ugrid cells for the entire UGrid using a cell stream. -/// \param[in] a_cellstream cells defined as follows: -/// Hexahedrons, polygons, quads, triangles etc: -/// Cell type (ElemTypeEnum), number of points, point numbers. -/// Generally 0-based, CCW, bottom, then top. Not true -/// for pixel or voxel. -/// Polyhedrons: -/// Cell type, number of faces, [num points in face, -/// point numbers (0-based, CCW when looking in)] repeated -/// for each face. -/// \return if successfully set -//------------------------------------------------------------------------------ -bool XmUGrid::Impl::SetCellstream(const VecInt& a_cellstream) -{ - if (IsValidCellstream(a_cellstream, (int)m_locations.size())) + m_pointIdxToPointsToCells.clear(); + m_pointIdxToPointsToCells.resize(m_locations.size(), 0); + m_pointsToCells.clear(); + + // get number of cells for each point + int numStreamItems = (int)m_cellstream.size(); + int cellIdx = 0; + int currIdx = 0; + VecInt pointLastUsedIdx(m_locations.size(), -1); + VecInt cellPoints; + while (currIdx < numStreamItems) { - m_cellstream = a_cellstream; - UpdateLinks(); - SetModified(); - return true; + // get cell type + int cellType = m_cellstream[currIdx++]; + if (currIdx >= numStreamItems) + return; + + // get the number of items (points or faces depending on cell type) + int numCellItems = m_cellstream[currIdx++]; + if (currIdx >= numStreamItems) + return; + + if (cellType == XMU_POLYHEDRON) + { + iGetUniquePointsFromPolyhedronCellstream(m_cellstream, numCellItems, currIdx, cellPoints, + pointLastUsedIdx); + + // Deemed to be slower than flat set-- Left only as a warning to others! + // std::stable_sort(cellPoints.begin(), cellPoints.end()); + // auto uniqueEnd = std::unique(cellPoints.begin(), cellPoints.end()); + + // flat set deemed to be slower than marking when a point was last used (pointLastUsedIdx) + + // for (auto pt = cellPoints.begin(); pt != uniqueEnd; ++pt) + for (auto cellPoint : cellPoints) + { + m_pointIdxToPointsToCells[cellPoint] += 1; + } + } + else + { + // iterate on points + int numPoints = numCellItems; + for (int ptIdx = 0; ptIdx < numPoints; ++ptIdx) + { + int pt = m_cellstream[currIdx++]; + m_pointIdxToPointsToCells[pt] += 1; + } + } + + ++cellIdx; } - else + m_pointIdxToPointsToCells.push_back(0); + + // change array of counts to array of offsets + int currCount = 0; + for (auto& item : m_pointIdxToPointsToCells) { - XM_LOG(xmlog::error, "Invalid cell stream data."); - return false; + int count = item + 1; + item = currCount; + currCount += count; } -} // XmUGrid::Impl::SetCellstream + + cellIdx = 0; + currIdx = 0; + m_pointsToCells.resize(currCount, 0); + std::fill(pointLastUsedIdx.begin(), pointLastUsedIdx.end(), -1); + while (currIdx < numStreamItems) + { + // get cell type + int cellType = m_cellstream[currIdx++]; + if (currIdx >= numStreamItems) + return; + + // get the number of items (points or faces depending on cell type) + int numCellItems = m_cellstream[currIdx++]; + if (currIdx >= numStreamItems) + return; + + if (cellType == XMU_POLYHEDRON) + { + iGetUniquePointsFromPolyhedronCellstream(m_cellstream, numCellItems, currIdx, cellPoints, + pointLastUsedIdx); + + // Deemed to be slower than flat set-- Left only as a warning to others! + // std::stable_sort(cellPoints.begin(), cellPoints.end()); + // auto uniqueEnd = std::unique(cellPoints.begin(), cellPoints.end()); + + // for (auto pt = cellPoints.begin(); pt != uniqueEnd; ++pt) + for (auto cellPoint : cellPoints) + { + int countIdx = m_pointIdxToPointsToCells[cellPoint]; + int& count = m_pointsToCells[countIdx]; // point's cell count + ++count; // incrementing point's cell count + int updateIdx = countIdx + count; + m_pointsToCells[updateIdx] = cellIdx; + } + } + else + { + // iterate on points + int numPoints = numCellItems; + for (int ptIdx = 0; ptIdx < numPoints; ++ptIdx) + { + int pt = m_cellstream[currIdx++]; + int countIdx = m_pointIdxToPointsToCells[pt]; + int& count = m_pointsToCells[countIdx]; // point's cell count + ++count; // incrementing point's cell count + int updateIdx = countIdx + count; + m_pointsToCells[updateIdx] = cellIdx; + } + } + + ++cellIdx; + } +} // CellStreamCache::UpdatePointLinks +} // namespace + +//////////////////////////////////////////////////////////////////////////////// +/// \class XmUGrid::Impl +/// \brief Implementation for XmUGrid which provides geometry for an +/// unstructured grid. +//////////////////////////////////////////////////////////////////////////////// //------------------------------------------------------------------------------ -/// \brief Get cell stream vector for a single cell. -/// \param[in] a_cellIdx the index of the cell -/// \param[in] a_cellstream The cellstream of the cell -/// @see SetCellstream for more detail on cell stream definitions. -/// \return whether it was successfull or not +/// \brief Returns the modified flag. Gets set when points or cells get changed. +/// \return the modified flag //------------------------------------------------------------------------------ -bool XmUGrid::Impl::GetCellCellstream(int a_cellIdx, VecInt& a_cellstream) const +bool XmUGrid::Impl::GetModified() const { - a_cellstream.clear(); - if (a_cellIdx < 0 || m_cellIdxToStreamIdx.size() < 2 || - a_cellIdx > m_cellIdxToStreamIdx.size() - 2) - { - return false; - } - int startIndex(m_cellIdxToStreamIdx[a_cellIdx]), endIndex(m_cellIdxToStreamIdx[a_cellIdx + 1]); - a_cellstream.assign(m_cellstream.begin() + startIndex, m_cellstream.begin() + endIndex); - return true; -} // XmUGrid::Impl::GetCellCellstream + return m_modified; +} // XmUGrid::Impl::GetModified //------------------------------------------------------------------------------ -/// \brief Get beginning index of cell in cell stream. -/// \param[in] a_cellIdx the index of the cell -/// \return The index of the cell in the cell stream. +/// \brief Resets the modified flag to false. //------------------------------------------------------------------------------ -int XmUGrid::Impl::GetCellCellstreamIndex(int a_cellIdx) const +void XmUGrid::Impl::SetUnmodified() { - if (a_cellIdx < 0 || a_cellIdx >= (int)m_cellIdxToStreamIdx.size()) + m_modified = false; +} // XmUGrid::Impl::SetUnmodified +//------------------------------------------------------------------------------ +/// \brief Sets the modified flag to true. +//------------------------------------------------------------------------------ +void XmUGrid::Impl::SetModified() +{ + ClearCacheValues(); + m_modified = true; +} // XmUGrid::Impl::SetModified +//------------------------------------------------------------------------------ +/// \brief Turn on or off use of caching. +/// \param a_useCache Flag to determine if caching will be used. +//------------------------------------------------------------------------------ +void XmUGrid::Impl::SetUseCache(bool a_useCache) +{ + m_useCache = a_useCache; +} // XmUGrid::Impl::SetUseCache +// Points +//------------------------------------------------------------------------------ +/// \brief Get the number of points. +/// \return the number of points +//------------------------------------------------------------------------------ +int XmUGrid::Impl::GetPointCount() const +{ + return m_cellStreamCache.GetPointCount(); +} // XmUGrid::Impl::GetPointCount +//------------------------------------------------------------------------------ +/// \brief Get vector of UGrid points. +/// \return constant reference to vector of points +//------------------------------------------------------------------------------ +const VecPt3d& XmUGrid::Impl::GetLocations() const +{ + return m_cellStreamCache.GetLocations(); +} // XmUGrid::Impl::GetLocations +//------------------------------------------------------------------------------ +/// \brief Set UGrid points. +/// \param[in] a_locations vector of point locations +//------------------------------------------------------------------------------ +void XmUGrid::Impl::SetLocations(const VecPt3d& a_locations) +{ + m_cellStreamCache.SetLocations(a_locations); + SetModified(); +} // XmUGrid::Impl::SetLocations +//------------------------------------------------------------------------------ +/// \brief Get the point +/// \param[in] a_pointIdx the index of the point +/// \return the point or an initialize point if the index is invalid +//------------------------------------------------------------------------------ +Pt3d XmUGrid::Impl::GetPointLocation(int a_pointIdx) const +{ + return m_cellStreamCache.GetPointLocation(a_pointIdx); +} // XmUGrid::Impl::GetPointLocation +//------------------------------------------------------------------------------ +/// \brief Set the point +/// \param[in] a_pointIdx the index of the point +/// \param[in] a_location The new location of the specified point +/// \return whether the point was successfully set +//------------------------------------------------------------------------------ +bool XmUGrid::Impl::SetPointLocation(int a_pointIdx, const Pt3d& a_location) +{ + bool modified = m_cellStreamCache.SetPointLocation(a_pointIdx, a_location); + if (modified) + SetModified(); + return modified; +} // XmUGrid::Impl::SetPointLocation + +//------------------------------------------------------------------------------ +/// \brief Get the X, Y location of a point. +/// \param[in] a_pointIdx The index of the point. +/// \return The location of the point with Z set to 0.0. +//------------------------------------------------------------------------------ +Pt3d XmUGrid::Impl::GetPointXy0(int a_pointIdx) const +{ + Pt3d pt = GetPointLocation(a_pointIdx); + pt.z = 0.0; + return pt; +} // XmUGrid::Impl::GetPointXy0 + +//------------------------------------------------------------------------------ +/// \brief Convert a vector of point indices into a vector of point 3d +/// \param[in] a_points a vector of point indices +/// \return vector of point 3d +//------------------------------------------------------------------------------ +VecPt3d XmUGrid::Impl::GetPointsLocations(const VecInt& a_points) const +{ + VecPt3d point3d; + for (auto point : a_points) { - return -1; + point3d.push_back(GetPointLocation(point)); } - return m_cellIdxToStreamIdx[a_cellIdx]; -} // XmUGrid::GetCellCellstreamIndex + return point3d; +} // XmUGrid::Impl::GetPointsLocations + //------------------------------------------------------------------------------ -/// \brief Get the cells neighboring a cell (cells associated with any of it's points) -/// \param[in] a_cellIdx the index of the cell -/// \return vector of cell indices +/// \brief Get extents of all points in UGrid +/// \param[out] a_min minimum extent of all points +/// \param[out] a_max maximum extent of all points +//------------------------------------------------------------------------------ +void XmUGrid::Impl::GetExtents(Pt3d& a_min, Pt3d& a_max) const +{ + m_cellStreamCache.GetExtents(a_min, a_max); +} // XmUGrid::Impl::GetExtents + +//------------------------------------------------------------------------------ +/// \brief Get the number of cells that use a point. +/// \param a_pointIdx The point to check. +/// \return The number of cells that use the point. +//------------------------------------------------------------------------------ +int XmUGrid::Impl::GetPointAdjacentCellCount(int a_pointIdx) const +{ + return m_cellStreamCache.GetPointAdjacentCellCount(a_pointIdx); +} // XmUGrid::Impl::GetPointAdjacentCellCount + //------------------------------------------------------------------------------ -VecInt XmUGrid::Impl::GetCellAdjacentCells(int a_cellIdx) const +/// \brief Get the cells that are associated with the specified point +/// \param[in] a_pointIdx the index of the point +/// \return a vector of the cell indices associated with this point +//------------------------------------------------------------------------------ +VecInt XmUGrid::Impl::GetPointAdjacentCells(int a_pointIdx) const { - VecInt neighbors; - GetCellAdjacentCells(a_cellIdx, neighbors); - return neighbors; -} // XmUGrid::Impl::GetCellAdjacentCells + VecInt cellsOfPoint; + GetPointAdjacentCells(a_pointIdx, cellsOfPoint); + return cellsOfPoint; +} // XmUGrid::Impl::GetPointAdjacentCells //------------------------------------------------------------------------------ -/// \brief Get the cells neighboring a cell (cells associated with any of it's points) -/// \param[in] a_cellIdx the index of the cell -/// \param[out] a_cellNeighbors vector of cell indices +/// \brief Get the cells that are associated with the specified point +/// \param[in] a_pointIdx the index of the point +/// \param[out] a_adjacentCells a vector of the adjacent cell indices +/// \return a vector of the cell indices associated with this point //------------------------------------------------------------------------------ -void XmUGrid::Impl::GetCellAdjacentCells(int a_cellIdx, VecInt& a_cellNeighbors) const +void XmUGrid::Impl::GetPointAdjacentCells(int a_pointIdx, VecInt& a_adjacentCells) const { - a_cellNeighbors.clear(); - VecInt pointsOfCell; - pointsOfCell = GetCellPoints(a_cellIdx); - for (int pointIdx : pointsOfCell) - { - VecInt associatedCells = GetPointAdjacentCells(pointIdx); - for (int associatedCell : associatedCells) - { - if (associatedCell != a_cellIdx) - { - bool found(false); - for (int a_cellNeighbor : a_cellNeighbors) - { - if (a_cellNeighbor == associatedCell) - { - found = true; - } - } - if (!found) - { - a_cellNeighbors.push_back(associatedCell); - } - } - } - } -} // XmUGrid::Impl::GetCellAdjacentCells + m_cellStreamCache.GetPointAdjacentCells(a_pointIdx, a_adjacentCells); +} // XmUGrid::Impl::GetPointAdjacentCells //------------------------------------------------------------------------------ -/// \brief Get a plan view polygon of a specified cell -/// \param[in] a_cellIdx The index of the cell. -/// \param[out] a_polygon Vector of Pt3d that is the plan view polygon. -/// \return False if the cell index does not exist or if the cell is not 2 or 3 -/// dimensional. +/// \brief Get the cells that are associated with the specified point +/// \param[in] a_pointIdx The index of the point. +/// \return The cell indices associated with this point. //------------------------------------------------------------------------------ -bool XmUGrid::Impl::GetCellPlanViewPolygon(int a_cellIdx, VecPt3d& a_polygon) const +void XmUGrid::Impl::GetPointAdjacentCells(int a_pointIdx, Span& a_adjacentCells) const { - a_polygon.clear(); - if (!IsValidCellIdx(a_cellIdx)) - return false; - int dimension = GetCellDimension(a_cellIdx); - if (dimension == 3) - return GetPlanViewPolygon3d(a_cellIdx, a_polygon); - else if (dimension == 2) - return GetPlanViewPolygon2d(a_cellIdx, a_polygon); - else - return false; -} // XmUGrid::Impl::GetCellPlanViewPolygon - + m_cellStreamCache.GetPointAdjacentCells(a_pointIdx, a_adjacentCells); +} // XmUGrid::Impl::GetPointAdjacentCells //------------------------------------------------------------------------------ -/// \brief Get the centroid location of a cell. -/// \param[in] a_cellIdx The index of the cell. -/// \param[out] a_centroid The location of the cell centroid. -/// \return False if the cell index does not exist. +/// \brief Gets the common cells from a vector of points +/// \param[in] a_points a vector of unique points +/// \return a vector of cell indices //------------------------------------------------------------------------------ -bool XmUGrid::Impl::GetCellCentroid(int a_cellIdx, Pt3d& a_centroid) const +VecInt XmUGrid::Impl::GetPointsAdjacentCells(const VecInt& a_points) const { - VecPt3d polyPoints; - bool retVal = true; - Pt3d centroid = Pt3d(0.0, 0.0, 0.0); - if (GetCellPlanViewPolygon(a_cellIdx, polyPoints)) - { - centroid = gmComputePolygonCentroid(polyPoints); - } - else if (!IsValidCellIdx(a_cellIdx)) - { - retVal = false; - } - else - { - VecPt3d cellPoints; - GetCellLocations(a_cellIdx, cellPoints); - for (Pt3d& pt : cellPoints) - { - centroid += pt; - } - centroid /= (double)cellPoints.size(); - } - a_centroid = centroid; - return retVal; -} // XmUGrid::Impl::GetCellCentroid + VecInt commonCells; + if (!a_points.empty()) + GetPointsAdjacentCells(&a_points[0], (int)a_points.size(), commonCells); + return commonCells; +} // XmUGrid::Impl::GetPointsAdjacentCells //------------------------------------------------------------------------------ -/// \brief Determine whether a cell is valid after a point is moved. -/// \param[in] a_cellIdx the index of the cell -/// \param[in] a_changedPtIdx index of the point to be changed -/// \param[in] a_newPosition location the point is to be moved to -/// \return whether the cell is valid +/// \brief Gets the common cells from a vector of points +/// \param[in] a_pointIdxs an array of point indices +/// \param[in] a_numPointIdxs number of points in array +/// \param[out] a_commonCellIdxs a vector of cell indices //------------------------------------------------------------------------------ -bool XmUGrid::Impl::IsCellValidWithPointChange(int a_cellIdx, - int a_changedPtIdx, - const Pt3d& a_newPosition) const +void XmUGrid::Impl::GetPointsAdjacentCells(const int* a_pointIdxs, + int a_numPointIdxs, + VecInt& a_commonCellIdxs) const { - if (GetCellDimension(a_cellIdx) == 2) - { - std::vector edges = GetCellEdges(a_cellIdx); - return !DoEdgesCrossWithPointChange(a_changedPtIdx, a_newPosition, edges); - } - else if (GetCellDimension(a_cellIdx) == 3) + a_commonCellIdxs.clear(); + if (a_numPointIdxs == 0) + return; + a_commonCellIdxs = GetPointAdjacentCells(a_pointIdxs[0]); + for (int i = 1; i < a_numPointIdxs; ++i) { - // Go through the affected faces and if they are not a side face - // check if the edges cross each other - VecInt2d faces = GetCell3dFacesPoints(a_cellIdx); - for (auto& face : faces) + VecInt tempAssociatedCells = GetPointAdjacentCells(a_pointIdxs[i]); + VecInt remove; + for (int j = 0; j < a_commonCellIdxs.size(); ++j) { - bool faceIsAffected = false; - for (int ptIdx : face) + bool found = false; + for (int tempAssociatedCell : tempAssociatedCells) { - if (ptIdx == a_changedPtIdx) + if (a_commonCellIdxs[j] == tempAssociatedCell) { - faceIsAffected = true; + found = true; + break; } } - if (faceIsAffected) - { - std::vector edges; - GetEdgesOfFace(face, edges); - if (DoEdgesCrossWithPointChange(a_changedPtIdx, a_newPosition, edges)) - return false; - } + if (!found) + remove.push_back(j); } + for (int j = (int)remove.size() - 1; j >= 0; --j) + { + a_commonCellIdxs.erase(a_commonCellIdxs.begin() + remove[j]); + } + if (a_commonCellIdxs.empty()) + break; } - return true; -} // XmUGrid::Impl::IsCellValidWithPointChange +} // XmUGrid::Impl::GetPointsAdjacentCells //------------------------------------------------------------------------------ -/// \brief Determine if a cell index is valid. -/// \param[in] a_cellIdx the index of the cell -/// \return True if the cell index is valid. +/// \brief Gets the cells adjacent to all of a vector of points. +/// \param[in] a_pointIdxs a vector of unique points +/// \param[out] a_adjacentCellIdxs a vector of cell indices //------------------------------------------------------------------------------ -bool XmUGrid::Impl::IsValidCellIdx(int a_cellIdx) const +void XmUGrid::Impl::GetPointsAdjacentCells(const VecInt& a_pointIdxs, + VecInt& a_adjacentCellIdxs) const { - return a_cellIdx >= 0 && a_cellIdx < GetCellCount(); -} // XmUGrid::Impl::IsValidCellIdx + if (!a_pointIdxs.empty()) + GetPointsAdjacentCells(&a_pointIdxs[0], (int)a_pointIdxs.size(), a_adjacentCellIdxs); +} // XmUGrid::Impl::GetPointsAdjacentCells //------------------------------------------------------------------------------ -/// \brief Determine whether adjacent cells are valid after a point is moved. -/// \param[in] a_changedPtIdx index of the point to be changed -/// \param[in] a_newPosition location the point is to be moved to -/// \return whether the change is valid +/// \brief Gets the cells adjacent to both of the two points. +/// \param[in] a_pointIdx1 first point index +/// \param[in] a_pointIdx2 second point index +/// \param[out] a_adjacentCellIdxs a vector of cell indices //------------------------------------------------------------------------------ -bool XmUGrid::Impl::IsValidPointChange(int a_changedPtIdx, const Pt3d& a_newPosition) const +void XmUGrid::Impl::GetPointsAdjacentCells(int a_pointIdx1, + int a_pointIdx2, + VecInt& a_adjacentCellIdxs) const { - bool validChange = false; - if (a_changedPtIdx >= 0 && a_changedPtIdx < m_locations.size()) - { - validChange = true; - VecInt affectedCells = GetPointAdjacentCells(a_changedPtIdx); - for (int affectedCell : affectedCells) - { - if (!IsCellValidWithPointChange(affectedCell, a_changedPtIdx, a_newPosition)) - { - validChange = false; - break; - } - } - } - return validChange; -} // XmUGrid::Impl::IsValidPointChange + int points[] = {a_pointIdx1, a_pointIdx2}; + GetPointsAdjacentCells(points, 2, a_adjacentCellIdxs); +} // XmUGrid::Impl::GetPointsAdjacentCells + //------------------------------------------------------------------------------ -/// \brief Get the edges of a cell given a face -/// \param[in] a_face a vector of point indices of a face -/// \param[out] a_edges a vector of point indices of an edge +/// \brief Given a point gets point indices attached to the point by an edge. +/// \param[in] a_pointIdx The point to get adjacent points from. +/// \param[out] a_edgePoints The indices of the adjacent points. //------------------------------------------------------------------------------ -void XmUGrid::Impl::GetEdgesOfFace(const VecInt& a_face, std::vector& a_edges) const +void XmUGrid::GetPointAdjacentPoints(int a_pointIdx, VecInt& a_edgePoints) const { - a_edges.reserve(a_face.size()); - for (int i = 1; i < a_face.size(); ++i) - { - a_edges.emplace_back(a_face[i - 1], a_face[i]); - } - if (a_edges.size() > 1) - a_edges.emplace_back(a_face[a_face.size() - 1], a_face[0]); -} // XmUGrid::Impl::GetEdgesOfFace + m_impl->GetPointAdjacentPoints(a_pointIdx, a_edgePoints); +} // XmUGrid::GetPointAdjacentPoints //------------------------------------------------------------------------------ -/// \brief Gets whether or not edges cross with a point change -/// \param[out] a_edges The edges to check -/// \param[in] a_changedPtIdx index of the point to be changed -/// \param[in] a_newPosition location the point is to be moved to -/// \return whether the edges cross +/// \brief Given a point gets point locations attached to the point by an edge. +/// \param[in] a_pointIdx The point to get attached point from. +/// \param[out] a_edgePoints A vector of points attached across edges. //------------------------------------------------------------------------------ -bool XmUGrid::Impl::DoEdgesCrossWithPointChange(int a_changedPtIdx, - const Pt3d& a_newPosition, - const std::vector& a_edges) const +void XmUGrid::GetPointAdjacentLocations(int a_pointIdx, VecPt3d& a_edgePoints) const +{ + m_impl->GetPointAdjacentLocations(a_pointIdx, a_edgePoints); +} // XmUGrid::GetPointAdjacentLocations +// Cells +//------------------------------------------------------------------------------ +/// \brief Get the number of cells. +/// \return the number of cells +//------------------------------------------------------------------------------ +int XmUGrid::Impl::GetCellCount() const +{ + return m_cellStreamCache.GetCellCount(); +} // XmUGrid::Impl::GetCellCount +//------------------------------------------------------------------------------ +/// \brief Get the number of cell points (including polyhedron). +/// \param[in] a_cellIdx the index of the cell +/// \return a vector of point indices +//------------------------------------------------------------------------------ +int XmUGrid::Impl::GetCellPointCount(int a_cellIdx) const { - std::vector> changedEdges; - std::vector> unChangedEdges; - for (const auto& edge : a_edges) - { - if (edge.GetFirst() == a_changedPtIdx) - { - changedEdges.emplace_back(a_newPosition, GetPointLocation(edge.GetSecond())); - } - else if (edge.GetSecond() == a_changedPtIdx) - { - changedEdges.emplace_back(GetPointLocation(edge.GetFirst()), a_newPosition); - } - else - { - unChangedEdges.emplace_back(GetPointLocation(edge.GetFirst()), - GetPointLocation(edge.GetSecond())); - } - } - for (const auto& changedEdge : changedEdges) + if (GetCellType(a_cellIdx) == XMU_POLYHEDRON) { - for (const auto& unChangedEdge : unChangedEdges) + boost::container::flat_set uniqueIdxs; + Span cellstream = GetCellCellstream(a_cellIdx); + if (!cellstream.empty()) { - if (gmLinesCross(changedEdge.first, changedEdge.second, unChangedEdge.first, - unChangedEdge.second)) - return true; + int currIdx = 1; + int numFaces = cellstream[currIdx++]; + for (int faceIdx = 0; faceIdx < numFaces; faceIdx++) + { + int numPoints = cellstream[currIdx++]; + for (int pointIdx = 0; pointIdx < numPoints; ++pointIdx) + { + int point = cellstream[currIdx++]; + uniqueIdxs.insert(point); + } + } } + return static_cast(uniqueIdxs.size()); } - return false; -} // XmUGrid::Impl::DoEdgesCrossWithPointChange -// Edges + return m_cellStreamCache.GetNumberOfItemsForCell(a_cellIdx); +} // XmUGrid::Impl::GetCellPointCount //------------------------------------------------------------------------------ -/// \brief Get the number of edges for a cell. +/// \brief Get the points of a cell (including polyhedron) /// \param[in] a_cellIdx the index of the cell -/// \return the count of cell edges +/// \return a vector of point indices //------------------------------------------------------------------------------ -int XmUGrid::Impl::GetCellEdgeCount(int a_cellIdx) const +VecInt XmUGrid::Impl::GetCellPoints(int a_cellIdx) const { - int cellType(GetCellType(a_cellIdx)); - switch (cellType) - { - // invalid - case XMU_INVALID_CELL_TYPE: - return -1; - break; - - case XMU_POLY_LINE: - return GetNumberOfItemsForCell(a_cellIdx) - 1; - break; - - case XMU_POLYGON: - return GetNumberOfItemsForCell(a_cellIdx); - break; - - case XMU_POLYHEDRON: - return GetNumberOfPolyhedronEdges(a_cellIdx); - break; - - default: - { - const VecEdge& edgeTable = iGetEdgeOffsetTable(cellType); - if (!edgeTable.empty()) - return (int)edgeTable.size(); - break; - } - } - - XM_ASSERT(0); - return -1; -} // XmUGrid::Impl::GetCellEdgeCount + VecInt pointsOfCell; + GetCellPoints(a_cellIdx, pointsOfCell); + return pointsOfCell; +} // XmUGrid::Impl::GetCellPoints //------------------------------------------------------------------------------ /// \brief Get the points of a cell. /// \param[in] a_cellIdx the index of the cell -/// \param[in] a_edgeIdx the index of the edge -/// \return a standard pair of point indices (which is an edge) +/// \param[out] a_cellPoints the points of the cell +/// \return if the cell index is valid //------------------------------------------------------------------------------ -XmEdge XmUGrid::Impl::GetCellEdge(int a_cellIdx, int a_edgeIdx) const +bool XmUGrid::Impl::GetCellPoints(int a_cellIdx, VecInt& a_cellPoints) const { - XmEdge edge; - if (a_edgeIdx < 0) - return edge; - - const int* cellstream = nullptr; - int streamLength; - GetCellCellstream(a_cellIdx, &cellstream, streamLength); - if (cellstream) + a_cellPoints.clear(); + VecInt cellstream; + if (GetCellCellstream(a_cellIdx, cellstream)) { - if (streamLength < 4) - return edge; - + XM_ENSURE_TRUE_NO_ASSERT(!cellstream.empty(), false); int cellType = cellstream[0]; - switch (cellType) - { - case XMU_POLY_LINE: - case XMU_POLYGON: + if (cellType == XMU_POLYHEDRON) { - int numPoints = cellstream[1]; - const int* cellPoints = cellstream + 2; - int idx1 = cellPoints[a_edgeIdx]; - int idx2 = cellPoints[(a_edgeIdx + 1) % numPoints]; - edge = XmEdge(idx1, idx2); - break; + iGetUniquePointsFromPolyhedronSingleCellstream(cellstream, a_cellPoints); + return true; } - - case XMU_POLYHEDRON: + else if (cellType == XMU_PIXEL) { - boost::container::flat_set cellEdges; - int currIdx = 2; - GetUniqueEdgesFromPolyhedronCellstream(cellstream, streamLength, cellEdges, currIdx); - - if (a_edgeIdx < (int)cellEdges.size()) - edge = *(cellEdges.begin() + a_edgeIdx); - break; + a_cellPoints.push_back(cellstream[2]); + a_cellPoints.push_back(cellstream[3]); + a_cellPoints.push_back(cellstream[5]); + a_cellPoints.push_back(cellstream[4]); + return true; } - - default: + else { - const VecEdge& edgeTable = iGetEdgeOffsetTable(cellType); - if (a_edgeIdx < (int)edgeTable.size()) - { - const XmEdge& edgeOffset = edgeTable[a_edgeIdx]; - const int* cellPoints = cellstream + 2; - int idx1 = cellPoints[edgeOffset.GetFirst()]; - int idx2 = cellPoints[edgeOffset.GetSecond()]; - edge = XmEdge(idx1, idx2); - } - break; - } + a_cellPoints.assign(cellstream.begin() + 2, cellstream.end()); + return true; } } + return false; +} // XmUGrid::Impl::GetCellPoints - return edge; -} // XmUGrid::Impl::GetCellEdge //------------------------------------------------------------------------------ -/// \brief Get the index of the adjacent cells (that shares the same cell edge) +/// \brief Get locations of cell points. /// \param[in] a_cellIdx the index of the cell -/// \param[in] a_edgeIdx the index of the edge -/// \return a vector of cell indices of the adjacent cells +/// \param[out] a_cellLocations The locations of the cell points //------------------------------------------------------------------------------ -VecInt XmUGrid::Impl::GetCellEdgeAdjacentCells(int a_cellIdx, int a_edgeIdx) const +void XmUGrid::Impl::GetCellLocations(int a_cellIdx, VecPt3d& a_cellLocations) const { - VecInt adjacentCells; - GetCellEdgeAdjacentCells(a_cellIdx, a_edgeIdx, adjacentCells); - return adjacentCells; -} // XmUGrid::Impl::GetCellEdgeAdjacentCells + VecInt ptIdxs = GetCellPoints(a_cellIdx); + a_cellLocations = GetPointsLocations(ptIdxs); +} // XmUGrid::Impl::GetCellPoints + //------------------------------------------------------------------------------ -/// \brief Get the index of the adjacent cells (that shares the same cell edge) +/// \brief Get the cell type of a specified cell. /// \param[in] a_cellIdx the index of the cell -/// \param[in] a_edgeIdx the index of the edge -/// \param[out] a_adjacentCellIdxs a vector of cell indices of the adjacent -/// cells +/// \return The type of the specified cell or -1 if invalid. //------------------------------------------------------------------------------ -void XmUGrid::Impl::GetCellEdgeAdjacentCells(int a_cellIdx, - int a_edgeIdx, - VecInt& a_adjacentCellIdxs) const +XmUGridCellType XmUGrid::Impl::GetCellType(int a_cellIdx) const { - a_adjacentCellIdxs.clear(); - VecInt adjacentCells; - VecInt cellNeighbors = GetCellAdjacentCells(a_cellIdx); - if (cellNeighbors.empty()) - { - return; - } - if (a_edgeIdx < 0 || a_edgeIdx >= GetCellEdgeCount(a_cellIdx)) - { - return; - } - - XmEdge currEdge, neighborEdge; - currEdge = GetCellEdge(a_cellIdx, a_edgeIdx); - for (int cellNeighbor : cellNeighbors) - { - for (int k(0); k < GetCellEdgeCount(cellNeighbor); k++) - { - neighborEdge = GetCellEdge(cellNeighbor, k); - if (currEdge.IsEquivalent(neighborEdge)) - { - a_adjacentCellIdxs.push_back(cellNeighbor); - } - } - } -} // XmUGrid::Impl::GetCellEdgeAdjacentCells + return m_cellStreamCache.GetCellType(a_cellIdx); +} // XmUGrid::Impl::GetCellType //------------------------------------------------------------------------------ -/// \brief Get the index of the adjacent cells (that shares the same cell edge) +/// \brief Count all number of the cells with each dimension (0, 1, 2, 3) +/// \return the count of dimensions of all of the cells of the ugrid +//------------------------------------------------------------------------------ +std::vector XmUGrid::Impl::GetDimensionCounts() const +{ + return m_cellStreamCache.GetDimensionCounts(); +} // XmUGrid::Impl::GetDimensionCounts +//------------------------------------------------------------------------------ +/// \brief Get the dimension of the specified cell. /// \param[in] a_cellIdx the index of the cell -/// \param[in] a_edgeIdx the index of the edge -/// \return index of the adjacent cell (or -1 if none). +/// \return the dimension of the cells or -1 if invalid index or invalid +/// dimension //------------------------------------------------------------------------------ -int XmUGrid::Impl::GetCell2dEdgeAdjacentCell(int a_cellIdx, int a_edgeIdx) const +int XmUGrid::Impl::GetCellDimension(int a_cellIdx) const { - if (GetCellDimension(a_cellIdx) != 2) - return -1; - - XmEdge edge = GetCellEdge(a_cellIdx, a_edgeIdx); - int pointIdx = edge.GetFirst(); - - int adjCellCount = GetPointAdjacentCellCount(pointIdx); - for (int cellIdx = 0; cellIdx < adjCellCount; cellIdx++) - { - int adjacentCell = m_pointsToCells[m_pointIdxToPointsToCells[pointIdx] + cellIdx + 1]; - if (adjacentCell == a_cellIdx) - continue; + return m_cellStreamCache.GetCellDimension(a_cellIdx); +} // XmUGrid::Impl::GetCellDimension - int edgeCount = GetCellEdgeCount(adjacentCell); - for (int adjacentEdgeIdx = 0; adjacentEdgeIdx < edgeCount; ++adjacentEdgeIdx) - { - XmEdge adjacentEdge = GetCellEdge(adjacentCell, adjacentEdgeIdx); - if (adjacentEdge.IsEquivalent(edge)) - return adjacentCell; - } +//------------------------------------------------------------------------------ +/// \brief Get the extents of the given cell. +/// \param[in] a_cellIdx The cell index to get the extents of. +/// \param[out] a_min The minimum location. +/// \param[out] a_max The maximum location. +//------------------------------------------------------------------------------ +void XmUGrid::Impl::GetCellExtents(int a_cellIdx, Pt3d& a_min, Pt3d& a_max) const +{ + VecPt3d pts; + GetCellLocations(a_cellIdx, pts); + iGetExtentsFromPoints(pts, a_min, a_max); +} // XmUGrid::Impl::GetCellExtents +//------------------------------------------------------------------------------ +/// \brief Get cell stream vector for the entire UGrid. +/// \return constant reference to the cell stream vector +//------------------------------------------------------------------------------ +const VecInt& XmUGrid::Impl::GetCellstream() const +{ + return m_cellStreamCache.GetCellstream(); +} // XmUGrid::Impl::GetCellstream +//------------------------------------------------------------------------------ +/// \brief Set the ugrid cells for the entire UGrid using a cell stream. +/// \param[in] a_cellstream cells defined as follows: +/// Hexahedrons, polygons, quads, triangles etc: +/// Cell type (ElemTypeEnum), number of points, point numbers. +/// Generally 0-based, CCW, bottom, then top. Not true +/// for pixel or voxel. +/// Polyhedrons: +/// Cell type, number of faces, [num points in face, +/// point numbers (0-based, CCW when looking in)] repeated +/// for each face. +/// \return if successfully set +//------------------------------------------------------------------------------ +bool XmUGrid::Impl::SetCellstream(const VecInt& a_cellstream) +{ + if (m_cellStreamCache.SetCellstream(a_cellstream)) + { + SetModified(); + return true; } - - return -1; -} // XmUGrid::Impl::GetCell2dEdgeAdjacentCell + return false; +} // XmUGrid::Impl::SetCellstream //------------------------------------------------------------------------------ -/// \brief Get the indices of the adjacent cells (that shares the same cell edge) -/// \param[in] a_edge the edge (a pair of point indices) -/// \param[out] a_adjacentCellIdxs a vector of cell indices of the adjacent cells +/// \brief Get cell stream vector for a single cell. +/// \param[in] a_cellIdx the index of the cell +/// \param[in] a_cellstream The cellstream of the cell +/// @see SetCellstream for more detail on cell stream definitions. +/// \return whether it was successfull or not //------------------------------------------------------------------------------ -void XmUGrid::Impl::GetEdgeAdjacentCells(const XmEdge& a_edge, VecInt& a_adjacentCellIdxs) const +bool XmUGrid::Impl::GetCellCellstream(int a_cellIdx, VecInt& a_cellstream) const { - a_adjacentCellIdxs.clear(); - GetPointsAdjacentCells(a_edge.GetFirst(), a_edge.GetSecond(), a_adjacentCellIdxs); -} // XmUGrid::Impl::GetEdgeAdjacentCells + return m_cellStreamCache.GetCellCellstream(a_cellIdx, a_cellstream); +} // XmUGrid::Impl::GetCellCellstream //------------------------------------------------------------------------------ -/// \brief Get the index of the adjacent cells (that shares the same cell edge) -/// \param[in] a_edge the edge (a pair of point indices) -/// \return a vector of cell indices of the adjacent cells +/// \brief Internal function to get start of cell stream for an individual cell. +/// Returns nullptr and zero length for invalid cell index. +/// \param[in] a_cellIdx the index of the cell. +/// \return The cell stream. //------------------------------------------------------------------------------ -VecInt XmUGrid::Impl::GetEdgeAdjacentCells(const XmEdge& a_edge) const +Span XmUGrid::Impl::GetCellCellstream(int a_cellIdx) const { - VecInt adjacentCellIdxs; - GetEdgeAdjacentCells(a_edge, adjacentCellIdxs); - return adjacentCellIdxs; -} // XmUGrid::Impl::GetEdgeAdjacentCells + return m_cellStreamCache.GetCellCellstream(a_cellIdx); +} // XmUGrid::Impl::GetCellCellstream //------------------------------------------------------------------------------ -/// \brief Get the Edges of a cell. -/// \param[in] a_cellIdx the cells to whom the edges belong -/// \return a vector of edges +/// \brief Get beginning index of cell in cell stream. +/// \param[in] a_cellIdx the index of the cell +/// \return The index of the cell in the cell stream. //------------------------------------------------------------------------------ -std::vector XmUGrid::Impl::GetCellEdges(int a_cellIdx) const +int XmUGrid::Impl::GetCellCellstreamIndex(int a_cellIdx) const { - std::vector edges; - GetCellEdges(a_cellIdx, edges); - return edges; -} // XmUGrid::Impl::GetCellEdges + return m_cellStreamCache.GetCellCellstreamIndex(a_cellIdx); +} // XmUGrid::GetCellCellstreamIndex //------------------------------------------------------------------------------ -/// \brief Get the Edges of a cell. -/// \param[in] a_cellIdx the cells to whom the edges belong -/// \param[out] a_edges a vector of edges (organized in std::pairs) +/// \brief Get the cells neighboring a cell (cells associated with any of it's points) +/// \param[in] a_cellIdx the index of the cell +/// \return vector of cell indices //------------------------------------------------------------------------------ -void XmUGrid::Impl::GetCellEdges(int a_cellIdx, std::vector& a_edges) const +VecInt XmUGrid::Impl::GetCellAdjacentCells(int a_cellIdx) const { - a_edges.clear(); - int numEdges = GetCellEdgeCount(a_cellIdx); - for (int i = 0; i < numEdges; ++i) - { - a_edges.push_back(GetCellEdge(a_cellIdx, i)); - } -} // XmUGrid::Impl::GetCellEdges - + VecInt neighbors; + GetCellAdjacentCells(a_cellIdx, neighbors); + return neighbors; +} // XmUGrid::Impl::GetCellAdjacentCells //------------------------------------------------------------------------------ -/// \brief Given a point gets point indices attached to the point by an edge. -/// \param[in] a_pointIdx The point to get adjacent points from. -/// \param[out] a_edgePoints The indices of the adjacent points. +/// \brief Get the cells neighboring a cell (cells associated with any of it's points) +/// \param[in] a_cellIdx the index of the cell +/// \param[out] a_cellNeighbors vector of cell indices //------------------------------------------------------------------------------ -void XmUGrid::Impl::GetPointAdjacentPoints(int a_pointIdx, VecInt& a_edgePoints) const +void XmUGrid::Impl::GetCellAdjacentCells(int a_cellIdx, VecInt& a_cellNeighbors) const { - a_edgePoints.clear(); - VecInt associatedCells = GetPointAdjacentCells(a_pointIdx); - if (associatedCells.empty()) - { - return; - } - for (int associatedCell : associatedCells) + a_cellNeighbors.clear(); + VecInt pointsOfCell; + pointsOfCell = GetCellPoints(a_cellIdx); + for (int pointIdx : pointsOfCell) { - for (int j = 0; j < GetCellEdgeCount(associatedCell); ++j) + VecInt associatedCells = GetPointAdjacentCells(pointIdx); + for (int associatedCell : associatedCells) { - XmEdge temp = GetCellEdge(associatedCell, j); - if (temp.GetFirst() == a_pointIdx) - { - a_edgePoints.push_back(temp.GetSecond()); - } - else if (temp.GetSecond() == a_pointIdx) + if (associatedCell != a_cellIdx) { - a_edgePoints.push_back(temp.GetFirst()); + bool found(false); + for (int a_cellNeighbor : a_cellNeighbors) + { + if (a_cellNeighbor == associatedCell) + { + found = true; + } + } + if (!found) + { + a_cellNeighbors.push_back(associatedCell); + } } } } - std::sort(a_edgePoints.begin(), a_edgePoints.end()); - auto it = std::unique(a_edgePoints.begin(), a_edgePoints.end()); - a_edgePoints.erase(it, a_edgePoints.end()); -} // XmUGrid::Impl::GetPointAdjacentPoints - -//------------------------------------------------------------------------------ -/// \brief Given a point gets point locations attached to the point by an edge. -/// \param[in] a_pointIdx The point to get attached point from. -/// \param[out] a_edgePoints A vector of points attached across edges. -//------------------------------------------------------------------------------ -void XmUGrid::Impl::GetPointAdjacentLocations(int a_pointIdx, VecPt3d& a_edgePoints) const -{ - VecInt edgePtIdxs; - GetPointAdjacentPoints(a_pointIdx, edgePtIdxs); - a_edgePoints = GetPointsLocations(edgePtIdxs); -} // XmUGrid::Impl::GetPointAdjacentLocations - -// Faces +} // XmUGrid::Impl::GetCellAdjacentCells //------------------------------------------------------------------------------ -/// \brief Get the number of cell faces for given cell. -/// \param[in] a_cellIdx the index of the cell -/// \return the count of cell faces +/// \brief Get a plan view polygon of a specified cell +/// \param[in] a_cellIdx The index of the cell. +/// \param[out] a_polygon Vector of Pt3d that is the plan view polygon. +/// \return False if the cell index does not exist or if the cell is not 2 or 3 +/// dimensional. //------------------------------------------------------------------------------ -int XmUGrid::Impl::GetCell3dFaceCount(int a_cellIdx) const +bool XmUGrid::Impl::GetCellPlanViewPolygon(int a_cellIdx, VecPt3d& a_polygon) const { - if (m_useCache) - { - CalculateCacheValues(); - int faceCount = IsValidCellIdx(a_cellIdx) ? m_numberOfFaces[a_cellIdx] : -1; - return faceCount; - } + a_polygon.clear(); + if (!IsValidCellIdx(a_cellIdx)) + return false; + int dimension = GetCellDimension(a_cellIdx); + if (dimension == 3) + return GetPlanViewPolygon3d(a_cellIdx, a_polygon); + else if (dimension == 2) + return GetPlanViewPolygon2d(a_cellIdx, a_polygon); else - { - return GetCell3dFaceCountNoCache(a_cellIdx); - } -} // XmUGrid::Impl::GetCell3dFaceCount + return false; +} // XmUGrid::Impl::GetCellPlanViewPolygon //------------------------------------------------------------------------------ -/// \brief Get the number of face points for a given cell and face. -/// \param[in] a_cellIdx The cell -/// \param[in] a_faceIdx The face -/// \return The number of face points or -1 if invalid face or cell index. +/// \brief Get the centroid location of a cell. +/// \param[in] a_cellIdx The index of the cell. +/// \param[out] a_centroid The location of the cell centroid. +/// \return False if the cell index does not exist. //------------------------------------------------------------------------------ -int XmUGrid::Impl::GetCell3dFacePointCount(int a_cellIdx, int a_faceIdx) const +bool XmUGrid::Impl::GetCellCentroid(int a_cellIdx, Pt3d& a_centroid) const { - if (!IsValidCellIdx(a_cellIdx)) + VecPt3d polyPoints; + bool retVal = true; + Pt3d centroid = Pt3d(0.0, 0.0, 0.0); + if (GetCellPlanViewPolygon(a_cellIdx, polyPoints)) { - return -1; + centroid = gmComputePolygonCentroid(polyPoints); } - - int cellType = GetCellType(a_cellIdx); - if (cellType == XMU_POLYHEDRON) + else if (!IsValidCellIdx(a_cellIdx)) { - const int* cellstream; - int length; - GetCellCellstream(a_cellIdx, &cellstream, length); - if (cellstream != nullptr) - { - auto currItem = cellstream; - currItem++; // skip cell type - int numFaces = *currItem++; - for (int faceIdx = 0; faceIdx < numFaces; ++faceIdx) - { - int numFacePoints = *currItem++; - if (a_faceIdx == faceIdx) - { - return numFacePoints; - } - currItem += numFacePoints; - } - } + retVal = false; } else { - const VecInt2d& faceTable = iGetFaceOffsetTable(cellType); - if (a_faceIdx >= 0 && a_faceIdx < (int)faceTable.size()) + VecPt3d cellPoints; + GetCellLocations(a_cellIdx, cellPoints); + for (Pt3d& pt : cellPoints) { - return (int)faceTable[a_faceIdx].size(); + centroid += pt; } + centroid /= (double)cellPoints.size(); } - - return 0; -} // XmUGrid::Impl::GetCell3dFacePointCount - -//------------------------------------------------------------------------------ -/// \brief Get the cell face for given cell and face index. -/// \param[in] a_cellIdx the index of the cell -/// \param[in] a_faceIdx the face index of the cell -/// \return a vector of point indices for the face index of the cell -//------------------------------------------------------------------------------ -VecInt XmUGrid::Impl::GetCell3dFacePoints(int a_cellIdx, int a_faceIdx) const -{ - VecInt facePoints; - GetCell3dFacePoints(a_cellIdx, a_faceIdx, facePoints); - return facePoints; -} // XmUGrid::Impl::GetCell3dFacePoints + a_centroid = centroid; + return retVal; +} // XmUGrid::Impl::GetCellCentroid //------------------------------------------------------------------------------ -/// \brief Get the cell face for given cell and face index. +/// \brief Determine whether a cell is valid after a point is moved. /// \param[in] a_cellIdx the index of the cell -/// \param[in] a_faceIdx the face index of the cell -/// \param[out] a_facePtIdxs a vector of point indices for the face +/// \param[in] a_changedPtIdx index of the point to be changed +/// \param[in] a_newPosition location the point is to be moved to +/// \return whether the cell is valid //------------------------------------------------------------------------------ -void XmUGrid::Impl::GetCell3dFacePoints(int a_cellIdx, int a_faceIdx, VecInt& a_facePtIdxs) const +bool XmUGrid::Impl::IsCellValidWithPointChange(int a_cellIdx, + int a_changedPtIdx, + const Pt3d& a_newPosition) const { - a_facePtIdxs.clear(); - if (a_faceIdx < 0) - { - return; - } - - const int* cellstream = nullptr; - int streamLength; - GetCellCellstream(a_cellIdx, &cellstream, streamLength); - if (!cellstream) - return; - - int cellType = *cellstream++; - if (cellType == XMU_POLYHEDRON) + if (GetCellDimension(a_cellIdx) == 2) { - int numFaces = *cellstream++; - for (int faceIdx = 0; faceIdx < numFaces; ++faceIdx) - { - int numFacePoints = *cellstream++; - if (faceIdx == a_faceIdx) - { - a_facePtIdxs.assign(cellstream, cellstream + numFacePoints); - break; - } - cellstream += numFacePoints; - } + std::vector edges = GetCellEdges(a_cellIdx); + return !DoEdgesCrossWithPointChange(a_changedPtIdx, a_newPosition, edges); } - else + else if (GetCellDimension(a_cellIdx) == 3) { - const VecInt2d& faceTable = iGetFaceOffsetTable(cellType); - if (a_faceIdx < (int)faceTable.size()) + // Go through the affected faces and if they are not a side face + // check if the edges cross each other + VecInt2d faces = GetCell3dFacesPoints(a_cellIdx); + for (auto& face : faces) { - if (cellstream) + bool faceIsAffected = false; + for (int ptIdx : face) { - const int* cellPoints = ++cellstream; - const VecInt& faceOffsets = faceTable[a_faceIdx]; - for (auto offset : faceOffsets) + if (ptIdx == a_changedPtIdx) { - int idx = cellPoints[offset]; - a_facePtIdxs.push_back(idx); + faceIsAffected = true; } } + if (faceIsAffected) + { + std::vector edges; + iGetEdgesOfFace(face, edges); + if (DoEdgesCrossWithPointChange(a_changedPtIdx, a_newPosition, edges)) + return false; + } } } -} // XmUGrid::Impl::GetCell3dFacePoints + return true; +} // XmUGrid::Impl::IsCellValidWithPointChange //------------------------------------------------------------------------------ -/// \brief Get the faces of a cell. -/// \param[in] a_cellIdx the cells to whom the faces belong -/// \return a vector of faces, which is a vector of point indices +/// \brief Determine if a cell index is valid. +/// \param[in] a_cellIdx the index of the cell +/// \return True if the cell index is valid. //------------------------------------------------------------------------------ -VecInt2d XmUGrid::Impl::GetCell3dFacesPoints(int a_cellIdx) const +bool XmUGrid::Impl::IsValidCellIdx(int a_cellIdx) const { - int numFaces = GetCell3dFaceCount(a_cellIdx); - - VecInt2d faces(numFaces); - for (int i(0); i < numFaces; i++) - { - faces[i] = GetCell3dFacePoints(a_cellIdx, i); - } - return faces; -} // XmUGrid::Impl::GetCell3dFacesPoints + return a_cellIdx >= 0 && a_cellIdx < GetCellCount(); +} // XmUGrid::Impl::IsValidCellIdx //------------------------------------------------------------------------------ -/// \brief Get the cell face neighbors for given cell and face index. -/// \param[in] a_cellIdx the index of the cell -/// \param[in] a_faceIdx the face index of the cell -/// \return a cell index of the neighbor +/// \brief Determine whether adjacent cells are valid after a point is moved. +/// \param[in] a_changedPtIdx index of the point to be changed +/// \param[in] a_newPosition location the point is to be moved to +/// \return whether the change is valid //------------------------------------------------------------------------------ -int XmUGrid::Impl::GetCell3dFaceAdjacentCell(int a_cellIdx, int a_faceIdx) const +bool XmUGrid::Impl::IsValidPointChange(int a_changedPtIdx, const Pt3d& a_newPosition) const { - if (m_useCache) + bool validChange = false; + if (a_changedPtIdx >= 0 && a_changedPtIdx < GetPointCount()) { - CalculateCacheValues(); - int faceNeighbor = -1; - int numFaces = GetCell3dFaceCount(a_cellIdx); - if (a_faceIdx >= 0 && a_faceIdx < numFaces) + validChange = true; + VecInt affectedCells = GetPointAdjacentCells(a_changedPtIdx); + for (int affectedCell : affectedCells) { - int faceOffset = m_cellFaceOffset[a_cellIdx]; - if (faceOffset >= 0) + if (!IsCellValidWithPointChange(affectedCell, a_changedPtIdx, a_newPosition)) { - faceOffset += a_faceIdx; - faceNeighbor = m_faceNeighbor[faceOffset]; - if (faceNeighbor == NEEDS_CALCULATION) - { - faceNeighbor = GetCell3dFaceAdjacentCellNoCache(a_cellIdx, a_faceIdx); - m_faceNeighbor[faceOffset] = faceNeighbor; - } - else - { - faceNeighbor = m_faceNeighbor[faceOffset]; - } + validChange = false; + break; } } - return faceNeighbor; - } - else - { - return GetCell3dFaceAdjacentCellNoCache(a_cellIdx, a_faceIdx); } -} // XmUGrid::Impl::GetCell3dFaceAdjacentCell + return validChange; +} // XmUGrid::Impl::IsValidPointChange + //------------------------------------------------------------------------------ -/// \brief Get the cell face neighbors for given cell and face index. -/// \param[in] a_cellIdx the index of the cell -/// \param[in] a_faceIdx the face index of the cell -/// \param[in] a_neighborCell the index of the neighboring cell -/// \param[in] a_neighborFace the face index of the neighboring cell adjacent -/// to the given face -/// \return a cell index of the neighbor +/// \brief Gets whether or not edges cross with a point change +/// \param[out] a_edges The edges to check +/// \param[in] a_changedPtIdx index of the point to be changed +/// \param[in] a_newPosition location the point is to be moved to +/// \return whether the edges cross //------------------------------------------------------------------------------ -bool XmUGrid::Impl::GetCell3dFaceAdjacentCell(int a_cellIdx, - int a_faceIdx, - int& a_neighborCell, - int& a_neighborFace) const +bool XmUGrid::Impl::DoEdgesCrossWithPointChange(int a_changedPtIdx, + const Pt3d& a_newPosition, + const std::vector& a_edges) const { - a_neighborCell = GetCell3dFaceAdjacentCell(a_cellIdx, a_faceIdx); - if (a_neighborCell < 0) + std::vector> changedEdges; + std::vector> unChangedEdges; + for (const auto& edge : a_edges) { - a_neighborCell = a_neighborFace = -1; - return false; + if (edge.GetFirst() == a_changedPtIdx) + { + changedEdges.emplace_back(a_newPosition, GetPointLocation(edge.GetSecond())); + } + else if (edge.GetSecond() == a_changedPtIdx) + { + changedEdges.emplace_back(GetPointLocation(edge.GetFirst()), a_newPosition); + } + else + { + unChangedEdges.emplace_back(GetPointLocation(edge.GetFirst()), + GetPointLocation(edge.GetSecond())); + } } - VecInt cellFacePts = GetCell3dFacePoints(a_cellIdx, a_faceIdx); - for (int faceIdx(0); faceIdx < GetCell3dFaceCount(a_neighborCell); faceIdx++) + for (const auto& changedEdge : changedEdges) { - VecInt curCellFacePts = GetCell3dFacePoints(a_neighborCell, faceIdx); - if (cellFacePts.size() == curCellFacePts.size()) + for (const auto& unChangedEdge : unChangedEdges) { - int cellFacePtsIdx(0); - for (; cellFacePtsIdx < (int)cellFacePts.size(); ++cellFacePtsIdx) - { - bool found = false; - for (int curCellFacePt : curCellFacePts) - { - if (cellFacePts[cellFacePtsIdx] == curCellFacePt) - { - found = true; - break; - } - } - if (!found) - break; // This face is not the face we are l. - } - if (cellFacePtsIdx == (int)cellFacePts.size()) // Every point is found - { - a_neighborFace = faceIdx; - a_neighborFace = faceIdx; + if (gmLinesCross(changedEdge.first, changedEdge.second, unChangedEdge.first, + unChangedEdge.second)) return true; - } } } - a_neighborCell = a_neighborFace = -1; return false; -} // XmUGrid::Impl::GetCell3dFaceAdjacentCell +} // XmUGrid::Impl::DoEdgesCrossWithPointChange +// Edges //------------------------------------------------------------------------------ -/// \brief Get the orientation of the face of a vertically prismatic cell. +/// \brief Get the number of edges for a cell. /// \param[in] a_cellIdx the index of the cell -/// \param[in] a_faceIdx the face index of the cell -/// \return The orientation of the face (TOP, BOTTOM, SIDE, UNKNOWN). +/// \return the count of cell edges //------------------------------------------------------------------------------ -XmUGridFaceOrientation XmUGrid::Impl::GetCell3dFaceOrientation(int a_cellIdx, int a_faceIdx) const +int XmUGrid::Impl::GetCellEdgeCount(int a_cellIdx) const { - if (m_useCache) - { - CalculateCacheValues(); - int faceOrientation = XMU_ORIENTATION_UNKNOWN; - int numFaces = GetCell3dFaceCount(a_cellIdx); - if (a_faceIdx >= 0 && a_faceIdx < numFaces) - { - int faceOffset = m_cellFaceOffset[a_cellIdx]; - if (faceOffset >= 0) - { - faceOffset += a_faceIdx; - faceOrientation = m_faceOrientation[faceOffset]; - if (faceOrientation == NEEDS_CALCULATION) - { - faceOrientation = (int)GetCell3dFaceOrientationNoCache(a_cellIdx, a_faceIdx); - m_faceOrientation[faceOffset] = faceOrientation; - } - else - { - faceOrientation = m_faceOrientation[faceOffset]; - } - } - } - return (XmUGridFaceOrientation)faceOrientation; - } - else + int cellType(GetCellType(a_cellIdx)); + switch (cellType) { - return GetCell3dFaceOrientationNoCache(a_cellIdx, a_faceIdx); + // invalid + case XMU_INVALID_CELL_TYPE: + return -1; + + case XMU_POLY_LINE: + return m_cellStreamCache.GetNumberOfItemsForCell(a_cellIdx) - 1; + + case XMU_POLYGON: + return m_cellStreamCache.GetNumberOfItemsForCell(a_cellIdx); + + case XMU_POLYHEDRON: + return GetNumberOfPolyhedronEdges(a_cellIdx); + + default: + // handled below + break; } -} // XmUGrid::Impl::GetCell3dFaceOrientation + const VecEdge& edgeTable = iGetEdgeOffsetTable(cellType); + if (!edgeTable.empty()) + return (int)edgeTable.size(); + return -1; +} // XmUGrid::Impl::GetCellEdgeCount //------------------------------------------------------------------------------ -/// \brief Find the orientation of a given 3D cell face. -/// \param[in] a_cellIdx The cell index. -/// \param[in] a_faceIdx The face index. -/// \return The orientation of the face (TOP, BOTTOM, SIDE, UNKNOWN). +/// \brief Get the points of a cell. +/// \param[in] a_cellIdx the index of the cell +/// \param[in] a_edgeIdx the index of the edge +/// \return a standard pair of point indices (which is an edge) //------------------------------------------------------------------------------ -XmUGridFaceOrientation XmUGrid::Impl::FaceOrientation(int a_cellIdx, int a_faceIdx) const +XmEdge XmUGrid::Impl::GetCellEdge(int a_cellIdx, int a_edgeIdx) const { - XmUGridFaceOrientation orientation = FaceOrientationWithFail(a_cellIdx, a_faceIdx); - if (orientation == XMU_ORIENTATION_UNKNOWN) + XmEdge edge; + if (a_edgeIdx < 0) + return edge; + + Span cellstream = GetCellCellstream(a_cellIdx); + if (!cellstream.empty()) { - // Is there any way to know if we reach this point? First face top and next - // bottom? What about multi-panel top with no cell above? - // XM_ASSERT(0); - orientation = VerticalOrientationFromOpposing(a_cellIdx, a_faceIdx); + if (cellstream.size() < 4) + return edge; + + int cellType = cellstream[0]; + switch (cellType) + { + case XMU_POLY_LINE: + case XMU_POLYGON: + { + int numPoints = cellstream[1]; + const int* cellPoints = cellstream.data() + 2; + int idx1 = cellPoints[a_edgeIdx]; + int idx2 = cellPoints[(a_edgeIdx + 1) % numPoints]; + edge = XmEdge(idx1, idx2); + break; + } + + case XMU_POLYHEDRON: + { + boost::container::flat_set cellEdges; + int currIdx = 2; + iGetUniqueEdgesFromPolyhedronCellstream(cellstream, cellEdges, currIdx); + + if (a_edgeIdx < (int)cellEdges.size()) + edge = *(cellEdges.begin() + a_edgeIdx); + break; + } + + default: + { + const VecEdge& edgeTable = iGetEdgeOffsetTable(cellType); + if (a_edgeIdx < (int)edgeTable.size()) + { + const XmEdge& edgeOffset = edgeTable[a_edgeIdx]; + const int* cellPoints = cellstream.data() + 2; + int idx1 = cellPoints[edgeOffset.GetFirst()]; + int idx2 = cellPoints[edgeOffset.GetSecond()]; + edge = XmEdge(idx1, idx2); + } + break; + } + } } - XM_ASSERT(orientation != XMU_ORIENTATION_UNKNOWN); - return orientation; -} // XmUGrid::Impl::FaceOrientation + return edge; +} // XmUGrid::Impl::GetCellEdge //------------------------------------------------------------------------------ -/// \brief Update internal links to navigate between associated points and -/// cells. +/// \brief Get the index of the adjacent cells (that shares the same cell edge) +/// \param[in] a_cellIdx the index of the cell +/// \param[in] a_edgeIdx the index of the edge +/// \return a vector of cell indices of the adjacent cells //------------------------------------------------------------------------------ -void XmUGrid::Impl::UpdateLinks() +VecInt XmUGrid::Impl::GetCellEdgeAdjacentCells(int a_cellIdx, int a_edgeIdx) const { - UpdateCellLinks(); - UpdatePointLinks(); -} // XmUGrid::Impl::UpdateLinks + VecInt adjacentCells; + GetCellEdgeAdjacentCells(a_cellIdx, a_edgeIdx, adjacentCells); + return adjacentCells; +} // XmUGrid::Impl::GetCellEdgeAdjacentCells //------------------------------------------------------------------------------ -/// \brief Update internal link from cells to cell stream index. +/// \brief Get the index of the adjacent cells (that shares the same cell edge) +/// \param[in] a_cellIdx the index of the cell +/// \param[in] a_edgeIdx the index of the edge +/// \param[out] a_adjacentCellIdxs a vector of cell indices of the adjacent +/// cells //------------------------------------------------------------------------------ -void XmUGrid::Impl::UpdateCellLinks() +void XmUGrid::Impl::GetCellEdgeAdjacentCells(int a_cellIdx, + int a_edgeIdx, + VecInt& a_adjacentCellIdxs) const { - m_cellIdxToStreamIdx.clear(); - - int currIdx = 0; - if (m_cellstream.empty()) + a_adjacentCellIdxs.clear(); + VecInt adjacentCells; + VecInt cellNeighbors = GetCellAdjacentCells(a_cellIdx); + if (cellNeighbors.empty()) { - m_cellIdxToStreamIdx.push_back(currIdx); return; } - - int numItems = (int)m_cellstream.size(); - while (currIdx < numItems) + if (a_edgeIdx < 0 || a_edgeIdx >= GetCellEdgeCount(a_cellIdx)) { - m_cellIdxToStreamIdx.push_back(currIdx); - - // get cell type - int cellType = m_cellstream[currIdx++]; - if (currIdx >= numItems) - return; - - // get the number of items - int numPoints = m_cellstream[currIdx++]; - if (currIdx >= numItems) - return; + return; + } - if (cellType == XMU_POLYHEDRON) + XmEdge currEdge, neighborEdge; + currEdge = GetCellEdge(a_cellIdx, a_edgeIdx); + for (int cellNeighbor : cellNeighbors) + { + for (int k(0); k < GetCellEdgeCount(cellNeighbor); k++) { - int numFaces = numPoints; - for (int faceIdx = 0; faceIdx < numFaces; ++faceIdx) + neighborEdge = GetCellEdge(cellNeighbor, k); + if (currEdge.IsEquivalent(neighborEdge)) { - int numFacePoints = m_cellstream[currIdx++]; - currIdx += numFacePoints; + a_adjacentCellIdxs.push_back(cellNeighbor); } } - else - { - currIdx += numPoints; - } } - - m_cellIdxToStreamIdx.push_back(currIdx); -} // XmUGrid::Impl::UpdateCellLinks +} // XmUGrid::Impl::GetCellEdgeAdjacentCells //------------------------------------------------------------------------------ -/// \brief Update internal links from points to associated cells. +/// \brief Get the index of the adjacent cells (that shares the same cell edge) +/// \param[in] a_cellIdx the index of the cell +/// \param[in] a_edgeIdx the index of the edge +/// \return index of the adjacent cell (or -1 if none). //------------------------------------------------------------------------------ -void XmUGrid::Impl::UpdatePointLinks() +int XmUGrid::Impl::GetCell2dEdgeAdjacentCell(int a_cellIdx, int a_edgeIdx) const { - m_pointIdxToPointsToCells.clear(); - m_pointIdxToPointsToCells.resize(m_locations.size(), 0); - m_pointsToCells.clear(); - - // get number of cells for each point - int numStreamItems = (int)m_cellstream.size(); - int cellIdx = 0; - int currIdx = 0; - VecInt pointLastUsedIdx(m_locations.size(), -1); - VecInt cellPoints; - while (currIdx < numStreamItems) - { - // get cell type - int cellType = m_cellstream[currIdx++]; - if (currIdx >= numStreamItems) - return; - - // get the number of items (points or faces depending on cell type) - int numCellItems = m_cellstream[currIdx++]; - if (currIdx >= numStreamItems) - return; - - if (cellType == XMU_POLYHEDRON) - { - GetUniquePointsFromPolyhedronCellstream(m_cellstream, numCellItems, currIdx, cellPoints, - pointLastUsedIdx); + if (GetCellDimension(a_cellIdx) != 2) + return -1; - // Deemed to be slower than flat set-- Left only as a warning to others! - // std::stable_sort(cellPoints.begin(), cellPoints.end()); - // auto uniqueEnd = std::unique(cellPoints.begin(), cellPoints.end()); + XmEdge edge = GetCellEdge(a_cellIdx, a_edgeIdx); + int pointIdx = edge.GetFirst(); - // flat set deemed to be slower than marking when a point was last used (pointLastUsedIdx) + int adjCellCount = GetPointAdjacentCellCount(pointIdx); + for (int cellIdx = 0; cellIdx < adjCellCount; cellIdx++) + { + int adjacentCell = m_cellStreamCache.GetPointAdjacentCell(pointIdx, cellIdx); + if (adjacentCell == a_cellIdx) + continue; - // for (auto pt = cellPoints.begin(); pt != uniqueEnd; ++pt) - for (auto cellPoint : cellPoints) - { - m_pointIdxToPointsToCells[cellPoint] += 1; - } - } - else + int edgeCount = GetCellEdgeCount(adjacentCell); + for (int adjacentEdgeIdx = 0; adjacentEdgeIdx < edgeCount; ++adjacentEdgeIdx) { - // iterate on points - int numPoints = numCellItems; - for (int ptIdx = 0; ptIdx < numPoints; ++ptIdx) - { - int pt = m_cellstream[currIdx++]; - m_pointIdxToPointsToCells[pt] += 1; - } + XmEdge adjacentEdge = GetCellEdge(adjacentCell, adjacentEdgeIdx); + if (adjacentEdge.IsEquivalent(edge)) + return adjacentCell; } - - ++cellIdx; } - m_pointIdxToPointsToCells.push_back(0); - // change array of counts to array of offsets - int currCount = 0; - for (auto& item : m_pointIdxToPointsToCells) + return -1; +} // XmUGrid::Impl::GetCell2dEdgeAdjacentCell +//------------------------------------------------------------------------------ +/// \brief Get the indices of the adjacent cells (that shares the same cell edge) +/// \param[in] a_edge the edge (a pair of point indices) +/// \param[out] a_adjacentCellIdxs a vector of cell indices of the adjacent cells +//------------------------------------------------------------------------------ +void XmUGrid::Impl::GetEdgeAdjacentCells(const XmEdge& a_edge, VecInt& a_adjacentCellIdxs) const +{ + a_adjacentCellIdxs.clear(); + GetPointsAdjacentCells(a_edge.GetFirst(), a_edge.GetSecond(), a_adjacentCellIdxs); +} // XmUGrid::Impl::GetEdgeAdjacentCells +//------------------------------------------------------------------------------ +/// \brief Get the index of the adjacent cells (that shares the same cell edge) +/// \param[in] a_edge the edge (a pair of point indices) +/// \return a vector of cell indices of the adjacent cells +//------------------------------------------------------------------------------ +VecInt XmUGrid::Impl::GetEdgeAdjacentCells(const XmEdge& a_edge) const +{ + VecInt adjacentCellIdxs; + GetEdgeAdjacentCells(a_edge, adjacentCellIdxs); + return adjacentCellIdxs; +} // XmUGrid::Impl::GetEdgeAdjacentCells +//------------------------------------------------------------------------------ +/// \brief Get the Edges of a cell. +/// \param[in] a_cellIdx the cells to whom the edges belong +/// \return a vector of edges +//------------------------------------------------------------------------------ +std::vector XmUGrid::Impl::GetCellEdges(int a_cellIdx) const +{ + std::vector edges; + GetCellEdges(a_cellIdx, edges); + return edges; +} // XmUGrid::Impl::GetCellEdges +//------------------------------------------------------------------------------ +/// \brief Get the Edges of a cell. +/// \param[in] a_cellIdx the cells to whom the edges belong +/// \param[out] a_edges a vector of edges (organized in std::pairs) +//------------------------------------------------------------------------------ +void XmUGrid::Impl::GetCellEdges(int a_cellIdx, std::vector& a_edges) const +{ + a_edges.clear(); + int numEdges = GetCellEdgeCount(a_cellIdx); + for (int i = 0; i < numEdges; ++i) { - int count = item + 1; - item = currCount; - currCount += count; + a_edges.push_back(GetCellEdge(a_cellIdx, i)); } +} // XmUGrid::Impl::GetCellEdges - cellIdx = 0; - currIdx = 0; - m_pointsToCells.resize(currCount, 0); - std::fill(pointLastUsedIdx.begin(), pointLastUsedIdx.end(), -1); - while (currIdx < numStreamItems) +//------------------------------------------------------------------------------ +/// \brief Given a point gets point indices attached to the point by an edge. +/// \param[in] a_pointIdx The point to get adjacent points from. +/// \param[out] a_edgePoints The indices of the adjacent points. +//------------------------------------------------------------------------------ +void XmUGrid::Impl::GetPointAdjacentPoints(int a_pointIdx, VecInt& a_edgePoints) const +{ + a_edgePoints.clear(); + VecInt associatedCells = GetPointAdjacentCells(a_pointIdx); + if (associatedCells.empty()) { - // get cell type - int cellType = m_cellstream[currIdx++]; - if (currIdx >= numStreamItems) - return; - - // get the number of items (points or faces depending on cell type) - int numCellItems = m_cellstream[currIdx++]; - if (currIdx >= numStreamItems) - return; - - if (cellType == XMU_POLYHEDRON) + return; + } + for (int associatedCell : associatedCells) + { + for (int j = 0; j < GetCellEdgeCount(associatedCell); ++j) { - GetUniquePointsFromPolyhedronCellstream(m_cellstream, numCellItems, currIdx, cellPoints, - pointLastUsedIdx); - - // Deemed to be slower than flat set-- Left only as a warning to others! - // std::stable_sort(cellPoints.begin(), cellPoints.end()); - // auto uniqueEnd = std::unique(cellPoints.begin(), cellPoints.end()); - - // for (auto pt = cellPoints.begin(); pt != uniqueEnd; ++pt) - for (auto cellPoint : cellPoints) + XmEdge temp = GetCellEdge(associatedCell, j); + if (temp.GetFirst() == a_pointIdx) { - int countIdx = m_pointIdxToPointsToCells[cellPoint]; - int& count = m_pointsToCells[countIdx]; // point's cell count - ++count; // incrementing point's cell count - int updateIdx = countIdx + count; - m_pointsToCells[updateIdx] = cellIdx; + a_edgePoints.push_back(temp.GetSecond()); } - } - else - { - // iterate on points - int numPoints = numCellItems; - for (int ptIdx = 0; ptIdx < numPoints; ++ptIdx) + else if (temp.GetSecond() == a_pointIdx) { - int pt = m_cellstream[currIdx++]; - int countIdx = m_pointIdxToPointsToCells[pt]; - int& count = m_pointsToCells[countIdx]; // point's cell count - ++count; // incrementing point's cell count - int updateIdx = countIdx + count; - m_pointsToCells[updateIdx] = cellIdx; + a_edgePoints.push_back(temp.GetFirst()); } } + } + std::sort(a_edgePoints.begin(), a_edgePoints.end()); + auto it = std::unique(a_edgePoints.begin(), a_edgePoints.end()); + a_edgePoints.erase(it, a_edgePoints.end()); +} // XmUGrid::Impl::GetPointAdjacentPoints + +//------------------------------------------------------------------------------ +/// \brief Given a point gets point locations attached to the point by an edge. +/// \param[in] a_pointIdx The point to get attached point from. +/// \param[out] a_edgePoints A vector of points attached across edges. +//------------------------------------------------------------------------------ +void XmUGrid::Impl::GetPointAdjacentLocations(int a_pointIdx, VecPt3d& a_edgePoints) const +{ + VecInt edgePtIdxs; + GetPointAdjacentPoints(a_pointIdx, edgePtIdxs); + a_edgePoints = GetPointsLocations(edgePtIdxs); +} // XmUGrid::Impl::GetPointAdjacentLocations - ++cellIdx; +// Faces +//------------------------------------------------------------------------------ +/// \brief Get the number of cell faces for given cell. +/// \param[in] a_cellIdx the index of the cell +/// \return the count of cell faces +//------------------------------------------------------------------------------ +int XmUGrid::Impl::GetCell3dFaceCount(int a_cellIdx) const +{ + if (m_useCache) + { + SetupFaceCacheValues(); + int faceCount = IsValidCellIdx(a_cellIdx) ? m_numberOfFaces[a_cellIdx] : -1; + return faceCount; + } + else + { + return GetCell3dFaceCountNoCache(a_cellIdx); } -} // XmUGrid::Impl::UpdatePointLinks +} // XmUGrid::Impl::GetCell3dFaceCount + //------------------------------------------------------------------------------ -/// \brief Get the dimension given the cell type (0d, 1d, 2d, or 3d). -/// \param[in] a_cellType the cell type -/// \return the dimension of the cell type +/// \brief Get the number of face points for a given cell and face. +/// \param[in] a_cellIdx The cell +/// \param[in] a_faceIdx The face +/// \return The number of face points or -1 if invalid face or cell index. //------------------------------------------------------------------------------ -int XmUGrid::Impl::DimensionFromCellType(XmUGridCellType a_cellType) +int XmUGrid::Impl::GetCell3dFacePointCount(int a_cellIdx, int a_faceIdx) const { - switch (a_cellType) + if (!IsValidCellIdx(a_cellIdx)) { - // invalid - case XMU_INVALID_CELL_TYPE: return -1; - break; - - // 0D - case XMU_EMPTY_CELL: - case XMU_VERTEX: - case XMU_POLY_VERTEX: - case XMU_CONVEX_POINT_SET: // Special class of cells formed by convex group of - // points - return 0; - break; - - // 1D - case XMU_LINE: - case XMU_POLY_LINE: - - case XMU_QUADRATIC_EDGE: - case XMU_PARAMETRIC_CURVE: - case XMU_HIGHER_ORDER_EDGE: - case XMU_CUBIC_LINE: // Cubic, isoparametric cell - return 1; - break; - // 2D - case XMU_TRIANGLE: - case XMU_TRIANGLE_STRIP: - case XMU_POLYGON: - case XMU_PIXEL: - case XMU_QUAD: - - case XMU_QUADRATIC_TRIANGLE: - case XMU_QUADRATIC_QUAD: - case XMU_BIQUADRATIC_QUAD: - case XMU_QUADRATIC_LINEAR_QUAD: - case XMU_BIQUADRATIC_TRIANGLE: - case XMU_QUADRATIC_POLYGON: - - case XMU_PARAMETRIC_SURFACE: - case XMU_PARAMETRIC_TRI_SURFACE: - case XMU_PARAMETRIC_QUAD_SURFACE: - - case XMU_HIGHER_ORDER_TRIANGLE: - case XMU_HIGHER_ORDER_QUAD: - case XMU_HIGHER_ORDER_POLYGON: - - return 2; - break; - - // assuming the rest are 3D - default: - return 3; - break; -#if 0 // Remaining definitions - XMU_TETRA = 10, - XMU_VOXEL = 11, - XMU_HEXAHEDRON = 12, - XMU_WEDGE = 13, - XMU_PYRAMID = 14, - XMU_PENTAGONAL_PRISM = 15, - XMU_HEXAGONAL_PRISM = 16, - - // Quadratic, isoparametric cells - XMU_QUADRATIC_TETRA = 24, - XMU_QUADRATIC_HEXAHEDRON = 25, - XMU_QUADRATIC_WEDGE = 26, - XMU_QUADRATIC_PYRAMID = 27, - XMU_TRIQUADRATIC_HEXAHEDRON = 29, - XMU_QUADRATIC_LINEAR_WEDGE = 31, - XMU_BIQUADRATIC_QUADRATIC_WEDGE = 32, - XMU_BIQUADRATIC_QUADRATIC_HEXAHEDRON = 33, - - // Polyhedron cell (consisting of polygonal faces) - XMU_POLYHEDRON = 42, - - // Higher order cells in parametric form - XMU_PARAMETRIC_TETRA_REGION = 55, - XMU_PARAMETRIC_HEX_REGION = 56, + } - // Higher order cells - XMU_HIGHER_ORDER_TETRAHEDRON = 64, - XMU_HIGHER_ORDER_WEDGE = 65, - XMU_HIGHER_ORDER_PYRAMID = 66, - XMU_HIGHER_ORDER_HEXAHEDRON = 67, -#endif + int cellType = GetCellType(a_cellIdx); + if (cellType == XMU_POLYHEDRON) + { + Span cellstream = GetCellCellstream(a_cellIdx); + if (!cellstream.empty()) + { + auto currItem = cellstream.data(); + currItem++; // skip cell type + int numFaces = *currItem++; + for (int faceIdx = 0; faceIdx < numFaces; ++faceIdx) + { + int numFacePoints = *currItem++; + if (a_faceIdx == faceIdx) + { + return numFacePoints; + } + currItem += numFacePoints; + } + } + } + else + { + const VecInt2d& faceTable = iGetFaceOffsetTable(cellType); + if (a_faceIdx >= 0 && a_faceIdx < (int)faceTable.size()) + { + return (int)faceTable[a_faceIdx].size(); + } } - return -1; -} // XmUGrid::Impl::DimensionFromCellType + return 0; +} // XmUGrid::Impl::GetCell3dFacePointCount //------------------------------------------------------------------------------ -/// \brief Get number of items given cell. For polyhedron number of items is -/// number of faces. For other cell types it is number of points. +/// \brief Get the cell face for given cell and face index. /// \param[in] a_cellIdx the index of the cell -/// \return the number of faces for polyhedron or number of points. +/// \param[in] a_faceIdx the face index of the cell +/// \return a vector of point indices for the face index of the cell //------------------------------------------------------------------------------ -int XmUGrid::Impl::GetNumberOfItemsForCell(int a_cellIdx) const +VecInt XmUGrid::Impl::GetCell3dFacePoints(int a_cellIdx, int a_faceIdx) const { - int cellType = GetCellType(a_cellIdx); - if (cellType == XMU_INVALID_CELL_TYPE) + VecInt facePoints; + GetCell3dFacePoints(a_cellIdx, a_faceIdx, facePoints); + return facePoints; +} // XmUGrid::Impl::GetCell3dFacePoints +//------------------------------------------------------------------------------ +/// \brief Get the cell face for given cell and face index. +/// \param[in] a_cellIdx the index of the cell +/// \param[in] a_faceIdx the face index of the cell +/// \param[out] a_facePtIdxs a vector of point indices for the face +//------------------------------------------------------------------------------ +void XmUGrid::Impl::GetCell3dFacePoints(int a_cellIdx, int a_faceIdx, VecInt& a_facePtIdxs) const +{ + if (a_cellIdx < 0 || a_faceIdx < 0 || a_cellIdx >= GetCellCount() || + a_faceIdx >= GetCell3dFaceCount(a_cellIdx)) { - return 0; + a_facePtIdxs.clear(); + return; } - else if (cellType != XMU_POLYHEDRON) + + int facePointCount = GetCell3dFacePointCount(a_cellIdx, a_faceIdx); + if (facePointCount <= 0) { - // Number of points is right after cell type - int startIndex = m_cellIdxToStreamIdx[a_cellIdx]; - return m_cellstream[startIndex + 1]; - } - else - { // Polyhedron case, same value is number of faces - int startIndex = m_cellIdxToStreamIdx[a_cellIdx]; - return m_cellstream[startIndex + 1]; + a_facePtIdxs.clear(); + return; } -} // XmUGrid::Impl::GetNumberOfItemsForCell + + a_facePtIdxs.resize(facePointCount); + CopyCell3dFacePoints(a_cellIdx, a_faceIdx, a_facePtIdxs.data()); +} // XmUGrid::Impl::GetCell3dFacePoints //------------------------------------------------------------------------------ -/// \brief Internal function to get start of cell stream for an individual cell. -/// Returns nullptr and zero length for invalid cell index. +/// \brief Get the cell face for given cell and face index. /// \param[in] a_cellIdx the index of the cell -/// \param[out] a_start pointer to the start of the stream for the cell -/// \param[out] a_length the length of the stream for the cell +/// \param[in] a_faceIdx the face index of the cell +/// \param[out] a_facePtIdxs location to write point indices for the face //------------------------------------------------------------------------------ -void XmUGrid::Impl::GetCellCellstream(int a_cellIdx, const int** a_start, int& a_length) const +void XmUGrid::Impl::CopyCell3dFacePoints(int a_cellIdx, int a_faceIdx, int* a_facePtIdxs) const { - int cellType = GetCellType(a_cellIdx); - if (cellType == XMU_INVALID_CELL_TYPE) + if (a_faceIdx < 0) + return; + + Span cellCellStream = GetCellCellstream(a_cellIdx); + if (cellCellStream.empty()) + return; + + auto cellstream = cellCellStream.data(); + int cellType = *cellstream++; + if (cellType == XMU_POLYHEDRON) { - *a_start = nullptr; - a_length = 0; + int numFaces = *cellstream++; + for (int faceIdx = 0; faceIdx < numFaces; ++faceIdx) + { + int numFacePoints = *cellstream++; + if (faceIdx == a_faceIdx) + { + std::copy(cellstream, cellstream + numFacePoints, a_facePtIdxs); + break; + } + cellstream += numFacePoints; + } } else { - int startIndex = m_cellIdxToStreamIdx[a_cellIdx]; - int nextCellIndex = m_cellIdxToStreamIdx[a_cellIdx + 1]; - *a_start = &m_cellstream[startIndex]; - a_length = nextCellIndex - startIndex; + const VecInt2d& faceTable = iGetFaceOffsetTable(cellType); + if (a_faceIdx < (int)faceTable.size()) + { + if (cellstream) + { + const int* cellPoints = ++cellstream; + const VecInt& faceOffsets = faceTable[a_faceIdx]; + int* copyTo = a_facePtIdxs; + for (auto offset : faceOffsets) + { + int idx = cellPoints[offset]; + *copyTo++ = idx; + } + } + } } -} // XmUGrid::Impl::GetCellCellstream +} // XmUGrid::Impl::CopyCell3dFacePoints //------------------------------------------------------------------------------ -/// \brief Get the number of edges for a polyhedron cell. +/// \brief Get the cell face for given cell and face index. Uses the cache. /// \param[in] a_cellIdx the index of the cell -/// \return the number of edges +/// \param[in] a_faceIdx the face index of the cell +/// \return The face points. //------------------------------------------------------------------------------ -int XmUGrid::Impl::GetNumberOfPolyhedronEdges(int a_cellIdx) const +void XmUGrid::Impl::GetCell3dFacePoints(int a_cellIdx, int a_faceIdx, Span& a_facePtIdxs) const +{ + SetupFaceCacheValues(); + int faceOffset = m_cellFaceOffset[a_cellIdx]; + int pointOffset = m_facePointOffset[faceOffset + a_faceIdx]; + int pointCount = m_facePointOffset[faceOffset + a_faceIdx + 1] - pointOffset; + a_facePtIdxs = Span(&m_facePoints[pointOffset], pointCount); +} // XmUGrid::Impl::GetCachedCell3dFacePoints +//------------------------------------------------------------------------------ +/// \brief Get the faces of a cell. +/// \param[in] a_cellIdx the cells to whom the faces belong +/// \return a vector of faces, which is a vector of point indices +//------------------------------------------------------------------------------ +VecInt2d XmUGrid::Impl::GetCell3dFacesPoints(int a_cellIdx) const { - const int* cellstream; - int streamLength; - GetCellCellstream(a_cellIdx, &cellstream, streamLength); - if (cellstream && streamLength > 0 && cellstream[0] == XMU_POLYHEDRON) + int numFaces = GetCell3dFaceCount(a_cellIdx); + + VecInt2d faces(numFaces); + for (int i(0); i < numFaces; i++) { - boost::container::flat_set edges; - int currItem = 2; - while (currItem < streamLength) - { - GetUniqueEdgesFromPolyhedronCellstream(cellstream, streamLength, edges, currItem); - } - return (int)edges.size(); + faces[i] = GetCell3dFacePoints(a_cellIdx, i); } - return 0; -} // XmUGrid::Impl::GetNumberOfPolyhedronEdges + return faces; +} // XmUGrid::Impl::GetCell3dFacesPoints //------------------------------------------------------------------------------ -/// \brief Get the unique points in a flat set -/// \param[in] a_cellstream the UGrid cell stream -/// \param[in] a_numCellItems the number of cell faces in the polyhedron -/// \param[in] a_currIdx the index of the cell stream; this should reference -/// the number of points in the first face. This variable will be updated -/// to the cell type of the next cell. -/// \param[out] a_uniqueGetCellPoints the unique points of the polyhedron -/// \param[out] a_pointLastUsedIdx a list of integers representing when a point -/// was last found while reading the stream. -/// \note: This function does NOT verify cellstream size!! This function -/// needs to be efficient! +/// \brief Get the cell face neighbors for given cell and face index. +/// \param[in] a_cellIdx the index of the cell +/// \param[in] a_faceIdx the face index of the cell +/// \param[in] a_neighborCell the index of the neighboring cell +/// \param[in] a_neighborFace the face index of the neighboring cell adjacent +/// to the given face +/// \return a cell index of the neighbor //------------------------------------------------------------------------------ -void XmUGrid::Impl::GetUniquePointsFromPolyhedronCellstream(const VecInt& a_cellstream, - int a_numCellItems, - int& a_currIdx, - VecInt& a_uniqueGetCellPoints, - VecInt& a_pointLastUsedIdx) +bool XmUGrid::Impl::GetCell3dFaceAdjacentCell(int a_cellIdx, + int a_faceIdx, + int& a_neighborCell, + int& a_neighborFace) const { - int stable = a_currIdx; - a_uniqueGetCellPoints.clear(); - int numFaces = a_numCellItems; - for (int faceIdx = 0; faceIdx < numFaces; ++faceIdx) + if (m_useCache) + return GetCell3dFaceAdjacentCellWithCache(a_cellIdx, a_faceIdx, a_neighborCell, a_neighborFace); + return GetCell3dFaceAdjacentCellNoCache(a_cellIdx, a_faceIdx, a_neighborCell,a_neighborFace); +} // XmUGrid::Impl::GetCell3dFaceAdjacentCell +//------------------------------------------------------------------------------ +/// \brief Get the cell face neighbors for given cell and face index. +/// \param[in] a_cellIdx the index of the cell +/// \param[in] a_faceIdx the face index of the cell +/// \param[in] a_neighborCell the index of the neighboring cell +/// \param[in] a_neighborFace the face index of the neighboring cell adjacent +/// to the given face +/// \return a cell index of the neighbor +//------------------------------------------------------------------------------ +bool XmUGrid::Impl::GetCell3dFaceAdjacentCellWithCache(int a_cellIdx, + int a_faceIdx, + int& a_neighborCell, + int& a_neighborFace) const +{ + a_neighborCell = GetCell3dFaceAdjacentCell(a_cellIdx, a_faceIdx); + if (a_neighborCell < 0) { - int numFacePoints = a_cellstream[a_currIdx++]; - for (int ptIdx = 0; ptIdx < numFacePoints; ++ptIdx) + a_neighborCell = a_neighborFace = -1; + return false; + } + + // find matching neighbor face + Span cellFacePts; + GetCell3dFacePoints(a_cellIdx, a_faceIdx, cellFacePts); + for (int faceIdx(0); faceIdx < GetCell3dFaceCount(a_neighborCell); faceIdx++) + { + Span curCellFacePts; + GetCell3dFacePoints(a_neighborCell, faceIdx, curCellFacePts); + if (cellFacePts.size() == curCellFacePts.size()) { - int pt = a_cellstream[a_currIdx++]; - if (a_pointLastUsedIdx[pt] < stable) + if (std::is_permutation(cellFacePts.begin(), cellFacePts.end(), curCellFacePts.begin())) { - a_pointLastUsedIdx[pt] = stable; - a_uniqueGetCellPoints.push_back(pt); + a_neighborFace = faceIdx; + a_neighborFace = faceIdx; + return true; } } } -} // XmUGrid::Impl::GetUniquePointsFromPolyhedronCellstream + a_neighborCell = -1; + a_neighborFace = -1; + return false; +} // XmUGrid::Impl::GetCell3dFaceAdjacentCell //------------------------------------------------------------------------------ -/// \brief Get the unique points for cell stream of a single polyhedron cell. -/// \param[in] a_cellstream a single cell stream that is a polyhedron type -/// \param[out] a_cellPoints the points of the cell -/// \return false if invalid stream +/// \brief Get the cell face neighbors for given cell and face index. +/// \param[in] a_cellIdx the index of the cell +/// \param[in] a_faceIdx the face index of the cell +/// \param[in] a_neighborCell the index of the neighboring cell +/// \param[in] a_neighborFace the face index of the neighboring cell adjacent +/// to the given face +/// \return a cell index of the neighbor //------------------------------------------------------------------------------ -bool XmUGrid::Impl::GetUniquePointsFromPolyhedronSingleCellstream(const VecInt& a_cellstream, - VecInt& a_cellPoints) +bool XmUGrid::Impl::GetCell3dFaceAdjacentCellNoCache(int a_cellIdx, + int a_faceIdx, + int& a_neighborCell, + int& a_neighborFace) const { - a_cellPoints.clear(); - if (a_cellstream.size() < 2) + a_neighborCell = GetCell3dFaceAdjacentCell(a_cellIdx, a_faceIdx); + if (a_neighborCell < 0) { + a_neighborCell = a_neighborFace = -1; return false; } - auto currItem = a_cellstream.begin(); - int cellType = *currItem++; - if (cellType != XMU_POLYHEDRON) + // find matching neighbor face + VecInt cellFacePts; + GetCell3dFacePoints(a_cellIdx, a_faceIdx, cellFacePts); + for (int faceIdx(0); faceIdx < GetCell3dFaceCount(a_neighborCell); faceIdx++) { - return false; + VecInt curCellFacePts; + GetCell3dFacePoints(a_neighborCell, faceIdx, curCellFacePts); + if (cellFacePts.size() == curCellFacePts.size()) + { + if (std::is_permutation(cellFacePts.begin(), cellFacePts.end(), curCellFacePts.begin())) + { + a_neighborFace = faceIdx; + a_neighborFace = faceIdx; + return true; + } + } } - - int numFaces = *currItem++; - for (int faceIdx = 0; faceIdx < numFaces; ++faceIdx) + a_neighborCell = -1; + a_neighborFace = -1; + return false; +} // XmUGrid::Impl::GetCell3dFaceAdjacentCell +//------------------------------------------------------------------------------ +/// \brief Get the orientation of the face of a vertically prismatic cell. +/// \param[in] a_cellIdx the index of the cell +/// \param[in] a_faceIdx the face index of the cell +/// \return The orientation of the face (TOP, BOTTOM, SIDE, UNKNOWN). +//------------------------------------------------------------------------------ +XmUGridFaceOrientation XmUGrid::Impl::GetCell3dFaceOrientation(int a_cellIdx, int a_faceIdx) const +{ + if (m_useCache) { - int numFacePoints = *currItem++; - for (int ptIdx = 0; ptIdx < numFacePoints; ++ptIdx) + SetupFaceCacheValues(); + int faceOrientation = XMU_ORIENTATION_UNKNOWN; + int numFaces = GetCell3dFaceCount(a_cellIdx); + if (a_faceIdx >= 0 && a_faceIdx < numFaces) { - int pt = *currItem++; - if (std::find(a_cellPoints.begin(), a_cellPoints.end(), pt) == a_cellPoints.end()) + int faceOffset = m_cellFaceOffset[a_cellIdx]; + if (faceOffset >= 0) { - a_cellPoints.push_back(pt); + faceOffset += a_faceIdx; + faceOrientation = m_faceOrientation[faceOffset]; + if (faceOrientation == NEEDS_CALCULATION) + { + faceOrientation = (int)GetCell3dFaceOrientationNoCache(a_cellIdx, a_faceIdx); + m_faceOrientation[faceOffset] = faceOrientation; + } + else + { + faceOrientation = m_faceOrientation[faceOffset]; + } } } + return (XmUGridFaceOrientation)faceOrientation; + } + else + { + return GetCell3dFaceOrientationNoCache(a_cellIdx, a_faceIdx); + } +} // XmUGrid::Impl::GetCell3dFaceOrientation +//------------------------------------------------------------------------------ +/// \brief Find the orientation of a given 3D cell face. +/// \param[in] a_cellIdx The cell index. +/// \param[in] a_faceIdx The face index. +/// \return The orientation of the face (TOP, BOTTOM, SIDE, UNKNOWN). +//------------------------------------------------------------------------------ +XmUGridFaceOrientation XmUGrid::Impl::FaceOrientation(int a_cellIdx, int a_faceIdx) const +{ + XmUGridFaceOrientation orientation = FaceOrientationWithFail(a_cellIdx, a_faceIdx); + if (orientation == XMU_ORIENTATION_UNKNOWN) + { + // Is there any way to know if we reach this point? First face top and next + // bottom? What about multi-panel top with no cell above? + // XM_ASSERT(0); + orientation = VerticalOrientationFromOpposing(a_cellIdx, a_faceIdx); } - return true; -} // XmUGrid::Impl::GetUniquePointsFromPolyhedronSingleCellstream + XM_ASSERT(orientation != XMU_ORIENTATION_UNKNOWN); + return orientation; +} // XmUGrid::Impl::FaceOrientation //------------------------------------------------------------------------------ -/// \brief Get the unique edges in a flat set for a given polyhedron -/// \param[in] a_start the UGrid cell stream (integer pointer) -/// \param[in] a_length the length of the cell stream -/// \param[in] a_cellEdges Unique cell edges of the polyhedron -/// \param[in] a_currIdx the index of the cell stream; this should reference -/// the number of points in the first face. This variable will be updated -/// to the cell type of the next cell. -/// \note: This function does NOT verify cellstream size!! This function -/// needs to be efficient! +/// \brief Get the number of edges for a polyhedron cell. +/// \param[in] a_cellIdx the index of the cell +/// \return the number of edges //------------------------------------------------------------------------------ -void XmUGrid::Impl::GetUniqueEdgesFromPolyhedronCellstream( - const int* a_start, - int& a_length, - boost::container::flat_set& a_cellEdges, - int& a_currIdx) +int XmUGrid::Impl::GetNumberOfPolyhedronEdges(int a_cellIdx) const { - int numFaces = a_start[1]; - for (int i(0); i < numFaces; i++) + Span cellstream = GetCellCellstream(a_cellIdx); + if (!cellstream.empty() && cellstream[0] == XMU_POLYHEDRON) { - int numPoints = a_start[a_currIdx++]; - for (int pointIdx = 0; pointIdx < numPoints; ++pointIdx) + boost::container::flat_set edges; + int currItem = 2; + while (currItem < cellstream.size()) { - int pt1Idx = pointIdx; - int pt2Idx = (pointIdx + 1) % numPoints; - // The % operator will reset the index back 0 so the final loop will - // provide the last and first point - int pt1 = a_start[a_currIdx + pt1Idx]; - int pt2 = a_start[a_currIdx + pt2Idx]; - // We want unique edges, so we add the lower point index first - if (pt1 < pt2) - a_cellEdges.insert(XmEdge(pt1, pt2)); - else - a_cellEdges.insert(XmEdge(pt2, pt1)); + iGetUniqueEdgesFromPolyhedronCellstream(cellstream, edges, currItem); } - a_currIdx += numPoints; + return static_cast(edges.size()); } -} // XmUGrid::Impl::GetUniqueEdgesFromPolyhedronCellstream + return 0; +} // XmUGrid::Impl::GetNumberOfPolyhedronEdges //------------------------------------------------------------------------------ /// \brief Get a plan view polygon of a specified 2D cell /// \param[in] a_cellIdx the index of the cell @@ -2551,6 +2938,9 @@ bool XmUGrid::Impl::GetPlanViewPolygon2d(int a_cellIdx, VecPt3d& a_polygon) cons //------------------------------------------------------------------------------ bool XmUGrid::Impl::GetPlanViewPolygon3d(int a_cellIdx, VecPt3d& a_polygon) const { + if (a_cellIdx < 0 || a_cellIdx >= GetCellCount()) + return false; + VecPt3d segments; if (GetCellXySegments(a_cellIdx, segments)) { @@ -2558,20 +2948,17 @@ bool XmUGrid::Impl::GetPlanViewPolygon3d(int a_cellIdx, VecPt3d& a_polygon) cons iMergeSegmentsToPoly(segments, a_polygon); return true; } - else + + // Non-prismatic cell + VecInt uniquePoints = GetCellPoints(a_cellIdx); + VecPt3d cellPoints(uniquePoints.size()); + for (int i(0); i < uniquePoints.size(); i++) { - // Non-prismatic cell - VecInt uniquePoints = GetCellPoints(a_cellIdx); - VecPt3d cellPoints(uniquePoints.size()); - for (int i(0); i < uniquePoints.size(); i++) - { - cellPoints[i] = GetPointLocation(uniquePoints[i]); - cellPoints[i].z = 0.0; // Insist that our z values are 0.0 for plan view - } - gmGetConvexHull(cellPoints, a_polygon, false); - return true; + cellPoints[i] = GetPointLocation(uniquePoints[i]); + cellPoints[i].z = 0.0; // Insist that our z values are 0.0 for plan view } - return false; + gmGetConvexHull(cellPoints, a_polygon, false); + return true; } // XmUGrid::Impl::GetPlanViewPolygon3d //------------------------------------------------------------------------------ /// \brief Get whether the cell face is of a side orientation. Only works for @@ -2621,22 +3008,6 @@ bool XmUGrid::Impl::GetCellXySegments(int a_cellIdx, VecPt3d& a_segments) const return foundSideFace; } // XmUGrid::Impl::GetCellXySegments -//------------------------------------------------------------------------------ -/// \brief Function to get the extents from a list of points. -/// \param[in] a_locations The point locations to get the extents of. -/// \param[out] a_min Minimum point location. -/// \param[out] a_max Maximum point location. -//------------------------------------------------------------------------------ -void XmUGrid::Impl::GetExtentsFromPoints(const VecPt3d& a_locations, Pt3d& a_min, Pt3d& a_max) const -{ - a_min.x = a_min.y = a_min.z = xms::XM_DBL_HIGHEST; - a_max.x = a_max.y = a_max.z = xms::XM_DBL_LOWEST; - for (const Pt3d& pt : a_locations) - { - gmAddToExtents(pt, a_min, a_max); - } -} // XmUGrid::Impl::GetExtentsFromPoints - //------------------------------------------------------------------------------ /// \brief Get the Xy locations of Face Points /// \param[in] a_cellIdx The index of the cells whose face points are wanted. @@ -2666,26 +3037,61 @@ bool XmUGrid::Impl::GetFaceXySegments(int a_cellIdx, int a_faceIdx, VecPt3d& a_s //------------------------------------------------------------------------------ /// \brief Calculate cached values for faster lookup. //------------------------------------------------------------------------------ -void XmUGrid::Impl::CalculateCacheValues() const +void XmUGrid::Impl::SetupFaceCacheValues() const { if (m_numberOfFaces.empty() && GetCellCount() != 0) { int cellCount = GetCellCount(); + // fill cache with number of faces and face offsets into other caches m_numberOfFaces.assign(cellCount, 0); m_cellFaceOffset.assign(cellCount + 1, 0); - int faceCount = 0; + int faceOffset = 0; for (int cellIdx = 0; cellIdx < cellCount; ++cellIdx) { int numberOfFaces = GetCell3dFaceCountNoCache(cellIdx); m_numberOfFaces[cellIdx] = numberOfFaces; - faceCount += numberOfFaces; - m_cellFaceOffset[cellIdx + 1] = faceCount; + faceOffset += numberOfFaces; + m_cellFaceOffset[cellIdx + 1] = faceOffset; + } + + // get initial array of face point offsets + m_facePointOffset.assign(faceOffset + 1, 0); + int totalFacePointCount = 0; + for (int cellIdx = 0; cellIdx < cellCount; ++cellIdx) + { + int numberOfFaces = m_numberOfFaces[cellIdx]; + for (int faceIdx = 0; faceIdx < numberOfFaces; ++faceIdx) + { + int faceCacheIndex = m_cellFaceOffset[cellIdx] + faceIdx; + int facePointCount = GetCell3dFacePointCount(cellIdx, faceIdx); + if (facePointCount < 0) + facePointCount = 0; + totalFacePointCount += facePointCount; + m_facePointOffset[faceCacheIndex + 1] = totalFacePointCount; + } + } + + // get initial array of face points + m_facePoints.assign(totalFacePointCount, 0); + for (int cellIdx = 0; cellIdx < cellCount; ++cellIdx) + { + int numberOfFaces = m_numberOfFaces[cellIdx]; + for (int faceIdx = 0; faceIdx < numberOfFaces; ++faceIdx) + { + int faceCacheIndex = m_cellFaceOffset[cellIdx] + faceIdx; + int facePointStart = m_facePointOffset[faceCacheIndex]; + int facePointCount = m_facePointOffset[faceCacheIndex + 1] - facePointStart; + if (facePointCount > 0) + { + CopyCell3dFacePoints(cellIdx, faceIdx, &m_facePoints[facePointStart]); + } + } } - m_faceOrientation.assign(faceCount, NEEDS_CALCULATION); - m_faceNeighbor.assign(faceCount, NEEDS_CALCULATION); + m_faceOrientation.assign(faceOffset, NEEDS_CALCULATION); + m_faceNeighborCell.assign(faceOffset, NEEDS_CALCULATION); } -} // XmUGrid::Impl::CalculateCacheValues +} // XmUGrid::Impl::SetupFaceCacheValues //------------------------------------------------------------------------------ /// \brief Clear cached so they will be recalculated. //------------------------------------------------------------------------------ @@ -2694,8 +3100,7 @@ void XmUGrid::Impl::ClearCacheValues() m_numberOfFaces.clear(); m_cellFaceOffset.clear(); m_faceOrientation.clear(); - m_faceNeighbor.clear(); - m_cellDimensionCounts.clear(); + m_faceNeighborCell.clear(); } // XmUGrid::Impl::ClearCacheValues //------------------------------------------------------------------------------ /// \brief Get the number of cell faces for given cell. @@ -2710,7 +3115,7 @@ int XmUGrid::Impl::GetCell3dFaceCountNoCache(int a_cellIdx) const int cellType(GetCellType(a_cellIdx)); if (cellType == XMU_POLYHEDRON) { - return GetNumberOfItemsForCell(a_cellIdx); + return m_cellStreamCache.GetNumberOfItemsForCell(a_cellIdx); } const VecInt2d& faceTable = iGetFaceOffsetTable(cellType); @@ -2722,6 +3127,18 @@ int XmUGrid::Impl::GetCell3dFaceCountNoCache(int a_cellIdx) const /// \param[in] a_faceIdx the face index of the cell /// \return a cell index of the neighbor //------------------------------------------------------------------------------ +int XmUGrid::Impl::GetCell3dFaceAdjacentCell(int a_cellIdx, int a_faceIdx) const +{ + if (m_useCache) + return GetCell3dFaceAdjacentCellWithCache(a_cellIdx, a_faceIdx); + return GetCell3dFaceAdjacentCellNoCache(a_cellIdx, a_faceIdx); +} // XmUGrid::Impl::GetCell3dFaceAdjacentCell +//------------------------------------------------------------------------------ +/// \brief Get the cell face neighbors for given cell and face index. +/// \param[in] a_cellIdx the index of the cell +/// \param[in] a_faceIdx the face index of the cell +/// \return a cell index of the neighbor +//------------------------------------------------------------------------------ int XmUGrid::Impl::GetCell3dFaceAdjacentCellNoCache(int a_cellIdx, int a_faceIdx) const { VecInt cellFace = GetCell3dFacePoints(a_cellIdx, a_faceIdx); @@ -2744,6 +3161,59 @@ int XmUGrid::Impl::GetCell3dFaceAdjacentCellNoCache(int a_cellIdx, int a_faceIdx return -1; } // XmUGrid::Impl::GetCell3dFaceAdjacentCellNoCache //------------------------------------------------------------------------------ +/// \brief Get the cell face neighbors for given cell and face index. +/// \param[in] a_cellIdx the index of the cell +/// \param[in] a_faceIdx the face index of the cell +/// \return a cell index of the neighbor +//------------------------------------------------------------------------------ +int XmUGrid::Impl::GetCell3dFaceAdjacentCellWithCache(int a_cellIdx, int a_faceIdx) const +{ + SetupFaceCacheValues(); + int faceNeighbor = -1; + int numFaces = GetCell3dFaceCount(a_cellIdx); + if (a_faceIdx >= 0 && a_faceIdx < numFaces) + { + int faceOffset = m_cellFaceOffset[a_cellIdx]; + if (faceOffset >= 0) + { + faceOffset += a_faceIdx; + faceNeighbor = m_faceNeighborCell[faceOffset]; + if (faceNeighbor != NEEDS_CALCULATION) + { + faceNeighbor = m_faceNeighborCell[faceOffset]; + return faceNeighbor; + } + } + + Span facePoints; + GetCell3dFacePoints(a_cellIdx, a_faceIdx, facePoints); + if (facePoints.empty()) + return -1; + + // find adjacent cell with a face that shares the same face points + Span pointCells; + GetPointAdjacentCells(facePoints[0], pointCells); + for (int checkCellIdx : pointCells) + { + if (checkCellIdx == a_cellIdx) + continue; + bool found = true; + for (int j = 1; found && j < facePoints.size(); ++j) + { + Span pointCells2; + GetPointAdjacentCells(facePoints[j], pointCells2); + found = std::find(pointCells2.begin(), pointCells2.end(), checkCellIdx) != pointCells2.end(); + } + if (found) + { + m_faceNeighborCell[faceOffset] = checkCellIdx; + return checkCellIdx; + } + } + } + return -1; +} // XmUGrid::Impl::GetCell3dFaceAdjacentCellWithCache +//------------------------------------------------------------------------------ /// \brief Get the orientation of the face of a vertically prismatic cell. /// \param[in] a_cellIdx the index of the cell /// \param[in] a_faceIdx the face index of the cell @@ -2883,8 +3353,6 @@ XmUGridFaceOrientation XmUGrid::Impl::FaceOrientationWithFail(int a_cellIdx, int XmUGridFaceOrientation XmUGrid::Impl::VerticalOrientationFromOpposing(int a_cellIdx, int a_faceIdx) const { - XmUGridFaceOrientation orientation = XMU_ORIENTATION_UNKNOWN; - // Assume only 1 top and bottom. Find other face that's top or bottom and // give answer as opposite. Otherwise first top other is bottom. bool firstUnknown = false; @@ -2892,27 +3360,25 @@ XmUGridFaceOrientation XmUGrid::Impl::VerticalOrientationFromOpposing(int a_cell bool foundBot = false; for (int face = 0; face < GetCell3dFaceCount(a_cellIdx); ++face) { - int orientation = FaceOrientationWithFail(a_cellIdx, face); - if (orientation == XMU_ORIENTATION_UNKNOWN) + int faceOrientation = FaceOrientationWithFail(a_cellIdx, face); + if (faceOrientation == XMU_ORIENTATION_UNKNOWN) { if (a_faceIdx == face) firstUnknown = true; } - else if (orientation == XMU_ORIENTATION_TOP) + else if (faceOrientation == XMU_ORIENTATION_TOP) foundTop = true; - else if (orientation == XMU_ORIENTATION_BOTTOM) + else if (faceOrientation == XMU_ORIENTATION_BOTTOM) foundBot = true; } if (foundTop) - orientation = XMU_ORIENTATION_BOTTOM; - else if (foundBot) - orientation = XMU_ORIENTATION_TOP; - else if (firstUnknown) - orientation = XMU_ORIENTATION_TOP; - else - orientation = XMU_ORIENTATION_BOTTOM; - return orientation; + return XMU_ORIENTATION_BOTTOM; + if (foundBot) + return XMU_ORIENTATION_TOP; + if (firstUnknown) + return XMU_ORIENTATION_TOP; + return XMU_ORIENTATION_BOTTOM; } // XmUGrid::Impl::VerticalOrientationFromOpposing //------------------------------------------------------------------------------ /// \brief Determines if a cell face is a vertical side face. @@ -2931,9 +3397,8 @@ bool XmUGrid::Impl::IsSideFace(int a_cellIdx, int a_faceIdx) const { int ptIdxLast = facePts.back(); Pt3d ptLast = GetPointXy0(ptIdxLast); - for (size_t facePtIdx = 0; facePtIdx < facePts.size(); ++facePtIdx) + for (int ptIdxCurr : facePts) { - int ptIdxCurr = facePts[facePtIdx]; Pt3d ptCurr = GetPointXy0(ptIdxCurr); // need to check for face that incorrectly has same point at beginning and ending // because some exist from old MODFLOW 6 models built from DISV package @@ -3011,7 +3476,7 @@ void XmUGrid::Impl::SetCellOrdering(XmUGridCellOrdering a_cellOrdering) m_cellOrdering = a_cellOrdering; if (m_useCache) { - m_faceOrientation.assign(m_faceNeighbor.size(), NEEDS_CALCULATION); + m_faceOrientation.assign(m_faceNeighborCell.size(), NEEDS_CALCULATION); } } // XmUGrid::Impl::SetCellOrdering //------------------------------------------------------------------------------ @@ -3249,26 +3714,7 @@ void XmUGrid::Swap(XmUGrid& a_xmUGrid) //------------------------------------------------------------------------------ bool XmUGrid::IsValidCellstream(const VecInt& a_cellstream, int a_points) { - if (a_cellstream.empty()) - { - return true; - } - - int currIdx = 0; - while (currIdx < a_cellstream.size()) - { - int step = iValidateCell(a_cellstream, a_points, currIdx); - if (step == BAD_CELL_STREAM) - { - return false; - } - else - { - currIdx += step; - } - } - - return true; + return iIsValidCellstream(a_cellstream, a_points); } // XmUGrid::IsValidCellstream //------------------------------------------------------------------------------ @@ -4026,7 +4472,7 @@ std::shared_ptr TEST_XmUBuildQuadUGrid(int a_rows, int a_cols, con { for (int j = 0; j < a_cols; ++j) { - points.push_back(Pt3d(j + a_origin.x, a_rows - i + a_origin.y)); + points.emplace_back(j + a_origin.x, a_rows - i + a_origin.y); } } @@ -4089,7 +4535,7 @@ std::shared_ptr TEST_XmUBuildDirectedHexahedronUgrid(int a_pointRo z = a_pointLayers - k - 1 + a_origin.z; else z = k + a_origin.z; - points.push_back(Pt3d(x, y, z)); + points.emplace_back(x, y, z); } } } @@ -4184,7 +4630,7 @@ std::shared_ptr TEST_XmUBuildPolyhedronUgrid(int a_rows, int a_col /// \param[in] a_rows number of rows in UGrid /// \param[in] a_cols number of columns in UGrid /// \param[in] a_lays number of layers in UGrid -/// \param[in] a_origin location wGetPointLocationd begins (min x, min y, min z) +/// \param[in] a_origin grid origin (min x, min y, min z) /// \return Returns the UGrid. //------------------------------------------------------------------------------ std::shared_ptr TEST_XmUBuildPolyhedronUgrid(int a_rows, @@ -4200,7 +4646,7 @@ std::shared_ptr TEST_XmUBuildPolyhedronUgrid(int a_rows, { for (int j = 0; j < a_cols; ++j) { - points.push_back(Pt3d(j + a_origin.x, a_rows - i + a_origin.y, a_lays - k + a_origin.z)); + points.emplace_back(j + a_origin.x, a_rows - i + a_origin.y, a_lays - k + a_origin.z); } } } @@ -5180,89 +5626,89 @@ void XmUGridUnitTests::testGetCell3dFacePoints() { // 2D Tests, include the bounds check VecInt emptyCellFaces = {}; - VecInt2d expectedCellFaces = {// XMU_QUAD - {}, - // XMU_PIXEL - {}, - // XMU_TRIANGLE - {}, - // XMU_POLYGON - {}, - // XMU_POLY_LINE - {}, - // XMU_LINE - {}}; - VecInt cellFaces; + VecInt2d expectedCellFacePoints = {// XMU_QUAD + {}, + // XMU_PIXEL + {}, + // XMU_TRIANGLE + {}, + // XMU_POLYGON + {}, + // XMU_POLY_LINE + {}, + // XMU_LINE + {}}; + VecInt cellFacePoints; // 2D Shapes std::shared_ptr ugrid2d = TEST_XmUGrid2dLinear(); TS_REQUIRE_NOT_NULL(ugrid2d); - cellFaces = ugrid2d->GetCell3dFacePoints(-1, 0); - TS_ASSERT_EQUALS(emptyCellFaces, cellFaces); - cellFaces = ugrid2d->GetCell3dFacePoints(0, -1); - TS_ASSERT_EQUALS(emptyCellFaces, cellFaces); + cellFacePoints = ugrid2d->GetCell3dFacePoints(-1, 0); + TS_ASSERT_EQUALS(emptyCellFaces, cellFacePoints); + cellFacePoints = ugrid2d->GetCell3dFacePoints(0, -1); + TS_ASSERT_EQUALS(emptyCellFaces, cellFacePoints); for (int i(0); i < ugrid2d->GetCellCount(); i++) { - cellFaces = ugrid2d->GetCell3dFacePoints(i, 0); - TS_ASSERT_EQUALS(expectedCellFaces[i], cellFaces); + cellFacePoints = ugrid2d->GetCell3dFacePoints(i, 0); + TS_ASSERT_EQUALS(expectedCellFacePoints[i], cellFacePoints); } - cellFaces = ugrid2d->GetCell3dFacePoints(ugrid2d->GetCellCount(), 0); - TS_ASSERT_EQUALS(emptyCellFaces, cellFaces); - cellFaces = ugrid2d->GetCell3dFacePoints(0, 1); - TS_ASSERT_EQUALS(emptyCellFaces, cellFaces); + cellFacePoints = ugrid2d->GetCell3dFacePoints(ugrid2d->GetCellCount(), 0); + TS_ASSERT_EQUALS(emptyCellFaces, cellFacePoints); + cellFacePoints = ugrid2d->GetCell3dFacePoints(0, 1); + TS_ASSERT_EQUALS(emptyCellFaces, cellFacePoints); // 3D Shapes std::shared_ptr ugrid3d = TEST_XmUGrid3dLinear(); TS_REQUIRE_NOT_NULL(ugrid3d); - expectedCellFaces = {// Tetra - {0, 1, 15}, - {1, 5, 15}, - {5, 0, 15}, - {0, 5, 1}, - // Voxel - {1, 16, 21, 6}, - {2, 7, 22, 17}, - {1, 2, 17, 16}, - {6, 21, 22, 7}, - {1, 6, 7, 2}, - {16, 17, 22, 21}, - // Hexahedron - {2, 17, 22, 7}, - {3, 8, 23, 18}, - {2, 3, 18, 17}, - {7, 22, 23, 8}, - {2, 7, 8, 3}, - {17, 18, 23, 22}, - // Polyhedron - {9, 8, 13, 14}, - {8, 9, 24, 23}, - {9, 14, 29, 24}, - {14, 13, 28, 29}, - {8, 13, 28, 23}, - {23, 24, 29, 28}, - // Wedge - {3, 4, 18}, - {8, 23, 9}, - {3, 8, 9, 4}, - {4, 9, 23, 18}, - {3, 18, 23, 8}, - // Pyramid - {5, 10, 11, 6}, - {5, 6, 20}, - {6, 11, 20}, - {11, 10, 20}, - {10, 5, 20}}; + expectedCellFacePoints = {// Tetra + {0, 1, 15}, + {1, 5, 15}, + {5, 0, 15}, + {0, 5, 1}, + // Voxel + {1, 16, 21, 6}, + {2, 7, 22, 17}, + {1, 2, 17, 16}, + {6, 21, 22, 7}, + {1, 6, 7, 2}, + {16, 17, 22, 21}, + // Hexahedron + {2, 17, 22, 7}, + {3, 8, 23, 18}, + {2, 3, 18, 17}, + {7, 22, 23, 8}, + {2, 7, 8, 3}, + {17, 18, 23, 22}, + // Polyhedron + {9, 8, 13, 14}, + {8, 9, 24, 23}, + {9, 14, 29, 24}, + {14, 13, 28, 29}, + {8, 13, 28, 23}, + {23, 24, 29, 28}, + // Wedge + {3, 4, 18}, + {8, 23, 9}, + {3, 8, 9, 4}, + {4, 9, 23, 18}, + {3, 18, 23, 8}, + // Pyramid + {5, 10, 11, 6}, + {5, 6, 20}, + {6, 11, 20}, + {11, 10, 20}, + {10, 5, 20}}; int currId(0); for (int i(0); i < ugrid3d->GetCellCount(); i++) { for (int j(0); j < ugrid3d->GetCell3dFaceCount(i); j++, currId++) { - cellFaces = ugrid3d->GetCell3dFacePoints(i, j); - TS_ASSERT_EQUALS(expectedCellFaces[currId], cellFaces); + cellFacePoints = ugrid3d->GetCell3dFacePoints(i, j); + TS_ASSERT_EQUALS(expectedCellFacePoints[currId], cellFacePoints); } } } // XmUGridUnitTests::testGetCell3dFacePoints