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

brep: one basis-derivative routine, bspline::basis_derivatives

Browse Source 
The ruled surface's far-curve solve kept its own copy of NURBS Book A2.3. The
geom crate now exposes bspline::basis_derivatives (the basis functions alive at
a parameter and their derivatives up to an order, each from its local knots via
local_basis, on the span tess::span_of picks); tess's second-derivative surface
evaluation and the ruled surface both use it, and the copy is gone. smlib's
mod list back in order.

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 08:14:38 -04:00
parent 9e2c1c57c8
commit 78937d1dbc
4 changed files with 29 additions and 87 deletions
Download Patch File Download Diff File Expand all files Collapse all files
crates/geom/src/bspline.rs
+16
View File
@@ -161,6 +161,22 @@ pub fn local_basis(knots: &[f64], u: f64, side: Side, n: usize) -> Vec<f64> {
out
}
/// The `degree + 1` basis functions of a B-spline with `knots` and `count` control points alive
/// at `t`, with their derivatives up to `order`: the index of the first of them and, for each,
/// its value and derivatives (element `k` the `k`-th). Each comes from its own local knots
/// ([`local_basis`]), on the span used at `t` (the last non-empty span starting at or before it,
/// the domain's last span at its end, from the left).
pub fn basis_derivatives(knots: &[f64], degree: usize, count: usize, t: f64, order: usize) -> Option<(usize, Vec<Vec<f64>>)> {
if degree == 0 || count <= degree || knots.len() != count + degree + 1 || !t.is_finite() {
return None;
}
let span = crate::tess::span_of(knots, degree, count, t);
let side = if t >= *knots.get(count)? { Side::Left } else { Side::Right };
let first = span.checked_sub(degree)?;
let ders = (first..=span).map(|i| Some(local_basis(knots.get(i..=i + degree + 1)?, t, side, order))).collect::<Option<Vec<_>>>()?;
Some((first, ders))
}
/// [`local_basis`] into a fixed array (no allocation): element `k` is the `k`-th derivative,
/// for `k ≤ min(n, MAX_SPLINE_DEGREE)`, the rest zero.
pub(crate) fn local_basis_array(knots: &[f64], u: f64, side: Side, n: usize) -> [f64; MAX_SPLINE_DEGREE + 1] {
crates/geom/src/tess.rs
+5 -12
View File
@@ -349,16 +349,10 @@ impl EvalSurface {
}
/// The values and first two derivatives of the `p + 1` basis functions of degree `p` alive at
/// `t` on `span`, each from its own local knots (`bspline::local_basis`).
fn basis_second(knots: &[f64], span: usize, p: usize, t: f64) -> Option<Vec<[f64; 3]>> {
let n = knots.len().checked_sub(p + 1)?;
let side = if t >= *knots.get(n)? { crate::bspline::Side::Left } else { crate::bspline::Side::Right };
(span.checked_sub(p)?..=span)
.map(|i| {
let d = crate::bspline::local_basis(knots.get(i..=i + p + 1)?, t, side, 2);
Some([*d.first()?, *d.get(1)?, *d.get(2)?])
})
.collect()
/// `t` ([`crate::bspline::basis_derivatives`]), and the span they live on.
fn basis_second(knots: &[f64], p: usize, n: usize, t: f64) -> Option<(usize, Vec<[f64; 3]>)> {
let (first, ders) = crate::bspline::basis_derivatives(knots, p, n, t, 2)?;
Some((first + p, ders.iter().map(|d| Some([*d.first()?, *d.get(1)?, *d.get(2)?])).collect::<Option<_>>()?))
}
fn nurbs_second(degree: [usize; 2], knots: &[Vec<f64>; 2], count: [usize; 2], ctrl: &[[f64; 4]], u: f64, v: f64) -> Option<(Vec3, Vec3, Vec3)> {
@@ -368,8 +362,7 @@ fn nurbs_second(degree: [usize; 2], knots: &[Vec<f64>; 2], count: [usize; 2], ct
if !ok(pu, nu, &knots[0]) || !ok(pv, nv, &knots[1]) || nu.checked_mul(nv)? != ctrl.len() {
return None;
}
let (su, sv) = (span_of(&knots[0], pu, nu, u), span_of(&knots[1], pv, nv, v));
let (bu, bv) = (basis_second(&knots[0], su, pu, u)?, basis_second(&knots[1], sv, pv, v)?);
let ((su, bu), (sv, bv)) = (basis_second(&knots[0], pu, nu, u)?, basis_second(&knots[1], pv, nv, v)?);
// Homogeneous sums: [A, A_u, A_v, A_uu, A_uv, A_vv].
let mut a = [[0.0f64; 4]; 6];
for (j, bj) in bv.iter().enumerate() {
crates/kernel/src/surface_ruled.rs
+7 -74
View File
@@ -550,82 +550,15 @@ fn end_directions(a: Vec3, b: Vec3) -> Option<(Vec3, Vec3)> {
Some((chord, chord))
}
/// The B-spline basis functions non-zero at `t` and their first two derivatives (The NURBS Book
/// A2.3): the span's first control point index and, per derivative, `degree + 1` values.
fn basis(knots: &[f64], degree: usize, count: usize, t: f64) -> Option<(usize, [Vec<f64>; 3])> {
let (lo, hi) = (*knots.get(degree)?, *knots.get(count)?);
let t = t.clamp(lo, hi);
let mut span = degree;
while span + 1 < count && *knots.get(span + 1)? <= t {
span += 1;
}
let p = degree;
let k = |i: usize| knots.get(i).copied().unwrap_or(0.0);
let mut ndu = vec![vec![0.0; p + 1]; p + 1];
let (mut left, mut right) = (vec![0.0; p + 1], vec![0.0; p + 1]);
ndu[0][0] = 1.0;
for j in 1..=p {
left[j] = t - k(span + 1 - j);
right[j] = k(span + j) - t;
let mut saved = 0.0;
for r in 0..j {
ndu[j][r] = right[r + 1] + left[j - r];
let tmp = if ndu[j][r] != 0.0 { ndu[r][j - 1] / ndu[j][r] } else { 0.0 };
ndu[r][j] = saved + right[r + 1] * tmp;
saved = left[j - r] * tmp;
}
ndu[j][j] = saved;
}
let mut ders = [vec![0.0; p + 1], vec![0.0; p + 1], vec![0.0; p + 1]];
for j in 0..=p {
ders[0][j] = ndu[j][p];
}
let mut a = vec![vec![0.0; p + 1]; 2];
for r in 0..=p {
let (mut s1, mut s2) = (0usize, 1usize);
a[0][0] = 1.0;
for kk in 1..=2usize.min(p) {
let mut d = 0.0;
let (rk, pk) = (r as isize - kk as isize, p - kk);
if r >= kk {
let ru = rk as usize;
a[s2][0] = a[s1][0] / if ndu[pk + 1][ru] != 0.0 { ndu[pk + 1][ru] } else { 1.0 };
d = a[s2][0] * ndu[ru][pk];
}
let j1 = if rk >= -1 { 1 } else { (-rk) as usize };
let j2 = if (r as isize - 1) <= pk as isize { kk - 1 } else { p - r };
for j in j1..=j2 {
let ru = (rk + j as isize) as usize;
let den = ndu[pk + 1][ru];
a[s2][j] = (a[s1][j] - a[s1][j - 1]) / if den != 0.0 { den } else { 1.0 };
d += a[s2][j] * ndu[ru][pk];
}
if r <= pk {
let den = ndu[pk + 1][r];
a[s2][kk] = -a[s1][kk - 1] / if den != 0.0 { den } else { 1.0 };
d += a[s2][kk] * ndu[r][pk];
}
ders[kk][r] = d;
std::mem::swap(&mut s1, &mut s2);
}
}
let mut factor = p as f64;
for (kk, row) in ders.iter_mut().enumerate().skip(1).take(2usize.min(p)) {
for v in row.iter_mut() {
*v *= factor;
}
factor *= (p - kk) as f64;
}
Some((span - p, ders))
}
/// The rational basis of a curve with `c`'s knots and weights at `t`, and its first two
/// derivatives: the coefficients of each control point in the curve's point, first and second
/// derivative there (by its parameter), as (control point index, [R, R', R'']).
fn rational_basis(c: &CurveN<4>, t: f64) -> Option<Vec<(usize, [f64; 3])>> {
let (first, n) = basis(&c.knots, c.degree, c.ctrl.len(), t)?;
let [t0, t1] = range(c);
let (first, n) = solidcraft_geom::bspline::basis_derivatives(&c.knots, c.degree, c.ctrl.len(), t.clamp(t0, t1), 2)?;
let w: Vec<f64> = (0..=c.degree).map(|j| c.ctrl.get(first + j).map_or(1.0, |h| h[3])).collect();
let weight = |k: usize| (0..=c.degree).map(|j| n[k][j] * w[j]).sum::<f64>();
let nd = |j: usize, k: usize| n.get(j).and_then(|d| d.get(k)).copied().unwrap_or(0.0);
let weight = |k: usize| (0..=c.degree).map(|j| nd(j, k) * w[j]).sum::<f64>();
let (w0, w1, w2) = (weight(0), weight(1), weight(2));
if !(w0.abs() > 1e-300) {
return None;
@@ -633,9 +566,9 @@ fn rational_basis(c: &CurveN<4>, t: f64) -> Option<Vec<(usize, [f64; 3])>> {
Some(
(0..=c.degree)
.map(|j| {
let r = n[0][j] * w[j] / w0;
let r1 = (n[1][j] * w[j] - r * w1) / w0;
let r2 = (n[2][j] * w[j] - 2.0 * r1 * w1 - r * w2) / w0;
let r = nd(j, 0) * w[j] / w0;
let r1 = (nd(j, 1) * w[j] - r * w1) / w0;
let r2 = (nd(j, 2) * w[j] - 2.0 * r1 * w1 - r * w2) / w0;
(first + j, [r, r1, r2])
})
.collect(),
crates/smlib/src/lib.rs
+1 -1
View File
@@ -23,8 +23,8 @@ pub mod history;
pub mod import;
mod measure;
mod ops;
mod split;
mod sheets;
mod split;
pub mod ssi;
pub mod sweeps;
pub mod trim;