International Journal of Antennas and Propagation
Volume 2018 (2018), Article ID 1612498, 6 pages
https://doi.org/10.1155/2018/1612498

Research Article

Efficient and Memory Saving Method Based on Pseudoskeleton Approximation for Analysis of Finite Periodic Structures

Chunbei Luo, Yong Zhang, and Hai Lin

State Key Laboratory of CAD&CG, Zhejiang University, Hangzhou 310058, China

Correspondence should be addressed to Hai Lin; lin@cad.zju.edu.cn

Received 5 April 2018; Revised 11 June 2018; Accepted 24 June 2018; Published 22 July 2018

Academic Editor: Paolo Baccarelli

Copyright © 2018 Chunbei Luo 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.

An efficient and memory saving method based on pseudoskeleton approximation (PSA) is presented for the effective and accurate analysis of finite periodic structures. Different from the macro basis function analysis model, our proposed method uses the formulations derived by the local Rao-Wilton-Glisson basis functions. PSA is not only used to accelerate the matrix-vector product (MVP) inside the single unit but also adopted to decrease the calculation burden of the coupling between the different cells. Moreover, the number of decomposed coupling matrices is minimized due to the displacement invariance of the periodic property. Consequently, even compared with the multilevel fast multipole algorithm (MLFMA), the new method saves much more memory resources and computation time, which is also demonstrated by the numerical examples.

1. Introduction

Periodic structures have recently found wide applications in the electromagnetic engineering such as antenna arrays and metamaterials with negative permittivity and negative permeability. Hence, the accurate and efficient analysis of periodic structures becomes quite essential. If the periodic structure is an infinite array, simple methods can be applied based on Floquet’s theorem [1] or periodic Green’s function [2], where only a single cell of the periodic structure is the domain of interest.

However, all periodic structures have finite size in real-life problems, although the size may be very large. Therefore, the numerical algorithms which accurately consider the mutual coupling between all cells should be used if the accurate results are required. The method of moments (MoM) [3] and its fast algorithms such as fast multipole method (FMM) [4], adaptive cross approximation (ACA) [5], and FFT-based methods [6] are flexible approaches to study the surface problems. However, the efficiency of numerical methods is still limited since the periodic property is not used in the algorithm framework. Recently, a hybrid method combining the accurate MoM and periodic method of moment (PMM) [7] has been proposed which can gain the balance between the two methods. Moreover, some physically based entire-domain basis functions [8] have been developed to reduce the number of unknowns. Further, the FMM and FFT techniques are integrated to accelerate the calculation [9].

Compared to the ACA method, pseudoskeleton approximation (PSA) [10] is also an efficient low-rank-based algebraic fast algorithm which makes it a really competitive alternative. In this paper, we propose an efficient method with low-memory requirement based on PSA to perform the analysis for finite periodic structures effectively and accurately. In consideration of the accuracy of the mutual interactions [8] and the simplicity of implementation, our proposed method uses the formulations derived from the local basis functions instead of macro basis function (MBF) [11, 12]. In this paper, PSA is not only used to accelerate the matrix-vector product (MVP) inside the single unit but also adopted to decrease the calculation burden of the coupling between the different cells. Moreover, the number of decomposed coupling matrices is minimized due to the displacement invariance of the periodic property. With these improvements, an efficient method with low-memory usage of finite periodic objects can be achieved. Several numerical examples are given to show the priority of the proposed method compared to the conventional multilevel fast multipole algorithm (MLFMA) [13] for periodic structures.

2. MoM and PSA Formulation

In this section, the basis principle of MoM and PSA is briefly introduced at first. Then, the choice of arguments in PSA is discussed.

2.1. MoM Equations and Its Fast Algorithms

Consider a time-harmonic electromagnetic wave scattering or radiation problem of an arbitrary perfect electrically conducting (PEC) object. The object is excited by an incident electric field , then the electric field integral equation (EFIE) associated with the surface equivalent currents can be expressed by

