Extended sampling method in inverse scattering

Let $u^{B_z}_{\infty}(\hat{x}; d), x\in \mathbb{S}$ be the far-field pattern of the sound-soft disc B_z centered at z with radius R for plane incident waves of all directions $d\in \mathbb{S}$. Define a far field operator $\mathcal{F}_z: L^2(\mathbb{S}) \to L^2(\mathbb{S})$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \label{FO} \mathcal{F}_z g(\hat{x}) = \int_{\mathbb{S}}u^{B_z}_{\infty}(\hat{x},d)g(d){\rm d} s(d), \quad \hat{x} \in \mathbb{S}. \nonumber \end{align} \tag{ 3.1 }$$

Let U^s(x) and $U_\infty(\hat{x})$ be the scattered field and far field pattern of the scatterer D due to one incident wave, respectively. Using the far field operator $\mathcal{F}_z$, we set up a far field equation

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \label{fe} \left(\mathcal{F}_z g\right)(\hat{x})=U_\infty(\hat{x}),\ \ \ \ \hat{x}\in \mathbb{S}. \nonumber \end{align} \tag{ 3.2 }$$

This integral equation is the main ingredient of the ESM. The advantage of using $\mathcal{F}_z$ is that it can be computed easily while the classical far field operator uses full aperture measured far field pattern.

We expect that the (approximate) solution of (3.2) for a sampling point z would provide useful information for the reconstruction of D. To this end, we first introduce some results that are useful in analyzing the far field equation (3.2) (see corollary 5.31, corollary 5.32 and theorem 3.22 of [4]).

Lemma 3.1. The Herglotz operator $\mathcal{H}: L^2(\mathbb{S})\rightarrow H^{{1}/{2}}(\partial B_z)$ defined by

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle (\mathcal{H}g)(x):=\int_{\mathbb{S}} {\rm e}^{{\rm i}kx\cdot d}g(d){\rm d}s(d),\ \ \ \ x\in \partial B_z \nonumber \end{align} \tag{ 3.3 }$$

is injective and has a dense range provided k² is not a Dirichlet eigenvalue for the negative Laplacian for B_z.

Lemma 3.2. The operator $A:H^{1/2}(\partial B_z)\rightarrow L^2(\mathbb{S})$ which maps the boundary values of radiating solutions $u\in H_{\rm loc}^1(\mathbb{R}^2\setminus\overline{B}_z)$ of the Helmholtz equation onto the far field pattern $u_\infty$ is bounded, injective and has a dense range.

Lemma 3.3. The far field operator $\mathcal{F}_z$ is injective and has a dense range provided that k² is not a Dirichlet eigenvalue for the negative Laplacian for B_z.

The following theorem is the main result for the far field equation (3.2).

Theorem 3.1. Let B_z be a sound-soft disc centered at z with radius R. Let D be an inhomogeneous medium or an obstacle with Dirichlet, Neumann, or the impedance boundary condition. Assume that kR does not coincide with any zero of the Bessel functions $J_n, n=0, 1, 2, \cdots$. Then the following results hold for the far field equation (3.2):

• 1.  
  If $D\subset B_z$, for a given $\varepsilon>0$, there exists a function $g_z^\varepsilon\in L^2(\mathbb{S})$ such that

  $$\begin{align} \newcommand{\bi}{\boldsymbol} \newcommand{\e}{{\rm e}} \displaystyle \label{fe2} \bigg\|\int_{\mathbb{S}}u^{B_z}_{\infty}(\hat{x},d)g_z^\varepsilon(d){\rm d} s(d)-U_\infty(\hat{x})\bigg\|_{L^2(\mathbb{S})}<\varepsilon \nonumber \end{align} \tag{ 3.4 }$$

  and the Herglotz wave function $v_{g_z^\varepsilon}(x):=\int_{\mathbb{S}}{\rm e}^{{\rm i}kx\cdot d}g_z^\varepsilon (d){\rm d}s(d), x\in B_z$ converges to the solution $w\in H^1(B_z)$ of the Helmholtz equation with w  =  −U^s on $\partial B_z$ as $\varepsilon\rightarrow 0$.
