[FEATURE] Proper diffraction edge detection in `TriangleMesh`

[rwydaegh]

Terms

• Checked the existing issues to see if my suggestion has not already been suggested;

Description

TriangleMesh.diffraction_edges_mask returns jnp.ones(...), marking every edge as a diffraction edge. The docstring says this is a placeholder. Without fixing this you can't identify which edges to diffract on, can't get wedge angles for UTD, and can't build path candidate graphs that include diffraction.

I'd add edge adjacency computation (ideally in Rust) and fix the placeholder:

@property
def edge_adjacency(self) -> tuple[Int[Array, "num_triangles 3"], Int[Array, "num_triangles 3"]]:
    """Adjacent triangle index and shared edge index per triangle edge.
    -1 for boundary edges. Edge i connects vertex[i] to vertex[(i+1) % 3].
    """

@property
def diffraction_edges_mask(self) -> Bool[Array, "num_triangles 3"]:
    """True for non-planar wedge edges: two adjacent triangles whose
    dihedral angle is not ~pi (not coplanar)."""

@property
def wedge_angles(self) -> Float[Array, "num_triangles 3"]:
    """Exterior wedge angle n*pi per edge. n=2 for flat, n=1.5 for 90-degree corner."""

@property
def diffraction_edge_vertices(self) -> Float[Array, "num_diffraction_edges 2 3"]:
    """Start and end vertices of each diffraction edge."""

The core operation is building an edge-to-face adjacency map: hash sorted vertex-index pairs per triangle, find pairs. Rust is the right place for this (differt-core/src/geometry/triangle_mesh.rs) since it's topology on large meshes (Munich has ~40K triangles), it's not JIT-traceable, and the Rust mesh infrastructure is already there. Cache the result as NumPy on the Python side.

For non-manifold edges (>2 adjacent triangles) I'd just mark them as non-diffraction and log a warning.

One nice property: when assume_quads is set, the shared diagonal of two triangles forming a quad is coplanar. The dihedral angle check handles that naturally.

For testing: a box scene should detect 12 exterior wedge edges at corners and 0 on coplanar diagonals. A flat plane should detect 0. Munich as a sanity check for building corners and rooftop edges.

This is a prerequisite for #429 (diffraction_coefficients), which needs wedge_angles (the n parameter for UTD). #431 (transition_matrices) also needs it for its diffraction branch. And fermat_path_on_linear_objects in differt.rt already knows how to trace paths on edges, it just doesn't know which edges to use yet.

[jeertmans]

Hi @rwydaegh, I agree this is a very important (currently missing) feature.

I'd add edge adjacency computation (ideally in Rust) and fix the placeholder:

At first, I wanted to keep diffraction edge computation in JAX, but I agree that it is not possible to differentiate the with respect to diffraction edge, since we don't know how many there are. The only downside of using Rust is if we want to compute diffraction edges from a custom scene, i.e., that was built from custom input data.

@property
def edge_adjacency(self) -> tuple[Int[Array, "num_triangles 3"], Int[Array, "num_triangles 3"]]:
    """Adjacent triangle index and shared edge index per triangle edge.
    -1 for boundary edges. Edge i connects vertex[i] to vertex[(i+1) % 3].
    """

@property
def diffraction_edges_mask(self) -> Bool[Array, "num_triangles 3"]:
    """True for non-planar wedge edges: two adjacent triangles whose
    dihedral angle is not ~pi (not coplanar)."""

@property
def wedge_angles(self) -> Float[Array, "num_triangles 3"]:
    """Exterior wedge angle n*pi per edge. n=2 for flat, n=1.5 for 90-degree corner."""

@property
def diffraction_edge_vertices(self) -> Float[Array, "num_diffraction_edges 2 3"]:
    """Start and end vertices of each diffraction edge."""

The core operation is building an edge-to-face adjacency map: hash sorted vertex-index pairs per triangle, find pairs. Rust is the right place for this (differt-core/src/geometry/triangle_mesh.rs) since it's topology on large meshes (Munich has ~40K triangles), it's not JIT-traceable, and the Rust mesh infrastructure is already there. Cache the result as NumPy on the Python side.

As mentioned earlier, we might want to verify how to deal with concatenated meshes. One assumption we could do is that concatenating two meshes does not create or delete any existing edge, but still allow the user to force the re-computation of the edges (if desired).

For non-manifold edges (>2 adjacent triangles) I'd just mark them as non-diffraction and log a warning.

This makes sense to do.

One nice property: when assume_quads is set, the shared diagonal of two triangles forming a quad is coplanar. The dihedral angle check handles that naturally.

You're right.

For testing: a box scene should detect 12 exterior wedge edges at corners and 0 on coplanar diagonals. A flat plane should detect 0. Munich as a sanity check for building corners and rooftop edges.

Yes, good idea, and we could very easily include an example 3D figures in the documentation showing the edges in various colors.

---

Another possibility, which might be slower than a pure Rust implementation, but would work in Python, with JAX, and support JIT, would be to proceed as follows:

• List all edges in the scenes;
• Count the number of triangles sharing a given edges (we can use jnp.unique with size=num_edges and a fill value to allow for JIT compilation):
  • if 1 triangle, mark this as a free-floating edge;
  • if 2 triangles, mark this as a candidate wedge for diffraction (which we will compute the interior/exterior wedge angle);
  • else, mark this a degenerate edge (possibly raising a warning);
• Compute the wedge interior/exterior angles.

We would thus have something like this:

def build_edges(self, *, update:bool=Fale) -> Self:
    """We use it to allow pre-computing triangle edges (this way, we can use the results when concatenating two meshes)."""

@property
def triangle_edges(self) -> Int[Array, "num_triangles 2"]:
    """Start and end indices of each triangle edges"""

@property
def triangle_edges_parents(self) -> tuple[Int[Array, " num_triangles"], Int[Array, "vnum_triangles"], Bool[Array, "num_triangles"]]:
    """Return first, second (optional) and whether it has more parents.
    This is very important information as we need it when computing diffraction coefficients (since we need to now the 0-face and n-face).
    """
    # note that the logic below might be actually computed inside `build_edges`
    # first is always != -1 and count == 1
    # when second == -1, count == 2
    # we may want to raise a warning here if (count > 2).any()
    return first, second, count

@property
def floating_edges_mask(self) -> Float[Array, "num_triangles"]:
    """Start and end vertices of each triangle edges"""
    _, _, count = self.wedges_triangle_indices
    return count == 1

@property
def wedges_mask(self) -> Bool[Array, "num_triangles"]:
    """True for any edge shared by exactly 2 triangles (excluding assume quads)."""
    _, _, count = self.wedges_triangle_indices
    return count == 2

@property
def wedge_angles(self) -> Float[Array, "num_triangles"]:
    """Exterior wedge angles n*pi per edge. n=2 for flat, n=1.5 for 90-degree corner."""
   # TODO: how do we determine the interior / exterior? It should really depend on the origin of the rays. It's actually a very challenging task (unless you have some add. knowledge on your mesh, like the fact that each sub mesh is a closed convex object).

While maybe not the fastest (on CPU at least), I believe my solution would be easier to implement. What do you think?