International Journal of Antennas and Propagation

Volume 2018, Article ID 8267504, 7 pages

https://doi.org/10.1155/2018/8267504

## The Application of Improved Spherical Harmonics Expansion-Based Multilevel Fast Multipole Algorithm in the Solution of Volume-Surface Integral Equation

^{1}School of Information Engineering, Communication University of China, Beijing 100024, China^{2}School of Information and Electronics, Beijing Institute of Technology, Beijing 100081, China

Correspondence should be addressed to Mang He; nc.ude.tib@gnameh

Received 27 November 2017; Accepted 1 March 2018; Published 11 April 2018

Academic Editor: Paolo Burghignoli

Copyright © 2018 Jinbo Liu et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

#### Abstract

During the solution of volume-surface integral equation (VSIE), to reduce the core memory requirement of the radiation patterns (RPs) of the basis functions, an improved spherical harmonics expansion-based multilevel fast multipole algorithm (SE-MLFMA) using the mixed-potential representation and the triangle-/tetrahedron-based grouping scheme is applied. Numerical results show that accompanying with a faster speed, the memory requirement of the RPs in the improved SE-MLFMA is several times less than that in the conventional MLFMA without compromising accuracy. A result employing the OpenMP parallelization and vector arithmetic logic unit (VALU) hardware acceleration technique is also shown to illustrate the robustness and scalability of the improved SE-MLFMA method.

#### 1. Introduction

In the computational electromagnetics, the method of moments (MoM) solution of volume-surface integral equation (VSIE) is an attractive approach, since it can be generally applied to the simulation of the arbitrary composite conductor-inhomogeneous dielectric objects [1–5]. As one of the most powerful fast solvers, the multilevel fast multipole algorithm (MLFMA) can make MoM be applied to the analysis of objects with large electrical size [4–8]. Based on the addition theorem of Green’s function and diagonalization of the translation operator, MLFMA drastically reduces the overall computational complexity from the order of *O*(*N*^{2}) to *O*(*N*log*N*) through three processes: aggregation, translation, and disaggregation. In the MLFMA, the far interactions between the well-separated groups are efficiently computed using the *k*-space integrations over the Ewald sphere. Before iterative solution, to enhance the computing efficiency, the *k*-space samples of the basis/testing functions, called radiation patterns (RPs), should be computed and stored prior. Because the quadrature sampling rate is determined by the bandwidth of the diagonalized translation matrix at the finest level, the RPs in the *k*-space are oversampled, which causes a significant amount of redundant memory usage.

To alleviate this problem, a spherical harmonics expansion-based MLFMA (SE-MLFMA) for the dyadic form of surface integral equation (SIE) was proposed [6]. With the use of symmetry, the memory requirement of RPs is significantly reduced from bytes in the conventional MLFMA to bytes in the SE-MLFMA. Here, *N* is the number of unknowns, while and are the order of multipole expansion and the SE degree of the RPs, respectively. The constant is equal to 1 in case of the electric surface field integral equation (ESFIE) or 2 in case of the combined surface field integral equation (CSFIE), and is 8 for single precision or 16 for double precision. In [6], for the RPs of the RWG [9] basis functions in the Cartesian coordinates is selected as . Further, based on the mixed-potential form of the SIE and the triangle-based grouping scheme, an improved SE-MLFMA with more efficiency was proposed [7]. In contrast to the traditional SE-MLFMA, the memory requirement of RPs is further reduced to bytes, which is proportional to the number of triangles . Besides, is eliminated because the memory cost is independent of the type of SIE. Moreover, since the RP computation is operated in the Cartesian coordinates, as the same as the aggregation, translation, and disaggregation processes, the improved SE-MLFMA totally obviates the repeated transformations between the spherical and the Cartesian coordinates, which are necessary for the traditional SE-MLFMA to eliminate the Gibbs phenomenon.

Nevertheless, the previous discussions were focused on SIE. Although the applicability of SE-MLFMA in SIE has been well verified, since the numerical characteristics and expression form of VSIE are very different from SIE, we still want to know whether this method can be effectively applied to the solution of VSIE. In this paper, when the RWG and SWG [10] basis functions are used to disperse the mixed-potential form of VSIE, how to efficiently apply the improved SE-MLFMA in the solution process is shown in detail. In contrast to the conventional MLFMA and the traditional SE-MLFMA, the improved SE-MLFMA shows the following improvements: firstly, the memory requirement of RPs is only proportional to the total number of triangles and tetrahedrons rather than that of the basis functions, leading to a considerable reduction of core memory requirement. Secondly, the memory requirement does not depend on the type of the SIE part in VSIE, which means that the application of CSFIE will not cause any memory increasing; Thirdly, since the Cartesian components of RPs are permanently used during the SE procedures, the repeated transformations between Cartesian and spherical coordinates are totally eliminated, resulting in a faster implementation.

