From f45e665ccf9e844eafee134a820e27b25079b189 Mon Sep 17 00:00:00 2001 From: Mattias Fagerlund Date: Fri, 30 Oct 2026 08:04:19 +0100 Subject: [PATCH 1/2] Exact geodesics: require vertex Gaussian curvatures once per propagate The saddle test in GeodesicAlgorithmExact::propagate required and unrequired the curvatures for every popped interval. Require them once before the loop and release them after it instead. Results are unchanged; add a test that the distance field does not depend on whether the caller holds the curvatures. --- src/surface/exact_geodesics.cpp | 4 ++-- test/CMakeLists.txt | 1 + test/src/exact_geodesics_test.cpp | 38 +++++++++++++++++++++++++++++++ 3 files changed, 41 insertions(+), 2 deletions(-) create mode 100644 test/src/exact_geodesics_test.cpp diff --git a/src/surface/exact_geodesics.cpp b/src/surface/exact_geodesics.cpp index 6b126962..cd76d1dd 100644 --- a/src/surface/exact_geodesics.cpp +++ b/src/surface/exact_geodesics.cpp @@ -444,6 +444,7 @@ void GeodesicAlgorithmExact::propagate(const std::vector& sources, IntervalWithStop candidates[2]; geom.requireCornerAngles(); + geom.requireVertexGaussianCurvatures(); while (!m_queue.empty()) { m_queue_max_size = std::max(m_queue.size(), m_queue_max_size); @@ -469,9 +470,7 @@ void GeodesicAlgorithmExact::propagate(const std::vector& sources, bool const last_interval = min_interval->next() == nullptr; auto saddleOrBoundary = [&](Vertex v) -> bool { - geom.requireVertexGaussianCurvatures(); bool saddle = geom.vertexGaussianCurvatures[v] < 0; - geom.unrequireVertexGaussianCurvatures(); return saddle || vertexIsBoundary[v]; }; @@ -607,6 +606,7 @@ void GeodesicAlgorithmExact::propagate(const std::vector& sources, m_time_consumed = (static_cast(stop) - static_cast(start)) / CLOCKS_PER_SEC; geom.unrequireCornerAngles(); + geom.unrequireVertexGaussianCurvatures(); /* for(unsigned i=0; i + +using namespace geometrycentral; +using namespace geometrycentral::surface; + +class ExactGeodesicSuite : public MeshAssetSuite {}; + +// The distance field does not depend on whether the caller holds the curvatures used by the saddle test +TEST_F(ExactGeodesicSuite, ResultIndependentOfRequiredQuantities) { + for (const MeshAsset& a : {getAsset("spot.ply", true), getAsset("bob_small.ply", true)}) { + a.printThyName(); + SurfaceMesh& mesh = *a.mesh; + VertexPositionGeometry& geometry = *a.geometry; + Vertex source = mesh.vertex(0); + + GeodesicAlgorithmExact mmp(mesh, geometry); + mmp.propagate(source); + VertexData dist = mmp.getDistanceFunction(); + + geometry.requireVertexGaussianCurvatures(); + GeodesicAlgorithmExact mmpHeld(mesh, geometry); + mmpHeld.propagate(source); + VertexData distHeld = mmpHeld.getDistanceFunction(); + geometry.unrequireVertexGaussianCurvatures(); + + for (Vertex v : mesh.vertices()) { + EXPECT_EQ(dist[v], distHeld[v]); + } + } +} From 4da494ec56c59a5d69bca341da7b56eed40a53df Mon Sep 17 00:00:00 2001 From: Mattias Fagerlund Date: Fri, 30 Oct 2026 08:04:21 +0100 Subject: [PATCH 2/2] Exact geodesics: honour max_propagation_distance without stop points check_stop_conditions returned false whenever there were no stop vertices, so propagate(source, max_propagation_distance) ignored the bound and covered the whole mesh. Stop once the queue passes the bound when no stop points are given; the behaviour with stop points is unchanged. Add tests with and without stop points. --- src/surface/exact_geodesics.cpp | 3 +- test/src/exact_geodesics_test.cpp | 77 +++++++++++++++++++++++++++++++ 2 files changed, 79 insertions(+), 1 deletion(-) diff --git a/src/surface/exact_geodesics.cpp b/src/surface/exact_geodesics.cpp index cd76d1dd..4ffd596d 100644 --- a/src/surface/exact_geodesics.cpp +++ b/src/surface/exact_geodesics.cpp @@ -631,7 +631,8 @@ bool GeodesicAlgorithmExact::check_stop_conditions(unsigned& index) const { } if (m_stop_vertices.empty()) { - return false; + // without stop vertices, the distance bound alone stops the propagation + return stop_dist < GEODESIC_INF; } while (index < m_stop_vertices.size()) { diff --git a/test/src/exact_geodesics_test.cpp b/test/src/exact_geodesics_test.cpp index 4772e56e..73d6f7bb 100644 --- a/test/src/exact_geodesics_test.cpp +++ b/test/src/exact_geodesics_test.cpp @@ -6,6 +6,7 @@ #include "gtest/gtest.h" +#include #include using namespace geometrycentral; @@ -13,6 +14,18 @@ using namespace geometrycentral::surface; class ExactGeodesicSuite : public MeshAssetSuite {}; +namespace { + +// A distance that about half of the vertices are closer than +double medianDistance(SurfaceMesh& mesh, const VertexData& dist) { + std::vector d; + for (Vertex v : mesh.vertices()) d.push_back(dist[v]); + std::nth_element(d.begin(), d.begin() + d.size() / 2, d.end()); + return d[d.size() / 2]; +} + +} // namespace + // The distance field does not depend on whether the caller holds the curvatures used by the saddle test TEST_F(ExactGeodesicSuite, ResultIndependentOfRequiredQuantities) { for (const MeshAsset& a : {getAsset("spot.ply", true), getAsset("bob_small.ply", true)}) { @@ -36,3 +49,67 @@ TEST_F(ExactGeodesicSuite, ResultIndependentOfRequiredQuantities) { } } } + +// With a distance bound and no stop points, propagation stops at the bound: vertices closer than the bound +// get their exact distance, and some vertices beyond it are left unreached +TEST_F(ExactGeodesicSuite, MaxPropagationDistanceWithoutStopPoints) { + for (const MeshAsset& a : {getAsset("spot.ply", true), getAsset("bob_small.ply", true)}) { + a.printThyName(); + SurfaceMesh& mesh = *a.mesh; + VertexPositionGeometry& geometry = *a.geometry; + Vertex source = mesh.vertex(0); + + GeodesicAlgorithmExact full(mesh, geometry); + full.propagate(source); + VertexData ref = full.getDistanceFunction(); + double bound = medianDistance(mesh, ref); + + GeodesicAlgorithmExact mmp(mesh, geometry); + mmp.propagate(source, bound); + VertexData dist = mmp.getDistanceFunction(); + + size_t unreached = 0; + for (Vertex v : mesh.vertices()) { + if (ref[v] < bound) { + EXPECT_NEAR(dist[v], ref[v], 1e-12 * bound); + } + if (dist[v] >= GEODESIC_INF) unreached++; + } + EXPECT_GT(unreached, 0u); + } +} + +// With stop points, propagation stops once they are covered; with a distance bound too, once both hold +TEST_F(ExactGeodesicSuite, MaxPropagationDistanceWithStopPoints) { + for (const MeshAsset& a : {getAsset("spot.ply", true), getAsset("bob_small.ply", true)}) { + a.printThyName(); + SurfaceMesh& mesh = *a.mesh; + VertexPositionGeometry& geometry = *a.geometry; + Vertex source = mesh.vertex(0); + + GeodesicAlgorithmExact full(mesh, geometry); + full.propagate(source); + VertexData ref = full.getDistanceFunction(); + double bound = medianDistance(mesh, ref); + + // A stop vertex well inside the bound + Vertex stop; + for (Vertex v : mesh.vertices()) { + if (v != source && ref[v] < 0.5 * bound && (stop == Vertex() || ref[v] > ref[stop])) stop = v; + } + ASSERT_TRUE(stop != Vertex()); + + GeodesicAlgorithmExact mmpStop(mesh, geometry); + mmpStop.propagate(source, GEODESIC_INF, {stop}); + EXPECT_NEAR(mmpStop.getDistance(stop), ref[stop], 1e-12 * bound); + + GeodesicAlgorithmExact mmpBoth(mesh, geometry); + mmpBoth.propagate(source, bound, {stop}); + VertexData dist = mmpBoth.getDistanceFunction(); + for (Vertex v : mesh.vertices()) { + if (ref[v] < bound) { + EXPECT_NEAR(dist[v], ref[v], 1e-12 * bound); + } + } + } +}