Application of Hierarchical Two-Level Spectral Preconditioning Method for Electromagnetic Scattering from the Rough Surface

The higher-order hierarchical Legendre basis functions combining the electrical field integral equations (EFIE) are developed to solve the scattering problems from the rough surface. The hierarchical two-level spectral preconditioning method is developed for the generalized minimal residual iterative method (GMRES). The hierarchical two-level spectral preconditioner is constructed by combining the spectral preconditioner and sparse approximate inverse (SAI) preconditioner to speed up the convergence rate of iterativemethods.Themultilevel fast multipolemethod (MLFMM) is employed to reducememory requirement and computational complexity of the method of moments (MoM) solution. The accuracy and efficiency are confirmed with a couple of numerical examples.


Introduction
The study of electromagnetic scattering from an object above a rough surface has a large number of applications, for example, in remote sensing, radar surveillance, and so on [1][2][3][4][5][6][7][8][9][10][11][12][13][14].Consequently, accurate and efficient investigation of the composite scattering problem is an important subject in computational electromagnetics and remote sensing society.Both the approximate and rigorous methods have been developed to tackle this problem.Among the approximate methods, some are based on small perturbation method (SPM) [12], Kirchhoff approximation (KA) [13], small slope approximation (SSA) [14] and so on.However, the height value must be very small compared to the electromagnetic wave length in the SPM, the radius of the surface must be larger than a wavelength in the KA, the slope must be small and the height must be moderate for the first order in the SSA.
For rigorous methods, the method of moments (MoM) [15][16][17][18][19] is one of the most widespread and generally accepted techniques for analyzing rough surface scattering.In MoM, the RWG functions [18] are often used for representing unknown current distributions due to the fact that it is convenient to model objects with arbitrary shape using triangular patches.When iterative solvers are used to solve the MoM matrix equation, the fast multipole method (FMM) or multilevel fast multipole method (MLFMM) [16,17] can be used to accelerate the calculation of matrix-vector multiplies.
For large rough surfaces, the RWG functions have a poor convergence and need a large number of unknowns for a desired accuracy.To circumvent this disadvantage, a remedy is to employ higher order basis functions [19][20][21][22][23][24][25].The high order functions can be categorized as interpolatory or hierarchical.The Hierarchical functions allow for much greater flexibility [22][23][24][25].The basis of order  is a subset of the basis of order  + 1, which allows mixing of different order bases in the same mesh.Andersen and Volakis proposed a class of hierarchical tangential vector finite elements (TVFE's) for FEM discretization [22].Webb also introduced a set of hierarchical vector basis functions of arbitrary order for triangular and tetrahedral finite elements [23].Jorgensen et al. proposed a new set of higher order hierarchical Legendre basis functions for a quadrilaterals element for integral equations solved with MoM [24].The higher order hierarchical Legendre basis functions can preserve almost perfect orthogonality while enforcing continuity of the normal current component across surface element, which can provide a much better condition number of the MoM matrix than other higher order hierarchical basis functions and can remarkably reduce the computational costs compared with the RWG basis function.Yang et al. successfully used higher order hierarchical Legendre basis functions to analyze the electromagnetic scattering from breaking waves at low grazing angles [25].
The system matrix in MoM from the divergenceconforming higher order basis functions is often an illconditioned matrix and result in the low convergence of the Krylov iterative method.Simple preconditioners like the Diagonal preconditioner (Diag) or symmetrical successive over-relaxation preconditioner (SSOR) can be effective only when the matrix has some degree of diagonal dominance [26].Incomplete LU (ILU) preconditioners have been successfully used on nonsysmmetric dense systems in [27,28].Pan and Sheng also developed the multilevel inverse-based ILU preconditioning technique based on MLFMA for effectively analysis of electromagnetic scattering from complex targets [29].However, the factors of the ILU preconditioner may become very ill-conditioned, and consequently the performance is very poor.An effective SAI preconditioner suitable for implementation in the FMM context has also been proposed [30], which is based on a Frobenius-norm minimization with a priori sparsity pattern selection strategy.Malas and Gürel proposed an efficient parallel SAI preconditioner combined with MLFMA in [31].Fan et al. proposed the preconditioning matrix interpolation technique for the fast analysis of the radar cross-section (RCS) over a broad frequency band [32].In [32], both the near field impedance and SAI preconditioning matrices are interpolated at intermediate frequencies over a relatively large frequency band with rational function interpolation technique.
Although the SAI preconditioner can improve matrix conditions by clustering most of the large eigenvalues close to one but still leave a few close to the origin, which can potentially slow down the convergence of Krylov methods.Inherited from the basic idea of the traditional two-grid cycle of the multigrid methods [33,34], the high frequency components of the iteration error belong to the subspace spanned by the eigenvectors associated with the large eigenvalues of the system matrix and can be represented on a fine grid defined by the higher order basis functions.The low frequency components left in the error are smooth enough and then can be represented on a coarse grid defined by the lower order basis functions.
In this paper, we apply similar ideas to improve the quality of a prescribed SAI preconditioner based on the higher order hierarchical basis functions.Due to the property of the hierarchical basis functions, the hierarchical two-level spectral preconditioning method is used to accelerate the convergence of the iterative method [35].In this method, the SAI preconditioner based on the higher order hierarchical basis functions is used to damp the high frequencies of the error, while the spectral preconditioner [36] in a two-level manner based on the lower order hierarchical basis functions is applied to smooth the low frequencies of the error.
The goal of this paper is to achieve accurate and fast analysis of the electromagnetic wave scattering from a rough surface.The remainder of this paper is organized as follows.In Section 2, theory and formulations are discussed.Numerical results are presented in Section 3 to demonstrate the accuracy and efficiency of the presented method.Section 4

Theory and Formulation
2.1.Rough Surface Modeling and EFIE.As shown in Figure 1, a perfectly electric conductor (PEC) cubic is located above a PEC rough surface.The rough surface is a Gaussian stationary stochastic process with a zero mean value, which can easily be generated by a spectral method.In order to eliminate the edge effects caused by truncation of the finite surface length, the tapered wave is usually employed [1].A TE-polarized tapered wave can be expressed as where  is the widith of the tapered wave,   and   are elevation angle and azimuth angle of the incident wave, respectively. is truncated length of rough surface.The width of the tapered wave should be large enough to illuminate upon the surface [2,5].When there is no object above the PEC rough surface, the induced current ⃗    ( ⃗ ) on the rough surface is where ⃗   ( ⃗ ) denotes the incident magnetic field.When there is no object above the rough surface, the scattering field ⃗    ( ⃗ ) above the rough surface can be written as where ⃗  is the position point on the rough surface  2 and When the object exists above the PEC rough surface, the scattering field ⃗    ( ⃗ ) above the rough surface can be written as where ⃗   ( ⃗ ) is the induced current on the surface  1 of object, ⃗ ⃗   ( ⃗ , ⃗   ) denotes the dyadic Green's function in rough surface [37] as the following: The first term at the right side of (7) represents the scattering field without object above the rough surface, the second term represents the coupling interactions between the rough surface and the object.

Higher Order Hierarchical Legendre Basis Functions.
Equations ( 7) can be solved by MoM.The factor that most limits the capability of the MoM is the number of unknowns.
As was pointed out in the introduction, this limitation can be relaxed by using higher order basis functions.For a given accuracy, the use of higher order basis functions allows us to use larger triangular patches to discrete the objects.In this paper, higher order hierarchical Legendre basis Functions defined on curvilinear quadrilateral surface [24,25] are chosen to represent the surface current ⃗   as where ⃗   and ⃗  V are the co-variant unitary vectors as ⃗   =  ⃗ / and ⃗  V =  ⃗ /V.
Without loss of generality, we consider only -directed currents where are unknown coefficients,   and  V are the basis orders along the -directed current flow direction and the transverse direction, respectively.  (V) are Legendre polynomials, and P () are the modified Legendre polynomial, which can be defined as C and   are the scaling factors and are defined as Finally, applying Galerkin's method to (8) in MoM, a matrix equation can be obtained as: In (13), the element of the impedance matrix A is evaluated with the presented hierarchical basis functions.The matrix equation can be solved by the restart GMRES [38] iterative method using MLFMM to accelerate the matrix-vector multiply.

Hierarchical Two-Level Spectral Preconditioner.
Though the use of high order basis functions reduces the number of unknowns, it still increases the condition number of the matrix system and decelerates the convergence of Krylov iterative methods.In order to speed up the convergence rate of Krylov methods, preconditioning techniques are employed to transform the EFIE matrix equations into an equivalent form [32] where A ℎ is the EFIE impedance matrix associated with the higher order hierarchical Legendre basis functions, and M ℎ is a improved SAI preconditioner based on near-field elements of A ℎ in MLFMA [30].
Although the SAI preconditioner can improve matrix conditions by clustering most of the large eigenvalues close to one but still leave a few close to the origin, which can potentially slow down the convergence of Krylov methods.In this paper, the hierarchical two-level spectral preconditioning method is used to improve the quality of a prescribed SAI preconditioner.It is based on the higher order hierarchical International Journal of Antennas and Propagation basis functions, which is used to damp the high frequencies of the error, while the low frequencies of the error is eliminated by a spectral preconditioner [36] in a two-level manner defined based on the lower order hierarchical basis functions.
Let U  be a set of eigenvectors of dimension   ×  associated with the smallest eigenvalues of the preconditioned matrix M ℎ A  , where A  denotes the EFIE impedance matrix associated with the lowerorder hierarchical basis functions.The second level spectral preconditioner can then be defined as where from U  by appending zeros and I ℎ and I  are unit matrix of dimension  ℎ and , respectively.The superscript  denotes the transpose and conjugate of a given complex matrix.
Combining this second spectral preconditioner with the prescribed preconditoner in the two-level manner, a hierarchical two-level preconditioner is derived and the linear system to be solved can be transformed into In order to evaluate approximations to smallest eigenvectors for constructing the second spectral preconditioner M  , the SAI preconditioned GMRES-DR algorithm [39] is used to solve the coarse grid system based on the lower order hierarchical basis functions, which generates approximations to eigenvectors as a byproduct.The approximate eigenvectors in GMRES-DR span a small Krylov subspace and so are generated in a compact form where U  is a   × matrix, whose columns span the subspace of approximate eigenvectors, U +1 is the same except an extra column and H  is a full  + 1 by  matrix.Note that (17) allows access to both the approximate eigenvectors and their products with SAI preconditioned matrix M ℎ A  while requiring storage of only  + 1 vectors of length   .Therefore, in order to construct the spectral preconditioner in a two-level manner, one can simply replace V  with U  by appending zeros, T  with H  , where H  is a full  by  matrix by discarding last row of H  .Therefore, a hierarchical two-level spectral preconditioned GMRES algorithm is presented for solving systems of linear equations, as follows: (i) Construct SAI preconditioning matrix M ℎ based on the higher order hierarchical basis functions.
(ii) Solve the coarse grid system based on the lower order hierarchical basis functions with the SAI preconditioned GMRES-DR algorithm, and use generated eigenvector information to construct the spectral preconditioning matrix M  .
(iii) Construct the two-level spectral preconditioning matrix M = M  M ℎ , and use it for GMRES algorithm to solve linear systems.

Numerical Experiments
In this section, the EFIE linear systems based on the presented hierarchical higher order Legendre basis functions are solved with MLFMM accelerated Krylov iterative methods.All numerical experiments are performed on a personal computer with Intel(R) Core(TM)2 Quad Q9500 2.99 GHz CPU and 4 GB RAM in single precision.The restarted version of GMRES () algorithm and the inner-outer Flexible GMRES algorithm (FGMRES) [40] are used as the iterative solvers.
is taken to be 30 in GMRES () algorithm.The inner and outer restart numbers of FGMRES are both 10.In the FGMRES algorithm, the required precision for the inner and outer iteration is 1.E-2, 1.E-3, respectively.Additional details and comments on the implementation are given below: (i) zero vector is taken as initial approximate solution for all examples and all systems in each example (ii) the iteration process is terminated when the normalized backward error is reduced by 10 −3 for all the examples.
The first example is a PEC orthogonal plane scatterer consisting of each plane of size (8 × 8) as shown in Figure 2, where  is the free-space wavelength.The scatterer is discretized with curvilinear quadrilateral patches for 2and 3-order hierarchical Legendre basis functions with 11502 and 9072 unknowns, respectively.The average mesh sizes are chosen as 0.3 for 2 order basis functions, and 0.5 for 3order basis functions.The incident angles of plane wave are   = 90 ∘ ,   = 45 ∘ .The bi-static RCS of the orthogonal plane scatterer with vertical polarization is given in Figure 1 for 2-and 3-order hierarchical Legendre basis functions.The results from FEKO software are also given in Figure 2 in order to verify the accuracy of the proposed the high order hierarchical Legendre basis functions.Good agreement is achieved between these two kinds of results.It can be seen that the use of the higher order hierarchical Legendre basis functions allows the number of unknowns to be reduced by a factor of 1.27.
In order to accelerate the convergence rate of GMRES algorithm for the EFIE system based on the hierarchical Legendre basis functions, the hierarchical two-level spectral preconditioning techniques are developed.Some numerical results are given to demonstrate the effectiveness of the hierarchical two-level spectral preconditioner for the solution of EFIE system.Figures 3 and 4 show the convergence histories of the hierarchical two-level spectral preconditioned GMRES, FGMRES and GMRES algorithms for the orthogonal plane scatterer in bistatic RCS computation.In the hierarchical two-level spectral preconditioned GMRES,  = 20 approximate smallest eigenvectors, extracted by GMRES-DR algorithm during solving EFIE system, are employed to build the two-level spectral preconditioner M  based on 1-order hierarchical Legendre basis functions for the first example.Table 1 further gives the construction time and solution time for different algorithms for orthogonal plane scatterer in bistatic RCS computation.In Table 1, "−" means that constructing time of the preconditioner is very few and   shown in Table 1, though the additional construct time of the hierarchical two-level spectral preconditioner is needed, the hierarchical two-level spectral preconditioned GMRES decreases the solution time by a factor of 3.0 compared to  FGMRES, a factor of 25 compared to GMRES for 2-order hierarchical Legendre basis functions.It is also found that the hierarchical two-level spectral preconditioned GMRES decreases the solution time by a factor of 4.0 compared to FGMRES, a factor of 53.0 compared to GMRES for 3-order ones.It can be concluded that the hierarchical two-level spectral preconditioning method improves the condition number of linear system based on high order hierarchical Legendre basis functions.Next, a PEC paraboloid antenna is considered in Figure 5.The size of aperture for the paraboloid antenna is 4.8 and the focus is 4.813.The paraboloid antenna is discretized with curvilinear quadrilateral patches for 2-and 3-order hierarchical Legendre basis functions with 19992 and 17652 unknowns, respectively.The average mesh sizes are also chosen as 0.3 for 2 order basis functions, and 0.5 for 3order basis functions.The incident angles of plane wave are   = 0 ∘ ,   = 0 ∘ .Figure 5 gives the bi-static RCS with vertical polarization for the paraboloid antenna.It is very clear that the results from 2-and 3-order hierarchical Legendre basis functions agree well with the results from FEKO.Figures 6 and 7 give the convergence histories of the hierarchical two-level spectral preconditioned GMRES, FGMRES and GMRES algorithms for the paraboloid antenna in bistatic RCS computation.Table 2 further lists the construction time and solution time for different algorithms for the paraboloid antenna.In these computations,  = 20 approximate smallest eigenvectors are employed to build the two-level spectral preconditioner M  .It can be observed that the hierarchical twolevel spectral preconditioning method can greatly improve the convergence.Compared to the GMRES, the hierarchical two-level spectral preconditioned GMRES decreases the total time by a factor of 7.5 for 2-order hierarchical Legendre basis functions, a factor of 14.0 for 3-order ones.
Finally, a PEC missile target above a rough surface in Figure 8 is considered to further demonstrate the hierarchical two-level spectral preconditioning method, in which a Gaussian PEC rough surface with the following Gaussian spectrum [41] is used.Here,   and   are the correlation lengths in and -directions, respectively, and ℎ is the rms height of the rough surface.As shown in Figure 8, the dimensions of the rough surface are 16 × 16, with ℎ = 0.1 and  the hierarchical two-level spectral preconditioning method based on high order hierarchical Legendre basis functions is very efficient in analyzing electromagnetic scattering from rough surfaces.It can also be observed that the hierarchical two-level spectral preconditioning method based on 3-order hierarchical Legendre basis functions improves the condition number of linear system better than that based on 2-order ones.

Conclusions
In this paper, the hierarchical two-level spectral preconditioning technique based on higher order hierarchical Legendre basis functions is presented for solving EFIE for scattering from conducting objects.The use of higher order basis functions leads to a reduction in the number of unknowns without compromising the accuracy of geometry modeling.Numerical experiments are performed and the comparison is made with the GMRES and FGMRES.It can be found that the proposed hierarchical two-level spectral preconditioning technique is more efficient and can significantly reduce the overall computational cost.It can be found that the hierarchical two-level spectral preconditioning technique is more efficient and can significantly reduce the overall computational cost.

Figure 1 :
Figure 1: Geometric model of EM scattering from the rough surface.

Figure 2 :Figure 3 :
Figure 2: (a) Geometry of a PEC orthogonal plane scatterer, (b) The bistatic RCS with VV polarization for the orthogonal plane scatterer.

Figure 4 :
Figure 4: Convergence history of the hierarchical two-level spectral preconditioned GMRES, FGMRES and GMRES algorithms for 3order hierarchical Legendre basis functions for the orthogonal plane scatterer.

FocusFigure 5 :
Figure 5: (a) Geometry of a PEC paraboloid antenna.(b) The bistatic RCS with VV polarization for the paraboloid antenna.

Figure 7 :
Figure 7: Convergence history of the hierarchical two-level spectral preconditioned GMRES, FGMRES and GMRES algorithms for 3-order hierarchical Legendre basis functions for the paraboloid antenna.

InternationalFigure 8 :Figure 9 :Figure 10 :
Figure 8: (a) Geometric model of EM scattering from the rough surface (b) The bistatic RCS with VV polarization for the missile target above a rough surface.

Table 1 :
Construction time and solution time (in seconds) for the orthogonal plane scatterer.

Table 2 :
Construction time and solution time (in seconds) for the paraboloid antenna.

Table 3 :
Construction time and solution time (in seconds) for the missile target above a rough surface.The tapered wave is incident at the elevation angle   = 30 0 and the azimuth angle   = 0 0 .The width of the tapered wave is 4.0.The dimension of the missile is 1.96×10.5×1.96.The height of the missile is 10 above the rough surface.The PEC missile target and rough surface are both discretized with curvilinear quadrilateral patches for 2and 3-order hierarchical Legendre basis functions with 55744 and 31836 unknowns, respectively.Figure8also gives the bistatic RCS for VV polarization of different methods, which shows a reasonably good agreement.Figures9 and 10give the convergence histories of the hierarchical two-level spectral preconditioned GMRES, FGMRES and GMRES algorithms for the missile target above a rough surface.Table3further shows the construction time and solution time for different algorithms for missile target above a rough surface.In these computations,  = 20 approximate smallest eigenvectors are also chosen to build the two-level spectral preconditioner M  .It can be observed that the hierarchical two-level spectral preconditioned GMRES decreases the solution time by a factor of 6.6 compared to FGMRES, a factor of 19.4 compared to GMRES for 2-order hierarchical Legendre basis functions.
=   = 1.5.It is also found that the hierarchical two-level spectral preconditioned GMRES decreases the solution time by a factor of 7.3 compared to FGMRES, a factor of 64.0 compared to GMRES for 3-order ones.It can be concluded that