diff --git a/crates/mesh/src/project/mod.rs b/crates/mesh/src/project/mod.rs
index 29a888d..4e7c92d 100644
--- a/crates/mesh/src/project/mod.rs
+++ b/crates/mesh/src/project/mod.rs
@@ -136,6 +136,63 @@ pub struct ProjectedMesh {
     pub report: ProjectReport,
 }
 
+/// Give each target face no boundary quad's centre landed on (a face narrower than the elements:
+/// a fillet's sliver between a cut and the walls) one row of the quads it runs through (of a
+/// target holding the face, centre within half the quad's shortest side of it): one quad of
+/// those beside each other across the face (projecting onto it at nearly the same point), all
+/// along it. A quad keeps its face when it is that face's
+/// last, so no face loses its only quad to a neighbour.
+fn claim_narrow_faces(boundary: &mut [BoundaryFace], centres: &[(P, f64, usize)], targets: &[Vec<u32>], g: &dyn Geometry) {
+    let mut count: BTreeMap<u32, usize> = BTreeMap::new();
+    for b in boundary.iter() {
+        *count.entry(b.brep_face).or_insert(0) += 1;
+    }
+    let targeted: BTreeSet<u32> = targets.iter().flatten().copied().collect();
+    for f in targeted {
+        if count.contains_key(&f) {
+            continue;
+        }
+        // (distance to the face, quad, its centre, the centre's foot on the face, shortest side)
+        let mut near: Vec<(f64, usize, P, P, f64)> = centres
+            .iter()
+            .enumerate()
+            .filter(|(_, (_, _, tg))| targets.get(*tg).is_some_and(|set| set.contains(&f)))
+            .filter_map(|(i, (centre, shortest, _))| {
+                let foot = g.project(f, *centre)?;
+                let d = dist(foot, *centre);
+                (d <= 0.5 * shortest).then_some((d, i, *centre, foot, *shortest))
+            })
+            .collect();
+        near.sort_by(|a, b| a.0.total_cmp(&b.0));
+        // The rows across the face: quads beside each other that project onto it at nearly the
+        // same point.
+        let mut rows: Vec<Vec<(usize, P, P, f64)>> = Vec::new();
+        for (_, i, centre, foot, shortest) in near {
+            let across = |(_, c, q, s): &(usize, P, P, f64)| {
+                let apart = dist(*c, centre);
+                apart < 1.5 * s.max(shortest) && dist(*q, foot) < 0.5 * apart
+            };
+            match rows.iter_mut().find(|row| row.iter().any(across)) {
+                Some(row) => row.push((i, centre, foot, shortest)),
+                None => rows.push(vec![(i, centre, foot, shortest)]),
+            }
+        }
+        // One quad per row: from the neighbour with the most quads left (the nearest of equals),
+        // so neighbours on both sides give alike.
+        for row in rows {
+            let left = |i: usize| boundary.get(i).and_then(|b| count.get(&b.brep_face)).copied().unwrap_or(0);
+            let Some(&(i, ..)) = row.iter().max_by_key(|(i, ..)| left(*i)) else { continue };
+            let Some(b) = boundary.get_mut(i) else { continue };
+            let Some(left) = count.get_mut(&b.brep_face) else { continue };
+            if *left > 1 {
+                *left -= 1;
+                b.brep_face = f;
+                *count.entry(f).or_insert(0) += 1;
+            }
+        }
+    }
+}
+
 /// The face groups joined across smooth B-rep edges.
 fn smooth_groups(surface: &BoundarySurface, faces: usize) -> Vec<usize> {
     let mut parent: Vec<usize> = (0..faces).collect();
@@ -626,6 +683,7 @@ pub fn project_mesh(t: &Topology, surface: &BoundarySurface, g: &dyn Geometry, o
         }
     }
     let mut boundary = Vec::with_capacity(quads.len());
