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

kernel: the last truck wire, face and loft helpers leave the construction files

Browse Source 
The base no longer links truck, so build.rs and sweep_path.rs drop what was
kept only for other modules: chain and profile wires, the profile face,
loft wires and their ruled and smooth faces, across/section_parameters, the
cubic interpolations and NetFit, the uniform wires, the six sphere patch
faces, and sweep_path's ruled_sheet and smooth_sheet with the station helpers
only they used. Their callers switch to planar_face_record,
smlib_sweeps::{chain_curves, body_from_record, sphere}, sketch_record_curve,
the SMLib sweeps with open ends, and bspline::CurveN::interpolate.
offset_backend and polyhedron pick faces through provenance::closest_on_tri,
boolean_guard being gone.

Callers left to switch (other ports): surfaces.rs, loft_rails.rs,
surface_offset.rs, surface_ruled.rs, step_in/brep.rs, trimmed_face_mesh_tests.rs.
(not built: base in progress; no error in these files)

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_017itGh5fpMEousnZ91cYcVt
This commit is contained in:
jyviko
2026-10-22 07:20:11 -04:00
parent a6b41dd825
commit 453d6395e6
2 changed files with 10 additions and 746 deletions
Download Patch File Download Diff File Expand all files Collapse all files
crates/kernel/src/build.rs
+5 -498
View File
@@ -1,18 +1,13 @@
//! Building bodies: extrude and revolve planar regions, primitives, lofts and sweeps, each
//! through SMLib ([`crate::smlib_sweeps`]). The wire, face and curve helpers below (chain and
//! profile wires, ruled and smooth loft faces, cubic interpolation, sphere patch faces) are kept
//! only for the modules that still build on them (surfaces, blends, loft rails, the STEP
//! reader); they go when those are ported.
//! through SMLib ([`crate::smlib_sweeps`]); a profile face as its schema 4 record
//! ([`planar_face_record`]).
use mt::builder;
use solidcraft_geom::{Grid, Loop2, Plane, Region2, Seg2, Vec2, Vec3};
use truck_modeling as mt;
use solidcraft_geom::{Loop2, Plane, Region2, Seg2, Vec2, Vec3};
use crate::body::{Body, from_p3, p3, v3};
use crate::body::Body;
use crate::numbered::{self, NumberedLoop, NumberedRegion};
use crate::record::{self, Side};
use crate::sweep_path::RailFit;
use crate::{KernelError, Result, guard};
use crate::{KernelError, Result};
const MAX_SEGS: usize = 20_000;
@@ -20,121 +15,6 @@ fn finite(x: f64, what: &str) -> Result<f64> {
if x.is_finite() && x.abs() < 1e7 { Ok(x) } else { Err(KernelError::Invalid(format!("{what} = {x}"))) }
}
/// A closed wire in 3D from a sketch loop, with each edge's segment and that segment's
/// position in the loop (a big arc makes two edges).
fn wire_sources(plane: &Plane, lp: &Loop2) -> Result<(mt::Wire, Vec<(mt::Edge, Seg2, usize)>)> {
let n = lp.segs.len();
if n == 0 || n > MAX_SEGS {
return Err(KernelError::Invalid(format!("loop with {n} segments")));
}
if !lp.is_closed(1e-6) {
return Err(KernelError::Invalid("profile loop is not closed".into()));
}
chain_sources(plane, &lp.segs, true)
}
/// An open wire in 3D through a chain of sketch segments, each starting where the one before
/// ends (a closed chain makes a closed wire).
pub(crate) fn chain_wire(plane: &Plane, segs: &[Seg2]) -> Result<mt::Wire> {
chain_wire_seeds(plane, segs).map(|(wire, _)| wire)
}
/// [`chain_wire`], with each edge's id and the position in the chain of the segment it was
/// built from (a big arc makes two edges), in wire order.
pub(crate) fn chain_wire_seeds(plane: &Plane, segs: &[Seg2]) -> Result<(mt::Wire, Vec<(mt::EdgeID, usize)>)> {
let n = segs.len();
if n == 0 || n > MAX_SEGS {
return Err(KernelError::Invalid(format!("curve chain with {n} segments")));
}
if !segs.iter().all(seg_ok) {
return Err(KernelError::Invalid("non-finite curve".into()));
}
for (s, t) in segs.iter().zip(segs.iter().skip(1)) {
if s.end().dist(t.start()) > 1e-6 {
return Err(KernelError::Invalid("the curves are not connected end to end".into()));
}
}
let closed = segs.first().zip(segs.last()).is_some_and(|(f, l)| l.end().dist(f.start()) <= 1e-6);
let (wire, sources) = chain_sources(plane, segs, closed)?;
Ok((wire, sources.into_iter().map(|(e, _, i)| (e.id(), i)).collect()))
}
/// The edges of a chain of segments (closed: the last one ends on the first's start vertex),
/// with each edge's segment and that segment's position in the chain (a big arc makes two
/// edges).
fn chain_sources(plane: &Plane, segs: &[Seg2], closed: bool) -> Result<(mt::Wire, Vec<(mt::Edge, Seg2, usize)>)> {
let n = segs.len();
let mut verts: Vec<mt::Vertex> = segs.iter().map(|s| builder::vertex(p3(plane.to_world(s.start())))).collect();
if !closed && let Some(last) = segs.last() {
verts.push(builder::vertex(p3(plane.to_world(last.end()))));
}
let mut edges = Vec::with_capacity(n);
let mut sources = Vec::with_capacity(n);
for (i, s) in segs.iter().enumerate() {
let (Some(a), Some(b)) = (verts.get(i), verts.get(if closed { (i + 1) % n } else { i + 1 })) else { continue };
for e in seg_edges(plane, s, a, b)? {
sources.push((e.clone(), *s, i));
edges.push(e);
}
}
Ok((edges.into(), sources))
}
/// The edges of one segment from `a` to `b` (a big arc makes two).
fn seg_edges(plane: &Plane, s: &Seg2, a: &mt::Vertex, b: &mt::Vertex) -> Result<Vec<mt::Edge>> {
Ok(match *s {
Seg2::Line { .. } => vec![builder::line(a, b)],
Seg2::Arc { sweep, .. } => {
if sweep.abs() < 1e-9 {
return Err(KernelError::Invalid("zero-length arc".into()));
}
// Split big arcs so every piece is well conditioned.
if sweep.abs() > std::f64::consts::PI * 1.5 {
let vm = builder::vertex(p3(plane.to_world(s.mid())));
vec![
builder::circle_arc(a, &vm, p3(plane.to_world(s.point_at(0.25)))),
builder::circle_arc(&vm, b, p3(plane.to_world(s.point_at(0.75)))),
]
} else {
vec![builder::circle_arc(a, b, p3(plane.to_world(s.mid())))]
}
}
Seg2::Cubic { .. } | Seg2::Conic { .. } | Seg2::Bezier { .. } => {
vec![crate::freeform::edge(plane, s, a, b).ok_or_else(|| KernelError::Invalid("profile curve".into()))?]
}
})
}
/// An edge of a profile wire: the segment it was built from and the profile segments under it.
pub(crate) struct ProfileEdge {
pub edge: mt::Edge,
pub segment: Seg2,
pub numbers: Vec<usize>,
}
/// A wire built from a profile loop, with its edges in wire order.
pub(crate) struct ProfileWire {
pub wire: mt::Wire,
pub edges: Vec<ProfileEdge>,
}
/// The wires of a region oriented as the builders build it ([`NumberedRegion::oriented`]), each
/// with what its edges were built from.
fn region_wires(plane: &Plane, r: &NumberedRegion) -> Result<Vec<ProfileWire>> {
r.loops()
.map(|lp| {
let (wire, s) = wire_sources(plane, &lp.lp)?;
let edges = s.into_iter().map(|(edge, segment, i)| ProfileEdge { edge, segment, numbers: lp.numbers(i).to_vec() }).collect();
Ok(ProfileWire { wire, edges })
})
.collect()
}
/// The wires alone.
fn wires_of(profile: &[ProfileWire]) -> Vec<mt::Wire> {
profile.iter().map(|p| p.wire.clone()).collect()
}
/// A region numbered on its own, as given.
fn numbered_region(r: &Region2) -> Result<NumberedRegion> {
numbered::number_regions(std::slice::from_ref(r)).pop().ok_or_else(|| KernelError::Invalid("profile".into()))
@@ -145,15 +25,6 @@ fn oriented(r: &Region2) -> Result<NumberedRegion> {
numbered_region(r).map(|n| n.oriented())
}
/// A planar face from a region (outer loop with holes).
pub(crate) fn face(plane: &Plane, r: &Region2) -> Result<mt::Face> {
attach_profile(&wires_of(&region_wires(plane, &oriented(r)?)?))
}
fn attach_profile(wires: &[mt::Wire]) -> Result<mt::Face> {
guard("face", || builder::try_attach_plane(wires).map_err(|e| KernelError::Invalid(format!("profile: {e}"))))
}
pub(crate) fn region_ok(r: &Region2) -> Result<()> {
if r.area().abs() < 1e-9 {
return Err(KernelError::Invalid("profile has no area".into()));
@@ -1133,16 +1004,6 @@ fn smlib_loft_regions(outlines: &[(Plane, NumberedLoop)], planes: &[Plane], numb
crate::smlib_sweeps::loft(&sections, first)
}
/// The wires of a loft's closed sections, matched: oriented alike, resampled to the same
/// segment count and turned to start alike; with the profile segments under each wire edge, and
/// whether the sections were turned round into the canonical order.
pub(crate) fn loft_wires(sections: &[(Plane, Loop2)]) -> Result<(Vec<mt::Wire>, Vec<Vec<Vec<usize>>>, bool)> {
let reversed = loft_reversed(sections)?;
let ordered = sections.iter().map(|(p, _)| *p).zip(numbered::number_loops(sections.iter().map(|(_, l)| l))).collect::<Vec<_>>();
let (wires, sources) = guard("loft sections", || loft_numbered_wires(&ordered, reversed))?;
Ok((wires, sources, reversed))
}
pub(crate) fn loft_reversed(sections: &[(Plane, Loop2)]) -> Result<bool> {
if sections.len() < 2 || sections.len() > 50 {
return Err(KernelError::Invalid("a loft needs 2…50 sections".into()));
@@ -1175,17 +1036,6 @@ pub(crate) fn loft_reversed(sections: &[(Plane, Loop2)]) -> Result<bool> {
Ok(big < 0.0)
}
fn loft_numbered_wires(sections: &[(Plane, NumberedLoop)], reversed: bool) -> Result<(Vec<mt::Wire>, Vec<Vec<Vec<usize>>>)> {
let mut wires: Vec<mt::Wire> = Vec::new();
// Per section, the profile segments under each of its wire's edges.
let mut edge_sources: Vec<Vec<Vec<usize>>> = Vec::new();
for (p, r) in loft_numbered_loops(sections, reversed)? {
wires.push(uniform_wire(&p, &r.lp)?);
edge_sources.push(r.from);
}
Ok((wires, edge_sources))
}
/// The loft's sections as the walls are built through them: in canonical order (`reversed`
/// turns them round), each oriented the same way, re-split to the same number of segments and
/// started at the segment that lines up with the lead section's first.
@@ -1243,134 +1093,6 @@ pub(crate) fn loft_numbered_loops(sections: &[(Plane, NumberedLoop)], reversed:
Ok(out)
}
/// Ruled faces between each two neighbouring section wires; each comes from the section edges
/// it spans.
pub(crate) fn ruled_loft_faces(wires: &[mt::Wire], edge_sources: &[Vec<Vec<usize>>]) -> Result<Vec<Side>> {
let mut sides: Vec<Side> = Vec::new();
for (k, w) in wires.windows(2).enumerate() {
let (Some(a), Some(b)) = (w.first(), w.get(1)) else { continue };
let shell = builder::try_wire_homotopy(a, b).map_err(|e| KernelError::Failed(format!("loft: {e}")))?;
match (edge_sources.get(k), edge_sources.get(k + 1)) {
(Some(na), Some(nb)) => sides.extend(record::ruled(&shell, a, b, na, nb)),
_ => sides.extend(shell.face_iter().map(|f| (f.clone(), None))),
}
}
Ok(sides)
}
/// Points across the sections (one per section, in order) as a smooth B-spline through them at the
/// parameters `u` (rising from 0 to 1, one per section): cubic through four or more, quadratic
/// through three.
pub(crate) fn across(pts: &[Vec3], u: &[f64]) -> Option<mt::BSplineCurve<mt::Point3>> {
match pts {
[a, b, c] => {
// The quadratic Bezier through `b` at t = u[1].
let t = *u.get(1)?;
if !(t > 0.0 && t < 1.0) {
return None;
}
let m = (*b - *a * ((1.0 - t) * (1.0 - t)) - *c * (t * t)) * (1.0 / (2.0 * t * (1.0 - t)));
mt::BSplineCurve::try_new(mt::KnotVec::bezier_knot(2), vec![p3(*a), p3(m), p3(*c)]).ok()
}
_ => interpolate_cubic_at(pts, u),
}
}
/// The parameters of the sections along a loft: in proportion to the distance between them (the
/// mean over the wires' vertices), so that unevenly spaced sections get a surface that follows
/// them. At a uniform parameter a section close to its neighbour makes the surface overshoot
/// between the others (by several times the section's size when the spacing differs tenfold).
/// Evenly spaced sections give the uniform parameters.
pub(crate) fn section_parameters(edges: &[Vec<mt::Edge>]) -> Vec<f64> {
let m = edges.len();
let uniform = || (0..m).map(|i| i as f64 / (m.max(2) - 1) as f64).collect::<Vec<f64>>();
let corners: Vec<Vec<Vec3>> = edges.iter().map(|es| es.iter().map(|e| from_p3(e.front().point())).collect()).collect();
let mut u = vec![0.0];
for w in corners.windows(2) {
let (Some(a), Some(b)) = (w.first(), w.get(1)) else { return uniform() };
let n = a.len().min(b.len());
if n == 0 {
return uniform();
}
let step = a.iter().zip(b).map(|(p, q)| p.dist(*q)).sum::<f64>() / n as f64;
u.push(u.last().copied().unwrap_or(0.0) + step);
}
let Some(total) = u.last().copied().filter(|t| t.is_finite() && *t > 0.0) else { return uniform() };
let u: Vec<f64> = u.iter().map(|x| x / total).collect();
// A section on top of its neighbour: no parameter for it, so the uniform ones.
if u.windows(2).any(|w| w.get(1).zip(w.first()).is_none_or(|(b, a)| b - a < 1e-6)) {
return uniform();
}
u
}
/// A loft whose side faces are B-spline surfaces through every section (smooth across the
/// middle sections, which leave no edges). Every section wire comes from `uniform_wire`, so the
/// k-th edges of all sections share one knot vector and their control points pair up. Face k
/// comes from the profile segments under the k-th edge of every section (`numbers`: by
/// section, by edge).
/// One smooth face through every section for each wire edge (cubic across four sections or
/// more, quadratic across three), neighbours sharing their seams. `closed` wires join their
/// last edge to their first; open ones have a seam at each end.
pub(crate) fn smooth_loft_faces(wires: &[mt::Wire], numbers: &[Vec<Vec<usize>>], closed: bool) -> Result<Vec<Side>> {
let fail = |m: &str| KernelError::Failed(format!("smooth loft: {m}"));
let edges: Vec<Vec<mt::Edge>> = wires.iter().map(|w| w.edge_iter().cloned().collect()).collect();
let n = edges.first().map(Vec::len).ok_or_else(|| fail("no sections"))?;
if edges.iter().any(|e| e.len() != n) {
return Err(fail("sections with different edge counts"));
}
let curve_of = |e: &mt::Edge| -> Option<mt::BSplineCurve<mt::Point3>> {
match e.curve() {
mt::Curve::BSplineCurve(c) if e.orientation() => Some(c),
_ => None,
}
};
let (Some(first), Some(last)) = (edges.first(), edges.last()) else { return Err(fail("sections")) };
let along = section_parameters(&edges);
// Seams: through the k-th vertex of every section (an open wire's last one at the end of
// its last edge).
let vertex = |e: &mt::Edge, end: bool| if end { e.back().clone() } else { e.front().clone() };
let seam_count = if closed { n } else { n + 1 };
let mut seams: Vec<mt::Edge> = Vec::with_capacity(seam_count);
for k in 0..seam_count {
let (at, end) = if k < n { (k, false) } else { (n - 1, true) };
let pts: Vec<Vec3> = edges.iter().filter_map(|es| es.get(at).map(|e| from_p3(vertex(e, end).point()))).collect();
let c = across(&pts, &along).ok_or_else(|| fail("seam"))?;
let (Some(a), Some(b)) = (first.get(at), last.get(at)) else { return Err(fail("seam")) };
seams.push(mt::Edge::new(&vertex(a, end), &vertex(b, end), mt::Curve::BSplineCurve(c)));
}
let mut faces = Vec::with_capacity(n);
for k in 0..n {
let curves: Vec<mt::BSplineCurve<mt::Point3>> =
edges.iter().map(|es| es.get(k).and_then(curve_of)).collect::<Option<_>>().ok_or_else(|| fail("section curves"))?;
let c0 = curves.first().ok_or_else(|| fail("section curves"))?;
let uk = c0.knot_vec().clone();
let m = c0.control_points().len();
if curves.iter().any(|c| c.control_points().len() != m || c.knot_vec() != &uk) {
return Err(fail("sections parameterised differently"));
}
let mut ctrl: Vec<Vec<mt::Point3>> = Vec::with_capacity(m);
let mut vk = None;
for j in 0..m {
let pts: Vec<Vec3> = curves.iter().filter_map(|c| c.control_points().get(j).map(|p| from_p3(*p))).collect();
let c = across(&pts, &along).ok_or_else(|| fail("surface"))?;
vk = Some(c.knot_vec().clone());
ctrl.push(c.control_points().clone());
}
let vk = vk.ok_or_else(|| fail("surface"))?;
let surface = mt::BSplineSurface::try_new((uk, vk), ctrl).map_err(|e| fail(&e.to_string()))?;
let (Some(e0), Some(e1), Some(s0), Some(s1)) = (first.get(k), last.get(k), seams.get(k), seams.get((k + 1) % seam_count)) else {
return Err(fail("face"));
};
// u along the section edge, v across the sections: counter-clockwise in (u, v).
let w: mt::Wire = vec![e0.clone(), s1.clone(), e1.inverse(), s0.inverse()].into();
let face = mt::Face::try_new(vec![w], mt::Surface::BSplineSurface(surface)).map_err(|e| fail(&e.to_string()))?;
let under: Vec<&usize> = numbers.iter().filter_map(|s| s.get(k)).flatten().collect();
faces.push((face, (under.len() >= numbers.len().max(1)).then(|| record::profile(under))));
}
Ok(faces)
}
/// A 3D path segment for sweeps.
#[derive(Clone, Copy, Debug, PartialEq)]
pub enum PathSeg {
@@ -1558,218 +1280,3 @@ pub fn sweep_with(plane: &Plane, region: &Region2, path: &[PathSeg], ctl: &crate
pub fn sweep(plane: &Plane, region: &Region2, path: &[PathSeg]) -> Result<Body> {
smlib_sweep(plane, region, path, &crate::SweepControl::default())
}
/// Cox–de Boor basis values N_{j,3}(u) for all j.
fn cubic_basis(knots: &[f64], n_ctrl: usize, u: f64) -> Vec<f64> {
let p = 3;
let m = knots.len();
let mut n = vec![0.0; m.saturating_sub(1)];
let last = n_ctrl.saturating_sub(1);
for (j, v) in n.iter_mut().enumerate() {
let (Some(&a), Some(&b)) = (knots.get(j), knots.get(j + 1)) else { continue };
// Half-open spans, with the last span closed at the end.
if (a <= u && u < b) || (u >= b && (b - a) > 0.0 && j == last && (u - b).abs() < 1e-12) {
*v = 1.0;
}
}
for k in 1..=p {
for j in 0..m.saturating_sub(1 + k) {
let (Some(&tj), Some(&tjk), Some(&tj1), Some(&tjk1)) = (knots.get(j), knots.get(j + k), knots.get(j + 1), knots.get(j + k + 1)) else {
continue;
};
let left = if tjk - tj > 0.0 { (u - tj) / (tjk - tj) * n.get(j).copied().unwrap_or(0.0) } else { 0.0 };
let right = if tjk1 - tj1 > 0.0 { (tjk1 - u) / (tjk1 - tj1) * n.get(j + 1).copied().unwrap_or(0.0) } else { 0.0 };
if let Some(v) = n.get_mut(j) {
*v = left + right;
}
}
}
n.truncate(n_ctrl);
n
}
/// Cubic B-spline through points at uniform parameters (clamped, knots by averaging).
pub(crate) fn interpolate_cubic(pts: &[Vec3]) -> Option<mt::BSplineCurve<mt::Point3>> {
let m = pts.len();
let u: Vec<f64> = (0..m).map(|i| i as f64 / (m.max(2) - 1) as f64).collect();
interpolate_cubic_at(pts, &u)
}
/// A cubic B-spline surface through a grid of samples at uniform parameters: cubics across v
/// through each u index's samples, then cubics across u through their control points.
pub(crate) struct NetFit {
pub surface: mt::BSplineSurface<mt::Point3>,
/// The surface's edges at the first and the last u (curves across v).
pub u_ends: [mt::BSplineCurve<mt::Point3>; 2],
/// The surface's edges at the first and the last v (curves across u).
pub v_ends: [mt::BSplineCurve<mt::Point3>; 2],
}
/// [`NetFit`] through `samples`; `None` with fewer than four samples across either direction.
pub(crate) fn interpolate_net(samples: &Grid<Vec3>) -> Option<NetFit> {
let across_v: Vec<mt::BSplineCurve<mt::Point3>> = samples.transposed().lines().map(interpolate_cubic).collect::<Option<_>>()?;
let across_u: Vec<mt::BSplineCurve<mt::Point3>> =
control_net(&across_v)?.lines().map(|run| interpolate_cubic(&run.iter().map(|&p| from_p3(p)).collect::<Vec<_>>())).collect::<Option<_>>()?;
let net = control_net(&across_u)?.transposed().into_u_rows();
let knots = (across_u.first()?.knot_vec().clone(), across_v.first()?.knot_vec().clone());
let surface = mt::BSplineSurface::try_new(knots, net).ok()?;
let u_ends = [across_v.first()?.clone(), across_v.last()?.clone()];
let v_ends = [across_u.first()?.clone(), across_u.last()?.clone()];
Some(NetFit { surface, u_ends, v_ends })
}
/// The control points of `curves` as a grid: value `(i, j)` is control point `j` of curve `i`.
/// `None` unless every curve has as many control points.
fn control_net(curves: &[mt::BSplineCurve<mt::Point3>]) -> Option<Grid<mt::Point3>> {
let rows: Vec<&[mt::Point3]> = curves.iter().map(|c| c.control_points().as_slice()).collect();
Grid::from_u_rows(&rows)
}
/// [`interpolate_cubic`] by chord length: the parameters (0 to 1) in proportion to the distance
/// along the points. Points closer together toward an end (where a curve leaves another
/// tangentially) then don't slow the curve to a stop there, where a uniform parameter lets it
/// double back on itself.
pub(crate) fn interpolate_cubic_chordal(pts: &[Vec3]) -> Option<(mt::BSplineCurve<mt::Point3>, Vec<f64>)> {
let mut u = vec![0.0];
for w in pts.windows(2) {
let (Some(a), Some(b)) = (w.first(), w.get(1)) else { return None };
u.push(u.last().copied().unwrap_or(0.0) + a.dist(*b));
}
let total = u.last().copied().filter(|t| t.is_finite() && *t > 0.0)?;
let u: Vec<f64> = u.iter().map(|x| x / total).collect();
if u.windows(2).any(|w| w.first() >= w.get(1)) {
return None;
}
interpolate_cubic_at(pts, &u).map(|c| (c, u))
}
/// The cubic B-spline through `pts` at parameters `u` (rising from 0 to 1), its knots averaged.
fn interpolate_cubic_at(pts: &[Vec3], u: &[f64]) -> Option<mt::BSplineCurve<mt::Point3>> {
let m = pts.len();
if m < 4 || u.len() != m {
return None;
}
let mut knots = vec![0.0; 4];
for j in 1..(m - 3) {
knots.push((j..j + 3).filter_map(|i| u.get(i)).sum::<f64>() / 3.0);
}
knots.extend([1.0; 4]);
// Dense solve of N · P = Q for each coordinate.
let mut a: Vec<Vec<f64>> = u.iter().map(|ui| cubic_basis(&knots, m, *ui)).collect();
let mut q: Vec<[f64; 3]> = pts.iter().map(|p| [p.x, p.y, p.z]).collect();
for col in 0..m {
let piv = (col..m).max_by(|x, y| a[*x][col].abs().total_cmp(&a[*y][col].abs()))?;
if a[piv][col].abs() < 1e-14 {
return None;
}
a.swap(col, piv);
q.swap(col, piv);
for r in 0..m {
if r == col {
continue;
}
let f = a[r][col] / a[col][col];
if f == 0.0 {
continue;
}
let prow = a[col].clone();
for (x, v) in a[r].iter_mut().zip(&prow).skip(col) {
*x -= f * v;
}
let qc = q[col];
for k in 0..3 {
q[r][k] -= f * qc[k];
}
}
}
let ctrl: Vec<mt::Point3> = (0..m).map(|i| mt::Point3::new(q[i][0] / a[i][i], q[i][1] / a[i][i], q[i][2] / a[i][i])).collect();
mt::BSplineCurve::try_new(mt::KnotVec::from(knots), ctrl).ok()
}
/// A closed wire whose edges are cubic B-splines parameterised (nearly) by arc length, so that
/// ruled surfaces between sections pair points proportionally along each edge.
pub(crate) fn uniform_wire(plane: &Plane, lp: &Loop2) -> Result<mt::Wire> {
uniform_chain_wire(plane, &lp.segs, true)
}
/// [`uniform_wire`] of a chain of segments, open or `closed`.
pub(crate) fn uniform_chain_wire(plane: &Plane, segs: &[Seg2], closed: bool) -> Result<mt::Wire> {
const SAMPLES: usize = 24;
let n = segs.len();
let mut verts: Vec<mt::Vertex> = segs.iter().map(|s| builder::vertex(p3(plane.to_world(s.start())))).collect();
if !closed && let Some(last) = segs.last() {
verts.push(builder::vertex(p3(plane.to_world(last.end()))));
}
let mut edges = Vec::with_capacity(n);
for (i, s) in segs.iter().enumerate() {
let (Some(a), Some(b)) = (verts.get(i), verts.get(if closed { (i + 1) % n } else { i + 1 })) else { continue };
let pts: Vec<Vec3> = (0..SAMPLES).map(|k| plane.to_world(s.point_at(k as f64 / (SAMPLES - 1) as f64))).collect();
let curve = interpolate_cubic(&pts).ok_or_else(|| KernelError::Failed("loft curve".into()))?;
edges.push(mt::Edge::new(a, b, mt::Curve::BSplineCurve(curve)));
}
Ok(edges.into())
}
/// Six outward patch faces of a sphere (the cube's faces projected onto it) on revolved surfaces
/// whose poles and seams lie outside each patch, for the STEP reader's spheres; callers guard them.
pub(crate) fn sphere_patch_faces(center: Vec3, radius: f64, rot: [Vec3; 3]) -> Result<Vec<mt::Face>> {
let radius = finite(radius, "radius")?;
let turn = |v: Vec3| rot[0] * v.x + rot[1] * v.y + rot[2] * v.z;
if radius <= 1e-6 || !center.is_finite() {
return Err(KernelError::Invalid("radius must be positive".into()));
}
let k = radius / 3f64.sqrt();
let corner = |sx: f64, sy: f64, sz: f64| center + turn(Vec3::new(sx * k, sy * k, sz * k));
// Cube faces: outward direction, rotation axis for the patch's surface (perpendicular to it),
// and the four corner sign triples counter-clockwise seen from outside.
let faces: [(Vec3, Vec3, [[f64; 3]; 4]); 6] = [
(Vec3::X, Vec3::Z, [[1.0, -1.0, -1.0], [1.0, 1.0, -1.0], [1.0, 1.0, 1.0], [1.0, -1.0, 1.0]]),
(-Vec3::X, Vec3::Z, [[-1.0, 1.0, -1.0], [-1.0, -1.0, -1.0], [-1.0, -1.0, 1.0], [-1.0, 1.0, 1.0]]),
(Vec3::Y, Vec3::Z, [[1.0, 1.0, -1.0], [-1.0, 1.0, -1.0], [-1.0, 1.0, 1.0], [1.0, 1.0, 1.0]]),
(-Vec3::Y, Vec3::Z, [[-1.0, -1.0, -1.0], [1.0, -1.0, -1.0], [1.0, -1.0, 1.0], [-1.0, -1.0, 1.0]]),
(Vec3::Z, Vec3::X, [[-1.0, -1.0, 1.0], [1.0, -1.0, 1.0], [1.0, 1.0, 1.0], [-1.0, 1.0, 1.0]]),
(-Vec3::Z, Vec3::X, [[-1.0, 1.0, -1.0], [1.0, 1.0, -1.0], [1.0, -1.0, -1.0], [-1.0, -1.0, -1.0]]),
];
let key = |s: [f64; 3]| ((s[0] > 0.0) as usize) | (((s[1] > 0.0) as usize) << 1) | (((s[2] > 0.0) as usize) << 2);
let mut verts: std::collections::HashMap<usize, mt::Vertex> = std::collections::HashMap::new();
let mut vert = |s: [f64; 3]| verts.entry(key(s)).or_insert_with(|| builder::vertex(p3(corner(s[0], s[1], s[2])))).clone();
let mut edges: std::collections::HashMap<(usize, usize), mt::Edge> = std::collections::HashMap::new();
let mut out: Vec<mt::Face> = Vec::new();
for (dir, axis, cs) in faces {
let (dir, axis) = (turn(dir), turn(axis));
let mut wire: Vec<mt::Edge> = Vec::new();
for i in 0..4 {
let (a, b) = (cs[i], cs[(i + 1) % 4]);
let (ka, kb) = (key(a), key(b));
let e = match edges.get(&(ka.min(kb), ka.max(kb))) {
Some(e) => e.clone(),
None => {
let (lo, hi) = if ka < kb { (a, b) } else { (b, a) };
let (va, vb) = (vert(lo), vert(hi));
let mid = (corner(lo[0], lo[1], lo[2]) + corner(hi[0], hi[1], hi[2])) * 0.5 - center;
let transit = center + mid.normalized().unwrap_or(dir) * radius;
let e = builder::circle_arc(&va, &vb, p3(transit));
edges.insert((ka.min(kb), ka.max(kb)), e.clone());
e
}
};
wire.push(if ka < kb { e } else { e.inverse() });
}
// The meridian runs pole to pole through the side away from the patch.
let (np, sp, back) = (center + axis * radius, center - axis * radius, center - dir * radius);
let meridian = builder::circle_arc(&builder::vertex(p3(np)), &builder::vertex(p3(sp)), p3(back));
let mut surface =
mt::Surface::RevolutedCurve(mt::Processor::new(mt::RevolutedCurve::by_revolution(meridian.oriented_curve(), p3(center), v3(axis))));
// Outward normal at the patch centre.
let probe = p3(center + dir * radius);
use mt::{ParametricSurface3D, SearchParameter};
if let Some((u, v)) = surface.search_parameter(probe, mt::SPHint2D::None, 100) {
let n = surface.normal(u, v);
if mt::InnerSpace::dot(n, v3(dir)) < 0.0 {
mt::Invertible::invert(&mut surface);
}
}
out.push(mt::Face::new(vec![wire.into()], surface));
}
Ok(out)
}
crates/kernel/src/sweep_path.rs
+5 -248
View File
@@ -2,39 +2,20 @@
//! and the laws that place a profile along a path: the rotation-minimising frame and the sweep
//! options (orientation, twist, end scale, guide rails; see [`SweepControl`]) give the stations
//! SMLib's controlled sweeps skin through ([`crate::smlib_sweeps`]). A solid along points is
//! SMLib's sweep along the lines between them ([`sweep_path_with`]). The open walls of chains
//! along a path ([`ruled_sheet`], [`smooth_sheet`]) are still built here for the surface sweeps.
//! SMLib's sweep along the lines between them ([`sweep_path_with`]).
use mt::builder;
use solidcraft_geom::{Grid, Plane, Region2, Seg2, Vec2, Vec3};
use truck_modeling as mt;
use solidcraft_geom::{Plane, Region2, Seg2, Vec2, Vec3};
use crate::body::{Body, p3};
use crate::build::{NetFit, interpolate_net};
use crate::body::Body;
use crate::numbered::NumberedLoop;
use crate::record::{self, Side};
use crate::{KernelError, Result, guard};
use crate::{KernelError, Result};
/// Most path points accepted.
const MAX_PATH_POINTS: usize = 100_000;
/// Stations of a smooth sweep (the B-spline walls interpolate them).
const MAX_SMOOTH_STATIONS: usize = 120;
/// Stations of a ruled sweep.
const MAX_RULED_STATIONS: usize = 2000;
/// A turn sharper than this between path segments is a corner (radians).
const CORNER: f64 = 0.35;
/// Samples per profile segment of the smooth walls.
/// Samples per profile segment where a closed sweep checks that its last section lands on its first.
const SEG_SAMPLES: usize = 12;
/// Stations of a smooth sweep that twists or scales, at the least.
const MIN_CONTROLLED_STATIONS: usize = 16;
/// Stations of a smooth sweep that twists or scales, at the most.
const MAX_CONTROLLED_STATIONS: usize = 160;
/// Stations of a ruled sweep that twists or scales, at the least.
const MIN_CONTROLLED_RULED: usize = 64;
/// The most twist one station of a smooth sweep may carry (radians).
const SMOOTH_TWIST_STEP: f64 = std::f64::consts::PI / 18.0;
/// The most twist one station of a ruled sweep may carry (radians).
const RULED_TWIST_STEP: f64 = std::f64::consts::PI / 36.0;
/// How the profile turns as it follows a sweep path.
#[derive(Clone, Copy, Debug, Default, PartialEq)]
@@ -211,28 +192,6 @@ fn turn(a: Vec3, b: Vec3) -> f64 {
}
}
/// Points evenly spaced by arc length along the polyline (ends kept).
fn resample(p: &[Vec3], n: usize) -> Vec<Vec3> {
let mut cum = vec![0.0];
for w in p.windows(2) {
let l = cum.last().copied().unwrap_or(0.0) + w[0].dist(w[1]);
cum.push(l);
}
let total = cum.last().copied().unwrap_or(0.0);
let mut out = Vec::with_capacity(n + 1);
let mut k = 0;
for i in 0..=n {
let s = total * i as f64 / n as f64;
while k + 2 < cum.len() && cum.get(k + 1).is_some_and(|c| *c < s) {
k += 1;
}
let (Some(a), Some(b), Some(ca), Some(cb)) = (p.get(k), p.get(k + 1), cum.get(k), cum.get(k + 1)) else { break };
let t = if cb - ca > 1e-15 { ((s - ca) / (cb - ca)).clamp(0.0, 1.0) } else { 0.0 };
out.push(*a + (*b - *a) * t);
}
out
}
/// Rotation-minimising frames (double reflection) at the points, for the given unit tangents;
/// each station maps the first frame onto its own.
fn frames(pts: &[Vec3], tangents: &[Vec3]) -> Option<Vec<[Vec3; 3]>> {
@@ -472,40 +431,6 @@ pub(crate) fn rail_frame(along: Vec3, across: Vec3, along_to: Vec3, across_to: V
[map(Vec3::X), map(Vec3::Y), map(Vec3::Z)]
}
/// The stations of a smooth sweep: as many as the path has points (4…[`MAX_SMOOTH_STATIONS`]),
/// more when the profile twists or scales.
fn smooth_count(points: usize, ctl: &SweepControl) -> usize {
let n = points.saturating_sub(1).clamp(4, MAX_SMOOTH_STATIONS);
if !ctl.varies() {
return n;
}
let twist = (ctl.twist.abs() / SMOOTH_TWIST_STEP).ceil() as usize;
n.max(MIN_CONTROLLED_STATIONS).max(twist).min(MAX_CONTROLLED_STATIONS)
}
/// The path with points added along its straight runs (its corners kept) when the profile
/// twists or scales, so ruled walls follow the change.
fn densify(pts: &[Vec3], ctl: &SweepControl) -> Vec<Vec3> {
if !ctl.varies() {
return pts.to_vec();
}
let total: f64 = pts.windows(2).map(|w| w[0].dist(w[1])).sum();
let want = MIN_CONTROLLED_RULED.max((ctl.twist.abs() / RULED_TWIST_STEP).ceil() as usize).min(MAX_RULED_STATIONS / 2);
let mut out = Vec::with_capacity(want + pts.len());
out.extend(pts.first().copied());
for w in pts.windows(2) {
let pieces = if total > 0.0 { ((w[0].dist(w[1]) / total) * want as f64).ceil() as usize } else { 1 };
let pieces = pieces.clamp(1, want.max(1));
out.extend((1..=pieces).map(|j| w[0].lerp(w[1], j as f64 / pieces as f64)));
}
out
}
/// The profile plane carried to a station.
fn station_plane(plane: &Plane, s0: Vec3, st: &Station) -> Plane {
Plane { origin: st.at + st.apply(plane.origin - s0), x: st.apply(plane.x), y: st.apply(plane.y) }
}
/// Sweep a planar region along a 3D polyline path. The profile keeps its relation to the path
/// (it need not touch it); the path should start near the profile.
pub fn sweep_path(plane: &Plane, region: &Region2, path: &[Vec3]) -> Result<Body> {
@@ -531,174 +456,6 @@ pub(crate) fn smlib_sweep_points(plane: &Plane, region: &Region2, path: &crate::
crate::sweep_with(plane, region, &segs, ctl)
}
fn central_tangents(p: &[Vec3]) -> Result<Vec<Vec3>> {
let n = p.len();
let closed = closed_path(p);
(0..n)
.map(|i| {
let a = p.get(if closed && (i == 0 || i + 1 == n) { n - 2 } else { i.saturating_sub(1) }).copied().unwrap_or_default();
let b = p.get(if closed && (i == 0 || i + 1 == n) { 1 } else { (i + 1).min(n - 1) }).copied().unwrap_or_default();
(b - a).normalized().ok_or_else(|| KernelError::Invalid("the sweep path has no direction".into()))
})
.collect()
}
/// Stations at the path's points (at most [`MAX_RULED_STATIONS`]), mitred at corners, with the
/// points they stand at.
fn ruled_stations(plane: &Plane, pts: &[Vec3], ctl: &SweepControl) -> Result<(Vec<Vec3>, Vec<Station>)> {
let pts = densify(pts, ctl);
let pts: Vec<Vec3> = if pts.len() > MAX_RULED_STATIONS { resample(&pts, MAX_RULED_STATIONS) } else { pts };
let n = pts.len();
let closed = closed_path(&pts);
let seg_dir = |i: usize| -> Option<Vec3> { (*pts.get(i + 1)? - *pts.get(i)?).normalized() };
let mut tangents = Vec::with_capacity(n);
let mut mitre = Vec::with_capacity(n);
for i in 0..n {
let d_in = if i > 0 {
seg_dir(i - 1)
} else if closed {
seg_dir(n - 2)
} else {
None
};
let d_out = if i + 1 < n {
seg_dir(i)
} else if closed {
seg_dir(0)
} else {
None
};
let (t, m) = match (d_in, d_out) {
(Some(a), Some(b)) => {
let t = (a + b).normalized().ok_or_else(|| KernelError::Invalid("the sweep path turns back on itself".into()))?;
let half = turn(a, b) / 2.0;
// Across the bend the section stretches by 1 / cos(half the turn).
let m = (b - a).normalized().map(|bend| (bend, 1.0 / half.cos().max(0.2)));
(t, m)
}
(Some(a), None) => (a, None),
(None, Some(b)) => (b, None),
(None, None) => return Err(KernelError::Invalid("sweep path".into())),
};
tangents.push(t);
mitre.push(m);
}
let st = if ctl.is_default() { stations(&pts, &tangents, &mitre)? } else { controlled_stations(plane, &pts, &tangents, &mitre, ctl)? };
Ok((pts, st))
}
/// Walls of chains of sketch segments on `plane` (open, or closed) carried along a 3D polyline
/// path that starts near them: ruled between stations at the path's points, mitred at corners,
/// with no caps. Neighbouring walls share their edges.
pub(crate) fn ruled_sheet(plane: &Plane, chains: &[Vec<Seg2>], path: &[Vec3], ctl: &SweepControl) -> Result<Vec<Side>> {
ctl.check()?;
let (pts, st) = ruled_stations(plane, &clean(path)?, ctl)?;
let edges: usize = chains.iter().map(Vec::len).sum();
if st.len().saturating_mul(edges) > 400_000 {
return Err(KernelError::Invalid("the sweep is too detailed (curves × path points)".into()));
}
let s0 = pts.first().copied().unwrap_or_default();
// Each chain's segments are numbered on from the last chain's; every station's wire of a
// chain has the same edges, each over the same segment.
let offsets: Vec<usize> = chains.iter().scan(0, |o, c| Some(std::mem::replace(o, *o + c.len()))).collect();
let mut wires: Vec<Vec<mt::Wire>> = Vec::with_capacity(st.len());
let mut numbers: Vec<Vec<Vec<usize>>> = Vec::new();
for s in &st {
let pl = station_plane(plane, s0, s);
let mut at = Vec::with_capacity(chains.len());
for (c, offset) in chains.iter().zip(&offsets) {
let (wire, seeds) = crate::build::chain_wire_seeds(&pl, c)?;
if numbers.len() < chains.len() {
numbers.push(seeds.iter().map(|(_, i)| vec![offset + i]).collect());
}
at.push(wire);
}
wires.push(at);
}
guard("path sweep", || {
let mut faces: Vec<Side> = Vec::new();
for w in wires.windows(2) {
let (Some(w0), Some(w1)) = (w.first(), w.get(1)) else { continue };
for ((a, b), n) in w0.iter().zip(w1).zip(&numbers) {
let sides: mt::Shell = builder::try_wire_homotopy(a, b).map_err(|e| KernelError::Failed(format!("path sweep: {e}")))?;
faces.extend(record::ruled(&sides, a, b, n, n));
}
}
Ok(faces)
})
}
/// Walls of chains of sketch segments on `plane` carried along a smooth 3D polyline path (no
/// turn sharper than [`CORNER`]): one B-spline wall per segment through evenly spaced,
/// rotation-minimising stations, neighbouring walls sharing their seam edges; no caps. `None`
/// for a path with corners (see [`ruled_sheet`]).
pub(crate) fn smooth_sheet(plane: &Plane, chains: &[Vec<Seg2>], path: &[Vec3], ctl: &SweepControl) -> Result<Option<Vec<Side>>> {
ctl.check()?;
let pts = clean(path)?;
if pts.windows(3).any(|w| turn(w[1] - w[0], w[2] - w[1]) > CORNER) {
return Ok(None);
}
let fail = |w: &str| KernelError::Failed(format!("path sweep: {w}"));
if chains.iter().map(Vec::len).sum::<usize>() > 4096 {
return Err(fail("too many curves for a smooth sweep"));
}
let n = smooth_count(pts.len(), ctl);
let sp = resample(&pts, n);
let tangents = central_tangents(&sp)?;
let st = if ctl.is_default() { stations(&sp, &tangents, &[])? } else { controlled_stations(plane, &sp, &tangents, &[], ctl)? };
let s0 = st.first().map(|s| s.at).unwrap_or_default();
let at = |i: usize, q: Vec2| st.get(i).map(|s| s.at + s.apply(plane.to_world(q) - s0)).unwrap_or_default();
let n = st.len() - 1;
let mut faces: Vec<Side> = Vec::new();
let mut offset = 0;
for chain in chains {
let closed = chain.first().zip(chain.last()).is_some_and(|(f, l)| l.end().dist(f.start()) <= 1e-6);
// Each segment: its curve at the first and last station, the columns through its start
// and end, and the wall surface.
let mut data = Vec::with_capacity(chain.len());
for s in chain {
let samples =
Grid::from_fn([n + 1, SEG_SAMPLES], |i, j| at(i, s.point_at(j as f64 / (SEG_SAMPLES - 1) as f64))).ok_or_else(|| fail("samples"))?;
data.push(interpolate_net(&samples).ok_or_else(|| fail("wall surface"))?);
}
let walls: Vec<mt::Face> = guard("path sweep", || {
let m = chain.len();
// Vertices along the chain (an open chain has one more than it has segments).
let corners: Vec<Vec2> = chain.iter().map(Seg2::start).chain(chain.last().map(Seg2::end).filter(|_| !closed)).collect();
let start: Vec<mt::Vertex> = corners.iter().map(|q| builder::vertex(p3(at(0, *q)))).collect();
let end: Vec<mt::Vertex> = corners.iter().map(|q| builder::vertex(p3(at(n, *q)))).collect();
let next = |j: usize| if closed { (j + 1) % m } else { j + 1 };
let mut seams = Vec::with_capacity(corners.len());
for (k, (a, b)) in start.iter().zip(&end).enumerate() {
// The column through each corner: the start column of the segment leaving it, or
// the end column of the last segment.
let col = match data.get(k) {
Some(NetFit { v_ends: [first, _], .. }) => first,
None => data.last().map(|NetFit { v_ends: [_, last], .. }| last).ok_or_else(|| fail("seam"))?,
};
seams.push(mt::Edge::try_new(a, b, mt::Curve::BSplineCurve(col.clone())).map_err(|e| fail(&e.to_string()))?);
}
let mut walls = Vec::with_capacity(m);
for (j, NetFit { surface: surf, u_ends: [c0, c1], .. }) in data.iter().enumerate() {
let (Some(a0), Some(b0), Some(a1), Some(b1), Some(sa), Some(sb)) =
(start.get(j), start.get(next(j)), end.get(j), end.get(next(j)), seams.get(j), seams.get(next(j)))
else {
return Err(fail("edges"));
};
let e0 = mt::Edge::try_new(a0, b0, mt::Curve::BSplineCurve(c0.clone())).map_err(|e| fail(&e.to_string()))?;
let e1 = mt::Edge::try_new(a1, b1, mt::Curve::BSplineCurve(c1.clone())).map_err(|e| fail(&e.to_string()))?;
let w: mt::Wire = vec![sa.clone(), e1, sb.inverse(), e0.inverse()].into();
walls.push(mt::Face::try_new(vec![w], mt::Surface::BSplineSurface(surf.clone())).map_err(|e| fail(&e.to_string()))?);
}
Ok(walls)
})?;
// Wall j comes from segment j of the chain.
faces.extend(walls.into_iter().enumerate().map(|(j, w)| (w, Some(record::profile(&[offset + j])))));
offset += chain.len();
}
Ok(Some(faces))
}
#[cfg(test)]
mod tests {
use std::f64::consts::{FRAC_PI_2, FRAC_PI_4, PI};