An inverse medium problem using Stekloff eigenvalues and a Bayesian approach

TL;DR

This paper introduces a novel approach combining integral equations, Bayesian estimation, and spectral indicators to reconstruct Stekloff eigenvalues and the index of refraction of inhomogeneous media from boundary data, validated by numerical experiments.

Contribution

It presents a new integral equation for the reciprocity gap method and integrates Bayesian and spectral indicator techniques for inverse spectral problems.

Findings

Effective reconstruction of Stekloff eigenvalues demonstrated

Bayesian method accurately estimates the index of refraction

Numerical experiments confirm the approach's robustness

Abstract

This paper studies the reconstruction of Stekloff eigenvalues and the index of refraction of an inhomogeneous medium from Cauchy data. The inverse spectrum problem to reconstruct Stekloff eigenvalues is investigated using a new integral equation for the reciprocity gap method. Given reconstructed eigenvalues, a Bayesian approach is proposed to estimate the index of refraction. Moreover, since it is impossible to know the multiplicities of the reconstructed eigenvalues and since the eigenvalues can be complex, we employ the recently developed spectral indicator method to compute Stekloff eigenvalues. Numerical experiments validate the effectiveness of the proposed methods.

Tables

Table 1. Table 1: The exact Stekloff eigenvalues and their reconstructions (in the parentheses) for n ​ ( x ) = 5 𝑛 𝑥 5 n(x)=5 .

disc	$1.2937$	$-0.4763$	$-0.5839$	$-1.2301$
	($1.30$)	($-0.48$)	($-0.58$)	($-1.23$)
square	$0.3792$	$-0.5418$	$-0.6148$
	($0.38$)	($-0.54$)	($-0.62)$
L-shaped	$3.0218$	$0.6882$	$-0.5266$	$-0.5834$	$-1.2188$
	($3.04$)	($0.70$)	($-0.52$)	($-0.58$)	($-1.20$)

Table 2. Table 2: The exact Stekloff eigenvalues and their reconstructions (in the parentheses) for n ​ ( x ) = 4 + 2 ​ | x | 𝑛 𝑥 4 2 𝑥 n(x)=4+2|x| .

disc	$2.0856$	$-0.4898$	$-0.5714$	$-1.2285$
	($2.10$)	($-0.48$)	($-0.54$)	($-1.20$)
square	$0.4137$	$-0.5758$	$-1.2348$
	($0.42$)	($-0.58$)	($-1.22$)
L-shaped	$0.9825$	$-0.4956$	$-0.6018$	$-1.2116$
	($1.00$)	($-0.50$)	($-0.60$)	($-1.20$)

Table 3. Table 3: The exact Stekloff eigenvalues and their reconstructions (in the parentheses) for n ​ ( x ) = 2 + 4 ​ i 𝑛 𝑥 2 4 𝑖 n(x)=2+4i .

disc	$-0.0549+0.4854i$	$-0.0211+0.2652i$	$-0.6361+0.0390i$
	($-0.06+0.46i$)	-0.02+0.26i	($-0.64+0.04i$)
square	$0.1303+0.5812i$	$0.0779+0.1415i$	$-0.6338+0.0170i$
	($0.12+0.56i$)	($0.08+0.14i$)	($-0.64+0.02i$)
L-shaped	$-0.0902+0.5468i$	$-0.1029+0.3600i$	$0.0364+0.2228i$
	($-0.10+0.52i$)	($-0.10+0.36i$)	($0.04+0.22i$)
	$-0.6385+0.0764i$	$-0.6346+0.0412i$
	($-0.64+0.08i$)	($-0.64+0.04i$)

Table 4. Table 4: The posterior means for three domains ( n ​ ( x ) = 5 𝑛 𝑥 5 n(x)=5 ) using one Stekloff eigenvalue. The value in the parentheses for the L-shaped domain is obtained using two Stekloff eigenvalues.

${\sigma }^2$	domain	$n_{CM}$
$0.05$	circle	4.9953
	square	5.0205
	L-shaped domain	6.2135 (5.0074)

Table 5. Table 5: The posterior mean of n ​ ( x ) 𝑛 𝑥 n(x) for three domains ( n ​ ( x ) = 4 + 2 ​ | x | 𝑛 𝑥 4 2 𝑥 n(x)=4+2|x| ).

${\sigma }^2$	domain	mean of $n_0$	${\beta }_{CM}^1+{\beta }_{CM}^2\left|x\right|$
$0.05$	circle	$4.9804$	$4.3916+1.9333\left|x\right|$
	square	$6.7791$	$3.8873+2.3409\left|x\right|$
	L-shaped	$6.2197$	$4.2953+1.7263\left|x\right|$

Table 6. Table 6: The posterior means of n ​ ( x ) 𝑛 𝑥 n(x) for three domains ( n ​ ( x ) = 2 + 4 ​ i 𝑛 𝑥 2 4 𝑖 n(x)=2+4i ).

${\sigma }^2$	domain	$n_{CM}$
$0.05$	circle	$1.8511+3.9849i$
	square	$2.2515+4.1935i$
	L-shaped	$2.1204+4.1978i$

