Monotonicity in inverse obstacle scattering on unbounded domains

We consider an inverse obstacle scattering problem for the Helmholtz equation with obstacles that carry mixed Dirichlet and Neumann boundary conditions. We discuss far field operators that map superpositions of plane wave incident fields to far field patterns of scattered waves, and we derive monotonicity relations for the eigenvalues of suitable modifications of these operators. These monotonicity relations are then used to establish a novel characterization of the support of mixed obstacles in terms of the corresponding far field operators. We apply this characterization in reconstruction schemes for shape detection and object classification, and we present numerical results to illustrate our theoretical findings.


Introduction
We discuss an inverse obstacle scattering problem for time-harmonic scalar waves governed by the Helmholtz equation. The goal is to recover the position and the shape of a collection of compactly supported scattering objects from far field observations of scattered waves. We consider impenetrable obstacles with mixed Dirichlet and Neumann boundary conditions, i.e., we assume that the scatterers D = D 1 ∪ D 2 consist of two components such that D 1 ∩ D 2 = ∅, where ∂D 1 carries a Dirichlet boundary condition while ∂D 2 carries a Neumann boundary condition. The Dirichlet part D 1 and the Neumann part D 2 of the scattering objects might consist of 2 Author to whom any correspondence should be addressed.
Original content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI. several connected components, and we do neither assume that the number of connected components nor whether they carry Dirichlet or Neumann boundary conditions are known a priori. Accordingly, qualitative reconstruction schemes (see, e.g., [3,6,8,9,30,31,35]), which do not make use of topological or physical properties of the scattering objects, are a natural choice. In addition to shape reconstruction, we will also show that the type of boundary condition on each connected component of the obstacle can be classified from scattering data, i.e., we show that the Dirichlet part D 1 and the Neumann part D 2 can be recovered separately.
Among qualitative methods for shape reconstruction, the linear sampling method has been successfully applied to inverse mixed obstacle scattering problems (see, e.g., [3][4][5]). The factorization method has been justified under the additional assumption that the Dirichlet part D 1 and the Neumann part D 2 of the scattering object can be separated a priori (see [18,19,31]). We build on and extend ideas from these works to develop a monotonicity based qualitative shape reconstruction technique. This monotonicity method is formulated in terms of far field operators that map superpositions of incident plane waves, which are being scattered at the unknown scattering objects, to the far field patterns of the corresponding scattered waves. It exploits monotonicity properties of the eigenvalues of suitable modifications of these operators. The main result of this work is a rigorous characterization of the support of mixed scattering obstacles in terms of the corresponding far field operators without any additional a priori information. This is a significant extension of the results in [18,19,31].
The monotonicity based approach to shape reconstruction has originally been developed for the inverse conductivity problem in [14,26], extending an earlier monotonicity based reconstruction scheme developed in [37]. The method is related to monotonicity principles for the Laplace equation established in [28,29]. It has been further developed in [22,23,27], its numerical implementation has been studied in [11][12][13], and recently an extension to impenetrable conductivity inclusions has been established in [7].
The analysis in [26] has been extended to inverse coefficient problems for the Helmholtz equation on bounded domains in [24,25], and in [16] the approach has been generalized to the inverse medium scattering problem on unbounded domains with plane wave incident fields and far field observations of scattered waves. An application of the monotonicity method to an inverse crack detection problem for the Helmholtz equation has recently been considered in [10]. For further recent contributions on monotonicity based reconstruction methods for various inverse problems for partial differential equations we refer to [1,2,20,21,32,36,38].
The main idea of the monotonicity method for inverse mixed obstacle scattering that we discuss in this work, is to compare the real part of the given (or observed) far field operator corresponding to the unknown scattering obstacles to various virtual (or simulated) probing operators corresponding to certain probing domains B. We show that suitable linear combinations of these operators are positive definite up to a finite dimensional subspace if and only if the probing domains are contained inside the support of the scattering objects (i.e., B ⊆ D), or if and only if the probing domains contain the unknown scattering object (i.e., D ⊆ B). Sampling the region of interest using sufficiently many suitable probing domains these observations can be translated into criteria and algorithms for shape reconstruction. A major difference to the linear sampling method or the factorization method is that the sampling strategy involves probing domains instead of test points, and that the shape characterization is formulated in terms of probing operators instead of test functions.
This article is organized as follows. In section 2 we briefly recall the mathematical formulation of the mixed obstacle scattering problem, and in section 3 we discuss a factorization of the corresponding far field operator from [18,19,31]. In section 4 we establish the existence of localized wave functions for the mixed obstacle scattering problem, and in section 5 we use these localized wave functions to prove a rigorous characterization of the support of scattering obstacles in terms of the far field operator. We discuss numerical algorithms based on these theoretical results in section 6, and we close with some concluding remarks.

