graphic design tool
git clone https://git.lucas.co/cce-designer.git
feat: detangle's Surface method, and a measure of what crosses
Method: Surface tests each point against the triangles near it and shares
the move between the point and the triangle's corners; Points, the first
version's solve, stays the default and what a node without the row runs.
self_intersections counts every edge through a triangle, and the node's
Tangled Group writes the points of what is still crossed after the solve.
Co-Authored-By: Claude Fable 5.1 <noreply@anthropic.com>
CLAUDE.md | 40 ++++++
nodes/detangle.json | 4 +-
src/detangle.rs | 399 ++++++++++++++++++++++++++++++++++++++++++++++++++--
src/main.rs | 206 +++++++++++++++++++++++++++
src/spatial.rs | 76 ++++++++++
5 files changed, 715 insertions(+), 10 deletions(-)
diff --git a/CLAUDE.md b/CLAUDE.md
index 0cf4f56..f63326e 100644
--- a/CLAUDE.md
+++ b/CLAUDE.md
@@ -1157,6 +1157,46 @@ to frame 240 from 113 to 71: once the whole surface is within a
thickness of itself every pass runs and the pairs themselves are the
work, and no bookkeeping saves that.
+**The Surface method and the measure (2026-09-29).** Everything above is
+the node's `Points` method, which is what a node without a `Method` row
+runs and what the template defaults to, so a save from before solves as it
+did. `Method: Surface` (`detangle::solve_surface`) tests each point against
+the TRIANGLES near it: a point over the middle of a triangle is near no
+corner of it, so where triangles are larger than the thickness the point
+test sees nothing at all. A contact is resolved along the line from the
+closest point on the triangle to the point (the triangle's normal where the
+point lies on it), and the move is SHARED — the point one way, the corners
+the other by how much of the closest point each is
+(`spatial::closest_weights_on_triangle`), with what is outside the Group
+taking none and the rest all of it, so a contact with a fixed triangle is
+resolved whole where Points resolves half. What a point receives from
+several contacts is their average weighted by depth, not their sum: a point
+over a shared edge touches both triangles and must move once. A triangle
+with a corner inside the point's excluded rings is not a contact. It does
+NOT know which side a point belongs on — one already through is pushed
+further through; that needs the positions the step began from and is the
+next piece of work, as are a cap on movement per substep and edge-edge
+contact.
+
+`detangle::self_intersections` is the MEASURE: every edge passing through a
+triangle (`spatial::segment_crosses_triangle`, tolerance relative to the
+lengths, so scale does not change the answer), no thickness and no rings.
+`crossings_beyond(geom, rings)` counts only those the solve is meant to see
+at a ring count, which is what separates a miss of the method from a fold
+inside the excluded neighbourhood. The node's **Tangled Group** row, when
+it names one, writes the points of what is STILL crossed after the solve
+(empty when nothing is); it costs a second search of the mesh and is off
+by default. `detangle_methods_compared` (ignored; release, `--ignored
+--nocapture`) pushes an icosphere's cap down into its own bowl a fifth of
+an edge a step. At 2562 points, Thickness 1, Rings 2: no detangle 1117
+crossings, all beyond the rings; Points 1699 (892 beyond), 5.1 ms a step;
+Surface 408, NONE beyond the rings, and none left at the end, 17 ms a step.
+At Thickness 0.5 Surface let 72 through beyond the rings at that size and
+none at 162 and 642 points. What Surface leaves is the fold at the cap's
+rim, inside the rings. Do not measure by pressing a sphere flat by the sign
+of y: that carries the equator's points past their own neighbours, which no
+setting is meant to see, and both methods look equally bad.
+
What still costs is the solver's, not the node's: an edit inside a simnet
re-solves from the seed, so a change at frame 120 is 120 steps. A
backward scrub no longer does — the next section.
diff --git a/nodes/detangle.json b/nodes/detangle.json
index af323fe..57ea4ae 100644
--- a/nodes/detangle.json
+++ b/nodes/detangle.json
@@ -5,9 +5,11 @@
"outputs": 1,
"params": [
{ "name": "Input", "type": "node", "default": "" },
+ { "name": "Method", "type": "choice:Points,Surface", "default": "Points" },
{ "name": "Thickness", "type": "slider", "default": "1.00", "min": 0.0, "max": 8.0, "step": 0.05 },
{ "name": "Rings", "type": "spinbox", "default": "2", "min": 0.0, "max": 6.0, "step": 1.0 },
{ "name": "Iterations", "type": "spinbox", "default": "4", "min": 1.0, "max": 32.0, "step": 1.0 },
- { "name": "Group", "type": "group", "default": "" }
+ { "name": "Group", "type": "group", "default": "" },
+ { "name": "Tangled Group", "type": "group", "default": "" }
]
}
diff --git a/src/detangle.rs b/src/detangle.rs
index 84fe2d3..c4dbec9 100644
--- a/src/detangle.rs
+++ b/src/detangle.rs
@@ -24,6 +24,18 @@
//! The results are the first version's BIT FOR BIT, which is what lets these
//! be optimizations and not changes: the same pairs, summed in the same
//! order. `the_detangle_solve_matches_its_reference` holds them together.
+//!
+//! That is the node's `Points` method. Two things here are NOT the first
+//! version's, and are asked for by name:
+//!
+//! - **The `Surface` method** ([`solve_surface`]) tests each point against
+//! the TRIANGLES near it, where Points tests it against points. A point
+//! over the middle of a triangle is near no corner of it, so on a mesh
+//! whose triangles are larger than the thickness the point test sees
+//! nothing at all.
+//! - **The measure** ([`self_intersections`]): every edge that passes
+//! through a triangle. It is what says whether a change to the solve
+//! helped, and what the node's `Tangled Group` is written from.
use crate::app::FsNode;
use crate::detail::Detail;
@@ -37,6 +49,10 @@ struct Topo {
key: u64,
rings: usize,
edges: Vec<[u32; 2]>,
+ /// The primitives as triangles, fanned as `Detail::triangulate` fans
+ /// them, and the primitive each came from.
+ tris: Vec<[u32; 3]>,
+ tri_prims: Vec<u32>,
/// Each point's excluded neighbourhood, itself included, ascending:
/// point `p`'s is `excluded[starts[p]..starts[p + 1]]`.
starts: Vec<u32>,
@@ -80,7 +96,16 @@ impl Topo {
excluded.extend_from_slice(&seen);
}
starts.push(excluded.len() as u32);
- Topo { key, rings, edges: geom.edges().to_vec(), starts, excluded }
+ let mut tris = Vec::new();
+ let mut tri_prims = Vec::new();
+ for prim in 0..geom.num_prims() {
+ let pts = geom.prim_points(prim);
+ for i in 1..pts.len().saturating_sub(1) {
+ tris.push([pts[0], pts[i], pts[i + 1]]);
+ tri_prims.push(prim as u32);
+ }
+ }
+ Topo { key, rings, edges: geom.edges().to_vec(), tris, tri_prims, starts, excluded }
}
fn excludes(&self, p: usize, q: u32) -> bool {
@@ -244,6 +269,11 @@ pub struct Work {
pub passes: usize,
pub grids: usize,
pub searched: usize,
+ /// Point-triangle contacts the Surface method resolved, over every pass.
+ pub contacts: usize,
+ /// Edges passing through a triangle when the solve was done — counted
+ /// only when the node names a Tangled Group to write them to.
+ pub crossings: usize,
}
pub fn apply(geom: &mut Detail, target: &FsNode) -> Work {
@@ -257,17 +287,48 @@ pub fn apply(geom: &mut Detail, target: &FsNode) -> Work {
if topo.edges.is_empty() {
return work;
}
- // Thickness in EDGE LENGTHS, so the setting means the same thing before
- // and after a remesh.
- let mean_edge = topo
- .edges
+ // A node without the row is one from before it, and solves as it did.
+ if node_param_str(target, "Method", "Points").trim().eq_ignore_ascii_case("Surface") {
+ solve_surface(geom, target, &topo, &mut work);
+ } else {
+ solve_points(geom, target, &topo, &mut work);
+ }
+ // What is STILL crossed once the solve is done, which is the part worth
+ // looking at. Only when asked for: it is a second search of the mesh.
+ let mark = node_param_str(target, "Tangled Group", "");
+ let mark = mark.trim();
+ if !mark.is_empty() {
+ let found = intersections_of(geom, &topo, false);
+ work.crossings = found.crossings;
+ geom.points_mut().create_group(mark);
+ for p in found.points {
+ geom.points_mut().add_to_group(mark, p as usize);
+ }
+ }
+ work
+}
+
+/// The node's thickness as a length: the setting is in EDGE LENGTHS, so it
+/// means the same thing before and after a remesh.
+fn thickness_of(geom: &Detail, target: &FsNode, topo: &Topo) -> f32 {
+ let mean_edge = mean_edge(geom, topo);
+ node_param_f32(target, "Thickness", 1.0).max(0.0) * mean_edge
+}
+
+fn mean_edge(geom: &Detail, topo: &Topo) -> f32 {
+ topo.edges
.iter()
.map(|e| (geom.pos(e[1] as usize) - geom.pos(e[0] as usize)).length())
.sum::<f32>()
- / topo.edges.len() as f32;
- let thickness = node_param_f32(target, "Thickness", 1.0).max(0.0) * mean_edge;
+ / topo.edges.len() as f32
+}
+
+/// The Points method: the first version's solve.
+fn solve_points(geom: &mut Detail, target: &FsNode, topo: &Topo, work: &mut Work) {
+ let n = geom.num_points();
+ let thickness = thickness_of(geom, target, topo);
if thickness <= 0.0 {
- return work;
+ return;
}
let iterations = node_param_f32(target, "Iterations", 4.0).clamp(1.0, 32.0) as usize;
let group = node_param_str(target, "Group", "");
@@ -340,5 +401,325 @@ pub fn apply(geom: &mut Detail, target: &FsNode) -> Work {
geom.set_pos(p, *v);
}
}
- work
+}
+
+/// Triangles filed by cell, flat, as [`FlatGrid`] files points: a triangle
+/// is in every cell its bounding box touched when it was filed.
+struct TriCells {
+ min: Vec3,
+ cell: f32,
+ dims: [i32; 3],
+ starts: Vec<u32>,
+ ids: Vec<u32>,
+ /// Where each POINT was when the triangles were filed.
+ filed: Vec<Vec3>,
+}
+
+impl TriCells {
+ fn build(points: &[Vec3], tris: &[[u32; 3]], cell: f32) -> TriCells {
+ let (min, max) = points.iter().fold(
+ (Vec3::splat(f32::MAX), Vec3::splat(f32::MIN)),
+ |(lo, hi), &p| (lo.min(p), hi.max(p)),
+ );
+ let (min, max) = if points.is_empty() { (Vec3::ZERO, Vec3::ZERO) } else { (min, max) };
+ let span = (max - min).max(Vec3::splat(1e-6));
+ let cell = cell.max(span.max_element() / 128.0).max(1e-6);
+ let dims = [
+ ((span.x / cell).ceil() as i32 + 1).clamp(1, 256),
+ ((span.y / cell).ceil() as i32 + 1).clamp(1, 256),
+ ((span.z / cell).ceil() as i32 + 1).clamp(1, 256),
+ ];
+ let mut grid = TriCells { min, cell, dims, starts: Vec::new(), ids: Vec::new(), filed: points.to_vec() };
+ let cells = (dims[0] * dims[1] * dims[2]) as usize;
+ let boxes: Vec<([i32; 3], [i32; 3])> = tris
+ .iter()
+ .map(|t| {
+ let [a, b, c] = t.map(|i| points[i as usize]);
+ (grid.coord(a.min(b).min(c)), grid.coord(a.max(b).max(c)))
+ })
+ .collect();
+ // The counting sort again, a triangle counted once per cell it is
+ // in. Placing in triangle order leaves each cell's ids ascending.
+ let mut starts = vec![0u32; cells + 1];
+ for (a, b) in &boxes {
+ for z in a[2]..=b[2] {
+ for y in a[1]..=b[1] {
+ for x in a[0]..=b[0] {
+ starts[grid.index([x, y, z]) + 1] += 1;
+ }
+ }
+ }
+ }
+ for c in 0..cells {
+ starts[c + 1] += starts[c];
+ }
+ let mut next = starts.clone();
+ let mut ids = vec![0u32; starts[cells] as usize];
+ for (i, (a, b)) in boxes.iter().enumerate() {
+ for z in a[2]..=b[2] {
+ for y in a[1]..=b[1] {
+ for x in a[0]..=b[0] {
+ let c = grid.index([x, y, z]);
+ ids[next[c] as usize] = i as u32;
+ next[c] += 1;
+ }
+ }
+ }
+ }
+ grid.starts = starts;
+ grid.ids = ids;
+ grid
+ }
+
+ fn coord(&self, p: Vec3) -> [i32; 3] {
+ let rel = (p - self.min) / self.cell;
+ [
+ (rel.x.floor() as i32).clamp(0, self.dims[0] - 1),
+ (rel.y.floor() as i32).clamp(0, self.dims[1] - 1),
+ (rel.z.floor() as i32).clamp(0, self.dims[2] - 1),
+ ]
+ }
+
+ fn index(&self, c: [i32; 3]) -> usize {
+ ((c[2] * self.dims[1] + c[1]) * self.dims[0] + c[0]) as usize
+ }
+
+ /// Every triangle filed in a cell the box touches, once, in the order
+ /// the cells are walked. `seen` is one mark per triangle and `stamp`
+ /// this gather's, which is cheaper than sorting what came back to find
+ /// the triangles that came back twice.
+ fn gather(&self, lo: Vec3, hi: Vec3, seen: &mut [u32], stamp: u32, out: &mut Vec<u32>) {
+ out.clear();
+ let (a, b) = (self.coord(lo), self.coord(hi));
+ for z in a[2]..=b[2] {
+ for y in a[1]..=b[1] {
+ for x in a[0]..=b[0] {
+ let c = self.index([x, y, z]);
+ for &t in &self.ids[self.starts[c] as usize..self.starts[c + 1] as usize] {
+ if seen[t as usize] != stamp {
+ seen[t as usize] = stamp;
+ out.push(t);
+ }
+ }
+ }
+ }
+ }
+ }
+
+ fn drift(&self, points: &[Vec3]) -> f32 {
+ points.iter().zip(&self.filed).map(|(p, f)| (*p - *f).length_squared()).fold(0.0f32, f32::max).sqrt()
+ }
+}
+
+/// The Surface method: each point against the triangles near it.
+///
+/// A contact is a point closer than the thickness to a triangle none of
+/// whose corners is in the point's excluded rings. It is resolved along the
+/// line from the closest point on the triangle to the point — the triangle's
+/// own normal where the point lies ON it, which is a direction, where two
+/// coincident points have none — and the move is SHARED: the point takes
+/// its part one way and the triangle's corners theirs the other, each corner
+/// by how much of the closest point it is. What cannot move (outside the
+/// Group) takes none and the rest take all of it, so a contact with a fixed
+/// triangle is resolved whole, where Points resolves half of it.
+///
+/// A pass is gathered against the positions at its start and applied at its
+/// end, as Points is. What a point receives from several contacts is their
+/// AVERAGE, weighted by how deep each is: a point over a shared edge is in
+/// contact with both triangles and must move once, not twice, and a sum is
+/// what makes a dense contact overshoot and ring.
+///
+/// What it does not know is which SIDE a point belongs on. A point already
+/// through a triangle is pushed further through; that takes the positions
+/// the step started from, and is not here.
+fn solve_surface(geom: &mut Detail, target: &FsNode, topo: &Topo, work: &mut Work) {
+ let n = geom.num_points();
+ let thickness = thickness_of(geom, target, topo);
+ if thickness <= 0.0 || topo.tris.is_empty() {
+ return;
+ }
+ let iterations = node_param_f32(target, "Iterations", 4.0).clamp(1.0, 32.0) as usize;
+ let group = node_param_str(target, "Group", "");
+ let group = group.trim().to_string();
+ let movable: Vec<bool> = (0..n).map(|p| group.is_empty() || geom.points().in_group(&group, p)).collect();
+ let free = |p: usize| if movable[p] { 1.0f32 } else { 0.0 };
+
+ let mut pos: Vec<Vec3> = (0..n).map(|p| geom.pos(p)).collect();
+ // Cells no smaller than a triangle, or each is filed in dozens.
+ let cell = thickness.max(mean_edge(geom, topo));
+ let mut grid = TriCells::build(&pos, &topo.tris, cell);
+ work.grids += 1;
+ let mut drift = 0.0f32;
+ let mut near = Vec::new();
+ let mut seen = vec![u32::MAX; topo.tris.len()];
+ let mut stamp = 0u32;
+ let mut push = vec![Vec3::ZERO; n];
+ let mut weight = vec![0.0f32; n];
+ let mut moved_at_all = false;
+ for _ in 0..iterations {
+ work.passes += 1;
+ if drift > DRIFT_CELLS * grid.cell {
+ grid = TriCells::build(&pos, &topo.tris, cell);
+ work.grids += 1;
+ drift = 0.0;
+ }
+ push.iter_mut().for_each(|v| *v = Vec3::ZERO);
+ weight.iter_mut().for_each(|w| *w = 0.0);
+ let mut any = false;
+ for p in 0..n {
+ work.searched += 1;
+ // A triangle within a thickness of here now had its box within
+ // a thickness and a drift of here when it was filed.
+ let reach = Vec3::splat(thickness + drift);
+ stamp = stamp.wrapping_add(1);
+ if stamp == u32::MAX {
+ seen.iter_mut().for_each(|m| *m = u32::MAX);
+ stamp = 0;
+ }
+ grid.gather(pos[p] - reach, pos[p] + reach, &mut seen, stamp, &mut near);
+ for &t in &near {
+ let corners = topo.tris[t as usize];
+ let [ia, ib, ic] = corners.map(|c| c as usize);
+ let (a, b, c) = (pos[ia], pos[ib], pos[ic]);
+ // Most of what a cell holds is nowhere near: a triangle
+ // whose box is a thickness away on any axis is further
+ // than that, and is turned away before it costs a search
+ // of the rings or a closest point.
+ let (lo, hi) = (a.min(b).min(c) - thickness, a.max(b).max(c) + thickness);
+ if pos[p].cmplt(lo).any() || pos[p].cmpgt(hi).any() {
+ continue;
+ }
+ if corners.iter().any(|&c| topo.excludes(p, c)) {
+ continue;
+ }
+ let w = crate::spatial::closest_weights_on_triangle(pos[p], a, b, c);
+ let d = pos[p] - (a * w[0] + b * w[1] + c * w[2]);
+ let len = d.length();
+ if len >= thickness {
+ continue;
+ }
+ // Inverse masses of one or none: the point's, and each
+ // corner's by the square of its share.
+ let give = free(p) + free(ia) * w[0] * w[0] + free(ib) * w[1] * w[1] + free(ic) * w[2] * w[2];
+ if give <= 0.0 {
+ continue;
+ }
+ let dir = if len < 1e-9 {
+ let normal = (b - a).cross(c - a).normalize_or_zero();
+ if normal == Vec3::ZERO {
+ continue;
+ }
+ normal
+ } else {
+ d / len
+ };
+ let deep = thickness - len;
+ let step = dir * (deep / give);
+ push[p] += step * (free(p) * deep);
+ weight[p] += free(p) * deep;
+ for (i, share) in [(ia, w[0]), (ib, w[1]), (ic, w[2])] {
+ // The corner moves by its share; its say in the average
+ // is its share too, so a corner that is barely part of
+ // one contact does not water down another it carries.
+ push[i] -= step * (share * free(i) * share * deep);
+ weight[i] += free(i) * share * deep;
+ }
+ work.contacts += 1;
+ any = true;
+ }
+ }
+ if !any {
+ break;
+ }
+ for p in 0..n {
+ if weight[p] > 0.0 {
+ pos[p] += push[p] / weight[p];
+ }
+ }
+ moved_at_all = true;
+ drift = grid.drift(&pos);
+ }
+ if moved_at_all {
+ for (p, v) in pos.iter().enumerate() {
+ geom.set_pos(p, *v);
+ }
+ }
+}
+
+/// Where a surface passes through itself.
+#[derive(Debug, Default, Clone, PartialEq)]
+pub struct Tangles {
+ /// Edge-triangle pairs where the edge passes through the triangle.
+ pub crossings: usize,
+ /// The points of those edges and triangles, ascending, each once.
+ pub points: Vec<u32>,
+}
+
+/// Every edge of `geom` that passes through one of its triangles.
+///
+/// This is the MEASURE, and it is geometric: no thickness and no rings. An
+/// edge is not tested against a triangle it shares a point with — they meet
+/// there by construction — nor against one fanned from a primitive the edge
+/// is a side of.
+pub fn self_intersections(geom: &Detail) -> Tangles {
+ if geom.num_points() == 0 || geom.num_prims() == 0 {
+ return Tangles::default();
+ }
+ let topo = topo_for(geom, 0);
+ if topo.edges.is_empty() {
+ return Tangles::default();
+ }
+ intersections_of(geom, &topo, false)
+}
+
+/// [`self_intersections`] counting only what the solve is MEANT to see at
+/// this ring count: an edge through a triangle none of whose corners is
+/// within `rings` of either end of it. The rest are folds inside the
+/// excluded neighbourhood, which the node leaves alone by design, and
+/// telling the two apart is what says whether a crossing is the method's
+/// miss or the setting's.
+pub fn crossings_beyond(geom: &Detail, rings: usize) -> usize {
+ if geom.num_points() == 0 || geom.num_prims() == 0 {
+ return 0;
+ }
+ let topo = topo_for(geom, rings);
+ if topo.edges.is_empty() {
+ return 0;
+ }
+ intersections_of(geom, &topo, true).crossings
+}
+
+fn intersections_of(geom: &Detail, topo: &Topo, beyond_rings: bool) -> Tangles {
+ let pos: Vec<Vec3> = (0..geom.num_points()).map(|p| geom.pos(p)).collect();
+ let grid = TriCells::build(&pos, &topo.tris, mean_edge(geom, topo));
+ let mut found = Tangles::default();
+ let mut marked = vec![false; pos.len()];
+ let mut near = Vec::new();
+ let mut seen = vec![u32::MAX; topo.tris.len()];
+ for (stamp, e) in topo.edges.iter().enumerate() {
+ let (a, b) = (pos[e[0] as usize], pos[e[1] as usize]);
+ grid.gather(a.min(b), a.max(b), &mut seen, stamp as u32, &mut near);
+ for &t in &near {
+ let tri = topo.tris[t as usize];
+ if tri.contains(&e[0]) || tri.contains(&e[1]) {
+ continue;
+ }
+ if beyond_rings && tri.iter().any(|&c| topo.excludes(e[0] as usize, c) || topo.excludes(e[1] as usize, c)) {
+ continue;
+ }
+ let of = geom.prim_points(topo.tri_prims[t as usize] as usize);
+ if of.contains(&e[0]) && of.contains(&e[1]) {
+ continue;
+ }
+ let [v0, v1, v2] = tri.map(|i| pos[i as usize]);
+ if crate::spatial::segment_crosses_triangle(a, b, v0, v1, v2) {
+ found.crossings += 1;
+ for i in e.iter().chain(tri.iter()) {
+ marked[*i as usize] = true;
+ }
+ }
+ }
+ }
+ found.points = (0..pos.len() as u32).filter(|&p| marked[p as usize]).collect();
+ found
}
diff --git a/src/main.rs b/src/main.rs
index 8ecb1c8..292d8f3 100644
--- a/src/main.rs
+++ b/src/main.rs
@@ -10134,6 +10134,212 @@ mod tests {
assert_eq!(a.positions(), b.positions());
}
+ /// A sheet of two large triangles with a fine patch over the middle of
+ /// one of them: the patch's points are nowhere near any corner of the
+ /// sheet, which is the case a point-to-point test cannot see. The patch
+ /// is tilted by `tilt` so its edges straddle the sheet when it is
+ /// lowered through it, and is the group "patch".
+ fn sheet_and_patch(height: f32, tilt: f32) -> Detail {
+ let mut d = Detail::new();
+ let s: Vec<u32> = [(-2.0, -2.0), (2.0, -2.0), (2.0, 2.0), (-2.0, 2.0)]
+ .iter()
+ .map(|&(x, z)| d.add_point(Vec3::new(x, 0.0, z)))
+ .collect();
+ d.add_prim(&[s[0], s[2], s[1]]);
+ d.add_prim(&[s[0], s[3], s[2]]);
+ let mut patch = Vec::new();
+ for i in 0..5 {
+ for j in 0..5 {
+ let (x, z) = (i as f32 * 0.1, j as f32 * 0.1);
+ patch.push(d.add_point(Vec3::new(0.9 + x, height + tilt * x, -1.1 + z)));
+ }
+ }
+ for i in 0..4 {
+ for j in 0..4 {
+ let at = |a: usize, b: usize| patch[a * 5 + b];
+ d.add_prim(&[at(i, j), at(i + 1, j + 1), at(i + 1, j)]);
+ d.add_prim(&[at(i, j), at(i, j + 1), at(i + 1, j + 1)]);
+ }
+ }
+ for &p in &patch {
+ d.points_mut().add_to_group("patch", p as usize);
+ }
+ d
+ }
+
+ /// How far the patch's nearest point is from the sheet's surface.
+ fn patch_clearance(d: &Detail) -> f32 {
+ let tris = [[0usize, 2, 1], [0, 3, 2]];
+ (4..d.num_points())
+ .flat_map(|p| tris.iter().map(move |t| (p, *t)))
+ .map(|(p, t)| {
+ let q = crate::spatial::closest_point_on_triangle(d.pos(p), d.pos(t[0]), d.pos(t[1]), d.pos(t[2]));
+ (d.pos(p) - q).length()
+ })
+ .fold(f32::MAX, f32::min)
+ }
+
+ /// The measure: every edge that passes through a triangle, and the
+ /// points of both. A surface that does not cross itself has none,
+ /// quads and shared corners included.
+ #[test]
+ fn the_tangle_measure_counts_edges_through_triangles() {
+ use crate::detangle::self_intersections;
+ assert_eq!(self_intersections(&Detail::new()).crossings, 0);
+ assert_eq!(self_intersections(&sphere_detail(Vec3::ZERO, 0.5, 10, 14)).crossings, 0);
+ assert_eq!(self_intersections(&box_detail(Vec3::ZERO, Vec3::ONE, 0.2)).crossings, 0, "quads fan into triangles that share a side");
+ // The patch above the sheet, then tilted through it.
+ assert_eq!(self_intersections(&sheet_and_patch(0.05, 0.0)).crossings, 0);
+ let through = self_intersections(&sheet_and_patch(-0.06, 0.3));
+ assert!(through.crossings > 0, "{through:?}");
+ // The points named are the crossing edges' and the crossed
+ // triangle's: some of the patch, not all of it, and the sheet's.
+ let of_patch = through.points.iter().filter(|&&p| p >= 4).count();
+ assert!(of_patch > 0 && of_patch < 25, "{through:?}");
+ assert!(through.points.iter().any(|&p| p < 4), "{through:?}");
+ assert!(through.points.windows(2).all(|w| w[0] < w[1]), "ascending, each once");
+ // The same answer at a thousandth of the size: the tolerance is
+ // relative.
+ let mut small = sheet_and_patch(-0.06, 0.3);
+ for p in 0..small.num_points() {
+ let v = small.pos(p);
+ small.set_pos(p, v * 0.001);
+ }
+ assert_eq!(self_intersections(&small), through);
+ }
+
+ /// What the Surface method is for: a point over the middle of a large
+ /// triangle is near none of its corners, so Points sees nothing and
+ /// Surface parts them.
+ #[test]
+ fn the_surface_method_sees_a_point_over_the_middle_of_a_triangle() {
+ let start = sheet_and_patch(0.03, 0.0);
+ let settings = |method: &'static str| {
+ phase3_node("detangle", &[("Method", method), ("Thickness", "0.25"), ("Rings", "2"), ("Iterations", "8")])
+ };
+ let mut by_points = start.clone();
+ let work = crate::detangle::apply(&mut by_points, &settings("Points"));
+ assert_eq!(by_points.positions(), start.positions(), "nothing is within a thickness of a corner: {work:?}");
+
+ let mut by_surface = start.clone();
+ let work = crate::detangle::apply(&mut by_surface, &settings("Surface"));
+ assert!(work.contacts > 0, "{work:?}");
+ let (before, after) = (patch_clearance(&start), patch_clearance(&by_surface));
+ assert!((before - 0.03).abs() < 1e-5);
+ assert!(after > 0.1, "the patch is {after} from the sheet");
+ // The patch moved as one: it was pushed off the sheet, not apart.
+ let edge = |d: &Detail| (d.pos(5) - d.pos(4)).length();
+ assert!((edge(&by_surface) / edge(&start) - 1.0).abs() < 0.05, "{} to {}", edge(&start), edge(&by_surface));
+ // The move is shared, so the sheet gave way too, downward.
+ assert!((0..4).all(|p| by_surface.pos(p).y <= 0.0) && (0..4).any(|p| by_surface.pos(p).y < 0.0));
+
+ // With the sheet outside the Group it stays where it is and the
+ // patch takes the whole move, not half of it.
+ let node = phase3_node(
+ "detangle",
+ &[("Method", "Surface"), ("Thickness", "0.25"), ("Rings", "2"), ("Iterations", "8"), ("Group", "patch")],
+ );
+ let mut held = start.clone();
+ crate::detangle::apply(&mut held, &node);
+ assert_eq!(held.positions()[..4], start.positions()[..4]);
+ assert!(patch_clearance(&held) > 0.1, "{}", patch_clearance(&held));
+
+ // A surface that touches itself nowhere is left alone, in one pass.
+ let round = sphere_detail(Vec3::ZERO, 0.5, 10, 14);
+ let mut d = round.clone();
+ let node = phase3_node("detangle", &[("Method", "Surface"), ("Thickness", "1.00"), ("Rings", "2"), ("Iterations", "8")]);
+ let work = crate::detangle::apply(&mut d, &node);
+ assert_eq!((work.passes, work.contacts), (1, 0), "{work:?}");
+ assert_eq!(d.positions(), round.positions());
+ }
+
+ /// The two methods as a simulation runs them: the patch is lowered a
+ /// little each step, by less than the thickness, toward a sheet that
+ /// does not move. Points lets it through; Surface holds it off. The
+ /// measure is what says so, and the node's Tangled Group is the measure
+ /// written onto the geometry.
+ #[test]
+ fn the_surface_method_holds_a_patch_off_a_sheet_step_after_step() {
+ let run = |method: &'static str| {
+ let node = phase3_node(
+ "detangle",
+ &[("Method", method), ("Thickness", "0.25"), ("Rings", "2"), ("Iterations", "8"), ("Group", "patch"), ("Tangled Group", "tangled")],
+ );
+ let mut d = sheet_and_patch(0.2, 0.3);
+ let (mut worst, mut marked) = (0, 0);
+ for _ in 0..25 {
+ for p in 4..d.num_points() {
+ let v = d.pos(p);
+ d.set_pos(p, v - Vec3::new(0.0, 0.02, 0.0));
+ }
+ let work = crate::detangle::apply(&mut d, &node);
+ assert_eq!(work.crossings, crate::detangle::self_intersections(&d).crossings);
+ assert_eq!(work.crossings > 0, d.points().group_len("tangled") > 0);
+ worst = worst.max(work.crossings);
+ marked = marked.max(d.points().group_len("tangled"));
+ }
+ (d, worst, marked)
+ };
+ let (through, worst, marked) = run("Points");
+ assert!(worst > 0 && marked > 0, "Points should have let it through: {worst}");
+ assert!((4..through.num_points()).all(|p| through.pos(p).y < 0.0), "and out the other side");
+
+ let (held, worst, marked) = run("Surface");
+ assert_eq!((worst, marked), (0, 0), "no edge went through at any step");
+ assert!((4..held.num_points()).all(|p| held.pos(p).y > 0.0), "the patch is still above the sheet");
+ assert!(held.points().has_group("tangled"), "the group is written, empty");
+ }
+
+ /// The two methods side by side, by the measure and by the clock: a
+ /// sphere's cap pushed down into its own bowl a little each step, until
+ /// it would have come out underneath. Run in release with `--ignored
+ /// --nocapture`. (Not pressed flat by the sign of y: that carries the
+ /// equator's points past their own neighbours, inside the excluded
+ /// rings, which no setting of the node is meant to see.)
+ #[test]
+ #[ignore]
+ fn detangle_methods_compared() {
+ for frequency in ["4", "8", "16"] {
+ let sphere = crate::shapes::sphere_node_detail(
+ &phase3_node("sphere", &[("Method", "Icosphere"), ("Frequency", frequency), ("Radius", "0.5")]),
+ Some(Vec3::ZERO),
+ );
+ let edges = sphere.edges();
+ let edge = edges.iter().map(|e| (sphere.pos(e[1] as usize) - sphere.pos(e[0] as usize)).length()).sum::<f32>() / edges.len() as f32;
+ // A fifth of an edge a step, until the pole has travelled the
+ // diameter and a little more.
+ let rate = edge * 0.2;
+ let steps = (1.1 / rate).ceil() as usize;
+ let cap: Vec<usize> = (0..sphere.num_points()).filter(|&p| sphere.pos(p).y > 0.2).collect();
+ for thickness in ["0.50", "1.00"] {
+ for method in ["None", "Points", "Surface"] {
+ let node = phase3_node("detangle", &[("Method", method), ("Thickness", thickness), ("Rings", "2"), ("Iterations", "4")]);
+ let mut d = sphere.clone();
+ let (mut worst, mut far, mut spent) = (0, 0, std::time::Duration::ZERO);
+ for _ in 0..steps {
+ for &p in &cap {
+ let v = d.pos(p);
+ d.set_pos(p, v - Vec3::new(0.0, rate, 0.0));
+ }
+ if method != "None" {
+ let t = std::time::Instant::now();
+ crate::detangle::apply(&mut d, &node);
+ spent += t.elapsed();
+ }
+ worst = worst.max(crate::detangle::self_intersections(&d).crossings);
+ far = far.max(crate::detangle::crossings_beyond(&d, 2));
+ }
+ let last = crate::detangle::self_intersections(&d).crossings;
+ println!(
+ "{:>5} points, {steps:>3} steps, thickness {thickness}, {method:>7}: worst {worst:>5} crossings ({far:>5} beyond the rings), last {last:>5}, {:.2} ms a step",
+ d.num_points(),
+ spent.as_secs_f64() * 1000.0 / steps as f64
+ );
+ }
+ }
+ }
+ }
+
#[test]
fn test_suture_counts_sustained_contact_before_it_fuses() {
// A grid sitting just above a collider it is in contact with.
diff --git a/src/spatial.rs b/src/spatial.rs
index c166b0e..2b913aa 100644
--- a/src/spatial.rs
+++ b/src/spatial.rs
@@ -59,6 +59,82 @@ pub fn closest_point_on_triangle(p: Vec3, a: Vec3, b: Vec3, c: Vec3) -> Vec3 {
a + ab * (vb / denom) + ac * (vc / denom)
}
+/// [`closest_point_on_triangle`] as WEIGHTS: how much of the closest point
+/// each corner is, summing to one. The point is `a * w[0] + b * w[1] +
+/// c * w[2]`, and the weights are what lets a caller hand a push on that
+/// point back to the corners that carry it.
+///
+/// The same Voronoi-region walk, region for region.
+pub fn closest_weights_on_triangle(p: Vec3, a: Vec3, b: Vec3, c: Vec3) -> [f32; 3] {
+ let (ab, ac, ap) = (b - a, c - a, p - a);
+ let (d1, d2) = (ab.dot(ap), ac.dot(ap));
+ if d1 <= 0.0 && d2 <= 0.0 {
+ return [1.0, 0.0, 0.0];
+ }
+ let bp = p - b;
+ let (d3, d4) = (ab.dot(bp), ac.dot(bp));
+ if d3 >= 0.0 && d4 <= d3 {
+ return [0.0, 1.0, 0.0];
+ }
+ let vc = d1 * d4 - d3 * d2;
+ if vc <= 0.0 && d1 >= 0.0 && d3 <= 0.0 {
+ let denom = d1 - d3;
+ let v = if denom.abs() < 1e-20 { 0.0 } else { d1 / denom };
+ return [1.0 - v, v, 0.0];
+ }
+ let cp = p - c;
+ let (d5, d6) = (ab.dot(cp), ac.dot(cp));
+ if d6 >= 0.0 && d5 <= d6 {
+ return [0.0, 0.0, 1.0];
+ }
+ let vb = d5 * d2 - d1 * d6;
+ if vb <= 0.0 && d2 >= 0.0 && d6 <= 0.0 {
+ let denom = d2 - d6;
+ let w = if denom.abs() < 1e-20 { 0.0 } else { d2 / denom };
+ return [1.0 - w, 0.0, w];
+ }
+ let va = d3 * d6 - d5 * d4;
+ if va <= 0.0 && (d4 - d3) >= 0.0 && (d5 - d6) >= 0.0 {
+ let denom = (d4 - d3) + (d5 - d6);
+ let w = if denom.abs() < 1e-20 { 0.0 } else { (d4 - d3) / denom };
+ return [0.0, 1.0 - w, w];
+ }
+ let denom = va + vb + vc;
+ if denom.abs() < 1e-20 {
+ return [1.0, 0.0, 0.0];
+ }
+ let (v, w) = (vb / denom, vc / denom);
+ [1.0 - v - w, v, w]
+}
+
+/// Whether the segment `a`-`b` passes through triangle `(v0, v1, v2)`,
+/// strictly between its ends.
+///
+/// Möller–Trumbore with the segment as the ray, and a tolerance RELATIVE to
+/// the lengths involved, so the answer does not change with the model's
+/// scale. A segment lying in the triangle's plane does not cross it.
+pub fn segment_crosses_triangle(a: Vec3, b: Vec3, v0: Vec3, v1: Vec3, v2: Vec3) -> bool {
+ let (d, e1, e2) = (b - a, v1 - v0, v2 - v0);
+ let h = d.cross(e2);
+ let det = e1.dot(h);
+ if det.abs() <= 1e-7 * e1.length() * e2.length() * d.length() {
+ return false;
+ }
+ let f = 1.0 / det;
+ let s = a - v0;
+ let u = f * s.dot(h);
+ if !(0.0..=1.0).contains(&u) {
+ return false;
+ }
+ let q = s.cross(e1);
+ let v = f * d.dot(q);
+ if v < 0.0 || u + v > 1.0 {
+ return false;
+ }
+ let t = f * e2.dot(q);
+ t > 0.0 && t < 1.0
+}
+
/// Where a ray meets a triangle, as a distance along the ray.
///
/// Möller–Trumbore. Lives here beside the other spatial queries because three