7. Command-line interface

The implemented command-line interface (CLI) can be powerful to use emg3d from any programming environment. It enables, for instance, to use emg3d as a forward modelling kernel in your inversion code running on a server. The capabilities of the CLI are more limited than of the Python interface, but they are sufficient to compute the responses of a given survey for a provided model, and to compute the gradient of the misfit function.

The biggest difficulty to work with the CLI is, at the moment, the file formats, or better the lack of their documentation. Surveys, models, and grids can be stored with emg3d to either h5, npz, or json files. The status of the documentation will improve towards the v1.0.0 release of emg3d. The “easiest” way at the moment is to generate such a file in Python, and reproduce its structure in your language or program of choice. Please feel free to open an issue on GitHub if there are questions in this regard. The idea for the future is to go through https://github.com/softwareunderground/subsurface to interface to various file formats.

import os
import pooch
import subprocess

# Adjust this path to a folder of your choice.
data_path = os.path.join('..', 'download', '')

Note that everything shown in this example is meant to be executed in a terminal, nothing is executed in Python. However, as this example gallery is generated in Python we have to use a work-around to show the functionality, the function bash(command). You would execute the provided command in your terminal.

def bash(command):
    "Prints the `command`, executes it in bash, and prints the output."""
    # Print command
    print(f"$ {command}")
    # Move to data_path
    command = f"cd {data_path}; " + command
    # Carry out command
    msg = subprocess.run(command, shell=True, stdout=subprocess.PIPE,
                         stderr=subprocess.STDOUT)
    # Print output
    print(msg.stdout.decode())

We have to fetch the example files. You can also download these manually or provide your own survey, model, and parameter files.

model = "GemPy-II.h5"
pooch.retrieve(
    'https://raw.github.com/emsig/data/2021-05-21/emg3d/models/'+model,
    'ea8c23be80522d3ca8f36742c93758370df89188816f50cb4e1b2a6a3012d659',
    fname=model,
    path=data_path,
)
survey = 'GemPy-II-survey-A.h5'
pooch.retrieve(
    'https://raw.github.com/emsig/data/2021-05-21/emg3d/surveys/'+survey,
    '5f2ed0b959a4f80f5378a071e6f729c6b7446898be7689ddc9bbd100a8f5bce7',
    fname=survey,
    path=data_path,
)

Out:

'/home/dtr/Codes/emsig/emg3d-gallery/examples/download/GemPy-II-survey-A.h5'

Getting started - help

The best way to get started is, as with any bash command, to have a look at its in-built help.

bash("emg3d --help")

Out:

$ emg3d --help
usage: emg3d [-h] [-n NPROC] [-f | -m | -g] [--path PATH] [--survey SURVEY]
             [--model MODEL] [--output OUTPUT]
             [--verbosity {-1,0,1,2} | -v | -q] [-d] [--report] [--version]
             [config]

Multigrid solver for 3D electromagnetic diffusion.

positional arguments:
  config                name of config file; default is 'emg3d.cfg'; consult
                        https://emg3d.emsig.xyz/en/stable/cli.html for its
                        format

optional arguments:
  -h, --help            show this help message and exit
  -n NPROC, --nproc NPROC
                        number of processors
  -f, --forward         compute synthetic data (default)
  -m, --misfit          compute synthetic data and their misfit
  -g, --gradient        compute synthetic data, misfit, and its gradient
  --path PATH           path (abs or rel); file names are relative to path
  --survey SURVEY       input survey file name; default is 'survey.h5'
  --model MODEL         input model file name; default is 'model.h5'
  --output OUTPUT       output files base name; default is 'emg3d_out'
  --verbosity {-1,0,1,2}
                        set verbosity; default is 0
  -v, --verbose         increase verbosity; can be used multiple times
  -q, --quiet           decrease verbosity
  -d, --dry-run         only display what would have been done
  --report              only display emg3d report
  --version             only display emg3d version

Configuration file

The CLI is driven by a configuration file, that is by default named emg3d.cfg and must be in the same directory as you execute the command. If it has a different name or a different location you have to provide the full or relative path and filename as the first argument to emg3d.

The configuration parameters are described in the documentation, consult Manual -> CLI.

