simpeg.regularization.SparseSmoothness#

class simpeg.regularization.SparseSmoothness(mesh, orientation='x', gradient_type='total', **kwargs)[source]#

Bases: BaseSparse, SmoothnessFirstOrder

Sparse smoothness (blockiness) regularization.

SparseSmoothness is used to recover models comprised of blocky structures. The level of blockiness is controlled by the choice in norm within the regularization function; with more blocky structures being recovered when a smaller norm is used. Optionally, custom cell weights can be included to control the degree of blockiness being enforced throughout different regions the model.

See the Notes section below for a comprehensive description.

Parameters:

meshregularization.RegularizationMesh

Mesh on which the regularization is discretized. Not the mesh used to define the simulation.

orientation{‘x’,’y’,’z’}

The direction along which sparse smoothness is applied.

normfloat, array_like

The norm defining sparseness thoughout the regularization function. Must be within the interval [0,2]. There are several options:

float: constant sparse norm throughout the domain.
(n_faces, ) array_like: define the sparse norm independently at each face set by orientation (e.g. x-faces).
(n_cells, ) array_like: define the sparse norm independently for each cell. Will be averaged to faces specified by orientation (e.g. x-faces).

active_cellsNone, (n_cells, ) numpy.ndarray of bool

Boolean array defining the set of RegularizationMesh cells that are active in the inversion. If None, all cells are active.

mappingNone, simpeg.maps.BaseMap

The mapping from the model parameters to the active cells in the inversion. If None, the mapping is the identity map.

reference_modelNone, (n_param, ) numpy.ndarray

Reference model. If None, the reference model in the inversion is set to the starting model. To include the reference model in the regularization, the reference_model_in_smooth property must be set to True.

reference_model_in_smoothbool, optional

Whether to include the reference model in the smoothness terms.

unitsNone, str

Units for the model parameters. Some regularization classes behave differently depending on the units; e.g. ‘radian’.

weightsNone, dict

Custom weights for the least-squares function. Each key points to a numpy.ndarray that is defined on the regularization.RegularizationMesh. A (n_cells, ) numpy.ndarray is used to define weights at cell centers, which are averaged to the appropriate faces internally when weighting is applied. A (n_faces, ) numpy.ndarray is used to define weights directly on the faces specified by the orientation input argument.

irls_scaledbool

If True, scale the IRLS weights to preserve magnitude of the regularization function. If False, do not scale.

irls_thresholdfloat

Constant added to IRLS weights to ensures stability in the algorithm.

gradient_type{“total”, “component”}

Gradient measure used in the IRLS re-weighting. Whether to re-weight using the total gradient or components of the gradient.

Attributes

`W`	Weighting matrix.
`active_cells`	Active cells defined on the regularization mesh.
`cell_gradient`	Partial cell gradient operator.
`cell_weights`	Deprecated property for 'volume' and user defined weights.
`gradientType`	gradient_type.gradientType has been deprecated.
`gradient_type`	Gradient measure used to update IRLS weights for sparse smoothness.
`indActive`	active_cells.indActive has been deprecated.
`irls_scaled`	Scale IRLS weights.
`irls_threshold`	Stability constant for computing IRLS weights.
`mapping`	Mapping from the inversion model parameters to the regularization mesh.
`model`	The model parameters.
`mref`	reference_model.mref has been deprecated.
`nP`	Number of model parameters.
`norm`	Norm for the sparse regularization.
`orientation`	Direction along which smoothness is enforced.
`parent`	The parent objective function
`reference_model`	Reference model.
`regmesh`	regularization_mesh.regmesh has been deprecated.
`regularization_mesh`	Regularization mesh.
`units`	Units for the model parameters.
`weights_keys`	Return the keys for the existing cell weights

reference_model_in_smooth

Methods

`__call__`(m)	Evaluate the regularization function for the model provided.
`deriv`(m)	Gradient of the regularization function evaluated for the model provided.
`deriv2`(m[, v])	Hessian of the regularization function evaluated for the model provided.
`f_m`(m)	Evaluate the regularization kernel function.
`f_m_deriv`(m)	Derivative of the regularization kernel function.
`get_lp_weights`(f_m)	Compute and return iteratively re-weighted least-squares (IRLS) weights.
`get_weights`(key)	Cell weights for a given key.
`map_class`	alias of `IdentityMap`
`remove_weights`(key)	Removes the weights for the key provided.
`set_weights`(**weights)	Adds (or updates) the specified weights to the regularization.
`test`([x, num, random_seed])	Run a convergence test on both the first and second derivatives.
`update_weights`(m)	Update the IRLS weights for sparse smoothness regularization.

Notes

The regularization function (objective function) for sparse smoothness (blockiness) along the x-direction as:

\[\phi (m) = \int_\Omega \, w(r) \, \Bigg | \, \frac{\partial m}{\partial x} \, \Bigg |^{p(r)} \, dv\]

where \(m(r)\) is the model, \(w(r)\) is a user-defined weighting function and \(p(r) \in [0,2]\) is a parameter which imposes sparseness throughout the recovered model. Sharper boundaries are recovered in regions where \(p(r)\) is small. If the same level of sparseness is being imposed everywhere, the exponent becomes a constant.

