Shape Reconstruction of Thin Electromagnetic Inclusions via Boundary Measurements: Level-Set Method Combined with the Topological Derivative

Abstract

An inverse problem for reconstructing arbitrary-shaped thin penetrable electromagnetic inclusions concealed in a homogeneous material is considered in this paper. For this purpose, the level-set evolution method is adopted. The topological derivative concept is incorporated in order to evaluate the evolution speed of the level-set functions. The results of the corresponding numerical simulations with and without noise are presented in this paper.

1. Introduction

In general, one of the main purposes of the inverse scattering problem is to determine the characteristics of an unknown target, for example, the shape, location, and internal constitution, based on the measured scattered field data. Therefore, there exist a considerable number of significant inverse scattering problems for identifying the information of arbitrary-shaped thin penetrable electromagnetic inhomogeneities or perfectly conducting cracks concealed in structures such as concrete walls of buildings and bridges (refer to [1] and the references therein).

In order to solve related inverse problems, various algorithms have been developed. Among them, an iterative algorithm based on level-set evolution has been established in [2] and many generalizations have been investigated (refer to [3] and the references therein). However, most research studies have focused on the shape reconstruction of volumetric targets. Therefore, in order to reconstruct crack-like targets via level-set evolution, a new concept is required.

In recent work [4], arbitrary-shaped thin penetrable conductivity inclusions have been successfully reconstructed via the level-set method by adopting two different level-set functions. Shape reconstruction of thin electromagnetic inclusions is considered in [5]. However, in [5], only the case of the purely permittivity contrast or permeability contrast with respect to the background material is considered. So, extension to both contrast cases is still remaining as a research topic.

Motivated by the mentioned above, the main objective of this paper is to reconstruct the shape of thin electromagnetic inclusions by adopting the level-set method in the case of both the dielectric permittivity and the magnetic permeability contrasts with respect to the background material. Note that, in general research papers, the evolution speed is evaluated via a complex calculation of the Fréchet derivatives using an adjoint technique. Instead of such evaluations, we incorporate a topological derivative to derive the evolution speed of two level-set functions introduced in [6–10].

This paper is organized as follows. In Section 2, we introduce two-dimensional direct scattering from an arbitrary-shaped thin penetrable electromagnetic inclusion and the level-set functions for describing such inclusions. In Section 3, the evolution speed of the level-set functions is derived via the topological derivative. We present some numerical simulation results in Section 4 to show the effectiveness of the level-set method. Finally, we provide short conclusions in Section 5.

2. Direct Scattering Problem and Two Level-Set Functions: A Brief Survey

2.1. Two-Dimensional Direct Scattering Problem

Let Ω be a two-dimensional homogeneous domain with a smooth boundary 𝜕Ω. Throughout this paper, we assume that Ω conceals an electromagnetic inclusion Γ with thickness 2ℎ, which is characterized in the neighborhood of the smooth curve 𝜎: Γ∶=𝐫+𝜂𝐧(𝐫)∶𝐫∈𝜎,𝜂∈(−ℎ,ℎ), (1) where 𝐧(𝐫) is a unit normal to 𝜎 at 𝐫. Throughout this paper, we assume that Γ does not touch 𝜕Ω; that is, we assume that Γ∩𝜕Ω=∅.

In this paper, the electric permittivities and magnetic permeabilities of Γ and Ω at a given nonzero frequency are assumed to be known; they are finite valued and they differ (either one or both of them) from the ones of the homogeneous embedding medium. For the sake of simplicity, let 𝜀₀ and 𝜇₀ denote the permittivity and permeability of Ω, and let 𝜀 and 𝜇 denote those of Γ, respectively. Then, we can define the following piecewise constants: 𝜀(𝐫)=𝜀₀for𝐫∈Ω⧵𝜇(𝐫)=𝜇Γ,𝜀for𝐫∈Γ,₀for𝐫∈Ω⧵Γ,𝜇for𝐫∈Γ. (2)

