I’ve noticed that in PF mag/gravity, SimPEG is using the integral formulation instead of the differential formulation used for almost all the other types of the SimPEGs problems.
My questions are:
- Why using integral formulation in this case ?
- What are the +/- of the integral VS differential formulation
- Are they leading to the same results in the end ?
Thank you very much for your time !
Hi @1apsus ,
Great question. Technically SimPEG has both the integral and PDE simulations, but we have been using the integral mostly for efficiency and legacy reasons. If you look at papers on grav/mag you will likely find that most of them use the integral formulation as well.
The main reason is that grav and mag are (in most cases) linear problems. So you only have to compute the forward operator once - no need to solve the forward again on model updates. And since linear, the forward operator is also the sensitivity matrix for the inverse. Once you have the sensitivities you get a few nice advantages: solving the CG iterations with simple matrix-vector operations, Jacobi pre-conditioning to speed up convergence, regularization weighting etc. The integral formulation also doesn’t require worrying about boundary conditions, so the mesh can be smaller.
That being said, PDE also has advantages. The size of the problem is technically not dependent on the number of receivers (only mesh size). It can also deal with non-linear magnetization effects (self-demag).
I might be missing a couple of points here, but I will let others chime in.
Thank you @domfournier for the insights !
By any chance, do you know from which paper the SimPEG integral formulation of the gravity or magnetics come from ? Actually, I try to understand how the forward operator is built from the integral form (I’ve looked to the “evaluate_integral” function of the simulation, but that’s not so clear for me yet)
Hey sorry for the delay:
Nagy, D., 1966, The gravitational attraction of a right rectangular prism: Geophysics, 31, 361–371.
Basically, it’s this in vector form:
And for the tensor components, you can take the partial derivatives…