375458-3-T FAST SPECTRAL DOMAIN METHOD FOR HYBRID FINIE ELEMENTBOUNDARY INTEGRAL MODELING OF DOUBLY PERIODIC STRUCTURES Hewlett-Packard EEsof Division 1400 Fountaingrove Pkwy Santa Rosa, CA Thomas F. Eibert John L. Volakis August 1998 375458-3-T = RL-2518

PROJECT INFORMATION PROJECT TITLE: REPORT TITLE: U-M REPORT No.: CONTRACT START DATE: END DATE: DATE: Algorithms for Phased Arrays Adaptive Integral Method for Hybrid FE-BI Modeling of 3D Doubly Periodic Structures 375458-3-T June 1, 1996 May 30, 1999 August 17, 1998 SPONSOR: Dave Wilson M/S 2US-H Hewlett-Packard Co. EEsof Division 1400 Fountaingrove Parkway Santa Rosa, CA 95403-1799 SPONSOR CONTRACT No.: U-M PRINCIPAL INVESTIGATOR: CONTRIBUTORS TO THIS REPORT: John L. Volakis EECS Dept. University of Michigan 1301 Beal Ave Ann Arbor, MI 48109-2122 Phone: (313) 764-0500 FAX: (313) 747-2106 volakis @umich.edu http://www-personal.engin.umich.edu/-volakis/ T. Eibert, J. Volakis

Fast Spectral Domain Algorithm for Hybrid Finite Element/Boundary Integral Modeling of Doubly Periodic Structures Thomas F. Eibert and John L. Volakis Radiation Laboratory, EECS Department The University of Michigan, Ann Arbor, MI 48109-2122 Abstract A new fast integral equation algorithm is introduced for an efficient evaluation of the boundary integral (BI) termination of hybrid finite element (FE)/BI methods applied to threedimensional doubly periodic structures. The proposed method is referred to as fast spectral domain algorithm (FSDA) since it makes use of the spectral Green's function representation to evaluate the matrix-vector products carried out within an iterative solver. FSDA avoids explicit generation of the usual fully populated method of moments (MoM) matrix. Instead, at each iteration the actual current distribution is summed up in the spectral domain and the spectral Floquet mode series for evaluation of the BI is carried out only once per testing function. Thus, given a fixed set of modes, FSDA leads to a central processing unit (CPU) time complexity of O(N) (N: number of BI unknowns) for the evaluation of the matrix-vector products in the iterative solver. In this case, the memory requirement is also of 0(N) since only the precomputed values of the spectral Floquet mode series for each basis function are stored. Validation and timing results of FSDA are given in this report and compared with results obtained by a conventional BI formulation as well as a BI formulation based on the 1

Adaptive Integral Method (AIM). In contrast- to AIM and other fast integral methods, FSDA does not require use of the original near-zone MoM matrix elements. Thus, it leads to considerable speed-ups of about two orders of magnitude even for relatively small numbers of unknowns. 1 Introduction Application of hybrid FE/BI methods to infinite periodic structures (antennas or frequency selective surfaces) [1, 2, 3, 4] provides full 3D modeling flexibility. Basically, the FE method is employed to model a unit cell representing the array, whereas the BI provides for a rigorous mesh truncation at the upper and/or lower surfaces of the discretized unit cell. Practical arrays can usually be analyzed by employing the appropriate half-space periodic Green's functions for the boundary integral. Because of the dense nature of the BI MoM matrix, the CPU requirements of the hybrid FE/BI method are dominated by the BI's 0(n2) complexity. Therefore, considerable speed-up of the method can only be achieved by accelerating the BI termination. Over the past few years, attention has been given to developing fast integral algorithms such as the fast multipole method (FMM) [5, 6, 7, 8, 9] or AIM [10, 11], both are aimed at accelerating the execution of the pertinent matrix-vector products within iterative solvers. The FMM achieves its speed-up by employing a multipole and plane wave expansion of the free-space Green's function to calculate the far-zone interactions efficiently for clusters of basis functions. AIM achieves its speed-up by recasting the basis functions on a uniform grid making use of the Toeplitz properties for the resulting BI matrix to calculate the far-zone interactions in the discrete Fourier domain. That is, in both techniques, the far-zone matrix elements are not explicitly generated. This is possible because in the discrete Fourier domain, the source 2

