diff --git a/crates/engine/src/sample.rs b/crates/engine/src/sample.rs index 1306875..bd72685 100644 --- a/crates/engine/src/sample.rs +++ b/crates/engine/src/sample.rs @@ -741,6 +741,21 @@ mod tests { assert!(r["mass"]["relative_change"].as_f64().unwrap().abs() < 0.01, "{r}"); } + /// The heated bracket's 4 mm tet mesh stays near its 23 000 elements and opens within the + /// control channel's 60 s request timeout. #303: refinement piled face vertices against a + /// feature curve, and the improvement then grew the mesh to 29 257 elements in about two + /// minutes. + #[test] + #[ignore = "heavy: cargo xtask test --heavy (CI heavy job)"] + fn the_heated_bracket_opens_quickly_with_its_usual_mesh() { + let start = std::time::Instant::now(); + let s = open("heated-bracket"); + let elapsed = start.elapsed(); + let elements = s.model.elements().len(); + assert!((20_000..=26_000).contains(&elements), "{elements} elements"); + assert!(elapsed < std::time::Duration::from_secs(60), "opened in {elapsed:?}"); + } + /// The heated bracket's temperatures lie between the air's and the heater's, falling along /// the bracket away from the heater. Its thermal expansion: with the air at the heater's /// temperature the bracket warms evenly by ΔT, and held at the foot's hole alone (at x = 40) diff --git a/crates/mesh/src/tet/refine.rs b/crates/mesh/src/tet/refine.rs index 6f6d383..5d77dcc 100644 --- a/crates/mesh/src/tet/refine.rs +++ b/crates/mesh/src/tet/refine.rs @@ -208,6 +208,19 @@ pub fn split_edge(cx: &mut TetComplex, ctx: &Context, u: u32, w: u32, tag: u32, None if tag != NONE => VertexKind::Face(tag), None => return Ok(None), }; + // A midpoint inside the diametral sphere of a nearby feature edge encroaches it: split that + // edge instead (Shewchuk 1998), or face vertices pile up against the curve, each split + // halving their distance to it. + if matches!(kind, VertexKind::Face(_)) { + let mut near = vec![u, w]; + if let Some(ring) = cx.edge_ring(u, w)? { + near.extend(ring.verts.first().copied()); + near.extend(ring.verts.last().copied()); + } + if let Some((a, b)) = encroached_curve_edge(cx, &near, m)? { + return split_edge(cx, ctx, a, b, tag, min_face); + } + } let Some(p) = ctx.place(kind, m) else { return Ok(None) }; // The projection may only move the midpoint a little (a quarter of the edge). if measure::dist(p.coords(), m) > 0.25 * len { @@ -216,6 +229,21 @@ pub fn split_edge(cx: &mut TetComplex, ctx: &Context, u: u32, w: u32, tag: u32, ops::insert_delaunay_apart(cx, p, kind, Site::Edge(u, w), MAX_CAVITY, MIN_GAP * ctx.h(m, len)) } +/// A feature-curve edge at one of the vertices `near` whose diametral sphere holds `p`. +fn encroached_curve_edge(cx: &mut TetComplex, near: &[u32], p: P) -> Result> { + for &c in near { + for n in cx.neighbours(c)? { + if cx.curve(c, n).is_some() { + let (a, b) = (cx.xyz(c)?, cx.xyz(n)?); + if measure::dist(p, measure::mid(a, b)) < 0.5 * measure::dist(a, b) { + return Ok(Some((c, n))); + } + } + } + } + Ok(None) +} + /// Whether a tetrahedron is bad, and its circumcentre. fn badness(x: &[P; 4], ctx: &Context, o: &RefineOptions) -> Option