Let’s write a simple example file.

with open(f"{data_path}emg3d.cfg", "w") as f:
    f.write("""[files]
survey = GemPy-II-survey-A.h5
model = GemPy-II.h5
output = emg3d_out.json

[simulation]
name = CLI Test
gridding = single

[solver_opts]
# nothing specified

[gridding_opts]
# nothing specified

[data]
sources = TxED-1
receivers = RxEP-08, RxEP-38
frequencies = f-1""")

The file defines the location of the survey and model files, the name and format of the output file, and selects a source, two receivers, and a frequency to compute.

Compute the forward responses

We just compute the forward response here, but the usage for the misfit or the gradient is very similar. We use the most verbose version here, to see what it does internally.

bash("emg3d --forward -vv")

Out:

$ emg3d --forward -vv
:: emg3d CLI forward START :: Sat Jul 17 18:00:59 2021 :: v1.1.0

--------------------------------------------------------------------------------
  Date: Sat Jul 17 18:00:59 2021 CEST

                OS : Linux
            CPU(s) : 4
           Machine : x86_64
      Architecture : 64bit
               RAM : 15.5 GiB
       Environment : Python

  Python 3.9.4 | packaged by conda-forge | (default, May 10 2021, 22:13:33)
  [GCC 9.3.0]

             numpy : 1.21.0
             scipy : 1.7.0
             numba : 0.53.1
             emg3d : 1.1.0
           empymod : 2.1.2
            xarray : 0.18.2
        discretize : 0.7.0
              h5py : 3.3.0
        matplotlib : 3.4.2
              tqdm : 4.61.2

  Intel(R) oneAPI Math Kernel Library Version 2021.2-Product Build 20210312
  for Intel(R) 64 architecture applications
--------------------------------------------------------------------------------

    :: CONFIGURATION ::

/home/dtr/Codes/emsig/emg3d-gallery/examples/download/emg3d.cfg
{
    "data": {
        "frequencies": [
            "f-1"
        ],
        "receivers": [
            "RxEP-08",
            "RxEP-38"
        ],
        "sources": [
            "TxED-1"
        ]
    },
    "files": {
        "log": "/home/dtr/Codes/emsig/emg3d-gallery/examples/download/emg3d_out.log",
        "model": "/home/dtr/Codes/emsig/emg3d-gallery/examples/download/GemPy-II.h5",
        "output": "/home/dtr/Codes/emsig/emg3d-gallery/examples/download/emg3d_out.json",
        "store_simulation": false,
        "survey": "/home/dtr/Codes/emsig/emg3d-gallery/examples/download/GemPy-II-survey-A.h5"
    },
    "simulation_options": {
        "gridding": "single",
        "name": "CLI Test"
    }
}

    :: LOAD SURVEY AND MODEL ::

Data loaded from «/home/dtr/Codes/emsig/emg3d-gallery/examples/download/GemPy-II-survey-A.h5»
[emg3d v0.17.1.dev22+ge23c468 (format 1.0) on 2021-04-14T22:24:03.229285].
Data loaded from «/home/dtr/Codes/emsig/emg3d-gallery/examples/download/GemPy-II.h5»
[emg3d v1.0.0rc3.dev5+g0cd9e09 (format 1.0) on 2021-05-21T18:40:16.721968].
/home/dtr/anaconda3/envs/emg3d-gallery/lib/python3.9/site-packages/emg3d/simulations.py:1501: RuntimeWarning: divide by zero encountered in double_scalars
  if get_y and xdiff/ydiff > 3:


    :: SIMULATION ::

:: Simulation «CLI Test» ::

- Survey «GemPy-II Survey A»: 1 sources; 2 receivers; 1 frequencies
- Model: resistivity; isotropic; 100 x 100 x 102 (1,020,000)
- Gridding: Single grid for all sources and frequencies; 128 x 64 x 64 (524,288)

    :: MESHES ::

= Source: all; Frequency: all =

         == GRIDDING IN X ==
