The scattering problem of obliquely incident electromagnetic waves by an inhomogeneous infinitely long cylinder

We study the scattering of a time-harmonic electromagnetic wave incident at an oblique angle on a dielectric, isotropic, and fully inhomogeneous cylinder. This infinitely long cylinder is oriented parallel to the z − axis in three dimensions. Our main focus is on the inverse problem, aiming to reconstruct the refractive index of the cylinder. The material parameters’ z − independence allows us to consider the scattering problem only on its horizontal cross-sections. We examine the well-posedness of the direct problem and proceed to numerically solve the corresponding inverse problem for one incident field. To achieve accurate reconstructions of the refractive index, we introduce a novel numerical scheme that yields satisfactory results.


Introduction
One of the most challenging problems in the field of inverse problems is the reconstruction of a material's refractive index, commonly known as the inverse medium problem.This problem captures the interest of researchers from diverse areas, including mathematics, physics, and mechanics, due to its wide applicability in fields such as non-destructive testing of materials, quantitative medical imaging, and radar imaging.
However, it poses significant difficulties, presenting both theoretical and numerical challenges.Many of these problems are ill-posed, which means that even if we establish the uniqueness of a solution for a medium with specific geometry, the solution remains unstable.From a numerical perspective, the task of recovering a three-dimensional functional from multiple two-dimensional data (far-fields) requires a considerable amount of time and memory, making the process computationally expensive.Consequently, iterative schemes are predominantly employed to tackle this issue.
Motivated by the initial works of Wang and Nakamura [8,9,14], two of the authors extended those results and investigated the well-posedness of the direct scattering problem involving a dielectric homogeneous cylinder.Additionally, they addressed the corresponding inverse problem, aiming to reconstruct the scatterer's shape [2,3].In [6], the direct problem was examined for a doubly-connected cylinder.This year, the authors examined the direct problem for a cylinder with piecewise constant material parameters (electric permittivity and magnetic permeability) [15].This marked the first attempt to consider some form of inhomogeneity inside the scatterer.Classical tools like the boundary integral equation method and the Shapiro-Lopatinskij condition were utilized to establish the well-posedness of the direct problem.
Although the general scattering problem is formulated in 3D, the z − independence of the material parameters allows us to reduce it to 2D.Previous works considered bounded domains, resulting in boundary conditions with both normal and tangential derivatives of the fields.The novelty of this work lies in formulating the problem in the entire space 3   with compactly supported electric permittivity while assuming a constant magnetic permeability.We achieve this by transforming the direct problem into a Lippmann-Schwinger type integral equation stated in 2  .
We establish the unique solvability of the direct problem by employing Green's formulas and the Rellich Lemma.To demonstrate the existence of a solution, we utilize an equivalent system of Lippmann-Schwinger integral equations.The properties of the corresponding integral operators allows us to apply the Riesz-Fredholm theory.For fundamental understanding of these theoretical tools, we refer readers to the books [16][17][18].
To numerically solve the inverse problem of reconstructing the refractive index (or equivalently, the electric permittivity) of the scatterer using only one incident field, we employ an iterative scheme.This iterative approach allows us to recover the unknown function point-by-point, taking into account both the equations in 2  and the knowledge of the far-field patterns of the scattered fields.
The volume integrals present in the system of equations are non-linear with respect to the refractive index.To address this challenging problem, we propose two distinct approaches.The first one, based on a locally homogeneous assumption, simplifies the final system of integral equations.The second approach involves an iterative scheme that effectively decouples the linear and non-linear terms.
We solve both systems in two distinct steps.Firstly, we focus on solving for the unknown density functions.Inspired by the well-known Jacobi method, we determine the unknown material parameter point-wise.This approach of decoupling the problem into two steps is inspired by similar methods used in various inverse problems for boundary integral equations, as seen in [19][20][21].
For the numerical implementation, we utilize quadrature rules to compute the volume integrals and approximate the partial derivatives.However, the discretized linear system poses numerical challenges due to its ill-conditioned nature.To address this, we employ Tikhonov regularization.The presented numerical results demonstrate the feasibility of the proposed scheme.
The paper is organized as follows: In section 2 we formulate the direct problem and we prove its unique solvability.In section 3 the corresponding inverse problem is stated.We present the equivalent system of integral equations to be to solved for the refractive index.The numerical examples and the reconstructions are presented in the last section.