where is the tangential unit vector on the surface of the object and and stand for the angular frequency and permeability, respectively. is 3D scalar Green’s function. The linear system of MoM is obtained by discretizing the unknown vector with Rao-Wilton-Glisson (RWG) [14] basis functions and applying Galerkin’s testing method. Let represent the above EFIE matrix system. The MVP process can be accelerated by fast algorithms which can be written as follows:

where is the matrix of near field interactions which are directly computed and stored and stands for couplings between the far-field interactions which will be accelerated together with in the iterative solving process.

2.2. Basic PSA Frameworks

According to the low-rank decomposition, the far-field interaction matrix with rank-deficient property can be approximated by a product of two much smaller submatrices and :

where is the effective rank and satisfies . While in the skeleton approximation (SA) theory [10], there is a nonsingular submatrix in . Denote the rectangular matrices as and which contain the columns and rows of , respectively, then is expressed as

where and have dimensions with and , respectively. In the PSA method, is reevaluated as

where is the pseudoinverse of with dimensions of , then columns and rows are chosen from to get the and . Three aspects need to be noted here: (i) is the inverse of and the computation of inverse is very expensive; (ii) the determination of the value of is a balance between accuracy and efficiency, where is a number large enough so that most important bases will be embedded; (iii) how to choose those working columns and rows of and .

For the (ii) and (iii) problems, we will discuss them in the next subsection. For the first problem, assuming that can be decomposed via singular value decomposition (SVD) as

where and are unitary matrices, is a diagonal matrix with nonnegative real numbers, and is the complex conjugate transpose of . In the actual implementation, the dimension of , , and can be further decreased by a preset threshold [15]. Let , , and represent the reduced submatrix of , , and , respectively. Then, calculation of the pseudoinverse of (i.e., ) is straightforward:

By combining (5) and (7), the original far-field interaction matrix can be decomposed as

2.3. Choice of Arguments in PSA

As mentioned in the previous subsection, the choice of the value of and the specific sampling rows and columns determine the performance of PSA. In the randomized PSA (RPSA) [15], is equal to . Then, the problem of (ii) transforms into how to estimate the rank of the original matrix. In [16], Chai and Jiao give the approximation of rank of the 3D EM problem by , where is the studying wave number. In this paper, the rank is approximated by

where is the diameter of the bounding box corresponding to the octree structure and is a preset positive parameter. The larger the is, the more accurate the matrix decomposition is. In this paper, when is set as 3, the satisfactory accuracy can be guaranteed.

For the problem of (iii), instead of using random numbers in RPSA, we use a strategy analogous to ACA in this paper. Firstly, initialize from the 0th row as the first row index. Then, find the largest entry in this row, and the corresponding column value in which this entry is located is chosen as the next column index. Similarly, find the largest entry in the current column and get the next row index which should be different from all previous row indexes. This process is carried out iteratively until rows and columns are found and stored. Please refer to [17] for more information.

3. Proposed Method for Periodic Structures

We consider a case of arbitrarily shaped PEC patch, for example, refer in Figure 1. The total impedance matrix contains two parts: self-coupling blocks and mutual coupling blocks. Therefore, both the blocks are analyzed and decomposed by PSA to gain better efficiency and low-memory usage.

Figure 1: Periodic structure with a finite size of .

3.1. PSA for Decomposition of Self-Coupling Matrix

In the previous mentioned periodic algorithms, the self-coupling block elements are calculated by a direct MoM method which is not efficient for a large array unit. Therefore, PSA is used to decompose the self-coupling matrix. As the same to all fast algorithms of MoM, the low-rank decompositions based on PSA are performed on the far-field groups while the near-field elements are directly computed and saved. Moreover, some operations are taken out to further improve the efficiency of PSA. Firstly, the pseudoinverse of the (i.e., ) will not be directly calculated since the complexity of SVD is very high. Instead, the LU decomposition will be performed when the dimension grows up. Different from the previously proposed PSA, , , and the LU decomposition of will be stored. Moreover, and in (8) are further decomposed by ACA technique. Finally, the original far-field coupling matrix can be decomposed into six submatrices which is showed by the following formula:

where (or ) is the ACA decomposition of (or ) and is the LU decomposition of . Note that the inverse of LU matrix of is not directly computed. Actually, the LU back substitute is performed after the MVP based on the matrix in each MVP process.

3.2. PSA for Decomposition of Mutual Coupling Matrices

In the traditional fast algorithms, the matrix decompositions are only performed on the far-field groups especially for the physically related methods such as MLFMA. However, the mutual interactions between two different cells (off-diagonal blocks) are all computed by PSA in our implementation. This treatment has been also used in [18] and demonstrated to be more efficient but losing little accuracy. Moreover, the displacement invariance property is explored in our implementation. Since the mutual coupling remains the same when the distance between two groups’ centers is not changed, the corresponding matrix will be calculated and stored only once under the index of relative distance. In the MVP process, the stored coupling matrices may be used several times when the relative distances are the same. By using this scheme, the number of mutual coupling matrices can be decreased dramatically. For example, the original number of off-diagonal matrix blocks of Figure 1 is while the reduced number will be 62 if the displacement invariance is used.

3.3. Preconditioner Considerations

The preconditioning stage should be also considered carefully while dealing with complex structures. In this paper, the preconditioning matrix for impedance matrix is built based on the self-coupling matrix .

Moreover, the inner-outer iteration scheme can be also applied when facing ill-conditioned matrix systems such as EFIE for radiation analysis. In our implementation, the inner iteration is the solution of the self-coupling matrices. Certainly, the inner solution area can be also extended to improve the convergence rate of outer iteration.

4. Numerical Results

In this section, several numerical experiments are presented to demonstrate the efficiency and low-memory requirement of the proposed method. For all the simulations, the mesh sizes are no less than 0.1λ (λ represents the wave length). The biconjugate gradient stabilized method (BiCGSTAB) is adopted as the iterative solver for the matrix equations, and the threshold of the iterative stopping criterion for residuum is set to 0.001. The sparse approximate inverse preconditioner is used in all simulations. All the computations were carried out on a workstation with four Xeon E5-4620 CPU and 256 GB of RAM with OpenMP technique, and the digits were stored in double precision.

4.1. Scattering Simulation

In order to verify the accuracy and efficiency of the proposed method, we first consider a periodic structure consisting of element cells (also shown in Figure 1), where the unit cell is a 0.5 m PEC sphere. The working frequency is set to be 1 GHz which leads to 20286 unknowns (degrees of freedom) in each patch. Hence, the total number of unknowns is . The gap between two unit cells is set to be 0.5 m. The radar cross sections (RCS) computed by MLFMA, the proposed method (periodic PSA), and the commercial software (FEKO) are illustrated in Figure 2(a). It could be clearly seen that the numerical results from the proposed method have excellent agreement with both the conventional MLFMA and FEKO. While as in Table 1, both the computation time and peak memory usage by periodic PSA are less than the MLFMA method. What is more, another case is considered when the gap between two unit cells is decreased to be 0.1 m. The RCS results and the use of computer resources are shown in Figure 2(b) and Table 1, respectively. Although the rank of mutual coupling matrix gets bigger due to the tight coupling interactions, the proposed method still shows the good efficiency and low-memory usage with losing little accuracy.

(a)

(b)

Figure 2: Radar cross sections (V-V polarization) of the sphere array with (a) 0.5 m and (b) 0.1 m gap computed by FEKO, MLFMA, and periodic PSA. The incident plane wave angles satisfy ° and °. The scattering angles work with varies from 0° to 360° and °.

Table 1: Computation time and memory usage of MLFMA and periodic PSA for scattering of elements.

Gap
Methods
Computation time
Memory usage

0.5 m
MLFMA
210 s
8726 MB

Periodic PSA
109 s
2145 MB

0.1 m
MLFMA
222 s
8801 MB

Periodic PSA
131 s
2305 MB