• 2.  
  If $\newcommand{\e}{{\rm e}} D\cap B_z=\emptyset$, every $g_z^\varepsilon\in L^2(\mathbb{S})$ that satisfies (3.4) for a given $\varepsilon>0$ is such that

  $$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \lim_{\varepsilon\rightarrow 0}\|v_{g_z^\varepsilon}\|_{H^1(B_z)}=\infty. \nonumber \end{align} \tag{ 3.5 }$$
• 3.  
  Let $\newcommand{\e}{{\rm e}} D\cap B_z\neq\emptyset$ and $D \not\subset B_z$. When the scattered solution of D can be extended from $\mathbb{R}^2\setminus D$ into $\mathbb{R}^2\setminus (D\cap B_z)$, then the conclusion is the same with case 1; otherwise, the conclusion is the same with case 2.

Proof. Since kR does not coincide with any zero of the Bessel functions $J_n, n=0, 1, 2, \cdots$, then k² is not a Dirichlet eigenvalue for the following eigenvalue problem for B_z [22]:

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \label{dep} \left\{\begin{array}{@{}ll@{}} \Delta u+k^2 u=0,\ \ \ {\rm in}\ B_z, \nonumber \\ u=0,\ \ \ {\rm on} \ \partial B_z. \end{array}\right. \nonumber \end{align} \tag{ 3.6 }$$

• Case 1: let $D\subset B_z$. From lemma 3.1, for any ε, there exists $g_z^\varepsilon$ such that

  $$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \label{equ2.1} \|\mathcal{H}g_z^\varepsilon+U^s\|_{L^2(\partial B_z)}\leqslant \frac{\varepsilon}{\|A\|}, \nonumber \end{align} \tag{ 3.7 }$$

  where $A:H^{1/2}(\partial B_z)\rightarrow L^2(\mathbb{S})$ is defined in lemma 3.2. Then

  $$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \label{equ2.2} A(-\mathcal{H}g_z^\varepsilon)=\int_\mathbb{S}u^{B_z}_{\infty}(\cdot,d)g_z^\varepsilon(d){\rm d} s(d). \nonumber \end{align} \tag{ 3.8 }$$

  Since $D\subset B_z$, we also have

  $$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \label{equ2.3} A(U^s)=U_\infty. \nonumber \end{align} \tag{ 3.9 }$$

  Taking the difference of (3.8) and (3.9), and from lemma 3.2 and (3.7), for any ε, we have $g_z^\varepsilon$ such that

  $$\begin{align} \newcommand{\bi}{\boldsymbol} \newcommand{\e}{{\rm e}} \displaystyle &\bigg\|\int_{\mathbb{S}}u^{B_z}_{\infty}(\hat{x},d)g_z^\varepsilon(d){\rm d} s(d)-U_\infty(\hat{x})\bigg\|_{L_2(\mathbb{S})}\nonumber \\ &=\|A(-\mathcal{H}g_z^\varepsilon)-A(U^s)\|_{L^2(\mathbb{S})}\nonumber \\ &\leqslant \|A\|\cdot \|\mathcal{H}g_z^\varepsilon+U^s\|_{L^2(\mathbb{S})}\nonumber \\ &\leqslant \varepsilon. \nonumber \end{align} \tag{ 3.10 }$$

  By (3.7), as $\varepsilon\rightarrow 0$, the Herglotz wave function $v_{g_z^\varepsilon}(x):=\int_{\mathbb{S}}{\rm e}^{{\rm i}kx\cdot d}g_z^\varepsilon (d){\rm d}s(d), x\in B_z$ converges to the unique solution $w\in H^1(B_z)$ of the following well-posed inner Dirichlet problem

  $$\begin{align*} \newcommand{\e}{{\rm e}} \displaystyle \left\{\begin{array}{@{}ll@{}} \Delta w+k^2 w=0,\ \ \ {\rm in}\ B_z, \nonumber \\ w=-U^s,\ \ \ {\rm on} \ \partial B_z. \end{array}\right. \nonumber \end{align*}$$
• Case 2: let $\newcommand{\e}{{\rm e}} D\cap B_z=\emptyset$. From lemma 3.3, for every ε, there exists $g_z^\varepsilon\in L^2(\mathbb{S})$ such that

  $$\begin{align} \newcommand{\bi}{\boldsymbol} \newcommand{\e}{{\rm e}} \displaystyle \label{equ2.4} \bigg\|\int_\mathbb{S} u^{B_z}_{\infty}(\hat{x},d)g_z^\varepsilon(d){\rm d}s(d)-U_\infty(\hat{x})\bigg\|_{L_2(\mathbb{S})}<\varepsilon. \nonumber \end{align} \tag{ 3.11 }$$

  Assume to the contrary that there exists a sequence $v_{g_z^{\varepsilon_n}}$ such that $\|v_{g_z^{\varepsilon_n}}\|_{H^1(B_z)}$ remains bounded as $\varepsilon_n \to 0, n \to \infty$. Without loss of generality, we assume that $v_{g_z^{\varepsilon_n}}$ converges to $v_{g_z}\in H^1(B_z)$ weakly as $n\rightarrow \infty$, where $v_{g_z}(x)=\int_{\mathbb{S}}{\rm e}^{{\rm i}kx\cdot d}g_z (d){\rm d}s(d), x\in B_z$. Let $V^s\in H_{\rm loc}^1(\mathbb{R}^2\setminus \overline{B}_z)$ be the unique solution of the following exterior Dirichlet problem

  $$\begin{align*} \newcommand{\e}{{\rm e}} \displaystyle \left\{\begin{array}{@{}lll@{}} \Delta V^s+k^2V^s=0\ \ \ {\rm in}\ \mathbb{R}^2\setminus \overline{B}_z, \nonumber \\ V^s=-v_{g_z} \ \ \ {\rm on}\ \partial B_z, \nonumber \\ \lim\limits_{r\rightarrow\infty}\sqrt{r}(\partial V^s/\partial r-{\rm i}kV^s)=0. \end{array}\right. \nonumber \end{align*}$$

  Its far field pattern is $V_\infty=\int_\mathbb{S} u^{B_z}_{\infty}(x, d)g_z(d){\rm d}s(d)$. While from (3.11), as $\varepsilon\rightarrow 0$, we have $\int_\mathbb{S} u^{B_z}_{\infty}(x, d)g_z(d){\rm d}s(d)=U_\infty$. Consequently,

  $$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle V_\infty=U_\infty. \nonumber \end{align} \tag{ 3.12 }$$

  Then by Rellich's lemma (see lemma 2.12 in [4]), the scattered fields coincide in $\mathbb{R}^2\setminus (\overline{B}_{z}\cup \overline{D})$ and we can identify $\tilde{U}^s:=V^s=U^s$ in $\mathbb{R}^2\setminus (\overline{B}_{z}\cup \overline{D})$. We have that $V^s$ has an extension into $\mathbb{R}^2\setminus B_z$ and U^s has an extension into $\mathbb{R}^2\setminus D$. Since $\newcommand{\e}{{\rm e}} D\cap B_z=\emptyset$, $\tilde{U}^s$ can be extended from $\mathbb{R}^2\setminus (\overline{B}_{z}\cup \overline{D})$ into all of $\mathbb{R}^2$, that is, $\tilde{U}^s$ is an entire solution to the Helmholtz equation. Since $\tilde{U}^s$ also satisfies the radiation condition, it must vanish identically in all of $\mathbb{R}^2$. This leads to a contradiction since $\tilde{U}^s$ is not a null function.
• Case 3: if the scattered field of D can be extended from $\mathbb{R}^2\setminus D$ into $\mathbb{R}^2\setminus (D\cap B_z)$, then $\tilde{D}:=D\cap B_z$ has the far field pattern $U_\infty$. One can simply replace D with $\tilde{D}$ in case 1 to arrive the same conclusion. Otherwise, if the scattered solution of D cannot be extended from $\mathbb{R}^2\setminus D$ into $\mathbb{R}^2\setminus (D\cap B_z)$, a similar proof to case 2 works. □

Remark 3.2. Note that the radius R of the disc B_z can be chosen such that kR is not a Dirichlet eigenvalue for the negative Laplacian for B_z. Hence the ESM is independent of wavenumbers, i.e. the ESM works for all wavenumbers. In contrast, the LSM needs to exclude wavenumbers which are certain eigenvalues related to the scatterer.

Now we are ready to present the ESM to reconstruct the location and approximate the support of a scatterer D using the far field pattern due to a single incident wave.

Recall that $u^{B_z}_{\infty}(\cdot, d)$ is the far field pattern of a sound-soft disc B_z centered at z with radius R and $U_\infty(\hat{x})$ is the measured far field pattern of an unknown scatterer D due to an incident plane wave. Let Ω be a domain containing D. For a sampling point $z \in \Omega$, we consider the linear ill-posed integral equation

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \label{ffe} \int_{\mathbb{S}}u^{B_z}_{\infty}(\hat{x},d)g_z(d){\rm d} s(d)=U_\infty(\hat{x}),\ \ \ \ \hat{x}\in \mathbb{S}. \nonumber \end{align} \tag{ 3.13 }$$

Suppose that (3.13) is solved by some regularization scheme, say Tikhonov regularization with parameter α. By theorem 3.1, one expects that the approximate solution $\|g^\alpha_z\|_{L^2(\mathbb{S})}$ is relatively large when D is not inside B_z and relatively small when D is inside B_z. Consequently, an approximation of the location and support of the scatterer D can be obtained by plotting $\|g^\alpha_z\|_{L^2(\mathbb{S})}$ for all sampling points $z \in \Omega$.

The algorithm of the ESM is as follows (see figure 1).

Zoom In Zoom Out Reset image size

3.2.1. The ESM.

An important step of the ESM is to choose the radius R of B_z. It would be ideal to choose B_z to be slightly larger than the scatterer D. However, this is not possible if no a priori information about the scatterer is available. To resolve the difficulty, we propose a multilevel technique to choose R.

3.3.1. Multilevel ESM.

The far field equation (3.2) of the ESM looks similar to the LSM proposed by Colton and Kirsch [3]. The far field equation for the LSM is

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \label{FELSM} \int_{\mathbb{S}}U_\infty(\hat{x}; d)g_z(d){\rm d}s(d)=\Phi_\infty(\hat{x}; z), \ \ \ \ \hat{x}\in \mathbb{S}. \nonumber \end{align} \tag{ 3.14 }$$

In (3.14), $U_\infty(\hat{x};d), \hat{x}, d\in \mathbb{S}$ is the full-aperture far field pattern of the scatterer D due to the incident plane wave $\newcommand{\e}{{\rm e}} u^i=\exp({\rm i}kx\cdot d)$, and $\Phi_\infty(\hat{x}, z)$ is the far field pattern of the fundamental solution of the Helmholtz equation $\Phi(x, z):=\frac{\rm i}{4}H^{(1)}_0(k|x-z|)$, where $H^{(1)}_0$ is the Hankel function of the first kind of order zero. For each sampling point z in the sampling region Ω which contains D, using some appropriate regularization scheme, one obtains an approximate solution $g^\alpha_z$. In general, the norm $\|g^\alpha_z\|$ is larger for $z\not \in D$ and smaller for $z\in D$ such that the support of D can be reconstructed accordingly.

Since $U_\infty(\hat{x}; d)$ is used as the kernel of the integral operator in (3.14), full aperture far field pattern is necessary. Namely, one needs $U_\infty(\hat{x}; d)$ for all $\hat{x}\in \mathbb{S}$ and $d \in \mathbb{S}$ to set up (3.14). As a consequence, the LSM cannot directly process limited aperture far field pattern, e.g. $U_\infty(\hat{x}; d_0)$ for a single incident direction d₀. This also applies to other classical sampling methods such as the Factorization Method and the Reciprocity Gap Method.

In contrast, the kernel of the integral operator $\mathcal{F}_z$ for the ESM is the full aperture far field pattern $u_{\infty}^{B_z}(\hat{x};d)$ of a sound-soft disc B_z. The measured far field pattern $U_{\infty}(\hat{x})$ is the right hand side of the far field equation (3.2). Hence, in principle, the ESM can take limited aperture data of any type as input.

Another related method is the range test due to Potthast, Sylvester and Kusiak [19]. Similar to the ESM, it processes $U_\infty(\hat{x})$ for a single incident wave. Let G be a convex test domain. The integral equation in [19] is given by

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \label{RTIE} \int_{\partial G} {\rm e}^{{\rm i} k \hat{x} \cdot{d}} g(d) {\rm d}s(d) = U_\infty(\hat{x}), \quad \hat{x} \in \mathbb{S}. \nonumber \end{align} \tag{ 3.15 }$$

It has a solution in $L^2(\partial G)$ if and only if the scattered field can be analytically extended up to the boundary $\partial G$. Then a convex support can be obtained numerically and the intersection of these supports provides information to reconstruct D.

The ESM uses a slightly different integral equation. More importantly, other than finding a convex domain, the ESM computes an indicator for a sampling point, which makes it possible to process data of multiple directions or even multiple frequencies.

Let $\Omega=[-10, 10]\times [-10, 10]$ and choose the sampling points to be

$$\begin{align*} \newcommand{\e}{{\rm e}} \displaystyle T:=\{(-10+0.1m, -10+0.1n), \quad m,n=0, 1, \cdots, 200\}. \nonumber \end{align*}$$

The radius of the sampling discs is set to R  =  1.

We discretize (4.2) to obtain a linear system $A^zg_z=F^z$, where A^z is the matrix given by

$$\begin{align*} \newcommand{\e}{{\rm e}} \displaystyle A^z_{l,j}={\rm e}^{{\rm i} kz\cdot d_j}u^B_{\infty}(d_l,d_j), \quad l, j=1,2,\cdots,52. \nonumber \end{align*}$$

Then the regularized solution is given by

$$\begin{align*} \newcommand{\bi}{\boldsymbol} \newcommand{\e}{{\rm e}} \displaystyle g^\alpha_z\approx \big((A^z)^*A^z+\alpha A^z\big)^{-1}(A^z)^*F^z. \nonumber \end{align*}$$

We plot the contours for the indicator function

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \label{Iz} I_z = \frac{\|g^\alpha_z\|_{l^2}}{\max\limits_{z\in T}\|g^\alpha_z\|_{l^2}} \nonumber \end{align} \tag{ 4.3 }$$

for all the sampling points $z\in T$.

Figure 2(a) is the triangular scatterer. Figures 2(b)–(d) are the contour plots of I_z for obstacles with Dirichlet, Neumann, and the impedance boundary conditions, respectively. The asterisks are the minimum points of I_z.

Zoom In Zoom Out Reset image size

We plot discs whose centers are z's minimizing I_z in figure 3 for different boundary conditions. These minimization points provide the correct locations of the obstacles.

Zoom In Zoom Out Reset image size

Figure 4 shows similar results of the EMS for the kite with Dirichlet, Neumann and the impedance boundary conditions. Figure 5 shows the reconstructions of the kite of different boundary conditions. Figure 6 show the results for the L-shape medium.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Next we show the reconstructions of the sound-soft triangle in figure 7 to demonstrate that the ESM works for different wavenumbers. The wavenumbers are k  =  0.2, k  =  3 and k  =  2.4048. In particular, k  =  2.4048 is a Dirichlet eigenvalue of the unit disc, which is not covered by theorem 3.1. Nonetheless, ESM reconstructs the target effectively. In fact, the ESM can avoid a particular wavenumber, e.g. a Dirichlet eigenvalue, by choosing a slightly different sampling disc.

Zoom In Zoom Out Reset image size

If the size of the scatterer is not known, one needs to decide the proper radius R of the sampling discs. We start the multilevel ESM using the sampling discs with a large radius R  =  2.4. Then we decrease the radius until a suitable R is found.

In the following figures, the solid lines are the reconstructions, red dashed lines are the exact scatterers, asterisks are the minimum points of I_z. Figure 8 shows the results of the multilevel ESM for the triangle object with Dirichlet boundary condition. The radius of the sampling discs reduces by half from 2.4 to 0.3. When the radius of the sampling discs is R  =  0.3, the minimum point is not inside the reconstruction for R  =  0.6. Hence the proper radius of the sampling discs should be R  =  0.6. Consequently, we obtain the estimation of the size of the target.

Zoom In Zoom Out Reset image size

Using sampling discs with radius 0.6, figure 9(a) shows the contours of I_z on the sampling points $(-10+0.1m, -10+0.1n), 0\leqslant m, n\leqslant 201$, and figure 9(b) shows the reconstruction result of the multilevel ESM for the triangle with Dirichlet boundary condition.

Zoom In Zoom Out Reset image size

Using the multilevel ESM, the proper radius R is 0.6 for the triangle and kite with different boundary conditions, and the proper radius R is computed to be 0.3 for the L-shaped medium scatterer. Figure 10 shows the reconstructions of the multilevel ESM for the triangle with Neumann and the impedance boundary condition ($\lambda=2$), the kite with Dirichlet, Neumann and the impedance boundary condition ($\lambda=2$), and the L-shape medium.

Zoom In Zoom Out Reset image size