The direct scattering problem
 We consider a dielectric, isotropic and inhomogeneous cylinder in ,

3
 that is infinitely long and parallel to the z-axis.It is characterized by the spatial-dependent electric permittivity ò(x) = ò(x) > 0 and the constant magnetic permeability μ(x) = μ 0 .We assume that the electric permittivity is constant (ò = ò 0 ) outside a ball containing the inhomogeneity.Then, the electromagnetic wave (E, H) satisfies the Maxwell equations Let θ ä (0, π) denote the incident angle with respect to the negative z − axis and f ä [0, 2π], the polar angle.Then, a time-harmonic transverse magnetic (TM) obliquely incident wave takes the form [2] for the polarization vector q f q f q = p cos cos , cos sin , sin , ( ) and the incident direction q f q f q = d sin cos , sin sin , cos . where ω > 0 is the frequency.In the following, we consider time-harmonic fields and we define the scattered wave by Considering the cylindrical symmetry of the scatterer, which is parallel to the z-axis and has z-independent material parameters, we proceed to examine fields of the form [2, 9] The Maxwell equations (1) associated with a wave of the form (2) are reduced in the following system: where e and h are the z − components of the fields e and h, respectively.Here k 2 (x) = ω 2 μ 0 ò(x) − β 2 and = -J 0 1 1 0 .

( )
We define the refractive index = n x where we have used that Similarly, the z − components of the incident fields, using (2) and (5), take the form With the presentation of the necessary conditions, we are now in position to precisely formulate the direct scattering problem.In this context, the main goal is to find the solution (e, h) of (4) in , 2  given an incident wave of the form (5), such that the scattered wave satisfies the Sommerfeld radiation condition uniformly for all directions x/r, where r = |x|.
Proof.To prove uniqueness, it is sufficient to show that the corresponding homogeneous problem, meaning for  denote a disk with radius r and boundary G , r which contains the inhomogeneity.We apply Green's first theorem in W r and considering (3) we obtain where ν denotes the unit normal vector.Similarly Green's theorem for the magnetic field implies The application of Green's divergence theorem in W r results in The left-hand side vanishes under the assumption that k remains constant at the boundary G .r In addition, the right-hand side admits the decomposition • ( ) for every scalar-valued function f .Combining the last two equations, we conclude that We take the imaginary parts of (8) and (9) and we add the two equations: The last equation together with (10) From the radiation condition (7) it follows that [14] ò ò ò Then, from (11) and the application of Rellich's lemma [17]  To establish the existence, we will transform the scattering problem into a set of integral equations of Lippmann-Schwinger type.To achieve this, we initially define the fundamental solution where H 0 The given operators depend on the function inside the brackets and operate on the function f.
 ( )is a solution of the scattering problem ( 4) - (7), then e h , solve where I is the 2 × 2 identity matrix and  ( )is a solution of the integral equation ( 14), then e h , is a solution of ( 4) - (7).
Proof.This proof follows the same lines of the proof of [17], theorem 8.3.Let e h , be the solution of (4)- (7).As in the proof of theorem 2.1, we consider an arbitrary point Î x Here we have used that x y e x y y ds y x y e e x y dy , , , , 0 , ⧹ since e sc satisfies the radiation condition (7).Remark that in the volume integral, we can integrate over all , 2  since m has support in W r and ò and k 2 are constant outside W .
r Similarly, for the magnetic field, we have

y y ds y k m y h y k y k y h y k y J k y e y x y dy
, , . 16 Combining (15) with (16) and using the definition (13), we conclude that (14) is satisfied.The second component of u inc is zero due to the special form of the incident field, see (5).
 ( )be a solution of (14), then since m has compact support, from [17], theorem 8.1, we get that Î e h C , .

