# Forward Simulation¶

## Simulation Class¶

The problem is a partial differential equation of the form:

$c(m, u) = 0$

Here, $$m$$ is the model and u is the field (or fields). Given the model, $$m$$, we can calculate the fields $$u(m)$$, however, the data we collect is a subset of the fields, and can be defined by a linear projection, $$P$$.

$d_\text{pred} = P u(m)$

For the inverse problem, we are interested in how changing the model transforms the data, as such we can take write the Taylor expansion:

$Pu(m + hv) = Pu(m) + hP\frac{\partial u(m)}{\partial m} v + \mathcal{O}(h^2 \left\| v \right\| )$

We can linearize and define the sensitivity matrix as:

$J = P\frac{\partial u}{\partial m}$

The sensitivity matrix, and it’s transpose will be used in the inverse problem to (locally) find how model parameters change the data, and optimize!

Working with the general PDE, $$c(m, u) = 0$$, where m is the model and u is the field, the sensitivity is defined as:

$J = P\frac{\partial u}{\partial m}$

We can take the derivative of the PDE:

$\nabla_m c(m, u) \partial m + \nabla_u c(m, u) \partial u = 0$

If the forward problem is invertible, then we can rearrange for $$\frac{\partial u}{\partial m}$$:

$J = - P \left( \nabla_u c(m, u) \right)^{-1} \nabla_m c(m, u)$

This can often be computed given a vector (i.e. $$J(v)$$) rather than stored, as $$J$$ is a large dense matrix.

## The API¶

### Simulation¶

class SimPEG.simulation.BaseSimulation(*args, **kwargs)[source]

BaseSimulation is the base class for all geophysical forward simulations in SimPEG.

Required Properties:

• counter (Counter): A SimPEG.utils.Counter object, an instance of Counter

• mesh (BaseMesh): a discretize mesh instance, an instance of BaseMesh

• sensitivity_path (String): path to store the sensitivty, a unicode string, Default: ./sensitivity/

• None

• solver_opts (Dictionary): solver options as a kwarg dict, a dictionary

• survey (BaseSurvey): a survey object, an instance of BaseSurvey

Optional Properties:

• model (Model): Inversion model., a numpy array of <class ‘float’>, <class ‘int’> with shape (*, *) or (*)

property mesh

mesh (BaseMesh): a discretize mesh instance, an instance of BaseMesh

property survey

survey (BaseSurvey): a survey object, an instance of BaseSurvey

property counter

counter (Counter): A SimPEG.utils.Counter object, an instance of Counter

property sensitivity_path

sensitivity_path (String): path to store the sensitivty, a unicode string, Default: ./sensitivity/

property solver_opts

solver_opts (Dictionary): solver options as a kwarg dict, a dictionary

deleteTheseOnModelUpdate = []
clean_on_model_update = []

List of matrix names to have their factors cleared on a model update

property Solver

Solver has been deprecated. See simulation.solver for documentation

property solverOpts

solverOpts has been deprecated. See solver_opts for documentation

property solver

None

fields(m=None)[source]

u = fields(m) The field given the model. :param numpy.ndarray m: model :rtype: numpy.ndarray :return: u, the fields

dpred(m, f=None)[source]

Create the projected data from a model. The fields, f, (if provided) will be used for the predicted data instead of recalculating the fields (which may be expensive!).

$d_\text{pred} = P(f(m))$

Where P is a projection of the fields onto the data space.

Jvec(m, v, f=None)[source]

Jv = Jvec(m, v, f=None) Effect of J(m) on a vector v. :param numpy.ndarray m: model :param numpy.ndarray v: vector to multiply :param Fields f: fields :rtype: numpy.ndarray :return: Jv

Jtvec(m, v, f=None)[source]

Jtv = Jtvec(m, v, f=None) Effect of transpose of J(m) on a vector v. :param numpy.ndarray m: model :param numpy.ndarray v: vector to multiply :param Fields f: fields :rtype: numpy.ndarray :return: JTv

Jvec_approx(m, v, f=None)[source]