Let 𝑢_Γ^(𝑛)(𝐫) denote a time-harmonic total field that satisfies the following Helmholtz equation: 1∇⋅𝜇(𝐫)∇𝑢_Γ^(𝑛)(𝐫)+𝜔²𝜀(𝐫)𝑢_Γ^(𝑛)(𝐫)=0inΩ (3) with boundary condition 1𝜇₀𝜕𝑢_Γ^(𝑛)(𝐫)(𝐫)=1𝜕𝜈𝜇₀𝜕exp𝑖𝑘₀𝐝_𝑛⋅𝐫(𝐫)(𝐫)𝜕𝜈=𝑔on𝜕Ω, (4) where 𝑘₀=𝜔√𝜀₀𝜇₀, 𝜈(𝐫) denotes the unit outward normal to 𝐫∈𝜕Ω, and 𝐝_𝑛 denotes a vector on the unit circle 𝕊¹. In the same manner, we denote 𝑢^(𝑛)(𝐫)=exp(𝑖𝑘₀𝐝_𝑛⋅𝐫) as a background solution of (3) with boundary condition (4).

2.2. Description of Thin Inclusion Using Two Level-Set Functions

Now, we describe the thin inclusion Γ by two different level-set functions (see [4, 5] for a detailed description). From now on, we assume that the level-set function 𝜑 is a 𝒞¹ function. With this assumption, we introduce the zero level-set function ℒ₀ and its boundary 𝜕ℒ₀ as ℒ₀∶=𝐫∈ℝ²∶𝜑(𝐫)≤0,𝜕ℒ₀∶=𝐫∈ℝ²∶𝜑(𝐫)=0, (5) respectively. Then, a thin region ℒ_𝑇 of thickness 2ℎ and inner boundary 𝜕ℒ₁ can be written as follows: ℒ_𝑇∶=ℒ₀∩𝐫−𝜂𝐧(𝐫)∶𝐫∈𝜕ℒ₀,0≤𝜂≤2ℎ,𝜕ℒ₁∶=𝜕ℒ_𝑇⧵𝜕ℒ₀, (6) respectively.

In order to describe a thin region of finite length (which is connected or disconnected) in Ω, we assume that another level-set function 𝜓 is also a continuously differentiable function. Let us define a band structure ℒ_𝐵⊂ℝ² by ℒ_𝐵∶=𝐫=𝑟₁,𝑟₂∈ℝ²∶𝜓𝑟₁<0. (7) Then, we can describe Γ as Γ=ℒ_𝑇∩ℒ_𝐵∩Ω (8) and the boundary 𝜕ℒ₂ of Γ for evolving at inclusion tips as 𝜕ℒ₂∶=𝐫=𝑟₁,𝑟₂∈𝜕Γ∶𝜓𝑟₁=0. (9)

Figure 1 shows the illustration of the thin inclusion Γ using two level-set functions.

(a)

(b)

Figure 1 

Illustration of the thin inclusion Γ using two level-set functions.

3. Evolution Speed of the Level-Set Functions and the Topological Derivative

The evolution of the level-set function 𝜙 satisfies the Hamilton-Jacobi type equation: 𝜕𝜙(𝐫)(𝐫(𝐫𝜕𝑡=𝔽)∇𝜙)=0, (10) where 𝜙 is 𝜑 or 𝜓 introduced in Section 2.2 and 𝔽(𝐫) is given by (see, e.g., [4, 5]) the following. (i)For level-set function 𝜑(𝐫) and 𝐫∈𝜕ℒ₀, 𝔽(𝐫)=Re𝕊(𝐫)𝛿𝐫;𝜕ℒ₀−𝛿𝐫+2ℎ𝐧(𝐫);𝜕ℒ₀. (11) (ii)For level-set function 𝜓(𝐫) and 𝐫=(𝑟₁,𝑟₂)∈𝜕ℒ₂, 𝔽(𝐫)=Re_𝜕ℒ2𝕊(𝐫)𝛿𝐫;𝜕ℒ₂𝜒(𝐫;ℒ)𝑑𝑟₂. (12)Here, 𝕊(𝐫) denotes the speed function of the level sets that we have to evaluate. In a recent work [5], 𝕊(𝐫) has been derived for the purely permittivity contrast (𝜀≠𝜀₀ and 𝜇=𝜇₀) or permeability contrast (𝜀=𝜀₀ and 𝜇≠𝜇₀) cases by evaluating the Fréchet derivative. In order to evaluate 𝕊(𝐫) for both the permittivity and permeability contrast cases, we introduce the following discrepancy functional: 𝔼(Ω)1∶=2^𝑁_𝑛=1𝑢_Γ^(𝑛)𝐫^′−𝑢^(𝑛)𝐫^′2_{𝐿2()𝜕Ω}=12^𝑁_𝑛=1𝜕Ω𝑢_Γ^(𝑛)𝐫^′−𝑢^(𝑛)𝐫^′2𝑑𝑆𝐫^′. (13) Then, we adopt the following relationship between the evolution speed of level sets and the topological derivative.

