InvProblem

class SimPEG.InvProblem.BaseInvProblem(dmisfit, reg, opt)[source]

Bases: SimPEG.Props.BaseSimPEG

Optional Properties:

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

Trade-off parameter

debug = False

Print debugging information

counter = None

Set this to a SimPEG.Utils.Counter() if you want to count things

deleteTheseOnModelUpdate = []

List of strings, e.g. [‘_MeSigma’, ‘_MeSigmaI’]

model

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

dmisfit = None

DataMisfit

reg = None

Regularization

opt = None

Optimization program

startup(m0)[source]

Called when inversion is first starting.

If you have things that also need to run in the method startup, you can create a method:

def _startup*(self, ... ):
    pass

Where the * can be any string. If present, _startup* will be called at the start of the default startup call. You may also completely overwrite this function.

warmstart
getFields(m, store=False, deleteWarmstart=True)[source]
get_dpred(m, f)[source]
evalFunction(m, return_g=True, return_H=True)[source]

Inversion

class SimPEG.Inversion.BaseInversion(invProb, directiveList=None, **kwargs)[source]

Bases: object

Inversion Class

name = 'BaseInversion'
debug = False

Print debugging information

counter = None

Set this to a SimPEG.Utils.Counter() if you want to count things

directiveList
run(m0)[source]

Runs the inversion!

Directives

class SimPEG.Directives.InversionDirective(**kwargs)[source]

Bases: properties.base.base.HasProperties

debug = False

Print debugging information

inversion

This is the inversion of the InversionDirective instance.

invProb
opt
reg
dmisfit
survey

Assuming that dmisfit is always a ComboObjectiveFunction, return a list of surveys for each dmisfit [survey1, survey2, ... ]

prob

Assuming that dmisfit is always a ComboObjectiveFunction, return a list of problems for each dmisfit [prob1, prob2, ...]

initialize()[source]
endIter()[source]
finish()[source]
validate(directiveList=None)[source]
class SimPEG.Directives.DirectiveList(*directives, **kwargs)[source]

Bases: object

dList = None

The list of Directives

debug
inversion

This is the inversion of the InversionDirective instance.

call(ruleType)[source]
validate()[source]
class SimPEG.Directives.BetaEstimate_ByEig(**kwargs)[source]

Bases: SimPEG.Directives.InversionDirective

BetaEstimate

beta0 = None

The initial Beta (regularization parameter)

beta0_ratio = 100.0

estimateBeta0 is used with this ratio

initialize()[source]

The initial beta is calculated by comparing the estimated eigenvalues of JtJ and WtW. To estimate the eigenvector of A, we will use one iteration of the Power Method:

\[\mathbf{x_1 = A x_0}\]

Given this (very course) approximation of the eigenvector, we can use the Rayleigh quotient to approximate the largest eigenvalue.

\[\lambda_0 = \frac{\mathbf{x^\top A x}}{\mathbf{x^\top x}}\]

We will approximate the largest eigenvalue for both JtJ and WtW, and use some ratio of the quotient to estimate beta0.

\[\beta_0 = \gamma \frac{\mathbf{x^\top J^\top J x}}{\mathbf{x^\top W^\top W x}}\]
Return type:float
Returns:beta0
class SimPEG.Directives.BetaSchedule(**kwargs)[source]

Bases: SimPEG.Directives.InversionDirective

coolingFactor = 8.0
coolingRate = 3
endIter()[source]
class SimPEG.Directives.TargetMisfit(**kwargs)[source]

Bases: SimPEG.Directives.InversionDirective

... note:: Currently the target misfit is not set up for joint inversions. Get in touch if you would like to help with the upgrade!

chifact = 1.0
phi_d_star = None
target
endIter()[source]
class SimPEG.Directives.SaveEveryIteration(**kwargs)[source]

Bases: SimPEG.Directives.InversionDirective

This directive saves an array at each iteration. The default direcroty is the current directoy and the models are saved as InversionModel-YYYY-MM-DD-HH-MM-iter.npy

Required Properties:

  • directory (String): directory to save results in, a unicode string, Default: .
  • name (String): root of the filename to be saved, a unicode string, Default: InversionModel
directory

directory (String): directory to save results in, a unicode string, Default: .

name

name (String): root of the filename to be saved, a unicode string, Default: InversionModel