Scattering by impenetrable obstacles
We consider the scattering of time-harmonic scalar waves in an unbounded homogeneous background medium by a collection of impenetrable obstacles carrying Dirichlet and Neumann boundary conditions. Suppose that 3, is open and Lipschitz bounded with connected complement R d \D such that D 1 ∩ D 2 = ∅. The subsets D 1 and D 2 may consist of finitely many connected components. Below we will impose Dirichlet boundary conditions on D 1 and Neumann boundary conditions on D 2 , and thus we refer to D 1 and D 2 as the Dirichlet and Neumann obstacles, respectively.
We assume that the wave motion is caused by an incident field u i ∈ H 1 loc (R d ) satisfying the Helmholtz equation with wave number k > 0 that is being scattered at the obstacle D. The scattered field u s ∈ H 1 loc (R d ) satisfies and the boundary conditions together with the Sommerfeld radiation condition uniformly with respect to all directions x : = x/|x| ∈ S d−1 . Throughout, the Helmholtz equation is to be understood in weak sense, but standard interior regularity results yield smoothness of u s in R d \D. In particular the Sommerfeld radiation condition (2.2c) is well defined. As usual, we call a (weak) solution to the Helmholtz equation on an unbounded domain that satisfies the Sommerfeld radiation condition uniformly with respect to all directions a radiating solution.
Furthermore, the solution has the asymptotic behavior