#### 2. Formulation

##### 2.1. Derivation of the Proposed Method

Using Garlerkin’s type of MoM, the VSIE is discretized into a matrix equation as where is a impedance matrix, and and are the unknown current coefficient and generalized voltage vectors, respectively. In the MLFMA, is decomposed into two parts as where and denote the near- and far-interaction matrices between the basis and test functions, respectively. As we know, in the VSIE, can be further decomposed as

In (3), and are the number of surface and volume unknowns, respectively, is a real constant, and is the free-space intrinsic impedance. The superscript and represent the type of interaction ( for conducting surface and for dielectric volume), while the secondary subscript *E* or *H* denotes that the related entry is generated by the electric field integral equation, or that by the magnetic surface field integral equation (MSFIE). The entry in (3) can be expressed in the mixed potential form as
where is either or , is for or 1 for , and are the *j*th RWG and SWG basis functions, respectively, and is the Green function.

Since RWG/SWG basis function is defined on two adjacent triangles/tetrahedrons, each double integral in (4) consists of four subintegrals of the same type, which represent the interactions between the double two triangles/tetrahedrons involved in the basis and test functions. The final result can be efficiently obtained by adding the contributions of all related triangle/tetrahedron-to-triangle/tetrahedron integrals to the basis-to-test function entry. Thus, in the following, only the subintegral of two “positive” triangles/tetrahedrons is considered, while other integrals can be obtained by simply changing the signs. If we assume that the and are the two “positive” triangle/tetrahedron containing the *j*th and *i*th basis functions, respectively, the contribution of the interaction between and to the far-interaction matrix entry in (4) will be expressed as
and

In (5), if as well as the *j*th face is on the boundary of dielectric bodies, in which case it is assumed that is defined only over the interior tetrahedron and that the exterior tetrahedron is fictitious, or other else. In this paper, the number of the faces of tetrahedrons constituting the boundary of dielectric bodies is counted as a part of the number of triangles.

Using the identical relation (5) can be expressed as where

In the MLFMA, the far-field interactions are done through groups. It is worthy to point out that our grouping scheme is based on the triangles/tetrahedrons but not the edges/faces. In other words, the index of the leaf box to which a given triangle/tetrahedron belongs is determined by comparing its center coordinates to those of the triangle/tetrahedron’s barycenter. Obviously, compared to the edge-/face-based grouping method, the triangle-/tetrahedron-based grouping scheme needs slightly more preprocessing time to collect the unknowns in various boxes at the finest MLFMA level. However, due to the utilization of binary-tree searching algorithm, whose complexity is the order of *O*(*N*log*N*) for *N* basis functions, this increasing time is very limited.

If and are grouped into the and leaf boxes, and the barycenter coordinates of the two boxes are and , respectively, then via the addition theorem, the Green function is transformed as with the translation operator [7]. Thus, (9) and (6) become and respectively, where the superscript “” denotes the complex conjugate and

If we substitute the definition and the divergence equation of the RWG and SWG basis functions into (13), we will obtain with and

In (16), and are the length of the *ς*th edge and the area of the *ς*th face, and and are the area of the *β*th triangle and the volume of the *β*th tetrahedron, respectively. It can be seen that the scalar and vector contain all information needed by the aggregation and disaggregation processes in the MLFMA, which means that they can be regarded as the RPs defined for given triangles/tetrahedrons. Further, the RPs can be expressed as a series of the spherical harmonics as
where is the orthonormalized spherical harmonics [3]. The expansion coefficients are computed by

Instead of the RPs, the expansion coefficients can be computed and stored in the setup of the MLFMA, which can reduce the memory requirement significantly. The remaining steps, such as the computation of SE coefficients as well as aggregation and disaggregation at the finest level using the spherical harmonics series, are done in the same way as presented in [6]. However, it is important to note that if MSFIE is involved into the simulation of the objects that contained closed conducting surfaces, when all incoming waves are collected in a certain group on the finest level, the disaggregation process needs to be operated elaborately. Let

(12) is translated into

If we further introduce and utilize the orthonormality of SE, (20) is simplified to

##### 2.2. Advantages of the Proposed Method