Theorem 1. Let 𝐷_𝑇𝔼(𝐫) denote the topological derivative of (13); then 𝐷_𝑇𝔼(𝐫) satisfies 𝐷_𝑇𝔼(𝐫)=(𝐫)𝕊(𝐫)−𝕊𝑓𝑜𝑟𝐫∈Γ𝑓𝑜𝑟𝐫∈Ω⧵Γ. (14)

Proof. A rigorous derivation can be found in [6, 8, 9].

Hence, we can evaluate the evolution speed 𝕊(𝐫) by computing the topological derivative 𝐷_𝑇𝔼(𝐫). The result is as follows.

Theorem 2. In both the permittivity and permeability contrast cases, 𝐷_𝑇𝔼(𝐫) can be written as 𝐷_𝑇𝔼(𝐫)=ℜ𝔢^𝑁_𝑛=1𝜔²𝜀−𝜀₀|𝒟|𝑣^(𝑛)(𝐫)𝑢^(𝑛)(𝐫)−2|𝒟|𝜇−𝜇₀𝜇+𝜇₀∇𝑣^(𝑛)(𝐫)⋅∇𝑢^(𝑛)(𝐫), (15) where 𝑣 satisfies the following adjoint problem: 1𝜇₀Δ𝑣^(𝑛)(𝐫)+𝜔²𝜀₀𝑣^(𝑛)(𝐫)1=0𝑖𝑛Ω,𝜇₀𝜕𝑣^(𝑛)(𝐫)𝜕𝜈=𝑢_Γ^(𝑛)(𝐫)−𝑢^(𝑛)(𝐫)𝑜𝑛𝜕Ω. (16)

Proof. See the appendix.

Note that, if either 𝜇=𝜇₀ or 𝜀=𝜀₀, we can obtain the Fréchet derivatives based on [5, Theorem 3.1] or [5, Theorem 3.2], respectively. This result is identical to that of [11, Theorem 3.2], but the proof is considerably easier.

4. Numerical Simulations

In this section, we choose the background material Ω as a unit circle centered at origin with 𝜀₀=𝜇₀=1. For thin inclusions Γ_𝑗, the thickness is set to ℎ=0.02 and the parameters of are chosen as 𝜀_𝑗=𝜇_𝑗=5, where 𝜀_𝑗 and 𝜇_𝑗 are the permittivity and permeability contrasts of Γ_𝑗, respectively. The supporting curve 𝜎_𝑗 is selected as 𝜎₁=𝑥,0.5𝑥²[,𝜎+0.1∶𝑥∈−0.5,0.5]₂=𝑥+0.2,𝑥³+𝑥²−0.6∶𝑥∈[−0.5,0.5]. (17)

For a given inclusion Γ_𝑗, let 𝑢_true(𝐫) be the solution to (3) with boundary condition (4). In the same manner, for an inclusion Γ_𝑗^(𝑠) at the iteration step 𝑠 (if 𝑠=0, it is an initial guess), let 𝑢^(𝑠)_comp denote the solution to (3) with boundary condition (4). For the least-square functional at a given frequency, we introduce a normalized error function: ℝ(𝑠)=𝑢_true𝐫_𝑛−𝑢^(𝑠)_comp𝐫_𝑛²_ℓ2𝑢^(𝑠)_comp𝐫_𝑛²_ℓ2for𝐫_𝑛∈𝜕Ω, (18) where 𝑓𝐫_𝑛²_ℓ2∶=^𝑁_𝑛=1𝑓𝐫_𝑛². (19) In this section, we adopt 𝑁=16 and the locations 𝐫_𝑛∈𝜕Ω are selected as 𝐫_𝑛=cos2𝑛𝜋𝑁,sin2𝑛𝜋𝑁,𝑛=1,2,…,𝑁. (20)