and observation points are decoupled whereas in the FMM the coupling between source and observation points is calculated via the centers of the expansion function clusters. The same type of decoupling observed in the discrete Fourier domain is inherent to the conventional spectral domain formulation. In the case of periodic arrays, the spectral FSDA simply takes advantage of mode decoupling to realize a fast alternative evaluation of the BI. FSDA starts with the conventional Floquet mode representation of the BI termination; however, the BI matrix elements are not explicitly calculated. Instead, at each iteration the Fourier transforms of the basis functions multiplied with their actual expansion coefficients are summed up. The spectral integral (Floquet mode series) is then calculated only once for every testing function to evaluate the matrix-vector products in the iterative solver. In contrast to AIM or FMM, this procedure does not even require an explicit calculation of the near-coupling BI matrix elements. That is, FSDA is completely free of any system matrix. Therefore, memory demand is mostly determined by the storage of the Fourier transforms of the basis functions which are calculated before the iteration process is started. For fixed numbers of Floquet modes, memory demand is O(N) and CPU complexity for the calculation of the matrix-vector products is also O(N). Typically, FSDA results in considerable speed-ups of about two orders of magnitude or more, even for relatively small BI systems. The report is organized as follows. After presenting the formulation for FSDA [12] within the context of hybrid FE/BI methods for infinite periodic structures, the convergence of the Floquet mode series is discussed and validation results are given. Subsequently, timing comparisons of FSDA versus conventional and AIM accelerated BI implementations are provided. 3

2 Formulation 2.1 Basic Hybrid FE/BI Formulation The conventional implementation of the hybrid FE/BI method for doubly periodic arrays leads to a linear algebraic system of the form [1, 2, 3, 4] At ACRO A S Eint 1 0 0 Eint 1 f int 1 ACRoSS Abound 0 End + 0 Ztop 0 Ebond = fbound (1) Aos 0 Aund Ebound 0 0 Z E d fbond L bot 't Abo t ' bottom where only a unit cell of the array is discretized (see Fig. 1). The A-matrices are sparse (20 to 40 non-zero elements per row) and are associated with the FE portion of the hybrid method. The BI Z-matrices are fully populated and are associated with the edges on the top and bottom boundaries of the discretized unit cell. As is usually the case, the right-hand side vector elements f represent excitations in the FE volume or the BI apertures. At the side walls, the unit cell mesh is terminated by imposing phase boundary conditions (PBCs) in accordance with Floquet's theorem. In our implementation, volume tessellation is based on triangular prismatic finite elements [13, 14]. This results in triangular surface meshes with Rao-Wilton-Glisson basis functions for the magnetic currents on the top and bottom BI surfaces of the unit cell. We do note that for arbitrary scan angles of the array, matrices in (1) can be non-symmetric due to the PBCs and the periodic Green's function. Thus, nonsymmetric biconjugate gradient (BiCG) or generalized minimal residual (GMRES) solvers are used for solution. The major computational burden in both of these solvers is the matrixvector product represented by the left-hand side of (1). Since the total A-matrix is sparse, its complexity (per matrix-vector product) is O(Nv), where NV denotes the number of volume edges. However, the complexity of the matrix-vector products [Ztop/bot] {EI t} is 0(N2,I) 4

