Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Feature/gridedit 766 improve mesh deletion speed #250

Merged
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
33 changes: 33 additions & 0 deletions libs/MeshKernel/include/MeshKernel/BoundingBox.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -106,6 +106,11 @@ namespace meshkernel
point.y >= m_lowerLeft.y && point.y <= m_upperRight.y;
}

/// @brief Checks if two bounding boxes overlaps
/// @param[in] boundingBox The input bounding box
/// @return True if the point if the this bounding box overlaps with another, false otherwise
bool Overlaps(const BoundingBox& boundingBox) const;

/// @brief Returns the lower left corner of the bounding box
/// @return The lower left corner of the bounding box
[[nodiscard]] auto& lowerLeft() const { return m_lowerLeft; }
Expand Down Expand Up @@ -141,6 +146,18 @@ namespace meshkernel
/// @brief Return the delta of the bounding box.
Vector Delta() const;

/// @brief Create a bounding box from two points
template <std::derived_from<Point> T>
static BoundingBox CreateBoundingBox(const T& first, const T& second)
{
const auto lowerLeftX = std::min(first.x, second.x);
const auto lowerLeftY = std::min(first.y, second.y);
const auto upperRightX = std::max(first.x, second.x);
const auto upperRightY = std::max(first.y, second.y);

return BoundingBox({lowerLeftX, lowerLeftY}, {upperRightX, upperRightY});
}

private:
Point m_lowerLeft; ///< The lower left corner of the bounding box
Point m_upperRight; ///< The upper right corner of the bounding box
Expand Down Expand Up @@ -203,3 +220,19 @@ inline meshkernel::Vector meshkernel::BoundingBox::Delta() const
{
return Vector(m_upperRight.x - m_lowerLeft.x, m_upperRight.y - m_lowerLeft.y);
}

inline bool meshkernel::BoundingBox::Overlaps(const BoundingBox& other) const
{

const auto& otherLowerleft = other.lowerLeft();
const auto& otherUpperRight = other.upperRight();
if (m_upperRight.x < otherLowerleft.x ||
otherUpperRight.x < m_lowerLeft.x ||
m_upperRight.y < otherLowerleft.y ||
otherUpperRight.y < m_lowerLeft.y)
{
return false;
}

return true;
}
10 changes: 10 additions & 0 deletions libs/MeshKernel/include/MeshKernel/Mesh2D.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -308,6 +308,16 @@ namespace meshkernel
/// it may be required to call Administrate after merging
static Mesh2D Merge(const Mesh2D& mesh1, const Mesh2D& mesh2);

/// @brief Get the mesh bounding box
///
/// @return The mesh bounding box
[[nodiscard]] BoundingBox GetBoundingBox() const;

/// @brief Get the bounding boxes of the mesh edges
///
/// @return The mesh edges bounding boxes
[[nodiscard]] std::vector<BoundingBox> GetEdgesBoundingBoxes() const;

private:
// orthogonalization
static constexpr double m_minimumEdgeLength = 1e-4; ///< Minimum edge length
Expand Down
5 changes: 5 additions & 0 deletions libs/MeshKernel/include/MeshKernel/Mesh2DIntersections.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -107,11 +107,13 @@ namespace meshkernel
/// @returns The intersection seed
std::tuple<UInt, UInt> GetIntersectionSeed(const Mesh2D& mesh,
const std::vector<Point>& polyLine,
const std::vector<BoundingBox>& polyLineBoundingBoxes,
const std::vector<bool>& vistedEdges) const;