A finite-difference time discretization form of (10) is 𝜙^(𝑠+1)(𝐫)=𝜙^(𝑠)(𝐫)+Δ𝑡(𝑠)𝔽^(𝑠)(𝐫)∇𝜙^(𝑠)(𝐫). (21) For the evolving level sets 𝜑 and 𝜓, the time step is chosen as Δ𝑡(𝑠)=3×10^−𝑠 (see [5, 12] for a detailed discussion). With this configuration, the forward and adjoint problems are solved via the finite element method (FEM) to avoid an inverse crime, and a 20 dB white Gaussian random noise is added to the boundary measurements data in order to present the robustness.

Remark 3 (generating a good initial guess). It is worth mentioning that, for a successful performance of the iterative based algorithm, one must start the iteration procedure with a good initial guess that is close to the unknown target. Without this guess, one will encounter unexpected scenarios such as nonconvergence (refer to Figures  7 and  14 in [5]). There are various noniterative shape reconstruction algorithms such as the multiple signal classification (MUSIC)-type algorithm [13–17], the end-point identification algorithm [18–20], the multifrequency based subspace migration imaging technique [21–27], the topological derivative strategy [10, 28–32], and the gradient-for-the-initiation (GFI) method [33].

First, we consider the shape reconstruction of Γ₁. In this example, we assume that the location of the end points of Γ₁ has been identified based on the location search algorithm in [18–20]. Thus, by connecting them by a straight line segment, a good initial guess can be obtained.

Figure 2 shows the shape reconstruction of Γ₁ at the operated frequency 𝜔=2𝜋/0.6. Based on this result, we can confirm that the shape of Γ₁ is successfully reconstructed and the residual considerably decreases accordingly.

(a) Initial guess

(b) Shape of Γ1(20)

(c) Shape of Γ1(60)

(d) Shape of Γ1(100)

(e) Shape of Γ1(140)

(f) Shape of Γ1(160)

(g) Shape of Γ1(200)

(h) Shape of Γ1

(i) Value of ℝ(𝑠)

Figure 2 

Reconstruction of Γ1 with noisy data.

Next, we consider the reconstruction of Γ₂. In contrast to the previous example, we apply the topological derivative strategy introduced in [10, 30] to obtain a good initial guess. In Figure 3, the obtained shape of Γ₂ via the topological derivative with 𝑁=16 and 𝜔=2𝜋/0.4 is illustrated. In this result, we can construct a supporting curve via well-identified points as 𝜎₂⁽⁰⁾=−0.4308𝑇₀(𝑥)+0.4160𝑇₁(𝑥)+0.1482𝑇₂(𝑥)+0.1895𝑇₃(𝑥)=0.7580𝑥³+0.2964𝑥²−0.1525𝑥−0.5790 (22) for 𝑥∈[−0.2451,0.6646] (see Figure  9 in [10]). Here, 𝑇_𝑗(𝑥) denotes the Chebyshev polynomials of the first kind. Considering this curve as the supporting curve of the initial guess, the reconstruction of the complete shape of Γ₂ is illustrated in Figure 3. Note that, because of the considerably close initial guess, only 31 iterations are sufficient to obtain a good result.

(a) Identified points

(b) Initial guess

(c) 31 iterations

(d) Value of ℝ(𝑠)

(e) Shape of Γ2(0) and Γ2

(f) Shape of Γ2(31) and Γ2

Figure 3 

Reconstruction of Γ2 with noisy data.

5. Conclusion Remarks

In this paper, a level-set method to perform the shape reconstruction of thin electromagnetic inclusion using both the dielectric permittivity and magnetic permeability contrasts with respect to the homogeneous domain was considered. The relationship between the topological derivative and the evolution speed of level-set functions allowed us to reconstruct the complete shape of the thin inclusion with a good initial guess. It is discovered from the numerical experiments that the proposed technique is stable with respect to random noise. Further, incorporating the topological derivative guarantees a successful reconstruction procedure.