(Ntop/bot: number of BI unknowns in top/bottom BI surfaces) and storage requirements are of the same order (using direct solution methods). For a more efficient implementation of the method, it is therefore crucial to reduce the complexity of the BI matrix-vector products. Also, for array problems it is essential to reduce the number of the explicitly calculated BI matrix elements as far as possible since the numerical cost for generating the periodic Green's function is relatively high. That is, if a fast integral method like AIM or FMM is applied to the periodic problem, the total solution time for most practical problems will mostly be determined by the evaluation of the near-zone coupling elements needed to overcome the approximations of these fast algorithms. Therefore, FSDA being completely free of generating any BI matrix elements, is especially attractive for reducing the complexity of the matrix-vector products in array modeling. 2.2 Fast Spectral Domain Algorithm The [Ztop/bot{Eb dt} portion of the equation system in (1) is a result of discretizing the boundary integral (ejwt time dependence is assumed and suppressed throughout) H(r) = -2jo J p(r, r,). (E(rs) x n) dss + He(r) (2) S via Galerkin's method. Here, HexC(r) is a possible external excitation in the presence of a metallic interface in the periodic aperture S (incident and reflected waves). The spatial domain periodic Green's function Gp(r, rs) is given by - - 1 -o oo p e-jkoRpq Gp(r, rs) = (I + 2VV) Z: e-jkOO.pq Rp (3) "0 p=-oo q=-oo 4rRpq with Rpq = Ir - r - Ppql (4) 5

and Ppq = PPa + q Pb (5) in which Pa' Pb are the lattice vectors parallel to the xy-plane (see Fig. 1). Also, I denotes the unit dyad and r and rs specify observation and source points. As usual, ko and Zo are the free-space wavenumber and characteristic impedance, respectively. Furthermore, ktoo in (3) is given by ktoo = kxoo x + kyoo 0 = ~(ko sin t9o cos po x + ko sin 9o sin oo y), (6) where to9, po are the spherical coordinates corresponding to the scan angles of a phased array (positive sign) or the arrival angles of an incident plane wave (negative sign). In the spectral domain, (2) can be written as H(r) = -2j J Gp(k~x ky). (E(kx, ky) x n) dkxdky + He r), (7) kxky where Gp(kx, ky) denotes the spectral representation of the periodic Green's function, given by [1] 0 00 k2 - k2 kxky Gp(k,) E E k k ](kt -ktpq). (8) k xkoy o kx o - kx Here, A = p. x pbl is the cross sectional area of the unit cell's top/bottom bounding surface (see Fig. 1), kt = kTx + kyy, (9) ktpq = ktoo + - [P(Pb x z) + q(z x Pa)] (10) is the so-called reciprocal lattice vector, and kz = /ko2-kt kt, (11) 6

