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

Research Article

Versatile Solver of Nonconformal Volume Integral Equation Based on SWG Basis Function

Chunbei Luo, Mingjie Pang, 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 14 June 2018; Accepted 24 September 2018; Published 9 December 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.

In this paper, a versatile solver of a nonconformal volume integral equation based on the Schaubert-Wilton-Glisson (SWG) basis function is presented. Instead of using a piecewise constant function, the robust conventional SWG basis function is chosen and used directly for discontinuous boundaries. A new map method technique is proposed for constructing SWG pairs, which reduces the complexity from to compared with a brute-force method. The integral equation is solved by the method of moments (MoM) and further accelerated by the multilevel fast multipole algorithm (MLFMA). What’s more, the hybrid scheme of MLFMA and adaptive cross approximation (ACA) is developed to resolve the low-frequency (LF) breakdown when dealing with over-dense mesh objects. Numerical results show that when in analysis of radiation or scattering problems from inhomogeneous dielectric objects or in LF conditions, the proposed solver shows high efficiency without loss of accuracy, which demonstrates the versatile performance of the proposed method.

1. Introduction

Electromagnetic (EM) scattering and radiation analysis of dielectric materials are attracting increasing attention for their potential for vast applications, such as in designing a microstrip dielectric antenna and scattering reduction of the coated layer in a stealth plane, dielectric radome, biological media, and plasma sheath. Generally, these problems can be solved by an integral equation using method of moments (MoM) [1] because it has lesser degrees of freedom than differential equation methods. Compared with a surface integral equation (SIE) such as a PMCHWT equation [2], a volume integral equation (VIE) is more flexible, robust, and accurate [3]; therefore, VIE is usually preferred or even the only option for complex dielectric anisotropic objects.

In [4], the discontinuous Galerkin (DG) VIE using the Schaubert-Wilton-Glisson (SWG) basis function is proposed and a hybrid discretization scheme that uses a mix of nonconformal and conformal meshes is adopted. In [4, 5], Zhang et al. found that the explicit enforcement of the continuity condition at the interface between two neighboring elements is not required in the DGVIE. This is because of the inherent quality of the Fredholm integral equation of the second kind, where the boundary conditions are naturally imbedded. Therefore, the conventional SWG basis function is a better choice to apply in nonconformal VIE than the piecewise constant function [6], since the SWG basis is more robust and accurate.

In this paper, we exploit a versatile solver of a nonconformal volume integral equation based on the SWG basis function. SWG pairs are used in the two neighboring elements sharing the common face, while half-SWG basis functions are adopted in the discontinuous boundary elements. Hence, before the simulation, it is important to find all the neighboring two elements (i.e., two neighboring tetrahedrons which share the common face) efficiently. In this paper a new map method technique is proposed to find SWG pairs faster than the brute-force method, which reduces the complexity of the computation time from to . To improve the efficiency of the proposed solver, the multilevel fast multipole algorithm (MLFMA) [7] is employed to accelerate the matrix vector multiplication (MVM) in the iterative solution process. In the analysis of dielectric objects containing different media, different scales of discrete elements are used for nonconformal VIE and this usually results in a dense mesh. However, when dealing with over-meshed objects, in other words, when the minimum box size is too small compared to the wavelength, the MLFMA suffers from a low-frequency (LF) breakdown problem [8, 9] limited by the addition theorem to calculate Green’s function. Another fast low-rank decomposition solver, namely adaptive cross approximation (ACA) [10], is less efficient compared to MLFMA, since it is a purely algebraic solver and is independent of Green’s function. Therefore, a hybrid scheme of MLFMA and ACA is developed to resolve the LF breakdown in this paper. For the different levels of an octree structure, different algorithms are used. This is because ACA does not exist in an LF breakdown and it is applied to accelerate MVM operations in the lowest levels, while MLFMA is more efficient, has lower complexity, and can be applied in higher levels. In summary, the contributions of this work can be listed as follows:

