Boundary Integral Equations Method for the Time-Harmonic Electromagnetic Scattering

with prescribed Dirichlet or Neumann boundary condition. For the uniqueness of this problem, usually the asymptotic behavior of the solution at infinity will also be given. Using the classical boundary layer approach (see [1,2]), this kind of problem can be solved uniquely using the double layer potential for the Dirichlet boundary condition, the single layer potential for the Neumann boundary condition, respectively, provided k is not an interior eigenvalue of the −∆ operator (see [3]). Since physically the uniqueness of the exterior boundary value problem is not in question, the cause of the non-uniqueness must have come from our mathematical method. Remaining in the framework of boundary integral equation method, various modifications have been developed to over-come this non-uniqueness for all k≠0 ([4] by adding a volume potential, [5-7] by using the combined boundary layers).


Introduction
Consider the exterior boundary value problem for the scalar Helmholtz equation in R n , n = 2, 3 ∆u + k 2 u = 0, k ∈ R, k≠0, in R n \D (1) with prescribed Dirichlet or Neumann boundary condition. For the uniqueness of this problem, usually the asymptotic behavior of the solution at infinity will also be given. Using the classical boundary layer approach (see [1,2]), this kind of problem can be solved uniquely using the double layer potential for the Dirichlet boundary condition, the single layer potential for the Neumann boundary condition, respectively, provided k is not an interior eigenvalue of the −∆ operator (see [3]). Since physically the uniqueness of the exterior boundary value problem is not in question, the cause of the non-uniqueness must have come from our mathematical method. Remaining in the framework of boundary integral equation method, various modifications have been developed to over-come this non-uniqueness for all k≠0 ( [4] by adding a volume potential, [5][6][7] by using the combined boundary layers).
In [8,9], the method of combined layers has been modified slightly and extended to the case of a vector Helmholtz equation. In this paper, we compare the numerical results of the classical boundary integral equation methods with those of the modified boundary integral equation methods with combined boundary layers for the vector Helmholtz equation. The plan of this paper is as follows. In the second section, we deduce the exterior boundary value problem from the time harmonic electromagnetic scattering problem from an ideal conductor. In the third section, we use the modified boundary layers approach to win a system of Fredholm integral equation of the second kind which is uniquely solvable. A section on numerical method is then followed by another giving some numerical examples. At the end of the paper, we give some properties of the system of the boundary integral equations from section 3 in the appendix.

Boundary Value Problem
Consider the electromagnetic scattering from a homogeneous, isotropic obstacle in R 3 . The electric Field ε and the magnetic Field Ή satisfy the Maxwell equations Where ε is the dielectricity, µ is the permeability and σ is the conductivity of the medium. For the time-harmonic electromagnetic wave of the form  with a frequency ω>0, the space-dependent complex-valued functions E and H solve the reduced Maxwell equations where the wave number k is give by with Im k≥0.
According to Stratton x x Now consider the case of an infinitely long cylinder which runs in the x 3 direction and has a cross section D parallel to the x 1 x 2 plane. We have thus the following two dimensional boundary value problems for the vector Helmholtz equation. D Denote ϑ(x) the unit tangent vector at x∈∂D which is given by ϑ(x):=e 3 ×ν(x). Given ξ, ρ∈C0,α (Γ,C),0<α< 1, find a vector field E: The last condition in the boundary value problem (H) is the Sommerfeld radiation condition. It describes an outgoing wave and is essential for the uniqueness of the problem. Knauff and Kress [8] proved the unique solvability of the 3-dimensional problem 1. With slight modifications, their proof can also be used here to show the wellposedness of the problem (H) (See sec. 3.5 in [3] for details).

Boundary Integral Equations
Motivated by the representation theorem from stratton and Chu [10] (see also [8]), we define the solution ansatz for our problem H.
For the numerical treatment of this problem, we need to parametrize this boundary integral equation. Write for a 2π-periodic, two times continuously differentiable function x : R→R 2 .
After a lengthy mathematical computation, which will be given in some details in the appendix, the boundary integral equation (4) will be brought into the following form which is an integral equation of the second kind with an bounded invertible L and a compact A.

