Merge pull request #1369 from hannobraun/triangulate

Fix triangulation of sharp, concave faces

crates/fj-kernel/src/algorithms/approx/face.rs

@@ -34,6 +34,7 @@ impl Approx for &FaceSet {
  let min_distance = ValidationConfig::default().distinct_min_distance;
  let mut all_points: BTreeSet<ApproxPoint<2>> = BTreeSet::new();
 
+ // Run some validation code on the approximation.
  for approx in &approx {
  let approx: &FaceApprox = approx;
 

crates/fj-kernel/src/algorithms/triangulate/delaunay.rs

@@ -1,17 +1,48 @@
+use std::collections::BTreeMap;
+
 use fj_math::{Point, Scalar, Triangle, Winding};
 use spade::HasPosition;
 
-use crate::objects::Handedness;
+use crate::{algorithms::approx::cycle::CycleApprox, objects::Handedness};
 
 /// Create a Delaunay triangulation of all points
 pub fn triangulate(
- points: Vec<TriangulationPoint>,
+ cycles: impl IntoIterator<Item = CycleApprox>,
  coord_handedness: Handedness,
 ) -> Vec<[TriangulationPoint; 3]> {
  use spade::Triangulation as _;
 
- let triangulation = spade::DelaunayTriangulation::<_>::bulk_load(points)
- .expect("Inserted invalid values into triangulation");
+ let mut triangulation = spade::ConstrainedDelaunayTriangulation::<_>::new();
+
+ let mut points = BTreeMap::new();
+
+ for cycle_approx in cycles {
+ let mut handle_prev = None;
+
+ for point in cycle_approx.points() {
+ let handle = match points.get(&point) {
+ Some(handle) => *handle,
+ None => {
+ let handle = triangulation
+ .insert(TriangulationPoint {
+ point_surface: point.local_form,
+ point_global: point.global_form,
+ })
+ .expect("Inserted invalid point into triangulation");
+
+ points.insert(point, handle);
+
+ handle
+ }
+ };
+
+ if let Some(handle_prev) = handle_prev {
+ triangulation.add_constraint(handle_prev, handle);
+ }
+
+ handle_prev = Some(handle);
+ }
+ }
 
  let mut triangles = Vec::new();
  for triangle in triangulation.inner_faces() {

crates/fj-kernel/src/algorithms/triangulate/mod.rs

@@ -6,7 +6,7 @@ mod polygon;
 use fj_interop::mesh::Mesh;
 use fj_math::Point;
 
-use self::{delaunay::TriangulationPoint, polygon::Polygon};
+use self::polygon::Polygon;
 
 use super::approx::{face::FaceApprox, Approx, Tolerance};
 
@@ -44,27 +44,20 @@ where
 
 impl Triangulate for FaceApprox {
  fn triangulate_into_mesh(self, mesh: &mut Mesh<Point<3>>) {
- let points: Vec<_> = self
- .points()
- .into_iter()
- .map(|point| TriangulationPoint {
- point_surface: point.local_form,
- point_global: point.global_form,
- })
- .collect();
  let face_as_polygon = Polygon::new()
  .with_exterior(
  self.exterior
  .points()
  .into_iter()
  .map(|point| point.local_form),
  )
- .with_interiors(self.interiors.into_iter().map(|interior| {
+ .with_interiors(self.interiors.iter().map(|interior| {
  interior.points().into_iter().map(|point| point.local_form)
  }));
 
+ let cycles = [self.exterior].into_iter().chain(self.interiors);
  let mut triangles =
- delaunay::triangulate(points, self.coord_handedness);
+ delaunay::triangulate(cycles, self.coord_handedness);
  triangles.retain(|triangle| {
  face_as_polygon
  .contains_triangle(triangle.map(|point| point.point_surface))
@@ -176,29 +169,26 @@ mod tests {
  Ok(())
  }
 
- #[ignore]
  #[test]
  fn sharp_concave_shape() -> anyhow::Result<()> {
  let objects = Objects::new();
 
- //
- // c
- // /|
- // e / |
- // |\ / |
- // | | / |
- // | \ / |
- // | \ / |
- // | d |
- // a ---------- b
- //
+ // e c
+ // |\ /|
+ // \ \ / b
+ // \ \ / /
+ // \ d /
+ // \a/
 
- let a = [0., 0.];
- let b = [0.4, 0.];
- //let b = [0.5, 0.]; // test passes with this change
- let c = [0.4, 1.0];
+ // Naive Delaunay triangulation will create a triangle (c, d, e), which
+ // is not part of the polygon. The higher-level triangulation will
+ // filter that out, but it will result in missing triangles.
+
+ let a = [0.1, 0.0];
+ let b = [0.2, 0.9];
+ let c = [0.2, 1.0];
  let d = [0.1, 0.1];
- let e = [0., 0.8];
+ let e = [0.0, 1.0];
 
  let surface = objects.surfaces.xy_plane();
  let face = Face::partial()
@@ -209,17 +199,15 @@ mod tests {
 
  let triangles = triangulate(face)?;
 
- let a3 = surface.geometry().point_from_surface_coords(a);
- let b3 = surface.geometry().point_from_surface_coords(b);
- let c3 = surface.geometry().point_from_surface_coords(c);
- let d3 = surface.geometry().point_from_surface_coords(d);
- let e3 = surface.geometry().point_from_surface_coords(e);
-
- assert!(triangles.contains_triangle([a3, b3, d3]));
- assert!(triangles.contains_triangle([b3, c3, d3]));
- assert!(triangles.contains_triangle([a3, d3, e3]));
+ let a = surface.geometry().point_from_surface_coords(a);
+ let b = surface.geometry().point_from_surface_coords(b);
+ let c = surface.geometry().point_from_surface_coords(c);
+ let d = surface.geometry().point_from_surface_coords(d);
+ let e = surface.geometry().point_from_surface_coords(e);
 
- assert!(!triangles.contains_triangle([b3, e3, d3]));
+ assert!(triangles.contains_triangle([a, b, d]));
+ assert!(triangles.contains_triangle([a, d, e]));
+ assert!(triangles.contains_triangle([b, c, d]));
 
  Ok(())
  }