TECHNICAL REPORT 031173-1-T A comparative study of an ABC and an artificial absorber for truncating finite element meshes T. Ozdemir and J. L. Volakis National Aeronautics and Space Administration Ames Research Center Moffett Field CA 94035 Naval Weapons Center China Lake CA 93555-6001 December 1993 31173-1-T = RL-2438

Report #031173-1-T NASA Interchange No. NCA2-839 PROGRESS REPORT Grant Interchange Title: Analysis of Cylindrical Wrap-Around and Doubly Conformal Patch Antennas Via the Finite ElementArtificial Absorber Method A comparative study of an ABC and an artificial absorber for truncating finite element meshes Report Title: Report Authors: Primary University Collaborator: Primary NASA-Ames Collaborator: T. Ozdemir and J. L. Volakis John L. Volakis Volakis@um.cc.umich.edu Telephone: (313) 764-0500 Alex Woo woo@ra-next.arc.nasa.gov Telephone: (415) 604-6010 Radiation Laboratory Department of Electrical Engineering and Computer Science The University of Michigan Ann Arbor MI 48109-2122 University Address: Date: December 1993 Funds for the support of this study have been allocated in part by the NASA-Ames Research Center, Moffett Field, California, under interchange No. NCA2-839.

A Comparative Study of an ABC and an Artificial Absorber for Truncating Finite Element Meshes T. Ozdemir and J. L. Volakis Radiation Laboratory University of Michigan Ann Arbor, Michigan 48109 Abstract The type of mesh termination used in the context of finite element formulations plays a major role on the efficiency and accuracy of the field solution. In this work, we evaluate the performance of an absorbing boundary condition (ABC) and an artificial absorber (a new concept) for terminating the finite element mesh. This analysis is done in connection with the problem of scattering by a finite slot array in a thick ground plane. The two approximate mesh truncation schemes are compared with the exact finite element-boundary integral (FEMBI) method in terms of accuracy and efficiency. It is demonstrated that both approximate truncation schemes yield reasonably accurate results even when the mesh is extended only 0.3 wavelengths away from the array aperture. However, the artificial absorber termination method leads to a substantially more efficient solution. Moreover, it is shown that the FEM-BI method remains quite competitive with the FEM-artificial absorber method when the FFT is used for computing the matrix-vector products in the iterative solution algorithm. These conclusions are indeed surprising and of major importance in electromagnetic simulations based on the finite element method. 1

1 Introduction The finite element method (FEM) is quite attractive for scattering and radiation analysis of complex, large scale structures because of its geometrical adaptability and low memory requirements. As is the case with any partial differential equation method, a difficulty with the FEM stems from the mesh termination condition which must be imposed at some distance from the scatterer or radiator. Typically, an absorbing boundary condition is employed to terminate the mesh, which is simply an approximate relation between the tangential magnetic (electric) field and the tangential derivatives of the electric (magnetic) field at the mesh termination boundary. A second order ABC retains up to second order tangential field derivatives, the highest to be accommodated in an FEM implementation with linear elements. The ABC mesh truncation scheme is quite attractive because it leads to sparse matrices but, unfortunately, it enlarges the computational domain since the ABC must be enforced at some distance away from the surface of the structure. The accuracy of the FEM solution improves as the ABC is enforced further away from the structure's surface but, to date, no hard recommendations are available as to the optimum/acceptable mesh truncation distance. This is because a strong correlation exists among the mesh truncation distance, the solution accuracy and the structure's geometry and composition. Among the different 3D FEM implementations, Mayergoyz and D'Angelo [1991] employed the Engquist and Majda [1977] ABCs and truncated the mesh about one wavelength away from the scatterer. Recently, Chatterjee, etc. [1993] showed that the ABC given by Webb and Kanellopoulos [1989] can be enforced 0.3 wavelengths away from the scatterer with acceptable solution accuracy. This ABC also leads to symmetric FE matrices and will be employed in this study. Alternatively, the boundary integral method can be used for truncating the mesh leading to the finite element - boundary integral method (FEMBI). The latter is a numerically exact method but has the disadvantage of rendering a partly full, partly sparse system whose solution may require greater CPU time and storage. However, if the BI is enforced on special mesh termination surfaces (i.e., cylindrical or flat), the fast Fourier transform (FFT) can be used in conjunction with an iterative solution of the FEMBI system to maintain a low O(N) storage and minimize the CPU time requirements. This has been demonstrated by Jin and Volakis [1991a,b] for 2

