masi/solvecraft
1
0
Fork 0
mirror of https://github.com/bherbruck/solidcraft.git synced 2026-11-03 02:25:55 +00:00
Code Issues Packages Projects Releases 1 Wiki Activity

geom: a face's box reaches its boundary's extremes, not its samples'

Browse Source 
form::extent took the box of 32 points per edge. A face's NURBS form is built
over that box, so where the boundary peaks between two samples the form
stopped short of its edge: the re-imported rod's bore ended at v = 9.99977
while its intersection with the rod reaches 10, and SMLib found the edge
1.9e-4 off the bore (SmBrep #4, SmEdge #13). Each sample higher or lower
than its neighbours now brackets the edge's extreme there, found by golden
section.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_0113gDVmStT5tCGKW9Fxgc5K
This commit is contained in:
jyviko
2026-10-25 19:51:02 -04:00
parent 098fac037e
commit 903e371fa7
2 changed files with 81 additions and 5 deletions
Download Patch File Download Diff File Expand all files Collapse all files
crates/geom/src/brep/form.rs
+67 -5
View File
@@ -317,6 +317,32 @@ pub fn edge<'a>(
/// Points per edge along a loop taken for a face's extent.
const EXTENT_POINTS: usize = 32;
/// The largest value of `f` on `[a, b]`, where it rises to one peak and falls (golden section, to
/// rounding; the ends count).
fn golden_max(f: &dyn Fn(f64) -> f64, a: f64, b: f64) -> f64 {
const R: f64 = 0.618_033_988_749_894_8;
let (mut a, mut b) = (a, b);
let best = f(a).max(f(b));
let (mut c, mut d) = (b - R * (b - a), a + R * (b - a));
let (mut fc, mut fd) = (f(c), f(d));
while b - a > 1e-15 * (1.0 + a.abs().max(b.abs())) {
if fc >= fd {
b = d;
d = c;
fd = fc;
c = b - R * (b - a);
fc = f(c);
} else {
a = c;
c = d;
fc = fd;
d = a + R * (b - a);
fd = f(d);
}
}
best.max(fc).max(fd)
}
/// The parameter box of a face on `label` over its boundary: each loop's points (`loops`, in
/// order along each loop) projected onto the label, a periodic parameter unwrapped from one point
/// to the next, and each loop after the first a whole number of periods along where that keeps
@@ -334,17 +360,23 @@ pub fn extent(label: &Surface, loops: &[Vec<(&Curve, [f64; 2], bool)>]) -> Optio
let mut prev: Option<[f64; 2]> = None;
for (c, [t0, t1], forward) in l {
let curve = super::eval::CurveEval::new(c).ok()?;
for k in 0..=EXTENT_POINTS {
let f = k as f64 / EXTENT_POINTS as f64;
let t = if *forward { t0 * (1.0 - f) + t1 * f } else { t1 * (1.0 - f) + t0 * f };
let mut q = eval.project(curve.point(t)?, prev, f64::INFINITY)?;
if let Some(p) = prev {
// The point at `t` on the label, a periodic parameter unwrapped next to `near`.
let at = |t: f64, near: Option<[f64; 2]>| -> Option<[f64; 2]> {
let mut q = eval.project(curve.point(t)?, near, f64::INFINITY)?;
if let Some(p) = near {
for d in 0..2 {
if let Some(Some(per)) = period.get(d) {
q[d] += ((p[d] - q[d]) / per).round() * per;
}
}
}
Some(q)
};
let mut samples: Vec<(f64, [f64; 2], [bool; 2])> = Vec::with_capacity(EXTENT_POINTS + 1);
for k in 0..=EXTENT_POINTS {
let f = k as f64 / EXTENT_POINTS as f64;
let t = if *forward { t0 * (1.0 - f) + t1 * f } else { t1 * (1.0 - f) + t0 * f };
let q = at(t, prev)?;
// At a pole the collapsed parameter is anything: it bounds nothing.
let pole = eval.pole(q);
for d in 0..2 {
@@ -356,6 +388,36 @@ pub fn extent(label: &Surface, loops: &[Vec<(&Curve, [f64; 2], bool)>]) -> Optio
if !pole.iter().any(|p| *p) {
prev = Some(q);
}
samples.push((t, q, pole));
}
// A sample is a point of the boundary, not its extreme: each sample higher (or lower)
// than its neighbours brackets an extreme of the edge, which is found there.
for k in 0..samples.len() {
let (_, q, pole) = samples[k];
let (a, b) = (samples[k.saturating_sub(1)].0, samples[(k + 1).min(samples.len() - 1)].0);
for d in 0..2 {
// Next to a pole the parameter it collapses is no function of the point: a
// meridian's angle there is rounding, not an extreme.
let near_pole = |j: Option<usize>| j.and_then(|j| samples.get(j)).is_some_and(|s| s.2[d]);
if pole[d] || near_pole(k.checked_sub(1)) || near_pole(Some(k + 1)) {
continue;
}
let side = |j: usize| samples.get(j).map(|s| s.1[d]);
let (before, after) = (if k > 0 { side(k - 1) } else { None }, side(k + 1));
for sign in [1.0, -1.0] {
// Not lower than either neighbour and higher than one (a level run, a
// seam's constant angle, has no extreme between its samples).
let higher = |x: Option<f64>| x.is_none_or(|x| sign * q[d] >= sign * x);
let above = |x: Option<f64>| x.is_some_and(|x| sign * q[d] > sign * x);
if !(higher(before) && higher(after) && (above(before) || above(after))) {
continue;
}
let f = |t: f64| at(t, Some(q)).map_or(f64::NEG_INFINITY, |p| sign * p[d]);
let x = sign * golden_max(&f, a.min(b), a.max(b));
l_lo[d] = l_lo[d].min(x);
l_hi[d] = l_hi[d].max(x);
}
}
}
}
for d in 0..2 {
crates/geom/src/brep/tests.rs
+14
View File
@@ -829,3 +829,17 @@ fn a_hole_beside_the_seam_lies_in_its_face_turn() {
assert!(hole_u > u[0] && hole_u + (a1 - a0) < u[1], "the hole in the turn: {u:?}");
assert!(v[0].abs() < 1e-9 && (v[1] - h).abs() < 1e-9, "{v:?}");
}
/// A face's box reaches its boundary's extremes, not its samples' nearest ones: an ellipse cut
/// from a cylinder of radius 20 by a plane at 45 degrees rises 20 above its centre, at an angle
/// between the points the box is sampled at. Short of it, the face's NURBS form stopped below the
/// edge (#4598: the bore of a re-imported rod ended 2.3e-4 short of where it meets the rod).
#[test]
fn a_face_box_reaches_its_boundary_extremes() {
let r = 20.0f64;
let cylinder = Surface::Cylinder { frame: Frame { origin: [0.; 3], x_axis: [1., 0., 0.], y_axis: [0., 1., 0.] }, radius: r };
let c = std::f64::consts::FRAC_1_SQRT_2;
let cut = Curve::Ellipse { frame: Frame { origin: [0., 0., 30.], x_axis: [c, 0., c], y_axis: [0., 1., 0.] }, radii: [r * 2f64.sqrt(), r] };
let [_, v] = form::extent(&cylinder, &[vec![(&cut, [0.1, 0.1 + TAU], true)]]).expect("an extent");
assert!((v[0] - 10.0).abs() < 1e-9 && (v[1] - 50.0).abs() < 1e-9, "{v:?}");
}