/// @brief Gets the next edge intersection
/// @returns The intersection seed
std::tuple<bool, UInt, UInt, double, double, double> GetNextEdgeIntersection(const std::vector<Point>& polyLine,
const std::vector<BoundingBox>& polyLineBoundingBoxes,
UInt edgeIndex,
UInt firstIndex,
UInt secondIndex,
Expand All @@ -120,6 +122,7 @@ namespace meshkernel
/// @brief Gets the next edge intersection
/// @returns The intersection seed
void IntersectFaceEdges(const std::vector<Point>& polyLine,
const std::vector<BoundingBox>& polyLineBoundingBoxes,
const std::vector<double>& cumulativeLength,
UInt currentCrossingEdge,
UInt currentFaceIndex,
Expand Down Expand Up @@ -164,6 +167,8 @@ namespace meshkernel
std::vector<FaceMeshPolyLineIntersection> m_facesIntersectionsCache; ///< A cache for saving the local face intersections of one inner or outer
std::vector<EdgeMeshPolyLineIntersection> m_edgesIntersections; ///< A vector collecting all edge intersection results
std::vector<FaceMeshPolyLineIntersection> m_faceIntersections; ///< A vector collecting all face intersection results
BoundingBox m_meshBoundingBox; ///< The mesh bounding box
std::vector<BoundingBox> m_meshEdgesBoundingBoxes; ///< The mesh edges bounding boxes
static constexpr UInt maxSearchSegments = 1000; ///< max number of steps in polyline intersection algorithm
};

Expand Down
30 changes: 30 additions & 0 deletions libs/MeshKernel/src/Mesh2D.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -2019,3 +2019,33 @@ meshkernel::Mesh2D Mesh2D::Merge(const Mesh2D& mesh1, const Mesh2D& mesh2)

return mergedMesh;
}

meshkernel::BoundingBox Mesh2D::GetBoundingBox() const
{

Point lowerLeft(std::numeric_limits<double>::max(), std::numeric_limits<double>::max());
Point upperRight = -lowerLeft;

const auto numNodes = GetNumNodes();
lucacarniato marked this conversation as resolved.
Show resolved Hide resolved
for (UInt e = 0; e < numNodes; ++e)
{
lowerLeft.x = std::min(lowerLeft.x, m_nodes[e].x);
lowerLeft.y = std::min(lowerLeft.y, m_nodes[e].y);
upperRight.x = std::max(upperRight.x, m_nodes[e].x);
upperRight.y = std::max(upperRight.y, m_nodes[e].y);
}

return BoundingBox(lowerLeft, upperRight);
}

std::vector<meshkernel::BoundingBox> Mesh2D::GetEdgesBoundingBoxes() const
{
std::vector<BoundingBox> result;
result.reserve(GetNumEdges());
for (const auto& e : m_edges)
{
result.emplace_back(BoundingBox::CreateBoundingBox(m_nodes[e.first], m_nodes[e.second]));
}

return result;
}
37 changes: 35 additions & 2 deletions libs/MeshKernel/src/Mesh2DIntersections.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -27,10 +27,14 @@ Mesh2DIntersections::Mesh2DIntersections(Mesh2D& mesh) : m_mesh(mesh)
// Declare local caches
m_edgesIntersectionsCache.resize(mesh.GetNumEdges());
m_facesIntersectionsCache.resize(mesh.GetNumFaces());

m_meshBoundingBox = m_mesh.GetBoundingBox();
m_meshEdgesBoundingBoxes = m_mesh.GetEdgesBoundingBoxes();
}

