From d76ed36e3ad657423f68662714987e6b6d47fc4a Mon Sep 17 00:00:00 2001 From: Akhil Velagapudi Date: Fri, 9 Oct 2026 17:35:07 +0000 Subject: [PATCH 1/4] Bound NURBS projection refinement near collapsed edges When the incumbent is the pole of a collapsed edge, every patch along that edge reaches it, so certifying one projection examined millions of patches: a 4x3 NURBS face on an FMC connector took minutes. Settle for the incumbent after 65536 patches, twice the most any ordinary query in the board corpus examines. Amp-Thread-ID: https://ampcode.com/threads/T-01a12192-0f95-7657-bde4-754232566545 Co-authored-by: Amp --- nurbs/src/sampled_surface.rs | 46 ++++++++++++++++++++++++++++++++++++ 1 file changed, 46 insertions(+) diff --git a/nurbs/src/sampled_surface.rs b/nurbs/src/sampled_surface.rs index f3e02e4..c35da05 100644 --- a/nurbs/src/sampled_surface.rs +++ b/nurbs/src/sampled_surface.rs @@ -26,6 +26,13 @@ pub struct ProjectionScratch { } const PROJECTION_TOL: f64 = 64. * f64::EPSILON; +/// Patches one projection may examine before it settles for its incumbent. +/// Ordinary queries examine at most a few hundred, and tens of thousands on +/// the worst surfaces seen. When the incumbent is the pole of a collapsed +/// edge, every patch along that edge reaches it: each subdivision halves +/// their gap to the incumbent but doubles their number, so certifying it at +/// f64 resolution would examine millions of patches for a single point. +const MAX_PROJECTION_PATCHES: usize = 1 << 16; #[derive(Debug, Clone)] struct SurfaceCell { @@ -544,7 +551,14 @@ where cell.spans)); } } + let mut examined = 0; while let Some((patch, lo, hi, spans)) = queue.pop_front() { + // The incumbent is a feasible point no farther than any patch + // evaluated so far; keep it rather than refine without end. + examined += 1; + if examined > MAX_PROJECTION_PATCHES { + break; + } let mid = (lo + hi) * 0.5; let d = self.surf.derivs_in_span::<2>(mid, spans, p); let normal = d[1][0].cross(&d[0][1]); @@ -682,6 +696,38 @@ mod tests { close((sampled.surf.point(uv)-target).norm_squared(), 0.02, 1e-12); } + #[test] + fn projection_near_a_collapsed_edge_terminates() { + // Face #62825 of a Samtec FMC connector (ASP-134486-01): the u = max + // row of control points collapses to one pole. Every patch along that + // edge reaches the pole, so certifying this projection at f64 + // resolution examined millions of patches. + let rows = [ + [(-5.589503327713219782, 3.795967870219560059, 24.88456541861770077, 1.), + (-5.632056464503211402, 3.715554243323587080, 24.90841966077580594, 0.8294535963452530636), + (-5.582626279250799683, 3.666787751982639953, 24.97187345677104986, 1.)], + [(-5.569374495857790208, 3.760876034061380047, 24.84803886198849909, 1.), + (-5.599499045440100353, 3.721231129362655210, 24.88421881863200369, 0.9379770243229920146), + (-5.582626266339880239, 3.666787755634880241, 24.90864600195275003, 1.)], + [(-5.570182848018389699, 3.711200168230984531, 24.80680506253315087, 1.), + (-5.574716636290549587, 3.689228228810422117, 24.82235407503383229, 1.000927775394292008), + (-5.582626295171389863, 3.666787741024605030, 24.83698483065439788, 1.)], + [(-5.582626325353579944, 3.666787764895184765, 24.76551111168974728, 1.); 3], + ]; + let surf = NDBSplineSurface::new(true, true, + KnotVector::from_multiplicities(3, &[0.001418343725273965962, 0.1860978202366415091], &[4, 4]), + KnotVector::from_multiplicities(2, &[0., 0.09081581946402909888], &[3, 3]), + rows.iter().map(|row| row.iter().map(|&(x, y, z, w)| + nalgebra_glm::DVec4::new(x * w, y * w, z * w, w)).collect()).collect()); + let sampled = SampledSurface::new(surf); + let target = DVec3::new(-5.58262627960619, 3.666787751487475, 24.765511111689747); + // The pole is the nearest point: the target lies just beyond it. + let uv = sampled.uv_from_point(target).unwrap(); + let pole = rows[3][0]; + let pole = DVec3::new(pole.0, pole.1, pole.2); + close((sampled.surf.point(uv) - target).norm(), (pole - target).norm(), 1e-12); + } + #[test] fn single_cell_projection_resolves_close_branches_in_parameter_space() { let surf = NDBSplineSurface::new(true, true, From b33fb0e49e736d6f44770e0b25da844140632d5a Mon Sep 17 00:00:00 2001 From: Akhil Velagapudi Date: Fri, 9 Oct 2026 17:35:07 +0000 Subject: [PATCH 2/4] Keep the first items of over-long bounded lists Some exporters write four-component directions, which panicked the parser and lost the whole model. Keep the first items, as OCCT does. Amp-Thread-ID: https://ampcode.com/threads/T-01a12192-0f95-7657-bde4-754232566545 Co-authored-by: Amp --- step/src/parse.rs | 12 +++++++++++- 1 file changed, 11 insertions(+), 1 deletion(-) diff --git a/step/src/parse.rs b/step/src/parse.rs index 325bc4c..912ac29 100644 --- a/step/src/parse.rs +++ b/step/src/parse.rs @@ -119,7 +119,10 @@ impl<'a, T: Parse<'a>, const CAP: usize> Parse<'a> for ArrayVec { e => e?, }; s = s_; - out.push(o); + // Some exporters write more items than the bound allows, such + // as four-component directions. Keep the first ones, as OCCT + // does, rather than reject the file. + let _ = out.try_push(o); } let (s, _) = char(')')(s)?; Ok((s, out)) @@ -367,6 +370,13 @@ mod tests { assert!(<&str>::parse("'unterminated").is_err()); } + #[test] + fn bounded_lists_keep_their_first_items() { + // Some exporters write four-component directions. + let (rest, ratios) = ArrayVec::::parse("(0.,0.,1.,0.),next").unwrap(); + assert_eq!((rest, ratios.as_slice()), (",next", &[0., 0., 1.][..])); + } + #[test] fn test_parse_entity_decl() { parse_entity_decl(b"#3=SHAPE_DEFINITION_REPRESENTATION(#4,#10);").unwrap(); From 7b516e6fc721d6d77290464cca6e8aa45fabb366 Mon Sep 17 00:00:00 2001 From: Akhil Velagapudi Date: Fri, 9 Oct 2026 17:35:07 +0000 Subject: [PATCH 3/4] Mesh faces in parallel and shapes in a fixed order Sample every edge of a shell up front, so faces sharing an edge still share its points, then mesh the faces on the rayon pool and merge them in face order. Shapes and their copies were meshed in hash map order, and colour groups of equal size were emitted in hash map order, so output changed from run to run. Sort shapes by entity id, copies by transform, and colour groups by size, then colour. Amp-Thread-ID: https://ampcode.com/threads/T-01a12192-0f95-7657-bde4-754232566545 Co-authored-by: Amp --- triangulate/src/colored_mesh.rs | 35 ++++++- triangulate/src/triangulate.rs | 156 ++++++++++++++++++++++++++------ 2 files changed, 157 insertions(+), 34 deletions(-) diff --git a/triangulate/src/colored_mesh.rs b/triangulate/src/colored_mesh.rs index d1396b6..6e4a0f4 100644 --- a/triangulate/src/colored_mesh.rs +++ b/triangulate/src/colored_mesh.rs @@ -124,10 +124,12 @@ pub fn group_mesh_by_color(mesh: &Mesh) -> Result { bucket.indices.extend([a, b, c]); } - let mut submeshes = buckets.into_values().collect::>(); - // Largest groups first for deterministic ordering. - submeshes.sort_by(|a, b| a.indices.len().cmp(&b.indices.len()).reverse()); - Ok(TessellatedMesh { submeshes }) + let mut submeshes = buckets.into_iter().collect::>(); + // Largest groups first. Map order is random, so colour breaks ties. + submeshes.sort_by_key(|(key, submesh)| (std::cmp::Reverse(submesh.indices.len()), *key)); + Ok(TessellatedMesh { + submeshes: submeshes.into_iter().map(|(_, submesh)| submesh).collect(), + }) } // --------------------------------------------------------------------------- @@ -135,7 +137,7 @@ pub fn group_mesh_by_color(mesh: &Mesh) -> Result { // --------------------------------------------------------------------------- /// Quantised RGBA colour key for bucketing. -#[derive(Clone, Copy, Debug, PartialEq, Eq, Hash)] +#[derive(Clone, Copy, Debug, PartialEq, Eq, PartialOrd, Ord, Hash)] struct ColorKey(u8, u8, u8, u8); impl ColorKey { @@ -261,6 +263,29 @@ mod tests { assert_eq!(green.indices, [0, 1, 2]); } + #[test] + fn orders_equal_sized_color_groups_the_same_every_time() { + let mut mesh = Mesh { + verts: (0..6) + .map(|i| vertex(DVec3::new(i as f64, (i % 3) as f64, 0.0), DVec3::new(0.0, 0.0, 1.0))) + .collect(), + triangles: vec![ + Triangle { verts: U32Vec3::new(0, 1, 2) }, + Triangle { verts: U32Vec3::new(3, 4, 5) }, + ], + }; + for (i, v) in mesh.verts.iter_mut().enumerate() { + v.color = if i < 3 { DVec3::x() } else { DVec3::y() }; + } + // Each call buckets colours in a new map, with new random keys. + let colors = || group_mesh_by_color(&mesh).unwrap().submeshes.iter().map(|s| s.color).collect::>(); + let first = colors(); + assert_eq!(first.len(), 2); + for _ in 0..16 { + assert_eq!(colors(), first); + } + } + #[test] fn keeps_coincident_vertices_separate_across_hard_edges() { let mesh = Mesh { diff --git a/triangulate/src/triangulate.rs b/triangulate/src/triangulate.rs index ea41292..739e7bd 100644 --- a/triangulate/src/triangulate.rs +++ b/triangulate/src/triangulate.rs @@ -62,7 +62,7 @@ fn collect_shape_instances<'a>( s: &'a StepFile, rep_instances: &HashMap, Vec>, shape_rep_relationship: &HashMap, Vec>>, -) -> HashMap, ShapeGeometry> { +) -> Vec<(RepresentationItem<'a>, ShapeGeometry)> { let mut todo: Vec<_> = rep_instances .iter() .flat_map(|(rep, mats)| mats.iter().copied().map(move |mat| (*rep, mat))) @@ -126,6 +126,15 @@ fn collect_shape_instances<'a>( }); } + // The maps above are visited in a random order. Mesh shapes, and their + // copies, in a fixed one so the output does not change from run to run. + let mut to_mesh: Vec<_> = to_mesh.into_iter().collect(); + to_mesh.sort_unstable_by_key(|(id, _)| id.0); + for (_, shape) in &mut to_mesh { + shape.instances.sort_unstable_by(|a, b| { + a.iter().map(|v| v.to_bits()).cmp(b.iter().map(|v| v.to_bits())) + }); + } to_mesh } @@ -839,13 +848,12 @@ fn shell( }; let v_start = mesh.verts.len(); let t_start = mesh.triangles.len(); - let mut scratch = Mesh::default(); - let mut edge_samples = HashMap::new(); - for face in faces { - scratch.verts.clear(); - scratch.triangles.clear(); - stats.num_faces += 1; - match advanced_face( + // Sample every edge once up front, so faces sharing an edge share its + // points exactly and the faces can then be meshed independently. + let edge_samples = sample_edges(s, faces, tolerance); + let mesh_face = |face: &Face| { + let mut scratch = Mesh::default(); + let result = advanced_face( s, *face, &mut scratch, @@ -853,8 +861,19 @@ fn shell( default_color, uncertainty, tolerance, - &mut edge_samples, - ) { + &edge_samples, + ); + (result, scratch) + }; + #[cfg(feature = "rayon")] + let meshed: Vec<_> = faces.par_iter().with_max_len(1).map(mesh_face).collect(); + #[cfg(not(feature = "rayon"))] + let meshed: Vec<_> = faces.iter().map(mesh_face).collect(); + // Results are merged in face order, so output does not depend on + // scheduling. + for (face, (result, mut scratch)) in faces.iter().zip(meshed) { + stats.num_faces += 1; + match result { Ok(()) => mesh.append(&mut scratch), Err(err) => { let surface_id = match &s[*face] { @@ -906,7 +925,7 @@ fn advanced_face( default_color: DVec3, uncertainty: f64, tolerance: f64, - edge_samples: &mut HashMap>, + edge_samples: &HashMap>, ) -> Result<(), Error> { // Closed shells may legally reference either ADVANCED_FACE or the more // general FACE_SURFACE; OCCT/KiCad translate both through FaceSurface. @@ -1605,7 +1624,7 @@ fn face_bound( b: FaceBound, edge_uses: &mut HashMap, tolerance: f64, - edge_samples: &mut HashMap>, + edge_samples: &HashMap>, ) -> Result<(Vec, usize), Error> { let (bound, orientation) = match &s[b] { Entity::FaceBound(b) => (b.bound, b.orientation), @@ -1636,8 +1655,46 @@ fn face_bound( } } +/// The 3D samples of every edge bounding `faces`, by edge curve id. A shell +/// has one STEP file and tolerance, so an edge's samples are the same for +/// every face that uses it; traversal and surface charts remain +/// face-specific. Edges that cannot be sampled are left out. +fn sample_edges(s: &StepFile, faces: &[Face], tolerance: f64) -> HashMap> { + let mut edges: Vec = faces + .iter() + .filter_map(|f| match &s[*f] { + Entity::AdvancedFace(face) => Some(&face.bounds[..]), + Entity::FaceSurface(face) => Some(&face.bounds[..]), + _ => None, + }) + .flatten() + .filter_map(|b| match &s[*b] { + Entity::FaceBound(b) => Some(b.bound), + Entity::FaceOuterBound(b) => Some(b.bound), + _ => None, + }) + .filter_map(|bound| match &s[bound] { + Entity::EdgeLoop(e) => Some(&e.edge_list[..]), + _ => None, + }) + .flatten() + .filter_map(|e| s.entity(*e).map(|edge| edge.edge_element.0)) + .collect(); + edges.sort_unstable(); + edges.dedup(); + let sample = |&id: &usize| { + edge_curve(s, EdgeCurve::new(id), true, tolerance) + .ok() + .map(|points| (id, points)) + }; + #[cfg(feature = "rayon")] + return edges.par_iter().filter_map(sample).collect(); + #[cfg(not(feature = "rayon"))] + return edges.iter().filter_map(sample).collect(); +} + fn edge_loop(s: &StepFile, edge_list: &[OrientedEdge], tolerance: f64, - edge_samples: &mut HashMap>, + edge_samples: &HashMap>, ) -> Result, Error> { let mut out = Vec::new(); for (i, e) in edge_list.iter().enumerate() { @@ -1649,12 +1706,15 @@ fn edge_loop(s: &StepFile, edge_list: &[OrientedEdge], tolerance: f64, let edge = s .entity(*e) .ok_or(Error::InvalidStepEntity("OrientedEdge"))?; - // A shell has one STEP file and tolerance. Cache only edge-owned 3D - // samples; traversal and surface charts remain face-specific. - let points = match edge_samples.entry(edge.edge_element.0) { - std::collections::hash_map::Entry::Occupied(entry) => entry.into_mut(), - std::collections::hash_map::Entry::Vacant(entry) => - entry.insert(edge_curve(s, edge.edge_element.cast(), true, tolerance)?), + // Edges that could not be sampled up front fail here, for each face + // that uses them. + let sampled; + let points = match edge_samples.get(&edge.edge_element.0) { + Some(points) => points, + None => { + sampled = edge_curve(s, edge.edge_element.cast(), true, tolerance)?; + &sampled + } }; if edge.orientation { out.extend_from_slice(points); @@ -2089,7 +2149,7 @@ mod tests { let step = StepFile::parse(&flat).unwrap(); let mut mesh = Mesh::default(); let color = DVec3::new(0.2, 0.4, 0.6); - advanced_face(&step, Id::new(72), &mut mesh, &HashMap::new(), color, 0., 0.01, &mut HashMap::new()).unwrap(); + advanced_face(&step, Id::new(72), &mut mesh, &HashMap::new(), color, 0., 0.01, &HashMap::new()).unwrap(); assert!(mesh.verts.iter().any(|v| v.pos == DVec3::new(1., 1., 0.))); assert!(mesh .verts @@ -2419,7 +2479,7 @@ mod tests { let mut previous = 0; for tolerance in [0.04,0.01] { let mut mesh = Mesh::default(); - advanced_face(&step, Id::new(13), &mut mesh, &HashMap::new(), DVec3::zeros(), 0., tolerance, &mut HashMap::new()).unwrap(); + advanced_face(&step, Id::new(13), &mut mesh, &HashMap::new(), DVec3::zeros(), 0., tolerance, &HashMap::new()).unwrap(); assert!(mesh.triangles.len() > previous); previous = mesh.triangles.len(); for t in &mesh.triangles { @@ -2488,12 +2548,14 @@ mod tests { let mut backward = edge_curve(&step, Id::new(id), false, 0.01).unwrap(); backward.reverse(); assert_eq!(forward, backward); - let mut samples = HashMap::new(); - let mut cached_backward = edge_loop(&step, &[Id::new(backward_id)], 0.01, &mut samples).unwrap(); - cached_backward.reverse(); - assert_eq!(forward, cached_backward); - assert_eq!(forward, edge_loop(&step, &[Id::new(forward_id)], 0.01, &mut samples).unwrap()); - assert_eq!(samples.len(), 1); + // Shared samples and sampling on demand agree. + let shared = HashMap::from([(id, forward.clone())]); + for samples in [&shared, &HashMap::new()] { + let mut cached_backward = edge_loop(&step, &[Id::new(backward_id)], 0.01, samples).unwrap(); + cached_backward.reverse(); + assert_eq!(forward, cached_backward); + assert_eq!(forward, edge_loop(&step, &[Id::new(forward_id)], 0.01, samples).unwrap()); + } } } } @@ -2528,8 +2590,44 @@ mod tests { .cloned() .collect(); let shapes = collect_shape_instances(&step, &roots, &HashMap::new()); - assert_eq!(shapes[&Id::new(10)].uncertainty, 2e-7); - assert!((shapes[&Id::new(11)].uncertainty - 2e-10).abs() < 1e-25); + // Shapes come in entity order. + assert_eq!(shapes.iter().map(|(id, _)| id.0).collect::>(), [10, 11]); + assert_eq!(shapes[0].1.uncertainty, 2e-7); + assert!((shapes[1].1.uncertainty - 2e-10).abs() < 1e-25); + } + + #[test] + fn shapes_and_their_copies_come_in_a_fixed_order() { + // Shape #10 is used by both representations, so it collects copies + // from each in whatever order the representation map is visited. + let text = b"ISO-10303-21;HEADER;ENDSEC;DATA; + #4=REPRESENTATION_CONTEXT('',''); + #7=SHAPE_REPRESENTATION('',(#10),#4); + #8=SHAPE_REPRESENTATION('',(#11,#10),#4); + #9=CLOSED_SHELL('',()); + #10=MANIFOLD_SOLID_BREP('',#9); + #11=MANIFOLD_SOLID_BREP('',#9); + ENDSEC;END-ISO-10303-21;"; + let flat = StepFile::strip_flatten(text).unwrap(); + let step = StepFile::parse(&flat).unwrap(); + let shift = |x: f64| glm::translation(&DVec3::new(x, 0., 0.)); + let order = || { + // Each new map hashes with new random keys. + let roots = HashMap::from([ + (Id::new(7), vec![shift(1.), shift(2.)]), + (Id::new(8), vec![shift(3.)]), + ]); + collect_shape_instances(&step, &roots, &HashMap::new()) + .iter() + .map(|(id, shape)| (id.0, shape.instances.iter().map(|m| m[(0, 3)]).collect::>())) + .collect::>() + }; + let first = order(); + assert_eq!(first.len(), 2); + assert_eq!(first.iter().map(|(_, copies)| copies.len()).sum::(), 4); + for _ in 0..16 { + assert_eq!(order(), first); + } } #[test] From bebd889b96c23d1a05e4d741158d6dce3abe8391 Mon Sep 17 00:00:00 2001 From: Akhil Velagapudi Date: Fri, 9 Oct 2026 17:35:07 +0000 Subject: [PATCH 4/4] Mesh whole spheres bounded by a vertex loop BGA solder balls are whole spheres bounded only by a vertex loop at the pole, which no chart covers, so every ball failed to tessellate. Mesh them as rings of latitude between two poles within the tolerance. Amp-Thread-ID: https://ampcode.com/threads/T-01a12192-0f95-7657-bde4-754232566545 Co-authored-by: Amp --- triangulate/src/surface.rs | 54 ++++++++++++++++++++++++++++++++++ triangulate/src/triangulate.rs | 48 ++++++++++++++++++++++++++++++ 2 files changed, 102 insertions(+) diff --git a/triangulate/src/surface.rs b/triangulate/src/surface.rs index 68e4fac..ac215d6 100644 --- a/triangulate/src/surface.rs +++ b/triangulate/src/surface.rs @@ -328,6 +328,9 @@ impl Surface { /// Tessellate a compact surface with no physical trim. Identify seam /// vertices by wrapped grid indices instead of cutting a planar polygon. pub fn untrimmed_mesh(&self, color: DVec3, same_sense: bool, tolerance: f64) -> Option { + if let Self::Sphere { location, radius } = self { + return Self::sphere_mesh(*location, *radius, color, same_sense, tolerance); + } let Self::Torus { mat, major_radius, @@ -376,6 +379,57 @@ impl Surface { Some(mesh) } + /// A whole sphere, as a face bounded only by a vertex loop at a pole is: + /// rings of latitude between two pole vertices. + fn sphere_mesh( + center: DVec3, + radius: f64, + color: DVec3, + same_sense: bool, + tolerance: f64, + ) -> Option { + if !(radius > 0.) { + return None; + } + // As for the torus: a chord of angle a sags r*a²/8 off the sphere; + // allocate half the error budget to each angular direction. + let nu = (2. * PI * (radius / (2. * tolerance)).sqrt()).ceil().max(6.) as u32; + let nv = (nu / 2).max(3); + nu.checked_mul(nv)?; + let sense = if same_sense { 1. } else { -1. }; + // The south pole, nv - 1 rings of nu vertices each, the north pole. + let mut normals = vec![DVec3::new(0., 0., -1.)]; + for v in 1..nv { + let (sv, cv) = (PI * v as f64 / nv as f64).sin_cos(); + normals.extend((0..nu).map(|u| { + let (su, cu) = (2. * PI * u as f64 / nu as f64).sin_cos(); + DVec3::new(sv * cu, sv * su, -cv) + })); + } + normals.push(DVec3::new(0., 0., 1.)); + let mut mesh = Mesh::default(); + mesh.verts = normals + .iter() + .map(|&n| Vertex { pos: center + n * radius, norm: n * sense, color }) + .collect(); + let (south, north) = (0, mesh.verts.len() as u32 - 1); + let ring = |v: u32, u: u32| 1 + (v - 1) * nu + u % nu; + // Counter-clockwise seen from outside. + let mut push = |a: u32, b: u32, c: u32| { + let verts = if same_sense { [a, b, c] } else { [a, c, b] }; + mesh.triangles.push(Triangle { verts: verts.into() }); + }; + for u in 0..nu { + push(south, ring(1, u + 1), ring(1, u)); + for v in 1..nv - 1 { + push(ring(v, u), ring(v, u + 1), ring(v + 1, u + 1)); + push(ring(v, u), ring(v + 1, u + 1), ring(v + 1, u)); + } + push(north, ring(nv - 1, u), ring(nv - 1, u + 1)); + } + Some(mesh) + } + pub fn new_plane(axis: DVec3) -> Result { let dropped = axis.iamax(); let scale = axis[dropped].abs(); diff --git a/triangulate/src/triangulate.rs b/triangulate/src/triangulate.rs index 739e7bd..fd504cf 100644 --- a/triangulate/src/triangulate.rs +++ b/triangulate/src/triangulate.rs @@ -2457,6 +2457,54 @@ mod tests { assert!(edges.is_empty()); } + #[test] + fn whole_sphere_bounded_by_a_pole_is_meshed() { + // A BGA solder ball: the whole sphere, bounded only by a vertex + // loop at the pole where it meets the package. + for (sense, outward) in [(".T.", 1.), (".F.", -1.)] { + let text = format!("ISO-10303-21;HEADER;ENDSEC;DATA; + #1=CARTESIAN_POINT('',(1.,2.,0.115)); + #2=DIRECTION('',(0.,0.,1.)); + #3=DIRECTION('',(1.,0.,0.)); + #4=AXIS2_PLACEMENT_3D('',#1,#2,#3); + #5=SPHERICAL_SURFACE('',#4,0.115); + #6=CARTESIAN_POINT('',(1.,2.,0.)); + #7=VERTEX_POINT('',#6); + #8=VERTEX_LOOP('',#7); + #9=FACE_BOUND('',#8,.T.); + #10=ADVANCED_FACE('',(#9),#5,{sense}); + ENDSEC;END-ISO-10303-21;"); + let flat = StepFile::strip_flatten(text.as_bytes()).unwrap(); + let step = StepFile::parse(&flat).unwrap(); + let tolerance = 0.01; + let mut mesh = Mesh::default(); + advanced_face(&step, Id::new(10), &mut mesh, &HashMap::new(), DVec3::zeros(), 0., tolerance, &HashMap::new()).unwrap(); + let center = DVec3::new(1., 2., 0.115); + for v in &mesh.verts { + assert!(((v.pos - center).norm() - 0.115).abs() < 1e-12); + assert!(((v.pos - center) / 0.115 * outward - v.norm).norm() < 1e-9); + } + // Closed: every edge is used once in each direction. + let mut edges = HashMap::<_, i32>::new(); + let mut area = 0.; + for t in &mesh.triangles { + let [a, b, c] = [t.verts.x, t.verts.y, t.verts.z].map(|i| mesh.verts[i as usize].pos); + let p = (a + b + c) / 3.; + assert!(0.115 - (p - center).norm() < tolerance); + let n = (b - a).cross(&(c - a)); + assert!(n.dot(&(p - center)) * outward > 0.); + area += n.norm() / 2.; + for (u, v) in [(t.verts.x, t.verts.y), (t.verts.y, t.verts.z), (t.verts.z, t.verts.x)] { + *edges.entry((u.min(v), u.max(v))).or_default() += if u < v { 1 } else { -1 }; + } + } + assert_eq!(edges.len() * 2, mesh.triangles.len() * 3); + assert!(edges.values().all(|&d| d == 0)); + let sphere = 4. * std::f64::consts::PI * 0.115f64.powi(2); + assert!((area - sphere).abs() < 0.05 * sphere, "{area} vs {sphere}"); + } + } + #[test] fn hemisphere_refines_interior_curvature_to_the_physical_budget() { let text = b"ISO-10303-21;HEADER;ENDSEC;DATA;