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

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 $O\left(N^2\right)$ to $O\left(N\log N\right)$ 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.

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.

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 $O\left(N^2\right)$ to $O\left(N\log N\right)$. 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 $O\left(N^2\right)$ to $O\left(N\log N\right)$. 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 ${\epsilon }_r\left(\mathbf{r}\right)$ and a permeability of ${\mu }_0$ residing in free space with EM parameters ${\epsilon }_0$ and ${\mu }_0$. The object is excited by an incident electric field $E^{inc}\left(\mathbf{r}\right)$, and the volume electric field integral equation (VEFIE) associated with the vector potential $A\left(\mathbf{r}\right)$ and scalar potential $\mathrm{\Phi }\left(\mathbf{r}\right)$ can be expressed by [11] $$\begin{array}{}\text{(1)}& \frac{D\left(\mathbf{r}\right)}{\epsilon \left(\mathbf{r}\right)}+j\omega A\left(\mathbf{r}\right)+\nabla \mathrm{\Phi }\left(\mathbf{r}\right)=E^{inc}\left(\mathbf{r}\right),\text{(2)}& A\left(\mathbf{r}\right)={\mu }_0{\int }_{V}J\left({\mathbf{r}}^{\prime }\right)\frac{e^{-{jk}_0\left|\mathbf{r}-{\mathbf{r}}^{\prime }\right|}}{4\pi \left|\mathbf{r}-{\mathbf{r}}^{\prime }\right|}d{v}^{\prime },\text{(3)}& \mathrm{\Phi }\left(\mathbf{r}\right)=-\frac{1}{j{\omega \epsilon }_0}{\int }_V{\nabla }^{\prime }\cdot J\left({\mathbf{r}}^{\prime }\right)\frac{e^{-{jk}_0\left|\mathbf{r}-{\mathbf{r}}^{\prime }\right|}}{4\pi \left|\mathbf{r}-{\mathbf{r}}^{\prime }\right|}d{v}^{\prime },\end{array}$$where $\omega$ is the angular frequency and $k_0$ is the wavenumber of free space. $J\left(\mathbf{r}\right)$ is the volume equivalent currents and is related to the electric flux density $D\left(\mathbf{r}\right)$ by $$\begin{array}{c}\text{(4)}& J\left(\mathbf{r}\right)=j\omega \kappa \left(\mathbf{r}\right)D\left(\mathbf{r}\right),\kappa \left(\mathbf{r}\right)=\frac{\epsilon \left(\mathbf{r}\right)-{\epsilon }_0}{\epsilon \left(\mathbf{r}\right)}.\end{array}$$

Discretizing the unknown vector $D\left(\mathbf{r}\right)$ by $$\begin{array}{}\text{(5)}& D\left(\mathbf{r}\right)=\sum _{n=1}^N{I}_n{f}_n\left(\mathbf{r}\right),\end{array}$$where $f_n\left(\mathbf{r}\right)$ is the SWG basis function [12] associated with the $n$th face, $I_n$ is the corresponding unknown coefficient, and $N$ 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: $$\begin{array}{}\text{(6)}& \left\{\left[D\right]+\left[A\right]+\left[\mathrm{\Phi }\right]\right\}\cdot \left[I\right]=\left[Z\right]\cdot \left[I\right]=\left[V\right],\end{array}$$where $\left[D\right]$, $\left[A\right]$, and $\left[\mathrm{\Phi }\right]$ are the matrices corresponding to the first, second, and third terms of the left side of (1), the unknown vector $\left[I\right]$ relates with the unknown coefficient $I_n$, and $\left[Z\right]$ and $\left[V\right]$ are the impedance matrix and excitation vector, respectively.

[figures omitted; refer to PDF]

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 $N_P$ SWG pairs defined in tetrahedron pairs and $N_S$ single (half) SWGs defined in boundary (including the nonconformal boundaries) tetrahedrons; thus, the impedance matrix in (6) can be expressed as $$\begin{array}{}\text{(7)}& \left[Z\right]=\left[\begin{array}{cc}{\left[Z_{PP}\right]}_{N_P\times N_P}& {\left[Z_{PS}\right]}_{N_P\times N_S}\\ {\left[Z_{SP}\right]}_{N_S\times N_P}& {\left[Z_{SS}\right]}_{N_S\times N_S}\end{array}\right].\end{array}$$

