Python API

The Python package under utils/ (importable as petgem) handles pre- and post-processing: it assembles the input bundle from the mesh, model, source, and receiver files, and reads the bundle and the kernel responses back into NumPy arrays. These are the functions used by utils/preprocess.py and by the per-case postprocess.py scripts. The file formats are documented in Data formats.

Preprocessing

petgem.runPreprocessing(*, mode, order, case_dir, mesh_filename, receiver_filename, source_filename=None, sigma_x, sigma_y, sigma_z, fixed_materials=(), input_filename='input.h5', params_filename='params.txt', im_source_filename=None, observed_filename=None, error_level=None, output_vtk=None, dm_view=False)

Shared preprocessing pipeline for the PETGEM forward / inverse kernels.

Pure library function - no CLI parsing. Builds the DMPlex, assigns the per-cell conductivity from sigma_{x,y,z} indexed by the gmsh:physical tag, parses the receivers and sources text files, and bundles everything into a single HDF5 input file (<case_dir>/<input_filename>). Also emits the matching params file.

Parameters

mode{‘forward’, ‘inverse’}

Selects which params-file template to emit.

orderint

Polynomial order, written into the params file.

case_dirstr

Directory containing the input data and where the bundle is written.

mesh_filename, receiver_filename, source_filenamestr

Filenames relative to case_dir.

sigma_x, sigma_y, sigma_zarray_like

Per-material conductivity. Index = 0-based material id (gmsh:physical - 1). Lengths must match.

input_filenamestr

Bundle filename inside case_dir. Default: input.h5

params_filenamestr

Filename of the PETGEM params file emitted inside case_dir. Default: params.txt

output_vtkstr or None

Optional VTU filename for the conductivity + materials_id view.

dm_viewbool

If True, dumps DMPlex info to stdout via PETSc options.

Readers

petgem.readBundle(filename)

Read the case-independent payload of a PETGEM input bundle HDF5.

Returns a dict with keys:

receivers : (N_recv, 3) ndarray of receiver positions (real-valued)
order      : int polynomial order
frequency : float - first transmitter frequency (Hz); for a single-
            frequency forward case this is the operating frequency
sources   : (N_src, 8) ndarray, columns =
            [freq, x, y, z, current, length, dipAngle, azimuthAngle]

The DMPlex / model_data fields inside the bundle are not returned - those are consumed by the C kernel via loadCsemInputs. This loader is for the Python postprocessing side and exposes only the parameters a case-specific validation script actually needs.

petgem.readResponses(filename, source=1)

Read one source’s responses from a PETGEM single-file output HDF5.

The fm.csem postprocessing writes a single file containing every transmitter under /sources/src{k}/ groups (k is 1-based). This helper returns the per-source dict for the requested transmitter, keeping the per-source result shape that pre-refactor callers used.

Parameters

filenamestr

Path to the responses HDF5 file.

sourceint, optional

1-based source index (default 1).

Returns

dict

Keys:

Ex, Ey, Ez, Hx, Hy, Hz : ndarrays (complex in PETSc complex builds)
source      : dict of /sources/src{k} attributes
              (frequency, x_pos, y_pos, z_pos, current, length,
               dip_angle, azimuth_angle)
provenance  : dict of root-level attributes
              (petgem_version, input_filename, date, order,
               mpi_tasks, num_sources, frequency)

Notes

Attributes are read via h5py (petsc4py’s attribute API is awkward for this read pattern). Vec components live under /sources/src{k}/fields/ in the new single-file layout.

petgem.readAllResponses(filename)

Read every source’s responses from a single PETGEM responses file.

Parameters

filenamestr

Path to the responses HDF5 file written by fm.csem postprocessing.

Returns

dict

Keys:

provenance : dict of root-level attributes (see ``readResponses``)
num_sources : int (mirrors provenance['num_sources'] for convenience)
sources    : dict keyed by 1-based source index; each value has the
             same shape as ``readResponses(filename, source=k)``
             (Ex..Hz arrays + 'source' attrs + 'provenance').

Notes

Source indices are discovered from root attribute num_sources; the underlying groups are named /sources/src1 … /sources/srcN.

petgem.readSigmaTable(path)

Read the per-material conductivity table (sigmas.txt).

Whitespace-delimited, one row per material (0-based row index = material id = gmsh:physical - 1), # comments and blank lines allowed:

# sigma_x sigma_y sigma_z [fixed]
0.1 0.1 0.1 1     # fixed in inversion (e.g. air, ocean)
1.0 1.0 1.0 0     # invertable
2.0 2.0 2.0       # 'fixed' column omitted -> defaults to 0

The 4th column (fixed) is optional and only meaningful for inverse modeling: a non-zero entry marks the material as held fixed during inversion (gradient zeroed, smoother self-only). Forward modeling ignores this column.

Returns (sigma_x, sigma_y, sigma_z, fixed_ids) where fixed_ids is a sorted list of 0-based material IDs flagged as fixed (empty list when the table has no 4th column).

petgem.readSourcesText(path)

Parse a transmitter file into a unified (N, 8) array with columns freq x y z current length dip azimuth (one row per transmitter).

Two on-disk layouts are accepted (auto-detected), so a single canonical representation feeds both forward and inverse preprocessing:

1. 8-field (canonical) - frequency on every row:

   freq x y z current length dip azimuth
   ...
   

2. Legacy forward - a lone frequency line followed by 7-field rows:

   <freq>
   x y z current length dip azimuth
   ...
   

   parsed into the same (N, 8) shape with that frequency on each row.

Forward modeling uses a single (repeated) frequency with one or more transmitters; inverse modeling carries one row per (frequency, dipole).

petgem.readObservedDataH5(path)

Read an observed-data HDF5 file produced by generate_observed_data.py or convert_invex_to_hdf5.py.

Returns (Ex, error_level) where Ex is the [N_freq, N_recv] complex128 ndarray stored at /Ex and error_level is the amplitude-relative noise level used to synthesize it (None if the source file has no error_level attribute). The compound HDF5 type {r: float64, i: float64} is materialized as complex128 by h5py.

petgem.readInvExDat(path)

Parse a MATLAB-style invEx.dat observed-data text file.

Format (one row per frequency):

freq_label  Re(Ex_1) Im(Ex_1)  Re(Ex_2) Im(Ex_2)  ...  Re(Ex_N) Im(Ex_N)

The leading frequency label is dropped (the bundle’s frequencies come from the inversion sources file); the remaining Re/Im pairs build the [N_freq, N_recv] complex array. Returns (Ex, None) to mirror readObservedDataH5() - a raw .dat carries no error-level attribute, so the noise level is supplied separately (-error_level / -im_error_level) or left to the kernel default.

Folding this in lets preprocess.py consume invEx.dat directly, so no separate convert-to-HDF5 step is required (the standalone utils/convert_invex_to_hdf5.py remains available for producing a reusable .h5).

Postprocessing

petgem.compareMagnitude(computed, reference)

Magnitude-domain error metrics between a computed and a reference field.

Shared postprocessing helper so every case reports the same three metrics the same way, instead of re-deriving them inline in each per-case script. Both inputs are 1-D field arrays (complex on a PETSc complex build, or real); the comparison is done on |field|.

Parameters

computed, referencearray_like

Field samples to compare (e.g. Ex at each receiver). Must have the same number of elements.

Returns

dict

nrmsd - RMS deviation of |computed| vs |reference|, normalized by the reference peak-to-peak range. rel_l2 - relative L2 norm of the magnitude difference. mape - mean absolute percentage error (%).