Furthermore, the scattering of a large-scale periodic structure consisting of element cells is considered. The working frequency is set to be 1 GHz. Different radii of the sphere are considered which lead to about 0.81, 3.38, and 5.46 million numbers of unknowns corresponding to 0.1 m, 0.2 m, and 0.25 m models, respectively. The distance between two centers of neighboring sphere is set to be a fixed value of 1.0 m. Table 2 shows the information of models and the numerical performance between MLFMA and the periodic PSA. It can be seen that for a large-scale problem, the proposed PSA method still needs less computation time and much less memory requirement. It should be noted that in the 3rd example in Table 2 when the sphere radius is set to 0.25 m, the finest box in MLFMA has to be set to 0.20λ to avoid low-frequency breakdown. Figure 3 shows that the RCS results of 0.25 m sphere models computed from the periodic PSA still agree well with FEKO and MLFMA.

Table 2: Computation time and memory usage of MLFMA and periodic PSA for scattering of elements.

Radius
Gap
Methods
Computation time
Memory usage

0.1 m
0.8 m
MLFMA
369 s
18229 MB

Periodic PSA
113 s
707 MB

0.2 m
0.6 m
MLFMA
1551 s
71455 MB

Periodic PSA
895 s
3740 MB

0.25 m
0.5 m
MLFMA
2341 s
96906 MB

Periodic PSA
2067 s
6136 MB

Figure 3: Radar cross sections (V-V polarization) of the sphere array with 0.5 m gap computed by FEKO, MLFMA, and periodic PSA. The incident plane wave angles satisfy ° and °. The scattering angles work with varies from 0° to 360° and °.

4.2. Radiation Simulation

Lastly, we consider a radiation problem, which contains antenna cells with delta-gap excitations. Each cell has 1145 unknowns, and the total number is . The gap between two cells in this model is 0.005 m which is smaller than 0.1λ. Since EFIE is the only choice for radiation problems, the inner-outer iteration scheme is adopted to accelerate the calculation process. Figure 4 shows the current density distribution of the antenna array. The far-field gain patterns computed by the proposed method and MLFMA are shown in Figure 5.

Figure 4: Current density distribution of patch array.

Figure 5: Radiation patterns of the antenna array with 0.005 m gap computed by MLFMA and periodic PSA. The far-field angles work with varies from 0° to 360° and °.

Obviously, the two methods give nearly the same results. However, as in Table 3, the computation time has been reduced from 419 seconds to 172 seconds in solving matrix equations for the radiation problem. Moreover, the memory usage for MLFMA is 8841 MB while the value is only 358 MB for the proposed method. Again, the proposed method performs much better than MLFMA. The first reason is that only 360 off-diagonal blocks are need to be calculated and stored in our method while the original value is . Secondly, the mesh size is about 0.01λ for the antenna cell (dense mesh). As it is known, the lowest-level box size is limited to be no less than 0.20λ since MLFMA suffers from the low-frequency breakdown, which results in a heavy near-field matrix consumption.

Table 3: Computation time and memory usage of MLFMA and periodic PSA for radiation of antenna array.

Methods
Computation time
Memory usage

MLFMA
419 s
8841 MB

Periodic PSA
172 s
358 MB

5. Conclusion

In this paper, an efficient and memory saving method based on pseudoskeleton approximation (PSA) for the effective and accurate analysis of finite periodic structures is presented. The appropriate choice of sampling rows and columns as well as sampling number in PSA is discussed and confirmed. Based on the algebraic fast algorithm PSA, we have carefully handled the self-coupling and mutual-coupling blocks of the original impedance matrix. Moreover, the number of decomposed coupling matrices is minimized by the employment of the displacement invariance of the periodic property. The simulation results of large-scale scattering and radiation examples show that the proposed periodic PSA method needs less computation time and much less memory compared to conventional MLFMA. Hence, numerical results demonstrate the efficiency and superiority of the proposed method.

Data Availability

No additional data are available.

Conflicts of Interest

The authors declare that there is no conflict of interest regarding the publication of this paper.

Acknowledgments

This work was supported in part by the National Natural Science Foundation of China under Grant 61171035.

References

F. Xu, Y. Zhang, W. Hong, K. Wu, and T. J. Cui, “Finite-difference frequency-domain algorithm for modeling guided-wave properties of substrate integrated waveguide,” IEEE Transactions on Microwave Theory and Techniques, vol. 51, no. 11, pp. 2221–2227, 2003.