with the branch of the square root taken so that Re(k.) > 0, Im(k.^) < 0. Also, E(kI, ky) is the Fourier transform of the electric field intensity on the top/bottom apertures of the unit cell. Throughout the report, the "-" denotes two-dimensional Fourier transforms. In proceeding with the description of FSDA, let us first introduce the expansion Ntop/bot E(r)= E EnuWn(r), (12) n=1 where wn(r) are Whitney edge elements [13, 14] and En represent the expansion coefficients. Introducing (12) into (7) and applying MoM testing with the weighting functions bm(r) = n x wm(r), (13) which are the well-known Rao-Wilton-Glisson basis functions, we get (bm(r,H(r)) = (kpq - kpq kxpqkypq koZoAkzpq kxpqkypq k - kpq Ntop/Nbot E Enbn(kxqp, kypq) + (bm(r),He(r)), (14) n=1 where "*" denotes complex conjugation. Apart from the latter term, (14) is the mth row-vector product of [Ztop/bot]{EJbop/t) } To generate all row-vector products, we let m = 1, 2,..., Ntop/bot For a conventional implementation of (14), the summation over n along with the coefficients En is moved in front of the spectral sums (p and q). The spectral series is then evaluated for each (m, n) combination to compute all matrix elements of the BI submatrices [Ztop/bot]. That is, the spectral series is evaluated N2 + N2 times. Similarly, in an iterative solver, the evaluation of the pertinent matrix-vector products requires N2 + N2bt multiplications. To evaluate the matrix-vector products [Ztop/bot]{Etbopnt} using O(Ntop/bot) CPU time complexity (fixed set of Floquet modes assumed), the following steps are instead executed: 7

Step 0: Precompute the Fourier transforms of the basis functions bn (k,, ky) and of -2j kGp(k;, ky) for the discrete wavenumbers k- = kxp and ky = kyp. Step 1: Compute the iteration field vector Ntop/bot Vpq = Enbn(kxp kyq) (15) n=l in the spectral domain for all (p, q) modes. If the totality of the (p, q) modes is assumed to be constant, this step requires CINtop/bot operations, where the constant C1 depends on the number of included modes. Step 2: Compute Wp, =-2jtGp(kxp, kyq) Vpq (16) To for all (p, q) modes. For a constant number of Floquet modes, this step requires a constant number of C2 operations. Step 3: Finally, compute bm (kxpq, kypq). Wpq m = 1,2,.., Ntop/Nbot (17) p q to complete the evaluation of [Ztop/bot]{Ebodt}. Again, for a constant number of the (p, q) modes, the evaluation of (17) for all m = 1, 2,..., Ntopbot requires C3Ntop/bot operations, with the constant C3 dependent on the number of included modes. Summing up all operations needed for the evaluation of [Ztop/bot] EbOndtt} within the iterative solver, we obtain CiNtop/bot + C2 + C3Ntop/bot operations, i. e. a CPU time complexity of O(Ntop/bot). Also, C4Ntop/bot operations are needed in Step 0 before the iteration process is started. Therefore, compared to other fast integral methods like FMM or AIM, where the near-coupling elements must be calculated explicitly, this overhead with CPU time complexity O(Ntop/bot) is extremely low. As seen later, the overhead reduction for matrix fill is the major 8

reason for the speed advantage of FSDA. The memory demand of FSDA is mostly determined by the pre-computed Fourier transforms of the basis function and is therefore CsNtop/bot, i. e. of O(Ntplbot) As noted above, the constants C1, C2, C3, C4 and C5 depend on the numbers of Floquet modes included to achieve convergence in the spectral series. Typically, for increasing unit cell sizes, increasing numbers of modes are needed. However, this is also the case for conventional BI implementations even when acceleration techniques are used like the Ewald transformation [16, 17]. Also, it should be noted that the convergence of the spectral series is determined by the total field distribution in the BI surfaces rather than the spectrum of the individual basis functions or the Green's function. Based on this observation, convergence of the spectral series can often be dramatically improved by shifting the BI surface a small distance from the inhomogeneities of the FE volume domain (see also [2]). 3 Results Example 1: Strip dipole array As a first example in demonstrating the efficiency of FSDA, we calculated the power reflection coefficient of the strip-dipole array shown in Fig. 2 which was proposed in [19]. We applied FSDA with different numbers of Floquet modes and the results are compared with data obtained by a MoM code involving a multilayered media Green's function [20]. In the figure, FSDA 0 refers to using only 1 Floquet mode and gives unacceptable errors; however, FSDA 1 (implying p and q range from -1 to 1, i. e. 9 terms of the modal representation) is already converged. The small frequency shift between FE/BI data with FSDA acceleration and the reference MoM results [20] was also observed for the conventional FE/BI implementation. 9