Note that $Np+Ns=N$. The entries of the matrix blocks in (7) are given by $$\begin{array}{c}\text{(8)}& {\left[Z_{PP}\right]}_{mn}={\left[D\right]}_{mn}^{T_m^{+},T_n^{+}}\cdot \delta \left(T_m^{+},T_n^{+}\right)+{\left[D\right]}_{mn}^{T_m^{+},T_n^{-}}\cdot \delta \left(T_m^{+},T_n^{-}\right)+{\left[D\right]}_{mn}^{T_m^{-},T_n^{+}}\cdot \delta \left(T_m^{-},T_n^{+}\right)+{\left[D\right]}_{mn}^{T_m^{-},T_n^{-}}\cdot \delta \left(T_m^{-},T_n^{-}\right)+{\left[A\right]}_{mn}^{T_m^{+},T_n^{+}}+{\left[A\right]}_{mn}^{T_m^{+},T_n^{-}}+{\left[A\right]}_{mn}^{T_m^{-},T_n^{+}}+{\left[A\right]}_{mn}^{T_m^{-},T_n^{-}}+{\left[\mathrm{\Phi }\right]}_{mn}^{T_m^{+},T_n^{+}}+{\left[\mathrm{\Phi }\right]}_{mn}^{T_m^{+},T_n^{-}}+{\left[\mathrm{\Phi }\right]}_{mn}^{T_m^{-},T_n^{+}}+{\left[\mathrm{\Phi }\right]}_{mn}^{T_m^{-},T_n^{-}},{\left[Z_{PS}\right]}_{mn}={\left[A\right]}_{mn}^{T_m^{+},T_n^{+}}+{\left[A\right]}_{mn}^{T_m^{-},T_n^{+}}+{\left[\mathrm{\Phi }\right]}_{mn}^{T_m^{+},T_n^{+}}+{\left[\mathrm{\Phi }\right]}_{mn}^{T_m^{-},T_n^{+}},\\ {\left[Z_{SP}\right]}_{mn}={\left[A\right]}_{mn}^{T_m^{+},T_n^{+}}+{\left[A\right]}_{mn}^{T_m^{+},T_n^{-}}+{\left[\mathrm{\Phi }\right]}_{mn}^{T_m^{+},T_n^{+}}+{\left[\mathrm{\Phi }\right]}_{mn}^{T_m^{+},T_n^{-}},\\ {\left[Z_{SS}\right]}_{mn}={\left[D\right]}_{mn}^{T_m^{+},T_n^{+}}\cdot \delta \left(T_m^{+},T_n^{+}\right)+{\left[A\right]}_{mn}^{T_m^{+},T_n^{+}}+{\left[\mathrm{\Phi }\right]}_{mn}^{T_m^{+},T_n^{+}},\end{array}$$where $m,n=1,2,3,\dots ,N$. $\delta \left(i,j\right)$ equals to 1 only if $i=j$; otherwise, it is equal to 0. ${\left[X\right]}_{mn}^{T_m^{\pm },T_n^{\pm }}\left(X=D,A\text{or}\mathrm{\Phi }\right)$ is the $\left(m,n\right)$-th entry in $\left[Z\right]$ of (6) corresponding to the $n$th basis function and $m$th test function, and $T_m^{\pm }$ and $T_n^{\pm }$ are the tetrahedron integral domains. Note that ${\left[D\right]}_{mn}^{T_m^{\pm },T_n^{\pm }}$ only appears in the $\left[Z_{PP}\right]$ and $\left[Z_{SS}\right]$, and half-SWGs are all defined in $T_m^{+}\left(T_n^{+}\right)$. The calculation of the integral formulation of ${\left[X\right]}_{mn}^{T_m^{\pm },T_n^{\pm }}$ 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 $T$ 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 $T$ 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 $T_n^{+}$; then, it searches the rest of the tetrahedrons as $T_n^{-}$. Note that a tetrahedron has 4 faces. Hence, there are twice the loops of searching that needs to be done with $4\sum \left(1+2+\cdots +T-1\right)=4T\left(T-1\right)/2$ times. Then, the complexity of the computation time is $O\left(N^2\right)$, where $N=cT$ and $c$ is a constant. The pseudocode for the brute-force method is shown in Figure 2(a).