e x k m y e y k y y y k y e y k y y J k y h y x y dy h x k m y h y k y k y h y k y J k y e y x y dy
, , 1 , .17 Since Φ satisfies the Sommerfeld radiation condition (7) uniformly with respect to y, and m is compactly supported, it is easily proved that e h , sc sc satisfy also the radiation condition.In addition, we obtain that the above set of equations results in (4).This completes the proof. , We would like to point out that the above theorem guarantees the equivalence between the scattering problem (4) -(7) and the system of integral equations (14).Next we prove existence of solution.
Proof.By theorem 2.2 it is sufficient to prove existence of a solution of (14).We remark that the operator ( ) have a weakly singular kernel and hence are compact [17].Thus, the operator ( ( )) is also compact and the Riesz- Fredholm theory applies.We have to show that the homogeneous equation ( and hence by theorem 2.1 it follows that = u 0. ,

The inverse scattering problem
We are interested in solving the inverse problem of reconstructing the refractive index n(x), for Î x 2  from the knowledge of the far-field patterns of the scattered wave To compute the far-fields, we utilize the representations given in (17) and proceed with the following system: where the superscript (∞) in the integral operators means that we deal with the far-field operators, where the Green function Φ is replaced by x y e e x S y , 8 , , .
We define the operators Then, the inverse problem is equivalent to the system of integral equations: In the above system of integral equations, we have interchanged the roles of the density and material functions, as our primary focus now lies on the function m (equivalently n).Notably, we have two equations, (19a) and (19b), defined in 2   , while the far-field equations (19c) and (19d) are defined on the unit circle S.
However, all three unknown functions (the two density functions and m) exist in 2  .This results in an underdetermined system of equations, motivating us to consider a point-wise solution.Thus, the consideration of the far-field equations independently, combined, or even just one of them does not result in any significant differences.In the following, we consider as far-field equation the summation of (19c) and (19d); however, the other two cases yield similar outcomes.We observe that the integral operators M and N act non-linearly on the unknown function m.To address this challenge, we propose a local and iterative approach for solving the inverse problem.The numerical solution will be obtained through a point-wise approach, which will be specified later.Scheme 1. Considering our intention to solve the system of equations point-wise, we make the assumption that the unknown function n is locally homogeneous.This assumption allows us to neglect the operators involving derivatives of n.Under this assumption, the system (19) is reduced to Positively, the obtained system of equations contains only integral operators that exhibit linear behavior when applied to both density and material functions, thereby eliminating the need for linearization, a common requirement in most inverse problems.
Scheme 2. Another way to handle the non-linearity of the full system (19) is to apply an approximation method.
We propose an iterative scheme where the terms involving gradients of n will be treated as known values from the previous iteration step.We propose to solve (19) by the iterative scheme: The initial step of this scheme coincides with the system (20).Afterwards, starting from the second step and onwards, the contributions of the terms with the gradients are incorporated.
For both schemes, we propose a double iteration method where the density functions are recovered 'globally' and the material function m 'point-by-point'.The outer loop, to split the system (20) into two sub-problems is well studied and it has been applied successfully to many inverse problems [19,23,24].The inner loop, corresponding to the point-wise reconstruction of the function m is based on the Jacobi iterative method.The idea comes from linear algebra and the block splitting methods.There, the linear equation for a matrix-valued unknown function is split into a sequence of equations for its columns and then the Jacobi method is applied, see for example [25,26].The sketch of the algorithm for the first scheme is presented in algorithm 1 and in figure 1 we illustrate how m is reconstructed point-by-point.The double iteration method is applied also to the system (21) of the second scheme.
Input: Far-fields ¥ e x ( ˆ) and ¥ h x , ( ˆ) for all Î x S ˆand initial guess m x , 0 ( ) ) solve the well-posed sub-system of equations to recover the density functions e k ( ) and h .

S e h m x e x h x S e h m
x , 2 2 to be solved for m .
( ) Figure 1.Graphic illustration of the Jacobi-type inner loop of algorithm 1.At every iteration step, only one coefficient (red) is considered as unknown.The recovered value (blue) is used as input for the next step.
is satisfied, where m * is the exact solution.