(1)The nonconformal VIE expanded by the SWG basis function in discontinued boundaries is explained in detail, and the solver is accelerated by MLFMA(2)A new technique for constructing SWG pairs is introduced which reduces the complexity of the computation time from to . To the best of our knowledge, this part has never been published before(3)A hybrid scheme of MLFMA and ACA is developed which resolves the LF breakdown problem

2. Formulation

First, the matrix system discretized by MoM based on the SWG basis function will be introduced, and the treatment of discontinued boundaries and the corresponding impedance matrix will be described in detail. Then, the brute-force method and map method for constructing SWG pairs are introduced and the complexities are analyzed. Finally, the hybrid MLFMA-ACA scheme is presented to settle the LF breakdown problem.

2.1. Volume Electric Field Integral Equation

Consider an inhomogeneous dielectric object with a permittivity of and a permeability of residing in free space with EM parameters and . The object is excited by an incident electric field , and the volume electric field integral equation (VEFIE) associated with the vector potential and scalar potential can be expressed by [11]

where is the angular frequency and is the wavenumber of free space. is the volume equivalent currents and is related to the electric flux density by

Discretizing the unknown vector by

where is the SWG basis function [12] associated with the th face, is the corresponding unknown coefficient, and is the total number of unknowns. Figure 1 shows the definition of the SWG pair and the half-SWG. Applying Galerkin’s test, then (1) can be rewritten in the matrix equations form as follows:

where , , and are the matrices corresponding to the first, second, and third terms of the left side of (1), the unknown vector relates with the unknown coefficient , and and are the impedance matrix and excitation vector, respectively.

(a)

(b)

Figure 1: (a) SWG pair and (b) half-SWG associated with the th element. and are the tetrahedron pairs and and are the position vectors, respectively. For the boundary element, the basis functions are all defined with half-SWG in in this paper.

Since (1) is the Fredholm integral equation of the second kind, it is naturally imbedded with the third type of boundary condition [4]. Therefore, the enforcement of the continuity condition is not required in this paper. Then, in the nonconformal elements of an object, the half-SWG is defined. Suppose that there are totally SWG pairs defined in tetrahedron pairs and single (half) SWGs defined in boundary (including the nonconformal boundaries) tetrahedrons; thus, the impedance matrix in (6) can be expressed as

Note that . The entries of the matrix blocks in (7) are given by

where . equals to 1 only if ; otherwise, it is equal to 0. is the -th entry in of (6) corresponding to the th basis function and th test function, and and are the tetrahedron integral domains. Note that only appears in the and , and half-SWGs are all defined in . The calculation of the integral formulation of has been developed in [12], and for the singularity treatment one can refer to [13].

2.2. Constructing SWG Pairs

Since SWG is defined in pairs for the conformal element which shares the same face, it is important and necessary to find the tetrahedron pairs efficiently before the simulation process. To the best of our knowledge, this part has never been published yet. Consider a CAD model which is discretized with numbers of independent tetrahedrons totally. In a tetrahedron pair (i.e., in Figure 1(a)), there are 3 vertices that makes the common face and the left two are SWG free vertices. By traversing the total tetrahedrons, the goal of constructing SWG pairs is to sort all the conformal elements into one group (hereinafter denoted as SWG-pair-group) and the remaining nonconformal elements (including physical boundaries of the model and the nonconformal boundaries due to the different media) into the other.

The brute-force method first searches for all the tetrahedrons as the ; then, it searches the rest of the tetrahedrons as . Note that a tetrahedron has 4 faces. Hence, there are twice the loops of searching that needs to be done with times. Then, the complexity of the computation time is , where and is a constant. The pseudocode for the brute-force method is shown in Figure 2(a).

(a)

(b)

Figure 2: Pseudocode for constructing SWG pairs by (a) the brute-force method and (b) map method. The special function of the map is marked with blue. All pairs of conformal SWG functions are stored in the SWG-pair-group.