Equations

$$\left\{\begin{array}[]{lll}&\Delta u+k^{2}n(x)u=0,\ \ \ \textrm{in}\ \mathbb{R}^{2}\setminus\{x_{0}\},\\ &u=u^{s}+u^{i},\\ &\lim\limits_{{r}\rightarrow\infty}r^{\frac{1}{2}}(\partial u^{s}/\partial r-{i}ku^{s})=0,\ \ \ r=|x|,\end{array}\right.$$

$$u^i:=\Phi (\cdot ,x_0)x_0\in \mathbb{R}\setminus \overline{D}$$

$$\left\{\begin{array}[]{ll}&\Delta w+k^{2}n(x)w=0,\ \ \ \textrm{in}\ B,\\ &\partial w/\partial\nu+\lambda w=0,\ \ \ \textrm{on}\ \Gamma,\end{array}\right.$$

$$\left(\nabla w,\nabla v\right)-k^2\left(nw,v\right)=-\lambda ⟨w,v⟩\forall v\in H^1(B),$$

$$\left\{\begin{array}[]{lll}&\Delta u+k^{2}nu=0,\ \ \ \textrm{in}\ B,\\ &u=0,\ \ \ \textrm{on}\ \Gamma.\end{array}\right.$$

$$n(x):=n_{c}=\left\{\begin{array}[]{lll}&1,\ \ \ \textrm{in}\ B\setminus\overline{D},\\ &c,\ \ \ \textrm{in}\ D.\end{array}\right.$$

$$\delta n_{c}:=\left\{\begin{array}[]{lll}&0,\ \ \ \textrm{in}\ B\setminus\overline{D},\\ &\delta c,\ \ \ \textrm{in}\ D,\end{array}\right.$$

$$\big{(}\nabla(w+\delta w),\nabla v\big{)}-k^{2}\big{(}(n_{c}+\delta n_{c})(w+\delta w),v\big{)}=-(\lambda^{-}_{1}+\delta\lambda^{-}_{1})\big{\langle}w+\delta w,v\big{\rangle}\quad\forall\,v\in H^{1}(B).$$

$$\displaystyle\big{(}\nabla\delta w,\nabla v\big{)}-k^{2}\big{(}\delta n_{c}(w+\delta w),v\big{)}-k^{2}\big{(}n_{c}\delta w,v\big{)}$$

$$\displaystyle\qquad\qquad\qquad\qquad=-\delta\lambda^{-}_{1}\big{\langle}w+\delta w,v\big{\rangle}-\lambda^{-}_{1}\big{\langle}\delta w,v\big{\rangle}\quad\forall\,v\in H^{1}(B).$$

$$k^{2}\big{(}\delta n_{c}(w+\delta w),w\big{)}=\delta\lambda^{-}_{1}\langle w+\delta w,w\rangle,$$

$$\delta {\lambda }_1^{-}$$

$$|(\delta w,w{)}_D|<\frac{1}{2}(w,w{)}_D\text{and}|⟨\delta w,w⟩|<\frac{1}{2}⟨w,w⟩.$$

$$0<\frac{k^{2}\delta c\big{(}w,w\big{)}_{D}}{3\langle w,w\rangle}\leq\delta\lambda^{-}_{1}\leq\frac{3k^{2}\delta c\big{(}w,w\big{)}_{D}}{\langle w,w\rangle}.$$

$$\left\{\begin{array}[]{lll}&\Delta u+k^{2}n_{c}u=0,\ \ \ \textrm{in}\ B,\\ &\frac{\partial u}{\partial\nu}=0,\ \ \ \textrm{on}\ \Gamma.\end{array}\right.$$

$$\left\{\begin{array}[]{lll}&\Delta u_{\lambda}+k^{2}u_{\lambda}=0,\ \ \ \textrm{in}\ \mathbb{R}^{2}\setminus\{\overline{B}\cup\{x_{0}\}\},\\ &\partial_{\nu}u_{\lambda}+\lambda u_{\lambda}=0,\ \ \ \textrm{on}\ \Gamma,\\ &\lim\limits_{{r}\rightarrow\infty}r^{\frac{1}{2}}(\partial u_{\lambda}^{s}/\partial r-{i}ku_{\lambda}^{s})=0,\ \ \ r=|x|,\end{array}\right.$$

$$R(v_1,v_2)={\int }_{\Gamma }(v_1{\partial }_{\nu }{v}_2-v_2{\partial }_{\nu }{v}_1)ds,$$

$$R(u_{\lambda }(\cdot ,x_0)-u(\cdot ,x_0),v_g(\cdot ))=R(u_{\lambda }(\cdot ,x_0),{\Phi }_z(\cdot ))\forall x_0\in C,$$

$$v_g(x):={\int }_{\mathbb{S}}{e}^{ikx\cdot d}g(d)ds(d).$$

$$\left\{\begin{array}[]{lll}&\Delta\tilde{u}^{s}+k^{2}\tilde{u}^{s}=0,\ \ \ \textrm{in}\ \mathbb{R}^{2}\setminus\overline{B},\\ &\partial_{\nu}\tilde{u}^{s}+\lambda\tilde{u}^{s}=f,\ \ \ \textrm{on }\Gamma,\\ &\lim\limits_{{r}\rightarrow\infty}r^{\frac{1}{2}}(\partial\tilde{u}^{s}/\partial r-{i}k\tilde{u}^{s})=0,\ \ \ r=|x|.\end{array}\right.$$

$${\tilde{u}}^s(x_0)$$

$$\mathcal{R}(g):=R\big{(}u_{\lambda}(\cdot,x_{0})-u(\cdot,x_{0}),v_{g}(\cdot)\big{)},\quad x_{0}\in C$$

$$\left\{\begin{array}[]{ll}&\Delta w^{s}+k^{2}nw^{s}=k^{2}(1-n)v_{g},\ \ \ \textrm{in}\ \mathbb{R}^{2},\\ &\lim\limits_{{r}\rightarrow\infty}r^{\frac{1}{2}}(\partial w^{s}/\partial\nu-{i}kw^{s})=0,\ \ \ r=|x|.\end{array}\right.$$

$$\displaystyle\int_{\Gamma}u_{\lambda}\big{(}\partial_{\nu}(w^{s}+v_{g})+\lambda(w^{s}+v_{g})\big{)}ds$$

$$\left\{\begin{array}[]{ll}&\Delta(w^{s}+v_{g})+k^{2}n(w^{s}+v_{g})=0,\ \ \ \textrm{in}\ B,\\ &\partial_{\nu}(w^{s}+v_{g})+\lambda(w^{s}+v_{g})=0\ \ \ \textrm{on}\ \Gamma.\end{array}\right.$$

$$\lim\limits_{n\rightarrow\infty}R\big{(}u_{\lambda}-u,v_{g_{n}}\big{)}=R(u_{\lambda},\Phi_{z}),\ \ \ \ u_{\lambda}\in U_{\lambda},\ u\in U$$

$$\lim\limits_{n\rightarrow\infty}R\big{(}u_{\lambda}-u,v_{g_{n}^{z}}\big{)}=R(u_{\lambda},\Phi_{z}),\ \ \ \ u_{\lambda}\in U_{\lambda},\ u\in U,$$

$$\left\{\begin{array}[]{ll}&\Delta w_{z}+k^{2}nw_{z}=0,\ \ \ \textrm{in}\ B,\\ &\partial_{\nu}w_{z}+\lambda w_{z}=\partial_{\nu}\Phi(\cdot,z)+\lambda\Phi(\cdot,z)\ \ \ \textrm{on}\ \Gamma.\end{array}\right.$$

$$w_z=w_z^i+w_z^s,$$

$$\left\{\begin{array}[]{ll}&\Delta w^{s}_{z}+k^{2}nw^{s}_{z}=k^{2}(1-n)w_{z}^{i},\ \ \ \textrm{in}\ \mathbb{R}^{2},\\ &\lim\limits_{{r}\rightarrow\infty}r^{\frac{1}{2}}(\partial w^{s}_{z}/\partial r-{i}kw^{s}_{z})=0,\ \ \ r=|x|.\end{array}\right.$$

$$v_{g_n}+w_z^s\to w_z^i+w_z^s=w_z,n\to \infty .$$

Full text

Key words: inverse medium problem, inverse spectrum problem, Stekloff eigenvalues, reciprocity gap, Bayesian approach, spectral indicator method

1 Introduction

Inverse scattering problems for inhomogeneous media have many applications such as medical imaging and nondestructive testing. In this paper, the inverse spectrum problem to reconstruct the Stekloff eigenvalues from Cauchy data is investigated first using a new integral equation for the reciprocity gap method. Then these eigenvalues are used to estimate the index of refraction of the inhomogeneous medium. Due to the lack of knowledge of the relation between Stekloff eigenvalues and the index of refraction, we propose a Bayesian approach. Since the eigenvalues are complex for absorbing media and the multiplicities are not known, the recently developed spectral indicator method is employed to compute the Stekloff eigenvalues [17, 18].

The reconstruction of certain eigenvalues from the scattering data has been studied by many researchers. In the context of qualitative methods in inverse scattering, it has been shown that interior eigenvalues such as Dirichlet eigenvalues and transmission eigenvalues can be determined from the scattering data [4, 31, 26] (see also the special issue edited by Lechleiter and Sun [23]). A related method, which can be used to compute interior eigenvalues using the scattering data, is the inside-outside duality [20, 21, 22, 29].

Given reconstructed eigenvalues, a legitimate question is what information about the obstacle can be obtained. For inhomogeneous non-absorbing media, transmission eigenvalues have been used to reconstruct the shape of the obstacle [32] and obtain useful information of the index of refraction [6, 31, 1, 14, 3, 24]. However, the use of transmission eigenvalues has two drawbacks: 1) multi-frequency data are necessary; and 2) only real transmission eigenvalues can be determined from the scattering data so far.

It has been shown that Stekloff eigenvalues associated with the scattering problem can be determined from far field data of a single frequency [5, 1]. Unlike transmission eigenvalues, Stekloff eigenvalues exist for absorbing media as well. Hence the use of Stekloff eigenvalues avoids the above two drawbacks and has the potential to work for a wider class of problems. In this paper, a new integral equation for the reciprocity gap (RG) method [8, 10, 28] is introduced to determine Stekloff eigenvalues from Cauchy data. Then a Bayesian approach is proposed to estimate the index of refraction. The Metropolis-Hastings (M-H) Algorithm is used to explore the posterior distribution. Numerical examples show that the proposed methods are effective. We refer the readers to [19, 30] and references therein on the Bayesian framework for inverse problems and [2, 11, 34, 15] for the Bayesian methods for some inverse scattering problems.

The rest of the paper is organized as follows. In Section 2, the forward scattering problem and the associated Stekloff eigenvalue problem are introduced. In Section 3, a new integral equation for the reciprocity gap method is proposed to reconstruct Stekloff eigenvalues using Cauchy data. In Section 4, a Bayesian approach and the MCMC method are proposed to estimate the index of refraction. Finally, numerical examples are provided in Section 5.

2 Scattering Problem and Stekloff Eigenvalues

In this section, we introduce the direct scattering problem, the Stekloff eigenvalue problem, and the inverse scattering problems using Cauchy data. Then a monotonicity of the largest negative Stekloff eigenvalue is proved.

Let $D$ be a bounded domain in $\mathbb{R}^{2}$ with boundary $\partial D$ of class $C^{2}$. Let $k$ be the wavenumber and $n(x)$ be the index of refraction such that $n(x)\in L^{\infty}(\mathbb{R}^{2})$. Assume that $n(x)=1$ for $\mathbb{R}^{2}\setminus\overline{D}$ and $\Re(n(x))>0,\Im(n(x))\geq 0$, where $\Re(\cdot)$ and $\Im(\cdot)$ denote the real and imaginary parts, respectively. The direct scattering problem is to find the total field $u$ such that

[TABLE]

where

[TABLE]

is the incident wave generated by a point source. Here $\Phi$ is the fundamental solution of the Helmholtz equation.

The associated Stekloff eigenvalue problem is defined as follows [5]. Find $\lambda\in\mathbb{C}$ and a non-trivial function $w$ such that

[TABLE]

where $B$ be a bounded domain in $\mathbb{R}^{2}$ and $\Gamma:=\partial B$ such that $D\subset B$.

Assume that the Cauchy data $u$ and $\partial_{\nu}u:=\partial u/\partial\nu$ are known on $\Gamma$ for each the incident wave $u^{i}:=\Phi(\cdot,x_{0}),x_{0}\in C$, where $C$ is a simple closed curve containing $B$ (see 1). The inverse scattering problems considered in this paper are:

• IP1

Reconstruct Stekloff eigenvalues from Cauchy data;

• IP2

Estimate the index of refraction $n(x)$ using Stekloff eigenvalues.

The weak formulation for the Stekloff eigenvalue problem (2.2) is to find $(\lambda,u)\in\mathbb{C}\times H^{1}(B)$ such that

[TABLE]

where $\left(f,g\right)=\int_{B}f\overline{g}dx\ \ \textrm{and}\ \ \left\langle f,g\right\rangle=\int_{\Gamma}f\overline{g}ds$.

When $n(x)$ is real, all Stekloff eigenvalues are real and they form an infinite discrete set [5]. We call $k^{2}$ a modified Dirichlet eigenvalue of $B$ if there exists a nontrivial $u\in H^{1}(B)$ such that

[TABLE]

Remark 2.1

Note that a standard Dirichlet eigenvalue problem is such that $n(x)\equiv 1$ in (2.4). For simplicity, in the rest of the paper, we call $k^{2}$ in (2.4) a Dirichlet eigenvalue.

It is shown in [1] that Stekloff eigenvalues accumulate at $-\infty$ if $k^{2}$ is not a Dirichlet eigenvalue. Next, we prove a property of the largest negative Stekloff eigenvalue $\lambda^{-}_{1}$ when $n(x)$ is given by

[TABLE]

Suppose $n_{c}$ is perturbed by

[TABLE]

where $\delta c$ is also a real constant. The perturbation $\delta n_{c}$ leads to $\delta w$ and $\delta\lambda^{-}_{1}$ of the eigenpair. From (2.3), $\delta w\in H^{1}(B)$ and $\delta\lambda^{-}_{1}$ satisfies

[TABLE]

Using the fact that $(w,\lambda^{-}_{1})$ is a real eigenpair, we have that

[TABLE]

Letting $v=w$ and noting that $n_{c}$ is real, we have that

[TABLE]

which implies that

[TABLE]

where $\left(f,g\right)_{D}=\int_{D}f\overline{g}dx$. If $\delta c>0$ is small enough, one has that

[TABLE]

From (2.6) and (2.7), we have

[TABLE]

This implies that $\lambda^{-}_{1}$ is monotonically increasing with respect to $n_{c}$. This breaks until $k^{2}$ becomes a (modified) Neumann eigenvalue, i.e., there exists a non-trivial $u$ such that

[TABLE]

Note that a standard eigenvalues is $k^{2}$ satisfying (2.9) for $n_{c}\equiv 1$. Again, we $k^{2}$ a Neumann eigenvalue for simplicity. Excluding this case, we actually proved the following theorem.

Theorem 2.2

Let the index of refraction be defined in (2.5) and $[a,b]$ be an interval that $k^{2}$ is not a Neumann eigenvalue of (2.9) for any $c\in[a,b]$. Then the largest negative Stekloff eigenvalue $\lambda^{-}_{1}$ is monotonically increasing on $[a,b]$.

Assume that the largest negative Stekloff eigenvalue $\lambda^{-}_{n}$ is obtained. If the shape of $D$ is known, $\lambda^{-}_{1}$ uniquely determines $n_{c}$ on some suitable interval $[a,b]$ by 2.2. However, it is not true on $\mathbb{R}$ as it is known that different $n_{c}$’s can give the same $\lambda^{-}_{1}$.

3 Reconstruction of Stekloff Eigenvalues

Now we consider IP1 to reconstruct Stekloff eigenvalues from Cauchy data. The main ingredient is the reciprocity gap method using Cauchy data [8, 10, 28]. Assume that $u$ and $\partial_{\nu}u:=\partial u/\partial\nu$ are known on $\Gamma$ for each point source incident wave $u^{i}:=\Phi(\cdot,x_{0}),x_{0}\in C$ (see 1). The following auxiliary scattering problem will be useful in the subsequent analysis. Find $u_{\lambda}(\cdot,x_{0}):=u_{\lambda}^{s}(\cdot,x_{0})+\Phi(\cdot,x_{0})$ such that

[TABLE]

where $\nu$ is the unit outward normal to $\Gamma$ and $\lambda$ is a constant such that $\Im(\lambda)\geq 0$. It is shown in [5] that (3.1) has a unique solution.

Denote by $U$ and $U_{\lambda}$ the sets of solutions $u(x,x_{0})$ to (2.1) and $u_{\lambda}(x,x_{0})$ to (3.1) for the incident wave $\Phi(\cdot,x_{0}),x_{0}\in C$, respectively. Define the reciprocity gap functional by

[TABLE]

where $v_{1}$ and $v_{2}$ are solutions of the Helmholtz equation. Let $\mathbb{S}:=\{d\in\mathbb{R}^{2};|d|=1\}$ and consider the integral equation of finding $g\in L^{2}(\mathbb{S})$ to

[TABLE]

where $\Phi_{z}(\cdot):=\Phi(\cdot,z)$ for some $z\in B$ and $v_{g}$ is the Herglotz wave function defined by

[TABLE]

Lemma 3.1

If $\int_{\Gamma}u_{\lambda}(x,x_{0})f(x)ds(x)=0$ for all $u_{\lambda}\in U_{\lambda}$, then $f(x)=0$ on $\Gamma$.

Proof:

Assume that $f(x)$ satisfies $\int_{\Gamma}u_{\lambda}(x,x_{0})f(x)ds(x)=0$ for all $x_{0}\in C$. Let $\tilde{u}^{s}$ be the solution of the following problem

[TABLE]

Using Green’s representation theorem [9], Green’s second theorem and the boundary condition $\partial_{\nu}u_{\lambda}+\lambda u_{\lambda}=0$ on $\Gamma$ for all $x_{0}\in C$, we have that

[TABLE]

The unique continuity principle implies that $\tilde{u}^{s}(x)=0$ in $\mathbb{R}^{2}\setminus\bar{B}$. By the trace theorem we have $f=0$ on $\Gamma$. $\square$

Theorem 3.1

If $\lambda$ is not a Stekloff eigenvalue of (2.2), then for $u\in U$ and $u_{\lambda}\in U_{\lambda}$, the operator $\mathcal{R}:L^{2}(\mathbb{S})\rightarrow L^{2}(C)$ defined by

[TABLE]

is injective.

Proof:

Let $g$ satisfy $R(u_{\lambda}(\cdot,x_{0})-u(\cdot,x_{0}),v_{g}(\cdot))=0$ for all $x_{0}\in C$. If $g\neq 0$, let $w^{s}$ solve

[TABLE]

Using the boundary condition $\partial_{\nu}u_{\lambda}+\lambda u_{\lambda}=0$ on $\Gamma$ and Green’s second theorem twice, the following holds

[TABLE]

From (3.7) and 3.1, we have $\partial_{\nu}(w^{s}+v_{g})+\lambda(w^{s}+v_{g})=0$ on $\Gamma$. Together with (3.6), $w^{s}+v_{g}$ satisfies

[TABLE]

Since $\lambda$ is not a Stekloff eigenvalue, (3.8) only has the trivial solution $w^{s}+v_{g}=0$ in $B$. From (3.6) and the unique continuity principle, $w^{s}+v_{g}=0$ in $\mathbb{R}^{2}$, i.e., the Herglotz wave function $v_{g}=-w^{s}$ satisfies the radiation condition. This is a contradiction. $\square$

The following theorem is the main result on the reconstruction of Stekloff eigenvalues from Cauchy data.

Theorem 3.2

If $\lambda$ is not a Stekloff eigenvalue of (2.2) and $z\in B$, then there exists a sequence $\{g_{n}\},g_{n}\in L^{2}(\mathbb{S})$, such that

[TABLE]

and $v_{g_{n}}$ converges in $L^{2}(B)$.

• 2.

If $\lambda$ is a Stekloff eigenvalue, then for every sequence $\{g_{n}^{z}\},g_{n}^{z}\in L^{2}(\mathbb{S})$ satisfying

[TABLE]

$\lim\limits_{n\rightarrow\infty}\|v_{g_{n}^{z}}\|_{H^{1}(B)}=\infty$* for almost every $z\in B$.*

Proof:

1. Let $w_{z}$ be the solution of the following problem

[TABLE]

From Lemma 3.1 of [5], we have

[TABLE]

where $w_{z}^{i}$ satisfies the Helmholtz equation in $B$ and $w_{z}^{s}\in H_{loc}^{2}(\mathbb{R}^{2})$ is a radiation solution to

[TABLE]

Due to the denseness property (Theorem 5.21 of [9]), there exists a sequence of Herglotz wave functions $\{v_{g_{n}}\}$ such that

[TABLE]

Next we show that $\{v_{g_{n}}\}$ satisfies \lim\limits_{n\rightarrow\infty}R\big{(}u_{\lambda}-u,v_{g_{n}}\big{)}=R(u_{\lambda},\Phi_{z}). Using Green’s second theorem twice, one has that

[TABLE]

where the last step is due to (3.11) and (3.13).

2. Assume on the contrary that for $z\in B_{\rho}$, where $B_{\rho}\subset B$ is a small ball of radius $\rho$, $\|v_{g_{n}^{z}}\|_{H^{1}(B)}$ is bounded as $n\rightarrow\infty$. Then there exists a subsequence of $v_{g_{n}^{z}}$, still denoted by $v_{g_{n}^{z}}$, converging weakly to a function $v^{i}\in H^{1}(B)$. Then

[TABLE]

Let $w^{s}\in H_{loc}^{2}(\mathbb{R}^{2})$ be a radiating solution to

[TABLE]

From the Green’s second theorem and (3.15), $w:=v^{i}+w^{s}$ satisfies

[TABLE]

From (3.17) and 3.1,

[TABLE]

which, together with (3.16), implies that $w$ satisfies

[TABLE]

From 3.1 and the proof of Theorem 3.3 of [5], (3.18) is solvable if and only if

[TABLE]

for each Stekloff eigenfunction ${\omega}_{\lambda}\in H^{1}(B)$. Since ${\omega}_{\lambda}$ satisfies $\partial_{\nu}{\omega}_{\lambda}+\lambda{\omega}_{\lambda}=0$ on $\Gamma$, (3.19) becomes

[TABLE]

Green’s representation theorem implies that $\omega_{\lambda}(z)=0$ for $z\in B_{\rho}$. The unique continuation principle now implies that the Stekloff eigenfunction $\omega_{\lambda}=0$ in $B$, which is a contradiction. $\square$

Based on 3.2, the following reciprocity gap algorithm can be used to reconstruct (several) Stekloff eigenvalues from Cauchy data.

The RG Algorithm

For a region of interests (e.g., an interval on $\mathbb{R}$ for real Stekloff eigenvalues or a rectangular region on $\mathbb{C}$ for complex Stekloff eigenvalues), generate a grid $T$.

• 2.

For each $\lambda\in T$, solve the scattering problem (3.1) to compute the auxiliary Cauchy data $u_{\lambda}(\cdot,x_{0})$ and $\partial_{\nu}u_{\lambda}(\cdot,x_{0})$ on $\Gamma$.

• 3.

Fix a point $z\in B$, use the Tikhonov regularization to compute an approximate solution $g_{\lambda}\in L^{2}(\mathbb{S})$ to the integral equation

[TABLE]

where $u(\cdot,x_{0})$ is the solution to (2.1) for $\Phi(\cdot,x_{0}),x_{0}\in C$.

• 4.

Choose $\lambda$ as a Stekloff eigenvalue of (2.2) if the norm of $g_{\lambda}$ is significantly larger (see Section 5.1).

Remark 3.3

The constructed solutions to the reciprocity gap equation may not form a divergent Herglotz wave function series. Hence the above numerical algorithm might not be able to construct all the eigenvalues.

4 Reconstruction of the Index of Refraction

The algorithm in the previous section can reconstruct Stekloff eigenvalues using Cauchy data. Given these reconstructed eigenvalues, in this section, we turn to IP2 to estimate the index of refraction. The relation between the index of refraction and Stekloff eigenvalues is complicated and, to a large extend, unknown. Even when $n(x)$ is constant, a single Stekloff eigenvalue cannot uniquely determine it. Note that 2.2 only holds on an appropriate interval.

To this end, we resort to the Bayesian approach, which has been popular for solving inverse problems in recent years [19, 30]. Firstly, the inverse problem is reformulated as a statistical inference for the index of refraction using a Bayes formula. Then the Metropolis-Hastings algorithm is employed to explore the posterior distribution of $n(x)$.

4.1 Bayesian Formulation

Denote by $\mathcal{N}$ the normal distribution and $\mathcal{U}$ the uniform distribution. IP2 can be written as the statistical inference of $n(x)$ such that

[TABLE]

where ${\boldsymbol{\lambda}}\in\mathbb{C}^{m}$ is a vector of (reconstructed) Stekloff eigenvalues, $n(x)$ is a random function, $\mathcal{G}$ is the operator mapping $n(x)$ to ${\boldsymbol{\lambda}}$ based on (2.2), and $E$ is the random noise. The noise $E\sim\mathcal{N}(0,\,\sigma^{2})$, which is modeled as additive and mutually independent of $n(x)$. In the Bayesian framework, the prior information can be coded into the prior density $\pi_{pr}(n)$. For example, if $n$ is known to be a real constant $n_{0}$ such that $a<n_{0}<b$, one may take the prior as the continuous uniform distribution, i.e., $n\sim\mathcal{U}(a,b)$.

Given Stekloff eigenvalues ${\boldsymbol{\lambda}}$, the goal of the Bayesian inverse problem is to seek statistical information of $n(x)$ by exploring the conditional probability distribution $\pi_{post}(n|{\boldsymbol{\lambda}})$, called the posterior distribution of $n$. An important quantity is the conditional mean (CM) of $n$ defined as

[TABLE]

which is an constant estimation of $n(x)$. If $n\sim\mathcal{U}(a,b)$, by the Bayes formula, the posterior distribution satisfies

[TABLE]

i.e.,

[TABLE]

where $I$ is the density function for $\mathcal{U}(a,b)$.

4.2 Markov Chain Monte Carlo Method

To explore $\pi_{post}(n|{\boldsymbol{\lambda}})$ given in (4.4), we employ the popular MCMC (Markov Chain Monte Carlo). MCMC to estimate CM is as follows: design a Markov Chain $\{X_{j}\}_{j=0}^{\infty}$ from required distribution and take the mean of the chain to approximate expectation. In particular, one could estimate ${\mathbb{E}}\{n|{\boldsymbol{\lambda}}\}$ by a sample mean using Monte Carlo integration:

[TABLE]

where $n_{j}$, $j=1,\cdots,m$, are samples drawing from $\pi_{post}(n|{\boldsymbol{\lambda}})$. Two popular methods are Metropolis-Hasting (M-H) algorithm [27, 16] and Gibbs sampler [12]. In this paper, we choose a delayed rejection adaptive Metropolis-Hasting algorithm [13].

The MH Algorithm

Choose initial value $n_{1}\in\mathbb{R}$ and set $j=1$; 2. 2.

Draw a sample $w$ from a proposal distribution

[TABLE]

and compute

[TABLE] 3. 3.

Draw $t\sim\mathcal{U}(0,1)$; 4. 4.

If $\alpha(n_{j},w)\geq t$, set $n_{j}=w$, else $n_{j+1}=n_{j}$. 5. 5.

When $j=K$, the maximum sample size, stop; else, $j\leftarrow j+1$ and go to 2.

4.3 Spectral Indicator Method

In the above algorithm, for each sample $n_{j}$, one needs to compute Stekloff eigenvalues, which are done by the finite element method for (2.2) [33, 25]. Note that the Stekloff eigenvalues are complex if $n(x)$ is a complex function. In addition, the reconstructed eigenvalues only approximate the exact ones and the multiplicities are not known in general. Thus the numerical methods need to compute Stekloff eigenvalues of (2.2) in a region on $\mathbb{C}$ close to the origin.

A recently developed spectral indicator method (SIM) is a good fit for this case [17, 18]. Given reconstructed eigenvalues ${\boldsymbol{\lambda}}$, a rectangular region containing these eigenvalues is chosen. Then SIM is used to compute all eigenvalues inside the region.

5 Numerical Examples

In this section, we present some numerical examples to use the RG method to reconstruct Stekloff eigenvalues from the Cauchy data and estimate the index of refraction using the Bayesian approach. Three scatterers are considered: a disc with radius $1$ centered at the origin, a square with vertices given by

[TABLE]

and an L-shaped domain given by

[TABLE]

Three different indices of refraction are chosen: i) $n(x)=5$, ii) $n(x)=4+2|x|$, and iii) $n(x)=2+4i$.

The synthetic scattering data is simulated by a finite element method with a perfectly matched layer (PML) for (2.1) [7]. The wavenumber is $k=1$. There are $100$ source points uniformly distributed on the curve $C$, a circle with radius $3$. We compute the Cauchy data at $100$ points uniformly distributed on $\Gamma:=\partial B$ (a circle with radius $2$) and add $3\%$ noise.

5.1 Reconstruction of Stekloff eigenvalues

Given Cauchy data, we show that several Stekloff eigenvalues close to the origin can be reconstructed effectively using the RG Algorithm presented in Section 3. For real $n(x)$, all the eigenvalues are real. We choose an interval $[-5,5]$ and use a uniform grid $T$ given by

[TABLE]

For complex $n(x)$, the Skeloff eigenvalues are complex. We choose a domain $[-1,0.5]\times[-0.5,1]$ and

[TABLE]

Using a fixed point $z=(0.2,0.6)$ in $B$, for each $\lambda\in T$, a discretization of (3.20) leads to a linear system

[TABLE]

where $A^{\lambda}$ is a matrix given by

[TABLE]

and $f^{\lambda}$ is a vector given by

[TABLE]

The Tikhonov regularization with the parameter $\alpha=10^{-5}$ is used to compute an approximate solution $g_{\lambda}$ to (5.3):

[TABLE]

In the following examples, we show the plots of $|g_{\lambda}|$. In all the figures, the crosses are the exact Stekloff eigenvalues computed using a finite element method [33, 25]. The results indicate that Stekloff eigenvalues close to the origin can be determined effectively.

Example 1 Real index of refraction $n(x)=5$. The exact and reconstructed eigenvalues are shown in 1. The plots of $|g_{\lambda}|$ for three domains are shown in 2.

Example 2 Real index of refraction $n(x)=4+2|x|$. The exact and reconstructed eigenvalues are shown in 2. The plots of $|g_{\lambda}|$ for three domains are shown in 3.

Example 3 Complex index of refraction $n(x)=2+4i$. The exact and reconstructed eigenvalues are shown in 3. The plots of $|g_{\lambda}|$ for three domains are shown in 4.

5.2 Estimation of the index of refraction

Given the reconstructed Stekloff eigenvalues, we present some numerical examples for the estimation of the index of refraction by the Bayesian approach. The examples are rather naive. Nonetheless, the results show the potential of statistical approaches for inverse scattering problems. Since the main goal is to show the effectiveness of the Bayesian approach, we assume that the shape of the scatterer is known in the following examples.

Example 4 Real constant index of refraction $n(x)=5$. Assume one Stekloff eigenvalue is reconstructed from Cauchy data: $-0.48$ for the disk, $-0.54$ for the square, and $-0.52$ for the L-shaped domain. Since $n$ is a real constant, we take a uniform prior $\mathcal{U}(0,8)$. The posterior density is given by

[TABLE]

We generate 3000 samples for each domain. The initial sample is chosen to be $n_{1}=2$. The rest samples $\{n_{i}\}_{i=2}^{3000}$ are drawn from the symmetric proposal distribution

[TABLE]

where $\gamma^{2}=2.4^{2}/2$.

4 shows $n_{CM}$ for three domains. The Markov chains are shown in 5. The samples concentrate around $n=5$ for the unit circle domain and square domain. However, for the L-shaped domain, the samples are accumulated around two values, $5$ and $7$. In fact, this implies that one Stekloff eigenvalue cannot uniquely determine the constant index of refraction. If two Stekloff eigenvalues, $0.70$ and $-0.52$, are used (see 6), we obtain $n_{CM}=5.0074$, which is a good approximation of $5$.

Example 5 Real function index of refraction $n(x)=4+2|x|$. Assume that two Stekloff eigenvalues are given and $n(x)$ is of the form $\beta^{1}+\beta^{2}|x|$. We first obtain a constant approximation $n_{0}$ for $n(x)$ as the above example. This provides some ideas of how to choose priors for $\beta^{1}$ and $\beta^{2}$. For the second step, two Stekloff eigenvalues are used. The posterior distribution is given by

[TABLE]

Two reconstructed Stekloff eigenvalues from 2 are used for each domain: $2.10,-0.48$ for disc, $0.42,-0.58$ for the square, and $1.00,-0.50$ for the L-shaped domain. 5 shows the reconstruction results.

Example 6 Complex index of refraction $n(x)=2+4i$. Assume that $\Re(n)\sim\mathcal{U}(0,8)$ and $\Im(n)\sim\mathcal{U}(0,8)$. The same proposal distribution $q(n_{j},n_{j-1})$ are used to sample both $\Re(n)$ and $\Im(n)$. We use Stekloff eigenvalues, $-0.02+0.26i$ for the circle, $-0.64+0.02i,0.12+0.56i$ for the square, and $-0.1+0.52i,0.04+0.22i$ for the L-shaped domain. 6 shows the reconstruction results and 7 shows the Markov chains.

6 Conclusions

In this paper, we show that the Cauchy data of a medium scattering problem can be used to reconstruct the related Stekloff eigenvalues, which is demonstrated by numerical examples. A Bayesian approach is proposed to estimate the index of refraction using the reconstructed Stkeloff eigenvalues. The method is particularly meaningful when there is a lack of understanding of the relation between the known data and the unknown quantities. In future, the authors plan to extend the preliminary study here to more challenging inverse scattering problems.

Acknowledgement

The research was supported in part by Guangdong Natural Science Foundation of China (2016A030313074) and NNSF of China under grants 11801218 and 11771068.