# TM Electromagnetic Scattering from PEC Polygonal Cross-Section Cylinders: A New Analytical Approach for the Efficient Evaluation of Improper Integrals Involving Oscillating and Slowly Decaying Functions

Authors: Mario Lucido; Chiara Santomassimo; Fulvio Schettino; Marco Donald Migliore; Daniele Pinchera; Gaetano Panariello
Journal: Advances in Mathematical Physics
 (2019)

---

## Abstract

The analysis of the TM electromagnetic scattering from perfectly electrically conducting polygonal cross-section cylinders is successfully carried out by means of an electric field integral equation formulation in the spectral domain and the method of analytical preconditioning which leads to a matrix equation at which Fredholm’s theory can be applied. Hence, the convergence of the discretization scheme is guaranteed. Unfortunately, the matrix coefficients are improper integrals involving oscillating and, in the worst cases, slowly decaying functions. Moreover, the classical analytical asymptotic acceleration technique leads to faster decaying integrands without overcoming the most important problem of their oscillating nature. Thus, the computation time rapidly increases as higher is the accuracy required for the solution. The aim of this paper is to show a new analytical technique for the efficient evaluation of such kind of integrals even when high accuracy is required for the solution.


---

## Body

## 1. Introduction

Spectral domain formulations are particularly suitable for the analysis of a wide class of electromagnetic problems ranging from the propagation in planar guides and waveguides or the radiation by planar antennas to the scattering from cylindrical structures or planar surfaces involving homogeneous or stratified media, just to give some examples. In general, the obtained integral equation in the spectral domain does not admit a closed form solution; hence, numerical schemes have to be adopted. The fast convergence of such methods is a key point. When dealing with polygonal cross-section cylindrical structures or canonical shape planar surfaces, just for examples, a well-posed matrix operator equation can be obtained by means of the method of analytical preconditioning [1]. It consists of the discretization of the integral equation by means of Galerkin’s method with a suitable set of expansion functions leading to a matrix equation at which Fredholm’s or Steinberg’s theorems can be applied [2, 3]. In the literature, it has been widely shown that this goal can be fully reached by selecting expansion functions reconstructing the physical behaviour of the fields on the involved objects with a closed-form spectral domain counterpart [4–15]. With such a choice, few expansion functions are needed to achieve highly accurate results and the convolution integrals are reduced to algebraic products. However, the obtained matrix coefficients are improper integrals of oscillating and, in the worst cases, slowly decaying functions to be numerically evaluated. The classical analytical asymptotic acceleration technique (CAAAT), consisting of the extraction from the kernels of such kind of integrals of their asymptotic behaviour while the slowly converging integrals of the extracted parts are expressed in closed form, allows us to obtain faster decaying integrands without overcoming the most important problem of their oscillating nature. Consequently, the convergence of the accelerated integrals becomes slower and slower as higher is the accuracy required for the solution.In order to overcome this problem, a novel technique has been proposed for the analysis of the propagation in single and multiple coupled microstrip lines in planarly layered media [16], the propagation in multilayered single and coplanar coupled striplines [17, 18], the scattering form a rectangular plate in homogeneous medium or buried in a lossy half-space [19, 20], the complex resonances of a rectangular patch in a multi-layered medium [21], the scattering from a tilted strip buried in a lossy half-space at oblique incidence [22], and recently the scattering from a circular plate in homogeneous medium [23]. By means of algebraic manipulations and a suitable integration procedure in the complex plane, the matrix coefficients, which are single/double integrals involving products of Bessel functions of the first kind, are expressed as linear combinations of proper integrals and/or improper integrals of nonoscillating functions which can be quickly evaluated. It is interesting to observe that, despite the same line of reasoning, the procedures developed in the quoted papers are, in general, different from problem to problem.In the papers [13, 15], the analysis of the electromagnetic scattering from perfectly electrically conducting (PEC) polygonal cross-section cylinders when a TM polarized plane wave orthogonally impinges onto the scatterer surface has been successfully approached. The problem has been formulated in terms of electric field integral equation (EFIE) in the spectral domain and discretized by means of Galerkin’s method. Jacobi polynomials multiplied by their orthogonality weight have been used as expansion functions. By means of a suitable choice of the polynomials’ parameters, the physical behaviour of the surface current density on each side of the polygonal cross-section and even on the adjacent wedges has been reconstructed. Moreover, it has been shown that the spectral domain counterpart of the selected expansion functions can be expressed in closed form in terms of confluent hypergeometric functions of the first kind. Due to the reciprocity theorem, it has been demonstrated that it is always possible to reduce the convolution integrals to algebraic products; i.e., the matrix coefficients are single improper integrals involving products of confluent hypergeometric functions of the first kind and complicated arguments numerically evaluated by means of CAAAT.The aim of this paper is the introduction of a new analytical technique for the efficient evaluation of such kind of integrals. Algebraic manipulations and a suitable integration procedure in the complex plane allow us to reduce each integral to a linear combination of proper integrals and improper integrals of nonoscillating functions involving confluent hypergeometric functions of the first kind and the second kind. As will be shown in the following, the proposed technique outperforms the CAAAT especially when a high accuracy is required for the solution.


## 2. Formulation and Solution of the Problem

### 2.1. Background

