diff --git a/cdt/src/contour.rs b/cdt/src/contour.rs index c970c6a..66c0155 100644 --- a/cdt/src/contour.rs +++ b/cdt/src/contour.rs @@ -161,7 +161,7 @@ impl Contour { ContourData::None => (), ContourData::Hull(hull_index, sign) => { t.hull.update(hull_index, e_ca); - t.half.set_sign(e_bc, sign); + t.half.set_sign(e_ca, sign); }, ContourData::Buddy(b) => t.half.link_new(b, e_ca), }; diff --git a/cdt/tests/concave_erase.rs b/cdt/tests/concave_erase.rs new file mode 100644 index 0000000..a2f559a --- /dev/null +++ b/cdt/tests/concave_erase.rs @@ -0,0 +1,54 @@ +// A deeply concave (but simple -- no self-intersections) heptagon whose constrained +// triangulation must be the 5-triangle interior partition. On the unfixed code the +// flood erase removes nothing and the full 6-triangle convex-hull triangulation is +// returned instead: constraint (3,4) sits on the partial hull when constraint (2,3)'s +// cavity walk consumes it, and `Contour::try_clip`'s positive-branch `a` arm restores +// the preserved sign onto the wrong rebuilt edge (`e_bc` instead of `e_ca`), so the +// constraint comes back unfixed and the odd-even erase parity leaks through it. +#[test] +fn concave_heptagon_erases_its_pocket() { + let pts: Vec<(f64, f64)> = vec![ + (41.3, 9.1), + (22.1, 63.1), + (14.2, 64.4), + (33.6, 8.0), + (1.6, 12.8), + (0.2, 5.9), + (40.3, 0.1), + ]; + let edges: Vec<(usize, usize)> = (0..7).map(|i| (i, (i + 1) % 7)).collect(); + let tris = cdt::triangulate_with_edges(&pts, &edges).expect("triangulation"); + + // n-2 triangles for a simple n-gon's interior partition. + assert_eq!( + tris.len(), + 5, + "expected the interior partition, got {:?}", + tris + ); + + // Unsigned triangle area must equal the polygon's net shoelace area (a partition -- + // the convex-hull triangulation covers ~2.2x this polygon's net area instead). + let tri_area: f64 = tris + .iter() + .map(|&(a, b, c)| { + let (ax, ay) = pts[a]; + let (bx, by) = pts[b]; + let (cx, cy) = pts[c]; + ((bx - ax) * (cy - ay) - (cx - ax) * (by - ay)).abs() * 0.5 + }) + .sum(); + let mut net = 0.0f64; + for i in 0..pts.len() { + let (x0, y0) = pts[i]; + let (x1, y1) = pts[(i + 1) % pts.len()]; + net += x0 * y1 - x1 * y0; + } + let net = (net * 0.5).abs(); + assert!( + (tri_area - net).abs() <= net * 1e-9, + "covered area {} != polygon net area {}", + tri_area, + net + ); +}