[figures omitted; refer to PDF]

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 $O\left(N\log N\right)$ [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 $ti$) of a tetrahedron (denoted as $i$) is set as a key; (3) find if the key ($ti$) 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 ($i$). It is obvious that two loops for searching all the faces need $4T$ times and searching a member in T-map in step (3) needs $\log 4T$ times at most. Hence, the total computation time is $4T\log 4T$; then, the complexity of the computation time by map is $O\left(N\log N\right)$.

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 $O\left(N^2\right)$ ($N$ 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 $O\left(N\log N\right)$. For the VEFIE matrix equations in (6), the acceleration of MVM by MLFMA can be expressed by [17] $$\begin{array}{}\text{(9)}& \left[Z\right]\cdot \left[I\right]=\left[Z_{\mathrm{n}\mathrm{e}\mathrm{a}\mathrm{r}}\right]\cdot \left[I\right]+\left[D\right]\cdot \left[T\right]\cdot \left[A\right]\cdot \left[I\right],\end{array}$$where $\left[Z_{\mathrm{n}\mathrm{e}\mathrm{a}\mathrm{r}}\right]$ is the matrix of near field interactions and $\left[D\right]$, $\left[T\right]$, and $\left[A\right]$ 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 $\left[Z_{\mathrm{f}\mathrm{a}\mathrm{r}}^{m\times n}\right]$ is approximated by $\left[{\stackrel{\text{~}}{Z}}^{m\times n}\right]$ with the product of two low-rank matrices $\left[U^{m\times r}\right]$ and $\left[V^{r\times n}\right]$ as [18] $$\begin{array}{}\text{(10)}& \left[Z_{\mathrm{f}\mathrm{a}\mathrm{r}}^{m\times n}\right]\approx \left[{\stackrel{\text{~}}{Z}}^{m\times n}\right]=\left[U^{m\times r}\right]\left[V^{r\times n}\right]=\sum _{i=1}^r{u}_i^{m\times 1}{v}_i^{1\times n},\end{array}$$where $m\times n$ is the matrix dimension and $r$ is the effective rank of the matrix and satisfies $r<\min \left(m,n\right)$. Combining (9) and (10), the hybrid MLFMA-ACA scheme for the MVM calculation can then be expressed as $$\begin{array}{}\text{(11)}& \left[Z\right]\cdot \left[I\right]=\left[Z_{\mathrm{n}\mathrm{e}\mathrm{a}\mathrm{r}}\right]\cdot \left[I\right]+\sum _{l=2}^k\left[Z_{\mathrm{M}\mathrm{L}\mathrm{F}\mathrm{M}\mathrm{A}}\right]\cdot \left[I\right]+\sum _{l=k+1}^L\left[Z_{\mathrm{A}\mathrm{C}\mathrm{A}}\right]\cdot \left[I\right],\end{array}$$where $\left[Z_{\mathrm{M}\mathrm{L}\mathrm{F}\mathrm{M}\mathrm{A}}\right]$ and $\left[Z_{\mathrm{A}\mathrm{C}\mathrm{A}}\right]$ denotes the MLFMA and ACA work matrices, respectively. $l=0,1,2,\dots ,L$ represent the different levels of an octree structure. Since the first two levels $\left(l=0,1\right)$ are all near field interaction groups, the MVM acceleration starts at $l=2$. 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 $k+1$ to $L$ in the lower levels without LF breakdown and MLFMA acceleration from 2 to $k\left(2\le k\le L\right)$ in higher levels to maintain high efficiency. If $k=L$, it means that only MLFMA is used. Generally $k$ is determined according to the box size of the $k$th level which is less than a preset value, e.g., $0.1{\lambda }_0$ (${\lambda }_0$ 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 $O\left(N^2\right)$, while it is $O\left(N\log N\right)$ 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 omitted; refer to PDF]

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 ${\epsilon }_{r1}=1.5$–$1.5j$ and ${\epsilon }_{r2}=1.0$–$1.0j$. The layered sphere is illuminated at $\theta =0^{°}$ by an incident plane wave at 500 MHz. For consistency, the standard mesh size is set as $0.1{\lambda }_0/\sqrt{{\epsilon }_r}$, 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.

[figures omitted; refer to PDF]

[figure omitted; refer to PDF]

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 $R$ 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–$0.005j$, and those of the middle layer are 4 mm and 1.5–$0.035j$, respectively. The array is symmetrically located in the sphere center $O$ along the $x$-axis, and the corresponding 8 dipoles placed with equal intervals ($∆=10\mathrm{m}\mathrm{m}$) are fed with an equal excitation amplitude, equal orientation ($\theta ={90}^{°}$, $\Phi =0^{°}$), and $0.25\pi$ 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 omitted; refer to PDF] [figure omitted; refer to PDF] [figure omitted; refer to PDF]

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 $r_1=0.25\mathrm{m}$, an outer radius of $r_2=0.5\mathrm{m}$, and a relative permittivity of ${\epsilon }_{r1}=2.0-0.1j$, ${\epsilon }_{r2}=4.0-0.2j$. The two-layer sphere is illuminated at $\theta =0^{°}$ by an incident plane wave at 300 MHz. The mesh sizes of the inner and outer parts are $0.027{\lambda }_0$ and $0.025{\lambda }_0$, 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 $0.025{\lambda }_0$, thereby a 6-level octree ($L=5$) is constructed. For the hybrid MLFMA-ACA method, the lowest levels where $l$ is from 4 to 5 are accelerated by ACA, while the higher levels where $l$ 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 $O\left(N^2\right)$ to $O\left(N\log N\right)$. 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.

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