In Figure1, a polygonal cross-section PEC cylinder is sketched. A coordinate system x
,
y
,
z is introduced such that the z axis coincides with the cylinder axis. The L sides of the polygonal cross-section are numbered clockwise and a local coordinate system x
i
,
y
i
,
z is introduced on the i-t
h side with the origin at the centre of the side itself in the position x
=
x
-
i
,
y
=
y
-
i and the y
i axis oriented in the outward direction. The angle ϕ
i
∈ ]−π, π] denotes the orientation of the x
i axis with respect to the x axis and 2
a
i denotes the length of the i-th side. A TM polarized plane wave impinges onto the cylinder surface with an angle ϕ with respect to the x axis and orthogonally with respect to the cylinder axis, namely,(1)
E
_
i
n
c
x
,
y
=
E
0
z
^
e
-
j
k
x
x
+
k
y
y
,where k
x
=
-
k
cos
⁡
ϕ and k
y
=
-
k
sin
⁡
ϕ, k
=
2
π
/
λ
=
ω
ε
μ is the wavenumber, λ is the wavelength, ε is the dielectric permittivity and μ the magnetic permeability of the external medium, ω is the angular frequency. With such a choice, only TM solutions can be obtained; i.e., the induced current is longitudinal and the electromagnetic field is invariant along the z axis.Figure 1
Geometry of the problem.By means of Meixner’s theory [24], the following edge behaviour can be established for the longitudinal current density on the i-th surface:(2)
J
i
z
x
i
=
1
-
x
i
a
i
t
i
+
1
+
x
i
a
i
t
i
-
J
-
i
z
x
iwith i
=
1,2
,
…
L, where (3)
t
i
±
=
ψ
i
±
-
π
2
π
-
ψ
i
±
ψ
i
±
≤
3
π
2
1
ψ
i
±
≥
3
π
2
,
ψ
i
± being the angles of the wedges at the abscissas x
i
=
±
a
i and the well-behaved function J
-
i
z
x
i belongs to the weighted Hilbert space L
i
2
-
a
i
,
a
i with inner product f
,
g
=
∫
-
a
i
a
i
1
-
x
i
/
a
i
t
i
+
1
+
x
i
/
a
i
t
i
-
f
x
i
g
∗
x
i
d
x
i. The edge behaviour in (2) is strictly related to the asymptotic behaviour of the Fourier transform of the current itself with respect to x
i (J
~
i
z
u). Indeed, by means of Watson’s lemma [25] it can be stated that(4)
J
~
i
z
u
~
u
→
+
∞
J
~
i
z
∞
u
=
η
i
-
e
-
j
u
a
i
u
t
i
-
+
1
+
η
i
+
e
j
u
a
i
u
t
i
+
+
1
,where η
i
± are suitable parameters depending on the geometry of the problem and the incident field.Hence, as shown in [13], the following spectral domain representation for the longitudinal electric field can be obtained by invoking the superposition principle:(5)
E
z
x
,
y
=
∑
i
=
1
L
∫
-
∞
+
∞
J
~
i
z
u
g
i
u
,
x
,
y
d
u
,where (6)
g
i
u
,
x
,
y
=
ω
μ
2
j
e
-
j
u
x
i
x
,
y
-
y
i
x
,
y
u
2
-
k
2
u
2
-
k
2
.By imposing the total electric field to be vanishing on the cylinder surface, an EFIE is obtained(7)
E
z
x
,
y
y
j
=
0
=
-
E
0
e
-
j
k
x
x
+
k
y
y
y
j
=
0with x
j
≤
a
j and j
=
1,2
,
…
L.The obtained integral equation is discretized by means of Galerkin’s method. In order to achieve fast convergence, the functionJ
-
i
z
x
i is expanded in series of Jacobi polynomials P
n
t
i
+
,
t
i
-
x
i
/
a
i which form an orthogonal basis in the weighted Hilbert space L
i
2
-
a
i
,
a
i. Hence, the following expansion can be established for the longitudinal component of the surface current on the i-th side:(8)
J
i
z
x
i
=
∑
n
=
0
+
∞
J
n
i
φ
n
t
i
+
,
t
i
-
x
i
a
i
,where (9a)
φ
n
α
,
β
x
a
=
1
-
x
a
α
1
+
x
a
β
P
n
α
,
β
x
/
a
ξ
n
α
,
β
Π
x
a
,
(9b)
ξ
n
α
,
β
=
∫
-
a
a
1
-
x
a
α
1
+
x
a
β
P
n
α
,
β
x
a
2
d
x
=
a
2
α
+
β
+
1
Γ
n
+
α
+
1
Γ
n
+
β
+
1
n
!
2
n
+
α
+
β
+
1
Γ
n
+
α
+
β
+
1are suitable normalization quantities, Π
· is the unitary rectangular window, and Γ
· denotes the Gamma function [26]. Moreover, the Fourier transform of the selected expansion functions can be expressed in closed-form in terms of confluent hypergeometric functions of the first kind F
1
1
·
;
·
;
· [26], i.e.,(10a)
φ
~
n
α
,
β
a
u
=
1
2
π
∫
-
∞
∞
φ
n
α
,
β
x
a
e
j
u
x
d
x
=
ξ
-
n
α
,
β
2
j
u
a
n
e
-
j
u
a
F
1
1
n
+
β
+
1
;
2
n
+
α
+
β
+
2
;
2
j
u
a
,
(10b)
ξ
-
n
α
,
β
=
a
2
α
+
β
B
n
+
α
+
1
,
n
+
β
+
1
π
n
!
ξ
n
α
,
β
.As shown in [15], by means of reciprocity, it is always possible to individuate a representation of the matrix coefficients such that the convolution integrals can be always interpreted as the Fourier transform in the complex plane of the expansion functions. Hence, the general element of the scattering matrix can be written as(11)
I
n
,
m
i
,
j
=
∫
-
∞
+
∞
f
^
n
,
m
i
,
j
u
d
uwhere(12a)
f
^
n
,
m
i
,
j
u
=
g
i
u
,
x
-
j
,
y
-
j
φ
~
n
t
i
+
,
t
i
-
a
i
u
φ
~
m
t
j
+
,
t
j
-
-
a
j
G
i
j
u
∗
,
(12b)
G
i
j
u
=
u
cos
⁡
φ
i
-
φ
j
+
j
u
2
-
k
2
sgn
⁡
y
i
x
-
j
,
y
-
j
sin
⁡
φ
i
-
φ
j
.Remembering the asymptotic behaviour of the confluent hypergeometric function of the first kind [26], (13)
F
1
1
a
;
b
;
z
~
z
→
+
∞
e
±
j
π
a
z
-
a
Γ
b
Γ
b
-
a
+
e
z
z
a
-
b
Γ
b
Γ
a
,where z
=
u
+
j
v is a complex number, the upper sign has to be chosen when -
π
/
2
≤
arg
⁡
(
z
)
<
3
π
/
2 and the lower sign in the case -
3
π
/
2
<
arg
⁡
(
z
)
≤
π
/
2, and it is simple to obtain the following oscillating asymptotic behaviour for the integrand in (11):(14)
f
^
n
,
m
i
,
j
u
~
u
→
+
∞
∑
t
=
0
1
e
-
β
t
i
,
j
u
∑
s
=
0
1
c
n
,
m
,
s
,
t
i
,
j
e
j
α
s
,
t
i
,
j
u
u
t
i
s
+
t
j
t
+
3
,where c
n
,
m
,
s
,
t
i
,
j are suitable parameters, t
l
0
=
t
l
+, with l
∈
i
,
j,(15a)
α
s
,
t
i
,
j
=
-
-
1
s
a
i
+
-
1
t
a
j
cos
⁡
φ
i
-
φ
j
-
x
i
x
-
i
,
y
-
i
,
(15b)
β
t
i
,
j
=
-
1
t
a
j
sgn
⁡
y
i
x
-
j
,
y
-
j
sin
⁡
φ
i
-
φ
j
+
y
i
x
-
j
,
y
-
j
≥
0
,where s
g
n
· is the signum function.It is interesting to observe that(16)
min
⁡
β
t
i
,
j
=
-
a
j
sin
⁡
φ
i
-
φ
j
+
y
i
x
-
j
,
y
-
jis the distance between the ordinate of the vertex of the j-th side closer to the i-th side in the coordinate system of the i-th side. It is not difficult to understand that when min
⁡
β
t
i
,
j
~
0, e.g., the i-th and the j-th sides are adjacent/almost adjacent or coplanar/almost coplanar sides, the integral in (11) is slowly convergent. In order to speed up the convergence of such kind of integrals, CAAAT can be applied; i.e., it is possible to extract the asymptotic behaviour of their integrands and express the integrals of the extracted contribution in closed form. Unfortunately, such a technique allows us to obtain faster decaying integrands without overcoming the most important problem of their oscillating nature.


