Contents of /branches/doubleplusgood/doc/inversion/ForwardGravity.tex

Revision 4345 - (show annotations)
Fri Mar 29 07:09:41 2013 UTC (5 years, 11 months ago) by jfenwick
File MIME type: application/x-tex
File size: 7183 byte(s)
Spelling fixes

 1 2 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 3 % Copyright (c) 2003-2013 by University of Queensland 4 5 % 6 % Primary Business: Queensland, Australia 7 % Licensed under the Open Software License version 3.0 8 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=1e-8} 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)}_j-g_{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)}_j-g_{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