For this example, the unit cell was meshed using triangular prismatic elements resulting in a total of 55370 volume unknowns and 4880 BI unknowns in each of the top and bottom BI surfaces. The corresponding CPU time requirements on a Pentium II PC/266 MHz were 244 sec, 281 sec, and 330 sec per frequency for p and q ranging from 0 to 0, -1 to 1, and -2 to 2, respectively. Example 2: Slot array The second example is the slot array depicted in Fig. 3 and also analyzed in [19]. Again, we present FSDA calculations for different numbers of Floquet modes and contrast them to MoM results [20]. For FSDA, convergence is now achieved when p and q range from -5 to 5 (FSDA 5). However, data using FSDA 3 and FSDA 4 (not shown in the diagram) are very close to this curve, as well. In comparison to the previous example, the larger number of Floquet modes is necessary due to the "more complicated" field distribution in the top BI surface which is exactly in the slot plane. As noted earlier, the convergence behavior can be considerably improved if the top BI surface is slightly shifted away from the slot. This is shown by an additional curve in Fig. 3 obtained when an air layer of thickness 0.5 mm and modeled by a single vertical prism is included on top of the slot. Then, FSDA 1 (p and q ranging from -1 to 1, named FSDA 1+) agrees almost exactly with the MoM data. For this example, the first unit cell mesh contained a total number of 3361 volume edges and the second mesh with the additional air layer involved 4961 volume edge. For the first mesh, the number of BI edges was 1240 at the bottom surface and 161 over the slot plane on the top BI surface. For the second mesh, with the air layer added at the top BI, both apertures consisted of 1240 unknowns. The corresponding CPU times on a Pentium II PC/266 MHz were 7.5 sec, 11.0 sec, 34.0 sec, and 11.8 sec per frequency for FSDA 1, FSDA 2, FSDA 5, and FSDS 1+, respectively. That is, by a slight shift of the top BI surface away from the slot, not only was the accuracy of 10

the results improved, but the CPU time was also reduced by nearly a factor of 3. A better understanding of the convergence behavior of FSDA can be obtained by looking at the spectra of the search vectors of the iterative solver, i. e. Vpq = Enl/bot Enbn(kp, kyq). Fig. 4 shows the normalized magnitude of the spectral distribution bn(kx, ky) of one cartesian component of a Rao-Wilton-Glisson basis function. In contrast to this relatively broad spectrum, the spectral distribution of the search vectors at the bottom BI surface (shown in Fig. 5) is much narrower. This type of spectrum results in a very fast convergence of the Floquet mode series used in FSDA. If the BI surface is in the plane of the slot, the spectral distribution of the search vectors is determined by the shape of the slot, especially in the direction along the width of the slot, and the spectra are relatively broad (see Fig. 6). However, on moving the BI surface away from the slot the spectrum of the search vectors is narrower and consequently the convergence of the corresponding Floquet mode series is faster. Example 3: Frequency selective surface (FSS) of slot-coupled microstrip patches In Fig. 7, the geometry of an FSS array consisting of slot-coupled microstrip patches is depicted. This array was investigated in [21] and acts as a strongly resonant band-pass structure. Because of its resonant characteristics, modeling of the fields between the patches with a local method like FEM is a difficult task and the mesh density must be relatively high. Also, it should be noted that the incidence angle Vo varies with frequency according to the waveguide measurement setup in [21]. Fig. 8 compares the transmission coefficients obtained by FE/BI-AIM and FE/BI-FSDA 1+ solutions to measured data presented in [21]. The agreement between the three curves is quite good and we note that the FSDA results are even closer to the measured data than those obtained by the FE/BI-AIM approach. With respect to complexity, the FE/BI-AIM solution used 98636 volume and 5308 BI unknowns in each of the top and bottom BI surfaces. Also, in the FE/BI-AIM approach, the BI was implemented 11