### 2.2. Proposed Solution

A new analytical technique is now introduced in order to speed up the numerical evaluation of the slowly converging scattering matrix coefficients.In the complex planez
=
u
+
j
v, the following relation can be established [26]:(17)
F
1
1
a
;
b
;
z
=
Γ
b
e
±
j
π
a
Γ
b
-
a
U
a
;
b
;
z
+
Γ
b
e
z
e
∓
j
π
b
-
a
Γ
a
U
b
-
a
;
b
;
-
z
,where U
·
;
·
;
· is the confluent hypergeometric function of the second kind [26] and the upper sign has to be chosen when -
π
/
2
≤
arg
⁡
(
z
)
<
3
π
/
2 while the lower sign in the case -
3
π
/
2
<
arg
⁡
(
z
)
≤
π
/
2.Therefore, by posing(18)
ψ
~
n
α
,
β
z
=
ξ
-
n
α
,
β
Γ
2
n
+
α
+
β
+
1
Γ
n
+
α
+
1
e
±
j
π
n
+
β
+
1
2
j
z
n
e
-
j
z
U
n
+
β
+
1
;
2
n
+
α
+
β
+
2
;
2
j
z
,where the upper sign has to be chosen when 0
≤
arg
⁡
(
z
)
<
2
π while the lower sign in the case -
π
<
arg
⁡
(
z
)
≤
π, it is possible to write(19)
φ
~
n
α
,
β
z
=
ψ
~
n
α
,
β
z
+
-
1
n
ψ
~
n
β
,
α
-
z
.It is interesting to observe thatψ
~
n
α
,
β
z is a discontinuous function for R
z
=
0. Moreover, since [26](20)
U
a
,
b
,
z
=
π
sin
⁡
π
b
F
1
1
a
,
b
,
z
Γ
1
+
a
-
b
Γ
b
-
z
1
-
b
F
1
1
1
+
a
-
b
,
2
-
b
,
z
Γ
a
Γ
2
-
b
,it is simple to state that ψ
~
n
α
,
β
z
~
z
→
0
d
n
α
,
β
/
z
2
n
+
α
+
β
+
1, where d
n
α
,
β is a suitable parameter. Then, ψ
~
n
α
,
β
z is singular for z
=
0 (except when α
=
β
=
-
1
/
2 and n
=
0 simultaneously) and it is a many-valued function with a branch-cut for I
z
>
0 if α
+
β is a noninteger number.Therefore, the general matrix coefficient can be rewritten as(21)
I
n
,
m
i
,
j
=
∫
-
r
+
r
f
^
n
,
m
i
,
j
u
d
u
+
∑
w
,
s
,
t
=
0
1
∫
r
+
∞
f
n
,
m
,
w
,
s
,
t
i
,
j
u
d
uwhere (22)
f
n
,
m
,
w
,
s
,
t
i
,
j
u
=
-
1
n
s
+
m
t
g
i
-
1
w
u
,
x
-
j
,
y
-
j
ψ
~
n
t
i
s
,
t
i
s
-
1
-
1
w
+
s
a
i
u
ψ
~
m
t
j
t
,
t
j
t
-
1
-
-
1
t
a
j
G
i
j
-
1
w
u
∗
,it is simple to show that the integration path of the improper integrals in (21) does not intersect the branch-cuts of the integrands, and it has been chosen r
>
k in order to avoid the possible singularities of f
n
,
m
,
w
,
s
,
t
i
,
j
z which can only be located on the real axis for -
k
≤
u
≤
k (see Figure 2, just for an example).Figure 2
Integration path.On the other hand, since [26](23)
U
a
;
b
;
z
~
z
→
+
∞
z
-
a
,for -
3
π
/
2
<
arg
⁡
(
z
)
<
3
π
/
2, the many-valued function f
n
,
m
,
w
,
s
,
t
i
,
j
z can be rewritten on its principal sheet as(24)
f
n
,
m
,
w
,
s
,
t
i
,
j
z
=
f
-
n
,
m
,
w
,
s
,
t
i
,
j
z
e
j
α
s
,
t
i
,
j
-
1
w
z
-
β
t
i
,
j
z
2
-
k
2
,where f
-
n
,
m
,
w
,
s
,
t
i
,
j
z
~
z
→
∞
c
-
n
,
m
,
w
,
s
,
t
i
,
j
/
z
γ
n
,
m
,
w
,
s
,
t
i
,
j with γ
n
,
m
,
w
,
s
,
t
i
,
j
≥
2 and c
-
n
,
m
,
w
,
s
,
t
i
,
j is a suitable parameter.Hence, the path in the complex plane along with the exponential function in (24) does not oscillate and, simultaneously, f
n
,
m
,
w
,
s
,
t
i
,
j
z is an integrable function can be obtained by posing (25)
α
u
-
β
I
R
z
=
0
α
v
+
β
R
R
z
≥
0
,where the positions α
=
α
s
,
t
i
,
j
-
1
w and β
=
β
t
i
,
j have been introduced in order to simplify the notation, while R
· and I
· denote the real part and the imaginary part of a complex number, respectively.Forα
=
β
=
0, formulas (25) are always verified. For α
≠
0 and β
=
0 the integration path defined by (25) is u
=
0, sgn
⁡
α
v
≥
0, while for α
=
0 and β
>
0 the integration path is given by u
≥
k, v
=
0.Now, let us consider the more general case forα
≠
0 and β
>
0. By means of algebraic manipulations and supposing to consider u
≥
0, (25) can be reduced to(26)
u
2
β
2
-
v
2
α
2
-
k
2
α
2
+
β
2
=
0
sgn
⁡
α
v
≥
0
,which leads to the following expressions of v and u as a function of u and v, respectively:(27a)
v
=
α
u
2
β
2
-
k
2
α
2
+
β
2
for
u
≥
β
k
α
2
+
β
2
,
(27b)
u
=
β
v
2
α
2
+
k
2
α
2
+
β
2
for
sgn
⁡
α
v
≥
0
.It is worth noting that, from a numerical point of view, (27a) should be preferred for α
/
β
~
0 while (27b) is more suitable for α
/
β
~
∞. Hence, for α
,
β
>
0, just for an example, a possible alternative integration path is the one in Figure 3.Figure 3
Alternative integration path.Now, let us consider the closed contourC sketched in Figure 4(a) nonintersecting any of the branch-cuts of f
n
,
m
,
w
,
s
,
t
i
,
j
z.Figure 4
Integration contours in the complex plane. (a) Integration contour nonintersecting the branch-cuts of the integrand. (b) Integration contour intersecting one of the branch-cuts of the integrand.
(a)
(b)By means of Cauchy’s integral theorem and Jordan’s lemma, it is possible to write(28)
lim
R
→
∞
⁡
∮
C
f
n
,
m
,
w
,
s
,
t
i
,
j
z
d
z
=
0
,and then, the general improper integral at the right-hand side of (21) can be rewritten as (29a)
∫
r
+
∞
f
n
,
m
,
w
,
s
,
t
i
,
j
u
d
u
=
j
∫
0
v
-
f
n
,
m
,
w
,
s
,
t
i
,
j
r
+
j
v
d
v
+
∫
r
u
-
f
n
,
m
,
w
,
s
,
t
i
,
j
u
+
j
v
-
d
u
+
I
-
n
,
m
,
w
,
s
,
t
i
,
j
(29b)
I
-
n
,
m
,
w
,
s
,
t
i
,
j
=
∫
u
-
+
∞
f
n
,
m
,
w
,
s
,
t
i
,
j
u
+
j
α
u
2
β
2
-
k
2
α
2
+
β
2
1
+
j
α
/
β
2
u
u
2
/
β
2
-
k
2
/
α
2
+
β
2
d
u
for
α
β
≤
1
∫
v
-
+
∞
f
n
,
m
,
w
,
s
,
t
i
,
j
β
v
2
α
2
+
k
2
α
2
+
β
2
+
j
v
β
/
α
2
v
v
2
/
α
2
+
k
2
/
α
2
+
β
2
+
j
d
v
for
α
β
≥
1where sgn
⁡
α
v
-
>
0 and u
-
=
β
v
-
2
/
α
2
+
k
2
/
α
2
+
β
2, which is a linear combination of proper integrals and an improper integral of an asymptotically nonoscillating function.Unfortunately, the integration contour can intersect one of the branch-cuts off
n
,
m
,
w
,
s
,
t
i
,
j
z along the horizontal/vertical path of the alternative integration path (see Figure 4(b)). In such a case, one of the confluent hypergeometric functions of the second kind of the integrand of the integral in (29b) takes values on its upper/lower secondary sheet. However [26],(30)
U
a
,
b
,
z
e
j
2
π
M
=
1
-
e
-
j
2
π
M
b
Γ
1
-
b
Γ
1
+
a
-
b
F
1
1
a
,
b
,
z
+
e
-
j
2
π
M
b
U
a
,
b
,
zfor -
π
<
arg
⁡
z
≤
π, where M
=
-
1,0
,
1 identifies the lower secondary sheet, the principal sheet, and the upper secondary sheet of U
·
;
·
;
·, respectively. Hence, it is not difficult to conclude that the integrand of the integral in (29b) is an asymptotically oscillating function for M
=
±
1. Such a problem is completely overcome by introducing the following generalization of the representation (19)(31)
φ
~
n
α
,
β
z
=
χ
~
n
,
M
+
,
M
-
α
,
β
z
+
-
1
n
χ
~
n
,
M
-
,
M
+
β
,
α
-
z
,being(32)
χ
~
n
,
M
+
,
M
-
α
,
β
z
=
e
j
2
π
M
-
n
+
β
+
1
ψ
~
n
α
,
β
z
+
1
-
e
j
2
π
M
+
n
+
α
+
1
-
1
n
ψ
~
n
β
,
α
-
z
,where the parameters M
±
∈
-
1,0
,
+
1 take into account the possible intersection of the integration contour in Figure 4 with one of the branch-cuts of f
n
,
m
,
w
,
s
,
t
i
,
j
z.


