Inversion of magnetization amplitude with different direction than external field

Potential Fields
1

Hi,
I would like to know if this is possible to do an inversion of magnetization or pseudo-susceptibility assuming the magnetization direction is different from the external field. I have tested the mvi inversion by imposing a starting and reference model with a fixed direction of magnetization different from the external field and i see the recovered magnetization remain close to the starting magnetization direction which is what i wanted but if possible i would like to keep the magnetization direction completely fixed, in order to limit the number of free parameters in the inversion. Is this possible with simpeg and if so how could I do that ?
Thank you,
Claire.

2

Hi Claire,

This is not something built in but I wrote a mapping class that can do this. I attached a snippet here. You would use simulation_scalar to run your inversion, and pass scalar_map to your regularization. If you have any questions on it let me know.

Cheers,

John

from simpeg.maps import IdentityMap
from simpeg import utils
from simpeg.utils import mkvc,sdiag
import scipy.sparse as sp


class MagnetizationMap(IdentityMap):

    def __init__(self, mesh=None, nP=None, inducing_field = None, **kwargs):
        super(MagnetizationMap, self).__init__(mesh=mesh, nP=nP, **kwargs)
        self.inducing_field = inducing_field
        num = np.ones(nP)
        Mx = sp.diags(num*inducing_field[0])
        My = sp.diags(num*inducing_field[1])
        Mz = sp.diags(num*inducing_field[2])
        self.M = sp.vstack([Mx, My, Mz])

    @property
    def shape(self):

        return (self.nP * 3, self.nP)


    def _transform(self, m):

        return self.M * m

    def deriv(self, m, v=None):
        if v is not None:
            return self.M * v
        return self.M * sp.eye(self.nP)

actv = active_from_xyz(mesh, topo)
nC = int(actv.sum())

dip = 45.0
azm = 90.0

dip_model = dip*np.ones(nC)
azm_model = azm*np.ones(nC)

#model amp is the scalar amplitude throughout the mesh
model_vector = sdiag(model_amp) * utils.mat_utils.dip_azimuth2cartesian(
    dip_model, azm_model
)

# Map you use to invert for vector inversion
vectorMap = maps.IdentityMap(nP=nC * 3)

# Pass this map to your regularization for scalar inversion
scalarMap = maps.IdentityMap(nP=nC)

mag_map_indfield = 1 * utils.mat_utils.dip_azimuth2cartesian(
    dip, azm
)[0]

# pass this map to scalar  forward simulation
mag_map = MagnetizationMap(nP=nC,inducing_field=mag_map_indfield)


#setup vector simulation
simulation_vector = magnetics.simulation.Simulation3DIntegral(
    survey=survey, mesh=mesh, chiMap=vectorMap, active_cells=actv, model_type="vector",
)

#setup scalar simulation
simulation_scalar = magnetics.simulation.Simulation3DIntegral(
    survey=survey, mesh=mesh, chiMap=mag_map, active_cells=actv, model_type="vector",
)

#check that they give the same result
d1 = simulation_vector.dpred(mkvc(model_vector))
d2 = simulation_scalar.dpred(model_amp)

assert(np.allclose(d1, d2))
3

Dear John,
Thank you very much for your prompt answer.
If we want to continue your script with the inversion step, can we then simply follow the instructions in the tutorial for susceptibility inversion at 3D Inversion of TMI Data to Recover a Susceptibility Model - SimPEG User Tutorials with the following definitions :

  • dmis = data_misfit.L2DataMisfit(data=data_object, simulation=simulation_scalar)
  • reference_tensor_model and starting_tensor_model defined as magnetization amplitude 1D array similar to the model_amp array

Is that all or are there additional changes that i am missing ?
Thank you,
Claire.

4

Hi Claire,

Yes that is right for all the other pieces you can just run a normal susceptibility style inversion. The only thing that changes is the simulation which uses that mapping.

John