as a conventional spatial domain formulation with Ewald acceleration of the Green's function modal series [16, 17] and AIM was employed to speed-up the BI portion of the code [10, 11]. The required CPU time for the FE/BI-AIM implementation on a Pentium II PC/266 MHz was between 3 and 4 hours per frequency. On using FSDA, the CPU times were considerably decreased even though the number of unknowns was increased by putting additional air layers between the microstrip patches and the BI termination. Specifically, the number of volume unknowns grew to a value of 137804 and the number on BI unknowns was increased to 14828 in each of the top and bottom BI surfaces. The required CPU time on the PC for FSDA 3+ (p and q ranging from -3 to 3.) was about 50 minutes per frequency. CPU timings for FSDA A final look at the CPU time dependence of FSDA as a function of BI unknowns is given in Figs. 9 and 10. These CPU times correspond to a dipole array with metallic backing (BI only on top of the mesh) and were obtained at a frequency very close to resonance with plane wave excitation. To obtain CPU times as a function on the unknown count, several elements were used to form a unit cell mesh. The largest investigated unit cell contained a 3 x 3 array with one prism layer used across the height of the mesh. Again, the calculations were performed on a Pentium II PC/266 MHz and the FSDA CPU curves are compared to those obtained by a conventional and AIM accelerated BI implementation [10, 11]. Further, the Green's function Floquet mode series for the conventional spatial domain (mixed potential) BI formulation were accelerated by the Ewald transformation [16, 17]. The number of terms in the FSDA Floquet mode series and in the Ewald series were chosen appropriately for the largest analyzed unit cell (p and q ranging from -9 to 9 in FSDA). In all cases, a non-symmetric BiCG solver was used which requires two matrix-vector products per iteration. Fig. 9 displays the CPU time per iteration of the BiCG solver whereas Fig. 10 illustrates the total solution time including 12

