/[escript]/trunk/doc/inversion/ForwardMagnetic.tex
ViewVC logotype

Contents of /trunk/doc/inversion/ForwardMagnetic.tex

Parent Directory Parent Directory | Revision Log Revision Log


Revision 4125 - (show annotations)
Wed Jan 2 06:15:00 2013 UTC (6 years, 7 months ago) by gross
File MIME type: application/x-tex
File size: 9042 byte(s)
some fixes in the inversion set-up
1
2 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3 % Copyright (c) 2003-2012 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/osl-3.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
16 \section{Magnetic Inversion}\label{sec:forward magnetic}
17 For the magnetic inversion we use the anomaly of the magnetic flux
18 density~\index{magnetic flux density} of the Earth.
19 The controlling material parameter is the susceptibility~\index{susceptibility}
20 $k$ of the rock.
21 With magnetization $M$ and inducing magnetic field anomaly $H^s$, the magnetic
22 flux density anomaly $B^s$ is given as
23 \begin{equation}\label{ref:MAG:EQU:1}
24 B_i = \mu_0 \cdot ( H^s_i + M_i )
25 \end{equation}
26 where $\mu_0 = 4 \pi \cdot 10^{-7} \frac{Vs}{Am}$.
27 In this forward model we make the simplifying assumption that the magnetization
28 is proportional to the known geomagnetic flux density $B^b$:
29 \begin{equation}\label{ref:MAG:EQU:4}
30 \mu_0 \cdot M_i = k \cdot B^b_i \;.
31 \end{equation}
32 Values for the magnetic flux density can be obtained by the International
33 Geomagnetic Reference Field (IGRF)~\cite{IGRF}
34 (or the Australian Geomagnetic Reference Field (AGRF)~\cite{AGRF}).
35 A rough approximation at latitude $\theta$ is given by
36 \begin{equation}\label{ref:MAG:EQU:5}
37 \begin{array}{rcl}
38 B^b_{\theta} & = & \displaystyle{ \frac{ \mu_0 \cdot m_{earth}}{4 \pi \cdot R_{earth}^3} sin(\theta) } \\
39 B^b_r & = & \displaystyle{ \frac{\mu_0 \cdot m_{earth}}{1 \pi \cdot R_{earth}^3} cos(\theta) }
40 \end{array}
41 \end{equation}
42 with the vacuum permeability\index{vacuum permeability} $\mu_0 = 4 \pi \cdot 10^{-7} \frac{Vs}{Am}$,
43 the magnetic dipole moment of Earth $m_{earth}=8.22 \cdot 10^{22} Am^2$ and
44 earth radius $R_{earth}= 6378137m$.
45 $B^b_r$ and $b^b_{\theta}$ denote the radial and latitudinal component of the
46 geomagnetic flux density.
47 Notice that convention~(\ref{REF:EQU:INTRO 9}) applies if Cartesian
48 coordinates\index{Cartesian coordinates} are used.
49 In most cases it is reasonable to assume that that the background field is
50 constant across the domain.
51
52 The magnetic field anomaly $H^s$ can be represented by the gradient of a
53 magnetic scalar potential\index{scalar potential!magnetic} $\psi$.
54 We use the form
55 \begin{equation}\label{ref:MAG:EQU:6}
56 \mu_0 \cdot H^s_i = - \psi_{,i}
57 \end{equation}
58 With this notation one gets from Equations~(\ref{ref:MAG:EQU:1}) and~(\ref{ref:MAG:EQU:4}):
59 \begin{equation}\label{ref:MAG:EQU:7}
60 B_i = - \psi_{,i} + k \cdot B^b_i
61 \end{equation}
62 As the $B^s$ magnetic flux density anomaly we obtain the PDE
63 \begin{equation}\label{ref:MAG:EQU:8}
64 - \psi_{,ii} = - (k B^b_i)_{,i}
65 \end{equation}
66 which needs to be solved for a given susceptibility $k$.
67 The magnetic scalar potential is set to zero at the top of the domain
68 $\Gamma_{0}$.
69 On all other faces the normal component of the magnetic flux density anomaly
70 $B_i$ is set to zero, i.e. $n_i \psi_{,i} = k \cdot n_i B^b_i$ with outer
71 normal field $n_i$.
72
73 From the magnetic scalar potential we can calculate the magnetic flux density
74 anomaly via Equation~(\ref{ref:MAG:EQU:8}) to calculate the defect to the given
75 data.
76 If $B^{(s)}_i$ is a measurement of the magnetic flux density anomaly for
77 survey $s$ and $\omega^{(s)}_i$ is a weighting factor the data defect
78 $J^{mag}(k)$ in the notation of Chapter~\ref{Chp:ref:introduction} is given as
79 \begin{equation}\label{ref:MAG:EQU:9}
80 J^{mag}(k) = \frac{1}{2}\sum_{s} \int_{\Omega} ( \omega^{(s)}_i \cdot (B_{i}- B^{(s)}_i) ) ^2 dx
81 \end{equation}
82 Summation over $i$ is performed.
83 The cost function kernel\index{cost function!kernel} is given as
84 \begin{equation}\label{ref:MAG:EQU:10}
85 K^{mag}(\psi_{,i},k) = \frac{1}{2}\sum_{s} ( \omega^{(s)}_i \cdot (k \cdot B^b_i - \psi_{,i} - B^{(s)}_i) ) ^2
86 \end{equation}
87 Notice that if magnetic flux density is measured in air one can ignore the
88 $k\cdot B^b_i$ as the susceptibility is zero.
89
90 In practice the magnetic flux density $b^{(s)}$ is measured along a certain
91 direction $d^{(s)}_i$ with a standard error deviation $\sigma^{(s)}$ at
92 certain locations in the domain.
93 In this case one sets $B^{(s)}_i=b^{(s)} \cdot d^{(s)}_i$ and the weighting
94 factors $\omega^{(s)}$ as
95 \begin{equation}\label{ref:MAG:EQU:11}
96 \omega^{(s)}_i
97 = \left\{
98 \begin{array}{lcl}
99 f \cdot \frac{d^{(s)}_i}{\sigma^{(s)}} & & \mbox{data are available} \\
100 & \mbox{ where } & \\
101 0 & & \mbox{ otherwise } \\
102 \end{array}
103 \right.
104 \end{equation}
105 where it is assumed that $d^{(s)}_i \cdot d^{(s)}_i =1$. With the objective to control the
106 gradient of the cost function the scaling factor $f$ is chosen in the way that
107 \begin{equation}\label{ref:MAG:EQU:12}
108 \sum_{s} \int_{\Omega} ( \omega^{(s)}_i B^{(s)}_i )
109 \cdot ( \omega^{(s)}_j \frac{1}{L_j} ) \cdot L^2 \cdot
110 ( B^b_n \frac{1}{L_n} )
111 \cdot k' \;
112 dx =\alpha
113 \end{equation}
114 where $\alpha$ defines a scaling factor which is typically set to one and $L$ is defined by equation~(\ref{ref:EQU:REG:6b}).
115 $k'$ is considering the
116 derivative of the density with respect to the level set function.
117
118 \subsection{Usage}
119
120 \LG{Add example}
121
122 \begin{classdesc}{MagneticModel}{domain, w, B, background_field,
123 \optional{, useSphericalCoordinates=False}
124 \optional{, tol=1e-8}}
125 opens a magnetic forward model over the \Domain \member{domain} with
126 weighting factors \member{w} ($=\omega^{(s)}$) and measured magnetic flux
127 density anomalies \member{B} ($=B^{(s)}$).
128 The weighting factors and the measured magnetic flux density anomalies must be vectors.
129 \member{background_field} defines the background magnetic flux density $B^b$
130 as a vector.
131 If \member{useSphericalCoordinates} is \True spherical coordinates are used.
132 \member{tol} sets the tolerance for the solution of the PDE~(\ref{ref:MAG:EQU:8}).
133 \end{classdesc}
134
135 \begin{methoddesc}[MagneticModel]{rescaleWeights}{
136 \optional{scale=1.}
137 \optional{k_scale=1.}}
138 rescale the weighting factors such condition~(\ref{ref:MAG:EQU:12}) holds where
139 \member{scale} sets the scale $\alpha$
140 and \member{k_scale} sets $k'$. This method should be called before any inversion is started
141 in order to make sure that all components of the cost function are appropriately scaled.
142 \end{methoddesc}
143
144
145 \subsection{Gradient Calculation}
146 This section briefly explains how the gradient
147 $\frac{\partial J^{mag}}{\partial k}$ of the cost function $J^{mag}$ with
148 respect to the susceptibility $k$ is calculated.
149 The magnetic potential $\psi$ from PDE~(\ref{ref:MAG:EQU:8}) is solved in weak form:
150 \begin{equation}\label{ref:MAG:EQU:201}
151 \int_{\Omega} q_{,i} \psi_{,i} \; dx = \int_{\Omega} k \cdot q_{,i} B^b_i \; dx
152 \end{equation}
153 for all $q$ with $q=0$ on $\Gamma_{0}$.
154 In the following we set $\Psi[k]=\psi$ for a given susceptibility $k$ as
155 solution of the variational problem~(\ref{ref:MAG:EQU:201}).
156 If $\Gamma_{k}$ denotes the region of the domain where the susceptibility is
157 known and for a given direction $p$ with $p=0$ on $\Gamma_{k}$ one has
158 \begin{equation}\label{ref:MAG:EQU:201aa}
159 \int_{\Omega} \frac{\partial J^{mag}}{\partial k} \cdot p \; dx = \int_{\Omega}
160 \sum_{s} (\omega^{(s)}_j
161 ( B^{(s)}_j-B_{j})) \cdot ( \omega^{(s)}_i ( \Psi[p]_{,i} - p \cdot B^b_i ) ) \; dx
162 \end{equation}
163 with
164 \begin{equation}\label{ref:MAG:EQU:202c}
165 Y_i[\psi]= \sum_{s} (\omega^{(s)}_j
166 (B^{(s)}_j - B_{j}) ) \cdot \omega^{(s)}_i
167 \end{equation}
168 This is written as
169 \begin{equation}\label{ref:MAG:EQU:202cc}
170 \int_{\Omega} \frac{\partial J^{mag}}{\partial k} \cdot p \; dx = \int_{\Omega}
171 Y_i[\psi] \Psi[p]_{,i} - p \cdot Y_i[\psi]B^b_i \; dx
172 \end{equation}
173 We then set $Y^*[\psi]$ as the solution of the equation
174 \begin{equation}\label{ref:MAG:EQU:202d}
175 \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_{0}
176 \end{equation}
177 with $Y^*[\psi]=0$ on $\Gamma_{0}$. With $r=\Psi[p]$ we get
178 \begin{equation}\label{ref:MAG:EQU:202dd}
179 \int_{\Omega} \Psi[p]_{,i} Y^*[\psi]_{,i} \; dx = \int_{\Omega} \Psi[p]_{,i} ,Y_i[\psi] \; dx
180 \end{equation}
181 and from Equation~(\ref{ref:MAG:EQU:201}) with $q=Y^*[\psi]$ we get
182 \begin{equation}\label{ref:MAG:EQU:20e}
183 \int_{\Omega} Y^*[\psi]_{,i} \Psi[p]_{,i} \; dx = \int_{\Omega} p \cdot Y^*[\psi]_{,i} B^b_i \; dx
184 \end{equation}
185 which leads to
186 \begin{equation}\label{ref:MAG:EQU:20ee}
187 \int_{\Omega} \Psi[p]_{,i} ,Y_i[\psi] \; dx = \int_{\Omega} p \cdot Y^*[\psi]_{,i} B^b_i \; dx
188 \end{equation}
189 and finally
190 \begin{equation}\label{ref:MAG:EQU:201a}
191 \int_{\Omega} \frac{\partial J^{mag}}{\partial k} \cdot p \; dx = \int_{\Omega}
192 p \cdot (Y_i[\psi] - Y^*[\psi]_{,i}) B^b_i \; dx
193 \end{equation}
194 or
195 \begin{equation}\label{ref:MAG:EQU:201b}
196 \frac{\partial J^{mag}}{\partial k} = (Y^*[\psi]_{,i}-Y_i[\psi]) B^b_i
197 \end{equation}
198

  ViewVC Help
Powered by ViewVC 1.1.26