Field#

class emg3d.fields.Field(grid, data=None, frequency=None, dtype=None, electric=True)[source]#

Bases: object

A Field contains the x-, y-, and z- directed electromagnetic fields.

A Field is a simple container that has a 1D array Field.field containing the x-, y-, and z-directed fields one after the other. The field can be any field, such as an electric field, a magnetic field, or a source field (which is an electric field).

The particular fields can be accessed via the Field.f{x;y;z} attributes, which are 3D arrays corresponding to the shape of the edges (electric fields) or the faces (magnetic fields) in this direction; sort-order is Fortran-like (‘F’).

Parameters:
gridTensorMesh

The grid; a emg3d.meshes.TensorMesh instance.

datandarray, default: None

The actual data, a ndarray of size grid.n_edges (electric fields) or grid.n_faces (magnetic fields). If None, it is initiated with zeros.

frequency{float, None}, default: None

Field frequency (Hz), used to compute the Laplace parameter s. Either positive or negative:

  • frequency > 0: Frequency domain, hence \(s = \mathrm{i}\omega = 2\mathrm{i}\pi f\) (complex);

  • frequency < 0: Laplace domain, hence \(s = f\) (real).

This is primarily important for source fields. However, frequency information is also required to obtain the magnetic field from an electric field.

dtypedtype, default: complex

Data type of the initiated field; only used if both frequency and data are None.

electricbool, default: True

If electric, the properties live on the edges of the grid, if magnetic, they live on the faces.

Attributes Summary

field

Entire field as 1D array [fx, fy, fz].

frequency

Frequency (Hz).

fx

Field in x direction.

fy

Field in y direction.

fz

Field in z direction.

smu0

s*mu_0; mu_0 = magn permeability of free space [H/m].

sval

Laplace parameter; s=iw in f-domain and s=f in Laplace-domain.

Methods Summary

copy()

Return a copy of the Field.

from_dict(inp)

Convert dictionary into emg3d.fields.Field instance.

get_receiver(receiver[, method])

Return the field (response) at receiver coordinates.

interpolate_to_grid(grid, **interpolate_opts)

Interpolate the field to a new grid.

to_dict([copy])

Store the necessary information of the Field in a dict.

Attributes Documentation

field#

Entire field as 1D array [fx, fy, fz].

frequency#

Frequency (Hz).

fx#

Field in x direction.

Shape:

  • electric: (grid.cell_centers_x, grid.nodes_y, grid.nodes_z)

  • magnetic: (grid.nodes_x, grid.cell_centers_y, grid.cell_centers_z)

fy#

Field in y direction.

Shape:

  • electric: (grid.nodes_x, grid.cell_centers_y, grid.nodes_z)

  • magnetic: (grid.cell_centers_x, grid.nodes_y, grid.cell_centers_z)

fz#

Field in z direction.

Shape:

  • electric: (grid.nodes_x, grid.nodes_y, grid.cell_centers_z)

  • magnetic: (grid.cell_centers_x, grid.cell_centers_y, grid.nodes_z)

smu0#

s*mu_0; mu_0 = magn permeability of free space [H/m].

sval#

Laplace parameter; s=iw in f-domain and s=f in Laplace-domain.

Methods Documentation

copy()[source]#

Return a copy of the Field.

classmethod from_dict(inp)[source]#

Convert dictionary into emg3d.fields.Field instance.

Parameters:
inpdict

Dictionary as obtained from emg3d.fields.Field.to_dict. The dictionary needs the keys field, frequency, grid, and electric; grid itself is also a dict which needs the keys hx, hy, hz, and origin.

Returns:
fieldField

A emg3d.fields.Field instance.

get_receiver(receiver, method='cubic')[source]#

Return the field (response) at receiver coordinates.

Note that in order to avoid boundary effects from the PEC boundary the outermost cells are neglected. Field values for coordinates outside of the grid are set to NaN’s. However, take into account that for good results all receivers should be far away from the boundary.

Parameters:
receiver{Rx*, list, tuple}

Receiver coordinates. The following formats are accepted:

  • Rx* instance, any receiver object from emg3d.electrodes.

  • list: A list of Rx* instances.

  • tuple: (x, y, z, azimuth, elevation); receiver coordinates and angles (m, °). All values can either be a scalar or having the same length as number of receivers.

Note that the actual receiver type of the Rx* instances has no effect here, it just takes the coordinates from the receiver instances.

methodstr, default: ‘cubic’

Interpolation method to obtain the response at receiver location; ‘cubic’ or ‘linear’.

Returns:
responsesEMArray

Responses at receiver.

interpolate_to_grid(grid, **interpolate_opts)[source]#

Interpolate the field to a new grid.

If the provided grid is identical to the grid of the field, it returns the actual field (not a copy).

Parameters:
gridTensorMesh

New grid; a emg3d.meshes.TensorMesh instance.

interpolate_optsdict

Passed through to emg3d.maps.interpolate. Defaults are method='cubic', log=False, and extrapolate=False.

Returns:
fieldField

A new emg3d.fields.Field instance on grid.

to_dict(copy=False)[source]#

Store the necessary information of the Field in a dict.

Parameters:
copybool, default: False

If True, returns a deep copy of the dict.

Returns:
outdict

Dictionary containing all information to re-create the Field.