Approximate effect of J(m) on a vector v :param numpy.ndarray m: model :param numpy.ndarray v: vector to multiply :param Fields f: fields :rtype: numpy.ndarray :return: approxJv

Jtvec_approx(m, v, f=None)[source]

Approximate effect of transpose of J(m) on a vector v. :param numpy.ndarray m: model :param numpy.ndarray v: vector to multiply :param Fields f: fields :rtype: numpy.ndarray :return: JTv

residual(m, dobs, f=None)[source]

The data residual:

$\mu_\text{data} = \mathbf{d}_\text{pred} - \mathbf{d}_\text{obs}$
Parameters
Return type

numpy.ndarray

Returns

data residual

make_synthetic_data(m, relative_error=0.05, noise_floor=0.0, f=None, add_noise=False, **kwargs)[source]

Make synthetic data given a model, and a standard deviation. :param numpy.ndarray m: geophysical model :param numpy.ndarray relative_error: standard deviation :param numpy.ndarray noise_floor: noise floor :param numpy.ndarray f: fields for the given model (if pre-calculated)

pair(survey)[source]

class SimPEG.simulation.BaseTimeSimulation(*args, **kwargs)[source]

Base class for a time domain simulation

Required Properties:

• counter (Counter): A SimPEG.utils.Counter object, an instance of Counter

• mesh (BaseMesh): a discretize mesh instance, an instance of BaseMesh

• sensitivity_path (String): path to store the sensitivty, a unicode string, Default: ./sensitivity/

• None

• solver_opts (Dictionary): solver options as a kwarg dict, a dictionary

• survey (BaseSurvey): a survey object, an instance of BaseSurvey

• t0 (Float): Origin of the time discretization, a float, Default: 0.0

• time_steps (TimeStepArray):

Sets/gets the time steps for the time domain simulation. You can set as an array of dt’s or as a list of tuples/floats. Tuples must be length two with […, (dt, repeat), …] For example, the following setters are the same:

sim.time_steps = [(1e-6, 3), 1e-5, (1e-4, 2)]
sim.time_steps = np.r_[1e-6,1e-6,1e-6,1e-5,1e-4,1e-4]

, an array or list of tuples specifying the mesh tensor of <class ‘float’> with shape (*)

Optional Properties:

• model (Model): Inversion model., a numpy array of <class ‘float’>, <class ‘int’> with shape (*, *) or (*)

property time_steps

time_steps (TimeStepArray): Sets/gets the time steps for the time domain simulation. You can set as an array of dt’s or as a list of tuples/floats. Tuples must be length two with […, (dt, repeat), …] For example, the following setters are the same:

sim.time_steps = [(1e-6, 3), 1e-5, (1e-4, 2)]
sim.time_steps = np.r_[1e-6,1e-6,1e-6,1e-5,1e-4,1e-4]

, an array or list of tuples specifying the mesh tensor of <class ‘float’> with shape (*)

property t0

t0 (Float): Origin of the time discretization, a float, Default: 0.0

property time_mesh
property nT
property times

Modeling times

property timeSteps

timeSteps has been deprecated. See time_steps for documentation

property timeMesh

time_mesh.timeMesh has been deprecated. See time_mesh for documentation

dpred(m, f=None)[source]

Create the projected data from a model. The fields, f, (if provided) will be used for the predicted data instead of recalculating the fields (which may be expensive!).

$d_\text{pred} = P(f(m))$

Where P is a projection of the fields onto the data space.

class SimPEG.simulation.LinearSimulation(*args, **kwargs)[source]

Class for a linear simulation of the form

$d = Gm$

where $$d$$ is a vector of the data, G is the simulation matrix and $$m$$ is the model. Inherit this class to build a linear simulation.

Required Properties:

• counter (Counter): A SimPEG.utils.Counter object, an instance of Counter

• mesh (BaseMesh): a discretize mesh instance, an instance of BaseMesh

• sensitivity_path (String): path to store the sensitivty, a unicode string, Default: ./sensitivity/

• None

• solver_opts (Dictionary): solver options as a kwarg dict, a dictionary

