diff --git a/.gitignore b/.gitignore index f5c85c3..d67dbc0 100644 --- a/.gitignore +++ b/.gitignore @@ -2,3 +2,5 @@ build/ffi/ .DS_Store *.o *.a +target/ +Cargo.lock diff --git a/Cargo.toml b/Cargo.toml new file mode 100644 index 0000000..252bee4 --- /dev/null +++ b/Cargo.toml @@ -0,0 +1,17 @@ +[workspace] +resolver = "2" +members = ["crates/*"] + +[workspace.package] +edition = "2024" +rust-version = "1.85" +publish = false + +[workspace.dependencies] +femm-sys = { path = "crates/femm-sys" } +femm-doc = { path = "crates/femm-doc" } + +[profile.release] +opt-level = 3 +lto = "thin" +codegen-units = 1 diff --git a/ROADMAP.md b/ROADMAP.md new file mode 100644 index 0000000..704d2d4 --- /dev/null +++ b/ROADMAP.md @@ -0,0 +1,210 @@ +# FEMM port roadmap + +## status snapshot + +The four C++ solvers (`fkn`, `belasolv`, `csolv`, `hsolv`) are linked into Rust via per-engine static archives, wrapped behind C ABI under `ffi/`, surfaced through `crates/femm-sys`. Document model is mirrored in Rust across four discipline crates (`femm-doc` for magnetostatic, `femm-doc-elec`, `femm-doc-curr`, `femm-doc-heat`); parsers and writers round-trip the on-disk text formats. An iced 0.14 canvas in `crates/femm-app` opens .fem, pans, zooms, and supports Select / AddNode / AddBlockLabel as click-tools against `femm-doc::edit`. Mag editing primitives exist (`add_node`, `add_segment`, `add_arc_segment`, `add_block_label`, `delete_selected_*`, closest-entity queries) but PSLG enforcement (intersection splitting, duplicate collapse, on-line break) is not yet implemented. There is no preprocessor for the three non-mag disciplines wired into the GUI, no post-processor at all, no property editing, no problem-definition dialog, no mesh launch, and no .ans loader. + +## pre-processor + +The original draws a planar straight-line graph (PSLG) of nodes, segments, arc segments, and block labels. Each discipline has its own doc class (`CFemmeDoc`, `CbeladrawDoc`, `CcdrawDoc`, `ChdrawDoc`) sharing identical geometry storage but differing in `MaterialProp`, `BoundaryProp`, `PointProp`, `Circuit` fields. The View class for each discipline (`CFemmeView`, `CbeladrawView`, `CcdrawView`, `ChdrawView`) is structurally identical — the discipline-specific behavior is confined to which property dialogs open from `OnOpenSelected`. + +### edit modes and mouse gestures + +`EditAction` in the View is an `int`: `0`=Node, `1`=Segment, `2`=BlockLabel, `3`=ArcSegment, `4`=Group. Each mode rebinds the meaning of left-click, right-click, right-double-click, and TAB. Source: `FemmeView.cpp` (mag), `beladrawView.cpp` (elec), `cdrawView.cpp` (current), `hdrawView.cpp` (heat). + +- node mode (`EditAction==0`) — left-click places a node at the snapped mouse position; right-click toggles selection of the closest node; right-double-click pops an info message-box for that node; TAB opens the CEnterPt dialog for typed coordinates. Tier S. +- segment mode (`EditAction==1`) — left-click picks the closest node as first endpoint; a second left-click picks the second endpoint and emits an `AddSegment`; right-click toggles selection of the closest segment; right-double-click reports length, BC, mesh size; ESC clears the half-built segment. Tier S. +- arc segment mode (`EditAction==3`) — same two-click flow as segment, but between the two clicks a `CArcDlg` opens for arc angle, max segment length, and boundary marker. Tier M. +- block-label mode (`EditAction==2`) — left-click places a label; right-click toggles selection; right-double-click reports material, circuit, mesh size, group, magnetization direction. Tier S. +- group mode (`EditAction==4`) — right-click toggles every entity sharing the closest entity's group number; TAB pops a `CGroupNumber` dialog to type a group number and toggle that whole group. Tier S. +- box select (`SelectWndFlag`) — drag a rectangle; on release selects all entities of the current mode (or all four kinds when in group mode) wholly inside. XOR-pixel rubber-band drawn in `OnMouseMove`. Tier M. +- circle select (`SelectCircFlag`, mag only via `OnFDSelectCirc`) — drag from center to edge; selects entities inside the resulting disc. Tier M. +- zoom-window (`ZoomWndFlag`) — drag a rectangle that becomes the new view extents. Tier S. +- pan via arrow keys; PgUp/PgDn zoom by step; Home resets zoom; `Ctrl+Z` undoes; DEL cuts selection; SPACE opens the property dialog for current selection; F3 halves mesh density on selected blocks; F4 doubles. Tier S. +- create-radius (`CreateRadiusFlag`) — next click chooses a node; a prompt asks for radius; if `CanCreateRadius` succeeds the corner is filleted in. Tier M. + +The XOR-pixel rubber-band approach in the original is a Win32 trick (read pixels, XOR-write, remember to restore). In iced this becomes a transient overlay shape on the canvas — simpler and tier-S to implement. + +### dialogs + +Every dialog below corresponds to one `Xxx.h`/`Xxx.cpp` pair under `femm/`. Discipline-prefixed variants share a layout but differ in fields. The .rc resource definitions are gone, so layout is redesigned; field set, validation, and ok-button effect come from the .h and the OnInitDialog/OnOK in the .cpp. + +operation dialogs opened by SPACE / Edit > Properties + +- `COpNodeDlg` / `bd_OpNodeDlg` / `cd_OpNodeDlg` / `hd_OpNodeDlg` — `OpNodeDlg.h` and prefixed variants. Fields: nodal-property name (dropdown over `nodeproplist`), in-group. Tier S. +- `COpSegDlg` / `bd_OpSegDlg` / `cd_OpSegDlg` / `hd_OpSegDlg` — `OpSegDlg.h`. Fields: boundary-property name, line mesh size (with auto-mesh toggle), hide-in-postprocessor, in-group. Tier S. +- `COpArcSegDlg` / `bd_OpArcSegDlg` / `cd_OpArcSegDlg` / `hd_OpArcSegDlg` — `OpArcSegDlg.h`. Fields: boundary marker, max segment length (degrees), hide, in-group. Tier S. +- `COpBlkDlg` / `bd_OpBlkDlg` / `cd_OpBlkDlg` / `hd_OpBlkDlg` — `OpBlkDlg.h`. Fields: block material, circuit (mag only), is-external (axisymmetric only), is-default, side length / area limit with auto-mesh toggle, mag direction (and a Lua "magdirfctn" textbox for functional MagDir), turns, in-group. Tier M. +- `COpGrp` — `OpGrp.h`. Bulk-edit the group number across a heterogeneous selection. Tier S. + +property-list dialogs opened from menu Properties > Materials / Boundaries / Points / Circuits + +- `CPtProp` — `PtProp.h`. List with Add / Modify / Delete. Generic shell parameterized by which underlying property list it is editing. Tier S. + +property-editor dialogs reached from the Add/Modify buttons of `CPtProp` (one per discipline per property kind) + +- `CMatDlg` — `MatDlg.h`. Magnetostatic material: name, mu_x, mu_y, H_c, J_re+J_im, conductivity, lam direction, lam thickness, lam fill, hysteresis lag angles `Theta_h*`, BH-curve edit button (opens `CBHData`), wire diameter, num strands, and a nonlinearity-mode combo. Tier L. +- `bdCMatDlg` — `bd_MatDlg.h`. Electrostatic material: name, eps_x, eps_y, volume charge density qv. Tier S. +- `cdCMatDlg` — `cd_MatDlg.h`. Current-flow material: name, conductivity ox/oy (real+imag), and three-curve nonlinearity arrays. Tier M. +- `hdCMatDlg` — `hd_MatDlg.h`. Heat material: name, Kx, Ky, Kt (volumetric heat capacity), qv, and a nonlinear KH-curve combo with edit button. Tier M. +- `CBdryDlg` — `BdryDlg.h`. Magnetostatic boundary: BdryFormat combo (Prescribed-A, small-skin-depth, mixed, strategic dual, periodic, antiperiodic, periodic-air-gap, antiperiodic-air-gap), A0/A1/A2 + phi, c0/c1, mu, sigma, inner/outer angle. Tier M. +- `bdCBdryDlg` — `bd_BdryDlg.h`. Electrostatic boundary: format combo, v (potential), qs (surface charge), c0/c1. Tier S. +- `cdCBdryDlg` — `cd_BdryDlg.h`. Current-flow boundary: format combo, vs (complex voltage), qs (complex surface current), c0/c1. Tier S. +- `hdCBdryDlg` — `hd_BdryDlg.h`. Heat boundary: format combo, Tset (fixed T), q (heat flux), beta (radiation), htc (convection), To1/To2 (ambient temps). Tier S. +- `CNodeProp` — `NodeProp.h`. Mag point property: A (potential, complex), J (point current, complex). Tier S. +- equivalent for elec/curr/heat (`bd_NodeProp.h`, `cd_NodeProp.h`, `hd_NodeProp.h`) — small variants. Tier S. +- `CCircProp` — `CircProp.h`. Mag circuit: name, total current (complex), series/parallel radio. Tier S. +- `bdCCircProp`, `cdCCircProp`, `hdCCircProp` — `bd_CircProp.h`, `cd_CircProp.h`, `hd_CircProp.h`. V/Q name plus voltage/charge or voltage/heat-source radio. Tier S. + +geometry-input and modal helpers (shared) + +- `CEnterPt` — `EnterPt.h`. TAB-coordinate entry (rectangular or polar). Tier S. +- `CArcDlg` — `ArcDlg.h`. Arc creation: angle, max-side-length, boundary. Tier S. +- `CGroupNumber` — `GroupNumber.h`. Single int. Tier S. +- `CCopyDlg` — `CopyDlg.h`. Rotate or translate, with about-point or delta, n-copies, IsMove flag. Tier S. +- `CMirrorDlg` — `MirrorDlg.h`. Two-point axis pa/pb. Tier S. +- `CScaleDlg` — `ScaleDlg.h`. Base point, scale factor. Tier S. +- `CGridMod` / `GRIDDLG` — `GridMod.h`, `GRIDDLG.H`. Grid size and coords (cartesian/polar). Tier S. +- `CPromptBox` — typed-radius prompt used by `OnCreateRadius`. Tier S. +- `CMakeABCDlg` — `MakeABCDlg.h`. Open-boundary construction: number of layers, radius, center x/y, edge-type combo. Tier M. +- `CExteriorProps` — `ExteriorProps.h`. Axisymmetric exterior region: Ri, Ro, Zo. Tier S. +- `CDXFImport` — `DXFImport.h`. DXF import tolerance. Tier S. + +problem-definition dialog (one per discipline) + +- `probdlg` — `probdlg.h`. Magnetostatic: problem type (planar / axisymmetric), length units, smart-mesh, frequency, precision, min-angle, depth, problem note, solver type (succ/newton), previous-solution path and prev-type. Tier M. +- `bdCProbDlg` — `bd_probdlg.h`. Electrostatic. Same shape minus frequency/solver/prev. Tier S. +- `cdCProbDlg` — `cd_probdlg.h`. Current flow: adds frequency. Tier S. +- `hdCProbDlg` — `hd_probdlg.h`. Heat: adds dt (time step) and previous-solution path. Tier S. + +material-library editor + +- `fe_CLibDlg` — `fe_libdlg.h` and `bd_libdlg.h` / `cd_libdlg.h` / `hd_libdlg.h`. Drag-and-drop tree (library on left, model on right). Backed by `LibFolderInfo`. Tier L. The .lib file format is read-line-by-line in `fe_libdlg.cpp::ParseLine`. The drag/drop behavior depends on Win32 CTreeCtrl, so the iced version becomes two side-by-side scrollable lists with arrow buttons or drag handles. + +preferences + +- `CPref` / `bdCPref` / `cdCPref` / `hdCPref` — `Pref.h` and prefixed variants. Per-discipline-preprocessor defaults: edit action, coords (cart/polar), problem type, length units, solver, default freq, default grid size, default pixels-per-unit, depth, precision, min-angle, show-grid, snap-grid, show-origin, show-names, and a per-color picker driven by a combo. Tier M. +- `CGeneralPrefs` — `GeneralPrefs.h`. App-wide: default-doc type, default-Lua-console (drop), default-XY-plot, default-output-window, default-smartmesh. Tier S. +- the post-processor variant: `CViewPref`, `bvCPref`, etc. — see post-processor section. + +BH curve editor (mag-only material substructure) + +- `CBHData` — `BHData.h`. Two text columns (B, H) plus name, with linear/log plot preview. Reads BH datafiles via `BHDatafile.cpp`. Tier M. + +### geometry processing + +- PSLG enforcement (`MOVECOPY.CPP::EnforcePSLG`, `FancyEnforcePSLG`, and per-discipline `bd_movecopy.cpp` / `cd_movecopy.cpp` / `hd_movecopy.cpp`) — the rebuild-list approach: empty the four lists, re-`AddNode` / `AddSegment` / `AddArcSegment` / `AddBlockLabel` everything, and the `Add*` functions check tolerance, split crossings, and reject overlaps. The fancy version takes an explicit tolerance for DXF import. Tier L. +- intersection math — `CFemmeDoc::GetIntersection` (line-line), `GetLineArcIntersection`, `GetArcArcIntersection`, `ShortestDistanceFromArc`, `GetCircle` (recover center/radius from two endpoints + arc angle). Pure geometry; port directly. Tier M. +- `CanCreateRadius` / `CreateRadius` — fillet two adjacent segments at a node with a given radius. Tier M. +- mirror, rotate, translate, scale — `RotateMove`, `TranslateMove`, `ScaleMove`, `MirrorSelected`, `RotateCopy`, `TranslateCopy`. The move variants act on the selection; copy variants accept n-copies and emit fresh entities. Each ends with a call to `EnforcePSLG`. Tier M. +- DXF import — `DXFImport.cpp` plus per-discipline overrides. Reads ASCII DXF r12: LINE, ARC, CIRCLE, LWPOLYLINE entities translate into segments and arcs, then `FancyEnforcePSLG(tol)` runs. Tier L. +- DXF export — `WriteDXF`. Walks lists and writes LINE/ARC entities. Tier M. +- mesh launch — `OnWritePoly` in `writepoly.cpp` (mag) + discipline variants. Writes a `.poly` file: nodes (with index, x, y, marker), segments (with marker referencing boundary id), holes (block labels marked ``), region attributes (each non-hole block label gets a region id and area target). Then `CreateProcess` invokes `triangle.exe -p -P -q -e -A -a -z -Q -I `. `FunnyOnWritePoly` handles periodic-BC enforcement (nodes on the periodic boundary must come in matched pairs). After Triangle returns, `LoadMesh` reads back `.node` / `.ele` / `.edge` to populate `meshnode` / `meshline` / `greymeshline`. Tier L. We exec the Triangle binary the same way the original does — only difference is `posix_spawn` instead of `CreateProcess`. + +## post-processor + +The post-processor opens an `.ans` file (for elec `.res`, current `.res`, heat `.anh`) produced by the solver. Each discipline has its own pair of files. The view classes are again structurally identical; the differences are in plot quantities and integral kinds. + +### loading + +- `CFemmviewDoc::OnOpenDocument` — `FemmviewDoc.cpp:226`. Reopens the matching .fem (to get nodes/lines/labels/properties), then reads the `.ans` body: solver settings header, mesh nodes (each with A, NumList), mesh elements (each with 3 node indices, lab, b1, b2, mu1, mu2), point/segment/block/arc property arrays, optional air-gap-element block. Allocates the `NumList`/`ConList` neighbor-list adjacency. Tier L. +- `belaviewDoc.cpp`, `cviewDoc.cpp`/`CVIEWDOC.CPP`, `hviewDoc.cpp` — same shape, different per-element quantities (V/D, V/J, T/F respectively). Tier M each once mag is done. + +### plot modes + +The view stores `DensityPlot` (an int code), `NumContours`, `ShowAr`, `ShowAi`, `ShowMask`, plus `LegendFlag`, `GreyContours`, `Smooth`, `PtsFlag`, `MeshFlag`, `VectorPlot`, `VectorScaleFactor`, `PlotBounds[9][2]`. `OnDraw` (FemmviewView.cpp:777) renders, in order: density fill, mesh wireframe (if MeshFlag), contour lines (if ShowAr/ShowAi), the PSLG geometry, vectors (if VectorPlot), the user-defined contour, the mask, names, mesh points, and the legend. + +mag density-plot options (`cvCDPlotDlg2` — yes, the mag view borrows the current-flow dialog class; see "notes and surprises") cover, depending on DC vs AC, the quantities |B|, Bx/Br, By/Bz, |H|, Hx/Hr, Hy/Hz, |J|, J_re, J_im. Sources: `FemmviewView.cpp::OnDplot`, `cv_DPlotDlg2.h`. + +- mesh wireframe — `OnShowMesh`. Tier S. +- input geometry overlay (always on). Tier S. +- equipotential / flux-line contours (mag) — uses real-A and optionally imag-A; banded by `NumContours`; `DoContours` is the line walker in elements. Tier M. +- density plot (mag) — fills each element with a color-mapped shade of B/H/J magnitude or component, optionally smoothed using nodal values. `PlotFluxDensity` per-element fill. Tier L. +- vector plot (mag) — `OnVplot` -> `CVPlotDlg`. Choices: none, B, H (DC) or B_re, H_re, B_im, H_im, B_re&B_im, H_re&H_im (AC). Tier M. +- point-info readout — `OnRButtonDblClk` in EditAction==0 mode shows mesh node (x,y,A,B,H,mu) at the click. `DisplayPointProperties` formats. Tier S. +- mesh-point markers — `PtsFlag`. Tier S. +- legend overlay — `LegendFlag` plus the 20-bin colormap (`Color00..Color19`, `Grey00..Grey19`). Tier S. +- mask display — show the Henrotte weighting-function mask used by the "Henrotte" force/torque integrals. Built by `MakeMask` / `bv_makemask.cpp` / `cv_makemask.cpp` / `hv_makemask.cpp`. Tier M. + +post-processor selection modes (`EditAction` again, but different meanings) + +- `0` Point — click reports values; right-double-click shows closest input node. Tier S. +- `1` Contour — left-click extends a user-defined contour (`pDoc->contour`) by snapping to closest node or to the closest segment/arc; SHIFT pops `CBendContourDlg` to bend the last contour leg into an arc; DEL removes the last point; ESC clears. Tier M. +- `2` Block — left-click toggles selection of the block under the cursor (uses `InTriangle` to find which mesh element, then floods to the connected region via `ConList`). Tier M. + +### integrations and contour ops + +The line integral dialog `CLIntDlg` exposes (mag): normal flux (B.n), MMF along H.t, contour length / surface area, force-via-stress-tensor, torque-via-stress-tensor, and (B.n)^2. The block integral dialog `CBlockInt` exposes 17 entries enumerated in `FemmviewView.cpp::OnMenuIntegrate` cases 0..16: A.J (coil energy), integral of A, energy in B.H/2, hysteresis losses, total losses, area, ohmic losses, total current, integral of B, volume, Lorentz force (DC + 2X), Lorentz torque, energy density, Maxwell-stress (Henrotte) force, Maxwell-stress torque, integral of R^2, R/g (specific reluctance). Air-gap integrals (`GapIntegral`) compute force/torque/energy across an Air Gap Element boundary using harmonics. Sources: `FemmviewView.cpp:2952-3515`, `FemmviewDoc.cpp::BlockIntegral`, `FemmviewDoc.cpp::LineIntegral`, `FemmviewDoc.cpp::lua_gapintegral`. Tier XL — these are the algorithmically heaviest part of the post-processor. + +elec / current / heat integrals are smaller subsets (force and torque only, for `cv_BlockInt.cpp` and `cv_LIntDlg.cpp`). + +### plot dialogs + +- `CCPlotDlg` / `CCplotDlg2` — `CPlotDlg.h`, `CplotDlg2.h`. Mag contour-plot setup: num contours, A-bounds (lower/upper), show real, show imag, show mask. Two variants because AC has imag part. Tier S. +- `cvCDPlotDlg2` — `cv_DPlotDlg2.h`. Density-plot setup: which scalar to plot (combo), show-it toggle, show legend, grayscale toggle, lower/upper bound. Tier S. +- `CVPlotDlg` / `bv_VPlotDlg` / `cv_VPlotDlg` / `hv_VPlotDlg` — vector-plot kind combo and scale factor. Tier S. +- `CXYPlotDlg` / `bv_XYPlotDlg` / `cv_XYPlotDlg` / `hv_XYPlotDlg` — XY plot setup for the user-contour: quantity combo (potential, |B|, B.n, B.t, |H|, H.n, H.t, plus J_eddy, Js+J_eddy in AC; for materials also Mu_x/Mu_y per block), N points, to-file plus file-format combo (text / matlab / mathematica / gnuplot / csv). The plot itself renders to a metafile via `Xyplot.cpp`. Tier M. +- `CGapPlotDlg` — `GapPlotDlg.h`. Air-gap plot: pick AGE, quantity, N points, file. Tier S. +- `CLIntDlg`, `cv_LIntDlg`, `hv_LIntDlg` — line integral type combo. Tier S. +- `CBlockInt`, `cv_BlockInt`, `hv_BlockInt` — block integral type combo. Tier S. +- `bv_CircDlg`, `cv_CircDlg`, `hv_CircDlg` — circuit-results readout. Tier S. + +### post-processor preferences + +- `CViewPref` / `bvCPref` / `cvCPref` / `hvCPref` — `viewpref.h` and prefixed variants. Edit-action default, default density-plot kind, default vector-plot kind, weighting scheme (Henrotte / uniform / extended), num contours default, grey-contours, legend, grid, show-Ar, show-Ai, show-mask, snap, smooth, shift-H, line-integral N points, plot N points, show-mesh, show-points, show-names, color overrides. Tier M. + +## shared infrastructure + +- grid (size, snap, show, polar vs cartesian) — same code path on the canvas in both pre and post. Tier S. +- screen↔world transforms — `ScreenToDwg` / `DwgToScreen`, currently implemented in `doc_canvas.rs`. Tier already done. +- zoom-to-fit (`OnZoomNatural`) — bounding-box over all nodes/labels/arcs with margin. Tier S. +- keyboard zoom (`CKbdZoom`) — typed window extents. Tier S. +- preferences file IO — `femme.cfg`, `belasolv.cfg`, `csolv.cfg`, `hsolv.cfg` are text key/value pairs (``, ``, ``, ...). `ScanPreferences` / `WritePreferences` per view+doc pair. Tier S. +- undo stack — `UpdateUndo` snapshots all four lists into `undonodelist` / `undolinelist` / `undoarclist` / `undoblocklist` before any mutation; `Undo` swaps. One level deep. Tier S. +- group numbering — every entity carries an `InGroup` int; group-mode selection and modal group editing are the only consumers. Tier S. +- material library on-disk — text in a custom block-keyed format; one library file per discipline, with vendor URL/name metadata. Tier M. +- BH curve evaluation — `BHData` and `BHDatafile.cpp` implement cubic-spline interpolation for nonlinear materials; the result is consumed by the solver, not the GUI, so the GUI side is just edit + serialize. Tier M. + +## suggested order + +each milestone names the smallest visible end-to-end result, not just an internal step. + +1. **finish mag PSLG enforcement and segment / arc click-tools.** Port `EnforcePSLG`, `GetIntersection`, `GetLineArcIntersection`, `GetArcArcIntersection`, `CanCreateRadius`/`CreateRadius`. Add AddSegment/AddArc canvas tools that pick the first node on click 1 and finalize on click 2. Demo: draw two crossing segments and watch them get split at the intersection. + +2. **selection model + DEL/move/copy.** Add per-entity selected flag, right-click toggle, box-select drag, group right-click, and the `CCopyDlg`/`CMirrorDlg`/`CScaleDlg` modals. Demo: build a quarter geometry, mirror-copy it to full. + +3. **mag problem-definition dialog + property list dialogs.** Implement `probdlg` and `CPtProp` (the four-discipline shell), plus `CMatDlg`, `CBdryDlg`, `CNodeProp`, `CCircProp`. Demo: open the bundled `demo.fem`, edit a material, save, reopen, see the change. + +4. **per-selection property dialogs.** `COpNodeDlg`, `COpSegDlg`, `COpArcSegDlg`, `COpBlkDlg`. Demo: click a label, hit space, change its block material from the dropdown, see the canvas redraw with the new label text. + +5. **mesh launch end-to-end (mag).** Port `OnWritePoly` (skip periodic for now), call out to the existing `triangle` binary, parse `.node`/`.ele` back into a mesh layer, draw it. Demo: hit F2 (or a "Mesh" button), see triangles overlay. + +6. **mag solver launch.** Wire the existing FFI to call `fkn` on the in-memory doc (or via the saved .fem/.poly path), then move on to step 7. Demo: hit "Analyze", solver runs, .ans appears on disk. + +7. **mag post-processor loading + mesh display + density plot.** Port `CFemmviewDoc::OnOpenDocument`, the neighbor-list builder, `PlotFluxDensity` for |B| only, plus the colormap. Demo: open a solved .ans, see a colored flux density plot. + +8. **mag post-processor contour plot + point-info readout.** `DoContours` for real-A, the equipotential drawing, plus mouse-click point query. Demo: click anywhere, get A/B/H/mu at that point. + +9. **mag line integrals.** Contour edit mode (`EditAction=1` in the postprocessor), `CLIntDlg`, `LineIntegral` for cases 0 (flux), 1 (MMF), 2 (length/area). Demo: draw a contour across an airgap, read off normal flux. + +10. **mag block integrals.** Block selection mode, `CBlockInt`, all 17 `BlockIntegral` cases. Demo: select coil, read coil energy and current; select rotor, read Lorentz torque. + +11. **port the other three disciplines.** Reuse the now-proven shell. The doc parsers exist already; the work is dialog families and (smaller) integral sets per discipline. + +12. **DXF, material library, BH editor, open-boundary builder.** Power-user features. Defer until the four disciplines work end-to-end. + +## notes and surprises + +- mag's density-plot dialog class is named `cvCDPlotDlg2` and lives in `cv_DPlotDlg2.h` (the current-flow viewer's file). The mag viewer doesn't have its own — `FemmviewView.cpp::OnDplot` instantiates the current-flow class directly. Not a bug, just an MFC class-rename that didn't get propagated. New iced port should have one `DensityPlotDialog` parameterized by an enum of plot kinds, not four near-copies. +- `EditAction` has different meanings in pre and post processor. In pre: 0=node, 1=segment, 2=label, 3=arc, 4=group. In post: 0=point, 1=contour, 2=block. The two views share the field name but not the semantics. +- the post-processor opens a *separate document* (`CFemmviewDoc`) rather than augmenting the pre-processor doc. The `.ans` loader rereads the `.fem` first to grab the original geometry and property lists, then reads the solved-element block on top. In the Rust port this should be the same: one `PostDoc` that holds a `FemmDoc` plus the solved mesh and per-element fields, not a giant union. +- `OnMenuAnalyze` shells out to `fkn.exe` via `CreateProcess` with the document path. In the port, the same model is fine — we already build per-engine archives, and the `analyze(path)` Rust call replaces the process spawn. Air-gap-element support (the `agelist` array) is a recent addition driving the gap-integral and gap-plot dialogs; many demo .fem files won't have any AGEs, so the GUI must handle `agelist.size==0` everywhere. +- the BH-curve dialog stores B and H as multi-line text strings (`CString m_Bdata; CString m_Hdata;`) and parses them on OK via `StripBHData`. Idiomatic Rust replacement: two synchronized Vec. +- `bd_nosebl.cpp` / `cd_nosebl.cpp` / `hd_nosebl.cpp` plus `NOSEBL.CPP` carry per-discipline copies of the geometry parser, even though the geometry sections are identical across disciplines. The Rust parser crates have already collapsed these into shared types where appropriate — keep doing that for any new shared code. +- the original references a "Lua Console" prominently in menus and the `lnu_*` view functions, and many doc methods are split into a non-Lua and a Lua entry point (`AddNode` and `lua_addnode`). Per project rules, the Lua console is out of scope; the `lua_*` static methods can be ignored entirely. The Lua interpreter linked through the FFI is used by the solver for material formulas (functional MagDir, BH expressions) only. +- the `triangle/triangle.c` mesher is left as a separate executable invocation. Don't try to inline it. +- `cFemmeView::OnMakeABC` produces an open-boundary "asymptotic boundary condition" approximation by laying down N concentric circular shells, each with a calculated mu and sigma — it builds geometry and material/boundary properties together. The dialog has a combo for the edge type (Dirichlet / Neumann), defaults sourced from `MakeABCDlg.cpp`. Useful but not on the critical path. +- two grid dialogs exist: `CGridMod` (`GridMod.h`) and `GRIDDLG` (`GRIDDLG.H`). They differ only in the coord-system combo wiring. Pick one for the port. +- `KbdZoom`, `MyMsgBox`, `PromptBox`, `OutBox`, `MaskProgress`, `LuaEdit`, `MyTabCtrl`, `MDITabs`, `MyRecentFileList`, `MyCommandLineInfo` are MFC convenience widgets that have no analogue and no port. The native iced widgets cover their roles. +- `ActiveFEMM.cpp` / `mathlink.h` are the ActiveX/Mathematica IPC shim — explicitly out of scope per the brief. +- `femmplotDoc.cpp` + `femmplotView.cpp` are a tiny separate document type that pops up to display XY-plot metafiles from the post-processor. It is not a "post processor" itself, just a viewer. The Rust port can render the same plot inside the same iced window — no separate doc class needed. diff --git a/crates/femm-app/Cargo.toml b/crates/femm-app/Cargo.toml new file mode 100644 index 0000000..688c3c8 --- /dev/null +++ b/crates/femm-app/Cargo.toml @@ -0,0 +1,17 @@ +[package] +name = "femm-app" +version = "0.0.1" +edition.workspace = true +rust-version.workspace = true +publish.workspace = true +description = "iced shell for the FEMM 4.2 port" + +[[bin]] +name = "femm" +path = "src/main.rs" + +[dependencies] +femm-sys = { workspace = true } +femm-doc = { workspace = true } +iced = { version = "0.14", features = ["canvas"] } +rfd = "0.17" diff --git a/crates/femm-app/assets/demo.fem b/crates/femm-app/assets/demo.fem new file mode 100644 index 0000000..2fba2a8 --- /dev/null +++ b/crates/femm-app/assets/demo.fem @@ -0,0 +1,107 @@ +[Format] = 4.0 +[Frequency] = 60 +[Precision] = 1e-08 +[MinAngle] = 30 +[DoSmartMesh] = 1 +[Depth] = 25.4 +[LengthUnits] = millimeters +[ProblemType] = planar +[Coordinates] = cartesian +[ACSolver] = 0 +[PrevType] = 0 +[PrevSoln] = "" +[Comment] = "demo: square air domain with a copper coil block in the middle" +[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] = 2 + + = "Air" + = 1 + = 1 + = 0 + = 0 + = 0 + = 0 + = 0 + = 0 + = 0 + = 0 + = 0 + = 0 + = 1 + = 0 + = 0 + = 0 + + + = "Copper" + = 1 + = 1 + = 0 + = 0 + = 0 + = 0 + = 58 + = 0 + = 0 + = 0 + = 0 + = 0 + = 1 + = 0 + = 0 + = 0 + +[CircuitProps] = 1 + + = "Coil" + = 100 + = 0 + = 1 + +[NumPoints] = 8 +-50 -50 1 0 +50 -50 1 0 +50 50 1 0 +-50 50 1 0 +-10 -10 0 0 +10 -10 0 0 +10 10 0 0 +-10 10 0 0 +[NumSegments] = 8 +0 1 -1 1 0 0 +1 2 -1 1 0 0 +2 3 -1 1 0 0 +3 0 -1 1 0 0 +4 5 -1 0 0 0 +5 6 -1 0 0 0 +6 7 -1 0 0 0 +7 4 -1 0 0 0 +[NumArcSegments] = 0 +[NumHoles] = 0 +[NumBlockLabels] = 2 +0 0 2 1 1 0 0 100 0 +30 30 1 1 0 0 0 1 0 diff --git a/crates/femm-app/src/doc_canvas.rs b/crates/femm-app/src/doc_canvas.rs new file mode 100644 index 0000000..294bc89 --- /dev/null +++ b/crates/femm-app/src/doc_canvas.rs @@ -0,0 +1,338 @@ +//! draws a FemmDoc on an iced canvas: nodes, segments, arcs, block labels, with pan/zoom and click-to-add. + +use femm_doc::FemmDoc; +use iced::widget::canvas::{ + self, Action, Canvas, Event, Frame, Geometry, Path, Stroke, Text, path::Builder, +}; +use iced::{Color, Element, Length, Point, Radians, Rectangle, Renderer, Theme, Vector, mouse}; + +const PADDING_PX: f32 = 24.0; +const NODE_RADIUS: f32 = 3.0; +const STROKE_WIDTH: f32 = 1.2; +const LABEL_TICK_PX: f32 = 6.0; +const ZOOM_STEP: f32 = 1.1; +const ZOOM_MIN: f32 = 0.05; +const ZOOM_MAX: f32 = 200.0; +const CLICK_DRAG_THRESHOLD_PX: f32 = 4.0; +const BG: Color = Color::WHITE; +const GEOM: Color = Color::BLACK; +const LABEL_COLOR: Color = Color::from_rgb(0.25, 0.45, 0.85); + +/// active editing mode on the canvas. +#[derive(Debug, Default, Clone, Copy, PartialEq, Eq)] +pub enum Tool { + #[default] + Select, + AddNode, + AddBlockLabel, +} + +/// messages emitted by the canvas back to the app. +#[derive(Debug, Clone)] +pub enum CanvasMessage { + /// left-click at the given doc-world coordinate, intent depends on the active tool. + Click { world: (f64, f64), tool: Tool }, +} + +/// pan offset and zoom factor applied on top of fit-to-view. +#[derive(Debug, Default, Clone, Copy)] +pub struct ViewState { + pan: Vector, + zoom: f32, + drag_origin: Option, + press_origin: Option, + dragged: bool, +} + +/// constructs the canvas widget for a doc reference. +pub fn view<'a>(doc: &'a FemmDoc, tool: Tool) -> Element<'a, CanvasMessage> { + Canvas::new(DocCanvas { doc, tool }) + .width(Length::Fill) + .height(Length::Fill) + .into() +} + +struct DocCanvas<'a> { + doc: &'a FemmDoc, + tool: Tool, +} + +impl<'a> canvas::Program for DocCanvas<'a> { + type State = ViewState; + + fn update( + &self, + state: &mut Self::State, + event: &Event, + bounds: Rectangle, + cursor: mouse::Cursor, + ) -> Option> { + match event { + Event::Mouse(mouse::Event::ButtonPressed(mouse::Button::Left)) => { + if let Some(p) = cursor.position_in(bounds) { + state.press_origin = Some(p); + state.dragged = false; + if self.tool == Tool::Select { + state.drag_origin = Some(p); + } + return Some(Action::capture()); + } + } + Event::Mouse(mouse::Event::ButtonReleased(mouse::Button::Left)) => { + let press = state.press_origin.take(); + state.drag_origin = None; + let was_dragged = std::mem::take(&mut state.dragged); + if !was_dragged && self.tool != Tool::Select { + if let (Some(start), Some(now)) = (press, cursor.position_in(bounds)) { + if (now.x - start.x).abs() < CLICK_DRAG_THRESHOLD_PX + && (now.y - start.y).abs() < CLICK_DRAG_THRESHOLD_PX + { + let view = ViewTransform::fit(self.doc, bounds, state); + let world = view.inverse_map(now); + return Some(Action::publish(CanvasMessage::Click { + world, + tool: self.tool, + })); + } + } + } + if press.is_some() { + return Some(Action::capture()); + } + } + Event::Mouse(mouse::Event::CursorMoved { .. }) => { + if let (Some(prev), Some(now)) = (state.drag_origin, cursor.position_in(bounds)) { + state.pan = state.pan + Vector::new(now.x - prev.x, now.y - prev.y); + state.drag_origin = Some(now); + state.dragged = true; + return Some(Action::request_redraw().and_capture()); + } + if let (Some(start), Some(now)) = (state.press_origin, cursor.position_in(bounds)) { + if (now.x - start.x).abs() > CLICK_DRAG_THRESHOLD_PX + || (now.y - start.y).abs() > CLICK_DRAG_THRESHOLD_PX + { + state.dragged = true; + } + } + } + Event::Mouse(mouse::Event::WheelScrolled { delta }) => { + if let Some(focus) = cursor.position_in(bounds) { + let lines = match delta { + mouse::ScrollDelta::Lines { y, .. } => *y, + mouse::ScrollDelta::Pixels { y, .. } => *y / 40.0, + }; + if state.zoom == 0.0 { state.zoom = 1.0; } + let prev = state.zoom; + let next = (prev * ZOOM_STEP.powf(lines)).clamp(ZOOM_MIN, ZOOM_MAX); + let factor = next / prev; + state.pan.x = focus.x - factor * (focus.x - state.pan.x); + state.pan.y = focus.y - factor * (focus.y - state.pan.y); + state.zoom = next; + return Some(Action::request_redraw().and_capture()); + } + } + Event::Keyboard(iced::keyboard::Event::KeyPressed { key, .. }) => { + if let iced::keyboard::Key::Character(c) = key { + if c.as_str().eq_ignore_ascii_case("r") { + *state = ViewState::default(); + return Some(Action::request_redraw()); + } + } + } + _ => {} + } + None + } + + fn draw( + &self, + state: &Self::State, + renderer: &Renderer, + _theme: &Theme, + bounds: Rectangle, + _cursor: mouse::Cursor, + ) -> Vec> { + let mut frame = Frame::new(renderer, bounds.size()); + frame.fill_rectangle(Point::ORIGIN, bounds.size(), BG); + + let view = ViewTransform::fit(self.doc, bounds, state); + + for s in &self.doc.segments { + if let (Some(p0), Some(p1)) = + (self.doc.nodes.get(s.n0 as usize), self.doc.nodes.get(s.n1 as usize)) + { + let a = view.map(p0.x, p0.y); + let b = view.map(p1.x, p1.y); + frame.stroke(&Path::line(a, b), + Stroke::default().with_width(STROKE_WIDTH).with_color(GEOM)); + } + } + + for a in &self.doc.arcs { + if let (Some(p0), Some(p1)) = + (self.doc.nodes.get(a.n0 as usize), self.doc.nodes.get(a.n1 as usize)) + { + if let Some((center, radius, start_angle, end_angle)) = + arc_geometry(p0.x, p0.y, p1.x, p1.y, a.arc_length, a.normal_direction) + { + let cx = view.map(center.0, center.1); + let r_px = (radius * view.scale) as f32; + let mut b = Builder::new(); + b.arc(canvas::path::Arc { + center: cx, + radius: r_px, + start_angle, + end_angle, + }); + frame.stroke(&b.build(), + Stroke::default().with_width(STROKE_WIDTH).with_color(GEOM)); + } else { + let a_px = view.map(p0.x, p0.y); + let b_px = view.map(p1.x, p1.y); + frame.stroke(&Path::line(a_px, b_px), + Stroke::default().with_width(STROKE_WIDTH).with_color(GEOM)); + } + } + } + + for n in &self.doc.nodes { + let p = view.map(n.x, n.y); + frame.fill(&Path::circle(p, NODE_RADIUS), GEOM); + } + + for label in &self.doc.block_labels { + let p = view.map(label.x, label.y); + let cross = Path::new(|b| { + b.move_to(Point::new(p.x - LABEL_TICK_PX, p.y)); + b.line_to(Point::new(p.x + LABEL_TICK_PX, p.y)); + b.move_to(Point::new(p.x, p.y - LABEL_TICK_PX)); + b.line_to(Point::new(p.x, p.y + LABEL_TICK_PX)); + }); + frame.stroke(&cross, + Stroke::default().with_width(STROKE_WIDTH).with_color(LABEL_COLOR)); + frame.fill_text(Text { + content: label.block_type.clone(), + position: Point::new(p.x + LABEL_TICK_PX + 4.0, p.y - 8.0), + color: LABEL_COLOR, + size: 12.0.into(), + ..Text::default() + }); + } + + vec![frame.into_geometry()] + } + + fn mouse_interaction( + &self, + state: &Self::State, + bounds: Rectangle, + cursor: mouse::Cursor, + ) -> mouse::Interaction { + if state.drag_origin.is_some() { + return mouse::Interaction::Grabbing; + } + if cursor.position_in(bounds).is_some() { + return match self.tool { + Tool::Select => mouse::Interaction::Grab, + Tool::AddNode | Tool::AddBlockLabel => mouse::Interaction::Crosshair, + }; + } + mouse::Interaction::default() + } +} + +/// affine map from doc coordinates to canvas pixels, composing fit-to-view with user pan/zoom. +struct ViewTransform { + scale: f64, + offset: Vector, + y_max: f64, +} + +impl ViewTransform { + fn fit(doc: &FemmDoc, bounds: Rectangle, view: &ViewState) -> Self { + let (xmin, xmax, ymin, ymax) = doc_bounds(doc); + let dx = (xmax - xmin).max(1e-9); + let dy = (ymax - ymin).max(1e-9); + let avail_w = (bounds.width as f64 - 2.0 * PADDING_PX as f64).max(1.0); + let avail_h = (bounds.height as f64 - 2.0 * PADDING_PX as f64).max(1.0); + let base_scale = (avail_w / dx).min(avail_h / dy); + let user_zoom = if view.zoom <= 0.0 { 1.0 } else { view.zoom }; + let scale = base_scale * user_zoom as f64; + + let drawn_w = dx * scale; + let drawn_h = dy * scale; + let pad_x = (bounds.width as f64 - drawn_w) / 2.0 - xmin * scale + view.pan.x as f64; + let pad_y = (bounds.height as f64 - drawn_h) / 2.0 + view.pan.y as f64; + + let _ = ymin; + Self { + scale, + offset: Vector::new(pad_x as f32, pad_y as f32), + y_max: ymax, + } + } + + fn map(&self, x: f64, y: f64) -> Point { + let px = (x * self.scale) as f32 + self.offset.x; + let py = ((self.y_max - y) * self.scale) as f32 + self.offset.y; + Point::new(px, py) + } + + /// inverse of [`map`]: converts a canvas pixel back to a doc-world coordinate. + fn inverse_map(&self, p: Point) -> (f64, f64) { + let x = (p.x - self.offset.x) as f64 / self.scale; + let y = self.y_max - (p.y - self.offset.y) as f64 / self.scale; + (x, y) + } +} + +fn doc_bounds(doc: &FemmDoc) -> (f64, f64, f64, f64) { + let mut xmin = f64::INFINITY; + let mut xmax = f64::NEG_INFINITY; + let mut ymin = f64::INFINITY; + let mut ymax = f64::NEG_INFINITY; + let mut had_point = false; + for n in &doc.nodes { + xmin = xmin.min(n.x); xmax = xmax.max(n.x); + ymin = ymin.min(n.y); ymax = ymax.max(n.y); + had_point = true; + } + for l in &doc.block_labels { + xmin = xmin.min(l.x); xmax = xmax.max(l.x); + ymin = ymin.min(l.y); ymax = ymax.max(l.y); + had_point = true; + } + if !had_point { return (-1.0, 1.0, -1.0, 1.0); } + if (xmax - xmin).abs() < 1e-9 { xmin -= 0.5; xmax += 0.5; } + if (ymax - ymin).abs() < 1e-9 { ymin -= 0.5; ymax += 0.5; } + (xmin, xmax, ymin, ymax) +} + +/// derives center, radius, and angular extent of an arc segment from its endpoints +/// and degree sweep. returns None when endpoints coincide. +fn arc_geometry( + x0: f64, y0: f64, + x1: f64, y1: f64, + arc_length_deg: f64, + _normal_direction: bool, +) -> Option<((f64, f64), f64, Radians, Radians)> { + let theta = arc_length_deg.to_radians(); + if theta.abs() < 1e-9 { return None; } + let dx = x1 - x0; + let dy = y1 - y0; + let chord = (dx * dx + dy * dy).sqrt(); + if chord < 1e-12 { return None; } + let radius = chord / (2.0 * (theta / 2.0).sin().abs()); + let mid_x = (x0 + x1) / 2.0; + let mid_y = (y0 + y1) / 2.0; + let perp_x = -dy / chord; + let perp_y = dx / chord; + let h = radius * (theta / 2.0).cos(); + let sign = if theta > 0.0 { 1.0 } else { -1.0 }; + let cx = mid_x + sign * perp_x * h; + let cy = mid_y + sign * perp_y * h; + + let a0 = (-(y0 - cy)).atan2(x0 - cx); + let a1 = (-(y1 - cy)).atan2(x1 - cx); + Some(((cx, cy), radius, Radians(a0 as f32), Radians(a1 as f32))) +} diff --git a/crates/femm-app/src/main.rs b/crates/femm-app/src/main.rs new file mode 100644 index 0000000..2352a1f --- /dev/null +++ b/crates/femm-app/src/main.rs @@ -0,0 +1,137 @@ +//! iced shell entry point for the FEMM 4.2 port. + +mod doc_canvas; + +use doc_canvas::{CanvasMessage, Tool}; +use femm_doc::FemmDoc; +use iced::widget::{button, column, container, row, text}; +use iced::{Element, Length, Task}; + +const DEMO_FEM: &str = include_str!("../assets/demo.fem"); +const ADD_TOLERANCE: f64 = 0.5; + +#[derive(Debug, Clone)] +enum Message { + OpenFem, + SelectTool(Tool), + Canvas(CanvasMessage), +} + +struct App { + doc: FemmDoc, + source_label: String, + status: String, + tool: Tool, +} + +impl App { + fn new() -> (Self, Task) { + let doc = FemmDoc::parse(DEMO_FEM).unwrap_or_default(); + let app = App { + doc, + source_label: String::from("demo.fem (embedded)"), + status: String::new(), + tool: Tool::Select, + }; + (app, Task::none()) + } + + fn title(&self) -> String { + format!("femm42 - {}", self.source_label) + } + + fn update(&mut self, msg: Message) -> Task { + match msg { + Message::OpenFem => { + let picked = rfd::FileDialog::new() + .add_filter("FEMM magnetostatic", &["fem", "FEM"]) + .add_filter("All files", &["*"]) + .pick_file(); + + if let Some(path) = picked { + let label = path + .file_name() + .and_then(|n| n.to_str()) + .unwrap_or("?") + .to_string(); + match FemmDoc::open(&path) { + Ok(d) => { + self.doc = d; + self.source_label = label; + self.status = String::new(); + } + Err(e) => { + self.status = format!("failed to open {label}: {e}"); + } + } + } + } + Message::SelectTool(t) => { + self.tool = t; + } + Message::Canvas(CanvasMessage::Click { world, tool }) => { + match tool { + Tool::AddNode => { + let idx = self.doc.add_node(world.0, world.1, ADD_TOLERANCE); + 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.status = format!("block label {idx} at ({:.3}, {:.3})", world.0, world.1); + } + Tool::Select => {} + } + } + } + Task::none() + } + + fn view(&self) -> Element<'_, Message> { + let stats = text(format!( + "{} nodes {} segments {} arcs {} labels", + self.doc.nodes.len(), + self.doc.segments.len(), + self.doc.arcs.len(), + self.doc.block_labels.len(), + )) + .size(12); + + let toolbar = row![ + button("Open .fem...").on_press(Message::OpenFem), + tool_button("Select", Tool::Select, self.tool), + tool_button("Add Node", Tool::AddNode, self.tool), + tool_button("Add Label", Tool::AddBlockLabel, self.tool), + text(&self.source_label).size(13), + stats, + ] + .spacing(8); + + let canvas = doc_canvas::view(&self.doc, self.tool).map(Message::Canvas); + + let mut body = column![toolbar].spacing(8).padding(12); + body = body.push(canvas); + if !self.status.is_empty() { + body = body.push(text(&self.status).size(12)); + } + + container(body) + .width(Length::Fill) + .height(Length::Fill) + .into() + } +} + +fn tool_button(label: &str, this_tool: Tool, active: Tool) -> Element<'_, Message> { + let btn = button(text(label).size(13)).on_press(Message::SelectTool(this_tool)); + if this_tool == active { + btn.style(button::primary).into() + } else { + btn.style(button::secondary).into() + } +} + +fn main() -> iced::Result { + iced::application(App::new, App::update, App::view) + .title(App::title) + .run() +} diff --git a/crates/femm-doc-curr/Cargo.toml b/crates/femm-doc-curr/Cargo.toml new file mode 100644 index 0000000..d76aef3 --- /dev/null +++ b/crates/femm-doc-curr/Cargo.toml @@ -0,0 +1,14 @@ +[package] +name = "femm-doc-curr" +version = "0.0.1" +edition.workspace = true +rust-version.workspace = true +publish.workspace = true +description = "rust-native .fec document model: parser, writer, geometry editing" + +[lib] +path = "src/lib.rs" + +[dependencies] +num-complex = "0.4" +thiserror = "2" diff --git a/crates/femm-doc-curr/src/geom.rs b/crates/femm-doc-curr/src/geom.rs new file mode 100644 index 0000000..4fae6b4 --- /dev/null +++ b/crates/femm-doc-curr/src/geom.rs @@ -0,0 +1,51 @@ +//! geometric primitives mirroring the .fec control-point model. + +/// control point in the planar geometry. +#[derive(Debug, Clone, Default)] +pub struct Node { + pub x: f64, + pub y: f64, + pub boundary_marker: String, + pub in_group: i32, + pub selected: bool, +} + +/// straight segment joining two control points by index. +#[derive(Debug, Clone, Default)] +pub struct Segment { + pub n0: i32, + pub n1: i32, + pub max_side_length: f64, + pub boundary_marker: String, + pub hidden: bool, + pub in_group: i32, + pub selected: bool, +} + +/// arc segment from n0 to n1, sweeping the given degree arc. +#[derive(Debug, Clone, Default)] +pub struct ArcSegment { + pub n0: i32, + pub n1: i32, + pub arc_length: f64, + pub max_side_length: f64, + pub boundary_marker: String, + pub hidden: bool, + pub in_group: i32, + pub normal_direction: bool, + pub selected: bool, +} + +/// region label assigning a material, conductor, and mesh-area constraint to an enclosed area. +#[derive(Debug, Clone, Default)] +pub struct BlockLabel { + pub x: f64, + pub y: f64, + pub max_area: f64, + pub block_type: String, + pub in_conductor: String, + pub in_group: i32, + pub is_external: bool, + pub is_default: bool, + pub selected: bool, +} diff --git a/crates/femm-doc-curr/src/lib.rs b/crates/femm-doc-curr/src/lib.rs new file mode 100644 index 0000000..cf6f41f --- /dev/null +++ b/crates/femm-doc-curr/src/lib.rs @@ -0,0 +1,74 @@ +//! current-flow pre-processor document model: geometry, properties, parser, writer. + +pub mod geom; +pub mod props; +pub mod parser; +pub mod writer; + +use num_complex::Complex64; + +pub use geom::{ArcSegment, BlockLabel, Node, Segment}; +pub use props::{BoundaryProp, ConductorProp, MaterialProp, PointProp}; + +#[derive(Debug, Clone, Copy, PartialEq, Eq)] +pub enum ProblemType { + Planar, + Axisymmetric, +} + +impl Default for ProblemType { + fn default() -> Self { ProblemType::Planar } +} + +#[derive(Debug, Clone, Copy, PartialEq, Eq, Default)] +pub enum Coords { + #[default] + Cartesian, + Polar, +} + +#[derive(Debug, Clone, Copy, PartialEq, Eq)] +pub enum LengthUnit { + Inches = 0, + Millimeters = 1, + Centimeters = 2, + Meters = 3, + Mils = 4, + Microns = 5, +} + +impl Default for LengthUnit { + fn default() -> Self { LengthUnit::Inches } +} + +/// current-flow .fec document: geometry, property tables, problem attributes. +#[derive(Debug, Clone, Default)] +pub struct FemmDoc { + 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 problem_type: ProblemType, + pub coords: Coords, + pub comment: String, + + pub ext_ro: f64, + pub ext_ri: f64, + pub ext_zo: f64, + + pub nodes: Vec, + pub segments: Vec, + pub arcs: Vec, + pub block_labels: Vec, + + pub materials: Vec, + pub boundaries: Vec, + pub points: Vec, + pub conductors: Vec, +} + +/// f64 complex number alias matching the C++ CComplex memory layout. +pub type Complex = Complex64; diff --git a/crates/femm-doc-curr/src/parser.rs b/crates/femm-doc-curr/src/parser.rs new file mode 100644 index 0000000..a2ad2de --- /dev/null +++ b/crates/femm-doc-curr/src/parser.rs @@ -0,0 +1,402 @@ +//! .fec text parser: `[Section] = value` scalars and `...` property stanzas. + +use crate::{ + ArcSegment, BlockLabel, BoundaryProp, ConductorProp, Coords, FemmDoc, LengthUnit, + MaterialProp, Node, PointProp, ProblemType, Segment, +}; +use num_complex::Complex64; +use std::path::Path; +use thiserror::Error; + +#[derive(Debug, Error)] +pub enum ParseError { + #[error("io error: {0}")] + Io(#[from] std::io::Error), + #[error("legacy 3.2-format .fec files are not supported")] + LegacyFormat, + #[error("unexpected end of file while reading {context}")] + UnexpectedEof { context: &'static str }, + #[error("malformed value in section {section}: {reason}")] + MalformedValue { section: &'static str, reason: String }, +} + +impl FemmDoc { + /// parses a .fec text buffer into a document. + pub fn parse(src: &str) -> Result { + let mut doc = FemmDoc::default(); + doc.depth = -1.0; + doc.smart_mesh = true; + doc.format = 1.0; + + if src + .lines() + .next() + .map(|l| l.trim_start().starts_with("Frequency")) + .unwrap_or(false) + { + return Err(ParseError::LegacyFormat); + } + + let mut lines = src.lines(); + + let mut point = PointProp::default(); + let mut bdry = BoundaryProp::default(); + let mut mat = MaterialProp::default(); + let mut cond = ConductorProp::default(); + + while let Some(raw) = lines.next() { + let tok = first_token(raw); + if tok.is_empty() { continue; } + + // problem attributes + if tok.eq_ignore_ascii_case("[format]") { + doc.format = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[frequency]") { + doc.frequency = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[depth]") { + doc.depth = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[precision]") { + doc.precision = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[minangle]") { + doc.min_angle = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[dosmartmesh]") { + doc.smart_mesh = parse_i32(strip_key(raw)) != 0; + } else if tok.eq_ignore_ascii_case("[lengthunits]") { + doc.length_units = parse_length_units(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[problemtype]") { + doc.problem_type = parse_problem_type(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[coordinates]") { + doc.coords = parse_coords(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[comment]") { + doc.comment = unescape_comment(unquote(strip_key(raw))); + } else if tok.eq_ignore_ascii_case("[extzo]") { + doc.ext_zo = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[extro]") { + doc.ext_ro = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[extri]") { + doc.ext_ri = parse_f64(strip_key(raw)); + + // point property stanza + } else if tok.eq_ignore_ascii_case("") { + point = PointProp::default(); + point.name = String::from("New Point Property"); + } else if tok.eq_ignore_ascii_case("") { + point.name = unquote(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + point.vp.re = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + point.vp.im = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + point.qp.re = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + point.qp.im = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + doc.points.push(std::mem::take(&mut point)); + + // boundary property stanza + } else if tok.eq_ignore_ascii_case("") { + bdry = BoundaryProp::default(); + bdry.name = String::from("New Boundary"); + } else if tok.eq_ignore_ascii_case("") { + bdry.name = unquote(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.format = parse_i32(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.vs.re = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.vs.im = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.qs.re = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.qs.im = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.c0.re = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.c0.im = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.c1.re = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.c1.im = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + doc.boundaries.push(std::mem::take(&mut bdry)); + + // material property stanza + } else if tok.eq_ignore_ascii_case("") { + mat = MaterialProp::default(); + mat.name = String::from("New Material"); + mat.ex = 1.0; + mat.ey = 1.0; + } else if tok.eq_ignore_ascii_case("") { + mat.name = unquote(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + mat.ox = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + mat.oy = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + mat.ex = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + mat.ey = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + mat.ltx = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + mat.lty = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + doc.materials.push(std::mem::take(&mut mat)); + + // conductor stanza + } else if tok.eq_ignore_ascii_case("") { + cond = ConductorProp::default(); + cond.name = String::from("New Conductor"); + } else if tok.eq_ignore_ascii_case("") { + cond.name = unquote(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + cond.vc.re = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + cond.vc.im = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + cond.qc.re = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + cond.qc.im = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + cond.conductor_type = parse_i32(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + doc.conductors.push(std::mem::take(&mut cond)); + + // counted geometry sections + } else if tok.eq_ignore_ascii_case("[numpoints]") { + let n = parse_i32(strip_key(raw)); + for _ in 0..n { + let line = lines.next().ok_or(ParseError::UnexpectedEof { context: "[NumPoints]" })?; + let (x, rest) = take_f64(line); + let (y, rest) = take_f64(rest); + let (bi, rest) = take_i32(rest); + let (group, rest) = take_i32(rest); + let (_cond_idx, _) = take_i32(rest); + let bm = resolve_name(bi, doc.points.iter().map(|p| p.name.as_str())); + doc.nodes.push(Node { + x, y, + boundary_marker: bm, + in_group: group, + selected: false, + }); + } + } else if tok.eq_ignore_ascii_case("[numsegments]") { + let n = parse_i32(strip_key(raw)); + for _ in 0..n { + let line = lines.next().ok_or(ParseError::UnexpectedEof { context: "[NumSegments]" })?; + let (n0, rest) = take_i32(line); + let (n1, rest) = take_i32(rest); + let (msl, rest) = take_f64(rest); + let (bi, rest) = take_i32(rest); + let (hidden, rest) = take_i32(rest); + let (group, rest) = take_i32(rest); + let (_cond_idx, _) = take_i32(rest); + let bm = resolve_name(bi, doc.boundaries.iter().map(|b| b.name.as_str())); + doc.segments.push(Segment { + n0, n1, + max_side_length: msl, + boundary_marker: bm, + hidden: hidden != 0, + in_group: group, + selected: false, + }); + } + } else if tok.eq_ignore_ascii_case("[numarcsegments]") { + let n = parse_i32(strip_key(raw)); + for _ in 0..n { + let line = lines.next().ok_or(ParseError::UnexpectedEof { context: "[NumArcSegments]" })?; + let (n0, rest) = take_i32(line); + let (n1, rest) = take_i32(rest); + let (al, rest) = take_f64(rest); + let (msl, rest) = take_f64(rest); + let (bi, rest) = take_i32(rest); + let (hidden, rest) = take_i32(rest); + let (group, rest) = take_i32(rest); + let (_cond_idx, _rest) = take_i32(rest); + let bm = resolve_name(bi, doc.boundaries.iter().map(|b| b.name.as_str())); + doc.arcs.push(ArcSegment { + n0, n1, + arc_length: al, + max_side_length: msl, + boundary_marker: bm, + hidden: hidden != 0, + in_group: group, + normal_direction: true, + selected: false, + }); + } + } else if tok.eq_ignore_ascii_case("[numholes]") { + let n = parse_i32(strip_key(raw)); + for _ in 0..n { + let line = lines.next().ok_or(ParseError::UnexpectedEof { context: "[NumHoles]" })?; + let (x, rest) = take_f64(line); + let (y, rest) = take_f64(rest); + let (group, _) = take_i32(rest); + doc.block_labels.push(BlockLabel { + x, y, + max_area: 0.0, + block_type: String::from(""), + in_conductor: String::from(""), + in_group: group, + is_external: false, + is_default: false, + selected: false, + }); + } + } else if tok.eq_ignore_ascii_case("[numblocklabels]") { + let n = parse_i32(strip_key(raw)); + for _ in 0..n { + let line = lines.next().ok_or(ParseError::UnexpectedEof { context: "[NumBlockLabels]" })?; + let (x, rest) = take_f64(line); + let (y, rest) = take_f64(rest); + let (mi, rest) = take_i32(rest); + let (ma_diam, rest) = take_f64(rest); + let (group, rest) = take_i32(rest); + let (ext_flags, _rest) = take_i32(rest); + + let block_type = if mi == 0 { + String::from("") + } else { + resolve_name(mi, doc.materials.iter().map(|m| m.name.as_str())) + }; + let max_area = if ma_diam < 0.0 { + 0.0 + } else { + std::f64::consts::PI * ma_diam * ma_diam / 4.0 + }; + + doc.block_labels.push(BlockLabel { + x, y, + max_area, + block_type, + in_conductor: String::from(""), + in_group: group, + is_external: ext_flags & 1 != 0, + is_default: ext_flags & 2 != 0, + selected: false, + }); + } + } + // unknown tokens are skipped silently, matching the original parser. + } + + // 3.2-era files omitted [Depth]; fill the default in the current length unit. + if doc.depth == -1.0 { + doc.depth = match doc.length_units { + LengthUnit::Millimeters => 1000.0, + LengthUnit::Centimeters => 100.0, + LengthUnit::Meters => 1.0, + LengthUnit::Mils => 1000.0 / 0.0254, + LengthUnit::Microns => 1.0e6, + LengthUnit::Inches => 1.0 / 0.0254, + }; + } + + Ok(doc) + } + + /// loads a .fec file by path. + pub fn open(path: impl AsRef) -> Result { + let text = std::fs::read_to_string(path)?; + Self::parse(&text) + } +} + +fn first_token(line: &str) -> &str { + line.trim_start().split_whitespace().next().unwrap_or("") +} + +fn strip_key(line: &str) -> &str { + match line.find('=') { + Some(i) => &line[i + 1..], + None => "", + } +} + +fn unquote(s: &str) -> String { + let s = s.trim(); + let first = s.find('"'); + let last = s.rfind('"'); + match (first, last) { + (Some(a), Some(b)) if a < b => s[a + 1..b].to_string(), + _ => s.to_string(), + } +} + +fn unescape_comment(s: String) -> String { + let bytes = s.as_bytes(); + let mut out = String::with_capacity(s.len()); + let mut i = 0; + while i < bytes.len() { + if bytes[i] == b'\\' && i + 1 < bytes.len() && bytes[i + 1] == b'n' { + out.push('\r'); + out.push('\n'); + i += 2; + } else { + out.push(bytes[i] as char); + i += 1; + } + } + out +} + +fn parse_f64(s: &str) -> f64 { + let s = s.trim().trim_matches(|c: char| c == ',' || c.is_whitespace()); + let end = s.find(|c: char| c.is_whitespace() || c == ',').unwrap_or(s.len()); + s[..end].parse::().unwrap_or(0.0) +} + +fn parse_i32(s: &str) -> i32 { + let s = s.trim().trim_matches(|c: char| c == ',' || c.is_whitespace()); + let end = s.find(|c: char| c.is_whitespace() || c == ',').unwrap_or(s.len()); + s[..end].parse::().unwrap_or(0) +} + +/// peels one f64 off the head of a whitespace/comma-separated row, returning the remainder. +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/comma-separated row, returning the remainder. +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) +} + +fn parse_length_units(s: &str) -> LengthUnit { + let t = first_token(s).to_ascii_lowercase(); + if t.starts_with("millimeters") { LengthUnit::Millimeters } + else if t.starts_with("c") { LengthUnit::Centimeters } + else if t.starts_with("meters") { LengthUnit::Meters } + else if t.starts_with("mils") { LengthUnit::Mils } + else if t.starts_with("microns") { LengthUnit::Microns } + else { LengthUnit::Inches } +} + +fn parse_problem_type(s: &str) -> ProblemType { + let t = first_token(s).to_ascii_lowercase(); + if t.starts_with("axi") { ProblemType::Axisymmetric } else { ProblemType::Planar } +} + +fn parse_coords(s: &str) -> Coords { + let t = first_token(s).to_ascii_lowercase(); + if t.starts_with("polar") { Coords::Polar } else { Coords::Cartesian } +} + +/// resolves a 1-based property index into the matching property name, falling back to empty. +fn resolve_name<'a, I: Iterator>(idx: i32, names: I) -> String { + if idx <= 0 { return String::new(); } + let want = idx as usize; + for (i, name) in names.enumerate() { + if i + 1 == want { return name.to_string(); } + } + String::new() +} + +#[allow(dead_code)] +fn _force_complex_use(_c: Complex64) {} diff --git a/crates/femm-doc-curr/src/props.rs b/crates/femm-doc-curr/src/props.rs new file mode 100644 index 0000000..7144beb --- /dev/null +++ b/crates/femm-doc-curr/src/props.rs @@ -0,0 +1,43 @@ +//! current-flow material, boundary, point, and conductor properties. + +use crate::Complex; + +/// material property entry referenced by name from block labels. +#[derive(Debug, Clone, Default)] +pub struct MaterialProp { + pub name: String, + pub ox: f64, + pub oy: f64, + pub ex: f64, + pub ey: f64, + pub ltx: f64, + pub lty: f64, +} + +/// boundary condition entry referenced by name from segments and arcs. +#[derive(Debug, Clone, Default)] +pub struct BoundaryProp { + pub name: String, + pub format: i32, + pub vs: Complex, + pub qs: Complex, + pub c0: Complex, + pub c1: Complex, +} + +/// point property entry referenced by name from nodes. +#[derive(Debug, Clone, Default)] +pub struct PointProp { + pub name: String, + pub vp: Complex, + pub qp: Complex, +} + +/// conductor property entry referenced by name from block labels. +#[derive(Debug, Clone, Default)] +pub struct ConductorProp { + pub name: String, + pub vc: Complex, + pub qc: Complex, + pub conductor_type: i32, +} diff --git a/crates/femm-doc-curr/src/writer.rs b/crates/femm-doc-curr/src/writer.rs new file mode 100644 index 0000000..6222b04 --- /dev/null +++ b/crates/femm-doc-curr/src/writer.rs @@ -0,0 +1,189 @@ +//! .fec text writer, inverse of [`parser`](crate::parser). + +use crate::{Coords, FemmDoc, LengthUnit, ProblemType}; +use std::fmt::Write; +use std::path::Path; + +impl FemmDoc { + /// renders the document to .fec text. + pub fn write(&self) -> String { + let mut out = String::new(); + writeln!(out, "[Format] = 1").unwrap(); + writeln!(out, "[Precision] = {:.17e}", self.precision).unwrap(); + writeln!(out, "[Frequency] = {:.17e}", self.frequency).unwrap(); + writeln!(out, "[MinAngle] = {:.17e}", self.min_angle).unwrap(); + writeln!(out, "[DoSmartMesh] = {}", self.smart_mesh as i32).unwrap(); + writeln!(out, "[Depth] = {:.17e}", self.depth).unwrap(); + + let units = match self.length_units { + LengthUnit::Inches => "inches", + LengthUnit::Millimeters => "millimeters", + LengthUnit::Centimeters => "centimeters", + LengthUnit::Meters => "meters", + LengthUnit::Mils => "mils", + LengthUnit::Microns => "microns", + }; + writeln!(out, "[LengthUnits] = {units}").unwrap(); + + match self.problem_type { + ProblemType::Planar => writeln!(out, "[ProblemType] = planar").unwrap(), + ProblemType::Axisymmetric => { + writeln!(out, "[ProblemType] = axisymmetric").unwrap(); + if self.ext_ro != 0.0 && self.ext_ri != 0.0 { + writeln!(out, "[extZo] = {:.17e}", self.ext_zo).unwrap(); + writeln!(out, "[extRo] = {:.17e}", self.ext_ro).unwrap(); + writeln!(out, "[extRi] = {:.17e}", self.ext_ri).unwrap(); + } + } + } + + let coords = match self.coords { + Coords::Cartesian => "cartesian", + Coords::Polar => "polar", + }; + writeln!(out, "[Coordinates] = {coords}").unwrap(); + + writeln!(out, "[Comment] = \"{}\"", escape_comment(&self.comment)).unwrap(); + + // point properties + writeln!(out, "[PointProps] = {}", self.points.len()).unwrap(); + for p in &self.points { + writeln!(out, " ").unwrap(); + writeln!(out, " = \"{}\"", p.name).unwrap(); + writeln!(out, " = {:.17e}", p.vp.re).unwrap(); + writeln!(out, " = {:.17e}", p.vp.im).unwrap(); + writeln!(out, " = {:.17e}", p.qp.re).unwrap(); + writeln!(out, " = {:.17e}", p.qp.im).unwrap(); + writeln!(out, " ").unwrap(); + } + + // boundary properties + writeln!(out, "[BdryProps] = {}", self.boundaries.len()).unwrap(); + for b in &self.boundaries { + writeln!(out, " ").unwrap(); + writeln!(out, " = \"{}\"", b.name).unwrap(); + writeln!(out, " = {}", b.format).unwrap(); + writeln!(out, " = {:.17e}", b.vs.re).unwrap(); + writeln!(out, " = {:.17e}", b.vs.im).unwrap(); + writeln!(out, " = {:.17e}", b.qs.re).unwrap(); + writeln!(out, " = {:.17e}", b.qs.im).unwrap(); + writeln!(out, " = {:.17e}", b.c0.re).unwrap(); + writeln!(out, " = {:.17e}", b.c0.im).unwrap(); + writeln!(out, " = {:.17e}", b.c1.re).unwrap(); + writeln!(out, " = {:.17e}", b.c1.im).unwrap(); + writeln!(out, " ").unwrap(); + } + + // material properties + writeln!(out, "[BlockProps] = {}", self.materials.len()).unwrap(); + for m in &self.materials { + writeln!(out, " ").unwrap(); + writeln!(out, " = \"{}\"", m.name).unwrap(); + writeln!(out, " = {:.17e}", m.ox).unwrap(); + writeln!(out, " = {:.17e}", m.oy).unwrap(); + writeln!(out, " = {:.17e}", m.ex).unwrap(); + writeln!(out, " = {:.17e}", m.ey).unwrap(); + writeln!(out, " = {:.17e}", m.ltx).unwrap(); + writeln!(out, " = {:.17e}", m.lty).unwrap(); + writeln!(out, " ").unwrap(); + } + + // conductor properties + writeln!(out, "[ConductorProps] = {}", self.conductors.len()).unwrap(); + for c in &self.conductors { + writeln!(out, " ").unwrap(); + writeln!(out, " = \"{}\"", c.name).unwrap(); + writeln!(out, " = {:.17e}", c.vc.re).unwrap(); + writeln!(out, " = {:.17e}", c.vc.im).unwrap(); + writeln!(out, " = {:.17e}", c.qc.re).unwrap(); + writeln!(out, " = {:.17e}", c.qc.im).unwrap(); + writeln!(out, " = {}", c.conductor_type).unwrap(); + writeln!(out, " ").unwrap(); + } + + // nodes + writeln!(out, "[NumPoints] = {}", self.nodes.len()).unwrap(); + for n in &self.nodes { + let t = lookup_idx(&n.boundary_marker, self.points.iter().map(|p| p.name.as_str())); + writeln!(out, "{:.17e}\t{:.17e}\t{}\t{}\t0", n.x, n.y, t, n.in_group).unwrap(); + } + + // segments + writeln!(out, "[NumSegments] = {}", self.segments.len()).unwrap(); + for s in &self.segments { + let t = lookup_idx(&s.boundary_marker, self.boundaries.iter().map(|b| b.name.as_str())); + let msl = if s.max_side_length < 0.0 { + String::from("-1") + } else { + format!("{:.17e}", s.max_side_length) + }; + writeln!(out, "{}\t{}\t{}\t{}\t{}\t{}\t0", + s.n0, s.n1, msl, t, s.hidden as i32, s.in_group).unwrap(); + } + + // arc segments + writeln!(out, "[NumArcSegments] = {}", self.arcs.len()).unwrap(); + for a in &self.arcs { + let t = lookup_idx(&a.boundary_marker, self.boundaries.iter().map(|b| b.name.as_str())); + writeln!(out, "{}\t{}\t{:.17e}\t{:.17e}\t{}\t{}\t{}\t0\t{:.17e}", + a.n0, a.n1, a.arc_length, a.max_side_length, t, + a.hidden as i32, a.in_group, a.max_side_length).unwrap(); + } + + // holes (block labels typed "") + let holes: Vec<_> = self.block_labels.iter().filter(|b| b.block_type == "").collect(); + writeln!(out, "[NumHoles] = {}", holes.len()).unwrap(); + for h in &holes { + writeln!(out, "{:.17e}\t{:.17e}\t{}", h.x, h.y, h.in_group).unwrap(); + } + + // regional attributes + let labels: Vec<_> = self.block_labels.iter().filter(|b| b.block_type != "").collect(); + writeln!(out, "[NumBlockLabels] = {}", labels.len()).unwrap(); + for b in labels { + let mi = lookup_idx(&b.block_type, self.materials.iter().map(|m| m.name.as_str())); + let area_field = if b.max_area > 0.0 { + format!("{:.17e}", (4.0 * b.max_area / std::f64::consts::PI).sqrt()) + } else { + String::from("-1") + }; + let flags = (b.is_external as i32) + ((b.is_default as i32) << 1); + writeln!(out, "{:.17e}\t{:.17e}\t{}\t{}\t{}\t{}", + b.x, b.y, mi, area_field, b.in_group, flags).unwrap(); + } + + out + } + + /// writes the document to a .fec file by path. + pub fn save(&self, path: impl AsRef) -> std::io::Result<()> { + std::fs::write(path, self.write()) + } +} + +/// finds the 1-based index of `name` in the given name iterator, or 0 when absent or empty. +fn lookup_idx<'a, I: Iterator>(name: &str, iter: I) -> i32 { + if name.is_empty() { return 0; } + for (i, n) in iter.enumerate() { + if n == name { return (i as i32) + 1; } + } + 0 +} + +/// inverse of `unescape_comment`: turns CR/LF pairs back into `\n` escapes. +fn escape_comment(s: &str) -> String { + let mut out = String::with_capacity(s.len()); + let bytes = s.as_bytes(); + let mut i = 0; + while i < bytes.len() { + if i + 1 < bytes.len() && bytes[i] == b'\r' && bytes[i + 1] == b'\n' { + out.push('\\'); + out.push('n'); + i += 2; + } else { + out.push(bytes[i] as char); + i += 1; + } + } + out +} diff --git a/crates/femm-doc-curr/tests/roundtrip.rs b/crates/femm-doc-curr/tests/roundtrip.rs new file mode 100644 index 0000000..dd1fdd1 --- /dev/null +++ b/crates/femm-doc-curr/tests/roundtrip.rs @@ -0,0 +1,198 @@ +use femm_doc_curr::FemmDoc; + +const FIXTURE: &str = r#"[Format] = 1 +[Precision] = 1e-08 +[Frequency] = 1000 +[MinAngle] = 30 +[DoSmartMesh] = 1 +[Depth] = 1 +[LengthUnits] = millimeters +[ProblemType] = planar +[Coordinates] = cartesian +[Comment] = "two-electrode current flow" +[PointProps] = 1 + + = "V=0" + = 0.5 + = 0.25 + = 1.5 + = -0.75 + +[BdryProps] = 1 + + = "outer" + = 2 + = 12 + = 3 + = 0.1 + = -0.2 + = 4 + = 5 + = 6 + = 7 + +[BlockProps] = 2 + + = "Copper" + = 5.8e7 + = 5.8e7 + = 1 + = 1 + = 0 + = 0 + + + = "FR4" + = 0.01 + = 0.02 + = 4.3 + = 4.5 + = 0.02 + = 0.03 + +[ConductorProps] = 1 + + = "Pad" + = 10 + = 1 + = 0.5 + = -0.5 + = 1 + +[NumPoints] = 4 +0 0 1 0 0 +10 0 0 0 0 +10 10 0 0 0 +0 10 0 0 0 +[NumSegments] = 4 +0 1 -1 1 0 0 0 +1 2 -1 1 0 0 0 +2 3 -1 1 0 0 0 +3 0 -1 1 0 0 0 +[NumArcSegments] = 0 +[NumHoles] = 0 +[NumBlockLabels] = 1 +5 5 2 0.1 0 0 +"#; + +#[test] +fn parses_fixture_geometry() { + let doc = FemmDoc::parse(FIXTURE).expect("parse"); + assert_eq!(doc.frequency, 1000.0); + assert_eq!(doc.nodes.len(), 4); + assert_eq!(doc.segments.len(), 4); + assert_eq!(doc.arcs.len(), 0); + assert_eq!(doc.block_labels.len(), 1); + assert_eq!(doc.materials.len(), 2); + assert_eq!(doc.boundaries.len(), 1); + assert_eq!(doc.points.len(), 1); + assert_eq!(doc.conductors.len(), 1); + assert_eq!(doc.nodes[0].boundary_marker, "V=0"); + assert_eq!(doc.segments[0].boundary_marker, "outer"); + assert_eq!(doc.block_labels[0].block_type, "FR4"); + + let p = &doc.points[0]; + assert!((p.vp.re - 0.5).abs() < 1e-12); + assert!((p.vp.im - 0.25).abs() < 1e-12); + assert!((p.qp.re - 1.5).abs() < 1e-12); + assert!((p.qp.im + 0.75).abs() < 1e-12); + + let b = &doc.boundaries[0]; + assert_eq!(b.format, 2); + assert!((b.vs.re - 12.0).abs() < 1e-12); + assert!((b.vs.im - 3.0).abs() < 1e-12); + assert!((b.qs.re - 0.1).abs() < 1e-12); + assert!((b.qs.im + 0.2).abs() < 1e-12); + assert!((b.c0.re - 4.0).abs() < 1e-12); + assert!((b.c0.im - 5.0).abs() < 1e-12); + assert!((b.c1.re - 6.0).abs() < 1e-12); + assert!((b.c1.im - 7.0).abs() < 1e-12); + + let m1 = &doc.materials[1]; + assert_eq!(m1.name, "FR4"); + assert!((m1.ox - 0.01).abs() < 1e-12); + assert!((m1.oy - 0.02).abs() < 1e-12); + assert!((m1.ex - 4.3).abs() < 1e-12); + assert!((m1.ey - 4.5).abs() < 1e-12); + assert!((m1.ltx - 0.02).abs() < 1e-12); + assert!((m1.lty - 0.03).abs() < 1e-12); + + let c = &doc.conductors[0]; + assert_eq!(c.name, "Pad"); + assert_eq!(c.conductor_type, 1); + assert!((c.vc.re - 10.0).abs() < 1e-12); + assert!((c.vc.im - 1.0).abs() < 1e-12); + assert!((c.qc.re - 0.5).abs() < 1e-12); + assert!((c.qc.im + 0.5).abs() < 1e-12); +} + +#[test] +fn round_trips_parse_write_parse() { + let a = FemmDoc::parse(FIXTURE).expect("parse a"); + let text = a.write(); + let b = FemmDoc::parse(&text).expect("parse b"); + + assert_eq!(a.frequency, b.frequency); + assert_eq!(a.precision, b.precision); + assert_eq!(a.depth, b.depth); + assert_eq!(a.nodes.len(), b.nodes.len()); + assert_eq!(a.segments.len(), b.segments.len()); + assert_eq!(a.materials.len(), b.materials.len()); + assert_eq!(a.conductors.len(), b.conductors.len()); + assert_eq!(a.boundaries.len(), b.boundaries.len()); + assert_eq!(a.points.len(), b.points.len()); + + for (x, y) in a.nodes.iter().zip(b.nodes.iter()) { + assert!((x.x - y.x).abs() < 1e-12); + assert!((x.y - y.y).abs() < 1e-12); + assert_eq!(x.boundary_marker, y.boundary_marker); + } + for (x, y) in a.segments.iter().zip(b.segments.iter()) { + assert_eq!(x.n0, y.n0); + assert_eq!(x.n1, y.n1); + assert_eq!(x.boundary_marker, y.boundary_marker); + } + for (x, y) in a.materials.iter().zip(b.materials.iter()) { + assert_eq!(x.name, y.name); + assert!((x.ox - y.ox).abs() < 1e-12); + assert!((x.oy - y.oy).abs() < 1e-12); + assert!((x.ex - y.ex).abs() < 1e-12); + assert!((x.ey - y.ey).abs() < 1e-12); + assert!((x.ltx - y.ltx).abs() < 1e-12); + assert!((x.lty - y.lty).abs() < 1e-12); + } + for (x, y) in a.points.iter().zip(b.points.iter()) { + assert_eq!(x.name, y.name); + assert!((x.vp.re - y.vp.re).abs() < 1e-12); + assert!((x.vp.im - y.vp.im).abs() < 1e-12); + assert!((x.qp.re - y.qp.re).abs() < 1e-12); + assert!((x.qp.im - y.qp.im).abs() < 1e-12); + } + for (x, y) in a.boundaries.iter().zip(b.boundaries.iter()) { + assert_eq!(x.name, y.name); + assert_eq!(x.format, y.format); + assert!((x.vs.re - y.vs.re).abs() < 1e-12); + assert!((x.vs.im - y.vs.im).abs() < 1e-12); + assert!((x.qs.re - y.qs.re).abs() < 1e-12); + assert!((x.qs.im - y.qs.im).abs() < 1e-12); + assert!((x.c0.re - y.c0.re).abs() < 1e-12); + assert!((x.c0.im - y.c0.im).abs() < 1e-12); + assert!((x.c1.re - y.c1.re).abs() < 1e-12); + assert!((x.c1.im - y.c1.im).abs() < 1e-12); + } + for (x, y) in a.conductors.iter().zip(b.conductors.iter()) { + assert_eq!(x.name, y.name); + assert_eq!(x.conductor_type, y.conductor_type); + assert!((x.vc.re - y.vc.re).abs() < 1e-12); + assert!((x.vc.im - y.vc.im).abs() < 1e-12); + assert!((x.qc.re - y.qc.re).abs() < 1e-12); + assert!((x.qc.im - y.qc.im).abs() < 1e-12); + } + + assert_eq!(a.block_labels.len(), b.block_labels.len()); + let la = &a.block_labels[0]; + let lb = &b.block_labels[0]; + assert!((la.x - lb.x).abs() < 1e-12); + assert_eq!(la.block_type, lb.block_type); + assert!((la.max_area - lb.max_area).abs() < 1e-10); +} diff --git a/crates/femm-doc-elec/Cargo.toml b/crates/femm-doc-elec/Cargo.toml new file mode 100644 index 0000000..2277691 --- /dev/null +++ b/crates/femm-doc-elec/Cargo.toml @@ -0,0 +1,14 @@ +[package] +name = "femm-doc-elec" +version = "0.0.1" +edition.workspace = true +rust-version.workspace = true +publish.workspace = true +description = "rust-native .fee document model: parser, writer, geometry editing" + +[lib] +path = "src/lib.rs" + +[dependencies] +num-complex = "0.4" +thiserror = "2" diff --git a/crates/femm-doc-elec/src/geom.rs b/crates/femm-doc-elec/src/geom.rs new file mode 100644 index 0000000..444e000 --- /dev/null +++ b/crates/femm-doc-elec/src/geom.rs @@ -0,0 +1,54 @@ +//! geometric primitives mirroring the .fee control-point model. + +/// control point in the planar geometry. +#[derive(Debug, Clone, Default)] +pub struct Node { + pub x: f64, + pub y: f64, + pub boundary_marker: String, + pub in_conductor: String, + pub in_group: i32, + pub selected: bool, +} + +/// straight segment joining two control points by index. +#[derive(Debug, Clone, Default)] +pub struct Segment { + pub n0: i32, + pub n1: i32, + pub max_side_length: f64, + pub boundary_marker: String, + pub in_conductor: String, + pub hidden: bool, + pub in_group: i32, + pub selected: bool, +} + +/// arc segment from n0 to n1, sweeping the given degree arc. +#[derive(Debug, Clone, Default)] +pub struct ArcSegment { + pub n0: i32, + pub n1: i32, + pub arc_length: f64, + pub max_side_length: f64, + pub boundary_marker: String, + pub in_conductor: String, + pub hidden: bool, + pub in_group: i32, + pub normal_direction: bool, + pub selected: bool, +} + +/// region label assigning a material and mesh-area constraint to an enclosed area. +#[derive(Debug, Clone, Default)] +pub struct BlockLabel { + pub x: f64, + pub y: f64, + pub max_area: f64, + pub block_type: String, + pub in_conductor: String, + pub in_group: i32, + pub is_external: bool, + pub is_default: bool, + pub selected: bool, +} diff --git a/crates/femm-doc-elec/src/lib.rs b/crates/femm-doc-elec/src/lib.rs new file mode 100644 index 0000000..95276e4 --- /dev/null +++ b/crates/femm-doc-elec/src/lib.rs @@ -0,0 +1,68 @@ +//! electrostatic pre-processor document model: geometry, properties, parser, writer. + +pub mod geom; +pub mod props; +pub mod parser; +pub mod writer; + +pub use geom::{ArcSegment, BlockLabel, Node, Segment}; +pub use props::{BoundaryProp, ConductorProp, MaterialProp, PointProp}; + +#[derive(Debug, Clone, Copy, PartialEq, Eq)] +pub enum ProblemType { + Planar, + Axisymmetric, +} + +impl Default for ProblemType { + fn default() -> Self { ProblemType::Planar } +} + +#[derive(Debug, Clone, Copy, PartialEq, Eq, Default)] +pub enum Coords { + #[default] + Cartesian, + Polar, +} + +#[derive(Debug, Clone, Copy, PartialEq, Eq)] +pub enum LengthUnit { + Inches = 0, + Millimeters = 1, + Centimeters = 2, + Meters = 3, + Mils = 4, + Microns = 5, +} + +impl Default for LengthUnit { + fn default() -> Self { LengthUnit::Inches } +} + +/// electrostatic .fee document: geometry, property tables, problem attributes. +#[derive(Debug, Clone, Default)] +pub struct FemmDoc { + pub format: f64, + pub precision: f64, + pub min_angle: f64, + pub smart_mesh: bool, + pub depth: f64, + pub length_units: LengthUnit, + pub problem_type: ProblemType, + pub coords: Coords, + pub comment: String, + + pub ext_ro: f64, + pub ext_ri: f64, + pub ext_zo: f64, + + pub nodes: Vec, + pub segments: Vec, + pub arcs: Vec, + pub block_labels: Vec, + + pub materials: Vec, + pub boundaries: Vec, + pub points: Vec, + pub conductors: Vec, +} diff --git a/crates/femm-doc-elec/src/parser.rs b/crates/femm-doc-elec/src/parser.rs new file mode 100644 index 0000000..33b4f22 --- /dev/null +++ b/crates/femm-doc-elec/src/parser.rs @@ -0,0 +1,394 @@ +//! .fee text parser: `[Section] = value` scalars and `...` property stanzas. + +use crate::{ + ArcSegment, BlockLabel, BoundaryProp, ConductorProp, Coords, FemmDoc, LengthUnit, + MaterialProp, Node, PointProp, ProblemType, Segment, +}; +use std::path::Path; +use thiserror::Error; + +#[derive(Debug, Error)] +pub enum ParseError { + #[error("io error: {0}")] + Io(#[from] std::io::Error), + #[error("unexpected end of file while reading {context}")] + UnexpectedEof { context: &'static str }, + #[error("malformed value in section {section}: {reason}")] + MalformedValue { section: &'static str, reason: String }, +} + +impl FemmDoc { + /// parses a .fee text buffer into a document. + pub fn parse(src: &str) -> Result { + let mut doc = FemmDoc::default(); + doc.depth = -1.0; + doc.smart_mesh = true; + doc.format = 1.0; + + let mut lines = src.lines(); + + let mut point = PointProp::default(); + let mut bdry = BoundaryProp::default(); + let mut mat = MaterialProp::default(); + let mut cond = ConductorProp::default(); + + while let Some(raw) = lines.next() { + let tok = first_token(raw); + if tok.is_empty() { continue; } + + // problem attributes + if tok.eq_ignore_ascii_case("[format]") { + doc.format = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[depth]") { + doc.depth = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[precision]") { + doc.precision = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[minangle]") { + doc.min_angle = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[dosmartmesh]") { + doc.smart_mesh = parse_i32(strip_key(raw)) != 0; + } else if tok.eq_ignore_ascii_case("[lengthunits]") { + doc.length_units = parse_length_units(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[problemtype]") { + doc.problem_type = parse_problem_type(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[coordinates]") { + doc.coords = parse_coords(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[comment]") { + doc.comment = unescape_comment(unquote(strip_key(raw))); + } else if tok.eq_ignore_ascii_case("[extzo]") { + doc.ext_zo = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[extro]") { + doc.ext_ro = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[extri]") { + doc.ext_ri = parse_f64(strip_key(raw)); + + // point property stanza + } else if tok.eq_ignore_ascii_case("") { + point = PointProp::default(); + point.name = String::from("New Point Property"); + } else if tok.eq_ignore_ascii_case("") { + point.name = unquote(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + point.vp = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + point.qp = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + doc.points.push(std::mem::take(&mut point)); + + // boundary property stanza + } else if tok.eq_ignore_ascii_case("") { + bdry = BoundaryProp::default(); + bdry.name = String::from("New Boundary"); + } else if tok.eq_ignore_ascii_case("") { + bdry.name = unquote(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.format = parse_i32(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.vs = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.qs = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.c0 = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.c1 = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + doc.boundaries.push(std::mem::take(&mut bdry)); + + // material property stanza + } else if tok.eq_ignore_ascii_case("") { + mat = MaterialProp::default(); + mat.name = String::from("New Material"); + mat.ex = 1.0; + mat.ey = 1.0; + } else if tok.eq_ignore_ascii_case("") { + mat.name = unquote(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + mat.ex = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + mat.ey = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + mat.qv = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + doc.materials.push(std::mem::take(&mut mat)); + + // conductor stanza + } else if tok.eq_ignore_ascii_case("") { + cond = ConductorProp::default(); + cond.name = String::from("New Conductor"); + } else if tok.eq_ignore_ascii_case("") { + cond.name = unquote(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + cond.vc = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + cond.qc = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + cond.conductor_type = parse_i32(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + doc.conductors.push(std::mem::take(&mut cond)); + + // counted geometry sections + } else if tok.eq_ignore_ascii_case("[numpoints]") { + let n = parse_i32(strip_key(raw)); + for _ in 0..n { + let line = lines.next().ok_or(ParseError::UnexpectedEof { context: "[NumPoints]" })?; + let (x, rest) = take_f64(line); + let (y, rest) = take_f64(rest); + let (bi, rest) = take_i32(rest); + let (group, rest) = take_i32(rest); + let (ci, _) = take_i32(rest); + let bm = resolve_name(bi, doc.points.iter().map(|p| p.name.as_str())); + let cn = if ci == 0 { + String::from("") + } else { + resolve_name(ci, doc.conductors.iter().map(|c| c.name.as_str())) + }; + doc.nodes.push(Node { + x, y, + boundary_marker: bm, + in_conductor: cn, + in_group: group, + selected: false, + }); + } + } else if tok.eq_ignore_ascii_case("[numsegments]") { + let n = parse_i32(strip_key(raw)); + for _ in 0..n { + let line = lines.next().ok_or(ParseError::UnexpectedEof { context: "[NumSegments]" })?; + let (n0, rest) = take_i32(line); + let (n1, rest) = take_i32(rest); + let (msl, rest) = take_f64(rest); + let (bi, rest) = take_i32(rest); + let (hidden, rest) = take_i32(rest); + let (group, rest) = take_i32(rest); + let (ci, _) = take_i32(rest); + let bm = resolve_name(bi, doc.boundaries.iter().map(|b| b.name.as_str())); + let cn = if ci == 0 { + String::from("") + } else { + resolve_name(ci, doc.conductors.iter().map(|c| c.name.as_str())) + }; + doc.segments.push(Segment { + n0, n1, + max_side_length: msl, + boundary_marker: bm, + in_conductor: cn, + hidden: hidden != 0, + in_group: group, + selected: false, + }); + } + } else if tok.eq_ignore_ascii_case("[numarcsegments]") { + let n = parse_i32(strip_key(raw)); + for _ in 0..n { + let line = lines.next().ok_or(ParseError::UnexpectedEof { context: "[NumArcSegments]" })?; + let (n0, rest) = take_i32(line); + let (n1, rest) = take_i32(rest); + let (al, rest) = take_f64(rest); + let (msl, rest) = take_f64(rest); + let (bi, rest) = take_i32(rest); + let (hidden, rest) = take_i32(rest); + let (group, rest) = take_i32(rest); + let (ci, _) = take_i32(rest); + let bm = resolve_name(bi, doc.boundaries.iter().map(|b| b.name.as_str())); + let cn = if ci == 0 { + String::from("") + } else { + resolve_name(ci, doc.conductors.iter().map(|c| c.name.as_str())) + }; + doc.arcs.push(ArcSegment { + n0, n1, + arc_length: al, + max_side_length: msl, + boundary_marker: bm, + in_conductor: cn, + hidden: hidden != 0, + in_group: group, + normal_direction: true, + selected: false, + }); + } + } else if tok.eq_ignore_ascii_case("[numholes]") { + let n = parse_i32(strip_key(raw)); + for _ in 0..n { + let line = lines.next().ok_or(ParseError::UnexpectedEof { context: "[NumHoles]" })?; + let (x, rest) = take_f64(line); + let (y, rest) = take_f64(rest); + let (group, _) = take_i32(rest); + doc.block_labels.push(BlockLabel { + x, y, + max_area: 0.0, + block_type: String::from(""), + in_conductor: String::from(""), + in_group: group, + is_external: false, + is_default: false, + selected: false, + }); + } + } else if tok.eq_ignore_ascii_case("[numblocklabels]") { + let n = parse_i32(strip_key(raw)); + for _ in 0..n { + let line = lines.next().ok_or(ParseError::UnexpectedEof { context: "[NumBlockLabels]" })?; + let (x, rest) = take_f64(line); + let (y, rest) = take_f64(rest); + let (mi, rest) = take_i32(rest); + let (ma_diam, rest) = take_f64(rest); + let (group, rest) = take_i32(rest); + let (ext_flags, _) = take_i32(rest); + + let block_type = if mi == 0 { + String::from("") + } else { + resolve_name(mi, doc.materials.iter().map(|m| m.name.as_str())) + }; + let max_area = if ma_diam < 0.0 { + 0.0 + } else { + std::f64::consts::PI * ma_diam * ma_diam / 4.0 + }; + + doc.block_labels.push(BlockLabel { + x, y, + max_area, + block_type, + in_conductor: String::from(""), + in_group: group, + is_external: ext_flags & 1 != 0, + is_default: ext_flags & 2 != 0, + selected: false, + }); + } + } + // unknown tokens are skipped silently, matching the original parser. + } + + // 3.2-era files omitted [Depth]; fill the default in the current length unit. + if doc.depth == -1.0 { + doc.depth = match doc.length_units { + LengthUnit::Millimeters => 1000.0, + LengthUnit::Centimeters => 100.0, + LengthUnit::Meters => 1.0, + LengthUnit::Mils => 1000.0 / 0.0254, + LengthUnit::Microns => 1.0e6, + LengthUnit::Inches => 1.0 / 0.0254, + }; + } + + Ok(doc) + } + + /// loads a .fee file by path. + pub fn open(path: impl AsRef) -> Result { + let text = std::fs::read_to_string(path)?; + Self::parse(&text) + } +} + +fn first_token(line: &str) -> &str { + line.trim_start().split_whitespace().next().unwrap_or("") +} + +fn strip_key(line: &str) -> &str { + match line.find('=') { + Some(i) => &line[i + 1..], + None => "", + } +} + +fn unquote(s: &str) -> String { + let s = s.trim(); + let first = s.find('"'); + let last = s.rfind('"'); + match (first, last) { + (Some(a), Some(b)) if a < b => s[a + 1..b].to_string(), + _ => s.to_string(), + } +} + +fn unescape_comment(s: String) -> String { + let bytes = s.as_bytes(); + let mut out = String::with_capacity(s.len()); + let mut i = 0; + while i < bytes.len() { + if bytes[i] == b'\\' && i + 1 < bytes.len() && bytes[i + 1] == b'n' { + out.push('\r'); + out.push('\n'); + i += 2; + } else { + out.push(bytes[i] as char); + i += 1; + } + } + out +} + +fn parse_f64(s: &str) -> f64 { + let s = s.trim().trim_matches(|c: char| c == ',' || c.is_whitespace()); + let end = s.find(|c: char| c.is_whitespace() || c == ',').unwrap_or(s.len()); + s[..end].parse::().unwrap_or(0.0) +} + +fn parse_i32(s: &str) -> i32 { + let s = s.trim().trim_matches(|c: char| c == ',' || c.is_whitespace()); + let end = s.find(|c: char| c.is_whitespace() || c == ',').unwrap_or(s.len()); + s[..end].parse::().unwrap_or(0) +} + +/// peels one f64 off the head of a whitespace/comma-separated row, returning the remainder. +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/comma-separated row, returning the remainder. +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) +} + +/// extracts the contents of the first double-quoted span in `s`, or returns empty. +#[allow(dead_code)] +fn take_quoted(s: &str) -> String { + let first = s.find('"'); + if let Some(a) = first { + let after = &s[a + 1..]; + if let Some(b) = after.find('"') { + return after[..b].to_string(); + } + } + String::new() +} + +fn parse_length_units(s: &str) -> LengthUnit { + let t = first_token(s).to_ascii_lowercase(); + if t.starts_with("millimeters") { LengthUnit::Millimeters } + else if t.starts_with("c") { LengthUnit::Centimeters } + else if t.starts_with("meters") { LengthUnit::Meters } + else if t.starts_with("mils") { LengthUnit::Mils } + else if t.starts_with("microns") { LengthUnit::Microns } + else { LengthUnit::Inches } +} + +fn parse_problem_type(s: &str) -> ProblemType { + let t = first_token(s).to_ascii_lowercase(); + if t.starts_with("axi") { ProblemType::Axisymmetric } else { ProblemType::Planar } +} + +fn parse_coords(s: &str) -> Coords { + let t = first_token(s).to_ascii_lowercase(); + if t.starts_with("polar") { Coords::Polar } else { Coords::Cartesian } +} + +/// resolves a 1-based property index into the matching property name, falling back to empty. +fn resolve_name<'a, I: Iterator>(idx: i32, names: I) -> String { + if idx <= 0 { return String::new(); } + let want = idx as usize; + for (i, name) in names.enumerate() { + if i + 1 == want { return name.to_string(); } + } + String::new() +} diff --git a/crates/femm-doc-elec/src/props.rs b/crates/femm-doc-elec/src/props.rs new file mode 100644 index 0000000..1829e36 --- /dev/null +++ b/crates/femm-doc-elec/src/props.rs @@ -0,0 +1,38 @@ +//! electrostatic material, boundary, point, and conductor properties. + +/// material property entry referenced by name from block labels. +#[derive(Debug, Clone, Default)] +pub struct MaterialProp { + pub name: String, + pub ex: f64, + pub ey: f64, + pub qv: f64, +} + +/// boundary condition entry referenced by name from segments and arcs. +#[derive(Debug, Clone, Default)] +pub struct BoundaryProp { + pub name: String, + pub format: i32, + pub vs: f64, + pub qs: f64, + pub c0: f64, + pub c1: f64, +} + +/// point property entry referenced by name from nodes. +#[derive(Debug, Clone, Default)] +pub struct PointProp { + pub name: String, + pub vp: f64, + pub qp: f64, +} + +/// conductor property entry referenced by name from nodes, segments, and arcs. +#[derive(Debug, Clone, Default)] +pub struct ConductorProp { + pub name: String, + pub vc: f64, + pub qc: f64, + pub conductor_type: i32, +} diff --git a/crates/femm-doc-elec/src/writer.rs b/crates/femm-doc-elec/src/writer.rs new file mode 100644 index 0000000..a5c59a8 --- /dev/null +++ b/crates/femm-doc-elec/src/writer.rs @@ -0,0 +1,179 @@ +//! .fee text writer, inverse of [`parser`](crate::parser). + +use crate::{Coords, FemmDoc, LengthUnit, ProblemType}; +use std::fmt::Write; +use std::path::Path; + +impl FemmDoc { + /// renders the document to .fee text. + pub fn write(&self) -> String { + let mut out = String::new(); + writeln!(out, "[Format] = 1").unwrap(); + writeln!(out, "[Precision] = {:.17e}", self.precision).unwrap(); + writeln!(out, "[MinAngle] = {:.17e}", self.min_angle).unwrap(); + writeln!(out, "[DoSmartMesh] = {}", self.smart_mesh as i32).unwrap(); + writeln!(out, "[Depth] = {:.17e}", self.depth).unwrap(); + + let units = match self.length_units { + LengthUnit::Inches => "inches", + LengthUnit::Millimeters => "millimeters", + LengthUnit::Centimeters => "centimeters", + LengthUnit::Meters => "meters", + LengthUnit::Mils => "mils", + LengthUnit::Microns => "microns", + }; + writeln!(out, "[LengthUnits] = {units}").unwrap(); + + match self.problem_type { + ProblemType::Planar => writeln!(out, "[ProblemType] = planar").unwrap(), + ProblemType::Axisymmetric => { + writeln!(out, "[ProblemType] = axisymmetric").unwrap(); + if self.ext_ro != 0.0 && self.ext_ri != 0.0 { + writeln!(out, "[extZo] = {:.17e}", self.ext_zo).unwrap(); + writeln!(out, "[extRo] = {:.17e}", self.ext_ro).unwrap(); + writeln!(out, "[extRi] = {:.17e}", self.ext_ri).unwrap(); + } + } + } + + let coords = match self.coords { + Coords::Cartesian => "cartesian", + Coords::Polar => "polar", + }; + writeln!(out, "[Coordinates] = {coords}").unwrap(); + writeln!(out, "[Comment] = \"{}\"", escape_comment(&self.comment)).unwrap(); + + // point properties + writeln!(out, "[PointProps] = {}", self.points.len()).unwrap(); + for p in &self.points { + writeln!(out, " ").unwrap(); + writeln!(out, " = \"{}\"", p.name).unwrap(); + writeln!(out, " = {:.17e}", p.vp).unwrap(); + writeln!(out, " = {:.17e}", p.qp).unwrap(); + writeln!(out, " ").unwrap(); + } + + // boundary properties + writeln!(out, "[BdryProps] = {}", self.boundaries.len()).unwrap(); + for b in &self.boundaries { + writeln!(out, " ").unwrap(); + writeln!(out, " = \"{}\"", b.name).unwrap(); + writeln!(out, " = {}", b.format).unwrap(); + writeln!(out, " = {:.17e}", b.vs).unwrap(); + writeln!(out, " = {:.17e}", b.qs).unwrap(); + writeln!(out, " = {:.17e}", b.c0).unwrap(); + writeln!(out, " = {:.17e}", b.c1).unwrap(); + writeln!(out, " ").unwrap(); + } + + // material properties + writeln!(out, "[BlockProps] = {}", self.materials.len()).unwrap(); + for m in &self.materials { + writeln!(out, " ").unwrap(); + writeln!(out, " = \"{}\"", m.name).unwrap(); + writeln!(out, " = {:.17e}", m.ex).unwrap(); + writeln!(out, " = {:.17e}", m.ey).unwrap(); + writeln!(out, " = {:.17e}", m.qv).unwrap(); + writeln!(out, " ").unwrap(); + } + + // conductor properties + writeln!(out, "[ConductorProps] = {}", self.conductors.len()).unwrap(); + for c in &self.conductors { + writeln!(out, " ").unwrap(); + writeln!(out, " = \"{}\"", c.name).unwrap(); + writeln!(out, " = {:.17e}", c.vc).unwrap(); + writeln!(out, " = {:.17e}", c.qc).unwrap(); + writeln!(out, " = {}", c.conductor_type).unwrap(); + writeln!(out, " ").unwrap(); + } + + // nodes + writeln!(out, "[NumPoints] = {}", self.nodes.len()).unwrap(); + for n in &self.nodes { + let t = lookup_idx(&n.boundary_marker, self.points.iter().map(|p| p.name.as_str())); + let c = lookup_idx(&n.in_conductor, self.conductors.iter().map(|c| c.name.as_str())); + writeln!(out, "{:.17e}\t{:.17e}\t{}\t{}\t{}", n.x, n.y, t, n.in_group, c).unwrap(); + } + + // segments + writeln!(out, "[NumSegments] = {}", self.segments.len()).unwrap(); + for s in &self.segments { + let t = lookup_idx(&s.boundary_marker, self.boundaries.iter().map(|b| b.name.as_str())); + let c = lookup_idx(&s.in_conductor, self.conductors.iter().map(|c| c.name.as_str())); + let msl = if s.max_side_length < 0.0 { + String::from("-1") + } else { + format!("{:.17e}", s.max_side_length) + }; + writeln!(out, "{}\t{}\t{}\t{}\t{}\t{}\t{}", + s.n0, s.n1, msl, t, s.hidden as i32, s.in_group, c).unwrap(); + } + + // arc segments + writeln!(out, "[NumArcSegments] = {}", self.arcs.len()).unwrap(); + for a in &self.arcs { + let t = lookup_idx(&a.boundary_marker, self.boundaries.iter().map(|b| b.name.as_str())); + let c = lookup_idx(&a.in_conductor, self.conductors.iter().map(|c| c.name.as_str())); + writeln!(out, "{}\t{}\t{:.17e}\t{:.17e}\t{}\t{}\t{}\t{}\t{:.17e}", + a.n0, a.n1, a.arc_length, a.max_side_length, t, + a.hidden as i32, a.in_group, c, a.max_side_length).unwrap(); + } + + // holes (block labels typed "") + let holes: Vec<_> = self.block_labels.iter().filter(|b| b.block_type == "").collect(); + writeln!(out, "[NumHoles] = {}", holes.len()).unwrap(); + for h in &holes { + writeln!(out, "{:.17e}\t{:.17e}\t{}", h.x, h.y, h.in_group).unwrap(); + } + + // regional attributes + let labels: Vec<_> = self.block_labels.iter().filter(|b| b.block_type != "").collect(); + writeln!(out, "[NumBlockLabels] = {}", labels.len()).unwrap(); + for b in labels { + let mi = lookup_idx(&b.block_type, self.materials.iter().map(|m| m.name.as_str())); + let area_field = if b.max_area > 0.0 { + format!("{:.17e}", (4.0 * b.max_area / std::f64::consts::PI).sqrt()) + } else { + String::from("-1") + }; + let flags = (b.is_external as i32) + ((b.is_default as i32) << 1); + writeln!(out, "{:.17e}\t{:.17e}\t{}\t{}\t{}\t{}", + b.x, b.y, mi, area_field, b.in_group, flags).unwrap(); + } + + out + } + + /// writes the document to a .fee file by path. + pub fn save(&self, path: impl AsRef) -> std::io::Result<()> { + std::fs::write(path, self.write()) + } +} + +/// finds the 1-based index of `name` in the given name iterator, or 0 when absent or empty. +fn lookup_idx<'a, I: Iterator>(name: &str, iter: I) -> i32 { + if name.is_empty() || name == "" { return 0; } + for (i, n) in iter.enumerate() { + if n == name { return (i as i32) + 1; } + } + 0 +} + +/// inverse of `unescape_comment`: turns CR/LF pairs back into `\n` escapes. +fn escape_comment(s: &str) -> String { + let mut out = String::with_capacity(s.len()); + let bytes = s.as_bytes(); + let mut i = 0; + while i < bytes.len() { + if i + 1 < bytes.len() && bytes[i] == b'\r' && bytes[i + 1] == b'\n' { + out.push('\\'); + out.push('n'); + i += 2; + } else { + out.push(bytes[i] as char); + i += 1; + } + } + out +} diff --git a/crates/femm-doc-elec/tests/roundtrip.rs b/crates/femm-doc-elec/tests/roundtrip.rs new file mode 100644 index 0000000..d1ac174 --- /dev/null +++ b/crates/femm-doc-elec/tests/roundtrip.rs @@ -0,0 +1,137 @@ +use femm_doc_elec::FemmDoc; + +const FIXTURE: &str = r#"[Format] = 1 +[Precision] = 1e-08 +[MinAngle] = 30 +[DoSmartMesh] = 1 +[Depth] = 1 +[LengthUnits] = millimeters +[ProblemType] = planar +[Coordinates] = cartesian +[Comment] = "parallel-plate capacitor" +[PointProps] = 1 + + = "V=0" + = 0 + = 0 + +[BdryProps] = 1 + + = "outer" + = 0 + = 0 + = 0 + = 0 + = 0 + +[BlockProps] = 2 + + = "Air" + = 1 + = 1 + = 0 + + + = "Dielectric" + = 4.2 + = 4.2 + = 0 + +[ConductorProps] = 2 + + = "Plate+" + = 100 + = 0 + = 1 + + + = "Plate-" + = 0 + = 0 + = 1 + +[NumPoints] = 4 +0 0 1 0 1 +10 0 0 0 0 +10 10 0 0 2 +0 10 0 0 0 +[NumSegments] = 4 +0 1 -1 1 0 0 0 +1 2 -1 1 0 0 0 +2 3 -1 1 0 0 0 +3 0 -1 1 0 0 0 +[NumArcSegments] = 0 +[NumHoles] = 0 +[NumBlockLabels] = 1 +5 5 2 0.1 0 0 +"#; + +#[test] +fn parses_fixture_geometry() { + let doc = FemmDoc::parse(FIXTURE).expect("parse"); + assert_eq!(doc.nodes.len(), 4); + assert_eq!(doc.segments.len(), 4); + assert_eq!(doc.arcs.len(), 0); + assert_eq!(doc.block_labels.len(), 1); + assert_eq!(doc.materials.len(), 2); + assert_eq!(doc.boundaries.len(), 1); + assert_eq!(doc.points.len(), 1); + assert_eq!(doc.conductors.len(), 2); + + assert_eq!(doc.nodes[0].boundary_marker, "V=0"); + assert_eq!(doc.nodes[0].in_conductor, "Plate+"); + assert_eq!(doc.nodes[2].in_conductor, "Plate-"); + + assert_eq!(doc.segments[0].boundary_marker, "outer"); + assert_eq!(doc.block_labels[0].block_type, "Dielectric"); + assert_eq!(doc.block_labels[0].in_conductor, ""); + assert_eq!(doc.conductors[0].vc, 100.0); + assert_eq!(doc.conductors[0].conductor_type, 1); + assert_eq!(doc.materials[1].ex, 4.2); + assert_eq!(doc.comment, "parallel-plate capacitor"); +} + +#[test] +fn round_trips_parse_write_parse() { + let a = FemmDoc::parse(FIXTURE).expect("parse a"); + let text = a.write(); + let b = FemmDoc::parse(&text).expect("parse b"); + + assert_eq!(a.precision, b.precision); + assert_eq!(a.depth, b.depth); + assert_eq!(a.nodes.len(), b.nodes.len()); + assert_eq!(a.segments.len(), b.segments.len()); + assert_eq!(a.materials.len(), b.materials.len()); + assert_eq!(a.conductors.len(), b.conductors.len()); + assert_eq!(a.points.len(), b.points.len()); + assert_eq!(a.boundaries.len(), b.boundaries.len()); + + for (x, y) in a.nodes.iter().zip(b.nodes.iter()) { + assert!((x.x - y.x).abs() < 1e-12); + assert!((x.y - y.y).abs() < 1e-12); + assert_eq!(x.boundary_marker, y.boundary_marker); + assert_eq!(x.in_conductor, y.in_conductor); + } + for (x, y) in a.segments.iter().zip(b.segments.iter()) { + assert_eq!(x.n0, y.n0); + assert_eq!(x.n1, y.n1); + assert_eq!(x.boundary_marker, y.boundary_marker); + assert_eq!(x.in_conductor, y.in_conductor); + } + for (x, y) in a.materials.iter().zip(b.materials.iter()) { + assert_eq!(x.name, y.name); + assert!((x.ex - y.ex).abs() < 1e-12); + assert!((x.ey - y.ey).abs() < 1e-12); + } + for (x, y) in a.conductors.iter().zip(b.conductors.iter()) { + assert_eq!(x.name, y.name); + assert!((x.vc - y.vc).abs() < 1e-12); + assert_eq!(x.conductor_type, y.conductor_type); + } + assert_eq!(a.block_labels.len(), b.block_labels.len()); + let la = &a.block_labels[0]; + let lb = &b.block_labels[0]; + assert!((la.x - lb.x).abs() < 1e-12); + assert_eq!(la.block_type, lb.block_type); + assert!((la.max_area - lb.max_area).abs() < 1e-10); +} diff --git a/crates/femm-doc-heat/Cargo.toml b/crates/femm-doc-heat/Cargo.toml new file mode 100644 index 0000000..8baa235 --- /dev/null +++ b/crates/femm-doc-heat/Cargo.toml @@ -0,0 +1,14 @@ +[package] +name = "femm-doc-heat" +version = "0.0.1" +edition.workspace = true +rust-version.workspace = true +publish.workspace = true +description = "rust-native .feh document model: parser, writer, geometry editing" + +[lib] +path = "src/lib.rs" + +[dependencies] +num-complex = "0.4" +thiserror = "2" diff --git a/crates/femm-doc-heat/src/geom.rs b/crates/femm-doc-heat/src/geom.rs new file mode 100644 index 0000000..724b8e5 --- /dev/null +++ b/crates/femm-doc-heat/src/geom.rs @@ -0,0 +1,54 @@ +//! geometric primitives mirroring the .feh control-point model. + +/// control point in the planar geometry. +#[derive(Debug, Clone, Default)] +pub struct Node { + pub x: f64, + pub y: f64, + pub boundary_marker: String, + pub in_conductor: String, + pub in_group: i32, + pub selected: bool, +} + +/// straight segment joining two control points by index. +#[derive(Debug, Clone, Default)] +pub struct Segment { + pub n0: i32, + pub n1: i32, + pub max_side_length: f64, + pub boundary_marker: String, + pub in_conductor: String, + pub hidden: bool, + pub in_group: i32, + pub selected: bool, +} + +/// arc segment from n0 to n1, sweeping the given degree arc. +#[derive(Debug, Clone, Default)] +pub struct ArcSegment { + pub n0: i32, + pub n1: i32, + pub arc_length: f64, + pub max_side_length: f64, + pub boundary_marker: String, + pub in_conductor: String, + pub hidden: bool, + pub in_group: i32, + pub normal_direction: bool, + pub selected: bool, +} + +/// region label assigning a material and mesh-area constraint to an enclosed area. +#[derive(Debug, Clone, Default)] +pub struct BlockLabel { + pub x: f64, + pub y: f64, + pub max_area: f64, + pub block_type: String, + pub in_conductor: String, + pub in_group: i32, + pub is_external: bool, + pub is_default: bool, + pub selected: bool, +} diff --git a/crates/femm-doc-heat/src/lib.rs b/crates/femm-doc-heat/src/lib.rs new file mode 100644 index 0000000..2fa5ab2 --- /dev/null +++ b/crates/femm-doc-heat/src/lib.rs @@ -0,0 +1,75 @@ +//! heat-flow pre-processor document model: geometry, properties, parser, writer. + +pub mod geom; +pub mod props; +pub mod parser; +pub mod writer; + +use num_complex::Complex64; + +pub use geom::{ArcSegment, BlockLabel, Node, Segment}; +pub use props::{BoundaryProp, ConductorProp, MaterialProp, PointProp}; + +#[derive(Debug, Clone, Copy, PartialEq, Eq)] +pub enum ProblemType { + Planar, + Axisymmetric, +} + +impl Default for ProblemType { + fn default() -> Self { ProblemType::Planar } +} + +#[derive(Debug, Clone, Copy, PartialEq, Eq, Default)] +pub enum Coords { + #[default] + Cartesian, + Polar, +} + +#[derive(Debug, Clone, Copy, PartialEq, Eq)] +pub enum LengthUnit { + Inches = 0, + Millimeters = 1, + Centimeters = 2, + Meters = 3, + Mils = 4, + Microns = 5, +} + +impl Default for LengthUnit { + fn default() -> Self { LengthUnit::Inches } +} + +/// heat-flow .feh document: geometry, property tables, problem attributes. +#[derive(Debug, Clone, Default)] +pub struct FemmDoc { + pub format: f64, + pub precision: f64, + pub min_angle: f64, + pub smart_mesh: bool, + pub depth: f64, + pub length_units: LengthUnit, + pub problem_type: ProblemType, + pub coords: Coords, + pub comment: String, + pub prev_soln: String, + pub dt: f64, + + pub ext_ro: f64, + pub ext_ri: f64, + pub ext_zo: f64, + + pub nodes: Vec, + pub segments: Vec, + pub arcs: Vec, + pub block_labels: Vec, + + pub materials: Vec, + pub boundaries: Vec, + pub points: Vec, + pub conductors: Vec, +} + +/// f64 complex number alias matching the C++ CComplex memory layout. +pub type Complex = Complex64; diff --git a/crates/femm-doc-heat/src/parser.rs b/crates/femm-doc-heat/src/parser.rs new file mode 100644 index 0000000..4a0cf95 --- /dev/null +++ b/crates/femm-doc-heat/src/parser.rs @@ -0,0 +1,414 @@ +//! .feh text parser: `[Section] = value` scalars and `...` property stanzas. + +use crate::{ + ArcSegment, BlockLabel, BoundaryProp, ConductorProp, Coords, FemmDoc, LengthUnit, + MaterialProp, Node, PointProp, ProblemType, Segment, +}; +use std::path::Path; +use thiserror::Error; + +#[derive(Debug, Error)] +pub enum ParseError { + #[error("io error: {0}")] + Io(#[from] std::io::Error), + #[error("legacy 3.2-format .feh files are not supported")] + LegacyFormat, + #[error("unexpected end of file while reading {context}")] + UnexpectedEof { context: &'static str }, + #[error("malformed value in section {section}: {reason}")] + MalformedValue { section: &'static str, reason: String }, +} + +impl FemmDoc { + /// parses a .feh text buffer into a document. + pub fn parse(src: &str) -> Result { + let mut doc = FemmDoc::default(); + doc.depth = -1.0; + doc.smart_mesh = true; + doc.format = 1.0; + + if src + .lines() + .next() + .map(|l| l.trim_start().starts_with("Frequency")) + .unwrap_or(false) + { + return Err(ParseError::LegacyFormat); + } + + let mut lines = src.lines(); + + let mut point = PointProp::default(); + let mut bdry = BoundaryProp::default(); + let mut mat = MaterialProp::default(); + let mut cond = ConductorProp::default(); + + while let Some(raw) = lines.next() { + let tok = first_token(raw); + if tok.is_empty() { continue; } + + // problem attributes + if tok.eq_ignore_ascii_case("[format]") { + doc.format = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[depth]") { + doc.depth = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[precision]") { + doc.precision = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[minangle]") { + doc.min_angle = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[dosmartmesh]") { + doc.smart_mesh = parse_i32(strip_key(raw)) != 0; + } else if tok.eq_ignore_ascii_case("[lengthunits]") { + doc.length_units = parse_length_units(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[problemtype]") { + doc.problem_type = parse_problem_type(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[coordinates]") { + doc.coords = parse_coords(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[comment]") { + doc.comment = unescape_comment(unquote(strip_key(raw))); + } else if tok.eq_ignore_ascii_case("[prevsoln]") { + doc.prev_soln = unquote(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[dt]") { + doc.dt = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[extzo]") { + doc.ext_zo = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[extro]") { + doc.ext_ro = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[extri]") { + doc.ext_ri = parse_f64(strip_key(raw)); + + // point property stanza + } else if tok.eq_ignore_ascii_case("") { + point = PointProp::default(); + point.name = String::from("New Point Property"); + } else if tok.eq_ignore_ascii_case("") { + point.name = unquote(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + point.tp = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + point.qp = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + doc.points.push(std::mem::take(&mut point)); + + // boundary property stanza + } else if tok.eq_ignore_ascii_case("") { + bdry = BoundaryProp::default(); + bdry.name = String::from("New Boundary"); + } else if tok.eq_ignore_ascii_case("") { + bdry.name = unquote(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.format = parse_i32(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.tset = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.qs = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.beta = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.h = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.tinf = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + doc.boundaries.push(std::mem::take(&mut bdry)); + + // material property stanza + } else if tok.eq_ignore_ascii_case("") { + mat = MaterialProp::default(); + mat.name = String::from("New Material"); + mat.kx = 1.0; + mat.ky = 1.0; + mat.kt = 3.0; + mat.qv = 0.0; + } else if tok.eq_ignore_ascii_case("") { + mat.name = unquote(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + mat.kx = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + mat.ky = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + mat.kt = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + mat.qv = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + let count = parse_i32(strip_key(raw)); + if count > 0 { + mat.tk_curve.clear(); + mat.tk_curve.reserve(count as usize); + for _ in 0..count { + let line = lines.next().ok_or(ParseError::UnexpectedEof { + context: "", + })?; + let (t, rest) = take_f64(line); + let (k, _) = take_f64(rest); + mat.tk_curve.push((t, k)); + } + } + } else if tok.eq_ignore_ascii_case("") { + doc.materials.push(std::mem::take(&mut mat)); + + // conductor stanza + } else if tok.eq_ignore_ascii_case("") { + cond = ConductorProp::default(); + cond.name = String::from("New Conductor"); + } else if tok.eq_ignore_ascii_case("") { + cond.name = unquote(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + cond.tc = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + cond.qc = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + cond.conductor_type = parse_i32(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + doc.conductors.push(std::mem::take(&mut cond)); + + // counted geometry sections + } else if tok.eq_ignore_ascii_case("[numpoints]") { + let n = parse_i32(strip_key(raw)); + for _ in 0..n { + let line = lines.next().ok_or(ParseError::UnexpectedEof { context: "[NumPoints]" })?; + let (x, rest) = take_f64(line); + let (y, rest) = take_f64(rest); + let (bi, rest) = take_i32(rest); + let (group, rest) = take_i32(rest); + let (ci, _) = take_i32(rest); + let bm = resolve_name(bi, doc.points.iter().map(|p| p.name.as_str())); + let ic = resolve_conductor(ci, doc.conductors.iter().map(|c| c.name.as_str())); + doc.nodes.push(Node { + x, y, + boundary_marker: bm, + in_conductor: ic, + in_group: group, + selected: false, + }); + } + } else if tok.eq_ignore_ascii_case("[numsegments]") { + let n = parse_i32(strip_key(raw)); + for _ in 0..n { + let line = lines.next().ok_or(ParseError::UnexpectedEof { context: "[NumSegments]" })?; + let (n0, rest) = take_i32(line); + let (n1, rest) = take_i32(rest); + let (msl, rest) = take_f64(rest); + let (bi, rest) = take_i32(rest); + let (hidden, rest) = take_i32(rest); + let (group, rest) = take_i32(rest); + let (ci, _) = take_i32(rest); + let bm = resolve_name(bi, doc.boundaries.iter().map(|b| b.name.as_str())); + let ic = resolve_conductor(ci, doc.conductors.iter().map(|c| c.name.as_str())); + doc.segments.push(Segment { + n0, n1, + max_side_length: msl, + boundary_marker: bm, + in_conductor: ic, + hidden: hidden != 0, + in_group: group, + selected: false, + }); + } + } else if tok.eq_ignore_ascii_case("[numarcsegments]") { + let n = parse_i32(strip_key(raw)); + for _ in 0..n { + let line = lines.next().ok_or(ParseError::UnexpectedEof { context: "[NumArcSegments]" })?; + let (n0, rest) = take_i32(line); + let (n1, rest) = take_i32(rest); + let (al, rest) = take_f64(rest); + let (msl, rest) = take_f64(rest); + let (bi, rest) = take_i32(rest); + let (hidden, rest) = take_i32(rest); + let (group, rest) = take_i32(rest); + let (ci, _) = take_i32(rest); + let bm = resolve_name(bi, doc.boundaries.iter().map(|b| b.name.as_str())); + let ic = resolve_conductor(ci, doc.conductors.iter().map(|c| c.name.as_str())); + doc.arcs.push(ArcSegment { + n0, n1, + arc_length: al, + max_side_length: msl, + boundary_marker: bm, + in_conductor: ic, + hidden: hidden != 0, + in_group: group, + normal_direction: true, + selected: false, + }); + } + } else if tok.eq_ignore_ascii_case("[numholes]") { + let n = parse_i32(strip_key(raw)); + for _ in 0..n { + let line = lines.next().ok_or(ParseError::UnexpectedEof { context: "[NumHoles]" })?; + let (x, rest) = take_f64(line); + let (y, rest) = take_f64(rest); + let (group, _) = take_i32(rest); + doc.block_labels.push(BlockLabel { + x, y, + max_area: 0.0, + block_type: String::from(""), + in_conductor: String::from(""), + in_group: group, + is_external: false, + is_default: false, + selected: false, + }); + } + } else if tok.eq_ignore_ascii_case("[numblocklabels]") { + let n = parse_i32(strip_key(raw)); + for _ in 0..n { + let line = lines.next().ok_or(ParseError::UnexpectedEof { context: "[NumBlockLabels]" })?; + let (x, rest) = take_f64(line); + let (y, rest) = take_f64(rest); + let (mi, rest) = take_i32(rest); + let (ma_diam, rest) = take_f64(rest); + let (group, rest) = take_i32(rest); + let (ext_flags, _) = take_i32(rest); + + let block_type = if mi == 0 { + String::from("") + } else { + resolve_name(mi, doc.materials.iter().map(|m| m.name.as_str())) + }; + let max_area = if ma_diam < 0.0 { + 0.0 + } else { + std::f64::consts::PI * ma_diam * ma_diam / 4.0 + }; + + doc.block_labels.push(BlockLabel { + x, y, + max_area, + block_type, + in_conductor: String::from(""), + in_group: group, + is_external: ext_flags & 1 != 0, + is_default: ext_flags & 2 != 0, + selected: false, + }); + } + } + // unknown tokens are skipped silently, matching the original parser. + } + + // 3.2-era files omitted [Depth]; fill the default in the current length unit. + if doc.depth == -1.0 { + doc.depth = match doc.length_units { + LengthUnit::Millimeters => 1000.0, + LengthUnit::Centimeters => 100.0, + LengthUnit::Meters => 1.0, + LengthUnit::Mils => 1000.0 / 0.0254, + LengthUnit::Microns => 1.0e6, + LengthUnit::Inches => 1.0 / 0.0254, + }; + } + + Ok(doc) + } + + /// loads a .feh file by path. + pub fn open(path: impl AsRef) -> Result { + let text = std::fs::read_to_string(path)?; + Self::parse(&text) + } +} + +fn first_token(line: &str) -> &str { + line.trim_start().split_whitespace().next().unwrap_or("") +} + +fn strip_key(line: &str) -> &str { + match line.find('=') { + Some(i) => &line[i + 1..], + None => "", + } +} + +fn unquote(s: &str) -> String { + let s = s.trim(); + let first = s.find('"'); + let last = s.rfind('"'); + match (first, last) { + (Some(a), Some(b)) if a < b => s[a + 1..b].to_string(), + _ => s.to_string(), + } +} + +fn unescape_comment(s: String) -> String { + let bytes = s.as_bytes(); + let mut out = String::with_capacity(s.len()); + let mut i = 0; + while i < bytes.len() { + if bytes[i] == b'\\' && i + 1 < bytes.len() && bytes[i + 1] == b'n' { + out.push('\r'); + out.push('\n'); + i += 2; + } else { + out.push(bytes[i] as char); + i += 1; + } + } + out +} + +fn parse_f64(s: &str) -> f64 { + let s = s.trim().trim_matches(|c: char| c == ',' || c.is_whitespace()); + let end = s.find(|c: char| c.is_whitespace() || c == ',').unwrap_or(s.len()); + s[..end].parse::().unwrap_or(0.0) +} + +fn parse_i32(s: &str) -> i32 { + let s = s.trim().trim_matches(|c: char| c == ',' || c.is_whitespace()); + let end = s.find(|c: char| c.is_whitespace() || c == ',').unwrap_or(s.len()); + s[..end].parse::().unwrap_or(0) +} + +/// peels one f64 off the head of a whitespace/comma-separated row, returning the remainder. +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/comma-separated row, returning the remainder. +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) +} + +fn parse_length_units(s: &str) -> LengthUnit { + let t = first_token(s).to_ascii_lowercase(); + if t.starts_with("millimeters") { LengthUnit::Millimeters } + else if t.starts_with("c") { LengthUnit::Centimeters } + else if t.starts_with("meters") { LengthUnit::Meters } + else if t.starts_with("mils") { LengthUnit::Mils } + else if t.starts_with("microns") { LengthUnit::Microns } + else { LengthUnit::Inches } +} + +fn parse_problem_type(s: &str) -> ProblemType { + let t = first_token(s).to_ascii_lowercase(); + if t.starts_with("axi") { ProblemType::Axisymmetric } else { ProblemType::Planar } +} + +fn parse_coords(s: &str) -> Coords { + let t = first_token(s).to_ascii_lowercase(); + if t.starts_with("polar") { Coords::Polar } else { Coords::Cartesian } +} + +/// resolves a 1-based property index into the matching property name, falling back to empty. +fn resolve_name<'a, I: Iterator>(idx: i32, names: I) -> String { + if idx <= 0 { return String::new(); } + let want = idx as usize; + for (i, name) in names.enumerate() { + if i + 1 == want { return name.to_string(); } + } + String::new() +} + +/// resolves a 1-based conductor index into a conductor name, falling back to "". +fn resolve_conductor<'a, I: Iterator>(idx: i32, names: I) -> String { + if idx <= 0 { return String::from(""); } + let want = idx as usize; + for (i, name) in names.enumerate() { + if i + 1 == want { return name.to_string(); } + } + String::from("") +} diff --git a/crates/femm-doc-heat/src/props.rs b/crates/femm-doc-heat/src/props.rs new file mode 100644 index 0000000..2eb1bcc --- /dev/null +++ b/crates/femm-doc-heat/src/props.rs @@ -0,0 +1,41 @@ +//! heat-flow material, boundary, point, and conductor properties. + +/// material property entry referenced by name from block labels. +#[derive(Debug, Clone, Default)] +pub struct MaterialProp { + pub name: String, + pub kx: f64, + pub ky: f64, + pub kt: f64, + pub qv: f64, + pub tk_curve: Vec<(f64, f64)>, +} + +/// boundary condition entry referenced by name from segments and arcs. +#[derive(Debug, Clone, Default)] +pub struct BoundaryProp { + pub name: String, + pub format: i32, + pub tset: f64, + pub qs: f64, + pub beta: f64, + pub h: f64, + pub tinf: f64, +} + +/// point property entry referenced by name from nodes. +#[derive(Debug, Clone, Default)] +pub struct PointProp { + pub name: String, + pub tp: f64, + pub qp: f64, +} + +/// conductor property entry referenced by name from nodes, segments, and arcs. +#[derive(Debug, Clone, Default)] +pub struct ConductorProp { + pub name: String, + pub tc: f64, + pub qc: f64, + pub conductor_type: i32, +} diff --git a/crates/femm-doc-heat/src/writer.rs b/crates/femm-doc-heat/src/writer.rs new file mode 100644 index 0000000..69406b7 --- /dev/null +++ b/crates/femm-doc-heat/src/writer.rs @@ -0,0 +1,199 @@ +//! .feh text writer, inverse of [`parser`](crate::parser). + +use crate::{Coords, FemmDoc, LengthUnit, ProblemType}; +use std::fmt::Write; +use std::path::Path; + +impl FemmDoc { + /// renders the document to .feh text. + pub fn write(&self) -> String { + let mut out = String::new(); + writeln!(out, "[Format] = 1").unwrap(); + writeln!(out, "[Precision] = {:.17e}", self.precision).unwrap(); + writeln!(out, "[MinAngle] = {:.17e}", self.min_angle).unwrap(); + writeln!(out, "[DoSmartMesh] = {}", self.smart_mesh as i32).unwrap(); + writeln!(out, "[Depth] = {:.17e}", self.depth).unwrap(); + + let units = match self.length_units { + LengthUnit::Inches => "inches", + LengthUnit::Millimeters => "millimeters", + LengthUnit::Centimeters => "centimeters", + LengthUnit::Meters => "meters", + LengthUnit::Mils => "mils", + LengthUnit::Microns => "microns", + }; + writeln!(out, "[LengthUnits] = {units}").unwrap(); + + match self.problem_type { + ProblemType::Planar => writeln!(out, "[ProblemType] = planar").unwrap(), + ProblemType::Axisymmetric => { + writeln!(out, "[ProblemType] = axisymmetric").unwrap(); + if self.ext_ro != 0.0 && self.ext_ri != 0.0 { + writeln!(out, "[extZo] = {:.17e}", self.ext_zo).unwrap(); + writeln!(out, "[extRo] = {:.17e}", self.ext_ro).unwrap(); + writeln!(out, "[extRi] = {:.17e}", self.ext_ri).unwrap(); + } + } + } + + let coords = match self.coords { + Coords::Cartesian => "cartesian", + Coords::Polar => "polar", + }; + writeln!(out, "[Coordinates] = {coords}").unwrap(); + + writeln!(out, "[PrevSoln] = \"{}\"", self.prev_soln).unwrap(); + writeln!(out, "[dT] = {:.17e}", self.dt).unwrap(); + writeln!(out, "[Comment] = \"{}\"", escape_comment(&self.comment)).unwrap(); + + // point properties + writeln!(out, "[PointProps] = {}", self.points.len()).unwrap(); + for p in &self.points { + writeln!(out, " ").unwrap(); + writeln!(out, " = \"{}\"", p.name).unwrap(); + writeln!(out, " = {:.17e}", p.tp).unwrap(); + writeln!(out, " = {:.17e}", p.qp).unwrap(); + writeln!(out, " ").unwrap(); + } + + // boundary properties + writeln!(out, "[BdryProps] = {}", self.boundaries.len()).unwrap(); + for b in &self.boundaries { + writeln!(out, " ").unwrap(); + writeln!(out, " = \"{}\"", b.name).unwrap(); + writeln!(out, " = {}", b.format).unwrap(); + writeln!(out, " = {:.17e}", b.tset).unwrap(); + writeln!(out, " = {:.17e}", b.qs).unwrap(); + writeln!(out, " = {:.17e}", b.beta).unwrap(); + writeln!(out, " = {:.17e}", b.h).unwrap(); + writeln!(out, " = {:.17e}", b.tinf).unwrap(); + writeln!(out, " ").unwrap(); + } + + // material properties + writeln!(out, "[BlockProps] = {}", self.materials.len()).unwrap(); + for m in &self.materials { + writeln!(out, " ").unwrap(); + writeln!(out, " = \"{}\"", m.name).unwrap(); + writeln!(out, " = {:.17e}", m.kx).unwrap(); + writeln!(out, " = {:.17e}", m.ky).unwrap(); + writeln!(out, " = {:.17e}", m.kt).unwrap(); + writeln!(out, " = {:.17e}", m.qv).unwrap(); + if !m.tk_curve.is_empty() { + writeln!(out, " = {}", m.tk_curve.len()).unwrap(); + for (t, k) in &m.tk_curve { + writeln!(out, " {:.17e}\t{:.17e}", t, k).unwrap(); + } + } + writeln!(out, " ").unwrap(); + } + + // conductor properties + writeln!(out, "[ConductorProps] = {}", self.conductors.len()).unwrap(); + for c in &self.conductors { + writeln!(out, " ").unwrap(); + writeln!(out, " = \"{}\"", c.name).unwrap(); + writeln!(out, " = {:.17e}", c.tc).unwrap(); + writeln!(out, " = {:.17e}", c.qc).unwrap(); + writeln!(out, " = {}", c.conductor_type).unwrap(); + writeln!(out, " ").unwrap(); + } + + // nodes + writeln!(out, "[NumPoints] = {}", self.nodes.len()).unwrap(); + for n in &self.nodes { + let t = lookup_idx(&n.boundary_marker, self.points.iter().map(|p| p.name.as_str())); + let c = lookup_conductor(&n.in_conductor, self.conductors.iter().map(|c| c.name.as_str())); + writeln!(out, "{:.17e}\t{:.17e}\t{}\t{}\t{}", n.x, n.y, t, n.in_group, c).unwrap(); + } + + // segments + writeln!(out, "[NumSegments] = {}", self.segments.len()).unwrap(); + for s in &self.segments { + let t = lookup_idx(&s.boundary_marker, self.boundaries.iter().map(|b| b.name.as_str())); + let c = lookup_conductor(&s.in_conductor, self.conductors.iter().map(|c| c.name.as_str())); + let msl = if s.max_side_length < 0.0 { + String::from("-1") + } else { + format!("{:.17e}", s.max_side_length) + }; + writeln!(out, "{}\t{}\t{}\t{}\t{}\t{}\t{}", + s.n0, s.n1, msl, t, s.hidden as i32, s.in_group, c).unwrap(); + } + + // arc segments + writeln!(out, "[NumArcSegments] = {}", self.arcs.len()).unwrap(); + for a in &self.arcs { + let t = lookup_idx(&a.boundary_marker, self.boundaries.iter().map(|b| b.name.as_str())); + let c = lookup_conductor(&a.in_conductor, self.conductors.iter().map(|c| c.name.as_str())); + writeln!(out, "{}\t{}\t{:.17e}\t{:.17e}\t{}\t{}\t{}\t{}\t{:.17e}", + a.n0, a.n1, a.arc_length, a.max_side_length, t, + a.hidden as i32, a.in_group, c, a.max_side_length).unwrap(); + } + + // holes (block labels typed "") + let holes: Vec<_> = self.block_labels.iter().filter(|b| b.block_type == "").collect(); + writeln!(out, "[NumHoles] = {}", holes.len()).unwrap(); + for h in &holes { + writeln!(out, "{:.17e}\t{:.17e}\t{}", h.x, h.y, h.in_group).unwrap(); + } + + // regional attributes + let labels: Vec<_> = self.block_labels.iter().filter(|b| b.block_type != "").collect(); + writeln!(out, "[NumBlockLabels] = {}", labels.len()).unwrap(); + for b in labels { + let mi = lookup_idx(&b.block_type, self.materials.iter().map(|m| m.name.as_str())); + let area_field = if b.max_area > 0.0 { + format!("{:.17e}", (4.0 * b.max_area / std::f64::consts::PI).sqrt()) + } else { + String::from("-1") + }; + let flags = (b.is_external as i32) + ((b.is_default as i32) << 1); + writeln!(out, "{:.17e}\t{:.17e}\t{}\t{}\t{}\t{}", + b.x, b.y, mi, area_field, b.in_group, flags).unwrap(); + } + + out + } + + /// writes the document to a .feh file by path. + pub fn save(&self, path: impl AsRef) -> std::io::Result<()> { + std::fs::write(path, self.write()) + } +} + +/// finds the 1-based index of `name` in the given name iterator, or 0 when absent or empty. +fn lookup_idx<'a, I: Iterator>(name: &str, iter: I) -> i32 { + if name.is_empty() { return 0; } + for (i, n) in iter.enumerate() { + if n == name { return (i as i32) + 1; } + } + 0 +} + +/// finds the 1-based index of a conductor name, treating "" or empty as 0. +fn lookup_conductor<'a, I: Iterator>(name: &str, iter: I) -> i32 { + if name.is_empty() || name == "" { return 0; } + for (i, n) in iter.enumerate() { + if n == name { return (i as i32) + 1; } + } + 0 +} + +/// inverse of `unescape_comment`: turns CR/LF pairs back into `\n` escapes. +fn escape_comment(s: &str) -> String { + let mut out = String::with_capacity(s.len()); + let bytes = s.as_bytes(); + let mut i = 0; + while i < bytes.len() { + if i + 1 < bytes.len() && bytes[i] == b'\r' && bytes[i + 1] == b'\n' { + out.push('\\'); + out.push('n'); + i += 2; + } else { + out.push(bytes[i] as char); + i += 1; + } + } + out +} diff --git a/crates/femm-doc-heat/tests/roundtrip.rs b/crates/femm-doc-heat/tests/roundtrip.rs new file mode 100644 index 0000000..79d6392 --- /dev/null +++ b/crates/femm-doc-heat/tests/roundtrip.rs @@ -0,0 +1,149 @@ +use femm_doc_heat::FemmDoc; + +const FIXTURE: &str = r#"[Format] = 1 +[Precision] = 1e-08 +[MinAngle] = 30 +[DoSmartMesh] = 1 +[Depth] = 1 +[LengthUnits] = millimeters +[ProblemType] = planar +[Coordinates] = cartesian +[PrevSoln] = "" +[dT] = 0.01 +[Comment] = "heat-flow planar slab" +[PointProps] = 1 + + = "Tfix" + = 300 + = 0 + +[BdryProps] = 1 + + = "ambient" + = 2 + = 0 + = 0 + = 0 + = 10 + = 293 + +[BlockProps] = 2 + + = "Air" + = 0.026 + = 0.026 + = 1006 + = 0 + + + = "Copper" + = 401 + = 401 + = 385 + = 0 + = 3 + 250 410 + 300 401 + 400 390 + +[ConductorProps] = 1 + + = "Heater" + = 0 + = 50 + = 0 + +[NumPoints] = 4 +0 0 1 0 1 +10 0 0 0 0 +10 10 0 0 0 +0 10 0 0 0 +[NumSegments] = 4 +0 1 -1 1 0 0 0 +1 2 -1 1 0 0 0 +2 3 -1 1 0 0 0 +3 0 -1 1 0 0 0 +[NumArcSegments] = 0 +[NumHoles] = 0 +[NumBlockLabels] = 1 +5 5 2 0.1 0 0 +"#; + +#[test] +fn parses_fixture_geometry() { + let doc = FemmDoc::parse(FIXTURE).expect("parse"); + assert_eq!(doc.dt, 0.01); + assert_eq!(doc.nodes.len(), 4); + assert_eq!(doc.segments.len(), 4); + assert_eq!(doc.arcs.len(), 0); + assert_eq!(doc.block_labels.len(), 1); + assert_eq!(doc.materials.len(), 2); + assert_eq!(doc.boundaries.len(), 1); + assert_eq!(doc.points.len(), 1); + assert_eq!(doc.conductors.len(), 1); + assert_eq!(doc.nodes[0].boundary_marker, "Tfix"); + assert_eq!(doc.nodes[0].in_conductor, "Heater"); + assert_eq!(doc.nodes[1].in_conductor, ""); + assert_eq!(doc.segments[0].boundary_marker, "ambient"); + assert_eq!(doc.block_labels[0].block_type, "Copper"); + assert_eq!(doc.materials[1].name, "Copper"); + assert_eq!(doc.materials[1].tk_curve.len(), 3); + assert_eq!(doc.materials[1].tk_curve[1], (300.0, 401.0)); + assert_eq!(doc.boundaries[0].h, 10.0); + assert_eq!(doc.boundaries[0].tinf, 293.0); + assert_eq!(doc.points[0].tp, 300.0); + assert_eq!(doc.conductors[0].qc, 50.0); +} + +#[test] +fn round_trips_parse_write_parse() { + let a = FemmDoc::parse(FIXTURE).expect("parse a"); + let text = a.write(); + let b = FemmDoc::parse(&text).expect("parse b"); + + assert_eq!(a.dt, b.dt); + assert_eq!(a.precision, b.precision); + assert_eq!(a.depth, b.depth); + assert_eq!(a.nodes.len(), b.nodes.len()); + assert_eq!(a.segments.len(), b.segments.len()); + assert_eq!(a.materials.len(), b.materials.len()); + assert_eq!(a.conductors.len(), b.conductors.len()); + + for (x, y) in a.nodes.iter().zip(b.nodes.iter()) { + assert!((x.x - y.x).abs() < 1e-12); + assert!((x.y - y.y).abs() < 1e-12); + assert_eq!(x.boundary_marker, y.boundary_marker); + assert_eq!(x.in_conductor, y.in_conductor); + } + for (x, y) in a.segments.iter().zip(b.segments.iter()) { + assert_eq!(x.n0, y.n0); + assert_eq!(x.n1, y.n1); + assert_eq!(x.boundary_marker, y.boundary_marker); + } + for (x, y) in a.materials.iter().zip(b.materials.iter()) { + assert_eq!(x.name, y.name); + assert!((x.kx - y.kx).abs() < 1e-12); + assert!((x.kt - y.kt).abs() < 1e-12); + assert_eq!(x.tk_curve.len(), y.tk_curve.len()); + for ((t0, k0), (t1, k1)) in x.tk_curve.iter().zip(y.tk_curve.iter()) { + assert!((t0 - t1).abs() < 1e-12); + assert!((k0 - k1).abs() < 1e-12); + } + } + for (x, y) in a.boundaries.iter().zip(b.boundaries.iter()) { + assert_eq!(x.name, y.name); + assert!((x.h - y.h).abs() < 1e-12); + assert!((x.tinf - y.tinf).abs() < 1e-12); + } + for (x, y) in a.conductors.iter().zip(b.conductors.iter()) { + assert_eq!(x.name, y.name); + assert!((x.qc - y.qc).abs() < 1e-12); + assert_eq!(x.conductor_type, y.conductor_type); + } + assert_eq!(a.block_labels.len(), b.block_labels.len()); + let la = &a.block_labels[0]; + let lb = &b.block_labels[0]; + assert!((la.x - lb.x).abs() < 1e-12); + assert_eq!(la.block_type, lb.block_type); + assert!((la.max_area - lb.max_area).abs() < 1e-10); +} diff --git a/crates/femm-doc/Cargo.toml b/crates/femm-doc/Cargo.toml new file mode 100644 index 0000000..1e2225e --- /dev/null +++ b/crates/femm-doc/Cargo.toml @@ -0,0 +1,14 @@ +[package] +name = "femm-doc" +version = "0.0.1" +edition.workspace = true +rust-version.workspace = true +publish.workspace = true +description = "rust-native .fem document model: parser, writer, geometry editing" + +[lib] +path = "src/lib.rs" + +[dependencies] +num-complex = "0.4" +thiserror = "2" diff --git a/crates/femm-doc/src/edit.rs b/crates/femm-doc/src/edit.rs new file mode 100644 index 0000000..59174a8 --- /dev/null +++ b/crates/femm-doc/src/edit.rs @@ -0,0 +1,212 @@ +//! geometry editing primitives on [`FemmDoc`]: add, delete, closest-point queries. + +use crate::{ArcSegment, BlockLabel, FemmDoc, Node, Segment}; + +impl FemmDoc { + /// adds a node at (x, y), returning the index of an existing node within `tol` distance when present. + pub fn add_node(&mut self, x: f64, y: f64, tol: f64) -> usize { + if tol > 0.0 { + for (i, n) in self.nodes.iter().enumerate() { + let dx = n.x - x; + let dy = n.y - y; + if (dx * dx + dy * dy).sqrt() <= tol { + return i; + } + } + } + self.nodes.push(Node { + x, y, + boundary_marker: String::new(), + in_group: 0, + selected: false, + }); + self.nodes.len() - 1 + } + + /// adds a block label at (x, y), returning the index of an existing label within `tol` when present. + pub fn add_block_label(&mut self, x: f64, y: f64, tol: f64) -> usize { + if tol > 0.0 { + for (i, b) in self.block_labels.iter().enumerate() { + let dx = b.x - x; + let dy = b.y - y; + if (dx * dx + dy * dy).sqrt() <= tol { + return i; + } + } + } + self.block_labels.push(BlockLabel { + x, y, + max_area: 0.0, + mag_dir: 0.0, + mag_dir_fctn: String::new(), + turns: 1, + block_type: String::from(""), + in_circuit: String::from(""), + in_group: 0, + is_external: false, + is_default: false, + selected: false, + }); + self.block_labels.len() - 1 + } + + /// adds a segment between two node indices, or returns the index of an existing duplicate. + pub fn add_segment(&mut self, n0: i32, n1: i32) -> Option { + if n0 == n1 { return None; } + let (a, b) = if n0 < n1 { (n0, n1) } else { (n1, n0) }; + for (i, s) in self.segments.iter().enumerate() { + let (sa, sb) = if s.n0 < s.n1 { (s.n0, s.n1) } else { (s.n1, s.n0) }; + if sa == a && sb == b { return Some(i); } + } + self.segments.push(Segment { + n0, n1, + max_side_length: -1.0, + boundary_marker: String::new(), + hidden: false, + in_group: 0, + selected: false, + }); + Some(self.segments.len() - 1) + } + + /// adds an arc segment between two node indices, sweeping `arc_length_deg` degrees. + pub fn add_arc_segment(&mut self, n0: i32, n1: i32, arc_length_deg: f64) -> Option { + if n0 == n1 { return None; } + self.arcs.push(ArcSegment { + n0, n1, + arc_length: arc_length_deg, + max_side_length: 10.0, + boundary_marker: String::new(), + hidden: false, + in_group: 0, + normal_direction: true, + selected: false, + }); + Some(self.arcs.len() - 1) + } + + /// removes selected nodes and rewrites segment/arc endpoint indices to drop references. + pub fn delete_selected_nodes(&mut self) -> usize { + let keep: Vec = self.nodes.iter().map(|n| !n.selected).collect(); + let mut remap: Vec = Vec::with_capacity(keep.len()); + let mut next: i32 = 0; + for &k in &keep { + remap.push(if k { let r = next; next += 1; r } else { -1 }); + } + let removed = self.nodes.len() - next as usize; + + let mut new_nodes = Vec::with_capacity(next as usize); + for (i, n) in self.nodes.drain(..).enumerate() { + if keep[i] { new_nodes.push(n); } + } + self.nodes = new_nodes; + + self.segments.retain(|s| { + let a = s.n0 as usize; + let b = s.n1 as usize; + a < keep.len() && b < keep.len() && keep[a] && keep[b] + }); + for s in &mut self.segments { + s.n0 = remap[s.n0 as usize]; + s.n1 = remap[s.n1 as usize]; + } + + self.arcs.retain(|a| { + let i = a.n0 as usize; + let j = a.n1 as usize; + i < keep.len() && j < keep.len() && keep[i] && keep[j] + }); + for a in &mut self.arcs { + a.n0 = remap[a.n0 as usize]; + a.n1 = remap[a.n1 as usize]; + } + + removed + } + + /// removes selected segments and returns the count. + pub fn delete_selected_segments(&mut self) -> usize { + let before = self.segments.len(); + self.segments.retain(|s| !s.selected); + before - self.segments.len() + } + + /// removes selected arc segments and returns the count. + pub fn delete_selected_arcs(&mut self) -> usize { + let before = self.arcs.len(); + self.arcs.retain(|a| !a.selected); + before - self.arcs.len() + } + + /// removes selected block labels and returns the count. + pub fn delete_selected_block_labels(&mut self) -> usize { + let before = self.block_labels.len(); + self.block_labels.retain(|b| !b.selected); + before - self.block_labels.len() + } + + /// returns the index of the closest node to (x, y), or None when the node list is empty. + pub fn closest_node(&self, x: f64, y: f64) -> Option { + let mut best: Option<(usize, f64)> = None; + for (i, n) in self.nodes.iter().enumerate() { + let d = (n.x - x).hypot(n.y - y); + match best { + None => best = Some((i, d)), + Some((_, bd)) if d < bd => best = Some((i, d)), + _ => {} + } + } + best.map(|(i, _)| i) + } + + /// returns the index of the closest block label to (x, y), or None when none exist. + pub fn closest_block_label(&self, x: f64, y: f64) -> Option { + let mut best: Option<(usize, f64)> = None; + for (i, b) in self.block_labels.iter().enumerate() { + let d = (b.x - x).hypot(b.y - y); + match best { + None => best = Some((i, d)), + Some((_, bd)) if d < bd => best = Some((i, d)), + _ => {} + } + } + best.map(|(i, _)| i) + } + + /// returns the index of the segment whose nearest point is closest to (x, y). + pub fn closest_segment(&self, x: f64, y: f64) -> Option { + let mut best: Option<(usize, f64)> = None; + for (i, s) in self.segments.iter().enumerate() { + let (Some(p0), Some(p1)) = (self.nodes.get(s.n0 as usize), self.nodes.get(s.n1 as usize)) else { continue }; + let d = point_to_segment_distance(x, y, p0.x, p0.y, p1.x, p1.y); + match best { + None => best = Some((i, d)), + Some((_, bd)) if d < bd => best = Some((i, d)), + _ => {} + } + } + best.map(|(i, _)| i) + } + + /// clears the selection flag on every geometric entity in the doc. + pub fn clear_selection(&mut self) { + for n in &mut self.nodes { n.selected = false; } + for s in &mut self.segments { s.selected = false; } + for a in &mut self.arcs { a.selected = false; } + for b in &mut self.block_labels { b.selected = false; } + } +} + +/// Euclidean distance from (px, py) to the segment between (ax, ay) and (bx, by). +fn point_to_segment_distance(px: f64, py: f64, ax: f64, ay: f64, bx: f64, by: f64) -> f64 { + let dx = bx - ax; + let dy = by - ay; + let len2 = dx * dx + dy * dy; + if len2 < 1e-18 { + return (px - ax).hypot(py - ay); + } + let t = (((px - ax) * dx + (py - ay) * dy) / len2).clamp(0.0, 1.0); + let cx = ax + t * dx; + let cy = ay + t * dy; + (px - cx).hypot(py - cy) +} diff --git a/crates/femm-doc/src/geom.rs b/crates/femm-doc/src/geom.rs new file mode 100644 index 0000000..375e2f9 --- /dev/null +++ b/crates/femm-doc/src/geom.rs @@ -0,0 +1,54 @@ +//! geometric primitives mirroring the .fem control-point model. + +/// control point in the planar geometry. +#[derive(Debug, Clone, Default)] +pub struct Node { + pub x: f64, + pub y: f64, + pub boundary_marker: String, + pub in_group: i32, + pub selected: bool, +} + +/// straight segment joining two control points by index. +#[derive(Debug, Clone, Default)] +pub struct Segment { + pub n0: i32, + pub n1: i32, + pub max_side_length: f64, + pub boundary_marker: String, + pub hidden: bool, + pub in_group: i32, + pub selected: bool, +} + +/// arc segment from n0 to n1, sweeping the given degree arc. +#[derive(Debug, Clone, Default)] +pub struct ArcSegment { + pub n0: i32, + pub n1: i32, + pub arc_length: f64, + pub max_side_length: f64, + pub boundary_marker: String, + pub hidden: bool, + pub in_group: i32, + pub normal_direction: bool, + pub selected: bool, +} + +/// region label assigning a material, circuit, and mesh-area constraint to an enclosed area. +#[derive(Debug, Clone, Default)] +pub struct BlockLabel { + pub x: f64, + pub y: f64, + pub max_area: f64, + pub mag_dir: f64, + pub mag_dir_fctn: String, + pub turns: i32, + pub block_type: String, + pub in_circuit: String, + pub in_group: i32, + pub is_external: bool, + pub is_default: bool, + pub selected: bool, +} diff --git a/crates/femm-doc/src/geom_math.rs b/crates/femm-doc/src/geom_math.rs new file mode 100644 index 0000000..37eeb54 --- /dev/null +++ b/crates/femm-doc/src/geom_math.rs @@ -0,0 +1,322 @@ +//! 2D intersection math and distance queries over (x, y) pairs and degree-sweep arcs. + +use num_complex::Complex64; + +/// closeness factor for common-endpoint and small-interval rejection. +const EPS_FRAC: f64 = 1.0e-8; + +/// recovers the center and radius of the circle through an arc's endpoints with the given sweep. +pub fn circle_from_arc(p0: (f64, f64), p1: (f64, f64), sweep_deg: f64) -> (f64, f64, f64) { + let a0 = Complex64::new(p0.0, p0.1); + let a1 = Complex64::new(p1.0, p1.1); + let chord = (a1 - a0).norm(); + let t = (a1 - a0) / chord; + let tta = sweep_deg.to_radians(); + let r = chord / (2.0 * (tta / 2.0).sin()); + let c = a0 + (Complex64::new(chord / 2.0, (r * r - chord * chord / 4.0).sqrt())) * t; + (c.re, c.im, r.abs()) +} + +/// shortest distance from a point to a segment. +pub fn shortest_distance_from_segment(p: (f64, f64), a: (f64, f64), b: (f64, f64)) -> f64 { + let dx = b.0 - a.0; + let dy = b.1 - a.1; + let len2 = dx * dx + dy * dy; + if len2 < 1e-18 { + return ((p.0 - a.0).powi(2) + (p.1 - a.1).powi(2)).sqrt(); + } + let t = (((p.0 - a.0) * dx + (p.1 - a.1) * dy) / len2).clamp(0.0, 1.0); + let cx = a.0 + t * dx; + let cy = a.1 + t * dy; + ((p.0 - cx).powi(2) + (p.1 - cy).powi(2)).sqrt() +} + +/// shortest distance from a point to a sweep-defined arc. +pub fn shortest_distance_from_arc( + p: (f64, f64), + a0: (f64, f64), + a1: (f64, f64), + sweep_deg: f64, +) -> f64 { + let (cx, cy, r) = circle_from_arc(a0, a1, sweep_deg); + let pc = Complex64::new(p.0 - cx, p.1 - cy); + let d = pc.norm(); + if d == 0.0 { + return r; + } + let t = pc / d; + let foot = pc - Complex64::new(r * t.re, r * t.im); + let l = foot.norm(); + let a0_dir = Complex64::new(a0.0 - cx, a0.1 - cy); + let z = ((t / a0_dir).arg() * 180.0 / std::f64::consts::PI + 360.0) % 360.0; + let sweep_abs = sweep_deg.abs(); + if z > 0.0 && z < sweep_abs { + return l; + } + let e0 = ((p.0 - a0.0).powi(2) + (p.1 - a0.1).powi(2)).sqrt(); + let e1 = ((p.0 - a1.0).powi(2) + (p.1 - a1.1).powi(2)).sqrt(); + e0.min(e1) +} + +/// intersection of the open segments p0->p1 and q0->q1, or None when none exists. +pub fn line_line_intersection( + p0: (f64, f64), + p1: (f64, f64), + q0: (f64, f64), + q1: (f64, f64), +) -> Option<(f64, f64)> { + let p0c = Complex64::new(p0.0, p0.1); + let p1c = Complex64::new(p1.0, p1.1); + let q0c = Complex64::new(q0.0, q0.1); + let q1c = Complex64::new(q1.0, q1.1); + + // shared endpoint -> no other intersection + if (p0c - q0c).norm() < f64::EPSILON + || (p0c - q1c).norm() < f64::EPSILON + || (p1c - q0c).norm() < f64::EPSILON + || (p1c - q1c).norm() < f64::EPSILON + { + return None; + } + + let ee = ((p1c - p0c).norm()).min((q1c - q0c).norm()) * EPS_FRAC; + + let denom = p1c - p0c; + if denom.norm() < f64::EPSILON { + return None; + } + let r0 = (q0c - p0c) / denom; + let r1 = (q1c - p0c) / denom; + + if r0.re <= 0.0 && r1.re <= 0.0 { return None; } + if r0.re >= 1.0 && r1.re >= 1.0 { return None; } + if r0.im <= 0.0 && r1.im <= 0.0 { return None; } + if r0.im >= 0.0 && r1.im >= 0.0 { return None; } + + let denom_im = r0.im - r1.im; + if denom_im.abs() < f64::EPSILON { return None; } + let z = r0.im / denom_im; + let x = ((1.0 - z) * r0 + z * r1).re; + if x < ee || x > 1.0 - ee { + return None; + } + + let hit = Complex64::new(p0.0, p0.1) * (1.0 - z) + Complex64::new(q0.0, q0.1) * 0.0; + let _ = hit; + let result = Complex64::new(q0.0, q0.1) * (1.0 - z) + Complex64::new(q1.0, q1.1) * z; + Some((result.re, result.im)) +} + +/// intersection points of an open line segment with an open arc, returning 0, 1, or 2 hits. +pub fn line_arc_intersection( + p0: (f64, f64), + p1: (f64, f64), + a0: (f64, f64), + a1: (f64, f64), + sweep_deg: f64, +) -> Vec<(f64, f64)> { + let p0c = Complex64::new(p0.0, p0.1); + let p1c = Complex64::new(p1.0, p1.1); + let a0c = Complex64::new(a0.0, a0.1); + let (cx, cy, r) = circle_from_arc(a0, a1, sweep_deg); + let c = Complex64::new(cx, cy); + + let d = (p1c - p0c).norm(); + if d < f64::EPSILON { + return Vec::new(); + } + let t = (p1c - p0c) / d; + let v = (c - p0c) / t; + + if v.im.abs() > r { + return Vec::new(); + } + let l = (r * r - v.im * v.im).sqrt(); + let tta = sweep_deg.to_radians().abs(); + + let mut out = Vec::with_capacity(2); + + if (l / r) < 1.0e-5 { + let hit = p0c + Complex64::new(v.re, 0.0) * t; + let r_param = ((hit - p0c) / t).re; + let z = ((hit - c) / (a0c - c)).arg(); + if r_param > 0.0 && r_param < d && z > 0.0 && z < tta { + out.push((hit.re, hit.im)); + } + return out; + } + + for sign in [1.0, -1.0] { + let hit = p0c + Complex64::new(v.re + sign * l, 0.0) * t; + let r_param = ((hit - p0c) / t).re; + let z = ((hit - c) / (a0c - c)).arg(); + if r_param > 0.0 && r_param < d && z > 0.0 && z < tta { + out.push((hit.re, hit.im)); + } + } + out +} + +/// intersection points of two open arcs, returning 0, 1, or 2 hits. +pub fn arc_arc_intersection( + a0_p0: (f64, f64), + a0_p1: (f64, f64), + a0_sweep_deg: f64, + a1_p0: (f64, f64), + a1_p1: (f64, f64), + a1_sweep_deg: f64, +) -> Vec<(f64, f64)> { + let (c0x, c0y, r0) = circle_from_arc(a0_p0, a0_p1, a0_sweep_deg); + let (c1x, c1y, r1) = circle_from_arc(a1_p0, a1_p1, a1_sweep_deg); + let c0 = Complex64::new(c0x, c0y); + let c1 = Complex64::new(c1x, c1y); + let d = (c1 - c0).norm(); + + if d > r0 + r1 || d < 1.0e-8 { + return Vec::new(); + } + + let l = ((r0 + r1 - d) * (d + r0 - r1) * (d - r0 + r1) * (d + r0 + r1)).sqrt() / (2.0 * d); + let c = 1.0 + (r0 / d) * (r0 / d) - (r1 / d) * (r1 / d); + let t = (c1 - c0) / d; + let tta0 = a0_sweep_deg.to_radians().abs(); + let tta1 = a1_sweep_deg.to_radians().abs(); + let a0c = Complex64::new(a0_p0.0, a0_p0.1); + let a1c = Complex64::new(a1_p0.0, a1_p0.1); + + let mut out = Vec::with_capacity(2); + + let first = c0 + Complex64::new(c * d / 2.0, l) * t; + let z0 = ((first - c0) / (a0c - c0)).arg(); + let z1 = ((first - c1) / (a1c - c1)).arg(); + if z0 > 0.0 && z0 < tta0 && z1 > 0.0 && z1 < tta1 { + out.push((first.re, first.im)); + } + + // tangent-touch case: only one intersection + if (d - r0 + r1).abs() / (r0 + r1) < 1.0e-5 { + return out; + } + + let second = c0 + Complex64::new(c * d / 2.0, -l) * t; + let z0 = ((second - c0) / (a0c - c0)).arg(); + let z1 = ((second - c1) / (a1c - c1)).arg(); + if z0 > 0.0 && z0 < tta0 && z1 > 0.0 && z1 < tta1 { + out.push((second.re, second.im)); + } + + out +} + +#[cfg(test)] +mod tests { + use super::*; + + fn close(a: (f64, f64), b: (f64, f64), tol: f64) -> bool { + ((a.0 - b.0).powi(2) + (a.1 - b.1).powi(2)).sqrt() < tol + } + + #[test] + fn line_line_crosses_at_origin() { + // (-1,-1)->(1,1) crosses (-1,1)->(1,-1) at origin. + let hit = line_line_intersection((-1.0, -1.0), (1.0, 1.0), (-1.0, 1.0), (1.0, -1.0)); + assert!(hit.is_some()); + assert!(close(hit.unwrap(), (0.0, 0.0), 1e-9)); + } + + #[test] + fn line_line_parallel_no_intersection() { + let hit = line_line_intersection((0.0, 0.0), (1.0, 0.0), (0.0, 1.0), (1.0, 1.0)); + assert!(hit.is_none()); + } + + #[test] + fn line_line_t_junction_at_endpoint_rejects() { + // the second segment's q0 sits ON the first segment's interior + // but the prospective endpoint check should reject (intersection at q0). + let hit = line_line_intersection((0.0, 0.0), (2.0, 0.0), (1.0, 0.0), (1.0, 1.0)); + // shape: one endpoint sits on the other segment, no clean crossing. + // GetIntersection returns FALSE when an intersection lies within ee of an endpoint of the + // prospective line — here the meeting point IS the prospective q0, so no split is wanted. + assert!(hit.is_none()); + } + + #[test] + fn line_line_shared_endpoint_no_intersection() { + let hit = line_line_intersection((0.0, 0.0), (1.0, 0.0), (1.0, 0.0), (2.0, 1.0)); + assert!(hit.is_none()); + } + + #[test] + fn circle_from_quarter_arc_unit() { + // a quarter circle (90 deg sweep) from (1, 0) to (0, 1) sits on the unit circle centered at origin. + let (cx, cy, r) = circle_from_arc((1.0, 0.0), (0.0, 1.0), 90.0); + assert!((cx - 0.0).abs() < 1e-9); + assert!((cy - 0.0).abs() < 1e-9); + assert!((r - 1.0).abs() < 1e-9); + } + + #[test] + fn line_through_quarter_arc_crosses_once() { + // diagonal line from (-2, 0.5) to (2, 0.5) crosses the quarter circle (1,0)->(0,1) at one point. + let hits = line_arc_intersection((-2.0, 0.5), (2.0, 0.5), (1.0, 0.0), (0.0, 1.0), 90.0); + assert_eq!(hits.len(), 1); + let h = hits[0]; + assert!((h.0 * h.0 + h.1 * h.1 - 1.0).abs() < 1e-9, "point should lie on unit circle: {h:?}"); + assert!((h.1 - 0.5).abs() < 1e-9); + assert!(h.0 > 0.0); + } + + #[test] + fn line_outside_arc_does_not_cross() { + // line clearly outside the arc swept region. + let hits = line_arc_intersection((-2.0, 2.5), (2.0, 2.5), (1.0, 0.0), (0.0, 1.0), 90.0); + assert!(hits.is_empty()); + } + + #[test] + fn arc_arc_two_intersections() { + // two unit circles offset by 1.0 along x, both swept 360 degrees... but we use partial arcs + // shaped to clip both crossings in their swept range. + // circle A: center (0,0), arcs covering full top half (180 deg from (1,0) to (-1,0)). + // circle B: center (1,0), arcs covering full top half (180 deg from (2,0) to (0,0)). + let hits = arc_arc_intersection( + (1.0, 0.0), (-1.0, 0.0), 180.0, + (2.0, 0.0), (0.0, 0.0), 180.0, + ); + assert_eq!(hits.len(), 1, "two top-half arcs of overlapping unit circles meet at one upper crossing: {hits:?}"); + let h = hits[0]; + assert!((h.0 - 0.5).abs() < 1e-9, "x = 0.5 by symmetry, got {h:?}"); + assert!(h.1 > 0.0, "intersection should be above the x-axis: {h:?}"); + } + + #[test] + fn arc_arc_disjoint() { + let hits = arc_arc_intersection( + (1.0, 0.0), (-1.0, 0.0), 180.0, + (11.0, 0.0), (9.0, 0.0), 180.0, + ); + assert!(hits.is_empty()); + } + + #[test] + fn shortest_distance_from_segment_perpendicular() { + // point (0, 1) to segment (-1, 0)->(1, 0) is distance 1. + let d = shortest_distance_from_segment((0.0, 1.0), (-1.0, 0.0), (1.0, 0.0)); + assert!((d - 1.0).abs() < 1e-9); + } + + #[test] + fn shortest_distance_from_segment_clamps_to_endpoint() { + // point (-2, 0) to segment (-1, 0)->(1, 0) is distance 1 (clamps to left endpoint). + let d = shortest_distance_from_segment((-2.0, 0.0), (-1.0, 0.0), (1.0, 0.0)); + assert!((d - 1.0).abs() < 1e-9); + } + + #[test] + fn shortest_distance_from_arc_radial() { + // point at origin distance to unit-circle quarter arc is 1. + let d = shortest_distance_from_arc((0.0, 0.0), (1.0, 0.0), (0.0, 1.0), 90.0); + assert!((d - 1.0).abs() < 1e-9); + } +} diff --git a/crates/femm-doc/src/lib.rs b/crates/femm-doc/src/lib.rs new file mode 100644 index 0000000..9143cd3 --- /dev/null +++ b/crates/femm-doc/src/lib.rs @@ -0,0 +1,94 @@ +//! magnetostatic pre-processor document model: geometry, properties, parser, writer. + +pub mod geom; +pub mod geom_math; +pub mod props; +pub mod parser; +pub mod writer; +pub mod edit; + +use num_complex::Complex64; + +pub use geom::{ArcSegment, BlockLabel, Node, Segment}; +pub use props::{BoundaryProp, CircuitProp, MaterialProp, PointProp}; + +#[derive(Debug, Clone, Copy, PartialEq, Eq)] +pub enum ProblemType { + Planar, + Axisymmetric, +} + +impl Default for ProblemType { + fn default() -> Self { ProblemType::Planar } +} + +#[derive(Debug, Clone, Copy, PartialEq, Eq, Default)] +pub enum Coords { + #[default] + Cartesian, + Polar, +} + +#[derive(Debug, Clone, Copy, PartialEq, Eq)] +pub enum LengthUnit { + Inches = 0, + Millimeters = 1, + Centimeters = 2, + Meters = 3, + Mils = 4, + Microns = 5, +} + +impl Default for LengthUnit { + fn default() -> Self { LengthUnit::Inches } +} + +#[derive(Debug, Clone, Copy, PartialEq, Eq, Default)] +pub enum PrevType { + #[default] + None = 0, + Incremental = 1, + Frozen = 2, +} + +#[derive(Debug, Clone, Copy, PartialEq, Eq, Default)] +pub enum ACSolver { + #[default] + SuccessiveApprox = 0, + Newton = 1, +} + +/// magnetostatic .fem document: geometry, property tables, problem attributes. +#[derive(Debug, Clone, Default)] +pub struct FemmDoc { + 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 prev_soln: String, + pub prev_type: PrevType, + + pub ext_ro: f64, + pub ext_ri: f64, + pub ext_zo: f64, + + pub nodes: Vec, + pub segments: Vec, + pub arcs: Vec, + pub block_labels: Vec, + + pub materials: Vec, + pub boundaries: Vec, + pub points: Vec, + pub circuits: Vec, +} + +/// f64 complex number alias matching the C++ CComplex memory layout. +pub type Complex = Complex64; diff --git a/crates/femm-doc/src/parser.rs b/crates/femm-doc/src/parser.rs new file mode 100644 index 0000000..4c2807e --- /dev/null +++ b/crates/femm-doc/src/parser.rs @@ -0,0 +1,486 @@ +//! .fem text parser: `[Section] = value` scalars and `...` property stanzas. + +use crate::{ + ACSolver, ArcSegment, BlockLabel, BoundaryProp, CircuitProp, Coords, FemmDoc, LengthUnit, + MaterialProp, Node, PointProp, PrevType, ProblemType, Segment, +}; +use num_complex::Complex64; +use std::path::Path; +use thiserror::Error; + +#[derive(Debug, Error)] +pub enum ParseError { + #[error("io error: {0}")] + Io(#[from] std::io::Error), + #[error("legacy 3.2-format .fem files are not supported")] + LegacyFormat, + #[error("unexpected end of file while reading {context}")] + UnexpectedEof { context: &'static str }, + #[error("malformed value in section {section}: {reason}")] + MalformedValue { section: &'static str, reason: String }, +} + +impl FemmDoc { + /// parses a .fem text buffer into a document. + pub fn parse(src: &str) -> Result { + let mut doc = FemmDoc::default(); + doc.depth = -1.0; + doc.smart_mesh = true; + doc.format = 4.0; + + if src + .lines() + .next() + .map(|l| l.trim_start().starts_with("Frequency")) + .unwrap_or(false) + { + return Err(ParseError::LegacyFormat); + } + + let mut lines = src.lines(); + let mut version: i32 = 0; + + let mut point = PointProp::default(); + let mut bdry = BoundaryProp::default(); + let mut mat = MaterialProp::default(); + let mut circ = CircuitProp::default(); + + while let Some(raw) = lines.next() { + let tok = first_token(raw); + if tok.is_empty() { continue; } + + // problem attributes + if tok.eq_ignore_ascii_case("[format]") { + let v = strip_key(raw); + let f = parse_f64(v); + doc.format = f; + version = (10.0 * f + 0.5) as i32; + } else if tok.eq_ignore_ascii_case("[frequency]") { + doc.frequency = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[depth]") { + doc.depth = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[precision]") { + doc.precision = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[acsolver]") { + doc.ac_solver = match parse_i32(strip_key(raw)) { + 1 => ACSolver::Newton, + _ => ACSolver::SuccessiveApprox, + }; + } else if tok.eq_ignore_ascii_case("[minangle]") { + doc.min_angle = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[dosmartmesh]") { + doc.smart_mesh = parse_i32(strip_key(raw)) != 0; + } else if tok.eq_ignore_ascii_case("[lengthunits]") { + doc.length_units = parse_length_units(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[problemtype]") { + doc.problem_type = parse_problem_type(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[coordinates]") { + doc.coords = parse_coords(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[comment]") { + doc.comment = unescape_comment(unquote(strip_key(raw))); + } else if tok.eq_ignore_ascii_case("[prevsoln]") { + doc.prev_soln = unquote(strip_key(raw)); + if doc.prev_soln.is_empty() { doc.prev_type = PrevType::None; } + } else if tok.eq_ignore_ascii_case("[prevtype]") { + doc.prev_type = match parse_i32(strip_key(raw)) { + 1 => PrevType::Incremental, + 2 => PrevType::Frozen, + _ => PrevType::None, + }; + } else if tok.eq_ignore_ascii_case("[extzo]") { + doc.ext_zo = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[extro]") { + doc.ext_ro = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("[extri]") { + doc.ext_ri = parse_f64(strip_key(raw)); + + // point property stanza + } else if tok.eq_ignore_ascii_case("") { + point = PointProp::default(); + point.name = String::from("New Point Property"); + } else if tok.eq_ignore_ascii_case("") { + point.name = unquote(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + point.ap.re = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + point.ap.im = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + point.jp.re = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + point.jp.im = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + doc.points.push(std::mem::take(&mut point)); + + // boundary property stanza + } else if tok.eq_ignore_ascii_case("") { + bdry = BoundaryProp::default(); + bdry.name = String::from("New Boundary"); + } else if tok.eq_ignore_ascii_case("") { + bdry.name = unquote(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.format = parse_i32(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.mu = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.sig = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.a0 = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.a1 = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.a2 = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.phi = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.c0.re = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.c0.im = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.c1.re = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.c1.im = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.inner_angle = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + bdry.outer_angle = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + doc.boundaries.push(std::mem::take(&mut bdry)); + + // material property stanza + } else if tok.eq_ignore_ascii_case("") { + mat = MaterialProp::default(); + mat.name = String::from("New Material"); + mat.mu_x = 1.0; + mat.mu_y = 1.0; + mat.lam_fill = 1.0; + } else if tok.eq_ignore_ascii_case("") { + mat.name = unquote(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + mat.mu_x = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + mat.mu_y = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + mat.h_c = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + mat.theta_m = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + mat.j_src.re = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + mat.j_src.im = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + mat.cduct = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + mat.theta_hn = parse_f64(strip_key(raw)); + if version == 30 { + mat.theta_hx = mat.theta_hn; + mat.theta_hy = mat.theta_hn; + } + } else if tok.eq_ignore_ascii_case("") { + mat.theta_hx = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + mat.theta_hy = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + mat.lam_d = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + mat.lam_fill = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + mat.lam_type = parse_i32(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + mat.n_strands = parse_i32(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + mat.wire_d = parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + let count = parse_i32(strip_key(raw)); + if count > 0 { + mat.bh_curve.clear(); + mat.bh_curve.reserve(count as usize); + for _ in 0..count { + let line = lines.next().ok_or(ParseError::UnexpectedEof { + context: "", + })?; + let (b, rest) = take_f64(line); + let (h, _) = take_f64(rest); + mat.bh_curve.push((b, h)); + } + } + } else if tok.eq_ignore_ascii_case("") { + doc.materials.push(std::mem::take(&mut mat)); + + // circuit stanza + } else if tok.eq_ignore_ascii_case("") { + circ = CircuitProp::default(); + circ.name = String::from("New Circuit"); + } else if tok.eq_ignore_ascii_case("") { + circ.name = unquote(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + circ.amps.re += parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + circ.amps.im += parse_f64(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + circ.circ_type = parse_i32(strip_key(raw)); + } else if tok.eq_ignore_ascii_case("") { + doc.circuits.push(std::mem::take(&mut circ)); + + // counted geometry sections + } else if tok.eq_ignore_ascii_case("[numpoints]") { + let n = parse_i32(strip_key(raw)); + for _ in 0..n { + let line = lines.next().ok_or(ParseError::UnexpectedEof { context: "[NumPoints]" })?; + let (x, rest) = take_f64(line); + let (y, rest) = take_f64(rest); + let (bi, rest) = take_i32(rest); + let (group, _) = take_i32(rest); + let bm = resolve_name(bi, doc.points.iter().map(|p| p.name.as_str())); + doc.nodes.push(Node { + x, y, + boundary_marker: bm, + in_group: group, + selected: false, + }); + } + } else if tok.eq_ignore_ascii_case("[numsegments]") { + let n = parse_i32(strip_key(raw)); + for _ in 0..n { + let line = lines.next().ok_or(ParseError::UnexpectedEof { context: "[NumSegments]" })?; + let (n0, rest) = take_i32(line); + let (n1, rest) = take_i32(rest); + let (msl, rest) = take_f64(rest); + let (bi, rest) = take_i32(rest); + let (hidden, rest) = take_i32(rest); + let (group, _) = take_i32(rest); + let bm = resolve_name(bi, doc.boundaries.iter().map(|b| b.name.as_str())); + doc.segments.push(Segment { + n0, n1, + max_side_length: msl, + boundary_marker: bm, + hidden: hidden != 0, + in_group: group, + selected: false, + }); + } + } else if tok.eq_ignore_ascii_case("[numarcsegments]") { + let n = parse_i32(strip_key(raw)); + for _ in 0..n { + let line = lines.next().ok_or(ParseError::UnexpectedEof { context: "[NumArcSegments]" })?; + let (n0, rest) = take_i32(line); + let (n1, rest) = take_i32(rest); + let (al, rest) = take_f64(rest); + let (msl, rest) = take_f64(rest); + let (bi, rest) = take_i32(rest); + let (hidden, rest) = take_i32(rest); + let (group, _) = take_i32(rest); + let bm = resolve_name(bi, doc.boundaries.iter().map(|b| b.name.as_str())); + doc.arcs.push(ArcSegment { + n0, n1, + arc_length: al, + max_side_length: msl, + boundary_marker: bm, + hidden: hidden != 0, + in_group: group, + normal_direction: true, + selected: false, + }); + } + } else if tok.eq_ignore_ascii_case("[numholes]") { + let n = parse_i32(strip_key(raw)); + for _ in 0..n { + let line = lines.next().ok_or(ParseError::UnexpectedEof { context: "[NumHoles]" })?; + let (x, rest) = take_f64(line); + let (y, rest) = take_f64(rest); + let (group, _) = take_i32(rest); + doc.block_labels.push(BlockLabel { + x, y, + max_area: 0.0, + mag_dir: 0.0, + mag_dir_fctn: String::new(), + turns: 1, + block_type: String::from(""), + in_circuit: String::from(""), + in_group: group, + is_external: false, + is_default: false, + selected: false, + }); + } + } else if tok.eq_ignore_ascii_case("[numblocklabels]") { + let n = parse_i32(strip_key(raw)); + for _ in 0..n { + let line = lines.next().ok_or(ParseError::UnexpectedEof { context: "[NumBlockLabels]" })?; + let (x, rest) = take_f64(line); + let (y, rest) = take_f64(rest); + let (mi, rest) = take_i32(rest); + let (ma_diam, rest) = take_f64(rest); + let (ci, rest) = take_i32(rest); + let (mag_dir, rest) = take_f64(rest); + let (group, rest) = take_i32(rest); + let (turns, rest) = take_i32(rest); + let (ext_flags, rest) = take_i32(rest); + let mdf = take_quoted(rest); + + let block_type = if mi == 0 { + String::from("") + } else { + resolve_name(mi, doc.materials.iter().map(|m| m.name.as_str())) + }; + let in_circuit = if ci == 0 { + String::from("") + } else { + resolve_name(ci, doc.circuits.iter().map(|c| c.name.as_str())) + }; + let max_area = if ma_diam < 0.0 { + 0.0 + } else { + std::f64::consts::PI * ma_diam * ma_diam / 4.0 + }; + + doc.block_labels.push(BlockLabel { + x, y, + max_area, + mag_dir, + mag_dir_fctn: mdf, + turns, + block_type, + in_circuit, + in_group: group, + is_external: ext_flags & 1 != 0, + is_default: ext_flags & 2 != 0, + selected: false, + }); + } + } + // unknown tokens are skipped silently, matching the original parser. + } + + // 3.2-era files omitted [Depth]; fill the default in the current length unit. + if doc.depth == -1.0 { + doc.depth = match doc.length_units { + LengthUnit::Millimeters => 1000.0, + LengthUnit::Centimeters => 100.0, + LengthUnit::Meters => 1.0, + LengthUnit::Mils => 1000.0 / 0.0254, + LengthUnit::Microns => 1.0e6, + LengthUnit::Inches => 1.0 / 0.0254, + }; + } + + Ok(doc) + } + + /// loads a .fem file by path. + pub fn open(path: impl AsRef) -> Result { + let text = std::fs::read_to_string(path)?; + Self::parse(&text) + } +} + +fn first_token(line: &str) -> &str { + line.trim_start().split_whitespace().next().unwrap_or("") +} + +fn strip_key(line: &str) -> &str { + match line.find('=') { + Some(i) => &line[i + 1..], + None => "", + } +} + +fn unquote(s: &str) -> String { + let s = s.trim(); + let first = s.find('"'); + let last = s.rfind('"'); + match (first, last) { + (Some(a), Some(b)) if a < b => s[a + 1..b].to_string(), + _ => s.to_string(), + } +} + +fn unescape_comment(s: String) -> String { + let bytes = s.as_bytes(); + let mut out = String::with_capacity(s.len()); + let mut i = 0; + while i < bytes.len() { + if bytes[i] == b'\\' && i + 1 < bytes.len() && bytes[i + 1] == b'n' { + out.push('\r'); + out.push('\n'); + i += 2; + } else { + out.push(bytes[i] as char); + i += 1; + } + } + out +} + +fn parse_f64(s: &str) -> f64 { + let s = s.trim().trim_matches(|c: char| c == ',' || c.is_whitespace()); + let end = s.find(|c: char| c.is_whitespace() || c == ',').unwrap_or(s.len()); + s[..end].parse::().unwrap_or(0.0) +} + +fn parse_i32(s: &str) -> i32 { + let s = s.trim().trim_matches(|c: char| c == ',' || c.is_whitespace()); + let end = s.find(|c: char| c.is_whitespace() || c == ',').unwrap_or(s.len()); + s[..end].parse::().unwrap_or(0) +} + +/// peels one f64 off the head of a whitespace/comma-separated row, returning the remainder. +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/comma-separated row, returning the remainder. +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) +} + +/// extracts the contents of the first double-quoted span in `s`, or returns empty. +fn take_quoted(s: &str) -> String { + let first = s.find('"'); + if let Some(a) = first { + let after = &s[a + 1..]; + if let Some(b) = after.find('"') { + return after[..b].to_string(); + } + } + String::new() +} + +fn parse_length_units(s: &str) -> LengthUnit { + let t = first_token(s).to_ascii_lowercase(); + if t.starts_with("millimeters") { LengthUnit::Millimeters } + else if t.starts_with("c") { LengthUnit::Centimeters } // matches original's 1-char check + else if t.starts_with("meters") { LengthUnit::Meters } + else if t.starts_with("mils") { LengthUnit::Mils } + else if t.starts_with("microns") { LengthUnit::Microns } + else { LengthUnit::Inches } +} + +fn parse_problem_type(s: &str) -> ProblemType { + let t = first_token(s).to_ascii_lowercase(); + if t.starts_with("axi") { ProblemType::Axisymmetric } else { ProblemType::Planar } +} + +fn parse_coords(s: &str) -> Coords { + let t = first_token(s).to_ascii_lowercase(); + if t.starts_with("polar") { Coords::Polar } else { Coords::Cartesian } +} + +/// resolves a 1-based property index into the matching property name, falling back to empty. +fn resolve_name<'a, I: Iterator>(idx: i32, names: I) -> String { + if idx <= 0 { return String::new(); } + let want = idx as usize; + for (i, name) in names.enumerate() { + if i + 1 == want { return name.to_string(); } + } + String::new() +} + +// silence dead-code warning in lib.rs until the writer module consumes it. +#[allow(dead_code)] +fn _force_complex_use(_c: Complex64) {} diff --git a/crates/femm-doc/src/props.rs b/crates/femm-doc/src/props.rs new file mode 100644 index 0000000..e3fddb9 --- /dev/null +++ b/crates/femm-doc/src/props.rs @@ -0,0 +1,57 @@ +//! magnetostatic material, boundary, point, and circuit properties. + +use crate::Complex; + +/// material property entry referenced by name from block labels. +#[derive(Debug, Clone, Default)] +pub struct MaterialProp { + pub name: String, + pub mu_x: f64, + pub mu_y: f64, + pub h_c: f64, + pub theta_m: f64, + pub j_src: Complex, + pub cduct: f64, + pub lam_d: f64, + pub theta_hn: f64, + pub theta_hx: f64, + pub theta_hy: f64, + pub lam_type: i32, + pub lam_fill: f64, + pub n_strands: i32, + pub wire_d: f64, + pub bh_curve: Vec<(f64, f64)>, +} + +/// boundary condition entry referenced by name from segments and arcs. +#[derive(Debug, Clone, Default)] +pub struct BoundaryProp { + pub name: String, + pub format: i32, + pub a0: f64, + pub a1: f64, + pub a2: f64, + pub phi: f64, + pub mu: f64, + pub sig: f64, + pub c0: Complex, + pub c1: Complex, + pub inner_angle: f64, + pub outer_angle: f64, +} + +/// point property entry referenced by name from nodes. +#[derive(Debug, Clone, Default)] +pub struct PointProp { + pub name: String, + pub jp: Complex, + pub ap: Complex, +} + +/// circuit property entry referenced by name from block labels. +#[derive(Debug, Clone, Default)] +pub struct CircuitProp { + pub name: String, + pub amps: Complex, + pub circ_type: i32, +} diff --git a/crates/femm-doc/src/writer.rs b/crates/femm-doc/src/writer.rs new file mode 100644 index 0000000..cdba9f8 --- /dev/null +++ b/crates/femm-doc/src/writer.rs @@ -0,0 +1,222 @@ +//! .fem text writer, inverse of [`parser`](crate::parser). + +use crate::{ACSolver, Coords, FemmDoc, LengthUnit, PrevType, ProblemType}; +use std::fmt::Write; +use std::path::Path; + +impl FemmDoc { + /// renders the document to .fem text. + pub fn write(&self) -> String { + let mut out = String::new(); + writeln!(out, "[Format] = {:.17e}", 4.0_f64).unwrap(); + writeln!(out, "[Frequency] = {:.17e}", self.frequency).unwrap(); + writeln!(out, "[Precision] = {:.17e}", self.precision).unwrap(); + writeln!(out, "[MinAngle] = {:.17e}", self.min_angle).unwrap(); + writeln!(out, "[DoSmartMesh] = {}", self.smart_mesh as i32).unwrap(); + writeln!(out, "[Depth] = {:.17e}", self.depth).unwrap(); + + let units = match self.length_units { + LengthUnit::Inches => "inches", + LengthUnit::Millimeters => "millimeters", + LengthUnit::Centimeters => "centimeters", + LengthUnit::Meters => "meters", + LengthUnit::Mils => "mils", + LengthUnit::Microns => "microns", + }; + writeln!(out, "[LengthUnits] = {units}").unwrap(); + + match self.problem_type { + ProblemType::Planar => writeln!(out, "[ProblemType] = planar").unwrap(), + ProblemType::Axisymmetric => { + writeln!(out, "[ProblemType] = axisymmetric").unwrap(); + if self.ext_ro != 0.0 && self.ext_ri != 0.0 { + writeln!(out, "[extZo] = {:.17e}", self.ext_zo).unwrap(); + writeln!(out, "[extRo] = {:.17e}", self.ext_ro).unwrap(); + writeln!(out, "[extRi] = {:.17e}", self.ext_ri).unwrap(); + } + } + } + + let coords = match self.coords { + Coords::Cartesian => "cartesian", + Coords::Polar => "polar", + }; + writeln!(out, "[Coordinates] = {coords}").unwrap(); + + let ac = match self.ac_solver { + ACSolver::SuccessiveApprox => 0, + ACSolver::Newton => 1, + }; + writeln!(out, "[ACSolver] = {ac}").unwrap(); + let prev_idx = match self.prev_type { + PrevType::None => 0, + PrevType::Incremental => 1, + PrevType::Frozen => 2, + }; + writeln!(out, "[PrevType] = {prev_idx}").unwrap(); + writeln!(out, "[PrevSoln] = \"{}\"", self.prev_soln).unwrap(); + writeln!(out, "[Comment] = \"{}\"", escape_comment(&self.comment)).unwrap(); + + // point properties + writeln!(out, "[PointProps] = {}", self.points.len()).unwrap(); + for p in &self.points { + writeln!(out, " ").unwrap(); + writeln!(out, " = \"{}\"", p.name).unwrap(); + writeln!(out, " = {:.17e}", p.jp.re).unwrap(); + writeln!(out, " = {:.17e}", p.jp.im).unwrap(); + writeln!(out, " = {:.17e}", p.ap.re).unwrap(); + writeln!(out, " = {:.17e}", p.ap.im).unwrap(); + writeln!(out, " ").unwrap(); + } + + // boundary properties + writeln!(out, "[BdryProps] = {}", self.boundaries.len()).unwrap(); + for b in &self.boundaries { + writeln!(out, " ").unwrap(); + writeln!(out, " = \"{}\"", b.name).unwrap(); + writeln!(out, " = {}", b.format).unwrap(); + writeln!(out, " = {:.17e}", b.a0).unwrap(); + writeln!(out, " = {:.17e}", b.a1).unwrap(); + writeln!(out, " = {:.17e}", b.a2).unwrap(); + writeln!(out, " = {:.17e}", b.phi).unwrap(); + writeln!(out, " = {:.17e}", b.c0.re).unwrap(); + writeln!(out, " = {:.17e}", b.c0.im).unwrap(); + writeln!(out, " = {:.17e}", b.c1.re).unwrap(); + writeln!(out, " = {:.17e}", b.c1.im).unwrap(); + writeln!(out, " = {:.17e}", b.mu).unwrap(); + writeln!(out, " = {:.17e}", b.sig).unwrap(); + writeln!(out, " = {:.17e}", b.inner_angle).unwrap(); + writeln!(out, " = {:.17e}", b.outer_angle).unwrap(); + writeln!(out, " ").unwrap(); + } + + // material properties + writeln!(out, "[BlockProps] = {}", self.materials.len()).unwrap(); + for m in &self.materials { + writeln!(out, " ").unwrap(); + writeln!(out, " = \"{}\"", m.name).unwrap(); + writeln!(out, " = {:.17e}", m.mu_x).unwrap(); + writeln!(out, " = {:.17e}", m.mu_y).unwrap(); + writeln!(out, " = {:.17e}", m.h_c).unwrap(); + writeln!(out, " = {:.17e}", m.theta_m).unwrap(); + writeln!(out, " = {:.17e}", m.j_src.re).unwrap(); + writeln!(out, " = {:.17e}", m.j_src.im).unwrap(); + writeln!(out, " = {:.17e}", m.cduct).unwrap(); + writeln!(out, " = {:.17e}", m.lam_d).unwrap(); + writeln!(out, " = {:.17e}", m.theta_hn).unwrap(); + writeln!(out, " = {:.17e}", m.theta_hx).unwrap(); + writeln!(out, " = {:.17e}", m.theta_hy).unwrap(); + writeln!(out, " = {}", m.lam_type).unwrap(); + writeln!(out, " = {:.17e}", m.lam_fill).unwrap(); + writeln!(out, " = {}", m.n_strands).unwrap(); + writeln!(out, " = {:.17e}", m.wire_d).unwrap(); + writeln!(out, " = {}", m.bh_curve.len()).unwrap(); + for (b, h) in &m.bh_curve { + writeln!(out, " {:.17e}\t{:.17e}", b, h).unwrap(); + } + writeln!(out, " ").unwrap(); + } + + // circuit properties + writeln!(out, "[CircuitProps] = {}", self.circuits.len()).unwrap(); + for c in &self.circuits { + writeln!(out, " ").unwrap(); + writeln!(out, " = \"{}\"", c.name).unwrap(); + writeln!(out, " = {:.17e}", c.amps.re).unwrap(); + writeln!(out, " = {:.17e}", c.amps.im).unwrap(); + writeln!(out, " = {}", c.circ_type).unwrap(); + writeln!(out, " ").unwrap(); + } + + // nodes + writeln!(out, "[NumPoints] = {}", self.nodes.len()).unwrap(); + for n in &self.nodes { + let t = lookup_idx(&n.boundary_marker, self.points.iter().map(|p| p.name.as_str())); + writeln!(out, "{:.17e}\t{:.17e}\t{}\t{}", n.x, n.y, t, n.in_group).unwrap(); + } + + // segments + writeln!(out, "[NumSegments] = {}", self.segments.len()).unwrap(); + for s in &self.segments { + let t = lookup_idx(&s.boundary_marker, self.boundaries.iter().map(|b| b.name.as_str())); + let msl = if s.max_side_length < 0.0 { + String::from("-1") + } else { + format!("{:.17e}", s.max_side_length) + }; + writeln!(out, "{}\t{}\t{}\t{}\t{}\t{}", + s.n0, s.n1, msl, t, s.hidden as i32, s.in_group).unwrap(); + } + + // arc segments + writeln!(out, "[NumArcSegments] = {}", self.arcs.len()).unwrap(); + for a in &self.arcs { + let t = lookup_idx(&a.boundary_marker, self.boundaries.iter().map(|b| b.name.as_str())); + writeln!(out, "{}\t{}\t{:.17e}\t{:.17e}\t{}\t{}\t{}\t{:.17e}", + a.n0, a.n1, a.arc_length, a.max_side_length, t, + a.hidden as i32, a.in_group, a.max_side_length).unwrap(); + } + + // holes (block labels typed "") + let holes: Vec<_> = self.block_labels.iter().filter(|b| b.block_type == "").collect(); + writeln!(out, "[NumHoles] = {}", holes.len()).unwrap(); + for h in &holes { + writeln!(out, "{:.17e}\t{:.17e}\t{}", h.x, h.y, h.in_group).unwrap(); + } + + // regional attributes + let labels: Vec<_> = self.block_labels.iter().filter(|b| b.block_type != "").collect(); + writeln!(out, "[NumBlockLabels] = {}", labels.len()).unwrap(); + for b in labels { + let mi = lookup_idx(&b.block_type, self.materials.iter().map(|m| m.name.as_str())); + let ci = lookup_idx(&b.in_circuit, self.circuits.iter().map(|c| c.name.as_str())); + let area_field = if b.max_area > 0.0 { + format!("{:.17e}", (4.0 * b.max_area / std::f64::consts::PI).sqrt()) + } else { + String::from("-1") + }; + let flags = (b.is_external as i32) + ((b.is_default as i32) << 1); + let trailing = if b.mag_dir_fctn.is_empty() { + String::new() + } else { + format!("\t\"{}\"", b.mag_dir_fctn) + }; + writeln!(out, "{:.17e}\t{:.17e}\t{}\t{}\t{}\t{:.17e}\t{}\t{}\t{}{}", + b.x, b.y, mi, area_field, ci, b.mag_dir, b.in_group, b.turns, flags, trailing).unwrap(); + } + + out + } + + /// writes the document to a .fem file by path. + pub fn save(&self, path: impl AsRef) -> std::io::Result<()> { + std::fs::write(path, self.write()) + } +} + +/// finds the 1-based index of `name` in the given name iterator, or 0 when absent or empty. +fn lookup_idx<'a, I: Iterator>(name: &str, iter: I) -> i32 { + if name.is_empty() { return 0; } + for (i, n) in iter.enumerate() { + if n == name { return (i as i32) + 1; } + } + 0 +} + +/// inverse of `unescape_comment`: turns CR/LF pairs back into `\n` escapes. +fn escape_comment(s: &str) -> String { + let mut out = String::with_capacity(s.len()); + let bytes = s.as_bytes(); + let mut i = 0; + while i < bytes.len() { + if i + 1 < bytes.len() && bytes[i] == b'\r' && bytes[i + 1] == b'\n' { + out.push('\\'); + out.push('n'); + i += 2; + } else { + out.push(bytes[i] as char); + i += 1; + } + } + out +} diff --git a/crates/femm-doc/tests/roundtrip.rs b/crates/femm-doc/tests/roundtrip.rs new file mode 100644 index 0000000..e2e29d1 --- /dev/null +++ b/crates/femm-doc/tests/roundtrip.rs @@ -0,0 +1,157 @@ +use femm_doc::FemmDoc; + +const FIXTURE: &str = r#"[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] = "two-coil planar magnetics" +[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] = 2 + + = "Air" + = 1 + = 1 + = 0 + = 0 + = 0 + = 0 + = 0 + = 0 + = 0 + = 0 + = 0 + = 0 + = 1 + = 0 + = 0 + = 0 + + + = "Copper" + = 1 + = 1 + = 0 + = 0 + = 0 + = 0 + = 58 + = 0 + = 0 + = 0 + = 0 + = 0 + = 1 + = 0 + = 0 + = 0 + +[CircuitProps] = 1 + + = "Coil" + = 100 + = 0 + = 1 + +[NumPoints] = 4 +0 0 1 0 +10 0 0 0 +10 10 0 0 +0 10 0 0 +[NumSegments] = 4 +0 1 -1 1 0 0 +1 2 -1 1 0 0 +2 3 -1 1 0 0 +3 0 -1 1 0 0 +[NumArcSegments] = 0 +[NumHoles] = 0 +[NumBlockLabels] = 1 +5 5 2 0.1 1 0 0 100 0 +"#; + +#[test] +fn parses_fixture_geometry() { + let doc = FemmDoc::parse(FIXTURE).expect("parse"); + assert_eq!(doc.frequency, 60.0); + assert_eq!(doc.nodes.len(), 4); + assert_eq!(doc.segments.len(), 4); + assert_eq!(doc.arcs.len(), 0); + assert_eq!(doc.block_labels.len(), 1); + assert_eq!(doc.materials.len(), 2); + assert_eq!(doc.boundaries.len(), 1); + assert_eq!(doc.points.len(), 1); + assert_eq!(doc.circuits.len(), 1); + assert_eq!(doc.nodes[0].boundary_marker, "A=0"); + assert_eq!(doc.segments[0].boundary_marker, "outer"); + assert_eq!(doc.block_labels[0].block_type, "Copper"); + assert_eq!(doc.block_labels[0].in_circuit, "Coil"); +} + +#[test] +fn round_trips_parse_write_parse() { + let a = FemmDoc::parse(FIXTURE).expect("parse a"); + let text = a.write(); + let b = FemmDoc::parse(&text).expect("parse b"); + + assert_eq!(a.frequency, b.frequency); + assert_eq!(a.precision, b.precision); + assert_eq!(a.depth, b.depth); + assert_eq!(a.nodes.len(), b.nodes.len()); + assert_eq!(a.segments.len(), b.segments.len()); + assert_eq!(a.materials.len(), b.materials.len()); + + for (x, y) in a.nodes.iter().zip(b.nodes.iter()) { + assert!((x.x - y.x).abs() < 1e-12); + assert!((x.y - y.y).abs() < 1e-12); + assert_eq!(x.boundary_marker, y.boundary_marker); + } + for (x, y) in a.segments.iter().zip(b.segments.iter()) { + assert_eq!(x.n0, y.n0); + assert_eq!(x.n1, y.n1); + assert_eq!(x.boundary_marker, y.boundary_marker); + } + for (x, y) in a.materials.iter().zip(b.materials.iter()) { + assert_eq!(x.name, y.name); + assert!((x.mu_x - y.mu_x).abs() < 1e-12); + assert!((x.cduct - y.cduct).abs() < 1e-12); + } + assert_eq!(a.block_labels.len(), b.block_labels.len()); + let la = &a.block_labels[0]; + let lb = &b.block_labels[0]; + assert!((la.x - lb.x).abs() < 1e-12); + assert_eq!(la.block_type, lb.block_type); + assert_eq!(la.in_circuit, lb.in_circuit); + assert_eq!(la.turns, lb.turns); + assert!((la.max_area - lb.max_area).abs() < 1e-10); +} diff --git a/crates/femm-sys/Cargo.toml b/crates/femm-sys/Cargo.toml new file mode 100644 index 0000000..5cb0024 --- /dev/null +++ b/crates/femm-sys/Cargo.toml @@ -0,0 +1,14 @@ +[package] +name = "femm-sys" +version = "0.0.1" +edition.workspace = true +rust-version.workspace = true +publish.workspace = true +description = "raw FFI bindings to the four FEMM solver engines (mag/elec/heat/curr)" +links = "femm_engines" + +[lib] +path = "src/lib.rs" + +[build-dependencies] +bindgen = "0.72" diff --git a/crates/femm-sys/build.rs b/crates/femm-sys/build.rs new file mode 100644 index 0000000..1ce31b3 --- /dev/null +++ b/crates/femm-sys/build.rs @@ -0,0 +1,55 @@ +//! builds the engine archives and generates FFI bindings for femm-sys. + +use std::env; +use std::path::PathBuf; +use std::process::Command; + +fn main() { + let manifest = PathBuf::from(env::var("CARGO_MANIFEST_DIR").unwrap()); + let repo_root = manifest.parent().unwrap().parent().unwrap().to_path_buf(); + let ffi_dir = repo_root.join("ffi"); + let build_dir = repo_root.join("build").join("ffi"); + + // re-run on changes to any FFI header, build script, or solver source dir. + println!("cargo:rerun-if-changed={}", ffi_dir.display()); + println!("cargo:rerun-if-changed={}", repo_root.join("scripts/macos/build_ffi.sh").display()); + for engine in ["fkn", "belasolv", "csolv", "hsolv", "liblua", "compat"] { + println!("cargo:rerun-if-changed={}", repo_root.join(engine).display()); + } + + let target_os = env::var("CARGO_CFG_TARGET_OS").unwrap(); + let script = match target_os.as_str() { + "macos" => repo_root.join("scripts/macos/build_ffi.sh"), + other => panic!("femm-sys: no engine build script wired up for target_os={other}"), + }; + + let status = Command::new("bash") + .arg(&script) + .status() + .expect("failed to invoke bash for engine build script"); + if !status.success() { + panic!("engine build script exited with {status}"); + } + + println!("cargo:rustc-link-search=native={}", build_dir.display()); + for lib in ["femm_mag", "femm_elec", "femm_heat", "femm_curr"] { + println!("cargo:rustc-link-lib=static={lib}"); + } + println!("cargo:rustc-link-lib=dylib=c++"); + + let bindings = bindgen::Builder::default() + .header(manifest.join("wrapper.h").to_str().unwrap()) + .clang_arg(format!("-I{}", ffi_dir.display())) + .allowlist_function("femm_.*") + .allowlist_type("Femm.*") + .allowlist_var("FEMM_.*") + .derive_default(true) + .derive_copy(true) + .derive_debug(true) + .layout_tests(false) + .generate() + .expect("bindgen failed to generate FFI bindings"); + + let out = PathBuf::from(env::var("OUT_DIR").unwrap()).join("bindings.rs"); + bindings.write_to_file(&out).expect("failed to write bindings.rs"); +} diff --git a/crates/femm-sys/src/lib.rs b/crates/femm-sys/src/lib.rs new file mode 100644 index 0000000..e7ac224 --- /dev/null +++ b/crates/femm-sys/src/lib.rs @@ -0,0 +1,5 @@ +//! raw bindgen output for the four FEMM solver C ABIs. + +#![allow(non_upper_case_globals, non_camel_case_types, non_snake_case, dead_code)] + +include!(concat!(env!("OUT_DIR"), "/bindings.rs")); diff --git a/crates/femm-sys/wrapper.h b/crates/femm-sys/wrapper.h new file mode 100644 index 0000000..184f52c --- /dev/null +++ b/crates/femm-sys/wrapper.h @@ -0,0 +1,5 @@ +// umbrella header for bindgen. +#include "femm_mag.h" +#include "femm_elec.h" +#include "femm_heat.h" +#include "femm_curr.h"