030601-5-T Simulation of Patch and Slot Antennas Using FEM with Prismatic Elements and Investigations of Artificial Absorber Mesh Termination Schemes J. Gong, T. Ozdemir, J. Volakis, and M. Nurnberger National Aeronautics and Space Administration Langley Research Center Hampton, VA 23681-0001 July 1995 30601-5-T = RL-2424

Report #030601-5-T NASA Langley Grant NAG 1-1478 Grant Title: Simulation of Conformal Spiral Slot Antennas on Composite Platforms Report Title: Simulation of Patch and Slot AntenR nas Using FEM with Prismatic Elements and Investigations of Artificial Absorber Mesh Termination Schemes Report Authors: Primary University Collaborator: Primary NASA-Langley Collaborator: University Address: Date: J. Gong, T. Ozdemir, J. Volakis and M. Nurnberger John L. Volakis Fred Beck Telephone: (804) 864-1829 Radiation Laboratory Department of Electrical Engineering and Computer Science Ann Arbor, Michigan 48109-2122 Email: volakis@umich.edu Phone: (313)764-0500 July 1995 1

NASA Grant No. NAG-1-1478 Simulation of Spiral Slot Antennas on Composite Platforms Project Goals * Develop robust finite element codes for narrow slot and other antennas on doubly curved platforms. * Investigate feeding techniques for conformal slot antennas and develop improved feed models in the context of FEM simulations. * Develop antenna miniaturization techniques using high contrast low loss materials and new anisotropic material * Deliver codes capable of analyzing a wide class of conformal antenna antennas without a need to use sophisticated meshing packages. Initial Project Plan Year 1 * Develop the necessary codes for the analysis of slot antennas, including spirals * Fabricate and measure slot spiral and examine new feeding techniques. * Validate new codes * Develop mesh truncation schemes using new ABCs and artificial absorbers Year 2 * Examine miniaturization techniques for slot spirals and other printed antennas using high contrast and anisotropic material * Measure performance of miniaturization procedures * Develop slot antenna analysis techniques for doubly curved platforms * Examine effects of curvature on antenna resonance characteristics Year 3 * Integrate high frequency codes and FEM codes for pattern prediction on complex platforms * New antenna development with multi-frequency and multi-function characteristics * Further code development with interactive features Year 1 Progress Summary Our year 1 progress can be characterized with four major achievement which are crucial toward the development of robust, easy to use antenna analysis code on doubly conformal platforms. Specifically, * A new FEM code was developed using prismatic meshes. This code is based on a new edgebased distorted prism and is particularly attractive for growing meshes associated with printed slot and patch antennas on doubly conformal platforms. It is anticipated that this technology will lead to interactive, simple to use codes for a large class of antenna geometries. Moreover, the codes can be expanded to include modeling of the circuit characteristics. The attached i

report describes the theory and validation of the new prismatic code using reference calculations and measured data collected at the NASA Langley facilities. The agreement between the measured and calculated data is indeed impressive even for the coated patch configuration. * A scheme was developed for improved feed modeling in the context of FEM. A new approach based on the voltage continuity condition was devised and successfully tested in modeling coax cables and aperture fed antennas. An important aspect of this new feed modeling approach is the ability to completely separate the feed and antenna mesh regions. In this manner, different elements can be used in each of the regions leading to substantially more improved accuracy and meshing simplicity. * A most important development this year has been the introduction of the perfectly matched interface (PMI) layer for truncating finite element meshes. So far the robust boundary integral method has been used for truncating the finite element meshes. However, this approach is not suitable for antennas on non-planar platforms. The PMI layer is a lossy anisotropic absorber with zero reflection at its interface. The Univ. of Michigan was first to make use of this absorber in FEM simulations with remarkable success. We expect widespread use of this absorber. * Quite unexpectedly, during year 1 (instead of year 3 as originally planned) we were able to interface our antenna code FEMACYL (for antennas on cylindrical platforms) with a standard high frequency code. This interface was achieved by first generating equivalent magnetic currents across the antenna aperture using the FEM code. These currents were employed as the sources in the high frequency code. As shown in the attached figures, the calculations are in very good agreement with measured data. As planned, during year-1 we also fabricated a slot spiral antenna for operation between 1.1 Ghz and 2.5 Ghz. Of most importance is the stripline feed used in the design. Three models of the feed were investigated before arriving at a successful feed design. As shown in the preliminary measurements, the slot spiral performed very well, having a return loss of -30dB over the entire design bandwidth. Papers Submitted for Publication 1. T. Ozdemir and J.L. Volakis, "Triangular prisms for edge-based vector finite element analysis," Submitted to IEEE Trans. on MTT 2. J. Gong and J.L. Volakis, "An efficient and accurate model of the coax cable feeding structure for FEM simulations," accepted for publication in IEEE Trans. on AP 3. J. Gong, et. al., "Performance of an anisotropic artificial absorber for truncating finite element meshes" submitted to IEEE Trans. on AP. 4. J. Gong and J.L. Volakis," Optimal Selection of a Uniaxial Artificial Absorber Layer for Truncating Finite Element Meshes" submitted to IEEE Microwave and Guided Wave Letters Papers in Preparation 1. J. Gong and J.L. Volakis, " FEM modeling of slot spirals using prismatic elements" 2. J. Gong and J.L. Volakis, "Mixed element FEM modeling of aperture coupled stacked patch antennas of arbitrary shape" 3. T. Ozdemir, J.L. Volakis and M. Gilreath, "Resonance characteristics of coated patch antennas on cylindrical platforms." 4. M. Numberger and J.L. Volakis, "A hybrid FEM/High frequency approach for antenna pattern calculations on complex platforms" In all, a total of eight papers resulted out of the Year-1 and related work. Our Year-2/3 work will proceed as planned. Basically, we will now concentrate on antenna design and miniaturization techniques for developing small scale UHF and VHF antennas. ii

Table of Contents Page Triangular Prisms for Edge-based Vector Finite Element 1 Analysis of Conformal Antennas 1. Introduction 2 2. Vector Edge-Based Basis Functions 3 3. Eigenvalue Computation 8 Rectangular Cavity 11 Circular Cavity 11 Pie-shell Cavity 11 4. Radiation and Scattering from Cavity Backed Structures 13 Radiation from dipoles in a cavity on planar platform 13 Plane wave scattering from a cavity on planar platform 14 Radiation from a patch antenna on cylindrical platform 14 5. Input Impedance Computations 16 Rectangular patch on planar platform 19 Circular patch on planar platform 19 Rectangular and circular patches on circular-cylindrical platform 22 6. Conclusions 24 7. Appendix 24 Efficient Finite Element Simulation of Spiral Slot Antennas Using 28 Triangular Prism Elements 1. Introduction 28 2. Edge-Based Triangular Prism Elements 29 Shape Functions: Type One 29 Shape Functions: Type Two 33 3. Mesh Generation and Code Interface 36 4. Preliminary Results 37 5. Appendix 38 Performance of an Anisotropic Artificial Absorber for Truncating 41 Finite Element Meshes 1. Introduction 41 2. Formulation 42 3. Results 44 Rectangular Waveguide 44 Microstrip Line 45 Band Eliminator 45 Sphere Scattering 45 4. Conclusions 45 Optimal Selection of a Uniaxial Artifical Absorber Layer for Truncating 56 Finite Element Meshes 1. Introduction 56 2. Uniaxial Absorbing Layer 57 3. Application to Shielded Microstrip Line 57 4. Conclusion 58 iii

Table of Contents Page A Hybridization of Finite Element and High Frequency Methods for 62 Pattern Prediction of Antennas on Aircraft Structures 1. Introduction 63 2. Finite Element Code Description 64 FEMA-CYL 64 FEMA-TETRA 66 FEMA-PRISM 67 3. High Frequency Code Descriptions 68 GTD/UTD Code 68 APATCH Code 69 4. Comparison of Measurments and Calculations 70 5. Conclusions 72 iv

Triangular prismns for edge-based vector finite element analysis of conformrnal antennas T. Ozdemir and J. L. Volakis Radiation Laboratory Department of Electrical Engineering and Computer Science University of Michigan Ann Arbor, Michigan 48109-2122 June 5, 1995 Abstract This paper deals with the derivation and validation of edge-based shape functions for the distorted triangular prism. Although the tetrahedron is often the element of choice for volume tessellation, mesh generation using tetrahedra is cumbersome and CPU intensive. On the other hand, the distorted triangular prism allows for meshes which are unstructured in two dimensions and structured in the third dimension. This leads to substantial simplifications in the meshing algorithm and many printed antenna and microwave circuit geometries can be easily tessellated using such a mesh. The new edge-based shape functions presented in this paper are validated by eigenvalue, radiation and scattering, and input impedance computations involving conformal antennas on planar and cylindrical platforms. Also presented are results showing the effect of platform curvature on resonance behavior of such antennas. 1

1 Introduction The brick and tetrahedron are popular elements in finite element analysis of electromagnetic problems. The first is indeed attractive because of its simplicity in constructing volume meshes whereas the tetrahedron is a highly adaptable, fail-safe element. It is often the element of choice for three dimensional (3D) meshing but requires sophisticated and CPU-intensive meshing packages. The distorted prism (see Figure 1) is another volume element which provides a compromise between the adaptibility of the tetrahedron and the simplicity of the brick. Basically, the distorted prism allows for unstructured meshing (free-meshing) on a surface and structured meshing in the third dimension. An approach for growing prismatic meshes is illustrated in Figure 2 and most volumetric regions in antenna and microwave circuit analysis can be readily tessellated using such a mesh. As seen, once the surface mesh at the different surface levels is constructed, the prismatic mesh can be generated by simply connecting the nodes between the adjacent mesh surfaces. This avoids use of CPU-intensive volumetric meshing packages and in many cases the user/analyst can construct the mesh without even resorting to a surface meshing package. Examples include printed circuits and antennas on planar platforms. Moreover, because of their triangular cross-section, the prisms overcome modeling difficulties associated with bricks at corners formed by planes or edges intersecting at small angles. A special case of the distorted prism is the right prism which is characterized by the right angles formed between the vertical arms and the triangular faces [1]. The top and bottom faces of the right prism are necessarily parallel and equal, restricting them to a limited range of applications, namely, geometries with planar surfaces. In contrast, the distorted triangular prism is almost as adaptable as the tetrahedron with the exception of cone-tips which are not likely to occur in printed antenna and microwave circuit configurations. In this paper, we introduce edge-based basis functions for the most general distorted prisms. These prisms have non-parallel triangular faces and each of their three vertical edges can be arbitrarily oriented. In the following, we first present the edge-based shape functions and then proceed with the derivation of the finite element matrix. Eigenvalue, radiation and scattering, and input impedance computations are used to validate the new basis functions for prismatic elements. Also included are data showing the effect of platform 2

ZA x Figure 1: The distorted triangular prism shown with the directions of the edge vectors. curvature on the resonance behavior of conformal patch antennas. 2 Vector Edge-based Basis functions Consider the distorted prism shown in Figure 3. The prism's top and bottom triangular faces are not necessarily parallel to each other and the three vertical arms are not perpendicular to the triangular faces. On our way to specifying a set of shape functions for the prism, we proceed with the identification of a triangular cross-section of the prism which can be uniquely defined given a point (x, y, z) interior to the prism. The cross-section contains this point as illustrated in Figure 3. A way to specify the nodes of such a triangular cross-section is to use the parametric representation r b + s(rt - rib), i 1,2,3, (1) where ri are the nodal position vectors of the triangular cross-section (see Figure 4(a)) whereas rrt and rib are associated with the top and bottom triangular faces, respectively, and 0 < s < 1. Thus, for s = 0 r specify the 3

4 L Surface normal-E.. '. S urface geometric fidelity Distorted triangular prisms stacked up Figure 2: Illustration of an approach for constructing prismatic meshes (unstructured surface mesh and structured mesh along the third dimension) Perspective view Top view I', \'7, I' \ I" ' <.v s^^ X, / I' \; /-: I; Interior point (x,y,z) Triangular cross-section containing the interior point Figure 3: Variation of the triangular cross-section along the length of the prism 4