• survey (BaseSurvey): a survey object, an instance of BaseSurvey

Optional Properties:

• linear_model (PhysicalProperty): The model for a linear problem, a physical property

• model (Model): Inversion model., a numpy array of <class ‘float’>, <class ‘int’> with shape (*, *) or (*)

• model_map (Mapping): Mapping of The model for a linear problem to the inversion model., a SimPEG Map

Other Properties:

• model_deriv (Derivative): Derivative of The model for a linear problem wrt the model.

property linear_model

The model for a linear problem

property model_map

Mapping of The model for a linear problem to the inversion model.

property model_deriv

Derivative of The model for a linear problem wrt the model.

property mesh

mesh (BaseMesh): a discretize mesh instance, an instance of BaseMesh

property survey

survey (BaseSurvey): a survey object, an instance of BaseSurvey

property G
fields(m)[source]

u = fields(m) The field given the model. :param numpy.ndarray m: model :rtype: numpy.ndarray :return: u, the fields

dpred(m, f=None)[source]

Create the projected data from a model. The fields, f, (if provided) will be used for the predicted data instead of recalculating the fields (which may be expensive!).

$d_\text{pred} = P(f(m))$

Where P is a projection of the fields onto the data space.

getJ(m, f=None)[source]
Jvec(m, v, f=None)[source]

Jv = Jvec(m, v, f=None) Effect of J(m) on a vector v. :param numpy.ndarray m: model :param numpy.ndarray v: vector to multiply :param Fields f: fields :rtype: numpy.ndarray :return: Jv

Jtvec(m, v, f=None)[source]

Jtv = Jtvec(m, v, f=None) Effect of transpose of J(m) on a vector v. :param numpy.ndarray m: model :param numpy.ndarray v: vector to multiply :param Fields f: fields :rtype: numpy.ndarray :return: JTv

class SimPEG.simulation.TimeStepArray(*args, **kwargs)[source]
class_info = 'an array or list of tuples specifying the mesh tensor'
validate(instance, value)[source]

Determine if array is valid based on shape and dtype

### Fields¶

class SimPEG.fields.Fields(*args, **kwargs)[source]
Fancy Field Storage

Required Properties:

• aliasFields (Dictionary):

a dictionary of the aliased fields with [alias, location, function], e.g. {“b”:[“e”,”F”,lambda(F,e,ind)]} , a dictionary

• knownFields (Dictionary):

a dictionary with the names of the know fields and their location on a mesh e.g. {“e”: “E”, “phi”: “CC”} , a dictionary

• simulation (BaseSimulation): a SimPEG simulation, an instance of BaseSimulation

property knownFields

knownFields (Dictionary): a dictionary with the names of the know fields and their location on a mesh e.g. {“e”: “E”, “phi”: “CC”} , a dictionary

property aliasFields

aliasFields (Dictionary): a dictionary of the aliased fields with [alias, location, function], e.g. {“b”:[“e”,”F”,lambda(F,e,ind)]} , a dictionary

dtype

alias of builtins.float

property simulation

simulation (BaseSimulation): a SimPEG simulation, an instance of BaseSimulation

property mesh
property survey
startup()[source]
property approxSize

The approximate cost to storing all of the known fields.

class SimPEG.fields.TimeFields(*args, **kwargs)[source]
Fancy Field Storage for time domain problems
fields = TimeFields(simulation=simulation, knownFields={'phi':'CC'})
fields[:,'phi', timeInd] = phi
print(fields[src0,'phi'])

Required Properties:

• aliasFields (Dictionary):

a dictionary of the aliased fields with [alias, location, function], e.g. {“b”:[“e”,”F”,lambda(F,e,ind)]} , a dictionary

• knownFields (Dictionary):

a dictionary with the names of the know fields and their location on a mesh e.g. {“e”: “E”, “phi”: “CC”} , a dictionary

• simulation (BaseTimeSimulation): a SimPEG time simulation, an instance of BaseTimeSimulation

property simulation