[1] R. F. Harrington, Field Computation by Moment Methods,DOI: 10.1109/9780470544631, 1993.

[2] A. J. Poggio, E. K. Miller, "Integral equation solutions of three-dimensional scattering problems," Computer Techniques for Electromagnetics: International Series of Monographs in Electrical Engineering,DOI: 10.1016/B978-0-08-016888-3.50008-8, 1970.

[3] A. F. Peterson, S. L. Ray, R. Mittra, Institute of Electrical, Electronics Engineers, Computational Methods for Electromagnetics, Vol. 2, 1998.

[4] L.-M. Zhang, 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, DOI: 10.1109/TAP.2016.2647686, 2017.

[5] L.-M. Zhang, 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, DOI: 10.1109/TAP.2015.2490254, 2015.

[6] N. A. Ozdemir, 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, DOI: 10.1109/TMAG.2007.891031, 2007.

[7] C.-C. Lu, 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, DOI: 10.1002/mop.4650071013, 1994.

[8] E. Darve, 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, DOI: 10.1098/rsta.2003.1337, 2004.

[9] H. Wallen, J. Sarvas, "Translation procedures for broadband MLFMA," Progress In Electromagnetics Research, vol. 55, pp. 47-78, DOI: 10.2528/PIER05021001, 2005.

[10] M. Bebendorf, "Approximation of boundary element matrices," Numerische Mathematik, vol. 86 no. 4, pp. 565-589, DOI: 10.1007/PL00005410, 2000.

[11] Q.-M. Cai, Y.-W. Zhao, Y.-T. Zheng, M.-M. Jia, Z. Zhao, 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, DOI: 10.1109/TAP.2015.2481925, 2015.

[12] D. Schaubert, D. Wilton, 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, DOI: 10.1109/TAP.1984.1143193, 1984.

[13] D. Wilton, S. Rao, A. Glisson, D. Schaubert, O. Al-Bundak, 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, DOI: 10.1109/TAP.1984.1143304, 1984.

[14] E. P. Becker, "Working draft, standard for programming language c++," 2006.

[15] N. M. Josuttis, The C++ Standard Library: A Tutorial and Reference, 2012.

[16] T. H. Cormen, C. E. Leiserson, R. L. Rivest, C. Stein, Introduction to Algorithms Second Edition, 2001.

[17] W. C. Chew, E. Michielssen, J. Song, J.-M. Jin, Fast and Efficient Algorithms in Computational Electromagnetics, 2001.

[18] K. Zhao, M. N. Vouvakis, 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, DOI: 10.1109/TEMC.2005.857898, 2005.