## 2.1. Background

In Figure1, a polygonal cross-section PEC cylinder is sketched. A coordinate system x
,
y
,
z is introduced such that the z axis coincides with the cylinder axis. The L sides of the polygonal cross-section are numbered clockwise and a local coordinate system x
i
,
y
i
,
z is introduced on the i-t
h side with the origin at the centre of the side itself in the position x
=
x
-
i
,
y
=
y
-
i and the y
i axis oriented in the outward direction. The angle ϕ
i
∈ ]−π, π] denotes the orientation of the x
i axis with respect to the x axis and 2
a
i denotes the length of the i-th side. A TM polarized plane wave impinges onto the cylinder surface with an angle ϕ with respect to the x axis and orthogonally with respect to the cylinder axis, namely,(1)
E
_
i
n
c
x
,
y
=
E
0
z
^
e
-
j
k
x
x
+
k
y
y
,where k
x
=
-
k
cos
⁡
ϕ and k
y
=
-
k
sin
⁡
ϕ, k
=
2
π
/
λ
=
ω
ε
μ is the wavenumber, λ is the wavelength, ε is the dielectric permittivity and μ the magnetic permeability of the external medium, ω is the angular frequency. With such a choice, only TM solutions can be obtained; i.e., the induced current is longitudinal and the electromagnetic field is invariant along the z axis.Figure 1
Geometry of the problem.By means of Meixner’s theory [24], the following edge behaviour can be established for the longitudinal current density on the i-th surface:(2)
J
i
z
x
i
=
1
-
x
i
a
i
t
i
+
1
+
x
i
a
i
t
i
-
J
-
i
z
x
iwith i
=
1,2
,
…
L, where (3)
t
i
±
=
ψ
i
±
-
π
2
π
-
ψ
i
±
ψ
i
±
≤
3
π
2
1
ψ
i
±
≥
3
π
2
,
ψ
i
± being the angles of the wedges at the abscissas x
i
=
±
a
i and the well-behaved function J
-
i
z
x
i belongs to the weighted Hilbert space L
i
2
-
a
i
,
a
i with inner product f
,
g
=
∫
-
a
i
a
i
1
-
x
i
/
a
i
t
i
+
1
+
x
i
/
a
i
t
i
-
f
x
i
g
∗
x
i
d
x
i. The edge behaviour in (2) is strictly related to the asymptotic behaviour of the Fourier transform of the current itself with respect to x
i (J
~
i
z
u). Indeed, by means of Watson’s lemma [25] it can be stated that(4)
J
~
i
z
u
~
u
→
+
∞
J
~
i
z
∞
u
=
η
i
-
e
-
j
u
a
i
u
t
i
-
+
1
+
η
i
+
e
j
u
a
i
u
t
i
+
+
1
,where η
i
± are suitable parameters depending on the geometry of the problem and the incident field.Hence, as shown in [13], the following spectral domain representation for the longitudinal electric field can be obtained by invoking the superposition principle:(5)
E
z
x
,
y
=
∑
i
=
1
L
∫
-
∞
+
∞
J
~
i
z
u
g
i
u
,
x
,
y
d
u
,where (6)
g
i
u
,
x
,
y
=
ω
μ
2
j
e
-
j
u
x
i
x
,
y
-
y
i
x
,
y
u
2
-
k
2
u
2
-
k
2
.By imposing the total electric field to be vanishing on the cylinder surface, an EFIE is obtained(7)
E
z
x
,
y
y
j
=
0
=
-
E
0
e
-
j
k
x
x
+
k
y
y
y
j
=
0with x
j
≤
a
j and j
=
1,2
,
…
L.The obtained integral equation is discretized by means of Galerkin’s method. In order to achieve fast convergence, the functionJ
-
i
z
x
i is expanded in series of Jacobi polynomials P
n
t
i
+
,
t
i
-
x
i
/
a
i which form an orthogonal basis in the weighted Hilbert space L
i
2
-
a
i
,
a
i. Hence, the following expansion can be established for the longitudinal component of the surface current on the i-th side:(8)
J
i
z
x
i
=
∑
n
=
0
+
∞
J
n
i
φ
n
t
i
+
,
t
i
-
x
i
a
i
,where (9a)
φ
n
α
,
β
x
a
=
1
-
x
a
α
1
+
x
a
β
P
n
α
,
β
x
/
a
ξ
n
α
,
β
Π
x
a
,
(9b)
ξ
n
α
,
β
=
∫
-
a
a
1
-
x
a
α
1
+
x
a
β
P
n
α
,
β
x
a
2
d
x
=
a
2
α
+
β
+
1
Γ
n
+
α
+
1
Γ
n
+
β
+
1
n
!
2
n
+
α
+
β
+
1
Γ
n
+
α
+
β
+
1are suitable normalization quantities, Π
· is the unitary rectangular window, and Γ
· denotes the Gamma function [26]. Moreover, the Fourier transform of the selected expansion functions can be expressed in closed-form in terms of confluent hypergeometric functions of the first kind F
1
1
·
;
·
;
· [26], i.e.,(10a)
φ
~
n
α
,
β
a
u
=
1
2
π
∫
-
∞
∞
φ
n
α
,
β
x
a
e
j
u
x
d
x
=
ξ
-
n
α
,
β
2
j
u
a
n
e
-
j
u
a
F
1
1
n
+
β
+
1
;
2
n
+
α
+
β
+
2
;
2
j
u
a
,
(10b)
ξ
-
n
α
,
β
=
a
2
α
+
β
B
n
+
α
+
1
,
n
+
β
+
1
π
n
!
ξ
n
α
,
β
.As shown in [15], by means of reciprocity, it is always possible to individuate a representation of the matrix coefficients such that the convolution integrals can be always interpreted as the Fourier transform in the complex plane of the expansion functions. Hence, the general element of the scattering matrix can be written as(11)
I
n
,
m
i
,
j
=
∫
-
∞
+
∞
f
^
n
,
m
i
,
j
u
d
uwhere(12a)
f
^
n
,
m
i
,
j
u
=
g
i
u
,
x
-
j
,
y
-
j
φ
~
n
t
i
+
,
t
i
-
a
i
u
φ
~
m
t
j
+
,
t
j
-
-
a
j
G
i
j
u
∗
,
(12b)
G
i
j
u
=
u
cos
⁡
φ
i
-
φ
j
+
j
u
2
-
k
2
sgn
⁡
y
i
x
-
j
,
y
-
j
sin
⁡
φ
i
-
φ
j
.Remembering the asymptotic behaviour of the confluent hypergeometric function of the first kind [26], (13)
F
1
1
a
;
b
;
z
~
z
→
+
∞
e
±
j
π
a
z
-
a
Γ
b
Γ
b
-
a
+
e
z
z
a
-
b
Γ
b
Γ
a
,where z
=
u
+
j
v is a complex number, the upper sign has to be chosen when -
π
/
2
≤
arg
⁡
(
z
)
<
3
π
/
2 and the lower sign in the case -
3
π
/
2
<
arg
⁡
(
z
)
≤
π
/
2, and it is simple to obtain the following oscillating asymptotic behaviour for the integrand in (11):(14)
f
^
n
,
m
i
,
j
u
~
u
→
+
∞
∑
t
=
0
1
e
-
β
t
i
,
j
u
∑
s
=
0
1
c
n
,
m
,
s
,
t
i
,
j
e
j
α
s
,
t
i
,
j
u
u
t
i
s
+
t
j
t
+
3
,where c
n
,
m
,
s
,
t
i
,
j are suitable parameters, t
l
0
=
t
l
+, with l
∈
i
,
j,(15a)
α
s
,
t
i
,
j
=
-
-
1
s
a
i
+
-
1
t
a
j
cos
⁡
φ
i
-
φ
j
-
x
i
x
-
i
,
y
-
i
,
(15b)
β
t
i
,
j
=
-
1
t
a
j
sgn
⁡
y
i
x
-
j
,
y
-
j
sin
⁡
φ
i
-
φ
j
+
y
i
x
-
j
,
y
-
j
≥
0
,where s
g
n
· is the signum function.It is interesting to observe that(16)
min
⁡
β
t
i
,
j
=
-
a
j
sin
⁡
φ
i
-
φ
j
+
y
i
x
-
j
,
y
-
jis the distance between the ordinate of the vertex of the j-th side closer to the i-th side in the coordinate system of the i-th side. It is not difficult to understand that when min
⁡
β
t
i
,
j
~
0, e.g., the i-th and the j-th sides are adjacent/almost adjacent or coplanar/almost coplanar sides, the integral in (11) is slowly convergent. In order to speed up the convergence of such kind of integrals, CAAAT can be applied; i.e., it is possible to extract the asymptotic behaviour of their integrands and express the integrals of the extracted contribution in closed form. Unfortunately, such a technique allows us to obtain faster decaying integrands without overcoming the most important problem of their oscillating nature.


