From 520188a82a0f4139a918fff96862734097cd7caa Mon Sep 17 00:00:00 2001 From: jess Date: Wed, 13 May 2026 01:14:06 -0700 Subject: [PATCH] demo --- assets/belasolv.svg | 47 +++- assets/idr_femmplottype.svg | 18 +- crates/femm-app/src/main.rs | 99 +++++++- crates/femm-doc-curr/src/lib.rs | 2 + crates/femm-doc-curr/src/mesh.rs | 161 ++++++++++++ crates/femm-doc-curr/src/poly.rs | 166 ++++++++++++ crates/femm-doc-elec/src/lib.rs | 2 + crates/femm-doc-elec/src/mesh.rs | 161 ++++++++++++ crates/femm-doc-elec/src/poly.rs | 166 ++++++++++++ crates/femm-doc-heat/src/lib.rs | 2 + crates/femm-doc-heat/src/mesh.rs | 161 ++++++++++++ crates/femm-doc-heat/src/poly.rs | 166 ++++++++++++ crates/femm-doc-mag/src/ans.rs | 422 +++++++++++++++++++++++++++++++ crates/femm-doc-mag/src/lib.rs | 1 + crates/femm-doc-mag/src/mesh.rs | 6 +- scripts/macos/build.sh | 10 + 16 files changed, 1569 insertions(+), 21 deletions(-) create mode 100644 crates/femm-doc-curr/src/mesh.rs create mode 100644 crates/femm-doc-curr/src/poly.rs create mode 100644 crates/femm-doc-elec/src/mesh.rs create mode 100644 crates/femm-doc-elec/src/poly.rs create mode 100644 crates/femm-doc-heat/src/mesh.rs create mode 100644 crates/femm-doc-heat/src/poly.rs create mode 100644 crates/femm-doc-mag/src/ans.rs diff --git a/assets/belasolv.svg b/assets/belasolv.svg index 5a8ad0c..bb59796 100644 --- a/assets/belasolv.svg +++ b/assets/belasolv.svg @@ -1,9 +1,38 @@ - - - - - - - - - \ No newline at end of file + + + + + + + + + + diff --git a/assets/idr_femmplottype.svg b/assets/idr_femmplottype.svg index c86a47f..e935b38 100644 --- a/assets/idr_femmplottype.svg +++ b/assets/idr_femmplottype.svg @@ -1,7 +1,21 @@ - - + + + + + + + + + + + + + + + + \ No newline at end of file diff --git a/crates/femm-app/src/main.rs b/crates/femm-app/src/main.rs index 0185351..5af1b63 100644 --- a/crates/femm-app/src/main.rs +++ b/crates/femm-app/src/main.rs @@ -4,9 +4,12 @@ mod doc_canvas; use doc_canvas::{CanvasMessage, Tool}; use femm_doc_mag::FemmDoc; +use femm_doc_mag::ans::MagSolution; use femm_doc_mag::mesh::Mesh; use iced::widget::{button, column, container, row, text}; use iced::{Element, Length, Task}; +use std::ffi::CString; +use std::path::Path; const DEMO_FEM: &str = include_str!("../assets/demo.fem"); const ADD_TOLERANCE: f64 = 0.5; @@ -18,11 +21,13 @@ enum Message { SelectTool(Tool), Canvas(CanvasMessage), RunMesh, + RunAnalyze, } struct App { doc: FemmDoc, mesh: Option, + solution: Option, source_label: String, status: String, tool: Tool, @@ -34,6 +39,7 @@ impl App { let app = App { doc, mesh: None, + solution: None, source_label: String::from("demo.fem (embedded)"), status: String::new(), tool: Tool::Select, @@ -63,6 +69,7 @@ impl App { Ok(d) => { self.doc = d; self.mesh = None; + self.solution = None; self.source_label = label; self.status = String::new(); } @@ -83,20 +90,52 @@ impl App { } Err(e) => { self.mesh = None; + self.solution = None; self.status = format!("mesh failed: {e}"); } } } + Message::RunAnalyze => { + self.solution = None; + match run_mesh(&self.doc) { + Ok(m) => self.mesh = Some(m), + Err(e) => { + self.mesh = None; + self.status = format!("mesh failed: {e}"); + return Task::none(); + } + } + let stem = active_stem(); + if let Err(e) = run_solve(&stem) { + self.status = format!("solve failed: {e}"); + return Task::none(); + } + let ans_path = stem.with_extension("ans"); + match MagSolution::open(&ans_path) { + Ok(sol) => { + self.status = format!( + "solved: {} mesh nodes, {} elements", + sol.mesh_nodes.len(), sol.mesh_elements.len(), + ); + self.solution = Some(sol); + } + Err(e) => { + self.status = format!("read .ans failed: {e}"); + } + } + } Message::Canvas(CanvasMessage::Click { world, tool }) => { match tool { Tool::AddNode => { let idx = self.doc.add_node(world.0, world.1, ADD_TOLERANCE); self.mesh = None; + self.solution = None; self.status = format!("node {idx} at ({:.3}, {:.3})", world.0, world.1); } Tool::AddBlockLabel => { let idx = self.doc.add_block_label(world.0, world.1, ADD_TOLERANCE); self.mesh = None; + self.solution = None; self.status = format!("block label {idx} at ({:.3}, {:.3})", world.0, world.1); } Tool::Select | Tool::AddSegment => {} @@ -153,6 +192,7 @@ impl App { tool_button("Add Segment", Tool::AddSegment, self.tool), tool_button("Add Label", Tool::AddBlockLabel, self.tool), button("Mesh").on_press(Message::RunMesh), + button("Analyze").on_press(Message::RunAnalyze), text(&self.source_label).size(13), stats, ] @@ -174,11 +214,16 @@ impl App { } } +/// stem path used by all mesh and solve operations for the active doc. +fn active_stem() -> std::path::PathBuf { + let dir = std::env::temp_dir().join("femm42-mesh"); + std::fs::create_dir_all(&dir).ok(); + dir.join("active") +} + /// saves the doc to a temp .fem and matching .poly, invokes Triangle, and reads back the mesh. fn run_mesh(doc: &FemmDoc) -> Result { - let dir = std::env::temp_dir().join("femm42-mesh"); - std::fs::create_dir_all(&dir).map_err(|e| format!("temp dir: {e}"))?; - let stem = dir.join("active"); + let stem = active_stem(); let fem_path = stem.with_extension("fem"); let poly_path = stem.with_extension("poly"); @@ -204,11 +249,51 @@ fn run_mesh(doc: &FemmDoc) -> Result { Mesh::load(&stem).map_err(|e| format!("read mesh: {e}")) } -/// resolves the Triangle binary path, preferring the build dir then $PATH. +/// runs the magnetostatic solver FFI pipeline against an on-disk .fem + .node + .ele rooted at `stem`. +fn run_solve(stem: &Path) -> Result<(), String> { + let stem_str = stem.to_str().ok_or_else(|| String::from("non-utf8 path"))?; + let cstem = CString::new(stem_str).map_err(|e| format!("path: {e}"))?; + + unsafe { + let doc = femm_sys::femm_mag_doc_new(); + if doc.is_null() { + return Err(String::from("femm_mag_doc_new returned null")); + } + let load = femm_sys::femm_mag_doc_load_fem(doc, cstem.as_ptr()); + if load == 0 { + femm_sys::femm_mag_doc_free(doc); + return Err(String::from("load_fem returned 0")); + } + if femm_sys::femm_mag_doc_load_mesh(doc) == 0 { + femm_sys::femm_mag_doc_free(doc); + return Err(String::from("load_mesh returned 0")); + } + if femm_sys::femm_mag_doc_renumber(doc) == 0 { + femm_sys::femm_mag_doc_free(doc); + return Err(String::from("renumber returned 0")); + } + if femm_sys::femm_mag_doc_solve(doc) == 0 { + femm_sys::femm_mag_doc_free(doc); + return Err(String::from("solve returned 0")); + } + femm_sys::femm_mag_doc_free(doc); + } + Ok(()) +} + +/// resolves the Triangle binary path. +/// search order: sibling of the running exe (bundled .app), repo build dir under target/.., FEMM_TRIANGLE. fn locate_triangle() -> Option { - let here = std::env::current_dir().ok()?; - let built = here.join("build/triangle/triangle"); - if built.is_file() { return Some(built); } + if let Ok(exe) = std::env::current_exe() { + if let Some(parent) = exe.parent() { + let bundled = parent.join("triangle"); + if bundled.is_file() { return Some(bundled); } + let dev = parent.join("../../build/triangle/triangle"); + if let Ok(canon) = dev.canonicalize() { + if canon.is_file() { return Some(canon); } + } + } + } if let Ok(env_path) = std::env::var("FEMM_TRIANGLE") { let p = std::path::PathBuf::from(env_path); if p.is_file() { return Some(p); } diff --git a/crates/femm-doc-curr/src/lib.rs b/crates/femm-doc-curr/src/lib.rs index 8f61f72..ecd5fe0 100644 --- a/crates/femm-doc-curr/src/lib.rs +++ b/crates/femm-doc-curr/src/lib.rs @@ -6,6 +6,8 @@ pub mod props; pub mod parser; pub mod writer; pub mod edit; +pub mod poly; +pub mod mesh; use num_complex::Complex64; diff --git a/crates/femm-doc-curr/src/mesh.rs b/crates/femm-doc-curr/src/mesh.rs new file mode 100644 index 0000000..c85669a --- /dev/null +++ b/crates/femm-doc-curr/src/mesh.rs @@ -0,0 +1,161 @@ +//! readers for Triangle's .node and .ele output files. + +use std::path::Path; +use thiserror::Error; + +/// 2D triangular mesh consumed by the canvas overlay and by the post-processor. +#[derive(Debug, Default, Clone)] +pub struct Mesh { + pub nodes: Vec, + pub elements: Vec, +} + +/// position and optional boundary marker for a single mesh node. +#[derive(Debug, Clone, Copy)] +pub struct MeshNode { + pub x: f64, + pub y: f64, + pub marker: i32, +} + +/// triangle element with three zero-based node indices and an optional region attribute. +#[derive(Debug, Clone, Copy)] +pub struct MeshElement { + pub v0: u32, + pub v1: u32, + pub v2: u32, + pub attribute: i32, +} + +#[derive(Debug, Error)] +pub enum MeshLoadError { + #[error("io error: {0}")] + Io(#[from] std::io::Error), + #[error("malformed header in {file}: {reason}")] + BadHeader { file: &'static str, reason: String }, + #[error("malformed row in {file}: {reason}")] + BadRow { file: &'static str, reason: String }, +} + +impl Mesh { + /// loads .node and .ele next to the given path stem. + pub fn load(stem: impl AsRef) -> Result { + let stem = stem.as_ref(); + let node_path = with_extension(stem, "node"); + let ele_path = with_extension(stem, "ele"); + let nodes = parse_node(&std::fs::read_to_string(&node_path)?)?; + let elements = parse_ele(&std::fs::read_to_string(&ele_path)?)?; + Ok(Mesh { nodes, elements }) + } +} + +/// parses a .node file: header `count 2 0 nbm`, then count rows of `idx x y [marker]`. +pub fn parse_node(src: &str) -> Result, MeshLoadError> { + let mut lines = src.lines().filter(|l| !l.trim_start().starts_with('#') && !l.trim().is_empty()); + let header = lines.next().ok_or(MeshLoadError::BadHeader { + file: ".node", + reason: "empty file".into(), + })?; + let mut hparts = header.split_whitespace(); + let count: usize = hparts.next() + .and_then(|s| s.parse().ok()) + .ok_or(MeshLoadError::BadHeader { file: ".node", reason: format!("expected count, got {header:?}") })?; + let _dim = hparts.next(); + let _attrs = hparts.next(); + let has_marker = hparts.next().map(|s| s.trim() == "1").unwrap_or(false); + + let mut nodes = Vec::with_capacity(count); + for line in lines.take(count) { + let mut parts = line.split_whitespace(); + let _idx: u32 = parts.next().and_then(|s| s.parse().ok()) + .ok_or(MeshLoadError::BadRow { file: ".node", reason: line.into() })?; + let x: f64 = parts.next().and_then(|s| s.parse().ok()) + .ok_or(MeshLoadError::BadRow { file: ".node", reason: line.into() })?; + let y: f64 = parts.next().and_then(|s| s.parse().ok()) + .ok_or(MeshLoadError::BadRow { file: ".node", reason: line.into() })?; + let marker = if has_marker { + parts.next().and_then(|s| s.parse().ok()).unwrap_or(0) + } else { + 0 + }; + nodes.push(MeshNode { x, y, marker }); + } + Ok(nodes) +} + +/// parses a .ele file: header `count 3 natt`, then count rows of `idx v0 v1 v2 [attr0 ...]`. +pub fn parse_ele(src: &str) -> Result, MeshLoadError> { + let mut lines = src.lines().filter(|l| !l.trim_start().starts_with('#') && !l.trim().is_empty()); + let header = lines.next().ok_or(MeshLoadError::BadHeader { + file: ".ele", + reason: "empty file".into(), + })?; + let mut hparts = header.split_whitespace(); + let count: usize = hparts.next() + .and_then(|s| s.parse().ok()) + .ok_or(MeshLoadError::BadHeader { file: ".ele", reason: format!("expected count, got {header:?}") })?; + let _nodes_per = hparts.next(); + let natt: usize = hparts.next().and_then(|s| s.parse().ok()).unwrap_or(0); + + let mut elements = Vec::with_capacity(count); + for line in lines.take(count) { + let mut parts = line.split_whitespace(); + let _idx: u32 = parts.next().and_then(|s| s.parse().ok()) + .ok_or(MeshLoadError::BadRow { file: ".ele", reason: line.into() })?; + let v0: u32 = parts.next().and_then(|s| s.parse().ok()) + .ok_or(MeshLoadError::BadRow { file: ".ele", reason: line.into() })?; + let v1: u32 = parts.next().and_then(|s| s.parse().ok()) + .ok_or(MeshLoadError::BadRow { file: ".ele", reason: line.into() })?; + let v2: u32 = parts.next().and_then(|s| s.parse().ok()) + .ok_or(MeshLoadError::BadRow { file: ".ele", reason: line.into() })?; + let attribute = if natt > 0 { + parts.next().and_then(|s| s.parse::().ok()).map(|f| f as i32).unwrap_or(0) + } else { + 0 + }; + elements.push(MeshElement { v0, v1, v2, attribute }); + } + Ok(elements) +} + +fn with_extension(stem: &Path, ext: &str) -> std::path::PathBuf { + let mut s = stem.to_path_buf(); + s.set_extension(ext); + s +} + +#[cfg(test)] +mod tests { + use super::*; + + #[test] + fn parses_node_with_markers() { + // 3 nodes, 2D, no attributes, with markers; comment line at the start is ignored. + let src = "# triangle .node\n3 2 0 1\n0 0.0 0.0 2\n1 1.0 0.0 0\n2 0.0 1.0 0\n"; + let nodes = parse_node(src).unwrap(); + assert_eq!(nodes.len(), 3); + assert_eq!(nodes[0].marker, 2); + assert!((nodes[2].y - 1.0).abs() < 1e-12); + } + + #[test] + fn parses_ele_with_region_attribute() { + // 2 triangles, 3 nodes per element, 1 region attribute column. + let src = "2 3 1\n0 0 1 2 1.0\n1 1 2 3 2.0\n"; + let els = parse_ele(src).unwrap(); + assert_eq!(els.len(), 2); + assert_eq!(els[0].v0, 0); + assert_eq!(els[0].v1, 1); + assert_eq!(els[0].v2, 2); + assert_eq!(els[0].attribute, 1); + assert_eq!(els[1].attribute, 2); + } + + #[test] + fn parses_node_without_marker_column() { + let src = "2 2 0 0\n0 0.0 0.0\n1 1.0 0.0\n"; + let nodes = parse_node(src).unwrap(); + assert_eq!(nodes.len(), 2); + assert_eq!(nodes[0].marker, 0); + } +} diff --git a/crates/femm-doc-curr/src/poly.rs b/crates/femm-doc-curr/src/poly.rs new file mode 100644 index 0000000..ae5949b --- /dev/null +++ b/crates/femm-doc-curr/src/poly.rs @@ -0,0 +1,166 @@ +//! .poly emitter consumed by the Triangle mesher. + +use crate::FemmDoc; +use std::fmt::Write; + +impl FemmDoc { + /// renders the doc geometry to the Triangle .poly text format, returning the file contents. + pub fn write_poly(&self) -> String { + let mut out = String::new(); + let bbox_diag = self.bbox_diagonal(); + let default_area = if bbox_diag > 0.0 { bbox_diag } else { -1.0 }; + + writeln!(out, "{}\t2\t0\t1", self.nodes.len()).unwrap(); + for (i, n) in self.nodes.iter().enumerate() { + let marker = point_marker_index(self, &n.boundary_marker); + writeln!(out, "{i}\t{:.17e}\t{:.17e}\t{marker}", n.x, n.y).unwrap(); + } + + writeln!(out, "{}\t1", self.segments.len()).unwrap(); + for (i, s) in self.segments.iter().enumerate() { + let marker = -boundary_marker_index(self, &s.boundary_marker); + writeln!(out, "{i}\t{}\t{}\t{marker}", s.n0, s.n1).unwrap(); + } + + let holes: Vec<&_> = self.block_labels.iter().filter(|b| b.block_type == "").collect(); + writeln!(out, "{}", holes.len()).unwrap(); + for (i, h) in holes.iter().enumerate() { + writeln!(out, "{i}\t{:.17e}\t{:.17e}", h.x, h.y).unwrap(); + } + + let regions: Vec<&_> = self.block_labels.iter().filter(|b| b.block_type != "").collect(); + writeln!(out, "{}", regions.len()).unwrap(); + for (i, r) in regions.iter().enumerate() { + let attr = (i as i32) + 1; + let max_area = if r.max_area > 0.0 && r.max_area < default_area { + r.max_area + } else { + default_area + }; + writeln!(out, "{i}\t{:.17e}\t{:.17e}\t{attr}\t{:.17e}", r.x, r.y, max_area).unwrap(); + } + + out + } + + /// writes the .poly text to disk at the given path. + pub fn save_poly(&self, path: impl AsRef) -> std::io::Result<()> { + std::fs::write(path, self.write_poly()) + } + + fn bbox_diagonal(&self) -> f64 { + if self.nodes.len() < 2 { return 0.0; } + let mut xmin = f64::INFINITY; + let mut xmax = f64::NEG_INFINITY; + let mut ymin = f64::INFINITY; + let mut ymax = f64::NEG_INFINITY; + for n in &self.nodes { + xmin = xmin.min(n.x); xmax = xmax.max(n.x); + ymin = ymin.min(n.y); ymax = ymax.max(n.y); + } + let dx = xmax - xmin; + let dy = ymax - ymin; + (dx * dx + dy * dy).sqrt() + } +} + +/// resolves a node's boundary-marker name to Triangle's marker integer (0 for unmarked). +fn point_marker_index(doc: &FemmDoc, name: &str) -> i32 { + if name.is_empty() { return 0; } + for (i, p) in doc.points.iter().enumerate() { + if p.name == name { return (i as i32) + 2; } + } + 0 +} + +/// resolves a segment's boundary-marker name to Triangle's marker integer, returned positive. +/// the caller negates for Triangle's segment-marker convention. +fn boundary_marker_index(doc: &FemmDoc, name: &str) -> i32 { + if name.is_empty() { return 0; } + for (i, b) in doc.boundaries.iter().enumerate() { + if b.name == name { return (i as i32) + 2; } + } + 0 +} + +#[cfg(test)] +mod tests { + use super::*; + use crate::{BoundaryProp, PointProp, BlockLabel}; + + #[test] + fn write_poly_header_lines_count_correctly() { + // four-corner square with one block label inside; verify the section headers and counts. + let mut d = FemmDoc::default(); + d.points.push(PointProp { name: "V=0".into(), ..PointProp::default() }); + d.boundaries.push(BoundaryProp { name: "outer".into(), ..BoundaryProp::default() }); + d.add_node(0.0, 0.0, 0.0); + d.add_node(1.0, 0.0, 0.0); + d.add_node(1.0, 1.0, 0.0); + d.add_node(0.0, 1.0, 0.0); + d.nodes[0].boundary_marker = "V=0".into(); + d.add_segment_with_marker(0, 1, "outer"); + d.add_segment_with_marker(1, 2, "outer"); + d.add_segment_with_marker(2, 3, "outer"); + d.add_segment_with_marker(3, 0, "outer"); + d.block_labels.push(BlockLabel { + x: 0.5, y: 0.5, + max_area: 0.0, + block_type: "Air".into(), + ..BlockLabel::default() + }); + + let poly = d.write_poly(); + let lines: Vec<&str> = poly.lines().collect(); + + assert!(lines[0].starts_with("4\t2\t0\t1"), "node header: {}", lines[0]); + assert!(lines[5].starts_with("4\t1"), "segment header: {}", lines[5]); + assert_eq!(lines[10], "0", "no holes section: {}", lines[10]); + assert_eq!(lines[11], "1", "one region: {}", lines[11]); + } + + #[test] + fn write_poly_marks_nodes_and_segments_by_property_index() { + // first point property listed -> marker 2; first boundary property listed -> -2 on segment. + let mut d = FemmDoc::default(); + d.points.push(PointProp { name: "V=0".into(), ..PointProp::default() }); + d.boundaries.push(BoundaryProp { name: "outer".into(), ..BoundaryProp::default() }); + d.add_node(0.0, 0.0, 0.0); + d.add_node(1.0, 0.0, 0.0); + d.nodes[0].boundary_marker = "V=0".into(); + d.add_segment_with_marker(0, 1, "outer"); + + let poly = d.write_poly(); + let lines: Vec<&str> = poly.lines().collect(); + + assert!(lines[1].ends_with("\t2"), "node 0 marker = 2: {}", lines[1]); + assert!(lines[2].ends_with("\t0"), "node 1 marker = 0: {}", lines[2]); + assert!(lines[4].ends_with("\t-2"), "segment 0 marker = -2: {}", lines[4]); + } + + #[test] + fn write_poly_holes_partitioned_from_regions() { + // one "" label counts as a hole; everything else is a region. + let mut d = FemmDoc::default(); + d.add_node(0.0, 0.0, 0.0); + d.add_node(1.0, 0.0, 0.0); + d.block_labels.push(BlockLabel { + x: 0.5, y: 0.5, + block_type: "".into(), + ..BlockLabel::default() + }); + d.block_labels.push(BlockLabel { + x: 0.2, y: 0.2, + block_type: "Air".into(), + ..BlockLabel::default() + }); + + let poly = d.write_poly(); + let lines: Vec<&str> = poly.lines().collect(); + + // node header (0), 2 node rows (1,2), segment header "0\t1" (3), + // hole count (4), 1 hole row (5), region count (6), 1 region row (7). + assert_eq!(lines[4], "1", "one hole: {}", lines[4]); + assert_eq!(lines[6], "1", "one region: {}", lines[6]); + } +} diff --git a/crates/femm-doc-elec/src/lib.rs b/crates/femm-doc-elec/src/lib.rs index 4e1b98c..085cc7b 100644 --- a/crates/femm-doc-elec/src/lib.rs +++ b/crates/femm-doc-elec/src/lib.rs @@ -6,6 +6,8 @@ pub mod props; pub mod parser; pub mod writer; pub mod edit; +pub mod poly; +pub mod mesh; pub use geom::{ArcSegment, BlockLabel, Node, Segment}; pub use props::{BoundaryProp, ConductorProp, MaterialProp, PointProp}; diff --git a/crates/femm-doc-elec/src/mesh.rs b/crates/femm-doc-elec/src/mesh.rs new file mode 100644 index 0000000..c85669a --- /dev/null +++ b/crates/femm-doc-elec/src/mesh.rs @@ -0,0 +1,161 @@ +//! readers for Triangle's .node and .ele output files. + +use std::path::Path; +use thiserror::Error; + +/// 2D triangular mesh consumed by the canvas overlay and by the post-processor. +#[derive(Debug, Default, Clone)] +pub struct Mesh { + pub nodes: Vec, + pub elements: Vec, +} + +/// position and optional boundary marker for a single mesh node. +#[derive(Debug, Clone, Copy)] +pub struct MeshNode { + pub x: f64, + pub y: f64, + pub marker: i32, +} + +/// triangle element with three zero-based node indices and an optional region attribute. +#[derive(Debug, Clone, Copy)] +pub struct MeshElement { + pub v0: u32, + pub v1: u32, + pub v2: u32, + pub attribute: i32, +} + +#[derive(Debug, Error)] +pub enum MeshLoadError { + #[error("io error: {0}")] + Io(#[from] std::io::Error), + #[error("malformed header in {file}: {reason}")] + BadHeader { file: &'static str, reason: String }, + #[error("malformed row in {file}: {reason}")] + BadRow { file: &'static str, reason: String }, +} + +impl Mesh { + /// loads .node and .ele next to the given path stem. + pub fn load(stem: impl AsRef) -> Result { + let stem = stem.as_ref(); + let node_path = with_extension(stem, "node"); + let ele_path = with_extension(stem, "ele"); + let nodes = parse_node(&std::fs::read_to_string(&node_path)?)?; + let elements = parse_ele(&std::fs::read_to_string(&ele_path)?)?; + Ok(Mesh { nodes, elements }) + } +} + +/// parses a .node file: header `count 2 0 nbm`, then count rows of `idx x y [marker]`. +pub fn parse_node(src: &str) -> Result, MeshLoadError> { + let mut lines = src.lines().filter(|l| !l.trim_start().starts_with('#') && !l.trim().is_empty()); + let header = lines.next().ok_or(MeshLoadError::BadHeader { + file: ".node", + reason: "empty file".into(), + })?; + let mut hparts = header.split_whitespace(); + let count: usize = hparts.next() + .and_then(|s| s.parse().ok()) + .ok_or(MeshLoadError::BadHeader { file: ".node", reason: format!("expected count, got {header:?}") })?; + let _dim = hparts.next(); + let _attrs = hparts.next(); + let has_marker = hparts.next().map(|s| s.trim() == "1").unwrap_or(false); + + let mut nodes = Vec::with_capacity(count); + for line in lines.take(count) { + let mut parts = line.split_whitespace(); + let _idx: u32 = parts.next().and_then(|s| s.parse().ok()) + .ok_or(MeshLoadError::BadRow { file: ".node", reason: line.into() })?; + let x: f64 = parts.next().and_then(|s| s.parse().ok()) + .ok_or(MeshLoadError::BadRow { file: ".node", reason: line.into() })?; + let y: f64 = parts.next().and_then(|s| s.parse().ok()) + .ok_or(MeshLoadError::BadRow { file: ".node", reason: line.into() })?; + let marker = if has_marker { + parts.next().and_then(|s| s.parse().ok()).unwrap_or(0) + } else { + 0 + }; + nodes.push(MeshNode { x, y, marker }); + } + Ok(nodes) +} + +/// parses a .ele file: header `count 3 natt`, then count rows of `idx v0 v1 v2 [attr0 ...]`. +pub fn parse_ele(src: &str) -> Result, MeshLoadError> { + let mut lines = src.lines().filter(|l| !l.trim_start().starts_with('#') && !l.trim().is_empty()); + let header = lines.next().ok_or(MeshLoadError::BadHeader { + file: ".ele", + reason: "empty file".into(), + })?; + let mut hparts = header.split_whitespace(); + let count: usize = hparts.next() + .and_then(|s| s.parse().ok()) + .ok_or(MeshLoadError::BadHeader { file: ".ele", reason: format!("expected count, got {header:?}") })?; + let _nodes_per = hparts.next(); + let natt: usize = hparts.next().and_then(|s| s.parse().ok()).unwrap_or(0); + + let mut elements = Vec::with_capacity(count); + for line in lines.take(count) { + let mut parts = line.split_whitespace(); + let _idx: u32 = parts.next().and_then(|s| s.parse().ok()) + .ok_or(MeshLoadError::BadRow { file: ".ele", reason: line.into() })?; + let v0: u32 = parts.next().and_then(|s| s.parse().ok()) + .ok_or(MeshLoadError::BadRow { file: ".ele", reason: line.into() })?; + let v1: u32 = parts.next().and_then(|s| s.parse().ok()) + .ok_or(MeshLoadError::BadRow { file: ".ele", reason: line.into() })?; + let v2: u32 = parts.next().and_then(|s| s.parse().ok()) + .ok_or(MeshLoadError::BadRow { file: ".ele", reason: line.into() })?; + let attribute = if natt > 0 { + parts.next().and_then(|s| s.parse::().ok()).map(|f| f as i32).unwrap_or(0) + } else { + 0 + }; + elements.push(MeshElement { v0, v1, v2, attribute }); + } + Ok(elements) +} + +fn with_extension(stem: &Path, ext: &str) -> std::path::PathBuf { + let mut s = stem.to_path_buf(); + s.set_extension(ext); + s +} + +#[cfg(test)] +mod tests { + use super::*; + + #[test] + fn parses_node_with_markers() { + // 3 nodes, 2D, no attributes, with markers; comment line at the start is ignored. + let src = "# triangle .node\n3 2 0 1\n0 0.0 0.0 2\n1 1.0 0.0 0\n2 0.0 1.0 0\n"; + let nodes = parse_node(src).unwrap(); + assert_eq!(nodes.len(), 3); + assert_eq!(nodes[0].marker, 2); + assert!((nodes[2].y - 1.0).abs() < 1e-12); + } + + #[test] + fn parses_ele_with_region_attribute() { + // 2 triangles, 3 nodes per element, 1 region attribute column. + let src = "2 3 1\n0 0 1 2 1.0\n1 1 2 3 2.0\n"; + let els = parse_ele(src).unwrap(); + assert_eq!(els.len(), 2); + assert_eq!(els[0].v0, 0); + assert_eq!(els[0].v1, 1); + assert_eq!(els[0].v2, 2); + assert_eq!(els[0].attribute, 1); + assert_eq!(els[1].attribute, 2); + } + + #[test] + fn parses_node_without_marker_column() { + let src = "2 2 0 0\n0 0.0 0.0\n1 1.0 0.0\n"; + let nodes = parse_node(src).unwrap(); + assert_eq!(nodes.len(), 2); + assert_eq!(nodes[0].marker, 0); + } +} diff --git a/crates/femm-doc-elec/src/poly.rs b/crates/femm-doc-elec/src/poly.rs new file mode 100644 index 0000000..f752688 --- /dev/null +++ b/crates/femm-doc-elec/src/poly.rs @@ -0,0 +1,166 @@ +//! .poly emitter consumed by the Triangle mesher. + +use crate::FemmDoc; +use std::fmt::Write; + +impl FemmDoc { + /// renders the doc geometry to the Triangle .poly text format, returning the file contents. + pub fn write_poly(&self) -> String { + let mut out = String::new(); + let bbox_diag = self.bbox_diagonal(); + let default_area = if bbox_diag > 0.0 { bbox_diag } else { -1.0 }; + + writeln!(out, "{}\t2\t0\t1", self.nodes.len()).unwrap(); + for (i, n) in self.nodes.iter().enumerate() { + let marker = point_marker_index(self, &n.boundary_marker); + writeln!(out, "{i}\t{:.17e}\t{:.17e}\t{marker}", n.x, n.y).unwrap(); + } + + writeln!(out, "{}\t1", self.segments.len()).unwrap(); + for (i, s) in self.segments.iter().enumerate() { + let marker = -boundary_marker_index(self, &s.boundary_marker); + writeln!(out, "{i}\t{}\t{}\t{marker}", s.n0, s.n1).unwrap(); + } + + let holes: Vec<&_> = self.block_labels.iter().filter(|b| b.block_type == "").collect(); + writeln!(out, "{}", holes.len()).unwrap(); + for (i, h) in holes.iter().enumerate() { + writeln!(out, "{i}\t{:.17e}\t{:.17e}", h.x, h.y).unwrap(); + } + + let regions: Vec<&_> = self.block_labels.iter().filter(|b| b.block_type != "").collect(); + writeln!(out, "{}", regions.len()).unwrap(); + for (i, r) in regions.iter().enumerate() { + let attr = (i as i32) + 1; + let max_area = if r.max_area > 0.0 && r.max_area < default_area { + r.max_area + } else { + default_area + }; + writeln!(out, "{i}\t{:.17e}\t{:.17e}\t{attr}\t{:.17e}", r.x, r.y, max_area).unwrap(); + } + + out + } + + /// writes the .poly text to disk at the given path. + pub fn save_poly(&self, path: impl AsRef) -> std::io::Result<()> { + std::fs::write(path, self.write_poly()) + } + + fn bbox_diagonal(&self) -> f64 { + if self.nodes.len() < 2 { return 0.0; } + let mut xmin = f64::INFINITY; + let mut xmax = f64::NEG_INFINITY; + let mut ymin = f64::INFINITY; + let mut ymax = f64::NEG_INFINITY; + for n in &self.nodes { + xmin = xmin.min(n.x); xmax = xmax.max(n.x); + ymin = ymin.min(n.y); ymax = ymax.max(n.y); + } + let dx = xmax - xmin; + let dy = ymax - ymin; + (dx * dx + dy * dy).sqrt() + } +} + +/// resolves a node's boundary-marker name to Triangle's marker integer (0 for unmarked). +fn point_marker_index(doc: &FemmDoc, name: &str) -> i32 { + if name.is_empty() { return 0; } + for (i, p) in doc.points.iter().enumerate() { + if p.name == name { return (i as i32) + 2; } + } + 0 +} + +/// resolves a segment's boundary-marker name to Triangle's marker integer, returned positive. +/// the caller negates for Triangle's segment-marker convention. +fn boundary_marker_index(doc: &FemmDoc, name: &str) -> i32 { + if name.is_empty() { return 0; } + for (i, b) in doc.boundaries.iter().enumerate() { + if b.name == name { return (i as i32) + 2; } + } + 0 +} + +#[cfg(test)] +mod tests { + use super::*; + use crate::{BoundaryProp, PointProp, BlockLabel}; + + #[test] + fn write_poly_header_lines_count_correctly() { + // four-corner square with one block label inside; verify the section headers and counts. + let mut d = FemmDoc::default(); + d.points.push(PointProp { name: "V=0".into(), ..PointProp::default() }); + d.boundaries.push(BoundaryProp { name: "outer".into(), ..BoundaryProp::default() }); + d.add_node(0.0, 0.0, 0.0); + d.add_node(1.0, 0.0, 0.0); + d.add_node(1.0, 1.0, 0.0); + d.add_node(0.0, 1.0, 0.0); + d.nodes[0].boundary_marker = "V=0".into(); + d.add_segment_with_marker(0, 1, "outer"); + d.add_segment_with_marker(1, 2, "outer"); + d.add_segment_with_marker(2, 3, "outer"); + d.add_segment_with_marker(3, 0, "outer"); + d.block_labels.push(BlockLabel { + x: 0.5, y: 0.5, + max_area: 0.0, + block_type: "Air".into(), + ..BlockLabel::default() + }); + + let poly = d.write_poly(); + let lines: Vec<&str> = poly.lines().collect(); + + assert!(lines[0].starts_with("4\t2\t0\t1"), "node header: {}", lines[0]); + assert!(lines[5].starts_with("4\t1"), "segment header: {}", lines[5]); + assert_eq!(lines[10], "0", "no holes section: {}", lines[10]); + assert_eq!(lines[11], "1", "one region: {}", lines[11]); + } + + #[test] + fn write_poly_marks_nodes_and_segments_by_property_index() { + // first point property listed -> marker 2; first boundary property listed -> -2 on segment. + let mut d = FemmDoc::default(); + d.points.push(PointProp { name: "V=0".into(), ..PointProp::default() }); + d.boundaries.push(BoundaryProp { name: "outer".into(), ..BoundaryProp::default() }); + d.add_node(0.0, 0.0, 0.0); + d.add_node(1.0, 0.0, 0.0); + d.nodes[0].boundary_marker = "V=0".into(); + d.add_segment_with_marker(0, 1, "outer"); + + let poly = d.write_poly(); + let lines: Vec<&str> = poly.lines().collect(); + + assert!(lines[1].ends_with("\t2"), "node 0 marker = 2: {}", lines[1]); + assert!(lines[2].ends_with("\t0"), "node 1 marker = 0: {}", lines[2]); + assert!(lines[4].ends_with("\t-2"), "segment 0 marker = -2: {}", lines[4]); + } + + #[test] + fn write_poly_holes_partitioned_from_regions() { + // one "" label counts as a hole; everything else falls into the regions section. + let mut d = FemmDoc::default(); + d.add_node(0.0, 0.0, 0.0); + d.add_node(1.0, 0.0, 0.0); + d.block_labels.push(BlockLabel { + x: 0.5, y: 0.5, + block_type: "".into(), + ..BlockLabel::default() + }); + d.block_labels.push(BlockLabel { + x: 0.2, y: 0.2, + block_type: "Air".into(), + ..BlockLabel::default() + }); + + let poly = d.write_poly(); + let lines: Vec<&str> = poly.lines().collect(); + + // node header (0), 2 node rows (1,2), segment header "0\t1" (3), + // hole count (4), 1 hole row (5), region count (6), 1 region row (7). + assert_eq!(lines[4], "1", "one hole: {}", lines[4]); + assert_eq!(lines[6], "1", "one region: {}", lines[6]); + } +} diff --git a/crates/femm-doc-heat/src/lib.rs b/crates/femm-doc-heat/src/lib.rs index cfbead4..019510b 100644 --- a/crates/femm-doc-heat/src/lib.rs +++ b/crates/femm-doc-heat/src/lib.rs @@ -6,6 +6,8 @@ pub mod props; pub mod parser; pub mod writer; pub mod edit; +pub mod poly; +pub mod mesh; use num_complex::Complex64; diff --git a/crates/femm-doc-heat/src/mesh.rs b/crates/femm-doc-heat/src/mesh.rs new file mode 100644 index 0000000..c85669a --- /dev/null +++ b/crates/femm-doc-heat/src/mesh.rs @@ -0,0 +1,161 @@ +//! readers for Triangle's .node and .ele output files. + +use std::path::Path; +use thiserror::Error; + +/// 2D triangular mesh consumed by the canvas overlay and by the post-processor. +#[derive(Debug, Default, Clone)] +pub struct Mesh { + pub nodes: Vec, + pub elements: Vec, +} + +/// position and optional boundary marker for a single mesh node. +#[derive(Debug, Clone, Copy)] +pub struct MeshNode { + pub x: f64, + pub y: f64, + pub marker: i32, +} + +/// triangle element with three zero-based node indices and an optional region attribute. +#[derive(Debug, Clone, Copy)] +pub struct MeshElement { + pub v0: u32, + pub v1: u32, + pub v2: u32, + pub attribute: i32, +} + +#[derive(Debug, Error)] +pub enum MeshLoadError { + #[error("io error: {0}")] + Io(#[from] std::io::Error), + #[error("malformed header in {file}: {reason}")] + BadHeader { file: &'static str, reason: String }, + #[error("malformed row in {file}: {reason}")] + BadRow { file: &'static str, reason: String }, +} + +impl Mesh { + /// loads .node and .ele next to the given path stem. + pub fn load(stem: impl AsRef) -> Result { + let stem = stem.as_ref(); + let node_path = with_extension(stem, "node"); + let ele_path = with_extension(stem, "ele"); + let nodes = parse_node(&std::fs::read_to_string(&node_path)?)?; + let elements = parse_ele(&std::fs::read_to_string(&ele_path)?)?; + Ok(Mesh { nodes, elements }) + } +} + +/// parses a .node file: header `count 2 0 nbm`, then count rows of `idx x y [marker]`. +pub fn parse_node(src: &str) -> Result, MeshLoadError> { + let mut lines = src.lines().filter(|l| !l.trim_start().starts_with('#') && !l.trim().is_empty()); + let header = lines.next().ok_or(MeshLoadError::BadHeader { + file: ".node", + reason: "empty file".into(), + })?; + let mut hparts = header.split_whitespace(); + let count: usize = hparts.next() + .and_then(|s| s.parse().ok()) + .ok_or(MeshLoadError::BadHeader { file: ".node", reason: format!("expected count, got {header:?}") })?; + let _dim = hparts.next(); + let _attrs = hparts.next(); + let has_marker = hparts.next().map(|s| s.trim() == "1").unwrap_or(false); + + let mut nodes = Vec::with_capacity(count); + for line in lines.take(count) { + let mut parts = line.split_whitespace(); + let _idx: u32 = parts.next().and_then(|s| s.parse().ok()) + .ok_or(MeshLoadError::BadRow { file: ".node", reason: line.into() })?; + let x: f64 = parts.next().and_then(|s| s.parse().ok()) + .ok_or(MeshLoadError::BadRow { file: ".node", reason: line.into() })?; + let y: f64 = parts.next().and_then(|s| s.parse().ok()) + .ok_or(MeshLoadError::BadRow { file: ".node", reason: line.into() })?; + let marker = if has_marker { + parts.next().and_then(|s| s.parse().ok()).unwrap_or(0) + } else { + 0 + }; + nodes.push(MeshNode { x, y, marker }); + } + Ok(nodes) +} + +/// parses a .ele file: header `count 3 natt`, then count rows of `idx v0 v1 v2 [attr0 ...]`. +pub fn parse_ele(src: &str) -> Result, MeshLoadError> { + let mut lines = src.lines().filter(|l| !l.trim_start().starts_with('#') && !l.trim().is_empty()); + let header = lines.next().ok_or(MeshLoadError::BadHeader { + file: ".ele", + reason: "empty file".into(), + })?; + let mut hparts = header.split_whitespace(); + let count: usize = hparts.next() + .and_then(|s| s.parse().ok()) + .ok_or(MeshLoadError::BadHeader { file: ".ele", reason: format!("expected count, got {header:?}") })?; + let _nodes_per = hparts.next(); + let natt: usize = hparts.next().and_then(|s| s.parse().ok()).unwrap_or(0); + + let mut elements = Vec::with_capacity(count); + for line in lines.take(count) { + let mut parts = line.split_whitespace(); + let _idx: u32 = parts.next().and_then(|s| s.parse().ok()) + .ok_or(MeshLoadError::BadRow { file: ".ele", reason: line.into() })?; + let v0: u32 = parts.next().and_then(|s| s.parse().ok()) + .ok_or(MeshLoadError::BadRow { file: ".ele", reason: line.into() })?; + let v1: u32 = parts.next().and_then(|s| s.parse().ok()) + .ok_or(MeshLoadError::BadRow { file: ".ele", reason: line.into() })?; + let v2: u32 = parts.next().and_then(|s| s.parse().ok()) + .ok_or(MeshLoadError::BadRow { file: ".ele", reason: line.into() })?; + let attribute = if natt > 0 { + parts.next().and_then(|s| s.parse::().ok()).map(|f| f as i32).unwrap_or(0) + } else { + 0 + }; + elements.push(MeshElement { v0, v1, v2, attribute }); + } + Ok(elements) +} + +fn with_extension(stem: &Path, ext: &str) -> std::path::PathBuf { + let mut s = stem.to_path_buf(); + s.set_extension(ext); + s +} + +#[cfg(test)] +mod tests { + use super::*; + + #[test] + fn parses_node_with_markers() { + // 3 nodes, 2D, no attributes, with markers; comment line at the start is ignored. + let src = "# triangle .node\n3 2 0 1\n0 0.0 0.0 2\n1 1.0 0.0 0\n2 0.0 1.0 0\n"; + let nodes = parse_node(src).unwrap(); + assert_eq!(nodes.len(), 3); + assert_eq!(nodes[0].marker, 2); + assert!((nodes[2].y - 1.0).abs() < 1e-12); + } + + #[test] + fn parses_ele_with_region_attribute() { + // 2 triangles, 3 nodes per element, 1 region attribute column. + let src = "2 3 1\n0 0 1 2 1.0\n1 1 2 3 2.0\n"; + let els = parse_ele(src).unwrap(); + assert_eq!(els.len(), 2); + assert_eq!(els[0].v0, 0); + assert_eq!(els[0].v1, 1); + assert_eq!(els[0].v2, 2); + assert_eq!(els[0].attribute, 1); + assert_eq!(els[1].attribute, 2); + } + + #[test] + fn parses_node_without_marker_column() { + let src = "2 2 0 0\n0 0.0 0.0\n1 1.0 0.0\n"; + let nodes = parse_node(src).unwrap(); + assert_eq!(nodes.len(), 2); + assert_eq!(nodes[0].marker, 0); + } +} diff --git a/crates/femm-doc-heat/src/poly.rs b/crates/femm-doc-heat/src/poly.rs new file mode 100644 index 0000000..7fe667f --- /dev/null +++ b/crates/femm-doc-heat/src/poly.rs @@ -0,0 +1,166 @@ +//! .poly emitter consumed by the Triangle mesher. + +use crate::FemmDoc; +use std::fmt::Write; + +impl FemmDoc { + /// renders the doc geometry to the Triangle .poly text format, returning the file contents. + pub fn write_poly(&self) -> String { + let mut out = String::new(); + let bbox_diag = self.bbox_diagonal(); + let default_area = if bbox_diag > 0.0 { bbox_diag } else { -1.0 }; + + writeln!(out, "{}\t2\t0\t1", self.nodes.len()).unwrap(); + for (i, n) in self.nodes.iter().enumerate() { + let marker = point_marker_index(self, &n.boundary_marker); + writeln!(out, "{i}\t{:.17e}\t{:.17e}\t{marker}", n.x, n.y).unwrap(); + } + + writeln!(out, "{}\t1", self.segments.len()).unwrap(); + for (i, s) in self.segments.iter().enumerate() { + let marker = -boundary_marker_index(self, &s.boundary_marker); + writeln!(out, "{i}\t{}\t{}\t{marker}", s.n0, s.n1).unwrap(); + } + + let holes: Vec<&_> = self.block_labels.iter().filter(|b| b.block_type == "").collect(); + writeln!(out, "{}", holes.len()).unwrap(); + for (i, h) in holes.iter().enumerate() { + writeln!(out, "{i}\t{:.17e}\t{:.17e}", h.x, h.y).unwrap(); + } + + let regions: Vec<&_> = self.block_labels.iter().filter(|b| b.block_type != "").collect(); + writeln!(out, "{}", regions.len()).unwrap(); + for (i, r) in regions.iter().enumerate() { + let attr = (i as i32) + 1; + let max_area = if r.max_area > 0.0 && r.max_area < default_area { + r.max_area + } else { + default_area + }; + writeln!(out, "{i}\t{:.17e}\t{:.17e}\t{attr}\t{:.17e}", r.x, r.y, max_area).unwrap(); + } + + out + } + + /// writes the .poly text to disk at the given path. + pub fn save_poly(&self, path: impl AsRef) -> std::io::Result<()> { + std::fs::write(path, self.write_poly()) + } + + fn bbox_diagonal(&self) -> f64 { + if self.nodes.len() < 2 { return 0.0; } + let mut xmin = f64::INFINITY; + let mut xmax = f64::NEG_INFINITY; + let mut ymin = f64::INFINITY; + let mut ymax = f64::NEG_INFINITY; + for n in &self.nodes { + xmin = xmin.min(n.x); xmax = xmax.max(n.x); + ymin = ymin.min(n.y); ymax = ymax.max(n.y); + } + let dx = xmax - xmin; + let dy = ymax - ymin; + (dx * dx + dy * dy).sqrt() + } +} + +/// resolves a node's boundary-marker name to Triangle's marker integer (0 for unmarked). +fn point_marker_index(doc: &FemmDoc, name: &str) -> i32 { + if name.is_empty() { return 0; } + for (i, p) in doc.points.iter().enumerate() { + if p.name == name { return (i as i32) + 2; } + } + 0 +} + +/// resolves a segment's boundary-marker name to Triangle's marker integer, returned positive. +/// the caller negates for Triangle's segment-marker convention. +fn boundary_marker_index(doc: &FemmDoc, name: &str) -> i32 { + if name.is_empty() { return 0; } + for (i, b) in doc.boundaries.iter().enumerate() { + if b.name == name { return (i as i32) + 2; } + } + 0 +} + +#[cfg(test)] +mod tests { + use super::*; + use crate::{BoundaryProp, PointProp, BlockLabel}; + + #[test] + fn write_poly_header_lines_count_correctly() { + // four-corner square with one block label inside; verify the section headers and counts. + let mut d = FemmDoc::default(); + d.points.push(PointProp { name: "T=0".into(), ..PointProp::default() }); + d.boundaries.push(BoundaryProp { name: "outer".into(), ..BoundaryProp::default() }); + d.add_node(0.0, 0.0, 0.0); + d.add_node(1.0, 0.0, 0.0); + d.add_node(1.0, 1.0, 0.0); + d.add_node(0.0, 1.0, 0.0); + d.nodes[0].boundary_marker = "T=0".into(); + d.add_segment_with_marker(0, 1, "outer"); + d.add_segment_with_marker(1, 2, "outer"); + d.add_segment_with_marker(2, 3, "outer"); + d.add_segment_with_marker(3, 0, "outer"); + d.block_labels.push(BlockLabel { + x: 0.5, y: 0.5, + max_area: 0.0, + block_type: "Air".into(), + ..BlockLabel::default() + }); + + let poly = d.write_poly(); + let lines: Vec<&str> = poly.lines().collect(); + + assert!(lines[0].starts_with("4\t2\t0\t1"), "node header: {}", lines[0]); + assert!(lines[5].starts_with("4\t1"), "segment header: {}", lines[5]); + assert_eq!(lines[10], "0", "no holes section: {}", lines[10]); + assert_eq!(lines[11], "1", "one region: {}", lines[11]); + } + + #[test] + fn write_poly_marks_nodes_and_segments_by_property_index() { + // first point property listed -> marker 2; first boundary property listed -> -2 on segment. + let mut d = FemmDoc::default(); + d.points.push(PointProp { name: "T=0".into(), ..PointProp::default() }); + d.boundaries.push(BoundaryProp { name: "outer".into(), ..BoundaryProp::default() }); + d.add_node(0.0, 0.0, 0.0); + d.add_node(1.0, 0.0, 0.0); + d.nodes[0].boundary_marker = "T=0".into(); + d.add_segment_with_marker(0, 1, "outer"); + + let poly = d.write_poly(); + let lines: Vec<&str> = poly.lines().collect(); + + assert!(lines[1].ends_with("\t2"), "node 0 marker = 2: {}", lines[1]); + assert!(lines[2].ends_with("\t0"), "node 1 marker = 0: {}", lines[2]); + assert!(lines[4].ends_with("\t-2"), "segment 0 marker = -2: {}", lines[4]); + } + + #[test] + fn write_poly_holes_partitioned_from_regions() { + // one "" label counts as a hole; everything else is a region. + let mut d = FemmDoc::default(); + d.add_node(0.0, 0.0, 0.0); + d.add_node(1.0, 0.0, 0.0); + d.block_labels.push(BlockLabel { + x: 0.5, y: 0.5, + block_type: "".into(), + ..BlockLabel::default() + }); + d.block_labels.push(BlockLabel { + x: 0.2, y: 0.2, + block_type: "Air".into(), + ..BlockLabel::default() + }); + + let poly = d.write_poly(); + let lines: Vec<&str> = poly.lines().collect(); + + // node header (0), 2 node rows (1,2), segment header "0\t1" (3), + // hole count (4), 1 hole row (5), region count (6), 1 region row (7). + assert_eq!(lines[4], "1", "one hole: {}", lines[4]); + assert_eq!(lines[6], "1", "one region: {}", lines[6]); + } +} diff --git a/crates/femm-doc-mag/src/ans.rs b/crates/femm-doc-mag/src/ans.rs new file mode 100644 index 0000000..6ea3c0d --- /dev/null +++ b/crates/femm-doc-mag/src/ans.rs @@ -0,0 +1,422 @@ +//! .ans solution-file parser: reuses the .fem header grammar and adds the post-[Solution] mesh body. + +use crate::{ + ACSolver, BoundaryProp, CircuitProp, Complex, Coords, FemmDoc, LengthUnit, MaterialProp, + PointProp, ProblemType, +}; +use std::path::Path; +use thiserror::Error; + +#[derive(Debug, Error)] +pub enum AnsParseError { + #[error("io error: {0}")] + Io(#[from] std::io::Error), + #[error("header parse error: {0}")] + Header(#[from] crate::parser::ParseError), + #[error("missing [Solution] marker")] + MissingSolution, + #[error("unexpected end of file while reading {context}")] + UnexpectedEof { context: &'static str }, + #[error("malformed row in {section}: {reason}")] + MalformedRow { section: &'static str, reason: String }, +} + +/// solver header attributes shared by every magnetostatic .ans file. +#[derive(Debug, Clone, Default)] +pub struct MagHeader { + pub format: f64, + pub frequency: f64, + pub precision: f64, + pub min_angle: f64, + pub smart_mesh: bool, + pub depth: f64, + pub length_units: LengthUnit, + pub ac_solver: ACSolver, + pub problem_type: ProblemType, + pub coords: Coords, + pub comment: String, + pub ext_ro: f64, + pub ext_ri: f64, + pub ext_zo: f64, +} + +/// one finite-element mesh node with the solved magnetic vector potential. +#[derive(Debug, Clone, Copy, Default)] +pub struct SolMeshNode { + pub x: f64, + pub y: f64, + pub a: Complex, +} + +/// one mesh triangle with three node indices, block label, and post-load flux/permeability slots. +#[derive(Debug, Clone, Copy, Default)] +pub struct SolMeshElement { + pub n: [u32; 3], + pub lab: i32, + pub b1: Complex, + pub b2: Complex, + pub mu1: Complex, + pub mu2: Complex, +} + +/// per-block-label circuit assignment value sourced from the trailing .ans stanza. +#[derive(Debug, Clone, Copy, Default)] +pub struct BlockCircuitValue { + pub case: i32, + pub value: Complex, +} + +/// parsed magnetostatic solution file: header, property tables, mesh body, and circuit body. +#[derive(Debug, Clone, Default)] +pub struct MagSolution { + pub header: MagHeader, + pub points: Vec, + pub boundaries: Vec, + pub materials: Vec, + pub circuits: Vec, + pub mesh_nodes: Vec, + pub mesh_elements: Vec, + pub block_circuits: Vec, +} + +impl MagSolution { + /// parses a .ans text buffer into a magnetostatic solution. + pub fn parse(src: &str) -> Result { + let (head, body) = split_at_solution(src).ok_or(AnsParseError::MissingSolution)?; + let doc = FemmDoc::parse(head)?; + let header = MagHeader { + format: doc.format, + frequency: doc.frequency.abs(), + precision: doc.precision, + min_angle: doc.min_angle, + smart_mesh: doc.smart_mesh, + depth: doc.depth, + length_units: doc.length_units, + ac_solver: doc.ac_solver, + problem_type: doc.problem_type, + coords: doc.coords, + comment: doc.comment, + ext_ro: doc.ext_ro, + ext_ri: doc.ext_ri, + ext_zo: doc.ext_zo, + }; + let ac = header.frequency != 0.0; + let mut lines = body.lines(); + + let nnodes = read_count(&mut lines, "mesh nodes count")?; + let mut mesh_nodes = Vec::with_capacity(nnodes); + for _ in 0..nnodes { + let row = lines.next().ok_or(AnsParseError::UnexpectedEof { context: "mesh nodes" })?; + mesh_nodes.push(parse_node_row(row, ac)?); + } + + let nels = read_count(&mut lines, "mesh elements count")?; + let mut mesh_elements = Vec::with_capacity(nels); + for _ in 0..nels { + let row = lines.next().ok_or(AnsParseError::UnexpectedEof { context: "mesh elements" })?; + mesh_elements.push(parse_element_row(row)?); + } + + let ncirc = read_count(&mut lines, "block-circuit rows count")?; + let mut block_circuits = Vec::with_capacity(ncirc); + for _ in 0..ncirc { + let row = lines.next().ok_or(AnsParseError::UnexpectedEof { context: "block-circuit rows" })?; + block_circuits.push(parse_block_circuit_row(row, ac)?); + } + + Ok(MagSolution { + header, + points: doc.points, + boundaries: doc.boundaries, + materials: doc.materials, + circuits: doc.circuits, + mesh_nodes, + mesh_elements, + block_circuits, + }) + } + + /// loads a .ans file from disk. + pub fn open(path: impl AsRef) -> Result { + let text = std::fs::read_to_string(path)?; + Self::parse(&text) + } +} + +/// splits the raw text on the `[Solution]` marker into a header chunk and a body chunk. +fn split_at_solution(src: &str) -> Option<(&str, &str)> { + let lower_target = "[solution]"; + let mut offset = 0usize; + for line in src.lines() { + let line_len = line.len(); + let trimmed = line.trim_start(); + if trimmed.to_ascii_lowercase().starts_with(lower_target) { + let head_end = offset; + let body_start = offset + line_len; + let body_start = skip_newline(src, body_start); + return Some((&src[..head_end], &src[body_start..])); + } + offset += line_len; + offset = skip_newline(src, offset); + } + None +} + +fn skip_newline(src: &str, mut i: usize) -> usize { + let b = src.as_bytes(); + if i < b.len() && b[i] == b'\r' { i += 1; } + if i < b.len() && b[i] == b'\n' { i += 1; } + i +} + +/// reads the next non-blank line as a decimal integer count. +fn read_count<'a, I: Iterator>( + lines: &mut I, + ctx: &'static str, +) -> Result { + for line in lines.by_ref() { + let t = line.trim(); + if t.is_empty() { continue; } + return t.parse::().map_err(|e| AnsParseError::MalformedRow { + section: ctx, + reason: format!("{e}: {t:?}"), + }); + } + Err(AnsParseError::UnexpectedEof { context: ctx }) +} + +/// parses one mesh-node row, picking the column layout from the AC/DC flag. +fn parse_node_row(line: &str, ac: bool) -> Result { + let (x, rest) = take_f64(line); + let (y, rest) = take_f64(rest); + let (a_re, rest) = take_f64(rest); + let a_im = if ac { + let (v, _) = take_f64(rest); + v + } else { + 0.0 + }; + Ok(SolMeshNode { x, y, a: Complex::new(a_re, a_im) }) +} + +/// parses one mesh-element row into three node indices plus a block-label index. +fn parse_element_row(line: &str) -> Result { + let (p0, rest) = take_i32(line); + let (p1, rest) = take_i32(rest); + let (p2, rest) = take_i32(rest); + let (lab, _) = take_i32(rest); + Ok(SolMeshElement { + n: [p0 as u32, p1 as u32, p2 as u32], + lab, + ..Default::default() + }) +} + +/// parses one block-circuit row of `case value [imag]`. +fn parse_block_circuit_row(line: &str, ac: bool) -> Result { + let (case, rest) = take_i32(line); + let (re, rest) = take_f64(rest); + let im = if ac { + let (v, _) = take_f64(rest); + v + } else { + 0.0 + }; + Ok(BlockCircuitValue { case, value: Complex::new(re, im) }) +} + +/// peels one f64 off the head of a whitespace/tab/comma-separated row. +fn take_f64(s: &str) -> (f64, &str) { + let s = s.trim_start_matches(|c: char| c.is_whitespace() || c == ','); + let end = s.find(|c: char| c.is_whitespace() || c == ',').unwrap_or(s.len()); + let (head, tail) = s.split_at(end); + (head.parse::().unwrap_or(0.0), tail) +} + +/// peels one i32 off the head of a whitespace/tab/comma-separated row. +fn take_i32(s: &str) -> (i32, &str) { + let s = s.trim_start_matches(|c: char| c.is_whitespace() || c == ','); + let end = s.find(|c: char| c.is_whitespace() || c == ',').unwrap_or(s.len()); + let (head, tail) = s.split_at(end); + (head.parse::().unwrap_or(0), tail) +} + +#[cfg(test)] +mod tests { + use super::*; + + const DC_FIXTURE: &str = "\ +[Format] = 4.0 +[Frequency] = 0 +[Precision] = 1e-08 +[MinAngle] = 30 +[DoSmartMesh] = 1 +[Depth] = 1 +[LengthUnits] = millimeters +[ProblemType] = planar +[Coordinates] = cartesian +[ACSolver] = 0 +[PrevType] = 0 +[PrevSoln] = \"\" +[Comment] = \"dc fixture\" +[PointProps] = 1 + + = \"A=0\" + = 0 + = 0 + = 0 + = 0 + +[BdryProps] = 1 + + = \"outer\" + = 0 + = 0 + = 0 + = 0 + = 0 + = 0 + = 0 + = 0 + = 0 + = 0 + = 0 + = 0 + = 0 + +[BlockProps] = 1 + + = \"Air\" + = 1 + = 1 + = 0 + = 0 + = 0 + = 0 + = 0 + = 0 + = 0 + = 0 + = 0 + = 0 + = 1 + = 0 + = 0 + = 0 + +[CircuitProps] = 1 + + = \"coil\" + = 1 + = 0 + = 1 + +[NumPoints] = 0 +[NumSegments] = 0 +[NumArcSegments] = 0 +[NumHoles] = 0 +[NumBlockLabels] = 1 +0.5\t0.5\t1\t-1\t0\t0\t0\t1\t0\t\"\" +[Solution] +4 +0.0\t0.0\t0.0 +1.0\t0.0\t0.125 +0.0\t1.0\t0.25 +1.0\t1.0\t0.375 +2 +0\t1\t2\t0 +1\t3\t2\t0 +1 +1\t0.5 +"; + + const AC_FIXTURE: &str = "\ +[Format] = 4.0 +[Frequency] = 60 +[Precision] = 1e-08 +[MinAngle] = 30 +[DoSmartMesh] = 1 +[Depth] = 1 +[LengthUnits] = millimeters +[ProblemType] = planar +[Coordinates] = cartesian +[ACSolver] = 0 +[PrevType] = 0 +[PrevSoln] = \"\" +[Comment] = \"ac fixture\" +[PointProps] = 0 +[BdryProps] = 0 +[BlockProps] = 0 +[CircuitProps] = 0 +[NumPoints] = 0 +[NumSegments] = 0 +[NumArcSegments] = 0 +[NumHoles] = 0 +[NumBlockLabels] = 0 +[Solution] +3 +0.0\t0.0\t0.0\t0.0 +1.0\t0.0\t0.125\t0.0625 +0.0\t1.0\t0.25\t0.125 +1 +0\t1\t2\t0 +0 +"; + + #[test] + fn parses_dc_header_and_counts() { + let sol = MagSolution::parse(DC_FIXTURE).unwrap(); + assert_eq!(sol.header.frequency, 0.0); + assert_eq!(sol.header.problem_type, ProblemType::Planar); + assert_eq!(sol.header.length_units, LengthUnit::Millimeters); + assert!((sol.header.depth - 1.0).abs() < 1e-12); + assert_eq!(sol.points.len(), 1); + assert_eq!(sol.boundaries.len(), 1); + assert_eq!(sol.materials.len(), 1); + assert_eq!(sol.circuits.len(), 1); + assert_eq!(sol.mesh_nodes.len(), 4); + assert_eq!(sol.mesh_elements.len(), 2); + assert_eq!(sol.block_circuits.len(), 1); + } + + #[test] + fn parses_dc_node_a_values_real_only() { + let sol = MagSolution::parse(DC_FIXTURE).unwrap(); + assert!((sol.mesh_nodes[2].a.re - 0.25).abs() < 1e-12); + assert!(sol.mesh_nodes[2].a.im.abs() < 1e-12); + assert!((sol.mesh_nodes[3].x - 1.0).abs() < 1e-12); + assert!((sol.mesh_nodes[3].y - 1.0).abs() < 1e-12); + } + + #[test] + fn parses_element_indices_and_label_with_zero_b_field() { + let sol = MagSolution::parse(DC_FIXTURE).unwrap(); + let e0 = &sol.mesh_elements[0]; + assert_eq!(e0.n, [0, 1, 2]); + assert_eq!(e0.lab, 0); + assert_eq!(e0.b1, Complex::new(0.0, 0.0)); + assert_eq!(e0.b2, Complex::new(0.0, 0.0)); + assert_eq!(e0.mu1, Complex::new(0.0, 0.0)); + assert_eq!(e0.mu2, Complex::new(0.0, 0.0)); + } + + #[test] + fn parses_ac_node_complex_a_values() { + let sol = MagSolution::parse(AC_FIXTURE).unwrap(); + assert_eq!(sol.header.frequency, 60.0); + assert_eq!(sol.mesh_nodes.len(), 3); + assert!((sol.mesh_nodes[1].a.re - 0.125).abs() < 1e-12); + assert!((sol.mesh_nodes[1].a.im - 0.0625).abs() < 1e-12); + assert_eq!(sol.mesh_elements.len(), 1); + assert_eq!(sol.block_circuits.len(), 0); + } + + #[test] + fn rejects_input_missing_solution_marker() { + let bad = "[Format] = 4.0\n[Frequency] = 0\n"; + match MagSolution::parse(bad) { + Err(AnsParseError::MissingSolution) => {} + other => panic!("expected MissingSolution, got {other:?}"), + } + } +} diff --git a/crates/femm-doc-mag/src/lib.rs b/crates/femm-doc-mag/src/lib.rs index e2e54c5..b67707e 100644 --- a/crates/femm-doc-mag/src/lib.rs +++ b/crates/femm-doc-mag/src/lib.rs @@ -8,6 +8,7 @@ pub mod writer; pub mod edit; pub mod poly; pub mod mesh; +pub mod ans; use num_complex::Complex64; diff --git a/crates/femm-doc-mag/src/mesh.rs b/crates/femm-doc-mag/src/mesh.rs index 752c5a9..c85669a 100644 --- a/crates/femm-doc-mag/src/mesh.rs +++ b/crates/femm-doc-mag/src/mesh.rs @@ -38,11 +38,11 @@ pub enum MeshLoadError { } impl Mesh { - /// loads .1.node and .1.ele next to the given path stem (Triangle's default output). + /// loads .node and .ele next to the given path stem. pub fn load(stem: impl AsRef) -> Result { let stem = stem.as_ref(); - let node_path = with_extension(stem, "1.node"); - let ele_path = with_extension(stem, "1.ele"); + let node_path = with_extension(stem, "node"); + let ele_path = with_extension(stem, "ele"); let nodes = parse_node(&std::fs::read_to_string(&node_path)?)?; let elements = parse_ele(&std::fs::read_to_string(&ele_path)?)?; Ok(Mesh { nodes, elements }) diff --git a/scripts/macos/build.sh b/scripts/macos/build.sh index da5cde7..8314118 100755 --- a/scripts/macos/build.sh +++ b/scripts/macos/build.sh @@ -20,6 +20,15 @@ ICONSET="$BUILD/AppIcon.iconset" ICNS="$RESOURCES/AppIcon.icns" SVG="$ROOT/assets/femm.svg" +echo "Building Triangle..." +bash "$ROOT/scripts/macos/build_triangle.sh" + +TRI="$ROOT/build/triangle/triangle" +if [ ! -f "$TRI" ]; then + echo "ERROR: triangle binary missing at $TRI" >&2 + exit 1 +fi + echo "Building Rust workspace (release)..." cargo build --release -p femm-app @@ -32,6 +41,7 @@ fi rm -rf "$APP" mkdir -p "$MACOS" "$RESOURCES" cp "$BIN" "$MACOS/femm" +cp "$TRI" "$MACOS/triangle" if [ -f "$SVG" ]; then if ! command -v rsvg-convert >/dev/null; then