Skip to content

fix: correct FMM source propagation#250

Open
liamnwhite1 wants to merge 3 commits into
nmwsharp:masterfrom
liamnwhite1:fix/fmm-source-propagation
Open

fix: correct FMM source propagation#250
liamnwhite1 wants to merge 3 commits into
nmwsharp:masterfrom
liamnwhite1:fix/fmm-source-propagation

Conversation

@liamnwhite1

Copy link
Copy Markdown
Contributor

Summary

This PR fixes several Fast Marching Method distance initialization and propagation issues for vertex, edge, face, and curve sources.

Details

FMM source initialization

FMMDistance() now stores tentative distances for source-adjacent vertices when they are enqueued.

Previously, triangle updates could read infinity for vertices that had been queued from source geometry but had not yet been finalized. Recording the tentative distance at enqueue time ensures that adjacent triangle updates operate on the intended provisional source distances.

Face-source distance correction

Face-source initialization now uses the correct opposite edge lengths in the barycentric squared-distance cross terms.

This PR also clamps tiny negative roundoff error before calling sqrt(). This prevents invalid source-to-vertex distances and avoids numerical NaNs caused by near-zero negative values.

Source-branch shortcut prevention

The FMM priority queue now carries a source label with each candidate distance.

Triangle updates only combine distances from the same source branch unless the update crosses an explicit source support edge. This prevents unrelated source branches from being combined during propagation, which could otherwise create artificial shortcuts through the mesh.

Explicit source-base edges are recorded for:

  • edge sources,
  • face-source support edges, and
  • signed curve segments.

This preserves valid source-local initialization while preventing unrelated branches from collapsing into incorrect shorter paths.

This was the most significant error addressed by the PR. The images below show the issue and the corrected behavior:

  1. the underlying distance field,
  2. the extracted isocontour before this fix, and
  3. the extracted isocontour after this fix.

There is still some deviation along the diagonal in the post-fix result. That deviation appears to come from the marching triangles extraction step and the underlying distance field, rather than from the FMM propagation bug fixed here. I plan to address that separately in a follow-up PR.

distance-field pre-fix post-fix

Impact

  • Improves correctness of FMM geodesic distances for non-vertex sources.
  • Prevents triangle updates from combining unrelated source branches.
  • Avoids NaNs from small negative roundoff values in face-source initialization.
  • Preserves the existing public API.

Store tentative source-adjacent distances when queueing candidates so triangle updates do not read infinities before source vertices are finalized.
Use the opposite edge lengths in the barycentric squared-distance cross terms and clamp tiny negative roundoff before taking square roots.
Track source labels through FMM candidates so triangle updates do not combine distances from unrelated source branches. Allow updates across explicit source-base edges from edge/face sources and signed curve segments so valid source support edges still initialize adjacent triangles.
@nzfeng

nzfeng commented Jun 10, 2026

Copy link
Copy Markdown
Collaborator

Thanks for opening a PR!

I agree with the fixes around line 211 for the squared-distance-to-triangle-vertex calculations, and clamping before sqrt --- I just double-checked the former and agree that there were typos.

The new addCandidateDistance function looks good to me.

Just to clarify the main new things happening in this PR:

  • I'm not sure I'm fully understanding the implications of adding source labels --- what does it mean to "only combine distances from the source branch"? Before I merge, I'd like to better understand why we want to prevent "unrelated source branches from being combined during propagation, which could otherwise create artificial shortcuts through the mesh." This sounds a bit counter to the usual definition of geodesic distance, which aims to always to provide the shortest possible distance to any point in the set of source curves and points, regardless of their "identity". Of course the behavior of true distance may be different than the application you have in mind.
  • It seems like sourceLabels are only either 0 or 1, corresponding to a vertex source or an edge source, respectively; would it make more sense to make them a bool instead of an int? Or, if you envision potentially expanding this to more types of sources than two, we can make it an enum.
  • I assume the source in the above images is the curve comprising the boundary of the square mesh?
  • Lines 3-5 don't seem necessary (the includes for algorithm, cmath, and limits), perhaps I can remove them.

I just ran your code on a few examples, and it seems to work well on a few toy examples. Thanks again for contributing.

@liamnwhite1

Copy link
Copy Markdown
Contributor Author

Hi @nzfeng ! Thank you for the work you've put into this algorithm. It's been immensely useful for the work we're doing!

Regarding your comments:

  • I'm glad you like the addCandidateDistance() helper. I added it because the old initialization path pushed source-adjacent vertices onto the priority queue without also storing their tentative distance in distances[] and so later triangle updates could read distances[neighVert] before that queued vertex had been finalized. If the tentative value was never stored, the update saw infinity and computed bad candidates. Routing all candidate insertion through one helper fixes that by updating distances[v] at enqueue time and also centralizes the “does this candidate improve the current value?” logic.

  • On sourceLabels: The intent is to prevent a triangle update from combining distances that originated from unrelated source supports. In the motivating failure mode, two source vertices can both be marked as sources, and the triangle update across the edge between them can behave as if the whole edge were a source, producing an artificially short distance to the opposite vertex. Below is a diagram of the issue for a simple square mesh. Without the source labeling, the opposite vertex is assigned a distance of $1/\sqrt(2)$ when really it should be $1$.
    Screenshot from 2026-06-11 09-07-31

  • Yes, the source in the images was the curve comprising the boundary of the square mesh. For reference, here is the related code:

    ManifoldSurfaceMesh::VertexData generateContouringField(Island& island, const std::vector<ManifoldSurfaceMesh::Vertex>& boundary_vertices) {
        auto& mesh = island.mesh();
    
        std::vector<std::pair<ManifoldSurfaceMesh::Vertex, double>>initial_distances;
        initial_distances.reserve(boundary_vertices.size());
        for (const auto& vertex : boundary_vertices) {
            initial_distances.emplace_back(vertex, 0.0001);
        }
    
        return geometrycentral::surface::FMMDistance(mesh.geometry(), initial_distances, true);
    }
    
  • I added the extra includes for <algorithm><cmath>, and <limits> just because I'm a stickler for following the "Include What You Use" (IWYU) philosophy so as not to rely on transitive dependencies. It make the code more robust but since they are already available transitively here and are not central to the fix, they don't have to be kept if you don't want to.

Hopefully this helps explain the PR a little better but if you would like for us to set up a meeting for me to show you more visually, just let me know!

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

2 participants