the analysis of cavity backed apertures and slots. Recently, a third mesh termination scheme was proposed by Jin etc. [1992] for two-dimensional scattering. This termination, too, is not exact and relies on the use of a metal backed material absorber to terminate the mesh. Because this absorber can not, generally, be manufactured, it is often referred to as an artificial absorber. In this paper, we propose a metal backed artificial absorber (AA) for terminating 3D finite element meshes, and examine its accuracy. A comparative study on the computational efficiency of the FEM-ABC, FEM-BI and FEM-AA methods is also performed for scattering by 3D finite slot arrays. The results of this study are indeed surprising and of major importance to the electromagnetics community making use of the finite element methods for scattering and radiation. 2 Geometry Description The geometry which was used to perform the comparative study on the efficiency of the FEM-BI, FEM-ABC and FEM-AA analysis methods is shown in Figure 1. It is a finite slot array in an otherwise thick, metallic ground plane. The slots are of thickness d and, for this study, they are assumed to have a square aperture of width W. There are, of course, no restrictions placed on the dimensions and geometry of the slot aperture by any of the above methods. Such restrictions result only from the choice of the volume elements used in the discretization of the slot volume. For rectangular shaped volumes, the simple rectangular elements are sufficient to retain geometrical fidelity and will be used in connection with all three finite element formulations employed in this study. Specifically, the edge-based finite element formulation is employed. This formulation is described in [Jin and Volakis, 1991b] and explicit expressions for the matrix elements are given in [Volakis, etc., 1994]. We also remark that the associated computer codes were written such that, the boundary conditions on all metallic surfaces and sections separating the slots are enforced apriori before construction of the final FE system. Below, we describe other details specific to each of the three formulations. 3

3 Description of the Analysis Methods Figure 2 displays cross sectional cuts of the 3D slot array and illustrates the location of the mesh termination boundary used in connection with the application of the FEM-BI, FEM-ABC and FEM-AA implementations. Throughout this study, it will be assumed that the array is illuminated by a plane wave incoming from the upper half space. Given such an excitation, of interest is the evaluation of the reflected and transmitted fields. For the FEM-BI implementation, the boundary integral termination is enforced on the top and bottom apertures of each slot comprising the finite array. The specifics of the implementation are presented in [Jin and Volakis, 1991a,b], and the resulting system is solved via the biconjugate gradient method (BiCG). To retain an O(N) memory requirement, the fast Fourier transform (FFT) was employed in the BiCG algorithm to perform the matrix-vector products associated with the BI subsystem. The specifics of this implementation are outlined in [Jin and Volakis, 1992] and we note that use of the FFT avoids the generation of the full BI subsystem altogether, except for the non-redundant matrix elements. It also leads to substantial speed-up in the execution of the matrix-vector products since the operations are reduced from O(Nb2) down to O(Nb ln Nb)., where Nb is the number of unknowns on all slot apertures. For the implementation of the FEM-ABC formulation, the mesh at the bottom slot apertures was extended 0.3 wavelengths beyond the slot array surface as illustrated in Figure 2(b). The mesh was then truncated by imposing the second order ABC n x H = -Y Et + 2 V x [n (V x E)~] + 2 Vt (V * Et) (1) whose accuracy was already examined in [Chatterjee, etc., 1993]. Generally, three brick elements were used to fill the mesh between the slot array aperture and the ABC surface. As noted on Figure 2(b), the upper mesh surface was still terminated by the boundary integral method. This mixed truncation scheme was selected to minimize the differences among the three implementations, thus, allowing a somewhat more controlled comparison. The implementation of the ABC truncation scheme on the upper half space would require special considerations to remove the large specular field contributions from the ground plane. Since this study is concerned with the relative merits of the three different termination schemes, it was not deemed 4

necessary to apply the ABC truncation scheme on the upper half space as well. For the solution of the resulting system, the BiCG algorithm was again employed. As discussed earlier, the matrix-vector product for the BI subsystem pertaining to the upper aperture elements was executed using the FFT. Figure 2(c) illustrates the mesh termination at the lower half space using a thin 0.15Ao artificial absorber. The absorber dielectric constants were found by minimizing the sum of the TE and TM reflection coefficient magnitudes of the planar metal-backed layer over all real angles of incidence. The Simplex method was employed for this minimization. However, given that the reflection coefficients must be minimized over a path of complex angles to be effective in absorbing diffracted waves, the minimization only over real angles does not necessarily lead to an optimum absorber. Thus, it was necessary to consider many different metal-backed absorbers (including multiple layer absorbers). The effectiveness of each absorber candidate was tested by using it to terminate the mesh at the bottom aperture of a single slot. Also, of considerable importance was the selection of a metal-backed absorber whose refractive index was sufficiently small to avoid a requirement for higher density gridding in the material. Note that the top surface of the absorber layer was placed only 0.15Ao away from the lower slot array aperture. Typically, in our implementation, two brick elements were placed between the slot aperture and the layer's surface. Also, two brick elements were used to model the lossy absorber. Note that the impedance of the absorber material is equal to that of free space since the relative permittivity and permeability are equal. Thus, the absorber simply acts as an attenuator of all waves impinging upon its surface. 4 Comparison of the Mesh Termination Schemes Before performing a comparison on the efficiency of the aforementioned mesh termination schemes, it is first instructive to establish the accuracy of the FEM-ABC and FEM-AA solutions. For this purpose, numerous backscatter and bistatic patterns were computed for different slot and array configurations. In all cases, the truncation schemes illustrated in Figures 2(b) and 2(c) yielded results which were in good agreement with the rigorous FEM-BI solution. Some example bistatic patterns are given in Figures 3 and 4. Figure 3 5