Numerical Method
We will apply the Nyström method to solve the boundary integral equation (5) numerically. This means that we need some quadrature rules for our integrals appear in the operators L, A. To begin with, we define the space of trigonometric functions For the numerical solution of the boundary integral equation (5), we need to approximate the integrals. Because of the singular behavior of the Hankel functions, we have 3 different kinds of integration (see appendix). Therefore it is wise to use different quadrature rules. We use the following well-known convergent quadrature rules: The quadrature (6) is used by Garrick [11], (7) is used in Martensen [12] and in Kussmaul [13]. The quadrature rule (8) We call (9) semi-discrete because it is still a functional equation. To win a full discrete system, we solve this equation at the collocation points which are the same as the interpolation points used for the trigonometric interpolation.
Mathematically, we define a projection operator Pn: C 0,α → C n × C n . Then the full discrete system reads P n (L n + A n )φ n = P n f.
We note here that the convergence φ n →φ can be showed as a consequence of the convergent quadratures (6)-(8) and the interpolation. We omit the proofs which are not essential in this paper.

Numerical Results
In this section we will demonstrate the efficiency of our method through some examples. To test the advantage of the method, the parameter η is taken to be 1. The boundary value problem will be solved by the Nyström method with trigonometric polynomials as the underlying interpolatory space. From the asymptotic behavior of the Hankel function, we have the following representation for the far field pattern to the solution of our boundary value problem.

Example 1
In the first example, we choose the ellipse with the parametrization For the incident field Ei , we use the plane wave mulitplicated by a constant vector is calculated for the two directions d and −d. We deal with the case of a smaller wave number (k=1 , Table 1) and the case of a larger wave number (k=10, Table 2). We see that the convergence is very fast in the case of an ellipse.

Example 2
In the second example, we choose a domain which is non-convex and non-symmetric. The bean-shaped domain is parametrized by x(t)=(ρ(t) cos(t), ρ sin(t)), t ∈ [0, 2π] with the function
note here that the convergence is because of the shape of the boundary slower than that in the first example. However, it is still satisfactory.

Example 3
In this example, we' d like to compare the method η=1 with η=0 for the stability in the case where the explicit solution of the boundary value problem is known. For this purpose we take the unit circle Ω as the domain of interest: x(t)=(cos(t), sin(t)),t ∈ [0, 2π].
As mentioned earlier, the solution of the boundary value is not unique if the wave number k is an interior eigenvalue. In the case of a unit circle, it happens to be the zeros of the Bessel functions Jn and their derivatives The solution to the boundary value problem is The corresponding far field pattern is given by We' d like to compare the two methods on a smaller zero and a larger zero, namely z 1 = 2.4048255576957727686 . . . ,and z 5 = 14.9309177084877859477 . . . .
We will solve the boundary value problem with those wave number k's which are close to the numbers z 1 and z 5 , respectively. The results for η=1 and η=0 will then be compared. The results are listed in the Tables  5-10.
At this place, we want to draw some conclusions from our numerical results. The first thing to note is that in the case where k is not very close to the eigenvalue, both the modified method η=1 and the classical method η=0 converge fast (Table 5, Table 6, Table 8). In the case where k is close to the eigenvalue, while the modified method converges really fast regardless of the value of k, the classical method does not converge ( Table 7, Table 9, Table 10). Secondly, the modified approach is very robust in the sense that it is very stable even when the wave number k is numerical the same as the eigenvalue of the interior problem (Table 10).
From this we conclude that the modified method is stable for all wave numbers.

Appendix
Here we will explain the mathematical transition from (4) to (5). It is known that the method of boundary integral equations has the advantage that it reduces the dimension by one. On the other hand, it  Table 3: Bean, x(t) = ( (t) cos(t); _(t) sin(t)), k = 1.  Table 4: Bean, x(t) = (ρ(t) cos(t); ρ (t) sin(t)), k = 10.        is also known that one often has to struggle with stronger singularities. Therefore, one main point here is to split the singularities from the equations and to treat them in a unified profitable way. The singularities come from the Hankel functions (1) After some cumbersome calculations, which we' d like to omit, the functional equation (4) becomes (5) with   11  12  11  12  11  11  12  12   21  22  21  22  21  21  22  22   :  ,  : , where (L11 ψ)(s):= −ψ(s), Instead of writing down the details and all the functions, we make some remarks here: 1. The functions mij (s, t), r ij (s, t), i, j = 1, 2 are as smooth as the parametrization of Γ is (which is assumed to be of class C 2 ). η = 1 η = 0 n ( ) 2. It can be shown that both ln t and ln 2 (4sin ) 2 t has the same singularity as t → 0. Since a convergent quadrature for the latter is already at hand (7), we just split the log -singularity as we did here. The same reason is taken into account for the use of cot Summing up, we have splitted and grouped the singularities in such a way that it is beneficial both for the theoretical treatment (unique solvability using Fredholm's theory for the 2. kind equations) and the numerical computation (convergent quadrature rules).