From 8c262b7b697e26b0902b9b32be2793f924daa185 Mon Sep 17 00:00:00 2001 From: Chad <167274875+UnsignedChad@users.noreply.github.com> Date: Sun, 24 May 2026 22:23:26 -0400 Subject: [PATCH 1/2] emikit: openEMS cross-check scripts (WIP -- needs FDTD tuning) infrastructure for an apples-to-apples full-wave comparison against emikits closed form on the same single-trace fixture. scripts work, openEMS builds and runs, port voltage data comes out clean. but the default Gaussian + auto -40 dB end criterion terminates after ~2 ns which is too short for a 30 MHz spectrum. committing as research artifact + reproducible methodology. final numbers need more tuning -- either disable EndCriteria + extend NrTS, or run CW per frequency instead of broadband. the README documents where the experiment stands. the textbook calibration tests already pin emikits formulas to Ott Eq 11-2 and Hockanson 1996 at 0.05 dB. an FDTD cross-check would confirm those formulas are faithful reductions of the full Maxwell problem -- which is a known result in EMC literature, not the foundation of our validation chain. --- emikit/openems/README.md | 47 ++++++++++ emikit/openems/compare_openems_emikit.py | 64 ++++++++++++++ emikit/openems/openems_microstrip.py | 104 +++++++++++++++++++++++ 3 files changed, 215 insertions(+) create mode 100644 emikit/openems/README.md create mode 100755 emikit/openems/compare_openems_emikit.py create mode 100755 emikit/openems/openems_microstrip.py diff --git a/emikit/openems/README.md b/emikit/openems/README.md new file mode 100644 index 0000000..c888125 --- /dev/null +++ b/emikit/openems/README.md @@ -0,0 +1,47 @@ +# OpenEMS cross-check scripts + +Apples-to-apples full-wave FDTD comparison against emikit's closed +form for the **same** single-trace fixture. If emikit and OpenEMS +agree to a few dB on the geometry emikit models, the algorithm is +right and the chamber-vs-emikit gap on the TI EVM is purely scope +(multi-net, activity factor, cable common-mode). + +## Status: infrastructure ready, FDTD parameters need iteration + +* `openems_microstrip.py` -- builds the fixture (30 mm trace, 0.15 mm + wide, 0.2 mm above solid GND, lumped ports + NF2FF at 10 m) +* `compare_openems_emikit.py` -- consumes the FDTD output and runs + emikit's `loop_e_field` for the same I_port spectrum + +Run on a workstation with OpenEMS Python bindings installed: + + source ~/opt/openEMS/venv/bin/activate + python3 openems_microstrip.py # produces full_result.h5 + python3 compare_openems_emikit.py + +## Known issues to fix before getting useful numbers + +1. The Gaussian excitation auto-terminates via the -40 dB energy + criterion at ~2 ns, giving only ~17 port samples over the run. + That is below the 33 ns needed to resolve 30 MHz cleanly. Either: + - bump `NrTS` past the natural decay and disable EndCriteria, or + - run multiple CW simulations at fixed frequencies (more robust). + +2. Lumped port has shown numerical sensitivity to mesh near the port + terminals. Switching to `AddMSLPort` (microstrip-line port, + designed for this geometry) is worth trying first. + +3. Mesh budget kept under 2M cells with explicit fine lines near the + trace + max_res=2.0 mm out to absorbing boundaries. Going finer + blows memory past 30 GB. + +## Why this is still worth committing + +Even without final numbers, the scripts capture the methodology so the +experiment is reproducible by anyone with OpenEMS installed. The +closed-form algorithm has already been pinned to its textbook source +(Ott Eq 11-2, Hockanson 1996) at 0.05 dB; what an FDTD cross-check +would prove is that those textbook formulas are themselves a faithful +reduction of the full Maxwell problem at this scale. That is a known +result in the EMC literature -- the experiment is confirmation, not +the foundation of the validation chain. diff --git a/emikit/openems/compare_openems_emikit.py b/emikit/openems/compare_openems_emikit.py new file mode 100755 index 0000000..b44600b --- /dev/null +++ b/emikit/openems/compare_openems_emikit.py @@ -0,0 +1,64 @@ +#!/usr/bin/env python3 +""" +Feed openEMS's measured port-current spectrum into emikit's closed-form +loop_e_field and compare to openEMS's full-wave far-field. + +If the two land within a few dB, emikit's algorithm is correct for the +geometry it models -- and the entire 30-70 dB chamber gap on the TI +EVM is purely scope (multi-net, activity factor, cable common-mode +which is modeled separately). +""" +import math +import h5py +import numpy as np + +ETA0 = 376.730313668 +C0 = 2.99792458e8 + +# Geometry (must match openems_run.py). +TRACE_LEN_M = 30.0e-3 +LOOP_HEIGHT_M = 0.2e-3 +OBS_DIST_M = 10.0 +LOOP_AREA = TRACE_LEN_M * LOOP_HEIGHT_M + + +def emikit_loop_e_field(i_a, freq_hz, distance_m, area_m2): + """Ott Eq 11-2 -- same as emikit's LoopEmissions::loop_e_field.""" + return (ETA0 * math.pi * i_a * area_m2 * freq_hz ** 2) / (C0 ** 2 * distance_m) + + +def main(): + with h5py.File("/home/chad/openems-runs/emikit_microstrip_air_result.h5", + "r") as h: + freqs = h["freq_hz"][...] + i_port = h["i_port_a"][...] + e_full = h["e_v_m"][...] + e_full_dbuv = h["e_dbuv_m"][...] + + print(f"{'freq_MHz':>9} {'|I|mA':>8} {'E_FW_dBuV':>11} " + f"{'E_emikit_dBuV':>14} {'delta_dB':>10}") + print("-" * 60) + + e_emikit_all = [] + deltas = [] + for f, ip, efw_v, efw_db in zip(freqs, i_port, e_full, e_full_dbuv): + e_emikit_v = emikit_loop_e_field(abs(ip), f, OBS_DIST_M, LOOP_AREA) + if e_emikit_v <= 0: + e_emikit_db = -1000.0 + else: + e_emikit_db = 20.0 * math.log10(e_emikit_v * 1e6) + e_emikit_all.append(e_emikit_db) + delta = e_emikit_db - efw_db + deltas.append(delta) + print(f"{f / 1e6:9.1f} {abs(ip) * 1e3:8.4f} {efw_db:11.2f} " + f"{e_emikit_db:14.2f} {delta:+10.2f}") + + valid = [d for d in deltas if abs(d) < 200] + if valid: + print(f"\nmean |delta| = {np.mean(np.abs(valid)):.2f} dB") + print(f"max |delta| = {np.max (np.abs(valid)):.2f} dB") + print(f"rms = {np.sqrt(np.mean(np.array(valid) ** 2)):.2f} dB") + + +if __name__ == "__main__": + main() diff --git a/emikit/openems/openems_microstrip.py b/emikit/openems/openems_microstrip.py new file mode 100755 index 0000000..8458117 --- /dev/null +++ b/emikit/openems/openems_microstrip.py @@ -0,0 +1,104 @@ +#!/usr/bin/env python3 +"""Full-fixture FDTD with sane mesh grading.""" +import os +import numpy as np +import h5py +from openEMS import openEMS +from CSXCAD import ContinuousStructure + +sim = os.path.expanduser("~/openems-runs/full") +os.makedirs(sim, exist_ok=True) + +TRACE_LEN = 30.0 +TRACE_W = 0.15 +TRACE_T = 0.035 +SUB_H = 0.2 +GND_HALF = 30.0 +AIR = 80.0 + +FDTD = openEMS(EndCriteria=1e-4, NrTS=60000) +FDTD.SetGaussExcite(500e6, 500e6) +FDTD.SetBoundaryCond(['MUR'] * 6) +CSX = ContinuousStructure() +FDTD.SetCSX(CSX) +CSX.GetGrid().SetDeltaUnit(1e-3) + +gnd = CSX.AddMetal("gnd") +gnd.AddBox([-GND_HALF, -GND_HALF, 0.0], [GND_HALF, GND_HALF, 0.0]) + +trace = CSX.AddMetal("trace") +trace.AddBox([-TRACE_LEN / 2, -TRACE_W / 2, SUB_H], + [ TRACE_LEN / 2, TRACE_W / 2, SUB_H + TRACE_T]) + +FDTD.AddLumpedPort(1, 50.0, + [-TRACE_LEN / 2, -TRACE_W / 2, 0.0], + [-TRACE_LEN / 2, TRACE_W / 2, SUB_H], + 'z', excite=True) +FDTD.AddLumpedPort(2, 50.0, + [ TRACE_LEN / 2, -TRACE_W / 2, 0.0], + [ TRACE_LEN / 2, TRACE_W / 2, SUB_H], + 'z', excite=False) + +# Explicit mesh lines near the trace + coarse at boundary, then smooth +# with a generous max_res (2 mm) so the mesh GROWS away from the fine +# lines instead of being held at 0.05 mm everywhere. +mesh = CSX.GetGrid() + +# X: coarse end-of-domain + fine ends-of-trace + middle filler. +x_lines = [-GND_HALF - AIR, -GND_HALF, -TRACE_LEN/2 - 1.0, + -TRACE_LEN/2 - 0.2, -TRACE_LEN/2, + -TRACE_LEN/4, 0.0, TRACE_LEN/4, + TRACE_LEN/2, TRACE_LEN/2 + 0.2, + TRACE_LEN/2 + 1.0, GND_HALF, GND_HALF + AIR] +mesh.AddLine('x', x_lines) +mesh.SmoothMeshLines('x', 2.0, ratio=1.4) + +# Y: tight cluster around trace, then coarse out. +y_lines = [-GND_HALF - AIR, -GND_HALF, -2.0, -0.5, -TRACE_W, + -TRACE_W/2, -TRACE_W/4, 0.0, TRACE_W/4, TRACE_W/2, + TRACE_W, 0.5, 2.0, GND_HALF, GND_HALF + AIR] +mesh.AddLine('y', y_lines) +mesh.SmoothMeshLines('y', 2.0, ratio=1.4) + +# Z: fine through dielectric+trace, coarse out. +z_lines = [-AIR, -10.0, -1.0, 0.0, SUB_H/4, SUB_H/2, 3*SUB_H/4, + SUB_H, SUB_H + TRACE_T/2, SUB_H + TRACE_T, + SUB_H + 0.2, 1.0, 10.0, SUB_H + AIR] +mesh.AddLine('z', z_lines) +mesh.SmoothMeshLines('z', 2.0, ratio=1.4) + +nx, ny, nz = (len(mesh.GetLines(a)) for a in 'xyz') +cells = nx * ny * nz +print(f"Mesh: {nx} x {ny} x {nz} = {cells:,} cells", flush=True) +assert cells < 30_000_000, f"mesh too big ({cells:,}), tune ratios" + +nf2ff = FDTD.CreateNF2FFBox() + +FDTD.Run(sim, verbose=2, cleanup=False) + +freqs = np.linspace(30e6, 1.0e9, 40) +ff = nf2ff.CalcNF2FF(sim, freqs, np.array([90.0]), np.array([0.0]), + radius=10.0) +e_ff = np.abs(ff.E_norm[0])[:, 0, 0] + +ut = np.loadtxt(os.path.join(sim, 'port_ut_1')) +it = np.loadtxt(os.path.join(sim, 'port_it_1')) +dt = ut[1, 0] - ut[0, 0] +t = ut[:, 0] +i_t = it[:, 1] +i_f = np.array([np.sum(i_t * np.exp(-2j * np.pi * f * t)) * dt + for f in freqs]) +i_port = np.abs(i_f) + +out = os.path.expanduser("~/openems-runs/full_result.h5") +with h5py.File(out, "w") as h: + h.create_dataset("freq_hz", data=freqs) + h.create_dataset("i_port_a", data=i_port) + h.create_dataset("e_v_m", data=e_ff) + h.create_dataset("e_dbuv_m", data=20.0 * np.log10(e_ff * 1e6)) + +print(f"\n{'freq_MHz':>9} {'|I|mA':>9} {'E_dBuV/m':>10}", flush=True) +for k in range(len(freqs)): + print(f"{freqs[k] / 1e6:9.1f} {i_port[k] * 1e3:9.4f} " + f"{20.0 * np.log10(e_ff[k] * 1e6):10.2f}", flush=True) +print(f"\nwrote {out}", flush=True) From 1f352f38f9a2de2e3d2173807d2ad7f6ce9eaa2d Mon Sep 17 00:00:00 2001 From: Chad <167274875+UnsignedChad@users.noreply.github.com> Date: Sun, 24 May 2026 22:34:34 -0400 Subject: [PATCH 2/2] tighten readme and script docstrings --- emikit/openems/README.md | 63 ++++++++++-------------- emikit/openems/compare_openems_emikit.py | 11 +---- 2 files changed, 27 insertions(+), 47 deletions(-) diff --git a/emikit/openems/README.md b/emikit/openems/README.md index c888125..de1fa6c 100644 --- a/emikit/openems/README.md +++ b/emikit/openems/README.md @@ -1,47 +1,34 @@ -# OpenEMS cross-check scripts +# openEMS cross-check scripts -Apples-to-apples full-wave FDTD comparison against emikit's closed -form for the **same** single-trace fixture. If emikit and OpenEMS -agree to a few dB on the geometry emikit models, the algorithm is -right and the chamber-vs-emikit gap on the TI EVM is purely scope -(multi-net, activity factor, cable common-mode). +Full-wave FDTD against emikit's closed form on the same single-trace +fixture. If the two agree to a few dB on the geometry emikit models, +the chamber-vs-emikit gap on the TI EVM is scope (multi-net, activity +factor, cable CM), not math. -## Status: infrastructure ready, FDTD parameters need iteration +## status -* `openems_microstrip.py` -- builds the fixture (30 mm trace, 0.15 mm - wide, 0.2 mm above solid GND, lumped ports + NF2FF at 10 m) -* `compare_openems_emikit.py` -- consumes the FDTD output and runs - emikit's `loop_e_field` for the same I_port spectrum +scripts run, openEMS builds and produces port + near-field data, but +the broadband Gaussian + -40 dB auto-terminate stops at ~2 ns -- too +short for a clean 30 MHz spectrum. needs more iteration before the +numbers are usable. -Run on a workstation with OpenEMS Python bindings installed: +files: - source ~/opt/openEMS/venv/bin/activate - python3 openems_microstrip.py # produces full_result.h5 - python3 compare_openems_emikit.py - -## Known issues to fix before getting useful numbers +* `openems_microstrip.py` -- 30 mm trace, 0.15 mm wide, 0.2 mm above + solid GND, lumped ports + NF2FF at 10 m +* `compare_openems_emikit.py` -- runs emikit's `loop_e_field` over the + FDTD-measured I_port spectrum -1. The Gaussian excitation auto-terminates via the -40 dB energy - criterion at ~2 ns, giving only ~17 port samples over the run. - That is below the 33 ns needed to resolve 30 MHz cleanly. Either: - - bump `NrTS` past the natural decay and disable EndCriteria, or - - run multiple CW simulations at fixed frequencies (more robust). +run on a workstation with openEMS Python bindings: -2. Lumped port has shown numerical sensitivity to mesh near the port - terminals. Switching to `AddMSLPort` (microstrip-line port, - designed for this geometry) is worth trying first. - -3. Mesh budget kept under 2M cells with explicit fine lines near the - trace + max_res=2.0 mm out to absorbing boundaries. Going finer - blows memory past 30 GB. + source ~/opt/openEMS/venv/bin/activate + python3 openems_microstrip.py + python3 compare_openems_emikit.py -## Why this is still worth committing +## todo -Even without final numbers, the scripts capture the methodology so the -experiment is reproducible by anyone with OpenEMS installed. The -closed-form algorithm has already been pinned to its textbook source -(Ott Eq 11-2, Hockanson 1996) at 0.05 dB; what an FDTD cross-check -would prove is that those textbook formulas are themselves a faithful -reduction of the full Maxwell problem at this scale. That is a known -result in the EMC literature -- the experiment is confirmation, not -the foundation of the validation chain. +- bump `NrTS`, disable `EndCriteria`, or run CW per frequency instead + of broadband Gaussian +- try `AddMSLPort` instead of `AddLumpedPort` (less sensitivity near + the port terminals) +- mesh capped under 2M cells; tighter spacing blows past 30 GB diff --git a/emikit/openems/compare_openems_emikit.py b/emikit/openems/compare_openems_emikit.py index b44600b..5d6b9ea 100755 --- a/emikit/openems/compare_openems_emikit.py +++ b/emikit/openems/compare_openems_emikit.py @@ -1,13 +1,6 @@ #!/usr/bin/env python3 -""" -Feed openEMS's measured port-current spectrum into emikit's closed-form -loop_e_field and compare to openEMS's full-wave far-field. - -If the two land within a few dB, emikit's algorithm is correct for the -geometry it models -- and the entire 30-70 dB chamber gap on the TI -EVM is purely scope (multi-net, activity factor, cable common-mode -which is modeled separately). -""" +"""feeds openEMS's measured port-current spectrum into emikit's +loop_e_field and compares to the full-wave far-field.""" import math import h5py import numpy as np