Preface This report was also submitted by the first author (LCK) as a dissertation in partial fulfillment of the requirements for the Ph.D. degree at the University of Michigan. The earlier part of this work on the finite element-boundary method and its applications to antenna radiation was funded by the NASA Langley Research Center. The latter part of the work on scattering and the finite element-ABC method was funded by the Naval Weapons Center (China Lake) through NASA-Ames and the Rome Laboratory through Mission Research Corp. i

031173-1-T NASA Interchange No. NCA2-839 Grant Interchange Title: Analysis of Cylindrical Wrap-Around and Doubly Conformal Patch Antennas Via the Finite Element-Artificial Absorber Method Report Authors: Leo C. Kempel and John L. Volakis Primary University Collaborator: John L. Volakis Volakis@ um.cc.umich.edu Telephone: (313) 746-0500 Primary NASA-Ames Collaborator: Alex Woo woo@ra-next.arc.nasa.gov Telephone: (415) 604-6010 University Address: Radiation Laboratory Department of Electrical Engineering and Computer Science The University of Michigan Ann Arbor, MI 48109-2122 Date: April 1994 The support of NASA-Ames Research Center, Moffett Field, CA was received under interchange No. NCA2-839.

TABLE OF CONTENTS P R E FA C E..................................... i LIST OF FIGURES............................... v LIST OF TABLES................................ x LIST OF APPENDICES............................ xi CHAPTER I. Introduction.............................. 1 1.1 Background............................ 1 1.2 Overview............................. 4 II. Basic Electromagnetics and Numerical Techniques...... 8 2.1 Basis Electromagnetics..................... 8 2.2 Numerical Techniques...................... 14 2.2.1 The Rayleigh-Ritz Method.............. 15 2.2.2 The Method of Weighted Residuals........ 16 2.2.3 Matrix Assembly................... 17 III. Finite Element-Boundary Integral Method........... 19 3.1 FE-BI Formulation........................ 20 3.2 Vector Weight Functions..................... 24 3.3 Finite Element Matrix...................... 26 3.4 Boundary Integral Matrix.................... 29 3.5 Excitation: FE-BI........................ 33 3.6 System Solution......................... 38 3.7 Scattering and Radiation.................... 41 3.7.1 Far-field Evaluation.................. 41 3.7.2 Input Impedance....................... 43 IV. Finite Element-Absorbing Boundary Condition Method... 45 11iii

4.1 Form ulation....................... -.. 4.2 Conformal ABCs......................... 52 4.2.1 First-order ABC.................... 54 4.2.2 Second-order ABC................... 56 4.3 Excitation: FE-ABC....................... 61 4.3.1 Aperture Term..................... 61 4.3.2 Material Term........................ 62 V. Scattering................................ 68 5.1 FE-BI: Validation........................ 70 5.2 Scattering Studies........................ 75 5.3 FE-ABC: Validation...................... 79 V I. R adiation................................ 93 6.1 FE-BI: Validation...................... 96 6.2 Radiation Studies........................ 101 6.3 FE-ABC: Validation....................... 112 V II. Conclusion............................... 122 7.1 Sum m ary............................. 122 7.2 Future Work........................... 125 APPENDICES..................................128 BIBLIOGRAPHY............................. 161 iv

LIST OF FIGURES Figure 2.1 Material interface............................. 10 3.1 Illustration of the cavity geometry situated on a metallic cylinder along with the coordinate system used in this work.......... 21 3.2 Cylindrical shell element......................... 25 3.3 Geodesic paths on a circular cylinder.................. 31 3.4 Probe-fed conformal patch antenna............... 37 4.1 Typical coated cavity-backed patch antenna with ABC mesh termination................................... 47 5.1 Comparison of RCS for a planar patch (3.678 cm x 2.75 cm) residing on a 7.34 cm x 5.334 cm x 0.1448 cm cavity filled with or = 4 dielectric and a corresponding quasi-planar patch on a large radius (10Xo) cylinder.............................. 71 5.2 Comparison of 2-D moment method results and FE-BI results for a 45~ x 5Ao x 0.1Xo air-filled cavity which is recessed in a 1Ao cylinder. 73 5.3 Far-zone subtraction scheme for the simulation of an infinite cylinder. 74 5.4 Comparison of the RCS computed via the FE-BI method and a BOR code for a 3Ao x 0.1oA air-filled wraparound cavity excited by a Hpolarized (y = 90~) plane wave with oblique incidence (6i = 80~, Xi = 0~)..................................... 76 5.5 RCS frequency response for a 2 cm x 3 cm x patch residing in a 5 cm x 6 cm x 0.07874 cm cavity with cr = 2.17 as a function of curvature for E-polarization (- = 0~)............... 77 v

5.6 RCS frequency response for a 2 cm x 3 cm x patch residing in a 5 cm x 6 cm x 0.07874 cm cavity with r = 2.17 as a function of curvature for H-polarization (, = 90~)................. S 5.7 Illustration of two types of wraparound arrays: (a) continlous cavity array; (b) discrete cavity array..................... 80 5.8 Comparison of E-polarized monostatic RCS at 3 GHz for a four patch array placed on a wraparound collar or in four discrete cavities. The patches and cavities are identical to the ones used in Figure 5.5. The observation plane is 0 = 90~...................... 81 5.9 Comparison of H-polarized monostatic RCS at 3 GHz for a four patch array placed on a wraparound collar or in four discrete cavities. The patches and cavities are identical to the ones used in Figure 5.6. The observation plane is 0 = 90~....................... 82 5.10 Comparison of H-polarized monostatic scattering at 3 GHz by a four patch array placed on a wraparound collar or in four discrete cavities. The patches and cavities are identical to the ones used in Figure 5.6. The observation plane is ) = 0~..................... 83 5.11 Cavity-backed patch antenna with ABC mesh termination. 86 5.12 Convergence of the second-order ABC as a function of displacement from the cavity aperture for E-polarized bistatic scattering. The reference data is provided by a rigorous FE-BI formulation for the same cavity-backed antenna....................... 87 5.13 Magnified view of Figure 5.12...................... 88 5.14 Convergence of the second-order ABC as a function of displacement from the cavity aperture for H-polarized bistatic scattering. The reference data is provided by a rigorous FE-BI formulation for the same cavity-backed antenna....................... 89 5.15 Magnified view of Figure 5.14..................... 90 5.16 Monostatic scattering by a 2 cm x 3 cm x patch residing in a 5 cm x 6 cm x 0.07874 cm cavity with cr = 2.17 for E-polarization (7 = 0~). 91 5.17 Monostatic scattering by a 2 cm x 3 cm x patch residing in a 5 cm x 6 cm x 0.07874 cm cavity with cr = 2.17 for H-polarization ( = 90~)................................. 92 vi

6.1 Illustration of (a) a circumferentially polarized patch element: and (b) an axially polarized patch element. The radius of the cylinder is denoted by a.............................. 7 6.2 Comparison of measured [54] and computed data for a circumferentially polarized element (E-plane) and an axially polarized element (H-plane). The 3.5 cm x 3.5 cm antenna was printed on a cylinder with a radius of a = 14.95 cm and a 0.3175 cm coating having or = 2.32. The probe feed was place at (af., zf) = (-1.0,0.0) for the circumferentially polarized patch and at (a0f. zf) = (0.0,-1.0) for the axially polarized antenna.................................. 99 6.3 H-plane pattern for a four element patch array. Each of the patches is 2 cm x 3 cm and are placed symmetrically around the cylinder. Only the patch centered at 0~ is fed while the other patches are terminated with 50Q loads....................... 100 6.4 Effect of cavity size on the E-plane radiation pattern of a circumferentially polarized patch antenna..................... 102 6.5 Effect of cavity size on the H-plane radiation pattern of an axially polarized patch antenna........................ 103 6.6 Resonance behavior of a circumferentially polarized patch antenna for various cylinder radii....................... 104 6.7 Resonance behavior of an axially polarized patch antenna for various cylinder radii........................ 106 6.8 Variation of the radiation pattern shape with respect to curvature for a circumferentially polarized antenna............... 107 6.9 Variation of the radiation pattern shape with respect to curvature for an axially polarized antenna.................... 108 6.10 Input impedance of a circumferentially polarized patch antenna for various cylinder radii. The frequency range was 2.4 GHz to 2.7 GHz and the cavity size was 14 cm x 14 cm................ 109 6.11 Input impedance of an axially polarized patch antenna for various cylinder radii. The frequency range was 2.4 GHz to 2.7 GHz and the cavity size was 14 cm x 14 cm..................... 110 vii

6.12 Coverage of various uniformly fed discrete wraparound arrays. Each patch is 2 cm x 3 cm and the feed point is at (of = ocenter-f = -0.375 cm) where Ocentr denotes the azimuthal center of each patchi. The array is operated at 3 GHz...............1....... 11 6.13 Coverage of various uniformly fed continuous wraparound arrays. Each patch is 2 cm x 3 cm and the feed point is at (fy = center^, f = -0.375 cm) where,center denotes the azimuthal center of each patch. The array is operated at 3 GHz..................... 113 6.14 Coverage of various uniformly fed collar wraparound antenna. The metallic collar has a height of 3 cm and is fed at symmetric points around the cylinder with zf = -0.375 cm. The antenna is operated at 3 G H z................................. 114 6.15 Convergence of the second-order ABC as a function of displacement from the cavity aperture for axially polarized patch. The reference data is provided by a rigorous FE-BI formulation for the same cavitybacked antenna.............................. 116 6.16 Magnified view of Figure 6.15...................... 117 6.17 Convergence of the second-order ABC as a function of displacement from the cavity aperture for circumferentially polarized patch. The reference data is provided by a rigorous FE-BI formulation for the same cavity-backed antenna....................... 118 6.18 Magnified view of Figure 6.17...................... 119 6.19 Input impedance calculations as a function of ABC surface displacement from the cavity aperture. The frequency range was 3.0 to 3.2 GHz with data taken every 5 MHz................... 120 C.1 Flat resistive strip geometry....................... 139 C.2 Periodic sequence with even symmetry................ 144 C.3 Asymmetric periodic sequence...................... 146 C.4 Comparison of operation count for a full matrix-vector product verses a FFT-based matrix-vector product for a 2-D formulation...... 148 C.5 Circular arc geometry......................... 149 viii

C.6 3-D geometries: (a) flat impedance insert in a groundplane (or impedance plate in free-space) (b) impedance insert in a circular cylinder (2a A 0)........................... 1..' ) 2 C.7 Rectangular patch discretization: (a) u-directed edges (b) v-directed edges..................................... 155 C.8 Local edge numbering scheme..................... 156 C.9 Comparison of operation count for a full matrix-vector product verses a FFT-based matrix-vector product for a 3-D formulation. N is the number of nodes in each direction of the grid, M1 = 2(N-3) and A12 = 2(N -2)................................. 160 ix