+    let mut centres = Vec::with_capacity(quads.len());
     let mut enclosed = 0.0;
     for (q, e, f, tg) in &quads {
         let pts = q.map(|v| nodes.get(v as usize).copied().unwrap_or_default());
@@ -634,7 +692,10 @@ pub fn project_mesh(t: &Topology, surface: &BoundarySurface, g: &dyn Geometry, o
         let set = targets.get(*tg).map_or(&[][..], Vec::as_slice);
         let brep_face = on_target(g, set, centre).map_or_else(|| set.first().copied().unwrap_or(0), |h| h.2);
         boundary.push(BoundaryFace { element: *e, face: *f, brep_face });
+        let shortest = (0..4).map(|i| dist(pts[i], pts[(i + 1) % 4])).fold(f64::INFINITY, f64::min);
+        centres.push((centre, shortest, *tg));
     }
+    claim_narrow_faces(&mut boundary, &centres, &targets, g);
     let surface_volume = surface.volume();
     let report = ProjectReport {
         blocks: t.blocks.len(),
diff --git a/crates/mesh/src/tet/geometry.rs b/crates/mesh/src/tet/geometry.rs
index 03fcd3f..c12e86e 100644
--- a/crates/mesh/src/tet/geometry.rs
+++ b/crates/mesh/src/tet/geometry.rs
@@ -507,6 +507,19 @@ impl Tessellation {
         self.face_tree(face).nearest4(p, &accept, &dist).map(|r| ([r[0], r[1], r[2]], r[3]))
     }
 
+    /// The B-rep face whose triangles lie nearest to `p`.
+    pub fn nearest_face(&self, p: P) -> Option<u32> {
+        let dist = |t: u32| {
+            let tri = self.triangles.get(t as usize).copied().unwrap_or([0; 3]);
+            let d = sub(closest_on_triangle(p, self.pt(tri[0]), self.pt(tri[1]), self.pt(tri[2])), p);
+            (dot(d, d), [0.0, 0.0, 0.0, f64::from(t)])
+        };
+        let r = self.tri_tree.nearest4(p, &|_| true, &dist)?;
+        // The triangle index travels through the f64 exactly (triangle counts are capped far
+        // below 2^53).
+        self.tri_face.get(r[3] as usize).copied()
+    }
+
     /// The tree over `face`'s triangles (the whole boundary's for a face out of range, which
     /// the queries' filter then finds nothing in).
     fn face_tree(&self, face: u32) -> &Tree {
@@ -528,38 +541,51 @@ impl Tessellation {
     }
 }
 
-/// How far off a plane face's triangles, in triangle sizes, a point of the plane may lie and
-/// still be on the face (rounding only: the triangles lie on the plane).
+/// How far off a face's triangles, in triangle sizes, a point of its surface may lie and still
+/// be on the face, beyond the triangles' own gap to the surface (rounding only).
 const PLANE_SLACK: f64 = 1e-6;
 
 impl Tessellation {
     /// The closest point to `p` on the trimmed face `face`: the exact surface's closest point
-    /// when it lies on the face (within a fifth of the local triangle size of the face's
-    /// triangles, or on them for a plane, whose triangles can span the whole face), else the closest point of the face's triangles. Unlike
-    /// [`BoundaryGeometry::project_face`] it never returns a point of the surface's extension
-    /// beyond the face's boundary, whatever the distance from `p`.
+    /// when it lies on the face (no further off the face's nearest triangle than twice that
+    /// triangle's own chordal gap to the surface, which is rounding only for a plane), else the
+    /// closest point of the face's triangles. Unlike [`BoundaryGeometry::project_face`] it never
+    /// returns a point of the surface's extension beyond the face's boundary, whatever the
+    /// distance from `p` — nor of the same surface's other faces: the pieces of one fillet split
+    /// by a cut share their cylinder, and a point on one piece is off the others.
     pub fn project_trimmed(&self, face: u32, p: P) -> Option<P> {
         let accept = |t: u32| self.tri_face.get(t as usize) == Some(&face);
-        let nearest = |x: P| {
+        let corners = |t: u32| {
+            let tri = self.triangles.get(t as usize).copied().unwrap_or([0; 3]);
+            (self.pt(tri[0]), self.pt(tri[1]), self.pt(tri[2]))
+        };
+        let near = {
             let dist = |t: u32| {
-                let tri = self.triangles.get(t as usize).copied().unwrap_or([0; 3]);
-                let (a, b, c) = (self.pt(tri[0]), self.pt(tri[1]), self.pt(tri[2]));
-                let q = closest_on_triangle(x, a, b, c);
-                let d = sub(q, x);
-                let mut e = [measure::dist(a, b), measure::dist(b, c), measure::dist(c, a)];
-                e.sort_by(f64::total_cmp);
-                (dot(d, d), [q[0], q[1], q[2], e[1]])
+                let (a, b, c) = corners(t);
+                let q = closest_on_triangle(p, a, b, c);
+                let d = sub(q, p);
+                (dot(d, d), [q[0], q[1], q[2], 0.0])
             };
-            self.face_tree(face).nearest4(x, &accept, &dist).map(|r| ([r[0], r[1], r[2]], r[3]))
+            let r = self.face_tree(face).nearest4(p, &accept, &dist)?;
+            [r[0], r[1], r[2]]
         };
-        let (near, _) = nearest(p)?;
         let Some(surface) = self.analytic(face) else { return Some(near) };
         let Some(q) = surface.project(p) else { return Some(near) };
-        let (back, size) = nearest(q)?;
-        // A plane's triangles lie on it, so a point of the plane off them is off the face; a
-        // curved face's triangles are chords, a fraction of their size off the surface.
-        let slack = if matches!(surface, AnalyticSurface::Plane { .. }) { PLANE_SLACK } else { 0.2 };
-        Some(if measure::dist(back, q) <= slack * size { q } else { near })
+        // The nearest triangle to the surface's point, with its slack: the largest gap between
+        // the triangle and the surface (at its edges' midpoints and centre: a chord's deepest
+        // points), doubled, and a millionth of its size for rounding.
+        let gap = |x: P| surface.project(x).map_or(0.0, |y| measure::dist(x, y));
+        let dist = |t: u32| {
+            let (a, b, c) = corners(t);
+            let o = closest_on_triangle(q, a, b, c);
+            let d = sub(o, q);
+            let mid = |u: P, v: P| scale(add(u, v), 0.5);
+            let deepest = [mid(a, b), mid(b, c), mid(c, a), scale(add(add(a, b), c), 1.0 / 3.0)].into_iter().map(gap).fold(0.0, f64::max);
+            let size = measure::dist(a, b).max(measure::dist(b, c)).max(measure::dist(c, a));
+            (dot(d, d), [o[0], o[1], o[2], 2.0 * deepest + PLANE_SLACK * size])
+        };
+        let r = self.face_tree(face).nearest4(q, &accept, &dist)?;
+        Some(if measure::dist([r[0], r[1], r[2]], q) <= r[3] { q } else { near })
     }
 
     /// The outward unit normal of B-rep face `face` near `p`: the exact surface's (oriented like