A new map method technique is proposed for constructing SWG pairs which is much more efficient than the brute-force method. Map is one of the associative containers that refer to a group of class in the standard template library (STL) of the C++ programming language [14, 15]. It allows mapping from one data item (a key) to another (a value) by constructing a self-balancing binary search tree, typically a red-black tree, which guarantees search complexity in [16]. The pseudocode for map is shown in Figure 2(b). The implementation of map is described as follows: (1) a map container named T-map is created; (2) a face (denoted as ) of a tetrahedron (denoted as ) is set as a key; (3) find if the key () is already existing in the T-map; (4) if an existing key is found, then a SWG pair is constructed and this key must be deleted since a face is shared by two tetrahedrons at most; otherwise, save this key and the value (). It is obvious that two loops for searching all the faces need times and searching a member in T-map in step (3) needs times at most. Hence, the total computation time is ; then, the complexity of the computation time by map is .

Numerical results in Subsection 3.1 verified the above analysis. It should be noted that the sort of nonconformal elements is not shown in the pseudocode. In fact, after searching for all the SWG pairs, the remaining members in T-map are all of the half-SWG functions.

2.3. Hybrid MLFMA-ACA Scheme

Due to the ( is the number of unknowns) complexity of MoM in storage, impedance matrix fill-in, and MVM, the matrix system is therefore solved with MLFMA to reduce the complexity to . For the VEFIE matrix equations in (6), the acceleration of MVM by MLFMA can be expressed by [17]

where is the matrix of near field interactions and , , and are the disaggregation, translation, and aggregation matrices, respectively. However, MLFMA breaks down in LF conditions. To tackle this problem, we propose a hybrid MLFMA-ACA method which takes advantage of ACA with the absence of the LF breakdown problem and MLFMA with high efficiency. In ACA, the far field interaction matrix is approximated by with the product of two low-rank matrices and as [18]

where is the matrix dimension and is the effective rank of the matrix and satisfies . Combining (9) and (10), the hybrid MLFMA-ACA scheme for the MVM calculation can then be expressed as

where and denotes the MLFMA and ACA work matrices, respectively. represent the different levels of an octree structure. Since the first two levels are all near field interaction groups, the MVM acceleration starts at . It should be noted that the larger l (denotes the lower level) corresponds to the smaller box size. By (11), the far field interactions of MVM are divided into ACA acceleration from to in the lower levels without LF breakdown and MLFMA acceleration from 2 to in higher levels to maintain high efficiency. If , it means that only MLFMA is used. Generally is determined according to the box size of the th level which is less than a preset value, e.g., ( is the wavelength in free space), when MLFMA gives poor accuracy.

3. Numerical Results

In this section, several numerical results are presented to demonstrate the versatility of the proposed nonconformal VIE solver. The example in Subsection 3.1 was carried out on a desktop with Intel® i5-4670 CPU @3.40 GHz and 32 GB RAM, others in Subsections 3.2 and 3.3 were done on a workstation with 32 cores of Intel® E5-4620 CPU @2.60 GHz and 256 GB RAM.

3.1. Performance of Constructing SWG Pairs

First, the efficiency of constructing SWG pairs by brute-force and map method in Subsection 2.2 is investigated. A sphere is modeled with different mesh sizes generating different numbers of tetrahedrons which vary from thousands to millions. Computation time of constructing SWG pairs by the two methods and their complexity are depicted in Figure 3. It can be seen that the computation complexity of the brute-force method is , while it is of the map method. Note that for over 1.4 million tetrahedrons, the computation times by map is still within seconds while by brute force it is unable to finish being limited by the computer source. Using the map method ensures the superior efficiency of the proposed solver in constructing SWG basis functions.

Figure 3: Computation time of constructing SWG pairs versus the number of tetrahedrons.

3.2. Examples of Scattering and Radiation

To demonstrate the accuracy and efficiency of the proposed nonconformal VIE solver, the bistatic RCS of a coated dielectric sphere in Figure 4 is studied firstly. The inner and outer radii of the two-layer sphere are 0.3 m and 0.5 m, respectively, with the relative permittivity of – and –. The layered sphere is illuminated at by an incident plane wave at 500 MHz. For consistency, the standard mesh size is set as , which results in 40,311 tetrahedrons for conformal VIE and 26,990 (8751 for the inner part and 18,239 for the outer) for nonconformal VIE. Figure 5 shows the bistatic RCS results from the two methods. It can be seen that both methods agree well with the Mie series. A comparison of the numerical performance between two methods is given in Table 1. It is obvious that nonconformal VIE has an advantage over conformal VIE in both memory usage and total computation time because of the reduction of the number of unknowns.

