Methods

be a smooth bounded domain in $$\mathbb{R}^{d},d = 2$$ or 3 and let ν x denote the outward normal to $$\partial \Omega $$ at x. Suppose that the conductivity of $$\Omega $$ is equal to 1. Let D denote a smooth anomaly inside $$\Omega $$ with conductivity 0 < k ≠ 1 < +. The voltage potential in the presence of the set D of conductivity anomalies is denoted by u. It is the solution to the conductivity problem



$$\displaystyle{ \ \left \{\begin{array}{@{}l@{}} \nabla \cdot \left (\chi (\Omega \;\setminus \;\overline{D}) + k\chi (D)\right )\nabla u = 0\quad \mbox{ in }\Omega, \\ \frac{\partial u} {\partial \nu } \bigg\vert _{\partial \Omega } = g\quad \left (g \in L^{2}(\partial \Omega ),\int _{ \partial \Omega }g\,d\sigma = 0\right ), \\ \int _{\partial \Omega }u\,d\sigma = 0, \end{array} \right. }$$

(1)
where χ(D) is the characteristic function of D.

The background voltage potential U in the absence of any anomaly satisfies


$$\displaystyle{ \ \left \{\begin{array}{@{}l@{}} \Delta U = 0\quad \mbox{ in}\Omega, \\ \frac{\partial U} {\partial \nu } \bigg\vert _{\partial \Omega } = g, \\ \int _{\partial \Omega }U\,d\sigma = 0. \end{array} \right. }$$

(2)

Let N(x, z) be the Neumann function for $$-\Delta \ \mbox{ in}\ \Omega $$ corresponding to a Dirac mass at z. That is, N is the solution to


$$\displaystyle{ \ \left \{\begin{array}{@{}l@{}} - \Delta _{x}N(x,z) =\delta _{z}\quad \mbox{ in }\Omega, \\ \frac{\partial N} {\partial \nu _{x}} \Big\vert _{\partial \Omega } = - \frac{1} {\vert \partial \Omega \vert }\;,\int _{\partial \Omega }N(x,z)\,d\sigma (x) = 0\quad \mbox{ for }z \in \Omega. \end{array} \right. }$$

(3)
Note that the Neumann function N(x, z) is defined as a function of $$x \in \overline{\Omega }$$ for each fixed $$z \in \Omega.$$

For B a smooth bounded domain in $$\mathbb{R}^{d}\ \mbox{ and}\ 0 < k\neq 1 < +\infty $$ a conductivity parameter, let $$\hat{v} =\hat{ v}(B,k)$$ be the solution to


$$\displaystyle{ \left \{\begin{array}{@{}l@{}} \Delta \hat{v} = 0\quad \mbox{ in }\mathbb{R}^{d}\;\setminus \;\overline{B}, \\ \Delta \hat{v} = 0\quad \mbox{ in }B, \\ \hat{v}\vert _{-}-\hat{ v}\vert _{+} = 0\quad \mbox{ on }\partial B, \\ k\frac{\partial \hat{v}} {\partial \nu } \bigg\vert _{-}-\frac{\partial \hat{v}} {\partial \nu } \bigg\vert _{+} = 0\quad \mbox{ on }\partial B, \\ \hat{v}(\xi )-\xi \rightarrow 0\quad \mbox{ as }\vert \xi \vert \rightarrow +\infty. \\ \end{array} \right. }$$

(4)
Here, one denotes


$$\displaystyle{v\vert _{\pm }(\xi ):=\lim _{t\rightarrow 0^{+}}v\left (\xi \pm t\nu _{\xi }\right ),\quad \xi \in \partial B,}$$
and


$$\displaystyle{\frac{\partial v} {\partial \nu _{\xi }} \bigg\vert _{\pm }(\xi ):=\lim _{t\rightarrow 0^{+}}\left \langle \nabla v\left (\xi \pm t\nu _{\xi }\right ),\nu _{\xi }\right \rangle,\quad \xi \in \partial B,}$$
if the limits exist, where ν ξ is the outward unit normal to ∂ B at ξ and $$\langle,\rangle$$ is the scalar product in $$\mathbb{R}^{d}$$. For ease of notation, the dot will be sometimes used for the scalar product in $$\mathbb{R}^{d}$$.

Recall that $$\hat{v}$$ plays the role of the first-order corrector in the theory of homogenization [79].



Asymptotic Analysis of the Voltage Perturbations



In this section, an asymptotic expansion of the voltage potentials in the presence of a diametrically small anomaly with conductivity different from the background conductivity is provided.

The following theorem gives asymptotic formulas for both boundary and internal perturbations of the voltage potential that are due to the presence of a conductivity anomaly.


Theorem 1 (Voltage perturbations).

Suppose that D = δB + z,δ being the characteristic size of D, and let u be the solution of (1), where 0 < k ≠ 1 < +∞. Denote by U the background solution, that is, the solution of (2).

(i)

The following asymptotic expansion of the voltage potential on $$\partial \Omega $$ holds for d = 2,3:


$$\displaystyle{ u(x) \approx U(x) -\delta ^{d}\nabla U(z)M(k,B)\nabla _{ z}N(x,z). }$$

(5)
Here, N(x,z) is the Neumann function, that is, the solution to (3), and $$M(k,B) = \left (m_{pq}\right )_{p,q=1}^{d}$$ is the polarization tensor given by


$$\displaystyle{ M(k,B):= (k - 1)\int _{B}\nabla \hat{v}(\xi )\;d\xi, }$$

(6)
where $$\hat{v}$$ is the solution to (1).

 

(ii)

Let w be a smooth harmonic function in $$\Omega $$. The weighted boundary measurements I w [U] satisfies


$$\displaystyle{ I_{w}[U]:=\int _{\partial \Omega }(u - U)(x)\frac{\partial w} {\partial \nu } (x)\,d\sigma (x) \approx -\delta ^{d}\nabla U(z) \cdot M(k,B)\nabla w(z). }$$

(7)

 

(iii)

The following inner asymptotic formula holds:


$$\displaystyle{ u(x) \approx U(z) +\delta \hat{ v}\left (\frac{x - z} {\delta } \right ) \cdot \nabla U(z)\quad \mbox{ for }x\mbox{ near }z. }$$

(8)

 

The inner asymptotic expansion (8) uniquely characterizes the shape and the conductivity of the anomaly. In fact, suppose for two Lipschitz domains B and B and two conductivities k and k that $$\hat{v}(B,k) =\hat{ v}(B^{{\prime}},k^{{\prime}})$$ in a domain englobing B and B , then using the jump conditions satisfied by $$\hat{v}(B,k)$$ and $$\hat{v}(B^{{\prime}},k^{{\prime}})$$, one can easily prove that B = B and k = k .

The asymptotic expansion (5) expresses the fact that the conductivity anomaly can be modeled by a dipole far away from z. It does not hold uniformly in $$\Omega $$. It shows that, from an imaging point of view, the location z and the polarization tensor M of the anomaly are the only quantities that can be determined from boundary measurements of the voltage potential, assuming that the noise level is of order δ d+1. It is then important to precisely characterize the polarization tensor and derive some of its properties, such as symmetry, positivity, and isoperimetric inequalities satisfied by its elements, in order to develop efficient algorithms for reconstructing conductivity anomalies of small volume.

Some important properties of the polarization tensor are listed in the next theorem.


Theorem 2 (Properties of the polarization tensor).

For 0 < k ≠ 1 < +∞, let $$M = M(k,B) = \left (m_{pq}\right )_{p,q=1}^{d}$$ be the polarization tensor associated with the bounded domain B in $$\mathbb{R}^{d}$$ and the conductivity k. Then

(i)

M is symmetric.

 

(ii)

If k > 1, then M is positive definite, and it is negative definite if 0 < k < 1.

 

(iii)

The following isoperimetric inequalities for the polarization tensor


$$\displaystyle{ \left \{\begin{array}{@{}l@{}} \frac{1} {k - 1}\,{\mathrm{trace}}(M) \leq \left (d - 1 + \frac{1} {k}\right )\vert B\vert, \\ (k - 1)\,{\mathrm{trace}}(M^{-1}) \leq \frac{d - 1 + k} {\vert B\vert }, \end{array} \right. }$$

(9)
hold, where trace denotes the trace of a matrix and |B| is the volume of B.

 

The polarization tensor M can be explicitly computed for disks and ellipses in the plane and balls and ellipsoids in three-dimensional space. See [25, pp. 81–89]. The formula of the polarization tensor for ellipses will be useful here. Let B be an ellipse whose semiaxes are on the x 1– and x 2-axes and of lengths a and b, respectively. Then, M(k, B) takes the form


$$\displaystyle{ M(k,B) = (k-1)\vert B\vert \left (\begin{array}{*{10}c} \frac{a + b} {a + kb} & 0 \\ 0 & \frac{a + b} {b + ka} \end{array} \right )\;. }$$

(10)
Formula (5) shows that from boundary measurements, one can always represent and visualize an arbitrary-shaped anomaly by means of an equivalent ellipse of center z with the same polarization tensor. Further, it is impossible to extract the conductivity from the polarization tensor. The information contained in the polarization tensor is a mixture of the conductivity and the volume. A small anomaly with high conductivity and a larger anomaly with lower conductivity can have the same polarization tensor.

The bounds (9) are known as the Hashin–Shtrikman bounds. By making use of these bounds, size and thickness estimations of B can be obtained. An inclusion whose trace of the associated polarization tensor is close to the upper bound must be infinitely thin [40].


Numerical Methods for Anomaly Detection



In this section, one applies the asymptotic formula (5) for the purpose of identifying the location and certain properties of the shape of the conductivity anomalies. Two simple fundamental algorithms that take advantage of the smallness of the anomalies are singled out: projection-type algorithms and multiple signal classification (MUSIC)-type algorithms. These algorithms are fast, stable, and efficient.


Detection of a Single Anomaly: A Projection-Type Algorithm



One briefly discusses a simple algorithm for detecting a single anomaly. The reader can refer to [31, 73] for further details. The projection-type location search algorithm makes use of constant current sources. One wants to apply a special type of current that makes ∇U constant in D. The injection current $$g = a\cdot \nu$$ for a fixed unit vector $$a \in \mathbb{R}^{d}$$ yields $$\nabla U = a$$ in $$\Omega $$.

Assume for the sake of simplicity that d = 2 and D is a disk. Set


$$\displaystyle{w(y) = -(1/2\uppi )\log \vert x - y\vert \quad \mbox{ for }x \in \mathbb{R}^{2}\setminus \overline{\Omega },y \in \Omega \;.}$$
Since w is harmonic in $$\Omega $$, then from (7) to (10), it follows that


$$\displaystyle{ I_{w}[U] \approx \frac{(k - 1)\vert D\vert } {\uppi (k + 1)} \frac{(x - z) \cdot a} {\vert x - z\vert ^{2}} \;,\quad x \in \mathbb{R}^{2}\setminus \overline{\Omega }\;. }$$

(11)

The first step for the reconstruction procedure is to locate the anomaly. The location search algorithm is as follows. Take two observation lines $$\Sigma _{1}$$ and $$\Sigma _{2}$$ contained in $$\mathbb{R}^{2}\setminus \overline{\Omega }$$ given by


$$\displaystyle\begin{array}{rcl} \Sigma _{1}&:=& \mbox{ a line parallel to }a, {}\\ \Sigma _{2}&:=& \mbox{ a line normal to }a. {}\\ \end{array}$$
Find two points $$P_{i} \in \Sigma _{i},i = 1,2,$$ so that


$$\displaystyle\begin{array}{rcl} I_{w}[U](P_{1}) = 0,\quad I_{w}[U](P_{2}) =\max _{x\in \Sigma _{2}}\vert I_{w}[U](x)\vert \;.& & {}\\ \end{array}$$
From (11), one can see that the intersecting point P of the two lines


$$\displaystyle\begin{array}{rcl} \Pi _{1}(P_{1})&:=& \{x\mid a \cdot (x - P_{1}) = 0\},{}\end{array}$$

(12)



$$\displaystyle\begin{array}{rcl} \Pi _{2}(P_{2})&:=& \{x\mid (x - P_{2})\text{ is parallel to }a\}{}\end{array}$$

(13)
is close to the center z of the anomaly D: | Pz | = O(δ 2).

Once one locates the anomaly, the factor | D | (k − 1)∕(k + 1) can be estimated. As it has been said before, this information is a mixture of the conductivity and the volume. A small anomaly with high conductivity and larger anomaly with lower conductivity can have the same polarization tensor.

An arbitrary-shaped anomaly can be represented and visualized by means of an ellipse or an ellipsoid with the same polarization tensor. See Fig. 1.

A183156_2_En_47_Fig1_HTML.gif


Fig. 1
Detection of the location and the polarization tensor of a small arbitrary-shaped anomaly by a projection-type algorithm. The shape of the anomaly is approximated by an ellipse with the same polarization tensor


We refer the reader to [57] for a discussion on the limits of the applicability of the projection-type location search algorithm and the derivation of a second efficient method, called the effective dipole method.


Detection of Multiple Anomalies: A MUSIC-Type Algorithm



Consider m well-separated anomalies $$D_{s} =\delta B_{s} + z_{s}$$ (these are a fixed distance apart), with conductivities k s , $$s = 1,\mathop{\ldots },m$$. Suppose for the sake of simplicity that all the domains B s are disks. Let $$y_{l} \in \mathbb{R}^{2}\setminus \Omega $$ for l = 1, , n denote the source points. Set


$$\displaystyle{U_{y_{l}} = w_{y_{l}}:= -(1/2\uppi )\log \vert x - y_{l}\vert \quad \mbox{ for }x \in \Omega,\quad l = 1,\ldots,n.}$$
The MUSIC-type location search algorithm for detecting multiple anomalies is as follows. For $$n \in \mathbb{N}$$ sufficiently large, define the response matrix $$A = (A_{ll^{{\prime}}})_{l,l^{{\prime}}=1}^{n}$$ by


$$\displaystyle{A_{ll^{{\prime}}} = I_{w_{y_{ l}}}[U_{y_{l^{{\prime}}}}]:=\int _{\partial \Omega }(u - U_{y_{l^{{\prime}}}})(x)\frac{\partial w_{y_{l}}} {\partial \nu } (x)\,d\sigma (x)\;.}$$
Expansion (7) yields


$$\displaystyle{A_{ll^{{\prime}}} \approx -\sum _{s=1}^{m}\frac{2(k_{s} - 1)\vert D_{s}\vert } {k_{s} + 1} \nabla U_{y_{l^{{\prime}}}}(z_{s})\nabla U_{y_{l}}(z_{s})\;.}$$

Introduce


$$\displaystyle{g(x) = \left (U_{y_{1}}(x),\ldots,U_{y_{n}}(x)\right )^{{\ast}}\;,}$$
where v denotes the transpose of the vector v.


Lemma 1 (MUSIC characterization of the range of the response matrix).

There exists n 0 > dm such that for any n > n 0, the following characterization of the location of the anomalies in terms of the range of the matrix A holds:


$$\displaystyle{ g(x) \in \mathop{\mathrm{{\mathrm{Range}}}}\nolimits (A)\mbox{ if and only if }x \in \{ z_{1},\ldots,z_{m}\}\;. }$$

(14)

The MUSIC-type algorithm to determine the location of the anomalies is as follows. Let P noise = IP, where P is the orthogonal projection onto the range of A. Given any point $$x \in \Omega $$, form the vector g(x). The point x coincides with the location of an anomaly if and only if P noise g(x) = 0. Thus, one can form an image of the anomalies by plotting, at each point x, the cost function


$$\displaystyle{W_{{\mathrm{MU}}}(x) = \frac{1} {\vert \vert P_{\text{noise}}g(x)\vert \vert }.}$$
The resulting plot will have large peaks at the locations of the anomalies.

Once one locates the anomalies, the factors $$\vert D_{s}\vert (k_{s} - 1)/(k_{s} + 1)$$ can be estimated from the significant singular values of A.


Bibliography and Open Questions



Part (i) in Theorem 1 was proven in [21, 44, 51]. The proof in [21] is based on a decomposition formula of the solution into a harmonic part and a refraction part first derived in [61]. Part (iii) is from [28]. The Hashin–Shtrikman bounds for the polarization tensor were proved in [43, 77]. The projection algorithm was introduced in [31, 73]. The MUSIC algorithm was originally developed for source separation in signal theory [94]. The MUSIC-type algorithm for locating small conductivity anomalies from the response matrix was first developed in [38]. The strong relation between MUSIC and linear sampling methods was clarified in [16]. The results of this section can be generalized to the detection of anisotropic anomalies [60].

As it has been said before, it is impossible to extract separately from the detected polarization tensor information about the material property and the size of the anomaly. However, if the measurement system is very sensitive, then making use of higher-order polarization tensors yields such information. See [25] for the notion of the higher-order polarization tensors.

One of the most challenging problems in anomaly detection using electrical impedance tomography is that in practical measurements, one usually lacks exact knowledge of the boundary of the background domain. Because of this, the numerical reconstruction from the measured data is done using a model domain that represents the best guess for the true domain. However, it has been noticed that an inaccurate model of the boundary causes severe errors for the reconstructions. An elegant and original solution toward eliminating the error caused by an incorrectly modeled boundary has been proposed and implemented numerically in [69]. As nicely shown in [67], another promising approach is to use multifrequency data. The anomaly can be detected from a weighted frequency difference of the measured boundary voltage perturbations. Moreover, this method eliminates the need for numerically simulated background measurements at the absence of the conductivity anomaly. See [58, 67].



3 Ultrasound Imaging for Anomaly Detection




Physical Principles


Ultrasound imaging is a noninvasive, easily portable, and relatively inexpensive diagnostic modality which finds extensive clinical use. The major applications of ultrasound include many aspects of obstetrics and gynecology involving the assessment of fetal health, intra-abdominal imaging of the liver and kidney, and the detection of compromised blood flow in veins and arteries.

Operating typically at frequencies between 1 and 10 MHz, ultrasound imaging produces images via the backscattering of mechanical energy from interfaces between tissues and small structures within tissue. It has high spatial resolution, particularly at high frequencies, and involves no ionizing radiation. The weaknesses of the technique include the relatively poor soft tissue contrast and the fact that gas and bone impede the passage of ultrasound waves, meaning that certain organs cannot easily be imaged. However, ultrasound imaging is a valuable technique for anomaly detection. It can be done in the time domain and the frequency domain.

Mathematical models for acoustical soundings of biological media involve the Helmholtz equation in the frequency domain and the scalar wave equation in the time domain.


Asymptotic Formulas in the Frequency Domain



Let k and ρ be positive constants. With the notation of section “Asymptotic Analysis of the Voltage Perturbations,” ρ is the compressibility of the anomaly D and k is its volumetric mass density. The scalar acoustic pressure u generated by the Neumann data g in the presence of the anomaly D is the solution to the Helmholtz equation:


$$\displaystyle{ \left \{\begin{array}{@{}ll@{}} \nabla \cdot \left (\chi \left (\Omega \setminus \overline{D}\right ) + k\chi (D)\right )\nabla u +\omega ^{2}(\chi \left (\Omega \setminus \overline{D}\right ) +\rho \chi (D))u = 0\quad \mbox{ in }\Omega, \\ \frac{\partial u} {\partial \nu } = g\quad \mbox{ on }\partial \Omega, \end{array} \right. }$$

(15)
while the background solution U satisfies


$$\displaystyle{ \left \{\begin{array}{@{}ll@{}} \Delta U +\omega ^{2}U = 0\quad \mbox{ in }\Omega, \\ \frac{\partial U} {\partial \nu } = g\quad \mbox{ on }\partial \Omega. \end{array} \right. }$$

(16)

Here, ω is the operating frequency. A relevant boundary data g is the normal derivative of a plane wave $$e^{i\omega x\cdot \theta }$$, with the wavelength λ: = 2πω, traveling in the direction of the unit vector θ.

Introduce the Neumann function for $$-\left (\Delta +\omega ^{2}\right )$$ in $$\Omega $$ corresponding to a Dirac mass at z. That is, N ω is the solution to


$$\displaystyle{ \ \left \{\begin{array}{@{}ll@{}} -\left (\Delta _{x} +\omega ^{2}\right )N_{\omega }(x,z) =\delta _{ z}\quad \mbox{ in }\Omega, \\ \frac{\partial N} {\partial \nu _{x}} \left \vert _{\partial \Omega } = 0\right.\quad \mbox{ on }\partial \Omega. \end{array} \right. }$$

(17)

Assuming that ω 2 is not an eigenvalue for the operator $$-\Delta \ \mbox{ in}\ L^{2}(\Omega )$$ with homogeneous Neumann boundary conditions, one can prove, using the theory of relatively compact operators, existence and uniqueness of a solution to (15) at least for δ small enough [95]. Moreover, the following asymptotic formula holds.


Theorem 3 (Pressure perturbations).

Let u be the solution of (15) and let U be the background solution. Suppose that $$D =\delta B + z,\ \mbox{ with}\ 0 < (k,\rho )\neq (1,1) < +\infty $$. Suppose that ωδ ≪ 1.

(i)

For any $$x \in \partial \Omega $$,


$$\displaystyle\begin{array}{rcl} u(x)& \approx & U(x) -\delta ^{d}\left (\nabla U(z) \cdot M(k,B)\nabla _{ z}N_{\omega }(x,z)\right. \\ & & \left.+\omega ^{2}(\rho -1)\vert B\vert U(z)N_{\omega }(x,z)\right ), {}\end{array}$$

(18)
where M(k,B) is the polarization tensor associated with B and k.

 

(ii)

Let w be a smooth function such that $$\left (\Delta +\omega ^{2}\right )w = 0$$ in $$\Omega $$. The weighted boundary measurements I w [U,ω] satisfy


$$\displaystyle\begin{array}{rcl} I_{w}[U,\omega ]&:=& \int _{\partial \Omega }(u - U)(x)\frac{\partial w} {\partial \nu } (x)\,d\sigma (x) \\ & \approx &-\delta ^{d}\left (\nabla U(z) \cdot M(k,B)\nabla w(z) +\omega ^{2}(\rho -1)\vert B\vert U(z)w(z)\right ). {}\end{array}$$

(19)

 

(iii)

The following inner asymptotic formula holds:


$$\displaystyle{ u(x) \approx U(z) +\delta \hat{ v}\left (\frac{x - z} {\delta } \right ) \cdot \nabla U(z)\quad \mbox{ for }x\mbox{ near }z, }$$

(20)
where $$\hat{v}$$ is the solution to (1).

 

Compared to the conductivity equation, the only extra difficulty in establishing asymptotic formulas for the Helmholtz equation (15) as the size of the acoustic anomaly goes to zero is that the equations inside and outside the anomaly are not the same.


Asymptotic Formulas in the Time Domain



Suppose that ρ = 1. Consider the initial boundary value problem for the (scalar) wave equation


$$\displaystyle{ \left \{\begin{array}{ll} &\partial _{t}^{2}u -\nabla \cdot \left (\chi \left (\Omega \setminus \overline{D}\right ) + k\chi (D)\right )\nabla u = 0\quad \text{in} \Omega _{T}, \\ &u(x,0) = u_{0}(x),\quad \partial _{t}u(x,0) = u_{1}(x)\quad \text{for}\ x \in \Omega, \\ &\frac{\partial u} {\partial \nu } = g\quad \text{on} \partial \Omega _{T},\end{array} \right. }$$

(21)
where T < + is a final observation time, $$\Omega _{T} = \Omega \times ]0,T[,\mbox{ and}\ \partial \Omega _{T} = \partial \Omega \times ]0,T[$$. The initial data $$u_{0},u_{1} \in \mathcal{C}^{\infty }\left (\overline{\Omega }\right )$$ and the Neumann boundary data $$g \in \mathcal{C}^{\infty }(0,T;\mathcal{C}^{\infty }(\partial \Omega ))$$ are subject to compatibility conditions.

Define the background solution U to be the solution of the wave equation in the absence of any anomalies. Thus, U satisfies


$$\displaystyle{ \left \{\begin{array}{ll} &\partial _{t}^{2}U - \Delta U = 0\quad \text{in} \Omega _{T}, \\ &U(x,0) = u_{0}(x),\quad \partial _{t}U(x,0) = u_{1}(x)\quad \text{for}\ x \in \Omega, \\ &\frac{\partial U} {\partial \nu } = g\quad \text{on}\ \partial \Omega _{T}.\end{array} \right. }$$

For ρ > 0, define the operator P ρ on tempered distributions by


$$\displaystyle{ P_{\rho }[\psi ](x,t) =\int _{\vert \omega \vert \leq \rho }e^{-\sqrt{-1}\omega t}\hat{\psi }(x,\omega )\;d\omega, }$$

(22)
where $$\hat{\psi }(x,\omega )$$ denotes the Fourier transform of ψ(x, t) in the t-variable. Clearly, the operator P ρ truncates the high-frequency component of ψ.

The following asymptotic expansion holds as $$\delta \rightarrow 0$$.


Theorem 4 (Perturbations of weighted boundary measurements).

Let $$w \in C^{\infty }\left (\overline{\Omega }_{T}\right )$$ satisfy $$\left (\partial _{t}^{2} - \Delta \right )w(x,t) = 0$$ in $$\Omega _{T}$$ with ∂ t w(x,T) = w(x,T) = 0 for $$x \in \Omega $$. Suppose that ρ ≪ 1∕δ. Define the weighted boundary measurements


$$\displaystyle{I_{w}[U,T]:=\int _{\partial \Omega _{T}}P_{\rho }[u - U](x,t)\frac{\partial w} {\partial \nu } (x,t)\,d\sigma (x)\,dt.}$$
Then, for any fixed $$T> \mbox{ diam}(\Omega )$$” src=”/wp-content/uploads/2016/04/A183156_2_En_47_Chapter_IEq59.gif”></SPAN>, <SPAN class=EmphasisTypeItalic>the following asymptotic expansion for I</SPAN> <SUB>w</SUB> <SPAN class=EmphasisTypeItalic>[U,T] holds as δ → 0:</SPAN><br />
<DIV id=Equ25 class=Equation><br />
<DIV class=EquationContent><br />
<DIV class=MediaObject><IMG alt=

(23)
where M(k,B) is defined by (6).

Expansion (23) is a weighted expansion. Pointwise expansions similar to those in Theorem 1 which is for the steady-state model can also be obtained.

Let $$y \in \mathbb{R}^{3}$$ be such that | yz | ≫ δ. Choose


$$\displaystyle{ U(x,t):= U_{y}(x,t):= \frac{\delta _{t=\vert x-y\vert }} {4\pi \vert x - y\vert }\quad \mbox{ for }x\neq y. }$$

(24)
It is easy to check that U y is the outgoing Green function to the wave equation:


$$\displaystyle{\left (\partial _{t}^{2} - \Delta \right )U_{ y}(x,t) =\delta _{x=y}\delta _{t=0}\quad \mbox{ in }\mathbb{R}^{3}\times ]0,+\infty [\;.}$$
Moreover, U y satisfies the initial conditions: $$U_{y}(x,0) = \partial _{t}U_{y}(x,0) = 0$$ for $$x\neq y$$. Consider now for the sake of simplicity the wave equation in the whole three-dimensional space with appropriate initial conditions:


$$\displaystyle{ \left \{\begin{array}{ll} &\partial _{t}^{2}u -\nabla \cdot \left (\chi \left (\mathbb{R}^{3}\setminus \overline{D}\right ) + k\chi (D)\right )\nabla u =\delta _{ x=y}\delta _{t=0}\quad \text{in}\ \mathbb{R}^{3}\times ]0,+\infty [\;, \\ &u(x,0) = 0,\quad \partial _{t}u(x,0) = 0\quad \text{for}\ x \in \mathbb{R}^{3},x\neq y.\end{array} \right. }$$

(25)

The following theorem holds.


Theorem 5 (Pointwise perturbations).

Let u be the solution to (25). Set U y to be the background solution. Suppose that ρ ≪ 1∕δ.

(i)

The following outer expansion holds


$$\displaystyle{ P_{\rho }[u - U_{y}](x,t) \approx -\delta ^{3}\int _{ \mathbb{R}}\nabla P_{\rho }[U_{z}](x,t-\tau ) \cdot M(k,B)\nabla P_{\rho }[U_{y}](z,\tau )\;d\tau \;, }$$

(26)
for x away from z, where M(k,B) is defined by (6) and U y and U z by (24).

 

(ii)

The following inner approximation holds:


$$\displaystyle{ P_{\rho }[u - U_{y}](x,t) \approx \delta \hat{ v}\left (\frac{x - z} {\delta } \right ) \cdot \nabla P_{\rho }[U_{y}](x,t)\quad \mbox{ for }x\mbox{ near }z, }$$

(27)
where $$\hat{v}$$ is given by (4) and U y by (24).

 

Formula (26) shows that the perturbation due to the anomaly is in the time domain a wave front emitted by a dipolar source located at point z.

Taking the Fourier transform of (26) in the time variable yields the expansions given in Theorem 3 for the perturbations resulting from the presence of a small anomaly for solutions to the Helmholtz equation at low frequencies (at large wavelengths compared to the size of the anomaly).


Numerical Methods




MUSIC-Type Imaging at a Single Frequency



Consider m well-separated anomalies $$D_{s} =\delta B_{s} + z_{s},s = 1,\mathop{\ldots },m$$. The compressibility and volumetric mass density of D s are denoted by ρ s  and k s , respectively. Suppose as before that all the domains B s are disks. Let $$(\theta _{1},\ldots,\theta _{n})$$ be n unit vectors in $$\mathbb{R}^{d}$$. For arbitrary $$\theta \in \{\theta _{1},\ldots,\theta _{n}\},$$ one assumes that one is in the possession of the boundary data u when the object $$\Omega $$ is illuminated with the plane wave $$U(x) = e^{i\omega \theta \cdot x}$$. Therefore, taking $$w(x) = e^{-i\omega \theta ^{{\prime}}\cdot x }$$ for $$\theta ^{{\prime}}\in \{\theta _{1},\ldots,\theta _{n}\}$$ shows that one is in possession of


$$\displaystyle{\sum _{s=1}^{m}\vert D_{ s}\vert \left (2\frac{(k_{s} - 1)} {k_{s} + 1} \theta \cdot \theta ^{{\prime}} + (\rho _{ s} - 1)\right )e^{i\omega (\theta -\theta ^{{\prime}})\cdot z_{ s}},}$$
for $$\theta,\theta ^{{\prime}}\in \{\theta _{1},\ldots,\theta _{n}\}.$$ Define the response matrix $$A = (A_{ll^{{\prime}}})_{l,l^{{\prime}}=1}^{n} \in \mathbb{C}^{n\times n}$$ by


$$\displaystyle{A_{ll^{{\prime}}} =\sum _{ s=1}^{m}\vert D_{ s}\vert \left (2\frac{(k_{s} - 1)} {k_{s} + 1} \theta _{l} \cdot \theta _{l^{{\prime}}} + (\rho _{s} - 1))\right )e^{i\omega (\theta _{l}-\theta _{l^{{\prime}}})\cdot z_{s} },\quad l,l^{{\prime}} = 1,\ldots,n\;.}$$
Introduce


$$\displaystyle{g(x) = \left ((1,\theta _{1})^{{\ast}}e^{i\omega \theta _{1}\cdot x},\ldots,(1,\theta _{ n})^{{\ast}}e^{i\omega \theta _{n}\cdot x}\right )^{{\ast}}.}$$

Analogously to Lemma 1, one has the following characterization of the location of the anomalies in terms of the range of the matrix A.


Lemma 2 (MUSIC characterization of the range of the response matrix).

There exists $$n_{0} \in \mathbb{N},n_{0}> (d + 1)m,$$” src=”/wp-content/uploads/2016/04/A183156_2_En_47_Chapter_IEq74.gif”></SPAN> <SPAN class=EmphasisTypeItalic>such that for any n ≥ n</SPAN> <SUB>0</SUB>, <SPAN class=EmphasisTypeItalic>the following statement holds:</SPAN><br />
<DIV id=Equo class=Equation><br />
<DIV class=EquationContent><br />
<DIV class=MediaObject><IMG alt=
where g j (x) is the jth column of g(x).

The MUSIC algorithm can now be used as before to determine the location of the anomalies. Let P noise = IP, where P is the orthogonal projection onto the range of A. The imaging functional


$$\displaystyle{W_{{\mathrm{MU}}}(x):= \frac{1} {\sum _{j=1}^{d+1}\vert \vert P_{\text{noise}}g^{j}(x)\vert \vert }}$$
has large peaks only at the locations of the anomalies. See Fig. 2.

A183156_2_En_47_Fig2_HTML.gif


Fig. 2
MUSIC-type reconstruction from the singular value decomposition of A represented in Fig. 3

The significant singular vectors of A can be computed by the singular value decomposition. The number of significant singular values determines the number of anomalies. If, for example, k s ≠ 1 and ρ s ≠ 1 for all s = 1, , m, then there are exactly (d + 1)m significant singular values of A and the rest are zero or close to zero. See Fig. 3. The significant singular values of A can be used to estimate $$\frac{(k_{s}-1)} {k_{s}+1} \vert D_{s}\vert $$ and (ρ s − 1) | D s |.

A183156_2_En_47_Fig3_HTML.gif


Fig. 3
Singular value decomposition of the response matrix A corresponding to two well-separated anomalies of general shape for n = 20, using a standard log scale. Six singular values emerge from the 14 others in the noise subspace



Backpropagation-Type Imaging at a Single Frequency



A backpropagation imaging functional at a single frequency ω is given by


$$\displaystyle{W_{{\mathrm{BP}}}(x):= \frac{1} {n}\sum _{l=1}^{n}e^{-2i\omega \theta _{l}\cdot x}I_{ w_{l}}[U_{l}],}$$
where $$U_{l}(x) = w_{l}(x) = e^{i\omega \theta _{l}\cdot x}$$ for $$\theta _{l} \in \{\theta _{1},\ldots,\theta _{n}\}\;.$$ Suppose that $$(\theta _{1},\ldots,\theta _{n})$$ are equidistant points on the unit sphere S d−1. For sufficiently large n, since


$$\displaystyle{ \frac{1} {n}\sum _{l=1}^{n}e^{i2\omega \theta _{l}\cdot x} \approx \left \{\begin{array}{@{}l@{}} j_{0}(2\omega \vert x\vert )\quad \mbox{ for }d = 3, \\ J_{0}(2\omega \vert x\vert )\quad \mbox{ for }d = 2, \end{array} \right.}$$
where j 0 is the spherical Bessel function of order zero and J 0 is the Bessel function of the first kind and of order zero, it follows that


$$\displaystyle{W_{{\mathrm{BP}}}(x) \approx \sum _{s=1}^{m}\vert D_{ s}\vert \left (2\frac{(k_{s} - 1)} {k_{s} + 1} + (\rho _{s} - 1))\right )\times \left \{\begin{array}{@{}l@{}} j_{0}(2\omega \vert x - z_{s}\vert )\quad \mbox{ for }d = 3, \\ J_{0}(2\omega \vert x - z_{s}\vert )\quad \mbox{ for }d = 2. \end{array} \right.}$$

An analogy between the backpropagation and MUSIC-type imaging can be established. Suppose that k s = 1 for s = 1, , m. One can see that


$$\displaystyle{W_{{\mathrm{MU}}}(x) \propto \frac{1} {\vert D_{s}\vert (\rho _{s} - 1) - W_{{\mathrm{BP}}}(x)}}$$
for x near z s [9].


Kirchhoff-Type Imaging Using a Broad Range of Frequencies



Let $$y_{l} \in \mathbb{R}^{2}\setminus \Omega $$ for l = 1, , n denote an array of source points. Set


$$\displaystyle{w_{y_{l}}(x) = \frac{i} {4}H_{0}^{(1)}(\omega \vert x - y_{ l}\vert )\quad \mbox{ and}\quad U_{y_{l^{{\prime}}}}(x) = \frac{i} {4}H_{0}^{(1)}(\omega \vert x - y_{ l^{{\prime}}}\vert ),}$$
where H 0 (1) is the Hankel function of first kind and order zero. Using the asymptotic form of the Hankel function, one finds that for ω | xy | ≫ 1,


$$\displaystyle{ \frac{i} {4}H_{0}^{(1)}(\omega \vert x - y\vert ) \approx \frac{1} {2\sqrt{2\pi }} \frac{e^{i\,\pi /4}} {\sqrt{\omega \vert x - y\vert }}e^{i\omega \vert x-y\vert },}$$
and


$$\displaystyle{ \frac{i} {4}\nabla H_{0}^{(1)}(\omega \vert x - y\vert ) \approx \frac{1} {2\sqrt{2\pi }}\left (\frac{i\omega (x - y)} {\vert x - y\vert } \right ) \frac{e^{i\,\pi /4}} {\sqrt{\omega \vert x - y\vert }}e^{i\omega \vert x-y\vert }.}$$

Assume a high-frequency regime with ω L ≫ 1 for L the distance from the array center point to the locations z s , s = 1, , m. It follows that


$$\displaystyle\begin{array}{rcl} & & I_{w_{l}}[U_{l^{{\prime}}},\omega ] \propto {}\\ & &\quad \sum _{s=1}^{m}\vert D_{ s}\vert \left (-2\frac{k_{s} - 1} {k_{s} + 1} \frac{(z_{s} - y_{l}) \cdot (z_{s} - y_{l^{{\prime}}})} {\vert z_{s} - y_{l}\vert \vert z_{s} - y_{l^{{\prime}}}\vert } + (\rho _{s} - 1)\right )e^{i\omega (\vert z_{s}-y_{l}\vert +\vert z_{s}-y_{l^{{\prime}}}\vert )}. {}\\ \end{array}$$
Introduce the response matrix $$A(\omega ) = (A_{ll^{{\prime}}}(\omega ))$$ by


$$\displaystyle{A_{ll^{{\prime}}}(\omega ):= I_{w_{l}}[U_{l^{{\prime}}},\omega ]}$$
and the illumination vector


$$\displaystyle{g(x,\omega ):= \left (\left (1, \frac{x - y_{1}} {\vert x - y_{1}\vert }\right )^{{\ast}}e^{i\omega \vert x-y_{1}\vert },\ldots,\left (1, \frac{x - y_{n}} {\vert x - y_{n}\vert }\right )^{{\ast}}e^{i\omega \vert x-y_{n}\vert }\right )^{{\ast}}.}$$
In the case of measurements at multiple frequencies (ω j ), we construct the weighted Kirchhoff imaging functional as


$$\displaystyle{W_{{\mathrm{KI}}}(x) = \frac{1} {J}\sum _{\omega _{j},j=1,\ldots,J}\sum _{l}\left (g(x,\omega _{j}),u_{l}(\omega _{j})\right )\left (g(x,\omega _{j}),\overline{v_{l}}(\omega _{j})\right ),}$$
where $$(a,b) = \overline{a} \cdot b$$, J is the number of frequencies, and u l and v l are, respectively, the left and right singular vectors of A. As for W MU, W KI is written in terms of the singular-value decompositions of the response matrices A(ω j ).


Time-Reversal Imaging



Unlike the three previous imaging methods, the one in this section is in time domain. It is based on time reversal.

The main idea of time reversal is to take advantage of the reversibility of the wave equation in a non-dissipative unknown medium in order to back propagate signals to the sources that emitted them. In the context of anomaly detection, one measures the perturbation of the wave on a closed surface surrounding the anomaly and retransmits it through the background medium in a time-reversed chronology. Then the perturbation will travel back to the location of the anomaly. One can show that the time-reversal perturbation focuses on the location z of the anomaly with a focal spot size limited to one-half the wavelength which is in agreement with the Rayleigh resolution limit.

In mathematical terms, suppose that one is able to measure the perturbation uU y and its normal derivative at any point x on a sphere S englobing the anomaly D and for a large time t 0. The time-reversal operation is described by the transform $$t\mapsto t_{0} - t$$. Both the perturbation and its normal derivative on S are time reversed and emitted from S. Then a time-reversed perturbation propagates inside the volume surrounded by S.

To detect the anomaly from measurements of the wavefield uU y away from the anomaly, one can use a time-reversal technique. Taking into account the definition of the outgoing fundamental solution (24) to the wave equation, spatial reciprocity, and time-reversal invariance of the wave equation, one defines the time-reversal imaging functional W TR by


$$\displaystyle\begin{array}{rcl} W_{{\mathrm{TR}}}(x,t)& =& \int _{\mathbb{R}}\int _{S}\left [U_{x}\left (x^{{\prime}},t - s\right )\frac{\partial P_{\rho }\left [u - U_{y}\right ]} {\partial \nu } \left (x^{{\prime}},t_{ 0} - s\right )\right. \\ & & \left.-\frac{\partial U_{x}} {\partial \nu } \left (x^{{\prime}},t - s\right )P_{\rho }\left [u - U_{ y}\right ]\left (x^{{\prime}},t_{ 0} - s\right )\right ]\,d\sigma \left (x^{{\prime}}\right )\,ds,{}\end{array}$$

(28)
where


$$\displaystyle{U_{x}\left (x^{{\prime}},t-\tau \right ) = \frac{\delta \left (t -\tau -\vert x - x^{{\prime}}\vert \right )} {4\pi \vert x - x^{{\prime}}\vert }.}$$
The imaging functional W TR corresponds to propagating inside the volume surrounded by S, the time-reversed perturbation P ρ [uU y ], and its normal derivative on S. Theorem 5 shows that


$$\displaystyle{P_{\rho }\left [u - U_{y}\right ](x,t) \approx -\delta ^{3}\int _{ \mathbb{R}}\nabla P_{\rho }[U_{z}](x,t-\tau ) \cdot m(z,\tau )\;d\tau,}$$
where


$$\displaystyle{ m(z,\tau ) = M(k,B)\nabla P_{\rho }[U_{y}](z,\tau )\;. }$$

(29)
Therefore, since


$$\displaystyle{ \begin{array}{l} \int _{\mathbb{R}}\int _{S}\left [U_{x}(x^{{\prime}},t - s)\frac{\partial P_{\rho }[U_{z}]} {\partial \nu } (x^{{\prime}},t_{ 0} - s-\tau )\right. \\ \left.\qquad \quad -\frac{\partial U_{x}} {\partial \nu } (x^{{\prime}},t - s)P_{\rho }[U_{ z}](x^{{\prime}},t_{ 0} - s-\tau )\right ]\,d\sigma (x^{{\prime}})\,ds \\ = P_{\rho }[U_{z}](x,t_{0} -\tau -t) - P_{\rho }[U_{z}](x,t - t_{0}+\tau ), \end{array} }$$

(30)
one obtains the approximation


$$\displaystyle\begin{array}{rcl} W_{{\mathrm{TR}}}(x,t) \approx -\delta ^{3}\int _{ \mathbb{R}}m(z,\tau ) \cdot \nabla _{z}\left [P_{\rho }[U_{z}](x,t_{0} -\tau -t) - P_{\rho }[U_{z}](x,t - t_{0}+\tau )\right ]d\tau,& & {}\\ \end{array}$$
which can be interpreted as the superposition of incoming and outgoing waves, centered on the location z of the anomaly. Since


$$\displaystyle{P_{\rho }[U_{y}](x,\tau ) = \frac{\sin \rho (\tau -\vert x - y\vert )} {2\pi (\tau -\vert x - y\vert )\vert x - y\vert },}$$
m(z, τ) is concentrated at the travel time τ = T = | zy |. It then follows that


$$\displaystyle{ W_{{\mathrm{TR}}}(x,t) \approx -\delta ^{3}m(z,T) \cdot \nabla _{ z}\left [P_{\rho }[U_{z}](x,t_{0} - T - t) - P_{\rho }[U_{z}](x,t - t_{0} + T)\right ]. }$$

(31)
The imaging functional W TR is clearly the sum of incoming and outgoing polarized spherical waves.

Approximation (31) has an important physical interpretation. By changing the origin of time, T can be set to 0 without loss of generality. Then by taking a Fourier transform of (31) over the time variable t, one obtains that


$$\displaystyle{\hat{W}_{{\mathrm{TR}}}(x,\omega ) \propto \delta ^{3}m(z,T) \cdot \nabla j_{ 0}(\omega \vert x - z\vert ),}$$
where ω is the wave number. This shows that the time-reversal perturbation W TR focuses on the location z of the anomaly with a focal spot size limited to one-half the wavelength.

An identity parallel to (30) can be derived in the frequency domain. In fact, one has


$$\displaystyle{ \begin{array}{ll} \int _{S}\left [\hat{U}_{x}\left (x^{{\prime}}\right )\frac{\partial \overline{\widehat{U_{z}}}} {\partial \nu } \left (x^{{\prime}}\right ) -\frac{\partial \hat{U}_{x}} {\partial \nu } \left (x^{{\prime}}\right )\overline{\widehat{U_{ z}}}\left (x^{{\prime}}\right )\right ]\,d\sigma \left (x^{{\prime}}\right )& = 2i\mathfrak{I}m\,\widehat{U_{ z}}(x) \\ & \propto j_{0}(\omega \vert x - z\vert ), \end{array} }$$

(32)
which shows that in the frequency domain, W TR coincides with W BP.


Bibliography and Open Questions



The initial boundary-value problems for the wave equation in the presence of anomalies of small volume have been considered in [5, 24]. Theorem 5 is from [12]. In [12], a time-reversal approach was also designed for locating the anomaly from the outer expansion (26). The physics literature on time reversal is quite rich. One refers, for instance, to [48] and the references therein. See [93] for clinical applications of time reversal. Many interesting mathematical works have dealt with different aspects of time-reversal phenomena: see, for instance, [33] for time reversal in the time domain, [4547, 82] for time reversal in the frequency domain, and [37, 50] for time reversal in random media.

The MUSIC-type algorithm for locating small acoustic or electromagnetic anomalies from the multi-static response matrix at a fixed frequency was developed in [18]. See also [1820], where a variety of numerical results was presented to highlight its potential and its limitation. It is worth mentioning that the MUSIC-type algorithm is related to time reversal [82, 87].

MUSIC and Kirchhoff imaging functionals can be extended to the time domain in order to detect the anomaly and its polarization tensor from (dynamical) boundary measurements [24].

The inner expansion in Theorem 5 can be used to design an efficient optimization algorithm for reconstructing the shape and the physical parameter of an anomaly from the near-field perturbations of the wavefield, which can be used in radiation force imaging.

In radiation force imaging, one uses the acoustic radiation force of an ultrasonic focused beam to remotely generate mechanical vibrations in organs. A spatiotemporal sequence of the propagation of the induced transient wave can be acquired, leading to a quantitative estimation of the physical parameters of the anomaly. See, for instance, [35, 36].

The proposed location search algorithms using transient wave or broad-range multifrequency boundary measurements can be extended to the case with limited-view measurements. Using the geometrical control method [34], one can still exploit those algorithms and perform imaging with essentially the same resolution using partial data as using complete data, provided that the geometric optics condition holds.

An identity similar to (32) can be derived in an inhomogeneous medium, which shows that the sharper the behavior of the imaginary part of the Green function around the location of the anomaly is, the higher is the resolution. It would be quite challenging to explicitly see how this behavior depends on the heterogeneity of the surrounding medium. This would yield super-resolved ultrasound imaging systems.


4 Infrared Thermal Imaging



Physical Principles



Infrared thermal imaging is becoming a common screening modality in the area of breast cancer. By carefully examining the aspects of temperature and blood vessels of the breasts in thermal images, signs of possible cancer or precancerous cell growth may be detected up to 10 years prior to being discovered using any other procedure. This provides the earliest detection of cancer possible.

Because of thermal imaging’s extreme sensitivity, these temperature variations and vascular changes may be among the earliest signs of breast cancer and/or a precancerous state of the breast. An abnormal infrared image of the breast is an important marker of high risk for developing breast cancer. See [3, 83].


Asymptotic Analysis of Temperature Perturbations



Suppose that the background $$\Omega $$ is homogeneous with thermal conductivity 1 and that the anomaly D = δ B + z has thermal conductivity 0 < k ≠ 1 < +. In this section, one considers the following transmission problem for the heat equation:


$$\displaystyle{ \left \{\begin{array}{ll} &\partial _{t}u -\nabla \cdot \left (\chi \left (\Omega \setminus \overline{D}\right ) + k\chi (D)\right )\nabla u = 0\quad \text{in}\ \Omega _{T}, \\ &u(x,0) = u_{0}(x)\quad \text{for}\ x \in \Omega, \\ &\frac{\partial u} {\partial \nu } = g\quad \text{on}\ \partial \Omega _{T},\end{array} \right. }$$

(33)
where the Neumann boundary data g and the initial data u 0 are subject to a compatibility condition. Let U be the background solution defined as the solution of


$$\displaystyle{ \left \{\begin{array}{ll} &\partial _{t}U - \Delta U = 0\quad \text{in}\ \Omega _{T}, \\ &U(x,0) = u_{0}(x)\quad \text{for}\ x \in \Omega, \\ &\frac{\partial U} {\partial \nu } = g\quad \text{on}\ \partial \Omega _{T}.\end{array} \right. }$$

The following asymptotic expansion holds as $$\delta \rightarrow 0$$.


Theorem 6 (Perturbations of weighted boundary measurements).

Let $$w \in C^{\infty }\left (\overline{\Omega }_{T}\right )$$ be a solution to the adjoint problem, namely, satisfy $$(\partial _{t} + \Delta )w(x,t) = 0$$ in $$\Omega _{T}$$ with w(x,T) = 0 for $$x \in \Omega $$. Define the weighted boundary measurements


$$\displaystyle{I_{w}[U,T]:=\int _{\partial \Omega _{T}}(u - U)(x,t)\frac{\partial w} {\partial \nu } (x,t)\,d\sigma (x)\,dt.}$$
Then, for any fixed T > 0, the following asymptotic expansion for I w [U,T] holds as δ → 0:


$$\displaystyle{ I_{w}[U,T] \approx -\delta ^{d}\int _{ 0}^{T}\nabla U(z,t) \cdot M(k,B)\nabla w(z,t)\,dt, }$$

(34)
where M(k,B) is defined by (6).

Note that (34) holds for any fixed positive final time T, while (23) holds only for $$T> \mbox{ diam}(\Omega ).$$” src=”/wp-content/uploads/2016/04/A183156_2_En_47_Chapter_IEq89.gif”></SPAN> This difference comes from the finite speed propagation property for the wave equation compared to the infinite one for the heat equation.</DIV><br />
<DIV class=Para>Consider now the background solution to be the Green function of the heat equation at <SPAN class=EmphasisTypeItalic>y</SPAN>:<br />
<DIV id=Equ39 class=Equation><br />
<DIV class=EquationContent><br />
<DIV class=MediaObject><IMG alt= 0, \\ 0\quad \mbox{ for }t < 0\;. \end{array} \right. }$$" src="/wp-content/uploads/2016/04/A183156_2_En_47_Chapter_Equ39.gif">

(35)
Let u be the solution to the following heat equation with an appropriate initial condition:


$$\displaystyle{ \left \{\begin{array}{ll} &\partial _{t}u -\nabla \cdot \left (\chi \left (\mathbb{R}^{d}\setminus \overline{D}\right ) + k\chi (D)\right )\nabla u = 0\quad \text{in}\ \mathbb{R}^{d}\times ]0,+\infty [\;, \\ &u(x,0) = U_{y}(x,0)\quad \text{for}\ x \in \mathbb{R}^{d}.\end{array} \right. }$$

(36)

Proceeding as in the derivation of (26), one can prove that δ u(x, t): = uU is approximated by


$$\displaystyle{ -(k - 1)\int _{0}^{t} \frac{1} {(4\pi (t-\tau ))^{d/2}}\int _{\partial D}e^{-\frac{\vert x-x^{{\prime}}\vert ^{2}} {4(t-\tau )} } \frac{\partial \hat{v}} {\partial \nu } \bigg\vert _{-}\left (\frac{x^{{\prime}}- z} {\delta } \right ) \cdot \nabla U_{y}\left (x^{{\prime}},\tau \right )\,d\sigma \left (x^{{\prime}}\right )\,d\tau, }$$

(37)
for x near z. Therefore, analogously to Theorem 5, the following pointwise expansion follows from the approximation (37).


Theorem 7 (Pointwise perturbations).

Let $$y \in \mathbb{R}^{d}$$ be such that |y − z|≫δ. Let u be the solution to (36). The following expansion holds

Only gold members can continue reading. Log In or Register to continue

Stay updated, free articles. Join our Telegram channel

Apr 9, 2016 | Posted by in GENERAL RADIOLOGY | Comments Off on Methods

Full access? Get Clinical Tree

Get Clinical Tree app for offline access