## 2.2. Proposed Solution

A new analytical technique is now introduced in order to speed up the numerical evaluation of the slowly converging scattering matrix coefficients.In the complex planez
=
u
+
j
v, the following relation can be established [26]:(17)
F
1
1
a
;
b
;
z
=
Γ
b
e
±
j
π
a
Γ
b
-
a
U
a
;
b
;
z
+
Γ
b
e
z
e
∓
j
π
b
-
a
Γ
a
U
b
-
a
;
b
;
-
z
,where U
·
;
·
;
· is the confluent hypergeometric function of the second kind [26] and the upper sign has to be chosen when -
π
/
2
≤
arg
⁡
(
z
)
<
3
π
/
2 while the lower sign in the case -
3
π
/
2
<
arg
⁡
(
z
)
≤
π
/
2.Therefore, by posing(18)
ψ
~
n
α
,
β
z
=
ξ
-
n
α
,
β
Γ
2
n
+
α
+
β
+
1
Γ
n
+
α
+
1
e
±
j
π
n
+
β
+
1
2
j
z
n
e
-
j
z
U
n
+
β
+
1
;
2
n
+
α
+
β
+
2
;
2
j
z
,where the upper sign has to be chosen when 0
≤
arg
⁡
(
z
)
<
2
π while the lower sign in the case -
π
<
arg
⁡
(
z
)
≤
π, it is possible to write(19)
φ
~
n
α
,
β
z
=
ψ
~
n
α
,
β
z
+
-
1
n
ψ
~
n
β
,
α
-
z
.It is interesting to observe thatψ
~
n
α
,
β
z is a discontinuous function for R
z
=
0. Moreover, since [26](20)
U
a
,
b
,
z
=
π
sin
⁡
π
b
F
1
1
a
,
b
,
z
Γ
1
+
a
-
b
Γ
b
-
z
1
-
b
F
1
1
1
+
a
-
b
,
2
-
b
,
z
Γ
a
Γ
2
-
b
,it is simple to state that ψ
~
n
α
,
β
z
~
z
→
0
d
n
α
,
β
/
z
2
n
+
α
+
β
+
1, where d
n
α
,
β is a suitable parameter. Then, ψ
~
n
α
,
β
z is singular for z
=
0 (except when α
=
β
=
-
1
/
2 and n
=
0 simultaneously) and it is a many-valued function with a branch-cut for I
z
>
0 if α
+
β is a noninteger number.Therefore, the general matrix coefficient can be rewritten as(21)
I
n
,
m
i
,
j
=
∫
-
r
+
r
f
^
n
,
m
i
,
j
u
d
u
+
∑
w
,
s
,
t
=
0
1
∫
r
+
∞
f
n
,
m
,
w
,
s
,
t
i
,
j
u
d
uwhere (22)
f
n
,
m
,
w
,
s
,
t
i
,
j
u
=
-
1
n
s
+
m
t
g
i
-
1
w
u
,
x
-
j
,
y
-
j
ψ
~
n
t
i
s
,
t
i
s
-
1
-
1
w
+
s
a
i
u
ψ
~
m
t
j
t
,
t
j
t
-
1
-
-
1
t
a
j
G
i
j
-
1
w
u
∗
,it is simple to show that the integration path of the improper integrals in (21) does not intersect the branch-cuts of the integrands, and it has been chosen r
>
k in order to avoid the possible singularities of f
n
,
m
,
w
,
s
,
t
i
,
j
z which can only be located on the real axis for -
k
≤
u
≤
k (see Figure 2, just for an example).Figure 2
Integration path.On the other hand, since [26](23)
U
a
;
b
;
z
~
z
→
+
∞
z
-
a
,for -
3
π
/
2
<
arg
⁡
(
z
)
<
3
π
/
2, the many-valued function f
n
,
m
,
w
,
s
,
t
i
,
j
z can be rewritten on its principal sheet as(24)
f
n
,
m
,
w
,
s
,
t
i
,
j
z
=
f
-
n
,
m
,
w
,
s
,
t
i
,
j
z
e
j
α
s
,
t
i
,
j
-
1
w
z
-
β
t
i
,
j
z
2
-
k
2
,where f
-
n
,
m
,
w
,
s
,
t
i
,
j
z
~
z
→
∞
c
-
n
,
m
,
w
,
s
,
t
i
,
j
/
z
γ
n
,
m
,
w
,
s
,
t
i
,
j with γ
n
,
m
,
w
,
s
,
t
i
,
j
≥
2 and c
-
n
,
m
,
w
,
s
,
t
i
,
j is a suitable parameter.Hence, the path in the complex plane along with the exponential function in (24) does not oscillate and, simultaneously, f
n
,
m
,
w
,
s
,
t
i
,
j
z is an integrable function can be obtained by posing (25)
α
u
-
β
I
R
z
=
0
α
v
+
β
R
R
z
≥
0
,where the positions α
=
α
s
,
t
i
,
j
-
1
w and β
=
β
t
i
,
j have been introduced in order to simplify the notation, while R
· and I
· denote the real part and the imaginary part of a complex number, respectively.Forα
=
β
=
0, formulas (25) are always verified. For α
≠
0 and β
=
0 the integration path defined by (25) is u
=
0, sgn
⁡
α
v
≥
0, while for α
=
0 and β
>
0 the integration path is given by u
≥
k, v
=
0.Now, let us consider the more general case forα
≠
0 and β
>
0. By means of algebraic manipulations and supposing to consider u
≥
0, (25) can be reduced to(26)
u
2
β
2
-
v
2
α
2
-
k
2
α
2
+
β
2
=
0
sgn
⁡
α
v
≥
0
,which leads to the following expressions of v and u as a function of u and v, respectively:(27a)
v
=
α
u
2
β
2
-
k
2
α
2
+
β
2
for
u
≥
β
k
α
2
+
β
2
,
(27b)
u
=
β
v
2
α
2
+
k
2
α
2
+
β
2
for
sgn
⁡
α
v
≥
0
.It is worth noting that, from a numerical point of view, (27a) should be preferred for α
/
β
~
0 while (27b) is more suitable for α
/
β
~
∞. Hence, for α
,
β
>
0, just for an example, a possible alternative integration path is the one in Figure 3.Figure 3
Alternative integration path.Now, let us consider the closed contourC sketched in Figure 4(a) nonintersecting any of the branch-cuts of f
n
,
m
,
w
,
s
,
t
i
,
j
z.Figure 4
Integration contours in the complex plane. (a) Integration contour nonintersecting the branch-cuts of the integrand. (b) Integration contour intersecting one of the branch-cuts of the integrand.
(a)
(b)By means of Cauchy’s integral theorem and Jordan’s lemma, it is possible to write(28)
lim
R
→
∞
⁡
∮
C
f
n
,
m
,
w
,
s
,
t
i
,
j
z
d
z
=
0
,and then, the general improper integral at the right-hand side of (21) can be rewritten as (29a)
∫
r
+
∞
f
n
,
m
,
w
,
s
,
t
i
,
j
u
d
u
=
j
∫
0
v
-
f
n
,
m
,
w
,
s
,
t
i
,
j
r
+
j
v
d
v
+
∫
r
u
-
f
n
,
m
,
w
,
s
,
t
i
,
j
u
+
j
v
-
d
u
+
I
-
n
,
m
,
w
,
s
,
t
i
,
j
(29b)
I
-
n
,
m
,
w
,
s
,
t
i
,
j
=
∫
u
-
+
∞
f
n
,
m
,
w
,
s
,
t
i
,
j
u
+
j
α
u
2
β
2
-
k
2
α
2
+
β
2
1
+
j
α
/
β
2
u
u
2
/
β
2
-
k
2
/
α
2
+
β
2
d
u
for
α
β
≤
1
∫
v
-
+
∞
f
n
,
m
,
w
,
s
,
t
i
,
j
β
v
2
α
2
+
k
2
α
2
+
β
2
+
j
v
β
/
α
2
v
v
2
/
α
2
+
k
2
/
α
2
+
β
2
+
j
d
v
for
α
β
≥
1where sgn
⁡
α
v
-
>
0 and u
-
=
β
v
-
2
/
α
2
+
k
2
/
α
2
+
β
2, which is a linear combination of proper integrals and an improper integral of an asymptotically nonoscillating function.Unfortunately, the integration contour can intersect one of the branch-cuts off
n
,
m
,
w
,
s
,
t
i
,
j
z along the horizontal/vertical path of the alternative integration path (see Figure 4(b)). In such a case, one of the confluent hypergeometric functions of the second kind of the integrand of the integral in (29b) takes values on its upper/lower secondary sheet. However [26],(30)
U
a
,
b
,
z
e
j
2
π
M
=
1
-
e
-
j
2
π
M
b
Γ
1
-
b
Γ
1
+
a
-
b
F
1
1
a
,
b
,
z
+
e
-
j
2
π
M
b
U
a
,
b
,
zfor -
π
<
arg
⁡
z
≤
π, where M
=
-
1,0
,
1 identifies the lower secondary sheet, the principal sheet, and the upper secondary sheet of U
·
;
·
;
·, respectively. Hence, it is not difficult to conclude that the integrand of the integral in (29b) is an asymptotically oscillating function for M
=
±
1. Such a problem is completely overcome by introducing the following generalization of the representation (19)(31)
φ
~
n
α
,
β
z
=
χ
~
n
,
M
+
,
M
-
α
,
β
z
+
-
1
n
χ
~
n
,
M
-
,
M
+
β
,
α
-
z
,being(32)
χ
~
n
,
M
+
,
M
-
α
,
β
z
=
e
j
2
π
M
-
n
+
β
+
1
ψ
~
n
α
,
β
z
+
1
-
e
j
2
π
M
+
n
+
α
+
1
-
1
n
ψ
~
n
β
,
α
-
z
,where the parameters M
±
∈
-
1,0
,
+
1 take into account the possible intersection of the integration contour in Figure 4 with one of the branch-cuts of f
n
,
m
,
w
,
s
,
t
i
,
j
z.