LIST OF TABLES Table A.1 Zeros of the Airy Function....................... 130 A.2 Constants for g(O)(J) and g()()...................... 132 A.3 Constants for f(O)(().......................... 133 x

LIST OF APPENDICES Appendix A. Fock Functions............................... 129 B. Iterative Solvers......................... 134 B.1 BiCG Algorithm......................... 134 B.2 CGS Algorithm.......................... 136 C. FFT-based Matrix-Vector Product.................... 138 C.1 2-D Integral Equations...................... 138 C.2 3-D Integral Equations...................... 151 xi

CHAPTER I Introduction Microstrip patch antennas offer considerable advantages in terms of weight, aerodynamic drag, cost, flexibility and observables over more conventional protruding antennas. These flat patch antennas were first proposed over thirty years ago by Deschamps [1] in the United States and Gutton and Baisinot [2] in France. Such antennas have been analyzed and developed for planar as well as curved platforms. However, the methods used in these designs employ gross approximations, suffer from extreme computational burden or require expensive physical experiments. The goal of this thesis is to develop accurate and efficient numerical modeling techniques which represents actual antenna structures mounted on curved surfaces with a high degree of fidelity. 1.1 Background Microstrip patch antennas are conformal metallic elements printed on a grounded dielectric substrate. Although a single patch element offers low gain as compared to other antennas, an array of such elements offers considerable gain and pattern shaping flexibility in a low profile package. Traditionally, such antennas were designed using either expensive measurements or approximate formulae such as the cavity model. 1

2 The cavity model [3] represents the radiating structure by magnetic currellt walls placed around the perimeter of the patch between the groundplane and the surface of the substrate. The original cavity model assumed a planar-rectangular patch but it has been extended to cylindrical-rectangular patches [4. 5]. The cavity model has proven quite accurate for radiation pattern calculations near resonance. However. it is not suitable for operation away from resonance since it is based upon a single mode assumption. Additionally, it is not especially accurate for input impedance calculations since probe feeds typically excite multiple modes within the cavity. The cavity model also ignores mutual and surface wave coupling between antenna elements which is often a concern to array designers since it results in parasitic energy loss. Finally, the cavity model cannot be used for scattering calculations. To overcome the deficiencies of the cavity model, rigorous integral equation formulations have been investigated [6]-[9]. Integral equations utilize an appropriate Green's function and equivalent currents to represent the fields. Such an approach offers a rigorous representation of the local fields which includes mutual coupling and surface wave effects. However, this approach suffers from two primary drawbacks: computational expense and memory demand. In particular, accurate computation of the substrate Green's function is both difficult and expensive since an accurate evaluation of a Sommerfeld integral is required [10]. Nevertheless, the integral equation technique has been extended to patch antennas printed on a coated circular cylinder [11, 12], which requires an extremely complex substrate dyadic Green's function. To reduce the memory demand and matrix fill time from 0(N2) to 0(N), the Biconjugate Gradient-Fast Fourier Transform (BiCG-FFT) may be used provided the substrate is either flat or circular and uniform discretization is used. Although integral equation formulations have proven quite useful for design and

3 analysis of patch antennas, they require an assumption which limits their fidelity for real-world antennas. Since a substrate Green's function is used in the formulation. the substrate and underlying metallic surface must be infinite. In practice, antenna designers place each radiator within a metallic cavity for suppression of unwanted coupling [13]. These structures have finite substrates which require a substantially more complex (and expensive) dyadic Green's function for their implementation. Thus, integral equation formulations cannot efficiently model such arrays and hence are of limited value in their design. To overcome the computational burden of integral equations, an alternative methodology was first proposed separately by Silvester and Hsieh [14] and McDonald and Wexler [15] which is now known as the Finite ElementBoundary Integral (FE-BI) method. Jin and Volakis [16] adapted the FE-BI technique for planar cavity-backed patch antennas. In this approach, planar-rectangular cavities were discretized using brick-shaped finite elements and the computational domain was closed with a boundary integral which provided an exact mesh closure condition. This method was found to be capable of modeling large complex arrays of cavity-backed patch antennas with low memory and computational burden. In addition to radiation by finite, possibly aperiodic, arrays of such structures, Jin and Volakis also investigated their scattering behavior [17]. The FE-BI method has been applied to two-dimensional [18], Body of Revolution [19] and planar three-dimensional problems [16, 20]. The primary reason for the lack of interest in applying the FE-BI approach to more general curved geometries is that the method is no longer efficient in terms of memory or computational demand since FFTs may not be used in conjunction with the boundary integral. However, this technique can be applied to cylindrical bodies as well as planar geometries without excessive computational cost. A formulation can be developed which yields compa

4 rable performance by utilizing special finite elements. dvadic Green's functions aInd the BiCG-FFT solver. The FE-BI technique is both accurate and efficient for cavities embedded in a ground plane or in a metallic cylinder. However, for any other platform, the computational burden of the FE-BI approach is excessive since the BiCG-FFT solver may not be used. In fact, the only dyadic Green's function available for arbitrary surfaces is the free-space Green's function which would require both electric and magnetic currents over the entire mesh boundary. Accordingly, Absorbing Boundary Conditions (ABCs) have been proposed to provide an approximate mesh closure condition which avoids the use of an expensive boundary integral. ABCs have been used extensively for large scattering applications [21] but not for antennas. Unfortunately, most ABCs assume a spherical mesh boundary [22, 23] which would result in excessive unknowns. Chatterjee and Volakis [24, 25] have proposed conformal vector ABCs which follow the shape of the underlying body and in this work they are applied to scattering and radiation by cavity-backed conformal antennas which are recessed in a metallic cylinder. Conformal ABCs have also recently been proposed by Stupfel [26]. Verification of the accuracy of the FE-ABC method for antenna applications is important since it is the only practical means of modeling large doubly curved and possibly coated conformal arrays. 1.2 Overview In this thesis, the FE method is extended to cavity-backed conformal antenna arrays embedded in a circular, metallic, infinite cylinder. Both the boundary integral and absorbing boundary mesh closure conditions will be used for terminating the mesh. These two approaches will be contrasted and used to study the scattering and

radiation behavior of several useful antenna configurations. An important feature of this study will be to examine the effect of curvature and cavity size on the scattering and radiation properties of wraparound conformal antenna arrays. Chapter 2 introduces the fundamental electromagnetic principles used in the development of the FE-BI and FE-ABC formulations. Included are Maxwell's equations, natural boundary conditions, the radiation condition, integral representation of fields in terms of the equivalence principal and the definition of radar cross section. In addition, the notation used throughout this thesis will be presented. An introduction to numerical formulations will be presented along with a discussion of matrix assembly. In particular, the integral equation and finite element methods will be compared with respect to matrix cost. In Chapter 3 the FE-BI formulation is developed for circular cylinders. Of particular interest are the newly introduced cylindrical shell elements and an efficient evaluation of the cylindrical dyadic Green's function. Appropriate plane wave, interior source, far-zone field and input impedance formulae are given. The emphasis of this formulation is on efficiency in terms of both memory and computational demand. In view of this goal, the formulation will assume uniform surface gridding so that the BiCG-FFT solver may be used. Also, large radius asymptotic formulae for the cylindrical dyadic Green's function will be presented along with the traditional modal formulae and criteria will be established for the appropriate use of these formulae. Chapter 4 introduces the conformal ABCs for mesh termination. The FE-ABC equations are developed which utilize a novel mixed total/scattered field formulation. This approach minimizes the potential error propagation usually associated with a total field formulation applied in the exterior region while retaining the lower computational burden of a total field formulation when used in the interior region.

6 New conformal vector ABCs which were developed by Chatterjee aild \olakis of the Radiation Laboratory are used for mesh closure. Both first and second order ABCs are presented for cylindrical-rectangular surfaces. The first order ABC always produces a symmetric system while the second order ABC is asymmetric. The second order ABC becomes symmetric for spherical boundaries and in the limit as the radius of the ABC surface becomes large. In Chapter 5 we investigate the scattering behavior of conformal antenna arrays. The FE-BI method is used to study the effect of curvature on the scattering pattern and the resonant frequency. Both principal polarizations are discussed. Also, the effect of cavity size is examined. One array of particular interest is the wraparound array. The difference in scattering behavior between discrete and continuous wraparound arrays is examined and related to known creeping wave and diffraction phenomena. The FE-ABC method will be compared to the FE-BI formulation and in particular, the minimum ABC displacement from the cavity aperture is determined. Such information is critical in minimizing the computational domain. Also, the FE-ABC method will be shown to recover the FE-BI data for monostatic scattering even when the transmitter/receiver is located behind the cylinder. Chapter 6 provides a corresponding study of radiation characteristics for conformal patch antennas. The effect of curvature and cavity size will be determined for both axially and circumferentially polarized antenna elements. Additionally, the effect of curvature on the input impedance of both types of patches will be examined. An example of pattern coverage will be given for discrete, continuous and collar-type wraparound antennas. The pattern behavior as the number of antenna elements increases is demonstrated. Finally, the FE-ABC method is compared with the FE-BI method in order to determine the minimum required ABC displacement for accurate

radiation pattern calculations. The performance of the FE-ABC method for input impedance computations is also discussed.

CHAPTER II Basic Electromagnetics and Numerical Techniques 2.1 Basis Electromagnetics James Clerk Maxwell presented a set of six, coupled, time-dependent equations which completely describe all electromagnetic phenomena. Although Maxwell's equations allow for very general field variations in both space and time, in this work only time-harmonic fields will be considered where a e+jwt time-dependency is assumed and suppressed. In addition, in this work only fields in a linear, isotropic medium will be considered. Accordingly, the following constitutive relations D = eE (2.1) B = iH (2.2) have been used to develop Maxwell's equations in differential or point form V x E = -Mi-jkZH (2.3) V xH = J + jkYE (2.4) where D is the electric flux density, E is the electric field, B is the magnetic flux density, H is the magnetic field, Ji is an impressed electric current and Mi is a fictitious magnetic current which is included so that (2.3) and (2.4) are symmetric. The medium is completely described by its intrinsic impedance (admittance) 8

9 Z = () = - / where c and p denote the medium's permittivity and permeability. respectively. The wavenumber is denoted by k = - =;/7. where A is the wavelength and wc is the corresponding angular frequency. An important medium is free-space where k = k = 2 = w/oco with Ao being the free-space wave length. 60 = 8.854 x 10'2F/m, yo = 47r x 10- H/m and the intrinsic impedance (admittance) is now denoted as Zo = ( -) = / 377Q. Faraday's (2.3) and Ampere-Maxwell's (2.4) Laws are not sufficient for the complete specification of a general field. One of four dependent scalar equations must also be chosen which may be derived from (2.3) and (2.4) V.E = pe/ (2.5) V -H = pm/l (2.6) V. J = -j'pe (2.7) V M, = j=pm (2.8) In these, pe is the electric charge density (c-) and likewise Pm denotes the fictitious magnetic charge density (s-). The latter is not a physical quantity but its introduction is important for mathematical constructions. Evidently, (2.5) and (2.6) are a statement of Gauss' Law in differential form whereas (2.7) and (2.8) are referred to as the continuity or charge conservation equations. Typically, (2.3), (2.4) and either (2.5) or (2.7) are used to determine the field quantities. Such a solution is, however, not unique and will require the specification of additional constraints. Application of Maxwell's equations (2.3) and (2.4) along with a derived equation (2.5)-(2.8) will result in an infinite number of solutions. A unique solution is obtained only after the imposition of various conditions such as boundary conditions, radiation conditions, edge conditions and transition conditions. The natural boundary

10 conditions satisfied by E and H at any dielectric interface are [27] x (El-E2) = 0 (2.9) n x (H- H2) = 0 (2.10) where the subscripts denote the different media which are joined at the interface and n is the normal unit vector directed into medium one as shown in Figure 2.1. Medium one is the external medium while medium two is the internal medium. medium 1 n ___ ElH1 ______ n_______ E:1,111 medium 2 Figure 2.1: Material interface. The first condition (2.9) states that the tangential components of the electric field are continuous across the dielectric boundary while the tangential magnetic field components are also continuous for a dielectric. Boundary conditions may also be written in terms of normal field components n (6, - ) - = 0 (2.11) n. ( H1- 2H2) = 0 (2.12)

11 which simply states that the normal components of both the electric and magnetic fields are discontinuous if there is a discontinuity in the permittivity and permeability, respectively. Another important interface involves a dielectric (medium 1) and perfect electric conductor (medium 2). The tangential electric field must now vanish at such an interface n x E =0 (2.13) since a perfect electric conductor cannot sustain an internal field. Note that a surface current may now flow along the interface n x H1 = J, (2.14) Similar conditions for the normal field components now may support a surface charge density. For perfect magnetic conductors, appropriate boundary conditions may be determined via duality. For this work, open scattering and radiation problems will be considered. Consequently, any valid and unique solution must also satisfy the Sommerfeld radiation condition [28] which describes the field behavior at infinity limnr Vx H +jk0rx H =0 (2.15) This simply states that the field is outgoing and of the form e-kr/r as r -. o Faraday's Law (2.3) and Ampere-Maxwell's Law (2.4) are coupled first order partial differential equations (PDEs). However, these may be combined to yield a single second order PDE in terms of E or H. Specifically, by taking the curl of one of these and making use of the other, the following second-order PDE is obtained E 2 E P oi one VxVx { } k0 } { } V x {'} (2.16) These are commonly referred to as the vector wave equations. Solution of the wave equation also requires sufficient boundary conditions for uniqueness.

12 Direct solution of NMaxwell's equations or the wave equation is only possible for a few special conditions. As a consequence, equivalent sources are often enlI)loyed in the solution of practical problems since often such an equivalent problem may be readily solved numerically. The surface equivalence principle states that any field exterior to a given (possibly fictitious) surface may be exactly represented by equivalent surface currents which are allowed to radiate into that external region. These equivalent currents are given in terms of a jump discontinuity between the total exterior (E1, H1) and interior (E2, H2) fields n x (H-H2) = J (2.17) n (2 - E) = M (2.18) The radiated fields due to these equivalent currents are given by the integral expressions E(r) = (s i(') x E(') V. x G(R) dS' -jkZ f i(') x H(r') G (R) dS' (2.19)!-/I H() = ff n(') x H(r') V x G(R) dS + jkY jf n(r) x E(r') G (R) dS' (2.20) where R = Ir - F', r denotes the observation point and r' is the integration point. The surface on which the equivalence theorem is applied is denoted by S. In (2.19) and (2.20), a dyadic Green's function is required which satisfies at least the radiation condition (2.15). Such a Green's function is the free-space dyadic Green's function Go I + ( 2 Go(R) (2.21) which utilizes the corresponding scalar Green's function e-jkoR Go(R) = (2.22) 4irR

13 Note that use of (2.21) requires both equivalent currents (2.18). Alternatively, w\lhen certain bodies such as an infinite metallic plane. cylinder or half plane is present. a dyadic Green's function can be chosen which satisfies an additional boundary condition such as the Dirchlet condition n x G1= 0 on S (2.23) so that (2.19) and (2.20) may be written in terms of only one current per field expression E(r = (f n(~) x E(r') V x Gi(R)dS' (2.24) H(r) = f n(i ) x H(r). V x Gi(R) dS (2.25) The Neumann boundary condition is also useful n V x G2 = 0 onS (2.26) with the resulting field integral expressions in terms of only one current E(F) = -jkZ n( f') x H(r') G2(R) dS' (2.27) H() = jkY nr (') x E(') G2(R) dS' (2.28) Such dyadic Green's functions are commonly referred to as a first kind dyadic Green's function and a second kind dyadic Green's function, respectively [29]. Often, engineering quantities such as the radar cross section (RCS) and the radiation pattern require evaluation of one of the previous integral expressions as the observation point (r) recedes to infinity. For example, the RCS is defined by 2 liM I ElE1s the = 47rr2 limr E( )I (2.29) ~IEi(rl where ES(3) is the scattered field and E(rF) is the incident field. The scattered field may be defined by either (2.19), (2.24) or (2.27). Analytical expressions for

14 these fields ma' sometimes be obtained: however. the required currents arc often determined via an appropriate numerical method. 2.2 Numerical Techniques Determination of electromagnetic fields within a region typically involve either integral or partial differential equations. Such formulations may be written in operator notation as Lu-f = 0 (2.30) subject to appropriate boundary or transition conditions B(u) = 0 (2.31) within the domain (Q) and on its boundary (&Q). In these, the operator L is based either on an integral representation of the fields such as (2.19)-(2.20) or upon the vector wave equation (2.16). The forcing function f is a known quantity while u, which is often a field or current, is unknown. Unfortunately, very few analytical solutions for (2.30) are available in electromagnetics. One such solution, the fields due to a magnetic dipole in the presence of a metallic cylinder, will be used in Chapter III to form the appropriate dyadic Green's function. However, most useful electromagnetic scattering and radiation problems cannot be solved using analytical methods. Rather, an approximate numerical solution is sought which in some way closely resembles the exact solution. Two methods of formulating such an approximate solution are: the Rayleigh-Ritz method and the method of weighted residuals. Both of these approaches are discussed by Zienkiewicz and Taylor in their classic finite element method text [30].

15 2.2.1 The Rayleigh-Ritz Method The Rayleigh-Ritz method seeks a stationary point of a variational functional. For operators which are self-adjoint and positive definite, the stationary point of the following functional F(i) = (u)-(u) (2.32) is an approximate solution of (2.30). In (2.32), the inner product over the domain Q of the two functions is defined as (a, b) = ab dQ (2.33) or for vector functions (a, b) = j a'bdQ (2.34) The choice of this inner product extends the validity of the variational expressions to complex fields. The trial function, u, is expanded in terms of N basis functions N u = ujw = {u}T {w} (2.35) j=i where wj are the basis functions and uj are the unknown expansion coefficients. In (2.35), column data vectors are denoted with {.} while row data vectors involve a transposition {})T. Substituting (2.35) into (2.32), the functional may now be written F(u) = 2 {u} {w} {w} dQi {u}- {u} I {w} f d (2.36) This functional is minimized by allowing all partial derivatives with respect to the coefficients, {u}, to vanish uF(u) = 2f [wiw}T dQ] {u}+ {u}T j{w}wId-n wif d = 0 (2.37)

16 which may be written as a matrix system [A] {u = {b} (2.38) The elements of the matrix [A] and excitation vector {b} are given by Aij = j wiCwj dQ b = J wif dQ (2.39) A word of caution, unlike other branches of engineering, no physical significance may be attached to the stationary point of the functional (2.32). In mechanical systems for example, minimizing this functional represents minimization of the total potential energy of the system. However, since electromagnetics involves complex quantities, such a statement may not be asserted. Additionally, the variational formulation given above requires a symmetric operator ~ which is not always available in computational electromagnetics (see for example the FE-ABC formulation in Chapter IV). However, the method of weighted residuals does not require a symmetric operator and is perhaps more familiar to electrical engineers since it does not involve variational principles. 2.2.2 The Method of Weighted Residuals The method of weighted residuals attempts to minimize the residual 1? = CLu-f >0 (2.40) where R = 0 would imply that i = u for all points within Q. Such a solution in general cannot be obtained and it is more realistic to find a solution which satisfies the minimal residual condition in some average or weighted sense / [ti {w}T {u} - tif ] d = 0 i= 1,2,3,...,N (2.41)

17 In general. any testing function (ti) may be used: hoNwever. since these functions modify the enforcement of the boundary conditions throughout the domain, the choice of testing functions effects the quality of the solution (2.41). One popular testing procedure is called collocation or point matching. In this, the testing function is a Dirac delta function, ti = 6(x - xi) which implies enforcement of the boundary conditions only at discrete points e.g. xi, i = 1,2,3,...,N. Another popular choice is termed Galerkin's procedure. In this, the testing function is identical to the expansion function used in (2.35) e.g. t2 = w2 and the weighted residual equation is given by [j wCiLw dQ]{ u} = j wifdQ (2.42) which is identical to the Ritz procedure given above. In Galerkin's method, the boundary conditions are satisfied throughout the computational domain in a least squares sense. 2.2.3 Matrix Assembly The construction of the system matrix (2.39) associated with either the Ritz or Galerkin's formulation represents a significant portion of the overall computational cost of either the moment method or the finite element method. However, these two popular techniques have vastly different cost characteristics which is the primary attraction of the finite element method for volume formulations. A general moment method matrix requires the computation of all matrix entries which is a O(N2) process where N is the number of unknowns. These matrices are full since each unknown couples with each other through the Green's function. On the other hand, FE matrices are sparse due to the locality of the PDE formulation and hence each degree of freedom interacts only with its neighbors. Such a matrix assembly process

18 is 0(C N) where C is a constant which depends on the type of finite elements used. For typical finite elements used in computational electromagnetics. this constant is between twenty and forty. Therefore, the matrix build time (and Ineinory requirement) for a FE matrix grows only linearly with the problem size while the build time (and memory requirement) for a moment method matrix grows as the square of the problem size. Since the number of unknowns increases roughly as the cube of the linear dimensions for volume formulations, an FE matrix grows as a cubic function of the linear dimension while the corresponding volumetric moment method system increases as the sixth power of the same linear dimensions.

CHAPTER III Finite Element-Boundary Integral Method The Finite Element (FE) Method is a computational technique which approximates a continuous boundary-value problem with a corresponding discrete problem. The latter is then solved numerically, thus obtaining a solution which approximates in some sense the original continuous problem. The FE method has been used in mathematical physics since the 1940's. It was first coupled with an exact Boundary Integral (BI) termination condition in an electromagnetics application by Silvester and Hsieh [14] and McDonald and Wexler [15] in 1971 and 1972, respectively. Recently, there has been renewed interest in the Finite Element-Boundary Integral (FE-BI) method of electromagnetics principally due to the work of Jin and Volakis [16], [17] and [20] who demonstrated the successful application of the method for cavity-backed antennas. Their major contribution was coupling the FE-BI approach with the Biconjugate Gradient-Fast Fourier Transform (BiCG-FFT) technique, thus allowing high fidelity simulations with a low 0(N) memory and computational demand. This chapter will extend the FE-BI method to cavities which are recessed in a metallic circular cylinder of infinite extent. The boundary integral will employ uniform zoning and hence the BiCG-FFT solver may be employed to retain low mem19

20 ory and 0(N log2:N) computational burden.'We will first present tile FE Ilmethod derivation starting with the vector wave equation and then proceeding with the introduction the boundary integral as an exact mesh truncation condition. An efficient computational strategy for the admittance matrix will be presented along with accurate expressions for interior and exterior source modeling. Finally, the computed electric fields will be used to evaluate the far-zone fields and, if applicable, the input impedance of patch antennas. 3.1 FE-BI Formulation Consider a cavity recessed in an infinite metallic cylinder, shown in Figure 3.1. The cavity walls are assumed to coincide with constant p-, 4- and z-surfaces and the cavity is filled with inhomogeneous material. Interior sources and lumped loads may be present as well as surface metallization patterns. The FE-BI system is developed directly from the inhomogeneous vector wave equation (2.16). To generate a system of equations from (2.16), the method of weighted residuals is applied which results in the symmetric inner product of a vectorvalued weight function and the vector wave equation as described in Chapter II. The integro-differential equation is given by -kJ, r(p, z)E(p, X, z) Wi(p,, Z) pddd = Ivix M(P[ z) Wi (p,, z)pdp dc dz -jkZo,J J (p, 4, z)e W(p,i, o z)p dp d dz (3.1) where Wi(p,, iz) is a subdomain vector-valued weight function to be specified and Vi is the ith volume element resulting from a discretization of the cavity. Given

21 l //////// // f/////'///////////^^ r,/, /////// ^ Figure 3.1: Illustration of the cavity geometry situated on a metallic cylinder along with the coordinate system used in this work.

99 the impressed sources (Jf. A1'). the right-hand side of (3.1) is the interior excitation function defined by lnt Aft (Pp 0. t J=,(P. O. 2) J + jZ J p Z) (p, O -) ( l(p, o )p dpdo dz(3.2) Upon application of the first vector Green's theorem [31], we recognize (3.1) as the weak form of the wave equation r V x E(p, q,z) V x Wi (p, z) IVi — r-(prPz)p dp do dz -ko j r(P p,, )E(p,, z) Wi(p,, z)pdp d dz -jkZo is (p,, z) x H(p,, z) W~(p, X, z)dS = ft (3.3) where n(p, (, z) indicates the outward pointing normal of the element surface associated with the ith unknown, Si is the surface area of that element, and H(p, X, z) is the total magnetic field. It can be shown that the surface integral of (3.3) vanishes for all elements which do not border the cavity aperture. Furthermore, their nonzero contribution is limited to the portion of their surface which coincides with the aperture. Unfortunately, (3.3) contains unknown electric and magnetic fields on the surface of the cylinder. To eliminate H from (3.3), we utilize the surface equivalence theorem introduced in Chapter II (2.28) with a second kind dyadic Green's function. This formula provides a relationship between the electric field and the magnetic field on the aperture of the cavity. We expand the total magnetic fields as H(i) = H i( ++Hr() + Hs(I Hcl() + Hs(F) (3.4) where Ht is the incident field, HT is the reflected field due to the infinite metallic cylinder and HS is the scattered field attributed to the cavity. The incident and

23 reflected fields will be specified later in this chapter while the scattered field is represented by (2.28). Simplifying (3.3) with (2.28) and (3.4). the resulting system is given by EX EV (p, ~; ) 7 X V i(p, z) -k2e,(p,, z)Wlj(p, O z)* p, z)pdpdd4dz +(koa)2 EjSa(j)6a(i) j [Wi(a, z), z) (a,, z)x j t G2(a, z, ) x p(a,, z) Wj(a,, z)] do' dz'd) dz = ft + f"' (3.5) where the unknown electric filed is expanded in terms of subdomain basis functions E = EjWj(p,, z) (3.6) 3 In (3.6), the subscript indicates the jth unknown and Wj(p, ), z) is the same edgebased expansion function as that used for testing in (3.1); e.g. Galerkin's procedure. The cylinder has radius a and the relative constitutive parameters of the material filled cavity are denoted by or and p,. The function ba(i)&a(j) is a product of two Kronecker functions and it simply indicates that the boundary integral only contributes when both the test and source unknowns are on the aperture. In addition to the interior excitation functional, fnt, an exterior excitation functional, f t, is present and will be discussed later in this chapter. The appropriate weight functions and the evaluation of the FE matrix entries which are represented by A.j in the FE-BI matrix system AlEl + [ [0] 1 {EP} {fter (3.7) [ {Ejt} } [0] [0] J {E7nt} {fjnt } will be presented next.

24 3.2 Vector Weight Functions An important factor in choosing the finite elements for gridding the cavity is the element's suitability for satisfying the mathematical requirements of the formulation as well as the physical features of the antenna system. Traditional node-based finite elements associate the degrees of freedom with the nodal fields but have proven unsatisfactory for three-dimensional electromagnetics applications since they do not correctly represent the null space of the curl operator and hence spurious modes are generated [32]-[34]. In contrast, edge-based elements correctly model the curl operator and therefore the electromagnetic fields. In addition, edge-based elements avoid explicit specification of the fields at corners where edge conditions may require a singularity. Jin and Volakis [20] presented edge-based brick elements which are convenient for rectangular-type structures and cavities, but for cavities residing in a circular cylinder, shell elements are the natural choice. Cylindrical shell elements possess both geometrical fidelity and simplicity for cylindrical-rectangular cavities. Figure 3.2 illustrates a typical shell element which has eight nodes connected by twelve edges: four edges aligned along each of the three orthogonal directions of the cylindrical coordinate system. Each element is associated with twelve vector shape functions given by W12(p, ), z) = Wp(p, (;, z;., rZt, +), W43(p),,) = Wp(p, (, Z;, 4l, Zt, -) W6(p, (, z) ) = Wp(p, (, z;', r, zb, -), W87(p, (, z) = Wp(p, ), Z; *,, Zb, +) W14(P, O, Z) = Wz(p, Sb, Z; Pb, r, Z, 4+), W263(p,, Z) = WO(p, O, Z; pa, r, Zt,-)

25 z Pb zt xr 2Cn Figure 3.2: Cylindrical shell element.

26 Ws48(p, za) = I (p o. -: pb. o.-), 1.(p., z) = li.(p, o. z: pa. *. +) (3.S) where 1'lk is associated with the edge which is delimited by local nodes (l.k) as shown in Figure 3.2. As seen from (3.8), three fundamental vector weight functions are required for the complete representation of the shell element. They are UW/(p,4 $,z; X ) = SPb ( - O)(Z - Z) W~,(p,, z; p,, z, s)' -P ah p W (p, z s) =,(p -)(- Z) WZ{p,^z, )= -I 0-(Z;p- 1)(q^-q()$ P(3.9) ta where the element parameters (Pa, pb, 0I, 4r, Zb, Zt) are shown in Figure 3.2, t = pb-pa, a =,r - 01 and h = Zt - Zb. Each local edge is distinguished by p, 0, z and s as given in (3.8). The k-term which appears in the definition of the'-directed weight (3.9) is essential in satisfying the divergence free requirement, i.e. so that V. Wj = 1. Note that as the radius of the cylinder becomes large, the curvature of these elements decreases, resulting in weight functions which are functionally similar to the bricks presented by Jin and Volakis [20]. Having specified the vector basis functions, we may proceed to develop the matrix entries for the system (3.7). 3.3 Finite Element Matrix The inherent locality of a partial differential equation formulation results in a highly sparse system, [A]. This FE matrix arises from the discretization of the first two integrals of (3.5) and owing to the finite support of the basis functions its entries are identically zero unless both the test and source edges are within an element. 1 Wj (p, T, z) will only satisfy this requirement within the volume of the element. These weighting functions introduce artificial charges on the faces of the element and are not divergenceless at element interfaces. This is allowable since these elements do not guarantee normal field continuity across the element faces.

27 These matrix entries can be represented as A- - = -Ist - ko 2 1: (3.10) Or where we have assumed constant material properties within an element (6r and Pr) and the subscripts (i,j) refer to the row and column of the FE matrix, respectively. The two auxiliary functions are defined by I~)'j = X, VX ~(p,, z; P, i,,is ). V x Wt(p,, z; pi, i, zi, si)p dp do dz I!2)iJ = /X Ws(p, q, z; tpSj, Zj'zJ' Wt(p, I, z;pi Xi zi si)pdp dodz (3.11) where (s,t) E {fp,,z} indicate the direction of the source and test edges, respectively. Since the fundamental weight vectors (3.9) are aligned along orthogonal directions, (3.11) is symmetric with respect to the source and test edges (e.g. /1)ij = If(l)i and thus, only six unique combinations of (s,t) for 1() need be determined and only three such combinations are required for I(2). They are I ( h) [2iIn (P ) r (- s)(- o )d + (ah)2 U \Pa2 4 (4 - i) t (z - is)(z - t)dz] iP)t( 1 Pb ) I(f) = -1h2 [ n (PA ) +~I (^ pa ) ]~ \ ~\ I = SSttPb ( O((1) 9,(ta h 1l 4 4> 1 13 3 I)2 2 I - = h2 4 b a) (3 s + t) Pa-b + t b f a) (2 (p2- p2 ) - 2t ((is + it) + pspt n ( Pz - )(- I(l) s- tS —t [ (pP - )(p — )dp zz -P- t_-,(2 + + n ) I(1) _ s",s~~th I~(2I \Ybf,- ~- P2P

28 (P(-P2) - (-o o-o ol pp() -:()Pa )J(o - ),(I)t OJ)( - o )do (2) = [ (Pb -P4) + (s + pt) (P3 -P ) + PP (p p2)] x i) =~.. n — ~ (o-o)(-otd ( -'s)(~ -'t)d~ zb /i(Z = t [1 (Pb -a) + (s -+ tt ) (P -Pb + 2Pst (P - P ) ] (ta)2 [4 3 (- Os —)( -t)d0 (3.12) where each of the unevaluated integrals are of the form 1 (1 (f-s) ((-t)do= (L2 U2)(s+ t+ (3- L3) + ( ( L) (3.13) The FE matrix is assembled by evaluating (3.10) for each element combination which contains both the ith and jth unknowns. Since the integrals (3.11) are symmetric, only the upper triangle or lower triangle of the FE matrix need be computed and stored. A lumped impedance post may be included in the formulation by adding a term to (3.5) and equivalently to (3.10); surface or sub-surface metallization layers may also be modeled. Lumped loads are approximated in the FE-BI formulation by a filamentary load [35] located at (4L, ZL) for radial, (PL, ZL) for azimuthal and (PL, AL) for axial posts, respectively. These posts have length 1, cross-sectional area s and impedance ZL. Such posts may be represented as t a(O - OL)a(Z - ZL) Zp(p,, z) = - L(- ZL) for radial posts ZL P Zp(p, 4, z) = — (P- PL)6(z - ZL) for azimuthal posts ZL h S(p- pL)&(b-) L) posts z,(p, O, z) = for axial posts (3.14) ZL P

29 assuming that each post spans only one element. The contribution to [A] is Pixvel b! A.j = jkoZo Zp(p., Oz)l.(po. ~)l11(poz)pdpdodz (3.15) which may be evaluated in closed form for the three post directions given il (3.14) as A^ = jkoZ sisj ( t2h2 [(P - (L - Lj) (ZL - i) (ZL -- j)] A = jkoZLo ij (t22) [(PL - ) (PL - P,)(L - Oi) (0L - j)] (3.16) where the superscript indicates the orientation of the post. In addition, infinitesimally thin metallization layers may be represented by simply fixing a priori the weight coefficients to zero for weights associated with edges which are tangential to the metal. This is a consequence of using a total electric field formulation. The symmetry and sparsity of the FE system [A] is maintained after the addition of these loads while the BI system [9] remains fully populated and symmetric. 3.4 Boundary Integral Matrix The boundary integral provides an exact boundary condition for mesh closure and its construction relies on a cylindrical dyadic Green's function. The entries of the boundary integral sub-matrix are Gj = (k a)2 Wt(a,,z; i, z i, zi, si) L[(a,, z) x G2(a,;), z) x p(a, (;(, z ) Ws(a, (, z'; pj, 4j, zj, sd )d' dz d dz (3.17)

30 where the weight functions are given by (3.9) and are evaluated at tle surface ) = (1. In (3.17), the dSadic Green's function (G2) satisfies both tile radiation condition alnd the Neumann boundary condition at p = a. This dyadic Green's function may be expressed exactly as in [31] 1 - 0 k 21 H,2) GZZ(a )=- ( H2 dk ( \' ) (27r)2 n=ooJ (-O k2H,)(-n) G ~~> n - -(2.. JJ o ['k^ H 2)1G~~~~a 6 ~~ ~0- - e - 0( \2 2 1 - Ia00 Q00 I H,(2)(r) nkz 2 H ()(2) ]jcne —z"dl Goo (a, Z - --- ~ C, - A G (a,'4)- W(27r)2 H)) nkoakp H'Q(y) e d (3.18) where y = kpa and kp = /k -k2. However, for large radius cylinders e.g. ka > 3, (3.18) is computationally prohibitive. In these cases, which are the main concern of this work, it is advantageous to employ an asymptotic expression for G2 [36]-[39]. Such expressions employ a creeping wave series expansion of which only the two direct path contributions, as shown in Figure 3.3, will be retained. The formula due to Pathak and Wang [36] has proven quite accurate and is given by GZZ(a,, z) j e-jkos {(cos0+ q(1 -q)(2-3cs20)) v(/3)}, -3-~eGZ (a,z) Jk qekossincos{( -3q(1 (jko -0ko 22 G (a,,) ~ -Z- qe-ks (sin20 + q(1 - q)(2 - 3sin 0)) v(/3) 27t +q [sec20 (u(i)- v())] (3.19) where = ks [r o0 3 and q =. In the definition of 3, s is the usual geodesic path length (s = (a)2 + z2) and 0 is the direction of the trajectory ( = tan-1 [aj). Depending on which of the two direct paths as shown in Figure 3.3 is used, b = ( (short path) or b = 27r - (long path). The soft and hard Fock functions, u(f)

31 geodesic paths,.................... Figure 3.3: Geodesic paths on a circular cylinder.

32 and v(3) respectively, are characteristic of on-surface creeping wave interactions anld have been extensively investigated by Logan [40]. Care must be taken in evaluating (3.17) so that the overall storage requirement remains 0(N) and the singular integrals of (3.17) are accurately computed. If uniform zoning is used, the resulting sub-matrix ([9]) is block Toeplitz and hence amenable to solution using the BiCG-FFT method. For the non-selfcell contributions, mid-point integration may be used while a regularization procedure must be employed for the self-cell along with a more accurate numerical integration scheme such as four point gaussian quadrature. Bird [39] noted that (3.19) recovers the metallic screen Green's function when 3 = 0 within the available approximation order. This suggests that (3.17) may be regularized by adding and subtracting metallic screen Green's function 2Go(a,ql, z)= I + 2; R= - (3.20):2 7rR; R — from (3.19). The free-space dyadic Green's function is given in closed form by (2.21). The resulting regularized Green's function (curvature contribution only) is given by GZZ(aj)Z,) ~ jk qe-ko { (cos20 + q(1 -q)(2 -3cos20)) [() - 1] GOZ (a, ) 2 qe-kossin cos0 (1-3g(1-q)) [v(3)G4Z(a, 27z) J G(a, z) 2 - qeo (sin 20+ q(1 - q)(2 - 3sin2))[v(n)-1] +q [sec20 (u(3) - v(3))] } (3.21) and since it is not singular it may be numerically calculated. When the regularized integrand is used, the contribution of the planar Green's function G i = 2(koa)2 f J Wt(a, X, z; pi,, i,Zi, Si)

33 (a. O. ) x G'(a..-) x p(a. o'.z) *~ I (a. o. ij.. j,~j Z Sj)do' dz'do dz (3.22) is added to (3.21) Gij = Gj + GP (3.23) and this is used instead of (3.17) with R = s. Upon use of a common vector identity and the divergence theorem [16], we obtain from (3.22) P= (ka)2 r,,-' 27 Js t(a,' z; WS(a, Z) z; pji, Oj,,jZ, sj)d z'd dz d)dz 2 -2 J/ V ['(a, (, z) x Wt(a, Oz; Pi, izir Si)] ejkR V [(a,') x Wsa, (a z; pjaj~. Zj,sj)] e-R Pi, O',,]R do) dz'dr) dz (3.24) where I = xx + yy + zz. and this form of the boundary integral may be readily evaluated even as R vanishes by employing the standard regularization procedure (see for example [16] or [41]). This procedure amounts to subdividing the surface into smaller patches and using an analytical evaluation of an additional regularized integrand. We note that v(3) -+ 1 as a -- oo and hence the regularized integrand (3.21) vanishes leaving only the planar contribution (3.24). 3.5 Excitation: FE-BI Two sources of electromagnetic fields are considered: external sources (plane waves) for scattering analysis and internal sources (probe-feed) for antenna parameter calculations. The use of the exact boundary condition in (3.5) allows the coupling

34 of an exterior excitation field into the cavity while the FE fornulation itself readily permits modeling of interior sources. In this section we describe the form of tlie source functionals ftet and fint along with their numerical implementation. The forcing function, due to exterior sources (fet) is given by ~i o o A fs lJ^ (a W 1z ) p(a xz ) x H'(a pz )do' dz (3.25) 0 III cy'(a,6','l. (3. 2) where Pri is the permittivity of the interior element associated with the ith aperture unknown, Wi(p,,,z) is the corresponding testing weight and Hcy' represents the magnetic field on the cylinder's surface in the absence of the cavity. A plane wave Ei = -ie-jko(k'F) [p (sin sin sin + cos 7os cos c i) + (siny cosi - cos 7cosOi sin oi)z cos 7 sin 0i] ejk~[asinO' cs c +cos ] Hi = Yo [p (sin cos Oi cos i - cos 7 sisin cs) - (ssin + cos sin + cos cos)z sin y sin oi] ek[in osi+zcos i ] (3.26) is assumed to be incident on the cylinder from the direction (0i,Oi) where y is the polarization angle and e = BOcosy + Oqsin7y is the electric field polarization. In these, the difference between the observation and incidence angles is denoted by fi = - qi. The total surface field is given by the sum of the incident and the corresponding scattered field from the infinite metallic cylinder [42]. Specifically, Hcyl(a, z) = Hi(a, ), z) + Hyl(a, 4, z) = Hyl +:zHCl (3.27) From traditional modal analysis, we have H (a,, z) -2Yo2, *C, z [E o^s(+ ~~b~ _ 1 ~ sn ~ = __t.....asi~i

35 n sin ) cos, 1 n( +~-) koa sin 0, H'(2)(koa sin 0i) cin 3 oc r, ^-' i -K'-C'i) Cyfy'(, s 9 ~ Cos 0, ean " HI (a,,z) = j2,, k E ['(2) a (3.28).ka, =-o H 2)(koa sin i) These expressions may be approximated by retaining only a few terms of the series if koa sin 0, is small. However, as this parameter becomes large (e.g. for large a and Qi - 90~), (3.28) may be replaced with equivalent asymptotic representations similar to those considered earlier. Utilizing Watson's transformation and Fock theory [42] in connection with (3.28), we find the following direct path creeping wave terms 2 H? 1 -yo siny sin 0iejko c~SZ > e- jkasin ip [g(O)(mq)p)]* p=l o2 2 HCY l j2Yo os y km ekocosiz y e-jkoasin [f( ( ka sin 02; 1: -Y, sin 7 cos Oie0kocoS" (_l)Pe-jkoPini p g(o)(mp) p=1 1' -jm in gw (mn p)] (3.29) k. a sin 0, p in which 1 = 3 ( - i), 42 = (- i) -, m = [k L 3, and the asterisk indicates complex conjugation. The appropriate Fock functions are [40]2 g9')() = _ Wl(tdt iJ r eJet f()(1) = J J e dt (3.30) V/7rrwi(t) where wi(t) denotes the Airy function of the first kind, w (t) is its derivative with respect to the argument and the integration contour (F) is given by Logan [40]. The asymptotic formulas (3.29) are quite accurate except when 4 h 4i. In this region, Goriainov's [43] expressions Hzc l -Yo sin a sin ieJkoosz { e-jkoasini1 [g(~)(ml )]* 2Logan [40] uses the e-jt time dependency in the definition of these functions requiring the complex conjugation in (3.29)

36 +_jkoa sin, cos(o-o,) [G(-mn cO (o - 0o))]' } H~ ~ j2iY Cos iQ- ejko cosi~ -tjkoasin6,l ( 1(0 l) HYl N2 ocos k0a sin O c' Le j +e.kosino, cos(0-~i) [F(-m cos (O - (,))] } +Yo sin a cos ieJkocos z {e-Jkoasinl9 l[g(~)(n(Il) -ka sin i (1)( )] __jkoasinO, cos(k-ki) [G(-_ COS (O - i)) -Jka sin O G(1)(-mcos (0$-E))]J (3.31) with [40] G(')() g() F()() - f(1)()eJ3 (3.32) have been found to be more accurate and can be used instead of (3.29). These surface field expressions may now be used to calculate the entries of the column vector {fiet} via a numerical evaluation of (3.25). In particular, the modal series (3.28) is used when koa sin i < 10 and either (3.29) or (3.31) when kca sin Oi > 10 as appropriate. Having now specified the external source, we present the more compact internal source used for this work. Conformal antenna patches are typically fed by a microstrip line printed along with the radiator on the surface of the substrate. The microstrip lines are in turn fed by a coaxial probe which originates behind the cavity as shown in Figure 3.4. The patch may also be fed directly by the coax feed or through some form of aperture coupling. Nevertheless, the principal excitation for the system is given by (3.2). The

37 Patch or microstrip line.~-, ~ [, ~ / Cavity V 1r. 9^*'^ Metallic cylinder FigureCoax feedconformalpatch antenna. Figure 3.4: Probe-fed conformal patch antenna.

38 impressed excitation current for an infinitesimally thin probe is given as t(p,, z) = 1o6( -of)6( ) for radial probe P Jint(p, z) = (Io6(p - pf)6(z - Zf) for azimuthal probe fnt(,,Z) = o -P - ) for axial probe (3.33) p which results in an excitation function (3.2) =- -jAkZj lhn O) [(f - Pi) (f - Zi)] for radial probe fi~nt - -jkoZolo [(pf - pi) (Zf - Zi)] for azimuthal probe -h [(pf 1i (O-, ) th fint = - j ko ZoI [(pf - P ) (-Xi )] for axial probe (3.34) With both types of excitation and the FE-BI matrix now specified, the BiCG method may be used to determine the unknown electric fields within the cavity. 3.6 System Solution The FE-BI system (3.7) may be solved using one of the popular direct methods such as Gauss-Jordan elimination, Gaussian elimination and LU decomposition [44]. Alternatively an iterative method such as the Gauss-Seidel method [44], the conjugate gradient method (CG) [44] or the biconjugate gradient (BiCG) method [44]-[46] may be used. The symmetric form of the BiCG solver has been chosen since it requires much less memory than a direct method and it can be implemented in a manner which is computationally efficient; i.e. utilizing only one matrix-vector product per iteration. Although use of an iterative method such as the BiCG method can require more wall clock time than a direct solver when multiple right-hand sides are considered, in-core memory is deemed the most critical and expensive resource. The BiCG algorithm is given in Appendix B and the efficient FFT-based calculation

39 of the boundary integral matrix-vector product is discussed in Appendix C. In this section, we will present only details of this iterative solver specific to this application. The BiCG algorithm requires one matrix-vector product per iteration (see the algorithm in Appendix B). This operation represents the bulk of the computational demand of the method and requires 0(N2) complex operations per iteration for fully populated matrices. The matrix-vector product is carried out by executing the sum N y[n] = [A] {x} = E A[n,n']x[n'] n = 1,2,3,...,N (3.35) n/=l and if the matrix is sparse, a storage scheme such as the Compressed Sparse Row (CSR) format may be used to reduce both the memory demand and computational load. The FE matrix [A] in (3.7) is such a sparse matrix. The CSR scheme retains only the non-zero entries of the matrix in one long data vector along with another data vector, known as the offset vector, which contains the number of non-zero entries per row of the matrix. Using the CSR scheme, (3.35) can be rewritten as r[n] y[n] = [A] {x}= A [e(n, n')] x[n'] n = 1,2,3,..., N (3.36) n'=1 where r[n] is the number of non-zero entries per row of the matrix and e(n,n') indicates which entry of the long data vector is associated with the matrix entry A[n, n']. Typically, the rows of the matrix are stored consecutively which results in a counter being used for the mapping function e(n, n'). Although additional memory saving is possible since the FE matrix is symmetric (thus only the upper triangle of the matrix need be stored), experience has shown that use of a symmetric matrixvector product leads to a severe performance degradation on vector computers (such as a CRAY) due to the resulting short vector lengths. Therefore, the entire sparse FEM matrix is stored and used in the product so that the code's performance is maximized when executed on vector architectures.

40 The boundary integral matrix-vector product involves the fully populated matrix. [g]. If uniform surface elements are used in the discretizationi. tis matrix-vector product may be expressed as a truncated discrete linear convolution and is tlhus amenable to efficient calculation using the fast Fourier transform (FFT). Although uniform zoning imposes restrictions on the geometries which can be analyzed by this FE-BI technique, the resulting memory and computational efficiencies have proven to be well worth this sacrifice. The boundary integral product is implemented as described in Appendix C with the following problem specific exception. The cross-term arrays do not possess the property: Gz[m - m, n - n = G,[m - m, n -n] and hence the periodic replication rule used by Jin and Volakis [47] cannot be used here. Replication refers to the padding process necessary in the computation of a linear convolution via a circular convolution. A different replication rule was developed using traditional techniques to perform the discrete convolution which is described in Appendix C. A general discussion of replication is provided in that appendix. Cylindrical-rectangular cavities are bounded by metal on all four lateral sides while wraparound cavities encircle the cylinder as one unbroken groove. These two types of cavities are shown in Figure 5.7. The cross-term replication rule must be changed for wraparound arrays since the cavity discretization is fundamentally different than is the case for cylindricalrectangular cavities. We note that these replication rules are not unique and are implementation dependent. In particular, since the columns of the discretization wraparound on a circular cylinder forming closed loops, the replication rule must account for this periodicity to reduce the computational burden.

41 3.7 Scattering and Radiation Once the cavity aperture and volume electric fields have been determined for either an external excitation (scattering) or an internal excitation (radiation), seseral engineering quantities may be calculated. The aperture fields may be used to determine the radar cross section (RCS) for scattering or the radiated fields for antennas. This entails the integration of the aperture fields with an appropriate Green's function. In addition, the input impedance may be calculated by using the interior cavity fields. In this section we will present the relevant formulae or calculating the far-zone radiated/scattered fields and the input impedance from the electric fields. 3.7.1 Far-field Evaluation Two of the most important applications of the presented formulation deal with the calculation of the cavity's RCS and the radiation pattern due to sources placed within or on the aperture of the cavity. We begin with the integral representation of the scattered magnetic field in terms of the aperture fields. We have, H(r,,) = jYokoa jG2(r, 0,; a,4, z ) [p(a, z) x E(a, X', z)] d(f dz (3.37) with (r,0,q) indicating the observation point in spherical coordinates. When the observation point is very far from the cylinder, the dyadic Green's function in (3.37) can be replaced by its far-zone representation G2(r, a, z; a,') - e-k~ [CG~' + GGz + Gt ) (3.38) where the unprimed unit vectors are functions of the observation position and primed ones are functions of the integration point in (3.37). The components of this far-zone

42 Green's function are determined by a mode matching procedure givilg Go,. j 2 ^jkCOS 6 jkocosS:' 1: ":jn( +(n-o')) (2r)2 (ka sin 0)2 - H'(2 (ka sin 0) Gfz 2' ejkocossl 1 ejn(-+(-o-)) (27 )2 a n=- H2) (ka sin 0) (27r)2 asin 6n=- H2)((kasin0) ( As one might expect, these series converge rather slowly for large ka sin0. They must therefore be recast in another form by employing Watson's transformation and Fock theory as described before. We have, G e47r kcs 1 ( ka sin 0 ( k^ cs0 e8kocosz Z)eioasin P [(0)(mi)sJ p=1 47r Pp=l,. ecos 0 O2 e -koa sin [g()(m0) ( G O O, 42a 7r s in 0 C^ L J-. f ( 0) ( r n ID p) ( 3.4 0 ) ekoa. cos(s [G(L)(-mcos ( - )) -kasin ] kGsin eijkoos Oz' {e -jkoasin9 [g(0)(rml i)p] * + 247r k sin i 0 where theappropriatcoejko s are giv(e) [G(~.(0) cos (e, - he) )]' 47r I si 2 eiko cos Oz' {-jkoa sinl [(~)(m4)11 * + 2a7r sin^ 0 eI c ~jkoasinOcos(s-c') [F(O)(-mcos ( )'))] } (3.41) are more useful. The far-zone scattered or radiated field can be computed numerically by using either (3.39),(3.40), or (3.41) in (3.37).

43 For the scattering problem. the RCS is most often the quantity of interest NwNhich is given by (2.29). Alternatively, the antenna gain may be computed from the far-zone fields as GdB(0,) = 10 log) [47r( )E(0 ] +10 loglo [ZR (3.42) where Acm is the wavelength in centimeters, Rin is the input resistance which is given in the next section and Er is the radiated electric field as r -x oo. In the far-zone, the electric field is related to the magnetic field by EO = -ZoHo E = ZoHe (3.43) 3.7.2 Input Impedance In addition to the antenna pattern and gain, designers are concerned with the input impedance of an antenna for feedline matching purposes. The FE approach allows the calculation of the input impedance of the radiating structure in a rather elegant manner. The input impedance is comprised of two contributions [48] Zin = ZP+ZD (3.44) where the first term is the probe's self-impedance which is the impedance of the probe in the absence of the patch and the second term is the contribution of the patch current to the total input impedance. The self-impedance of a radial probe is typically approximated as [49] Zp = 60 [(ko)2 r+ jk1s- (1-) + r 21- +] (3.45) where 6r is the relative permittivity of the substrate, re is the radius of the probe and I is the length of the probe. For very thin substrates and thin probe wires, this

44 contribution is negligible. Ignoring the probe-feed's self impedance, we have [-1S] Zln = - E(p). (.-) Jt(p. o.:)pdpdodz (3.46) where the impressed current is given by (3.33), VI refers to the volume elements containing the probe-feed, the electric field is the interior field associated with the feed edge and Io is the constant current impressed on the probe. Utilizing (3.9) and (3.33) into (3.46) yields n - - 1' — E(i) -hl -l ()pi (zs - zi)] for radial probes Z = i [(ps - i) ( - for azimuthal probes Z = - E(i)t h [(ps - i) (ZS - ()] for axial probes (3.47) which must be summed over the four parallel edges of the element which contains the feed. In (3.47), the probe location is indicated by (.5, z5), (p, zs) or (ps, Os) for radial, azimuthal or axial probes, respectively.

CHAPTER IV Finite Element-Absorbing Boundary Condition Method In a previous chapter, the Finite Element-Boundary Integral (FE-BI) Method suitable for cavity-backed antenna analysis was introduced for cavities which are recessed in an infinite, metallic cylinder. The interior fields of the cavity are determined via the FE method while the aperture fields are determined using a boundary integral in conjunction with the FE method. The boundary integral provided an exact mesh closure condition; however, use of such a global condition partially destroys the sparsity of the FE-BI system. In fact, as the surface area of the aperture increases, the poor scaling property of the boundary integral dominates the computational and memory demand of this methodology. In addition, the FE-BI method is ill-suited for use with coated conformal antennas or protruding antennas. Therefore, it is advantageous to seek an alternative mesh closure condition which preserves the inherent sparsity of the FE system. One approximate closure condition utilizes an absorbing boundary condition (ABC). This condition provides an approximate relationship between the tangential electric and magnetic fields on the boundary of the computational domain. Wilcox [50] introduced a spherical wave expansion of these fields which is the typical ap45

46 proach taken for deriving ABCs. A major disadvantage of two popular ABCs for the vector wave equation, which were proposed by Peterson [22] and \Webb and KaInnelopolous [23], is the requirement that the mesh boundary be spherical. For most geometries, including cylindrical-rectangular cavities, the use of a spherical mesh closure results in an excessive number of unknowns. Consequently, conformal ABCs which minimize the number of unknowns expended in the region between the geometry and the outer boundary have been actively sought. One conformal ABC which has proven remarkably accurate for scattering applications was recently proposed by Chatterjee and Volakis [24]. This ABC uses a second-order expansion of the local fields on the boundary to relate the electric and magnetic fields and it explicitly includes the curvature of the boundary in its formulation. In this chapter, the FE method will be married with these new ABCs to form a FE-ABC method suitable for cavity-backed antennas which are recessed in a circular cylinder. The complete formulation requires development of the FE equations along with the specification of the ABC for a cylindrical-rectangular box boundary. An important aspect of this development is the use of a domain decomposition which employs a scattered field formulation in the exterior region while utilizing a total field formulation within the cavity. Such an approach minimizes the number of unknowns within the cavity while removing the incident and reflected fields from the exterior region unknowns. The ABC presented by Chatterjee and Volakis [24] absorbs only the scattered field in the exterior region (see Figure 4.1) and due to known error propagation effects associated with a total field formulation, a scattered field formulation is preferred in the exterior region. For radiation analysis, such a decomposition is unnecessary since the radiated field is identical to the total field. Previously, the FE-ABC method had only been used for scattering analysis of finite

47 bodies (see for example Chatterjee et al. [21]). A primary goal of this chapter is to present the first application of this method to infinite bodies as well as the first use of ABCs for antenna parameter calculations. We begin with the formulation of the FE-ABC equations. 4.1 Formulation Consider the computational domain shown in Figure 4.1. There are two volume --- -.Sat Patches Region I Composite skin *... ~.... I Composite skin Composite skin;... ":.. /. *. \\>..............................................,'.C o m p o s ite.sk i n Metal Figure 4.1: Typical coated cavity-backed patch antenna with ABC mesh termination. regions: an exterior region, V', which includes any radome overlay and an interior region, V". Both regions may be inhomogeneous and are separated by the aperture surface, Sap, and the surface metallization surface, S1m both of which lie on the surface of the metallic cylinder (p = a). Thus, the exterior region is defined by p > a while the interior region has p < a. The computational domain is bounded by the union of the metallic surface, Smetal = Ssm + Scm where SCm is the metallic walls of the cavity and the ABC surface, Sabc. Within the computational volume, the total electric fields may be written as E(r) = EI(j?+)EC^Y(j? rEV'

48 = E"() 7-rC (4.1) where Ecy(r) = E~i(r) + Er(r') as before. The total magnetic fields are likewise written H() = HI(r+HCl(r^) FrcV -H1( r^ V "r (4.2) where HCY'(r) = Hi(fr + H(r'). The boundary conditions are readily written in terms of the electric and magnetic fields. Within the cavity, the tangential electric field vanishes on the metallic walls n x E"(r') = 0'E Scm (4.3) while on the aperture, the total tangential fields are continuous fix ( x E() E= ixE e sap n x H(r) = fn H(r)-n x H ('(F) re Sp (4.4) On metallic surfaces, all tangential electric fields vanish, i.e. n x ECl( = -) - n x I( E() = 0 r E Ssm (4.5) while n x Ec1 (r) also vanishes over the aperture nLx ECYl(r) = 0 E Sap (4.6) since it contains both the incident and reflected fields. Thus, the only non-zero electric fields on the surface of the metallic cylinder correspond to the unknown fields within each region which are continuous across the surface aperture as implied by (4.4).

49 The FE equations may be developed by considering the inhomogeneous vector wave equation (2.16). Employing the method of weighted residuals and Green's first vector identity, the weak form of the vector wave equation (3.3) is obtained where the interior source functional, fint, is once again given by (3.2). A domain decomposition is accomplished by substituting the total field relationships (4.1) and (4.2) into (3.3) producing the FE-ABC equation jX Vx > Vx _ EE- W~i dV + X [V x EClU.V x Wi 2 dV (dS jkZ c (" x Hcl) W dS + jkoZo [x (HI + Hc)' - jkoZof (P x Hi") i dS = Z/nt # (4.7) However, due to the boundary conditions at the surface of the cylinder, (4.4) and (4.5), the two surface integrals over the aperture cancel one another and with some further manipulation, (4.7) can be written as [7 X El V] dV + I;' [VxP1l Vx W _ ^2/ErE/l Wi] dV + J (n x V x E') Wi dS = [f + jkZ JSb (n x Hc' ) Wi dS - [ V x EY.V x - _ k2 X c W] d (4.8) the ABC surface in terms of the curl of the electric field. The second term on the right-hand side of (4.8) serves the same role as the exterior source functional, fet,

50 given previously for the FE-BI formulation in (3.25): however, as will be shown. thils term will vanish when combined with the succeeding volume contribution. The last term on the right-hand side is a source functional which is present due to the use of a scattered field formulation in the exterior region, V1'. As written in (4.8). this final term must be computed over the entire exterior region. However, by applying Green's first vector identity, this final term may be written as fV J W [V (7 X E) 2 - EcYl] dV + fI Wi x E cy] dS (4.9) The quantity within the brackets in the volume integral of (4.9) is recognized as the vector wave equation and since ECYl already satisfies the free-space wave equation, this contribution is non-zero only for material (cr 57 1 and/or Pur:7 1) regions. The surface integral is the boundary term which is non-zero only over the ABC surface, the aperture surface and the interface between two materials in the exterior region which have different permeabilities. Thus, (4.9) may be rewritten as = I- i> [V (V X EYl) - -. dV + JY i / Wi' V x -k-kEc'l Wi dV + koZo (n x H ). w ds + x H ). d + quantity f o r thrlre i Tri ILre]SId xn * [Wi XV XECY] dS (4.10) where Stri is associated with the interior material layer and Pre is the corresponding quantity for the exterior material layer. The third surface integral need only be computed for each dissimilar material interface in the exterior region. Utilizing the source-free vector wave equation and assuming that material properties are constant

51 within an element. (4.10) may be included into (4.8) to get the final FE-ABC equatiol /,, -— k~ CrE 1V + x. [V x~ Xl dX~ ~ /vi~i[ x E V x 11; _ k,2CF" l 1-I dl + ^/r. X -6. X EI IV -^,.-t - /y) k, x dS kJsab ( v Ecyl) Wi.dVS =- fkoZo1 ] ( s x i') Wi'x'] ddS[- - W dV re] J[d [i (4.11) The source contribution associated with the ABC surface has been canceled by the addition of (4.10) and the remaining terms of the right-hand side of (4.11) are consistent with scattered field formulations. The aperture contribution is analogous to the exterior source functional used for the FE-BI formulation with the exception that the permeability is now associated with the finite element on the exterior side of the aperture. Implementation of (4.11) will give superior results to (4.8) due to numerical error propagation which is inevitable when the volume source term is evaluated throughout the exterior region as in (4.8). In addition, since (4.8) need be computed over the entire exterior region rather than over the presumably smaller material region, it would require substantially more computational effort than is required for (4.11). This later set of FE-ABC equations may be written as a linear system of equations e-abc] { } (4.12) { Ent} {f int where the FE matrix [Afe-abc] may be written as a sum of the FE matrix used in the FE-BI formulation (3.10) and a second term attributed to the ABC surface [Afe-abc] = [A4] + [AabC4] (4.13)

52 The new FE-ABC equation (4.11) is comparable to (3.3) except that the latter utilizes a total field formulation throughout the computational domlainI. However. previously we utilized an integral expression for the total magnetic field across the aperture which resulted in the FE-BI equation (3.5). Such an integral expression provides an exact relationship between the total tangential electric and magnetic fields over the aperture surface which also formed the computational domain boundary. Alternatively, we may employ an approximate relationship between these two fields with the goal of retaining the sparsity of the resulting linear system. Additionally, as shown in (4.11), this FE-ABC formulation may be used for coated as well as uncoated geometries. In the next section, we will develop an approximate relationship suitable for mesh closure. 4.2 Conformal ABCs Traditional three-dimensional vector ABCs [22, 23] require a spherical outer boundary which results in an excessive number of unknowns. New conformal ABCs have recently been proposed by Chatterjee and Volakis [24] which have an outer boundary that follows the contour of the enclosed geometry resulting in a minimal number of unknowns. In this section, the specific expressions required by this new ABC for a cylindrical-rectangular box boundary will be derived. A definition of ABC order will be given and subsequently the first- and second-order ABC expressions will be presented. For the purposes of discussion, we define a secondary field as the field which is a consequence of equivalent currents that are supported by some primary source which is either external or internal to the computational domain. Thus for scattering problems, the scattered field is the secondary field while the incident and reflected

53 fields are considered primary fields. Likewise. for a radiation problem, the radiated field is the secondary field whereas the source field due to an impressed current is the primary field. In (4.11), we recognize that an ABC must supply a relationship between the tangential components of magnetic and electric secondary fields on the absorbing boundary, Sabc. The secondary field may be expressed as a Wilcox expansion Es(nt,t2) = 1/ lim E (4.14) 47rfipU Pi/ P where u = v/RR2, Ri = Pi + n and Pi is a principal radius of curvature. In this form, the curvature of the non-spherical wavefront is explicitly used. The point of observation is given in Dupin coordinates [31] as x = n + o(t1, t2) (4.15) where n is the unit normal and xo(t1, t2) denotes the surface of the reference phase front and therefore, tl and t2 denote tangential coordinates on that surface. Absorbing boundary conditions annihilate outward propagating waves up to a certain order. A zeroth-order (P = 0) ABC represents the usual Sommerfeld radiation condition (2.15). A first-order ABC (P = 1) annihilates all fields with up to a u-1 dependency while all higher order fields are reflected back into the computational domain. For a cylindrical surface, u = Vpf, thus the zeroth-order ABC is simply the geometrical optics spread factor while the first order ABC annihilates fields up to 0(p-1). Evidently, as the ABC order increases, the reflected fields have an increasingly higher attenuation factor and hence the boundary may be placed closer to the geometry without inducing erroneous reflections. In the following two sections, we will present the first- and second-order conformal ABCs attributed to Chatterjee and Volakis[24]. In particular, the appropriate

54 expressions for a cylindrical-rectangular box boundary will be given. The matrix cIltries which enforce these ABCs will be presented assuming the use of the edge-based weight functions (3.9). 4.2.1 First-order ABC The first-order ABC is similar to the impedance boundary condition attributed to Rytov [51]. Specification of the ABC requires the relationship between the tangential components of the electric field and its curl. For the first-order ABC, we have x V x Es = jko +m m-7] E (4.16) where Et = tlE1 + t2Et2, rm =, = = Iiltlti + K2t2t2 and h1,2 are the principal curvatures of the surface. These curvatures are given by [31]' = -h (4.17) where a denotes a normal derivative and hi are the appropriate metric coefficients. The ABC formulae are going to differ for each surface of the cylindrical-rectangular box. There is a singly curved surface which has n = p with p being constant and a total of four flat surfaces having either n = ~-f with q being constant or n = ~z with z being constant. Since the curvature and tangential edge orientations are different for these three cases, each surface will be examined separately. For singly curved surfaces, the normal unit vector is radially directed and the tangential edges which form the surface are >- and i-directed. Accordingly, (4.16) may be written as x V x E = jko E= +- jko-2 \E:z (4.18) since t1 = X, t2 = z, 1 = -- and "S2 = 0. The matrix entry for the ABC surface is

55 given by Ab = \ab x v x I] * I I dS (4.19) (4.20) where the subscripts denote the orientation of the test and source edges, respectively, and the edge-based weights (3.9) are used for Wi and Wj in (4.19). For the flat surfaces which have fn = ~+, the ABC (4.16) is given by x V x E = jko [pE +zEz] (4.21) since the tangential edges are oriented in the p- and i-directions. The corresponding matrix entries (4.19) are abc = jk0 (h) [ - ] [3 ( 3 - z) - (Zi + ) (z - z) + zZh'k o v/ + SPa - - abc = 3 s (t) [3 (p -p) - (Pi + P,) (pb -- ) + pipjt] (4.22) For the flat surfaces which have n = ~i, the ABC (4.16) is given by x Vx xEs = jk [pE + E] (4.23) since the tangential edges are oriented in the p- and +-directions. The corresponding matrix entries (4.19) are abc = jk (-) n [3] 3 Z2.) i ~h2 _( 4243 2 A~ i s () [ab (p -P) - (Pi + Pi) (P -_ P) + \P p _ p) (4.24)

56 These ABC matrix contributions may be assembled il tile usual fashlion alnd added to the FE matrix (3.10) previously used with the FE-BI formulation. Since the preceding ABC has only first-order accuracy, the mesh closure must often be placed far from the cavity aperture which results in a large number of unknowns. To reduce the computational burden of the method, second-order conformal ABCs need to be explored as presented in the next section. 4.2.2 Second-order ABC The conformal second-order ABC is necessarily more complex than the first-order ABC. The tangential component of the curl may be approximated as n x Vx Es = E +' Vx [i (. V )] + Vt (n * E) (4.25) where Vt denotes the tangential surface gradient operator [31]. Unfortunately, use of (4.25) would result in an asymmetric system [Aabc] due to the last term which possesses only one differential operator. An asymmetric system requires an iterative solver which utilizes two vector-matrix products such as the conjugate gradient squared (CSG) solver presented in Appendix B which is based upon the work of Vorst [52]. A symmetric system, as seen previously in Chapter III, requires only one matrix-vector product if the BiCG solver is used. Additionally, for symmetric systems, only the upper or lower triangle of the matrix need be computed and stored. Chatterjee and Volakis [24] realized that the gradient in (4.25) may be approximated by Vt(V E) = jkoV (n E + (n-5) (4.26) With both a gradient and a divergence operator present, one operator can be transferred to the test vector while the other may remain with the source vector. Hence,

the resulting matrix may be symmetric since both the test and source fields are differentiated. With (4.26). (4.25) may be written n x V x = E+ 3 V7 x [n (n V x E )]+ t (p. E) (4.27) For the basis vectors (3.9) considered in this thesis, V. Es is always zero on Sabc and hence the third term of (4.27) will not contribute to this form of the ABC. For surfaces with a common constant curvature for both tangential directions on a surface, this new vector ABC (4.27) will lead to a symmetric FE system [Aabc] as shown by Chatterjee [25]. However, if the principal curvatures on a surface are unequal, the system will be asymmetric. For either (4.25) or (4.27), the three coefficient dyads are given by x4 a D [ - + D (jk -i) + mK2 } = 2{ D-Ai-2t }i D - AnK - 2K;i i= 2 1 D- -2i=l{( \<2) [iko + 3Um - -2K,] tti} (4.28) where Kg = 1/K2, A/N = /1 - /K2 and D = jko + 51i - $-. In the case of (4.27), - must be divided by jko due to (4.26). As was the case for the first-order ABC, it is advantageous to consider the secondorder ABC for singly curved and flat surfaces separately. For a singly curved surface, the unit normal direction and curvature parameters are n = P, I1 = - and K2 = 0 as before. After some manipulation, we find that (4.28) becomes 2kop + D (jko i + l [1 -o[+ +jkoD ] = ^2p 2p j2kop+ 1 +2kj2kopp9 5 = [~>~t>+zz] (4.29)

58 where D = jk. - - Note that /3 is not symmetric unless p -- x For the second-order ABC, it is advantageous to segment the matrix entry (4.19) into three parts. After some vector manipulation [25]. these contributions are given by i''abc~ = f,.ab [n. Vx x] [W n* V x (B T,)] dS; I(3)abc = I [ [* vt (n.W dS (4.30) surface, we find (2-)abc [2Pb 1 ] (PbA A z ~ [j2kz p 3] s- ) (2)abc [j2kopb + 1] zz [ j2kopb - 3] cP/ I(2)abc - j2opb + 1 Sisj h I(l)abc __3 3 - hi 2 2 ah () = ( ) [1 (Z-3) - (Zi +2 Z) (2 - ) + ZiZ I (3)abc= ij ( ) [b (q -SS ) - 2 (. +,3) (q - 4) + ~iA2] (4.31) where the first subscript denotes the test edge direction while the second subscript indicates the source edge direction. Note that the cross terms of the second contribution has I(2)abc I()bc which results in an asymmetric FE system regardless of which form of the ABC, (4.25) or (4.27), is used. This asymmetry is a direct consequence of the different principal radii of curvature for a cylindrical surface. However,

59 these terms are asymptotically identical as the radius of the ABC surface beconles large since the surface will then be approximately planar. A symmetric ABIC ina be obtained by dividing both the numerator and the denominator of (4.31) by the ABC radius, pb which results in 1(2)abc - - [i0 I(2)abc _ 2 l (4.32) 124> ~ ko S (ko s2sj To determine whether the symmetric evaluations given above may be used in place of the asymmetric terms, the RCS was computed using all three forms of the ABC for a 0.6Ao x 0.6Ao cavity recessed in a 0.2Ao metallic cylinder. The cavity if filled with a material which has Or = 2. The ABC surface had a radius of 0.5Ao and the lateral walls were place 0.3A0 away from the cavity sides. Thus, the ABC surface is placed at least 0.3A0 away from the cavity aperture which is shown to be sufficient in Chapter V. The computed RCS for all three ABCs were within 0.01 dB of each other and 0.5 dB of the value computed by the FE-BI method. This comparison shows that the symmetric form of the ABC may be used in place of the asymmetric ABC in order to reduce the computational and storage demand of this FE-ABC method even for a highly curved surface. For the flat surfaces (n = or n = ~i), we find that the first contribution to the ABC entry, I(1)abc, is identical to the first-order case, (4.21) or (4.23), since a = jko [ttl1 + t22] (4.33) where tl =-, t2 = p if = ~fo or t= -, t2 = ( if n = ~i. However, for the second-order ABC, there are additional contributions due to = — [iltl + it2t2] - = jko3 (4.34)

60 If n = ~0, these additional contributions are given by 2)abc _- 1- (Pb [Pb1 I(2)abc - S( ) - - PP ~(ko) (k h Pa (2)abc -_ Pb n(Pb) 3kot) Pa I(2)abc - 1iS (h zP' t= ] [)o ZZ jk- 7 I(X) = Si^i (Sj) [3 -Pb Pa) - 2 (pi + Pj) (Pb - Pa) + Pi t] (4.35) However, if n = ~i, these additional contributions are given by 3(2)abc _ i5J [() l.., rc, (z3lnfz)- f+-i r(2)ai[c In () [1 1 ) I(3)abc Pj ( 3 3 l )- - )] P th t i j t (2)abc = ((^ 2 (Pa)-Pi()l-n1) ] 1 kt2) [2 (pP PPJ)a) 12 PJn 1)b = ^ ( ) [ ( - <) -) (4 + i) (r -s b) + i j] 4?3)ab = c (-) f[ (P- P) - (Pi + Pi) (P2 - Pa) + - t] (4.36) Although it appears that 1(,,)abc 7 I42)bc for n = ~i, in fact these terms will be identical when evaluated on the ABC surface since pi = Pj = pa. We have now fully specified these adittance matributions are given both first- and secondorder ABC mesh closure conditions. When coupled with the previously given FE matrix, the resulting system retains (N) sparsity regardless of the problem size or composition. Furthermore, coated geometries may be readily modeled with this PPb I I(a)abc ~ 1 f a _2 21 P I 22 j + +t Pjt ] (4.36 identical when evaluated on the ABC surface since pi - pj 3 p,. order ABC mesh closure conditions. When coupled with the previously given FE

61 formulation since the FE method permits inhomogeneous materials. However. the source terms given previously as the right-hand side of (4.11) must still be specified. 4.3 Excitation: FE-ABC In the previous section, vector ABCs were presented along with the associated matrix entries required to implement the FE-ABC formulation (4.11). AWe can noN generate the FE matrix, [Af e-abc], which is fully sparse. In order to find the scattered or radiated field, we need to develop the excitation function (or right-hand side) of (4.11). The first term, fnt, is the interior current probe source given previously by (3.2). The second term is analogous to the external soure functional, fe't, used for the FE-BI formulation(3.25); however, this term now contains the permeability of the material exterior to the cavity aperture rather than that of the cavity fill. The final term is a consequence of the scattered field formulation in the exterior region. It vanishes for free-space since the incident and refected fields satisfy the vector wave equation (2.16). For coated geometries, this term will not vanish and must be retained for that case as discussed later. 4.3.1 Aperture Term Previously, the FE-BI formulation required the external excitation function given by (3.25). In this, the magnetic field attributed to the incident and reflected fields is evaluated on the mesh boundary. Since the mesh boundary coincided with the cylinder surface, considerable simplification was possible in generating these magnetic fields (see (3.28),(3.29) and (3.31)). However, for an FE-ABC formulation, the mesh boundary does not correspond to the cylinder surface although as shown above, a source term may be attributed to the aperture. Therefore, in the absence of any radome covering, the exterior source term is similar to the one used for the FE-BI

62 formulation (3.25). The forcing function due to exterior sources is given by fext(abc) koaZo s -'''''' d~' f ext(abc) = *k(aZ0? fc W2(a, ( ) q z (fa, ), ) x H Yl(a, o' ) do' dz (4.37) where the integration surface, St, is associated with the ith unknown which lies on the cavity aperture and HcY1 is the sum of the incident and reflected fields. Note that (4.37) contains the permeability of the exterior material (tre) while the corresponding function for the FE-BI formulation (3.25) contains the interior material (or cavity fill) permeability (/Lri). This is a consequence of applying Green's first vector identity to either the exterior or interior volume regions, respectively. Thus, the same aperture excitation term may be used for the FE-ABC and FE-BI formulations. 4.3.2 Material Term As previously mentioned, if a material is present in the exterior region, an additional volume and surface source terms need be included in (4.12). For example, these terms are necessary for a printed antenna with a radome covering. Such terms are a consequence of using a scattered field formulation in the exterior region. Typically, a radome is a complex, inhomogeneous layered structure which exhibit desired transmission and reflection characteristics within a design bandwidth. Traditional method of moments formulations can be used to model such complex layered overlays; however, such an approach requires an extremely complex and costly dyadic Green's function or alternatively a full volume formulation, either of which is associated with a fully populated matrix. Since the exterior region in this FE-ABC formulation is modeled using finite elements, radomes or any material covering may have arbitrary composition without the necessity of a complex Green's function. In fact, the same algorithm may be used for both coated and uncoated bodies by simply

63 including the material source terms which are given next. The volume contribution may be recognized as the third termn on the right-hand side of (4.11) fV = -ko [1 -r] lId 4i E Y' dV (4.38) This contribution is associated with the volumetric electric field and using the cylindrical shell elements (3.9) presented previously, this term is given by fip = -kl [ r - h J J (j U- ) (z - i) Ep (p, 6, z) dpdcdz fr.OhJb J4t Jko vr ceh / Ir _ In Jzb pa, oIrrth zb JI a - o rth ] ] (P -pi (" - i) ) Eq (P. 4, z)pdpd~dz f[7 = _k2c - Pa (p — p,) EZ' (p, ), z)pdpdqpdz _ko2 [ P.r r at, J I a'' (4.39) where the subscripts indicate the orientation of the ith edge. It remains to determine the electric field anywhere in the exterior region which will require traditional modal analysis. An incident plane wave (also presented in (3.26)) is given by Ei eie-jko(k'. = [ (sin 7 sin i + co os O c os cosi) + (sin y cosi - cos cos Osin i)z cos 7 sin ] ej[sin cos$+z cos s] H. = Yo [p (sinc cos Oi cos - cos y sin -i) - (sin 7 cos Oi sin <i + cos 7 cos i) - z sin 7 sin Oi] ejko[, sin cos4i+zcos&] which is assumed to be incident on the cylinder from the direction (Oi,qi) and y is the polarization angle where e = cos yO + sin yq is the electric field polarization.

64 In these expressions, the difference between the observation and incidence angles is denoted by o,- = - i. The electric field associated with the incident and reflected fields is determined by traditional modal analysis which enforces a vanishing tangential electric field on the surface of the cylinder (p = a). This field is given by 00 ECYI - cosa sinOieiko^sGiz co2ose0n+ski) (4 +40) Ec? e= J (k cosi - cos (k+ sin H 2(kosin Oip i( = Jn (ksin )) H2) (k3 sin p) noo - (k" sin Ecyl = - cos sin iejkOc~z C Q,,en(C+') (4.40) where the modal coefficients are given by (k Sin - - Jn (k sin 0a) H(2) (ko sin Oip) Hn( (ko sin Oia) In these, J () and YJn () are Bessel functions of the first and( second kind, respecJ'm (ko sin a (a)) sod not fields (are gvn yina) - n = Jn (Yoe sin J [ s (incos 0 ) H(2)co (kosinO p) Ye= k seZ [sin ycosH k sin a jcosYo = J(k (, sin - p) (4.41) Hn' (c, sin 0 a) In these, Jn (*) and Yn (-) are Bessel functions of the first and second kind, respectively. Primed expressions in (4.41) denote differentiation with respect to the argument (e.g. J' (() - 7, (()). Note that the polarization angle (7) should not be confused with the modal coefficient (77). For reference, the corresponding magnetic fields are given by Hfyl -Yo0 c j sin 7co cos 0 - co si O1p — n e2 n n-= - 00 sn p HCl = Yoe siny sin iek cos7iz O yinC+s k o) sin p e(4 ) H' = -Yosinsiniek ~8~ neJn(+n)o (4.42) n= — oo

65 If (4.42) is evaluated on the surface of the cylinder (p = a). a, 3, = 0 as required since the normal component of the magnetic field (or the tangential compoilent of the electric field) must vanish on a metallic surface and (4.42) reduces to (3.28). These doubly infinite series may be converted into singly infinite series which are more amenable to numerical implementation Ecy = _ekoO { cos 7 cos 0osi + E = _-ejkocos~Z{jsin Ty/o - 2 Z jn+l c Co7 c CoOi k7 a-n sin (n4) - sin n73 cos (ni n=1 [ ksi 0ip Ji ) EY = cos y sin iejkocOSz o + 2 > j'na cos (nq) } (4.43) Hcyl = _yoejkcos0Z{j sin 7cos i3o + 2 jn+l [sin 7 cos 0n cos (ni) - cos 7 k i sin sin(n i)] } HY = yoejkOcOs, {jcosz 7jcos + 2 jn+l1 [sin 7 cos 0i- i sn sin (ni)) +cos cos ) Hl = -Yo sin 7 sin ieok~ cosZ 7 + 2 Z 7 j cos (nqi) (4.44) n=1 All of the formulae given above require numerical evaluation of the integrals which may be performed with a mid-point rule. At each integration point, the magnetic field is determined by truncating the infinite series in (4.43) so that the field value satisfies the convergence criterion lE (n) - Ey(n - 0.00001 (4.45)< n P~ ~ ~ ~ ~ ~ ~ ~ ~ ~~~~~~~~~~r,

66 where u C {pO, z}. Since these series exhibit poor convergence for large radius cylinders, the modal formulae (4.43) are only suitable for relatively small radius cylinders (i.e. k osin Oa < 50). For larger radius cylinders. asymptotic expressions similar to (3.29) and (3.31) are more practical. When there is a discontinuity in permeability between material lavers, a surface term must also be included on the right-had side of (4.11). That term is given by fsd = jZZk [- - n-1 I. [ x Hyl] dS (4.46) [Mri Mre ] S' L 4 and may be evaluated in a similar manner as the ABC surface term (4.37). After some manipulation, (4.46) is given by Sd S i 1 zt J k =- -jkopbZo -( - - -) H' (b,, z, z) d>dz I I p,. zt r fSd = jkopbZo- [1 - i H-' (pb, H, z) dqdz (4.47) O [ire tri JJ assuming that the curved integration surface has n = p and that the cylindrical shell shape functions (3.9) are used for field representation. For an interface surface which has = -n, the functional (4.46) is evaluated as Sd _ iPb 1 1 zt Pb Z MY ( ip = jko Zo-,- -l - I4 - H (p, 4),z)dpdz ire Pri IJ Zb Pa P fis = -jkoZo ] jl (P - p i)HY' (p, ~I, z) dpdz (4.48) while if n = (, (4.46) is given by f 3kZ h [lre / 1 I A HcY1 (p, Or, z) dpdz n h _re pri J zb a P fid = -jk3Zo i: (P - p i)HpYH (p, )zr, z) dpdz (4.49) For an interface surface which has n = -, the functional (4.46) is evaluated as fsd jkZ0 /Z-I (4)- ) Hf) (p, 4) zb) dpd1 Oa [lire lri J'P J Pa Sd _ Si 1 ii tr [Pb fi = jkoZo - - ] P( - p-3i)H Y' (p,, Zb) pdpd4 (4.50) t Li/re PLri] J~kl I Pa

67 while if n = 5. (4.46) is given by fSd _SPb [i up = -jkZ o Ho (p, o[. ot) dpdo a -orf' [ rn. ]/ a rSd 1i i - = jkZ I - I - Z(p- P)H (p. i.z) pdpdo (4.51)'t tire r i I a P H The magnetic fields are computed using (4.44) -where the series are truncated according to (4.45). With the specification of the surface and volume source functionals, the FE-ABC methodology is now fully developed for a full range of typical applications. For radiation problems, since the total field is equal to the radiated field, the FE-ABC formulation is quite similar to the FE-BI formulation with the exception of replacing the exact mesh closure condition on the surface of the cylinder with an approximate condition at some distance above to the cylinder's surface. However, for scattering analysis, in addition to the new mesh closure surface/condition, the excitation term may include radome contributions. Thus for coated antennas, scattering analysis requires additional volume and surface source terms which are a consequence of using a scattered field formulation in the exterior region.

CHAPTER V Scattering One of the principal goals of this research was to conduct studies of the scattering by cavity-backed antennas recessed in a circular cylinder. Over the course of some forty years of radar cross section (RCS) investigations, the electromagnetics community has learned techniques which suppress traditionally dominant scattering contributions such as reflected and creeping waves associated with curved surfaces. Having reduced these bulk scattering features through the use of exotic materials and platform shaping, previously insignificant scattering contributions have become increasingly important. Scattering by an array of cavity-backed antenna elements is one such secondary contribution. In the past, most conformal antenna array designs were performed using either rigorous numerical formulations involving an infinite ground plane or approximate techniques such as the cavity model. Planar solutions are suitable for near broadside scattering provided the cylinder radius is large. However, since no circulating waves are supported by a metallic plane, bistatic scattering analysis in the penumbra and shadow regions of a cylinder cannot be obtained from planar results. Additionally, creeping wave coupling creates a significant scattered field on the rear-side of the cylinder even for backscatter measurements. Modern military aircraft demand low 68

69 observables over a full range of aspect angles necessitating the inclusion of creeping wave effects in any useful analysis technique. Approximate cavity models are useful for radiation analysis at resonance: however, they are unsuitable for scattering applications. The two hybrid finite element techniques presented in this thesis provide a means of performing, in the case of the FE-BI method, a rigorous analysis of scattering by cavity-backed antenna elements which are recessed in an infinite, metallic cylinder. The FE-ABC method provides an approximate but more flexible technique which may include the presence of a radome. In this chapter, the FE-BI method will be validated by comparison with data generated by different numerical models. This method will be shown to be highly accurate and, as implemented via a BiCG-FFT solver, highly efficient in terms of computational and memory demand. The FEBI method will then be used to study the scattering behavior of typical conformal antenna elements. In particular, the effect of cavity size and curvature will be investigated. The results will suggest useful design concepts which will aid RCS engineers in developing low observable antenna arrays. The FE-ABC formulation presented in this thesis will be compared to the results generated by the FE-BI formulation. Their agreement will confirm, for the first time, that ABC termination schemes may be used above infinite, metallic structures. Previously, all FE-ABC applications involved closed, finite bodies. The FE-ABC approach offers considerable advantages over the FE-BI method since it may be efficiently used for computing the scattering by conformal antennas recessed in doubly curved platforms. Thus, it is essential that this technique be validated with a rigorous solution, before its application to the necessarily more complex doubly curved structures. Additionally, the FE-ABC method provides a means of readily including

70 a complex inhomogeneous radome which is oinevitably placed above a conforial element. This approach also permits the modeling of a radar absorbing material (1RAM) coating which is typically placed around a transparent radome. Such a RAM layer is used to suppress creeping waves and, therefore, it is an important feature in low observable designs and it is also useful in suppressing Electromagnetic Interference (EMI). 5.1 FE-BI: Validation Having solved for the electric fields induced by an incident plane wave, the computed RCS data must be validated with known results. Available measured or computed data are rather scarce and as a consequence, we are forced to rely on limiting cases in order to validate this work. For example, as the radius of curvature decreases, a cylindrical-rectangular cavity will approximate a planar-rectangular cavity. Another limiting case involves comparison of an elongated 3-D cavity with a corresponding 2-D cavity for normal incidence (Oi = 90~). Finally, the infinite cylinder results are compared with a finite Body of Revolution (BOR) model for oblique angles of incidence. The first validation effort for scattering by cavity-backed patch antennas relies on the fact that a small patch on a very large radius cylinder is quasi-planar and approximates rather well a planar-rectangular patch. For our test we chose as a reference a planar 3.678 cm x 2.75 cm patch residing on a 7.34 cm x 5.334 cm x 0.1448 cm cavity filled with a dielectric having cr = 4. The equivalent patch on a 32.6 cm cylinder is 6.46~x 2.75 cm residing on a 12.90~x 5.334 cm x 0.1448 cm cavity. At the operating frequency of 9.2 GHz, the cylinder has an electrical radius of 10o. Figure 5.1 shows the results for the patch on a large radius cylinder with correspond

71 2 0.0.....I.....I......I........ 10.0 -. -10.0,, -0.0 M3: /';*' ---- Cylinder: E-pol 0 " J.^ \ / Q~C0 Planar: E-pol 13 -40-0, - o13 Planar: H-pol -40.0 - 0.0 10.0 20.0 30.0 40.0 50.0 60.0 70.0 80.0 90.0 Angle (0) [deg] Figure 5.1: Comparison of RCS for a planar patch (3.678 cm x 2.75 cm) residing on a 7.34 cm x 5.334 cm x 0.1448 cm cavity filled with cr = 4 dielectric and a corresponding quasi-planar patch on a large radius (lOo) cylinder. a corresponding quasi-planar patch on a large radius (10X,o) cylinder.

72 ing data for the planar cavity-backed patch. Clearly, the two RCS patterns are ill excellent agreement. and although Figure 5.1 illustrates only monostatic scattering in the 4 = 00 plane, additional runs for normally incident monostatic scattering and various bistatic situations yield similar agreement. Comparisons may also be made for elongated cavities and 2-D moment method results. Long narrow cavities have very little axial interaction for principal plane (0 = 900) incidence and therefore results based on this formulation should compare well with corresponding 2-D data. It is well known that the RCS of the 3-D scattering body of length L > AO is related to the corresponding 2-D scattering of the same physical cross section via the relation U3D = 2 (-) a2D (5.1) Such a comparison is shown in Figure 5.2 for monostatic scattering by a 45~ x 5Ao x 0.1A0 cavity recessed in a cylinder with a radius of 1Ao for both principal polarizations. Once again the agreement between the two results is excellent, thus providing a partial validation of the formulation for highly curved geometries. Similar agreement has been observed for bistatic scattering in the 0 = 90~ plane. The planar approximation eliminates the effects of curvature, which is of primary interest in this work, and the 2-D comparisons done above are only valid for normal incidence. To consider oblique incidence on a highly curved structure, comparisons with a Body of Revolution (BOR) code are made for wraparound cavities. Since the BOR code can only model finite structures, an infinite cylinder is simulated by coherently subtracting the far-zone fields of the finite structure without a cavity from similar data which includes the cavity. Such an procedure mimics common measurement practices and was found suitable for near normal incidence in the case of H-polarization (y = 90~). This procedure is illustrated in Figure 5.3. An example

73 30.0 | | | I 20.0 Ti - -10.0 Il -' -'3n I -30.0 i -10.0 I- ~ DE. -, c -20.0. l -30.0 3. ~,-40.0.. FEM-BI: E-pol |~~ ~~~-O -5. FEM-BI: H-pol -50.0 2 C ai of 2-D MoM: E-pol -60.0 in a 2-D MoM: H-pol -70.0 o -80.0 0.0 30.0 60.0 90.0 120.0 150.0 180.0 Angle (<) [deg] Figure 5.2: Comparison of 2-D moment method results and FE-BI results for a 450 x 5Ao x O.1Ao air-filled cavity which is recessed in a 1Ao cylinder.

74 Far-Zone Subtraction Procedure: Wraparound Cavity Ogive End-Cap Ogive End-Cap Ogive End-Cap Ogive End-Cap minus approximately Figure 5.3: Far-zone subtraction scheme for the simulation of an infinite cylinder. Figure 5.3: Far-zone subtraction scheme for the simulation of an infinite cylinder.

calculation for the latter case is given in Figure 5.4 where a bistatic scattering pattern is presented in the t = 0~ plane due to a plane wave incident at (0, = SOc)o = 0~). The BOR model was a 10Ao long cylinder with 360 ogival end-caps. The wraparouind groove had a height of 3Ao and a depth of 0.1kA. Clearly, there is good agreement between the FE-BI results and data based on the BOR formulation. 5.2 Scattering Studies The previous comparisons serve to validate the formulation. Having done so, it is instructive to examine the effect that curvature has on the scattering properties of cavity-backed patch antennas. Consider a 2 cm x 3 cm patch residing on a 5.0 cm x 6.0 cm x 0.07874 cm cavity which is filled with a dielectric having Cr = 2.17. Figures 5.5 and 5.6 illustrate the behavior of this geometry as a function of frequency and curvature for E- and H-polarization, respectively. Evidently, the resonance behavior of this patch is sensitive to curvature for both principal polarizations. The frequency response for E-polarization is more sensitive to curvature since the radiating surface field component is parallel to the long side of the patch and cavity. If the patch and cavity were oriented so that the long side is parallel to the +-direction, the response to H-polarization would exhibit greater sensitivity. Such an effect is important to low observable antenna designers since they want to operate the antenna in the region of lowest RCS. This low return region is a consequence of delicate cancellations due to the physical layout of the aperture. Such cancellations are not as complete for highly curved structures as they are for planar cavities. Conformal antenna designers often use wraparound antenna arrays to achieve omnidirectional coverage. Two different configurations are typically used: a wraparound cavity where the cavity is filled with a single continuous collar of dielectric and dis

76 20.0. FE-BI o BOR f 10.0 0.0- 0 L -10.0 0 0 -20.0 0.0 30.0 60.0 90.0 120.0 150.0 180.0 Observation Angle (0) [deg] Figure 5.4: Comparison of the RCS computed via the FE-BI method and a BOR code for a 3Ao x 0.1Ao air-filled wraparound cavity excited by a H-polarized (y = 90~) plane wave with oblique incidence (Oi = 80~, i = 0~).

10.0.0 7 -10.0 hi -20.0 (^C3~~~~ -20.0~~~ -^a= 5 cm Be' *, t —------ a= 10cm P -30.0 --- a= 15 cm - I'\,- - - - - - a = 20 cm -40.0 -- - - = 100 cm l * Planar -50.0'' {' 3.0 3.5 4.0 4.5 5.0 Frequency [GHz] Figure 5.5: RCS frequency response for a 2 cm x 3 cm x patch residing in a 5 cm x 6 cm x 0.07874 cm cavity with cr = 2.17 as a function of curvature for E-polarization (y = 0~).

78 20.0 I 10.0.0 "0' -10.0.=5..... —------- a= 10cm -20.0 ------ a 15 cm'7 a = 20 cm. / -30.0 - ---- a= 100cm 0 * Planar -40.0: I. I 4.0 4.5 5.0 5.5 6.0 Frequency [GHz] Figure 5.6: RCS frequency response for a 2 cm x 3 cm x patch residing in a 5 cm x 6 cm x 0.07874 cm cavity with er = 2.17 as a function of curvature for H-polarization (y = 90~).

79 crete cavities symmetrically placed around the circumference of the cylinder. These two configurations are shown in Figure 5.7. Since near resonance, the radiation properties of these two types of antennas are similar if all elements are uniformly driven, any RCS advantage which one might possess could govern the appropriate choice of arrays. Figure 5.8 compares the E-polarized monostatic scattering at 3 GHz in the 0 = 90~ plane for a continuous wraparound cavity and four discrete cavities; where the patches and cavities are identical to those used in the previous example. Not surprisingly, the continuous wraparound structure has a higher return due to coupling within the substrate. However, since in this case the scattered field is due to the z-component of the surface field (4-directed magnetic currents), both cavities yield large scattered fields in the four directional lobes. Figure 5.9 is the corresponding comparison for H-polarization. In this case, the scattered field is attributed to the (-component of the surface fields (z-directed magnetic currents). Therefore, substrate modes diffract near the patch resulting in discrete lobes for the discrete array while creeping waves shed isotropically for the continuous wraparound cavity. Low observable designs will favor discrete cavity arrays over wraparound cavities since the scattering may be channeled in preferred directions and the overall scattering level is consistently lower. A final example is shown in Figure 5.10 where we observe that other than the expected higher scattering from the wraparound cavity, the scattering behavior of the two arrays in the 0 = 00 plane for H-polarization is very similar. 5.3 FE-ABC: Validation In the previous section, the FE-BI formulation was used to investigate the role that both curvature and cavity dimensions play in determining the overall scattering properties of cavity-backed conformal antennas. The FE-BI method is extremely

. ~ ~ ~ ~.............................................................................. —.................... -................-........i..i. -...;.:.... —. i- -: f- EE i E. E! iEi.E - i;: _ - E.:. -... -..:... -...,.. i. i:..:.:- _ _.........-............... _................ -.. E.....................I...-................... -.-......................-..... —............- SS..i.-...:........:..-.. _.............-. _...............E............................., _..... -. -................-_ _...........................-._.......... E.. i...................;...E:.::-................................................................. -.-1-...... 1.........E......-E. - -... -... I............ I.:............................-........I..............I............ -.................................................::::::::::::::::::::.:::::::::::::::::.::........................:.:..::..::..:::.....:::::::-::.........::::::::::::::::::::::::::::::::::::.............:::::::::::.:::::::::......:::.::::.:.:.:.:::.:::::::::::::::........:...............:::::::::.......:::::::::::::::::::::::.::::::::::::..::::::::::::::::::::::::::..:::::::::::::::::::::::....::::::::::.::::::.::::::::::::.::::::::::::.........::::::::::::::::::::::::::::........:: i:.......:.:::.:.:::::.: i::::E EE:::::::E:::::::::::~:E E:::::::::....::::::::::..... -..:::::::::::::::::::.::::::::::..I......::::::::::::::::..............::::.........:::::::::........::::::.::::::................:::...::::::::::::::::::::::iE i:::::::::::::::::::::::::::1:.::::::::.::::::::::::::::::.:::::...:::.::::::::::::::.:::::::.:::.:::::::.:::::::................::::.:::::::::.:.:::::::::::.........::...............:::: i iE::: iE::: i::::::::::i:::E::::EE:E i:;ES:: i:::D::.:::::::::Ei:::..:::.::.::::.:::::::......::::::::::::..:::..::.::::.::I:............::::::...............:::::::::::::R::::i:: E:~:E::E:::..:.:::::::::::::::......:::::::::::.::..::.::::::.::::...-:::..........:::..............(b).. Figure5.7: Illustration of two types of wraparound arrays: (a) continuous cavity......... array; ~ ~ ~..........et ca it.array...

81 2 0.0............. I...........I................ 0 Wraparound 15.0 5.0 10.0 -5.0 -5.0 -', v - 0.0 10.0 20.0 30.0 40.0 50.0 60.0 70.0 80.0 90.0 Angle (~) [deg] Figure 5.8: Comparison of E-polarized monostatic RCS at 3 GHz for a four patch array placed on a wraparound collar or in four discrete cavities. The patches and cavities are identical to the ones used in Figure 5.5. The observation plane is 0 = 90~.

82 2 0.0............... I............... 15Wraparound 15.0 -------- Discrete mq 10.0 V IO. 5.0 o t,, 3 o 9 \ i'#g/s\ -%:... I, i l''' -10.0 0.0 10.0 20.0 30.0 40.0 50.0 60.0 70.0 80.0 90.0 Angle (4) [deg] Figure 5.9: Comparison of H-polarized monostatic RCS at 3 GHz for a four patch array placed on a wraparound collar or in four discrete cavities. The patches and cavities are identical to the ones used in Figure 5.6. The observation plane is 0 = 90~.

83 10.0..............I.............. I.... Wraparound 5.0 -------- Discrete 0.0 Q 0d -5.0 I — 10.0 S'''' \ /I -15.0' i -20.0................... 30.0 40.0 50.0 60.0 70.0 80.0 90.0 Angle (0) [deg] Figure 5.10: Comparison of H-polarized monostatic scattering at 3 GHz by a four patch array placed on a wraparound collar or in four discrete cavities. The patches and cavities are identical to the ones used in Figure 5.6. The observation plane is q = 0~.

84 powerful and efficient when implemented with a BiCG-FFT solver. However, as previously mentioned, it is not very flexible due to the uniform surface griddilng'requirement and limited number of available dvadic Green's functions. Rather. the FE-ABC method holds greater promise for impact since it can readily be used to model arbitrary patch and cavity shapes, coated and doubly curved conformal antennas. To date, the FE-ABC method has only been used to study the scattering by finite bodies. In this section, the FE-ABC method will be compared with rigorous FE-BI results for typical uncoated cavity-backed antennas which are recessed in an infinite metallic cylinder. Since the ABC surface must be displaced a distance from the cavity aperture in order to recover accurate results, the minimum displacement will be determined for computational domain minimization. Since the ABC absorbs outgoing waves, each principal polarization and typical incident angle directions will be examined. Consider a 2 cm x 3 cm patch antenna which is printed atop a metallic cavity encasing a 5 cm x 6 cm x 0.07874 cm substrate which has a dielectric constant of. = 2.17 and scattering calculations are made at 3 GHz. The second-order ABC (4.25) is placed rA, from the cavity aperture while the lateral walls (e.g. n = ~f or n = ~2) of the ABC are placed 0.5Ao from the cavity aperture. Experience suggests that 0.5Ao lateral displacement is sufficiently large so that the only surface that need be considered is the normal (A = p) surface which is placed at different distances, T, from the cavity aperture as shown in Figure 5.11. The bistatic scattering for normal incidence (qi = 0~,0i = 90~) and observation in the 0 = 900 plane due to an E-polarized plane wave is shown in Figure 5.12. To determine of the effectiveness of the ABC, Figure 5.13 illustrates a ten degree section of Figure 5.12 near the

~S 85 backscatter region. Evidently, if the ABC surface is placed at least 0.3A. fromT the surface, the FE-ABC method recovers the FE-BI data to an acceptable accuracy. The corresponding patterns for H-polarization are shown in Figures 5.11 and 5.15,. respectively. For this polarization, the FE-BI results are recovered when the ABC surface is only 0.1Xo from the cavity aperture. Such discrepancy in performance between the two principal polarizations illustrates the main drawback of the FEABC method. Since the ABC absorbs only the propagating waves, it must be placed beyond the region where near-fields are significant. However, the near-field/farfield boundary is problem and excitation dependent. Nevertheless, experience has shown that an ABC displacement of 0.3A0 is a good rule-of-thumb for scattering computations. In the previous examples, the incident field was normal to the patch where the ABC is seen to work quite well. However, the performance of the ABC must be shown for grazing and rear-side incidence as well as normal. To do so, the same antenna element used previously, with the ABC surface displaced 0.3A0 from the cavity aperture, is compared in Figures 5.16 and 5.17 with the FE-BI results for observation in the 0 = 90~ plane. The FE-ABC method is seen to recover the FE-BI results even for incidence from the rear side of the cylinder.

86 Patch Sabc ^^^^<~~~~ \*~ ~Region I ^ ^.Probe feed 0.5 0.5' Figure 5.11: Cavity-backed patch antenna with ABC mesh termination.

87 0.0' I' -10.0 o -20.0 I5 -30.0 -50.0 - - = 0.3X - t- --- - =o 0.4)X -60.0 F o FE-BI -70.0.. 0.0 30.0 60.0 90.0 120.0 150.0 180.0 Angle (~) [deg] Figure 5.12: Convergence of the second-order ABC as a function of displacement from the cavity aperture for E-polarized bistatic scattering. The reference data is provided by a rigorous FE-BI formulation for the same cavity-backed antenna.

88 0.0 -1.0 - = O. I x -2.0 ------- = 0.2 o -3.0 ------ t= 0.3k - - - - - - = 0.4X -4.0 - 4. o I E FE-BI c/s -5.0 o -6.0 -7.0. —. -g.0 ~- O o o 0 0 o o o 0 -9.0 -10.0 - -10.0 L, I. I I. i I. i, 0.0 1.0 2.0 3.0 4.0 5.0 6.0 7.0 8.0 9.0 10.0 Angle (p) [deg] Figure 5.13: Magnified view of Figure 5.12.

89 -10.0......... -20.0? -30.0 ---- T=O.1X I: ------ x = 0.3X -: ------ -=0.4 -50.0 o FE-BI -60.0..... 0.0 30.0 60.0 90.0 120.0 150.0 180.0 Angle (Q) [deg] Figure 5.14: Convergence of the second-order ABC as a function of displacement from the cavity aperture for H-polarized bistatic scattering. The reference data is provided by a rigorous FE-BI formulation for the same cavity-backed antenna.

90 -10.0 I' I' I' I' -' -11.0 __.l -12.0 - - = 0.2X:m -13.0 -- - T= 0.4X < -14.0'~ o FE-BI C) -15.0 -16.0 -17.0 -18.0 -19.0 -20.0 0.0 1.0 2.0 3.0 4.0 5.0 6.0 7.0 8.0 9.0 10.0 Angle (p) [deg] Figure 5.15: Magnified view of Figure 5.14.

91 0.0 -10.0 FE-ABC -20.0 > ~~-260.0 -~o FE-BI -s -30.0 N -40.0 -50.0 Q -60.0 -70.0 C -80.0 O -90.0 -100.0 -110.0 -120.0..... 0.0 30.0 60.0 90.0 120.0 150.0 180.0 Angle (~) [deg] Figure 5.16: Monostatic scattering by a 2 cm x 3 cm x patch residing in a 5 cm x 6 cm x 0.07874 cm cavity with cr = 2.17 for E-polarization (y = 0~).

92 -10.0: -20.0: -*H~~~~ —- FE-ABC O ~~-20.0 F<o FE-BI -30.0 U -40.0 -50.0 -60.0 -70.0 0.0 30.0 60.0 90.0 120.0 150.0 180.0 Angle (Q) [deg] Figure 5.17: Monostatic scattering by a 2 cm x 3 cm x patch residing in a 5 cm x 6 cm x 0.07874 cm cavity with ~r = 2.17 for H-polarization (y = 900).

CHAPTER VI Radiation In the previous chapter, the scattering behavior of cavity-backed antennas which are recessed in a metallic cylinder were investigated. Both the FE-BI and FE-ABC methods were used and found to be accurate as compared to available reference data. In this chapter, the radiation parameters of conformal patch antennas mounted on singly curved surfaces are examined. In particular, the antenna pattern and input impedance of such patch antennas are investigated. Approximate analysis methods such as the cavity model have been used to design conformal antennas. Initially, planar solutions [3] were obtained which have been found to match measured antenna patterns for operation at resonant frequencies. Since the cavity model assumes single mode fields below the patch, such agreement is understandable. However, in contrast, input impedance calculations using the cavity model are not accurate with respect to either input resistance or resonant frequency localization. The reason for this inaccuracy is the fact that a typical probe feed generally excites multiple modes and the fact that, in reality, the field is not confined to the region immediately below the patch. This latter bias is partially compensated by utilizing effective permittivities rather than measured permittivities. Typically, with this empirical modification, the resonant frequency localization is acceptable; 93

94 however, the input resistance calculations remain inaccurate in comparison wit i measured data. Although the cavity model is limited to resonant patches and its input impedance predictions are inaccurate, the cavity model has been successfully extended to rectangularcylindrical patch antennas [4, 5, 53]. Such extensions to the cavity model have been used to explore the effect of curvature on both the radiation pattern and the input impedance of curved patches. Once again, the calculated patterns compare favorably with measurements [5, 54] while the input impedance calculations are rather poor. Use of an effective permittivity improves the resonant frequency localization; however, the input resistance is still unacceptably inaccurate for design use. Patch antennas fed by a probe have a very small bandwidth. In an effort to analyze printed broadband antennas such as spirals, rigorous numerical techniques based on integral equation formulations have been presented in the literature [6]-[9]. Since the metallization pattern is assumed to be printed atop a dielectric coated groundplane, all of these formulations employ a complex dyadic Green's function. Typically, a rather simplistic feed model for a current probe is utilized and although the radiation pattern of such a model agrees with measurement, the input impedance once again has unacceptable error. Recently, Aberle, Pozar and Birtcher [55] have investigated improved feed models which include a finite probe diameter and an impressed voltage source rather than a current source. Although most numerical formulations require a planar structure, these integral equation techniques have also been extended to patches printed on a coated cylinder [11, 12]. Even though such formulations have proven useful for antenna analysis and design, the assumption of a coated, infinite substrate is not realistic. In practice, spurious coupling between elements within an array will occur on such a structure and as a consequence, each

95 antenna element is typicallS placed within a substrate filled metallic cavity. Aberle recently developed an integral equation formulation for planar cavity-backed antenna elements [13]. However, this formulation suffers, as do most integral equation techniques, from a very large computational and memory burden due to the fully populated matrix and complex dyadic Green's function. If uniform zoning is used on the patch surface, the resulting matrix will be block Toeplitz and accordingly it may be efficiently solved using the BiCG-FFT solver employed in this work. However, such an implementation would still require a highly complex dyadic Green's function. In addition, inhomogeneous substrates or radomes may not be efficiently solved using a surface integral equation. In response to the inefficiencies and limitations of the integral equation formulations, the FE-BI method was introduced by Jin and Volakis [17]. This approach has proven quite successful for planar antenna radiation pattern calculations. As discussed earlier, the inaccurate input impedance data computed by that code may be remedied by using an accurate feed model which is currently under investigation at the Radiation Laboratory. In this thesis, the FE-BI formulation is extended to circular cylinder platforms. In this chapter, the two finite element formulations, the rigorous FE-BI method and the approximate FE-ABC method, will be utilized to study cylindrical- rectangular patch antenna radiation. First, the FE-BI method will be validated for antenna pattern applications by comparison with measured data. The input impedance calculations will be seen to suffer from the same inaccuracies as the planar implementation due to the use of a similar simplistic probe. After validation, we will examine the effect of curvature on the pattern and the input impedance of the cylindricalrectangular patches. An important class of antennas which wrap around the cir

96 cumference of the platform cannot be studied using planar models. Suchl anteIlias are used for communications and guidance control since omnidirectional coverage is essential. The proposed FE-BI has the capability of examining three typical types of wraparound antennas: discrete cavity, continuous cavity and single collar antennas. Since all three types of antennas may be analyzed by the FE-BI method, comparisons may be made between the three in terms of coverage and mutual coupling between elements. Finally, as was done in Chapter V for scattering calculations, the FE-ABC method will be compared to the FE-BI method and a minimal ABC displacement will be established. 6.1 FE-BI: Validation Two types of patch antenna elements are investigated and these are shown in Figure 6.1. Each patch is aa~ x b in size where a denotes the radius of the cylinder and a is the subtended angle of the patch. A patch whose radiating side walls form constant 4-surfaces is termed an axially polarized patch and is fed at Of =. Circumferentially (or azimuthally) polarized patches have radiating walls forming constant z-surfaces and are typically fed at z/ =-. Observation in the 0 = 90~ plane is the E-plane for circumferentially polarized patches and the H-plane for axially polarized elements. The terminology originates with the cavity model for patch antennas. We will now characterize a typical cavity-backed patch antenna. Several computed and measured antenna patterns have been published for patches printed on a coated cylinder. One such patch, which is 3.5 cm x 3.5 cm, was used by Sohtell [54] to compare the accuracy of the cavity model [5] to a surface current integral equation [12]. The measured data was taken at 2.615 GHz for a metallic cylinder which was 63.5 cm long and had a radius of 14.95 cm. The cylinder was

97 T b/2 probe feed + radiating sides b/2 I| a:_____ (a) radiating sides probe feed -t byFig: In of () a cy p d p h e; ad (b) an axially polarized patch element. The radius of the cylinder is denoted by a.

98 coated with a 0.3175 cm uniform dielectric having relative permittivity of r =. 2. Data was taken for -180~ < Q < 180~ in the 0 = 90~ plane corresponding to thie iplane for circumferentially polarized elements and the 1-plane for axially polarized patches. Figure 6.2 compares these measured patterns with data generated using the FE-BI formulation for an identical patch placed within a 360~ x 7 cm cavity. This wraparound cavity best simulated the measured coated cavity. Note that the H-plane patterns are symmetric due to the symmetric placement of the feed, whereas the E-plane patterns are not symmetric. The placement of the feed was not specified in [54]; however, the agreement for the E-plane pattern shown in Figure 6.2 indicates that the position used in the FE-BI model (abf = -1 cm) is reasonable. The feed was placed at Zf = -1 cm for the axially polarized (H-plane) case. In addition to single patches, the FE-BI formulation may be used to design microstrip arrays. Such an approach includes mutual coupling between elements which is ignored by the cavity model. Furthermore, the FE-BI formulation consumes less computational resources than a comparable integral equation formulation due to the sparsity of the FE matrix. The H-plane pattern of a four element array was measured to gauge the accuracy of the FE-BI approach. Each element is 2 cm x 3 cm and placed within a 5 cm x 6 cm x 0.07874 cm cavity which is filled with a dielectric having ET = 2.17. The cylinder is 91.44 cm long and has a radius of 15.24 cm. The cavities are placed symmetrically around the cylinder (e.g. a patch is centered at 0~, 90~, 1800 and 270~). Only the patch centered at 0~ was excited while the remaining patches were terminated with a 50Q load. The driving patch is axially polarized and the feed is located at zf = -0.375 cm. Figure 6.3 illustrates the excellent agreement between the FE-BI formulation and the measured H-plane data.

99 0.0' I' -10.0 "E =: 3'/ e _: -20.0 E 0 _~ -30.0 -, FE-BI: E-plane E 0 a ------—. FE-BI: H-plane'. Eo Measured: E-plane -40.0 El Measured: H-plane -50.0................. -180.0 -90.0 0.0 90.0 180.0 Angle (p) [deg] Figure 6.2: Comparison of measured [54] and computed data for a circumferentially polarized element (E-plane) and an axially polarized element (H-plane). The 3.5 cm x 3.5 cm antenna was printed on a cylinder with a radius of a = 14.95 cm and a 0.3175 cm coating having cr = 2.32. The probe feed was place at (a4q, zf) = (-1.0,0.0) for the circumferentially polarized patch and at (a/f, zf) = (0.0,-1.0) for the axially polarized antenna.

100 0.0 -10.0- - I~ -20.0 0~ -30.0-,@~:. / ---- FE-BI: Z\'.....Measured 1 -40.0 iI -50.0._ -180.0 -90.0 0.0 90.0 180.0 Angle (~) [deg] Figure 6.3: H-plane pattern for a four element patch array. Each of the patches is 2 cm x 3 cm and are placed symmetrically around the cylinder. Only the patch centered at 0~ is fed while the other patches are terminated with 50Q loads.

101 6.2 Radiation Studies In Chapter V., discrete cavity arrays were found to have a significantly! loower radar cross section compared to a continuous wraparound array. Thus, the size of the cavity had a significant effect on the scattering properties of the array. The axially and circumferentially polarized antennas used by Sohtell [54] were placed within cavities which were 7 cm high and approximately 300, 50~, 90~, 180~, 270~ or 3600 in angular extent. Figure 6.4 illustrates that azimuthal cavity size has little effect on the radiation pattern for a circumferentially polarized element. A similar comparison for the axially polarized patch is shown in Figure 6.5. The back lobe of the antenna (near 4 = 1800) is very small for cavities less than 180~ in extent but increases for larger cavities. For cavities which lie on the forward face of the cylinder, the substrate modes diffract off the cavity walls; an effect which has little influence on the main lobe of the pattern. However, for wraparound cavities and cavities which extend into the back side of the cylinder, the substrate modes shed like creeping waves giving rise to the back lobe. Having established the effect of cavity size on the antenna patterns, it is instructive to gauge the effect that curvature has on the resonance behavior (or gain) of patch antennas. The two antennas were placed in 14 cm x 14 cm cavities recessed in cylinders with increasing radius. The frequency was allowed to vary from 2.4 GHz to 2.7 GHz and the peak radiated power (IEe|,2) was recorded at each frequency. Figure 6.6 illustrates that the resonance frequency increases with increasing curvature for a circumferentially polarized antenna, but the maximum radiated power is similar regardless of element curvature. Note in the cavity model, the radiating edges for a circumferentially polarized patch are the azimuthal walls of the cavity

102 0.0.. -10.0 LI 300 Cavity *0 -20.0 -.. N 5-0 —-- 50~ Cavity...........90~ Cavity Z0:0 1800 Cavity Z -30.0 - - - -. 2700 Cavity - —. 3600 Cavity -40.0I..... -180.0 -90.0 0.0 90.0 180.0 Angle (<) [deg] Figure 6.4: Effect of cavity size on the E-plane radiation pattern of a circumferentially polarized patch antenna.

103 0.0.. -10.0 -20.0 ~r3" ~!I' ~~I" 4 ~-~0, 30~ Cavity\ 1/ i "0 -30.0'A " [ dg* F:''!,' - - 500 Cavity'; l id - 900 Cavity -40.0! -40.0 i',;.... 180~ Cavity Z' 270~ Cavity -5- - - -- 3600 Cavity -60.0 -180.0 -90.0 0.0 90.0 180.0 Angle (Q) [deg] Figure 6.5: Effect of cavity size on the H-plane radiation pattern of an axially polarized patch antenna.

104 35.0... |.... 5'* 0 /, \- \ 1 SE 50 // / /, \\\ 30.0 - / // \ \. vau c d r // rasi!1, /,,, /', " \,,10cm a =14.95 cm ------ a =20 cm L, / / i ~ / ~. / ------ a 200 cm 2.40 2.45 2.50 2.55 2.60 2.65 2.70 Frequency (f) [GHz] Figure 6.6: Resonance behavior of a circumferentially polarized patch antenna for various cylinder radii. ~'~ ~ ", various cylinder radii.

105 (see Figure 6.1) which have a constant separation regardless of the cylinder radius. However, the axially polarized patch has decreasing resonant gain Kwith increasing curvature as shown in Figure 6.7. For this patch. radiation is attributed to the axial magnetic walls of the cavity model which have increasing angular separation with decreasing curvature. These walls radiate strongly away from the pattern peak (4> = 0~). Accordingly, the gain of an axially polarized antenna decreases with increasing curvature. The radiation pattern of a circumferentially polarized antenna is largely unaffected by curvature as shown in Figure 6.8 when excited at a resonant frequency. However, the radiation pattern of the axially polarized antenna broadens as the curvature increases which is illustrated in Figure 6.9. Once again, both relationships are readily explained by considering the effect curvature has on the orientation of the cavity model radiating walls. In addition to the gain and pattern of an antenna, designers require the input impedance for matching purposes. For the antenna examined above (in a 14 cm x 14 cm cavity), the input impedance was calculated from 2.4 GHz to 2.7 GHz for various cylinder radii. Figure 6.10 illustrates that the input impedance of a circumferentially polarized patch antenna is not affected by curvature while Figure 6.11 shows that increased curvature reduces the input impedance of an axially polarized patch. This observation agrees the the results reported by Luk et. al. [53]. An important application of conformal antennas on drones is for guidance control. An omnidirectional pattern is required for this application since the presence of a null would have a disastrous effect on the drone's guidance capabilities. Omnidirectional coverage may be obtained by placing a wraparound array on the cylinder such as the ones shown in Figure 5.7. Consider the axially fed 2 cm x 3 cm patch used previously. Figure 6.12 illustrates the coverage obtained by various discrete wraparound

106 35.0 | 0.0............ II -30.0 a =,, cm 2.40 2.45 2.50 2.55 2.60 2.65 2.70 Frequency (f) [GHz] Figure6.7: Resonance behavior of an axially polarized patch antenna for various cylinder radii.

107 4 0. L * * *....r |,. i i i i |........,......... 30.0- - 20.0 10.0 I'*~.~ ^ f _'_ L — a= 10 cm ". 0 - ------- a =14.95 cm.: ------ a=20 cm -10.0 -; ------ a = 30 cm -180.0 -90.0 0.0 90.0 180.0 Angle (0) [deg] Figure 6.8: Variation of the radiation pattern shape with respect to curvature for a circumferentially polarized antenna.

108 40.0 30.0 I ". ------ a =20 cm -20.0 - -10.0 r0 // -= \\\\3m -180.0 -90.0 0.0 90.0 180.0 Angle (4) [deg] Figure 6.9: Variation of the radiation pattern shape with respect to curvature for an axially polarized antenna.

109 a= 10 cm —.... a= 14.95 cm ---- a = 20 cm - - - a = 200 cm o a-*o Figure 6.10: Input impedance of a circumferentially polarized patch antenna for various cylinder radii. The frequency range was 2.4 GHz to 2.7 GHz and the cavity size was 14 cm x 14 cm.

110 0; a -f *j =- a=10 cm ---- a = 200cm, a!o Figure 6.11: Input impedance of an axially polarized patch antenna for various cylinder radii. The frequency range was 2.4 GHz to 2.7 GHz and the cavity size was 14 cm x 14 cm.

111 10.0 0.0. I....., -180.0 -90.0 0.0 " 90.0 180.0,, ",,,:. Angle (p) [deg] One patch -------- Two patches ------ Four patches ------ Eight patches - -- - -Sixteen patches Figure 6.12: Coverage of various uniformly fed discrete wraparound arrays. Each patch is 2 cm x 3 cm and the feed point is at (f ='center Zf = -0375 cm) where,center denotes the azimuthal center of each patch. The array is operated at 3 GHz. I::'I! /i?1 O'5 s 3~~~~~~~~~~~~~~~~~~~~ Z(~~~~~~~;~'"'':','",,:',,~'. -20.0~~~~~~~~~: "~ -]o.oL ~,',',.:':',l.',:: 0.:,;'! q:'''~~:'~ ~ ~~nl'~,,I'I,', I r,,c,,,,m,, ",n,,e,,e,,in,,,t,~f,,,c,z', —~ 7 cm -'hr 7,,, eoe h zmta etro ahpth h ra is~~~~~~~~!'""ate':',':'

112 arrays (see Figure 5.7a) operated at 3 GHz which corresponds to a cylinder radius of 1.52789Ao. Evidently, as elements are added to the array, phase interference cat1uses side lobes to develop which are seen as ripples in the pattern. The geodesic distance between each element center is: 4.8A0, 2.4Ao, 1.2Ao and 0.AG for two, four, eight and sixteen patches, respectively. When sixteen elements are used, the inter-element spacing is sufficiently small as to suppress side lobes and the desired omnidirectional pattern is obtained. Generally, the required number of radiators will scale with the radius of the cylinder. Figure 6.13 shows similar results for a continuous wraparound array (see Figure 5.7b). Another common type of wraparound antenna utilizes a single, continuous metallic collar as a radiating element. Figure 6.14 illustrates the pattern for such an antenna when uniformly excited probe feeds are place symmetrically around the collar with Zf = -0.375 cm. The inter-feed spacing is identical to the center element spacing used in the previous examples. Evidently, as the number of feeds increases, the currents supported by the collar become more uniform and for the case of sixteen feeds, the pattern is omnidirectional. However, notice that for this case, even when only eight feeds are used, a ~1 dB ripple is obtained. 6.3 FE-ABC: Validation As previously mentioned in Chapter V, the FE-ABC method is considerably more flexible than the FE-BI approach. However, the FE-ABC formulation requires a computational domain which extends some distance beyond the cavity aperture. In Chapter V, a ABC surface displacement of 0.3A0 was shown to be sufficient for most calculations. In this section, similar comparisons are made for both radiation pattern and input impedance calculations.

113 10.0 I,,,,., I,,,,,,'.,,,,, 10.0 i. I:',^ i /m,' 1i,,'^' ~,,. -180.0 -90.0 0.0 90.0 180.0 Angle (() [deg] ---- One patch -------- Two patches ------ Four patches ----- - - Eight patches, I,',I I, Angle () [ldeg] is operated at 3 GHz. Four patches.... Sixteen patches is operated at 3 GHz.

114 10.0......................... M:: \ / \ "0 -10.0 -20.0 -30.0.... -180.0 -90.0 0.0 90.0 180.0 Angle (~) [deg] One feed -------- Two feeds ------ Four feeds ------ Eight feeds.- - - - Sixteen feeds Figure 6.14: Coverage of various uniformly fed collar wraparound antenna. The metallic collar has a height of 3 cm and is fed at symmetric points around the cylinder with Zf = -0.375 cm. The antenna is operated at 3 GHz.

115 Consider a 2 cm x 3 cm patch antenna which was printed atop a metallic cavit! encasing a 5 cm x 6 cm x 0.0i7874 cm substrate having a dielectric constant of - = 2.17. This is the same antenna used for the FE-ABC validation in Chapter V for scattering calculations. The second order ABC was place rTA from the cavity aperture while the lateral walls (e.g. n = ~O or n = H~) of the ABC were placed 0.5AO from the cavity aperture. The H-plane antenna pattern for an axially polarized patch is shown in Figure 6.15 where the probe feed is placed at Of = 00 and Zf = -0.375 cm and the operating frequency is 3 GHz. Figure 6.16 shows the main lobe peak value of this antenna for various ABC surface displacements. In this example, the radiation pattern calculated using the FE-ABC formulation was seen to be in good agreement with the pattern computed using the more rigorous FEBI approach. Similar comparisons for a circumferentially polarized antenna element (e.g. 3 cm x 2 cm patch within a 6 cm x 5 cm x 0.07874 cm cavity with a radial feed at of = -1.41~,zf = 0 cm) are given in Figures 6.17 and 6.18. Evidently, placement of the ABC 0.3AO away from the cylinder's surface is sufficiently accurate for patterns obtained from circumferentially polarized as well as an axially polarized patch antennas. However, as illustrated in Figure 6.19, the FE-ABC method is not sufficiently accurate for input impedance calculations. For the antenna used in this comparison, the FE-ABC solution produces an erroneous value even when the ABC surface is 0.4AX from the aperture. Evidently, the FE-ABC method does not accurately model the fields below the patch radiator unless the ABC surface is placed prohibitively far from the aperture. Even a very large ABC displacement will be inaccurate since the increasing computational error associated with a very large mesh will dominate the decreasing approximation error on the ABC surface. However, as shown in the previous example, far-field calculations are more forgiving

116 20.0.......... |... 10.0 - 0.0 X -10.0 -20.0 = ~^~~~ /:~~~~. =0.2X.. -30.0 - -f - - - -r=0.3X -T -- ---- T = 0.4X -0.0 -50.0:,,.....I...... I,,, 1 -180.0 -90.0 0.0 90.0 180.0 Angle (~) [deg] Figure 6.15: Convergence of the second-order ABC as a function of displacement from the cavity aperture for axially polarizedpatch. The reference data is provided by a rigorous FE-BI formulation for the same cavity-backed antenna.

117 2 0.0... | I.. |... | [.. |... |... |... | I.. | I.. |... 19.0 -9.-____ I = o. 1 18.0 ------ =0.2X 17.0 - T=0.3X ------ T=O.4X )E 16.0 o FE-BI 15.0 ~ 14.0 c.7"13.0 13.0 ----------------------------------------—... 12.0 00 11.0 10.0...... I...... I -10.0 -8.0 -6.0 -4.0 -2.0 0.0 2.0 4.0 6.0 8.0 10.0 Angle (4) [deg] Figure 6.16: Magnified view of Figure 6.15.

118 20.0 L....... |... 10.0 0.0 o -20.0W | =.1 =0.2X -30.0 ------ -=0.3X - -= 0.4X -40.0 - o FE-BI -50.0. -180.0 -90.0 0.0 90.0 180.0 Angle (p) [deg] Figure 6.17: Convergence of the second-order ABC as a function of displacement from the cavity aperture for circumferentially polarized patch. The reference data is provided by a rigorous FE-BI formulation for the same cavity-backed antenna.

119 20.0,. i I.' I... I..! I. [. I. 19.0 __ 18.0 ---- T=0.2X 17.0 - T= 0. --— x- T= O.4X ~) 16.0 0o FE-BI:,,.: 15.0 ~ 14.0 13.0.- - ---—.. -.. - 12.0 o o o o o 0 0 o 0 0 0 11.0 10.0 -10.0 -8.0 -6.0 -4.0 -2.0 0.0 2.0 4.0 6.0 8.0 10.0 Angle (p) [deg] Figure 6.18: Magnified view of Figure 6.17.

120 - =0.1 -.- t=0.2k = —-- 0.3.... — t =0.4k - o- FE-BI Figure 6.19: Input impedance calculations as a function of ABC surface displacement from the cavity aperture. The frequency range was 3.0 to 3.2 GHz with data taken every 5 MHz.

1'21 and hence the FE-ABC method is useful for radiation pattern calculations.

CHAPTER VII Conclusion 7.1 Summary In this thesis, a new approach for analyzing cavity-backed cylindrically conformal antenna arrays was presented. Previously, such antennas were designed using approximate planar cavity models or rigorous integral equation formulations. Although both of these techniques have been extended to coated cylinders, planar solutions were typically used even for curved platform designs due to their simplicity. Both techniques work reasonably well for resonant antenna pattern calculations. However, for scattering or input impedance calculations, the integral equation method has high computational demands whereas the cavity model is too approximate. Additionally, both methods typically involve metallic patches printed atop a coated groundplane. A cavity-backed integral equation formulation is available; however, the calculation of the internal dyadic Green's function limits its utility. Due to our past favorable experience with the application of FE-BI formulations to planar cavity-backed planar arrays, in this work the FE-BI method was applied to curved conformal arrays. A critical feature of the formulation is the use of the iterative BiCG-FFT solver which greatly reduces the memory and computational burden as compared to a direct implementation. Thus by using uniform surface gridding, 122

123 impressive computation speeds have been obtained at the expense of grid flexibility. Another important aspect of this work is the use of creeping rwave expressions for the calculation of the cylindrical dyadic Green's function. The traditional modal series expressions of this Green's function imposes a prohibitive computational burden on the boundary integral calculation while the creeping wave series expressions are both efficient and accurate when used for a large radius cylinder. The FE-BI formulation has been successfully used to study the scattering behavior of conformal antennas. Discrete cavity arrays have been shown to have a lower overall RCS as compared to equivalent continuous wraparound arrays. The primary reason is that the cavities suppress creeping and substrate wave coupling. Such waves diffract into preferred directional lobes which allow the designers to channel the scattered field away from the source. Conversely, wraparound collar arrays radiate isotropically resulting in a more uniform pattern. Therefore, low observable antenna designers would prefer discrete arrays over continuous arrays due to the increased control available with the former conformal elements. The effect of curvature on scattering resonance was also examined. Evidently, for either principal polarization, the resonant scattering is not affected by curvature. However, depending on the orientation of the patch, the null which succeeds the resonant peak is curvature dependent. Highly curved patches have a more shallow null than planar or nearly planar patches. RCS engineers would like to place such a null in an expected threat frequency and accordingly, they would prefer a planar patches in such circumstances. In addition to scattering studies, the FE-BI formulation was used for radiation analysis. Although this implementation is sufficient for radiation pattern calculations, as was the case with the planar formulation, the infinitesimally thin probe

124 feed was found to be too simplistic for accurate input impedance computations. Nevertheless, since this feed is sufficient for resonant frequency localization, the FE-131 method was used to study the effect of curvature on the resonant frequency of curved patch antennas. Evidently, the resonant frequency of an axially polarized element is sensitive to curvature while a circumferentially polarized element is insensitive to curvature. Likewise, the co-polarized pattern associated with the axially polarized antenna exhibited a marked increase in the back lobe while the pattern due to the corresponding circumferentially polarized element was essentially unchanged by curvature. Also, the axially polarized element's co-polarized pattern is dramatically effected by cavity size while the corresponding circumferentially polarized pattern is barely altered by increasing curvature. Likewise, the input impedance of an axially polarized element is dependent on the curvature of the substrate while the circumferentially polarized element has a input impedance which is insensitive to curvature. Another important aspect of this thesis was the application of the FE-ABC method for radiation and scattering by cavity-backed antennas on an infinite structure. In this, a mixed total/scattered field formulation was introduced which substantially reduced the likelihood of error propagation and decreased the computational effort of the implementation. New second-order conformal ABCs [25] were used in this FE-ABC method which greatly reduced the computational domain. An ABC displacement of only 0.3A is usually sufficient for both scattering and radiation pattern calculations. However, no reasonable displacement was found which yielded accurate input impedance data. The FE-ABC method was found to be accurate for both bistatic and backscatter calculations even for incidence in the rear side of the cylinder. Thus, the presence of the ABC does not appear to alter the creeping wave interactions on the surface of the cylinder. Both principal polarizations were used

125 for the scattering study as well as axial and circumferentially polarized radiating aintenna elements. We expect that the FE-ABC method will be accurate for arlbitrarincidence and polarization. 7.2 Future Work Within the framework of this thesis, two obvious avenues of research are suggested: development of an improved feed model and application of the FE-ABC method to doubly curved and coated conformal antennas. Historically, the development of feed models for finite element applications has paralleled the similar development of feeds for integral equation formulations. Both methods initially employed a simple infinitesimally thin wire probe. Integral equations then used increasingly more advanced feeds until the current state-of-the-art, the base voltage gap. Such a voltage gap source was used by Aberle [55] where the current supported by the metal wire is determined via an integral equation. The finite thickness of the wire is simulated by placing the current expansion functions at the center of the wire while testing at the outer radius. Unfortunately, such an approach is not readily implemented in a finite element formulation though advanced FE feed models are desirable. For example, Gong and Volakis meshed a coaxial cable feed using tetrahedral elements to a point somewhat below the cavity where the lowest order coaxial mode was used as a coaxial aperture source. This technique, although rigorous, yielded poorly conditioned matrices and imposed a prohibitively expensive computational and memory burden. A different approach is to apply the lowest order coaxial mode across the aperture of the cable at the base of the cavity itself. Such a model will work quite well at the resonant frequency since in a matched condition, the only coaxial mode present at this aperture will be of the lowest order. However, for operation away from resonance,

126 the mismatch at the cable aperture would require higher-order coaxial modes. A key aspect of this feed model is to use a voltage continuity condition to couple the cable aperture fields to the interior fields. Thus. the challenge is in matching the voltage due to higher order modes with the potential within the cavity. If however a field continuity condition is used, mixed elements must be combined perhaps through the use of a boundary integral. Although this work is useful to the design community in and of itself, its main purpose was to provide a key step in the evolution of conformal antenna simulation from planar structures to doubly curved and coated platforms. The FE-BI technique was used to validate the FE-ABC method for singly curved antennas. Additionally, the FE-ABC method was seen to be accurate for scattering and radiation pattern calculations by cavity-backed structures which are recessed in an infinite metallic surface. With confidence, similar formulations may be developed for doubly curved and coated antenna elements. To do so, new distorted brick elements should be developed. Such elements will not be divergence free and hence the method will require a penalty term; however, use of distorted bricks will allow custom mesh generators to be used which would alleviate the need to utilize a sophisticated (and expensive) general purpose mesh generator such as SDRC-IDEAS. The FE-ABC method should be used rather than the FE-BI approach since a doubly curved surface does not allow use of a BiCG-FFT solver and hence the computational and memory burden of a FE-BI method is prohibitive. Additionally, the FE-ABC approach allows the introduction of complex, inhomogeneous material radomes. For such coated geometries, the equivalence principle will require two vector currents on the computational domain border for far-zone field computations. The reflected fields will need be determined for a doubly curved surface via high frequency techniques and efficient

127 on-surface and far-zone formulae will have to be developed.

APPENDICES 128

129 APPENDIX A Fock Functions The asymptotic form of the dyadic Green's function with observation both on the surface of the cylinder and in the far field involves Fock functions. These have been extensively studied and tabulated by Logan [40]. In this appendix, the numerical evaluation of these functions is performed for both small arguments and large arguments. The on-surface Fock functions used in this work are rv(() w ej ()w2(r) -e' dT 2 V7Jo-J27r/3 w(7r) 33,7r/4~2 W002(T)er d U(g) = ej3/4 ei -jd (A.1) N/7r Jooe-Jz2/3 W2(r) where w2(r) denotes Airy functions of the second kind and w(Tr) is its derivative with respect to T. For small arguments (( < 0.6), the asymptotic expansion of (A.1) is given by V7r3 ~ 7 - j+.. 9 V(g) - 1.0-+2 + 3 + 7 - e 4 g +... 4 60 512 1.0 ___ 5 5 9 U() 1.0- ej-+ 3 + - 3. (A.2) 2 12 64 while a rapidly converging residue series is used for n > 0.6 v(g) ~ = e- J ~- ( - en=l

130 Zeros of the w2(r) and z2(r) Tn -= rn I-j 3 and rT = I| rJ-j |n | Tn l | | Tn 1 2.33811 1.011879 2 4.08795 3.24819 3 5.52056 4.82010 4 6.78661 6.16331 5 7.94413 7.37218 6 9.02265 8.48849 7 10.0402 9.53545 8 11.0085 10.5277 9 11.9300 11.4751 10 12.8288 12.3848 Table A.1: Zeros of the Airy Function 10 U(s) 2e4/ 4 y ()n,)-1 e-jTn (A.3) n=l where rT and Tr are zeros of w2(r) and w2(r), respectively. Those zeros are given in the Table A.1 The far-zone Fock functions are given by g(l)(() - j r f _ J/~ r w,(r) G(')(() = g(')()ej3 F()() = ((e 3 (A.4)

131 where wL'(T) denotes Airy functions of the first kind. U(Tr) is its derivative with respect to T and the appropriate integration contour (F) is given by Logan [401. T wo of the functions required in this work. g(~)(() and g()(F), may be calculated using g(0)(0 = 2.06e-J < -1.3 1.39937 + 6 c(m) =1.39937 + M m! (K)m -1.3< <0.5 m=1 10 ["'( )i( 0.5 < ~ < 4.0 m=1 a'(m)Ai(m) 05 < <40 = 1.8325e [-(o.s823-jO. 594)-j] > 4.0 (A.5) g)() = -j2.0 (2 + j - 025e3 <-2.8 10 e[ka(m)~] m-1 = K Ai m) 0.5 < 0 ~4.0 m=, Ai(m) -1.8325 (0.8823- jO.5094 + e(2).[883o 34 ] > > 4.0 (A.6) with constant c = e-j5"/6 and the coefficients for (A.5) and (A.6) are given in Table A.2. In a similar manner, f(o)(~) may be computed using f(o)() _ j2,(1 0.25 0.5) _j. ) = J2( ~05~3 )6 e<-1'1 =0.77582 + e-j7/3 y, c (,>)m -1.1 < < 0.5 m=1 m=10 [A(m)] =0.0 > 4.0 (A.7) which utilizes the constants in Table A.3

132 Constants for (A.5) and (A.6) m c(m) a'(m) Ai(m) 1 0.7473831 1.01879297 0.5356566 2 -0.6862081 3.2481975 -0.41901548 3 -2.9495325 4.82009921 0.38040647 4 -3.4827075 6.16330736 -0.35790794 5 8.9378967 7.37217726 0.34230124 6 56.1946214 8.48848673 -0.33047623 7 9.53544905 0.32102229 8 10.52766040 -0.31318539 9 11.47505663 0.30651729 10 12.38478837 -0.30073083 Table A.2: Constants for g(O)(() and g(l)(.

133 Constants for (A.7) m c(m) a(m) Ai'(m) 1 1.146730417 2.33810741 0.70121082 2 0.86284558 4.08794944 -0.80311137 3 -2.0192636 5.52055983 0.86520403 4 -9.977776 6.78670809 -0.91085074 5 -14.59904 7.94413359 0.94733571 6 49.0751 9.02265085 -0.97792281 7 10.04017434 1.00437012 8 11.00852430 -1.02773869 9 11.93601556 1.04872065 10 12.82877675 -1.06779386 Table A.3: Constants for f(O)(~).

134 APPENDIX B Iterative Solvers Iterative solvers are useful for large sparse linear systems since they typically require considerably less memory as compared to direct solvers such as LU decomposition. Two iterative solvers are used in this work: the biconjugate gradient and conjugate gradient squared algorithms. The former is used for the symmetric systems which are a result of using the FE-BI method or the FE-ABC method in conjunction with the symmetric ABC (4.32). The latter is used for the FE-ABC formulation which employs an asymmetric ABC such as (4.25) or (4.27) without any artificial symmetry. This appendix will present both algorithms beginning with the biconjugate gradient. Note that the algorithms presented are do not employ any preconditioning. Preconditioned versions are readily found in the literature (see for example [52]). B.1 BiCG Algorithm The biconjugate gradient (BiCG) algorithm is one member of a family of iterative solvers which have proven useful in computational electromagnetics [45]. The BiCG unlike the conjugate gradient (CG) method has the advantage of utilizing only one matrix-vector product in its symmetric implementation. Although the convergence characteristics of the BiCG algorithm are erratic (see for example [45]), it often

135 converges in fewer iterations than the CG algorithm. This section presents the BiCG algorithm appropriate for use with symmetric matrices [46]. Consider the system [A]{x} = {b} (B.1) where [A] is a symmetric matrix, {x} is the unknown data vector and {b} is the excitation data vector. The BiCG pseudocode follows (assuming an initial guess {x}1 which may be {0}): Initialize: {r}- = {b}-[A] {x} {P}1 = {b}-[A] {x} Iterate: an < {r}n X {r}n > < [A] {P}n {P) > {}n+ = {}n + an {p}n {r}n+1 = {r}n-an[A] {P} < {rn+l, {r}n+l > Cn - < {r},I {r}* > {} = r}n+l + Cn {P}n (B.2) where the Euclidean norm is given by M < {x},{b} > = x[m]b*[m] (B.3) m=l In (B.2) the subscripts refer to the iteration and the asterisk denotes complex conjugation. It is necessary to terminate the iteration when the error present in the

136 solution is suitably small. Although the error cannot be directly measured (otlierwise there would be no need to solve the system), the error does control the correctionI factor between each iteration. When the error is large. the correction applied to the new guess is large while a small error results in only a small change in the result from one iteration to the next. Many termination criteria have been used in the past. One of the most popular is < {r}n+l {r}+l > K E (B.4) < {b, {b} > where E < 0.01 is a typical acceptable tolerance. As can be expected, the tolerance may need to be tightened or relaxed depending on the problem at hand and the desired accuracy of the result. Note that e should be kept small for antenna impedance computations but can be relaxed for scattering and radiation pattern calculations. B.2 CGS Algorithm The conjugate gradient squared (CGS) solver, like the BiCG, does not guarantee convergence. Unlike the symmetric version of the BiCG algorithm, the CGS requires two vector-matrix multiplies; however, the CSG solver can be used for either symmetric or asymmetric matrices. As is the case with the BiCG solver, the convergence of the CGS algorithm is erratic; however, it often requires fewer iterations than the either the BiCG or the CG algorithm. This section presents the CGS algorithm appropriate for use with either symmetric or asymmetric matrices [52]. Consider the system [A]{x} = {b} (B.5) where [A] is a symmetric matrix, {x} is the unknown data vector and {b} is the excitation data vector. The CGS pseudocode follows (assuming an initial guess {x}1

137 which may be {0}): Initialize: {rl = {b} -[A] {x} {r = {r}l (B.6) {p}1 = {1} (B.7) {q} = {1} (B.8) pi = 0 (B.9) (B.10) Iterate: Pn = < {)1, {}n-1 > Pn /3 = pn-1 {u} = r}n-1 + S q}n-1 {P}n = u + ({q}n-_ + 1Pn-l) {v = [A] {p}, = ^Pn <{r}l,{ V}> {q},n = {u}- {v {w} {u) + {q} {X}n = {X}nl + a { {r} = )l_ - a [A] {w} (B.11) The same termination criterion as was used for the BiCG algorithm is used for the CGS algorithm.

138 APPENDIX C FFT-based Matrix-Vector Product The numerical solution of integral equations (IE) (or boundary integrals) is often performed by converting the IE into a linear system of equations using the moment method (MM) procedure. The MM solution often requires the generation and storage of 0(N2) matrix entries and, if a direct solver such as LU decomposition is utilized, O(N3) operations are required for its solution, where N is the number of unknowns. If an iterative solver such as the biconjugate gradient (BiCG) method is used, the solution can be found in O(r i N2) operations where r is the number of right-hand sides and i is the number of iterations per right-hand side. However, for circulant or block circulant matrices, the solution may be reached in 0(r -i' Nlog2N) operations. This requires the use of FFTs for computing the matrix-vector products in the BiCG algorithm and consequently the resulting solution scheme is often referred to as the BiCG-FFT method. In this appendix, we will present specific examples of this efficient technique for 2-D and 3-D geometries. We will first look at the relatively straightforward 2-D problem followed by the necessarily more involved 3-D case. C.1 2-D Integral Equations Suppose a flat resistive strip having width w centered at the origin of the y=O plane is excited by an E-polarized (TMz) plane wave as shown in Figure C.1. The

139 y p -w/2 w/2 Figure C.1: Flat resistive strip geometry.

140 appropriate Electric Field Integral Equation (EFIE) may be formed by enforcing thle resistive transition condition on the strip giving EZ(x) = R(x)Jz(x) + 4 Jz( )H2) (kol - x') d (C.1) 2 where Et(x) is the incident electric field, R(x) is the normalized resistivity as a function of lateral position, Jz(x) is the equivalent electric current on the strip induced by the incident field and ko = 2 — is the free-space wavenumber. The time dependency ejCt is assumed and suppressed. In (C.1), the primed coordinate denotes the source point while the unprimed one indicates the test point. Once the current has been determined by solving (C.1), the scattered field in the far-zone is given by,-j(kop-f) E - - p — [ -Zo Jz( x)ejkox cos] (C.2) where Zo is the free-space intrinsic impedance. We proceed with the numerical solution of (C.1) by expanding the unknown current in terms of subdomain basis functions as N-1 Jz(x) = J[n]Wn(x) (C.3) n=O where W t Wn(x) = W(x) nAx — < (n + )Ax — 2 - - 2 = 0 else (C.4) and Ax = N. The weight functions W(x) may assume various forms, one of the simplest being a pulse (W(x) = 1). Substituting (C.3) into (C.1) and performing Galerkin's testing, we obtain ]2 Wn(x)Ei(x)dx = J[n}]{ 2 W2(x)R(x)dx,+ -n=o n= 2 n=n k; 4 Wn(d)Wx('dx)H}')(k,\x -( ) dxdx} (C.5)

141 which is the discrete form of the integral equation. \Nlhen pulse basis functions are used. (C.5) becomes (n+l)A x — -1 f /(n+l)Ar- u ]/ 2 E,(x)dx = J[n'] 68[n -n] R(x)dx + 2 n =0 -k J(n+Al)zA- J(n +l)Ax — 4 2 2H(2) (ko I x x- |) d 4 lJntx-^ Wx- 2 2 (C.6) where the Kronecker delta function [n-n - ] = 1 n= n 0 n n' (C.7) has been utilized to indicate the self-cell. Assuming that the resistivity is constant within a segment (e.g. R[n] = R(x) for nAx < x < (n + l)Ax) and making the change of variables g = x +, = +, we have j(n+ E (- ) dg E J[n']{R[n]Ax[n - r'] + n =0 k j(n+l)A J(n+1)Ax H)(ko - )dd} 4 Ax Axl) (C.8) We now observe that the double integral is in convolutional form and since the segments are of uniform length (Ax), we may introduce the discrete function k, /(n+l)Ax /(nI'+l)Ax g[n - n] - I H(2) (k, x - x'I)d(kdg (C.9) T 4,, JAx AX, and rewrite (C.8) as jn) E2()d = \ J[n']R[n][n - n] +, J[n']g[n - ri] / (n+l )a. ^ ^' An- ] N -1 dnAx n =0 n =0 (C.10)

142 The first sum in (C.10) is recognized as the product of a diagonal matrix (AxR[7l]l6[I1n']) and the unknown data vector (J[n],n7' = 0. 1,2...... - 1) whereas the secolld sum is a truncated, discrete, linear convolution which may be efficiently calculated using FFTs. The discrete form of the IE (C.10) may be written as a matrix equation [Z]{J} = {f} (C.11) where the excitation vector entries are given by I (n+l)A2) f[n] (n+l d= (C.12) and the impedance matrix is given by Z[n, n'] = AxR[n]6[n- n'] + g[n -n] (C.13) Alternatively, in preparing to take advantage of the convolution property, we may rewrite (C.10) as [AxR[n]6[n- n]] {J} + [g[n - n {J} = {f} (C.14) Obviously, the first matrix-vector product can be trivially computed in N operations y[n] = JP[n]AxR[n] n = 0,1,2,..., N - 1 (C.15) where the superscript p denotes the iteration. The second matrix-vector product involves a fully populated matrix and would thus normally require 0(N2) operations per iteration for its execution. However, since g[n- n'] is a discrete convolution operator, the product may be computed by invoking the discrete convolution theorem {J} [g[n- n']] = H[N]Fb1 {FD {j } *FD {g}} (C.16)

143 where * denotes a Hadamard product and the discrete Fourier transforml (DFT) pair is given by M-1 FD{fi[ ]} = Ef[m]e-Ijm m=O 1 M -1.D {F[q]} = F[q]c fmq (C.17) q=O In practice, a radix-2 FFT algorithm is used which requires only O(NAlog2 0N) operations as compared to the DFT given in (C.17) which requires 0(N2) operations. The pulse function H[N] indicates that only N of a possible M values of the discrete convolution are retained. Also, the DFT and its inverse operate on periodic sequences of period M (all sequences with a tilde are periodic in this appendix). For example, the unknown iterate (JP) is given by JP[n] = J[n] n=O,1,2,...,N-1 0 n = N, N + 1, N + 2,..., M - 1 (C.18) It is instructive at this point to examine a couple of fundamental properties of DFTs and discrete convolutions. DFTs always assume periodic sequences as shown for example in Figure C.2. In this, a single period of the data sequence to be transformed is indicated by the standard sequence window. This particular data segment is the one usually provided to the DFT. However, since the sequence is periodic, any full period sequence window is permitted as shown in Figure C.2. Accordingly, in the BiCG-FFT method, te user may select which sequence window to use provided an appropriate mapping functions between the data sequence and the unknowns are specified. This flexibility may ease the programming task as is discussed below. Computation of a discrete linear convolution via a discrete circular convolution requires that M > 2N - 1 values of the matrix sequence be provided. In particular, these

144 A 4 A f A A 4 41 i A standard sequence window alternative sequence window Figure C.2: Periodic sequence with even symmetry.

143 values must span all possible combinations of source and test unknown interactions. The utility of the BiCG-FFT method rests upon the fact that for certain problems. a particular redundancy exists within the matrix. For example, suppose the strip shown in Figure C.1 were divided into N segments numbered consecutively from left to right. The resulting matrix will be Toeplitz and there will be N distinct matrix entries corresponding to the first row of the matrix. The first row can be loaded into the left half of the standard sequence window shown in Figure C.2. The Nth row of the matrix forms the right half of that sequence and since the matrix is Toeplitz, the sequence is even symmetric. Thus, the left half of the sequence corresponds to interactions between a test cell and source cells which are physically to the right of the test cell. Likewise, the left half of the sequence corresponds to interactions when the source cells are to the left of the test cell. The first entry in the sequence is the self-cell contribution. Once again, the data sequence contains the full span of possible physical alignment between the test and source cells. This example is a case where the even symmetry replication rule may be applied g[n] = g[n] n- =0,1,2,..., N -1 M = 0 n = N,N+1,M+2,..., 1 MM M M = g[M-n] n= 2 +' + 2,..., M-1 (C.19) 2J 2 2 2 and we have assumed that g[n-n'] = g[n -n]. As shown in (C.19), if M > 2N- 1 it is necessary to pad g[n] with zeros in the middle of the sequence. Use of a replication rule reduces the computational burden to the method since only one row of the matrix need be computed. For situations where first and last rows of the matrix are not mirror images of each other (e.g. it matters whether the source cell is to the right or left of the test cell), an asymmetric sequence is requires as shown in Figure C.3. In

146 A A A A A A! A A A A ^ AA A ^ A A' I I ii I _ standard sequence window alternative sequence window Figure C.3: Asymmetric periodic sequence.

147 this case. all 2N - 1 interactions need be explicitly computed. The 2-D examples given in this appendix all possess even symmetry. However. the cross-terms presented in the next section for 3-D problems are not even symnmetric. As it turns out however, these terms do have some redundancy which may be exploited by using a non-standard sequence window in order to reduce the matrix fill time. The details of such an operation are beyond the scope of this appendix. Combining (C.15) and (C.16), the matrix vector product in the iteration cycle of the BiCG-FFT method is efficiently computed as [Z]{JP} = Axz JP[n]R[n]}+ H[N]' {FD {JP} FDi}} (C.20) This computation requires O(N + (M) log2(M)) operations, provided radix-2 FFTs are used in (C.20), and should be compared to the standard matrix vector product N-1 [Z]JP} = JP[n']Z[n, n] n=0,1,2,...,N-1 (C.21) n =0 Figure C.4 illustrates the comparison between (C.21) and the discrete convolution (C.20) computed using radix-2 FFTs. Clearly, for N > 30, (C.16) is more efficient than (C.21). Another 2-D geometry which yields systems which may be converted into circulant matrices is the circular strip such as the one shown in Figure C.5. The EFIE for a resistive circular strip is given by E'(4) R() Jz() + -2 Jz( ) ( 2koa sin ( 2') ) do'do z\V 4 i- 0 2 (C.22) where a is the radius of the arc and a is the subtended angle. Note that if the arc is closed, c = 360~. A discrete form of (C.22) may be obtained by using pulse basis

148 10000.0.......... 9000.0 I ---- 0(N2 / 2I 8000.0 O(N) --------- O(N+(2N- 1)log2(2N- 1)) 7000.0 /. 6000.0 5000.0 4000.0 3000.0 2000.0 1000.0 0.0........................... 1........ 0.0 10.0 20.0 30.0 40.0 50.0 60.0 70.0 80.0 90.0 100.0 Number of Unknowns (N) Figure C.4: Comparison of operation count for a full matrix-vector product verses a FFT-based matrix-vector product for a 2-D formulation.

149 y a/' F I' Figure C.5: Circular arc geometry.

150 and uniform zoning (Ao =.). We have (n+l)AO aE Q - 1 " z( - -)d J[71 ] R[n]Ao6[n-n + JnA Z L ko r(n+l)A1= (n'+l)Ao -- 4 J / H(2 (2koa sin ( ) do'do (C.23) which is similar to (C.8). The entries of the Green's function matrix are given by [o] 4 Jna f Jn/~ (A ( 2 ) 2) s[n - n'] o 2ka sindod (C.24) Indeed, (C.23) may be solved in exactly the same manner as (C.8) using FFTs. The only difference between the flat strip and the circular arc is a possible additional symmetry present in the Green's function sequence. If the arc is closed (a = 360~), we find that N g[N-n] = g[n] n=1,2,3,..., 2- 1 (C.25) which indicates that only N + 1 entries of the sequence need be computed rather than N. Therefore, the solution of (C.23) for a = 360~ requires roughly half the number of matrix entry evaluations as compared to an equal length flat strip. For a < 360~, a similar symmetry exists to a lesser extent as long as a > 180~. In this case, although the Green's function sequence (C.24) is periodic, it is incorrect to assume that (C.23) now involves a circular convolution. The convolution is still a truncated, linear, discrete convolution and in practice the additional symmetry of (C.25) is not exploited unless matrix build time is excessive.

151 C.2 3-D Integral Equations As was the case for 2-D geometries, there are certain 3-D geometries which admit efficient solution methodologies by- making use of the FFT-based implementation of matrix-vector products. Three of these geometries are the flat plate, an impedance insert in a ground plane and an impedance insert in an infinite metallic circular cylinder. These geometries are shown in Figure C.6. We shall now proceed with the formulation for the later two problems (e.g. the planar and cylindrical inserts) using a general coordinate system. In the following, the general coordinates (u,v) may be considered as (u = x,v = y) for the planar case and (u = aopv = z) for the cylindrical case. An integral equation may be obtained by enforcing the standard impedance boundary condition (SIBC) n x n x E(u,v) = -r1Zon x H(u,v) (C.26) where n is an outward directed unit normal and the total magnetic field is given by H(u, v) = H9(u, v) + jkoYo j G2(u -u, vv). [n x E(u, v )]du'dv (C.27) where S denotes the surface of the insert. In (C.27), the surface field H9O(u, v) is comprised of the incident field H (u, v) and the field reflected by the cylinder when no insert is present (Hr(u, v)). For the planar geometry, the second kind Green's function is given by =,, =,, VV1 e-jkoW(u-^)2+v-^ G(u-u,v-v ) = 2Go(u-u, v-v ) = 2 [I +(- 4r+() k2 47r2(u-u')2 + (v-')2 (C.28)

152 2b:. A.....,........... 2a- A.2b + S {>":!SS:>SS X ~~~,.~CC~ / ~ ~~~ ~~~~~~~~~~ I ~- ~~~~~~~;;,ua;Ij:

153 where Go is the free-space dyadic Green's function and I = xx 4 + yy - z is the ideni factor. For a cylinder. the appropriate Green's function is available as al eigenvalue series [29] for small radius cylinders or as a creeping wave series [36] for large radius cylinders. Irrespective of the form of G2, upon making use of (C.27) into (C.26). we have A -., n x n x E(u,v) -Zon x H9, (u, v) + 71 jkon x J G2(u -, -v ) [ n x E(u',v )] du'dv' (C.29) We note that for planar geometries, an appropriate right-handed system is (uv,, n) whereas for cylindrical problems, the right-handed system is given by (n, u, v). Upon decomposing the fields (E, H90) and the dyadic Green's function (G2) in terms of tangential components (u, v), (C.29) yields the following coupled set of integral equations ~~A 7 rgo( Eu(u, v) U: ZoHv9~(u,v) = - ( + 17 jko [E(u', v )GU(u - u, v - v ) - Ev(u', v)Guv (u -u',v')]du'dv' Eu(u,v) v:-ZoHu0(u, v) = + jko j [E(u',v')Gu((u- u', v -) EV(u,vI)GVV(U -- UV - )]du'dv' (C.30) To convert (C.30) into a discrete set of equations, we expand the unknown electric field components in terms of subdomain basis functions NV-3 NU-2 E,(u,v) = E E Eu[t,s]W}W (u,v) t=O s=0

154 C,.-2 AN'u-3 E,,(u,t) = v A E, [t ] s'(U. v) IC.) t=O s=O where NAr denotes the number of nodes in the ui-direction and AX. is the number of nodes in the'-direction. As shown in Figure C. 7, t denotes the row number of the edge discretization and s is the column number. For this example. the unknoxwns are associated with the free-edges and the associated basis functions may be represented as (V - Vb) W-(u,, V) ( if local edge = 1 (Vt - Vb) (vt - V) Vt -Vb) if local edge = 2 (Vt - Vb) O else Ws(u,v)) = Ur if local edge = 3 (u-u) if local edge (Ur - l) (Ur - u) if local edge = 4 O else (C.32) where the local edge numbering is illustrated in Figure C.8. Note that a free edge is one that is not tangential to any metallic walls and that use of (C.32) requires a finite element type assembly. For the discretization shown in Figure C.8, there are (Nu - 1)(Nv - 1) elements and a total of 2NuNV - 3(N, + Nv) unknowns. Of these, (Nu - 1)(Nv - 2) unknowns are associated with the ^-directed field while (Nu - 2)(Nv - 1) unknowns represent the v-directed field. Substituting the field expansions (C.31) into the coupled integral equations (C.30) and employing the Galerkin's testing procedure, we have Nv-3 Nu-2 F[t,s] = > E[ts']gu[t —ts-s]+ t'=O s'=O Nv-2 Nu-3 S 5 E,[t', s]g,,[t-t, s-s'] t = 0, 1,2,..., N-3 t'=0 s =0

155 1 1 ~ V N - 3 t=. v = 1 2 3 N - 2 (a) NV - 2 —---------------------.................................................................................................... 3 i...................................................... s: 1 2 3 N- 3 (b) Figure C.7: Rectangular patch discretization: (a) u-directed edges (b) v-directed edges.

156 Vb -— 1 b - l — I I I1 I I I I I I I I I I Ul Ur U Figure C.8: Local edge numbering scheme.

157 s = 0.1,.2......, -'2',. -3 NU -2 F,[t.s] = E E it'. s']g,[t - t'..-s.] t'=O s'=O Nv-2 Nu-3 > >3 E,[ts']gs4[t- t',s - s'] t =0 1.2,....,\-2 t'=0 s'=0 s =0,1,2,..., A-3 (C.33) where g"''[t - t, s - s(u dudv -jko J j WU(u', v')Wu(u, v)Guv(u - uv - v)du dv' dudv gvv[t-t,s-s] = - [-[]] j We(u, v)W(u, v),_ dudv +jko j J W v(u, v)Wu(u,v)Gvu(u - u, v -v)du'dv dudv jko Wv(u', )Wv(u, v) ( ( - \, v -. d dv g[t- t', s-s] jko j J Wv(u, v)Wu(u, v)G (u- u, v- v)du dv dudv vu-t t', s-] = jkJo /j Wu(u', v')Wv(u, v)GVv(u - u', v - v')du'dv'dudv guv [t, s] == ko \ W ~v'() u - vd vd Fu[t,s] = Zo W(u, v)H9~(u,v)dudv Fv[t,s = -Z Js W( v(u, v)H0(u,v)dudv (C.34) and e refers to the test element while e' denotes the source element. We note that each of the double sums in (C.33) is a truncated, discrete, linear convolution and hence amenable to the BiCG-FFT method. To proceed, we define the 2-D DFT pair (analogous to (C.17)) Ml-1 M2-1 {F2D {f[t, S} = E A f t, Z]e -M 2(sp+t ) t=O s=O {F[ p]} = 1 M -1 M2-1 (C.35)'T'2z}'[q'P]} MiM2 >3 y i[qp]eJMiM2(sP+t) 1V\112 t=O.Q=O

158 Using (C.35), the convolutions in (C.33) can then be rewritten as Af1-1 AM2-1 E Z E[t'. s']g[t - t'. s - ] = 7D {2D {E[t s]} * F2D {g[t.s]}} t'=O sl=O ((C.36) The order of the relevant DFTs must be M1 > 2(number of rows) - 1 and Al12 > 2(number of columns) - 1 where the number of rows and columns of the discretization may vary with each convolution (see Figure C.7). For example, the first convolution in (C.33) is associated with iu-testing and u-source edges and hence the number of rows and columns is (N, - 2) and (N, - 1), respectively. The field sequences are loaded into a M1 x M2 array in row/column order of the field discretization and the remaining entries form a zero pad. The Green's function sequence must be loaded into a similar array (in the same manner) and periodic replication must be performed to provide the necessary "negative lags". If the sequence has the property, g[t - t', s - SI] = g[t - t, sI - s], then this replication process takes the form g[t,s] = g[t,s] 0 t< 1 < s M2 1 - -2 - -2 = g[M1 +2-t,s] < t < M-1 0 < s < 1 2 - 2 = g[t, M2+2-s] 0<t< -1 M2 < < M2- 1 - - 2 2 - - = g[M +2-t,M2+2-s] t < t < M -1M21 2 2 (C.37) For the presented example, guu[t, s] and g,, [t, s] possess this property while the crossterms, guv[t, s] and gv,[t, s], may not. In the case of a planar insert, the cross-terms do possess even symmetry due to the use of a symmetric closed form dyadic Green's function. On the other hand, the cylindrical insert does not have such symmetry since the equidistant interactions may differ by a sign on account of the form of the dyadic

159 Green's function. For example, relative to the first quadrant of the array, interactions in the second and fourth quadrants have the opposite sign while the third quadrant have the like sign. These signs may be accounted for by using a complex replication rule and an non-standard sequence window. Otherwise, all possible lags must be computed requiring longer matrix build time since more than the first u-directed and vi-directed edges need be used as sources. Once the periodic arrays are loaded, the required matrix-vector product for the uu-interactions may be performed in 0((MilogM1)(M2log M2)) operations rather than the 0(((Nu - 2)(Nv - 3))2) operations required for a standard matrix-vector product. The comparison is shown in Figure C.9 with M1 = 2(Nv - 3), M2 = 2(N, - 2) and Nv = NU = N. Clearly, when the number of nodes per side exceeds 10-15, the FFT-based matrix-vector product is more efficient than a conventional matrix-vector product. In practice, the FFT-based product is more efficient than a standard product in terms of wall clock time for N < 10 since in order to exploit the memory savings afforded by uniform zoning of a convolutional kernel without using FFTs, additional overhead is incurred to match the appropriate matrix entry with the correct vector entries. Similar results are obtained for the other convolutions in (C.33).

160.60E+06.55E+06.50E+06.45E+06 —-- 0o( ((N-2)(N-3))2).45E+06 --------........ ( (Mlog2M)(M210g2M2)).40E+06.35E+06 3,.30E+06 o.25E+06.20E+06.15E+06.1OE+06- /.50E+05.....-.OOE+00 I.... O 0.0 5.0 10.0 15.0 20.0 25.0 30.0 Number of Rows/Columns in Grid (N) Figure C.9: Comparison of operation count for a full matrix-vector product verses a FFT-based matrix-vector product for a 3-D formulation. N is the number of nodes in each direction of the grid, M1 = 2(N-3) and M2 = 2(N-2).

BIBLIOGRAPHY 161

162 BIBLIOGRAPHY [1] G.A. Deschamps, "Microstrip microwave antennas" presented at the 3rd USAF Symp. on Antennas, 1953. [2] H. Gutton and G. Baissinot, "Flat aerial for ultra high frequencies," French Patent No. 703113, 1955. [3] Y.T. Lo, D. Solomon, F.R. Ore and W.F. Richards, "Theory and experiment on microstrip antennas," IEEE Trans. Antennas Propagat., Vol. 27, pp. 137-145, Mar. 1979. [4] C.M. Krowne, "Cylindrical-rectangular microstrip antenna," IEEE Antennas Propagat., Vol. 31, pp. 194-199, 1983. [5] J.S. Dahele, R.J. Mitchell, K.M. Luk and K.F. Lee, "Effect of curvature on characteristics of rectangular patch antenna," Electronics Letters,Vol. 23, No. 14, pp. 748-749, 1987. [6] E.H. Newman and P.Tulyathan, "Analysis of microstrip antennas using moment methods," IEEE Antennas Propagat., Vol. 29, pp. 47-53, Jan. 1981. [7] J.R. Mosig and F.E. Gardiol, "General integral equation formulation for microstrip antennas and scatterers," Proc. Inst. Elec. Eng., Vol. 132, Pt. H, pp. 424-432, Dec. 1985. [8] D.M. Pozar, "Input impedance and mutual coupling of rectangular microstrip antennas," IEEE Antennas Propagat., Vol. 30, pp. 1191-1196, Nov. 1982. [9] D.M. Pozar and S.M. Voda, "A rigorous analysis of a microstripline fed patch antenna," IEEE Antennas Propagat., Vol. 35, pp. 1343-1350, Dec. 1987. [10] P.B. Katehi and N.G. Alexopoulos, "Real axis integration of Sommerfeld integrals with applications to printed circuit antennas," J. Math. Phys., Vol. 14, pp. 527-533, 1983. [11] S.B. Fonseca and A.J. Giarola, "Analysis of microstrip wraparound antennas using dyadic Green's functions," IEEE Antennas Propagat., Vol. 31, pp. 248253, 1983.

163 [12] J. Ashkenazv. S. Shtrikman and D. Treves, -'Electric surface current model for the analysis of microstrip antennas on cylindrical bodies." IEEE Antenna. Propagat.. Vol. 33. pp. 295-300, Mar. 1985. [13] J. Aberle, "On the use of metallized cavities backing microstrip antennas", 1991 IEEE Antennas and Propagat. Soc. Int. Symp., Vol. 1, pp. 60-63. June 1991. [14] P. Silvester and M.S. Hsieh, "Finite-element solution of 2-dimensional exterior field problems," Proc. Inst. Elec. Eng., Vol. 118, pp. 1743-1747, Dec. 1971. [15] B.H. McDonald and A. Wexler, "Finite-element solution of unbounded field problems," IEEE Trans. Microwave Theory Tech., Vol. 20, pp. 841-847, Dec. 1972. [16] 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, No. 1, pp. 97-104, Jan. 1991. [17] J-M Jin and J.L. Volakis, "A hybrid finite element method for scattering and radiation by microstrip patch antennas and arrays residing in a cavity," IEEE Trans. Antennas and Propagat., Vol. 39, No. 11, pp. 1598-1604, Nov. 1991. [18] J-M Jin and V.V. Liepa, "Application of hybrid finite element method to electromagnetic scattering by coated cylinders," IEEE Trans. Antennas Propagat., Vol. 36, pp. 50-54, Jan. 1988. [19] J.D. Collins, "A finite element - boundary integral method for electromagnetic scattering," Ph.D. dissertation, Radiation Laboratory, The University of Michigan, 1992. [20] J-M Jin and J.L. Volakis, "Electromagnetic scattering by and transmission through a three-dimensional slot in a thick conducting plane," IEEE Trans. Antennas and Propagat., Vol. 39, No. 4, pp. 543-550, Apr. 1991. [21] A. Chatterjee J.M. Jin and J.L. Volakis, "Edge-based finite elements and vector ABC's applied to 3-D scattering," IEEE Trans. Antennas Propagat., Vol. 41, pp. 221-226, Feb. 1993. [22] A.F. Peterson, "Absorbing boundary conditions for the vector wave equation," Microwave Opt. Technol. Lett., Vol. 1, pp. 62-64, April 1988. [23] J.P. Webb and V.N. Kanellopoulos, "Absorbing boundary conditions for finite element solution of the vector wave equation," Microwave Opt. Technol. Lett., Vol. 2, pp. 370-372, Oct. 1989. [24] A. Chatterjee and J.L. Volakis, "Conformal absorbing boundary conditions for the vector wave equation," Microwave Opt. Technol. Lett., Vol. 6, pp. 886-889, Dec. 1993.

164 [25] A. Chatterjee. Scattering from large three-dimensional struct ures. PhI.D. disse'rtation, Radiation Laboratory. The University of Michigan. 1994. [26] B. Stupfel "Absorbing boundary conditions on arbitrary boundaries for the scalar and vector wave equations." 1993 IEEE.Antennas and Propagat. Soc. Int. Symp., Vol. 1, pp. 296-299, June 1993. [27] J. van Bladel, Electromagnetic Fields, New York: Hemisphere Publishing Co.. 1985. [28] A. Sommerfeld, Electrodynamics, New York: Academic Press, 1952. [29] C.T. Tai, Dyadic Green functions in electromagnetic theory, New York: IEEE Press, 1994. [30] O.C. Zienkiewicz and R.L. Taylor, The finite element method, New York: McGraw-Hill Book Co., 1989 [31] C.T. Tai, Generalized vector and dyadic analysis. New York: IEEE Press, 1992. [32] Z.J. Cendes, "Vector finite elements for electromagnetic field computation," IEEE Trans Magnetics, Vol. 27, No. 5, pp. 3958-3966, Sept. 1991. [33] A. Bossavit, "A rationale for edge-elements in 3D fields computations," IEEE Trans. Magnetics Vol. 24, No. 1, pp. 74-79, Jan 1988. [34] A. Chatterjee, 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. [35] J.L. Volakis, A. Chatterjee and J. Gong, "A class of hybrid finite element methods for electromagnetics: A review," JEWA, 1993. [36] P.H. Pathak and N.N. Wang, "An analysis of the mutual coupling between antennas on a smooth convex surface," Ohio State Univ. ElectroScience Lab., Report 784583-7, Oct. 1978. [37] J. Boersma and S.W. Lee, "Surface field due to a magnetic dipole on a cylinder: Asymptotic expansions of exact solution," Univ. Illinois Electromagnetics Lab., Report 78-17, 1978. [38] T.S. Bird, "Comparison of asymptotic solutions for the surface field excited by a magnetic dipole on a cylinder," IEEE Trans. Antennas and Propagat., Vol. 32, No. 11, pp. 1237-1244, Nov. 1984. [39] T.S. Bird, "Accurate asymptotic solution for the surface field due to apertures in a conducting cylinder," IEEE Trans. Antennas and Propagat., Vol. 33, No. 10, pp. 1108-1117, Oct. 1985. [40] N. A. Logan, "General research in diffraction theory," Lockheed Aircraft Corp., Missiles and Space Div., Vol. 1 and 2, Report LMSD-288088, Dec. 1959.

165 [41] K. Barkeshli and J.L. Volakis.'On the implementation of the conjugate gradient Fourier transform method for scattering by planar plates." IEEE Ante7lnas and Propagat. Mag.. Vol. 32. No. 2. pp. 20-26. Apr. 1990. [42] 0. Einarsson, R.E. Kleinman. P. Laurin. and P.L.E. Uslenghi. "Studies ill radar cross sections L - Diffraction and scattering by regular bodies IV: the circular cylinder," University of Michigan Technical Report No. 7133-3-T, 1966. [43] A.S. Goriainov, "An asymptotic solution of the problem of diffraction of a plane electromagnetic wave by a conducting cylinder," Radio Engr. and Electr. Phys.. Vol. 3, pp. 23-39, 1958. (English translation of Radiotekhnica i Elektronica, Vol. 3) [44] Matrix Computation, Baltimore: Johns Hopkins University Press, 1989. [45] T.K. Sarkar, "On the application of the generalized biconjugate gradient method," JEWA, Vol. 1, No. 3, pp. 223-242, 1987. [46] C.F. Smith, A.F. Peterson and R. Mittra, "The biconjugate gradient method for electromagnetic scattering," IEEE Trans. Antennas and Propagat., Vol. 38, No. 6, pp. 938-940, June 1990. [47] J.M. Jin and J.L. Volakis, "Biconjugate gradient FFT solution for scattering by planar plates," Electromagnetics, Vol. 12, pp. 105-119, 1992. [48] T.M. Hashaby, S.M. Ali and J.A. Kong, "Input impedance and radiation pattern of cylindrical-rectangular and wraparound microstrip antennas," IEEE Trans. Antennas and Propagat., Vol. 38, No. 5, pp. 722-731, May 1990. [49] J.R. Mosig, R.C. Hall and F.E. Gardiol, "Numerical analysis of microstrip patch antennas" in Handbook of microstrip antennas, ed. J.R. James and P.S. Hall, London, UK: Peter Peregrinus, 1989. [50] C.H. Wilcox, "An expansion theorem for electromagnetic fields," Comm. Pure Appl. Math., Vol. 9, pp. 115-134, May 1956. [51] S.M. Rytov, "Computation of the skin effect by the perturbation method," J. Exp. Theor. Phys., Vol. 10, pp. 180, 1940 (translation by V. Kerdemelidis and K.M. Mitzner). [52] H.A. van der Vorst, "BI-CGSTAB: A fast and smoothly converging variant of BI-CG for the solution of nonsymmetric linear systems," SIAM J. Sci. Stat. Comput., Vol. 13, No. 2, pp. 631-644, Mar. 1992. [53] K-M Luk, K-F Lee and J.S. Dahele, "Analysis of the cylindrical-rectangular patch antenna," IEEE Trans. Antennas and Propagat., Vol. 37, No. 2, pp. 143147, Feb. 1989. [54] E.V. Sohtell, "Microstrip antennas on a cylindrical surface," in Handbook of microstrip antennas, Ed. J.R. James and P.S. Hall, Peregrinus: London, pp. 1227-1255, 1989.

166 [55] J.T. Aberle, D.M. Pozar and C.R. Birtcher. "Evaluation of input impedance aI(l radar cross section of probe-fed microstrip patch elements using an accurate feed model," IEEE Trans. Antennas and Propagat.. Vol. 39, No. 12. Ip. 1691-169)i. Dec. 1991.

UNIVERSITY OF MICHIGAN 3 9015 03126 3620