In this paper, we performed shape reconstruction with a priori information of unknown targets, for example, thickness, permittivity, and permeability values. Performing simultaneous reconstruction without a priori information will be a challenge (refer to [3, Section 10.3]).

Appendix

Proof of Theorem 2. In this appendix, we prove Theorem 2. In order to derive 𝐷_𝑇𝔼(𝐫) of (13), we construct a small electromagnetic inclusion Σ=𝐫+𝜌𝒟 in Ω with the same permittivity and permeability of Γ and denote this domain by Ω∨Σ. Then, by virtue in [10], the following asymptotic expansion holds: 𝔼(Ω∨Σ)=𝔼(Ω)+𝑓𝜌𝐷_𝑇𝔼(𝐫)+𝑜𝑓𝜌, (A.1) where 𝑓(𝜌)→0 as 𝜌→0+. Let 𝑢_Σ^(𝑛)(𝐫) be a solution to (3) with boundary condition (4) in the existence of Σ. Then, the following asymptotic expansion formula holds (see [13]): 𝑢_Σ^(𝑛)𝐫^′−𝑢^(𝑛)𝐫^′=𝜌²𝜇−𝜇₀𝜇+𝜇₀|𝒟|∇𝑢^(𝑛)(𝐫)⋅∇𝒩𝐫,𝐫^′+𝜔²𝜀−𝜀₀|𝒟|𝑢^(𝑛)(𝐫)𝒩𝐫,𝐫^′+𝑜𝜌², (A.2) where 𝒩(𝐫,𝐫^′) denotes the Neumann function that is the solution to 1𝜇₀Δ𝒩𝐫,𝐫^′+𝜔²𝜀₀𝒩𝐫,𝐫^′=−𝛿𝐫,𝐫^′1,𝐫∈Ω,𝜇₀𝜕𝒩𝐫,𝐫^′𝜕𝜈(𝐫)=0,𝐫∈𝜕Ω. (A.3)
Then, applying (A.2) to (13), we can compute 𝔼()=1Ω∨Σ2^𝑁_𝑛=1𝑢_Γ^(𝑛)𝐫^′−𝑢_Σ^(𝑛)𝐫^′2_{𝐿2()𝜕Ω}=𝔼(Ω)+^𝑁_𝑛=1𝔼_Σ(𝐫)+𝑜𝜌², (A.4) where 𝔼_Σ(𝐫)=_𝜕Ω𝑢_Γ^(𝑛)𝐫^′−𝑢^(𝑛)𝐫^′×𝑢^(𝑛)𝐫^′−𝑢_Σ^(𝑛)𝐫^′𝑑𝑆𝐫^′. (A.5) Then, applying boundary condition (16) and asymptotic formula (A.2) yields 𝔼_Σ(𝐫)=_𝜕Ω1𝜇₀𝜕𝑣^(𝑛)𝐫𝜕𝜈^′×𝜌²𝜇−𝜇₀𝜇+𝜇₀|𝒟|∇𝑢^(𝑛)(𝐫)⋅∇𝒩𝐫,𝐫^′+𝜌²𝜔²𝜀−𝜀₀|𝒟|𝑢^(𝑛)(𝐫)𝒩𝐫,𝐫^′𝑑𝑆𝐫^′=𝜌²𝜔²𝜀−𝜀₀|𝒟|𝔼_𝜀(𝐫)+𝜌²𝜇−𝜇₀𝜇+𝜇₀|𝒟|𝔼_𝜇(𝐫). (A.6)
First, applying integration by parts and (A.3), we can compute 𝔼_𝜀(𝐫) as 𝔼_𝜀(𝐫)=_𝜕Ω1𝜇₀𝜕𝑣^(𝑛)𝐫𝜕𝜈^′𝑢^(𝑛)(𝐫)𝒩𝐫,𝐫^′𝑑𝑆𝐫^′=_Ω1∇⋅𝜇₀𝜕𝑣^(𝑛)𝐫𝜕𝜈^′𝑢^(𝑛)(𝐫)𝒩𝐫,𝐫^′𝑑𝐫^′+_Ω∇𝑣^(𝑛)𝐫^′⋅𝑢^(𝑛)(𝐫)1𝜇₀∇𝒩𝐫,𝐫^′𝑑𝐫^′=_Ω1𝜇₀Δ𝑣^(𝑛)𝐫^′𝑢^(𝑛)(𝐫)𝒩𝐫,𝐫^′𝑑𝐫^′+_Ω𝜔²𝜀₀𝑣^(𝑛)𝐫^′𝑢^(𝑛)(𝐫)𝒩𝐫,𝐫^′𝑑𝐫^′−_Ω𝑣^(𝑛)𝐫^′𝑢^(𝑛)(𝐫)1𝜇₀Δ𝒩𝐫,𝐫^′𝑑𝐫^′+_Ω𝜔²𝜀₀𝑣^(𝑛)𝐫^′𝑢^(𝑛)(𝐫)𝒩𝐫,𝐫^′𝑑𝐫^′=𝑣^(𝑛)(𝐫)𝑢^(𝑛)(𝐫). (A.7)
Similarly, 𝔼_𝜇(𝐫) becomes 𝔼_𝜇(𝐫)=_𝜕Ω1𝜇₀𝜕𝑣^(𝑛)𝐫𝜕𝜈^′∇𝑢^(𝑛)(𝐫)⋅∇𝒩𝐫,𝐫^′𝑑𝑆𝐫^′=_Ω1𝜇₀Δ𝑣^(𝑛)𝐫^′∇𝑢^(𝑛)(𝐫)∇𝒩𝐫,𝐫^′𝑑𝐫^′+_Ω∇𝑣^(𝑛)𝐫^′⋅∇𝑢^(𝑛)(𝐫)1𝜇₀Δ𝒩𝐫,𝐫^′𝑑𝐫^′=−_Ω𝜔²𝜀₀𝑣^(𝑛)𝐫^′∇𝑢^(𝑛)(𝐫)∇𝒩𝐫,𝐫^′𝑑𝐫^′−_Ω𝑣^(𝑛)𝐫^′∇𝑢^(𝑛)(𝐫)∇1𝜇₀Δ𝒩𝐫,𝐫^′𝑑𝐫^′=_Ω−𝑣^(𝑛)𝐫^′∇𝑢^(𝑛)(𝐫)×∇1𝜇₀Δ𝒩𝐫,𝐫^′+𝜔²𝜀₀𝒩𝐫,𝐫^′𝑑𝐫^′=_Ω∇𝑣^(𝑛)𝐫^′⋅∇𝑢^(𝑛)(𝐫)×1𝜇₀Δ𝒩𝐫,𝐫^′+𝜔²𝜀₀𝒩𝐫,𝐫^′𝑑𝐫^′=−∇𝑣^(𝑛)(𝐫)∇𝑢^(𝑛)(𝐫). (A.8) Hence, 𝔼_Σ(𝐫)=𝜌²𝜔²𝜀−𝜀₀𝑣^(𝑛)(𝐫)𝑢^(𝑛)(𝐫)−2|𝒟|𝜇−𝜇₀𝜇+𝜇₀∇𝑣^(𝑛)(𝐫)⋅∇𝑢^(𝑛)(𝐫), (A.9) and by considering the real part (A.9) and comparing (A.1) and (A.4), we can observe that 𝑓(𝜌)=𝜌² and 𝐷_𝑇𝔼(𝐫)=ℜ𝔢^𝑁_𝑛=1𝜔²𝜀−𝜀₀|𝒟|𝑣^(𝑛)(𝐫)𝑢^(𝑛)(𝐫)−2|𝒟|𝜇−𝜇₀𝜇+𝜇₀∇𝑣^(𝑛)(𝐫)⋅∇𝑢^(𝑛)(𝐫). (A.10) This completes the proof.

Acknowledgments

The author would like to acknowledge the anonymous reviewer for his valuable comments. This work was supported by the Basic Science Research Program through the National Research Foundation of Korea (NRF), funded by the Ministry of Education, Science, and Technology (no. 2011-0007705), and the Research Program of the Kookmin University in Korea.