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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 4095 - (show annotations)
Wed Dec 5 05:32:22 2012 UTC (7 years, 1 month ago) by caltinay
File MIME type: application/x-tex
File size: 10044 byte(s)
A bit of doco cleanup.

1 \chapter{Regularization}\label{Chp:ref:regularization}
2
3 The general cost function $J_{total}$ to be minimized has some of the cost
4 function $J_{forward}$ measuring the defect of the result from the
5 forward model with the data, and the cost function $J_{reg}$ introducing the
6 regularization into the problem and makes sure that a unique answer exists.
7 The regularization term is a function of, possibly vector-valued, level set
8 function $m$ which represents the physical properties to be represented and is,
9 from a mathematical point of view, the unknown of the inversion problem.
10 It is the intention that the values of $m$ are between zero and one and that
11 actual physical values are created from a mapping before being fed into a
12 forward model. In general the cost function $J_{reg}$ is defined as
13 \begin{equation}\label{EQU:REG:1}
14 J_{reg}(m) = \frac{1}{2} \int_{\Omega}
15 \sum_{k} \mu^{(0)}_k \cdot s^{(0)}_k \cdot m_k^2 + \mu^{(1)}_{ki} \cdot s^{(1)}_{ki} \cdot L_i^2 \cdot m_{k,i} \cdot m_{k,i}
16 + \sum_{l<k} \mu^{(c)}_{lk} \cdot s^{(c)}_{lk} \cdot L^4 \cdot \sigma(m_l,m_k) dx
17 \end{equation}
18 where $s^{(0)}_k$, $s^{(1)}_{ki}$ and $s^{(c)}_{lk}$ are scaling factors with
19 values between zero and one (limits including).
20 They may vary with their location within the domain $\Omega$.
21 $L_i$ is the width of the domain in $x_i$ direction and $L^2=L_i \cdot L_i$.
22 $\sigma$ is a given symmetric, non-negative cross correlation function.
23 We use
24 \begin{equation}\label{EQU:REG:4}
25 \sigma(a,b) = ( a_{,i} \cdot a_{,i}) \cdot ( b_{,i} \cdot b_{,i}) - ( a_{,i} \cdot b_{,i})^2
26 \end{equation}
27 Notice that the cross correlation function is measuring the angle between the
28 surface normals of contours of level set functions.
29 So minimizing the cost function will align the surface normals of the contours.
30 The additional weight factors $\mu^{(0)}_k$, $ \mu^{(1)}_{ki}$ and
31 $\mu^{(c)}_{lk}$ are between zero and one and constant across the domain.
32 They are potentially modified during the inversion in order to improve the
33 balance between the different terms in the cost function.
34 Notice that values for $\mu^{(0)}_k$, $ \mu^{(1)}_{ki}$ and $\mu^{(c)}_{lk}$
35 are relevant for which the values of the corresponding entries in scaling
36 factors $s^{(0)}_k$, $s^{(1)}_{ki}$ and $s^{(c)}_{lk}$ are non-zero.
37 Also notice that the factors $L^4$ and $L_i^2$ take care of any length scale changes.
38 With the notation
39 \begin{equation}\label{EQU:REG:2}
40 \begin{array}{rcl}
41 \widehat{s}^{(0)}_k & = & \mu^{(0)}_k \cdot w^{(0)}_k \\
42 \widehat{s}^{(1)}_{ki} & = & \mu^{(1)}_{ki} \cdot w^{(1)}_{ki} \cdot L_i^2 \\
43 \widehat{s}^{(c)}_{lk} & = & \mu^{c}_{lk} \cdot w^{(c)}_{lk} \cdot L^4 \\
44 \end{array}
45 \end{equation}
46 with $k<l$ we can write
47 \begin{equation}\label{EQU:REG:1b}
48 J_{reg}(m) = \frac{1}{2} \int_{\Omega} \left(
49 \sum_{k} \widehat{s}^{(0)}_k \cdot m_k^2 + \widehat{s}^{(1)}_{ki} \cdot m_{k,i} \cdot m_{k,i}
50 + \sum_{l<k} \widehat{s}^{(c)}_{lk} \cdot \sigma(m_l,m_k) \right) dx
51 \end{equation}
52 We need to provide the derivative $\frac{ \partial J_{reg}}{\partial q}$ of
53 the cost function $J_{reg}$ with respect to a given direction $q$ which equals
54 zero at locations where $m$ is assumed to be zero.
55 The derivative is given as
56 \begin{equation}\label{EQU:REG:3}
57 \frac{ \partial J_{reg}}{\partial q}(m) =
58 \int_{\Omega} Y_k \cdot q_k + X_{k,i} q_{k,i} dx
59 \end{equation}
60 where
61
62
63
64 For a single-valued level set function this takes the form
65 \begin{equation}\label{EQU:REG:3}
66 \frac{ \partial J_{reg}}{\partial q}(m) =
67 \mu^{reg} \int_{\Omega} \omega \cdot m \cdot q + \omega_i \cdot L_i^2 \cdot m_{,i} \cdot q_{,i} dx
68 \end{equation}
69 So we can represent the gradient $\nabla J_{reg}$ of the cost function
70 $J_{reg}$ by the pair of values $(Y,X)$ where we set
71 \begin{equation}\label{EQU:REG:3b}
72 Y=\mu^{reg} \cdot \omega \cdot m \mbox{ and } X_i = \mu^{reg} \cdot \omega_i \cdot L_i^2 \cdot m_{,i}
73 \end{equation}
74 and
75 \begin{equation}\label{EQU:REG:3c}
76 \frac{ \partial J_{reg}}{\partial q}(m) = [ \nabla J_{reg}(m), q ] =
77 \int_{\Omega} Y \cdot q + X_i \cdot q_{,i} dx
78 \end{equation}
79 where $[.,.]$ is called the dual product.
80
81
82
83 ==============================================================
84
85 where $<.,.>$ is an inner product. If the level set function $m$ has several components $m_j$ the inner product $<.>$ is given
86 in the form
87 $\omega^{(k)}$ and $\omega^{(k)}_i$ are fixed non-negative weighting factors and $\mu^{reg}_k$ are weighting factors
88 which may be modified during the inversion. $L_i$ is the length of the domain in $x_i$ direction. In the special case that
89 the level set function $m$ has a single component the inner product takes the form
90 \begin{equation}\label{EQU:REG:2b}
91 <p,q> =
92 \mu^{reg} \int_{\Omega} \omega \cdot p \cdot q + \omega_i \cdot L_i^2 \cdot p_{,i} \cdot q_{,i} dx
93 \end{equation}
94 In practice it is assumed that the level set function is known to be zero in certain regions in the domain. Typically these regions
95 corresponds to region above the surface or regions explored by drilling.
96
97 We need to provide the derivative of the cost function $J_{reg}$ with respect to a given direction $q$ which equals zero at locations
98 where $m$ is assumed to be zero. For a single-valued
99 level set function th is takes the form
100 \begin{equation}\label{EQU:REG:3}
101 \frac{ \partial J_{reg}}{\partial q}(m) =
102 \mu^{reg} \int_{\Omega} \omega \cdot m \cdot q + \omega_i \cdot L_i^2 \cdot m_{,i} \cdot q_{,i} dx
103 \end{equation}
104 So we can represent the gradient $\nabla J_{reg}$ of the cost function $J_{reg}$ by the pair of values $(Y,X)$ where we set
105 \begin{equation}\label{EQU:REG:3b}
106 Y=\mu^{reg} \cdot \omega \cdot m \mbox{ and } X_i = \mu^{reg} \cdot \omega_i \cdot L_i^2 \cdot m_{,i}
107 \end{equation}
108 and
109 \begin{equation}\label{EQU:REG:3c}
110 \frac{ \partial J_{reg}}{\partial q}(m) = [ \nabla J_{reg}(m), q ] =
111 \int_{\Omega} Y \cdot q + X_i \cdot q_{,i} dx
112 \end{equation}
113 where $[.,.]$ is called the dual product.
114
115 For a multi-valued level set function an additional correlation term is introduced into the cost function $J_{total}$:
116 \begin{equation}\label{EQU:REG:1c}
117 J_{reg}(m) = \frac{1}{2} < m, m > + \frac{1}{2} \sum_{k,l} \mu_{kl}^{sec} \cdot \int_{\Omega} \sigma(m_k,m_l) dx
118 \end{equation}
119 where $sigma$ is a given symmetric, non-negative correlation function, and $\mu_{kl}^{sec}$ are symmetric, weighting factors
120 ($\mu_{kl}^{sec} = \mu_{lk}^{sec}$, $\mu_{kk}^{sec}=0$) which may
121 be altered during the inversion. We use the correlation function
122 \begin{equation}\label{EQU:REG:4}
123 \sigma(a,b) = \frac{L^2}{2} \cdot ( ( a_{,i} \cdot a_{,i}) \cdot ( b_{,i} \cdot b_{,i}) - ( a_{,i} \cdot b_{,i})^2 )
124 \end{equation}
125 with $L=L_i \cdot L_i$. Minimizing $J_{reg}(m)$ is minimizing the angle between the surface normals of the contours formed by
126 two level set function. the derivative of the cost function $J_{reg}$ with respect to a given direction $q$ which equals zero at locations
127 where $m$ is assumed to be zero:
128 \begin{equation}\label{EQU:REG:5}
129 \begin{array}{ll}
130 \displaystyle{\frac{ \partial J_{reg}}{\partial q}(m)} =
131 \displaystyle{\sum_{k} \mu^{reg}_k \int_{\Omega} \omega^{(k)} \cdot m_k \cdot q_k + \omega^{(k)}_i \cdot L_i^2 \cdot m_{k,i} \cdot q_{k,i} dx } \\
132 + \displaystyle{\sum_{k,l} \mu_{kl}^{sec} \cdot {L^2} \int_{\Omega} ( m_{k,i} \cdot q_{k,i}) \cdot ( m_{l,j} \cdot m_{l,j}) - ( m_{k,j} \cdot m_{l,j}) \cdot ( q_{l,i} \cdot m_{k,i}) } dx
133 \end{array}
134 \end{equation}
135 Similar to the single-case we can represent
136 the gradient $\nabla J_{reg}$ of the cost function $J_{reg}$ by the pair of values $(Y,X)$ where we set
137 \begin{equation}\label{EQU:REG:6}
138 Y_k= \mu^{reg}_k \cdot \omega^{(k)} \cdot m_k
139 \end{equation}
140 and
141 \begin{equation}\label{EQU:REG:6b}
142 X_{ki} = \mu^{reg}_k \cdot \omega^{(k)} \cdot L_i^2 \cdot m_{k,i} +
143 \sum_{l} \mu_{kl}^{sec} \cdot {L^2} \cdot ( ( m_{l,j} \cdot m_{l,j}) \cdot m_{k,i} - ( m_{l,j} \cdot m_{k,j}) \cdot m_{l,i} )
144 \end{equation}
145 and
146 \begin{equation}\label{EQU:REG:7}
147 \frac{ \partial J_{reg}}{\partial q}(m) = [ \nabla J_{reg}(m), q ] =
148 \int_{\Omega} Y_j \cdot q_j + X_{ki} \cdot q_{k,i} dx
149 \end{equation}
150 where $[.,.]$ is the dual product.
151
152 We also need to provide an approximation of the inverse of the Hessian operator which provides a
153 level set function $h$ for a given value $r$ represented by the pair of values $(Y,X)$. If one ignores the correlation function
154 the inner product defines the Hessian operator of the cost function. In this approach we set
155 \begin{equation}\label{EQU:REG:8}
156 < p, h > = [p, r]
157 \end{equation}
158 for all $p$. This problem can be solved using \escript \class{LinearPDE} class by setting
159 \begin{equation}\label{EQU:REG:8b}
160 \begin{array}{rcl}
161 A_{ij} & =& (\omega_i \cdot L_i^2) \cdot \delta_{ij} \\
162 D & = & \mu^{reg} \cdot \omega
163 \end{array}
164 \end{equation}
165 and $X$ and $Y$ as defined by $r$ for the case of a single-valued level set function.
166 For a vector-valued level-set function one sets:
167 \begin{equation}\label{EQU:REG:8c}
168 \begin{array}{rcl}
169 A_{kilj} & = & (\mu^{reg}_l \omega^{(l)}_i L_i^2) \cdot \delta_{kl} \cdot \delta_{ij} \\
170 D_{kl} & = & \mu^{reg}_l \cdot \omega^{(l)} \delta_{kl} \omega
171 \end{array}
172 \end{equation}
173 ====================================================
174 \begin{classdesc}{Regularization}{domain
175 \optional{, s0=\None}
176 \optional{, s1=\None}
177 \optional{, sc=\None}
178 \optional{, location_of_set_m=Data()}
179 \optional{, numLevelSets=1}
180 \optional{, useDiagonalHessianApproximation=\True}
181 \optional{, tol=1e-8}}
182 opens a linear, steady, second order PDE on the \member{domain}.
183 The parameters \member{numEquations} and \member{numSolutions} give the number
184 of equations and the number of solution components.
185 If \member{numEquations} and \member{numSolutions} are non-positive, then the
186 number of equations and the number of solutions, respectively, stay undefined
187 until a coefficient is defined.
188 \end{classdesc}
189
190
191
192 \section{The general regularization class}
193 \begin{classdesc}{RegularizationBase}{}
194
195 \end{classdesc}

  ViewVC Help
Powered by ViewVC 1.1.26