4)
and w ∞ ∈ L 2 (S d−1 ) is called the far field pattern of w.
Choosing f = −u i | ∂D 1 and g = −∂u i /∂ν| ∂D 2 in lemma 2.1 proves the existence and uniqueness of solutions to the scattering problem (2.2). For the special case of a plane wave incident field u i (x; θ) := e ikx·θ , x ∈ R d , we explicitly indicate the dependence on the incident direction θ ∈ S d−1 by a second argument, and accordingly we write u s (·; θ), and u ∞ (·; θ) for the corresponding scattered field and its far field pattern, respectively.
We define the far field operator and we note that F mix D is compact and normal (see, e.g., [31, theorem 3.3]). Moreover, the scattering operator is defined by where C d is again the constant from (2.4). The operator S mix D is unitary, and consequently the eigenvalues of F mix D lie on the circle of radius 1/(2k|C d | 2 ) centered in i/(2k|C d | 2 ) in the complex plane (cf, e.g., [31, theorem 3.3]

Factorizations of the far field operator
Next we briefly recall three factorizations of the far field operators F mix D , F dir D 1 , and F neu D 2 , which have been used in the traditional factorization method, and that will be applied to develop the monotonicity based shape characterization in section 5 below. As usual, the single layer operator is defined by (3.1) and the normal derivative of the double layer potential is given by Here, Φ k denotes the fundamental solution to the Helmholtz equation in R d .

Dirichlet or Neumann obstacles
The first result describes the factorization of the far field operator for Dirichlet obstacles. (a) The far field operator F dir Proof. This is shown in [31, lemmas 1.13-1.14 and theorem 1.15].

Remark 3.3.
An immediate consequence of theorem 3.2 is that the real part 3 is a compact perturbation of a self-adjoint and coercive operator. This implies that Re(F dir D 1 ) has only finitely many positive eigenvalues. In theorem 5.3 below we will significantly refine and extend this observation. ♦ Next we consider the factorization of the far field operator for Neumann obstacles.

Theorem 3.4.
(a) The far field operator F neu to the far field pattern w ∞ of the unique radiating solution to the exterior Neumann boundary value problem (2.3a) and (2.3c). (b) G neu D 2 is compact and one-to-one with dense range in L 2 (S d−1 ). (c) If k 2 is not a Neumann eigenvalue of −Δ in D 2 , then N D 2 is an isomorphism. 3 As usual the real part of a linear operator A : X → X on a Hilbert space X is the self-adjoint operator given by Re(A) := Proof. This is shown in [31, theorem 1.26].
Remark 3.5. An immediate consequence of theorem 3.4 is that the real part is a compact perturbation of a self-adjoint and coercive operator. This implies that Re(F neu D 2 ) has only finitely many negative eigenvalues. In theorem 5.4 below we will significantly refine and extend this observation. ♦

Mixed obstacles
In the mixed case the obstacle D = D 1 ∪ D 2 consists of two bounded components and carries Dirichlet boundary conditions on ∂D 1 and Neumann boundary conditions on ∂D 2 .
(a) The far field operator F mix D : to the far field pattern w ∞ of the unique radiating solution to the exterior mixed boundary value problem (2.3a)-(2.3c). The operator T mix D : D is compact and one-to-one with dense range in L 2 (S d−1 ). Proof. This is shown in [ (3.4) and From the asymptotic behavior of the fundamental solution Φ k we obtain that H * B φ, φ ∈ H − 1 2 (∂B), is just the far field pattern of the single layer potential We will use the relation (see [31, p 18]).
The operators H * B H B corresponding to certain probing domains B will be utilized as probing operators in the shape characterization results that are developed in section 5 below.

Localized wave functions
In this section we establish the existence of localized wave functions. These are pairs of certain wave functions such that one component has arbitrarily large norm on some prescribed boundary while the other component has arbitrarily small norm on some different boundary. These localized wave functions will be essential in the proof of the monotonicity based shape characterization in section 5 below.

Dirichlet or Neumann obstacles
To start with, we consider the case when either only Dirichlet or Neumann obstacles are present.
Let B ⊆ R d be open and Lipschitz bounded, and let Γ ⊆ ∂B be relatively open. We define the restriction operator and we note that the adjoint operator satisfies Here, H − 1 2 (Γ) denotes the dual space of H 1 2 (Γ) (see, e.g., [33, p 99]). Accordingly, we introduce Since (4.1) remains true if we modify ∂B away from Γ, we can w.l.o.g. assume that k 2 is not a Dirichlet eigenvalue of −Δ in B. Then S B and G dir B are injective (cf theorem 3.2(b) and (c)), and since R(R * Γ ) has infinite dimensional range, this shows that R(H * Γ ) is infinite dimensional as well.
The proof of theorem 4.1 relies on the following lemmas.
respectively. Rellich's lemma (cf, e.g., [9, theorem 2.14]) and unique continuation guarantee Then w is an entire radiating solution to the Helmholtz equation, and thus w = 0 in R d . This shows that h = w ∞ 1 = 0. In the next lemma we quote a special case of lemma 2.5 in [14].

Lemma 4.3. Let X, Y and Z be Hilbert spaces, and let A : X → Y and B : X → Z be bounded linear operators. Then,
Now we give the proof of theorem 4.1.
Proof of theorem 4.1. Let D 2 = ∅, and let B, D 1 ⊆ R d be open and Lipschitz bounded such that R d \D 1 is connected, and suppose that B D 1 . Let V ⊆ L 2 (S d−1 ) be a finite dimensional subspace. We denote by P V : and we have seen before that R(H * Γ ) is infinite dimensional. Using a simple dimensionality argument (see [25, lemma 4.7]) it follows that Accordingly, lemma 4.3 implies that there is no constant C > 0 such that , which ends the proof.
Using similar arguments the following result for Neumann obstacles can be shown.

Mixed obstacles
For the general mixed case, i.e., when both Dirichlet and Neumann obstacles are present, we require a refined version of the theorems 4.1 and 4.4, which we call simultaneously localized wave functions (see also [21], where a similar construction has been used).
To begin with, we define additional restriction operators. Let Then the adjoint operators satisfy

Furthermore, given an open and Lipschitz bounded
We note that the adjoint operator satisfies Here, H The proof of theorem 4.5 relies on the following lemma. Accordingly, are radiating solutions to Rellich's lemma and unique continuation guarantee that . Therefore, Since Γ is C 1 smooth, 4 this and the smoothness of SL B f B away from ∂B imply that f Γ ∈ C 1 (Γ). Without loss of generality we assume that Denoting by X ⊆ H 1 2 (Γ) the subspace of piecewise linear continuous functions on Γ that vanish on ∂Γ, 5 Since X is infinite dimensional and G mix D R * D 1 R Γ * is one-to-one (see theorem 3.6(b)), we find that Z is infinite dimensional as well.
Now we give the proof of theorem 4.5.
Since D 1 B and ∂D 1 is piecewise C 1 smooth, there exists Γ ⊆ ∂D 1 \B relatively open such that Γ is C 1 smooth. Combining lemma 4.6 with a simple dimensionality argument (see [25, lemma 4.7]) we find that denotes the subspace in lemma 4.6, and thus Accordingly, lemma 4.3 implies that there is no constant C > 0 such that Therefore, there exists a sequence ( ψ m ) m∈N ⊆ L 2 (S d−1 ) satisfying as m → ∞. We define ψ m := ψ m − P V ψ m ⊆ V ⊥ for any m ∈ N to obtain 5 The latter condition ensures that X does not contain any smooth functions except for zero.
This ends the proof.
The following result can be shown proceeding similarly to the proof of theorem 4.5.

Monotonicity based shape reconstruction
Using the localized wave functions developed in the previous section we will establish monotonicity relations for far field operators in terms of the following extension of the Loewner order to compact self-adjoint operators, which has been introduced in [25]. Let A 1 , A 2 : X → X be compact self-adjoint operators on a Hilbert space X, and let r ∈ N. We write

Dirichlet or Neumann obstacles
In the following we consider the case when either only Dirichlet or Neumann obstacles are present. We discuss criteria to characterize the support of an unknown scattering obstacle D in terms of the corresponding far field operator and a probing operator H * B H B corresponding to a probing domain B. These criteria established in theorems 5.3 and 5.4 below describe whether the probing domain B is contained in the support D of the scattering obstacles or not. To begin with, we discuss the case when only Dirichlet obstacles are present. Proof.
(a) Let B ⊆ D 1 . We define P dir is the unique radiating solution to the exterior Dirichlet boundary value problem (2.3a) and (2.3b) with D 1 replaced by B (and D 2 = ∅). Then P dir is a compact linear operator by standard interior regularity results, and the uniqueness of solutions to the exterior Dirichlet boundary value problem (see lemma 2.1) implies that G dir

Recalling (3.5), this shows that
Substituting the factorization (3.3) gives Using theorem 3.2(d) and (e) we find that 1 2 (S D 1 + S * D 1 ) is a compact perturbation of the self-adjoint and coercive operator S D 1 ,i , i.e., with some compact self-adjoint operator K. Accordingly, , where c 1 denotes the coercivity constant of S D 1 ,i (see theorem 3.2(d)). We define the subspace associated to an eigenvalue larger than c 1 .
The spectral theorem for compact self-adjoint operators shows that V is finite dimensional. Accordingly, (b) Let B D 1 . We suppose that there exists a finite dimensional subspace V ⊆ L 2 (S d−1 ) such that Again, using the factorization (3.3) we find that, for all ψ ⊆ V ⊥ , for some C > 0, and on the other hand .
Substituting this into (5.1) and applying theorem 4.1 gives a contradiction.
The following result for Neumann obstacles can be shown using similar arguments as in the proof of theorem 5.3.

Mixed obstacles
Next we consider the general mixed case, i.e., when both Dirichlet and Neumann obstacles are present. While the criteria developed in theorems 5.3 and 5.4 determine whether a certain probing domain B is contained in the support D of the scattering obstacles or not, the criterion for the mixed case established in theorem 5.5 below characterizes whether a certain probing domain B contains the support D of the scattering obstacles or not.
Remark 5.6. The results in theorem 5.5 remain true in the special case, when D 2 = ∅ and F mix D = F dir D 1 , and also in the special case, when D 1 = ∅ and F mix D = F neu D 2 . The corresponding shape characterizations complement the results established in theorems 5.3 and 5.4. ♦

Proof.
(a) Let D 1 ⊆ B. It has been shown in [31, lemma 3.5] that with some compact self-adjoint operator K, i.e., Re(F mix D ) + H * B H B is a compact perturbation of a self-adjoint and coercive operator. Similar to the proof of theorem 5.3(a) this implies (a). (b) Let R d \(B ∪ D 1 ∪ D 2 ) be connected and let ∂D 1 be piecewise C 1 smooth. We suppose that there exists a finite dimensional subspace Combining theorem 3.6(a) with theorems 3.2 and 3.4(d) and (e), we find that with some compact self-adjoint operator K. Accordingly, we define the subspace associated to an eigenvalue with absolute value larger than c 1 2 , where c 1 denotes the coercivity constant of S D 1 ,i (see theorem 3.2(d)). The spectral theorem for compact self-adjoint operators shows that V 2 ⊆ L 2 (S d−1 ) is finite dimensional, and thus Using the restriction operators R D 1 and R D 2 from (4.2), we find that, for all ψ ∈ (V 1 + V 2 ) ⊥ , 0 Re (d) This can be shown proceeding similarly to the proof of part (b), and using theorem 4.7 to obtain a contradiction.

Numerical examples
We now work toward numerical implementations of the shape characterizations developed in section 5. The main issue here is that numerical approximations of the operators F dir D 1 , F neu D 2 , F mix D , and H B are necessarily finite dimensional. Accordingly, the question, whether suitable combinations of these operators as considered in theorems 5.3-5.5 are positive definite up to some finite dimensional subspace, needs to be carefully relaxed to obtain reliable numerical algorithms. We present some preliminary ideas in this direction, restricting the discussion to the two-dimensional case.

An explicit radially symmetric example
We illustrate the shape characterization results from theorems 5.3 and 5.5 for the special case of a single radially symmetric Dirichlet obstacle by an explicit example. Let D 2 = ∅, and let D 1 = B r (0) ⊆ R 2 be the disk of radius r > 0 centered at the origin.
We first derive series expansions for the incident and scattered fields and use them to compute the eigenvalue value decomposition of the far field operator F dir D 1 . The Jacobi-Anger expansion (see, e.g., [9, (3.89)]) shows that for each incident direction θ = (cost, sint) ∈ S 1 the incident field satisfies A short calculation yields that the scattered field is given by Substituting the asymptotic behavior of the Hankel functions (see, e.g., [9, (3.82)]) into this expansion we find that the far field pattern of u s is Let g ∈ L 2 (S 1 ) with Fourier expansion g(θ) = m∈Z g m e imt , θ = (cost, sint) ∈ S 1 . Then the far field operator F dir D 1 : L 2 (S 1 ) → L 2 (S 1 ) from (2.5) satisfies x = (cos φ x , sin φ x ) ∈ S 1 . Accordingly, the eigenvalues and eigenvectors of F dir D 1 are given by (λ (r) Here we used the integral representation of J 0 (see, e.g., [34, 10.9.2]). Writing φ = (cosτ , sinτ ) and substituting the Fourier expansion of g we find that (see, e.g., [9, (3.88)]). Accordingly, the eigenvalues and eigenvectors of the operator H * B H B are given by (μ (R) n , v n ) n∈Z with From (6.1) and (6.3) we conclude that in the special case, when D 1 = B r (0) and B = B R (0), the eigenvalues and eigenvectors of Re(F dir Using the criteria established in theorems 5.3 and 5.5(a) and (b) we obtain that (a) If R < r, then Re(F dir D 1 ) + H * B H B has only finitely many positive but infinitely many negative eigenvalues, and (b) If R > r, then Re(F dir D 1 ) + H * B H B has only finitely many negative but infinitely many positive eigenvalues.
We illustrate how this can be utilized to reconstruct the radius of the scatterer D 1 = B r (0) from observations of F dir D 1 by a numerical example. We evaluate the eigenvalues Re(λ (r) n ), μ (R) n , and Re(λ (r) n ) + μ (R) n with wave number k = 1, radius of the obstacle r = 4, and n = 0, . . . , 1000 for different values of R ∈ [0, 40] in Matlab using the explicit formulas given in (6.1), (6.  Accordingly, both sequences (Re(λ (r) n )) n∈Z and (μ (R) n ) n∈Z decay rapidly once the value of |n| is sufficiently large. Since eigenvalues below some threshold are rounded to zero in Matlab, and since the eigenvalues μ (R) n are on average increasing with respect to R, this explains the increasing but somewhat low numbers of positive eigenvalues of H * B H B in the left plot in figure 1. A similar reasoning explains the seemingly low numbers of negative eigenvalues of Re(F dir D 1 ) in the right plot in figure 1.
In practice the far field data will usually be corrupted by measurement errors, and it will not be possible to compute the eigenvalues of Re(F dir D 1 ) + H * B H B with very high precision, as done in this example so far. To see how this influences the numerical results, we repeat the previous computations but consider only those eigenvalues with absolute value larger than a threshold δ = 0.01. For comparison, we note that the eigenvalue of largest magnitude of Re(F dir D 1 ) in this example is Re(λ (r) 2 ) ≈ 11.03. In figure 2 we show plots of the number of positive eigenvalues Re(λ (r) n ) (dotted), μ (R) n (dashed), and Re(λ (r) n ) + μ (R) n (solid) within the range n = 0, . . . , 1000 that are larger than δ (left plot) and of the number of negative eigenvalues that are smaller than −δ (right plot) as a function of R. The transition in the behavior of the eigenvalues of Re(F dir D 1 ) + H * B H B at R = r = 4 is not nearly as pronounced as before. The reason for this behavior is the rapid decay of the sequences of eigenvalues for larger values of |n|, which implies that only few eigenvalues remain above the threshold δ. However, a rough estimate of r would still be possible by visual inspection of these plots, in particular from the plot on the right-hand side of figure 2.

A sampling strategy for Dirichlet or Neumann obstacles
In the special case, when only Dirichlet obstacles are present, the number of positive eigenvalues of Re(F dir D 1 ) + H * B H B can be utilized to decide whether a probing domain B ⊆ R 2 is contained in the support of the scatterer D 1 or not. We discuss this approach in the following and comment on the special case, when only Neumann obstacles are present at the end of this section.
Let D 2 = ∅, and let D 1 ⊆ R 2 be open and Lipschitz bounded. We assume that far field observations u ∞ ( x l ; θ m ) are available for N equidistant observation and incident directions x l , θ m ∈ {(cos φ n , sin φ n ) ∈ S 1 | φ n = (n − 1)2π/N , n = 0, . . . , N − 1}, (6.5) 1 l, m N. Accordingly, applying the trapezoid rule to (2.5), we find that the matrix approximates the far field operator F dir D 1 . Assuming that the support of the scatterer D 1 is contained in the disk B R (0) for some R > 0, we require N 2kR, (6.7) where as before k denotes the wave number, to fully resolve the relevant information contained in the far field patterns (see, e.g., [17]). We consider an equidistant grid of points with step size h = R/J in the region of interest [−R, R] 2 . For each z ij ∈ Δ we consider a probing operator H * B i j H B i j with B ij = B h/2 (z ij ). Applying the trapezoid rule to (6.2) we find that for each z ij ∈ Δ this operator is approximated by the matrix Therewith, we compute the eigenvalues λ (i j) 1 , . . . , λ (i j) N ∈ R of the self-adjoint matrix For numerical stabilization, we discard eigenvalues with absolute values smaller than some threshold. This threshold depends on the quality of the given far field data. If there are good reasons to believe that A dir B i j is known up to a perturbation of size δ > 0 with respect to the spectral norm, then we can only trust in those eigenvalues with magnitude larger than δ (see, e.g., [15, theorem 7.2.2]). To obtain a reasonable estimate for δ, we use the magnitude of the non-normal part of F dir D 1 , i.e., we take δ = (F dir since this quantity should be zero for exact data and be of the order of the data error, otherwise. The characterization of the support of Dirichlet obstacles in theorem 5.3 suggests that we count for each sampling point z ij ∈ Δ the number of positive eigenvalues of A dir B i j . We define the indicator function I dir : Δ → N, Example 6.1. We consider a single Dirichlet obstacle that has the shape of a kite as sketched in figure 3 (left), and simulate the corresponding far field matrix F dir D 1 ∈ C N×N for N observation and incident directions as in (6.5) using a Nyström method for a boundary integral formulation of the scattering problem (2.2) for two different wave numbers k = 1 (with N = 32), and for k = 5 (with N = 128) in accordance with the sampling condition (6.7).
In figure 3 we show color coded plots of the indicator function I dir from (6.11) with threshold parameter δ = 10 −14 (i.e., the number of positive eigenvalues of the matrix A dir B i j from (6.10) that are larger than δ = 10 −14 evaluated at each grid point z ij ∈ Δ) in the region of interest [−10, 10] 2 ⊆ R 2 . The sampling grid Δ from (6.8) consists of 2J + 1 = 201 grid points in each direction.
The number of positive eigenvalues of the matrix A dir B i j increases with increasing wave number, and it is larger at test points z ij sufficiently far away from the support of the scatterers than at test points z ij inside, as suggested by theorem 5.3. The lower value always coincides with the number of positive eigenvalues of the real part Re(F dir D 1 ) of the far field matrix from (6.6) that are larger than the threshold δ. The total number of eigenvalues of A dir B i j , j = 1, . . . , J, whose absolute values are larger than δ is approximately (on average over all grid points) 24 (for k = 1) and 50 (for k = 5).
Depending on the wave number, the lowest level set of the indicator function I dir nicely approximates the support of the scatterer.
The reconstruction algorithm is rather sensitive to noise in the far field data. To see this, we repeat the previous computation but add 0.1% complex-valued uniformly distributed error to the far field matrix F dir D 1 that we simulate using the Nyström method. We estimate the nonnormality error of the corresponding scattering operator and accordingly we choose δ = 0.1 for the threshold in the reconstruction algorithm. In figure 4, we show color coded plots of the indicator function I dir from (6.11) for wave numbers k = 1 (with N = 32), and for k = 5 (with N = 128). The total number of eigenvalues of A dir B i j , j = 1, . . . , J, whose absolute values are larger than δ is approximately (on average over all grid points) 8 (for k = 1) and 22 (for k = 5).
The reconstruction for k = 5 is better than the reconstruction for k = 1 because more eigenvectors are stably propagated into the far field for larger wave numbers (the number of eigenvalues with absolute values above the threshold δ increases with k). However, despite the low noise level, the shape of the obstacle is not reconstructed very well. ♦ If only Neumann obstacles are present, i.e., D 1 = ∅ and D 2 ⊆ R 2 is open and Lipschitz bounded, then we use the corresponding far field matrix F neu D 2 ∈ C N×N as in (6.6) and the matrix T B i j ∈ C N×N from (6.9) to compute for each sampling point z ij ∈ Δ the eigenvalues The characterization of the support of Neumann obstacles in theorem 5.4 suggests that we count for each sampling point z ij ∈ Δ the number of positive eigenvalues of A neu B i j . We define the indicator function I neu : Δ → N, Theorem 5.4 suggests that I neu admits smaller values at test points z ij inside the obstacle than for test points outside of D 2 .
Example 6.2. In the second example, we consider a Neumann obstacle that has the shape of a peanut as sketched in figure 5. We simulate the corresponding far field matrix F neu D 2 ∈ C N×N for N observation and incident directions using a Nyström method for a boundary integral formulation of the scattering problem (2.2) for two different wave numbers k = 1 (with N = 32), and for k = 5 (with N = 128).
In figure 5 we show color coded plots of the indicator function I neu from (6.12) with threshold parameter δ = 10 −14 in the region of interest [−10, 10] 2 ⊆ R 2 . The equidistant rectangular sampling grid on the region of interest from (6.8) consists of 201 grid points in each direction.
Again, the number of positive eigenvalues of the matrix A B i j increases with increasing wave number, and it is larger at test points z ij sufficiently far away from the support of the scatterers than at test points z ij inside, in compliance with theorem 5.4. The lower value always coincides with the number of negative eigenvalues of the matrix Re(F neu D 2 ) that are smaller than the threshold −δ = −10 −14 . The number of eigenvalues of A neu B i j , j = 1, . . . , J, whose absolute values are larger than δ = 10 −14 is approximately (on average over all grid points) 25 (for k = 1) and 55 (for k = 5). Depending on the wave number the support of the indicator function I neu approximates the support of the scatterer rather well. ♦

Separating mixed obstacles
We return to the general mixed case, i.e., when both Dirichlet and Neumann obstacles are present. While the algorithm developed for Dirichlet or Neumann obstacles in the previous subsection determines whether a sufficiently small probing domain B is contained inside the support of the unknown scattering obstacle D or not, the shape characterization for mixed obstacles established in theorem 5.5 describes whether a sufficiently large probing domain B contains the support D of the scattering objects or not. A corresponding numerical algorithm that implements a similar criterion for the inverse conductivity problem has recently been proposed in [13]. However, since in contrast to the inverse conductivity problem, the monotonicity relations in theorem 5.5 only hold up to certain finite dimensional subspaces of unknown dimension, an extension of the reconstruction algorithm from [13] to the mixed inverse obstacle problem is not straightforward. The numerical implementation of the shape characterization for mixed obstacles from theorem 5.5 is much more delicate than the implementation for the pure Dirichlet or Neumann obstacles considered in the previous subsection. In the following we consider a reduced problem, and utilize theorem 5.5 to develop an algorithm to recover the convex hulls of the Dirichlet obstacle D 1 and of the Neumann obstacle D 2 separately. This is clearly less than full shape reconstruction, but on the other hand our results confirm that we can separate Dirichlet and Neumann obstacles from far field data, at least when their convex hulls do not overlap. We treat the Dirichlet part first, and comment on the Neumann part below. The idea is to consider a sufficiently large number of probing disks B = B R (z) ⊆ R 2 , where for each center z ∈ R 2 the radius R > 0 is chosen as small as possible but such that B still completely covers D 1 . Intersecting those disks then gives an approximation of the convex hull of D 1 . To determine the optimal radius R for each of these disks, we use theorem 5.5, which (under some additional assumptions) says that Example 6.3. We consider a kite-shaped Dirichlet obstacle and a peanut-shaped Neumann obstacle as shown in figure 6 (left). We simulate the corresponding far field matrix F mix D ∈ C N×N analogous to (6.6) for wave number k = 1 and N = 64 observation and incident directions using a Nyström method for a boundary integral formulation of the mixed scattering problem (2.2).
To begin with, we fix the center z = (15, 0) of a single probing disk B = B R (z) and evaluate the matrix on a whole interval of radii 0 < R < 40. Here the matrix T B R (z) ∈ C N×N is defined analogous to (6.9). As in our previous examples we choose a threshold parameter δ = 10 −14 , and in figure 6 (right) we show the number of negative eigenvalues of A mix,+ B R (z) that are smaller than −δ as a function of the radius R (solid).
We observe a similar behavior as for the concentric disks studied in section 6.1 (cf the plots on the right-hand side of figures 1 and 2). The number of negative eigenvalues of A mix,+ B R (z) decreases with increasing R until it becomes stationary up to small oscillations around R ≈ 18. Our theoretical results suggest that radius R, where this transition from decreasing to almost stationary appears, corresponds to the radius of the smallest disk that still completely covers the Dirichlet obstacle D 1 .
To evaluate this transition numerically, we fit a smoothing spline curve through the number of negative eigenvalues of A mix,+ B R (z) as shown on the right-hand side of figure 6 (dotted). We determine the point of maximum signed curvature of this smoothing spline and use the corresponding value of R as approximation of the radius of the smallest disk around z that still completely covers the Dirichlet obstacle D 1 . On the right-hand side of figure 6  we can proceed as in example 6.3 for the Dirichlet obstacle to determine minimal radii of probing disks B R (z) containing the Neumann obstacle.
Example 6.4. We continue with example 6.3 and pick 16 evenly spaced points z 1 , . . . , z 16 on a circle of radius 15 around the origin, which are shown as solid pluses in the two plots on the right-hand side of figure 7. The points z , = 1, . . . , 16, are the centers of 16 probing disks that are used to approximate the convex hulls of the Dirichlet obstacle and of the Neumann obstacle separately. For each center z we estimate the radii R dir and R neu of the smallest disks B R dir (z ) and B R neu (z ) centered at z that completely cover the Dirichlet obstacle and the Neumann obstacle, respectively. These estimates are obtained as described in example 6.3. Therewith we compute approximations Reasonable approximations of the convex hulls of the Dirichlet obstacle and of the Neumann obstacle are recovered by the algorithm. These reconstructions could, e.g., be used as a priori information that is required in the factorization method for mixed obstacles from [18,19,31] to obtain an improved shape reconstruction. Different configurations for the centers of the probing disks that are used to approximate the convex hulls (i.e., of the solid pluses in the two plots on the right-hand side of figure 7) are possible. However, arranging them on a circle around the region of interest worked well in all examples that we considered so far. ♦ Example 6.5. We discuss another example with two Dirichlet obstacles (kite-shaped and peanut-shaped) and one Neumann obstacle (an ellipse) as shown in figure 8 (left). We simulate the corresponding far field matrix F mix D ∈ C N×N for wave number k = 1 and N = 64 observation and incident directions using a Nyström method, and we apply the reconstruction scheme to approximate the convex hulls of the Dirichlet obstacles D 1 and of the Neumann obstacle D 2 with the same parameters as in the previous example. The reconstructions C dir and C neu are shown in figure 8 (middle and right).
Again the approximations of the convex hulls of the Dirichlet obstacles and of the Neumann obstacle are satisfactory in the sense that they do not overlap and allow to separate the two components of the scatterer. Accordingly, they could also be used as a priori information that is required in the factorization method for mixed obstacles from [18,19,31] to obtain an improved shape reconstruction. ♦ Further numerical experiments show that this algorithm is also very sensitive to noise in the data.

Conclusions
Locating and estimating the shape of scatterers based on far field observations is a basic problem in remote sensing. In this work we have established a monotonicity based shape characterization for a mixed inverse obstacle scattering problem. Using this technique we have shown that the shape of Dirichlet and Neumann obstacles are uniquely determined independently by the corresponding far field operator without additional a priori information. Numerical examples have been presented to illustrate the potential and limitations of applications of these theoretical results in reconstruction algorithms. However, the question of how to apply the novel monotonicity principles in an efficient and robust shape reconstruction algorithm for mixed inverse obstacle problems requires further research efforts.