{ let c = measure::circumcenter(x)?; diff --git a/crates/mesh/src/tet/tests.rs b/crates/mesh/src/tet/tests.rs index 1986b17..4950641 100644 --- a/crates/mesh/src/tet/tests.rs +++ b/crates/mesh/src/tet/tests.rs @@ -195,6 +195,53 @@ mod pipeline { assert_eq!(a, b); } + /// Refinement splitting a face edge whose midpoint lies inside a feature edge's diametral + /// sphere splits the feature edge instead (#303): otherwise each split halves the distance + /// of the face vertices to the curve and they pile up against it. + #[test] + fn a_face_edge_split_that_encroaches_a_curve_splits_the_curve() { + use super::super::complex::{NONE, VertexKind, face_toward}; + use super::super::geometry::Tessellation; + use super::super::ops::{Site, insert_delaunay}; + use super::super::refine::{Context, MAX_CAVITY, split_edge}; + use super::super::{complex_from_cdt, measure}; + use crate::predicates::Point3; + use crate::size::UniformSize; + // Whether `p` lies in triangle `f` of the z = 0 plane. + let contains = |f: [[f64; 3]; 3], p: [f64; 3]| { + let side = |a: [f64; 3], b: [f64; 3]| (b[0] - a[0]) * (p[1] - a[1]) - (b[1] - a[1]) * (p[0] - a[0]); + let s = [side(f[0], f[1]), side(f[1], f[2]), side(f[2], f[0])]; + s.iter().all(|&x| x >= 0.0) || s.iter().all(|&x| x <= 0.0) + }; + let s = box_surface(); + let cdt = tetrahedralize(&s, &CdtOptions::default()).unwrap(); + let (mut cx, curve_faces) = complex_from_cdt(&s, &cdt, &QualityOptions::default()).unwrap(); + let (geometry, size) = (Tessellation::new(&s, Vec::new()), UniformSize(5.0)); + let ctx = Context { size: &size, geometry: &geometry, curve_faces: &curve_faces }; + // A curve edge on the z = 0 side, and a face vertex beside it on that side. + let ((a, b), _) = cx.curve_edges().into_iter().find(|&((a, b), _)| [a, b].iter().all(|&v| cx.xyz(v).unwrap()[2] == 0.0)).unwrap(); + let (xa, xb) = (cx.xyz(a).unwrap(), cx.xyz(b).unwrap()); + let len = measure::dist(xa, xb); + let centre = [5.0, 5.0, 0.0]; + let mid = measure::mid(xa, xb); + let inward = measure::scale(measure::sub(centre, mid), 0.3 * len / measure::dist(centre, mid)); + let q = Point3::try_from(measure::add(mid, inward)).unwrap(); + let (t, k) = cx + .live_tets() + .flat_map(|(t, vs)| (0..4).map(move |k| (t, vs, k))) + .find(|&(t, vs, k)| { + cx.adj(t).unwrap()[k] == NONE && cx.tags(t).unwrap()[k] == 4 && contains(face_toward(vs, k).map(|v| cx.xyz(v).unwrap()), q.coords()) + }) + .map(|(t, _, k)| (t, k)) + .unwrap(); + let v = insert_delaunay(&mut cx, q, VertexKind::Face(4), Site::Face(t, k), MAX_CAVITY).unwrap().unwrap(); + assert!(cx.edge_tet(v, a).unwrap().is_some()); + let w = split_edge(&mut cx, &ctx, v, a, 4, 0.0).unwrap().unwrap(); + assert!(matches!(cx.kind(w), VertexKind::Curve(_)), "{:?}", cx.kind(w)); + assert!(measure::dist(cx.xyz(w).unwrap(), mid) < 1e-9); + cx.validate().unwrap(); + } + #[test] fn a_seeded_edge_gets_the_seeded_spacing_and_the_size_grows_from_it() { use crate::seeds::EdgeSeed; diff --git a/docs/examples.md b/docs/examples.md index 295d119..fe66083 100644 --- a/docs/examples.md +++ b/docs/examples.md @@ -76,8 +76,9 @@ The corpus filleted bracket (`tests/corpus/filleted-bracket`: an 8 mm L-bracket, hole through each leg) in steel, both hole bores clamped. The end of the upright is held at 120 °C by a heater; the rest of its outside loses heat to still air at 20 °C (h = 10 W/(m² K)); steel conducts with k = 50 W/(m K), expands with α = 12 × 10⁻⁶/K and is -stress free at 20 °C. One coupled temperature–displacement step (steady) on 23 641 tet10 -elements: about 8 s to open (meshing) and 13 s to solve. +stress free at 20 °C. One coupled temperature–displacement step (steady) on about 23 000 tet10 +elements (the count differs a little between platforms): about 5 s to open (meshing) and 13 s to +solve. **Look at:** Contours ▸ Temperature: the upright stays within a few degrees of the heater and the foot falls to 95 °C, a steady balance between conduction and the air. Contours ▸ Von Mises