For implementation within SimPEG, the regularization function and its variables must be discretized onto a mesh. The discrete approximation for the regularization function (objective function) is expressed in linear form as:

\[\phi (\mathbf{m}) = \sum_i \tilde{w}_i \, \Bigg | \, \frac{\partial m_i}{\partial x} \, \Bigg |^{p_i}\]

where \(m_i \in \mathbf{m}\) are the discrete model parameters defined on the mesh. \(\tilde{w}_i \in \mathbf{\tilde{w}}\) are amalgamated weighting constants that 1) account for cell dimensions in the discretization and 2) apply user-defined weighting. \(p_i \in \mathbf{p}\) define the norm for each face (set using norm).

It is impractical to work with the general form directly, as its derivatives with respect to the model are non-linear and discontinuous. Instead, the iteratively re-weighted least-squares (IRLS) approach is used to approximate the sparse norm by iteratively solving a set of convex least-squares problems. For IRLS iteration \(k\), we define:

\[\phi \big (\mathbf{m}^{(k)} \big ) = \sum_i \tilde{w}_i \, \Bigg | \, \frac{\partial m_i^{(k)}}{\partial x} \Bigg |^{p_i} \approx \sum_i \tilde{w}_i \, r_i^{(k)} \Bigg | \, \frac{\partial m_i^{(k)}}{\partial x} \Bigg |^2\]

where the IRLS weight \(r_i\) for iteration \(k\) is given by:

\[r_i^{(k)} = \Bigg [ \Bigg ( \frac{\partial m_i^{(k-1)}}{\partial x} \Bigg )^2 + \epsilon^2 \; \Bigg ]^{{p_i}/2 - 1}\]

and \(\epsilon\) is a small constant added for stability (set using irls_threshold). For the set of model parameters \(\mathbf{m}\) defined at cell centers, the objective function for IRLS iteration \(k\) can be expressed as follows:

\[\phi \big ( \mathbf{m}^{(k)} \big ) \approx \Big \| \, \mathbf{W}^{(k)} \, \mathbf{G_x} \, \mathbf{m}^{(k)} \Big \|^2\]

where

\(\mathbf{m}^{(k)}\) are the discrete model parameters at iteration \(k\),

\(\mathbf{G_x}\) is the partial cell-gradient operator along x (x-derivative),

\(\mathbf{W}^{(k)}\) is the weighting matrix for iteration \(k\). It applies the IRLS weights, user-defined weighting, and accounts for cell dimensions when the regularization function is discretized.

Note that since \(\mathbf{G_x}\) maps from cell centers to x-faces, the weighting matrix acts on variables living on x-faces.

Reference model in smoothness:

Gradients/interfaces within a discrete reference model \(\mathbf{m}^{(ref)}\) can be preserved by including the reference model the smoothness regularization. In this case, the least-squares problem for IRLS iteration \(k\) becomes:

\[\phi \big ( \mathbf{m}^{(k)} \big ) \approx \Big \| \, \mathbf{W}^{(k)} \mathbf{G_x} \big [ \mathbf{m}^{(k)} - \mathbf{m}^{(ref)} \big ] \Big \|^2\]

This functionality is used by setting \(\mathbf{m}^{(ref)}\) with the reference_model property, and by setting the reference_model_in_smooth parameter to True.

IRLS weights, user-defined weighting and the weighting matrix:

Let \(\mathbf{w_1, \; w_2, \; w_3, \; ...}\) each represent an optional set of custom weights defined on the faces specified by the orientation property; i.e. x-faces for smoothness along the x-direction. Each set of weights were either defined directly on the faces or have been averaged from cell centers. And let \(\mathbf{r_x}^{\!\! (k)}\) represent the IRLS weights for iteration \(k\). The net weighting applied within the objective function is given by:

\[\mathbf{w}^{(k)} = \mathbf{r_x}^{\!\! (k)} \odot \mathbf{v_x} \odot \prod_j \mathbf{w_j}\]

where \(\mathbf{v_x}\) are cell volumes projected to x-faces; i.e. where the x-derivative lives. For a description of how IRLS weights are updated at every iteration, see the documentation for update_weights().

The weighting matrix used to apply the weights is given by:

\[\mathbf{W}^{(k)} = \textrm{diag} \Big ( \sqrt{\mathbf{w}^{(k)} \, } \Big )\]

Each set of custom weights is stored within a dict as an numpy.ndarray. A (n_cells, ) numpy.ndarray is used to define weights at cell centers, which are averaged to the appropriate faces internally when weighting is applied. A (n_faces, ) numpy.ndarray is used to define weights directly on the faces specified by the orientation input argument. The weights can be set all at once during instantiation with the weights keyword argument as follows:

>>> array_1 = np.ones(mesh.n_cells)  # weights at cell centers
>>> array_2 = np.ones(mesh.n_faces_x)  # weights directly on x-faces
>>> reg = SparseSmoothness(
>>>     mesh, orientation='x', weights={'weights_1': array_1, 'weights_2': array_2}
>>> )

or set after instantiation using the set_weights method:

>>> reg.set_weights(weights_1=array_1, weights_2=array_2})