The advantages of the improved SE-MLFMA are stated as follows: (1)The memory requirement of RPs is further decreased. In the traditional SE-MLFMA [6], the basis functions are grouped according to the center of the edges/faces on which the basis functions are defined. As a consequence, the memory requirement of RPs is proportional to the number of basis functions . However, the improved SE-MLFMA groups the basis functions according to the barycenter of triangles/tetrahedrons containing the basis functions, which makes that the memory requirement of RPs is proportional to the total number of triangles and tetrahedrons . Since generally, is much less than , the memory requirement of RPs is drastically reduced. Further, according to (16), we need to store four SE coefficients (one for scalar and the other three for vector ) for each triangle/tetrahedron in the conventional MLFMA. With the symmetry of the RPs, the memory requirement of RPs for the conventional MLFMA in the mixed-potential form is bytes, while for the improved SE-MLFMA, according to (17), the memory requirement is bytes. If we define the memory cost ratio for RPs as , then from (23) and (24), we obtain Usually, is equal to [6, 7]. Compared to the conventional MLFMA, the improved SE-MLFMA shows a significant decrease in the memory requirement of RPs.(2)When MSFIE is involved in the computation, no extra memory requirement of RPs will be needed. In the traditional SE-MLFMA, when MSFIE is coupled with ESFIE to form the CSFIE to analyze the objects containing closed conducting surfaces, double times core memory for the RPs of RWG basis functions will be consumed [6]. However, in the improved SE-MLFMA, due to the symmetry of the triangle’s RPs, the complex conjugate of (16) can be reused as the receiving patterns in the disaggregation process. Therefore, the memory requirement of RPs is not dependent on the type of SIE.(3)The repeated translation between different coordinates is totally eliminated. For the dyadic representation of integral equations in spherical coordinates, when all incoming waves are collected in a certain finest box, due to the discontinuities in the spherical vector components, the Gibbs phenomenon will arise if the integral operation of the spherical harmonics representation is also evaluated in spherical coordinates [6]. Moreover, it is found that can be one degree smaller if the SE is performed for Cartesian vector components instead of spherical ones [6]. From the above, for the dyadic form of integral equations, two spherical components of the RPs tangential to the Ewald surface should be firstly transformed into three Cartesian components for the SE process, which need to be transformed back to the spherical coordinate system for upward aggregation. These repeated translations between different coordinates occupy a portion of the computation time of matrix-vector product (MVP). However, in the improved SE-MLFMA, not only the RPs but also the aggregation, translation, and disaggregation processes are all operated in the Cartesian coordinates. Therefore, the improved SE-MLFMA totally obviates the need for transforming the RPs between different coordinates.

#### 3. Numerical Results

To illustrate the effectiveness of the improved SE-MLFMA in VSIE, two composite closed conductor-dielectric objects are simulated. RWG and SWG basis functions are used to discretize the conducting surfaces and dielectric bodies, respectively. All matrix equations are iteratively solved by a restarted GMRES solver [11] with a simple diagonal preconditioner to reach the 0.001 convergence criterion, while the restarted number is fixed to 100. The finest box size in the MLFMA is 0.2*λ*. All experiments are carried out on a workstation with 2.4 GHz CPU in single precision ().

##### 3.1. Conductor-Dielectric Cylinder

The top half of the cylinder is perfectly conducting, while the bottom half is dielectric (). The object is 5*λ* in diameter and 10*λ* in height, with the conducting segment in the middle. ESFIE is adopted to model the conducting cylinder part. A moderate mesh size is chosen to generate totally 833,574 unknowns with respect to 33,490 triangles and 387,318 tetrahedrons. For the detailed comparison of the conventional MLFMA (conv) in the mixed-potential form and the improved SE-MLFMA (imp), Table 1 shows the influence of and on the memory cost of RPs (Mem), the number of iterations, and the computing time per MVP. The symbol denotes the various selections of with a certain . It is observed that due to the use of SE-MLFMA, a considerable core memory is saved, while the memory requirement of RPs matches with (23) or (24) that is only proportional to the total number of triangles and tetrahedrons but not that of the basis functions. Besides, the improved SE-MLFMA is even more efficient than the conventional MLFMA in the MVP implementation, since the aggregation of the outgoing waves and the disaggregation of the incoming waves at the finest level can also be computed in a faster way by summation of the spherical harmonics instead of the integrations. Moreover, since in the improved SE-MLFMA the whole process is totally executed in the Cartesian coordinates, the repeated transformations between Cartesian and spherical coordinates which are essential in the traditional SE-MLFMA are totally eliminated, resulting in more fast MVP implementation.