pimsr-forward

Forward modeling engine for the PIMSR (Physics-Informed Multi-modal Subsurface Reconstruction) project.

Turns stochastic geology models from pimsr-geogen into realistic geophysical observables:

• Magnetotellurics (MT), 2D TE+TM (mt2d.py) — SimPEG-based finite-volume forward for x-z sections invariant along strike y. The physical convention is TE = Ey/Hx = Zyx (Simulation2DMagneticField) and TM = Ex/Hy = Zxy (Simulation2DElectricField). The default survey has 12 stations x 8 frequencies. It is validated against the 1D analytic response to 0.13% (TE) / 0.4% (TM) median error on layered models.
• Magnetotellurics (MT), 1D analytic — exact impedance recursion for layered media (Wait recursion). Apparent resistivity + phase over a configurable period band.
• Gravity — vertical attraction of right rectangular prisms (Nagy 1966 closed form), evaluated on a surface profile.
• Sensor / noise model — multiplicative apparent-resistivity scatter, Gaussian phase noise, AR(1)-correlated galvanic distortion calibrated on real USArray residuals, per-station static shift, and per-mode severity. The stronger calibrated tail is applied to local TM/Zxy: the historical raw geographic Zyx curves from the approximately east-west field profiles rotate to this physical mode.
• Dataset builders — 1D per-station (dataset.py), 2D section (dataset2d.py, strict shard merging), and per-sample 3D volume (dataset3d.py) HDF5 training sets for pimsr-inversion. The explicit contracts are pimsr-mt-gravity-1d v2, pimsr-mt-2d v2 and pimsr-mt-3d v1. They record physical units, the modulo-180 phase convention, exact mode and scenario order, generation seeds/sample identities, source-artifact digests, forward/noise contracts and geometry. Files are fully validated under a .part name and published atomically; existing destinations are never overwritten. Legacy files without the complete physics and provenance contract are intentionally rejected by new consumers.

Physics validation

The 1D MT kernel is cross-checked against half-space and two-layer analytic limits in tests (tests/test_mt1d.py), and against SimPEG's 1D MT simulation in pimsr-benchmarks. The 2D TE/TM forwards are validated against the 1D analytic response on layered sections (tests/test_mt2d.py).

The gravity kernel is checked against the infinite-slab Bouguer limit and prism symmetry identities (tests/test_gravity.py).

Install

Python 3.11 or newer is required. SimPEG 0.25's discretize>=0.12 dependency and the resumable 3D worker lifecycle both require Python 3.11.

pip install -e .

Usage

from pimsr_forward import mt1d_response, prism_gz, SensorModel

rho_app, phase = mt1d_response(resistivities, thicknesses, periods)

Strict EMTF field assembly is shared by inversion and benchmarks. Use assemble_station_profile_modes for physical field comparisons: it rotates the complete tensor into the fitted profile frame and preserves the native projected station coordinates and station count. assemble_profile_modes additionally maps a profile onto a supplied synthetic station grid; its output is explicitly tagged publishable_physical_geometry=False and includes the horizontal compression factor.

from pimsr_forward import assemble_station_profile_modes

observations = assemble_station_profile_modes(
    "data/emtf", frequencies, profile_ids=["MTH15", "MTH16", "WYYS1"]
)

Dataset generation CLI:

# 1D per-station dataset
pimsr-forward-dataset --geology geology.h5 --out dataset.h5 --seed 42
# 2D TE+TM section dataset (requires the mt2d extra: pip install -e ".[mt2d]")
pimsr-forward-dataset2d --n 1000 --seed 7 --out ds2d.h5
# resumable 3D sample directory; matching completed samples are validated and skipped
pimsr-forward-dataset3d --out data/3d --start 0 --count 100 --seed 7 --workers 4

For parallel generation, write disjoint shards with the same --seed and non-overlapping --start-index ranges, then merge them with merge_shards. The merge rejects duplicate indices, mixed generation settings, incomplete shards and geometry mismatches before creating its output.

License

MIT