diff --git a/crates/kernel/src/build.rs b/crates/kernel/src/build.rs index 1b94a447..9078e4ab 100644 --- a/crates/kernel/src/build.rs +++ b/crates/kernel/src/build.rs @@ -2063,6 +2063,15 @@ pub fn sphere_patches_rotated(center: Vec3, radius: f64, rot: [Vec3; 3]) -> Resu /// [`sphere_patches_rotated`], every patch swept from profile segments `arcs`. fn sphere_solid(center: Vec3, radius: f64, rot: [Vec3; 3], arcs: &[usize]) -> Result { + guard("sphere", || { + let shell: mt::Shell = sphere_patch_faces(center, radius, rot)?.into(); + let solid = Solid::try_new(vec![shell]).map_err(|e| KernelError::Failed(format!("sphere: {e}")))?; + record::attach(Body::new(solid)?, &[], Some(vec![record::profile(arcs); 6]), Default::default()) + }) +} + +/// The six outward patch faces of [`sphere_patches_rotated`]; truck calls, so callers guard them. +pub(crate) fn sphere_patch_faces(center: Vec3, radius: f64, rot: [Vec3; 3]) -> Result> { 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() { @@ -2080,50 +2089,46 @@ fn sphere_solid(center: Vec3, radius: f64, rot: [Vec3; 3], arcs: &[usize]) -> Re (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]]), ]; - guard("sphere", || { - 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 = 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 = Vec::new(); - for (dir, axis, cs) in faces { - let (dir, axis) = (turn(dir), turn(axis)); - let mut wire: Vec = 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); + 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 = 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 = Vec::new(); + for (dir, axis, cs) in faces { + let (dir, axis) = (turn(dir), turn(axis)); + let mut wire: Vec = 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 } - } - out.push(mt::Face::new(vec![wire.into()], surface)); + }; + wire.push(if ka < kb { e } else { e.inverse() }); } - let shell: mt::Shell = out.into(); - let solid = Solid::try_new(vec![shell]).map_err(|e| KernelError::Failed(format!("sphere: {e}")))?; - record::attach(Body::new(solid)?, &[], Some(vec![record::profile(arcs); 6]), Default::default()) - }) + // 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) } diff --git a/crates/kernel/src/step_in/brep.rs b/crates/kernel/src/step_in/brep.rs index cd40bb3b..2ae1dbb1 100644 --- a/crates/kernel/src/step_in/brep.rs +++ b/crates/kernel/src/step_in/brep.rs @@ -339,11 +339,77 @@ impl Builder { return Err("no face could be imported".into()); } let read = out.len(); + let ball = match (out.as_slice(), origin.as_slice()) { + ([face], [file_face]) => whole_sphere(cx, *file_face, face).map(|patches| (*file_face, patches)), + _ => None, + }; + if let Some((file_face, patches)) = ball { + origin = vec![file_face; patches.len()]; + out = patches; + } let out = super::split::split_bands(out, &mut origin, failed); Ok((out.into(), closed, read, origin)) } } +/// A ball written as one face on a sphere bounded only by its seam, run out and back from pole +/// to pole. Truck meshes a face on a revolved surface in parameter space, where that loop +/// encloses nothing, so the ball is built as the kernel's own sphere instead: six patches with +/// no pole or seam inside them (`build::sphere_patch_faces`), turned inside out for a void. +/// `None` for any other face. +fn whole_sphere(cx: &Ctx, file_face: u64, face: &mt::Face) -> Option> { + use mt::{ParametricSurface3D, SearchParameter}; + let surface = entity(cx, file_face).ok()?.params().get(2).and_then(Param::as_ref_id)?; + if entity(cx, surface).ok()?.name() != "SPHERICAL_SURFACE" { + return None; + } + let wires = face.boundaries(); + let edges: Vec<&mt::Edge> = wires.iter().flat_map(|w| w.edge_iter()).collect(); + let mut uses: HashMap = HashMap::new(); + let mut poles: Vec = Vec::new(); + for e in &edges { + *uses.entry(e.id()).or_default() += 1; + for v in [e.front(), e.back()] { + if !poles.iter().any(|p| p.id() == v.id()) { + poles.push(v.clone()); + } + } + } + // Every edge a seam (used once each way), and the two vertices it joins are the poles. + if uses.is_empty() || uses.values().any(|n| *n != 2) { + return None; + } + let [a, b] = poles.as_slice() else { return None }; + let (a, b) = (a.point(), b.point()); + let center = a + (b - a) * 0.5; + let radius = a.distance(b) * 0.5; + let tol = cx.tol.max(radius * 1e-6); + // The seam is a meridian: all of it on the sphere whose diameter joins the poles. + let samples = |e: &&mt::Edge| { + let c = e.curve(); + let (t0, t1) = c.range_tuple(); + (0..=8).map(move |i| c.subs(t0 + (t1 - t0) * i as f64 / 8.0)) + }; + if !(radius.is_finite() && radius > tol) || !edges.iter().flat_map(samples).all(|p| (p.distance(center) - radius).abs() <= tol) { + return None; + } + // Ball or void: the face's normal at the middle of the seam, against the way out of the centre. + let probe = edges.first().and_then(|e| samples(e).nth(4))?; + let oriented = face.oriented_surface(); + let (u, v) = oriented.search_parameter(probe, mt::SPHint2D::None, 100)?; + let outward = oriented.normal(u, v).dot(probe - center) > 0.0; + let z = (crate::body::from_p3(a) - crate::body::from_p3(b)).normalized()?; + let x = z.any_perp(); + let c = crate::body::from_p3(center); + let mut patches = guard("step sphere", || crate::build::sphere_patch_faces(c, radius, [x, z.cross(x), z])).ok()?; + if !outward { + for f in &mut patches { + f.invert(); + } + } + Some(patches) +} + /// Read a solid-like representation item: MANIFOLD_SOLID_BREP, BREP_WITH_VOIDS, FACETED_BREP, /// or a SHELL_BASED_SURFACE_MODEL (closed shells give solids, open ones open bodies). pub(crate) fn read_item(cx: &Ctx, id: u64) -> R { diff --git a/crates/kernel/src/step_in/tests.rs b/crates/kernel/src/step_in/tests.rs index 8a195af6..41d841ba 100644 --- a/crates/kernel/src/step_in/tests.rs +++ b/crates/kernel/src/step_in/tests.rs @@ -7,6 +7,117 @@ fn rel(a: f64, b: f64) -> f64 { (a - b).abs() / b.abs().max(1e-9) } +/// A single spherical face bounded by one meridian edge used in both directions. +/// `end_angle` is measured from the north pole toward the seam; PI is the south pole. +fn seam_sphere_step(center: Vec3, radius: f64, axis: Vec3, reference: Vec3, azimuth: f64, reverse_edge: bool, end_angle: f64) -> String { + let seam = reference * azimuth.cos() + axis.cross(reference) * azimuth.sin(); + let end = center + (axis * end_angle.cos() + seam * end_angle.sin()) * radius; + let north = center + axis * radius; + let circle_axis = axis.cross(seam); + let point = |p: Vec3| format!("({:.16},{:.16},{:.16})", p.x, p.y, p.z); + let (start, finish, sense, first, second) = if reverse_edge { (9, 7, "F", 15, 14) } else { (7, 9, "T", 14, 15) }; + format!( + "{HEAD}#1=CARTESIAN_POINT('',{});\n#2=DIRECTION('',{});\n#3=DIRECTION('',{});\n#4=AXIS2_PLACEMENT_3D('',#1,#2,#3);\n#5=SPHERICAL_SURFACE('',#4,{radius:.16});\n#6=CARTESIAN_POINT('',{});\n#7=VERTEX_POINT('',#6);\n#8=CARTESIAN_POINT('',{});\n#9=VERTEX_POINT('',#8);\n#10=DIRECTION('',{});\n#11=AXIS2_PLACEMENT_3D('',#1,#10,#12);\n#12=DIRECTION('',{});\n#13=CIRCLE('',#11,{radius:.16});\n#16=EDGE_CURVE('',#{start},#{finish},#13,.{sense}.);\n#14=ORIENTED_EDGE('',*,*,#16,.T.);\n#15=ORIENTED_EDGE('',*,*,#16,.F.);\n#17=EDGE_LOOP('',(#{first},#{second}));\n#18=FACE_OUTER_BOUND('',#17,.T.);\n#19=ADVANCED_FACE('',(#18),#5,.T.);\n#20=CLOSED_SHELL('',(#19));\n#21=MANIFOLD_SOLID_BREP('Ball',#20);\n#22=SHAPE_REPRESENTATION('',(#21),$);\n{TAIL}", + point(center), + point(axis), + point(reference), + point(north), + point(end), + point(circle_axis), + point(seam), + ) +} + +fn assert_seam_sphere(text: &str, radius: f64) { + let imp = step_import(text).unwrap(); + assert_eq!(imp.bodies.len(), 1, "{:?}", imp.warnings); + assert!(imp.warnings.is_empty(), "{:?}", imp.warnings); + let body = &imp.bodies[0]; + assert!(body.closed); + let m = measure(&body.body).unwrap(); + let want_volume = 4.0 / 3.0 * std::f64::consts::PI * radius.powi(3); + let want_area = 4.0 * std::f64::consts::PI * radius.powi(2); + assert!(rel(m.volume, want_volume) < 1e-3, "volume {} vs {want_volume}", m.volume); + assert!(rel(m.area, want_area) < 1e-3, "area {} vs {want_area}", m.area); + let mesh = body.body.display_mesh(radius * 1e-3).unwrap(); + assert!(!mesh.positions.is_empty() && !mesh.triangles.is_empty()); +} + +#[test] +fn seam_only_sphere_imports_as_a_whole_ball() { + let text = seam_sphere_step(Vec3::ZERO, 2.0, Vec3::Z, Vec3::X, 0.0, false, std::f64::consts::PI); + assert_seam_sphere(&text, 2.0); +} + +#[test] +fn tilted_seam_only_sphere_imports_as_a_whole_ball() { + let axis = Vec3::new(0.0, std::f64::consts::FRAC_1_SQRT_2, std::f64::consts::FRAC_1_SQRT_2); + let text = seam_sphere_step(Vec3::new(30.0, -12.0, 7.0), 1.9845, axis, Vec3::X, 0.7, false, std::f64::consts::PI); + assert_seam_sphere(&text, 1.9845); +} + +#[test] +fn reversed_seam_only_sphere_imports_as_a_whole_ball() { + let text = seam_sphere_step(Vec3::ZERO, 2.0, Vec3::Z, Vec3::X, 0.0, true, std::f64::consts::PI); + assert_seam_sphere(&text, 2.0); +} + +#[test] +fn seam_only_sphere_whose_placement_axis_misses_its_poles_imports_as_a_whole_ball() { + // The seam runs between poles on the X axis, while the surface's placement points along Z. + let text = seam_sphere_step(Vec3::new(5.0, 6.0, 7.0), 3.0, Vec3::X, Vec3::Y, 0.4, false, std::f64::consts::PI); + let placement = |axis: &str, reference: &str| format!("#2=DIRECTION('',{axis});\n#3=DIRECTION('',{reference});\n"); + let (one, zero) = ("1.0000000000000000", "0.0000000000000000"); + let (x, y, z) = (format!("({one},{zero},{zero})"), format!("({zero},{one},{zero})"), format!("({zero},{zero},{one})")); + let text = text.replacen(&placement(&x, &y), &placement(&z, &x), 1); + assert!(text.contains(&placement(&z, &x)), "the placement was not replaced"); + assert_seam_sphere(&text, 3.0); +} + +#[test] +fn seam_only_sphere_void_is_turned_inside_out() { + // A 5 mm ball with a 2 mm ball-shaped void, written as voids are: the void's closed shell + // (faces facing out of it) used through an ORIENTED_CLOSED_SHELL with orientation .F. + let renumber = |text: &str, by: u64| { + let mut out = String::new(); + let mut chars = text.chars().peekable(); + while let Some(c) = chars.next() { + out.push(c); + if c == '#' { + let mut n = String::new(); + while let Some(d) = chars.next_if(char::is_ascii_digit) { + n.push(d); + } + out += &(n.parse::().unwrap() + by).to_string(); + } + } + out + }; + let inner = seam_sphere_step(Vec3::ZERO, 2.0, Vec3::Z, Vec3::X, 0.0, false, std::f64::consts::PI); + let inner = &inner[inner.find("#1=").unwrap()..inner.find("#21=").unwrap()]; + let outer = seam_sphere_step(Vec3::ZERO, 5.0, Vec3::Z, Vec3::X, 1.0, false, std::f64::consts::PI); + let brep = "#21=MANIFOLD_SOLID_BREP('Ball',#20);\n"; + assert!(outer.contains(brep)); + let void = format!("{}#200=ORIENTED_CLOSED_SHELL('',*,#120,.F.);\n#21=BREP_WITH_VOIDS('Hollow',#20,(#200));\n", renumber(inner, 100)); + let text = outer.replacen(brep, &void, 1); + let imp = step_import(&text).unwrap(); + assert_eq!(imp.bodies.len(), 1, "{:?}", imp.warnings); + assert!(imp.warnings.is_empty(), "{:?}", imp.warnings); + assert!(imp.bodies[0].closed); + let m = measure(&imp.bodies[0].body).unwrap(); + let want = 4.0 / 3.0 * std::f64::consts::PI * (125.0 - 8.0); + assert!(rel(m.volume, want) < 1e-3, "volume {} vs {want}", m.volume); + let area = 4.0 * std::f64::consts::PI * (25.0 + 4.0); + assert!(rel(m.area, area) < 1e-3, "area {} vs {area}", m.area); +} + +#[test] +fn non_antipodal_seam_sphere_does_not_panic() { + let text = seam_sphere_step(Vec3::ZERO, 2.0, Vec3::Z, Vec3::X, 0.0, false, 2.5); + let result = std::panic::catch_unwind(|| step_import(&text)); + assert!(result.is_ok()); +} + #[test] #[ignore = "manual STEP import timing; run without concurrent builds"] fn benchmark_independent_step_solids() {