diff --git a/ATTRIBUTION.md b/ATTRIBUTION.md index 2154986..f007f4b 100644 --- a/ATTRIBUTION.md +++ b/ATTRIBUTION.md @@ -92,6 +92,7 @@ is committed and listed here with its exports. | 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 | +| assets/icons/results.export.vtu.svg, assets/icons/results.export.vtu.vectorcraft | FieldCraft contributors | Original toolbar icon drawn with VectorCraft from `assets/icons/src/results.export.vtu.json` (`cargo xtask icons`) | MIT OR Apache-2.0 | | assets/icons/mesh.export.stl.svg, assets/icons/mesh.export.stl.vectorcraft | FieldCraft contributors | Original toolbar icon drawn with VectorCraft from `assets/icons/src/mesh.export.stl.json` (`cargo xtask icons`) | MIT OR Apache-2.0 | | assets/icons/mesh.region.remesh.svg, assets/icons/mesh.region.remesh.vectorcraft | FieldCraft contributors | Original toolbar icon drawn with VectorCraft from `assets/icons/src/mesh.region.remesh.json` (`cargo xtask icons`) | MIT OR Apache-2.0 | | assets/icons/mesh.size.adapt.svg, assets/icons/mesh.size.adapt.vectorcraft | FieldCraft contributors | Original toolbar icon drawn with VectorCraft from `assets/icons/src/mesh.size.adapt.json` (`cargo xtask icons`) | MIT OR Apache-2.0 | diff --git a/assets/icons/results.export.vtu.svg b/assets/icons/results.export.vtu.svg new file mode 100644 index 0000000..55f71b8 --- /dev/null +++ b/assets/icons/results.export.vtu.svg @@ -0,0 +1,10 @@ + + + results.export.vtu + + + + + + + diff --git a/assets/icons/results.export.vtu.vectorcraft b/assets/icons/results.export.vtu.vectorcraft new file mode 100644 index 0000000..140555c Binary files /dev/null and b/assets/icons/results.export.vtu.vectorcraft differ diff --git a/assets/icons/src/results.export.vtu.json b/assets/icons/src/results.export.vtu.json new file mode 100644 index 0000000..8752608 --- /dev/null +++ b/assets/icons/src/results.export.vtu.json @@ -0,0 +1 @@ +{"shapes": [{"d": "M4 3 L12 3 L16 7 L16 21 L4 21 Z", "fill": "soft", "stroke": "line"}, {"d": "M12 3 L12 7 L16 7", "stroke": "line"}, {"d": "M7 17 L7 14 M10 17 L10 11 M13 17 L13 13", "stroke": "accent", "width": 2}, {"d": "M15 15 L21 15 M18.19 16.53 L21 15 L18.19 13.47", "stroke": "accent", "width": 1.5}]} diff --git a/crates/engine/src/catalogue.rs b/crates/engine/src/catalogue.rs index 632015d..5c2d0fa 100644 --- a/crates/engine/src/catalogue.rs +++ b/crates/engine/src/catalogue.rs @@ -618,6 +618,7 @@ pub static CATALOGUE: &[Capability] = &[ soon("results.cut", "Cut", Results, "Plots", "M9", "Section cuts through the model"), soon("results.display_group", "Display Group", Results, "Plots", "M9", "Show and hide parts of the model"), soon("results.streamlines", "Streamlines", Results, "Plots", "M12", "Streamlines and flow fields"), + ready("results.export.vtu", "Export VTU", Results, "Files", "Write a solved step as VTK unstructured grids for ParaView"), soon("results.read", "Read Results", Results, "Files", "M9", "Read results in public formats"), soon("results.convergence", "Convergence", Results, "Files", "M9", "Mesh convergence studies"), // ---- Tools ---- diff --git a/crates/engine/src/cmd/forms.rs b/crates/engine/src/cmd/forms.rs index 8db144e..c035da9 100644 --- a/crates/engine/src/cmd/forms.rs +++ b/crates/engine/src/cmd/forms.rs @@ -1838,6 +1838,17 @@ pub const EXPORT_CSV: &[Field] = &[ int("mode", "Mode", 1, 1_000_000).help("Probe: the mode to read, for frequency and buckling steps."), int("frame", "Frame", 0, 1_000_000_000).help("Probe: the frame to read, for explicit and transient steps."), ]; +pub const EXPORT_VTU: &[Field] = &[ + path("path", "File", true).help("A .vtu file for one state, or a .pvd collection for every frame or mode.").req(), + SOLVED, + int("mode", "Mode", 1, 1_000_000).help("Frequency steps: the mode to write; every mode in a .pvd when empty.").hint("all (.pvd) or 1"), + choice("frames", "Frames", &[("all", "Every frame"), ("last", "The last frame")]) + .help("A .pvd series: every frame or only the last.") + .def("\"all\""), + choice("encoding", "Encoding", &[("ascii", "Text"), ("binary", "Binary (base64)")]) + .help("Text is readable; binary is smaller and exact.") + .def("\"ascii\""), +]; pub const REACTIONS: &[Field] = &[ SOLVED, references("sets", "Sets", C::Set) diff --git a/crates/engine/src/cmd/mod.rs b/crates/engine/src/cmd/mod.rs index e31668d..e4f91a3 100644 --- a/crates/engine/src/cmd/mod.rs +++ b/crates/engine/src/cmd/mod.rs @@ -182,6 +182,7 @@ pub static ALL: &[CommandSpec] = &[ results::FRAME, results::REACTIONS, results::EXPORT_CSV, + results::EXPORT_VTU, results::CLEAR, diagnostics::CRASH_REPORTS, diagnostics::CRASH_REPORT, diff --git a/crates/engine/src/cmd/results.rs b/crates/engine/src/cmd/results.rs index 1180924..0514717 100644 --- a/crates/engine/src/cmd/results.rs +++ b/crates/engine/src/cmd/results.rs @@ -70,6 +70,16 @@ pub const REACTIONS: CommandSpec = CommandSpec { run: reactions, }; +pub const EXPORT_VTU: CommandSpec = CommandSpec { + id: "results.export.vtu", + label: "Export VTU", + params: r#"{"step": "Load", "path": "out.vtu" | "out.pvd", "mode": 1, "frames": "all|last", "encoding": "ascii|binary"} — writes a solved step as VTK XML unstructured grids for ParaView: the mesh (points in node label order, cells in element label order, arrays `node_label` and `element_label`), the step's nodal fields as point data and per-element fields as cell data. Static: displacement, rotation (with shells or beams), reaction, stress (6 components in VTK order XX, YY, ZZ, XY, YZ, XZ), von_mises, temperature (coupled steps); nonlinear static: displacement per frame, the final fields on the last; frequency: mode_shape, with mode, frequency and eigenvalue as field data; explicit: displacement and velocity per frame, cell stress, plastic_strain, damage, eroded; heat: temperature per frame, heat_flux on the last. Point array `has_result` is 0 at nodes the step has no values for (written as 0, as are values that are not finite). A `.vtu` path writes one state (the last frame, or `mode`, default 1); a `.pvd` path writes every frame or mode (`frames: "last"` only the last) as `out_0000.vtu`, … beside it and the collection listing each with its time (load factor, mode number). `encoding` ascii (default) or binary (inline base64). Fluid steps are not supported yet. Replies {files, bytes, points, cells}"#, + form: Some(super::forms::EXPORT_VTU), + kind: CommandKind::Session, + requires: &[Requirement::Result], + run: export_vtu, +}; + pub const EXPORT_CSV: CommandSpec = CommandSpec { id: "results.export.csv", label: "Export Results (CSV)", @@ -379,6 +389,251 @@ fn reactions(s: &mut Session, p: &Value) -> Result { Ok(json!({"step": step, "stale": !current, "about": about, "sets": per_set, "total": total, "applied": r.applied_total})) } +/// Most files one `results.export.vtu` writes. +pub const MAX_VTU_FILES: usize = 10_000; +/// Most bytes one `results.export.vtu` writes, over all its files. +pub const MAX_VTU_BYTES: u64 = 8 << 30; + +/// One state of a result written as a `.vtu`: the end of a linear step, a mode (from 1) or a +/// frame (from 0). +#[derive(Clone, Copy, Debug, PartialEq)] +enum State { + Final, + Mode(usize), + Frame(usize), +} + +fn export_vtu(s: &mut Session, p: &Value) -> Result { + use fieldcraft_io::vtk::{self, Encoding}; + let step = step_for(s, p, Want::Any)?; + let mode: Option = opt(p, "mode")?; + let last_only = match opt::(p, "frames")?.as_deref() { + None | Some("all") => false, + Some("last") => true, + Some(_) => return Err(EngineError::invalid("frames", "must be \"all\" or \"last\"")), + }; + let encoding = match opt::(p, "encoding")?.as_deref() { + None | Some("ascii") => Encoding::Ascii, + Some("binary") => Encoding::Base64, + Some(_) => return Err(EngineError::invalid("encoding", "must be \"ascii\" or \"binary\"")), + }; + let path = super::file_param(s, p, "path", fieldcraft_io::roots::Access::Write)?; + let target = std::path::Path::new(&path); + let series = match target.extension().map(|e| e.to_string_lossy().to_ascii_lowercase()).as_deref() { + Some("vtu") => false, + Some("pvd") => true, + _ => return Err(EngineError::Params("`path` must end in .vtu (one file) or .pvd (a series of .vtu files)".into())), + }; + let r = match s.latest_result(&step) { + Some((_, false)) => return Err(EngineError::Params(format!("the model changed since step `{step}` was solved; solve it again"))), + Some((r, true)) => std::sync::Arc::clone(r), + None => return Err(EngineError::Params(format!("step `{step}` has no result; solve it first"))), + }; + let frames = |n: usize, what: &str| -> Result> { + let last = n.checked_sub(1).ok_or_else(|| EngineError::Params(format!("step `{step}` has no {what}")))?; + Ok(if series && !last_only { (0..n).map(State::Frame).collect() } else { vec![State::Frame(last)] }) + }; + if mode.is_some() && !matches!(r.kind(), Kind::Frequency(_)) { + return Err(EngineError::invalid("mode", "is for frequency and buckling steps")); + } + let states: Vec = match r.kind() { + Kind::Fluid(_) => return Err(EngineError::Params("fluid results are per cell; not supported yet".into())), + Kind::Static(st) => match &st.nonlinear { + Some(h) if !h.frames.is_empty() => frames(h.frames.len(), "frames")?, + _ => vec![State::Final], + }, + Kind::Frequency(f) => match mode { + Some(m) if m == 0 || m > f.modes.len() => { + return Err(EngineError::invalid("mode", format!("step `{step}` has modes 1 to {}", f.modes.len()))); + } + Some(m) => vec![State::Mode(m)], + None if series && !last_only => (1..=f.modes.len()).map(State::Mode).collect(), + None if f.modes.is_empty() => return Err(EngineError::Params(format!("step `{step}` has no modes"))), + None => vec![State::Mode(1)], + }, + Kind::Explicit(e) => frames(e.frames.len(), "field frames")?, + Kind::Heat(h) => frames(h.frames.len(), "temperatures")?, + }; + if states.len() > MAX_VTU_FILES { + return Err(EngineError::Params(format!("step `{step}` has {} frames; at most {MAX_VTU_FILES} files are written", states.len()))); + } + + let model = &s.model; + let (nodes, elements) = vtk::points_and_cells(model); + let build = |state: State| -> Result<(f64, String)> { + let (time, fields, data) = vtu_fields(&r, state, &nodes, &elements); + let data: Vec<(&str, f64)> = data.iter().map(|(n, v)| (*n, *v)).collect(); + Ok((time, vtk::write_vtu(model, &fields, &data, encoding)?)) + }; + // The last state has the most fields: its size times the count bounds the output before + // anything is written. + let Some(&last) = states.last() else { return Err(EngineError::Params(format!("step `{step}` has nothing to write"))) }; + let sample = build(last)?; + let estimate = (sample.1.len() as u64).saturating_mul(states.len() as u64); + if estimate > MAX_VTU_BYTES { + return Err(EngineError::Params(format!( + "the export would write about {estimate} bytes; at most {MAX_VTU_BYTES} (8 GiB) are written: choose `frames: \"last\"` or a single .vtu" + ))); + } + + let mut written: Vec = Vec::new(); + let mut bytes = 0usize; + let result = (|| -> Result> { + if !series { + fieldcraft_io::fs::write_atomic(target, sample.1.as_bytes(), false)?; + written.push(target.to_path_buf()); + bytes += sample.1.len(); + return Ok(vec![path.clone()]); + } + let stem = target.file_stem().map(|n| n.to_string_lossy().into_owned()).unwrap_or_default(); + let dir = target.parent().unwrap_or(std::path::Path::new("")); + let mut entries: Vec<(f64, String)> = Vec::with_capacity(states.len()); + let mut files = Vec::with_capacity(states.len() + 1); + for (k, &state) in states.iter().enumerate() { + let name = format!("{stem}_{k:04}.vtu"); + let file = dir.join(&name); + let built = if state == last { None } else { Some(build(state)?) }; + let (time, text) = built.as_ref().map_or((sample.0, sample.1.as_str()), |(t, x)| (*t, x.as_str())); + fieldcraft_io::fs::write_atomic(&file, text.as_bytes(), false)?; + written.push(file.clone()); + bytes += text.len(); + entries.push((time, name)); + files.push(file.to_string_lossy().into_owned()); + } + let list: Vec<(f64, &str)> = entries.iter().map(|(t, n)| (*t, n.as_str())).collect(); + let collection = vtk::pvd(&list); + fieldcraft_io::fs::write_atomic(target, collection.as_bytes(), false)?; + bytes += collection.len(); + files.push(path.clone()); + Ok(files) + })(); + match result { + Ok(files) => Ok(json!({"files": files, "bytes": bytes, "points": nodes.len(), "cells": elements.len()})), + Err(e) => { + // Leave nothing half written behind. + for f in &written { + let _ = std::fs::remove_file(f); + } + Err(e) + } + } +} + +/// A stress (xx, yy, zz, xy, xz, yz) in VTK's symmetric tensor order (XX, YY, ZZ, XY, YZ, XZ). +fn vtk_tensor(t: &[f64; 6]) -> [f64; 6] { + [t[0], t[1], t[2], t[3], t[5], t[4]] +} + +/// Values given per label (`labels` sorted, as a result's nodes and elements are) spread over +/// `order` (the points' or cells' labels): `N` components each, 0 where a label has none. +fn spread(order: &[u32], labels: &[u32], values: &[[f64; N]]) -> Vec { + let mut out = vec![0.0; order.len() * N]; + for (label, v) in labels.iter().zip(values) { + if let Ok(i) = order.binary_search(label) + && let Some(slot) = out.get_mut(i * N..(i + 1) * N) + { + slot.copy_from_slice(v); + } + } + out +} + +fn scalars(values: &[f64]) -> Vec<[f64; 1]> { + values.iter().map(|v| [*v]).collect() +} + +/// The fields of one state of a result, its time (load factor, mode number) and its field data. +fn vtu_fields( + r: &crate::StepResult, + state: State, + nodes: &[u32], + elements: &[u32], +) -> (f64, Vec, Vec<(&'static str, f64)>) { + use fieldcraft_io::vtk::{Field, Location}; + let mut fields = Vec::new(); + let mut data = Vec::new(); + let mut add = |name: &str, location: Location, components: usize, values: Vec, integer: bool| { + fields.push(Field { name: name.into(), location, components, values, integer }); + }; + let on = |labels: &[u32]| -> Vec { nodes.iter().map(|n| if labels.binary_search(n).is_ok() { 1.0 } else { 0.0 }).collect() }; + add("has_result", Location::Point, 1, on(r.nodes()), true); + let point = Location::Point; + let cell = Location::Cell; + let time = match (r.kind(), state) { + (Kind::Static(st), _) => { + let frame = match (state, &st.nonlinear) { + (State::Frame(k), Some(h)) => h.frames.get(k).map(|f| (f.load_factor, f.displacement.as_slice(), k + 1 == h.frames.len())), + _ => None, + }; + let (time, displacement, last) = frame.unwrap_or((0.0, st.displacement.as_slice(), true)); + add("displacement", point, 3, spread(nodes, &st.nodes, displacement), false); + if last { + if !st.shells.is_empty() || !st.beams.is_empty() { + add("rotation", point, 3, spread(nodes, &st.nodes, &st.rotation), false); + } + add("reaction", point, 3, spread(nodes, &st.nodes, &st.reaction), false); + if !st.stress.is_empty() { + let stress: Vec<[f64; 6]> = st.stress.iter().map(vtk_tensor).collect(); + add("stress", point, 6, spread(nodes, &st.nodes, &stress), false); + add("von_mises", point, 1, spread(nodes, &st.nodes, &scalars(&st.von_mises)), false); + } + if let Some(t) = &st.temperature { + add("temperature", point, 1, spread(nodes, &st.nodes, &scalars(t)), false); + } + } + if st.nonlinear.is_some() { + data.push(("load_factor", time)); + } + time + } + (Kind::Frequency(f), State::Mode(m)) => { + if let Some(mode) = m.checked_sub(1).and_then(|k| f.modes.get(k)) { + add("mode_shape", point, 3, spread(nodes, &f.nodes, &mode.shape), false); + data.extend([("mode", m as f64), ("frequency", mode.frequency), ("eigenvalue", mode.eigenvalue)]); + } + m as f64 + } + (Kind::Explicit(e), State::Frame(k)) => { + let Some(f) = e.frames.get(k) else { return (0.0, fields, data) }; + if !f.displacement.is_empty() { + add("displacement", point, 3, spread(nodes, &e.nodes, &f.displacement), false); + } + if !f.velocity.is_empty() { + add("velocity", point, 3, spread(nodes, &e.nodes, &f.velocity), false); + } + if !f.stress.is_empty() { + let stress: Vec<[f64; 6]> = f.stress.iter().map(vtk_tensor).collect(); + add("stress", cell, 6, spread(elements, &e.elements, &stress), false); + } + if !f.plastic_strain.is_empty() { + add("plastic_strain", cell, 1, spread(elements, &e.elements, &scalars(&f.plastic_strain)), false); + } + if !f.damage.is_empty() { + add("damage", cell, 1, spread(elements, &e.elements, &scalars(&f.damage)), false); + } + if !f.failed.is_empty() { + let eroded: Vec<[f64; 1]> = f.failed.iter().map(|&b| [if b { 1.0 } else { 0.0 }]).collect(); + add("eroded", cell, 1, spread(elements, &e.elements, &eroded), true); + } + data.push(("time", f.time)); + f.time + } + (Kind::Heat(h), State::Frame(k)) => { + let time = h.times.get(k).copied().unwrap_or(0.0); + if let Some(t) = h.frames.get(k) { + add("temperature", point, 1, spread(nodes, &h.nodes, &scalars(t)), false); + } + if k + 1 == h.frames.len() && !h.flux.is_empty() { + add("heat_flux", point, 3, spread(nodes, &h.nodes, &h.flux), false); + } + data.push(("time", time)); + time + } + _ => 0.0, + }; + (time, fields, data) +} + fn export_csv(s: &mut Session, p: &Value) -> Result { let what: String = super::req(p, "what")?; let path = super::file_param(s, p, "path", fieldcraft_io::roots::Access::Write)?; @@ -549,6 +804,164 @@ mod tests { use super::*; use crate::cmd::solve::tests::cantilever; + /// The numbers of the first `DataArray` called `name` in a written `.vtu`. + fn vtu_array(xml: &str, name: &str) -> Vec { + let at = xml.find(&format!("Name=\"{name}\"")).unwrap_or_else(|| panic!("no array {name}")); + let start = at + xml[at..].find('>').unwrap() + 1; + let end = start + xml[start..].find("").unwrap(); + xml[start..end].split_whitespace().map(|t| t.parse().unwrap()).collect() + } + + fn scratch(name: &str) -> std::path::PathBuf { + let dir = std::env::temp_dir().join(format!("fieldcraft-vtu-{name}-{}", std::process::id())); + let _ = std::fs::remove_dir_all(&dir); + std::fs::create_dir_all(&dir).unwrap(); + dir + } + + #[test] + fn export_vtu_writes_a_static_step_and_its_modes() { + let dir = scratch("static"); + let out = dir.join("out.vtu").to_string_lossy().into_owned(); + let mut s = Session::default(); + s.execute("file.sample", &json!({"name": "cantilever"})).unwrap(); + let e = s.execute("results.export.vtu", &json!({"step": "Load", "path": out})).unwrap_err().to_string(); + assert!(e.contains("Needs a solved step"), "{e}"); + s.execute("solve.static", &json!({"step": "Load"})).unwrap(); + let reply = s.execute("results.export.vtu", &json!({"step": "Load", "path": out})).unwrap(); + let xml = std::fs::read_to_string(&out).unwrap(); + assert_eq!(reply["points"], s.model.nodes().len()); + assert_eq!(reply["cells"], s.model.elements().len()); + assert_eq!(reply["bytes"], xml.len()); + assert!(xml.contains(&format!("NumberOfPoints=\"{}\"", s.model.nodes().len()))); + + // The tip node's displacement is the probe's, at the node's point. + let labels = vtu_array(&xml, "node_label"); + let u = vtu_array(&xml, "displacement"); + let r = s.result("Load").unwrap().clone(); + let Kind::Static(st) = r.kind() else { panic!() }; + let (k, tip) = st.displacement.iter().enumerate().max_by(|a, b| a.1[2].abs().total_cmp(&b.1[2].abs())).unwrap(); + let tip_label = st.nodes[k]; + let probe = s.execute("solve.probe", &json!({"step": "Load", "nodes": [tip_label]})).unwrap(); + let i = labels.iter().position(|&l| l == f64::from(tip_label)).unwrap(); + assert_eq!(&u[3 * i..3 * i + 3], tip.as_slice()); + assert_eq!( + probe["nodes"][0]["displacement"].as_array().map(|a| a.iter().filter_map(Value::as_f64).collect::>()), + Some(tip.to_vec()), + "{probe}" + ); + // Stress components in VTK order: XY at 3, YZ at 4, XZ at 5. + let stress = vtu_array(&xml, "stress"); + let t = st.stress[k]; + assert_eq!(&stress[6 * i..6 * i + 6], [t[0], t[1], t[2], t[3], t[5], t[4]].as_slice()); + assert_eq!(vtu_array(&xml, "has_result").iter().filter(|&&h| h == 1.0).count(), st.nodes.len()); + // Deterministic. + s.execute("results.export.vtu", &json!({"step": "Load", "path": out})).unwrap(); + assert_eq!(std::fs::read_to_string(&out).unwrap(), xml); + + // Modes: one file per mode with its frequency, and a collection. + s.execute("solve.frequency", &json!({"step": "Modes"})).unwrap(); + let pvd = dir.join("modes.pvd").to_string_lossy().into_owned(); + let reply = s.execute("results.export.vtu", &json!({"step": "Modes", "path": pvd})).unwrap(); + let files = reply["files"].as_array().unwrap(); + assert_eq!(files.len(), 7, "six modes and the collection: {reply}"); + let mode2 = std::fs::read_to_string(dir.join("modes_0001.vtu")).unwrap(); + assert_eq!(vtu_array(&mode2, "mode"), [2.0]); + let collection = std::fs::read_to_string(&pvd).unwrap(); + assert!(collection.contains("timestep=\"2.0\" part=\"0\" file=\"modes_0001.vtu\""), "{collection}"); + let one = dir.join("mode3.vtu").to_string_lossy().into_owned(); + s.execute("results.export.vtu", &json!({"step": "Modes", "path": one, "mode": 3, "encoding": "binary"})).unwrap(); + assert!(std::fs::read_to_string(&one).unwrap().contains("format=\"binary\"")); + + // Refusals write nothing. + let none = dir.join("none.vtu").to_string_lossy().into_owned(); + for (p, why) in [ + (json!({"step": "Modes", "path": none, "mode": 0}), "modes 1 to 6"), + (json!({"step": "Modes", "path": none, "mode": 1_000_000_000}), "modes 1 to 6"), + (json!({"step": "Load", "path": none, "mode": 1}), "frequency"), + (json!({"step": "Load", "path": dir.join("x.vtk").to_string_lossy()}), ".vtu"), + (json!({"step": "Load", "path": none, "frames": 3}), "frames"), + (json!({"step": "Load", "path": none, "encoding": "hex"}), "encoding"), + (json!({"step": "Nope", "path": none}), "has no result"), + ] { + let e = s.execute("results.export.vtu", &p).unwrap_err().to_string(); + assert!(e.contains(why), "{p}: {e}"); + } + s.execute("mesh.node.add", &json!({"label": 9999, "xyz": [50, 50, 50]})).unwrap(); + let e = s.execute("results.export.vtu", &json!({"step": "Load", "path": none})).unwrap_err().to_string(); + assert!(e.contains("the model changed since step `Load` was solved"), "{e}"); + assert!(!dir.join("none.vtu").exists() && !dir.join("x.vtk").exists()); + let _ = std::fs::remove_dir_all(&dir); + } + + #[test] + fn export_vtu_writes_nonlinear_frames_and_heat() { + let dir = scratch("nl-heat"); + let mut s = Session::default(); + cantilever(&mut s); + s.run_script(&json!([ + {"command": "model.section.edit", "params": {"name": "bar", "set": "all", "material": "Mat", "kind": {"type": "solid"}, "formulation": {}}}, + {"command": "model.step.create", "params": {"name": "Large", "procedure": {"type": "nonlinear_static", "initial_increment": 0.25}, + "boundary_conditions": ["clamp"], "loads": ["push"]}}, + {"command": "model.material.edit", "params": {"name": "Mat", "elastic": {"type": "isotropic", "young": 1000, "poisson": 0}, "density": 1, "thermal": {"conductivity": 2}}}, + {"command": "model.boundary_condition.create", "params": {"name": "Heater", "type": "temperature", "target": {"set": "root"}, "value": 100}}, + {"command": "model.set.create", "params": {"name": "skin", "faces": (1..=6).map(|f| json!({"element": 10, "face": f})).collect::>()}}, + {"command": "model.load.create", "params": {"name": "Air", "type": "convection", "set": "skin", "film": 0.5, "sink": 20}}, + {"command": "model.step.create", "params": {"name": "Therm", "procedure": {"type": "heat_transfer"}, "boundary_conditions": ["Heater"], "loads": ["Air"]}}, + ])) + .unwrap(); + s.execute("solve.nonlinear", &json!({"step": "Large"})).unwrap(); + s.execute("solve.heat", &json!({"step": "Therm"})).unwrap(); + + let pvd = dir.join("large.pvd").to_string_lossy().into_owned(); + let reply = s.execute("results.export.vtu", &json!({"step": "Large", "path": pvd})).unwrap(); + let Kind::Static(st) = s.result("Large").unwrap().kind() else { panic!() }; + let h = st.nonlinear.as_ref().unwrap(); + assert_eq!(reply["files"].as_array().unwrap().len(), h.frames.len() + 1); + let first = std::fs::read_to_string(dir.join("large_0000.vtu")).unwrap(); + let last = std::fs::read_to_string(dir.join(format!("large_{:04}.vtu", h.frames.len() - 1))).unwrap(); + assert!(first.contains("Name=\"displacement\"") && !first.contains("Name=\"stress\""), "final stresses on the last frame only"); + assert!(last.contains("Name=\"stress\"") && last.contains("Name=\"von_mises\"")); + assert_eq!(vtu_array(&last, "load_factor"), [h.frames.last().unwrap().load_factor]); + + let heat = dir.join("heat.vtu").to_string_lossy().into_owned(); + s.execute("results.export.vtu", &json!({"step": "Therm", "path": heat})).unwrap(); + let xml = std::fs::read_to_string(&heat).unwrap(); + let t = vtu_array(&xml, "temperature"); + assert!(t.contains(&100.0), "the heater's temperature"); + assert_eq!(vtu_array(&xml, "heat_flux").len(), 3 * s.model.nodes().len()); + let _ = std::fs::remove_dir_all(&dir); + } + + #[test] + fn export_vtu_writes_explicit_frames_as_a_series() { + let dir = scratch("impact"); + let mut s = Session::default(); + s.execute("file.sample", &json!({"name": "impact"})).unwrap(); + s.execute("solve.explicit", &json!({})).unwrap(); + let pvd = dir.join("impact.pvd").to_string_lossy().into_owned(); + let reply = s.execute("results.export.vtu", &json!({"path": pvd})).unwrap(); + let Kind::Explicit(e) = s.result("Impact").unwrap().kind() else { panic!() }; + let frames = e.frames.len(); + assert!(frames > 2); + assert_eq!(reply["files"].as_array().unwrap().len(), frames + 1); + let collection = std::fs::read_to_string(&pvd).unwrap(); + let times: Vec = + collection.match_indices("timestep=\"").map(|(i, m)| collection[i + m.len()..].split('"').next().unwrap().parse().unwrap()).collect(); + assert_eq!(times.len(), frames); + assert!(times.windows(2).all(|w| w[0] < w[1]), "{times:?}"); + let last = std::fs::read_to_string(dir.join(format!("impact_{:04}.vtu", frames - 1))).unwrap(); + let cell = &last[last.find("").unwrap()..last.find("").unwrap()]; + assert!(cell.contains("Name=\"stress\" NumberOfComponents=\"6\""), "{cell}"); + // A .vtu takes the last frame; `frames: "last"` makes a one-entry series. + let one = dir.join("last.vtu").to_string_lossy().into_owned(); + s.execute("results.export.vtu", &json!({"path": one})).unwrap(); + assert_eq!(std::fs::read_to_string(&one).unwrap(), last); + let lastp = dir.join("last.pvd").to_string_lossy().into_owned(); + assert_eq!(s.execute("results.export.vtu", &json!({"path": lastp, "frames": "last"})).unwrap()["files"].as_array().unwrap().len(), 2); + let _ = std::fs::remove_dir_all(&dir); + } + fn solved() -> Session { let mut s = Session::default(); cantilever(&mut s); diff --git a/crates/engine/src/lib.rs b/crates/engine/src/lib.rs index ab22c35..738f650 100644 --- a/crates/engine/src/lib.rs +++ b/crates/engine/src/lib.rs @@ -545,6 +545,12 @@ mod tests { json!({"format": "nastran", "path": "/no/such/dir/part.bdf", "steps": vec!["x".repeat(1000); 2000]}), json!({"format": "nastran", "path": "/no/such/dir/part.bdf", "steps": {}}), json!({"path": 7}), + json!({"path": "/no/such/dir/out.vtu", "step": "Nope"}), + json!({"path": "/no/such/dir/out.vtu", "mode": 0}), + json!({"path": "/no/such/dir/out.pvd", "mode": 1_000_000_000}), + json!({"path": "/no/such/dir/out.pvd", "frames": 3, "encoding": []}), + json!({"path": "x\u{0}.vtu"}), + json!({"path": format!("/no/such/{}.vtu", "x".repeat(100_000))}), json!({"path": "x\u{0}.stl"}), json!({"path": "/no/such/dir/skin.stl", "set": 7}), json!({"path": "/no/such/dir/skin.stl", "set": "x".repeat(100_000)}), diff --git a/crates/io/src/vtk.rs b/crates/io/src/vtk.rs index a704c88..eaf653a 100644 --- a/crates/io/src/vtk.rs +++ b/crates/io/src/vtk.rs @@ -9,9 +9,11 @@ //! //! Points are the model's nodes in label order and cells its elements in label order; the point //! array `node_label` and the cell array `element_label` keep the labels. Element node order is -//! the model's ([`ElementKind`]), which is VTK's for every kind but the wedge: a FieldCraft -//! wedge's triangle 1-2-3 has its right-hand normal towards 4-5-6, a VTK wedge's base 0-1-2 -//! points away from 3-4-5, so wedges are written as 0, 2, 1, 3, 5, 4. +//! the model's ([`ElementKind`]), which is VTK's for every kind. That holds for the wedge too: +//! VTK's wedge has the parametric coordinates (0,0,0), (1,0,0), (0,1,0) for points 0-2 and the +//! same at ζ = 1 for 3-5, so the base 0-1-2's right-hand normal points to 3-4-5, as a FieldCraft +//! wedge's 1-2-3 points to 4-5-6 (VTK 9 computes a positive volume for it and lists its outward +//! base face as 0-2-1). //! //! Data arrays are written as text (`ascii`, numbers in the shortest form that reads back to the //! same `f64`) or inline base64 (`binary`, each array prefixed with its byte count as a `UInt64`, @@ -41,15 +43,6 @@ pub fn cell_type(kind: ElementKind) -> u8 { } } -/// The positions in an element's node list in the order VTK lists the cell's points. -pub fn vtk_order(kind: ElementKind) -> &'static [usize] { - const IDENTITY: [usize; 20] = [0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19]; - match kind { - ElementKind::Wedge6 => &[0, 2, 1, 3, 5, 4], - k => IDENTITY.get(..k.node_count()).unwrap_or(&[]), - } -} - /// Where a field's values are: one tuple per point (node) or per cell (element). #[derive(Clone, Copy, Debug, PartialEq, Eq)] pub enum Location { @@ -106,8 +99,7 @@ pub fn write_vtu(model: &Model, fields: &[Field], field_data: &[(&str, f64)], en let mut offsets: Vec = Vec::with_capacity(elements.len()); let mut types: Vec = Vec::with_capacity(elements.len()); for e in model.elements().values() { - for &i in vtk_order(e.kind) { - let label = e.nodes.get(i).copied().ok_or_else(|| Error::Invalid(format!("a {} lists too few nodes", e.kind.name())))?; + for &label in &e.nodes { connectivity.push(index(label)? as i64); } offsets.push(connectivity.len() as i64); @@ -438,8 +430,7 @@ mod tests { let labels = numbers(array(&xml, "node_label").1); let written: Vec = numbers(array(&xml, "connectivity").1).iter().map(|&i| labels[i as usize] as u32).collect(); let e = m.element(7).unwrap(); - let expected: Vec = vtk_order(kind).iter().map(|&i| e.nodes[i]).collect(); - assert_eq!(written, expected, "{kind:?}"); + assert_eq!(written, e.nodes, "{kind:?}"); } let ids: Vec = ElementKind::ALL.iter().map(|&k| cell_type(k)).collect(); // rod2, rod3, beam2, tri3, tri6, quad4, quad8, tet4, tet10, wedge6, hex8, hex20 @@ -466,27 +457,25 @@ mod tests { [0, 1, 2].map(|k| p.iter().map(|q| q[k]).sum::() / n) } - /// The VTK document's convention for each 3D cell: the base face's right-hand normal points - /// towards the rest of the cell for the tetrahedron and the hexahedron, and away from it - /// (out of the cell) for the wedge. The signed measure is positive when the written cell - /// follows it. + /// VTK's convention for each 3D cell, from its parametric coordinates: the base face's + /// right-hand normal points to the rest of the cell (tetrahedron 0-1-2 to 3, hexahedron + /// 0-1-2-3 to 4-7, wedge 0-1-2 to 3-5). The signed measure is positive when the written cell + /// follows it; VTK 9.7's vtkCellSizeFilter and vtkMeshQuality agree (a unit wedge written + /// this way measures +0.5, reordered 0-2-1-3-5-4 it measures -0.5). #[test] fn solid_cells_have_positive_volume_in_vtk_convention() { - for (kind, base, rest, towards) in [ - (ElementKind::Tet4, 3, 3..4, true), - (ElementKind::Tet10, 3, 3..4, true), - (ElementKind::Hex8, 4, 4..8, true), - (ElementKind::Hex20, 4, 4..8, true), - (ElementKind::Wedge6, 3, 3..6, false), + for (kind, base, rest) in [ + (ElementKind::Tet4, 3, 3..4), + (ElementKind::Tet10, 3, 3..4), + (ElementKind::Hex8, 4, 4..8), + (ElementKind::Hex20, 4, 4..8), + (ElementKind::Wedge6, 3, 3..6), ] { let p = cell_points(&write_vtu(&one(kind), &[], &[], Encoding::Ascii).unwrap()); let b = &p[..base]; let n = cross(sub(b[1], b[0]), sub(b[2], b[0])); - let d = dot(n, sub(centre(&p[rest]), centre(b))); - assert!(if towards { d > 0.0 } else { d < 0.0 }, "{kind:?}: {d}"); + assert!(dot(n, sub(centre(&p[rest]), centre(b))) > 0.0, "{kind:?}"); } - // Without the wedge's reorder the base would point the wrong way. - assert_eq!(vtk_order(ElementKind::Wedge6), [0, 2, 1, 3, 5, 4]); } #[test] @@ -496,10 +485,7 @@ mod tests { (ElementKind::Tri6, &[(0, 1), (1, 2), (2, 0)]), (ElementKind::Quad8, &[(0, 1), (1, 2), (2, 3), (3, 0)]), (ElementKind::Tet10, &[(0, 1), (1, 2), (2, 0), (0, 3), (1, 3), (2, 3)]), - ( - ElementKind::Hex20, - &[(0, 1), (1, 2), (2, 3), (3, 0), (4, 5), (5, 6), (6, 7), (7, 4), (0, 4), (1, 5), (2, 6), (3, 7)], - ), + (ElementKind::Hex20, &[(0, 1), (1, 2), (2, 3), (3, 0), (4, 5), (5, 6), (6, 7), (7, 4), (0, 4), (1, 5), (2, 6), (3, 7)]), ]; for (kind, list) in edges { let p = cell_points(&write_vtu(&one(kind), &[], &[], Encoding::Ascii).unwrap()); diff --git a/crates/model/src/mesh.rs b/crates/model/src/mesh.rs index 104ec4a..fa0b576 100644 --- a/crates/model/src/mesh.rs +++ b/crates/model/src/mesh.rs @@ -20,9 +20,8 @@ pub struct Node { /// - hex8/hex20: face 1-2-3-4, then 5-6-7-8 above it; hex20's mid-side nodes 9–12 lie on the /// edges of 1-2-3-4, 13–16 on those of 5-6-7-8, and 17–20 on edges 1-5, 2-6, 3-7, 4-8. /// -/// This is the Abaqus order, and VTK's for every kind but wedge6: a VTK wedge's base 0-1-2 has its -/// right-hand normal pointing away from 3-4-5, so the VTK writer lists a wedge as 1-3-2-4-6-5. -/// Formats that order nodes differently, such as Nastran's +/// This is the Abaqus (and VTK) order; a VTK wedge too has its base 0-1-2 normal pointing to +/// 3-4-5 (its parametric coordinates put 0-1-2 at ζ = 0 counter-clockwise seen from ζ = 1). Formats that order nodes differently, such as Nastran's /// CHEXA with 20 nodes (vertical edges before the top face's edges), are reordered by their deck /// readers and writers; the model never stores another order. #[derive(Debug, Clone, Copy, PartialEq, Eq, PartialOrd, Ord, Hash, Serialize, Deserialize)] diff --git a/docs/features.md b/docs/features.md index 4e706bd..303c045 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. -74 ready, 39 coming soon. +75 ready, 39 coming soon. ## Geometry @@ -153,6 +153,7 @@ the control channel and AI agents over MCP (docs/mcp.md) run the same commands. | Plots | Cut — Section cuts through the model | coming soon (M9) | | | Plots | Display Group — Show and hide parts of the model | coming soon (M9) | | | Plots | Streamlines — Streamlines and flow fields | coming soon (M12) | | +| Files | Export VTU — Write a solved step as VTK unstructured grids for ParaView | ready | `results.export.vtu` | | Files | Read Results — Read results in public formats | coming soon (M9) | | | Files | Convergence — Mesh convergence studies | coming soon (M9) | |