matrix fill of the three implementations. For small BI matrices [Ztop/bot], the conventional BI needs the shortest time per iteration; however, due to the O(N log N) complexity of AIM and the (9(N) complexity of FSDA (constant number of Floquet modes assumed) this is quickly changed as the number of unknowns increases. For the investigated problem sizes, AIM was always faster per iteration than FSDA but the differences got less for increasing numbers of unknowns. In terms of total solution time, FSDA was almost two orders of magnitude faster than AIM and the conventional BI, even for small numbers of unknowns. Also, the speed advantage compared to the conventional BI got increasingly more pronounced for larger problems. The total solution time ratio between FSDA and AIM was almost constant for increasing problem sizes (i. e. about two orders of magnitude). The largest portion of the total solution time of the conventional BI and AIM is due to the CPU time needed for the explicit calculation of the BI matrix elements. Therefore, the major advantage of FSDA is that it is completely free of a need to form any elements of the BI matrix. 4 Conclusions We introduced a new method, the fast spectral domain algorithm (FSDA) for mesh termination in a hybrid FE/BI formulation as applied to doubly periodic structures. FSDA was derived starting from the well-known spectral domain formulation of the planar BI. It was shown that FSDA leads to large CPU time reductions as compared to conventional BI implementations and other fast integral methods such as AIM. The improved speed of FSDA is mainly due to that no matrix elements (including those in the near-zone) are explicitly calculated. Therefore, a speed-up of about 2 orders of magnitude is obtained even for small numbers of BI unknowns (as small as 500). For a fixed set of Floquet modes in the spetral series, CPU time complexity 13

per iteration is 0(N) and memory demand is of the same order. The presented results show that 9 Floquet modes are often sufficient for practical problems if the FE mesh is chosen appropriately. 14

References [1] D. T. McGrath and V. P. Pyati, "Phased Array Antenna Analysis with the Hybrid Finite Element Method," IEEE Trans. Antennas Propagat., vol. 42, no. 12, pp. 1625-1630, Dec. 1994. [2] E. W. Lucas and T. W. Fontana, "A 3-D Hybrid Finite Element/Boundary Element Method for the Unified Radiation and Scattering Analysis of General Infinite Periodic Arrays," IEEE Trans. Antennas Propagat., vol. 43, no. 2, pp. 145-153, Feb. 1995. [3] D. T. McGrath and V. P. Pyati, "Periodic Structure Analysis Using a Hybrid Finite Element Method," Radio Science, vol. 31, no. 5, pp. 1173-1179, Sep/Oct. 1996. [4] J. D'Angelo and I. Mayergoyz, "Pased Array Antenna Analysis," in Finite Element Software for Microwave Engineering, Edited by T. Itoh, G. Pelosi, and P. P. Silvester, John Wiley & Sons, Inc., 1996. [5] N. Engheta, W. D. Murphy, V. Rohklin, and M. S. Vassiliou, "The Fast Multipole Method (FMM) for Electromagnetic Scattering Problems," IEEE Trans. Antennas Propagat., vol. 40, no. 6, pp. 634-641, Jun. 1992. [6] R. Coifman, V. Rohklin, and S. Wandzura, "The Fast Multipole Method for the Wave Equation: A Pedestrian Prescription," IEEE Antennas Propagat. Mag., vol. 35, pp. 7-12, Jun. 1993. [7] N. Lu and J.-M. Jin, "Application of the Fast Multipole Method to Finite-Element Boundary-Integral Solution of Scattering Problems," IEEE Trans. Antennas Propagat., vol. 44, no. 6, pp. 781-786, Jun. 1996. 15

[8] S. S. Bindiganavale and J. L. Volakis, "Comparison of Three FMM Techniques for Solving Hybrid FE-BI Systems," IEEE Antennas Propagat. Mag., vol. 39, no. 4, pp. 47-60, Aug. 97. [9] J. M. Song, C. C. Lu, and W. C. Chew, "Multilevel Fast Multipole Algorithm for Electromagnetic Scattering by Large Complex Objects," IEEE Trans. Antennas Propagat., vol. 45, no. 10, pp. 1488-1493, Oct. 1997. [10] E. Bleszynski, M. Bleszynski, and T. Jaroszewicz, "AIM: Adaptive Integral Method Compression Algorithm for Solving Large-Scale Electromagnetic Scattering and Radiation Problems," Radio Science, vol. 31, no. 5, pp. 1225-1251, Sep./Oct. 1996. [11] H. T. Anastassiu, M. Smelyanskiy, S. Bindiganavale, and J. L. Volakis, "Scattering from Relatively Flat Surfaces Using the Adaptive Integral Method," Radio Science, vol. 33, no. 1, pp. 7-16, Jan.-Feb. 1998. [12] T. F. Eibert and J. L. Volakis, Fast Spectral Domain Algorithm for Rapid Solution of Integral Equations, Electronic Letters, vol. 34, no. 13, pp. 1297-1299, Jun. 25, 1998. [13] T. Ozdemir and J. L. Volakis, "Triangular Prisms for Edge-Based Vector Finite Element Analysis of Conformal Antennas," IEEE Trans. Antennas Propagat., vol. 45, no. 5, pp. 788-797, May 1997. [14] R. D. Graglia, D. R. Wilton, A. F. Peterson, and I.-L. Gheorma, "Higher Order Interpolatory Vector Bases on Prism Elements," IEEE Trans. Antennas Propagat., vol. 46, no. 3, Mar. 1998. [15] Y. Saad, Iterative Methods for Sparse Linear Systems, Boston: PWS Pub. Co., 1996. 16

[16] P. P. Ewald. Die Berechnung optischer und elektrostatischer Gitterpotentiale, Ann. Phys. 64, pp. 253-287, 1921. [17] K. E. Jordan, G. R. Richter, P. Sheng, An Efficient Numerical Evaluation of the Green's Function for the Helmholtz Operator on Periodic Structures, J. of Comp. Phy., 63, pp. 222-235, 1986. [18] B. Houshmand, W. C. Chew, and S.-W. Lee, "Fourier Transform of a Linear Distribution with Triangular Support and Its Applications in Electromagnetics," IEEE Trans. Antennas Propagat., vol. 39, pp. 252-254, Feb. 1991. [19] R. Mittra, C. H. Chan, and T. Cwik, "Techniques for Analyzing Frequency Selective Surfaces - A Review," Proc. IEEE, vol. 76, no. 12, pp. 1593-1615, Dec. 1988. [20] D. R. Wilton and D. R. Jackson, Personal Communications, 1997. [21] R. Pous, D. M. Pozar, A Frequency-Selective Surface Using Aperture-Coupled Microstrip Patches, IEEE Trans. Antennas Propagat., vol. 39, no. 12, pp. 1763-1769, Dec. 1991. 17

Figure Captions Figure 1: Doubly periodic array configuration and unit cell with BI termination on the top and bottom surfaces and periodic boundary conditions at the vertical side faces. Figure 2: TE power reflection coefficient of a strip dipole array as presented in [19]. For FSDA 0, 1, and 2, p and q range from 0 to 0, -1 to 1, and -2 to 2, respectively. The MoM results are from [20]. Figure 3: TM power reflection coefficient of a slot array as presented in [19]. For FSDA 1, 2, and 5, p and q range from -1 to 1, -2 to 2, and -5 to 5, respectively. For FSDA 1+, p and q range from -1 to 1 with the top BI surface placed one prism layer above the slot. MoM results are from [20]. Figure 4: Normalized magnitude of the x-component spectrum of a Rao-Wilton-Glisson basis function of the unit cell in Fig. 3 (p, q are according to (8)). Figure 5: Normalized magnitude of the x-component spectrum of the magnetic current search vector in the bottom BI surface of the unit cell in Fig. 3 (p, q are according to (8)). Figure 6: Normalized magnitude of the x-component spectrum of the magnetic current search vector in the top BI surface of the unit cell in Fig. 3 (BI imposed in plane of the slot, p, q are according to (8)). Figure 7: FSS unit cell of the aperture coupled microstrip patches as suggested in [21]. The parameters are or = 2.2, d = 1.6 mm, a = 36.07 mm, b = 34.04 mm, Ws = 2 mm, Ls = 8 mm, Wp = Lp = 28 mm. Figure 8: TE transmission coefficient of the FSS structure in Fig. 7 for different FE/BI models compared to measured results from [21]. (p = 0~, do varying from 57~ to 32~ for f varying from 2.5 GHz to 4.0 GHz, according to waveguide measurement setup in [21]. 18

Figure 9: CPU times per iteration as a function of BI unknowns for FSDA. conventional and AIM accelerated BI implementations. The CPU times were obtained on a Pentium II PC/266 MHz and refer to a microstrip dipole array with one prism layer. Figure 10: Total solution times as a function of BI unknowns for FSDA, conventional and AIM accelerated BI implementations. The CPU times were obtained on a Pentium II PC/266 MHz and refer to a microstrip dipole array with one prism layer. 19

- I --- - -7,- I W / 'I / I x I/ I-PBC PBC BI Figure 1: 20

1.0 5 5 5 5. - p c) C0 c) so 0.8 0.61 V — + r LW E X Er=2 4II AlOmml m -- ^^"0 1\ I I I I 'I I I I I I I I!! I I I IA I I ---— FSDA 0 ------ FSDA 1 FSDA 2 MoM Wilton et al. 0.4 0.2 I II II II I 0.0 2 a 2 m 16 m 10 14 18 Frequency in GHz Figure 2: 21

1.0 -I = 0.8 " / *3 \ \\ /, g * \1.875 mm U 0.6 \l, \\ l, E 9 0.4.FSDA1 I: -FSDA '' FS0.2 FS, DA1 ' ': — FSDA5 CIA'FSDA I + ~r=4.M O Wilton et al. \ - 10 mm -* 006 10 14 18 22 Frequency in GHz Figure 3: 22

15 10 Figure 4: 23

I I. I I I 1 I I I I I I I I -I1- - I I I I I __ -I- -- I I 1 I I I I I "I -I I - - - I - I >- I - I I —.I - -I I I -K I r- I I I - I N. I N. I 51 -15 -10 wo~ p Figure 5: 24

Figure 6: 25

Wp 7 mq x d d Figure 7: 26

1.0 10a Q 0 0~ 0.8 0.6 0.4 0.2 0.0 Frequency (GHz) Figure 8: 27

101 A 10' 10-1 I I I I I I I I i I I ~ I i i i i -i 1 t-t-4 —4-4 - I I I I I I FSDA 9 AIM r 7 conv. BI dr 103 104 Number of BI unknowns Figure 9: 28

1O -- _____ _ 1 0 2 ~ - _ _ _ _ _ 11O3 Figre10 104 29