simulation (BaseTimeSimulation): a SimPEG time simulation, an instance of BaseTimeSimulation

### Survey¶

class SimPEG.survey.BaseSurvey(*args, **kwargs)[source]

Survey holds the sources and receivers for a survey

Required Properties:

• counter (Counter): A SimPEG counter object, an instance of Counter

• source_list (a list of BaseSrc): A list of sources for the survey, a list (each item is an instance of BaseSrc)

property counter

counter (Counter): A SimPEG counter object, an instance of Counter

property source_list

source_list (a list of BaseSrc): A list of sources for the survey, a list (each item is an instance of BaseSrc)

getSourceIndex(sources)[source]
property nD

Number of data

property vnD

Vector number of data

property nSrc

Number of Sources

property srcList

srcList has been deprecated. See source_list for documentation

dpred(m=None, f=None)[source]
makeSyntheticData(m, std=None, f=None, force=False, **kwargs)[source]
pair(simulation)[source]
class SimPEG.survey.BaseTimeSurvey(*args, **kwargs)[source]

Required Properties:

• counter (Counter): A SimPEG counter object, an instance of Counter

• source_list (a list of BaseSrc): A list of sources for the survey, a list (each item is an instance of BaseSrc)

property unique_times
property times

unique_times.times has been deprecated. See unique_times for documentation

### Data¶

class SimPEG.data.Data(*args, **kwargs)[source]

Data storage. This class keeps track of observed data, relative error of those data and the noise floor.

data = Data(survey, dobs=dobs, relative_error=relative, noise_floor=floor)

or

data = Data(survey, dobs=dobs, standard_deviation=standard_deviation)

Required Properties:

• dobs (Array):

Vector of the observed data. The data can be set using the survey parameters:

data = Data(survey)
for src in survey.source_list:
data[src, rx] = datum

, a list or numpy array of <class ‘float’>, <class ‘int’> with shape (*)

• noise_floor (UncertaintyArray):

Noise floor of the data. This can be set using an array of the same size as the data (e.g. if you want to assign a different noise floor to each datum) or as a scalar if you would like to assign a the same noise floor to all data.

The standard_deviation is constructed as follows:

relative_error * np.abs(dobs) + noise_floor

For example, if you set

data = Data(survey, dobs=dobs)
data.noise_floor = 1e-10

then the contribution to the standard_deviation is equal to

data.noise_floor

, An array that can be set by a scalar value or numpy array of <class ‘float’>, <class ‘int’> with shape (*)

• relative_error (UncertaintyArray):

Relative error of the data. This can be set using an array of the same size as the data (e.g. if you want to assign a different relative error to each datum) or as a scalar if you would like to assign a the same relative error to all data.

The standard_deviation is constructed as follows:

relative_error * np.abs(dobs) + noise_floor

For example, if you set

data = Data(survey, dobs=dobs)
data.relative_error = 0.05

then the contribution to the standard_deviation is equal to

data.relative_error * np.abs(data.dobs)

, An array that can be set by a scalar value or numpy array of <class ‘float’>, <class ‘int’> with shape (*)

• survey (BaseSurvey): a SimPEG survey object, an instance of BaseSurvey

property survey

survey (BaseSurvey): a SimPEG survey object, an instance of BaseSurvey

property dobs

dobs (Array): Vector of the observed data. The data can be set using the survey parameters:

data = Data(survey)
for src in survey.source_list:
data[src, rx] = datum

, a list or numpy array of <class ‘float’>, <class ‘int’> with shape (*)

property relative_error

relative_error (UncertaintyArray): Relative error of the data. This can be set using an array of the same size as the data (e.g. if you want to assign a different relative error to each datum) or as a scalar if you would like to assign a the same relative error to all data.

The standard_deviation is constructed as follows:

relative_error * np.abs(dobs) + noise_floor

For example, if you set

data = Data(survey, dobs=dobs)
data.relative_error = 0.05

then the contribution to the standard_deviation is equal to

data.relative_error * np.abs(data.dobs)

, An array that can be set by a scalar value or numpy array of <class ‘float’>, <class ‘int’> with shape (*)

