diff --git a/ATTRIBUTION.md b/ATTRIBUTION.md index 24c62a8..f44c2d3 100644 --- a/ATTRIBUTION.md +++ b/ATTRIBUTION.md @@ -88,6 +88,7 @@ is committed and listed here with its exports. | assets/icons/mesh.layers.svg, assets/icons/mesh.layers.vectorcraft | FieldCraft contributors | Original toolbar icon drawn with VectorCraft from `assets/icons/src/mesh.layers.json` (`cargo xtask icons`) | MIT OR Apache-2.0 | | assets/icons/mesh.node.add.svg, assets/icons/mesh.node.add.vectorcraft | FieldCraft contributors | Original toolbar icon drawn with VectorCraft from `assets/icons/src/mesh.node.add.json` (`cargo xtask icons`) | MIT OR Apache-2.0 | | assets/icons/mesh.node.merge.svg, assets/icons/mesh.node.merge.vectorcraft | FieldCraft contributors | Original toolbar icon drawn with VectorCraft from `assets/icons/src/mesh.node.merge.json` (`cargo xtask icons`) | MIT OR Apache-2.0 | +| assets/icons/mesh.node.relax.svg, assets/icons/mesh.node.relax.vectorcraft | FieldCraft contributors | Original toolbar icon drawn with VectorCraft from `assets/icons/src/mesh.node.relax.json` (`cargo xtask icons`) | MIT OR Apache-2.0 | | assets/icons/mesh.order.quadratic.svg, assets/icons/mesh.order.quadratic.vectorcraft | FieldCraft contributors | Original toolbar icon drawn with VectorCraft from `assets/icons/src/mesh.order.quadratic.json` (`cargo xtask icons`) | MIT OR Apache-2.0 | | assets/icons/mesh.quad.faces.svg, assets/icons/mesh.quad.faces.vectorcraft | FieldCraft contributors | Original toolbar icon drawn with VectorCraft from `assets/icons/src/mesh.quad.faces.json` (`cargo xtask icons`) | MIT OR Apache-2.0 | | assets/icons/mesh.quality.svg, assets/icons/mesh.quality.vectorcraft | FieldCraft contributors | Original toolbar icon drawn with VectorCraft from `assets/icons/src/mesh.quality.json` (`cargo xtask icons`) | MIT OR Apache-2.0 | diff --git a/assets/icons/mesh.node.relax.svg b/assets/icons/mesh.node.relax.svg new file mode 100644 index 0000000..f95c15f --- /dev/null +++ b/assets/icons/mesh.node.relax.svg @@ -0,0 +1,10 @@ + + + mesh.node.relax + + + + + + + diff --git a/assets/icons/mesh.node.relax.vectorcraft b/assets/icons/mesh.node.relax.vectorcraft new file mode 100644 index 0000000..7261516 Binary files /dev/null and b/assets/icons/mesh.node.relax.vectorcraft differ diff --git a/assets/icons/src/mesh.node.relax.json b/assets/icons/src/mesh.node.relax.json new file mode 100644 index 0000000..4fdd68e --- /dev/null +++ b/assets/icons/src/mesh.node.relax.json @@ -0,0 +1 @@ +{"shapes": [{"circle": [12, 12, 9], "stroke": "line", "width": 1.25, "dash": [2, 2]}, {"d": "M5 12 L19 12 M12 5 L12 19", "stroke": "line", "width": 1.0}, {"d": "M7.5 8 L16.5 8 M7.5 16 L16.5 16 M8 7.5 L8 16.5 M16 7.5 L16 16.5", "stroke": "line", "width": 1.0}, {"circle": [12, 12, 2.4], "fill": "accent", "stroke": "accent"}]} diff --git a/crates/engine/src/catalogue.rs b/crates/engine/src/catalogue.rs index 5ca0c50..de65303 100644 --- a/crates/engine/src/catalogue.rs +++ b/crates/engine/src/catalogue.rs @@ -365,6 +365,14 @@ pub static CATALOGUE: &[Capability] = &[ "Merge coincident nodes, keeping sets, loads and boundary conditions on the survivors; check for duplicates and cracks, renumber", &["mesh.check.duplicates", "mesh.renumber"], ), + ready_more( + "mesh.node.relax", + "Relax Nodes", + Mesh, + "Improve", + "Relax the nodes inside a sphere to fix spikes and slivers; surface nodes slide on the geometry; kept as a field to apply again after a remesh", + &["mesh.node.relax.apply"], + ), soon("mesh.region.remesh", "Remesh Region", Mesh, "Edit", "M6", "Remesh a region in place; split and collapse elements"), soon("mesh.import.orphan", "Orphan Mesh", Mesh, "Edit", "M6", "Import a mesh without geometry"), ready("mesh.quality", "Quality", Mesh, "Verify", "Quality metrics and the certified validity check, failing elements listed"), diff --git a/crates/engine/src/cmd/forms.rs b/crates/engine/src/cmd/forms.rs index 3098d7a..beb3c57 100644 --- a/crates/engine/src/cmd/forms.rs +++ b/crates/engine/src/cmd/forms.rs @@ -1264,6 +1264,29 @@ pub const NODE_MERGE: &[Field] = &[ pub const CHECK_DUPLICATES: &[Field] = &[MERGE_TOLERANCE]; +pub const NODE_RELAX: &[Field] = &[ + xyz("center", "Centre", D::Length).help("The centre of the sphere of influence; click the mesh to place it.").hint("click the mesh"), + positive("radius", "Radius", D::Length).help("Nodes closer to the centre than this move; Alt + scroll over the view changes it."), + int("around_element", "Around element", 1, 99_999_999) + .help("Instead of a centre: a sphere round this element (radius default twice its longest edge).") + .hint("none"), + labels("nodes", "Nodes", LabelOf::Node).help("Instead of a sphere: these nodes.").hint("none"), + references("bodies", "Bodies", C::Body).help("Only the nodes of these bodies.").hint("every body"), + int("iterations", "Iterations", 1, 100).help("Sweeps over the nodes.").def("5"), + range("strength", "Strength", D::None, 0.0, 1.0).help("How far each node moves toward its relaxed position, 0 to 1.").def("1"), + choice("falloff", "Falloff", &[("smooth", "Smooth"), ("linear", "Linear"), ("none", "None")]) + .help("How the strength falls from the centre to the radius, so the edit blends into the mesh.") + .def("\"smooth\""), + text("field", "Keep as").help("The relaxation field to keep, to apply again after a remesh.").hint("Relax1, Relax2, …"), + flag("reapply", "Reapply after remesh").help("Apply the field again by itself after each remesh.").def("false"), +]; + +pub const NODE_RELAX_APPLY: &[Field] = &[ + text("field", "Field").req().help("The relaxation field kept by Relax Nodes."), + references("bodies", "Bodies", C::Body).help("Only the nodes of these bodies.").hint("the field's bodies"), + range("polish", "Polish", D::None, 0.0, 1.0).help("Strength of one guarded smoothing sweep after the field is applied; 0 for none.").def("0.3"), +]; + pub const RENUMBER: &[Field] = &[ choice("what", "Renumber", &[("both", "Nodes and elements"), ("nodes", "Nodes"), ("elements", "Elements")]) .help("Which labels change.") diff --git a/crates/engine/src/cmd/mod.rs b/crates/engine/src/cmd/mod.rs index 64d3c02..c73d471 100644 --- a/crates/engine/src/cmd/mod.rs +++ b/crates/engine/src/cmd/mod.rs @@ -18,10 +18,12 @@ pub use datum::edit_params as datum_edit_params; pub use datum::preview as datum_preview; pub use mesh::analytic_surfaces; pub use nodes::merge_preview as node_merge_preview; +pub use relax::{Preview as RelaxPreview, preview as relax_preview, sphere as relax_sphere}; mod mirror; mod model; pub(crate) mod nodes; pub(crate) mod octree; +pub(crate) mod relax; pub mod replace; mod results; pub(crate) mod seeds; @@ -86,6 +88,8 @@ pub static ALL: &[CommandSpec] = &[ nodes::MERGE, nodes::DUPLICATES, nodes::RENUMBER, + relax::RELAX, + relax::APPLY, tags::TAG, tags::TAG_CLEAR, model::UNITS, diff --git a/crates/engine/src/cmd/relax.rs b/crates/engine/src/cmd/relax.rs new file mode 100644 index 0000000..cc65039 --- /dev/null +++ b/crates/engine/src/cmd/relax.rs @@ -0,0 +1,540 @@ +//! Node relaxation (#398): `mesh.node.relax` moves the nodes inside a sphere of influence (or +//! listed nodes, elements or a set) to fix the worst elements there, nodes on the geometry +//! staying on it ([`fieldcraft_mesh::relax`]), and keeps what it did as a relaxation field of the +//! model; `mesh.node.relax.apply` does a field again on the current mesh (after a remesh), and a +//! field marked `reapply` is applied by itself after every remesh of its bodies. +//! +//! **Where a node may go** comes from the mesh's geometry tags: a node on one B-rep face slides +//! in that face's parameters, a node where two faces meet slides along their edge, a node on a +//! geometry vertex or on three or more faces stays. Faces of the same surface (the two halves a +//! kernel splits a cylinder into, or one face as cut into cells) count once. Partition cut faces +//! (`partition::a` and `:b`) are not geometry: a node on one moves freely, and the cells on +//! both sides share it, so the mesh stays conforming across the cut. A node on the skin whose +//! geometry is unknown (an imported mesh) stays. + +use std::collections::{BTreeMap, BTreeSet}; + +use fieldcraft_mesh::relax::{self, Constraints, Falloff, Freedom, RegionQuality, Relaxed}; +use fieldcraft_mesh::surface::{BoundarySurface, SurfaceOptions}; +use fieldcraft_mesh::tet::geometry::{AnalyticSurface, Tessellation}; +use fieldcraft_model::geometry::Entity; +use fieldcraft_model::{Model, Relaxation}; +use serde_json::{Value, json}; + +use crate::requirement::Requirement; +use crate::schema::ProblemKind; +use crate::{CommandKind, CommandSpec, EngineError, Result, Session}; + +pub const RELAX: CommandSpec = CommandSpec { + id: "mesh.node.relax", + label: "Relax Nodes", + params: r#"{"center": [x, y, z], "radius": mm (the sphere of influence) | "around_element": label (a sphere round that element; radius default twice its longest edge) | "nodes": [labels] | "elements": [labels] | "set": "Name", "bodies"?: ["Body1"] (only their nodes), "iterations"?: 5 (1..=100), "strength"?: 1 (0..1], "falloff"?: "smooth" (default: weight 1 at the centre to 0 at the radius) | "linear" | "none", "field"?: "Name" (the relaxation field to keep; default Relax1, Relax2, …), "reapply"?: false (apply the field again after each remesh of these bodies)} — moves the nodes to improve the worst elements around them (smart Laplacian, then optimization of the worst element; a move is kept only when the worst element around the node does not get worse, and never inverts one). Nodes on a B-rep face slide in its UV parameters and stay on it, nodes on a B-rep edge slide along it, nodes on vertices stay; nodes on partition cut faces move freely (shared by both cells); mid-side nodes follow, on the geometry. Reports the nodes moved, the largest displacement, the region's smallest scaled Jacobian, largest edge ratio and worst elements before and after, and the field kept. One undo step"#, + form: Some(super::forms::NODE_RELAX), + kind: CommandKind::Model, + requires: &[Requirement::Mesh], + run: relax_nodes, +}; + +pub const APPLY: CommandSpec = CommandSpec { + id: "mesh.node.relax.apply", + label: "Apply Relaxation", + params: r#"{"field": "Name", "bodies"?: ["Body1"] (default the field's), "polish"?: 0.3 (strength of one guarded smoothing sweep after; 0 for none)} — applies a relaxation field kept by mesh.node.relax to the current mesh (typically after a remesh): each node inside the relaxed region takes the displacement interpolated from the field at its parametric coordinates in the old cells, on its face or edge as mesh.node.relax moves it; a move that would invert an element is undone; then one guarded smoothing sweep. Reports the nodes moved and the largest displacement"#, + form: Some(super::forms::NODE_RELAX_APPLY), + kind: CommandKind::Model, + requires: &[Requirement::Mesh], + run: apply, +}; + +/// The commands after which `reapply` fields are applied again. +const REMESH: &[&str] = &[ + "mesh.tet.generate", + "mesh.hex.generate", + "mesh.hex.sweep", + "mesh.hex.mapped", + "mesh.hex.blocks", + "mesh.hex.octree", + "mesh.cells", + "mesh.block", +]; + +const DEFAULT_ITERATIONS: usize = 5; + +/// The region's quality as JSON. +fn quality_json(q: &RegionQuality) -> Value { + json!({ + "elements": q.elements, + "min_scaled_jacobian": q.min_quality, + "max_edge_ratio": if q.max_edge_ratio.is_finite() { json!(q.max_edge_ratio) } else { Value::Null }, + "invalid": q.invalid, + "worst": q.worst.iter().map(|(e, v)| json!({"element": e, "quality": v})).collect::>(), + }) +} + +/// Faces a node lies on (geometry references), from the element faces, shells and vertex tags. +struct NodeTags { + faces: BTreeMap>, + vertex: BTreeSet, + /// Nodes on a partition cut face. + cut: BTreeSet, +} + +fn node_tags(model: &Model, nodes: &BTreeSet) -> NodeTags { + let c = model.classification(); + let mut t = NodeTags { faces: BTreeMap::new(), vertex: BTreeSet::new(), cut: BTreeSet::new() }; + let note = |n: u32, tag: u32, t: &mut NodeTags| { + let Some(g) = c.entity(tag) else { return }; + if let Entity::Face(name) = &g.entity { + if name.starts_with("partition:") { + t.cut.insert(n); + } else { + t.faces.entry(n).or_default().insert((g.body.clone(), name.clone())); + } + } + }; + for (&label, tags) in c.elements() { + let Some(e) = model.element(label) else { continue }; + if !e.nodes.iter().any(|n| nodes.contains(n)) { + continue; + } + if e.kind.dimension() == 2 + && let Some(r) = tags.region + { + for &n in &e.nodes { + note(n, r, &mut t); + } + } + for (&face, &tag) in &tags.faces { + let Some(local) = e.kind.faces().get(usize::from(face).wrapping_sub(1)) else { continue }; + for &i in local.iter() { + if let Some(&n) = e.nodes.get(i) { + note(n, tag, &mut t); + } + } + } + } + for (&n, &tag) in c.nodes() { + if nodes.contains(&n) && c.entity(tag).is_some_and(|g| matches!(g.entity, Entity::Vertex(_))) { + t.vertex.insert(n); + } + } + t +} + +/// Whether two exact surfaces are the same surface. +fn same_surface(a: &AnalyticSurface, b: &AnalyticSurface, tol: f64) -> bool { + use AnalyticSurface as S; + let sub = |x: [f64; 3], y: [f64; 3]| [x[0] - y[0], x[1] - y[1], x[2] - y[2]]; + let dot = |x: [f64; 3], y: [f64; 3]| x[0] * y[0] + x[1] * y[1] + x[2] * y[2]; + let len = |x: [f64; 3]| dot(x, x).sqrt(); + let unit = |x: [f64; 3]| { + let l = len(x); + if l > 0.0 { [x[0] / l, x[1] / l, x[2] / l] } else { x } + }; + let parallel = |x: [f64; 3], y: [f64; 3]| (dot(unit(x), unit(y)).abs() - 1.0).abs() < 1e-9; + let on_axis = |p: [f64; 3], o: [f64; 3], a: [f64; 3]| { + let d = sub(p, o); + let a = unit(a); + let along = dot(d, a); + len([d[0] - a[0] * along, d[1] - a[1] * along, d[2] - a[2] * along]) <= tol + }; + match (*a, *b) { + (S::Plane { normal: n1, point: p1 }, S::Plane { normal: n2, point: p2 }) => parallel(n1, n2) && dot(sub(p2, p1), unit(n1)).abs() <= tol, + (S::Cylinder { origin: o1, axis: a1, radius: r1 }, S::Cylinder { origin: o2, axis: a2, radius: r2 }) => { + (r1 - r2).abs() <= tol && parallel(a1, a2) && on_axis(o2, o1, a1) + } + (S::Sphere { center: c1, radius: r1 }, S::Sphere { center: c2, radius: r2 }) => len(sub(c1, c2)) <= tol && (r1 - r2).abs() <= tol, + (S::Cone { apex: p1, axis: a1, half_angle: h1 }, S::Cone { apex: p2, axis: a2, half_angle: h2 }) => { + len(sub(p1, p2)) <= tol && dot(unit(a1), unit(a2)) > 1.0 - 1e-9 && (h1 - h2).abs() <= 1e-9 + } + (S::Torus { center: c1, axis: a1, major: m1, minor: n1 }, S::Torus { center: c2, axis: a2, major: m2, minor: n2 }) => { + len(sub(c1, c2)) <= tol && parallel(a1, a2) && (m1 - m2).abs() <= tol && (n1 - n2).abs() <= tol + } + _ => false, + } +} + +/// The design's geometry as relaxation constraints. +struct DesignGeometry { + /// Per body: its tessellation with exact surfaces. + bodies: Vec, + /// Face ids: (body, B-rep face). + faces: Vec<(usize, u32)>, + /// Curve ids: (body, the two B-rep faces). + curves: Vec<(usize, [u32; 2])>, + freedom: BTreeMap, +} + +impl Constraints for DesignGeometry { + fn freedom(&self, node: u32) -> Option { + self.freedom.get(&node).copied() + } + fn slide_face(&self, face: u32, p: [f64; 3], d: [f64; 3]) -> Option<[f64; 3]> { + let &(b, f) = self.faces.get(face as usize)?; + self.bodies.get(b)?.slide_face(f, p, d) + } + fn slide_curve(&self, curve: u32, p: [f64; 3], d: [f64; 3]) -> Option<[f64; 3]> { + let &(b, f) = self.curves.get(curve as usize)?; + self.bodies.get(b)?.slide_curve(f, p, d) + } +} + +/// The constraints of `nodes` from the mesh's tags and the open design. +fn design_geometry(s: &Session, nodes: &BTreeSet) -> Result { + let tags = node_tags(&s.model, nodes); + let mut g = DesignGeometry { bodies: Vec::new(), faces: Vec::new(), curves: Vec::new(), freedom: BTreeMap::new() }; + // Body name → (index into `g.bodies`, its exact surfaces), built when first needed. + let mut loaded: BTreeMap>)>> = BTreeMap::new(); + let mut load = |name: &str, g: &mut DesignGeometry| -> Result>)>> { + if let Some(l) = loaded.get(name) { + return Ok(l.clone()); + } + let shape = s.design.as_deref().and_then(|d| d.body(name)); + let entry = match shape { + None => None, + Some(shape) => { + let tol = (shape.size()? * 1e-3).max(1e-6); + let mesh = shape.tessellate(tol)?; + let surface = BoundarySurface::new( + &mesh.positions, + &mesh.triangles, + &mesh.tri_face, + &SurfaceOptions { weld_distance: mesh.weld_distance, ..SurfaceOptions::default() }, + )?; + let analytic = super::mesh::analytic_surfaces(shape); + g.bodies.push(Tessellation::new(&surface, analytic.clone())); + Some((g.bodies.len() - 1, analytic)) + } + }; + loaded.insert(name.to_string(), entry.clone()); + Ok(entry) + }; + let tol = 1e-6 * s.model.nodes().values().fold(0.0f64, |m, n| m.max(n.xyz.iter().fold(0.0f64, |a, x| a.max(x.abs())))).max(1.0); + let mut face_ids: BTreeMap<(usize, u32), u32> = BTreeMap::new(); + let mut curve_ids: BTreeMap<(usize, [u32; 2]), u32> = BTreeMap::new(); + for &n in nodes { + if tags.vertex.contains(&n) { + g.freedom.insert(n, Freedom::Fixed); + continue; + } + let Some(refs) = tags.faces.get(&n) else { + if tags.cut.contains(&n) { + g.freedom.insert(n, Freedom::Free); + } + continue; + }; + // Resolve each face; group those on one surface. + let mut groups: Vec)>> = Vec::new(); + let mut unknown = false; + for (body, face) in refs { + let Some((b, analytic)) = load(body, &mut g)? else { + unknown = true; + continue; + }; + let Some(f) = s.design.as_deref().and_then(|d| d.body(body)).and_then(|sh| sh.face_index(face)) else { + unknown = true; + continue; + }; + let surface = analytic.get(f as usize).copied().flatten(); + let group = groups.iter_mut().find(|grp| { + grp.iter().any(|(gb, gf, gs)| (*gb == b && *gf == f) || matches!((gs, surface), (Some(x), Some(y)) if same_surface(x, &y, tol))) + }); + match group { + Some(grp) => grp.push((b, f, surface)), + None => groups.push(vec![(b, f, surface)]), + } + } + let freedom = match groups.as_slice() { + _ if unknown => Freedom::Fixed, + [one] => { + let Some(&(b, f, _)) = one.first() else { continue }; + let next = g.faces.len() as u32; + let id = *face_ids.entry((b, f)).or_insert(next); + if id == next { + g.faces.push((b, f)); + } + Freedom::Face(id) + } + [x, y] => { + // A body with a face of each surface carries their edge. + let pair = x.iter().find_map(|&(bx, fx, _)| y.iter().find(|(by, _, _)| *by == bx).map(|&(_, fy, _)| (bx, [fx.min(fy), fx.max(fy)]))); + match pair { + Some(key) => { + let next = g.curves.len() as u32; + let id = *curve_ids.entry(key).or_insert(next); + if id == next { + g.curves.push(key); + } + Freedom::Curve(id) + } + None => Freedom::Fixed, + } + } + _ => Freedom::Fixed, + }; + g.freedom.insert(n, freedom); + } + Ok(g) +} + +/// The nodes of `bodies`' elements, when given. +fn body_nodes(model: &Model, p: &Value) -> Result, BTreeSet)>> { + let Some(bodies) = super::opt::>(p, "bodies")? else { return Ok(None) }; + let mut nodes = BTreeSet::new(); + for b in &bodies { + let elements = model.body_elements(b); + if elements.is_empty() { + return Err(EngineError::Params(format!("body `{}` has no mesh", b.chars().take(80).collect::()))); + } + for e in elements { + if let Some(el) = model.element(e) { + nodes.extend(el.nodes.iter().copied()); + } + } + } + Ok(Some((bodies, nodes))) +} + +fn point(p: &Value, key: &str) -> Result> { + match super::opt::<[f64; 3]>(p, key)? { + Some(c) if c.iter().all(|x| x.is_finite()) => Ok(Some(c)), + Some(_) => Err(EngineError::param(key, ProblemKind::NotANumber)), + None => Ok(None), + } +} + +/// The sphere around element `label`: its centroid, and twice its longest corner edge. +fn around_element(model: &Model, label: u32) -> Result<([f64; 3], f64)> { + let e = model.element(label).ok_or_else(|| EngineError::Params(format!("no element {label}")))?; + let x: Vec<[f64; 3]> = e.nodes.iter().take(e.kind.corner_count()).filter_map(|n| model.node(*n).map(|p| p.xyz)).collect(); + if x.is_empty() { + return Err(EngineError::Params(format!("element {label} has no nodes"))); + } + let c = x.iter().fold([0.0; 3], |s, p| [s[0] + p[0], s[1] + p[1], s[2] + p[2]]).map(|v| v / x.len() as f64); + let mut longest = 0.0f64; + for edge in e.kind.edges() { + let at = |i: usize| edge.get(i).and_then(|&k| x.get(k)).copied(); + if let (Some(a), Some(b)) = (at(0), at(1)) { + longest = longest.max(((a[0] - b[0]).powi(2) + (a[1] - b[1]).powi(2) + (a[2] - b[2]).powi(2)).sqrt()); + } + } + Ok((c, 2.0 * longest)) +} + +/// The nodes to move and their weights, and the sphere when there is one. +pub(crate) fn weights(model: &Model, p: &Value) -> Result<(BTreeMap, Option<([f64; 3], f64)>)> { + let strength = super::optional_number(p, "strength", 1.0)?; + if !(strength > 0.0 && strength <= 1.0) { + return Err(EngineError::param("strength", ProblemKind::NotInRange { min: 0.0, max: 1.0 })); + } + let falloff = match super::opt::(p, "falloff")? { + None => Falloff::Smooth, + Some(f) => Falloff::from_name(&f).ok_or_else(|| EngineError::param("falloff", ProblemKind::NotAChoice))?, + }; + let only = body_nodes(model, p)?.map(|(_, n)| n); + let around = match super::opt::(p, "around_element")? { + Some(e) => Some(around_element(model, e)?), + None => None, + }; + let center = point(p, "center")?.or(around.map(|a| a.0)); + let radius = match p.get("radius") { + Some(_) => Some(super::optional_number(p, "radius", 0.0)?) + .filter(|r| *r > 0.0 && r.is_finite()) + .ok_or_else(|| EngineError::param("radius", ProblemKind::NotPositive)) + .map(Some)?, + None => around.map(|a| a.1), + }; + let listed = ["nodes", "elements", "set"].iter().filter(|k| p.get(**k).is_some()).count(); + if let (Some(c), Some(r)) = (center, radius) { + if listed > 0 { + return Err(EngineError::Params( + "give a sphere (`center` and `radius`, or `around_element`) or `nodes`, `elements` or `set`, not both".into(), + )); + } + return Ok((relax::sphere_weights(model, c, r, strength, falloff, only.as_ref())?, Some((c, r)))); + } + if center.is_some() != radius.is_some() { + return Err(EngineError::Params("a sphere of influence needs both `center` and `radius`".into())); + } + let mut nodes: BTreeSet = BTreeSet::new(); + if p.get("nodes").is_some() { + nodes.extend(super::labels(p, "nodes", fieldcraft_model::MAX_NODES)?); + } + let mut elements: Vec = Vec::new(); + if p.get("elements").is_some() { + elements.extend(super::labels(p, "elements", fieldcraft_model::MAX_ELEMENTS)?); + } + if let Some(name) = super::opt::(p, "set")? { + let set = model.resolve_set(&name)?; + nodes.extend(set.nodes.iter().copied()); + elements.extend(set.elements.iter().copied()); + elements.extend(set.faces.iter().map(|f| f.element)); + } + for e in elements { + let el = model.element(e).ok_or_else(|| EngineError::Params(format!("no element {e}")))?; + nodes.extend(el.nodes.iter().copied()); + } + if listed == 0 { + return Err(EngineError::Params("say which nodes to relax: `center` and `radius`, `around_element`, `nodes`, `elements` or `set`".into())); + } + if let Some(only) = &only { + nodes.retain(|n| only.contains(n)); + } + Ok((relax::uniform_weights(model, &nodes, strength)?, None)) +} + +fn iterations(p: &Value) -> Result { + match p.get("iterations") { + None => Ok(DEFAULT_ITERATIONS), + Some(v) => v + .as_u64() + .and_then(|i| usize::try_from(i).ok()) + .filter(|i| (1..=relax::MAX_ITERATIONS).contains(i)) + .ok_or_else(|| EngineError::param("iterations", ProblemKind::NotInRange { min: 1.0, max: relax::MAX_ITERATIONS as f64 })), + } +} + +/// The nodes the relaxation may touch: the weighted ones and their elements' nodes. +fn touched(model: &Model, weights: &BTreeMap) -> BTreeSet { + let mut out: BTreeSet = weights.keys().copied().collect(); + for e in model.elements().values() { + if e.nodes.iter().any(|n| weights.contains_key(n)) { + out.extend(e.nodes.iter().copied()); + } + } + out +} + +/// What `mesh.node.relax` with `params` would do, without doing it: the open form's preview. +pub struct Preview { + pub sphere: Option<([f64; 3], f64)>, + pub relaxed: Relaxed, +} + +pub fn preview(s: &Session, params: &Value) -> Option { + let (w, sphere) = weights(&s.model, params).ok()?; + let g = design_geometry(s, &touched(&s.model, &w)).ok()?; + let relaxed = relax::relax(&s.model, &w, iterations(params).ok()?, &g).ok()?; + Some(Preview { sphere, relaxed }) +} + +/// The sphere `params` describe, for the form's gizmo before the rest is valid. +pub fn sphere(s: &Session, params: &Value) -> Option<([f64; 3], f64)> { + let around = params.get("around_element").and_then(Value::as_u64).and_then(|e| around_element(&s.model, u32::try_from(e).ok()?).ok()); + let c = point(params, "center").ok().flatten().or(around.map(|a| a.0))?; + let r = params.get("radius").and_then(Value::as_f64).filter(|r| *r > 0.0 && r.is_finite()).or(around.map(|a| a.1))?; + Some((c, r)) +} + +fn relax_nodes(s: &mut Session, p: &Value) -> Result { + let (w, sphere) = weights(&s.model, p)?; + if w.is_empty() { + return Err(EngineError::Params("no mesh nodes in the region to relax".into())); + } + let iterations = iterations(p)?; + let g = design_geometry(s, &touched(&s.model, &w))?; + let r = relax::relax(&s.model, &w, iterations, &g)?; + let field_name = match super::opt::(p, "field")? { + Some(n) => n, + None => (1..).map(|i| format!("Relax{i}")).find(|n| s.model.relaxation(n).is_none()).unwrap_or_else(|| "Relax".into()), + }; + let bodies = match body_nodes(&s.model, p)? { + Some((b, _)) => b, + None => Vec::new(), + }; + let reapply = super::opt::(p, "reapply")?.unwrap_or(false); + if !r.positions.is_empty() { + let cells = relax::field_cells(&s.model, &r); + s.model.set_relaxation(&field_name, Relaxation { bodies, cells, reapply })?; + } + for (n, x) in &r.positions { + s.model.set_node(*n, *x)?; + } + let summary = format!( + "moved {} node{} (largest {:.3e}); scaled Jacobian min {:.3} → {:.3}{}", + r.moved, + if r.moved == 1 { "" } else { "s" }, + r.max_displacement, + r.before.min_quality, + r.after.min_quality, + if r.positions.is_empty() { String::new() } else { format!("; kept as field {field_name}") } + ); + Ok(json!({ + "moved": r.moved, + "nodes_changed": r.positions.len(), + "max_displacement": r.max_displacement, + "iterations": r.iterations, + "sphere": sphere.map(|(c, r)| json!({"center": c, "radius": r})), + "before": quality_json(&r.before), + "after": quality_json(&r.after), + "field": (!r.positions.is_empty()).then_some(field_name), + "summary": summary, + })) +} + +/// Apply field `name` to the nodes of `bodies` (or of the whole mesh). +fn apply_field(s: &mut Session, name: &str, bodies: &[String], polish: f64) -> Result { + let field = s + .model + .relaxation(name) + .cloned() + .ok_or_else(|| EngineError::Params(format!("no relaxation field `{}`", name.chars().take(80).collect::())))?; + let mut nodes = BTreeSet::new(); + let bodies = if bodies.is_empty() { &field.bodies[..] } else { bodies }; + if bodies.is_empty() { + nodes.extend(s.model.nodes().keys().copied()); + } else { + for b in bodies { + for e in s.model.body_elements(b) { + if let Some(el) = s.model.element(e) { + nodes.extend(el.nodes.iter().copied()); + } + } + } + } + // Only nodes the field reaches need constraints. + nodes.retain(|n| s.model.node(*n).is_some_and(|p| relax::field_displacement(&field.cells, p.xyz).is_some())); + let g = design_geometry(s, &nodes)?; + let r = relax::apply_field(&s.model, &field.cells, Some(&nodes), &g, polish)?; + for (n, x) in &r.positions { + s.model.set_node(*n, *x)?; + } + Ok(r) +} + +fn apply(s: &mut Session, p: &Value) -> Result { + let name: String = super::req(p, "field")?; + let bodies = super::opt::>(p, "bodies")?.unwrap_or_default(); + let polish = super::optional_number(p, "polish", 0.3)?; + if !(0.0..=1.0).contains(&polish) { + return Err(EngineError::param("polish", ProblemKind::NotInRange { min: 0.0, max: 1.0 })); + } + let r = apply_field(s, &name, &bodies, polish)?; + Ok(json!({ + "moved": r.moved, + "nodes_changed": r.positions.len(), + "max_displacement": r.max_displacement, + "summary": format!("applied field {name}: moved {} node{} (largest {:.3e})", r.moved, if r.moved == 1 { "" } else { "s" }, r.max_displacement), + })) +} + +/// After a remesh (`id`), apply again every field marked `reapply`; its reply lists them. +pub(crate) fn after_remesh(s: &mut Session, id: &str, mut reply: Value) -> Result { + if !REMESH.contains(&id) { + return Ok(reply); + } + let fields: Vec = s.model.relaxations().iter().filter(|(_, f)| f.reapply).map(|(n, _)| n.clone()).collect(); + let mut applied = Vec::new(); + for name in fields { + let r = apply_field(s, &name, &[], 0.3)?; + applied.push(json!({"field": name, "moved": r.moved})); + } + if !applied.is_empty() + && let Some(o) = reply.as_object_mut() + { + o.insert("relaxations_applied".into(), Value::Array(applied)); + } + Ok(reply) +} + +#[cfg(test)] +mod tests; diff --git a/crates/engine/src/cmd/relax/tests.rs b/crates/engine/src/cmd/relax/tests.rs new file mode 100644 index 0000000..46dd1cc --- /dev/null +++ b/crates/engine/src/cmd/relax/tests.rs @@ -0,0 +1,194 @@ +use super::*; + +/// The plate (40 × 30 × 10) with its Ø10 hole round (15, 15): a point on the hole's wall. +const ON_HOLE: [f64; 3] = [20.0, 15.0, 5.0]; +const HOLE: ([f64; 2], f64) = ([15.0, 15.0], 5.0); + +fn plate() -> Session { + super::super::datum::tests::plate().0 +} + +fn hole_radius(p: [f64; 3]) -> f64 { + ((p[0] - HOLE.0[0]).powi(2) + (p[1] - HOLE.0[1]).powi(2)).sqrt() +} + +/// Which analytic face a node of the plate lies on, if any. +fn on_faces(p: [f64; 3]) -> Vec<&'static str> { + let mut f = Vec::new(); + if (hole_radius(p) - HOLE.1).abs() < 1e-6 { + f.push("hole"); + } + if p[2].abs() < 1e-6 || (p[2] - 10.0).abs() < 1e-6 { + f.push("cap"); + } + f +} + +/// Spikes near the hole: interior nodes pushed about; nodes on the hole turned round it and +/// moved along its axis (still on the cylinder), nodes on the caps moved in their plane. +fn perturb(s: &mut Session, centre: [f64; 3], radius: f64, amount: f64) { + let nodes: Vec<(u32, [f64; 3])> = s.model.nodes().iter().map(|(n, p)| (*n, p.xyz)).collect(); + for (k, (n, p)) in nodes.into_iter().enumerate() { + let d2: f64 = (0..3).map(|i| (p[i] - centre[i]).powi(2)).sum(); + if d2 > radius * radius { + continue; + } + let wiggle = |m: usize| amount * ((k * m % 7) as f64 - 3.0) / 3.0; + let faces = on_faces(p); + let q = match faces.as_slice() { + [] => [p[0] + wiggle(3), p[1] + wiggle(5), p[2] + wiggle(11)], + ["hole"] => { + let a = wiggle(3) / HOLE.1; + let (x, y) = (p[0] - HOLE.0[0], p[1] - HOLE.0[1]); + [HOLE.0[0] + x * a.cos() - y * a.sin(), HOLE.0[1] + x * a.sin() + y * a.cos(), p[2] + wiggle(5)] + } + ["cap"] if hole_radius(p) > HOLE.1 + 1e-6 => [p[0] + wiggle(3), p[1] + wiggle(5), p[2]], + _ => continue, + }; + s.model.set_node(n, q).unwrap(); + } +} + +fn check_on_geometry(before: &Model, after: &Model) { + for (n, p) in after.nodes() { + let old = before.node(*n).unwrap().xyz; + for f in on_faces(old) { + match f { + "hole" => assert!((hole_radius(p.xyz) - HOLE.1).abs() < 1e-6, "node {n} left the hole: r = {}", hole_radius(p.xyz)), + _ => assert!((p.xyz[2] - old[2]).abs() < 1e-6, "node {n} left its cap"), + } + } + } +} + +fn invalid(s: &mut Session) -> u64 { + s.execute("mesh.quality", &json!({})).unwrap()["invalid"]["count"].as_u64().unwrap_or(0) +} + +/// The issue's hex test on a design: spikes near the hole relax, every surface node stays on its +/// face (the hole's cylinder exactly), quality rises, nothing inverts, nothing outside the sphere +/// moves; one undo step; the field is kept and saved with the model. +#[test] +fn hex_spikes_near_a_hole_relax_on_the_geometry() { + let mut s = plate(); + s.execute("mesh.hex.generate", &json!({"size": 2.0})).unwrap(); + perturb(&mut s, ON_HOLE, 3.0, 0.35); + let perturbed = s.model.clone(); + let r = s.execute("mesh.node.relax", &json!({"center": ON_HOLE, "radius": 4.0, "iterations": 10})).unwrap(); + let (before, after) = (r["before"]["min_scaled_jacobian"].as_f64().unwrap(), r["after"]["min_scaled_jacobian"].as_f64().unwrap()); + assert!(after > before + 0.05, "{r}"); + assert_eq!(r["after"]["invalid"], 0, "{r}"); + assert!(r["moved"].as_u64().unwrap() > 0 && r["field"] == "Relax1", "{r}"); + check_on_geometry(&perturbed, &s.model); + for (n, p) in s.model.nodes() { + let old = perturbed.node(*n).unwrap().xyz; + if old != p.xyz { + let d2: f64 = (0..3).map(|i| (old[i] - ON_HOLE[i]).powi(2)).sum(); + assert!(d2 < 16.0, "node {n} outside the sphere moved"); + } + } + assert_eq!(invalid(&mut s), 0); + // The field is saved with the model. + let text = serde_json::to_string(&s.model).unwrap(); + let back: Model = serde_json::from_str(&text).unwrap(); + assert_eq!(back.relaxation("Relax1"), s.model.relaxation("Relax1")); + s.execute("edit.undo", &json!({})).unwrap(); + assert_eq!(s.model, perturbed); +} + +/// The same on tet10 with curved mid-side nodes on the hole: they stay on the cylinder. +#[test] +fn tet10_relaxes_with_curved_mid_side_nodes_on_the_hole() { + let mut s = plate(); + s.execute("mesh.tet.generate", &json!({"size": 2.5, "order": "quadratic"})).unwrap(); + perturb(&mut s, ON_HOLE, 3.0, 0.2); + // Put the mid-side nodes back on their edges' middles, on the hole where the edge is. + let mids: Vec<(u32, u32, u32)> = s + .model + .elements() + .values() + .flat_map(|e| e.kind.edges().iter().map(|ed| (e.nodes[ed[0]], e.nodes[ed[1]], e.nodes[ed[2]])).collect::>()) + .collect(); + let original = s.model.clone(); + for (a, b, m) in mids { + let (pa, pb, pm) = (original.node(a).unwrap().xyz, original.node(b).unwrap().xyz, original.node(m).unwrap().xyz); + let mut q = [(pa[0] + pb[0]) / 2.0, (pa[1] + pb[1]) / 2.0, (pa[2] + pb[2]) / 2.0]; + if on_faces(pm).contains(&"hole") { + let r = hole_radius(q); + q = [HOLE.0[0] + (q[0] - HOLE.0[0]) * HOLE.1 / r, HOLE.0[1] + (q[1] - HOLE.0[1]) * HOLE.1 / r, q[2]]; + } + s.model.set_node(m, q).unwrap(); + } + let perturbed = s.model.clone(); + let r = s.execute("mesh.node.relax", &json!({"center": ON_HOLE, "radius": 4.0, "iterations": 8})).unwrap(); + assert!(r["after"]["min_scaled_jacobian"].as_f64().unwrap() >= r["before"]["min_scaled_jacobian"].as_f64().unwrap(), "{r}"); + assert_eq!(r["after"]["invalid"], 0, "{r}"); + check_on_geometry(&perturbed, &s.model); +} + +/// A partition cut is not geometry: nodes on it move freely, shared by the cells on both sides. +#[test] +fn nodes_on_a_partition_cut_move_freely() { + let mut s = plate(); + s.execute("geometry.datum", &json!({"name": "Cut", "by": {"type": "offset", "plane": "yz", "offset": 27}})).unwrap(); + s.execute("geometry.partition", &json!({"name": "P", "datum": "Cut"})).unwrap(); + s.execute("mesh.cells", &json!({"size": 2.0, "default": "hex"})).unwrap(); + let centre = [27.0, 15.0, 5.0]; + perturb(&mut s, centre, 2.5, 0.3); + let before = s.model.clone(); + s.execute("mesh.node.relax", &json!({"center": centre, "radius": 3.5, "iterations": 10})).unwrap(); + let moved_off = s + .model + .nodes() + .iter() + .filter(|(n, p)| { + let old = before.node(**n).unwrap().xyz; + (old[0] - 27.0).abs() < 1e-6 && on_faces(old).is_empty() && (p.xyz[0] - 27.0).abs() > 1e-6 + }) + .count(); + assert!(moved_off > 0, "nodes on the cut are free"); + assert_eq!(invalid(&mut s), 0); +} + +/// A field marked `reapply` is applied again after a remesh; `mesh.node.relax.apply` applies it +/// on demand; around the worst element works. +#[test] +fn a_kept_field_is_applied_again_after_a_remesh() { + let mut s = plate(); + s.execute("mesh.hex.generate", &json!({"size": 3.0})).unwrap(); + perturb(&mut s, ON_HOLE, 4.0, 0.3); + let r = s.execute("mesh.node.relax", &json!({"center": ON_HOLE, "radius": 5.0, "field": "Spike", "reapply": true})).unwrap(); + assert_eq!(r["field"], "Spike"); + let r = s.execute("mesh.hex.generate", &json!({"size": 1.5})).unwrap(); + assert_eq!(r["relaxations_applied"][0]["field"], "Spike", "{r}"); + assert_eq!(invalid(&mut s), 0); + let r = s.execute("mesh.node.relax.apply", &json!({"field": "Spike"})).unwrap(); + assert!(r["summary"].as_str().unwrap().contains("Spike")); + let worst = s.execute("mesh.quality", &json!({})).unwrap()["metrics"]["scaled_jacobian"]["worst"].as_u64().unwrap(); + let r = s.execute("mesh.node.relax", &json!({"around_element": worst})).unwrap(); + assert!(r["after"]["min_scaled_jacobian"].as_f64() >= r["before"]["min_scaled_jacobian"].as_f64(), "{r}"); +} + +#[test] +fn hostile_parameters_never_panic() { + let mut s = plate(); + s.execute("mesh.hex.generate", &json!({"size": 4.0})).unwrap(); + for p in [ + json!({}), + json!({"center": [0, 0, 0]}), + json!({"center": [1e300, 0, 0], "radius": 1}), + json!({"center": [0, 0, 0], "radius": -1}), + json!({"center": [0, 0, 0], "radius": 1, "nodes": [1]}), + json!({"nodes": [1], "strength": 2}), + json!({"nodes": [1], "iterations": 0}), + json!({"nodes": [1], "iterations": 1000}), + json!({"nodes": [1], "falloff": "cubic"}), + json!({"around_element": 99999999}), + json!({"set": "nope"}), + json!({"nodes": [1], "bodies": ["nope"]}), + ] { + assert!(s.execute("mesh.node.relax", &p).is_err(), "{p}"); + } + assert!(s.execute("mesh.node.relax.apply", &json!({"field": "nope"})).is_err()); + assert!(s.execute("mesh.node.relax.apply", &json!({})).is_err()); +} diff --git a/crates/engine/src/lib.rs b/crates/engine/src/lib.rs index 8810584..e55c373 100644 --- a/crates/engine/src/lib.rs +++ b/crates/engine/src/lib.rs @@ -328,7 +328,8 @@ impl Session { let scope = crash::Scope::enter(); let journaled = crash::command_begin(&scope, self, id, params); let (revision, opened) = (self.revision, self.opened); - let outcome = std::panic::catch_unwind(std::panic::AssertUnwindSafe(|| (spec.run)(self, params))); + let outcome = + std::panic::catch_unwind(std::panic::AssertUnwindSafe(|| (spec.run)(self, params).and_then(|v| cmd::relax::after_remesh(self, id, v)))); let result = match outcome { Ok(Err(EngineError::Param(p))) => Err(EngineError::Param(match spec.form { Some(form) => schema::refine(form, params, p), diff --git a/crates/ui-egui/src/context.rs b/crates/ui-egui/src/context.rs index cce5c96..70d028b 100644 --- a/crates/ui-egui/src/context.rs +++ b/crates/ui-egui/src/context.rs @@ -126,6 +126,7 @@ fn form(app: &mut ModelApp, ui: &mut egui::Ui, t: &Tokens) { let selected = (!app.pick.selection.is_empty()).then(|| app.pick.selection.summary()); // While a job runs (this form's mesher, a solve) OK waits for it. let running = app.busy().map(str::to_string); + let relax_note = app.relax_note(); let units = app.session.model.units(); let session = &app.session; // The set Use selection will make on OK is offered with the others. @@ -252,6 +253,9 @@ fn form(app: &mut ModelApp, ui: &mut egui::Ui, t: &Tokens) { } else if d.attempted && !missing.is_empty() { ui.label(RichText::new(&missing).color(t.error)); } + if let Some(note) = &relax_note { + ui.label(RichText::new(note).small().color(t.text)); + } ui.collapsing(RichText::new("What it takes").size(12.0).color(t.text_dim), |ui| { ui.label(RichText::new(spec.params).small().color(t.text_dim)); }); diff --git a/crates/ui-egui/src/lib.rs b/crates/ui-egui/src/lib.rs index 1256493..476c23a 100644 --- a/crates/ui-egui/src/lib.rs +++ b/crates/ui-egui/src/lib.rs @@ -333,6 +333,8 @@ type DatumPreview = (Option<(u64, &'static str, String)>, Option, Arc>); +/// The open Relax Nodes form's preview, by model revision and form text. +type RelaxPreview = (Option<(u64, String)>, Option>); /// A screenshot asked for over the control channel. struct Shot { @@ -355,6 +357,7 @@ pub struct ModelApp { /// The datum the open form defines, and the revision, command and form text it is for. preview: std::cell::RefCell, merge_preview: std::cell::RefCell, + relax_preview: std::cell::RefCell, /// Forms opened so far (numbers each [`ParamsDialog::opening`]). forms_opened: u64, pub viewport: Viewport, @@ -493,6 +496,7 @@ impl ModelApp { dialog: None, preview: Default::default(), merge_preview: Default::default(), + relax_preview: Default::default(), forms_opened: 0, viewport: Viewport::default(), icons: icons::IconCache::default(), @@ -1157,6 +1161,81 @@ impl ModelApp { cache.1.clone() } + /// What the open Relax Nodes form would do (the nodes' new positions and the region's + /// quality), worked out once per form text and model revision. + pub fn relax_preview(&self) -> Option> { + let d = self.dialog.as_ref().filter(|d| d.command == "mesh.node.relax")?; + let key = (self.revision, d.text.clone()); + let mut cache = self.relax_preview.borrow_mut(); + if cache.0.as_ref() != Some(&key) { + let preview = serde_json::from_str::(&d.text).ok().and_then(|p| fieldcraft_engine::cmd::relax_preview(&self.session, &p)); + *cache = (Some(key), preview.map(Arc::new)); + } + cache.1.clone() + } + + /// The open Relax Nodes form's preview in a line: nodes moved and the region's quality. + pub fn relax_note(&self) -> Option { + let p = self.relax_preview()?; + let r = &p.relaxed; + Some(format!( + "Preview: {} nodes move (largest {:.3}); scaled Jacobian min {:.3} → {:.3}, edge ratio max {:.2} → {:.2}, {} → {} invalid", + r.moved, + r.max_displacement, + r.before.min_quality, + r.after.min_quality, + r.before.max_edge_ratio, + r.after.max_edge_ratio, + r.before.invalid, + r.after.invalid + )) + } + + /// Whether the Relax Nodes form is open. + pub fn relaxing(&self) -> bool { + self.dialog.as_ref().is_some_and(|d| d.command == "mesh.node.relax") + } + + /// Put the open Relax Nodes form's sphere at `at` (a click on the mesh); a sphere without a + /// radius yet gets a twentieth of the mesh's size. + pub fn relax_place(&mut self, at: Vec3) -> bool { + let mut bounds = Aabb3::EMPTY; + for n in self.session.model.nodes().values() { + bounds.add(Vec3::from(n.xyz)); + } + let size = bounds.diagonal(); + let Some(d) = self.dialog.as_mut().filter(|d| d.command == "mesh.node.relax") else { return false }; + let mut params = serde_json::from_str::(&d.text).ok().filter(Value::is_object).unwrap_or_else(|| json!({})); + params["center"] = json!([at.x, at.y, at.z]); + if let Some(o) = params.as_object_mut() { + o.remove("around_element"); + o.remove("nodes"); + } + if params.get("radius").and_then(Value::as_f64).is_none_or(|r| !(r > 0.0)) && size > 0.0 && size.is_finite() { + params["radius"] = json!(size / 20.0); + } + d.text = params.to_string(); + d.touched.insert("center".into()); + true + } + + /// Scale the open Relax Nodes form's radius by `factor` (Alt + scroll over the view). + pub fn relax_scale_radius(&mut self, factor: f64) -> bool { + let Some(sphere) = self + .dialog + .as_ref() + .filter(|d| d.command == "mesh.node.relax") + .and_then(|d| serde_json::from_str::(&d.text).ok().and_then(|p| fieldcraft_engine::cmd::relax_sphere(&self.session, &p))) + else { + return false; + }; + let Some(d) = self.dialog.as_mut() else { return false }; + let mut params = serde_json::from_str::(&d.text).ok().filter(Value::is_object).unwrap_or_else(|| json!({})); + params["radius"] = json!(sphere.1 * factor); + d.text = params.to_string(); + true + } + /// The datum the open datum form defines so far, drawn before OK. /// Worked out once per form text and model revision (a face centre costs a design query). fn datum_preview(&self) -> Option { @@ -1403,7 +1482,8 @@ impl ModelApp { let set = self.selected_set(); let cut = self.partition_preview(); let merge = self.merge_preview(); - let selected = hash(&|h| (scene, self.pick.generation, &set, &cut, &*merge).hash(h)); + let relax_key = self.dialog.as_ref().filter(|d| d.command == "mesh.node.relax").map(|d| d.text.clone()); + let selected = hash(&|h| (scene, self.pick.generation, &set, &cut, &*merge, &relax_key).hash(h)); let key = hash(&|h| (selected, self.pick.hover_generation).hash(h)); if self.pick.overlay.0 == Some(key) { return (key, self.pick.overlay.1.clone()); @@ -1439,6 +1519,10 @@ impl ModelApp { let m = self.with_picker(|p| p.items_overlay(&items, theme::rgb(t.pick_hover))); s.points.extend(m.points); } + // The Relax Nodes sphere and its preview (#398). + if relax_key.is_some() { + self.relax_overlay(&mut s, &t); + } self.pick.selected_overlay = (Some(selected), Arc::new(s)); } let mut s = (*self.pick.selected_overlay.1).clone(); @@ -1454,6 +1538,36 @@ impl ModelApp { (key, self.pick.overlay.1.clone()) } + /// The open Relax Nodes form drawn over the scene: its sphere of influence, translucent, and + /// the edges of the elements it changes at their relaxed positions. + fn relax_overlay(&self, s: &mut Scene, t: &theme::Tokens) { + let sphere = self + .dialog + .as_ref() + .and_then(|d| serde_json::from_str::(&d.text).ok()) + .and_then(|p| fieldcraft_engine::cmd::relax_sphere(&self.session, &p)); + if let Some((c, r)) = sphere { + scene::sphere(s, Vec3::from(c), r, theme::rgb(t.pick_hover)); + } + let Some(p) = self.relax_preview() else { return }; + let model = &self.session.model; + let at = |n: u32| p.relaxed.positions.get(&n).copied().or_else(|| model.node(n).map(|x| x.xyz)).map(Vec3::from); + let color = theme::rgb(t.pick_selected); + for e in model.elements().values().filter(|e| e.nodes.iter().any(|n| p.relaxed.positions.contains_key(n))) { + for edge in e.kind.edges() { + let ends = (edge.first().and_then(|&i| e.nodes.get(i)), edge.get(1).and_then(|&i| e.nodes.get(i))); + if let (Some(&a), Some(&b)) = ends + && let (Some(a), Some(b)) = (at(a), at(b)) + { + s.lines.push(fieldcraft_render::SceneLine { a, b, color, width: 1.5 }); + } + } + } + for x in p.relaxed.positions.values() { + s.points.push(fieldcraft_render::ScenePoint { p: Vec3::from(*x), color, size: 4.0 }); + } + } + /// The seed on `edge`, if it has one. pub fn seed(&self, edge: &fieldcraft_model::GeometryRef) -> Option<&fieldcraft_model::Seed> { self.session.model.seeds().iter().find(|s| s.edge == *edge) @@ -2274,6 +2388,35 @@ mod tests { assert!(app.pick.selection.is_empty()); } + /// Relax Nodes: a click places the sphere, Alt + scroll sizes it, the form previews the moved + /// nodes and the region's quality, and the sphere is drawn translucent (#398). + #[test] + fn relax_nodes_places_its_sphere_and_previews() { + let mut app = ModelApp::new(Session::default()); + let c = [[0., 0., 0.], [4., 0., 0.], [4., 4., 0.], [0., 4., 0.], [0., 0., 4.], [4., 0., 4.], [4., 4., 4.], [0., 4., 4.]]; + app.run("mesh.block", json!({"body": "A", "kind": "hex8", "corners": c, "divisions": [4, 4, 4]})).unwrap(); + // A spike inside. + let n = *app.session.model.nodes().iter().find(|(_, p)| p.xyz == [2.0, 2.0, 2.0]).unwrap().0; + app.session.model.set_node(n, [2.7, 2.6, 2.5]).unwrap(); + app.model_changed(); + assert!(!app.relaxing() && app.relax_preview().is_none()); + app.start("mesh.node.relax").unwrap(); + assert!(app.relaxing()); + assert!(app.relax_place(Vec3 { x: 2.0, y: 2.0, z: 2.0 })); + assert!(app.relax_scale_radius(10.0)); + let note = app.relax_note().unwrap(); + assert!(note.starts_with("Preview:"), "{note}"); + let p = app.relax_preview().unwrap(); + assert!(p.relaxed.positions.contains_key(&n) && p.relaxed.after.min_quality > p.relaxed.before.min_quality); + let (_, overlay) = app.overlay(); + assert!(!overlay.fills.is_empty() && !overlay.lines.is_empty(), "the sphere and the preview are drawn"); + // OK runs it as the form says. + let text = app.dialog.as_ref().unwrap().text.clone(); + app.dialog = None; + app.run("mesh.node.relax", serde_json::from_str(&text).unwrap()).unwrap(); + assert_ne!(app.session.model.node(n).unwrap().xyz, [2.7, 2.6, 2.5]); + } + /// Mesh Quality reads as a sentence in the status bar, not as its JSON report. #[test] fn mesh_quality_reports_in_a_line() { diff --git a/crates/ui-egui/src/reply.rs b/crates/ui-egui/src/reply.rs index 88ebf8e..070e622 100644 --- a/crates/ui-egui/src/reply.rs +++ b/crates/ui-egui/src/reply.rs @@ -214,6 +214,8 @@ enum Act { Pick(fieldcraft_model::GeometryRef), /// Select elements in the viewport and frame them (Mesh Quality's worst and invalid ones). Elements(Vec), + /// Open Relax Nodes round an element (Mesh Quality's worst). + Relax(u32), } /// The reply's section of the context panel. @@ -256,6 +258,10 @@ pub fn show(app: &mut ModelApp, ui: &mut egui::Ui, t: &Tokens) { } Act::Pick(g) => app.select([crate::selection::Item::Geometry(g)], crate::selection::Mode::Replace), Act::Elements(labels) => app.show_elements(&labels), + Act::Relax(element) => { + app.show_elements(&[element]); + app.open_form("mesh.node.relax", serde_json::json!({"around_element": element})); + } } } } @@ -413,6 +419,16 @@ fn quality(ui: &mut egui::Ui, t: &Tokens, v: &Value, acts: &mut Vec) { } } } + if let Some(worst) = v.pointer("/metrics/scaled_jacobian/worst").and_then(Value::as_u64).and_then(|w| u32::try_from(w).ok()) { + ui.add_space(6.0); + if ui + .button("Relax around worst element") + .on_hover_text("Open Relax Nodes with its sphere round the element of smallest scaled Jacobian; preview, then OK") + .clicked() + { + acts.push(Act::Relax(worst)); + } + } if let Some(skipped) = v.get("skipped").and_then(Value::as_object).filter(|s| !s.is_empty()) { ui.add_space(6.0); let text: Vec = skipped.iter().map(|(k, n)| format!("{k}: {}", n.as_u64().map_or_else(|| n.to_string(), count))).collect(); diff --git a/crates/ui-egui/src/scene.rs b/crates/ui-egui/src/scene.rs index b2496f0..0db5097 100644 --- a/crates/ui-egui/src/scene.rs +++ b/crates/ui-egui/src/scene.rs @@ -8,7 +8,7 @@ //! their own over the model, a constant size on screen, muted over a results plot. use fieldcraft_render::glyphs::{Definition, Drawn, GlyphShow, GlyphView, Glyphs}; -use fieldcraft_render::{Aabb3, Camera, DisplayMesh, Scene, SceneLabel, SceneLine, ScenePoint, SceneTri, Vec3}; +use fieldcraft_render::{Aabb3, Camera, DisplayMesh, Rgb, Scene, SceneFill, SceneLabel, SceneLine, ScenePoint, SceneTri, Vec3}; use serde::{Deserialize, Serialize}; use crate::theme::{Tokens, rgb}; @@ -99,6 +99,40 @@ pub fn datum_preview(s: &mut Scene, datum: fieldcraft_model::Datum, bounds: &Aab fieldcraft_render::glyphs::datums(&model, bounds, &|_| true, Some(NAME), &t.glyph_colors(), s) } +/// A sphere drawn translucent (a tool's region of influence, such as Relax Nodes'), with its +/// three great circles as lines. +pub fn sphere(s: &mut Scene, c: Vec3, r: f64, color: Rgb) { + if !(r > 0.0 && r.is_finite()) { + return; + } + const RINGS: usize = 12; + const SEGMENTS: usize = 24; + let at = |i: usize, j: usize| { + let theta = std::f64::consts::PI * i as f64 / RINGS as f64; + let phi = std::f64::consts::TAU * j as f64 / SEGMENTS as f64; + Vec3 { x: theta.sin() * phi.cos(), y: theta.sin() * phi.sin(), z: theta.cos() } + }; + for i in 0..RINGS { + for j in 0..SEGMENTS { + let (a, b, d, e) = (at(i, j), at(i + 1, j), at(i + 1, j + 1), at(i, j + 1)); + for tri in [[a, b, d], [a, d, e]] { + let n = (tri[0] + tri[1] + tri[2]).normalized().unwrap_or(tri[0]); + s.fills.push(SceneFill { tri: SceneTri { p: tri.map(|v| c + v * r), normal: n, color }, alpha: 0.18 }); + } + } + } + for k in 0..SEGMENTS { + let (p, q) = (std::f64::consts::TAU * k as f64 / SEGMENTS as f64, std::f64::consts::TAU * (k + 1) as f64 / SEGMENTS as f64); + for (u, v) in [ + (Vec3 { x: p.cos(), y: p.sin(), z: 0.0 }, Vec3 { x: q.cos(), y: q.sin(), z: 0.0 }), + (Vec3 { x: p.cos(), y: 0.0, z: p.sin() }, Vec3 { x: q.cos(), y: 0.0, z: q.sin() }), + (Vec3 { x: 0.0, y: p.cos(), z: p.sin() }, Vec3 { x: 0.0, y: q.cos(), z: q.sin() }), + ] { + s.lines.push(SceneLine { a: c + u * r, b: c + v * r, color, width: 1.0 }); + } + } +} + /// The scene for `display` drawn in `style`. pub fn build(display: &DisplayMesh, style: Style, show: Show, t: &Tokens) -> Scene { let mut s = Scene { background: Some((rgb(t.viewport_top), rgb(t.viewport_bottom))), ..Default::default() }; @@ -253,7 +287,7 @@ mod tests { .collect(); // An element edge across the seeded edge is not along it. d.edges.push(fieldcraft_render::DisplayEdge { - p: [Vec3::ZERO, Vec3::new(0.0, 10.0, 0.0)], + p: [Vec3::ZERO, Vec3 { x: 0.0, y: 10.0, z: 0.0 }], nodes: [0, 99], corners: [0, 99], feature: true, @@ -290,7 +324,7 @@ mod tests { let t = Tokens::dark(); let geometry = |s: &Scene| s.lines.iter().filter(|l| l.color == rgb(t.geometry_edge)).map(|l| [l.a, l.b]).collect::>(); let s = build(&d, Style::ShadedEdges, Show::DEFAULT, &t); - assert_eq!(geometry(&s), vec![[Vec3::new(0.0, 0.0, 0.0), Vec3::new(1.0, 0.0, 0.0)]]); + assert_eq!(geometry(&s), vec![[Vec3 { x: 0.0, y: 0.0, z: 0.0 }, Vec3 { x: 1.0, y: 0.0, z: 0.0 }]]); // Without the mesh the design's edges are drawn. let s = build(&d, Style::ShadedEdges, Show { mesh: false, ..Show::DEFAULT }, &t); assert_eq!(geometry(&s).len(), 3); diff --git a/crates/ui-egui/src/viewport.rs b/crates/ui-egui/src/viewport.rs index 1967151..fa36e77 100644 --- a/crates/ui-egui/src/viewport.rs +++ b/crates/ui-egui/src/viewport.rs @@ -47,8 +47,12 @@ pub fn show(app: &mut ModelApp, ui: &mut egui::Ui) { if response.clicked() && let Some(at) = response.interact_pointer_pos() { - let _ = app.select_at(at, mode(ui.input(|i| i.modifiers))); - probe_at(app, at, rect); + // With Relax Nodes open a click on the mesh places its sphere (#398). + let placed = app.relaxing() && app.surface_point(at).is_some_and(|p| app.relax_place(p)); + if !placed { + let _ = app.select_at(at, mode(ui.input(|i| i.modifiers))); + probe_at(app, at, rect); + } } if !app.capturing { hover(app, &response); @@ -383,7 +387,10 @@ fn navigate(app: &mut ModelApp, ui: &egui::Ui, response: &egui::Response, rect: let options = ui.ctx().options(|o| o.input_options); let (scroll, pinch, at) = ui.input(|i| (wheel_points(i, &options), i.zoom_delta(), i.pointer.hover_pos())); let factor = (-f64::from(scroll) * ZOOM_PER_POINT * zoom_sign(app)).exp() / f64::from(pinch).max(1e-3); - if (factor - 1.0).abs() > 1e-9 { + // Alt + scroll sizes the Relax Nodes sphere instead of zooming. + let sized = + app.relaxing() && ui.input(|i| i.modifiers.alt) && scroll != 0.0 && app.relax_scale_radius((f64::from(scroll) * ZOOM_PER_POINT).exp()); + if !sized && (factor - 1.0).abs() > 1e-9 { let anchor = at.and_then(|p| { app.surface_point(p).or_else(|| app.cam.target_plane_point(f64::from(p.x - rect.left()), f64::from(p.y - rect.top()), w, h)) }); diff --git a/docs/features.md b/docs/features.md index a62a9e8..55a908a 100644 --- a/docs/features.md +++ b/docs/features.md @@ -6,7 +6,7 @@ Every entry is a toolbar button. **Ready** entries run their command; **coming s shown greyed with the roadmap milestone that brings them (ROADMAP.md). Scripts, the command line, the control channel and AI agents over MCP (docs/mcp.md) run the same commands. -70 ready, 39 coming soon. +71 ready, 39 coming soon. ## Geometry @@ -50,6 +50,7 @@ the control channel and AI agents over MCP (docs/mcp.md) run the same commands. | Edit | Element — Add an element on existing nodes | ready | `mesh.element.add`, `mesh.element.remove`, `mesh.clear` | | Edit | Mirror Mesh — Mirror a body's mesh across a plane, joined on the plane (mesh half a symmetric part) | ready | `mesh.mirror` | | Edit | Merge Nodes — Merge coincident nodes, keeping sets, loads and boundary conditions on the survivors; check for duplicates and cracks, renumber | ready | `mesh.node.merge`, `mesh.check.duplicates`, `mesh.renumber` | +| Improve | Relax Nodes — Relax the nodes inside a sphere to fix spikes and slivers; surface nodes slide on the geometry; kept as a field to apply again after a remesh | ready | `mesh.node.relax`, `mesh.node.relax.apply` | | Edit | Remesh Region — Remesh a region in place; split and collapse elements | coming soon (M6) | | | Edit | Orphan Mesh — Import a mesh without geometry | coming soon (M6) | | | Verify | Quality — Quality metrics and the certified validity check, failing elements listed | ready | `mesh.quality` | diff --git a/docs/mcp.md b/docs/mcp.md index dc969ff..45e97a2 100644 --- a/docs/mcp.md +++ b/docs/mcp.md @@ -167,6 +167,24 @@ edges of partition cells matching counts. 4. `mesh.seed.list` lists the seeds and the edges they reach (`stale` when the design has no such edge); `mesh.seed.clear {edges?}` removes them. Seeds are saved with the model. +## Relaxing a few bad elements + +When `mesh.quality` shows a handful of spikes in an otherwise good mesh: + +1. `execute` `mesh.node.relax {"around_element": }` (or `{"center": + [x, y, z], "radius": r}`, or `nodes` / `elements` / `set`; `iterations`, `strength`, `falloff`). + Nodes there relax (smart Laplacian, then optimization of the worst element), each move kept + only when the worst element around the node does not get worse and nothing inverts. Nodes on a + B-rep face slide in its UV parameters and stay on it, nodes on a B-rep edge slide along it, + vertex nodes stay; nodes on a partition cut face move freely (both cells share them). The + reply gives the nodes moved, the largest displacement, and the region's smallest scaled + Jacobian, largest edge ratio and worst elements before and after. One undo step. +2. The relaxation is kept as a field of the model (`field`, default `Relax1`, …; saved with it). + After a remesh, `mesh.node.relax.apply {"field": "Relax1"}` puts the relaxed shape on the new + mesh (each node takes the displacement interpolated in the old cell it lies in, on its face, + then one guarded smoothing sweep); with `"reapply": true` on the relaxation every remesh does + it by itself (the mesher's reply lists `relaxations_applied`). + ## Limits The server treats its input as hostile. A message line may be at most 16 MiB (the control diff --git a/docs/meshing/literature.md b/docs/meshing/literature.md index 634ee48..531e69e 100644 --- a/docs/meshing/literature.md +++ b/docs/meshing/literature.md @@ -952,6 +952,25 @@ input points), Steiner points kept half the local spacing from the boundary. is exact: a planar bilinear quad's Jacobian is affine in each reference coordinate, so four exact `orient2d` corner tests decide it. +### B7a. Node relaxation of a finished mesh (implemented in `fieldcraft-mesh::relax`, `mesh.node.relax`) + +- **Method.** Smart Laplacian (S. Canann, J. Tristano, M. Staten, *An approach to combined Laplacian and + optimization-based smoothing for triangular, quadrilateral, and quad-dominant meshes*, IMR 1998; L. Freitag 1997 + above): a node's weighted step toward the centroid of its edge neighbours is kept only when the worst element + around it does not get worse and no element inverts, halved otherwise. Then one step up the gradient of the worst + element's quality (the active max–min element of Freitag and Ollivier-Gooch, IJNME 40, 1997; P. Knupp, IJNME 48, + 2000), kept only when it raises the worst. Quality is the scaled Jacobian of `crate::quality`; validity its + certified check. +- **Surface nodes in UV.** A node on a B-rep face moves in the face's own parameters (the displacement resolved on + the surface's coordinate directions, then the surface evaluated there), so it stays exactly on an analytic face; + free-form faces fall back to projection onto the face's tessellation (SolidCraft's kernel exposes the analytic + surfaces only). Nodes on a B-rep edge slide along it; vertex nodes are held; partition cut faces are not + geometry, so their nodes move freely. +- **Relaxation fields.** What a relaxation did is kept as the relaxed region's cells before and after; a later + mesh takes the displacement at each node's parametric coordinates in the old cell holding it (barycentric or + trilinear shape functions), slides it onto its face or edge, undoes any move that would invert an element, and + smooths once (`mesh.node.relax.apply`, or `reapply` after every remesh). + --- ### B8. Mapped and swept hexahedral meshing (implemented in `fieldcraft-mesh::hex`) diff --git a/skills/fieldcraft-meshing/SKILL.md b/skills/fieldcraft-meshing/SKILL.md index e7e869f..edb5b05 100644 --- a/skills/fieldcraft-meshing/SKILL.md +++ b/skills/fieldcraft-meshing/SKILL.md @@ -62,6 +62,11 @@ Before the first mesh command: 1. `mesh.quality`: all valid; scaled Jacobian min ≥ 0.3 (≥ 0.5 in the region of interest); edge ratio max ≤ 10; no element below these is acceptable without a stated reason. + **Isolated spikes** (a few bad elements in an otherwise good mesh) are fixed last, with + `mesh.node.relax {"around_element": }` (or `center` + `radius`): nodes there relax, surface + nodes sliding on their faces, nothing inverting. It is a final touch, never a substitute for a + better block plan: many bad elements mean the cuts are wrong (§2). Pass `"reapply": true` (or + run `mesh.node.relax.apply {field}`) so the relaxation follows a later remesh. 2. Look at it: screenshots of each body and a section. Fans, pinches, wedges, slivers → back to §2. 3. Volume: `mesh.quality` volume equals the geometry volume (a mismatch means orphan or missing cells, #385).