Numerical implementation
We consider the square domain We apply the differential quadrature method to numerically compute the partial derivatives appearing in the integrals.This method approximates the derivative of a function by a weighted linear sum at some nodal points [27].
Let x = (x 1 , x 2 ) ä D. We use the following non-uniformly distributed grid points The differential quadrature method combined with such discretization gives accurate results [28,29].Since all involved functions (density and material) are considered to be at least continuously differentiable, we apply the Clenshaw-Curtis quadrature method to compute numerically the double integrals in D [30][31][32].
In the chosen domain D, the Clenshaw-Curtis method is as accurate as the Gauss quadrature and its advantage is that the weights are calculated using the FFT [33].In the following examples, we consider W Í Ì D .
Here, λ 0 > 0 represents the initial regularization parameter.The value is determined through a trial-anderror process; however, the obtained reconstructions are sufficiently satisfactory, rendering any sophisticated selection method for λ 0 unnecessary at this stage.

Reconstructions
In the following examples, we consider Ω r to a the layered disk with center zero and radius one.We set ò 0 = μ 0 = 1.The far-field data are calculated for n = 32 and we add δ = 3% noise.The results after ten iterations are presented in figures 2 and 3.
In the second example, we consider a refractive index that violates the regularity properties of section 2 to demonstrate the feasibility of the proposed scheme.We set In figures 4 and 5 we present the true (left) and the recovered (right) function after eight iterations with N = 31 sampling points and initial guess n (0) =1 − m (0) = 2.8.The regularization parameter was set λ 0 = 0.001.Scheme 4. Here, we set w = 2, and the incident and polar angles are given by q p = 3 and f p = 8, respectively.In this example, the refractive index is given by n x e e x 1.9 , , The recovered refractive index is presented in figures 6 and 7 for noisy data after ten iterations and N = 41.The reconstructions demonstrate comparable results to those obtained from the simplified system, highlighting the applicability of the 'local' assumption.
We observe that we get satisfactory reconstructions after few iteration steps for both piecewise constant and spatial dependent refractive indices from noisy data.Both numerical schemes produce stable reconstructions.As the iterations progress, we see that the error decreases since the previously reconstructed values contribute to the next initial guesses, see also figure 1.However, since we solve point-by-point the increase from N = 31 (second example) to N = 41 (first and third examples) increase the computational time from few to several minutes.

Conclusions
In this study, we have addressed the scattering problem of a time-harmonic electromagnetic wave by an infinitely long, inhomogeneous, and penetrable cylinder.Considering the z − independence of the material parameters, we focused on analyzing the scattering problem on the horizontal cross-sections of the scatterer.Our investigation included examining the well-posedness of the direct problem and successfully solving the corresponding inverse problem through numerical methods.We have introduced two numerical algorithms, with a local and iterative nature, to reconstruct the refractive index satisfactorily while minimizing the impact of the non-linear terms and using only one incident field.
is satisfied uniformly for all directions x/|x|.We define b q the compactly supported function m ≔ 1 − n.Then, the equation (3) using the vector identity ∇ • ( f∇g) = fΔg + ∇f • ∇g, for scalar functions f and g, takes the form

1 ()
is the Hankel function of first kind and zero order and the volume integral operators

2 r
and we choose an open disk W , r with boundary G , r containing the support of m, such that Î W x .Using Green's theorem for e in W r and considering (4) we obtain Here ¥ S j refers to the discretized form of the far-field operator acting only on m .j The operator appearing in(22) is compact, and we apply Tikhonov regularization in order to solve it.Update = -

r 2 
Then, the integrals defined in (13) admit the discretized forms

Scheme 3 . 4 .⧹
Scheme 3. The frequency is w = 3, the incident angle q p = 3 5and the polar angle f p = 4.In the first example, we consider the function

Figure 2 .
Figure 2. The exact refractive index n of the first example (left).The reconstructed value after ten iterations from noisy data (right).

Figure 3 .
Figure 3.The refractive index n of the first example in x − y plane.Top left: true function.Top right: initial guess.Bottom left: reconstructed function.Bottom right: relative error.

Figure 4 .
Figure 4.The exact refractive index n of the second example (left).The reconstructed value after ten iterations from noisy data (right).

Figure 5 .
Figure 5.The refractive index n of the second example in x − y plane.Top left: true function.Top right: initial guess.Bottom left: reconstructed function.Bottom right: relative error.

Figure 6 .
Figure 6.The exact refractive index n of the third example (left) and its reconstructed value after ten iterations from data with 2% noise (right).