fileName
class SimPEG.Directives.SaveModelEveryIteration(**kwargs)[source]

Bases: SimPEG.Directives.SaveEveryIteration

This directive saves the model as a numpy array at each iteration. The default direcroty is the current directoy and the models are saved as InversionModel-YYYY-MM-DD-HH-MM-iter.npy

Required Properties:

  • directory (String): directory to save results in, a unicode string, Default: .
  • name (String): root of the filename to be saved, a unicode string, Default: InversionModel
initialize()[source]
endIter()[source]
class SimPEG.Directives.SaveOutputEveryIteration(**kwargs)[source]

Bases: SimPEG.Directives.SaveEveryIteration

Required Properties:

  • directory (String): directory to save results in, a unicode string, Default: .
  • name (String): root of the filename to be saved, a unicode string, Default: InversionModel
header = None
save_txt = True
beta = None
phi_d = None
phi_m = None
phi_m_small = None
phi_m_smooth_x = None
phi_m_smooth_y = None
phi_m_smooth_z = None
phi = None
initialize()[source]
endIter()[source]
load_results()[source]
plot_misfit_curves(fname=None, dpi=300, plot_small_smooth=False, plot_phi_m=True, plot_small=False, plot_smooth=False)[source]
plot_tikhonov_curves(fname=None, dpi=200)[source]
class SimPEG.Directives.SaveOutputDictEveryIteration(**kwargs)[source]

Bases: SimPEG.Directives.SaveEveryIteration

Saves inversion parameters at every iteraion.

Required Properties:

  • directory (String): directory to save results in, a unicode string, Default: .
  • name (String): root of the filename to be saved, a unicode string, Default: InversionModel
outDict = {}
saveOnDisk = False
initialize()[source]
endIter()[source]
class SimPEG.Directives.Update_IRLS(**kwargs)[source]

Bases: SimPEG.Directives.InversionDirective

f_old = 0
f_min_change = 0.01
beta_tol = 0.1
beta_ratio_l2 = None
prctile = 100
chifact_start = 1.0
chifact_target = 1.0
IRLSiter = 0
minGNiter = 1
maxIRLSiter = 20
iterStart = 0
sphericalDomain = False
updateBeta = True
betaSearch = True
coolingFactor = 2.0
coolingRate = 1
ComboObjFun = False
mode = 1
coolEpsOptimized = True
coolEps_p = True
coolEps_q = True
floorEps_p = 1e-08
floorEps_q = 1e-08
coolEpsFact = 1.2
silent = False
fix_Jmatrix = False
target
start
initialize()[source]
endIter()[source]
startIRLS()[source]
angleScale()[source]

Update the scales used by regularization for the different block of models

validate(directiveList)[source]
class SimPEG.Directives.UpdatePreconditioner(**kwargs)[source]

Bases: SimPEG.Directives.InversionDirective

Create a Jacobi preconditioner for the linear problem

update_every_iteration = True

Update every iterations if False

threshold = 1e-08
initialize()[source]
endIter()[source]
class SimPEG.Directives.Update_Wj(**kwargs)[source]

Bases: SimPEG.Directives.InversionDirective

Create approx-sensitivity base weighting using the probing method

k = None
itr = None
endIter()[source]
class SimPEG.Directives.UpdateSensitivityWeights(**kwargs)[source]

Bases: SimPEG.Directives.InversionDirective

Directive to take care of re-weighting the non-linear magnetic problems.

mapping = None
JtJdiag = None
everyIter = True
threshold = 1e-12
switch = True
initialize()[source]
endIter()[source]
update()[source]
getJtJdiag()[source]

Compute explicitely the main diagonal of JtJ Good for any problem where J is formed explicitely

getWr()[source]

Take the diagonal of JtJ and return a normalized sensitivty weighting vector

updateReg()[source]

Update the cell weights with the approximated sensitivity

updateOpt()[source]

Update a copy of JtJdiag to optimization for preconditioner

class SimPEG.Directives.ProjectSphericalBounds(**kwargs)[source]

Bases: SimPEG.Directives.InversionDirective

Trick for spherical coordinate system. Project heta and phi angles back to [-pi,pi] using back and forth conversion. spherical->cartesian->spherical

initialize()[source]
endIter()[source]