(x 3 y3, Z 3t~ 3t ( 3 (,.y.(X, YtZ 2? A (y" / /' (x' Y31z *.; (x,- y3 z, y, / d, (xy ' ~) ' - cross-section --- ( h1 d 3 3"b3 (a) (b) Figure 4: (a) Nodal coordinates, (b) triangular cross-section with the local coordinates ~ and q. nodes of the bottom triangle and for s = 1 ri reduce to the nodes of the top triangle. Since s assumes values between these two limits, the cross-section sweeps the entire volume of the prism. Clearly, each (x, y, z) point belongs to one of the triangles specified by ri and consequently there is a unique correspondence between s and a given point in the prism's volume. To determine s in terms of (x,y,z), we first note that (1) represents nine equations with ten unknowns (s and the coordinates of the three nodal points). Thus, additional equations are needed and these are the equations defining the plane coinciding with the triangle, viz. + + i = 0 = 1, 2 3 (2) a b c x+ + - + = 0 (3) a b c in which the constants a, b, and c correspond to the locations at which the plane intersects the x,y, and z axes, respectively. Substituting in (2) the values of xi, yi and zi obtained from (1) yields three equations involving four unknowns a, b, c and s. Solving these three linear equations for -, 1 and 1 in terms of s and substituting them in (3), we end up with the third order polynomial a3s3 + a2s2 + als + ao = 0 (4) whose real root is the desired value of s. The expressions for the coefficients ai can be given in closed form in terms of x, y, and z, but are rather lengthy 5

and bear no significance on the rest of the analysis. Therefore, they have been omitted but can easily be derived as outlined above. The appropriate solution of (4) in terms of ao,1,2,3 is [2] a2 s = (S31 + 2)-3 (5) 3a3 where S1 = [r +(q3+ r2)]2 (6) 82 = [r-(q3+r2)2]3 (7) 1 1 r = 6a (aa2- 3aoa3) - (a2/a3)3 (8) 1 1 2 q = al- a2 (9) and this completely specifies s in terms of x, y, and z. We next proceed with the derivation of the basis functions. We choose to represent the field variation across the triangular cross-section (defined by s) using the two-dimensional Whitney form [3]. A simple linear variation will be assumed along the length of the prism. Specifically, the vector basis functions for the top triangle edges can be expressed as N1 = d (L2VL3-L3VL2)s N2 d2(L3VL -L1VL3)s (10) N3 =d3(L1IVL2-L2VL1)s and correspondingly those for the bottom triangle edges will be M1 = d (L2VL3- L3VL2)(1 - s) M2 = d2(L3VL1-LiVL(3)( - s) (11) M3 = d3(LlVL2-L2VLi)(1 -s). The subscripts in these expressions identify the edge numbers as shown in Figure 1 and the distance parameters di are equal to the side lengths of the triangular cross-section containing the observation point (see Figure 4(b)). Also, 1 6

/ A '- (x,y,z) 2 1 'A2 (a) (b) (c) A Figure 5: (a) Vector map of N 2 or M2 (b) Vector map of K1 (c) Variation of v as a function of x, y, and z. COS &3 sinal /3 L2(Q,7) = ~h (12) cos a2 sin a2 L3(,7) = -, —^+ h 7 h3 h3 are the usual two dimensional scalar node-based basis functions [4] for the same triangle with ai denoting the interior vertex angles and hi being the node heights from the opposite side. The variables c and r represent the local coordinates and are illustrated in Figure 4(b). As required with all edge-based shape functions, N, * e^ and Mi * e^i have unity amplitude on the ith edge whereas N * Cj M * ej = 0 for i Z j. Their vector character is depicted in Figure 5(a) and as seen they simply "curl" around the node opposite to the edge on which their tangential components become unity. It remains to define the shape functions for the three vertical edges and we chose to express these by the linear representations KAs) = Vb^ )Lio (i ) K2( ) = v(, 7) 2, ) (13) K3(g, 7) = V(,71) L3 71) As before, Li are the node-based shape functions defined in (12) and a pic 7

torial description of K1 is found in Figure 5(b). Of particular importance in (13) is the unit vector v(^, ). It is a linear weighting of the unit vectors 1i, v2 and i3 associated with the vertical arms (see Figure 5(c)), and is given by 3 ) Li((, ) (14) i=l This particular choice of v minimizes tangential field discontinuity across the faces. Contrary to the rectangular brick, tetrahedron and right prism, the edge-based vector basis functions for the distorted prism do not ensure tangential field continuity across the faces, nor are they divergenceless. The same is true for the curvilinear bricks [5] and, in the case of the distorted prism, the field discontinuity and divergence increase with surface curvature (see Figure 2). This is certainly an undesirable feature but in most cases (particularly when sampling at 15 or so points per linear wavelength and the radii of curvature are greater than a wavelength), the angular deviation of the vertical arms is quite small. Consequently, for all practical purposes the field discontinuity and divergence are negligible. 3 Eigenvalue Computation In this section, we examine the validity of the presented edge-based functions. Specifically, we consider the eigenvalues of three different cavities using the edge-based distorted prism as the tessellation element. We begin by first deriving the matrix elements following Galerkin's testing. The weighted residuals of the vector wave equation are J N. (V xV xE- kE)dV =0, i=1,2,3, (15) JJJ " * (V xV xE - k2E)dV = 0, = 1,2,3, (16) f vKi(V x VxE- k2E)dV=O, i=1,2,3, (17) in which Ni, Mi and Ki comprise the nine edge-based vector basis functions defined in the previous section and E is the electric field vector. 8

The matrix equations are generated by introducing the representation 3 E = E [ENNi(r) (r)+ MM i(r+ EKKi(r)] (18) i=l where EiN, EiM and ElK are the expansion coefficients, and correspond to the average amplitudes of the field vector at the ith edge. Substituting (18) into (15)-(17), and invoking the divergence theorem yields the element equations 3 3 S E-N[NNCj - k2NNDij] + EjM[NMCij - k2NMDij] + j=l j=l 3 5 EjK[NKC j - kNKDij] = 0 (19) j=1 3 3 E E,[MNCj - k2MNDj] + E M[[MMC j - k2MMDij] + j=l j=l 3 EjK[MIKCj - k2MKDj] = 0 (20) j=1 3 3 EjN[IKNJC - k2KNDij] + EM[IKMCt - ko KMDii] + j=1 j=l 3 EJK[KKICij- k2KKDi] = 0, i = 1,2,3. (21) j=1 where NNCie (V x N) i)(V x N)dV (22) NMCie = f / (V x Ni) * (V x Me)dV (23) NKCig = jJ (V x Ni). (V x K)dV (24) MMCe = J (V x M.) (V x Me)dV (25) MKC7 = J /i/ (V x M.) * (V x Ke)dV (26) KKC = I I (V x Ki) * (V x Ke)dV (27) 9

Region 1 (O < z < z,) Region 2 ( z< < zz2) Region 3 (z2< z < zz) Triangular cross-section Quadrilateral cross-section Triangular cross-section fz z3 / 4\ / I \ / I \ / \ \ / \/ / X - " I _ Figure 6: The three regions of integration over the volume of the distorted prism. - - Pism volume, integration subregions, subregion cross-section NNDe = J Ni * Ned V (28) NMDi = I Ni.* MdV (29) NKDI = N* KAdV (30) i = I MMdV (31) KKDI = J K i* K d. I (33) quadrilateral cross-section. In the case of the right prism, the integrals have closed form expressions and are given in the Appendix. Upon assembly of (19)-(21) and boundary condition enforcement, we ob 10

cm Mode k,cm1 % Error 0.5cm (Exact) Prism Brick Tetra.::- S^ —: TE 101 5.236 -0.99 -1.36 0.44 - TM1io 7.025 -4.44 -2.23 0.70 0.75cm T -73. -28. \ - - TE01 7.531 0.07 -2.58 1.00 i:;?/ ii-'Sf ----!" TE201 7.531 -0.25 -3.13 -0.56 TMlll 8.179 -0.31 -2.09 2.29 (Actual mesh) (a) (b) Figure 7: (a) Rectangular cavity discretized using right triangular prisms, (b) eigenvalues for the air-filled cavity. tain the generalized eigenvalue system [A]{x}= ko2[B]{x} (34) in which A = k2 represent the eigenvalues of the problem. The matrices [A] and [B] are square, symmetric, real and positive definite. Example 1: Rectangular Cavity The first example is the rectangular cavity shown in Figure 7(a). Results based on brick, tetrahedron [6] and prism discretizations are given in Figure 7(b), where they are compared to the exact values. We observe that, overall, the data based on the prism are better than those based on the brick and, at least, as good as those associated with the tetrahedron. Example 2: Circular Cavity The tessellation of the circular cavity (drum) is shown in Figure 8(a). The mesh was constructed by first creating a surface grid of triangles at the base of the drum. The prisms were then generated by growing the mesh along the axis of the cylinder. In this case, the resulting volume element is the right prism and the eigenvalues are given in Figure 8(b) along with the exact. Except for the higher order mode TE211, the remaining eigenvalues were computed to within three percent of the exact. Example 3: Pie-shell Cavity The third example is a pie-shell sector as shown in Figure 9(a). It is obtained by bending the rectangular cavity considered earlier and the resulting 11

Radius = lcm (Actual mesh) 4i': ---<',0=-X II Mode TMoio k, cfiml (Exact) 2.405 % Error (Computed) 1.29 TE 111 3.640 2.17 TM lo 3.830 -2.90 TMoii 3.955 0.81 TE 211 4.380 -8.97 (a) (b) Figure 8: (a) Circular cavity discretized using right triangular prisms, (b) eigenvalues for the air-filled cavity (Actual mesh) / 0.75cm Mode k, cm-1 % Error (Exact) (Computed) TM 1l 4.693 1.52 TM210 6.009 0.95 TE 1, 6.640 4.96 TE 21 7.513 -1.11 TE 01 7.579 1.71..- \ - — ' \.... #- lcm 0.5cm (70a) (a) (b) Figure 9: (a) Pie-shell cavity discretized using distorted triangular prisms, (b) eigenvalues for the air-filled cavity 12

Side view:!z Top view:! x FEM mesh ABC surface - f l ABC surface —.3 ko Figure 10: The geometrical set up for finite element-ABC formulation of radiation and scattering from cavities in infinite ground planes. volume element is the distorted prism with a vertical arm angular deviation of five degrees. The computed and exact eigenvalues for the first five dominant modes are given in Figure 9(b), and these testify to the accuracy of the distorted prisms in modeling curved geometries. 4 Radiation and scattering from cavity backed structures As the next set of tests, we have looked at radiation and scattering from cavity backed radiators recessed in an infinite planar or cylindrical metallic platform. The general set up geometry is shown in Figure 10. As it is seen, a three dimensional cavity is discretized using the triangular prisms. The finite element mesh is extended outside the cavity in all directions about a fraction of the wave length and truncated using the socond order ABC reported in [7]. Example 1: Radiation from dipoles in a cavity on planar platform We have first looked at the radiation from a pair of current elements residing inside the cavity; one horizontally and one vertically oriented as shown in Figure 11. Radiation pattern is shown in Figures 12a and b. The 13

0.3,0 '. | x (Side view) y-directed element Actual gridding z-directed element 0.525o S tYA (Top view) 0.75 Xo Figure 11: Locations of current elements with respect to the gridding inside the cavity solid and dotted line results were obtained using brick elements and closing the finite element domain at the aperture of the cavity via boundary integral (BI) method, which is an exact formulation and is explained in [11]. The triangular prism results were obtained in conjunction with the absorbing boundary condition (ABC) to truncate the mesh 0.3Ao away from the cavity aperture in all directions. As the plots clearly demonstrate, the triangular prisms together with the ABC perform well. Example 2: Plane wave scattering from a cavity on planar platform Next, we have looked at the plane wave scattering from a larger cavity whose dimensions in x,y, and z directions are 1.2Ao, 0.75AX, and 0.3Ao, respectively. The geometrical set up for the finite element-ABC formulation is the same as in the previous radiation problem which is illustrated in Figure 10. Figures 13 and 14 show bistatic and backscattering radar cross sections (RCS), respectively, of the cavity. It is seen that the FEM-ABC method using triangular prism elements exhibit sufficient accuracy in computing plane wave scattering from cavities in a ground plane. Example 3: Radiation from a patch antenna on cylindrical platform As the last test, we have looked at the radiation from a cavity backed patch antenna on a circular-cylindrical platform shown in Figure 15(a). The distorted triangular prisms used to discretize the cavity in Figure 9 have been used here also to discretize the computational domain. Finite element mesh 14

(a) (b) 0 n I I -10 -o 0 - br" i " ----.... Brick-BI (o-pol) -- Brick-BI (o-pol) \ * Triangular prism-ABC \ I__________ I -10 0"x *N Brick-BI (0-pol) \ - - Brick-BI (0-pol) * Triangular prism-ABC -20 -20 -_ I Ii..,1 I. 0 30 60 90 -J)V 0 30 60 90 Observation angle ~ (deg) Observation angle 0 (deg) Figure 12: Radiation pattern of two current elements located inside the cavity; (a) )=45~; (b) 0=45~ (a) (b) ll*... 1(0 I I 0 -10 -20 F..... t u,t -,. XN. - Brick-BI (8-pol) l - - Brick-BI (*-pol) i * Triangle prism-ABC (o-pol) * Triangle prism-ABC (s-pol) -10 F. 10 0 -20 h ^^', -ar* — A -.a |- Brick-BI (s-pol) |- -- Brick-BI (*-pol) * Triangular prism-ABC (o-pol) Triangular prism-ABC (#-pol)j -30 -30 h -40 -40 1 _ f, I I Iin I 0 30 60 90 Observation angle e (deg) 0 30 60 Observation angle, (deg) 90 Figure 13: Bistatic RCS of a cavity of dimensions 1.2 Xo x 0.75 )o x 0.3 Xo for a plane wave with its electric field vector making an angle of 45 witl the unit vector in 0 direction incident on the cavity aperutre at 0 9n=45 0 inc=30 Q (a) 0=60Q (b) 0=45 ~ 15

20,, 20 -,o -10 -10 -20 Brick-BI (o-pol) -20 '~ Brick-BI (0-pol) -30 - Brick-BI (o-pol) -30 - - Brick-BI (#-pol) Triangle prism-ABC (9-pol) Triangle prism-ABC (8-pol) -40 * Triangle prism-ABC (*-pol) -40. Triangle prism-ABC (,-pol -50 -50 0 30 60 90 0 30 60 90 Observation angle 0 (deg) Observation angle g (deg) (a) (b) Figure 14: Backscatter RCS of a cavity of dimensions 1.2Xo x 0.75 Xo x 0.3;o for an incident plane wave with its electric field vector making an angle of 30~ with the unit vector in 0 direction; (a) (=30~; (b) 0=60~. has been terminated 0.3 wavelengths away from the cavity aperture by the first order cylindrical ABC derived in [7]. Figure 15(b) shows the radiation pattern computed as compared to the more rigorous finite element-boundary integral solutions [8]. Again the results are excellent. The reference solution (FE-BI) has been proven to be accurate by comparing against experimental data as shown in Figure 16(a) which displays a rectangular patch antenna backed by a continous wraparound cavity recessed in a circular ogive cylinder [9]. Figure 16(b) shows the same antenna but this time with a dielectric superstrate covering the entire surface of the test body. Here the experimental data is compared against finite element-ABC solution (FE-ABC) which terminates the finite element mesh 0.3 wavelengths away from the surface of the cylinder [10]. Agreement is again good. 5 Input impedance computations As the last set of tests, we computed input impedances of a rectangular and a circular patch antenna backed by a cavity and recessed in an infinite metallic platform. The results are compared to more accurate boundary integral solutions for planar platforms, and the effect of curvature on the resonant 16

A z PEC infinite cylinder (Actual mesh) Feed location Cavity wall Patch ----- 0.05 Xo I 0 0 I C)o o I r I J r (a) ABC at 0.3 X0 )=30 o 0=60 o clI E 0 40 30 20 10 0 -10 -20 -30 -40 - C'4 1-1 0 to gJO 0 0) 30 20 10 0 -10 -20 -30 -40 (b) 0 30 60 90 120 150 180 Observation angle 0 (deg) 0 90 180 270 360 Observation angle 4 (deg) FE-BI (0-Pol) * FE-ABC (0-Pol) --- FE-BI ((-Pol) * FE-ABC (0-Pol) Figure 15: Radiation from a rectangular patch on circular-cylindrical platform: (a) geometry, (b) radiation pattern. 17

(a) (b) 0 0. z g-4 2 10 0 -10 -20 zla "3 1.0 0t -c o N 0 Z 10 0 -10 -20 -30 ' -180 -30 ' -180 -90 0 90 180 Angle (4) [deg] (0=90~) -90 0 90 180 Angle (~) [deg] (0=90~) y I/- ABC mesh - termination '- surface ~r=2.33 BI mesh ~r =2.09 I.I x (z=0) termination surface 1211 _ 0. 12" -1- 12" - (c) z -e 12" ^f Frequency = 3.52 GHz X= 3.36" 0.69" 1.973" Figure 16: Comparison with experimental results, (a) rectangular patch on a wraparound cavity, (b) with a dielectric superstrate, (c) test body: circular ogive cylinder. 18

frequency is observed for cylindrical platforms. Example 1: Rectangular patch on planar platform As the first test, we have looked at the input impedance of a 2cm x 3cm rectangular patch on the aperture of a 5cm x 6cm x 0.079cm cavity recessed in a planar metallic ground plane (see Figure 17). The meshing of the computation volume is the same as in Figure 10 but this time the mesh is terminated using a metal- backed isotropic absorber instead of an absorbing boundary condition. This new termination technique has already been shown to work well in these types of geometries [12]. In Figure 17(a), the results refered to as FEMA-PRISM have been obtained using prisms as the discretization volume element in conjuction with the absorber to terminate the mesh. The results refered to as FEMA-BRICK, which are the accurate ones, have been obtained using bricks as the discretization volume elements in conjuction with the boundary integral technique to terminate the mesh at the cavity aperture [11]. The error in estimating the resonant frequency (the frequency at which the resistance is maximum and the reactance is zero) is about 20 MHz. This is less than 1% of the resonant frequency (3.17 GHz), which is well within the range of accuracy of the moment method techniques. It does not make sense to compare the input resistance values at resonance since both techniques use an overly simplistic modeling of the feeding structure and the level of input resistance is very much dependant on the accuracy of the model. For each frequency computation, the metal wall at the back of tie absorber has been kept at 3cm (0.3Ao at resonance) away from the edges of and above the cavity aperture as shown in Figure 17(b). The patch is fed 0.375cm off the center along the longer dimension so the shorter edges are the radiating edges as shown in Figure 17(c) and the radiated electric field is polarized in the direction of the longer edge. Example 2: Circular patch on planar platform As the second test, we have looked at the input impedance of a circular patch of radius 1.3cm on the aperture of a circular cavity of radius 3.47cm and depth 0.45cm recessed in a planar metallic ground plane (see Figure 18). In principle, the meshing of the computation volume is the same as in Figure 10 but the surface mesh is as in Figure 8(a). As in the preceeding example, the mesh is terminated using a metal-backed isotropic absorber. In Figure 18(a), 19

Input Impedance vs. Frequency lcm - - 5cm I IE 0 x o -.= = Rr = 1-j2.7 (a) At 7 41 2I (b) 2.0 2.5 3.0 3.5 4.0 4.5 5.0 Frequency (GHz) Finite Element-Boundary Integral (FEMA-BRICK) Finite Element-Metal Backed Absorber (FEMA-PRISM) 1.5cml -- tl.5cm 1.5cm 1.5cm l /0.08cm r =2.17 Resonant frequency = 3.17 GHz Error in resonant frequency estimation < 1% Vertically directed field amplitude under the patch r ' 71 2000 V/cm 1600 1200 800 400 0 (c) Figure 17: Rectangular patch on planar platform, (a) input impedance computations, (b) antenna geometry and absorber termination boundary, (c) radiating sides of the patch. 20

Input Impedance vs. Frequency r = 1.3cm (0.15 WL at resonance) r = 2.17cm (0.25 WL at resonance) r = 2.6cm (0.3 WL at resonance) r = 1.3cm (0.15 WL at resonance) E 150 ~ 100 0 0 a4 [.0 2.5 3.0 3.5 4.0 4.5 5.0 Frequency (GHz) (a) (b) ~- Finite Element - Boundary Integral (FEMA-TETRA) - - - - Finite Element - Metal Backed Absorber (FEMA-PRISM) Resonant frequency = 3.55 GHz Error in resonant frequency estimation = 1% r= 2.4 Vertically directed field magnitude under the patch 1....:..S. 500 V/cm 400 300 200 100 0 (c) Figure 18: Circular Patch on Planar Platform, (a) Input impedance computations, (b) Antenna geometry and the absorber termination boundary, (c) Radiating sides of the patch 21

the results refered to as FEMA-TETRA, which are the more accurate ones, have been obtained using tetrahedra as the discretization volume elements in conjuction with the boundary integral technique to terminate the mesh at the cavity aperture [13]. The error in estimating the resonant frequency is about 50 MHz. This is 1% of the resonant frequency (3.51 GHz), which is again within the range of accuracy of the moment method techniques. For each frequency computation, the metal wall at the back of the absorber has been kept at 2.6cm (0.3Ao at resonance) away from the edges of and above the cavity aperture as shown in Figure 18(b). The patch is fed 0.8cm off the center of the patch so, as shown in Figure 18(c), the glowing parts of the circular patch are two apposing sides symmetrically located with respect to the center of the patch at the same angular location as the feed. In case of both examples, the agreement between the absorber and boundary integral termination results improved as the distance of the absorber termination boundary from the cavity apereture was increased, as expected. Example 3: Rectangular and circular patches on circular-cylindrical platform As the last example, we looked at the effect of the curvature on the resonance behavior of the patch antennas using prisms. We took the rectangular and circular patch antennas examined in Example 1 and Example 2, respectively, and curved the platforms circularly first along x-axis with Rx =: 5cm, Ry = oo then along y-axis with Rx = oo, Ry = 5cm where Rx and Ry are the radii of curvature along x and y-axis, respectively. While curving the platform, linear dimensions of the cavity aperture have been kept the same as in the planar case. Figures 19(a) and (b) show the computed input resistance as a function of frequency for rectangular and circular patches, respectively. Each figure contains three input resistance plots showing the shift in the resonant frequency as a result of curving the patches. First observation is that the resonant frequency increases with the curvature. Secondly, the amount of increase depends on the way the platform is curved with respect to the radiating sides of the patch: When the platform is curved in a way to bring the radiating sides closer to each other (R, = oo, Ry = 5cm), the shift in the frequency is greater than (almost twice as much as) the case in which the radiating sides stay the same distance apart as the platform is curved 22

A ---w y y 70 60 50 40 30 (a) o._:4 Radiating sides x 20 - 10 - 0 -3.10, x. R: Radius of curvature along x I R: Radius of curvature along y y I, (Af)x= 12 MHz 3.18 3.20 (Af)y = 22 MHz 3.12 3.14 3.16 Frequency (GHz) (Af) (AfDX y Radiating sides x E cy) (b) 11 1-1.c i=4 (Af)x= 70 MHz (Af)y= 150 MHz 0 L 3.0 3.2 3.4 3.6 3.8 4.0 Frequency (GHz) Figure 19: Shift in the resonant frequency due to curving, (a) rectangular patch, (b) circular patch 23

(Rx = 5cm, Ry = oc). This observation is true for both the rectangular and circular patches. 6 Conclusions The distorted prism is an indispensable tool for the analysis of many doubly curved antenna and microwave geometries. It provides simple volume meshing without compromising accuracy (see the comparison with the tetrahedron in Figure 7(b)). In this paper, we introduced a set of edge-based shape functions for the distorted prism and used them to generate the element equations. For validation purposes, the eigenvalues of three different cavities (rectangular, cylindrical and pie-shell) were computed and compared to the exact solutions. In almost all cases, the error was less than three percent. Radiation, scattering and input impedance computation results are also excellent when compared to the finite element-boundary integral solutions. Having validated the new prismatic elements and the basis functions, we used them to look at the effect of curvature on the resonance berhavior of conformal patch antennas. The result is that the resonance frequency increases with the curvature irrespective of the type of the patch. Appendix For the right prism, closed form expressions of the integrals in (22)-(33) are possible. Referring to Figure Al, we have dide cos /kn COS/ jm COS 3km ENNC-: = ( Xjm + Xkn - Xjn - c hkhn hjhm hkhm cos/3 jn 2 c2 hid1 hh~ hXkm + 3 h h h h sin/jk sintmn ) hjhn 03 hhijkhmhn,Cr e dide ( cos/kn cos/33m cos/3km ENMCe = -- ~ Xjm - Xkn + - Xjn + c hkhn hj hm hk hm COS /jn I C2h di hjhn 3m hhkhmhn k Sinhhmn ) 3lfi 0 in3~ 24

ENKCZ~ EMMC~i EMKCie EKKCie ENND,-e ENMDje ENKDie EMMDe EMKD-2e EKK Die hi di Cos/3.e - 6 h-he - ENNC-je - -ENKCie h1d,1o COS/3i = C 2 h-h1 COS /3k -hkhe tdAd = C 3 COS flkn Xkh m+ Cos05/33mk COS/3Am 3n Cos /3jn X hjhm Xn hkhm h h~km 1 - -ENNDie 2 = 0 -ENNDie =0 = Cxie where I3rs {=~ if r = s otherwise Xrs = id2 {ww.S + -i[(COta3 - COta2)(rlsWr + 7lrws) + 2(~swr + &Ws)] + 2 3 -1-[3(COta3 - COta2)(q,~r + 71rc~s) + 2YrTlIs (Cot2aC2 - COtac2cota3 + 12 cot2aC3) + 6Gr~s]}, Tr, S L1,23, the set {i, J, k}I takes on the value { 1, 2, 31, {2, 3, 1} or { 3, 1, 2}1, and the same is true for the set {l, m, n1 25

(3 l'A N2, pI,3) 7\ (1 II -- KK K, m2,.//' A,.... N3 K2 c vT X/,/, /< h1 3 y>"O~1 rsj a2 2 d3 (b) 1 M3 (a) Figure Al: Right triangular prism: (a) Indexing of nodes and edges, (b) parameters defining the triangular cross-section. W2 = ~ CS2 = 2 sin ca3 w2=0 ~2 — c~ 1 -8 h2h2 W3 = 6 C2 3 sin ca2 0 3 h3 References [1] Sacks, Z., S. Mohan, N. Buris and J. F. Lee, "A prism finite element time domain method with automatic mesh generation for solving microwave cavities," IEEE APS Int. Symp. Digest, Vol. 3, pp. 2084-87, Seattle, Washington, June 19-24, 1994. [2] Abramowitz, M and I. A. Stegun, Handbook of Mathematical Functions, Dover Publications, Inc., New York, 1972. [3] Whitney, H., Geometric Integration Theory, Princeton University Press, Princeton, New Jersey, 1957. 26

[4] Zienkiewics, 0. C. and R. L. Taylor, The Finite Element Method, 4th ed., Vol. 1, McGraw-Hill, New York, 1989. [5] Antilla, G. E. and N. G. Alexopoulos, "Scattering from complex threedimensional geometries by a curvilinear hybrid finite-element-integral equation approach," J. Opt. Soc. Ame. A, Vol. 11, No. 4, pp. 1445-57, April 1994. [6] Chatterjee, A., J. M. Jin and J. L. Volakis, "Computation of Cavity Resonances Using Edge-Based Finite Elements," IEEE Trans. Microwave Theory Tech., Vol. 40, No. 11, pp. 2106-2108, Nov. 1992. [7] Chatterjee, A., "Investigation of finite element-ABC methods for electromagnetic field simulation", Ph.D. Thesis, Chap. VI, EECS Dept., Univ. of Michigan, Ann Arbor, MI 48109-2122, September 1994. [8] Kempel, L. C., J. L. Volakis and R.J. Silva, "Radiation by cavity-backed antennas on a circular cylinder," IEE Proceedings, 1994. [9] Gilreath, M., Personal communication, NASA Langley Research Center, September 1994. [10] Kempel, L. C., "Radiation and scattering from cylindrically conformal printed antennas," Ph.D. Thesis, Chap. IV, EECS Dept., Univ. of Michigan, Ann Arbor, MI 48109-2122, April 1994. [11] Jin, J. -M. and J. L. Volakis, "A finite element-boundary integral formulation for scattering by three-dimensional cavity-backed apertures," IEEE Trans. Antenn. Prop., vol. 39, No. 1, pp. 97-104, January 1991. [12] Ozdemir, T. and J. L. Volakis, "A comparative study of an absorbing boundary condition and an artificial absorber for truncating finite element meshes," Radio Science, Vo. 29, No. 5, pp. 1255-1263, Sept. - Oct. 1994. [13] Gong, J., J. L. Volakis, A. C. Woo and H. T. G. Wang, "A hybrid finite element-boundary integral method for the analysis of cavity- backed antennas of arbitrary shape," IEEE Trans. Antenn. Prop., vol. 42, No. 9, pp. 1233-1242, Sept. 1994. 27

Efficient Finite Element Simulation of Spiral Slot Antennas Using Triangular Prism Elements Jian Gong and John L. Volakis Radiation Laboratory University of Michigan 1301 Beal Avenue, EECS Building Ann Arbor, MI 48109-2122 I. Introduction It has been reported that a hybrid finite element-boundary integral technique can be employed for characterizing conformal antennas of arbitrary shape. Indeed, planar/non-planar, rectangular/non-rectangular designs, ring slot or spiral slot antennas with probe feed, coaxial cable feed, or microstrip line feed, can be simulated with this formulation. This is because of the geometrical adaptability of tetrahedral elements used in the implementation. However, in practice, certain configurations require extremely high sampling rate to accommodate the computational domain. Among them is a multilayered geometry with thin substrate(s) or bonding film(s) in the presence of a thick spacer is a typical design, which is found attractive for increasing the bandwidth of the microstrip antenna. Another example is the thin (spiral) slot antenna, where the slot width is much smaller than other dimensions (diameter or inter-distance between spiral turns). In both cases, the mesh is extremely dense (with over 50 or even 100 samples per wavelength depending on geometry) whereas typical discritizations involve only 10-20 elements per wavelength. This redundant sampling problem is especially severe for 3-D tetrahedral meshes, where the geometrical details usually distort the tetras. 28

The numerical system assembled from this type of mesh often leads to large system conditions due to a degraded mesh quality. Also, the mesh generation is difficult and the solution CPU time is large. In this report, we propose a finite element-boundary integral formulation using edge-based triangular prism elements. It can be shown that this element choice is ideally suited for planar antenna configurations, include spirals, circular and triangular slots. Among the many advantages of these elements, the most important is the simplicity of the mesh generation. Also, it requires much smaller number of unknowns and efforts for accurate modeling of complex geometries without significant increasing the system size. As is typically the case with finite element formulations, layered inhomogeneous or anisotropic material are readily modeled. This is particularly important since the purpose of our study is to investigate antenna miniaturization techniques and resonance control. II. Edge-Based Triangular Prism Elements Edge-based elements for vector field modeling are free of spurious modes and inherently satisfy the tangential continuity boundary conditions. In this section, we describe two types of shape functions for triangular prisms and briefly derive the corresponding FEM matrix elements. Shape Functions: Type One Consider a triangle shown in figure l(a), where the three edges are numbered in the counterclockwise direction. In this numbering fashion, the relations of the edges and the corresponding two nodes of each edge are uniquely determined and tabulated as follows edges nodel node2 1 Z2 2 (~D 3 @ ~ Note that the surface area may be written as S" = 2h, where li is the length of the iZh edge and h is the height from this edge to the corresponding vertex. 29

~ 1 a) (a) ~ ~ (b) Figure 1: (a) Triangle edge and node numbering system; (b) Geometrical interpretation of eq(2) and eq(3) for constructing shape functions. z A Z=Zc +AZ I - -- -- - - A 3 I I ' 9 I 7 I I I 5, - - " 0 4 -- - -0 4% 8 - -- Z=Z c 6 i=1,2,3 j=4,5,6 k=7,8,9 Figure 2: Illustration of a typical triangular prism element and the numbering system. 30

It is also noted that the height may be expressed as h = esi z x (r-r.), (1) where r is a position vector right on the ith edge, ri is the position vector at the vertex and ei is the unit vector along the edge. Hence, the triangle area becomes I. S le iz x (r-r). (2) 2 Referring to figure l(b), it is evident that if r is measured within the triangle, then (2) gives the area of another triangle (12'3') with the same length of the base edge shift upward so that r is located on it. In this case, if one writes a vector Si such that S. - - r.(3) S= 2Se x(r-r), (3) then e Si becomes the ratio of this triangle area to the area Se of triangle (123) and this is the design of the edge basis functions for (2-D) triangular elements. The corresponding volumetric basis functions associated with the triangular prism shown in figure 2 may then be written, by inspection, as V. =, S i= 1,2,3 V (zc+ Az-)5 j4,5,6 (4) Vk = Z(k k =7,8,9 where k is the triangle simplex coordinate associated with the vertex (Xk, Yk) to which the kth vertical edge is linked; zc and Az are the z coordinate of the base surface and the height of the prism, respectively. It should be remarked that the basis functions of this format are divergence-less. In order to simplify our formulation, we will rewrite the above basis functions in Cartesian coordinates, yielding, 2Se Az(Z - Az) [-(y - yi)x + (x - x)y] = 1,2,3 2SeAz Vj = 2S (Zc + AZ- Zy ) [-(y- yj)^ + (x - X)y] = 4, 5, (5) Vk = Z [(XklYk2 - Xk2Ykl) - X(Yk2 - Ykl) + Y(Xk2 - Xkl)] k = 7,8,9 31

where (xi, yi) and (xj, yj) are the vertices of the ith and jth edges, respectively; (xk, yk), (Xkl, Ykl) and (xk2, Yk2) are three nodes or vertices in counterclockwise rotation. It is straightforward to derive the curl of (5) and they are given by V x Vi =-2SA [(x -xi) + (y - - (z - z)] = 1, 2, 3 2SeAz VxVj =2Sa [(X-Xx4+(y-yY5+Z +AZ-Z =45 (6) V x Vk = 2S [(Xk2 - xk)X + (yk2 ( yk )y] k = 7,8,9 For FEM implementation, the following FEM elementary entries em = JV X Vm V X VndV (7) Qmn = V VmVndv (8) need to be explicitly expressed on the assumption of isotropic material property. To this end, we follow the notation defined in (5) and (6), where i, i'=1,2,3 represent the top triangle edges, j,j'=1,2,3 denote the bottom triangle edges and k, k'=1,2,3 stand for the vertical three edges. It is found that (7) and (8) can be analytically evaluated and we tabulate the results as follows Pii, = Cii [DiilaZ+ Se(AlZ)3] (9) Pj, = C, [Djj,,,z + S (AZ)3] (10) Pkk' (11) 4Se PI = Pji [D - Se(AZ)3] (12) Pik = Pki= - x *k(SX - Se) + t k(SY-yiSe)] (13) Pj = Pkj = (S)2 [X k (SX - XijSe) +. * k(SY- yjSe)] (14) Qii = (A3Ci i (15) 3 32

Qjj, = (3C Djj (16) Qkk' = AZSeTkk, (17) Qij Qj CijDij (18) 6 Q k - Qk = Qk = Qkj= (19) where Tkk' = 1/6 for k = k'; 1/12 for k # k' C.- (20) = 4(SA)2 (20) DAj = SXX - (x x + xj)SX + xxjSe + SYY - (y + yj)SY + yyjSe The related quantities in the above list of the expressions are defined as Se = J dxdy SXe SX = x dxdy SYSY = y dxdy SXX = x2dxdy SYY = J y2dxdy SXe SXY = xy dxdy Se and these integrals can be analytically evaluated and the results will be given in the appendix. Shape Functions: Type Two Edge-based shape functions for triangular elements can be derived and expressed in various ways and we demonstrate here another type of shape functions. Either of these two sets may be employed for FEM implementation and the solutions from both FEM implementations are shown to be identical. We start from the nodal basis functions for 3-D triangular prism elements. It is well known that the nodal basis functions for the prism shown in figure 2 may be simply expressed in terms of the simplex coordinates of triangles 33

Si C1, I2, and C3 -= 1-I1-2, where (i (i=1,2,3) are equal to the ratio - where Si is the area of the triangle formed by an inner point and the two corners opposite to the ith corner. This definition provides a geometrical explanation as to why (i is equal to 1 at node i and zero at the other two nodes,. The difference of the functions for the top and the bottom triangle is taken into account by introducing a linear dependence in z so that the node shape functions associated with the upper triangle vanishes at the lower triangle. for nodes on top surface for nodes on bottom surface Nl = - (z - zC)C1 N4 = (zc+ Az - )i N2 = (Z - N5 = (zC + Az-z} (21) Az Az N3 = A zc)3 N6 = A(zc + Az -z)3 The corresponding edge basis functions can be constructed using the above simplex coordinates. Specifically, the variations of functions along the z axis remain the same as in (21). However, in the x-y plane, we let S = V - 2V(i, i = 1,2,3 (22) where il, i2 represent node 1 and node 2 of the ith edge, whose direction points from node 1 to node 2. There are several interesting properties associated with this type of functions. First, since the simplex coordinates (mI m = 1, 2, 3 are linear functions with respect to x and y, the gradient of the coordinates become constants. Specifically, V(m = Z x -2Se (23) where 1m is the length of the mth edge opposite to node m and Se is the triangle area, as usual. Also, it can be shown that the relation r E edge i ( 0 r E other edges holds for r at any location of the edges, taking advantage of the features of the simplex coordinates. Clearly, (22) is also divergence-less, since V. S= 0. 34

It is again straightforward to write down the edge-based shape functions for triangular prisms. They are given, by inspection, as VI = - z - z)((iV( - (iV(O' Z = 1,2,3 Vj = z (z + A z - z)(C V(j2 - (J2 V(jl ) j = 4,56 (25) Vk = Z(k k - 7,8,9 It is noted that Vk (k=1,2,3) are the same as those associated with the type one shape functions. The curl of (25) is explicitly given by V X Vi = I {.i1 2 X VCi2 - z2 x Vi,1 + 2(z - zc)Vi x V2} (26) V x V j = / {jz x VXj2 - j2z x V(j, + 2(z - zc)Vj1 x V(j2}(27) V Vk = V(k x z (28) and these expressions can be used to evaluate the matrix elements (7) and (8). The integral relations needed for this purpose are given in the appendix. Below we simply list the analytical expressions for the FEM entries: C.-A + (29) Pii, = C [Dii, + Gii] (29) Pjj, = Cjj + G (30) Pkk' = A tk. tk' (31) Pi = Ci [- + G] (32) Pik = Pki 12S [(1 - 12) k] (33) P3k =Pkj - 12S [(I -2) ] (34) Qii' = B-ii, Di (35) Q = Bjj, Djj, (36) Qkk' A= ZSTkk (37) -Q = 1BD (38) Q 2 U U 35

Qik = Qki = Qk = Qkj = O (39) where Az — B-j = 12Se llj 1 A= 4AzS j Dij = nllT22 + II22T11 - I12T21 - I21T12 G — = f J) (1111122-H12H21) Gj \ /IImn = V(jim V(jn m,n=1,2 and for im n 12 Tm 1 for im # Jn Also, the length of ith edge, for instance, is denoted by li whereas 1k represents the length of the edge opposite to the vertex k. It can be seen that the two sets of shape functions have different geometrical meanings and both lead to closed form matrix entries. Also, the complexity of two formulations is about the same (careful observation shows that type-two is a bit more concise in form but type-one may be a bit easier to program) and the solutions rendered by their corresponding FEM implementations are identical. III. Mesh Generation And Code Interface Triangular prismatic elements are practically attractive for conformal antenna simulations because they simplify mesh generation procedures rather substantially. Basically, once the surface triangular mesh of the printed antenna aperture is completed, the volume mesh can be generated automatically by growing the prismatic elements vertically into the cavity. Inhomogeneities, layers and parasitic elements can be readily incorporated into the triangular prism mesh without difficulty. 36

Typically, upon completion of the mesh, the following parameters/tables must be compiled for use by the finite element solver. They are: * node numbers and coordinates on top surface; * table of element numbers and corresponding nodes on top surface; * table for element and corresponding edge numbers and node numbers; * group of nodes at the outer boundary of top surface; * group of elements in the slot(s)/aperture(s); * additional group(s) of nodes (on surface) for interior PEC's. Other information, such as number of layers along z axis, thickness and material properties (Er, ir) associated with each dielectric layer, as well as frequency scanning range, output information (scattering RCS/radiation pattern or input impedance) and iterative solution tolerance, should be inputted interactively. IV. Preliminary Results To validate our recently developed finite element-boundary integral code, we chose an annular ring slot antenna, whose tetrahedral mesh is available. Both RSC and radiation patterns are tested for the antenna. The configuration used for this testing is shown in figure 3, where the radius and the depth of the cavity are 12.35 cm and 3 cm, respectively; the ring has a radius 7.7 cm and the slot width is 0.75 cm; the dielectric constant is 1.35. In both scattering and radiation cases, we chose the operational frequency to 'be 0.8 GHz. Figure4 shows the comparisons of the RCS pattern with the computed results using the tetrahedral FE-BI code, whereas figure 5 demonstrates both the co-polarization and cross-polarization radiation patterns. They all are in good agreement. The mesh generation procedure for the prismatic element implementation is significantly simplified. Our near future work is to test various spiral antennas using this code and and model the input impedance for different feed structures. Also, the impact of the dielectric substrate/superstrate on the performance and geometry of antennas, such as polarization control and structural miniaturization will be especially focused. 37

t I I probe ' I I a = 12.35 cm b = 0.75 cm po= 7.7 cm 0.7 <f< 1 GHz d = 3 cm.. 135 Figure 3: The geometry of the annular ring slot antenna used for comparisons of the prismatic FE-BI code with the tetra FE-BI code. Appendix Certain integral relations have been employed in the report for analytical evaluation of the elementary FEM entries. These are often found in FEM books but they are usually in local coordinates. Since the expressions in global coordinates are also needed, we tabulate them in this appendix without proof. Assume the three nodes i,j and m of a triangle are in counterclockwise rotation, their global coordinates are given by (Xi, Yi), (Xj, Yj), (Xm, Ym). we have, Se = Js dxdy= 2 1 xj yj \ 1 Xm Ym / SX =/ x dxdy SY e SY =X y dxdy = Se S e (Xi + Xj + Xm) 3 (Yx + Yj + Ym) S3 38

-30 - 0 Prism Elements - Tetrahedral Elements CC -40 _____-________ m -50 - -60 - -70 - -80C -_so o --- o i ----o ---o ----o --- —o-i -o - ' -o - o -oo 0 10 20 30 40 50 60 70 80 90 100 theta (degree) Figure 4: Comparisons of the bistatic (co-pol) RCS patterns computed using the tetra FE-BI code and the prismatic FE-BI code. The normal incident plane wave is polarized along the 0b = 0 plane and the observation cut is perpendicular to that plane.. I I 0 -10 - — 20 m-T -30,-. C -40 n- -50 -60 - -70 80 - -100 -80 -60 -40 -2C ) 0 20 theta (degree) 40 60 80 100 Figure 5: Illustration of x-pol and co-pol radiation patterns from the ring slot antenna shown in figure 3. The solid lines are computed using the tetrahedral FE-BI code whereas the dotted lines are computed using the prismatic FE-BI code. The excitation probe is place at the point (y=0O) in the slot marked in figure 3. 39

SXX = x2 dxdy =2{(X+Xj+Xm)2+(X7+XJ2+X2?)J} SYy = + y2 dxdy =X{(Ym+Y )+Y2+(Y~ +Yj+Y )} sxY = J xydxdy {(X + Xj + Xm)(Yi + Y + Ym) + (XiY + XY + XmYm)} Also, if the simplex coordinates for these three nodes are denoted by (C, Cn (m), then 2 xy (C C dxdy - (p + q + r + 2)! s where p, q, r are all non-negative integer numbers. Since the elements we used for the FEM implementation are linear, the associated integrands involve only linear or quadratic functions of C. Therefore, 2 / i dxdy = S 2 J3 1j Id for i j (j dxdY = j for i j 12 / 40

Performance of an Anisotropic Artificial Absorber for Truncating Finite Element Meshes Jian Gong, John L. Volakis, David M. Kingsland1 and Jin-Fa Lee1 Radiation Lab, Dept. of Electrical Engin. and Comp. Sci.// University of Michigan// Ann Arbor, MI 48109-2122 Abstract A new artificial material absorber for truncating finite element meshes is investigated. The interface of the absorber is made reflectionless by choosing er and tar to be complex diagonal tensors. With some loss, a metal backed thin absorber layer is then sufficient for terminating the mesh. This scheme is simpler to implement than conventional absorbing boundary conditions and offers the potential for higher accuracy. In this paper, we investigate the effectiveness of this anisotropic absorber on the basis of results obtained for problems in propagation (waveguide and printed circuits) and scattering. Based on the results, specific recommendations are given for selecting the absorber loss parameters to optimize its performance for a given thickness. I. Introduction One of the most important aspects of finite difference and finite element implementations is the truncation of the computational volume. An ideal truncation scheme must ensure that outgoing waves are not reflected backwards at the mesh termination surface, i.e. the mesh truncation scheme must simulate a surface which actually does not exist. To date, a variety of non-reflecting or absorbing boundary conditions (ABCs) have been employed for truncating the computational volume at some distance from the radiating or scattering surface, and applications to microwave circuits and devices have also been reported. The ABCs are typically second or higher order boundary conditions and are applied at the mesh termination surface to truncate the computational volume as required by any PDE solution. Among them, a class of ABCs is based on the one-way wave equation method [1, 2] and another is derived starting with the Wilcox Expansion [3, 4]. Also, higher order ABCs using Higdon's [5, 6]formulation and problem specific numerical ABCs have been successfully used, 1EM CAD Lab, ECE Department, Worcester Polytechnic Institute, Worcester, MA 01609 41

particularly for truncating meshes in guided structures [7]. There are several difficulties with traditional ABCs. Among them are accuracy limitations and a need to enforce them at some distance from the perturbance (leading to a large computational volume). Also, although higher order ABCs are considered more accurate, they compromise the convergence of the solution and, furthermore, the resulting system is not as easily amenable to parallelization. An alternative to traditional ABCs is to employ an artificial absorber for mesh truncation. Basically, instead of an ABC, a thin layer of absorbing material is used to truncate the mesh, and the performance for a variety of such absorbers has been considered [8], [9]. Nevertheless, these lossy artificial absorbers (homogeneous or not) still exhibit a non-zero reflection at incidence angles away from normal. Recently, though, Berenger [10] introduced a new approach for modeling an artificial absorber that is reflectionless at its interface for all incidence angles. In two dimensions, his approach requires the splitting of the field components involving normal (to the boundary) derivatives and assigning to each component a different conductivity. In this manner an additional degree of freedom is introduced that is chosen to simulate a reflectionless medium at all incidence angles. Provided the medium is lossy, this property is maintained for a finite thickness layer. Berenger refers to the latter as a perfectly matched layer(PML) and generalizations of his idea to three dimensions have already been considered [11], [12]. Also, implementations of the absorber for truncating finite difference-time domain(FDTD) solutions have so far been found highly successful. Nevertheless, it should be noted that Berenger's PML does not satisfy Maxwell's equations and cannot be easily implemented in finite element (FEM) solution. In this paper we examine a new anisotropic(uniaxial) artificial absorber [13] for truncating FEM meshes. This artificial absorber is also reflectionless at all incidence angles. Basically, by making appropriate choices for the constitutive parameter tensors, the medium impedance can be made independent of frequency, polarization, and wave incidence angle. A PML layer can then be constructed by introducing sufficient loss in the material properties. The implementation of this artificial absorber for truncating finite element meshes is straightforward and, moreover, the absorber is Maxwellian. Below, we begin with a brief presentation of the proposed artificial absorber, and this is followed by an examination of the absorber's performance in terminating guided structures and volume meshes in scattering problems. Results are presented which show the absorber's performance as a function of thickness/frequency and for different loss factors. II. Formulation Consider the waveguide, shielded microstrip and scatterer shown in Figures 1 and 2. We are interested in modeling the wave propagation in these structures using the finite element method. For a general anisotropic medium, the functional to be 42

minimized is j V = XE. VV x E. (-i. V x E) EdV J- E x (r1. V x E) dS, (1) i->n +Sout in which =r and r denote the permeability and permittivity tensors whereas E is the total electric field in the medium. The surface integrals over Si, and S,,t must be evaluated by introducing an independent boundary condition and the ABC serves for this purpose but alternatively an absorbing layer may be used. An approach to evaluate the performance of an absorbing layer for terminating the FE mesh is to extract the reflection coefficient computed in the presence of the absorbing layer used to terminate the computational domain. In this study we consider the performance of a thin uniaxial layer for terminating the FE mesh for certain wave guided structures and scatterers. Such a uniaxial layer was proposed by Sacks et.al. [13] who considered the plane wave reflection from an anisotropic interface (see Figure 3 ). If =r and Er are the relative constitutive parameter tensors of the form a2 0 0 = E = 0 b2 0 (2) 0 0 C2 the TE and TM reflection coefficients at the interface (assuming free space as the background material) are cosei- -A cosOt RTE a2 cosOi + pcosOt V + 2 (3) TM a2cosf t - cosei cos0i + acosOt and by choosing a2 = b2 and c2 = it follows that RTE = RTM = 0 for all incidence angles, implying a perfectly matched material interface. If we set a2 = a - j/3, the reflected field for such a metal-backed uniaxial layer is IR(0,)l =e-2ptcs (4) where t is the thickness of the layer and Q0 is the plane wave incidence angle. The parameter ak is simply the wavenumber in the absorber. Basically, the proposed metal backed uniaxial layer has a reflectivity of -30 dB if /tcosOi = 0.275A or -55dB if 3tcosOi = 0.5A, where A is the wavelength of the background material. The reflection coefficient (4) can be reduced further by backing the layer with an ABC rather than a PEC. However, the PEC backing is more attractive because it eliminates altogether the integrals over the surfaces. Clearly, although the interface is reflectionless, the finite thickness layer is not, and this is also true for Berenger's PML absorber. 43

Below we present a number of results which show the performance of the proposed uniaxial absorbing layer as a function of the parameter /, the layer thickness t and frequency for the guided and scattering structures shown in Figures 1 and 2. We remark that for the microstrip line example it is necessary to set a2 = rb(( - j/3) for the permittivity tensor and a2 = rb(a - jp) for the permeability tensor, where orb and u1rb are the relative constitutive parameters of the background material (i.e. the substrate). Results Rectangular Waveguide Let us first consider the rectangular waveguide shown in Figure 1. The guide's crosssection has dimensions 4.755 cm x 2.215 cm and is chosen to propagate only the TE10 mode. It is excited by an electric probe at the left, and Figure 4 shows the mode field strength inside the waveguide which has been terminated by a perfectly matched uniaxial layer. As expected, the field decay inside the absorber is exponential and for 3 values less than unity the wave does not have sufficient decay to suppress reflections from the metal backing of this 5cm layer. Consequently, a VSWR of about 1.1 is observed for / = 0.5. However, as is increased to unity, the VSWR is nearly 1.0 and the wave decay is precisely given by e-~kdcosi, where d is the wave travel distance measured from the absorber interface and here Oi = 44.5~. Nevertheless, when /3 is increased to larger values, the rapid decay is seen to cause unacceptable VSWR's. One is therefore prompted to look for an optimum decay factor /3 for a given absorber thickness and Figure 5 provides a plot of the TE1o mode reflection coefficient as a function of /. Figure 5 is typical of the absorber performance and demonstrates its broadband nature and the existence of an optimum value of f3 for minimizing the reflection coefficient. Basically, the results suggest that P must be chosen to provide the slowest decay without causing reflections from the absorber backing. That is, the lowest reflection coefficient is achieved when the entire absorber width is used to reduce the wave amplitude before it reaches the absorber's backing. As expected, this optimum value of /3 changes with frequency but the broadband properties of the absorber are still maintained since acceptable low reflections can still be achieved for unoptimized / values. For example, in the case of f = 4.5Ghz the optimum value of 3 = 1 gives a reflection coefficient of -45dB whereas the value of / = 3 gives a reflection coefficient of -37 dB which is still acceptable. It should be noted though that setting 3 = 3 allows use of an absorber which is about 2cm or 0.3 free space wavelengths. Not surprisingly (see (4)), for this example, the value of a does not play an important role in the performance of the absorber and this is demonstrated in Figure 6. As seen, setting a = /3 gives the same performance as the case of a = 1 shown in Figure 5. Our tests also show that other choices of a give the same absorber performance. However, it is expected that a will play a role in the presence of attenuating modes and it is therefore recommended to choose a = 3 to ensure that all modes are 44

absorbed. Microstrip line The performance of the perfectly matched uniaxial layer in absorbing the shielded microstrip line mode is illustrated in Figure 7 where the reflection coefficient is plotted as a function of /3. In this case, the microstrip line is extended into the 5cm absorbing layer up to the metallic wall. Similarly to the waveguide, we again observe that an optimum h3 value exists and it was verified that (in the absorber) the quasi-TEM mode has the same attenuation behavior as shown in Figure 4. The reflection coefficient at the optimum 3 = 1 is now -28.5dB and if better performance is required, a thicker absorbing layer may be required. Again, as in the case of the waveguide example, the value of a is of small importance to the performance of the absorber and this is demonstrated in Figure 8. Band Eliminator Another guided structure that was considered is the stripline band eliminator shown in Figure 9. The band eliminator consists of a circular disk 1.64 cm in radius, placed on a substrate having a relative dielectric constant of 2.4. The disk is directly coupled to 50 Q striplines and to take advantage of symmetry, only half of the problem's domain is used in the FEM computation. The thickness of the absorber used in this calculation is 1.5cm, and in Figure 10 we present the insertion loss results from 2GHz to 4GHz. The values of a and d for the absorber were both set to 3 and, as seen, our calculations compare very well with measured data and reference calculations given in [14], thus, demonstrating the accuracy of the absorber termination. Sphere Scattering As another example we examined the scattering by a metallic sphere illuminated by a plane wave (see Figure 2). The unknown quantity in the FEM implementation was the electric field, and the excitation due to the incident plane wave was incorporated by enforcing the boundary condition Et = -E on the surface of the sphere. The FEM mesh was terminated by a cubical O.lAo thick metal backed absorber whose interface was placed only 0.05Ao away from the sphere's surface. Based on the attenuation curves shown in Figure 4, a value of / = 4 or slightly greater is needed to suppress reflections from the absorber's backing. Indeed, lower values of / did not perform well for this application. The bistatic scattering pattern computed with a mesh density of 14 edges per wavelength is shown in Figure 11 and is within IdB of the exact result. Better agreement would clearly require thicker absorber and a lower /3 value. Again, the value of a does not have appreciable effect on the absorber's performance. About 56,000 unknowns were used for this simulation and the system converged in about 2000 iterations. We also remark that higher sampling densities did not appreciably improve the accuracy of the result or the relative convergence rate. Conclusions A new method of mesh truncation based on an artificial absorbing layer with anisotropic material properties was implemented for terminating the FEM mesh. Explicit recom 45

mendations were given for choosing the absorbers thickness and attenuation properties to optimize performance. It was demonstrated that good absorption can be achieved with relatively thin absorbing layers (a small fraction of a wavelength) placed close to the computational domain. Also, because there is no surface integral or ABC at the mesh truncation boundary, the method is easily implemented in FEM solutions and is well suited for parallel computations. References [1] B. Engquist and A. Majda, "Absorbing boundary conditions for the numerical simulation of waves," Math. Comput. Vol. 31, pp. 629-651, 1977 [2] L.Halpern and L.N. Trefethen, "Wide-angle one-way wave equations," J. Acoust. Soc. Amer. Vol. 84, pp. 1397-1404, 1988 [3] J.P. Webb and V.N. Kanellopoulos, " Absorbing boundary conditions for the finite element solution of the vector wave equation," Microwave Opt. Tech. Lett. Vol. 2, pp. 370-372, 1989 [4] A. Chatterjee and J.L.Volakis, "Conformal absorbing boundary conditions for the vector wave equation," Microwave Opt. Tech. Lett. Vol. 6, pp. 886-889,1993 [5] R.L. Higdon, "Absorbing boundary conditions for acoustic and elastic waves in stratified media," J. Comp. Phys. Vol. 101, pp. 386-418, 1992 [6] J.Fang, "ABCs applied to model wave propagation in Microwave integratedcircuits," IEEE Trans. MTT, Vol. 42, No. 8, pp. 1506-1513, Aug. 1994 [7] J.-S. Wang and R. Mittra, "Finite element analysis of MMIC structures and electronic packages using absorbing boundary conditions," IEEE Trans. Microwave Theory and Techn., Vol. 42, pp. 441-449, March 1994. [8] T. Ozdemir and J.L.Volakis, "A comparative study of an absorbing boundary condition and an artificial absorber for truncating finite element meshes," Radio Science Vol. 29, No. 5, pp. 1255-1263, Sept. 1994 [9] C. Rappaport and L. Bahrmasel, "An absorbing boundary condition based on anechoic absorber for EM scattering computation," J. Electromagn. Waves Appl., Vol. 6, No. 12, pp. 1621-1634, Dec. 1992. [10] J.P. Berenger "A Perfectly Matched Layer for the Absorption of Electromagnetic Waves," J. Comp. Physics, Vol. 114, pp. 185-200, 1994. [11] D.S. Katz, E.T. Thiele, and A. Taflove "Validation and Extension to Three Dimensions of the Berenger PML Absorbing Boundary Condition for FD-TD Meshes," IEEE Microwave and Guided Wave Letters, pp.268-270, August 1994 46

[12] W.C. Chew and W.H. Weedon "A 3-D Perfectly Matched Medium from Modified Maxwell's Equations with Stretched Coordinates," Microwave and Optical Technology Letters, pp. 599-604. Sept. 1994. [13] Z.S. Sacks, D.M. Kingsland, R. Lee and J.F. Lee, "A perfectly matched anisotropic absorber for use as an absorbing boundary condition," submitted for publication, 1994 [14] R. R. Bonetti and P. Tissi, "Analysis of Planar Disk Networks", IEEE MTT-26, July 1978, pp. 471-477. 47

List of Figures 1 A rectangular waveguide (a) and a microstrip line (b) truncated using the perfectly matched uniaxial absorbing layer.................... 9 2 Geometry of sphere scattering example........................ 9 3 Plane wave incidence on an interface between two diagonally anisotropic half-spaces................................ 10 4 Field values of the TE10 mode inside a waveguide terminated by a perfectly matched uniaxial layer. The absorber is 10 elements thick from the 71st to the 80th element and each element was 0.5 cm which translates to about 13 samples per wavelength at 4.5 GHz. 10 5 Reflection coefficient vs /3 (a = 1) for the perfectly matched uniaxial layer used to terminate the waveguide shown in Figure l(a)....... 11 6 Reflection coefficient vs /3, with a = /3, for the perfectly matched uniaxial layer used to terminate the waveguide shown in Figure l(a). 11 7 Reflection coefficient vs /3 with a=1, for the shielded microstrip line terminated by the perfectly matched uniaxial layer.............. 12 8 Reflection coefficient vs /3 with a = f3, for the shielded microstrip line terminated by the perfectly matched uniaxial layer............... 12 9 Geometry of a band eliminator............................ 13 10 Insertion loss for the band eliminator shown in Figure 9............ 14 11 Bistatic scattering for a 0.4A diameter sphere for different values of /3. The FEM mesh density was approximately 14 elements per A and was terminated by a cubical thin absorber whose inner and outer diameter was 0.5A and 0.7A, respectively.................................. 15 48

Electric Probe \= = d (a). waveguide s=-= y z= 1-j d,=4y= = 1 -jp d+t=40cm t=2cm Stripline Absorbing layer Stripline er' r ' h — IH ~r r' _L W t::: Er L 3H - b -r=32 LL~ Er=3*2 d + t = 12.0 cm L = 2.38 cm h = 0.21 cm H = 1.06 cm w = 0.47 cm (b). Stripline Figure 1: A rectangular waveguide (a) and a microstrip line perfectly matched uniaxial absorbing layer. (b) truncated using the. Anisotropic Absorber jI 0.1 X -.,! Figure 2: Geometry of sphere scattering example. 49

Region 1 Reflected wave Oir ----— Incident wave Incident wave Region 2 Transmitted wave - - - - -O - - ------ a20 0 0 b20 0 0C2 V z Figure 3: Plane wave incidence on an interface between two diagonally anisotropic half-spaces. 40...., 35 - beta=0.5 0 + 3 + 3. 25-.. 5.. 20 15- - 01.. 60 62 64 66 68 70 72 74 76 78 80 Segments number along waveguide (alpha=1.0, f=4.5GHz) Figure 4: Field values of the TE1o mode inside a waveguide terminated by a perfectly matched uniaxial layer. The absorber is 10 elements thick from the 71St to the 80th element and each element was 0.5 cm which translates to about 13 samples per wavelength at 4.5 GHz. 50

'~-20 - ^^ - -. -o 3-25 - -30 - " f ~30 ^ ",,,, -35 - 1 -' ^ -40 -. / / -45 --50 -------------— ' --- —---- 0 1 2 3 4 5 6 7 8 9 10 Lossy parameter beta (alpha=1.0) Figure 5: Reflection coefficient vs 3 (a = 1) for the perfectly matched uniaxial layer used to terminate the waveguide shown in Figure l(a). 0 1 2 3 4 5 6 7 8 9 parameter beta=alpha Figure 6: Reflection coefficient vs /3, with a = /3, for the perfectly layer used to terminate the waveguide shown in Figure l(a). 10 matched uniaxial 51

( eI -5 -10 -15 a., ' -20 0 o -25 m - 30 n'" - f=4.0GHz f=4.5GHz f=5.0GHz 0 2 3 4 I 6 7I I -35 - -40 - -45 i-lll I.... 0 1 2 3 4 5 6 7 Lossy parameter beta (alpha=1.0) Figure 7: Reflection coefficient vs P with a=1, for the nated by the perfectly matched uniaxial layer. 8 9 1 D shielded microstrip line termi-.. v -5 -10 -15 - -- f=4.0GHz f=4.5GHz f=5.0GHz -a._0. 1 -20 o -25 -30 a: -35 - -40 - -45 F -Fr50 1 0 1 2 3 4 5 6 7 parameter beta=alpha Figure 8: Reflection coefficient vs 73 with a = /, for terminated by the perfectly matched uniaxial layer. 8 9 10 the shielded microstrip line 52

6.0 cm.4 __ _ _ _ __ _ _ _ __ _ _ _ __ _ _ _ __ _ _ 1.5 cm,I output port Top View 5mm 6.0 cm input port 6.0 cm.4 — O Cross Section View 6.4mm Figure 9: Geometry of a band eliminator. 53

0 E...... r^c -,[: CM c /cv. 0 CD (Ep) ssoa) uoiJSU| U) [] c6 6 6 d (gP) sso-] UO!1JeSUI Figure 10: Insertion loss for the band eliminator shown in Figure 9. 54

Bistatic RCS (14 elements per wave length) 0,00 -5,00 E x -10,00 -15,00 -20,00 V. I I Mie series solution ----— Beta=4.2, Alpha=,.0 ---- Beta=4.2, Alpha=4.2 0 50 100 150 200 Azimuthal angle Figure 11: Bistatic scattering for a 0.4A diameter sphere for different values of d. The FEM mesh density was approximately 14 elements per A and was terminated by a cubical thin absorber whose inner and outer diameter was 0.5A and 0.7A, respectively. 55

Optimal Selection of a Uniaxial Artificial Absorber Layer for Truncating Finite Element Meshes Jian Gong and John L. Volakis Radiation Lab, Dept. of Electrical Engin. and Comp. Sci. University of Michigan Ann Arbor, MI 48109-2212 Abstract We examine the performance and optimal selection of a new uniaxial artificial absorber layer for truncating finite element meshes. The absorber is reflectionless at its interface for all incidence angles but in practice a finite layer must be used, and as a result the absorber is no longer reflectionless. In this Letter we examine the absorber's optimal selection of its thickness and dielectric properties in the context of a finite element implementation for a shielded microstrip line. Introduction A most important aspect of a finite element implementation is the termination of the computational domain. This is carried out by either using an absorbing boundary condition (ABC) [1, 2, 4, 3] or a thin layer absorber [5, 6]. However, the ABCs and thin layer absorbers considered so far have non-zero reflections for waves impinging at oblique incidences and this limits their accuracy and leads to requirements for larger computational domains. Recently, Sack et al. [7] proposed a new anisotropic (uniaxial) absorber -for truncating finite element meshes. This absorber is indeed reflectionless (i.e. perfectly matched at its interface) for all incidence angles. However, in practice a finite thickness absorber must be employed, as shown in figure 1. As a result, although the absorber has a perfectly matched interface (PMI), the layer is no longer reflectionless due to reflection from the back of the layer. It is therefore of interest to optimize the layer's absorptivity in terms of its thickness and material properties. Of interest is to also examine the efficiency of the absorber in actual numerical simulations, since to date no application of the absorber's performance and optimization has yet been published. In the following, we 56

examine the absorber's performance and optimization when used to terminate the finite element mesh of a shielded microstrip line. Uniaxial Absorbing Layer Consider the uniaxial absorbing layer shown in Figure 1 whose relative permittivity and permeability are given by a 0 0 a 0 0 br = trb0 b 0 Cr = rb 0 b 0 (1) 0 0 c 0 0 c where (Erb, rb) refer to the relative permittivity and permeability of the background medium. The corresponding TE and TM reflection coefficients at the interface are given by TE - cosOi - -CosOt COSO + acos~t (2) cosOi + -cosOt fTM b-cosot - cosOi cosOi + a/cosOt and by choosing a = b and c = I/a, it follows that RTE = RTM = 0 for all incidence angles. This implies a perfectly matched interface (PMI) and the transmitted wave decays as e-wkdcosi, where d is the distance from the interface, k = w /,ocovIrcr is the wavenumber of the background material and we have set a = a - jf. Clearly, unless the transmitted wave is completely absorbed while traveling in the metal-backed absorber, the layer will have a nonzero reflected field. It is therefore of interest to optimize the absorbing layer performance through an appropriate selection of its parameters. Application to Shielded Microstrip Line Consider the microstrip line shown in Fig. 2. The line is excited by a vertical probe at its left and is terminated to the right by a thin PMI layer of the uniaxial absorber. The length of the line is d+t=12crm and a finite element implementation with the brick element [8] was used for modeling the shielded structure. We also note that for this implementation the microstrip line was extended into the 5cm absorbing layer. Figure 2 shows the vertical field magnitude below the microstrip line for different 3 values and a = 1. Only the fields across and near the absorbing layer are depicted since our interest is to examine the field behavior in that region. As expected the field decay inside the absorber is exponential and for 3 less than unity the wave does not have sufficient decay to suppress reflections from the metal backing of the 5cm PMI layer. As 3 increases to unity, the VSWR decreases 57

to about 1.05 but for / values above 3 we again observe deterioration of the VSWR. One is therefore prompted to seek an optimum decay factor P for a given absorber thickness. Our study led to the conclusion that P should be selected so that the entire absorber thickness is utilized in attenuating the field to less than 5 percent of its original value. For the chosen 5cm PMIJ layer, a value of / = 1 or so is appropriate at a frequency of 4.5GHz for a 30 dB reduction. Lower reflectivities can be achieved by increasing the layer thickness and changing P appropriately. It turns out that a plays very little role in the performance of the absorber since the curves in Figure 2 remain nearly unchanged as a is increased above unity. For large / values the decay becomes very rapid and the PMI interface begins to resemble a magnetic wall. This behavior is due to the numerical FEM implementation and stems from the inaccurate modeling of the rapidly decaying fields in the layer. It is basically more efficient to use a thicker absorber layer with a lower 3 value and use nominal sampling (10 to 15 per wavelength) rather than a thinner absorber with high 3 values and consequently excessive sampling requirements. Our study shows that the absorber fields with 1 < / < 2 can be accurately modeled using cell sizes on the order of 0.05 guided wavelengths. The accuracy of the IPMI layer for circuit parameter computations is illustrated in Figure 3. It is seen that use of the optimized 5cm PMI layer, with a = 1 and / = 1, yields very accurate input impedance values. The shown microstrip line impedances were computed by measuring the vertical field at the probe's location without a need to extract the VSWR which is often difficult with unstructured finite element meshes. A recent measurement technique has also permitted the validation of the computed tangential fields (to the microstrip) above the dielectric. As seen from Figure 4, the calculated and measured fields are nearly overlaid and this is of importance since inexact mesh truncation schemes are prone to inaccuracies in the near zone. Note that the shielded microstrip line dimensions for the data in Fig. 4 are different from those corresponding to Fig. 3. Conclusion We examined the performance and optimal use of a new uniaxial absorber layer with a perfectly matched interface(PMI) for terminating finite element meshes. It was found that the loss parameter of the PMI layer must be selected to have a nominal value and make full use of the entire layer thickness for field decay. Moreover, it was demonstrated that such a PMI layer is extremely accurate for truncating finite element meshes and rivals the Berenger's PML [10] layer which was recently introduced for truncating finite difference -time domain meshes. It is therefore expected that the PMI layer will find wide-spread use for truncating finite element meshes. 58

References [1] L.Halpern and L.N. Trefethen, "Wide-angle one-way wave equations," J. Acoust. Soc. Amer. Vol. 84, pp. 1397-1404, 1988 [2] R.L. Higdon, "Absorbing boundary conditions for acoustic and elastic waves in stratified media," J. Comp. Phys. Vol. 101, pp. 386-418, 1992 [3] J.Fang, "ABCs applied to model wave propagation in Microwave integrated- circuits," IEEE Trans. MTT, Vol. 42, No. 8, pp. 1506-1513, Aug. 1994 [4] A. Chatterjee and J.L.Volakis, "Conformal absorbing boundary conditions for the vector wave equation," Microwave Opt. Tech. Lett. Vol. 6, pp. 886-889,1993 [5] T. Ozdemir and J. L.Volakis, "A comparative study of an absorbing boundary condition and an artificial absorber for truncating finite element meshes," Radio Science Vol. 29, No. 5, pp. 1255-1263, Sept. 1994 [6] C. Rappaport and L. Bahrmasel, "An absorbing boundary condition based on anechoic absorber for EM scattering computation," J. Electromagn. Waves Appl., Vol. 6, No. 12, pp. 1621-1634, Dec. 1992. [7] Z.S. Sacks, D.M. IKingsland, R. Lee and J.F. Lee, "A perfectly matched anisotropic absorber for use as an absorbing boundary condition," submitted for publication [8] J.L. Volakis, J.Gong and A. Alexanian, "Electromagnetic Scattering From Microstrip Patch Antennas and Spirals Residing in a Cavity," Electromagnetics, Vol. 14, No. 1, pp 63-85, 1994 [9] B.C. Wadell, "Tramsmission Line Design Handbook," Artech House, Boston, 1991 [10] J.P. Berenger "A Perfectly Matched Layer for the Absorption of Electromagnetic Waves," J. Comp. Physics, Vol. 114, pp. 185-200, 1994. 59

t rb rrb Prz E rb rb rz Figure 1: Illustration of wave incidence upon a perfectly match interface (PMI) with and without metal backing. 62 64 66 68 70 72 74 76 Segments number along waveguide (alpha=1.0, f=4.5GHz) 60

70....... 65 60 55 Theory --- FEM \ \Absoring layer \ \ Microstrip line lar t = 12 0 \\: l., -' "....=hi 2 -\ d1t=1Oc, m hi, \ w L / s F - h: 021 m, H L=2.3oan w:O 0.47 cm 50 CAl E O 45 N 40 35 30 25 20 ' 1 I I I I I' 0.2 0.4 0.6 0.8 1 width of microstrip (cm) 1.2 1.4 1.6 Figure 2: 0.1 *E 2 0.01 0.001 I I I I I I I I I I I I I I — Measured Tangential Field --- -3 - FDTD Calculated -- - FEM\ Calculated I I I II I, I I I I I I I I 10 100 Position (microns) 1000 Figure 3: Comparison of the near field as a function of locations above the microstrip 61

A hybridization of finite element and high frequency methods for pattern prediction of antennas on aircraft structures T. Ozdenmir, M. W. Nurnberger, and J. L. Volakis Radiation Laboratory Department of Electrical Engineering and Computer Science University of Michigan Ann Arbor, Michigan 48109-2122 volakis@um. cc.umich.edu 1 Abstract This paper considers the hybridization of the finite element and high frequency methods for predicting the radiation pattern of printed antennas mounted on aircraft platforms. The finite element method is used to model the cavity-backed antennas whereas the interactions between the radiators and the substructures are treated via a high frequency technique such as the GTD, PO/PTD or SBR. We present comparisons between measurements and calculations along with a qualitative description of the employed finite element and high frequency codes. 62

2 Introduction Over the past few years the finite element method (FEM) has been applied to various scattering [1]-[3] and antenna [4]-[6] analysis problems with remarkable success. One of the main advantages of the Finite Element Method (FEM) is its inherent geometrical adaptability (see Figure 1) and ease in handling materials. Consequently rather complex configurations can be readily handled without a need for code rewriting or modification. Typically, the FEM solver is interfaced with sophisticated mesh generators [7]-[8], which can provide the necessary volume mesh using a variety of tessellation elements such as bricks, tetrahedrals or prisms (see Figure 2). This permits an accurate geometrical modeling of the structure which is essential for antenna related problems since fine geometrical features play an important role in the characterization of the antenna parameters. However, as is often the case with robust numerical simulations, the combined modeling of large and small geometrical features within a single mesh or computational domain leads to inefficiencies and unacceptable CPU requirements. An example of this situation is the problem of antenna radiation in the presence of a large platform such as an aircraft, a satellite or other airborne and land vehicles. In this case, the antenna requires accurate modeling of its fine geometrical features and material parameters (see Figure 3) using, for example, the FEM. However, modeling of the large aircraft platform using the FEM requires several million degrees of freedom, thus, making it unattractive for such an implementation. Instead, it is more appropriate to model the vehicle substructure using a high frequency method such as the Geometrical Theory of Diffraction(GTD) [9] or the Physical Optics(PO)/Physical Theory of Diffraction (PTD) [10] and Shooting and Bouncing Ray(SBR) method [11]. To overcome the inefficiencies of a single formulation, in this paper we present a simple hybridization of the finite element and high frequency methods. Basically, the finite element method is employed for the analysis of the antenna in the absence of substructures such as wings, fins, engines, or surface details. The pattern or antenna aperture currents generated by the FEM code are then used as the source(s) of the high frequency code. Two high frequency codes 63

were investigated in examining the effectiveness of this hybridization. One of the codes is based on the GTD and the other high frequency code combines the PO/PTD and SBR to compute the interactions of the antenna pattern with the surrounding substructures. Basically, the proposed hybridization takes advantage of the robust FEM formulation to obtain an accurate characterization of the antenna in isolation whereas the high frequency technique is used for modeling the interactions with the large features of the substructure(see Figure 4). The effectiveness of this hybridization is evaluated by comparing the results of the analysis with measurements corresponding to a patch antenna radiating on a finite cylinder in the presence of a flat plate attached to the cylinder as illustrated in Figure 5. In the following sections we first give a qualitative description of the finite element and high frequency codes, and then proceed with a more detailed description of the hybridization procedure followed by a comparison of measurements with calculations. 3 Finite Element Code Description Several finite element codes have been developed at the University of Michigan for the analysis and design of printed antenna configurations. Typically, the printed antenna configuration is assumed to be recessed in some metallic or coated platform and the codes differ in the element used for the tessellation of the antenna, the type of platform assumed in the analysis (planar, cylindrical or doubly curved) and the closure condition employed for terminating the finite element mesh as illustrated in Figure 6. The following codes are available for antenna radiation and scattering analysis. FEMA-CYL: This code is specialized to cavity-backed antennas recessed in a metallic cylindrical platform(see Figure 6). The FEM is used to model the cavity region and the boundary integral or absorbing boundary conditions (ABC) are employed for truncating the mesh. In order to simplify the user interface, cylindrical shell elements are used for the discretization of the cavity and the aperture (see Figure 2). The result is a structured mesh and the user interface is reduced to specifying the dimensions of the antenna. Multiple 64

patches and patch arrays (consisting of individual cavity-backed elements or elements on a single substrate) can be considered but their geometry is restricted by the generated grid. This situation is similar to finite difference codes, where the model's boundary is modified to fit the geometry. Consequently, although this code is user-friendly, it is not appropriate for computing the radiation parameters of circular patches or spirals. However, scattering parameters are not as sensitive to minor geometrical modifications and therefore FEMACYL can still be used for computing scattering from non-rectangular printed antenna configurations and cavities. The theoretical background of FEMA-CYL is described in [12] and [13], and the code has been validated for antenna and scattering applications. An example calculation using FEMA-CYL is illustrated in Figure 7. The measured data given in Figure 7 were obtained at NASA-Langley [14] with the patch antenna placed on the cylindrical section of the ogive+cylinder structure as illustrated in the Figure. For this calculation, the cylinder Green's function of the second kind was employed for truncating the finite element mesh to generate a combined finite element-boundary integral system. The latter is solved using the biconjugate gradient [15] or QMR [16],[17] iterative solver and the FFT is used to speed-up the matrix-vector product evaluation in the solver [18]. Consequently, in spite of the fact that the boundary integral matrix is fully populated, the code's computational and memory requirements remain at O(N). The system solution generates the fields within and on the surface of the antenna cavity. For the radiation and scattering pattern calculations only the surface electric fields are used to generate equivalent magnetic currents. These are subsequently substituted in the radiation integral with the appropriate platform Green's function to generate the antenna parameters. For antenna excitation, a probe or a coaxial cable feed model is employed and in the case of scattering the excitation becomes the aperture fields due to the incoming field. When a probe feed model is used, the input impedance is the line integral of the electric field along the length of the probe divided by the magnitude of the probe current. When a coaxial feed is employed, the excitation is introduced by setting the fields 65

coinciding with the edge elements bordering the opening of the coax cable equal to AV/AL, where AV is the given potential between the outer and inner surface of the cable and AL = b - a, where b and a correspond to the outer and inner radius of the cable, respectively. The center conductor is modeled by setting to zero the field unknowns that are associated with the edges coinciding with the conductor. The improved results of this feed modeling approach are shown in Figure 8 and the details of the modeling are described in [19]. FEMA-CYL can also be used to include the effect of superstrate materials (see Figure 6) which may extend over the antenna platform. To avoid use of complicated Green's functions, in this case the mesh is extended a fraction of a wavelength over the cavity's aperture and an absorbing boundary condition (ABC) is employed for truncating the FEM mesh [10]. With this type of mesh truncation, the entire system is sparse and thus the memory and CPU requirement are again 0(N) without a need to make use of the FFT. We have found that placement of the mesh at a distance of 0.3 wavelengths away from the coatings surface is sufficient to obtain the radiation pattern with reasonable accuracy as illustrated in Figure 7(b). The pattern is computed by introducing equivalent electric and magnetic currents placed at the surface of the dielectric coating and then propagating them using the free-space Green's function. FEMA-TETRA: This code employs the same FEM formulation as FEMA-CYL for computing the parameters of printed antennas recessed in a metallic platform. In contrast to FEMA-CYL, this code employs tetrahedral elements for modeling the radiating structure and bricks or shells for modeling the feed (see Figure 9). As a result, it incorporates maximum geometrical adaptability and has already been employed for the analysis of antenna+feed configurations beyond the capabilities of moment method codes. The FFT can also be invoked for carrying out the matrix-vector products associated with the boundary integral portion of the system in the iterative solver, thus maintaining an 0(N) memory requirement for the entire system as described in [4]. Its accuracy has already been 66

demonstrated for a variety of antennas, including stacked patches and finite arrays in ground planes, printed and slot spirals and aperture antennas. In the case of aperture fed patch configurations such as that shown in Figure 10, the mesh generation is divided into the regions above and below the feeding slot. The mesh in each region can employ different elements and this is important in avoiding meshing bottlenecks due to the narrowness of the slot. Two meshing regions are then mathematically connected by enforcing electric field continuity across the slot. However, since the mesh elements in each regions are different, this condition is implemented by first relating the fields in the slot to the potential existing across the slot. Using this procedure, there is no need to have coincidence of the edges bordering the slot and thus different elements can be used in the regions above and below the slot aperture. The FEMA-TETRA code is interfaced with commercial solid- modeling packages such as SDRC I-DEAS and PATRAN [8] which are used to enter the geometry and generate the volume mesh. Using such geometry modeling packages and owed to the adaptability of tetrahedrals, nearly any printed antenna configuration can be modeled and analyzed using FEMA-TETRA. As can be expected, because FEMA-TETRA is not specialized to any antenna configuration all geometrical data must be entered through the solid- modeling package and this can be a time consuming task for inexperienced users. FEMA-PRISM: To avoid the time consuming task of volume mesh generation without compromising geometrical adaptability, this code provides an attractive compromise by making use of prisms (see Figure 2). Prisms have triangular faces at their top and bottom faces and can thus be used to model any patch configuration with sufficient geometrical fidelity. Typically, since the substrate/superstrate is of constant thickness, once the surface grid is constructed, the prismatic volume mesh can be generated by growing the elements above and below the surface of the printed antenna as illustrated in Figure 11. Thus using prisms reduces the mesh generation process from a three-dimensional to a two-dimensional task. That is, the user needs to only generate the surface mesh in the aperture containing the printed antenna and the volume mesh 67

is then built automatically by specifying the substrate/superstrate thickness and information relating to the closure condition. FEMAPRISM makes use of an ABC or metal-backed absorber for truncating the mesh at some distance from the aperture surface. The accuracy of this approach is illustrated in Figure 12 where the input impedance of a circular patch computed using a metal-backed absorber closure is shown to be in good agreement with results based on the boundary integral closure. As in the other codes, the radiation and scattering patterns are computed by propagating the equivalent magnetic currents (existing at the cavity aperture) using the free-space Green's function. Also, since a non-integral closure condition is used for truncating the mesh, no restriction is placed upon the curvature of the platform. Results generated by using the FEM codes refered to above will not, of course, include effects due to interactions between the antenna and the substructures residing on the same platform. To include these, the FEM codes are interfaced with the high frequency analysis packages described next. 4 High Frequency Code Descriptions GTD/UTD Code: The hybridization of GTD and low frequency codes such as the moment method were considered in the early 1980s [20],[21] for scattering and antenna radiation analysis. Typically, the GTD was employed to account for diffraction contributions from edges and corners near the radiating elements and a review of some of this work was given more recently in [22]. In the case of scattering, hybrid GTD and moment method codes have been used to model small details on larger structures. Specific applications include uses of the moment method to model a small wire or a crack in the presence of a large platform which is associated with several scattering centers caused by geometrical optics reflections and diffractions from edges and corners. For the most part these hybrid codes have been restricted to specialized implementations where the scattering substructure is hard-wired into the code and is characterized by a few scattering centers. Also, for radiation analysis available GTD codes have been restricted to include a small class of antenna elements such as narrow slots, monopoles and infinitesimal current elements. 68

A well developed UTD code, referred to as NEC-BSC was developed by Marhefka and Burnside [23] at the Ohio State University and was used in conjunction with the FEM codes for printed antenna pattern calculations on simple aircraft platforms. The UTD code computes the interactions between current elements and the aircraft substructure using ray diffraction methods. In the employed UTD code the fuselage is modeled as a truncated cylinder or an ellipsoid and the attached wing and fins are modeled as combination of flat plates. The overall radiation pattern of the specified current element is then computed by adding the direct, geometrical optics, and diffraction contribution from all possible primary and secondary sources which arrive at the observer. Shadowing is also performed as part of the pattern calculation. In our implementation, the UTD code is supplied with the location and amplitude of the magnetic surface current elements as computed by the FEM code in the absence of the fuselage appendages. Consequently, these currents include all coupling effects internal to the antenna cavity and among the antenna elements. However, they do not include coupling effects due to substructure contributions which return back to the antenna aperture. One way to include these effects is by re-executing the FEM code with the secondary external diffraction contribution as part of the excitation in addition to the feed. Alternatively, reciprocity can be used to first compute the secondary fields arriving at the antenna aperture due to plane wave incidence before executing the FEM code. APATCH Code: This code was developed at DEMACO [24] and models the aircraft structure using facets. It is an outgrowth of the XPATCH code [25] and since the facets can be as small as 1/10 of a wavelength, it permits a good geometrical fidelity of realistic aircraft. The field calculations are performed using PO/PTD for the specular and first order diffraction contributions. An advantage of using the PO is the inclusion of contributions from coated surfaces using the plane wave reflection coefficients of the multilayer coating. Another advantage of the code is the inclusion of multiple specular interactions using the SBR method [26]. Shadowing is typically a 69

major task when hundreds of thousands of facet elements are used for modeling the structure, but can be rapidly executed by making use of the Z-buffer of the SGI workstation. The APATCH code accepts as input either the surface magnetic current elements computed by the FEM code over the entire nonmetallic aperture of the printed antenna or the corresponding volumetric radiation pattern. If the source information is given by the magnetic currents, the code generates ray bundles which are then interacted with the aircraft structure using the principles of ray and physical optics. Each element is treated individually and it is therefore appropriate to group nearby elements to reduce CPU time. The final pattern is computed by a linear addition of the contributions from all current elements or groups of elements. If the APATCH source is the overall antenna pattern in the absence of substructure interactions, the code introduces weighted ray fields which are subsequently interacted with the aircraft structure using ray and physical optics principles. There is enough flexibility into the code to compute the overall pattern by adding the primary and secondary patterns each provided by different procedures. For example, the final pattern may be computed as the linear addition of the primary pattern given by the FEM code and the secondary pattern due to the interference of the substructure as provided by APATCH. The APATCH code includes a well developed graphical interface for rendering the entire aircraft geometry using the open GL library of the SGI workstation. Ray interactions can also be displayed along with the surface field strengths and this is particularly useful in identifying the major disturbances due to the substructure. 5 Comparison of Measurements and Calculations To examine the effectiveness of the proposed hybridization of finite element and high frequency codes, the cylinder+wing structure shown in Figure 5 was constructed by Naval Air Warfare Center Weapons Division (China Lake, CA)1 and the University of Michi1R. Sliva and H. Wang provided the original cylinder model with the antenna cavities. 70

gan. The patch antenna geometry is shown in Figure 13 and was placed on the cylinder's surface at a = 28.7~, 45~, and 90~ (see Figure 5). The radiation patterns of the patch antenna in the presence of the wing were obtained at the facilities of Mission Research Corporation (Dayton, Ohio) and are shown in Figure 14. It is seen that the patterns at a = 28.7~ and a = 45~ are affected quite substantially by the presence of the plate since the plate is visible to the antenna. However, when the patch is on the cylinder's top (a = 90~), the plate's effect diminishes and the overall pattern is nearly identical to that in the absence of the plate. The computations for the configuration in Figure 5 were done by first using the FEMA-CYL code to obtain the surface electric fields over the non-metallic portion of the aperture housing the patch. These were then grouped in patches of 3x3 pixels in size and turned into equivalent magnetic current infinitesimal dipoles whose stength was set equal to (zMZ + qMO) a6q6&z, where Mz and MO are the surface current densities over the pixel and a 6X 6z is its area. The calculated radiation patterns for the case where a = 45~ are shown in Figure 15 and seen to be in good agreement with measurements. In particular, the pattern based on the UTD code is in better agreement with the measured data in the shadow region (below the plate) but both high frequency codes give nearly identical results in the lit region. This is primarily because APATCH employs the less accurate PO integration to obtain the fields in the shadow region below the plate. In contrast, the UTD code uses geometrical optics and diffraction theory which is known to be more accurate for flat plates with straight edges. This better performance of the UTD code is even more pronounced when the patch antenna is brought closer to the plate's surface as shown in Figure 16. Clearly, this is due to the proximity of the patch to the plate's surface, resulting in stronger secondary currents on plate's surface whose PO/PTD approximation is much less accurate than that provided by the UTD which accounts for ray curvature and surface diffraction effects. The UTD pattern is in good agreement with the measured data except near the shadow boundary associated with the plate. For practical aircraft configurations, the antenna will be placed at far distances from the 71

wing and in that case the APATCH code is more attractive because of its geometrical adaptability and capability to handle materials in the context of the PO approximation. 6 Conclusions We presented a rather simple hybridization of finite element and high frequency codes for antenna pattern calculations in the presence of a complex structure such as an aircraft. Basically, the finite element code was employed to generate the aperture fields on the antennas surface and these were then turned into equivalent magnetic currents. These currents were subsequently used as the sources to the high frequency codes and the antenna pattern was calculated as the sum of the direct antenna radiation pattern and that generated by ray interactions with the substructure. The accuracy of this procedure was verified by comparing the calculated results with measured data. Since the measured set-up consisted of canonical components, the hybridization with the UTD code leads to more accurate results in the shadow region. However, for general airframe configurations, the hybridization with the APATCH code which combines the PO, PTD and SBR methods is more attractive due to its geometrical adaptability and capability to handle material coatings. Although we did not consider antenna loading effects due to the substructure re-radiation, this can be easily incorporated into the analysis by executing the finite element code in the presence of the antenna feed excitation and the secondary excitation arriving at the antenna aperture after undergoing reflection and/or diffraction. References [1] J. M. Jin and J. L. Volakis,"A finite element- boundary integral formulation for scattering by three-dimensional cavity-backed apertures," IEEE Trans. Antennas and Propagat. Vol.39, pp.97 -104, January 1991. 72

[2] J. L. Volakis, A. Chatterjee and L. C. Kempel, "A review of the finite element method for three dimensioanal scattering," J. Optical Society of America A, pp.1422-1433, April 1994. [3] J. L. Volakis, A. Chetterjee, and J. Gong, "A class of hybrid finite element methods for electromagnetics: A review," J. Electromagnetic Wave and Applications, Vol.8, No.9/10, September 1994. [4] J. Gong, J. L. Volakis, A. C. Woo, and H. T. G. Wang, "A hybrid finite element-boundary integral method for the analysis of cavity-backed antennas of arbitrary shape," IEEE Trans. Anten. Prop., Vol.42, No.9, pp.1233-1242, September 1994. [5] J. M. Jin and J. L. Volakis,"A hybrid finite element method for scattering and radiation from microstrip patch antennas and arrays residing in a cavity," IEEE Trans. Antennas and Propagat., Vol.39, pp.1598-1604, November 1991. [6] J. L. Volakis, J. Gong, and A. Alexanian, "Electromagnetic scattering from microstrip patch antennas and spirals residing in a cavity," Electromagnetics, Vol.4, No.1, pp.63-85, 1994. [7] A. C. Woo, J. L. Volakis, G. Hulbert, and J. Case, "Survey of volumetric grid generators," IEEE Anten. Prop. Magaz. Vol.36, No.4, pp.72-75, August 1994. [8] A. C. Woo, T. Ozdemir, J. L. Volakis, G. Hulbert, and J. Case, "''A survey of volumetric grid generators, Part II: Commercial grid generators," IEEE Anten. Prop. Magaz. Vol.37, No.1, pp.96-98, February 1995. [9] R. S. Hansen, ed., Geometrical theory of diffraction, IEEE Press, 1981. [10] T. B. A. Senior and J. L. Volakis, Approximate Boundary Conditions in Electromagnetics, IEE Press, 1995. [11] H. Ling, R. Chou, and S. W. Lee,"Shooting and Bouncing Rays: Calculating the RCS of an arbitrarily shaped cavity," IEEE Trans. Antennas Propagat., Vol.37 pp.194-205, 1988. 73

[12] L. C. Kempel and J. L. Volakis,"Scattering by cavity-backed antennas on a circular cylinder," IEEE Trans. Antennas Propagat., Vol.42, No.9, pp.1268-1279, 1994. [13] L. C. Kempel, J. L. Volakis, and R. Sliva,"Radiation by cavitybacked antennas on a circular cylinder," IEE Proceedings-H, 1995. [14] M. Gilreath, Personal communication, NASA Langley Research Center, September 1994. [15] J. M. Jin and J. L. Volakis,"A biconjugate gradient FFT solution for scattering by planar plates," Electromagnetics, Vol.12, No.1, pp.105-119, 1992. [16] R. W. Freund," Conjugate gradient type methods for linear systems," SIAM J. Sci. Stat. Comput. Vol.13, pp.425-448, 1992. [17] R. W. Freund,"A transpose-free quasi-minimal residual algorithm for non-hermitian linear systems," SIAM J. Sci. Comput., pp.470-482, 1993. [18] L. C. Kempel, "Radiation and scattering from cylindrically conformal printed antennas," Ph.D. Thesis, Dept. of Elect. Engr. and Comp. Sci., University of Michigan, Ann Arbor, MI 48105, April 1]994. [19] J. Gong and J. L. Volakis, "An efficient and accurate model of the coax cable feeding structure for FEM simulations," to appear in IEEE Trans. Anten. Prop., 1995. [20] G. A. Thiele and T. H. Newhouse,"A hybrid technique for combining moment methods with the geometrical theory of diffraction," IEEE Trans. Antennas Propagat., Vol. AP-23, pp. 62-69, January 1975. [21] W. D. Burnside, C. L. Yu, and R. J. Marhefka,"A technique to combine the geometrical theory of diffraction and the moment method," IEEE Trans. Antennas Propagat., vol. AP-23, pp.551 -558, July 1975. [22] G. A. Thiele, "Overview of selected hybrid methods in radiating system analysis," IEEE Proc., Vol.80, No.1, January 1992. 74

[23] R. J. Marhefka and W. D. Burnside,"Numerical Electromagnetic Code - Basic Scattering Code: NEC-BSC (Version 2), Part I: User's Manual," Tech. Rep. # 712242-14 (713742), ElectroScience Laboratory, Department of Electrical Engineering, The Ohio State University, Columbus, Ohio 43212. [24] Users Manual for XPATCH, DEMACO, 100 Trade Centre Drive, Champaign IL, 1993. [25] D. J. Andersh, M. Hazlett, S. W. Lee, D. D. Reeves, D. P. Sullivan, and Y. Chu,"XPATCH: A high-frequency electromagneticscattering prediction code and environment for complex threedimensional objects," IEEE Antenna and Propagation Magazine, Vol.36, No.1, February 1994. [26] J. Baldauf, S. W. Lee, L. Lin, S. K. Jeng, S. M. Scarborough, and C. L. Yu,"High frequency scattering from trihedral corner reflectors and other benchmark targets: SBR versus experiment," IEEE Trans. Antennas Propagat., pp.1345-1351, September 1991. 75

Figure 1. Triangular surface mesh of ridge horn used for a volume tetrahedral mesh. Note that the ridge is not in the principal plane of the horn. Right Angled Brick Skewed Brick Curvilinear Brick i. Skewed Triangular Prism Cylindrical shell Cylindrical shell Tetrahedron Figure 2. Some commonly usl tessellation elements for generating volume mesh. 76

Conformal 'Antenna/Arrav Conformal Artificial Termination Surface Vehicle Platform Figure 3. Printed antenna configuration recessed in a treated metallic surface. / Interference Lobe Figure 4. Illustration of antenna action with platform. pattern illumination and inter 77

antenna aperture dia= 30.48 cm Cl cylinder wing Figure 5. Measurement model for evaluating the hybrid FEM-High Frequeincy implementation. (a) ------------ abc Patches Region I ^... s. '. ""../ Composite skin R er Rgi on (b) stripline Figure 6. (a) Antenna array on cylindrical platform, (b) FE-ABC modeling of an antenna element 78

(a) (b) - 1 0 - 3.- -io, V s -20._ -30 z -180 z a-:) 3 a, -3 C3 c: c 180 Z 10 0 -10 -20 -30 L-180 -90 0 90 Angle ()) [deg] (0=90") -90 0 90 Angle (0) [deg] (0=90") 180 A --- ABC mesh P termination,- surface (z=0) E =2.33 JfBS}f -' 7 - j 'i::::: r!.:.-e,'.~;g:,, i......... 0.4' - -— * i^S;::^(^ ' j"'! -^ BI mesh -_ termination surface ~r =2.09 (c) 1" 12". 12" 12" z, 12" I Frequency = 3.52 GHz X= 3.36" ix 0.69" 1.973" Figure 7. Comparison of measured and calculated radiation patterns for a rectangular patch on a cylindrical platform, (a) without dielectric covering, (b) with dielectric covering, (c) test body: circular ogive cylinder. 79

1 o v s J 4 0 c A / 0 x Measunr Probe FEM (Voltage Mnclf-h eAll, t) u r_ -1 I...O C) C-I E C.L. r_ 0-4 C; 4-0 ctI ZLc"I C: 'a ct E * —4 K0 0 I 2.5 2.5 2& ZtS 2.7 t.S 2.s 2.o 2.9 Frequency (GHz) Frequency (GHz) improved coax Figure 8. Input impedance calculations using an ial cable modeling. -v v v v v v V \/\/ V-\/ Figure 9. Modeling of patch antennas using tetrahedral elements. 80

FEM Mesh Term inated by Assum ing Dom inant Mode 0.18128 cm E=1.12 ~=3.27 //0.08128cm -7 -- * I 6764 cm,QJ6764 cm +.00 ~-1.00 Feed 0.05588 cm Figure 10. Illustration of an aperture fed dual circular patch antenna. The portion of the cavity above the slot is modeled using tetrahedrals whereas the lower microstrip feed is modeled using brick elements. 81

Surface normal geometric Distorted prism x Figure 11. Illustration of the prismatic mesh and element. 82

Input Impedance \vs. FrequencN 0() -50 2.0 2.5 3.0 3.5 4.0 4.5 5.0 Frequency (GHz) - Finite Element - Boundary Integral (FEMA-TETRA) - Finite Element - Metal Backed Absorber (FEMA-PRIS M) Resonant frequency = 3.55 GHz Error = 1.5%c r = 1.3cm ((). 15 \\L at resonaltnc r, = 2. 17cm (0.25 \VL at rcsonancc t r. = 2.6cm (0.3 V\L at rcsonance r = 1.3cm (0.15 \WL at resoIance) O. CIl A r 1-j 2.7 ~,.= 2.4 Vertically directed Field under the patch with surface mesh superimposed Feed 500 V CavityiI P b 400 Patch 'LiI F i 200 Figure 12. Accuracy of the artificial absorber truncation schemes for input impedance computations. Figre12 Acuayo ~eartfca bsre rnato cee ~foariptimeacecmuttos /cm 83

2.36" 7 A 0.79" Patch -n,,.IL:~ - -~ ~ ~ Probe feed 1.97" -TCavity 0.031" 0.031" Freq. = 3.3 GHz. — 1.18" ----- 0.44" _...- T -...... i..; r= 2.17 Figure 13. Patch antenna geometry and material properties. m -15 E-20 |-25 2. (, N \ -150 -100 -50 0 50 Observation Angle 100 150 A (deg) Figure 14. Measured pattern of the patch for the configuration shown in Fig. 5. The shown pattern in the absence of the plate was calculated using FEMA-CYL 84

10 0 -10 m -20 a -30 -40 -50 -60 a =45~ -150 -100 -50 0 50 100 150 Observation Angle e (deg) Figure 15. Comparison of measured and calculated patterns when the pattch antenna is placed at a = 45~ using the hybridization of the finite element and high freguency codes. 10 0 -10 m -20 v ~ -30 -40 -50 -60 \ a =28.7~ -150 -100 -50 0 50 100 Observation Angle 0 (deg) 150 Figure 16. Comparison of measured and calculated patterns when the patch antenna is placed at a = 28.7~ using the hybridization of the finite element and high frequency codes. 85