(a)

(b)

Figure 4: Geometrical model of a coated dielectric sphere. (a) Geometrical parameter definition and (b) cross section of the nonconformal meshes.

Figure 5: Bistatic RCS in VV polarization from a coated dielectric sphere with the relative permittivities of – and –, an inner radius of , and an outer radius of .

Table 1: Comparison of numerical performance between conformal VIE and nonconformal VIE.

Methods
Number of unknowns
Memory usage
Total time

Conformal VIE
82,325
11,523 MB
341.25 s

Nonconformal VIE
56,125
6122 MB
199.79 s

To further demonstrate the accuracy and flexibility of the proposed solver, a radiation pattern of an eight-dipole array in a three-layer hemispherical dielectric radome is analyzed. The geometric model is shown in Figure 6 which is the same as the Figure 4(c) in [11]. The base radius of the radome is 50 mm, the thicknesses and the relative permittivities of both the innermost and the outermost layers are 3 mm and 2.0–, and those of the middle layer are 4 mm and 1.5–, respectively. The array is symmetrically located in the sphere center along the -axis, and the corresponding 8 dipoles placed with equal intervals () are fed with an equal excitation amplitude, equal orientation (, ), and progressive phase. The working frequency is set at 3.0 GHz. Figure 7 shows the normalized far-field radiation pattern with or without a radome. It is observed that the proposed method agrees well with the result calculated by the commercial software FEKO, which is also consistent with the result of [11] (in Figure 8). This indicates the accuracy and reliability of the nonconformal VIE solver in this paper.

Figure 6: Geometrical model of a dipole array in a three-layer hemispherical dielectric radome.

Figure 7: Normalized far-field radiation pattern in the XOZ plane of a small 8-dipole array in a hemispherical radome.

Figure 8: Bistatic RCS in VV polarization from a coated dielectric sphere with the relative permittivities of and , an inner radius of , and an outer radius of .

3.3. Example of LF Breakdown Problem

Finally, a numerical example of an LF breakdown is analyzed. Consider a coated sphere in Figure 4(a) with an inner radius of , an outer radius of , and a relative permittivity of , . The two-layer sphere is illuminated at by an incident plane wave at 300 MHz. The mesh sizes of the inner and outer parts are and , which result in 31,481 and 278,902 tetrahedrons, respectively. The total number of unknowns is 630,564. Note that this is an over-meshed object and MLFMA faces the LF breakdown problem. The minimal box size in the lowest level is set as , thereby a 6-level octree () is constructed. For the hybrid MLFMA-ACA method, the lowest levels where is from 4 to 5 are accelerated by ACA, while the higher levels where is from 2 to 3 are accelerated by MLFMA. Figure 8 shows the bistatic RCS results from MLFMA, ACA, and the hybrid MLFMA-ACA methods. It can be observed that both ACA and hybrid MLFMA-ACA agree well with the Mie analytical solution, while MLFMA breaks down as expected. In this case, MLFMA fails to converge in 500 iteration steps and gets a result of poor accuracy. MLFMA-ACA and ACA each have a total computation time of 5.67 h and 7.55 h, and a memory usage of 113.2 GB and 162.8 GB, respectively. It can be seen that the hybrid method is immune to LF breakdown and achieves higher efficiency and less memory usage compared with the MLFMA and ACA methods.

4. Conclusion

