EM Solver Evaluation for Nagelfluh Integration

Date: 2026-06-13 — Updated with deep comparative research

Context: What Nagelfluh Already Has

Nagelfluh has full SimPEG support for the following workflows (all production-ready):

AEM (TDEM)

• System format: GEX (SkyTEM proprietary). Dual-moment waveform, tilt, gate factors all read from GEX.
• Inversion type: deterministic, gradient-based, fixed layer count, produces point estimates only.
• Regularisation: laterally constrained L2 (smooth) + optional IRLS (sparse/blocky). No uncertainty output.
• AEM systems: SkyTEM only (GEX-driven). No native support for Tempest, VTEM, RESOLVE, or other systems.

Magnetics

Equivalent source gridding + full 3D OcTree susceptibility/MVI inversion, all SimPEG-based.

Integration Model

Every Nagelfluh process type — regardless of the solver's implementation language — follows the same wrapper pattern:

class invert_tem_foo:
    @classmethod
    def run(cls, storage_context=None, **kwargs):
        with tempfile.TemporaryDirectory() as tmp:
            # 1. Download inputs from S3/MinIO via fsspec
            download_via_fsspec(kwargs['input_data'], tmp, storage_context)
            # 2. Convert to solver's native format
            write_input_files(tmp, kwargs)
            # 3. Run the solver as a subprocess
            subprocess.run(['solver_binary', '--input', tmp, '--output', tmp], check=True)
            # 4. Convert outputs back to Aarhus XYZ / msgpack
            result_xyz = read_output_files(tmp)
            # 5. Upload to storage
            dataset_url = write_dataset(result_xyz, storage_context)
        return {'status': 'success', 'outputs': {'result': dataset_url}}

This means implementation language is not a criterion. C++, Fortran, Julia, or any other compiled language is equally feasible as long as the binary can be included in the Docker image. The real integration cost drivers are:

• Format conversion effort: how different is the solver's native I/O from Aarhus XYZ/GEX?
• Docker build complexity: does it need a compiler, special libs (FFTW, MPI, FEniCS, Julia runtime)?
• Licence: does it permit use in a commercial SaaS context?

Gap Analysis

The evaluation question for each remaining solver is: what new capability does it bring to Nagelfluh that SimPEG does not already provide?

Solvers that fill none of these gaps — because SimPEG already covers their method well — are ranked at the bottom regardless of their standalone quality.

Ranked Evaluation

Tier 1 — High marginal value: fills genuine gaps

1. GeoBIPy

Gap filled: G1 (Bayesian inversion), G2 (dual-moment TDEM)

What it adds beyond SimPEG:

GeoBIPy uses reversible-jump Markov Chain Monte Carlo (rj-MCMC), the same Bayesian algorithm used by GA-AEM's stochastic mode and HiQGA's MCMC mode. It does not produce a single resistivity model — it samples the full posterior probability distribution P(model | data), with the number of layers as a free parameter. Outputs include:

• Resistivity vs. depth posterior PDFs: a probability image showing not just the best model but how certain the data are at each depth
• Interface-depth probability: reveals where layer boundaries are well-constrained vs. ambiguous
• Layer-count posterior: shows whether the data support 2, 3, or 5 layers
• Credible intervals at 68% and 95% level — equivalent to ±1σ and ±2σ but fully nonlinear

This is the uncertainty quantification differentiator for Nagelfluh. A user who runs SimPEG's LCI gets a smooth model. A user who runs GeoBIPy on the same data gets the same model plus a probability map that tells them which parts of the model are geologically reliable.

GeoBIPy processes soundings independently (no spatial coupling), which is a trade-off vs. SimPEG LCI. Spatial coherence comes from data density; lateral constraints are an active research feature. Computationally it requires MPI and scales via embarrassingly parallel MCMC chains — designed for HPC clusters. Not suitable for realtime turnaround; suited for post-survey detailed analysis of priority areas.

Used by USGS for national AEM programs (Alaska permafrost, Great Plains aquifers). The Python API makes it the easiest Bayesian AEM solver to wrap as a Nagelfluh process type.