K. Yasumoto and K. Yoshitomi, “Efficient calculation of lattice sums for free-space periodic Green’s function,” IEEE Transactions on Antennas and Propagation, vol. 47, no. 6, pp. 1050–1055, 1999.

R. F. Harrington, Field Computation by Moment Methods, Wiley-IEEE Press, 1993.

C.-C. Lu and W. C. Chew, “A multilevel algorithm for solving a boundary integral equation of wave scattering,” Microwave and Optical Technology Letters, vol. 7, no. 10, pp. 466–470, 1994.

M. Bebendorf, “Approximation of boundary element matrices,” Numerische Mathematik, vol. 86, no. 4, pp. 565–589, 2000.

E. Bleszynski, M. Bleszynski, and T. Jaroszewicz, “Aim: adaptive integral method for solving large-scale electromagnetic scattering and radiation problems,” Radio Science, vol. 31, no. 5, pp. 1225–1251, 1996.

J. Su, X. Xu, and B. Hu, “Hybrid PMM-MoM method for the analysis of finite periodic structures,” Journal of Electromagnetic Waves and Applications, vol. 25, no. 2-3, pp. 267–282, 2011.

W. B. Lu, T. J. Cui, Z. G. Qian, X. X. Yin, and W. Hong, “Accurate analysis of large-scale periodic structures using an efficient sub-entire-domain basis function method,” IEEE Transactions on Antennas and Propagation, vol. 52, no. 11, pp. 3078–3085, 2004.

W. B. Lu, T. J. Cui, and H. Zhao, “Acceleration of fast multipole method for large-scale periodic structures with finite sizes using sub-entire-domain basis functions,” IEEE Transactions on Antennas and Propagation, vol. 55, no. 2, pp. 414–421, 2007.

S. A. Goreinov, N. L. Zamarashkin, and E. E. Tyrtyshnikov, “Pseudo-skeleton approximations by matrices of maximal volume,” Mathematical Notes, vol. 62, no. 4, pp. 515–519, 1997.

E. Suter and J. R. Mosig, “A subdomain multilevel approach for the efficient MoM analysis of large planar antennas,” Microwave and Optical Technology Letters, vol. 26, no. 4, pp. 270–277, 2000.

D. Gonzalez-Ovejero, F. Mesa, and C. Craeye, “Accelerated macro basis functions analysis of finite printed antenna arrays through 2D and 3D multipole expansions,” IEEE Transactions on Antennas and Propagation, vol. 61, no. 2, pp. 707–717, 2013.

J. Song, C.-C. Lu, and W. C. Chew, “Multilevel fast multipole algorithm for electromagnetic scattering by large complex objects,” IEEE Transactions on Antennas and Propagation, vol. 45, no. 10, pp. 1488–1493, 1997.

S. Rao, D. Wilton, and A. Glisson, “Electromagnetic scattering by surfaces of arbitrary shape,” IEEE Transactions on Antennas and Propagation, vol. 30, no. 3, pp. 409–418, 1982.

X. Zhu and W. Lin, “Randomised pseudo-skeleton approximation and its application in electromagnetics,” Electronics Letters, vol. 47, no. 10, pp. 590–592, 2011.

W. Chai and D. Jiao, “Theoretical study on the rank of integral operators for broadband electromagnetic modeling from static to electrodynamic frequencies,” IEEE Transactions on Components, Packaging and Manufacturing Technology, vol. 3, no. 12, pp. 2113–2126, 2013.

Y. Zhang and H. Lin, “Localized pseudo-skeleton approximation method for electromagnetic analysis on electrically large objects,” Progress In Electromagnetics Research Letters, vol. 57, pp. 103–109, 2015.

A. Heldring, J. M. Tamayo, C. Simon, E. Ubeda, and J. M. Rius, “Sparsified adaptive cross approximation algorithm for accelerated method of moments computations,” IEEE Transactions on Antennas and Propagation, vol. 61, no. 1, pp. 240–246, 2013.