Skin depth     [m] : 390 / 7117625 / 7117625  [corr. to `properties`]
Survey dom. DS [m] : 4000 - 16000
Comp. dom. DC  [m] : -96000 - 116000
Final extent   [m] : -101116 - 160572
Cell widths    [m] : 130 / 130 / 39449  [min(DS) / max(DS) / max(DC)]
Number of cells    : 128 (93 / 34 / 1)  [Total (DS/DC/remain)]
Max stretching     : 1.000 (1.000) / 1.374  [DS (seasurface) / DC]

         == GRIDDING IN Y ==
Skin depth     [m] : 390 / 7117625 / 7117625  [corr. to `properties`]
Survey dom. DS [m] : 5500 - 9500
Comp. dom. DC  [m] : -94500 - 109500
Final extent   [m] : -97120 - 112120
Cell widths    [m] : 130 / 130 / 29638  [min(DS) / max(DS) / max(DC)]
Number of cells    : 64 (32 / 32 / 0)  [Total (DS/DC/remain)]
Max stretching     : 1.000 (1.000) / 1.404  [DS (seasurface) / DC]

         == GRIDDING IN Z ==
Skin depth     [m] : 390 / 2251 / 7117625  [corr. to `properties`]
Survey dom. DS [m] : -6725 - -1722
Comp. dom. DC  [m] : -20867 - 98278
Final extent   [m] : -27189 - 98732
Cell widths    [m] : 130 / 130 / 32910  [min(DS) / max(DS) / max(DC)]
Number of cells    : 64 (40 / 24 / 0)  [Total (DS/DC/remain)]
Max stretching     : 1.000 (1.000) / 1.485  [DS (seasurface) / DC]

  TensorMesh: 524,288 cells

                      MESH EXTENT             CELL WIDTH      FACTOR
  dir    nC        min           max         min       max      max
  ---   ---  ---------------------------  ------------------  ------
   x    128   -101,115.97    160,571.57    129.95 39,449.49    1.37
   y     64    -97,119.50    112,119.50    129.95 29,637.95    1.40
   z     64    -27,189.50     98,732.09    129.95 32,909.74    1.48


    :: FORWARD COMPUTATION ::


Compute efields:           0/1  [00:00]
Compute efields: ██████████1/1  [01:04]
Compute efields: ██████████1/1  [01:04]

    - SOLVER INFO <efield> -

= Source TxED-1; Frequency 0.5 Hz = 2.3e-03; 0(1); 0:00:13
:: emg3d :: 1.3e-04; 0(2); 0:00:26
:: emg3d :: 1.6e-05; 0(3); 0:00:38
:: emg3d :: 2.9e-06; 0(4); 0:00:51
:: emg3d :: 5.7e-07; 0(5); 0:01:04
:: emg3d :: 5.6e-07; 1(5); 0:01:04
:: emg3d :: 5.6e-07; 1(5); 0:01:04; CONVERGED

    :: SAVE RESULTS ::

Data saved to «/home/dtr/Codes/emsig/emg3d-gallery/examples/download/emg3d_out.json»
[emg3d v1.1.0 (format 1.0) on 2021-07-17T18:02:04.813529].

:: emg3d CLI forward END   :: Sat Jul 17 18:02:04 2021 :: runtime = 0:01:05

Generated files

The run created the files emg3d_out.json (name and file-type can be changed in the config file), which contains the responses, and the file emg3d_out.log with some information.

bash("ls emg3d* -1")

Out:

$ ls emg3d* -1
emg3d.cfg
emg3d_out.json
emg3d_out.log

Log

Let’s have a look at the log, which is mostly the same as was printed above with the verbosity flag:

bash("cat emg3d_out.log")

Out:

$ cat emg3d_out.log
:: emg3d CLI forward START :: Sat Jul 17 18:00:59 2021 :: v1.1.0

--------------------------------------------------------------------------------
  Date: Sat Jul 17 18:00:59 2021 CEST

                OS : Linux
            CPU(s) : 4
           Machine : x86_64
      Architecture : 64bit
               RAM : 15.5 GiB
       Environment : Python

  Python 3.9.4 | packaged by conda-forge | (default, May 10 2021, 22:13:33)
  [GCC 9.3.0]

             numpy : 1.21.0
             scipy : 1.7.0
             numba : 0.53.1
             emg3d : 1.1.0
           empymod : 2.1.2
            xarray : 0.18.2
        discretize : 0.7.0
              h5py : 3.3.0
        matplotlib : 3.4.2
              tqdm : 4.61.2

  Intel(R) oneAPI Math Kernel Library Version 2021.2-Product Build 20210312
  for Intel(R) 64 architecture applications
