1 

2 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
3 
% Copyright (c) 20032013 by University of Queensland 
4 
% http://www.uq.edu.au 
5 
% 
6 
% Primary Business: Queensland, Australia 
7 
% Licensed under the Open Software License version 3.0 
8 
% http://www.opensource.org/licenses/osl3.0.php 
9 
% 
10 
% Development until 2012 by Earth Systems Science Computational Center (ESSCC) 
11 
% Development since 2012 by School of Earth Sciences 
12 
% 
13 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
14 

15 
\section{Gravity Inversion}\label{sec:forward gravity} 
16 
For the magnetic inversion we use the anomaly of the gravity acceleration~\index{gravity acceleration} of the Earth. 
17 
The controlling material parameter is the density~\index{density} $\rho$ of 
18 
the rock. 
19 
If the density field $\rho$ is known the gravitational potential $\psi$ is 
20 
given as the solution of the PDE 
21 
\begin{equation}\label{ref:GRAV:EQU:100} 
22 
\psi_{,ii} = 4\pi G \cdot \rho 
23 
\end{equation} 
24 
where $G=6.6730 \cdot 10^{11} \frac{m^3}{kg \cdot s^2}$ is the gravitational 
25 
constant. 
26 
The gravitational potential is set to zero at the top of the 
27 
domain $\Gamma_0$. 
28 
On all other faces the normal component of the gravity acceleration anomaly 
29 
$g_i$ is set to zero, i.e. $n_i \psi_{,i} = 0$ with outer normal field $n_i$. 
30 
The gravity force $g_i$ is given as the negative of the gradient of the gravity 
31 
potential $\psi$: 
32 
\begin{equation}\label{ref:GRAV:EQU:101} 
33 
g_i =  \psi_{,i} 
34 
\end{equation} 
35 
From the gravitational potential we can calculate the gravity acceleration 
36 
anomaly via Equation~(\ref{ref:GRAV:EQU:101}) to obtain the defect to the 
37 
given data. 
38 
If $g^{(s)}_i$ is a measurement of the gravity acceleration anomaly for 
39 
survey $s$ and $\omega^{(s)}_i$ is a weighting factor the data defect 
40 
$J^{grav}(k)$ in the notation of Chapter~\ref{chapter:ref:inversion cost function} is given as 
41 
\begin{equation}\label{ref:GRAV:EQU:9} 
42 
J^{grav}(k) = \frac{1}{2}\sum_{s} \int_{\Omega} ( \omega^{(s)}_i \cdot (g_{i} g^{(s)}_i) ) ^2 dx 
43 
\end{equation} 
44 
Summation over $i$ is performed. 
45 
The cost function kernel\index{cost function!kernel} is given as 
46 
\begin{equation}\label{ref:GRAV:EQU:10} 
47 
K^{grav}(\psi_{,i},k) = \frac{1}{2}\sum_{s} ( \omega^{(s)}_i \cdot (\psi_{,i}+ g^{(s)}_i) ) ^2 
48 
\end{equation} 
49 
In practice the gravity acceleration $g^{(s)}$ is measured in vertical 
50 
direction $z$ with a standard error deviation $\sigma^{(s)}$ at certain 
51 
locations in the domain. 
52 
In this case one sets the weighting factors $\omega^{(s)}$ as 
53 
\begin{equation}\label{ref:GRAV:EQU:11} 
54 
\omega^{(s)}_i 
55 
= \left\{ 
56 
\begin{array}{lcl} 
57 
f \cdot \frac{\delta_{iz}}{\sigma^{(s)}} & & \mbox{data are available} \\ 
58 
& \mbox{ where } & \\ 
59 
0 & & \mbox{ otherwise } \\ 
60 
\end{array} 
61 
\right. 
62 
\end{equation} 
63 
With the objective to control the 
64 
gradient of the cost function 
65 
the scaling factor $f$ is chosen in the way that 
66 
\begin{equation}\label{ref:GRAV:EQU:12} 
67 
\sum_{s} \int_{\Omega} ( \omega^{(s)}_i g^{(s)}_i ) \cdot ( \omega^{(s)}_j \frac{1}{L_j} ) \cdot 4\pi G L^2 \cdot \rho' \; dx =\alpha 
68 
\end{equation} 
69 
where $\alpha$ defines a scaling factor which is typically set to one and $L$ is defined by equation~(\ref{ref:EQU:REG:6b}). $\rho'$ is considering the 
70 
derivative of the density with respect to the level set function. 
71 

72 

73 
\subsection{Usage} 
74 

75 
\LG{Add example} 
76 

77 
\begin{classdesc}{GravityModel}{domain, 
78 
w, g, 
79 
\optional{, useSphericalCoordinates=False} 
80 
\optional{, fixPotentialAtBottom=False}, 
81 
\optional{, tol=1e8} 
82 
} 
83 
opens a gravity forward model over the \Domain \member{domain} with 
84 
weighting factors \member{w} ($=\omega^{(s)}$) and measured gravity acceleration anomalies\member{g} ($=g^{(s)}$). 
85 
The weighting factors and the measured gravity acceleration anomalies\member must be vectors. 
86 
If \member{useSphericalCoordinates} is \True spherical coordinates are used. 
87 
\member{tol} set the tolerance for the solution of the PDE~(\ref{ref:GRAV:EQU:100}). 
88 
If \member{fixPotentialAtBottom} is set to \True, the gravitational potential 
89 
at the bottom is set to zero in addition to the potential on the top. 
90 
\end{classdesc} 
91 

92 
\begin{methoddesc}[GravityModel]{rescaleWeights}{ 
93 
\optional{scale=1.} 
94 
\optional{rho_scale=1.}} 
95 
rescale the weighting factors such condition~(\ref{ref:GRAV:EQU:12}) holds where 
96 
\member{scale} sets the scale $\alpha$ 
97 
and \member{rho_scale} sets $\rho'$. This method should be called before any inversion is started 
98 
in order to make sure that all components of the cost function are appropriately scaled. 
99 
\end{methoddesc} 
100 

101 

102 
\subsection{Gradient Calculation} 
103 
This section briefly explains how the gradient 
104 
$\frac{\partial J^{grav}}{\partial \rho}$ of the cost function $J^{grav}$ with 
105 
respect to the density $\rho$ is calculated. We follow the concept as outlined in section~\ref{chapter:ref:inversion cost function:gradient}. 
106 
The magnetic potential $\psi$ from PDE~(\ref{ref:GRAV:EQU:100}) is solved in 
107 
weak form: 
108 
\begin{equation}\label{ref:GRAV:EQU:201} 
109 
\int_{\Omega} q_{,i} \psi_{,i} \; dx =  \int_{\Omega} 4\pi G \cdot q \rho\; dx 
110 
\end{equation} 
111 
for all $q$ with $q=0$ on $\Gamma_0$. 
112 
In the following we set $\Psi[\cdot]=\psi$ for a given density $\cdot$ as 
113 
solution of the variational problem~(\ref{ref:GRAV:EQU:201}). 
114 
If $\Gamma_{\rho}$ denotes the region of the domain where the density is known 
115 
and for a given direction $p$ with $p=0$ on $\Gamma_{\rho}$ one has 
116 
\begin{equation}\label{ref:GRAV:EQU:201aa} 
117 
\int_{\Omega} \frac{\partial J^{grav}}{\partial \rho} \cdot p \; dx = \int_{\Omega} 
118 
\sum_{s} (\omega^{(s)}_j \cdot 
119 
(g^{(s)}_jg_{j}) ) \cdot ( \omega^{(s)}_i \Psi[p]_{,i}) \; dx 
120 
\end{equation} 
121 
with 
122 
\begin{equation}\label{ref:GRAV:EQU:202c} 
123 
Y_i[\psi]= \sum_{s} (\omega^{(s)}_j \cdot 
124 
(g^{(s)}_jg_{j}) ) \cdot \omega^{(s)}_i 
125 
\end{equation} 
126 
This is written as 
127 
\begin{equation}\label{ref:GRAV:EQU:202cc} 
128 
\int_{\Omega} \frac{\partial J^{grav}}{\partial \rho} \cdot p \; dx = \int_{\Omega} 
129 
Y_i[\psi] \Psi[p]_{,i} \; dx 
130 
\end{equation} 
131 
We then set $Y^*[\psi]$ as the solution of the equation 
132 
\begin{equation}\label{ref:GRAV:EQU:202d} 
133 
\int_{\Omega} r_{,i} Y^*[\psi]_{,i} \; dx = \int_{\Omega} r_{,i} ,Y_i[\psi] \; dx \mbox{ for all } p \mbox{ with } r=0 \mbox{ on } \Gamma_{top} 
134 
\end{equation} 
135 
with $Y^*[\psi]=0$ on $\Gamma_0$. With $r=\Psi[p]$ we get 
136 
\begin{equation}\label{ref:GRAV:EQU:202dd} 
137 
\int_{\Omega} \Psi[p]_{,i} Y^*[\psi]_{,i} \; dx = \int_{\Omega} \Psi[p]_{,i} ,Y_i[\psi] \; dx 
138 
\end{equation} 
139 
and from Equation~(\ref{ref:GRAV:EQU:201}) with $q=Y^*[\psi]$ we get 
140 
\begin{equation}\label{ref:GRAV:EQU:20e} 
141 
\int_{\Omega} Y^*[\psi]_{,i} \Psi[p]_{,i} \; dx =  \int_{\Omega} 4\pi G \cdot Y^*[\psi] \cdot p\; dx 
142 
\end{equation} 
143 
which leads to 
144 
\begin{equation}\label{ref:GRAV:EQU:20ee} 
145 
\int_{\Omega} \Psi[p]_{,i} ,Y_i[\psi] \; dx =  \int_{\Omega} 4\pi G \cdot Y^*[\psi] \cdot p \; dx 
146 
\end{equation} 
147 
and finally 
148 
\begin{equation}\label{ref:GRAV:EQU:201a} 
149 
\int_{\Omega} \frac{\partial J^{grav}}{\partial \rho} \cdot p \; dx =  \int_{\Omega} 
150 
4\pi G \cdot Y^*[\psi] \cdot p \; dx 
151 
\end{equation} 
152 
or 
153 
\begin{equation}\label{ref:GRAV:EQU:201b} 
154 
\frac{\partial J^{grav}}{\partial \rho} = 4\pi G \cdot Y^*[\psi] 
155 
\end{equation} 
156 