## 3. Results and Discussion

The aim of this section is to show the efficiency of the presented technique even by means of comparisons with the CAAAT.The following normalized truncation error is introduced:(33)
e
r
r
N
=
J
N
+
1
-
J
N
J
N
,where · is the usual Euclidean norm and J
M is the vector of all the expansion coefficients of the longitudinal currents on all the sides evaluated with M expansion functions on each side. The simulations are performed on a laptop equipped with an Intel Core 2 Duo CPU T9600 2.8-GHz 3-GB RAM, running Windows XP and the integrals evaluated by means of a Gauss-Legendre quadrature routine. In this way, for the proposed examples, about 70 integrals per second can be numerically evaluated.First of all, the analysis of the scattering from an irregular triangular cross-section cylinder of sidesA
B
¯
=
a
1
+
1
/
3, B
C
¯
=
2
a
/
3, C
A
¯
=
a
2 for a
=
λ
/
2
,
λ
,
2
λ when a TM polarized plane wave impinges with ϕ
=
π
/
6 and E
0
=
1
V
/
m is performed. Such a configuration involves only adjacent sides for which, as underlined above, the integrands of the integrals of the coefficients’ matrix are slowly decaying functions. In Figure 5(a), the normalized truncation error by varying the number of expansion functions used on each side of the cylinder cross-section (N) is plotted revealing a very fast convergence in all the examined cases. As a matter of fact, a normalized truncation error of 1
0
-
1
,
1
0
-
2
,
1
0
-
3 is achieved for N ranging from 4, 6, 13 to 11, 14, and 17, respectively. It is interesting to note that less than 35 seconds are needed to fill the coefficients’ matrix for N
=
20. For the sake of completeness, in Figure 5(b) the reconstructed surface current densities for all the examined cases are plotted. In Figure 5(c) the CPU time ratio, i.e., the ratio between the computation time needed to reconstruct the solution as obtained by using the CAAAT with respect to the presented technique, is reported for all the examined cases as a function of N. As clearly shown, the presented technique always outperforms the CAAAT. Moreover, the CPU time ratio quickly increases as higher is N, i.e., as higher is the accuracy required for the solution.Figure 5
TM scattering from an irregular triangular cross-section cylinder withA
B
¯
=
a
1
+
1
/
3, B
C
¯
=
2
a
/
3, C
A
¯
=
a
2, and for a
=
λ
/
2
,
λ
,
2
λ, ϕ
=
π
/
6, and E
0
=
1
V
/
m. (a) Normalized truncation error. (b) Surface current density. (c) CPU time ratio.
(a)
(b)
(c)In a second example, the scattering from two identical irregular triangular cross-section cylinders of sidesA
B
¯
=
D
E
¯
=
λ
1
+
1
/
3, B
C
¯
=
E
F
¯
=
2
λ
/
3, and C
A
¯
=
F
D
¯
=
λ
2 with A
B
¯ and D
E
¯ coplanar sides and B
D
¯
=
λ
/
4
,
λ
/
2
,
λ when a TM polarized plane wave impinges with ϕ
=
π
/
6 and E
0
=
1
V
/
m is examined. This is an even more troublesome configuration because it involves even coplanar sides for which the integrands of the corresponding integrals of the coefficients’ matrix are slowly decaying functions which oscillate faster and faster as greater is the distance between the coplanar sides (see formula (16)). In Figures 6(a) and 6(b), the normalized truncation error and the surface current density are plotted, respectively. The convergence is very fast and substantially independent of the distance between the cylinders. An accuracy for the solution of 1
0
-
1
,
1
0
-
2
,
1
0
-
3 is achieved for N
=
7,9
,
14, respectively. It is interesting to note that less than 100 seconds are needed to fill the coefficients’ matrix for N
=
20. In Figure 6(c), the CPU time ratio is reported as a function ofN. In all the examined cases, the presented technique outperforms the CAAAT especially as higher is the distance between the cylinders and as higher is the number of expansion functions used.Figure 6
TM scattering from two identical irregular triangular cross-section cylinders withA
B
¯
=
D
E
¯
=
λ
1
+
1
/
3, B
C
¯
=
E
F
¯
=
2
λ
/
3, C
A
¯
=
F
D
¯
=
λ
2, A
B
¯ and D
E
¯ coplanar sides, and for B
D
¯
=
λ
/
4
,
λ
/
2
,
λ, ϕ
=
π
/
6, and E
0
=
1
V
/
m. (a) Normalized truncation error. (b) Surface current density. (c) CPU time ratio.
(a)
(b)
(c)In the last example, the scattering from an irregular pentagonal cross-section cylinder of sidesA
B
¯
=
D
E
¯
=
E
A
¯
=
2
λ, B
C
¯
=
λ
3, and C
D
¯
=
λ when a TM polarized plane wave impinges with ϕ
=
π
/
6 and E
0
=
1
V
/
m is analysed. As shown in Figure 7(a), the convergence is very fast even in such a case, as N
=
7,9
,
11 is enough to achieve an accuracy for the solution of 1
0
-
1
,
1
0
-
2
,
1
0
-
3, respectively. It is interesting to note that less than 70 seconds are needed to fill the coefficients’ matrix for N
=
20. Even if the correctness of the integral formulation and the discretization scheme has been widely discussed by the authors in the papers [13, 15], for the sake of completeness, a comparison with the commercial software CST Microwave Studio (CST-MWS) is shown in Figure 7(b) in terms of surface current density revealing a very good agreement. To conclude, the CPU time ratio is reported in Figure 7(c). As can be seen, the presented technique always outperforms the CAAAT. Moreover, the CPU time ratio increases by increasingN.Figure 7
TM scattering from an irregular pentagonal cross-section cylinder withA
B
¯
=
D
E
¯
=
E
A
¯
=
2
λ, B
C
¯
=
λ
3, C
D
¯
=
λ, and for ϕ
=
π
/
6 and E
0
=
1
V
/
m. (a) Normalized truncation error. (b) Surface current density. (c) CPU time ratio.
(a)
(b)
(c)


## 4. Conclusions

In recent papers, the TM scattering from polygonal cross-section PEC cylinders has been addressed by means of an EFIE formulation in the spectral domain and the method of analytical preconditioning leading to the numerical evaluation of improper integrals of asymptotically oscillating and slowly decaying functions. In this paper, a new analytical technique has been introduced in order to represent such kind of integrals as linear combinations of quickly converging integrals. The proposed numerical results show that such a technique is very effective and outperforms the CAAAT.



---

*Source: 7902836-2019-01-10.xml*