--------------------------------------------------------------------------------

    :: CONFIGURATION ::

/home/dtr/Codes/emsig/emg3d-gallery/examples/download/emg3d.cfg
{
    "data": {
        "frequencies": [
            "f-1"
        ],
        "receivers": [
            "RxEP-08",
            "RxEP-38"
        ],
        "sources": [
            "TxED-1"
        ]
    },
    "files": {
        "log": "/home/dtr/Codes/emsig/emg3d-gallery/examples/download/emg3d_out.log",
        "model": "/home/dtr/Codes/emsig/emg3d-gallery/examples/download/GemPy-II.h5",
        "output": "/home/dtr/Codes/emsig/emg3d-gallery/examples/download/emg3d_out.json",
        "store_simulation": false,
        "survey": "/home/dtr/Codes/emsig/emg3d-gallery/examples/download/GemPy-II-survey-A.h5"
    },
    "simulation_options": {
        "gridding": "single",
        "name": "CLI Test"
    }
}

    :: LOAD SURVEY AND MODEL ::

Data loaded from «/home/dtr/Codes/emsig/emg3d-gallery/examples/download/GemPy-II-survey-A.h5»
[emg3d v0.17.1.dev22+ge23c468 (format 1.0) on 2021-04-14T22:24:03.229285].
Data loaded from «/home/dtr/Codes/emsig/emg3d-gallery/examples/download/GemPy-II.h5»
[emg3d v1.0.0rc3.dev5+g0cd9e09 (format 1.0) on 2021-05-21T18:40:16.721968].
/home/dtr/anaconda3/envs/emg3d-gallery/lib/python3.9/site-packages/emg3d/simulations.py:1501: RuntimeWarning: divide by zero encountered in double_scalars
  if get_y and xdiff/ydiff > 3:


    :: SIMULATION ::

:: Simulation «CLI Test» ::

- Survey «GemPy-II Survey A»: 1 sources; 2 receivers; 1 frequencies
- Model: resistivity; isotropic; 100 x 100 x 102 (1,020,000)
- Gridding: Single grid for all sources and frequencies; 128 x 64 x 64 (524,288)

    :: MESHES ::

= Source: all; Frequency: all =

         == GRIDDING IN X ==
Skin depth     [m] : 390 / 7117625 / 7117625  [corr. to `properties`]
Survey dom. DS [m] : 4000 - 16000
Comp. dom. DC  [m] : -96000 - 116000
Final extent   [m] : -101116 - 160572
Cell widths    [m] : 130 / 130 / 39449  [min(DS) / max(DS) / max(DC)]
Number of cells    : 128 (93 / 34 / 1)  [Total (DS/DC/remain)]
Max stretching     : 1.000 (1.000) / 1.374  [DS (seasurface) / DC]

         == GRIDDING IN Y ==
Skin depth     [m] : 390 / 7117625 / 7117625  [corr. to `properties`]
Survey dom. DS [m] : 5500 - 9500
Comp. dom. DC  [m] : -94500 - 109500
Final extent   [m] : -97120 - 112120
Cell widths    [m] : 130 / 130 / 29638  [min(DS) / max(DS) / max(DC)]
Number of cells    : 64 (32 / 32 / 0)  [Total (DS/DC/remain)]
Max stretching     : 1.000 (1.000) / 1.404  [DS (seasurface) / DC]

         == GRIDDING IN Z ==
Skin depth     [m] : 390 / 2251 / 7117625  [corr. to `properties`]
Survey dom. DS [m] : -6725 - -1722
Comp. dom. DC  [m] : -20867 - 98278
Final extent   [m] : -27189 - 98732
Cell widths    [m] : 130 / 130 / 32910  [min(DS) / max(DS) / max(DC)]
Number of cells    : 64 (40 / 24 / 0)  [Total (DS/DC/remain)]
Max stretching     : 1.000 (1.000) / 1.485  [DS (seasurface) / DC]

  TensorMesh: 524,288 cells

                      MESH EXTENT             CELL WIDTH      FACTOR
  dir    nC        min           max         min       max      max
  ---   ---  ---------------------------  ------------------  ------
   x    128   -101,115.97    160,571.57    129.95 39,449.49    1.37
   y     64    -97,119.50    112,119.50    129.95 29,637.95    1.40
   z     64    -27,189.50     98,732.09    129.95 32,909.74    1.48


    :: FORWARD COMPUTATION ::


    - SOLVER INFO <efield> -

