masi/solvecraft
1
0
Fork 0
mirror of https://github.com/bherbruck/solidcraft.git synced 2026-11-08 04:52:10 +00:00
Code Issues Packages Projects Releases 1 Wiki Activity

Merge pull request #1533 from tangency-labs/port/approximate-surface

Browse Source 
port: SmSurface::ApproximateSurface (surface.surface_approximate)
This commit is contained in:
bherbruck
2026-10-29 20:11:24 -04:00
committed by GitHub
parent 6a4d9b3f88 fa675aef3f
commit ddfb389578
8 changed files with 3724 additions and 2 deletions
Download Patch File Download Diff File Expand all files Collapse all files
crates/difftest/src/units/surface_surface_approximate.rs
+585
View File
@@ -0,0 +1,585 @@
//! Runner for unit `surface.surface_approximate`: `SmSurface::ApproximateSurface` on real
//! B-spline surfaces, planes, cylinders, cones, spheres, tori and surfaces of revolution,
//! bare or inside up to three offset, curve bounded, STEP and volume (rigid transform)
//! surfaces, so every `GetRootSurface` override is reached.
//!
//! C++: `solidcore_ffi::units::surface_surface_approximate` (shim
//! `harness/cpp/shim/surface_surface_approximate.cpp`, which documents the block).
//! Rust: [`solidcore::surface::surface::Surface::approximate_surface_ignoring_status`].
//!
//! A case is the shim's input block (`tol level just_copy n_wraps wraps.. kind` and the
//! kind's construction data, as `surface.surface_assert_valid`'s); both sides build their
//! own surface from it. The output is the status, the achieved tolerance and, when a
//! surface comes back, its type and its `SmBSplineSurface` part (degrees, control point
//! counts, knots with multiplicities, every control point and weight). Generated cases
//! cover every kind with and without wraps (offsets at positive, negative, zero and
//! nested distances and past a radius of curvature; curve bounded, STEP and volume
//! surfaces, alone and around or inside offsets), copies of B-spline kinds
//! (`bJustCopyBSplines`), starters of degree above 3 (the single-patch starter), simple
//! knots (the patch index check), loose, tight and default tolerances, subdivision
//! levels 0 to 4, at huge and tiny scales. Fixtures: `fixtures/surface.surface_approximate/`
//! (issue #1441's offset quarter cylinder among them).
// The generator indexes the nets it has just built (as other runners do).
#![allow(clippy::indexing_slicing)]
use std::borrow::Cow;
use serde::{Deserialize, Serialize};
use serde_json::Value;
use solidcore::curve::curve::Curve;
use solidcore::error::SolidError;
use solidcore::geom::extent2d::Extent2d;
use solidcore::geom::tol::approx_tol3d;
use solidcore::geom::transform::Transform;
use solidcore::geom::vector2d::Vec2;
use solidcore::geom::vector3d::Vec3;
use solidcore::kernel::Kernel;
use solidcore::nurbs::b_spline_curve::{BSplineCurve, ControlPointForm};
use solidcore::nurbs::b_spline_surface::{BSplineSurface, SurfParamType};
use solidcore::nurbs::crv_basic::Curve as NlCurve;
use solidcore::nurbs::math::{CPoint, NOW};
use solidcore::nurbs::srf_basic::Surface as NlSurface;
use solidcore::surface::cone::Cone;
use solidcore::surface::curve_bounded_surface::CurveBoundedSurface;
use solidcore::surface::offset_surface::OffsetSurface;
use solidcore::surface::plane::Plane;
use solidcore::surface::sphere::{Sphere, SphereOptions};
use solidcore::surface::srf_in_volume::SrfInVolume;
use solidcore::surface::surf_of_revolution::SurfOfRevolution;
use solidcore::surface::surface::{StepSurface, Surface, SurfaceKind};
use solidcore::surface::torus::{Torus, TorusOptions};
use solidcore::surface::volume::{Volume, VolumeKind};
use solidcore_ffi::units::surface_surface_approximate as ffi;
use crate::case::CaseContext;
use crate::runner::default_cases;
use crate::{Case, HarnessError, Output, Rng, UnitDiff};
/// Generated cases per run, whatever `--seeds` says (it can only raise this).
pub const MIN_GENERATED: usize = 400;
/// `SM_USE_DEFAULT`.
const USE_DEFAULT: f64 = 9_876_543.21;
/// `surface.surface_approximate` runner.
pub struct Runner;
/// The shim's input block.
#[derive(Debug, Clone, Serialize, Deserialize)]
struct Inputs {
block: Vec<f64>,
}
fn herr(e: impl std::fmt::Display) -> HarnessError {
HarnessError::new(e.to_string())
}
fn output(flat: Vec<f64>) -> Output {
Output::new().values("out", flat)
}
fn cpp(i: &Inputs) -> Result<Output, HarnessError> {
// each case in a fresh child: a native crash fails that case, not the unit run
match solidcore_ffi::isolate::fresh_process(|| ffi::eval(&i.block)) {
Ok(v) => Ok(output(v)),
Err(solidcore_ffi::FfiError::InvalidArgument(r)) => Ok(Output::rejected(r)),
Err(e) => Err(herr(e)),
}
}
// ---- Rust side ----------------------------------------------------------------------
struct Rd<'a> {
v: &'a [f64],
at: usize,
}
impl Rd<'_> {
fn num(&mut self) -> Result<f64, HarnessError> {
let x = *self.v.get(self.at).ok_or_else(|| herr("short block"))?;
self.at += 1;
Ok(x)
}
fn int(&mut self) -> Result<usize, HarnessError> {
Ok(self.num()? as usize)
}
fn flag(&mut self) -> Result<bool, HarnessError> {
Ok(self.num()? != 0.0)
}
fn vec(&mut self) -> Result<Vec3, HarnessError> {
Ok(Vec3::new(self.num()?, self.num()?, self.num()?))
}
fn dom(&mut self) -> Result<Extent2d, HarnessError> {
Ok(Extent2d::from_bounds(
self.num()?,
self.num()?,
self.num()?,
self.num()?,
))
}
fn point(&mut self, rational: bool) -> Result<CPoint, HarnessError> {
let (x, y, z, w) = (self.num()?, self.num()?, self.num()?, self.num()?);
Ok(CPoint::new(x, y, z, if rational { w } else { NOW }))
}
}
/// The surface of `kind` (the builders of `surface.surface_assert_valid`'s runner).
fn build(k: &Kernel, rd: &mut Rd<'_>, kind: usize) -> Result<Surface, HarnessError> {
if kind == 0 {
let (p, q, nu, nv) = (rd.int()?, rd.int()?, rd.int()?, rd.int()?);
let rational = rd.flag()?;
let ku = (0..=nu + p).map(|_| rd.num()).collect::<Result<Vec<_>, _>>()?;
let kv = (0..=nv + q).map(|_| rd.num()).collect::<Result<Vec<_>, _>>()?;
let mut net = Vec::new();
for _ in 0..nu {
net.push((0..nv).map(|_| rd.point(rational)).collect::<Result<Vec<_>, _>>()?);
}
let nurb = NlSurface::new(net, p, q, ku, kv).map_err(herr)?;
return Ok(Surface::from(BSplineSurface::from_nurb(nurb, false)));
}
let o = rd.vec()?;
let x = rd.vec()?;
let y = rd.vec()?;
Ok(match kind {
1 => {
let s = Vec2::new(rd.num()?, rd.num()?);
let d = rd.dom()?;
Surface::from(Plane::new(o, x, y, s, &d))
}
2 => {
let (r, s0, s1, h) = (rd.num()?, rd.num()?, rd.num()?, rd.num()?);
let (swap, inside, mngc) = (rd.flag()?, rd.flag()?, rd.flag()?);
Surface::from(Cone::new_cylinder_bounded(
o, x, y, r, s0, s1, h, swap, inside, mngc,
))
}
3 => {
let (bot, top, s0, s1, h) = (rd.num()?, rd.num()?, rd.num()?, rd.num()?, rd.num()?);
let (swap, inside, mngc) = (rd.flag()?, rd.flag()?, rd.flag()?);
Surface::from(Cone::new_bounded(
o, x, y, bot, top, s0, s1, h, swap, inside, mngc,
))
}
4 => {
let d = rd.dom()?;
let r = rd.num()?;
let (swap_uv, inside_out) = (rd.flag()?, rd.flag()?);
let opts = SphereOptions {
swap_uv,
inside_out,
circle_nurb: None,
};
Surface::from(Sphere::new(k, o, x, y, d, r, opts).map_err(herr)?)
}
5 => {
let d = rd.dom()?;
let (major, minor) = (rd.num()?, rd.num()?);
let (swap_uv, inside_out) = (rd.flag()?, rd.flag()?);
let opts = TorusOptions {
swap_uv,
inside_out,
circle_nurb: None,
};
Surface::from(Torus::new(k, o, x, y, d, major, minor, opts).map_err(herr)?)
}
_ => {
let d = rd.dom()?;
let swap = rd.flag()?;
let (p, n) = (rd.int()?, rd.int()?);
let rational = rd.flag()?;
let pts = (0..n).map(|_| rd.point(rational)).collect::<Result<Vec<_>, _>>()?;
let knots = (0..=n + p).map(|_| rd.num()).collect::<Result<Vec<_>, _>>()?;
let gen = Curve::from(BSplineCurve::from_nurb(
NlCurve::new(pts, p, knots).map_err(herr)?,
3,
false,
));
Surface::from(SurfOfRevolution::new(k, gen, o, x, y, d, swap).map_err(herr)?)
}
})
}
/// The `SmStatus` of a Rust error.
fn status(e: &SolidError) -> f64 {
match e {
SolidError::Upstream { code, .. } => f64::from(*code),
_ => -1.0,
}
}
/// The `SmBSplineSurface` part of a result (`SM_CAST_PTR` of the returned pointer: a
/// copy keeps its kind).
fn b_spline_part<'a>(
k: &Kernel,
s: &'a Surface,
) -> Result<Cow<'a, BSplineSurface>, SolidError> {
Ok(match s.kind_with_nurb(k)? {
SurfaceKind::BSpline(b) => Cow::Borrowed(b),
SurfaceKind::Plane(p) => p.base_with_nurb()?,
SurfaceKind::Cone(c) => Cow::Borrowed(c.rev().base()),
SurfaceKind::Sphere(c) => Cow::Borrowed(c.base().base()),
SurfaceKind::Torus(c) => Cow::Borrowed(c.base().base()),
SurfaceKind::SurfOfRevolution(c) => Cow::Borrowed(c.base()),
_ => return Err(SolidError::sm_err("difftest::surface_approximate::b_spline_part")),
})
}
/// The shim's `ap_put` record of `s`.
fn put(v: &mut Vec<f64>, s: &BSplineSurface) {
let nu = s.number_control_points(SurfParamType::U);
let nv = s.number_control_points(SurfParamType::V);
v.extend([
s.degree(SurfParamType::U) as f64,
s.degree(SurfParamType::V) as f64,
nu as f64,
nv as f64,
]);
for dir in [SurfParamType::U, SurfParamType::V] {
match s.knots(dir, None) {
Ok((knots, mults)) => {
v.push(knots.len() as f64);
v.extend(knots);
v.extend(mults.iter().map(|&m| m as f64));
}
Err(_) => v.push(-1.0),
}
}
for i in 0..nu {
for j in 0..nv {
match s.control_point(ControlPointForm::EuclideanRational, i, j) {
Ok((p, w)) => v.extend([p.x, p.y, p.z, w]),
Err(_) => v.extend([0.0, 0.0, 0.0, -1.0]),
}
}
}
}
fn rust(i: &Inputs) -> Result<Output, HarnessError> {
let k = Kernel::default();
let mut rd = Rd {
v: &i.block,
at: 0,
};
let tol = rd.num()?;
let level = rd.int()?;
let just_copy = rd.flag()?;
let n_wrap = rd.int()?;
let mut wraps = Vec::new();
for _ in 0..n_wrap {
let t = rd.int()?;
let m = match t {
0 | 1 => 1,
3 => 12,
_ => 0,
};
wraps.push((t, (0..m).map(|_| rd.num()).collect::<Result<Vec<_>, _>>()?));
}
let kind = rd.int()?;
let mut surface = build(&k, &mut rd, kind)?;
for (t, d) in wraps {
let at = |i: usize| d.get(i).copied().ok_or_else(|| herr("short wrap"));
surface = match t {
0 => Surface::from(OffsetSurface::new(&k, at(0)?, surface, true, 0).map_err(herr)?),
1 => Surface::from(
CurveBoundedSurface::new(surface, Vec::new(), Vec::new(), at(0)? != 0.0)
.map_err(herr)?,
),
2 => Surface::from(StepSurface::new(surface, true).map_err(herr)?),
_ => {
let v = |i: usize| -> Result<Vec3, HarnessError> {
Ok(Vec3::new(at(i)?, at(i + 1)?, at(i + 2)?))
};
let t = Transform::new(v(0)?, v(3)?, v(6)?, Some(v(9)?));
let vol = Volume::with_maps(VolumeKind::Transform(t), None, None).map_err(herr)?;
Surface::from(SrfInVolume::new(&k, surface, false, vol, 3, 0).map_err(herr)?)
}
};
}
let r = match surface.approximate_surface_ignoring_status(
&k,
approx_tol3d(&k, tol),
level,
just_copy,
) {
Ok(r) => r,
Err(e) => {
return Ok(Output::unported(&e).unwrap_or_else(|| Output::err(e.to_string())));
}
};
let mut v = vec![
r.failure.as_ref().map_or(0.0, status),
r.achieved_tol,
f64::from(u8::from(r.surface.is_some())),
];
if let Some(s) = &r.surface {
v.push(s.type_code() as f64);
match b_spline_part(&k, s) {
Ok(b) => put(&mut v, &b),
Err(e) => {
return Ok(Output::unported(&e).unwrap_or_else(|| Output::err(e.to_string())));
}
}
}
Ok(output(v))
}
// ---- generation ---------------------------------------------------------------------
fn scale(rng: &mut Rng) -> f64 {
match rng.below(8) {
0 => 10f64.powi(rng.range_i64(-6, -2) as i32),
1 => 10f64.powi(rng.range_i64(2, 6) as i32),
_ => rng.uniform(0.5, 4.0),
}
}
/// A unit-ish frame: the identity, a rotation about z or y, or a slightly skewed one.
fn frame(rng: &mut Rng, b: &mut Vec<f64>, s: f64) {
let o = [rng.uniform(-2.0, 2.0) * s, rng.uniform(-2.0, 2.0) * s, rng.uniform(-2.0, 2.0) * s];
let (x, y) = match rng.below(4) {
0 => ([1.0, 0.0, 0.0], [0.0, 1.0, 0.0]),
1 => {
let a = rng.uniform(0.0, 6.0);
([a.cos(), a.sin(), 0.0], [-a.sin(), a.cos(), 0.0])
}
2 => {
let a = rng.uniform(0.0, 6.0);
([a.cos(), 0.0, a.sin()], [0.0, 1.0, 0.0])
}
_ => ([1.0, 0.0, 0.0], [rng.uniform(-0.3, 0.3), 1.0, rng.uniform(-0.2, 0.2)]),
};
b.extend(o);
b.extend(x);
b.extend(y);
}
fn flag(rng: &mut Rng, p: f64) -> f64 {
f64::from(u8::from(rng.chance(p)))
}
fn angles(rng: &mut Rng) -> (f64, f64) {
match rng.below(5) {
0 => (0.0, 360.0),
1 => (-360.0, 0.0),
2 => {
let a = rng.uniform(-300.0, 200.0);
(a, a + rng.uniform(1.0, 160.0))
}
3 => (0.0, 90.0),
_ => (0.0, rng.uniform(10.0, 359.0)),
}
}
/// A clamped knot vector for `n` control points of degree `p`: simple interior knots,
/// or interior knots of full multiplicity `p` (piecewise Bezier).
fn clamped(rng: &mut Rng, n: usize, p: usize) -> Vec<f64> {
let interior = n - p - 1;
let mut mid: Vec<f64> = if interior > 0 && interior.is_multiple_of(p) && rng.chance(0.5) {
let spans = interior / p;
(0..spans)
.flat_map(|i| std::iter::repeat_n((i + 1) as f64 / (spans + 1) as f64, p))
.collect()
} else {
(0..interior)
.map(|i| (i + 1) as f64 / (interior + 1) as f64)
.collect()
};
mid.sort_by(f64::total_cmp);
let mut k = vec![0.0; p + 1];
k.extend(mid);
k.extend(vec![1.0; p + 1]);
k
}
fn gen_bspline(rng: &mut Rng, b: &mut Vec<f64>, s: f64) {
let deg = |rng: &mut Rng| match rng.below(8) {
0 => 4,
1 => 5,
_ => 1 + rng.below(3) as usize,
};
let p = deg(rng);
let q = deg(rng);
let nu = p + 1 + [0, p, 2 * p, 1, 2][rng.below(5) as usize];
let nv = q + 1 + [0, q, 2 * q, 1, 2][rng.below(5) as usize];
let rational = rng.chance(0.25);
let ku = clamped(rng, nu, p);
let kv = clamped(rng, nv, q);
let bump = [0.0, 0.3, 1.0][rng.below(3) as usize];
let pole = rng.chance(0.1);
b.extend([p as f64, q as f64, nu as f64, nv as f64, f64::from(u8::from(rational))]);
b.extend(ku);
b.extend(kv);
for i in 0..nu {
for j in 0..nv {
let u = i as f64 / (nu - 1) as f64;
let v = j as f64 / (nv - 1) as f64;
let w = if rational { rng.uniform(0.5, 2.0) } else { 1.0 };
let pt = if pole && i == 0 {
// a collapsed boundary row
[s, 0.5 * s, 0.0]
} else {
[(u + 1.0) * s, (v + 0.5) * s, bump * rng.uniform(-1.0, 1.0) * s]
};
if rational {
b.extend([pt[0] * w, pt[1] * w, pt[2] * w, w]);
} else {
b.extend([pt[0], pt[1], pt[2], 1.0]);
}
}
}
}
fn gen_case(rng: &mut Rng) -> Inputs {
let kind = rng.below(7) as usize;
let s = scale(rng);
let tol = match rng.below(8) {
0 => USE_DEFAULT,
1 => 1e-6 * s,
2 => 0.0,
3 => -1e-3 * s, // fabs
4 => 0.1 * s,
_ => 1e-3 * s,
};
let level = match rng.below(6) {
0 => 0,
1 => 1,
2 => 4,
_ => 3,
};
let mut b = vec![tol, f64::from(level), flag(rng, 0.15)];
let n_wrap = match rng.below(10) {
0..=2 => 0,
3 => 2,
4 => 3,
_ => 1,
};
b.push(f64::from(n_wrap));
for _ in 0..n_wrap {
match rng.below(10) {
0 => b.extend([1.0, flag(rng, 0.5)]), // curve bounded
1 => b.push(2.0), // STEP
2 => {
// in a volume: a rigid motion
let a = rng.uniform(0.0, 6.0);
let t = [rng.uniform(-1.0, 1.0) * s, rng.uniform(-1.0, 1.0) * s, 0.5 * s];
b.push(3.0);
b.extend(t);
b.extend([a.cos(), a.sin(), 0.0, -a.sin(), a.cos(), 0.0, 0.0, 0.0, 1.0]);
}
_ => {
b.push(0.0);
b.push(match rng.below(8) {
0 => 0.0,
1 => -rng.uniform(0.05, 0.5) * s,
2 => rng.uniform(1.0, 3.0) * s, // past a radius of curvature
_ => rng.uniform(0.05, 0.5) * s,
});
}
}
}
b.push(kind as f64);
if kind == 0 {
gen_bspline(rng, &mut b, s);
return Inputs { block: b };
}
frame(rng, &mut b, s);
match kind {
1 => {
let sx = if rng.chance(0.2) { -1.0 } else { rng.uniform(0.2, 3.0) };
let sy = rng.uniform(0.2, 3.0);
let u0 = rng.uniform(-2.0, 2.0) * s;
let v0 = rng.uniform(-2.0, 2.0) * s;
let du = rng.uniform(0.1, 3.0) * s;
let dv = rng.uniform(0.1, 3.0) * s;
b.extend([sx, sy, u0, v0, u0 + du, v0 + dv]);
}
2 => {
let (a0, a1) = angles(rng);
b.extend([s * rng.uniform(0.2, 2.0), a0, a1, s * rng.uniform(0.1, 3.0)]);
b.extend([flag(rng, 0.3), flag(rng, 0.3), flag(rng, 0.3)]);
}
3 => {
let (a0, a1) = angles(rng);
let (bot, top) = match rng.below(4) {
0 => (0.0, rng.uniform(0.5, 2.0)),
1 => (rng.uniform(0.5, 2.0), 0.0),
_ => (rng.uniform(0.2, 2.0), rng.uniform(0.2, 2.0)),
};
b.extend([bot * s, top * s, a0, a1, s * rng.uniform(0.1, 3.0)]);
b.extend([flag(rng, 0.3), flag(rng, 0.3), flag(rng, 0.3)]);
}
4 => {
let (a0, a1) = angles(rng);
let (v0, v1) = match rng.below(4) {
0 => (-90.0, 90.0),
1 => (-90.0, rng.uniform(-80.0, 89.0)),
2 => (rng.uniform(-89.0, 80.0), 90.0),
_ => {
let v = rng.uniform(-85.0, 60.0);
(v, v + rng.uniform(1.0, 85.0 - v.max(-85.0)).min(89.0 - v))
}
};
b.extend([a0, v0, a1, v1, s * rng.uniform(0.2, 3.0)]);
b.extend([flag(rng, 0.3), flag(rng, 0.3)]);
}
5 => {
let (a0, a1) = angles(rng);
let (v0, v1) = angles(rng);
let major = s * rng.uniform(0.5, 3.0);
let minor = match rng.below(4) {
0 => major * rng.uniform(1.1, 2.0), // spindle
_ => major * rng.uniform(0.05, 0.9),
};
b.extend([a0, v0, a1, v1, major, minor]);
b.extend([flag(rng, 0.3), flag(rng, 0.3)]);
}
_ => {
let (a0, a1) = angles(rng);
let p = 1 + rng.below(4) as usize;
let n = p + 1 + rng.below(3) as usize;
let interior = n - p - 1;
let mut knots = vec![0.0; p + 1];
knots.extend((0..interior).map(|i| (i + 1) as f64 / (interior + 1) as f64));
knots.extend(vec![1.0; p + 1]);
b.extend([a0, 0.0, a1, 1.0, flag(rng, 0.3)]);
let rational = rng.chance(0.2);
b.extend([p as f64, n as f64, f64::from(u8::from(rational))]);
for i in 0..n {
let t = i as f64 / (n - 1) as f64;
let (r, z) = (s * rng.uniform(0.3, 2.0), s * (t + rng.uniform(0.0, 0.1)));
let w = if rational { rng.uniform(0.5, 2.0) } else { 1.0 };
b.extend([r * w, 0.0, z * w, w]);
}
b.extend(knots);
}
}
Inputs { block: b }
}
impl UnitDiff for Runner {
fn id(&self) -> &'static str {
"surface.surface_approximate"
}
fn generate(&self, rng: &mut Rng) -> Option<Value> {
serde_json::to_value(gen_case(rng)).ok()
}
fn cases(&self, ctx: &CaseContext) -> Result<Vec<Case>, HarnessError> {
let ctx = CaseContext {
seeds: ctx.seeds.max(MIN_GENERATED),
..ctx.clone()
};
default_cases(self, &ctx)
}
fn run_cpp(&self, case: &Case) -> Result<Output, HarnessError> {
cpp(&case.parse()?)
}
fn run_rust(&self, case: &Case) -> Result<Output, HarnessError> {
rust(&case.parse()?)
}
}
crates/solidcore-ffi/src/units/surface_surface_approximate.rs
+266
View File
@@ -0,0 +1,266 @@
//! Safe wrapper for unit `surface.surface_approximate`: `SmSurface::ApproximateSurface`
//! on a real B-spline surface, plane, cylinder, cone, sphere, torus or surface of
//! revolution built from plain data, optionally inside up to three offset, curve bounded,
//! STEP or volume (rigid transform) surfaces (block layout in
//! `harness/cpp/shim/surface_surface_approximate.cpp`).
use crate::{call_with_buffer, FfiError};
mod sys {
extern "C" {
pub fn sc_surface_surface_approximate_eval(
input: *const f64,
n: i64,
out: *mut f64,
cap: i64,
len: *mut i64,
) -> i32;
}
}
/// Largest subdivision level accepted (each level can double the spans in each
/// direction, and every span is evaluated).
pub const MAX_LEVEL: usize = 4;
fn bad(what: &str) -> FfiError {
FfiError::InvalidArgument(format!("surface_approximate: {what}"))
}
struct Rd<'a> {
v: &'a [f64],
at: usize,
}
impl Rd<'_> {
fn num(&mut self) -> Result<f64, FfiError> {
let x = *self.v.get(self.at).ok_or_else(|| bad("short block"))?;
self.at += 1;
if !x.is_finite() || x.abs() > 1e12 {
return Err(bad("non-finite or huge value"));
}
Ok(x)
}
fn int(&mut self, lo: f64, hi: f64) -> Result<usize, FfiError> {
let x = self.num()?;
if x.fract() != 0.0 || x < lo || x > hi {
return Err(bad("bad count"));
}
Ok(x as usize)
}
fn flag(&mut self) -> Result<bool, FfiError> {
let x = self.num()?;
if x != 0.0 && x != 1.0 {
return Err(bad("bad flag"));
}
Ok(x == 1.0)
}
fn vec(&mut self) -> Result<[f64; 3], FfiError> {
Ok([self.num()?, self.num()?, self.num()?])
}
fn dom(&mut self) -> Result<[f64; 4], FfiError> {
Ok([self.num()?, self.num()?, self.num()?, self.num()?])
}
}
fn len(v: [f64; 3]) -> f64 {
(v[0] * v[0] + v[1] * v[1] + v[2] * v[2]).sqrt()
}
/// Two axes long enough and far from parallel (the constructors unitize them).
fn axes(x: [f64; 3], y: [f64; 3]) -> Result<(), FfiError> {
let c = [
x[1] * y[2] - x[2] * y[1],
x[2] * y[0] - x[0] * y[2],
x[0] * y[1] - x[1] * y[0],
];
if len(x) < 1e-6 || len(y) < 1e-6 || len(c) < 1e-3 * len(x) * len(y) {
return Err(bad("degenerate or parallel axes"));
}
Ok(())
}
/// An angular interval the analytic constructors accept.
fn angular(d0: f64, d1: f64) -> Result<(), FfiError> {
if d0 >= d1 || d0 < -360.0 || d1 > 360.0 || d1 - d0 > 360.0 {
return Err(bad("bad angular domain"));
}
Ok(())
}
/// Non-decreasing knots of positive span, multiplicities at most the order.
fn knots(k: &[f64], p: usize) -> Result<(), FfiError> {
if k.windows(2).any(|w| matches!(w, [a, b] if b < a)) || k.last() <= k.first() {
return Err(bad("bad knots"));
}
let mut run = 1;
for w in k.windows(2) {
run = if matches!(w, [a, b] if a == b) {
run + 1
} else {
1
};
if run > p + 1 {
return Err(bad("knot multiplicity above the order"));
}
}
Ok(())
}
/// `n` control points (`x y z w`), weights positive when rational.
fn points(rd: &mut Rd<'_>, n: usize, rational: bool) -> Result<(), FfiError> {
for _ in 0..n {
rd.vec()?;
let w = rd.num()?;
if rational && w <= 1e-6 {
return Err(bad("weight not positive"));
}
}
Ok(())
}
fn check(block: &[f64]) -> Result<(), FfiError> {
let mut rd = Rd { v: block, at: 0 };
rd.num()?; // tol
rd.int(0.0, MAX_LEVEL as f64)?;
rd.flag()?;
let n_wrap = rd.int(0.0, 3.0)?;
for _ in 0..n_wrap {
match rd.int(0.0, 3.0)? {
0 => {
rd.num()?;
}
1 => {
rd.flag()?;
}
2 => {}
_ => {
rd.vec()?;
let x = rd.vec()?;
let y = rd.vec()?;
let z = rd.vec()?;
axes(x, y)?;
axes(y, z)?;
axes(z, x)?;
}
}
}
let kind = rd.int(0.0, 6.0)?;
if kind == 0 {
let p = rd.int(1.0, 5.0)?;
let q = rd.int(1.0, 5.0)?;
let nu = rd.int(2.0, 32.0)?;
let nv = rd.int(2.0, 32.0)?;
if nu <= p || nv <= q {
return Err(bad("too few control points for the degree"));
}
let rational = rd.flag()?;
let ku: Vec<f64> = (0..=nu + p).map(|_| rd.num()).collect::<Result<_, _>>()?;
let kv: Vec<f64> = (0..=nv + q).map(|_| rd.num()).collect::<Result<_, _>>()?;
knots(&ku, p)?;
knots(&kv, q)?;
points(&mut rd, nu * nv, rational)?;
} else {
rd.vec()?;
let x = rd.vec()?;
let y = rd.vec()?;
axes(x, y)?;
match kind {
1 => {
let sx = rd.num()?;
let sy = rd.num()?;
let d = rd.dom()?;
if sx.abs() < 1e-6 || sy.abs() < 1e-6 || d[0] >= d[2] || d[1] >= d[3] {
return Err(bad("degenerate plane scale or domain"));
}
}
2 | 3 => {
let (r0, r1) = if kind == 2 {
let r = rd.num()?;
(r, r)
} else {
(rd.num()?, rd.num()?)
};
let s0 = rd.num()?;
let s1 = rd.num()?;
let h = rd.num()?;
angular(s0, s1)?;
if h <= 1e-6 || r0.abs().max(r1.abs()) < 1e-6 || (kind == 2 && r0 <= 1e-6) {
return Err(bad("degenerate height or radius"));
}
rd.flag()?;
rd.flag()?;
rd.flag()?;
}
4 => {
let d = rd.dom()?;
let r = rd.num()?;
angular(d[0], d[2])?;
if d[1] >= d[3] || d[1] < -90.0 || d[3] > 90.0 || r <= 1e-6 {
return Err(bad("sphere domain or radius"));
}
rd.flag()?;
rd.flag()?;
}
5 => {
let d = rd.dom()?;
let major = rd.num()?;
let minor = rd.num()?;
angular(d[0], d[2])?;
angular(d[1], d[3])?;
if major <= 1e-6 || minor <= 1e-6 {
return Err(bad("torus radii"));
}
rd.flag()?;
rd.flag()?;
}
_ => {
let d = rd.dom()?;
angular(d[0], d[2])?;
if d[1] >= d[3] {
return Err(bad("revolution domain"));
}
rd.flag()?;
let p = rd.int(1.0, 5.0)?;
let n = rd.int(2.0, 32.0)?;
if n <= p {
return Err(bad("too few generator points"));
}
let rational = rd.flag()?;
points(&mut rd, n, rational)?;
let k: Vec<f64> = (0..=n + p).map(|_| rd.num()).collect::<Result<_, _>>()?;
knots(&k, p)?;
}
}
}
if rd.at != block.len() {
return Err(bad("trailing values"));
}
Ok(())
}
/// Runs `SmSurface::ApproximateSurface(ctx, tol, achieved, level, result, just_copy)` on
/// the surface `block` describes and returns the shim's output record.
///
/// # Errors
///
/// [`FfiError::InvalidArgument`] for a block the constructors cannot take safely
/// (non-finite or huge values, degenerate axes, radii, heights, scales or domains,
/// angular domains outside the constructors' ranges, bad knots, non-positive weights,
/// degenerate transform axes, more than three wraps, a level above [`MAX_LEVEL`]); the
/// shim's status otherwise.
pub fn eval(block: &[f64]) -> Result<Vec<f64>, FfiError> {
check(block)?;
call_with_buffer(1 << 16, |out, len| {
// SAFETY: the validated block outlives the call; output and length pointers come
// from call_with_buffer.
unsafe {
sys::sc_surface_surface_approximate_eval(
block.as_ptr(),
block.len() as i64,
out.as_mut_ptr(),
out.len() as i64,
&raw mut *len,
)
}
})
}
crates/solidcore/src/surface/surface.rs
+1
View File
@@ -7,6 +7,7 @@
//! each lives in its own submodule and adds methods to [`Surface`] (or free
//! functions).
pub mod approximate;
pub mod assert_valid;
pub mod p1;
pub mod p10;
crates/solidcore/src/surface/surface/approximate.rs
+887
View File
@@ -0,0 +1,887 @@
//! Approximation of any surface by a B-spline surface within a tolerance
//! (`ApproximateSurface`): piecewise cubic Bezier patches from the surface's corner
//! data, split until every patch is close enough, then knots removed again.
//!
//! Upstream: `source/SMLib/src/SmSurface.cpp`, `source/SMLib/inc/SmSurface.h`,
//! `source/SMLib/src/SmBSplineSurfaceNLib.cpp`, `source/SMLib/inc/SmBSplineSurface.h`,
//! `source/SMLib/inc/SmCubicBezierSurface.h`, `source/SMLib/inc/SmOffsetSurface.h`,
//! `source/SMLib/inc/SmCurveBoundedSurface.h`, `source/SMLib/inc/SmSTEPSurface.h`,
//! `source/SMLib/inc/SmSrfInVolume.h`.
//!
//! | Rust | Upstream |
//! |-----------------------------------------|----------------------------------------------|
//! | [`Surface::approximate_surface`], [`ApproximatedSurface`] | `SmSurface::ApproximateSurface` |
//! | [`Surface::approximate_surface_ignoring_status`] | the same, for a caller that ignores its status |
//! | [`Surface::compute_cubic_bezier_approx`] | `SmSurface::ComputeCubicBezierApprox` |
//! | [`CubicBezierSurface`] | `SmCubicBezierSurface` (its points only) |
//! | `root_surface` (private) | `GetRootSurface` (virtual) |
//! | [`BSplineSurface::degree_elevate`] | `SmBSplineSurface::DegreeElevate` |
//! | [`BSplineSurface::insert_knots`] | `SmBSplineSurface::InsertKnots` (NLib form) |
//! | [`BSplineSurface::remove_knots`] | `SmBSplineSurface::RemoveKnots` |
//!
//! None of these is in `harness/callgraph/units.json`: they are listed in
//! `harness/callgraph/unused_functions.txt` (no upstream test or `SM_API` entry point
//! reaches them; `SmBrepData`'s `bWriteAsBSplines` writer and the graphics output call
//! `ApproximateSurface`).
//!
//! Deviations:
//! - **Outputs.** `rpNewBSplineSurface` and `dAchievedTol` are the fields of
//! [`ApproximatedSurface`]. The new surface is a [`Surface`]: with
//! `bJustCopyBSplines` the copy of a plane, cone, sphere, torus or swept surface is
//! that kind (upstream returns it through an `SmBSplineSurface*`); otherwise it is
//! always a B-spline. [`Surface::approximate_surface`] returns the status as `Err`
//! (the achieved tolerance is then lost); a caller that ignores the status, as every
//! upstream caller does, uses [`Surface::approximate_surface_ignoring_status`].
//! - **Patch index check.** When a starter surface's knots are not all of full
//! multiplicity, a patch's control points lie outside the net; upstream signals
//! `SM_ERR` but returns `SM_SUCCESS` without surface, as the port does (`Ok` with
//! `surface: None`).
//! - **Root surface.** `GetRootSurface` of an offset, curve bounded surface or surface
//! in a volume follows its base; the chain is walked in a loop bounded by
//! [`MAX_GEOMETRY_NESTING`] (`InvalidInput` beyond). The root of a plane, cone,
//! sphere, torus or swept surface is its `SmBSplineSurface` part, its NURB made as
//! `GetDegree` and the copy constructor make it (`GetOrCreateGwNurbPointer`).
//! - **Degree elevation below the degree.** [`BSplineSurface::degree_elevate`] returns
//! `InvalidInput` for a new degree below the current one, where upstream passes the
//! negative increment on to `N_SrfElevateDegree` (`ApproximateSurface` only elevates
//! degrees below 3, so it never asks for that).
//! - **Subdivision rounds.** `2 * lSubdivisionLevel + 1` wraps around in upstream's
//! `ULONG`; the port returns `InvalidInput` when it overflows `usize`.
//! - **Starter construction.** A failed `CreateCanonical` of the single-patch starter
//! (a degenerate natural domain) is returned; upstream dereferences the `NULL`
//! surface it leaves.
//! - **Ignored statuses.** As upstream: `Copy`, `GetKnots`, `GetControlPoint`,
//! `SetControlPoint`, `DegreeElevate` and `RemoveKnots` (a failed elevation or
//! removal leaves the surface unchanged); an unported dependency still stops the
//! port.
//! - **Memory and contexts.** `SmContext`, `SmNLibStackHandler` and `NL_STACKS` are
//! dropped; the debug dumps and drawing (`SM_DEBUG_CODE`) are omitted.
use std::borrow::Cow;
use crate::error::{ignored, Result, SolidError};
use crate::geom::extent1d::Extent1d;
use crate::geom::extent2d::Extent2d;
use crate::geom::vector2d::Point2;
use crate::geom::vector3d::{Point3, Vec3};
use crate::kernel::Kernel;
use crate::math::UNDEF_DOUBLE;
use crate::nurbs::b_spline_curve::{ControlPointForm, KnotType};
use crate::nurbs::b_spline_surface::{
BSplineSurface, BSplineSurfaceForm, CanonicalSurfaceData, SurfParamType,
};
use crate::nurbs::b_spline_surface_n_lib::p2::EvaluateFlags;
use crate::nurbs::knots::Side;
use crate::nurbs::srf_basic::Direction;
use crate::nurbs::srf_tool::p1::insert_knots;
use crate::nurbs::srf_tool::p2::elevate_degree;
use crate::nurbs::srf_tool::p3::remove_all_knots;
use crate::nurbs::srf_tool::p4::remove_knots_constrained;
use super::p2::{Surface, SurfaceKind, B_SPLINE_SURFACE_TYPE, MAX_GEOMETRY_NESTING};
const OP_APPROX: &str = "surface::surface::approximate::approximate_surface";
const OP_BEZIER: &str = "surface::surface::approximate::compute_cubic_bezier_approx";
const OP_ROOT: &str = "surface::surface::approximate::root_surface";
const OP_ELEVATE: &str = "surface::surface::approximate::degree_elevate";
const OP_INSERT: &str = "surface::surface::approximate::insert_knots";
const OP_REMOVE: &str = "surface::surface::approximate::remove_knots";
/// A bicubic Bezier patch (`SmCubicBezierSurface`): its 4x4 control points,
/// `points[v][u]` (`m_vP[v][u]`).
#[derive(Debug, Clone, Copy, PartialEq)]
pub struct CubicBezierSurface {
/// `m_vP`: row `v`, column `u`.
pub points: [[Point3; 4]; 4],
}
impl Default for CubicBezierSurface {
fn default() -> Self {
Self {
points: [[Point3::ZERO; 4]; 4],
}
}
}
/// The result of [`Surface::approximate_surface`].
#[derive(Debug, Clone)]
pub struct ApproximatedSurface {
/// `rpNewBSplineSurface`: the approximating surface with the fewest subdivisions
/// within the tolerance; `None` when the finest allowed one is still out of it
/// (and in the cases the module deviations list).
pub surface: Option<Surface>,
/// `dAchievedTol`: the largest sampled distance of the last test (`0` for a copied
/// B-spline, `-1` when nothing was tested).
pub achieved_tol: f64,
/// The status upstream returns and its caller ignores
/// ([`Surface::approximate_surface_ignoring_status`]); always `None` from
/// [`Surface::approximate_surface`], which returns it as `Err`.
pub failure: Option<SolidError>,
}
impl Surface {
/// A B-spline surface within `approx_tol_3d` (an `SmApproxTol3d` value:
/// [`crate::geom::tol::approx_tol3d`] converts a caller's `double`) of this one
/// (`ApproximateSurface`),
/// splitting the spans up to `subdivision_level` times in each direction. With
/// `just_copy_b_splines` a surface that is an `SmBSplineSurface` (a B-spline,
/// plane, cone, sphere, torus or swept surface) is copied instead.
///
/// # Errors
///
/// `SM_ERR` when no subdivision level is within the tolerance, or the starter
/// surface's degree is above 3; those of the evaluations, of the starter surface's
/// construction and of knot insertion; `InvalidInput` for a root surface chain
/// deeper than [`MAX_GEOMETRY_NESTING`].
pub fn approximate_surface(
&self,
k: &Kernel,
approx_tol_3d: f64,
subdivision_level: usize,
just_copy_b_splines: bool,
) -> Result<ApproximatedSurface> {
let mut achieved_tol = -1.0;
let surface = self.approximate_surface_impl(
k,
approx_tol_3d,
&mut achieved_tol,
subdivision_level,
just_copy_b_splines,
)?;
Ok(ApproximatedSurface {
surface,
achieved_tol,
failure: None,
})
}
/// [`Surface::approximate_surface`] for a caller that ignores its status, as
/// `SmBrepData`'s B-spline writer does: the outputs as the call leaves them (no
/// surface on failure; the achieved tolerance of the last test) and the failed
/// status in [`ApproximatedSurface::failure`].
///
/// # Errors
///
/// [`SolidError::Unported`] and [`SolidError::PreparationRequired`] only.
pub fn approximate_surface_ignoring_status(
&self,
k: &Kernel,
approx_tol_3d: f64,
subdivision_level: usize,
just_copy_b_splines: bool,
) -> Result<ApproximatedSurface> {
let mut achieved_tol = -1.0;
let (surface, failure) = match self.approximate_surface_impl(
k,
approx_tol_3d,
&mut achieved_tol,
subdivision_level,
just_copy_b_splines,
) {
Ok(surface) => (surface, None),
Err(e) if e.is_unported() || matches!(e, SolidError::PreparationRequired { .. }) => {
return Err(e)
}
Err(e) => (None, Some(e)),
};
Ok(ApproximatedSurface {
surface,
achieved_tol,
failure,
})
}
/// The body of `ApproximateSurface`, writing `dAchievedTol` as it goes.
fn approximate_surface_impl(
&self,
k: &Kernel,
approx_tol_3d: f64,
achieved_tol: &mut f64,
subdivision_level: usize,
just_copy_b_splines: bool,
) -> Result<Option<Surface>> {
const OP: &str = OP_APPROX;
// init output
*achieved_tol = -1.0;
// low work - when asked, just copy BSplines
if just_copy_b_splines && self.is_kind_of(B_SPLINE_SURFACE_TYPE)? {
let surface = ignored(self.copy(k))?;
*achieved_tol = 0.0;
return Ok(surface);
}
// arrive here when working with any surface not derived from SmBSplineSurface
// get RootSurface - not every surface type has one, and make the starter
// output surface
let mut approx_surf: Option<BSplineSurface> = None;
if let Some(root) = self.root_surface(k)? {
// let starter ApproxSurf equal a copy of the SmBSplineSurface part of the
// root surface, unless it has degree > 3 in either direction: our
// algorithm requires cubic.
let u_deg = root.degree(SurfParamType::U);
let v_deg = root.degree(SurfParamType::V);
if u_deg <= 3 && v_deg <= 3 {
approx_surf = Some(root.copy());
}
}
// If not a BSplineSurface of degree 3 or less:
let mut approx_surf = if let Some(s) = approx_surf {
s
} else {
// let starter ApproxSurf equal a degree 3x3 single patch with the same
// domain as the RootSurface
let uv_domain = self.natural_uv_domain_with_nurb(k)?;
let mut control_points = [Point3::ZERO; 16];
for (ii, row) in control_points.chunks_exact_mut(4).enumerate() {
for (jj, p) in row.iter_mut().enumerate() {
*p = Point3::new(ii as f64, jj as f64, 1.0);
}
}
let mults = [4, 4];
let u_starter_knots = [uv_domain.u_min(), uv_domain.u_max()];
let v_starter_knots = [uv_domain.v_min(), uv_domain.v_max()];
BSplineSurface::create_canonical(&CanonicalSurfaceData {
u_degree: 3,
v_degree: 3,
control_points: &control_points,
form: BSplineSurfaceForm::Unspecified,
u_multiplicities: &mults,
v_multiplicities: &mults,
u_knots: &u_starter_knots,
v_knots: &v_starter_knots,
knot_type: KnotType::UniformKnots,
weights: None,
uv_domain: None,
})?
};
// specify the control point form
let form = if approx_surf.is_rational() {
ControlPointForm::EuclideanRational
} else {
ControlPointForm::NonRational
};
// Nurb locals: upstream's arrays live across the iterations, so a failed
// `GetKnots` (no NURB) leaves the previous round's knots
let mut u_knots: Vec<f64> = Vec::new();
let mut v_knots: Vec<f64> = Vec::new();
let mut more_knots: Vec<f64> = Vec::new();
let mut testing_u_dir = false; // We'll toggle this each big iteration.
// Process u- and v-direction for each subD level, so iterate to 2x SubD level.
let last_sub_d = subdivision_level
.checked_mul(2)
.and_then(|n| n.checked_add(1))
.ok_or_else(|| SolidError::invalid_input(OP, "subdivision level too large"))?;
for _sub_d in 0..=last_sub_d {
testing_u_dir = !testing_u_dir; // Toggle
let mut u_deg = approx_surf.degree(SurfParamType::U);
let mut v_deg = approx_surf.degree(SurfParamType::V);
// no work - high degree surfaces
if u_deg > 3 || v_deg > 3 {
return Err(SolidError::sm_err(OP));
}
// elevate surface degrees to 3 (statuses ignored)
if u_deg < 3 {
ignored(approx_surf.degree_elevate(SurfParamType::U, 3))?;
}
if v_deg < 3 {
ignored(approx_surf.degree_elevate(SurfParamType::V, 3))?;
}
u_deg = approx_surf.degree(SurfParamType::U);
v_deg = approx_surf.degree(SurfParamType::V);
// Update pApproxSurf using cubic Bezier approximations to each patch.
// get knots of the subdivided surface (statuses ignored: without NURB the
// arrays keep their values)
if let Some((knots, _mult)) = ignored(approx_surf.knots(SurfParamType::U, None))? {
u_knots = knots;
}
if let Some((knots, _mult)) = ignored(approx_surf.knots(SurfParamType::V, None))? {
v_knots = knots;
}
// For each surface patch, compute the cubic Bezier approximation from the
// actual surface and copy those control Points into new approx surface.
for (iu, u_span) in u_knots.windows(2).enumerate() {
for (iv, v_span) in v_knots.windows(2).enumerate() {
let (&[u0, u1], &[v0, v1]) = (u_span, v_span) else {
continue;
};
let dom = Extent2d::from_points(Point2::new(u0, v0), Point2::new(u1, v1));
// build cubic bezier approximation to this curve
let cb = self.compute_cubic_bezier_approx(k, &dom)?;
let u_min = iu * 3;
let v_min = iv * 3;
let num_ctrl_u = approx_surf.number_control_points(SurfParamType::U);
let num_ctrl_v = approx_surf.number_control_points(SurfParamType::V);
// copy span approximation control points into new approx surface
for iub in 0..=3 {
for ivb in 0..=3 {
if u_min + iub >= num_ctrl_u || v_min + ivb >= num_ctrl_v {
// SE(SM_ERR): some problem here with knots; upstream
// returns SM_SUCCESS without surface
return Ok(None);
}
// (the point read is overwritten; statuses ignored)
let _ =
ignored(approx_surf.control_point(form, u_min + iub, v_min + ivb))?;
let weight = 1.0;
let p = cb
.points
.get(ivb)
.and_then(|row| row.get(iub))
.copied()
.ok_or_else(|| SolidError::sm_err(OP))?;
ignored(approx_surf.set_control_point(
form,
u_min + iub,
v_min + ivb,
p,
weight,
))?;
}
}
} // end iter every V patch
} // end iter every U patch
// Test the accuracy of this control polygon approx surface.
// If good enough then we can stop here.
// We'll test points at 1/4, 2/4, 3/4:
let num_between: usize = 3;
let fraction = (num_between + 1) as f64;
let mut failed = false;
if let Some((knots, _)) = ignored(approx_surf.knots(SurfParamType::U, None))? {
u_knots = knots;
}
if let Some((knots, _)) = ignored(approx_surf.knots(SurfParamType::V, None))? {
v_knots = knots;
}
let mut ivl_inner = Extent1d::new();
let mut ivl_outer = Extent1d::new();
// init max value
*achieved_tol = -1.0;
more_knots.clear();
let (outer_knots, inner_knots) = if testing_u_dir {
(&u_knots, &v_knots)
} else {
(&v_knots, &u_knots)
};
// for every surface patch
for outer_span in outer_knots.windows(2) {
failed = false; // Starting over on a new row.
let &[o0, o1] = outer_span else {
continue;
};
// (knots ascend: the status is ignored, as upstream)
let _ = ivl_outer.set_min_max(o0, o1);
// for every 'outer' dir sample point
for jj in 1..=num_between {
let outer_param = ivl_outer.evaluate(jj as f64 / fraction);
// for every 'inner' dir patch
for inner_span in inner_knots.windows(2) {
let &[i0, i1] = inner_span else {
continue;
};
let _ = ivl_inner.set_min_max(i0, i1);
// for every 'inner' dir sample point
for kk in 1..=num_between {
let inner_param = ivl_inner.evaluate(kk as f64 / fraction);
let uv = if testing_u_dir {
Point2::new(outer_param, inner_param)
} else {
Point2::new(inner_param, outer_param)
};
// Eval points on both surfaces and check distance.
let pnt = self.evaluate_point(k, uv)?;
let approx_pnt = approx_surf.evaluate_point(uv)?;
let dist = pnt.distance(approx_pnt);
// save max achieved tol
if dist > *achieved_tol {
*achieved_tol = dist;
}
// check gap for failure
if dist > approx_tol_3d {
// mark the approx surface as no good and quit
failed = true;
// Add a U or V knot between existing knots, multiplicity
// lUDeg or lVDeg, into the 'outer' direction knots.
let deg = if testing_u_dir { u_deg } else { v_deg };
for _ in 0..deg {
more_knots.push((o0 + o1) / 2.0);
}
break; // out of 'inner' dir samples into 'inner' dir patches
}
} // end iter every 'inner' dir sample point
if failed {
break; // out of 'inner' dir patches into 'outer' dir samples
}
} // end iter every 'inner' dir patch
if failed {
break; // out of 'outer' dir samples into 'outer' dir patches
}
} // end iter every 'outer' patch sample point
// ... continue on to the next 'outer' dir row/column of patches
} // end iter every 'outer' dir patch
// when the approx surface is good enough (upstream tests the flag of the
// last row or column only)
if !failed {
// Remove excess knots: the result does not have to be piecewise Bezier.
ignored(approx_surf.remove_knots(approx_tol_3d, true, true, None, None))?;
// all done
return Ok(Some(Surface::from(approx_surf)));
}
// arrive here when current subdivision level is not within tolerance -
// subdivide again
if more_knots.is_empty() {
// Can't happen? Infinite loop if it did.
return Err(SolidError::sm_err(OP));
}
let param = if testing_u_dir {
SurfParamType::U
} else {
SurfParamType::V
};
approx_surf.insert_knots(param, &more_knots)?;
} // end iter lSubD
// arrive here when approximation failed to converge
Err(SolidError::sm_err(OP))
}
/// The bicubic Bezier patch over `uv_domain` from this surface's positions, first
/// derivatives and twists at the domain's corners, scaled to the domain
/// (`ComputeCubicBezierApprox`).
///
/// # Errors
///
/// Those of the evaluations.
pub fn compute_cubic_bezier_approx(
&self,
k: &Kernel,
uv_domain: &Extent2d,
) -> Result<CubicBezierSurface> {
const OP: &str = OP_BEZIER;
// exact tangent values
let flags = EvaluateFlags {
u_side: Side::Left,
v_side: Side::Left,
only_upper_half: false,
non_zero_tangents: false,
};
let size = uv_domain.size();
// `SmVector3d sMat[2][2]` starts at SM_UNDEF_DOUBLE (its members' default
// initializers); an evaluation that writes only the upper half leaves the twist
// `sMat[1][1]` at that value for every corner
let writes_twist = evaluation_writes_twist(self)?;
let undef = Vec3::new(UNDEF_DOUBLE, UNDEF_DOUBLE, UNDEF_DOUBLE);
let corner = |uv: Point2| -> Result<[[Vec3; 2]; 2]> {
let m = self.evaluate(k, uv, 1, 1, flags)?;
let at = |i: usize, j: usize| -> Result<Vec3> {
m.get(i)
.and_then(|row| row.get(j))
.copied()
.ok_or_else(|| SolidError::sm_err(OP))
};
let twist = if writes_twist { at(1, 1)? } else { undef };
Ok([[at(0, 0)?, at(0, 1)?], [at(1, 0)?, twist]])
};
let mut cb = CubicBezierSurface::default();
let p = &mut cb.points;
// Now load control polygon at corners using hermite information.
{
let mat = corner(uv_domain.min())?;
let p00 = mat[0][0];
let du = mat[1][0] * size.x;
let dv = mat[0][1] * size.y;
let duv = mat[1][1] * size.x * size.y;
let p01 = p00 + du / 3.0;
let p10 = p00 + dv / 3.0;
let p11 = p10 + p01 - p00 + duv / 9.0;
p[0][0] = p00;
p[0][1] = p01;
p[1][0] = p10;
p[1][1] = p11;
}
{
let mat = corner(uv_domain.max())?;
let p00 = mat[0][0];
let du = mat[1][0] * size.x;
let dv = mat[0][1] * size.y;
let duv = mat[1][1] * size.x * size.y;
let p01 = p00 - du / 3.0;
let p10 = p00 - dv / 3.0;
let p11 = p10 + p01 - p00 + duv / 9.0;
p[3][3] = p00;
p[3][2] = p01;
p[2][3] = p10;
p[2][2] = p11;
}
{
let mat = corner(uv_domain.evaluate(1.0, 0.0))?;
let p00 = mat[0][0];
let du = mat[1][0] * size.x;
let dv = mat[0][1] * size.y;
let duv = mat[1][1] * size.x * size.y;
let p01 = p00 + dv / 3.0;
let p10 = p00 - du / 3.0;
let p11 = p10 + p01 - p00 - duv / 9.0;
p[0][3] = p00;
p[0][2] = p10;
p[1][3] = p01;
p[1][2] = p11;
}
{
let mat = corner(uv_domain.evaluate(0.0, 1.0))?;
let p00 = mat[0][0];
let du = mat[1][0] * size.x;
let dv = mat[0][1] * size.y;
let duv = mat[1][1] * size.x * size.y;
let p10 = p00 + du / 3.0;
let p01 = p00 - dv / 3.0;
let p11 = p10 + p01 - p00 - duv / 9.0;
p[3][0] = p00;
p[2][0] = p01;
p[3][1] = p10;
p[2][1] = p11;
}
Ok(cb)
}
/// `GetRootSurface`: the `SmBSplineSurface` part of this surface (a B-spline, plane,
/// cone, sphere, torus or swept surface, its NURB made when missing), followed
/// through offset, curve bounded and volume surfaces to their base; `None` for a
/// STEP surface.
///
/// # Errors
///
/// `InvalidInput` for a chain deeper than [`MAX_GEOMETRY_NESTING`];
/// [`SolidError::Unported`] from a NURB's make.
fn root_surface(&self, k: &Kernel) -> Result<Option<Cow<'_, BSplineSurface>>> {
let mut cur = self;
for _ in 0..=MAX_GEOMETRY_NESTING {
let next = match cur.kind() {
SurfaceKind::OffsetSurface(s) => s.base_surface(),
SurfaceKind::CurveBounded(s) => s.basis(),
SurfaceKind::SrfInVolume(s) => s.surface(),
SurfaceKind::Step(_) => return Ok(None),
SurfaceKind::BSpline(s) => return Ok(Some(Cow::Borrowed(s))),
SurfaceKind::Plane(s) => return Ok(Some(s.base_for_query())),
SurfaceKind::SurfOfExtrusion(_) => {
// a failed make leaves the surface without NURB (status ignored)
return match cur.extrusion_base_with_nurb(k) {
Ok(base) => Ok(Some(base)),
Err(e) if e.is_unported() => Err(e),
Err(_) => match cur.kind() {
SurfaceKind::SurfOfExtrusion(s) => Ok(Some(Cow::Borrowed(s.base()))),
_ => Ok(None),
},
};
}
SurfaceKind::Cone(_)
| SurfaceKind::Sphere(_)
| SurfaceKind::Torus(_)
| SurfaceKind::SurfOfRevolution(_) => {
let base = match cur.kind_with_nurb(k)? {
SurfaceKind::Cone(s) => s.rev().base(),
SurfaceKind::Sphere(s) => s.base().base(),
SurfaceKind::Torus(s) => s.base().base(),
SurfaceKind::SurfOfRevolution(s) => s.base(),
_ => return Err(SolidError::sm_err(OP_ROOT)),
};
return Ok(Some(Cow::Borrowed(base)));
}
};
cur = next;
}
Err(SolidError::invalid_input(
OP_ROOT,
"root surface chain deeper than MAX_GEOMETRY_NESTING",
))
}
}
/// Whether `Evaluate(uv, 1, 1, .., bOnlyUpperHalf = FALSE, ..)` of `s` writes the
/// twist `D[1][1]`. A curve bounded surface passes the call to its basis surface and a
/// STEP surface to its base's `EvaluateSTEP`; there a curve bounded surface passes it
/// along again, the cone, sphere, torus and surface of revolution overrides write the
/// upper half only (`D[1][1]` needs `min(u, v) >= 2`), the surface of extrusion
/// override zeroes the whole array first, and every other kind evaluates in full
/// (`EvaluateSTEP`'s default is `Evaluate`).
///
/// # Errors
///
/// `InvalidInput` for a chain deeper than [`MAX_GEOMETRY_NESTING`] or a STEP surface
/// without base surface.
fn evaluation_writes_twist(s: &Surface) -> Result<bool> {
const OP: &str = OP_BEZIER;
let mut cur = s;
let mut step = false;
for _ in 0..=MAX_GEOMETRY_NESTING {
cur = match cur.kind() {
SurfaceKind::CurveBounded(c) => c.basis(),
// (a STEP surface's own `EvaluateSTEP` is the default, `Evaluate`: its base's
// `EvaluateSTEP` again)
SurfaceKind::Step(st) => {
step = true;
st.surface()
.ok_or_else(|| SolidError::invalid_input(OP, "STEP surface without base"))?
}
SurfaceKind::Cone(_)
| SurfaceKind::Sphere(_)
| SurfaceKind::Torus(_)
| SurfaceKind::SurfOfRevolution(_)
if step =>
{
return Ok(false);
}
_ => return Ok(true),
};
}
Err(SolidError::invalid_input(
OP,
"surface chain deeper than MAX_GEOMETRY_NESTING",
))
}
impl BSplineSurface {
/// Raises the degree in `param`'s direction to `new_degree` (`DegreeElevate`).
///
/// # Errors
///
/// `NULL` pointer without NURB; `InvalidInput` for a new degree below the
/// current one; those of the NLib elevation (the surface is then unchanged).
pub fn degree_elevate(&mut self, param: SurfParamType, new_degree: usize) -> Result<()> {
const OP: &str = OP_ELEVATE;
let dir = if param == SurfParamType::V {
Direction::V
} else {
Direction::U
};
let increm = new_degree
.checked_sub(self.degree(param))
.ok_or_else(|| SolidError::invalid_input(OP, "new degree below the degree"))?;
let sur_p = self
.nurb
.as_ref()
.ok_or_else(|| SolidError::null_pointer(OP))?;
let sur = elevate_degree(sur_p, increm, dir).map_err(SolidError::nl_ser(OP))?;
// place output sur in a temporary and swap the NURBs
self.nurb = Self::from_gw_nurb(&sur).nurb;
Ok(())
}
/// Refines the knot vector in `param`'s direction with `new_knots`, an ascending
/// list (`InsertKnots`, the `SmTArray<double>&` form).
///
/// # Errors
///
/// `NULL` pointer without NURB; those of the NLib insertion (the surface is then
/// unchanged).
pub fn insert_knots(&mut self, param: SurfParamType, new_knots: &[f64]) -> Result<()> {
const OP: &str = OP_INSERT;
let dir = if param == SurfParamType::U {
Direction::U
} else {
Direction::V
};
let sur_p = self
.nurb
.as_ref()
.ok_or_else(|| SolidError::null_pointer(OP))?;
// Make sur with the union of surP target knot vector and the new knots
let sur = insert_knots(sur_p, new_knots, dir).map_err(SolidError::nl_ser(OP))?;
// build a surface from it and swap the NURBs
self.nurb = Self::from_gw_nurb(&sur).nurb;
Ok(())
}
/// Removes the knots whose removal keeps the surface within `approx_tol_3d`, in u
/// and/or v (both when neither is asked for), except the knots in `kept_u_knots` /
/// `kept_v_knots` (`RemoveKnots`).
///
/// # Errors
///
/// `NULL` pointer without NURB; those of the NLib removal (the surface is then
/// unchanged).
pub fn remove_knots(
&mut self,
approx_tol_3d: f64,
remove_u_knots: bool,
remove_v_knots: bool,
kept_u_knots: Option<&[f64]>,
kept_v_knots: Option<&[f64]>,
) -> Result<()> {
const OP: &str = OP_REMOVE;
let dir = if remove_u_knots && remove_v_knots {
Direction::UV
} else if remove_u_knots {
Direction::U
} else if remove_v_knots {
Direction::V
} else {
Direction::UV
};
let sur_p = self
.nurb
.as_ref()
.ok_or_else(|| SolidError::null_pointer(OP))?;
let sur = if kept_u_knots.is_some() || kept_v_knots.is_some() {
remove_knots_constrained(sur_p, kept_u_knots, kept_v_knots, approx_tol_3d, dir)
} else {
remove_all_knots(sur_p, approx_tol_3d, dir)
}
.map_err(SolidError::nl_ser(OP))?;
self.nurb = Self::from_gw_nurb(&sur).nurb;
Ok(())
}
}
#[cfg(test)]
mod tests {
use super::*;
use crate::surface::cone::Cone;
use crate::surface::offset_surface::OffsetSurface;
fn kernel() -> Kernel {
Kernel::default()
}
fn bilinear() -> BSplineSurface {
let pts = [
Point3::new(0.0, 0.0, 0.0),
Point3::new(0.0, 1.0, 0.2),
Point3::new(1.0, 0.0, 0.1),
Point3::new(1.0, 1.0, 0.0),
];
BSplineSurface::create_canonical(&CanonicalSurfaceData {
u_degree: 1,
v_degree: 1,
control_points: &pts,
form: BSplineSurfaceForm::Unspecified,
u_multiplicities: &[2, 2],
v_multiplicities: &[2, 2],
u_knots: &[0.0, 1.0],
v_knots: &[0.0, 1.0],
knot_type: KnotType::Unspecified,
weights: None,
uv_domain: None,
})
.unwrap()
}
#[test]
fn copies_b_splines_when_asked() {
let k = kernel();
let s = Surface::from(bilinear());
let r = s.approximate_surface(&k, 1e-3, 3, true).unwrap();
assert_eq!(r.achieved_tol, 0.0);
assert!(r.surface.is_some());
}
#[test]
fn approximates_a_b_spline_exactly() {
// a bilinear patch elevated to bicubic is reproduced by its Hermite data
let k = kernel();
let s = Surface::from(bilinear());
let r = s.approximate_surface(&k, 1e-6, 3, false).unwrap();
assert!(
r.achieved_tol >= 0.0 && r.achieved_tol < 1e-12,
"{}",
r.achieved_tol
);
assert!(r.surface.is_some());
}
#[test]
fn approximates_an_offset_quarter_cylinder() {
// issue #1441: the offset of a quarter cylinder of radius 5
let k = kernel();
let cyl = Surface::from(Cone::new_cylinder_bounded(
Vec3::ZERO,
Vec3::new(1.0, 0.0, 0.0),
Vec3::new(0.0, 1.0, 0.0),
5.0,
0.0,
90.0,
10.0,
false,
false,
false,
));
let off = Surface::from(OffsetSurface::new(&k, 2.0, cyl, true, 0).unwrap());
let tol = 1e-3;
let r = off.approximate_surface(&k, tol, 3, false).unwrap();
assert!(r.failure.is_none());
assert!(
r.achieved_tol >= 0.0 && r.achieved_tol <= tol,
"{}",
r.achieved_tol
);
let s = r.surface.unwrap();
let uv = off.natural_uv_domain_with_nurb(&k).unwrap().mid();
let p = s.evaluate_point(&k, uv).unwrap();
let q = off.evaluate_point(&k, uv).unwrap();
assert!(p.distance(q) <= 2.0 * tol, "{}", p.distance(q));
}
#[test]
fn degree_elevate_and_insert_keep_the_shape() {
let mut s = bilinear();
let before = s.evaluate_point(Point2::new(0.3, 0.7)).unwrap();
s.degree_elevate(SurfParamType::U, 3).unwrap();
s.degree_elevate(SurfParamType::V, 3).unwrap();
assert_eq!(s.degree(SurfParamType::U), 3);
s.insert_knots(SurfParamType::U, &[0.5, 0.5, 0.5]).unwrap();
assert_eq!(s.number_control_points(SurfParamType::U), 7);
let after = s.evaluate_point(Point2::new(0.3, 0.7)).unwrap();
assert!(before.distance(after) < 1e-12);
s.remove_knots(1e-9, true, true, None, None).unwrap();
assert_eq!(s.number_control_points(SurfParamType::U), 4);
}
}
crates/solidcore/src/surface/surface/p1.rs
+68 -2
View File
@@ -155,6 +155,7 @@ use crate::geom::vector3d::{Point3, Vec3};
use crate::intersect::surface_silhouette::SurfaceSilhouette;
use crate::intersect::surface_tracer::SurfaceTrace;
use crate::kernel::Kernel;
#[cfg(test)]
use crate::math::EFF_ZERO_DEG;
use crate::math::{max2, min2, BIG_DOUBLE, EFF_ZERO, EFF_ZERO_SQ, EFF_ZERO_SQRT, PI};
@@ -2119,13 +2120,66 @@ impl Surface {
/// surface's) the same way ([`OffsetSurface::natural_uv_domain_with_nurb`]), and so
/// does a lazy surface of extrusion (`SmSurfOfExtrusion` does not override
/// `GetNaturalUVDomain`: its made NURB's domain, whatever the generator's kind).
/// Every other kind answers as [`Surface::natural_uv_domain_if_built`].
/// A curve bounded surface and a surface in a volume pass the call to their basis /
/// base surface, and a STEP surface answers its base's STEP domain
/// ([`Surface::step_uv_domain_with_nurb`]), so a lazy analytic under them is made
/// as upstream's virtual calls make it. Every other kind answers as
/// [`Surface::natural_uv_domain_if_built`].
///
/// # Errors
///
/// [`SolidError::Unported`] from the make.
/// [`SolidError::Unported`] from the make; `InvalidInput` for a curve bounded, STEP
/// or volume chain deeper than [`MAX_GEOMETRY_NESTING`](crate::surface::surface::p2::MAX_GEOMETRY_NESTING).
pub fn natural_uv_domain_with_nurb(&self, k: &Kernel) -> Result<Extent2d> {
self.natural_uv_domain_with_nurb_at(k, 0)
}
/// [`Surface::step_uv_domain`] as upstream's `GetSTEPUVDomain`, whose base default
/// is the virtual `GetNaturalUVDomain` ([`Surface::natural_uv_domain_with_nurb`]):
/// the kinds with their own STEP domain answer it, a surface in a volume its base's,
/// the others their natural domain with a lazy NURB made.
///
/// # Errors
///
/// As [`Surface::natural_uv_domain_with_nurb`].
pub fn step_uv_domain_with_nurb(&self, k: &Kernel) -> Result<Extent2d> {
self.step_uv_domain_with_nurb_at(k, 0)
}
/// [`Surface::step_uv_domain_with_nurb`] at chain depth `depth`.
fn step_uv_domain_with_nurb_at(&self, k: &Kernel, depth: usize) -> Result<Extent2d> {
Ok(match self.kind() {
SurfaceKind::Plane(_)
| SurfaceKind::Cone(_)
| SurfaceKind::Sphere(_)
| SurfaceKind::Torus(_)
| SurfaceKind::SurfOfRevolution(_)
| SurfaceKind::SurfOfExtrusion(_) => self.step_uv_domain(),
SurfaceKind::SrfInVolume(s) => s.surface().step_uv_domain_with_nurb_at(
k,
nesting_depth(depth, "surface::surface::p1::step_uv_domain_with_nurb")?,
)?,
SurfaceKind::BSpline(_)
| SurfaceKind::CurveBounded(_)
| SurfaceKind::Step(_)
| SurfaceKind::OffsetSurface(_) => self.natural_uv_domain_with_nurb_at(k, depth)?,
})
}
/// [`Surface::natural_uv_domain_with_nurb`] at chain depth `depth`.
fn natural_uv_domain_with_nurb_at(&self, k: &Kernel, depth: usize) -> Result<Extent2d> {
const OP: &str = "surface::surface::p1::natural_uv_domain_with_nurb";
Ok(match self.kind_with_nurb(k)? {
SurfaceKind::CurveBounded(s) => s
.basis()
.natural_uv_domain_with_nurb_at(k, nesting_depth(depth, OP)?)?,
SurfaceKind::SrfInVolume(s) => s
.surface()
.natural_uv_domain_with_nurb_at(k, nesting_depth(depth, OP)?)?,
SurfaceKind::Step(s) => match s.surface() {
Some(base) => base.step_uv_domain_with_nurb_at(k, nesting_depth(depth, OP)?)?,
None => s.natural_uv_domain(),
},
SurfaceKind::Cone(s) => s.rev().base().natural_uv_domain(),
SurfaceKind::SurfOfRevolution(s) => s.base().natural_uv_domain(),
SurfaceKind::Sphere(s) => s.base().base().natural_uv_domain(),
@@ -4259,6 +4313,18 @@ impl<'a> SssIntersectEnfo<'a> {
}
}
/// The next depth of a curve bounded, STEP or volume chain walk, `InvalidInput` past
/// [`MAX_GEOMETRY_NESTING`](crate::surface::surface::p2::MAX_GEOMETRY_NESTING).
fn nesting_depth(depth: usize, op: &'static str) -> Result<usize> {
if depth >= crate::surface::surface::p2::MAX_GEOMETRY_NESTING {
return Err(SolidError::invalid_input(
op,
"surface chain deeper than MAX_GEOMETRY_NESTING",
));
}
Ok(depth + 1)
}
#[cfg(test)]
mod tests {
use super::*;
fixtures/surface.surface_approximate/cases.json
+1578
View File
File diff suppressed because it is too large Load Diff
harness/cpp/shim/surface_surface_approximate.cpp
+327
View File
@@ -0,0 +1,327 @@
// SPDX-License-Identifier: Apache-2.0
// Unit surface.surface_approximate. Builds a real SmBSplineSurface, SmPlane, SmCylinder,
// SmCone, SmSphere, SmTorus or SmSurfOfRevolution from plain data, optionally wraps it in
// up to three SmOffsetSurface, SmCurveBoundedSurface, SmSTEPSurface or SmSrfInVolume
// objects, and runs SmSurface::ApproximateSurface on it.
//
// Input block: `tol level just_copy n_wraps wraps... kind` then the kind's data, as in
// surface_surface_assert_valid.cpp:
// 0 B-spline: p q nu nv rational knots_u[nu+p+1] knots_v[nv+q+1] pts[nu*nv*4] (row i = u)
// 1 plane: O[3] X[3] Y[3] scale[2] domain[4]
// 2 cylinder: O X Y r start end h swap inside mngc
// 3 cone: O X Y bot top start end h swap inside mngc
// 4 sphere: O X Y domain[4] r swap inside
// 5 torus: O X Y domain[4] major minor swap inside
// 6 revolution: O X Y domain[4] swap, generator: p n rational pts[n*4] knots[n+p+1]
// Wraps are applied innermost first, each owning its base: `0 d` SmOffsetSurface(d, base,
// TRUE, 0); `1 implicit` SmCurveBoundedSurface without curves; `2` SmSTEPSurface(base,
// TRUE); `3 O[3] X[3] Y[3] Z[3]` SmSrfInVolume(base, FALSE, SmTransform(O, X, Y, &Z), 3, 0)
// (both parts copied; an original is freed only when it was copied: a curve bounded or
// STEP base has no Copy, so the object keeps and owns the original).
// `tol` is passed as an SmApproxTol3d (fabs, SM_USE_DEFAULT = half the zone tolerance).
//
// Output: `status achieved has` then, when a surface came back, its `type` and its
// SmBSplineSurface part: `deg_u deg_v n_u n_v`, per direction `count knots[count]
// mults[count]`, then the control points `x y z w` (SM_CP_EUCLIDIAN_RATIONAL, u major).
#include "sc_shim.h"
#include "StdAfx.h"
#include <SmAxis2Placement.h>
#include <SmBSplineCurve.h>
#include <SmBSplineSurface.h>
#include <SmCone.h>
#include <SmContext.h>
#include <SmCurve.h>
#include <SmCurveBoundedSurface.h>
#include <SmCylinder.h>
#include <SmOffsetSurface.h>
#include <SmPlane.h>
#include <SmSTEPSurface.h>
#include <SmSphere.h>
#include <SmSrfInVolume.h>
#include <SmSurfOfRevolution.h>
#include <SmTArray.h>
#include <SmTol.h>
#include <SmTorus.h>
#include <SmTransform.h>
#include <nurbs.h>
#include <vector>
namespace
{
struct ApConePeek : SmCone
{
ApConePeek(const SmPoint3d& o, const SmVector3d& x, const SmVector3d& y, double bot, double top,
double start, double end, double h, SmBoolean swap, SmBoolean inside, SmBoolean mngc)
: SmCone(o, x, y, bot, top, start, end, h, swap, inside, mngc)
{
}
};
struct ApRevolution : SmSurfOfRevolution
{
ApRevolution(SmBSplineCurve* gen, const SmPoint3d& o, const SmVector3d& x, const SmVector3d& y,
const SmExtent2d& d, SmBoolean swap)
: SmSurfOfRevolution(gen, o, x, y, d, swap)
{
}
};
struct ApRd
{
const double* in;
int64_t n;
int64_t at;
bool ok = true;
double num()
{
if (at >= n)
{
ok = false;
return 0.0;
}
return in[at++];
}
int64_t idx() { return (int64_t)num(); }
SmBoolean b() { return num() != 0.0 ? TRUE : FALSE; }
SmVector3d vec()
{
double x = num();
double y = num();
double z = num();
return SmVector3d(x, y, z);
}
SmExtent2d dom()
{
double a = num();
double b2 = num();
double c = num();
double d = num();
return SmExtent2d(a, b2, c, d);
}
};
gw_CPOINT ap_point(ApRd& rd, bool rational)
{
gw_CPOINT p;
p.x = rd.num();
p.y = rd.num();
p.z = rd.num();
double w = rd.num();
p.w = rational ? w : NL_NOW;
return p;
}
// The surface of kind `kind` read from `rd`; NULL on a bad block.
SmSurface* ap_build(ApRd& rd, SmContext& ctx, int64_t kind)
{
if (kind == 0)
{
int64_t p = rd.idx(), q = rd.idx(), nu = rd.idx(), nv = rd.idx();
bool rational = rd.num() != 0.0;
if (p < 1 || q < 1 || nu < 2 || nv < 2 || nu > 64 || nv > 64)
return NULL;
std::vector<double> ku(nu + p + 1), kv(nv + q + 1);
for (auto& v : ku)
v = rd.num();
for (auto& v : kv)
v = rd.num();
std::vector<gw_CPOINT> pts(nu * nv);
for (auto& c : pts)
c = ap_point(rd, rational);
if (!rd.ok)
return NULL;
std::vector<gw_CPOINT*> rows(nu);
for (int64_t i = 0; i < nu; i++)
rows[i] = &pts[i * nv];
return new (ctx) SmBSplineSurface((ULONG)(nu - 1), (ULONG)(nv - 1), (short)p, (short)q,
(ULONG)(nu + p), (ULONG)(nv + q), rows.data(), ku.data(),
kv.data());
}
if (kind < 1 || kind > 6)
return NULL;
SmPoint3d o = rd.vec();
SmVector3d x = rd.vec();
SmVector3d y = rd.vec();
if (kind == 1)
{
double sx = rd.num(), sy = rd.num();
SmExtent2d d = rd.dom();
if (!rd.ok)
return NULL;
return new (ctx) SmPlane(o, x, y, SmVector2d(sx, sy), d, &ctx);
}
if (kind == 2)
{
double r = rd.num(), s0 = rd.num(), s1 = rd.num(), h = rd.num();
SmBoolean swap = rd.b(), inside = rd.b(), mngc = rd.b();
if (!rd.ok)
return NULL;
return new (ctx) SmCylinder(o, x, y, r, s0, s1, h, swap, inside, mngc, &ctx);
}
if (kind == 3)
{
double bot = rd.num(), top = rd.num(), s0 = rd.num(), s1 = rd.num(), h = rd.num();
SmBoolean swap = rd.b(), inside = rd.b(), mngc = rd.b();
if (!rd.ok)
return NULL;
return new (ctx) ApConePeek(o, x, y, bot, top, s0, s1, h, swap, inside, mngc);
}
if (kind == 4)
{
SmExtent2d d = rd.dom();
double r = rd.num();
SmBoolean swap = rd.b(), inside = rd.b();
if (!rd.ok)
return NULL;
return new (ctx) SmSphere(o, x, y, d, r, swap, inside, NULL, &ctx);
}
if (kind == 5)
{
SmExtent2d d = rd.dom();
double major = rd.num(), minor = rd.num();
SmBoolean swap = rd.b(), inside = rd.b();
if (!rd.ok)
return NULL;
return new (ctx) SmTorus(o, x, y, d, major, minor, swap, inside, NULL, &ctx);
}
SmExtent2d d = rd.dom();
SmBoolean swap = rd.b();
int64_t p = rd.idx(), cn = rd.idx();
bool rational = rd.num() != 0.0;
if (p < 1 || cn < 2 || cn > 64)
return NULL;
std::vector<gw_CPOINT> pts(cn);
for (auto& c : pts)
c = ap_point(rd, rational);
std::vector<double> k(cn + p + 1);
for (auto& v : k)
v = rd.num();
if (!rd.ok)
return NULL;
SmBSplineCurve* gen = new (ctx)
SmBSplineCurve(3, (ULONG)(cn - 1), (ULONG)(cn + p), (short)p, pts.data(), k.data());
return new (ctx) ApRevolution(gen, o, x, y, d, swap);
}
// The SmBSplineSurface part of `s`: degrees, counts, knots and control points.
void ap_put(std::vector<double>& v, SmBSplineSurface* s)
{
ULONG nu = s->SmBSplineSurface::GetNumberControlPoints(SM_SP_U);
ULONG nv = s->SmBSplineSurface::GetNumberControlPoints(SM_SP_V);
v.push_back((double)s->SmBSplineSurface::GetDegree(SM_SP_U));
v.push_back((double)s->SmBSplineSurface::GetDegree(SM_SP_V));
v.push_back((double)nu);
v.push_back((double)nv);
for (int dir = 0; dir < 2; dir++)
{
SmTArray<double> k;
SmTArray<ULONG> m;
SmStatus st = s->SmBSplineSurface::GetKnots(dir == 0 ? SM_SP_U : SM_SP_V, k, &m);
if (st != SM_SUCCESS)
{
v.push_back(-1.0);
continue;
}
v.push_back((double)k.GetSize());
for (ULONG i = 0; i < k.GetSize(); i++)
v.push_back(k[i]);
for (ULONG i = 0; i < m.GetSize(); i++)
v.push_back((double)m[i]);
}
for (ULONG i = 0; i < nu; i++)
for (ULONG j = 0; j < nv; j++)
{
SmPoint3d p;
double w = 0.0;
SmStatus st = s->GetControlPoint(SM_CP_EUCLIDIAN_RATIONAL, i, j, p, w);
if (st != SM_SUCCESS)
p = SmPoint3d(0.0, 0.0, 0.0), w = -1.0;
v.insert(v.end(), {p.x, p.y, p.z, w});
}
}
} // namespace
SC_EXPORT int32_t sc_surface_surface_approximate_eval(const double* input, int64_t n, double* out,
int64_t cap, int64_t* out_len)
{
SC_TRY
SmContext ctx;
ApRd rd{input, n, 0};
double tol = rd.num();
int64_t level = rd.idx();
bool just_copy = rd.num() != 0.0;
int64_t n_wrap = rd.idx();
if (!rd.ok || level < 0 || n_wrap < 0 || n_wrap > 3)
return SC_ERR_ARG;
struct Wrap
{
int64_t type;
double d[12];
};
std::vector<Wrap> wraps((size_t)n_wrap);
for (auto& w : wraps)
{
w.type = rd.idx();
int64_t m = w.type == 0 || w.type == 1 ? 1 : w.type == 3 ? 12 : 0;
if (w.type < 0 || w.type > 3)
return SC_ERR_ARG;
for (int64_t i = 0; i < m; i++)
w.d[i] = rd.num();
}
int64_t kind = rd.idx();
SmSurface* surface = ap_build(rd, ctx, kind);
if (surface == NULL)
return SC_ERR_ARG;
if (!rd.ok || rd.at != n)
{
delete surface;
return SC_ERR_ARG;
}
for (const Wrap& w : wraps)
{
if (w.type == 0)
surface = new (ctx) SmOffsetSurface(w.d[0], *surface, TRUE, 0);
else if (w.type == 1)
{
SmTArray<SmCurve*> c3d, cuv;
SmBoolean implicit = w.d[0] != 0.0 ? TRUE : FALSE;
surface = new (ctx) SmCurveBoundedSurface(*surface, &c3d, &cuv, implicit);
}
else if (w.type == 2)
surface = new (ctx) SmSTEPSurface(*surface, TRUE);
else
{
SmVector3d o(w.d[0], w.d[1], w.d[2]), x(w.d[3], w.d[4], w.d[5]);
SmVector3d y(w.d[6], w.d[7], w.d[8]), z(w.d[9], w.d[10], w.d[11]);
SmTransform* t = new (ctx) SmTransform(&ctx, o, x, y, &z);
SmSurface* base = surface;
SmSrfInVolume* siv = new (ctx) SmSrfInVolume(*base, FALSE, *t, 3, 0, &ctx);
// Copy flag 3 makes the object own (and later delete) what it points to. A
// failed Copy (SmSurface's default for curve bounded and STEP surfaces) leaves
// it pointing to the original: free only the originals it replaced by copies.
if (siv->GetSurface() != base)
delete base;
if (siv->GetVolume() != t)
delete t;
surface = siv;
}
}
SmBSplineSurface* result = NULL;
double achieved = 0.0;
SmStatus st = surface->ApproximateSurface(ctx, SmApproxTol3d(tol), achieved, (ULONG)level,
result, just_copy ? TRUE : FALSE);
std::vector<double> v;
v.push_back(st == SM_SUCCESS ? 0.0 : (double)st);
v.push_back(achieved);
v.push_back(result != NULL ? 1.0 : 0.0);
if (result != NULL)
{
v.push_back((double)result->GetType());
ap_put(v, result);
delete result;
}
delete surface;
return sc_copy_out(v.data(), (int64_t)v.size(), out, cap, out_len);
SC_CATCH
}
plan/units/surface.surface_approximate.toml
+12
View File
@@ -0,0 +1,12 @@
id = "surface.surface_approximate"
upstream = ["source/SMLib/src/SmSurface.cpp", "source/SMLib/src/SmBSplineSurfaceNLib.cpp", "source/SMLib/inc/SmCubicBezierSurface.h"]
rust_module = "surface::surface::approximate"
depends_on = ["surface.surface_p1", "surface.surface_p2", "surface.offset_surface", "nurbs.b_spline_surface_p1", "nurbs.b_spline_surface_p2", "nurbs.b_spline_surface_n_lib_p1", "nurbs.b_spline_surface_n_lib_p2", "nlib.srf_tool_p1", "nlib.srf_tool_p2", "nlib.srf_tool_p3", "nlib.srf_tool_p4"]
soft_depends_on = []
est_lines = 640
layer = 32
tests = []
status = "green"
claimed_by = ""
claimed_at = ""
notes = "PR #1533, issue #1441 (offset surfaces cannot be exported; solidcore 1.0.1 parity). NOT IN callgraph/units.json: SmSurface::ApproximateSurface and its own callees are in harness/callgraph/unused_functions.txt (no upstream test or SM_API entry point reaches them). Owns SmSurface::ApproximateSurface (SmSurface.cpp:877), SmSurface::ComputeCubicBezierApprox (SmSurface.cpp:8600), SmCubicBezierSurface (its points only), the GetRootSurface overrides (SmSurface, SmBSplineSurface, SmOffsetSurface, SmCurveBoundedSurface, SmSTEPSurface, SmSrfInVolume), and SmBSplineSurface::DegreeElevate, InsertKnots (SmTArray form, SmBSplineSurfaceNLib.cpp:6490) and RemoveKnots. Rust: Surface::approximate_surface / approximate_surface_ignoring_status / compute_cubic_bezier_approx, BSplineSurface::degree_elevate / insert_knots / remove_knots. Not ported (follow-ups): the callers SmBrepData's bWriteAsBSplines writer (brep_data_write writes native formats only), SmFace::ApproximateWithBSpline, SmFace::ApproximateSurface, SmSurface::ApproximateSurfaceLeastSqs, the graphics output; SmCubicBezierSurface's other members. Status pending the gate run."