std::tuple<UInt, UInt> Mesh2DIntersections::GetIntersectionSeed(const Mesh2D& mesh,
const std::vector<Point>& polyLine,
const std::vector<BoundingBox>& polyLineBoundingBoxes,
const std::vector<bool>& vistedEdges) const
{
UInt crossedSegmentIndex = constants::missing::uintValue;
Expand All @@ -48,6 +52,16 @@ std::tuple<UInt, UInt> Mesh2DIntersections::GetIntersectionSeed(const Mesh2D& me
continue;
}

if (!m_meshBoundingBox.Overlaps(polyLineBoundingBoxes[segmentIndex]))
{
continue;
}

if (!m_meshEdgesBoundingBoxes[edgeIndex].Overlaps(polyLineBoundingBoxes[segmentIndex]))
{
continue;
}

const auto [isEdgeCrossed,
intersectionPoint,
crossProductValue,
Expand Down Expand Up @@ -79,6 +93,7 @@ std::tuple<UInt, UInt> Mesh2DIntersections::GetIntersectionSeed(const Mesh2D& me
}

std::tuple<bool, UInt, UInt, double, double, double> Mesh2DIntersections::GetNextEdgeIntersection(const std::vector<Point>& polyLine,
const std::vector<BoundingBox>& polyLineBoundingBoxes,
UInt edgeIndex,
UInt firstIndex,
UInt secondIndex,
Expand Down Expand Up @@ -122,6 +137,18 @@ std::tuple<bool, UInt, UInt, double, double, double> Mesh2DIntersections::GetNex
firstIndex = firstIndex - 1;
}

if (!m_meshBoundingBox.Overlaps(polyLineBoundingBoxes[firstIndex]))
{
numSteps++;
continue;
}

if (!m_meshEdgesBoundingBoxes[edgeIndex].Overlaps(polyLineBoundingBoxes[firstIndex]))
{
numSteps++;
continue;
}

auto intersection = AreSegmentsCrossing(polyLine[firstIndex],
polyLine[secondIndex],
m_mesh.m_nodes[m_mesh.m_edges[edgeIndex].first],
Expand All @@ -141,6 +168,7 @@ std::tuple<bool, UInt, UInt, double, double, double> Mesh2DIntersections::GetNex
}

void Mesh2DIntersections::IntersectFaceEdges(const std::vector<Point>& polyLine,
const std::vector<BoundingBox>& polyLineBoundingBoxes,
const std::vector<double>& cumulativeLength,
UInt currentCrossingEdge,
UInt currentFaceIndex,
Expand Down Expand Up @@ -179,7 +207,7 @@ void Mesh2DIntersections::IntersectFaceEdges(const std::vector<Point>& polyLine,
secondIndex,
crossProductValue,
adimensionalPolylineSegmentDistance,
adimensionalEdgeDistance) = GetNextEdgeIntersection(polyLine, edgeIndex, segmentIndex, segmentIndex + 1, Direction::Forward);
adimensionalEdgeDistance) = GetNextEdgeIntersection(polyLine, polyLineBoundingBoxes, edgeIndex, segmentIndex, segmentIndex + 1, Direction::Forward);
}

if (!intersectionFound)
Expand All @@ -189,7 +217,7 @@ void Mesh2DIntersections::IntersectFaceEdges(const std::vector<Point>& polyLine,
secondIndex,
crossProductValue,
adimensionalPolylineSegmentDistance,
adimensionalEdgeDistance) = GetNextEdgeIntersection(polyLine, edgeIndex, segmentIndex, segmentIndex + 1, Direction::Backward);
adimensionalEdgeDistance) = GetNextEdgeIntersection(polyLine, polyLineBoundingBoxes, edgeIndex, segmentIndex, segmentIndex + 1, Direction::Backward);
}

// none of the polyline intersect the current edge
Expand Down Expand Up @@ -228,9 +256,12 @@ void Mesh2DIntersections::Compute(const std::vector<Point>& polyLine)
const auto polyLineSize = static_cast<UInt>(polyLine.size());

std::vector<double> cumulativeLength(polyLine.size(), 0.0);
std::vector<BoundingBox> polyLineBoundingBoxes(polyLine.size() - 1);
for (UInt i = 1; i < polyLineSize; ++i)
{
cumulativeLength[i] = cumulativeLength[i - 1] + ComputeDistance(polyLine[i], polyLine[i - 1], m_mesh.m_projection);

polyLineBoundingBoxes[i - 1] = BoundingBox::CreateBoundingBox(polyLine[i - 1], polyLine[i]);
}

std::queue<std::array<UInt, 2>> crossingEdges;
Expand All @@ -244,6 +275,7 @@ void Mesh2DIntersections::Compute(const std::vector<Point>& polyLine)
const auto [crossedEdgeIndex, crossedSegmentIndex] = GetIntersectionSeed(
m_mesh,
polyLine,
polyLineBoundingBoxes,
vistedEdges);

// no valid seed found in the entire mesh, we are done
Expand All @@ -267,6 +299,7 @@ void Mesh2DIntersections::Compute(const std::vector<Point>& polyLine)
continue;
}
IntersectFaceEdges(polyLine,
polyLineBoundingBoxes,
cumulativeLength,
currentCrossingEdge,
currentFaceIndex,
Expand Down
Loading