= Source TxED-1; Frequency 0.5 Hz = 2.3e-03; 0(1); 0:00:13
:: emg3d :: 1.3e-04; 0(2); 0:00:26
:: emg3d :: 1.6e-05; 0(3); 0:00:38
:: emg3d :: 2.9e-06; 0(4); 0:00:51
:: emg3d :: 5.7e-07; 0(5); 0:01:04
:: emg3d :: 5.6e-07; 1(5); 0:01:04
:: emg3d :: 5.6e-07; 1(5); 0:01:04; CONVERGED

    :: SAVE RESULTS ::

Data saved to «/home/dtr/Codes/emsig/emg3d-gallery/examples/download/emg3d_out.json»
[emg3d v1.1.0 (format 1.0) on 2021-07-17T18:02:04.813529].

:: emg3d CLI forward END   :: Sat Jul 17 18:02:04 2021 :: runtime = 0:01:05

Result / output

We stored the output as a json-file, just so we can easily cat it. However, you might want to use h5 files, as they are compressed. Note that json cannot store complex numbers, so the responses are stored as [[Re], [Im]].

bash("cat emg3d_out.json")

Out:

$ cat emg3d_out.json
{
  "configuration": {
    "files": {
      "survey": "/home/dtr/Codes/emsig/emg3d-gallery/examples/download/GemPy-II-survey-A.h5",
      "model": "/home/dtr/Codes/emsig/emg3d-gallery/examples/download/GemPy-II.h5",
      "output": "/home/dtr/Codes/emsig/emg3d-gallery/examples/download/emg3d_out.json",
      "log": "/home/dtr/Codes/emsig/emg3d-gallery/examples/download/emg3d_out.log",
      "store_simulation": false
    },
    "simulation_options": {
      "gridding": "single",
      "name": "CLI Test"
    },
    "data": {
      "sources": [
        "TxED-1"
      ],
      "receivers": [
        "RxEP-08",
        "RxEP-38"
      ],
      "frequencies": [
        "f-1"
      ]
    }
  },
  "data__complex__array-float64": [
    [
      [
        [
          -2.3279584817131315e-08
        ],
        [
          NaN
        ]
      ]
    ],
    [
      [
        [
          -4.8500089399022605e-09
        ],
        [
          NaN
        ]
      ]
    ]
  ],
  "_date": "2021-07-17T18:02:04.813529",
  "_version": "emg3d v1.1.0",
  "_format": "1.0"
}

Report

The report-command prints an extensive list of your python environment of components which are important for emg3d. If you report an issue it helps a lot if you provide this info.

bash("emg3d --report")

Out:

$ emg3d --report

--------------------------------------------------------------------------------
  Date: Sat Jul 17 18:02:09 2021 CEST

                OS : Linux
            CPU(s) : 4
           Machine : x86_64
      Architecture : 64bit
               RAM : 15.5 GiB
       Environment : Python

  Python 3.9.4 | packaged by conda-forge | (default, May 10 2021, 22:13:33)
  [GCC 9.3.0]

             numpy : 1.21.0
             scipy : 1.7.0
             numba : 0.53.1
             emg3d : 1.1.0
           empymod : 2.1.2
            xarray : 0.18.2
        discretize : 0.7.0
              h5py : 3.3.0
        matplotlib : 3.4.2
              tqdm : 4.61.2

  Intel(R) oneAPI Math Kernel Library Version 2021.2-Product Build 20210312
  for Intel(R) 64 architecture applications
--------------------------------------------------------------------------------

Total running time of the script: ( 1 minutes 16.978 seconds)

Estimated memory usage: 9 MB

Gallery generated by Sphinx-Gallery