Integration path: Write a libaarhusxyz → GeoBIPy HDF5 adapter; wrap as invert_tem_bayesian process with MPI-enabled K8s pod (multi-container job). System parameters entered via TdemSystem constructor (similar to GEX-to-SimPEG mapping already done).

2. empymod (as standalone FDEM forward process)

Gap filled: G2 (FDEM AEM: RESOLVE, Dighem, HEM)

What it adds beyond SimPEG:

Nagelfluh uses SimPEG's TDEM stack, which has no configured process type for frequency-domain helicopter AEM systems: RESOLVE (CGG), Dighem (Aerodat), HEM (Geotech). These measure in-phase and quadrature response at multiple frequencies and are widely used for base metal exploration and groundwater. empymod models all RESOLVE-type coil geometries (HCP, coaxial) and frequencies directly and is the de-facto standard forward engine for 1D FDEM AEM.

The v2.6.0 addition of arbitrary waveform support also makes it directly competitive with SimPEG for TDEM forward modeling — and about 100–1000× faster per sounding, because it uses a semi-analytical DLF approach rather than sparse matrix solvers. This matters as a standalone "forward model" process type for synthetic data generation, sensitivity testing, and as the inner engine of a custom inversion workflow.

empymod is the lowest integration cost of any solver on this list and the highest utility-to-effort ratio. A standalone forward_1d_fdem process type serving RESOLVE-type users can be written in under 100 lines.

3. GA-AEM

Gap filled: G2 (multi-system: Tempest, VTEM, GeoTEM, SPECTREM + stochastic mode) — with significant licence caveat

What it adds beyond SimPEG:

System coverage: GA-AEM has validated .stm system files for Tempest, VTEM, GeoTEM, and SPECTREM — the commercial systems widely used outside the SkyTEM ecosystem (particularly in Australia, Canada, South Africa). These systems cannot be run directly through Nagelfluh's GEX-based invert_tem without a GEX equivalent. GA-AEM's .stm format is the natural bridge.

ASEGGDF2 format: The international AEM exchange standard used by Geoscience Australia, CSIRO, and most survey contractors. Users with ASEGGDF2 archives cannot currently feed data into Nagelfluh without conversion.

Stochastic mode (garjmcmctdem): Same rj-MCMC algorithm as GeoBIPy but the GA reference implementation for Tempest. Validated against decades of Australian national AEM surveys.

Critical limitation: GPL v2 licence. Building a commercial SaaS on GPL v2 code is legally problematic without releasing the wrapper. Options: (a) negotiate a commercial licence from Geoscience Australia, or (b) take legal advice on whether a subprocess-isolated container is a "separate program" under GPL v2.

Integration path: Build Docker image with GA-AEM compiled (FFTW + optional MPI). Write Aarhus XYZ ↔ ASEGGDF2 adapter. Standard subprocess wrapper. The real cost is the format adapter and the licence resolution, not the language. Unlocks every operator using Tempest, VTEM, or ASEGGDF2 archives.

3. HiQGA.jl

Gap filled: G1 (MCMC inversion), G2 (multi-system), with more efficient MCMC

What it adds beyond SimPEG and beyond GeoBIPy:

HiQGA is the most algorithmically advanced of the three stochastic solvers. Its gradient MCMC (NUTS — No U-Turn Sampler / Hamiltonian MC) uses gradient information (the Jacobian) to propose more efficient MCMC steps. This can reduce the number of forward calls needed to characterise the posterior by 10–100× compared to random-walk rj-MCMC (GeoBIPy / GA-AEM). For surveys where per-sounding MCMC is computationally expensive, HiQGA's NUTS approach may make Bayesian inversion practically feasible on moderate hardware rather than requiring an HPC cluster.

Additional value: HiQGA covers the widest range of inversion strategies in a single package (deterministic, fixed-D MCMC, trans-D MCMC, gradient MCMC), enabling comparative analysis.

The integration challenge is Docker image size (Julia runtime adds ~500 MB) and precompilation, but both are solvable at image build time. MIT licence removes the GA-AEM licence risk entirely.

4. P223 Suite (LeroiAir / ArjunAir)