In this paper, we present a versatile solver of nonconformal VIE based on the SWG basis function to analyze scattering and radiation problems. The discontinuous boundaries are dealt with half-SWG, and the impedance matrix of VEFIE are deduced in detail. A new technique based on the map method is first introduced to construct SWG pairs which reduces the complexity of the computation time from to . Furthermore, the VIE solver is accelerated by MLFMA and the LF breakdown problem is resolved by the hybrid MLFMA-ACA scheme. This hybrid method is quite simple and easy to implement which utilizes ACA in the lowest levels and MLFMA in the higher levels to accelerate the MVM process. Numerical results verify the complexity analysis of the map method, the high efficiency and lesser memory usage of nonconformal VIE, and the absence of LF breakdown using the MLFMA-ACA hybrid scheme, which demonstrates the versatility of the proposed solver.

Data Availability

No data were used to support this study.

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

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

A. J. Poggio and E. K. Miller, “Integral equation solutions of three-dimensional scattering problems,” in Computer Techniques for Electromagnetics: International Series of Monographs in Electrical Engineering, MB Associates, 1970.

A. F. Peterson, S. L. Ray, R. Mittra, Institute of Electrical, and Electronics Engineers, Computational Methods for Electromagnetics, Vol. 2, IEEE press, New York, NY, USA, 1998.

L.-M. Zhang and X.-Q. Sheng, “Discontinuous galerkin volume integral equation solution of scattering from inhomogeneous dielectric objects by using the SWG basis function,” IEEE Transactions on Antennas and Propagation, vol. 65, no. 3, pp. 1500–1504, 2017.

L.-M. Zhang and X.-Q. Sheng, “A discontinuous galerkin volume integral equation method for scattering from inhomogeneous objects,” IEEE Transactions on Antennas and Propagation, vol. 63, no. 12, pp. 5661–5667, 2015.

N. A. Ozdemir and J.-F. Lee, “A nonconformal volume integral equation for electromagnetic scattering from penetrable objects,” IEEE Transactions on Magnetics, vol. 43, no. 4, pp. 1369–1372, 2007.

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.

E. Darve and P. Have, “A fast multipole method for Maxwell equations stable at all frequencies,” Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences, vol. 362, no. 1816, pp. 603–628, 2004.

H. Wallen and J. Sarvas, “Translation procedures for broadband MLFMA,” Progress In Electromagnetics Research, vol. 55, pp. 47–78, 2005.

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

Q.-M. Cai, Y.-W. Zhao, Y.-T. Zheng, M.-M. Jia, Z. Zhao, and Z.-P. Nie, “Volume integral equation with higher order hierarchical basis functions for analysis of dielectric electromagnetic scattering,” IEEE Transactions on Antennas and Propagation, vol. 63, no. 11, pp. 4964–4975, 2015.

D. Schaubert, D. Wilton, and A. Glisson, “A tetrahedral modeling method for electromagnetic scattering by arbitrarily shaped inhomogeneous dielectric bodies,” IEEE Transactions on Antennas and Propagation, vol. 32, no. 1, pp. 77–85, 1984.

D. Wilton, S. Rao, A. Glisson, D. Schaubert, O. Al-Bundak, and C. Butler, “Potential integrals for uniform and linear source distributions on polygonal and polyhedral domains,” IEEE Transactions on Antennas and Propagation, vol. 32, no. 3, pp. 276–281, 1984.

E. P. Becker, “Working draft, standard for programming language c++,” Tech. Rep., Tech. Rep. N2009 = 06-0079, ISO/IEC JTC 1, Information Technology, Subcommittee SC 22, Programming Language C++, 2006.

N. M. Josuttis, The C++ Standard Library: A Tutorial and Reference, Addison-Wesley, 2012.

T. H. Cormen, C. E. Leiserson, R. L. Rivest, and C. Stein, Introduction to Algorithms Second Edition, The MIT Press, 2001.

W. C. Chew, E. Michielssen, J. Song, and J.-M. Jin, Fast and Efficient Algorithms in Computational Electromagnetics, Artech House, Inc., 2001.

K. Zhao, M. N. Vouvakis, and J.-F. Lee, “The adaptive cross approximation algorithm for accelerated method of moments computations of emc problems,” IEEE Transactions on Electromagnetic Compatibility, vol. 47, no. 4, pp. 763–773, 2005.