displays the e00s and acr bistatic radar cross-section (RCS) of a 2x2 array. The array is illuminated by a plane wave incident along the direction (binc = 0~, inc = 80~) for the a0e computation and along (qinc = 0~, 'nc = 300) for the a+0 pattern. The displayed pattern is taken in the upper xz plane and the angle a is measured from the x-axis. Figure 4 displays the ae00 bistatic pattern of a 15x2 slot array illuminated by a plane wave at (inc = 300, 0ic = 400). The pattern is taken in the upper b = 30~ plane and a is measured from 0 = 90~. For all calculations shown in Figures 3 and 4, the slot size was 1AoxlAoX0.225Ao as shown in Figure 1, and the period was 1.15Ao in the x-direction and 1.225Ao in the y-direction. We observe from Figures 3 and 4 that all three solutions (FEM-BI, FEMABC and FEM-AA) are in good agreement with each other. As expected, the FEM-AA formulation required less unknowns than the FEM-ABC formulation because all degrees of freedom on the metal backing the absorber are set to zero apriori. The number of unknowns used in connection with each formulation (and the CPU time required to compute the bistatic pattern) are given to the right of the legend on each figure. Of course, the FEMBI formulation required the least unknowns since the mesh is terminated at the surface of the bottom aperture. Note that a fourth curve (identified by the legend cavity) is included in Figures 3 and 4, and this corresponds to the scattering by the same geometry with the bottom aperture of the slots shorted. The inclusion of this curve serves to illustrate the level of power which must be absorbed by the artificial absorber or the ABC. The number of unknowns used in connection with each of the finite element formulations are given in Figure 5. Not surprising, the unknowns grow linearly, but the slope for each of the formulations is different. The FEM-BI method has the smallest slope and the FEM-ABC method has the largest. The corresponding actual CPU times on an HP9000/750 workstation for calculating the scattered field (single incidence angle) are given in Figure 6. It is seen that the FEM-ABC method demands the most CPU time primarily because it is associated with the most unknowns, thus, requiring more iterations for solution convergence. The substantially better performance of the FEM-AA and FEM-BI methods is indeed surprising. In view of the unknown comparisons given in Figure 5, the superior CPU time performance of the FEM-AA method over the FEM-ABC method stems from improved solution convergence and not only from the differences in the degrees of freedom between the two methods. 6

The CPU time performance of the FEM-BI is indeed impressive given the additional computational burden inherent with this method. Obviously, the incorporation of the FFT for computing the matrix-vector products in the BiCG algorithm has played a major role in reducing the CPU time. Most importantly, the improved performance of the FEM-BI method remains unaltered as the array size increases. For this application, both the FEM-BI and FEM-AA methods have comparable CPU performances. But the impressive performance of the FEM-BI method can only be achieved for planar mesh terminations and provided uniform gridding is employed. However, we can contend that the artificial absorber termination can be generalized to terminate meshes on non-planar surfaces making it quite attractive for a variety of applications. 5 Conclusions The presented efficiency comparison of the three mesh termination schemes has demonstrated that, although both ABC and Artificial Absorber terminations had good and comparable accuracies, the Artificial Absorber termination was the one challenging the BI termination in terms of CPU time efficiency. Therefore, the conclusion is that the Artificial Absorber termination is a viable option when the boundary integral termination can not be used. Given its rigor, the boundary integral termination is, of course, the best choice for the analysis of planar array. This is true, provided the FFT is employed for computing the matrix-vector products appearing in the iterative solution algorithm. 6 References Chatterjee, A., J. M. Jin and J. L. Volakis, "Application of edgebased finite elements and ABCs to 3-D scattering," IEEE Trans. Antennas Propagat., Vol. 41, No. 2, pp. 221-226, February 1993 Engquist, B. and A. Majda, "Absorbing boundary conditions for the numerical simulation of waves," Math. Comp., Vol. 31, pp. 629-651, 1977 7

Jin, J. M. and J. L. Volakis, "A finite element-boundary integral formulation for scattering by three-dimensional cavity-backed apertures," IEEE Trans. Antennas Propagat., Vol. 39, No. 1, pp. 97-104, January 1991a Jin, J. M. and J. L. Volakis, "Electromagnetic Scattering by and Transmission Through a Three-Dimensional Slot in a Thick Conducting Plane," IEEE Trans. Antennas Propagat., Vol. 39, No. 4, pp. 543-550, April 1991b Jin, J. M. and J. L. Volakis, "A biconjugate gradient FFT solution for scattering by planar plates," Electromagnetics, 12:105-119, 1992 Jin, J. M., J. L. Volakis and V. V. Liepa, "Fictitious absorber for truncating finite element meshes in scattering," IEE ProceedingsH, Vol. 139, No. 5, October 1992 Mayergoyz, I. D. and J. D'Angelo, "New finite element formulation for 3-D scattering problems," IEEE Trans. Magnetics, vol. 27, no. 5, pp. 3967-3970, September 1991. Volakis, J. L., J. Gong and A. Alexanian, "Electromagnetic scattering from microstrip patch antennas and spirals residing in a cavity," Electromagnetics, 1994 (to appear) Webb, J. P. and V. N. Kanellopoulos, "Absorbing boundary conditions for finite element solution of the vector wave equation," Microwave and Opt. Tech. Letters, Vol. 2, No. 10, pp. 370 -372, October 1989 8

0.15Xo l7 -. /00 d = 0.225,o - (a) Upper Aperture (b) Lower Aperture I - - - - A A finite rectangular slot array in a thick metallic ground plane (a) three dimensional perspective view (b) side view of an isolated array cell

Boundary Integral (BI #1) '. ------ PEC PEC PEC ~ I........ - - - '. * * I.,, l,*. PEC PEC PEC 3I #2 (a) BI,. i.....,........- I r I. -., ~ - - ~r - -.- I.!; %?, II I I:.... i 0.3 0.3k Yo Yo ABC surface on which ft x H = -YoE, + -V x [i(V x E),] + -VO(V I E,) 2k0 2k2 ko 2k BI - *, '..106eo"-,% 11 - - 9 I PEC L rPEC PE* C PEC PEC 0-O.15X 0 15X - - 0.15A 0K/ /, o.1 5 (c) Absorbing Material (cE = /, = 1 - j2.7) PEC Figure 2: Illustration of the mesh termination schemes used in connection with each of the three finite element formulations a) Boundary Integral Termination (FEM-BI) b) ABC Termination (FEM-ABC) c) Artificial Absorber Termination (FEM-AA)

"Cl) qb C.) 30 20 10 0 -10 -20 -30 -40 -50 (a) 0:30 60 90 120 150 180 Observation Angle a [deg] Co eb U Co 30 20 10 0 -*0 -20 -30 -40 -50 (5) 0 30 60 90 120 150 180 Observation Angle a [deg] Figure 3: Bistatic RCS of a 2x2 slot array whose element geometry and sepa:ration distance are given in Figure 1. The excitation is a plane wave and the pattern is taken in the xz plane (a) orGo RCS with (ifnc = 0~ Oi - 80~) (b) cra0 RCS with (0inc = 0~, inc = 300)

5.0 40 30 20 10 0 -10 ri3 — 1 q~ b 0 V -30 -40 -50 0 0 FEM-BI FEM-ABC FEM-AA CAVITY Unk./CPU 50,400/5,624s 155,112/18,260s 133,998/4,946s I m l I I I 1 I.. I I I I I I I I I I 0 30 60 90 120 150 180 Observation Angle a [deg] Figure 4: Bistatic a00 RCS of the array shown in Figure 1 in the $ = 30~ cut for a plane wave incidence at (qinc = 30~, 0inc = 400). The array has 15 elements in the x-direction and 2 elements along the y-direction

rA O 2) T3 0:D z 175 150 125 100 75 50 25 l I - FEM-BI " -- FEM-ABC - - - FEM-AA. —. - I J.0 0 11 101, J, J. I 11 "I / I" I, 000, I I I, oll d, 41 0 5 10 15 20 Number of Array Elements, N, Figure 5: Comparison of the number of unknowns for the three termination schemes vs slot array size (see Figure 1). For this comparison, the array has NAx2 elements

rC, 0 o o o o,so., 100000 10000 1000 -'.. -— l -- FEM-BI -- FEM-ABC --- FEM-AA 0 5 10 15 20 Number of Array Elements, NX Figure 6: Comparison of the CPU time required by each of the finite element formulations to compute the scattering for one incidence angle vs array size (see Figure 1).The array is shown in Figure 1 and has Nx2 elements