Gap filled: G3 (3D target body forward modeling)

What it adds beyond SimPEG:

SimPEG's Simulation1DLayeredStitched assumes the earth is horizontally layered and laterally slowly varying. This is an excellent model for sedimentary basins, aquifer systems, and permafrost. However, it fails fundamentally for discrete 3D conductors such as:

• Massive sulfide ore bodies (VMS, Ni-Cu deposits)
• Conductive graphite horizons that dip steeply
• Anthropogenic targets (buried tanks, pipes)

When a SkyTEM or VTEM survey flies over a dipping conductive plate, the 1D LCI model interprets it as a conductive layer smeared horizontally — wrong depth, wrong geometry, wrong conductance estimate. LeroiAir models the conductor as a parameterized plate (strike, dip, plunge, depth, conductance, along-strike and down-dip extent) and computes the correct 3D EM response of the plate in a layered background. This is the standard tool for target characterisation in mineral exploration AEM.

This is not an inversion tool but a forward modeling tool. The workflow is: run SimPEG LCI to get the background layered model → identify anomalies in the LCI residuals → use LeroiAir to forward-model the anomaly with various plate geometries → match the observed data to determine plate parameters.

Integration: Subprocess wrapper + ASCII control file generator from Nagelfluh parameters. The most natural process type would be forward_3d_plate that takes a plate parameterization (from the UI) + background LCI model + survey geometry and produces synthetic AEM data for comparison.

5. EMagPy

Gap filled: G4 (ground FDEM coil surveys — entirely different instrument class from AEM)

What it adds beyond SimPEG:

EMagPy targets ground-deployed frequency-domain EM coil sensors: GF Instruments CMD series, Dualem 1S/2S/21S/42S, Geonics EM31/EM38/EM34. These are dragged or hand-carried sensors operating at low induction number, widely used for:

• Precision agriculture (soil salinity/clay mapping)
• Shallow environmental contamination surveys
• Archaeological geophysics
• Near-surface geological mapping

SimPEG does have a 1D FDEM simulation, but it has no dedicated processing pipeline for these instruments, no CSV input format for the specific column naming conventions used by CMD/Dualem exporters, and no automated sensor-geometry handling (coil separation, height above ground, orientation). EMagPy packages all of this in a clean Python API.

This opens a different user community from AEM operators: environmental consultants, agronomists, archaeologists. Ground FDEM surveys are far more numerous than AEM surveys (much lower entry cost), and a cloud-based processing service for them is an underserved market.

License note: GPL-3.0. The GPL's network use provision (Affero GPL) does not apply here since it's plain GPL-3.0, but using GPL libraries in a SaaS product should be reviewed with counsel. The process type container runs GPL code in isolation (not linked to the Nagelfluh codebase), which substantially reduces the concern.

6. ModEM

Gap filled: G5 (3D MT with EDI format compatibility)

What it adds beyond SimPEG:

SimPEG does have a Natural Source EM (NSEM) module for 3D MT, but it has known boundary condition issues (GitHub #1767) and is much less actively maintained than SimPEG's AEM modules. More critically, SimPEG does not speak EDI. EDI (Electrical Data Interchange) is the universal MT data format produced by every MT acquisition system and processed by every MT processing software (BIRRP, EMTF, Metronix, Phoenix). A MT practitioner's workflow starts with EDI files. An integration that can't read EDI can't serve MT users.

ModEM is the de facto standard 3D MT inversion code in the research and exploration community. Its output format (plain-text model files) is understood by visualisation tools like ModEM-gui, MTPy, and GMT. GPU support (CUDA/HIP) makes it competitive for large 3D models.

The case for ModEM integration: AEM surveys are routinely combined with MT surveys for complementary depth coverage. AEM resolves conductivity structure to ~300–500 m depth; MT extends this to several km. An operator doing both surveys in one field campaign naturally wants one platform. If Nagelfluh handles AEM and refers the MT data elsewhere, that's friction. If it handles both, it captures the whole multi-method workflow.

Integration path: Subprocess wrapper around the ModEM executable; write Python utilities to translate from project data to EDI input, run Mod3DMT, and parse the text output. pyModEM helps with the file I/O. Not trivial but well-defined.

Tier 2 — Moderate marginal value

8. emg3d

Gap filled: G6 (3D CSEM forward modeling)

9. PyGimli

Gap filled: ERT/IP surveys + multi-method joint inversion

10. P223 ArjunAir

Gap filled: G3 (2.5D profiling — complement to LeroiAir)

Lower priority than LeroiAir. ArjunAir is a more computationally expensive 2.5D FEM and is used for profiling surveys where LeroiAir's plate model is insufficient (e.g., gradational 2D conductors). It uses the same P223 Fortran codebase and would be packaged alongside LeroiAir naturally. Separately, not worth the integration effort.

11. TEM1D

Marginal gap: independent 1D TDEM reference / AarhusInv-lineage benchmarking

Tier 3 — Low or no marginal value over SimPEG

Solvers in this tier are either: (a) covered by SimPEG's existing AEM/EM/magnetics support (all could in principle be wrapped via subprocess, but doing so adds no capability); (b) MATLAB-only (MATLAB is not in the Docker image and is not freely licensable, making these genuinely unusable); (c) inactive or restrictively licensed; or (d) address methods — marine CSEM, academic 3D MT — outside the Nagelfluh user base. Language alone is not a reason for low ranking.

Summary: Ranking by Marginal Value

All solvers are integrated via the same subprocess wrapper pattern (fsspec download → subprocess → fsspec upload). Integration effort reflects format conversion complexity and Docker build complexity, not programming language.

¹ GeoBIPy Python codebase: USGS public domain. One Fortran component: GPL-2. The GPL component may be replaceable with scipy equivalents.

Recommended Integration Sequence

Phase 1 — Low-hanging fruit (trivial effort, real gaps)

empymod (forward_fdem_1d, forward_tdem_1d) — pure Python, Apache-2.0, zero new dependencies. Write one process type class to expose 1D FDEM forward modeling for RESOLVE/Dighem/HEM users and fast synthetic TDEM data generation. Validates the "forward model without inversion" process pattern. Estimated: 1–2 days.

GeoBIPy (invert_tem_bayesian) — Python-native, USGS public domain. Write the libaarhusxyz → GeoBIPy HDF5 format adapter, deploy an MPI-enabled K8s pod. The uncertainty quantification output (posterior PDFs, credible interval maps) is the clearest premium-tier feature on this list — directly differentiates from "another AEM platform". Estimated: 1–2 weeks.

Phase 2 — Unlock non-SkyTEM operators (highest user growth potential)

GA-AEM (invert_tem_gaaem, invert_tem_stochastic) — resolve licence first: either negotiate a commercial licence from Geoscience Australia, or seek legal opinion on subprocess isolation under GPL v2. Once resolved: write ASEGGDF2 ↔ Aarhus XYZ adapter, build Docker image with FFTW + compiled GA-AEM binaries. Unlock Tempest, VTEM, GeoTEM operators. Estimated: 3–4 weeks after licence resolution.

Phase 3 — Mineral exploration workflows

P223 LeroiAir (forward_3d_plate) — subprocess Fortran binary, CFL control file generator from UI parameters. Targets hard-rock mineral exploration (sulfide detection, plate characterisation). Estimated: 1–2 weeks.

ModEM (invert_mt_3d) — subprocess C++ binary, EDI format adapter. Covers the MT leg of multi-method mineral exploration surveys. Estimated: 2–3 weeks.

Phase 4 — Ecosystem expansion

HiQGA.jl (invert_tem_hiqga) — MIT licence, most advanced MCMC algorithm. Build Julia runner base image. Requires data format adapter and juliacall or subprocess bridge. Estimated: 2–3 weeks.

EMagPy (invert_fdem_ground) — GPL-3.0 requires review; otherwise trivial to integrate. Opens agriculture/environment/archaeology market. Estimated: 3–5 days after licence clearance.

Data Format Compatibility Notes

A shared nagelfluh_aem_adapters Python package containing these format converters would avoid duplication across process types and is worth building before Phase 2.