property noise_floor

noise_floor (UncertaintyArray): Noise floor of the data. This can be set using an array of the same size as the data (e.g. if you want to assign a different noise floor to each datum) or as a scalar if you would like to assign a the same noise floor to all data.

The standard_deviation is constructed as follows:

relative_error * np.abs(dobs) + noise_floor

For example, if you set

data = Data(survey, dobs=dobs)
data.noise_floor = 1e-10

then the contribution to the standard_deviation is equal to

data.noise_floor

, An array that can be set by a scalar value or numpy array of <class ‘float’>, <class ‘int’> with shape (*)

property standard_deviation

Data standard deviations. If a relative error and noise floor are provided, the standard_deviation is

data.standard_deviation = (
data.relative_error*np.abs(data.dobs) +
data.noise_floor
)

otherwise, the standard_deviation can be set directly

data.standard_deviation = 0.05 * np.absolute(self.dobs) + 1e-12

Note that setting the standard_deviation directly will clear the relative_error and set the value to the noise_floor property.

property nD
property shape
property index_dictionary

Dictionary of data indices by sources and receivers. To set data using survey parameters:

data = Data(survey)
for src in survey.source_list:
index = data.index_dictionary[src][rx]
data.dobs[index] = datum
tovec()[source]
fromvec(v)[source]
property std

std has been deprecated. See relative_error for documentation

property eps

eps has been deprecated. See noise_floor for documentation

class SimPEG.data.SyntheticData(*args, **kwargs)[source]

Data class for synthetic data. It keeps track of observed and clean data

Required Properties:

• dclean (Array):

Vector of the clean synthetic data. The data can be set using the survey parameters:

data = Data(survey)
for src in survey.source_list:
index = data.index_dictionary(src, rx)
data.dclean[indices] = datum

, a list or numpy array of <class ‘float’>, <class ‘int’> with shape (*)

• dobs (Array):

Vector of the observed data. The data can be set using the survey parameters:

data = Data(survey)
for src in survey.source_list:
data[src, rx] = datum

, a list or numpy array of <class ‘float’>, <class ‘int’> with shape (*)

• noise_floor (UncertaintyArray):

Noise floor of the data. This can be set using an array of the same size as the data (e.g. if you want to assign a different noise floor to each datum) or as a scalar if you would like to assign a the same noise floor to all data.

The standard_deviation is constructed as follows:

relative_error * np.abs(dobs) + noise_floor

For example, if you set

data = Data(survey, dobs=dobs)
data.noise_floor = 1e-10

then the contribution to the standard_deviation is equal to

data.noise_floor

, An array that can be set by a scalar value or numpy array of <class ‘float’>, <class ‘int’> with shape (*)

• relative_error (UncertaintyArray):

Relative error of the data. This can be set using an array of the same size as the data (e.g. if you want to assign a different relative error to each datum) or as a scalar if you would like to assign a the same relative error to all data.

The standard_deviation is constructed as follows:

relative_error * np.abs(dobs) + noise_floor

For example, if you set

data = Data(survey, dobs=dobs)
data.relative_error = 0.05

then the contribution to the standard_deviation is equal to

data.relative_error * np.abs(data.dobs)

, An array that can be set by a scalar value or numpy array of <class ‘float’>, <class ‘int’> with shape (*)

• survey (BaseSurvey): a SimPEG survey object, an instance of BaseSurvey

property dclean

dclean (Array): Vector of the clean synthetic data. The data can be set using the survey parameters:

data = Data(survey)
for src in survey.source_list:
index = data.index_dictionary(src, rx)
data.dclean[indices] = datum

, a list or numpy array of <class ‘float’>, <class ‘int’> with shape (*)

class SimPEG.data.UncertaintyArray(*args, **kwargs)[source]
class_info = 'An array that can be set by a scalar value or numpy array'
validate(instance, value)[source]

Determine if array is valid based on shape and dtype

### Sources¶

class SimPEG.survey.BaseSrc(*args, **kwargs)[source]

SimPEG Source Object

Required Properties:

• receiver_list (a list of BaseRx): receiver list, a list (each item is an instance of BaseRx)

Optional Properties:

• location (SourceLocationArray): Location of the source [x, y, z] in 3D, a 1D array denoting the source location of <class ‘float’>, <class ‘int’> with shape (*)

property loc

loc has been deprecated. See location for documentation

property rxList

rxList has been deprecated. See receiver_list for documentation

property nD

Number of data

property vnD

Vector number of data

receiver_list (a list of BaseRx): receiver list, a list (each item is an instance of BaseRx)

property location

location (SourceLocationArray): Location of the source [x, y, z] in 3D, a 1D array denoting the source location of <class ‘float’>, <class ‘int’> with shape (*)

class SimPEG.survey.SourceLocationArray(*args, **kwargs)[source]
class_info = 'a 1D array denoting the source location'
validate(instance, value)[source]

Determine if array is valid based on shape and dtype

class SimPEG.survey.BaseRx(*args, **kwargs)[source]

Required Properties:

• locations (RxLocationArray): Locations of the receivers (nRx x nDim), an array of receiver locations of <class ‘float’>, <class ‘int’> with shape (*, *)

• projGLoc (StringChoice): Projection grid location, default is CC, any of “CC”, “Fx”, “Fy”, “Fz”, “Ex”, “Ey”, “Ez”, “N”, Default: CC

• storeProjections (Boolean): Store calls to getP (organized by mesh), a boolean, Default: True

property projGLoc

projGLoc (StringChoice): Projection grid location, default is CC, any of “CC”, “Fx”, “Fy”, “Fz”, “Ex”, “Ey”, “Ez”, “N”, Default: CC

property storeProjections

storeProjections (Boolean): Store calls to getP (organized by mesh), a boolean, Default: True

property locations

locations (RxLocationArray): Locations of the receivers (nRx x nDim), an array of receiver locations of <class ‘float’>, <class ‘int’> with shape (*, *)

property locs

locs has been deprecated. See locations for documentation

property nD

Number of data in the receiver.

getP(mesh, projGLoc=None)[source]

Returns the projection matrices as a list for all components collected by the receivers.

Note

Projection matrices are stored as a dictionary listed by meshes.

eval(**kwargs)[source]
evalDeriv(**kwargs)[source]
class SimPEG.survey.BaseTimeRx(*args, **kwargs)[source]

SimPEG Receiver Object for time-domain simulations

Required Properties:

• locations (RxLocationArray): Locations of the receivers (nRx x nDim), an array of receiver locations of <class ‘float’>, <class ‘int’> with shape (*, *)

• projGLoc (StringChoice): Projection grid location, default is CC, any of “CC”, “Fx”, “Fy”, “Fz”, “Ex”, “Ey”, “Ez”, “N”, Default: CC

• projTLoc (StringChoice): location on the time mesh where the data are projected from, either “N” or “CC”, Default: N

• storeProjections (Boolean): Store calls to getP (organized by mesh), a boolean, Default: True

• times (Array): times where the recievers measure data, a list or numpy array of <class ‘float’>, <class ‘int’> with shape (*)

property projTLoc

projTLoc (StringChoice): location on the time mesh where the data are projected from, either “N” or “CC”, Default: N

property times

times (Array): times where the recievers measure data, a list or numpy array of <class ‘float’>, <class ‘int’> with shape (*)

property nD

Number of data in the receiver.

getSpatialP(mesh)[source]

Returns the spatial projection matrix.

Note

This is not stored in memory, but is created on demand.

getTimeP(timeMesh)[source]

Returns the time projection matrix.

Note

This is not stored in memory, but is created on demand.

getP(mesh, timeMesh)[source]

Returns the projection matrices as a list for all components collected by the receivers.

Note

Projection matrices are stored as a dictionary (mesh, timeMesh) if storeProjections is True

class SimPEG.survey.RxLocationArray(*args, **kwargs)[source]
class_info = 'an array of receiver locations'
validate(instance, value)[source]

Determine if array is valid based on shape and dtype