kernel: fillet_g2_corners, set-back corners filled G2 to every face round them

A corner where three or more rounded edges meet, every edge there rounded, is
set back (1.5 r unless a setback names it) and left open by SMLib
(blend_smlib::fillet_setback_open: the sheet and each corner's hole, in the
sheet's orders); each hole is filled by fill::matched_patch, the height field
over its plane fitted to the position, tangent plane and normal curvature
across the boundary of every face beside it (exact second derivatives), rounds
and planes alike, each edge's pcurve its exact projection; SMLib sews the
sheet and the patches into one solid, which must close. The history is the
sewing's carried through the sheet's, each patch face made from the rounded
edges its neighbours come from.

corner_g2 is that construction now. The polar patch over three sections it
held (never sewn in) could not be: its boundary is a closed cubic through
stations, which rounds off a corner of the loop, and the setback hole has six.

The 60 x 40 x 20 box corner of r5 rounds (07/015): 10 faces, a valid closed
solid (record and SMLib's check), curvature jump across the patch's seams
6.3e-4 per mm clear of the hole's corners against the ball's 0.2; set back 8
between the hull and prism bounds; set back 11 to 15 inside the box but by
under 7e-4; all eight corners of a box filled; refusals named.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
This commit is contained in:
jyviko
2026-10-26 12:32:39 -04:00
committed by Irakli Kourtis
parent d63a50c47c
commit 704640375d
6 changed files with 441 additions and 641 deletions
+112 -6
View File
@@ -228,17 +228,123 @@ pub fn fillet_setback(body: &Body, edges: &[Vec3], radius: f64, setbacks: &[(Vec
crate::blend_smlib::fillet_setback(body, edges, radius, setbacks).and_then(|b| closed(body, b))
}
/// [`fillet`] whose corners where three blended edges meet are G2 patches. SMLib's fillet
/// corners are a sphere or a Coons patch that meets the rounds tangentially, so this is refused
/// rather than quietly given those corners; [`fillet`] (the default) is unchanged.
pub fn fillet_g2_corners(body: &Body, edges: &[Vec3], radius: f64) -> Result<Body> {
/// [`fillet`] whose corners where three or more rounded edges meet are G2 patches
/// ([`crate::corner_g2`]): each such corner is set back, `setbacks` naming a corner by a point at
/// its vertex and the distance its patch starts along each edge there (the others take
/// [`crate::DEFAULT_SETBACK_FACTOR`] times the radius), SMLib leaves it open, and a patch that
/// meets the rounds and the faces round it to the second order fills it. The faces are sewn
/// into one solid. Without such a corner this is [`fillet`]. Refused: a corner where an edge is
/// not rounded (SMLib sets back only corners where all are), a setback that does not reach past
/// where the rounds would meet, a hole that is no height field over its plane or a patch that
/// does not meet its neighbours within the tolerances.
pub fn fillet_g2_corners(body: &Body, edges: &[Vec3], radius: f64, setbacks: &[(Vec3, f64)]) -> Result<Body> {
let _span = solidcraft_profile::span("kernel", "fillet_g2_corners");
body.require_brep("fillet")?;
if !(radius.is_finite() && radius > 0.0) {
return Err(KernelError::Invalid("the blend radius must be finite and positive".into()));
}
let _ = edges;
Err(crate::blend_smlib::unsupported("a fillet with G2 corner patches"))
let (picked, _) = crate::blend_smlib::picked_edges(body, edges)?;
let rec = body.record();
let picked: BTreeSet<usize> = picked.into_iter().collect();
// Each vertex's edges (a closed edge's one vertex once).
let mut at: std::collections::BTreeMap<u32, BTreeSet<usize>> = std::collections::BTreeMap::new();
for (i, e) in rec.edges.iter().enumerate() {
if e.uses.is_empty() {
continue;
}
for v in e.vertices {
at.entry(v).or_default().insert(i);
}
}
let point = |v: u32| rec.vertices.get(v as usize).map(|x| v3(x.point));
let mut corners: Vec<Vec3> = Vec::new();
for (v, around) in &at {
let rounded = around.iter().filter(|e| picked.contains(e)).count();
if rounded < 3 {
continue;
}
let p = point(*v).ok_or_else(|| KernelError::Failed("fillet: a vertex reference".into()))?;
if rounded < around.len() {
return Err(unsupported(&format!(
"a G2 corner at ({:.3}, {:.3}, {:.3}) where an edge is not rounded (the corner is set back, and SMLib sets back only a corner where every edge is)",
p.x, p.y, p.z
)));
}
corners.push(p);
}
if corners.is_empty() {
return fillet(body, edges, radius);
}
// Each corner's setback: the one given at its vertex, else the default.
let reach = body.size() * 0.05 + 1e-3;
let mut distance = vec![crate::corner_g2::DEFAULT_SETBACK_FACTOR * radius; corners.len()];
for (p, d) in setbacks {
let k = corners
.iter()
.enumerate()
.map(|(k, c)| (c.dist(*p), k))
.filter(|(d, _)| *d <= reach)
.min_by(|a, b| a.0.total_cmp(&b.0))
.map(|(_, k)| k)
.ok_or_else(|| {
KernelError::Invalid(format!(
"the setback at ({:.3}, {:.3}, {:.3}) names no G2 corner (a vertex where three or more rounded edges meet)",
p.x, p.y, p.z
))
})?;
if let Some(slot) = distance.get_mut(k) {
*slot = *d;
}
}
let set_back: Vec<(Vec3, f64)> = corners.iter().copied().zip(distance).collect();
let open = crate::blend_smlib::fillet_setback_open(body, edges, radius, &set_back)?;
let mut patches = Vec::with_capacity(open.holes.len());
for hole in &open.holes {
let (patch, _) = crate::corner_g2::corner_patch(&open.sheet, hole, radius)?;
patches.push(patch);
}
let mut parts: Vec<&Body> = vec![&open.sheet];
parts.extend(patches.iter());
let sewn = crate::stitch(&parts, 0.0)?;
let open_shells = sewn.validity();
if !open_shells.is_empty() {
return Err(KernelError::Failed(format!("fillet: the G2 corners do not close the body: {}", open_shells.join("; "))));
}
let history = corner_history(&sewn, body, &open.sheet, &patches);
closed(body, sewn.with_checked_history(history, &crate::record::counts(&[body])))
}
/// The history of the sewn G2-cornered fillet over `[body]`: the sewing's history (over the
/// sheet and the patches) carried through the sheet's (over the body); each patch face made from
/// the edges whose rounds meet at its corner, as SMLib's own corner patch is, each of its edges
/// the sheet's.
fn corner_history(sewn: &Body, body: &Body, sheet: &Body, patches: &[Body]) -> Option<crate::history::History> {
use crate::history::{EdgeOrigin, FaceOrigin};
let outer = sewn.history()?.clone();
let inner = sheet.history()?;
let mut parts = vec![inner.clone()];
for patch in patches {
let own = patch.history()?;
// The patch's face comes from the rounded edges the sheet's faces beside it come from.
let mut origins: BTreeSet<FaceOrigin> = BTreeSet::new();
for e in &own.edges {
for o in &e.origins {
let EdgeOrigin::Edge { edge, .. } = o else { continue };
for f in inner.edges.get(*edge).map(|x| x.faces.as_slice()).unwrap_or_default() {
let beside = inner.faces.get(*f).map(|x| x.origins.as_slice()).unwrap_or_default();
origins.extend(beside.iter().copied().filter(|g| matches!(g, FaceOrigin::Edge { .. })));
}
}
}
let mut h = own.clone().substituted(0, inner).ok()?;
for f in &mut h.faces {
*f = crate::record::generated(origins.iter().copied().collect());
}
parts.push(h);
}
let h = outer.over(&parts.iter().collect::<Vec<_>>()).ok()?;
h.check_counts((sewn.face_count(), sewn.edge_count()), &crate::record::counts(&[body])).ok()?;
Some(h)
}
/// The step of the radius search: three significant digits of the requested radius.
+68 -4
View File
@@ -43,6 +43,22 @@ pub(crate) fn fillet_setback(body: &Body, edges: &[Vec3], radius: f64, setbacks:
imp::blend_setback(body, edges, &vec![radius; edges.len()], setbacks)
}
/// A set-back blend whose corners SMLib leaves open ([`fillet_setback_open`]): the sheet, its
/// history over the body it was made from, and per corner the edges round its hole in order,
/// each its index in the sheet's edge order with whether the face beside it runs it forwards
/// (counter-clockwise about the solid's outward normal).
pub(crate) struct OpenCorners {
pub(crate) sheet: Body,
pub(crate) holes: Vec<Vec<(usize, bool)>>,
}
/// [`fillet_setback`] with each set-back corner left open: SMLib builds the rounds up to the
/// corners' setback sections and the curves on the faces between them, and leaves out the patch
/// that would fill each corner, for the caller to fill.
pub(crate) fn fillet_setback_open(body: &Body, edges: &[Vec3], radius: f64, setbacks: &[(Vec3, f64)]) -> Result<OpenCorners> {
imp::blend_setback_open(body, edges, &vec![radius; edges.len()], setbacks)
}
/// Curvature-continuous rounds with the setback of a round of `radius` ([`crate::fillet_g2`]).
pub(crate) fn fillet_g2(body: &Body, edges: &[Vec3], radius: f64) -> Result<Body> {
imp::blend(body, "G2 fillet", Kind::Curvature, edges, &vec![radius; edges.len()])
@@ -204,6 +220,43 @@ mod imp {
run(body, what, &picked.edges, |b| b.blend_setback(Blend::Round, &picked.points, sizes, &corners, &distances, tol))
}
/// [`blend_setback`] with the set-back corners left open (`super::OpenCorners`).
pub(super) fn blend_setback_open(body: &Body, edges: &[Vec3], sizes: &[f64], setbacks: &[(Vec3, f64)]) -> Result<super::OpenCorners> {
let what = "fillet with open setback corners";
body.require_brep(what)?;
sizes_ok(what, sizes)?;
if setbacks.is_empty() || setbacks.len() > solidcraft_smlib::MAX_POINTS {
return Err(KernelError::Invalid(format!("{what}: give 1…{} corners", solidcraft_smlib::MAX_POINTS)));
}
if setbacks.iter().any(|(p, d)| !p.is_finite() || !(d.is_finite() && *d > 1e-6 && *d < 1e6)) {
return Err(KernelError::Invalid(format!("{what}: each corner is a finite point with a positive, finite setback")));
}
let picked = resolve(body, edges)?;
let corners = corner_vertices(body, setbacks.iter().map(|(p, _)| *p))?;
let distances: Vec<f64> = setbacks.iter().map(|(_, d)| *d).collect();
let tol = pick_tolerance(body);
let (sheet, holes, order) =
run_with(body, what, &picked.edges, |b| b.blend_setback_open(Blend::Round, &picked.points, sizes, &corners, &distances, tol))?;
// SMLib's indices in the sheet's own orders.
let place = |list: &[usize], i: usize| if list.is_empty() { Some(i) } else { list.get(i).copied() };
let renumber = |hole: &solidcraft_smlib::Hole| -> Result<Vec<(usize, bool)>> {
hole.edges
.iter()
.map(|h| {
let edge = match &order {
Some(o) => place(&o.edges, h.edge),
None => Some(h.edge),
};
edge.filter(|e| *e < sheet.edge_count())
.map(|e| (e, h.forward))
.ok_or_else(|| failed(what, "hole", "an edge SMLib names is not in the sheet"))
})
.collect()
};
let holes = holes.iter().map(renumber).collect::<Result<Vec<_>>>()?;
Ok(super::OpenCorners { sheet, holes })
}
/// The body's vertex nearest each point, within the pick tolerance: the ends of its edges.
fn corner_vertices(body: &Body, picks: impl Iterator<Item = Vec3>) -> Result<Vec<[f64; 3]>> {
let mesh = body.mesh(&crate::MeshTolerance::chord((body.size() * 2e-3).max(1e-3)))?;
@@ -249,10 +302,21 @@ mod imp {
/// `op` on `body` in SMLib, tagged as operand 0, checked, decoded, with its history.
fn run(body: &Body, what: &str, edges: &[usize], op: impl FnOnce(&Brep) -> solidcraft_smlib::Result<Brep>) -> Result<Body> {
run_with(body, what, edges, |b| op(b).map(|r| (r, ()))).map(|(out, (), _)| out)
}
/// [`run`] for an `op` that also gives something of its own, with where SMLib's faces and
/// edges of the result are in the body's orders.
fn run_with<T>(
body: &Body,
what: &str,
edges: &[usize],
op: impl FnOnce(&Brep) -> solidcraft_smlib::Result<(Brep, T)>,
) -> Result<(Body, T, Option<solidcraft_smlib::history::Order>)> {
guard(what, || {
let (sb, order) = smlib_body_ordered(body)?;
let tagged = sb.tagged(0).map_err(|e| failed(what, "tag", e))?;
let r = op(&tagged).map_err(|e| match e {
let (r, extra) = op(&tagged).map_err(|e| match e {
SmError::Unsupported => super::unsupported(&format!("this {what}, which SMLib cannot build exactly")),
SmError::NoEffect => KernelError::Failed(format!("SMLib made no change: the {what} reported success but built no face")),
SmError::NoRoom => KernelError::TooLarge(format!("the {what} is larger than the faces beside its edges allow")),
@@ -266,8 +330,8 @@ mod imp {
return Err(failed(what, step, check));
}
let (out, result_order) = body_from_smlib_ordered(&r)?;
let history = match (order, result_order) {
(Some(o), Some(ro)) => match r.history().map_err(|e| failed(what, "history", e))?.history(&[o], &ro) {
let history = match (order, &result_order) {
(Some(o), Some(ro)) => match r.history().map_err(|e| failed(what, "history", e))?.history(&[o], ro) {
Ok(mut h) => {
let operand = History::identity(body).edges;
let replaced = blended(body, edges, &h);
@@ -281,7 +345,7 @@ mod imp {
},
_ => None,
};
Ok(out.with_checked_history(history, &record::counts(&[body])).inherit_regularization_protection(body))
Ok((out.with_checked_history(history, &record::counts(&[body])).inherit_regularization_protection(body), extra, result_order))
})
}
+37 -488
View File
@@ -1,502 +1,51 @@
//! G2 corner patches: the patch that fills a vertex where filleted edges meet, built to repeat each
//! blend's surface up to the second order instead of the rolling ball's sphere.
//! G2 fillet corners (`corner: G2`): the corner where three or more rounded edges meet, filled
//! with a patch that meets every face around it to the second order, instead of the rolling
//! ball's sphere (or SMLib's Coons patch), which is only tangent to the rounds and keeps its own
//! 1 / r across each seam.
//!
//! [`crate::blend_corner`] and [`crate::curvecorner`] fill the corner of a vertex with a piece of
//! the ball that rolls along the three edges (a sphere for equal radii, a rational patch for
//! unequal ones). That patch is tangent to every blend along the sections where they stop, so the
//! corner is G1 — but across a seam the ball keeps its own 1/r where a fillet is flat along its
//! axis, so the two surfaces meet with a jump in curvature.
//! The corner is set back: SMLib runs each round whole to its cross-section a setback from the
//! vertex and leaves the corner open ([`crate::blend_smlib::fillet_setback_open`]), so the hole's
//! boundary is the rounds' sections and, on each face between two rounds, SMLib's curve from one
//! round's contact line to the next. Each hole is filled with one face
//! ([`crate::fill::matched_patch`]): a bicubic height field over the hole's plane fitted to every
//! face beside the hole — its position, its tangent plane and its normal curvature across the
//! boundary, read off the face's own surface — the rounds and the planes alike, so the patch is
//! G2 to all of them. The faces are sewn into one solid by SMLib ([`crate::stitch`]).
//!
//! Here the corner is one tensor patch whose boundary is the blends' section curves themselves and
//! whose cross-sections carry each blend's own derivatives: the position, the direction the patch
//! leaves the blend's surface in, and the blend's curvature vector along that direction are the
//! first three control points of every cross-section, one control point per continuity condition,
//! exactly as [`crate::bridge`]'s sections do. The other end of every cross-section gathers at
//! `apex` — a point inside the material the corner caps — where the sections arrive in a common
//! tangent plane with a common curvature, so the patch closes there smoothly like a sphere's pole.
//! Along a seam the patch therefore repeats the blend up to the second order (G2).
//!
//! The boundary is a closed loop (the three sections joined corner to corner), so the direction
//! along it is a *closed* cubic B-spline: an interpolating curve with free ends swings wide at the
//! corner where the loop joins and folds the patch there. Where two blends meet at that corner
//! their sections leave in directions 90° apart and ask for different curvature at that one point
//! (a fillet is flat across its own seam, the other is not), which no single surface can hold; the
//! corner is a still point of the construction and [`seam_error`] is measured clear of it. A patch
//! that folds, or a corner with no room, is refused rather than replaced by the G1 ball.
//! Without a setback the boundary would be the three sections alone, meeting at points where
//! two rounds 90° apart ask for different curvatures; with it, the only points where two faces
//! disagree are the hole's corners, where a section ends on a plane: the round's G1 seam with the
//! plane, whose 1 / r jump every round has along its whole length, ends there, and no surface is
//! G2 to both the round and the plane at that point. The patch passes from one to the other
//! within 15% of each side from the corner (`solidcraft_geom::fill_matched::FADE_SHARE`) and
//! is G2 along the rest of every seam.
use solidcraft_geom::Vec3;
use solidcraft_geom::brep::{self as b, eval::SurfaceEval};
use solidcraft_geom::bspline::{HPoint, NurbsSurface};
use solidcraft_geom::fill_matched::MatchedFill;
use crate::blend::unsupported;
use crate::bridge::Frame;
use crate::body::Body;
use crate::{KernelError, Result};
/// A corner's name in an error.
fn fail(m: impl std::fmt::Display) -> KernelError {
KernelError::Failed(format!("G2 corner: {m}"))
}
/// A corner's default setback: this times its round's radius along each edge. It must reach
/// past where the rounds meet without one (the radius, at a square corner), and SMLib's own
/// patch bulges past the faces beyond about twice the radius; one and a half leaves the patch a
/// side as long as the round's section is wide, with room either way.
pub const DEFAULT_SETBACK_FACTOR: f64 = 1.5;
/// Conditions each end of a cross-section carries: position, tangent and curvature.
const END_CONDITIONS: usize = 3;
/// Degree across the patch: three conditions at each end, less the control point they share.
const ACROSS: usize = 2 * END_CONDITIONS - 1;
/// Degree along the seam loop.
const ALONG: usize = 3;
/// Fewest stations a section curve may be sampled at (a cubic through them).
const MIN_STATIONS: usize = 4;
/// Most stations one section curve may carry.
const MAX_STATIONS: usize = 512;
/// How much of a station's `leave` has to lie across its seam for the patch to have a direction to
/// leave the blend in at all.
const ACROSS_ENOUGH: f64 = 0.2;
/// How far the patch's normal curvature across a seam may be from the face's, clear of the hole's
/// corners: this share of the corner's round's curvature 1 / r, the jump a rolling ball makes.
pub const CURVATURE_SHARE: f64 = 5e-3;
/// A G2 corner patch: its surface, and its boundary as the section curves themselves, each an
/// interval of the surface's own boundary curve (the `v = 0` line, `seam_curve`) from one corner of
/// the vertex to the next, so that it lies on the patch exactly.
#[derive(Clone, Debug, PartialEq)]
pub struct CornerPatch {
pub surface: b::Surface,
pub seam_curve: b::Curve,
/// Each side's interval of `seam_curve` (its parameter is the surface's u).
pub seams: Vec<[f64; 2]>,
}
/// The G2 corner patch over three blends' section curves: `sides[i]` is one blend's section curve,
/// sampled in order from one corner of the vertex to the next, each side's last station the next
/// side's first (the three close a loop). `apex` is a point inside the material the corner caps;
/// the patch gathers there. Its boundary is the section curves themselves (one interval of the
/// patch's own boundary curve per side); the caller orients the face where it sews it in.
pub fn corner_g2_patch(sides: &[Vec<Frame>], apex: Vec3) -> Result<CornerPatch> {
build(sides, apex)
}
/// The seam's own direction at station `i` of a side, from a four-point central difference through
/// the stations (`crate::bridge::frames` samples a chain evenly, and a corner's own callers do the
/// same): a two-point chord would only be second order in the step and would tilt the corner's
/// tangent plane by more than G2 is judged at. `None` within two stations of a side's end.
fn seam_tangent(side: &[Frame], i: usize) -> Option<Vec3> {
let at = |j: isize| usize::try_from(j).ok().and_then(|j| side.get(j)).map(|f| f.p);
let (a, b, c, d) = (at(i as isize - 2)?, at(i as isize - 1)?, at(i as isize + 1)?, at(i as isize + 2)?);
(b * 8.0 - a - c * 8.0 + d).normalized()
}
/// The seam's direction at every station of the closed loop (the same difference, wrapping round
/// its ends), falling back to the chord to a neighbour where the four points are not there.
fn loop_tangents(frames: &[Frame]) -> Vec<Option<Vec3>> {
let m = frames.len();
if m == 0 {
return Vec::new();
/// The patch filling one open corner of `sheet` (a sheet from
/// [`crate::blend_smlib::fillet_setback_open`]): `hole` its edges in order, each with whether the
/// face beside it runs it forwards; `radius` the corner's round's radius, which sets the
/// curvature tolerance. The patch is a one-face sheet on the sheet's own edges.
pub(crate) fn corner_patch(sheet: &Body, hole: &[(usize, bool)], radius: f64) -> Result<(Body, MatchedFill)> {
if !(radius.is_finite() && radius > 0.0) {
return Err(KernelError::Invalid("a G2 corner's radius must be finite and positive".into()));
}
let at = |j: isize| frames.get(j.rem_euclid(m as isize) as usize).map(|f| f.p);
(0..m)
.map(|i| {
let four = || {
let (a, b, c, d) = (at(i as isize - 2)?, at(i as isize - 1)?, at(i as isize + 1)?, at(i as isize + 2)?);
(b * 8.0 - a - c * 8.0 + d).normalized()
};
four().or_else(|| {
let prev = at(i as isize - 1)?;
let next = at(i as isize + 1)?;
(next - prev).normalized()
})
})
.collect()
}
/// The three section curves as one closed loop of stations, the loop index each side starts at,
/// and the corner's own size (twice the mean distance of a station from their centre).
fn corner_loop(sides: &[Vec<Frame>]) -> Result<(Vec<Frame>, Vec<usize>, f64)> {
if sides.len() != 3 {
return Err(unsupported("a G2 corner needs three blends meeting at one vertex (a corner of two, or one, still ends on the ball's patch)"));
}
let mut frames: Vec<Frame> = Vec::new();
let mut starts: Vec<usize> = Vec::new();
for side in sides {
if side.len() < MIN_STATIONS || side.len() > MAX_STATIONS {
return Err(KernelError::Invalid(format!("a corner section curve needs {MIN_STATIONS} to {MAX_STATIONS} sampled stations")));
}
for f in side {
if ![f.p, f.leave, f.bend, f.bend_change].into_iter().all(Vec3::is_finite) || !(f.leave.len() > 0.5) {
return Err(KernelError::Invalid("a corner station must carry finite data and a unit direction".into()));
}
}
starts.push(frames.len());
// A side runs from one corner to the next; its last station is the next side's first, so the
// loop takes it once.
frames.extend(side.iter().take(side.len() - 1).copied());
}
let centre = frames.iter().fold(Vec3::ZERO, |a, f| a + f.p) * (1.0 / frames.len().max(1) as f64);
let size = frames.iter().map(|f| f.p.dist(centre)).fold(0.0, f64::max) * 2.0;
if !(size > 0.0 && size.is_finite()) {
return Err(unsupported("a corner whose section curves have no extent"));
}
for (i, side) in sides.iter().enumerate() {
let Some(prev) = sides.get((i + sides.len() - 1) % sides.len()).and_then(|s| s.last()) else { continue };
let Some(&here) = side.first() else { continue };
if prev.p.dist(here.p) > size * 1e-3 {
return Err(unsupported("the section curves do not meet at a corner"));
}
}
Ok((frames, starts, size))
}
/// Every station's `leave` squared up with the seam it belongs to: the part of it across the seam,
/// which is the direction the cross-section has to leave the blend in. A `leave` along the seam
/// itself leaves the cross-section no room, and is refused.
fn transverse(frames: &[Frame], tangents: &[Option<Vec3>]) -> Result<Vec<Frame>> {
let mut out = Vec::with_capacity(frames.len());
for (f, t) in frames.iter().zip(tangents) {
let across = match t {
Some(t) => f.leave - *t * f.leave.dot(*t),
// Stations that offer no direction along the seam: the leave stands as it came.
None => f.leave,
};
if across.len() < ACROSS_ENOUGH * f.leave.len() {
return Err(unsupported("a section curve meets the corner patch along its own tangent, so the patch has no direction to leave it in"));
}
let leave = across.normalized().ok_or_else(|| unsupported("a corner station has no direction across its seam"))?;
out.push(Frame { leave, ..*f });
}
Ok(out)
}
/// One station's cross-section: `ACROSS + 1` control points of a Bézier from the seam point to
/// `apex`. The first three hold the blend: its position, the direction the patch leaves it in
/// (one control point along `leave`) and its curvature vector across the seam (one more). The last
/// three hold the apex: the point, a tangent in the apex's own plane, and the curvature of a dome
/// over the corner, so the sections arrive smoothly instead of leaving a point.
fn section(f: &Frame, apex: Vec3, n_apex: Vec3, kappa: f64) -> Result<[Vec3; ACROSS + 1]> {
let k = ACROSS as f64;
let chord = f.p.dist(apex);
if !(chord.is_finite() && chord > 0.0) {
return Err(unsupported("a blend too small for the corner patch to fit"));
}
let speed = chord / 3.0;
let to_apex = (f.p - apex).normalized().ok_or_else(|| unsupported("a section curve starts on the apex"))?;
let leave_apex =
(to_apex - n_apex * to_apex.dot(n_apex)).normalized().ok_or_else(|| unsupported("a section curve reaches the apex along the apex normal"))?;
let mut c = [Vec3::ZERO; ACROSS + 1];
c[0] = f.p;
c[1] = f.p + f.leave * (speed / k);
c[2] = f.bend * (speed * speed / (k * (k - 1.0))) + c[1] * 2.0 - c[0];
c[ACROSS] = apex;
c[ACROSS - 1] = apex + leave_apex * (speed / k);
c[ACROSS - 2] = (-n_apex * kappa) * (speed * speed / (k * (k - 1.0))) + c[ACROSS - 1] * 2.0 - c[ACROSS];
if c.iter().any(|p| !p.is_finite()) {
return Err(fail("a cross-section of the corner is not finite"));
}
Ok(c)
}
/// A closed uniform cubic B-spline through `points`, once round the loop: the curve at each station
/// is `(q[i-1] + 4 q[i] + q[i+1]) / 6`, so the stations fix the net by a cyclic tridiagonal system,
/// solved here with the Sherman–Morrison form of the Thomas algorithm (a corner's loop runs to
/// hundreds of stations). The net comes back in the open form the tensor surface wants — the first
/// three control points repeated past the end — so it still closes, with one uniform span per
/// station and no free end anywhere.
fn closed_layer(points: &[Vec3]) -> Result<Vec<Vec3>> {
let m = points.len();
if m < 2 * ALONG + 2 {
return Err(fail("a corner loop needs at least eight stations"));
}
// Sub-diagonal 1, diagonal 4, super-diagonal 1, and 1 in both corners.
let (sub, diag, sup, corner) = (1.0f64, 4.0f64, 1.0f64, 1.0f64);
let gamma = -diag;
let mut bb = vec![diag; m];
if let Some(first) = bb.first_mut() {
*first = diag - gamma;
}
if let Some(last) = bb.last_mut() {
*last = diag - corner * corner / gamma;
}
// Thomas, one right-hand side at a time (three components solved together).
let thomas = |rhs: &[Vec3]| -> Result<Vec<Vec3>> {
let mut c = vec![0.0f64; m];
let mut d = rhs.to_vec();
let mut denom = bb.first().copied().unwrap_or(0.0);
if !(denom.abs() > 1e-300) {
return Err(fail("the corner's stations are degenerate"));
}
let (Some(c0), Some(d0)) = (c.first_mut(), d.first_mut()) else { return Err(fail("stations")) };
*c0 = sup / denom;
*d0 = *d0 * (1.0 / denom);
for i in 1..m {
let (prev_c, prev_d) = (c.get(i - 1).copied().unwrap_or(0.0), d.get(i - 1).copied().unwrap_or(Vec3::ZERO));
let b = bb.get(i).copied().unwrap_or(diag);
denom = b - sub * prev_c;
if !(denom.abs() > 1e-300) {
return Err(fail("the corner's stations are degenerate"));
}
if let Some(slot) = c.get_mut(i) {
*slot = if i + 1 < m { sup / denom } else { 0.0 };
}
if let Some(slot) = d.get_mut(i) {
*slot = (*slot - prev_d * sub) * (1.0 / denom);
}
}
let mut x = vec![Vec3::ZERO; m];
let last = d.last().copied().unwrap_or(Vec3::ZERO);
if let Some(slot) = x.last_mut() {
*slot = last;
}
for i in (0..m - 1).rev() {
let (Some(ci), Some(xn)) = (c.get(i).copied(), x.get(i + 1).copied()) else { continue };
let di = d.get(i).copied().unwrap_or(Vec3::ZERO);
if let Some(slot) = x.get_mut(i) {
*slot = di - xn * ci;
}
}
Ok(x)
};
let rhs: Vec<Vec3> = points.iter().map(|p| *p * 6.0).collect();
let x = thomas(&rhs)?;
// The unit right-hand side the corner terms add: gamma at the first station, the corner at the
// last.
let mut unit = vec![Vec3::ZERO; m];
if let Some(slot) = unit.first_mut() {
*slot = Vec3::new(gamma, gamma, gamma);
}
if let Some(slot) = unit.last_mut() {
*slot = Vec3::new(corner, corner, corner);
}
let z = thomas(&unit)?;
let (x0, xn) = (x.first().copied().unwrap_or(Vec3::ZERO), x.last().copied().unwrap_or(Vec3::ZERO));
let (z0, zn) = (z.first().copied().unwrap_or(Vec3::ZERO), z.last().copied().unwrap_or(Vec3::ZERO));
// z answers a scalar right-hand side, so the correction is one scalar factor per component.
let factor = |x0: f64, xn: f64, z0: f64, zn: f64| (x0 + corner * xn / gamma) / (1.0 + z0 + corner * zn / gamma);
let f = Vec3::new(factor(x0.x, xn.x, z0.x, zn.x), factor(x0.y, xn.y, z0.y, zn.y), factor(x0.z, xn.z, z0.z, zn.z));
let mut out = Vec::with_capacity(m + ALONG);
for i in 0..m {
let (xi, zi) = (x.get(i).copied().unwrap_or(Vec3::ZERO), z.get(i).copied().unwrap_or(Vec3::ZERO));
out.push(xi - Vec3::new(f.x * zi.x, f.y * zi.y, f.z * zi.z));
}
if !out.iter().all(|p| p.is_finite()) {
return Err(fail("the corner's boundary is not finite"));
}
for i in 0..ALONG {
let p = out.get(i).copied().unwrap_or(Vec3::ZERO);
out.push(p);
}
Ok(out)
}
/// The stations of the loop as one tensor control net: across (v) each station's cross-section is
/// the Bézier above, along the seam loop (u) each of its control points is a closed cubic through
/// the stations, so the face passes through every section exactly and its parameter line along the
/// loop closes without an end to swing from. The loop's station count comes back with it: station
/// `i` sits at parameter `i + ALONG`, and the loop closes at `count + ALONG`.
fn tensor(layers: &[[Vec3; ACROSS + 1]]) -> Result<(NurbsSurface, f64)> {
let m = layers.len();
let nu = m + ALONG;
let knots_u: Vec<f64> = (0..=(nu + ALONG) as i64).map(|i| i as f64).collect();
let mut ctrl: Vec<HPoint> = Vec::with_capacity(nu * (ACROSS + 1));
for j in 0..=ACROSS {
let points: Vec<Vec3> = layers.iter().map(|c| c.get(j).copied().unwrap_or(Vec3::ZERO)).collect();
let net = closed_layer(&points)?;
for p in &net {
ctrl.push(HPoint::weighted(*p, 1.0));
}
}
let mut knots_v = vec![0.0; ACROSS + 1];
knots_v.extend(std::iter::repeat_n(1.0, ACROSS + 1));
let s = NurbsSurface::from_homogeneous([ALONG, ACROSS], [knots_u, knots_v], ctrl).map_err(|e| fail(format!("the corner patch: {e}")))?;
Ok((s, m as f64))
}
/// The surface across its own domain: the patch must keep facing one way, or a blend would be met
/// by a piece of the patch folded back on itself.
fn no_fold(surface: &SurfaceEval, loop_range: (f64, f64), (v0, v1): (f64, f64)) -> Result<()> {
// The u knots run past the loop at both ends (a closed curve has no clamped ends to name), so
// the loop's own range is brought in rather than read off the knot vector.
let (u0, u1) = loop_range;
let normal = |u: f64, v: f64| -> Option<Vec3> { surface.derivs(u, v).and_then(|(_, su, sv)| su.cross(sv).normalized()) };
let mid = normal((u0 + u1) * 0.5, v0 + (v1 - v0) * 0.5).ok_or_else(|| unsupported("the corner patch has no normal at its middle"))?;
// Clear of the apex (v1), where the patch's parameter line along the seam collapses to a point
// and the normal is read from the sections approaching it.
for i in 0..=8 {
for j in 0..9 {
let u = u0 + (u1 - u0) * i as f64 / 8.0;
let v = v0 + (v1 - v0) * j as f64 / 9.0;
let n = normal(u, v).ok_or_else(|| unsupported("the corner patch folds back on itself"))?;
if n.dot(mid) <= 0.0 {
return Err(unsupported("the corner patch folds back on itself with these section curves"));
}
}
}
Ok(())
}
/// Build the patch: sections, tensor surface, and the intervals of the section curves themselves.
fn build(sides: &[Vec<Frame>], apex: Vec3) -> Result<CornerPatch> {
if !apex.is_finite() {
return Err(KernelError::Invalid("the corner apex must be a finite point".into()));
}
let (frames, starts, size) = corner_loop(sides)?;
let frames = transverse(&frames, &loop_tangents(&frames))?;
// The apex's own plane: the mean direction from the apex to the boundary it caps.
let n_apex = frames
.iter()
.fold(Vec3::ZERO, |a, f| a + (f.p - apex))
.normalized()
.ok_or_else(|| unsupported("the corner's apex sits on its own boundary curves"))?;
// How hard the patch bends into the apex: as much as its blends do across their seams, and
// never flatter than a dome of the corner's own size (a flat cap would end in a point).
let mean_bend = frames.iter().map(|f| f.bend.len()).sum::<f64>() / frames.len().max(1) as f64;
let reach = (frames.iter().map(|f| f.p.dist(apex)).sum::<f64>() / frames.len().max(1) as f64).max(size * 1e-3);
let kappa = mean_bend.max(1.0 / (2.0 * reach));
if !(kappa.is_finite() && kappa > 0.0) {
return Err(unsupported("a corner the patch cannot span"));
}
let layers: Vec<[Vec3; ACROSS + 1]> = frames.iter().map(|f| section(f, apex, n_apex, kappa)).collect::<Result<_>>()?;
let (nurbs, m) = tensor(&layers)?;
// Every weight is one (the sections are polynomial): the patch is a plain B-spline.
let [count_u, count_v] = nurbs.count();
let net = nurbs.net();
let point = |i: usize, j: usize| net.get(i, j).and_then(|h| h.project()).map(|p| [p.x, p.y, p.z]).ok_or_else(|| fail("the corner patch's net"));
let points = (0..count_u).map(|i| (0..count_v).map(|j| point(i, j)).collect::<Result<Vec<_>>>()).collect::<Result<Vec<_>>>()?;
let [knots_u, knots_v] = nurbs.knots().clone();
let surface = b::Surface::Spline { u_knots: knots_u.clone(), v_knots: knots_v.clone(), points };
let eval = SurfaceEval::new(&surface).map_err(|e| fail(format!("the corner patch: {e}")))?;
let (v0, v1) = (knots_v.first().copied().unwrap_or(0.0), knots_v.last().copied().unwrap_or(1.0));
no_fold(&eval, (ALONG as f64, m + ALONG as f64), (v0, v1))?;
// The boundary is the loop of section curves, intervals of the patch's own `v = 0` line: its
// control points are the net's first column, on the same u knots.
let seam_points = (0..count_u).map(|i| point(i, 0)).collect::<Result<Vec<_>>>()?;
let seam_curve = b::Curve::Spline { knots: knots_u, points: seam_points };
let mut cuts: Vec<f64> = starts.iter().map(|i| *i as f64 + ALONG as f64).collect();
cuts.push(m + ALONG as f64);
let seams = cuts.windows(2).filter_map(|w| Some([*w.first()?, *w.get(1)?])).collect::<Vec<_>>();
if seams.len() != sides.len() {
return Err(fail("the corner's seams"));
}
Ok(CornerPatch { surface, seam_curve, seams })
}
/// What the patch achieves along the seams of `sides`: the largest distance from a seam station to
/// the patch, the largest angle (radians) between the patch's normal there and the blend's, and the
/// largest difference between the curvature across the seam and the blend's own, relative to
/// `reference` — the curvature a rolling-ball corner shows across the same seam, which is what G2
/// removes. Measured on the built surface alone, by the same finite differences of its own samples
/// that `crate::bridge`'s tests read a cross-section's shape with, so it is evidence about the
/// geometry and not about the construction.
pub fn seam_error(surface: &b::Surface, sides: &[Vec<Frame>], reference: f64) -> Result<(f64, f64, f64)> {
if !(reference.is_finite() && reference > 0.0) {
return Err(KernelError::Invalid("a reference curvature must be finite and positive".into()));
}
// The surface's parameter box: its natural one, a periodic direction one turn.
let range = |n: b::Natural| match n {
b::Natural::Periodic => Some((0.0, std::f64::consts::TAU)),
b::Natural::Range(a, z) => (a.is_finite() && z.is_finite()).then_some((a, z)),
};
let [nu, nv] = surface.natural();
let (Some((u0, u1)), Some((v0, v1))) = (range(nu), range(nv)) else {
return Err(fail("the corner patch has no parameter domain"));
};
let surface = SurfaceEval::new(surface).map_err(|e| fail(format!("the corner patch: {e}")))?;
let at = |u: f64, v: f64| surface.point([u, v]).ok_or_else(|| fail("the corner patch has no point there"));
// Samples across the section for a five-point difference, as `bridge`'s tests use.
const INWARDS: usize = 5;
let hv = (v1 - v0).abs() * 1e-3;
let hu = (u1 - u0).abs() * 1e-3;
if !(hv > 0.0 && hu > 0.0) {
return Err(fail("the corner patch has no extent to measure"));
}
let (mut gap, mut angle, mut curvature) = (0.0f64, 0.0f64, 0.0f64);
for side in sides {
if side.len() < 5 {
return Err(KernelError::Invalid("a seam needs at least five stations to measure its direction".into()));
}
let mut hint: Option<[f64; 2]> = None;
for (i, f) in side.iter().enumerate() {
let Some(t) = seam_tangent(side, i) else { continue };
let found = surface
.project(f.p, hint, f64::INFINITY)
.or_else(|| surface.project(f.p, None, f64::INFINITY))
.ok_or_else(|| fail("a seam station is not on the corner patch"))?;
hint = Some(found);
let [u, v] = found;
gap = gap.max(at(u, v)?.dist(f.p));
// The cross-section, sampled inwards from the seam (or from the other end of the patch,
// with the odd differences turned round).
let out = v <= (v0 + v1) * 0.5;
let step = if out { hv } else { -hv };
let q: Vec<Vec3> = (0..INWARDS).map(|k| at(u, v + step * k as f64)).collect::<Result<_>>()?;
let (Some(&f0), Some(&f1), Some(&f2), Some(&f3), Some(&f4)) = (q.first(), q.get(1), q.get(2), q.get(3), q.get(4)) else {
return Err(fail("the corner patch has no room for a sample"));
};
let sign = if out { 1.0 } else { -1.0 };
let d1 = (f0 * -25.0 + f1 * 48.0 + f2 * -36.0 + f3 * 16.0 + f4 * -3.0) * (sign / (12.0 * hv));
let d2 = (f0 * 35.0 + f1 * -104.0 + f2 * 114.0 + f3 * -56.0 + f4 * 11.0) * (1.0 / (12.0 * hv * hv));
// Along the seam, to stand the tangent plane up.
let du = if u - hu >= u0 && u + hu <= u1 {
(at(u + hu, v)? - at(u - hu, v)?) * (0.5 / hu)
} else {
let (a, b) = if u + hu <= u1 { (u, u + hu) } else { (u, u - hu) };
(at(b, v)? - at(a, v)?) * (1.0 / hu)
};
let n = du.cross(d1).normalized().ok_or_else(|| fail("the corner patch has no normal at a seam"))?;
let blend_n = t.cross(f.leave).normalized().ok_or_else(|| fail("a blend has no normal at a seam"))?;
angle = angle.max(n.dot(blend_n).abs().clamp(0.0, 1.0).acos());
// The curvature vector of the cross-seam curve: what the blend's own surface does along
// the direction the patch leaves it in.
let speed2 = d1.len2().max(1e-300);
let bend = (d2 - d1 * (d2.dot(d1) / speed2)) * (1.0 / speed2);
curvature = curvature.max((bend - f.bend).len() / reference.max(f.bend.len()).max(1e-30));
}
}
Ok((gap, angle, curvature))
crate::fill::matched_patch(sheet, hole, CURVATURE_SHARE / radius)
}
#[cfg(test)]
#[path = "corner_g2_tests.rs"]
mod tests;
/// The G2 corner through the fillet path: a box's three edges at one corner rounded, rolling ball
/// against G2 patch, measured on the built bodies with the continuity analysis.
#[cfg(test)]
#[allow(clippy::unwrap_used)]
mod fillet_corner_tests {
use super::*;
use crate::{ContinuityTolerances, SeamReport, measure, seam_diagnostics};
/// A 20 mm box with its three edges at one corner rounded by 4: the rolling ball fills that
/// corner with a piece of its sphere, the G2 route with a patch on the blends' own sections.
fn cornered(g2: bool) -> crate::Body {
let b = crate::box_solid(Vec3::ZERO, Vec3::new(20.0, 20.0, 20.0)).unwrap();
let edges = [Vec3::new(10.0, 0.0, 0.0), Vec3::new(0.0, 10.0, 0.0), Vec3::new(0.0, 0.0, 10.0)];
if g2 { crate::fillet_g2_corners(&b, &edges, 4.0).unwrap() } else { crate::fillet(&b, &edges, 4.0).unwrap() }
}
#[test]
fn the_ball_corner_is_measurably_g1_and_a_g2_corner_is_refused_until_it_sews() {
// The rolling ball: tangent to every blend, keeping its own 1/r across each corner seam.
let body = cornered(false);
let volume = measure(&body).unwrap().volume;
assert!(volume > 7000.0 && volume < 8000.0, "the box is still a solid: {volume}");
let report = seam_diagnostics(&body, &[], 33, ContinuityTolerances { gap_mm: 1e-3, angle_deg: 1.0, curvature_per_mm: 1e-6 }).unwrap();
// Every seam of a filleted box is tangent; the corner's three are the short ones (the
// blend-to-face seams run the edge's own 12 mm and more).
let tangent: Vec<&SeamReport> = report.seams.iter().filter(|s| s.normal_angle_deg.value < 0.01).collect();
let corner: Vec<&&SeamReport> = tangent.iter().filter(|s| s.length < 8.0).collect();
assert_eq!(
corner.len(),
3,
"three corner seams among {} tangent ones: {:?}",
tangent.len(),
tangent.iter().map(|s| s.length).collect::<Vec<_>>()
);
let worst = corner.iter().map(|s| s.curvature_jump_per_mm.value).fold(0.0f64, f64::max);
eprintln!("rolling ball: corner seams jump {worst:.6} per mm (volume {volume:.3})");
assert!((worst - 0.25).abs() < 1e-3, "the ball keeps its own 1/r (r = 4) across the seam: {worst}");
// The G2 corner patch is built on the blends' own sections, but sewing it into the shell
// is not done yet: the fillet with `corner: g2` refuses with a clear message rather than
// leaving the shell open, and leaves the ball the corner.
let b = crate::box_solid(Vec3::ZERO, Vec3::new(20.0, 20.0, 20.0)).unwrap();
let edges = [Vec3::new(10.0, 0.0, 0.0), Vec3::new(0.0, 10.0, 0.0), Vec3::new(0.0, 0.0, 10.0)];
let err = crate::fillet_g2_corners(&b, &edges, 4.0).unwrap_err().to_string();
eprintln!("G2 corner: {err}");
assert!(err.contains("not supported yet:") && err.contains("G2 corner"), "{err}");
}
}
+157 -140
View File
@@ -1,149 +1,166 @@
//! What the G2 corner patch has to prove on a corner whose answer is known analytically: the box.
use solidcraft_geom::Vec3;
use solidcraft_geom::brep as b;
use super::*;
use std::f64::consts::FRAC_PI_2;
use crate::*;
/// One seam of a box corner: a quarter circle about `centre` from direction `u` to direction `w`,
/// with `leave` the way the patch leaves the blend's surface. A fillet is flat along its own axis,
/// so its curvature across the seam is nothing.
fn seam(centre: Vec3, u: Vec3, w: Vec3, leave: Vec3, r: f64, n: usize) -> Vec<Frame> {
(0..n)
.map(|i| {
let theta = FRAC_PI_2 * i as f64 / (n - 1) as f64;
Frame { p: centre + (u * theta.cos() + w * theta.sin()) * r, leave, bend: Vec3::ZERO, bend_change: Vec3::ZERO }
})
.collect()
/// The three r5 rounds at the corner (60, 40, 20) of a 60 x 40 x 20 box (ladder 07/009 and
/// 07/015).
fn box_corner() -> (Body, [Vec3; 3]) {
let b = box_solid(Vec3::new(0.0, 0.0, 0.0), Vec3::new(60.0, 40.0, 20.0)).unwrap();
(b, [Vec3::new(30.0, 40.0, 20.0), Vec3::new(60.0, 20.0, 20.0), Vec3::new(60.0, 40.0, 10.0)])
}
/// A box corner at the origin (material x, y, z ≥ 0) with three fillets of radius `r`: the three
/// section curves where the fillets stop, meeting at (0,r,r), (r,0,r) and (r,r,0), and the middle of
/// the ball's octant a rolling ball would leave there.
fn box_corner(r: f64, n: usize) -> (Vec<Vec<Frame>>, Vec3, Vec3) {
let c = Vec3::new(r, r, r);
let (a, b, d) = (Vec3::new(0.0, r, r), Vec3::new(r, 0.0, r), Vec3::new(r, r, 0.0));
let dir = |p: Vec3| (p - c).normalized().unwrap_or(Vec3::X);
let sides = vec![
// A → B on the z edge's cylinder (flat along z), B → C on the x edge's, C → A on the y
// edge's; each patch leaves the blend towards the material the corner caps.
seam(c, dir(a), dir(b), Vec3::new(0.0, 0.0, -1.0), r, n),
seam(c, dir(b), dir(d), Vec3::new(-1.0, 0.0, 0.0), r, n),
seam(c, dir(d), dir(a), Vec3::new(0.0, -1.0, 0.0), r, n),
];
let apex = c - Vec3::new(1.0, 1.0, 1.0).normalized().unwrap_or(Vec3::Z) * r;
(sides, apex, c)
}
/// The rolling-ball corner the kernel builds today: the sphere through the same three seams, made
/// as a sphere record about the ball's centre, its axis along the corner's diagonal.
fn rolling_ball(r: f64, c: Vec3) -> b::Surface {
let axis = Vec3::new(1.0, 1.0, 1.0).normalized().unwrap_or(Vec3::Z);
let side = axis.any_perp().normalized().unwrap();
let y = axis.cross(side);
b::Surface::Sphere { frame: b::Frame { origin: [c.x, c.y, c.z], x_axis: [side.x, side.y, side.z], y_axis: [y.x, y.y, y.z] }, radius: r }
}
/// The stations clear of the corners. Two blends meet at a corner with their seams 90° apart and
/// ask for different curvature at that one point, which no single patch can hold, and the corner's
/// own spline rounds the seam there; both effects fade within a few stations, so the measurement
/// starts clear of them.
fn interior(sides: &[Vec<Frame>]) -> Vec<Vec<Frame>> {
sides.iter().map(|s| s.iter().skip(4).take(s.len().saturating_sub(8)).copied().collect()).collect()
}
#[test]
fn a_box_corner_is_curvature_continuous_with_every_blend() {
let r = 2.0;
let (sides, apex, centre) = box_corner(r, 25);
let surface = corner_g2_patch(&sides, apex).unwrap().surface;
let inner = interior(&sides);
let (gap, angle, curvature) = seam_error(&surface, &inner, 1.0 / r).unwrap();
println!("G2 corner: gap {gap:e}, normal angle {angle:e} rad, curvature {curvature:e} (relative to 1/r)");
assert!(gap <= r * 1e-5, "the patch is {gap} from the seam");
// The stated tolerances: the normal angle within 0.01° and the curvature across the seam within
// 1% of the blend's own (here the fillet is flat across the seam, so that is 1% of 1/r).
assert!(angle < 0.01f64.to_radians(), "normal angle {angle} rad");
assert!(curvature < 1e-2, "curvature across the seam is {curvature} of the blend's");
// The same measurement on the rolling-ball corner the kernel builds today: G1 (angle nought)
// but the ball's own 1/r across each seam, which is the whole of the blend's curvature.
let ball = rolling_ball(r, centre);
let (ball_gap, ball_angle, ball_curvature) = seam_error(&ball, &inner, 1.0 / r).unwrap();
println!("rolling ball: gap {ball_gap:e}, normal angle {ball_angle:e} rad, curvature {ball_curvature:e}");
assert!(ball_angle < 0.01f64.to_radians(), "the ball is tangent to each blend");
assert!(ball_curvature > 0.5, "the ball keeps its own curvature across the seam: {ball_curvature}");
assert!(curvature * 20.0 < ball_curvature, "the G2 patch has to be far closer than the ball");
}
#[test]
fn a_corner_meets_three_blends_of_different_shape_apart() {
// Blends that differ along their seams: a fillet of radius 2 next to rounder walls, each asking
// for its own curvature across the seam, at its own direction.
let r = 2.0;
let (mut sides, apex, _) = box_corner(r, 33);
let bends = [0.25, 0.4, 0.6];
for (i, bend) in bends.iter().enumerate() {
let Some(side) = sides.get_mut(i) else { continue };
for f in side.iter_mut() {
// Across the seam, perpendicular to both the leave direction and the seam.
f.bend = f.leave.any_perp() * *bend;
}
}
let surface = corner_g2_patch(&sides, apex).unwrap().surface;
let inner = interior(&sides);
let (gap, angle, curvature) = seam_error(&surface, &inner, bends.iter().copied().fold(0.0, f64::max)).unwrap();
println!("three different blends: gap {gap:e}, angle {angle:e} rad, curvature {curvature:e}");
assert!(gap <= r * 1e-6, "{gap}");
assert!(angle < 0.01f64.to_radians(), "{angle}");
assert!(curvature < 1e-2, "each blend's own curvature across the seam: {curvature}");
}
#[test]
fn a_g2_corner_refuses_what_it_cannot_hold() {
let r = 2.0;
let (sides, apex, _) = box_corner(r, 9);
assert!(corner_g2_patch(&sides, apex).is_ok());
// Two blends (the ball's patch still fills that corner), or one.
assert!(corner_g2_patch(&sides[..2], apex).is_err());
assert!(corner_g2_patch(&sides[..1], apex).is_err());
assert!(corner_g2_patch(&[], apex).is_err());
// Too few stations to interpolate a cubic through, and too many to work with.
let thin: Vec<Vec<Frame>> = sides.iter().map(|s| s[..3].to_vec()).collect();
assert!(corner_g2_patch(&thin, apex).is_err());
let fat: Vec<Vec<Frame>> = sides.iter().map(|s| s[..1].iter().cycle().take(MAX_STATIONS + 1).copied().collect()).collect();
assert!(corner_g2_patch(&fat, apex).is_err());
// A vertex the construction cannot hold: the section curves do not meet at their corners.
let mut apart = sides.clone();
if let Some(first) = apart.first_mut().and_then(|s| s.first_mut()) {
first.p += Vec3::new(0.5, 0.0, 0.0);
}
assert!(corner_g2_patch(&apart, apex).is_err());
// A section curve that runs along the seam, so the patch has no direction to leave it in.
let along: Vec<Vec<Frame>> = sides
/// The worst curvature jump across the seams of the faces that are none of the box's planes or
/// rounds (the corner's: a sphere, or a patch), and the same clear of the hole's corners: a station counts there when
/// it is at least `share` of its seam's side (its edge and the edge that continues it smoothly:
/// SMLib splits each side in two) from the end of the side at a corner, a vertex where the seam
/// turns.
fn corner_seams(body: &Body, share: f64) -> (f64, f64, usize) {
let rec = brep_record(body).unwrap();
let corner_faces: Vec<usize> = rec
.faces
.iter()
.map(|s| {
s.windows(2)
.map(|w| {
let d = w.first().zip(w.get(1)).map(|(a, b)| (b.p - a.p).normalized().unwrap_or(Vec3::X)).unwrap_or(Vec3::X);
Frame { leave: d, ..w.first().copied().unwrap_or(s[0]) }
})
.chain(s.last().copied())
.collect()
.enumerate()
.filter(|(_, f)| !matches!(f.surface, b::Surface::Plane { .. } | b::Surface::PlaneThrough { .. } | b::Surface::Cylinder { .. }))
.map(|(i, _)| i)
.collect();
let edges: Vec<usize> = (0..rec.edges.len()).filter(|e| rec.edges[*e].uses.iter().any(|u| corner_faces.contains(&(u.face as usize)))).collect();
let report = seam_diagnostics(body, &edges, 65, ContinuityTolerances { gap_mm: 1e-3, angle_deg: 1.0, curvature_per_mm: 1e-6 }).unwrap();
// The corner's vertices where its boundary turns: an end of two of its seams whose
// directions there differ.
let ends: Vec<(Vec3, Vec3, f64)> = report
.seams
.iter()
.flat_map(|s| {
let (p0, p1) = (s.points[0], *s.points.last().unwrap());
let d0 = (s.points[1] - p0).normalized().unwrap();
let d1 = (s.points[s.points.len() - 2] - p1).normalized().unwrap();
[(p0, d0, s.length), (p1, d1, s.length)]
})
.collect();
assert!(corner_g2_patch(&along, apex).is_err());
// Hostile numbers never panic.
let mut nan = sides.clone();
if let Some(f) = nan.first_mut().and_then(|s| s.get_mut(1)) {
f.p = Vec3::new(f64::NAN, 0.0, 0.0);
let turns =
|p: Vec3, d: Vec3| ends.iter().any(|(q, e, _)| q.dist(p) < 1e-6 && (d - *e).len() > 1e-6 && d.dot(*e) > -0.9998 && d.dot(*e) < 0.9998);
let (mut all, mut clear, mut n) = (0.0f64, 0.0f64, 0usize);
for s in &report.seams {
let (p0, p1) = (s.points[0], *s.points.last().unwrap());
let d0 = (s.points[1] - p0).normalized().unwrap();
let d1 = (s.points[s.points.len() - 2] - p1).normalized().unwrap();
let corner_end = if turns(p0, d0) {
Some(0.0)
} else if turns(p1, d1) {
Some(1.0)
} else {
None
};
for st in &s.stations {
all = all.max(st.curvature_jump_per_mm);
// Half a side per edge: `share` of the side is twice that of the edge.
if corner_end.is_none_or(|c| (st.edge_fraction - c).abs() >= 2.0 * share) {
clear = clear.max(st.curvature_jump_per_mm);
n += 1;
}
}
}
assert!(corner_g2_patch(&nan, apex).is_err());
assert!(corner_g2_patch(&sides, Vec3::new(f64::INFINITY, 0.0, 0.0)).is_err());
// A corner with no room: every station on the apex.
let on_apex: Vec<Vec<Frame>> = sides.iter().map(|s| s.iter().map(|f| Frame { p: apex, ..*f }).collect()).collect();
assert!(corner_g2_patch(&on_apex, apex).is_err());
(all, clear, n)
}
/// The G2 corner of 07/015, set back the default 7.5 (1.5 r): one closed valid solid of the box's
/// six faces, the three rounds and the patch (10 faces, as the sphere corner's 07/009), within
/// the box, whose volume lies between the corner cube's hull of the body's points on its faces
/// (the body is convex) and the rounds' prisms, as the setback's own does; and across every
/// seam of the patch the curvature jump, a fifth per mm on the rolling ball's, is under
/// `CURVATURE_SHARE / r` clear of the hole's corners.
#[test]
fn the_box_corner_g2_patch_is_a_closed_solid_curvature_continuous_to_its_neighbours() {
let (body, edges) = box_corner();
let ball = fillet(&body, &edges, 5.0).unwrap();
let g2 = fillet_g2_corners(&body, &edges, 5.0, &[]).unwrap();
assert!(g2.validity().is_empty(), "{:?}", g2.validity());
let v = crate::validity::validity(&g2);
assert!(v.ok(), "{:?}", v);
let m = measure(&g2).unwrap();
let mb = measure(&ball).unwrap();
eprintln!("G2 corner: faces {} edges {} vertices {}, volume {:.6} (ball {:.6})", m.faces, m.edges, m.vertices, m.volume, mb.volume);
assert_eq!(m.faces, 10);
assert!((m.bbox.min - Vec3::ZERO).len() < 1e-6 && (m.bbox.max - Vec3::new(60.0, 40.0, 20.0)).len() < 1e-6, "{:?}", m.bbox);
let (all, clear, n) = corner_seams(&g2, 0.15);
let (ball_all, ball_clear, _) = corner_seams(&ball, 0.15);
eprintln!(
"G2 corner seams: worst jump {all:.4e} per mm, {clear:.4e} clear of the corners ({n} stations); the ball's {ball_all:.4e} / {ball_clear:.4e}"
);
assert!(clear <= CURVATURE_SHARE / 5.0, "{clear}");
}
/// Set back 8, as SMLib's own setback test: outside the corner cube [52, 60] x [32, 40] x
/// [12, 20] the body is the sphere corner's, and inside it the convex body holds the hull of its
/// points on the cube's faces (299.0201, by a hull of them sampled) and lies within the three
/// rounds' prisms (`fillet_setback_tests`).
#[test]
fn set_back_8_the_g2_corner_lies_between_the_hull_and_the_prisms() {
let pi = std::f64::consts::PI;
let (body, edges) = box_corner();
let sphere = measure(&fillet(&body, &edges, 5.0).unwrap()).unwrap().volume;
let g2 = fillet_g2_corners(&body, &edges, 5.0, &[(Vec3::new(60.0, 40.0, 20.0), 8.0)]).unwrap();
let v = measure(&g2).unwrap().volume;
let outside = sphere - (512.0 - 225.0 * (1.0 - pi / 4.0) - 125.0 * (1.0 - pi / 6.0));
let prisms = 512.0 - 225.0 * (1.0 - pi / 4.0) - 125.0 * (2f64.sqrt() - 1.0);
eprintln!("G2 corner set back 8: volume {v:.6} in [{:.6}, {:.6}]", outside + 299.0201, outside + prisms);
assert!(v > outside + 299.0201 && v < outside + prisms);
let (_, clear, _) = corner_seams(&g2, 0.15);
assert!(clear <= CURVATURE_SHARE / 5.0, "{clear}");
}
/// Long setbacks, where SMLib's own patch bulges past the faces (beyond about 2 r, by 0.25 at
/// 3 r): the G2 patch is fitted to the faces' own planes along their curves and stays inside the
/// box but by under a thousandth (where it passes from a round's curvature to a plane's), G2
/// clear of the hole's corners.
#[test]
fn a_long_setback_keeps_the_g2_patch_inside_the_box() {
let (body, edges) = box_corner();
for s in [11.0, 12.0, 13.0, 14.0, 15.0] {
let g2 = fillet_g2_corners(&body, &edges, 5.0, &[(Vec3::new(60.0, 40.0, 20.0), s)]).unwrap();
let m = g2.mesh(&MeshTolerance::chord(1e-3)).unwrap();
let over = m.positions.iter().map(|p| (p.x - 60.0).max(p.y - 40.0).max(p.z - 20.0)).fold(f64::NEG_INFINITY, f64::max);
let (_, clear, _) = corner_seams(&g2, 0.15);
eprintln!("set back {s}: past the box by {over:.2e}, curvature jump clear of the corners {clear:.3e}");
assert!(over < 1e-3 && clear <= CURVATURE_SHARE / 5.0, "{s}: {over} {clear}");
}
}
/// Every edge of a box rounded: eight corners, each set back and filled, one closed solid.
#[test]
fn a_box_with_every_edge_rounded_has_eight_g2_corners() {
let body = box_solid(Vec3::ZERO, Vec3::new(30.0, 20.0, 10.0)).unwrap();
let mut edges = Vec::new();
for (x, y) in [(0.0, 0.0), (30.0, 0.0), (0.0, 20.0), (30.0, 20.0)] {
edges.push(Vec3::new(x, y, 5.0));
}
for (y, z) in [(0.0, 0.0), (20.0, 0.0), (0.0, 10.0), (20.0, 10.0)] {
edges.push(Vec3::new(15.0, y, z));
}
for (x, z) in [(0.0, 0.0), (30.0, 0.0), (0.0, 10.0), (30.0, 10.0)] {
edges.push(Vec3::new(x, 10.0, z));
}
let g2 = fillet_g2_corners(&body, &edges, 2.0, &[]).unwrap();
assert!(g2.validity().is_empty(), "{:?}", g2.validity());
let v = crate::validity::validity(&g2);
let m = measure(&g2).unwrap();
eprintln!("eight G2 corners: {} faces, volume {:.4}, valid {}", m.faces, m.volume, v.ok());
assert!(v.ok(), "{v:?}");
assert_eq!(m.faces, 6 + 12 + 8);
}
/// What the G2 corner cannot do is refused: a setback naming no corner, one that does not reach
/// past where the rounds meet; two rounded edges at a corner are no G2 corner, the plain fillet.
#[test]
fn what_the_g2_corner_cannot_build_is_refused() {
let (body, edges) = box_corner();
let off = fillet_g2_corners(&body, &edges, 5.0, &[(Vec3::new(30.0, 20.0, 10.0), 8.0)]).unwrap_err().to_string();
assert!(off.contains("names no G2 corner"), "{off}");
let short = fillet_g2_corners(&body, &edges, 5.0, &[(Vec3::new(60.0, 40.0, 20.0), 4.0)]);
eprintln!("set back 4 at r5: {short:?}");
assert!(short.is_err());
let two = fillet_g2_corners(&body, &edges[..2], 5.0, &[]).unwrap();
let plain = fillet(&body, &edges[..2], 5.0).unwrap();
assert_eq!(measure(&two).unwrap().faces, measure(&plain).unwrap().faces);
}
+64 -1
View File
@@ -15,8 +15,9 @@ use solidcraft_geom::Vec3;
use solidcraft_geom::brep::eval::SurfaceEval;
use solidcraft_geom::brep::{self as b, Curve2, eval::CurveEval};
use solidcraft_geom::bspline::{HPoint, SplineError};
use solidcraft_geom::fill_general::{Across, MAX_FILL_POINTS};
use solidcraft_geom::fill_general::{Across, MAX_FILL_POINTS, Station};
use solidcraft_geom::fill_iso::{IsoEdge, IsoFill, IsoSpec};
use solidcraft_geom::fill_matched::{MatchedFill, MatchedSpec};
use crate::body::Body;
use crate::bridge::Continuity;
@@ -186,6 +187,60 @@ fn general_fill(
Ok((made, report))
}
/// How far (degrees) the loop may turn at a vertex for it not to be a corner of the loop.
const CORNER_TURN_DEG: f64 = 1.0;
/// The one-face sheet on the loop `chain` of `body`'s own edges (each a free edge of the face
/// beside it), a height field over the loop's plane that meets every face beside the loop with
/// G2 (`solidcraft_geom::fill_matched`), within [`face_tolerance`] of the loop, its normal within
/// [`PATCH_ANGLE_TOL_DEG`] and its normal curvature across the loop within `curvature_tol` of
/// each face's, clear of the loop's corners (a vertex where the loop turns, where two faces
/// meet). Each edge's pcurve is its exact projection onto the plane. Its history is over
/// `[body]`: the face is new, each edge the body's own.
pub(crate) fn matched_patch(body: &Body, chain: &[(usize, bool)], curvature_tol: f64) -> Result<(Body, MatchedFill)> {
body.require_brep("patch")?;
if chain.len() < 2 {
return Err(KernelError::Invalid("a patch needs a loop of two or more edges".into()));
}
let record = crate::native::brep_record(body)?;
let size = body.size().max(solidcraft_geom::LINEAR_TOL);
let mut sides = Vec::with_capacity(chain.len());
let mut beside = Vec::with_capacity(chain.len());
for &(e, forward) in chain {
sides.push(Side::new(&record, e, forward)?);
beside.push(Beside::new(&record, e, forward, true, size)?);
}
// A vertex of the loop is a corner where the loop turns there: the edge arriving and the one
// leaving go different ways.
let n = chain.len();
let mut corner = vec![false; n];
for (k, slot) in corner.iter_mut().enumerate() {
let prev = (k + n - 1) % n;
let (Some(a), Some(b)) = (sides.get(prev), sides.get(k)) else { continue };
let (Some(arrive), Some(leave)) = (a.direction(1.0), b.direction(0.0)) else { continue };
*slot = arrive.dot(leave) < CORNER_TURN_DEG.to_radians().cos();
}
let mut stations = |k: usize, ts: &[f64]| -> std::result::Result<Vec<Station>, String> {
let (Some(side), Some(face)) = (sides.get_mut(k), beside.get(k)) else { return Err("no such edge".into()) };
let edge = record.edges.get(face.edge).ok_or("edge reference")?;
let [a, z] = edge.interval;
ts.iter()
.map(|&t| {
let p = side.point(t).ok_or("an edge has no point there")?;
let across = face.across(&record, a + (z - a) * t).map_err(|e| e.to_string())?.ok_or("an edge has no tangent there")?;
let at_corner = t == 0.0 && corner.get(k).copied().unwrap_or(false);
Ok(Station { p, up: across.normal, across: Some(across), corner: at_corner, clash: false })
})
.collect()
};
let spec = MatchedSpec { tol: face_tolerance(body), angle_tol: PATCH_ANGLE_TOL_DEG.to_radians(), curvature_tol };
let fill = solidcraft_geom::fill_matched::fill_matched(n, &mut stations, &spec).map_err(spline_error)?;
let projection = Projection { origin: fill.fill.origin, axes: fill.fill.axes };
let surface = crate::surfaces::spline_surface(&fill.fill.surface);
let made = patch_body(body, &record, chain, surface, fill.fill.counter_clockwise, &Pcurves::Projected(&projection), face_tolerance(body))?;
Ok((made, fill))
}
/// A loop edge as the patch takes it: a non-rational B-spline on the edge's interval, running
/// the loop's way, its ends on the edge's vertices. Conics and traced curves have no such form
/// in their own parameter (an arc's angle is no rational function of a NURBS parameter), and a
@@ -483,6 +538,14 @@ impl Side {
let s = t.min(1.0);
self.curve.point(if self.forward { a + (z - a) * s } else { z - (z - a) * s })
}
/// The unit direction the loop runs the edge in at `t`.
fn direction(&self, t: f64) -> Option<Vec3> {
let [a, z] = self.interval;
let s = t.clamp(0.0, 1.0);
let (_, d) = self.curve.derivs(if self.forward { a + (z - a) * s } else { z - (z - a) * s })?;
(if self.forward { d } else { d * -1.0 }).normalized()
}
}
/// The plane of the loop when its edges certainly lie in one within `tol`, `None` otherwise. The
+3 -2
View File
@@ -166,9 +166,10 @@ pub use build::{
loft_with_rails, path_polyline, planar_face_record, revolve, revolve_whole, sphere, sweep, sweep_with, torus,
};
pub use continuity::{
ContinuityReport, ContinuityTolerances, EdgePick, MAX_SAMPLES, MIN_SAMPLES, SeamReport, SeamVerdict, SeamWorst, pick_edges, seam_diagnostics,
ContinuityReport, ContinuityTolerances, EdgePick, MAX_SAMPLES, MIN_SAMPLES, SeamReport, SeamStation, SeamVerdict, SeamWorst, pick_edges,
seam_diagnostics,
};
pub use corner_g2::{CornerPatch, corner_g2_patch, seam_error};
pub use corner_g2::{CURVATURE_SHARE, DEFAULT_SETBACK_FACTOR};
pub use delete_face::delete_faces;
pub use distance::{BodyDistance, DistanceShape, body_distance};
pub use draft_corners::{CornerFill, CornerOption, CornerRefusal, CornerReport, DraftCorner, draft_corners};