SECURITY CLASSIFICATION OF THIS PAGE (When Data Entered) _____ REPORT DOCUMENTATION PAGE | _READ INSTRUCTIONS REPORT DOCUMENTATION PAGE BEFORE COMPLETING FORM 1. REPORT NUMBER 2. GOVT ACCESSION NO. 3. RECIPIENT'S CATALOG NUMBER 018803-1-T 4. TITLE (and Subtitle) 5. TYPE OF REPORT & PERIOD COVERED Scattering by Resistive Plates Technical 1 Oct 80/30 Sept 81 6. PERFORMING ORG. REPORT NUMBER 7. AUTHOR(s) 8. CONTRACT OR GRANT NUMBER(s) M. Naor and T.B.A. Senior F30602-78-C-0148 9. PERFORMING ORGANIZATION NAME AND ADDRESS 10. PROGRAM ELEMENT. PROJECT, TASK The University of Michigan AREA & WORK UNIT NUMBERS Radiation Laboratory, ECE Department Ann Arbor, MI 48109 1. CONTROLLING OFFICF hAM- ^" 12. REPORT DATE Rome Air Development Center/01 9 Griffiss Air Force Base 13. NUMBER OF PAGES New York 13441 112 14. MONITORING AGENCY NAME & ADDRESS(if different from Controlling Office) IS. SECURITY CLASS. (of this report) Southeastern Center for Electrical Engineering Education Unclassified 11th & Massachusetts, St. Cloud, Florida 32769 s.a. DECLASSIFICATION/DOWNGRADING SCHEDULE 16. DISTRIBUTION STATEMENT (of this Report) Approved for public release; distribution unlimited. 17. DISTRIBUTION STATEMENT (of the abstract entered in Block 20, if different from Report) 18. SUPPLEMENTARY NOTES The findings in this report are not to be construed as an official Department of the Air Force position, unless so designated by other authorized documents. 19. KEY WORDS (Continue on reverse side if necessary and identify by block number) Scattering Computed data Electromagnetic waves RCS reduction Metallic plates Resistive plates 20. ABSTRACT (Continue on reverse side If necessary and identify by block number) The problem of electromagnetic scattering by resistive plates has been solved numerically using two coupled E-field integral equations and the method of moments. The Stratton-Chu representation for the electromagnetic field has been shown to be applicable to this type of problem and integral equations are derived from this representation. The numerical procedure developed guarantees the continuity of the potentials on the plate, a property which is shown to be essential for the convergence of the solution. DD JAN 73 1473 EDITION OF 1 NOV 6S IS OBSOLETE SECURITY CLASSIFICATION OF THIS PAGE (When Data Entered)

20. Abstract The solution was tested by comparison with measured RCS data of resistive and perfectly conducting plates and the agreement for plates of about one square wavelength' in area was within 1 to 2 dB for all angles of incidence. The low frequency behavior of a resistive plate has been compared to that of a perfectly conducting one based on the leading terms of the Rayleigh series. For a resistive plate the electric dipole moment is the same as for a perfectly conducting one but the magnetic dipole moment is zero. The transition to perfect conductivity is therefore discontinuous at low frequencies. The effect of the resistivity appears inthe quadrupole term, and this can be expressed in terms of potentials similar to those required for the dipole contribution. The dipole moments for plates of several shapes have been computed by solving the integral equations for the static potentials.

SCATTERING BY RESISTIVE PLATES Technical Report by M. Naor and T.B.A. Senior Radiation Laboratory Department of Electrical and Computer Engineering The University of Michigan Ann Arbor, Michigan 48109 Prepared for Southeastern Center for Electrical Engineering Education 11th & Massachusetts Avenue St. Cloud, FL 32769 Contract F30602-78-C-0148 September 1981

TABLE OF CONTENTS Page ACKNOWLEDGEMENT i ii ABSTRACT iv LIST OF TABLES v LIST OF FIGURES vi LIST OF APPENDICES viii LIST OF SYMBOLS ix CHAPTER I. INTRODUCTION 1 1.1 Background 1 1.2 Outline of the Work 3 1.3 The Resistive Boundary Condition 6 CHAPTER II. THEORETICAL BACKGROUND 10 2.1 Representations of the Electromagnetic Field 10 2.2 Numerical Solution of Integral and Differential Equations by the Moment Method 15 2.2.1 General Procedure 15 2.2.2 Subsectional Bases and CO[fJ Continuity 16 2.2.3 Subsectional Bases and C1[o] Continuity 21 CHAPTER III. SCATTERING BY RESISTIVE PLATE-DYNAMIC FORMULATION 22 3.1 Integral Equations 22 3.2 Uniqueness of the Solution 24 3.3 Numerical Solution 26 3.4 The Scattered. Field and Radar Cross Section 31 3.5 Computer Results and Comparison with Experiment 35 CHAPTER IV. LOW FREQUENCY SCATTERING BY PLATES 50 4.1 The Rayleigh Series 50 4.2 The Scattered Far Field 51 4.3 Perfectly Conducting Platelets 54 4.4 Resistive Platelets 58 -i

Page CHAPTER V. NUMERICAL SOLUTION OF THE STATIC INTEGRAL EQUATIONS 65 5.1 The Electrostatic Potential 65 5.2 The Magnetostatic Potential 68 5.3 Selected Numerical Results 69 CHAPTER VI. CONCLUSIONS AND RECOMMENDATIONS 83 APPENDICES 87 BIBLIOGRAPHY 100 -ii

ACKNOWLEDGEMENT The material contained in this report was also used as a dissertation submitted by the first author to the University of Michigan, in partial fulfillment for the degree of Doctor of Philosophy, Department of Electrical and Computer Engineering. -iii

ABSTRACT The problem of electromagnetic scattering by resistive plates has been solved numerically using two coupled E-field integral equations and the method of moments. The Stratton-Chu representation for the electromagnetic field has been shown to be applicable to this type of problem and integral equations are derived from this representation. The numerical procedure developed guarantees the continuity of the potentials on the plate, a property which is shown to be essential for the convergence of the solution. The solution was tested by comparison with measured RCS data of resistive and perfectly conducting plates and the agreement for plates of about one square wavelength in area was within 1 to 3 dB for all angles of incidence. The low frequency behavior of a resistive plate has been compared to that of a perfectly conducting one based on the leading terms of the Rayleigh series. For a resistive plate the electric dipole moment is the same as for a perfectly conducting one but the magnetic dipole moment is zero. The transition to perfect conductivity is therefore discontinuous at low frequencies. The effect of the resistivity appears in the quadrupole term, and this can be expressed in terms of potentials similar to those required for the dipole contribution. The dipole moments for plates of several shapes have been computed by solving the integral equations for the static potentials. -iv

LIST OF TABLES Table Page 5.1 Normalized Dipole Moments for the Disc 73 5.2 Dipole Moments of Various Plates 78 -v

LIST OF FIGURES Figure Page 1.1 Geometry and Coordinate System for a Plate 4 2.1 Local Coordinate System Attached to the Edge 13 2.2 (a) The Nodes of a Rectangular Cell; (b) A Corresponding General Quadrilateral Cell 20 3.1 The Interwoven Cells Scheme 27 3.2 Normalized RCS; E-Polarization II Long Edge 36 3.3 Normalized RCS; H-Polarization II Long Edge 37 3.4 Normalized RCS; E-Polarization II Short Edge 38 3.5 Normalized RCS; H-Polarization II Short Edge 39 3.6 Bistatic Scattering; Grazing Incidence and E-Polarization 42 3.7 Bistatic Scattering; 45 Degrees Incidence and E-Polarization 43 3.8 Bistatic Scattering; Normal Incidence and E-Polarization 44 3.9 Current Distribution Along AA'; Grazing Incidence 45 3.10 Current Distribution Along AA'; Normal Incidence 46 3.11 Current Distributions Along BB' and CC'; Grazing Incidence. 47 3.12 Edge Currents at Grazing Incidence 49 5.1 Convergence of P 71 xx 5.2 Convergence of Mzz 74 5.3 Mzz as a Function of 6 for Fixed N 75 -vi

Figure Page 5.4 The Various Shapes Used in Table 5.2 77 5.5 Convergence of the Dipole Moments to the Static Limit 82 B.1 Notation for Two Adjacent Cells 93 C.1 Geometries for the Calculation of the Matrix El ements 97 -vii

LIST OF APPENDICES Appendix Page A Derivation of the Representation from the Definitions of the Dyadic Green's Functions 88 B Proof of C [Q] Continuity of the Subsectional Basis 92 C Evaluation of the Matrix Elements 96 -viii

LIST OF SYMBOLS a Unit vector in the direction of the incident electric field. A Area of the plate. A = /n^H G dS' Electric vector potential. -* ^ - A = /n^E G dS' Magnetic vector potential. b = k^a Unit vector in the direction of the incident magnetic field. c A constant. csch-lx Hyperbolic cosecant of x. C Edge of S. C~[0] Class of functions which are continuous in the domain a. C'[2] Class of functions which, together with first order partial derivatives, are continuous in Q. E Total electric field. Es Scattered electric field. Einc Incident electric field. f A function. Gm i - Magnetic dyadic Green's function. 1 elr'r G = 4 Ir - r Free space Green's function. G = 1 Static free space Green's function. 0 4 i r l' Ge Electric dyadic Green's function. -ix

Gm Magnetic dyadic Green's function. h Thickness of a resistive sheeto H Total magnetic field. HS Scattered magnetic field. Hinc Incident magnetic field. i = -T Imaginary unit. i An integer index. I Unity dyad. j An integer index. J = n^H Induced electric current. -* J = n^E Induced magnetic current. k = 2Tr/X = w/J Wave number. k A unit vector in the direction of propagation of the incident field. m An integer index. mo Magnetic dipole moment. n An integer index. n Outwards pointing unit normal to S. p An integer index. pO Electric dipole moment. q An integer index. r A spherical coordinate. r Field point. r Unit vector in the direction of the observation point. Source point. -x

rect( - ) = ( Ix-xI< Ax2 Rectangular pulse function. O otherwise R Resistivity of the plate. Rp p-dimensional real Euclidean space. sinc x = sin x/x The sinc function. sinh-1 x = ln(x+Tl+x2) Arc Hyperbolic sine of x. S Surface of the scattering body. u A vector function. V Volume included inside S. w A weighting function. x A Cartesian coordinate. x Unit vector in the x-direction. Xm x coordinate of a sampling point. X A function space. y A Cartesian coordinate. y Unit vector in the y-direction. Yn y coordinate of a sampling point. Y A function space. Y = 1/Z Admittance of the medium. z A Cartesian coordinate. z Unit vector in the z-direction. Z = /7/ Impedance of the medium. l, = n %n {,n m # n Kronecker delta function. mn O,m n 6(x-x ) Dirac's delta function. -xi

AX Size of a rectangular cell, in the x-direction. Ay Size of a rectangular cell, in the y-direction. c Permittivity of the medium. C A cylindrical coordinate. n Normalized impedance. n Y-Yn e A spherical coordinate. e The complement of the spherical coordinate of the direction of incidence. oS The complement of the spherical coordinate of the direction of observation. K An integer or semi-integer index. X Wavelength. X An integer or semi-integer index. Permeability of the medium. An integer or semi-integer index. v An integer or semi-integer index. X-X p A cylindrical coordinate. p = cn-E Induced electric charge. * ^ p = pn-H Induced magnetic charge. o~~~a ~Conductivity. Cylindrical coordinate. A spherical coordinate. ^ -: = f n.E G dS' Electric scalar potential. -xii

Spherical coordinate of the direction of incidence. (s Spherical coordinate of the direction of observation. * ^= fn.H G dS' Magnetic scalar potential. w Frequency of the time harmonic field. n~~Q ~A domain CRP. v2 = a2/aX2 + 32/ay2 + a2/az2 The Laplacian operator. v2 = a2/ax2 + a2/ay2 The surface Laplacian operator. 1+ Difference of the values on the top and bottom faces of a plate. -xiii -

CHAPTER I. INTRODUCTION 1.1 Background The resistive boundary condition was devised in order to simplify the solution of boundary value problems involving thin lossy dielectrics,and Harrington and Mautz (1975) have shown that it provides a good model for such materials. Practical examples of "resistive materials" range from man-made lossy composite materials used, for example, for the fin of a missile to reduce its radar cross section, to snowflakes and smoke particles. A precise definition of the resistive boundary condition will be given later in this chapter. A number of dielectric structures have been simulated by this boundary condition. Infinite structures such as the half plane, strip (Senior, 1979a,b) and wedge (LaHaie, 1981) have been treated analytically, using transform techniques. Closed shells (Harrington and Mautz, 1975; Senior and Weil, 1977) were treated by numerical methods. There have been no previous numerical analyses of finite, open resistive structures because of the inherent numerical difficulty associated with such structures. This difficulty is common to any finite, open, thin structure, resistive or perfectly conducting alike. It arises because only the highly singular E-field integral equation is applicable and the situation has been aggravated by the common belief that only E-field formulations derivable from the -1

-2Franz representation are valid (Stratton, 1941; Bouwkamp, 1947). The various representations and their numerical implications will be discussed in Chapter II. Because the source of numerical difficulty in the treatment of plates is the same for all boundary conditions, it is helpful to review the techniques used for treating perfectly conducting plates. Mittra et al. (1973) and Rahmat-Sammi and Mittra (1974) suggested a way of partially integrating an E-field equation, thus reducing the order of its singularity. Their method is a two-dimensional analog of Hallen's integral equation for wire antennas. Unlike the onedimensional case, the general homogeneous solution of the differential operator is not known. Mittra et al. proceeded by assuming that the unknown general solution is expandable in a Fourier series and determined the coefficients of the expansion by imposing an edge condition along the rim of the plate. Their method works for perfectly conducting plates, but does not extend to resistive ones. Recently, Glisson and Wilton (1980) came up with a different method. They used, but without justification, an E-field integral equation based on the Stratton-Chu representation. As will be seen in Chapter III, their numerical scheme guarantees both continuity of the vector potential and satisfaction of the edge condition. Their method is a breakthrough from a numerical viewpoint, but requires a sounder theoretical basis. The Babinet complement of a perfectly conducting plate is an aperture in a metallic screen. Because of its great importance, the problem of diffraction by an aperture in a plane screen has been

-3investigated by many authors. Most treatments involve some kind of an approximation, either physical optics (cf. Bouwkamp, 1953), the geometrical theory of diffraction (Keller, 1962) or low frequency expansions. The last method is an example of Rayleigh scattering (Lord Rayleigh, 1897). A systematic procedure for treating diffraction by small obstacles was developed by Stevenson (1953). In the many papers that followed, a variety of low frequency diffraction problems were solved. Van Bladel (1977) provides a comprehensive list of the work published through 1975. Integral equations for the fields of small apertures were developed by Bethe (1944). Van Bladel (1970) and his colleagues (DeMeulenaere and Van Bladel, 1977; De Smedt and Van Bladel, 1980) attempted to solve them numerically, but with only limited success. Butler (1974) integrated those integral equations in a manner analogous to that used by Mittra, and numerical results have been reported by Umashankar and Butler (1974). A different approach based on a finite element formulation was used by Okon and Harrington (1980). 1.2 Outline of the Work The task is to determine the fields scattered by a resistive plate S, lying in the x-y plane (Fig. 1.1), when illuminated by a plane wave: -inc a ikk. r -in Yb ikk (1.: e (1.1)

-4b a S y Fg1.1Goeran X Fig. 1.1 Geometry and Coordinate Systems for a Plate.

-5-iwt The time convention e is assumed throughout. The unit vector A A A.A k = -(x cos cos + y cos i sin i +z sin i) (1.1) is in the direction of propagation of the incident wave. a, b and k form a right-handed orthonormal set. The resistivity R of the plate is in general complex and a function of position. Its value can be allowed to go to zero (the perfectly conducting limit). The problem is solved numerically by applying the moment method to a suitable integral equation for the induced currents. Because of the finite capacity of any computer, the largest plate that can be treated at present is about a square wavelength in area. Although this is not large enough to permit direct comparison of the results with those obtained by asymptotic methods such as the geometrical theory of diffraction, plates of this size are of practical interest and display similar characteristics to larger plates. Special attention is paid to the low frequency limit, where semi-analytic results can be derived from the Rayleigh series. The first two terms of the expansion are determined, and efficient methods are developed for the solution of the associated integral equations. The resistive boundary condition and its relation to other types of boundary conditions are discussed in the remaining section of this chapter. The first half of Chapter II is devoted to a discussion of the representation of an electromagnetic field. A proof

-6of the adequacy of the Stratton-Chu formulation for a plate is provided there. The method of moments and some of its ramifications are then presented, and these two topics lay the foundation for the numerical solution of the problem in Chapter III. Several integral equations are derived and evaluated. One is selected for numerical treatment on the basis that it provides the fastest convergence. Results are compared with experimental data and their implications are discussed. Chapters IV and V deal with the low frequency limit. The analytical derivations are carried out in Chapter IV and the numerical solution in Chapter V. Comparison with other methods and some experimental and theoretical bounds for the dipole moments are also covered in the latter chapter. Concluding remarks and suggestions for the extension of this work are reserved for the last chapter. 1.3 The Resistive Boundary Condition A resistive sheet is characterized by three properties. It is infinitesimally thin, carries only electric currents and these are proportional to the tangential components of the electric field. Mathematical expressions for these conditions are: J n^E_ = 0 (1.3) and n xE = Rn^J = Rn^(n^H). (1.4)

-7The proportionality factor R is related to the material properties by 1 iZ R = lim h lim k(- l)h (1.5) a 0-o r -+ oo h o h+ o where h is the thickness and a the conductivity of the material. In practice, this boundary condition is applied to nonmagnetic (p = p0), highly conductive materials whose thickness h is very small compared with the wavelength x and the penetration depth 6. The material is thus partially transparent. The value R = 0 corresponds to perfect conductivity. Equation (1.5) does not determine the frequency dependence of R. Empirically, the resistivities of resistive sheets remain constant over a broad range of frequencies, and in fact are measured at dc. Thus it can be assumed that the resistivity R is independent of k, at low frequencies at least. The resistive boundary condition is related to two other types of boundary conditions: the conductive and the impedance boundary conditions. The conductive boundary condition is the dual of the resistive one. A conductive sheet supports only magnetic currents ^ i (J = n^HI = O), which are related to the tangential components of the magnetic field by

-8A m A m A A + n^H = R*n^J* = -R*nA(n^E)j and R* = lim k(Ur- l)h h + o This is a fictitious kind of material. The only practical importance of the boundary condition is that combined with a resistive sheet it produces an impedance sheet. An impedance sheet is an opaque (i.e., h >> 6), thin (h ~< x) layer supporting both electric and magnetic currents. These currents are related to each other, and on each face of the sheet n^(n^E) = -nZn^H. (1.6) For a planar sheet this relation reads E = +ZH, E = +nZH x y y x The upper and lower signs refer to the top and bottom faces respectively. By writing the resistive and conductive boundary conditions for planar sheets in the forms E + = -2R(Hy -, Ey + Ey = 2R(H - x) E E = 1 + + Hy) + x 2 (Hy +H,) Ey -E =R x

-9of magnetic and electr and assuming coexistence of (1978) obtained - = - R + 4R Hy + - 4R Hx+ | lL- H- x hese relatins turn into the impedance boundary condition if =2R /2R*d on an impedance Thus the electric and magnetic currents induce e the same as those on a resistive and a conductive sheet sheet are the same a s Petos ve with R = (1/2)nZ and R* = /(resistive plate the w oed forthe currents on a r Having soVd onductive plate ca be. magnetic currents on an con e cne obtained by duality. The impedance plate is then the superposition f both types of current sheets. of both t..Ives

CHAPTER II. THEORETICAL BACKGROUND 2.1 Representations of the Electromagnetic Field Integral representations are powerful tools in the theory of electromagnetism. They provide a means for the derivation of integral equations and of asymptotic low- and high-frequency approximations. Several representations are available. Although they will be shown to be equivalent, each possesses certain properties that make it advantageous for some applications and not for others. A couple of cases will demonstrate this point. The total field is the sum of the scattered and incident fields. Since the latter is given, only the scattered field will be considered. The most widely used representation is due to Stratton and Chu (Stratton, 1941; p. 464). It expresses the scattered field in terms of the surface values of both the electric and magnetic fields on the surface S of a scattering body: E = ikZ f n'H G dS G dS' -+ v n'-E G dS' (2.1) S S S and H = -ikY n'^E G dS' - n'-H G dS' + v^ n'^H G dS'. (2.2) S S S -10

-11^ - ^ ^ ^ - Identifying n x H, n x E, n-E and n-H with J, J*, p and p* respectively, the various integrals in (2.1) and (2.2) are recognized as vector and scalar potentials. Hence the fields can be written as: E = ikZA - + v^A* (2.3) and Hs = -ikYA* -V* +VA.A (2.4) These two formulations are entirely equivalent. The charges and currents are related and the relations are given by the continuity equations: V-J = ikY (2.5) and v-J* = -ikZ. (2.6) These equations supplement (2.1) and (2.2), or (2.3) and (2.4) to form the complete representation. Another formulation attributed to Franz (Jones, 1964) requires knowledge of the induced currents alone: E = (k2 nH+ vv.) fn n' H G dS', (2.7) S S s = -iY (k2 + vv-) n' E G dS' + v n'AH G dS'. (2.8) S S

-12The terms of this representation are electric and magnetic Hertz potentials (linearly related to the vector potentials). An unconventional derivation of the representation which avoids the use of Green's theorem is presented in Appendix A. It is well known that the fields represented by the Franz formulation are Maxwellian. Stratton-Chu fields, however, satisfy Maxwell's equations if and only if they are equivalent to those of Franz. By applying Stokes' theorem and Maxwell's equation v^H = -ikYE to / V'^[HG] -n dS' it is. found that (ibid) vv ~ n'^H G dS' = -ikYv n'iE G dS' + v'G H-dc. (2.9) S S C An analogous result holds for the dual case (E H, H -E, Y + Z). It is thus evident that the two representations are equivalent when the line integrals vanish. The task is therefore to show that the line integrals are indeed zero. The case of a smooth closed surface S is trivial. The closed path C degenerates into an arbitrary single point of S and the line integrals disappear. Investigation of the behavior of the fields in the vicinity of an edge establishes the same result for surfaces with edges. Attaching a local cylindrical coordinate system p,,c (Fig. 2.1) and imposing the condition that the energy stored in the neighborhood of the edge is finite, the fields are seen to be 0(pq-l) with q > 0.

-13FLryatp Fig. 2.1 Local Coordinate System Attached to the Edge.

-14The p and q components of Maxwell's equation v x H = -ikYE are: aH aH 1 ai a- = -ikYE (2.10) P P and aH aH d ~ = -ikYE (2.11) If all components of E and H are O(pq ), then the first term on the left hand side of (2.10) and the second term on the left hand side of (2.11) are O(pq-2) Since the other terms are (pq-1) it follows that Hz lim Z 0 (2.12) p - a0 3H lim 0o = ~. (2.13) P o The other possibility is that Hz is O(pq) in which case lim Hz = 0. p + o Equations (2.12) and (2.13) have two distinct solutions: one has q t 1, implying lim Hz = 0 (the other components can be singular), P - o and the other has Hz independent of both p and q, implying q = 1. In both cases, Hz is continuous at the edge, independently of the boundary condition that holds on S. In particular, if S = S U S 1 2 where S and S are (open) smooth surfaces that meet at an angle 1 2 along C

-15V' (HG) * n dS' = f( H - I )G - dc (2.14) 1 2 holds. H and H are the values of H on S and S respectively. 1'2 1 2 Since H is independent of q the integral on the right-hand side M of (1.14) vanishes. This result can be generalized when S = uI S m=l m and S has several edges. The continuity of H also guarantees the disappearance of the line integrals in the case of an aperture in a perfectly conducting screen. The integral Jv'GH *-d vanishes on the screen, hence it C vanishes on the aperture side of the edge also. In the next section it will be shown that the choice of representation has a profound effect on the convergence properties of a numerical solution. 2.2 Numerical Solution of Integral and Differential Equations by the Moment Method 2.2.1 General Procedure. The method of moments (cf. Harrington, 1968) was devised to solve numerically equations of the form: Lf = g; in C RP, (2.15) where L is an operator mapping the unknown function fC X (over Rq) into the known function g EY (over RP). An inner product is defined in Y as <u,v> = f uv dxP; Vu,v C Y. (2.16)

-16A set {w}M CY of weighting (or sampling) functions is selected. m1 Taking the inner product of the given equation with wm we get <wm,Lf> = <w,g>; 1 < m < M The function f is then expanded in terms of a basis {f }' of X n1 =f anfn (2.17) n=1 and the series is truncated at some finite N. The result is a set of M algebraic equations for the N unknown coefficients: N <WmLfn> a = <wmg> (2.18) n=1 or, in matrix notation, Wa = j, (2.19) where Wmn = <wm,Lfn> and gm = <wm,g>. The system of equations (with M = N) is solved by Gaussian elimination and backwards substitution. 2.2.2 Subsectional Bases and C~[Q] Continuity. The choice of the basis {fn} in which f is expanded is a crucial step in the solution. Obviously a set of basis functions which resemble the required solution, that is, satisfy the same boundary conditions and have the same continuity properties, will lead to rapid convergence,

-17and only a few terms of the expansion are then needed to describe f. Generally it is difficult to "guess" a good basis, and it is an almost impossible task for a domain a of complex geometry with complex boundary conditions. It is much easier to divide Q into subdomains QN such that N U QN (2.20) n-= and use basis functions which have support on an: n; r 62Qn fn = (2.21) nn 0;r n Such a basis is called a subsectional basis. Quite often the nature of the operator L dictates that f possess certain continuity properties. A subsectional basis may lead to a representation of f which does not possess the required continuity, and thus results in poor convergence. For a planar domain Q C R2 there is more than one choice of a subsectional basis that guarantees Co[Q] continuity. One possibility, which will be extensively used in this work, uses rectangular cells:

-18Ax AXAYA ~mn (xY) Ix -7 < X < x +, y < 2 < Yn + "mn - x'Y) xm 2' m 2 Y - n 2 (2.22) and basis functions that are the tensor product of two quadratic polynomials in i and n respectively, i - x - xm, n - Y - Yn: a + a + a n + a in +f 2 1 2 3 4 5 + a nl2 + a 22n + a En2 + a 2n2; (x,y) Qn 6 7 8 9 mn fmn. (2.23) 0; (x,y) ~ mn The term in E2n2 may eventually be dropped, but it is more convenient to retain it for the time being. The coefficients ai (which depend on m,n) are selected so that f attains its exact value at nine predetermined points of mn' It is shown in Appendix B that if these points are taken to be the four corners of Qmn' the four midpoints of its sides and its center, then N f Z f mnE Cc[ n=i

-19With this choice of nodes in &mn (Fig. 2.2a) f (.x + AX y f (x X A f(Xm+ )- f(xm - 2,) af f(xm=Yn + f( ) (m2- ('Yn - 2 ) a2fl = 2 (X + 2 ) + - / n) - 2f(xy,) y ( AX6 Ay2 (n,yn) (2.24) and the right-hand sides are finite difference representations of these derivatives, numerically equivalent to (2.17) with (2.23). The expressions given by (2.24) are particularly useful when f is a weighted integral of another unknown function to be determined. The choice of sampling points is crucial. The dramatic effect of loss of continuity as the points 1,2,3, and 4 are moved away from the boundary of nmn will be demonstrated when numerical results are discussed.

-20- I~ (a) ( Fig. 2.2 (a) The Nodes of a Rectangular Cell; (b) A Corresponding General Quadrilateral Cell.

-21Before concluding this section, it is appropriate to remark that the rectangular cell Qmn can be mapped into a general quadrilateral cell. A quadratic transformation, based on the same nine points, is a very convenient means to effect the mapping at least approximately (Fig. 2.2b). The applicability of the method is thus extended to more general geometries. 2.2.3 Subsectional Bases and Cl[Q] Continuity. Some operators, such as a2/axay, impose higher continuity constraints on f, and this particular operator requires C1[Q] continuity. It is much more difficult to achieve this kind of continuity with subsectional bases. Both f and its partial derivatives must be specified at the nodes of Qmn. Even then, no simple polynomial expression which ensures continuity of f and its first derivatives is possible (Zienkiewicz, 1977, Chapter 10). Many attempts have been made to overcome this difficulty, with different degrees of success [ibid, Chapters 10-12]. For mechanical problems it has been recognized that the best approach is to reformulate the problem so that only C [Q] continuity is required. The original equation (2.15) is replaced by a set of equivalent equations involving operators of lower order. It turns out that it is more economical to increase the number of equations than to use more parameters per cell.

CHAPTER III. SCATTERING BY RESISTIVE PLATE-DYNAMIC FORMULATION 3.1 Integral Equations Integral equations for the induced currents are derived from the field representations by taking the observation point to lie on the surface S of the scatterer and applying the boundary conditions there. Several E-field formulations are possible but the H-field equations becomes an identity for an infinitely thin open body. Selection of the preferred equation for numerical treatment is based on the conclusions of the preceding chapter. From a numerical standpoint, the equation that imposes the lowest order of continuity requirements is the best. One equation is deduced from the Franz representation using the boundary conditions (1.3) and (1.4): J ^ Ainc z.(k2 + vv.) J G dS'. (3.1) S The dyadic operator vv contains second order mixed derivatives such as 32/axay and requires Cl[S] continuity of the vector potential. This is hardly encouraging and the need for other formulations is evident. Operating on one Cartesian component of (3.1) with a2/ax2 + k2 and with a2/ay2 + k2 on the other, subtracting and simplifying -22

-23the incident field term using all three relations -inc inc inc v. E =, v x i ikZH and (V2 + k2)-inc (2 + k2)E = 0 leads to an equation for one component of J. Repeating the whole process with the two differential operators interchanged results in an equation for the second component. Combined together they provide the integral equation: ^ inc ikZ(v2 + k2) G dS' = ikZz + R[k2 + zvV(z.v^)]J. (3.2) This is valid only if R is constant, but the equation has even more serious shortcomings. C1 continuity of 3 itself is now required. Inc The equation fails at grazing incidence since aH //z introduces a factor sin e i into the excitation term. The transition to a perfectly conducting plate is also discontinuous and when R = 0 there is no coupling between the two components of the currents. Satisfaction of the boundary condition z.H = 0 (on S) is no longer guaranteed and must be imposed separately. The only advantage this formulation may offer is the possibility of partial integration. When R = 0 this equation was integrated by Mittra et al., (1973) by using Green's theorem in a plane and the two-dimensional "free space" Green's function. Their method is capable of removing the differential operator v7 + k2.

-24^ A However, it is impossible to integrate z^v(z-v^) at the same time, and thus, if R f 0 we are still left with the requirement for C1 continuity of J. A third set of equations can be derived from the Stratton-Chu representation. The boundary conditions (1.3) and (1.4), when substituted into (2.1), produce the integral equation: Rz J z En + ikZz^J G dS' - zAv P G dS' (3.3) subject to the constraint V ~ J C= * ( (3.4) z E The highest order continuity required is CO[S], and the burden of continuity is now shared between the scalar potential (in 3.3) and the current (in 3.4). The price paid is the increased number of equations that have to be solved simultaneously. In a later section it is shown how that difficulty is circumvented. Having a valid and convenient formulation, it is feasible now to discuss such matters as uniqueness of the solution and numerical procedures. 3.2 Uniqueness of the Solution Jones (1974) has examined the question of uniqueness of the solution in the case of a perfectly conducting solid body. He showed that E-field integral equations have more then one solution

-25at wave numbers which are eigenvalues corresponding to the TE modes of a cavity, identical in size and shape to the body and having perfectly conducting walls. Regarding the plate as a limit of a hollow right prism of altitude h, and repeating Jones' argument, we see that if the equations (3.3) and (3.4) have more than one solution, the differences J,p of these solutions satisfy Rn^J = ikZn^ J G dS' - n.v PG dS' and ik v = J =i Z c These are the equations satisfied by the eigenmodes of the cavity. All the eigenvalues of the equations are complex, but even so the accuracy and stability of the solutions for k real may be seriously affected by the presence of a nearby complex pole. Every mode of these equations can be put into a one-to-one correspondence with a TE mode of a perfectly conducting cavity of the same geometry. As the altitude h goes to zero all TE modes disappear. Hence the solution of (3.3) and (3.4) is unique for any value of R, including R = 0.

-263.3 Numerical Solution Equations (3.3) and (3.4) are solved by the moment method. Written in scalar form, these equations are: RJ = Exn + ikZ JG dSG dS' = n iZA X: x - ax (3.3a) RJ = En +ikZfJG dS' - f G dS' = Enc +ikZAy - y yJ y: y ay (3.3b) P - Z -~x + 1 (3.4) piZ x a3 ay If the finite difference expression (2.22) is to be applied to (3.4), it is clear that the sampling points of J and Jy must be a y shifted by Ax/2 and Ay/2, in the x and y directions, respectively, relative to the sampling point of p. Bearing in mind that it is inadvisable to sample p on the edges of the plate (where p is singular if R = 0), and that Jx and Jy are zero on edges parallel to the y- and x-axes, respectively, one arrives at Glisson and Wilton's (1980) scheme of interwoven meshes (Fig. 3.1) for a rectangular plate. A corresponding expansion of the unknown functions is M-1 N J: I ) relct ~2. -rect rect 1 (3.5) x X p +1/2,q,qA + 1.J

Yn ~ ~ ~ ~I+l 12 LIC~ -+ -~I) — Yn+1/2 _____ ____~~~~~~~____ t t I -t t t _ ~I t I I f I yl" __________ ^^~~ x m \ m+1/2 M Fig. 3.1 The Interwoven Cells Scheme.

-28M N-1 Y =Yp,q+/2 rectL Ax rect Ay ) ) p rect -p - rect ] 1q (3.7) p=l q=i L L where i J J =0 1 /2q XM+1/2,q Vp,1/2 Yp,N+1/2 from the edge condition. Using sampling functions Wx = 6(x - X+l/2( - ) Wy = 6(X - xM)6(y - Yntl/2) W = 6(x - Xm)6(y - yn) for (3.3a), (3.3b) and (3.4), respectively, and using (2.22), equations (3.3) and (3.4) are converted to i nc ax m+l,n m,n W~ ~~~ = 1 m 7~ n

-29inc Rmn+l/2Jy = +lkZA mn+l /2 Ym, n+l1/2,n+l/2 - 1 (,n - ) (3.9) Ay m,n+l m,n iZ 1 1 ^- (J -J ))+ (J+ = p0 k ax Xm+l/2,n Xm-1/2,n A Ym,n+l/2 m,n-1/2 (3.10) An expansion of the potential corresponding to (3.5) through (3.7) is M-1 N = J G dS'I I Xm+l/2,n x m+l/2,n Xp+/2q p=1 q=1 pfl/2,q Gm+1/2,n,p+1/2,q etc. where 11,v X=X y=yv and G,vK,x =S Gpv dS' K,X Using the translational symmetry of G: m+a,n+,p+a,q+: m,n,pq

-30substituting (3.10) into (3.8) and (3.9) and rearranging, the final form of the equations is obtained as: M-1 N ikYEnc J k2G + 1(G i m+l/2,n p+/2,q m,n,p,q Ax2 m+,n,p,q p=: q= X'~' ~q1/ ~y r+1Gm,n,p,q rn+1,n-1,ppq + G-np 2m-l m,n,p,q) + ikYR p+l/2,q np nq ] M N-i n " z 1 ( G - G p=l qi Pq+l/2 AxAy m+ln, p q m+ln-1, - Gm + G ) (3.11) m,n,p,q m,n-1,p,q M-1 N ikYE = n (G - G y ^ p+1/2 q ^Y mn+lp,q m-1,n+lpq M N-1 Gm,n,p,q + Gm-ln,p,q_ ) iy'1 n1 p, p=i q=i p1 q+l/2 k2G + (G +G - 2G ) m,n,p,q Ay2 m,n+l p,q m,n-l,p,q m,n,p,q + ikYRp /,ql/2 6mSn ]1. (3.12)

-31The calculation of the coefficients G is discussed in ^,B,Y,S Appendix C. The final step is to solve (3.10) and (3.11) by Gaussian elimination and backwards substitution. One point deserve further comment. Equation (3.4) was substituted into (3.3) only after it was converted into the form (3.10) by using (2.22). This way the continuity was preserved, but the number of coupled equations was reduced approximately by one third. As a check, numerical results of (3.10) and (3.11) were compared with the currents calculated from the corresponding equations with (3.10) kept as a separate equation. The results were identical to five significant decimal figures. Having determined the induced currents, it is now possible to calculate the far field. 3.4 The Scattered Field and Radar Cross Section The far field is derived from the currents using the Franz representation and the far field approximation: {r l - r' -r * r' and is ikr A A ^k ikr Es ekr r rArA JJ eikrr dS' ekr S(r) (3.13) with ikr ~s = z r 5 = -zZ r ^s(r) e (3.14) H = ZSr^ E = -Z k r^S(r) (.3.14)

-32From formulas (3.5) and (3.6), S can be expressed as: S(r) = iZ si (Fy cos - F sin c) 4) = s s y s " s - sin2 9s]+ y [cos2 e cos s (F sin s - Fy cos ) - Fy sin2 e s + z sin 9s cos es(Fx cos qs + Fy sin s) (3.15) where F e JJ ikr.r dS' = aXAy sinc x cos es cos s M-1 N sin* O sin c sin Xp+l/2,q L 2 z z X p+1/2,q exp [-ik cos es(xp+l/2 cos s + Yq sin s)] M N-1 + yt z" z Jyp^ q~l2 exp[-ik cos es(xp cos s + q+l/2 sin s) p=l q=l pq+l/2 (3.16) es is measured from grazing and ^ ^ ^ ^ r = x cos es cos (s +y cos es sin ( s + z sin O is a unit vector in the direction of the observation point (Fig. 1.1).

-33A normalized quantity corresponding to the scattered far field is the bistatic scattering cross section, defined as (Crispin and Siegel, 1968): a(r) = lim iL, (3.17) r -). ooE -inc E is measured on the surface S. For a plane wave incidence (3.17) can be written as r) = 1 ISr)12. (3.18) The value of a corresponding to r = -k is called the backscattered, or radar cross section. Equation (3.17) or (3.18) does not reveal explicitly the dependence of the scattered field and its polarization on the polarization of the incident field. There is no uniform convention as to how to describe this dependence in the general bistatic situation for an arbitrarily shaped body. For a plate, it is convenient to decompose the fields into components parallel and perpendicular to the z = 0 plane (Fig. 1.1), and to treat separately the cases of inc -inc E- and H-polarizations, i.e., En or Hi parallel to the plane and perpendicular to the direction of incidence. A general incidence is a superposition of these two cases. Further simplification ^ results from choice of the x-axis parallel to a (b), the direction of polarization of the incident electric (magnetic) field. With this

-34choice of the coordinate system, the left-hand sides of (3.11) and (3.12) are ikY exp(.ikyn cos ei) and 0, respectively. Restricting the point of observation to the plane of incidence, i.e., =.i = ~ 7T/2, the direct polarized component of the scattered field becomes actually parallel to the incident field (of a given polarization). An explicit expression for (3.18) is then = k4Z2 IFXl2 (3.19) X2 167r3 for direct polarization. The cross polarized current is directed parallel to ^ ^ - ^ ^ ^ u = r xx = y sin e - z cos a (3.20) and = - F IF 2 cos2. (3.21) x2 163 y cos2 The corresponding integrals are given by M-i N F I Jx exp(-ikyq cos es) (3.22) pi1 q= p+l/2,q and M M N-1 Fy S= c l exp(-ikyq/2 cos 0,) (3.23) p=r q=1pti,vq+ respectively.

-353.5 Computed Results and Comparison with Experiment Computations were carried out for rectangular plates illuminated by an E- or H-polarizaed plane wave, incident in one of the principal planes of the plate, i.e., the principal axes of the plate coincide with the x- and y-axes. The direction of observation was confined to the plane of incidence. Computed and available measured results for 1.16x x 0.85x perfectly conducting and resistive plates are compared in Figs. 3.2 through 3.5. The resistive plate is a composite of a 800 a/square resistive sheet sandwiched between two 3.3 x 10-3 x thick dielectric plates. The relative dielectric constant of these plates is Er = 2.53. Figure 3.2 displays the normalized radar cross section as a function of the angle of incidence e. for E-polarization with the field parallel to the long side of the plate. Figure 3.3 presents the corresponding results for H-polarization. Figures 3.4 and 3.5 give the results for E- and H-polarizations with the fields parallel to the short side. The experimental results measured over an angular range of 1800 were folded into a 900 display, revealing some asymmetries about normal incidence. Examination of the levels of the peaks at that incidence at different polarizations show inconsistencies. Nevertheless, since this was the only available set of measured data, it was used for making a comparison. With the exception of the results for the resistive plate in Fig. 5.5 and the deep nulls in

-36-,....E Measured 15 t t Computed 10 -2!0o_,-</ /1 0 - R= /1// -10 // / - R = (2.12+0.2i)Z - \/ 30 o~ 60 90 Fig. 3.2 Normalized RCS; E-Polarization II Long Edge.

-37H t I Measured 15 Computed 10 -20 / [dB]i -o. /' /'\E 111. I -10~ / Y~ -R (2.12+0.2i)Z -20 / 30 8~ 60 90 Fig. 3.3 Normalized RCS; H-Polarization I Long Edge.

'e6p3 oq4s ii uo LezLe LOd-3'S3o pazLLeLUoN t'E'6.L 06 09 oe OE 4 / Z I - II|, Z( LZ O+ZL'Z) ='d t~ill "oI \; /1 - i/ yI I i0 /~~~~/ // pans/ ep] -8/

-39Measured 15~I __1 Computed 15 ///?/. "?q R = (2.1~+0.2i)Z!, \ // -' -10- / -20 11 / / \ 1; -3,/ i 30 0~ 60 90 Fig. 3.5 Normalized RCS; H-Polarization I Short Edge.

-40all the other cases, computed results are within 2 dB of the experimental data and usually much closer. There are two reasons to question the validity of the experimental data of the particular case in Fig. 3.5. At normal incidence, the maximum level is several dB higher than the corresponding peaks at other polarizations, whereas the computed values are mutually consistent. The measured lobe near grazing is so untypical for H-polarization (theoretically the plate should be invisible to H-polarization at grazing incidence), that one should suspect a stray reflection from the experimental setup. As to the other discrepancies, all computational methods consistently predict deeper nulls than are actually measured. This is the effect of noise filling up the nulls. Otherwise, the agreement between computed and measured results is quite remarkable. It is even more so considering the fact that owing to the zero edge conditions, the numerical procedure uses an effectively smaller plate. A zero current is enforced on half cells along the two edges of the plate perpendicular to the current. For a plate divided into 14 x 10 cells it means a reduction factor of 13/14 or 9/10, depending on the direction of incidence and polarization. At normal incidence the radar cross section is roughly proportional to the square of the area. Thus, the corresponding reduction factors for the cross section are -0.65 dB and -0.92 dB, respectively. This brings the results into even closer agreement with experimental data. The effective resistivity of the composite plate was computed as (Senior, 1981):

-41Req = R - 2ikd[cos e Y- R( - 1) (3.24) where d is the thickness of the dielectric plates. The imaginary part of Req is quite small and the value of Req at normal incidence was used for all angles. There are other uncertainties associated with the values of er: its imaginary part, the thickness d and possible air gaps between the dielectric plates and the resistive sheet. It is therefore surprising that such a simple model produced such good results. At E-polarization and grazing incidence the dominant induced currents on a perfectly conducting plate are edge currents along the leading and trailing edges of the plate (Knott et al., 1971). A possible qualitative model for the direct polarized scattered field is an endfire, two element array of wire antennas along these edges. This model is supported by the data for bistatic scattering in the principal plane for grazing incidence, presented in Fig. 3.6. Results for 450 and normal incidences with the same polarization are shown in Fig. 3.7 and Fig. 3.8, respectively. In both cases, maximum scattering occurs in the specular directions. Typical current distributions for grazing and normal incidences are shown in Figs. 3.9 through 3.11. The large front and rear edge currents on the metallic plate are clearly visible in those figures. The currents on the resistive plate do not possess such singularities and their variation in the y-direction resembles that of the physical optics currents (except, of course, for an

-42R =_0 10 dB ] R = (2.12+0.2i -10 -R/ /0 -30 L / 269-20-10 20 60 90 oo 120 150 180 Fig. 3.6 Bistatic Scattering; Grazing Incidence and E-Polarization.

-4310 R= 0 -1 0 I -30 -30 60',,,,,,, 30 60 eo 90 120 150 180 Fig. 3.7 Bistatic Scattering; 45 Degrees Incidence and E-Polarization.

-44R =0 10 R = (2.12+0.2i )Z -10 / I -20 t. / -30 30 60 eo 90 120 150 180 Fig. 3.8 Bistatic Scattering; Normal Incidence and E-Polarization.

Perfectly Conducting'esistive 90 / arg J/ Deg.] / / / / F5/ / -90 x A A I I,,, I I I A 0.5 y/x 1.0 A' Fig. 3.9 Current Distribution Along AA'; Grazing Incidence.

* ouap.LouI LPUJON,VV UOLV uo.LnqLjJSLla:uaLn3 OL''.Lj,V O'L Y/~ S'O V f..... I' I, I I I I f,' ~iz VV 1 3 IV V1' x 06 -9^~~~~~~~X ~

-47BB'... arg Jx [Deg. _ CC' CC' x -90 ~'~ I B C BB' Jx 5 CC' c c x/x 1.0 Fig. 3.11 Current Distributions Along BB' and CC'; Grazing Incidence.

-48attenuation factor). The variation in the x-direction is of the same form in almost all cases; therefore only one typical distribution is presented. Two transverse cross sections are shown for the current on the perfectly conducting plate and one for the resistive plate (on which there is very little variation from one location to another along the y-axis). An important effect of the cross polarized currents is the coupling they provide between the front and rear edges. This coupling mechanism was discussed by Knott et al. (ibid). Figure 3.12 displays the behavior of the edge currents, at grazing incidence, as a function of distance along the edge. Two prominent features of this graph are the side edge current buildup and the nonsinusoidal nature of the rear edge current. These characteristics have been detected experimentally by the above-mentioned authors. On the resistive plate there is no such buildup of the side edge current, which is small compared with the currents on the leading and trailing edges. Results for scattering by plates which are small compared to the wavelength computed by the procedures described in this chapter are presented in Chapter V, after the topic of low frequency scattering has been discussed in Chapter IV.

E r~~il 1~ Metallic 10QI [' ~ ~Resistive a b dxI T-I~~~~~~~~~~~~~~~~~~~~~~~~~~J,~~~~~~~~~~~~~~~~~~~~~~~' x~~~~~~~~~~~~ y = 0.0414x \fx = 0.425x /y = 1.1186X i ~ ~ ~ ~~~ - i- - -. 4 1 -..,-... 4-4 _ ^ i' i " 7 o a b' Fig. 3.12 Edge Currents at Grazing Incidence.

CHAPTER IV. LOW FREQUENCY SCATTERING BY PLATES 4.1 The Rayleigh Series The advantage of using asymptotic expansions is that many analytic and semi-analytic results can be obtained for scatterers of arbitrary shape. The Rayleigh series is such an expansion, valid in the low frequency regime. The expansion parameter is kD, where D is a characteristic dimension of the body. The series 00 E = (ik)m E (4.1) m=o and the corresponding expansion for H are convergent in the vicinity of the body when kD << 1. The individual terms of the Rayleigh series (4.1) satisfy the reduced Maxwell's equations ^E = 0, V^Ho = O0 VEm ZH, VHm = -YEm; (m > 1) Zm -1i - m i - _ v-Em = 0, vH = 0; (m >0). (4.2) Stevenson (1953) has offered a systematic method for the determination of these partial near fields. The essence of his method is that the terms are obtained successively as solutions of -50

-51potential problems. The near field terms are then used to produce a multipole expansion for the far field. Later, Kleinman (1967) showed that in some cases the method fails to produce terms beyond the first order (m = 1) term and amended the procedure. In principle, the revised method can be used to determine the partial fields of any order. However, if analytical expressions are desired, there is a difficulty associated with the resistive (or perfectly conducting) boundary condition. For every m, Em is the sum of a known vector and a gradient of a potential function (Kleinman, 1977). The boundary condition yields equations for the tangential derivatives of that potential. In some cases, such as m = 0 or bodies of special geometry, these equations can be integrated analytically to yield a standard Dirichlet boundary condition. (Numerical integration involving an open line integral is always possible.) In some instances at least, it is possible to derive analytic expressions for the first order far field, in terms of zeroth order potentials, without finding the corresponding near field first. 4.2 The Scattered Far Field When determining the scattered far field at low frequencies the Franz representation displays a clear-cut advantage over that of Stratton. Franz's formulation is the "natural" representation for the scattered far field, since that field is generated by the induced currents only. The induced charges have no direct effect on the radiation field. One near field term is needed to determine the corresponding far field term if the Franz representation is used,

-52but two near field terms are required when the Stratton-Chu representation is employed (Kleinman, 1967b). Using the Franz formulation with electric currents only and the far field approximation, the scattered field is given by ikr A A A ~S eikr ^ ^ ik^ r-S E i z r.r.J(n.H) e- dSI 4r rr I(nH) er dS' eikr ^ 00 A A -ikZe r Fr^ (-1)P (ik (n^H)(r')P dS' kD<< p=o Substituting expansion (4.1) for H and collecting equal powers of ik, the field is expressed as: ikr A A s: -ikZe rAr^ ( H)dS' + ik (n.H)dS' - 4rr ko - ik f(nAH)r-r' dS' + (ik)2 (n^H) dS' - (ik)2 (n.H )(r-r') dS' + ()2 (n^Ho)(r-.')2 dS' + O(k4) Further simplification is obtained by the vector identity (Jones, 1964, p. 531). n^V dS' = (nv'^V)F' dS' (4.3)

-53Hence by (4.2) fn.H dS' = (n-v' )' dS'= AH dS' = -Y (n-E )i' dS' n.H dS' = -Y (n-E )r' dS' — 2 1 Finally ~s - E^ - 4mr -r1-r ^r^PO + kZr^m + ikr^) (4.4) 4irr rrr~p0 ZrA Olr4rA(r.J where from Kleinman (1973) Po = f(n-Eo)r' dS' (4.5) m = -fr'c (n^H )dS' (4.6) and F(r) = (n.E )r' dS' + Z (nH )r.r' dS' Z (HdS' (4.7) _ z f(a~H )r-.')2 dS' (4-7) 2 ^

-54An alternative and more symmetric form is Es = k2e ) + Zrr(p + ik) ) + O(k4) (4.8) 4tr re r^,,, (po + i kpl MO + ik (4.8) which is obtained by further applications of (4.3) to V = H (r-r') 1 and to V = H (r-r')2. However, (4.4) proves to be more convenient for actual calculations. The quadrupole moments in (4.8) are given by p = { (nEl)r' dS' -. (n-Ed (r.r')F' dS' 1 and. (4.9) 1 = J r^'(n^H)dS - 4 rl ^ (n^Ho)(rr')dS' Once the induced currents and charges are determined, the far field can be calculated by (4.4) or (4.8). It is now necessary to deal with the particular geometry and constitution of the scattering body. The following sections will focus on perfectly conducting and resistive plates. The perfectly conducting plate is the more basic structure and is discussed first. 4.3 Perfectly Conducting Platelets The first two terms of the expansion (4.1) of the incident fields (1.1) and (1.2) are: -inc ^ -inc: Ein = a, HC = Yb,(4.10) Ein = a(k r), Hinc = Yb(k-r). (4.11) 1 1

-55The currents and charges on a flat infinitesimally thin plate are given as differences across the plate A -1 J z= (H (4.12) and _ = z.E+. (4.13) Dipole Moments. The static fields Es and Rs are 0 0 irrotational (see 4.2) and can be obtained from potentials. In the case of, Ess = -vE where d is an exterior potential. The 0 0o corresponding incident potential is -a*r, and from the boundary -s A -inc condition z^E = -z^E on the plate it follows that:t = a*r + c on S, (4.14) where c is a constant chosen to satisfy the zero induced charge condition J dS =. (4.15) S Because of the form of (4.14) it is convenient to write 3 ^ A^ = Z(a'xi)o (4.16) i=1

-56where the xi, i=1,2,3, are Cartesian coordinates with x = z. The 3 potential (i can be represented as i i o ~-~-~- Go dS' S and by applying the boundary condition (4.14) we obtain r 3 +I+ -x. - C. = 0 G dS' (4.17) S i + which is an integral equation from which to determine a- - On substituting (4.16) into (4.5) we have p = cP-a (4.18) 0 (Keller et al., 1972), where P is the electric polarizability tensor with elements i + P j = -J -& x. dS, (4.19) 1Pij - z 3' S which are functions only of the geometry of S. It is evident that P.. = 0 if i,j = 3, and from (4.15), (4.17) and (4.19) 0 0 P = JJGO az at dS dS' (4.20) J31 The tensor is threfore symmetric, i.e., P P and has at m ost three independent elements.

-57-s The procedure for the magnetic field is similar: H= -Yvo where *o is again an exterior potential in terms of which m= Yf _i z dS (4.21) S (Kleinman, 1973). The boundary condition on the plate is zH = -inc -zH,n implying 0 ^ o = b-z on S (4.22) Therefore i is normalized as )- = b-z 3. (4.23) The representation JI 3 3 i aG 93r = j m z l 31 dS' S then leads to the integral equation - a 2 f GI 31 dS' (4.24) 3z2 S - from which i31+ can be found. In terms of this m Y(b-z)z I3jo dS (4.25)

-58implying (Keller et al., 1972) a magnetic polarizability tensor M with just the single non-zero element M 31- + dS (4.26) ^3 3'o'0 S The solution of equations (4.17) and (4.24) and the computation of the dipole moments are discussed in Chapter V. 4.4 Resistive Platelets The resistive plate is described mathematically by (1.3) and (1.4). At k = 0 however, the condition V = o must be explicitly applied. The plate is therefore invisible as regards the static magnetic field, i.e., Ho = 0 and (1.4) then implies that for all R, Eo is the same as if R = O. Thus, for a resistive plate, mo = 0 and the electric dipole moment po is the same as for perfectly conducting plates. To determine how the low frequency scattered field is dependent on R, it is necessary to examine the first order terms in (4.4). Since H = 0 the third term of (4.7) is zero. To evaluate the second term, z^H has to be determined first. From Maxwell's equations it is found that

-59z A A A ^4Z A a CzH) = z^VH + Yz(z^E) (4.27) and by application of (1.4) on S [ E a Ey Hz - U ax ay ik L ax y ik az Hence, 8H + Z 1 A __- H = -b-z on S az = R oz= R In addition H1+ = 0. The Cartesian components of H are themselves exterior potentials and therefore a 3H + A'-r H = - z G dS' - 1 b-z G dS' on S. (4.28) iz a z' o R J o S S Also, ^ r + aGo z H = z^H |I a dS' S - az/J Go dS' S o and hence, on the plate,

-60a ^ - Z ^^ ^ -,y -zH ) = - - JZHI -GO dS' 1z., H 1,,H G dS' (4.29) Using (4.27), (4.28) and the fact that z^Eo = -zxa on S, the left-hand side of (4.29) can be evaluated, and the solution of (4.29) can then be written as. A A A 1 ^ A z^H = -Yz(za)31, + - b-z u (4.30) where u is a tangential vector satisfying the integral equation J> GOd'= 2 f (r') GO dS' (4.31) C S with the line integral taken around the perimeter C of the plate, and 43J+ satisfies 0- = az2 f GodS' (4.32) 31+ is proportional to the magnetostatic potential for a perfectly conducting plate. The numerical solution of (4.31) requires only a simple modification of the program used to solve (4.32).

-61Knowing 3o and u, the third term of (4.7) is A A A + Z J zH z r-r dS' = -zA(zAa)J 31 r'r dS' S S o + Z u(r')rr' dS'. (4.33) S The next step is to compute the first term on the right-hand side of (4.7). Since H = 0, E is the gradient of an exterior 0 1 potential (Kleinman, 1967), i.e., E =-vs with + r -;' dS = 0 jtz S Hence 2 z |dS' d= - (x + c.) dS' S i=1 S and from (4.14), (4.16) and the reciprocity theorem for exterior potentials, 2 z-E 1+' dS' = - xi i - z.v'l dS'. (4.34) S i=1 S

-62The zeroth order potential (^ is associated with an incident field having a = x.. If the corresponding first order scattered magnetic field is i, Maxwell's equations imply 1$ 1ze v J = ZZ-VAH, (4.35) and when this is substituted into (4.34), the vector relation in Jones (1964, p. 531) can be used to give 2 AiI+ zE r' dS' = Z x. zf v' Hzi dS' S i=1 S From the boundary condition (1.4) A -inc A A ZV: = zEin - RzA(zAH +) ^V1 1^ 1 A.A. A. 31- (A.)] zAla(kr) - a u(b-z) and hence 2 f z.l _ +' dS' = - -z (b'-) - z). S i=1 S ~z^Hl dS'. (4.36) 1 -

-63Finally, from (4.30) -zo+ 3 1" z^H1 = -Yz^ )z^x 3 +i R_ z.(k^xi)u and when the vectors are combined to eliminate the summation, ZElr' dS' = {[(kr')z^(z^a) - (bZ)ZA(Zu) (a u)(z k)] + -R z^(za)o+ 31+ - R (z^k) [ (kr')a u - (b)u. (4.37) The complete first order contribution is therefore F(r) = f [(k - r) r' z.(zAa) - (b*z)z^(z^5) S ^ A A + R A + + + (a u)(z^k) ] z^ (z^a) o -o 3+ ZR[(r-r')(bz)u - (zAk) (kr')a.u - (b-z)u u dS' (4.38) It is evident that F has components only parallel to the plate. As expected, F is a function of R, and the presence of terms inversely proportional to R makes explicit the discontinuous behavior as R + 0. For any given R f 0, there is a frequency (which

-64depends on the plate dimensions) below which F no longer represents a valid correction to the zeroth order contribution. A case of some interest is that in which the incident magnetic vector is parallel to the plate, so that b*z = O. Then, z^k = ^ ^A -(a.z)b and F(r) = ZA(z^a) [(k- r)r' 31 31 dS' S +.b(a-^z)f k~ [ - RP31I]a-u dS' (4.39) ~SR Z O If, in addition, az = 0 implying normal incidence on the plate, the second integral in (4.39) vanishes. The first two order terms in the low frequency expansion then go over smoothly into the result for a perfectly conducting plate as R - 0. On the other hand, if k*z = 0 corresponding to grazing (edge-on) incidence, a is in the z direction (and perpendicular to u) and F(r) = 0. This is consistent with the fact that any electric current sheet is invisible to a plane wave at grazing incidence for perpendicular polarization.

CHAPTER V. NUMERICAL SOLUTION OF THE STATIC INTEGRAL EQUATIONS 5.1 The Electrostatic Potential The electrostatic potential is derived from the integral equation (4.17). Although the equation is not explicit because of the unknown constant c., the form of (4.17) suggests the following splitting: -o + Oci~+ (5.1) 0o i(1) i(o where oi() and i(2) satisfy the integral equations: 0 0 () + = f ~ I G dS' (5.2) and i(2)+ -1 GO dS' (5.3) respectively. With this definition and the total induced charge condition (4.15), ci can be determined as r (~1j) + ~ qi(2) + C. = - z dS, dS' (5.4) Integral equations (5.2) and (5.3) differ only in the excitation term and therefore can be treated alike. They can be summarized by the single equation -65

-66fJt a i (i) Go dS' = -. (5.5) S (-2 The weakly singular kernel Go does not impose continuity restrictions on the solutions, and the simple pulse expansions M N ^^^J+ rFXx - X p yzl = Eapq rect rt Yy (5.6) p=i q=i are adequate,where + stands for ( () or () Equation (5.5) ox 0y o is converted into a set of algebraic equations M XmN 1 < m < M ~ ~ apqG n p q apq G~m~n~plq y (5.7) p=i q,np,q 1 N P=1 q=1 1 1 < n < N when sampled by the impulse functions Wmn = (.x-x )6(y-y ). G is similar to G of Section 3.3, except that it mn q m,n,p,q' m,n,p,q involves the static Green's function Go instead of G. Explicit expressions for Go are given in Appendix C. m,n,p,q The amount of numerical labor can be considerably reduced if the plate is symmetric with respect to one or two axes. If r is the reflection of r with respect to the x. axis, then

-67ai(.) ii+ () + and (2) (2) 0 - ()) ) as is evident from (5.2) and (5.3). Then, (5.4) implies that c. = 0 for plates with at least one axis of symmetry. Furthermore, one component of (5.7) can be written in the form M/2 N Z7 ^~~~ v ~1 < m < M/2 a (Go - Go ) -X; p= ql Pq m,n,-p,q m,n,-p,q 1 < n < N when that axis of symmetry is the x-axis and the cells are numbered accordingly. Similar equations hold for plates that are symmetric with respect to the y axis or with respect to both axes. Convergence of the solution of (5.7) is enhanced if the cell size is allowed to decrease toward the edge to better accommodate -1/2 a variation of p (Meixner, 1972), where p is the distance from the edge. However, the possibility of treating both electrostatic and magnetostatic potentials by the same program makes that variation of cell size undesirable. Once Eqs. (5.7) are solved, it is possible to compute the elements of the electric polarizability tensor using (4.19)

-68together with (5.1) and (5.4). Numerical results are presented in Section 5.3. 5.2 The Magnetostatic Potential The integral equation (4.24) from which the magnetostatic potential is determined is much more singular than (4.17). Consequently, the more elaborate scheme of Section 2.2.2 must be used. The first step is to replace the second normal derivative in (4.24) by the surface Laplacian (a2/3x2)+(a2/ay2), using the definition of the Green's function: v2Go = -6. Then, the tangential derivatives are represented by (2.22), and 31J itself is represented by a pulse expansion N x-x -1 y- y 1 = bpq rect rect q (5.8) o - pqx A Ay p=1 q=l The resulting algebraic equations are M N 1 = 4 Z ^ bpq G +G 2Gm, q 4pq AX 2 mm+/2 p,q m-l/2,n,p,q m-n,p, + G +G -2G \ >y O m,n+l/2,p,q om,n-1/2,p,q,n,p,q 1 <m<M -; - e ( 5 )(5.9) 1 <n< N

-69Symmetry considerations similar to those presented in the previous section apply to (5.9). Mzz is computed using (4.26). 5.3 Selected Numerical Results Several factors determine the efficiency of the numerical solution of (5.7) and (5.9). The translational symmetry of G mnpq allows a very considerable saving of matrix filling time if equal size cells are used. Because the matrix elements of (5.9) are combinations of terms similar to the matrix elements of (5.7), it is possible to incorporate the solution of both equations into a single program, leading to a further reduction in filling time. Several versions of the computer program were written. Programs EPLT and APLT solve both (5.7) and (5.9) and employ uniform cell size. EPLT is for plates that are symmetric with respect to both the x and the y axes and APLT is for asymmetric plates. Program DPLT solves (5.7) for symmetric plates. DPLT allows several optional variations of the cell size. It seems that the quadratic taper: 3n2 N(N + l)(2N +1) (5.10) is the best compromise of decreasing the cell size near the edge without increasing the middle cells too much. N is the total number of subdivisions per half plate along a given axis and n = 1,2,...,N is a running index of the subdivisions, counting from the edge towards the middle.

-70Results of a convergence test for Px (=Py) for a square plate Xx yy are displayed in Fig. 5.1. Curve No. 1 presents results computed by EPLT, and curve No. 2 the corresponding results of DPLT. The convergence rate is compared with that reported by Okon and Harrington (1980, p. 76). They used triangular cells with area coordinates but also tried uniform cells (curve No. 3) and tapered cell sizes (curve No. 4). Since their program was intended to be as general as possible, they have avoided taking advantage of the symmetry. To make the comparison in terms of sampling points on the plate, their reported number of subdivisions was divided by 2. For the same N their matrices are four times as large. Other adaptations of their results were multiplication by a factor 2 and an interchange of the roles of the electric and magnetic dipole moments, since they were concerned with apertures. The vertical scale of Fig. 5.1 was greatly magnified to emphasize the differences in the convergence rates of the various methods. As expected, the methods using tapered cell size converge much faster than those employing uniform cells. The quadratic taper (5.10) is slightly more efficient than the linear taper used by Okon and Harrington (ibid). Their results computed with uniform cells have not converged yet, but seem to follow the same curve as the results computed by EPLT. There was a question as to how well our programs would perform when applied to geometries that do not conform to rectangular cells. The circular disc for which an analytical solution is available was

-711.05 1.00 PXX A3/2 0.95 0.90 / I I 1 0*9 5 N 10 Fig. 5.1 Convergence of P. XX

-72used to test the performance. Only cells whose center lay inside the circle were retained and the effective shape was polygonal with serrated edges. Since the potentials grow indefinitely towards the edge, this procedure could conceivably produce large errors. Normalized values of Pxx (=Pyy) obtained by EPLT are listed in Table 5.1 along with the corresponding errors. Results of Okon and Harrington (ibid) are listed as well for comparison. Although their grid matches the circular shape much more closely, their results are better only when tapered cells are used. Comparison of the results of EPLT for a rectangle and a disc shows that the accuracy was reduced from a small fraction of a percent to 3 to 4 percent when the same number of sampling points were used. Figure 5.2 shows the convergence of Mzz for a square plate. The curve marked by 6 = 0.5Ax was produced by EPLT. The one marked 6 = 0.35AX depicts the results of a similar program where the nodes were allowed to move away from the edge of the cell, as shown in the insert. A third line, marked a = 0, shows the results of evaluating the second normal derivative in (4.24) analytically under the integral sign, expanding 31+ by pulses and point matching. The rapidly decreasing rate of convergence is the result of the loss of continuity. The utmost importance of CO[S] continuity is highlighted in a different way in Fig. 5.3. Mzz is plotted as a function of 6 for a fixed N (=2). Compared to the value displayed in Fig. 5.2, the error in the computed value of Mzz has a sharp minimum at 6 = 0.5Ax, where the numerical representation is continuous. Most

-73Table 5.1: Normalized Dipole Moments for the Disc M zz/A2 Error P /A/2 Error A Method N xx (%) Error (%) 7 0.4604 3.9 0.9082 5.2 1.3 9 0.4614 3.7 0.9184 4.1 0.6 W 0 0.4608 3.8 0.9249 3.4 0.6 12 0.4620 3.5 0.9297 3.0 1.0 o 1 E 4 0.4506 5.9 0.8904 7.02 0.64 r=1' 5 0.4572 4.5 0.9038 5.64 0.41 s- Q I ^ 4 0.4754 0.7 0.9644 0.69 1.35 o 1 5 0.4766 0.46 0.9570 0.09 0.83 o > Exact 0.4788993 0.9577986

-740.60 \ \ ~~- ax 6^\5=0 AX 0.55 \ 6 = 0.35 Ax Mzz A3/2 0.50 ^^ 6z = 0.50 Ax 0.45 Okon and Harrington 5 N 10 Fig. 5.2 Convergence of Mzz.

-750.70 0.60 M zz 3/2 0.500.5 Ax Ax Fig. 5.3 M as a Function of 6 for Fixed N. ZZ

-76noticeable is the jump occurring between 6 = 0.4995Ax and 6 = 0.50Ax. The same behavior was observed for larger values of N, where the scale was compressed due to the convergence. The convergence rate with 6 = 0.5Ax is comparable to that of Okon and Harrington (ibid) with variable cell size. The opposite direction of convergence is attributed to their use of different basis functions which underestimate the potentials over most of each cell. The accuracy of the results for other geometries was again tested using a disc, and the results are summarized in Table 5.1. The accuracy is comparable to that of Okon and Harrington (ibid) with uniform cells. The accuracy, however, depends on the particular choice of the grid and care must be exercised in selecting this grid. Meshes that consistently overestimate or underestimate the area of the plate should be avoided. Numerical results for the dipole moments of several plates are presented in Table 5.2. The shapes of these plates are displayed in Fig. 5.4. For each shape the normalized dipole moments are given for several values of the length to width ratio (L/W). For a given value of L/W, the dipole moments of plates that are convex in the plane show very little sensitivity to the shape. Shapes that are concave in the plane deviate from this rule, but the largest deviation is associated with doubly connected shapes. These observations are supported by the data of De Smedt and Van Bladel (1980). For such shapes, normalization to (area)3/2 is no longer appropriate and, in fact, the polarizability tensor elements diverge as the ratio

-77L L rectangle rhombus L L w?_ __LW W/2 L/2 L/2 cross frame L/2 W L "L" Shape Fig. 5.4 The Various Shapes Used in Table 5.2.

-78Table 5.2: Dipole Moments of Various Plates Shape L/W N M /A3/2 p /A3/2 P /A3/2 rectangle 1 10 0.4479 0.9889 0.9889 12 0.4474 0.9972 0.9972 2 7 0.4165 0.6055 1.655 9 0.4148 0.6151 1.681 10 0.4144 0.6185 1.691 5 7 0.3129 0.3513 3.782 9 0.3117 0.3566 3.839 10 0.3113 0.3584 3.859 rhombus 2 10 0.4050 0.5788 1.862 12 0.4052 0.5841 1.884 5 10 0.3028 0.3289 5.011 12 0.3027 0.3317 5.081 cross 1 8 0.4290 0.9816 0.9816 10 0.4268 0.9966 0.9966 2 8 0.3981 0.5954 1.762 10 0.3961 0.6045 1.787 5 8 0.3007 0.3393 4.343 10 0.2994 0.3442 4.399 frame 1 8 0.2490 1.496 1.496 10 0.2445 1.516 1.516 2 10 0.2423 0.9477 2.591 5 10 0.2201 0.5490 5.914 ~T"L" 1 10 0.4027 1.095 1.095 0.3040 12 0.3981 1.116 1.116 0.3070

-79of the area of the plate to the area of a singly connected plate of the same exterior contour becomes smaller. In the case of a solid body a variety of isoperimetric and other bounds on the tensor elements have been developed (Kleinman and Senior, 1972), and may serve to adequately approximate the elements, For a plate the more useful estimates are empirical in nature (De Smedt and Van Bladel, 1980), but bounds can be rigorously established from the fact that for a solid body 3 3 3 3 P = 7 Pij(a-xi)(axj) and M = Mj(b x)(b-x) i=i j=1 i=l j=1 are increasing set functions (Schiffer and Szego, 1949). By regarding the plate as the limit of a right prism, we have P. < P < Pit 11- ii -1 - 1 (i = 1,2) M' < M < M" 33 33 - 33 where the single and double primed quantities refer to the largest inscribed and smallest circumscribed circular disks respectively. For a disk of radius a, P = P = 16a3/3 and M = 8a3/3. These 11 22 33 bounds are not applicable for doubly connected shapes. It is evident from Table 5.2 and other data published in the literature (De Smedt and Van Bladel, 1980; Okon and Harrington, 1980), that P,, and M,, decrease with increasing L/W, whereas Pyy

-80increases (y is the long axis of the plate). The bounds in terms of inscribed and circumscribed circular disks become useless for large values of L/W. Other bounds, in terms of circular disks of the same area, can be conjectured, but no mathematical proof is yet available for them, and there is always a risk that some shape exists that violates them. As an additional test, a comparison was made of the dynamic and the static programs. It is possible to compute (less efficiently) the dipole moments by the program for dynamic scattering. Making the plate very small compared to the wavelength and using E-polarization at grazing incidence so that b*z = a-x = 1, the dipole moment m is given by: i = r' A^J dS' = (x'J - Y'Jx)dS' S S and by using (3.5) and (3.6), M N-i M-1 N 0m1 2 \y E m Ym,n+1/2 n Xm+1/2, m=i n=i m=n n i (5.10) Similarly, M N Pxx = xI dS' -= xLy n M-i N i m n n+l/2 (5.11)

-81Computed results for square plates are shown in Fig. 5.5. The values of PXx (=Pyy) for both perfectly conducting and resistive plates approach the static limit as W decreases. The magnetic dipole moment of the metallic plate is close to the static limit and the magnetic moment of the resistive plate disappears. These results were predicted in Chapter IV, and are used as an additional verification of the computer program.

-82~ Metal 1 ic *.-. * Resistive Static Limit 1.0 Pxx A3/2 Static 0.5 / Limit M xx 3/2 - Static / imit p,,, I I I I I I I 0.05 W/x 0.10 Fig. 5.5 Convergence of the Dipole Moments to the Static Limit.

CHAPTER VI. CONCLUSIONS AND RECOMMENDATIONS The problem of electromagnetic scattering by resistive and perfectly conducting rectangular plates was solved numerically. The numerical techniques developed for this problem can be easily applied to other problems involving E-field integral equations, e.g., wire antennas, strips and non-planar open shells. Impedance plates can be treated as well by superimposing resistive and conductive plates. The numerical procedure was tested by comparison with measured data and with analytical and computed results for very small plates. The accuracy of computed radar cross sections for plates of about a square wavelength in area is 1 to 2 dB. The qualitative behavior of the induced currents, especially the edge currents, checks well with currents measured on slightly smaller plates. Numerical representations were developed for a variety of operators. Comparison of Eqs. (3.11) and (3.12) with the Cartesian components of Eq. (3.1) reveals, a posteriori, the following numerical representations for the differential operators 32/ax3y and 32/ayax: 32f x +' Y - f X1 +' Y xa y (Xm'Ym) axay [ (x -, 2 m 2 - f(xmYn ) + f(XmYn - 4 )] (6.1) -83

-84and ay_ 2xy -fx' Y + ayaI( x xy [f(XYmYn x2 ) m,' n 2/ m, n - f(xmYn) + f(xm -^ yn)] (6.2) that preserve the continuity of f given by (2.22) and (2.23), and with these representations it is possible to deal with E-field integral equations derived from either the Franz or Stratton-Chu formulations. At present, the computer program has two major limitations: it handles only rectangular plates and the maximum area of these plates is restricted to about a square wavelength. Both shortcomings can be removed or alleviated. There are several possible courses of action which enable other shapes to be handled. The most straightforward one is to try to approximate that shape using rectangular cells. This method was used in the static program and produced results for the dipole moments to an accuracy of a few percent. The radar cross section is a weighted integral of the induced currents and is thus expected to be somewhat stationary. Another approach is to map the rectangular mesh onto a grid with general quadrilateral cells, to fit more accurately the given shape. Experience with the static programs may discourage one from taking this course. The cost of running the program that employs

-85unequal cells, was higher than that of the other program, using equal cells, with matrices about twice as large, (.144 x 144) compared with (64 x 64). The reason is that the latter program takes advantage of the translational symmetry of the matrix elements to reduce the number of different elements to M x N, where M and N are the number of subdivisions along the x and y axes, respectively. This cannot be done with unequal cells and the number of matrix elements to be computed rises to (M x N)2. The fill time of the corresponding matrix exceeds the inversion time of the much larger matrix in the other program. A possible compromise is to use special shapes for edge cells alone. Since the number of edge cells is but a fraction of the total number, the matrix fill time will increase less sharply. Another feature that can be incorporated into a special treatment of the edge cells is to use singular elements (Zienkiewicz, 1977, p. 670; Jones, 1979, p. 506). This may somewhat increase the rate of convergence, and its effect will be more pronounced in the static programs, where the number of required sampling points does not depend on the size of the plate. Non-planar open shells can be treated by mapping the rectangular grid onto a curved mesh. A simple example is Glisson and Wilton's (1980) treatment of a bent plate. In general, such a mapping will result in unequal cells. Unlike the shape limitation that can be relatively easily removed by a modification of the existing program, an increase

-86in the maximum plate size that can be handled requires improvement in the matrix inversion techniques. The interaction matrix for the coupled integral equations is a dense matrix. When equal cells are used, the matrix is symmetric and there is an additional hidden symmetry because the kernel of the equation is of a convolution type. The interaction matrices of one-dimensional structures with convolution kernels are of Toeplitz form. However, there is no way in which the corresponding interaction matrix of a two-dimensional structure can be put into that form. Some progress in inversion techniques for matrices that are close to Toeplitz matrices has been reported by Friedlander et al. (1979), but the technique is not applicable to plates with variable resistivity because the elements on the main diagonal of the corresponding interaction matrix are all different. Matrices of this type may have to wait for hardware improvements that will permit the handling of much larger matrices.

APPENDICES -87

APPENDIX A. DERIVATION OF THE REPRESENTATIONS FROM THE DEFINITIONS OF THE DYADIC GREEN'S FUNCTION The electric and magnetic dyadic Green's functions are defined (Tai, 1973) by: v^Ge =Gm (A.1) and ^ = k2e + I6(r - r'). (A.2) Gm e They also satisfy v-G = - r')] (A.3) e k2 and V-G = 0. (A.4) m Equations (A.1) and (A.2) can be replaced by the uncoupled, second order, vector Helmholtz equations v.^^G - k2G = Is(r- r') (A.5) e e and v^v^m - k2Gm = v^[L6(r - r')]. (A.6) -88

-89Levine and Schwinger (.1950) concluded from (A.3), (A.4), (A.5), (A.6) and the Helmholtz equation for the scalar Green's function v2G + k2G = -6(r - r') that G = [ + VV G (A.7) e k2 and G = V,(IG). (A.8) Equations (A.5) and (A.6) are functional equations in a distribution space (and so are (A.1)-(A.4)). Before interpreting them in that way, it is necessary to specify the distribution space. Gelfand and Shilov (1964) noted, in the preface to the second volume of their book, that different distribution spaces should be used for solving different types of problems. It is equally true when solving scattering problems with different geometries. For scattering by finite smooth bodies, the most familiar distribution space can be used. This is the space of linear functionals over the space of infinitely differentiable vector functions J(r'), with finite support on the scatterer. A different space is required for finite scatterers with edges. The test functions are differentiable and bounded only in the interior of the support, and some of their components may diverge towards the edges as p (0 < t < 1) (Meixner, 1972). A countable set of norms for this space is

-90IJIIp = sup iPJr)l = 1,2... rl' S With these norms the proof that the test function space is a fundamental one and the construction of the dual space, are exactly analogous to those of Gelfand and Shilov in Vol. 2, Ch. II and Vol. 3, Ch. I. Other distribution spaces should be used for infinite and semi-infinite bodies with or without edges, but they will not be discussed here. Equation (A.5) should be understood as a functional over the test functions J(r'): <vV^Ge,J> - k2 <Ge,J> = <I5,J> Since Ge is a convolution type distribution, the differential operator can be brought outside to give = - 2iY (vAv^ - k2) <Ge J> - (vvv - k2)E The last equality is obtained from the vector Helmholtz equation for E. Ge satisifes the radiation condition and the same boundary conditions as E. The uniqueness of the solution of the Helmholtz equation therefore implies

-91E = ikZ <G,J> = ikZ ++ 1VV <G,J> k2 The distribution G corresponds to the ordinary function G(r,r'); hence <G,J> = G(r,r')J(r') dS' S As far as the scattered field alone is concerned, the body is replaced by the induced currents, and the free space Green's function should be used. Hence = ik[ I 1 ]V J G dS' (A.9) k2 and a similar interpretation of (A.6) leads to: s = V. GJ dS. (A.10) Equations (A.9) and (A.10) are identical with the Franz representation for the scattered field when only electric currents are present. Magnetic sources can be treated analogously and the complete representation is then obtained by superposition.

APPENDIX B. PROOF OF co[Q] CONTINUITY OF THE SUBSECTIONAL BASIS Using the notation of Fig. B.1, the coefficients a. expansion (2.21) for cell I are given by: c( = f(O) I a f(3)- f(1) f(2)- f(4) 2 ax'A3 Ay f(6) + f(8) - f(5) - f(7) 4 AxAy = 2 f(3) + f(1) - 2f(0) a = 2 f(2) + f(4) - 2f(0) 5 ^XAX2 6 Ay2 2:Xa = [f(6) - f(8) + f(5) - f(7) + 2f(4) - 2f(2)] 7 Ax2Ay a = 2 [f(6) - f(8) - f(5) + f(7) +2f(1) - 2f(3)] 8 AxAy2 a = X22 f(5) + f(6) + f(7) + f(8) - 2[f(1) + f(2) 9 Ax2Ay2 L + f(3) + f(4)]+ 4f(O)). (B.1) -92

-93-'13 11 5 2 6 14~0. 1, 3 8 4 7 Fig. B.1 Notation for Two Adjacent Cells.

-94A similar expression with coefficients a. and variables * * 5 = 5 = x - x, nn = -Yn+l holds for cell II. Taking the difference A= f,(,n) - fmn+ (,n ) along the common line of the two cells yields A =: (ta -c*) A + (a -2) + (a - a*) + (a + a Y i I 6 6 4 2 2 4 4 2 2 +(8 a8*) 45 + (a c*)+( + - +)A + (' a*) 2 + 8 8 4 5 5 " 7 7 2 779 9 4 =A + BE + C2. (B.2) Using (B.1) and the corresponding expressions for cell No. II in (B.2), the coefficients A,B and C are obtained: A = f(O) - f() + f(2) - f(4) + f(ll - f(2) + f(2) + f(4) - 2f(0) 2 2 2 f(ll) + f(2) - 2f(9) = 0 2 Ax B = f(3) - f(l) -f(12) - f(10)j+ f(6) + f(8) - f(5)- f(7) 2 f(4) + f(5) - f(6) - f(13) f(6) - f(8) - f(5) + f(7) + 2f(l) -2f(3) 2 2 f(14) - f(5) - f(13) + f(6) + 2f(10) - 2f(12) 2

-95Ax2-C = 2[f(3) + f(l) - 2f(0)] -2[f(12) + f(10) - 2f(9)] + If(6) - f(8) + f(5) - f(7) + 2f(4) - 2f(2)] + [f(14) - f(5) + f(13) - f(6) + 2f(2) - 2f(11)] + [f(5) + f(6) + f(7) + f(8)] - 2[f(1) + f(2) + f(3) + f(4)] + 4f(0) - [f(13) + f(14) + f(6) + f(5)] + 2[f(10) + f(ll) + f(12) + f(2)]- 4f(9) = 0 N Therefore f = fmn is continuous across the cell boundary. N=1

APPENDIX C. EVALUATION OF THE MATRIX ELEMENTS No closed form is available for the integrals x +Ax/2 yx+Ay/2 Gf F G(x,yvx y')dx' dy' xK-ax/2 yX-Ay/2 of Section 3.3. When the observation point is not in the integration cell, the integrand is smooth and regular. Numerical integration by 2x2 Gaussian quadrature is optimal (Zienkiewicz, 1977; p. 280). The self-cell term is approximated, when the cell is small enough, as fol 1 ows: r ikj,-~'j lr 1 + ikrm -r'-{ G1 d Gm =1 emn ~ X~ dS' -mn ~ ~~ dS' mmnmn J 41J mS mn 41 - mn I 4 S mn- mmnmn G + ikAxAy. (C.1) mnmn G is the self-cell term of the mtarix of the electrostatic m,n,m,n equation (5.7). This integral can be evaluated, using cylindrical coordinates and the notation of Fig. Cl(a): -96

-97e i ty 3 Ay/2 _____ Ay (xSy) (x Sy) yl 2 (Xm'Yn ) | | Ln ) Ay/2 AX AX (a) Self-Cell (b) Self-Cell; node on the edge (XmYn) ] (XpYg) -2 1 (c) Observation point on the same row (Xm'Yn) (d) General case Fig. C.1 Geometries for the Calculation of the Matrix Elements.

-98qo ax/2secp iT/2 Ay/2cscp Go =n, 4 d dp + dp + df dp mn,m,n O O 0 2 [ x sinh-1 y + Ay csch- (C.2) The other cases of Fig. C.1 are evalauted in a similar fashion: G = 2Ax ~ sinh 2 + Ay csch (C.3) Gm,n,m+l/2,n as illustrated in Fig. C.l(b), and G = a * sinh Y - 3- sinh ^ + y[csch csch m,n,K,n C.4) where (X xm- x+ AX I5 "Ax + = Im K- XI +A x'm XK 2 and =4yA Y: 2 as illustrated in Fig. C.l(c) and

-99G = 2 (c[sin' I - sinhKl - ] + B[sinhl ^ - sinh 1 ] + y[csch - - csch B] +, [csch - -c csch ]} (C.5) where a and B are as before and Y = Y - YX + 2 and: IYn - yl-2 ^n XI-'

BIBLIOGRAPHY -100

BIBLIOGRAPHY Bethe, H. A. (1944), "Theory of diffraction by small holes", Phys. Rev. 66, 163-182. Bouwkamp, C. J. (1947), Review of: Copson, E. T., "An Integral equation method of solving plane diffraction problems" in Math. Rev. 8, 179-180. Bouwkamp, C. J. (1954), "Diffraction theory", Rep. Prog. Phys. 17, 35-100. Butler, C. M. (1974), "Formulation of integral equations for an electrically small aperture in a conducting screen", Proc. 1974 Int. IEEE AP-S Symposium, 101-104. Crispin, J. W. and Siegel, K. M. (1968), Methods of Radar Cross Section Analysis, Academic Press. De Meulenaere, F., and J. Van Bladel (1977),"Polarizability of some small apertures," IEEE Trans. Antennas Propagat. AP-25 (2), 198-204. De Smedt, R., and J. Van Bladel (1980), "Magnetic polarizability of some small apertures," IEEE Trans. Antennas Propagat. AP-28 (5), 703-707. Friedlander,B., M. Morf, T. Kailath and L. Ljung, (1979), "New inversion formulas for matrices classified in terms of their distance from Toeplitz matrices," Lin. Algb. and Appl., 27, 31-60. -101

-102Gelfand, I. M. and Shilov, G. E. (1964), Generalized Functions, Academic Press. Glisson, A. W. and Wilton, D. R., (1980), "Simple and efficient numerical methods for problems of electromagnetic radiation and scattering from surfaces," IEEE Trans. Antennas and Propag. AP-28, 593-603. Harrington, R. F. (1968), "Field computation by the Moment Method," MacMillan Co., New York. Harrington, R. F. and J. R. Mautz, (1975), "An impedance sheet approximation for thin dielectric shells," IEEE Trans. Antennas Propagat. AP-23, 531-534. Jones, D. S. (1964), "The Theory of Electromagnetism," MacMillan Co., New York; p. 531. Jones, D. S. (1974), "Numerical methods for antenna problems," Proc. IEE 121, 573-582. Jones, D. S. (1979), Methods in Electromagnetic Wave Propagation, Clarendon Press, Oxford, England. Jones, D. S. (1980), "The scattering of long electromagnetic waves," Quart. J. Mech. Appl. Math. 33 (Pt. 1), 105-122. Keller, J. B. (1962), "Geometric theory of diffraction," J. Opt. Soc. Am. 52, 116-130. Keller, J. B., R. E. Kleinman and T.B.A. Senior (1972), "Dipole moments in Rayleigh scattering," J. Inst. Math. Applics. 9, 14-22. Kleinman, R. E. (1967), "Low frequency solution of electromagnetic scattering problems," in "Electromagnetic Wave Theory" (ed. J. Brown), Pergamon, Press, New York.

-103Kleinman, R. E. (1967b), "Far field scattering at low frequencies," Appl. Sci. Res. 16, 1-8. Kleinman, R. E. (1973), "Dipole moments and near field potentials," Appl. Sci. Res. 27, 335-340. Kleinman, R. E., and T.B.A. Senior (1972), "Rayleigh scattering cross sections," Radio Sci., 7 (10), 937-942. Knott, E. F., T.B.A. Senior and V. V. Liepa (1971), "Plates and edges," IEEE Trans. Antennas and Propag. AP-19, 788-789. LaHaie, I. (1981), Function-Theoretic Techniques for the Electromagnetic Scattering by a Resistive Wedge, Ph.D. Dissertation, The University of Michigan. Levine, H. and Schwinger, J. (1950), "On the theory of electromagnetic wave diffraction by an aperture in an infinte plane conducting screen," Comm. Pure and Appl. Math. 3, 355-391. Lord Rayleigh (1897), "On the incidence of aerial and electric waves on small obstacles in the form of ellipsoids or elliptic cylinders, and on the passage of electric waves through circular aperture in a conducting screen," Phil Mag., 44, 28. Mittra, R, Y. Rahmat-Samii, D. V. Jamnejad and W. A. Davis, "A new look at the thin plate scattering problem," Radio Sci. 8, 869-875. Meixner, J. (1972), "The behavior of electromagnetic fields at edges," IEEE Trans. Antennas and Propag. AP-20, 442-446. Okon, E. E. and R. F. Harrington (1980), "The polarizability of apertures of arbitrary shape," Technical Report No. 12, Syracuse University Department of Electrical and Computer Engineering, Syracuse, NY 13210.

-104Polya, G. and Szego, G. (1951), "Isoperimetric Inequalities in Mathematical Physics," Princeton University Press. Rahmat-Samii Y. and R. Mittra (1974), "Integral equation solution and RCS computation of a thin rectangular plate," Trans. IEEE Antennas and Propag. AP-22, 608-610. Schiffer, M., and G. Szego (1949), "Virtual mass and polarization," Trans. Amer. Math. Soc. 67, 130-205. Senior, T.B.A. (1975), "Diffraction tensors for imperfectly conducting edges," Radio Sci. 10, 911-919. Senior, T.B.A. and Weil, H. (1977), "Electromagnetic scattering and absorption by thin walled dielectric cylinders with applications to ice crystals," Appl. Opt. 16, 2979-2985. Senior, T.B.A. (1978), "Some problems involving imperfect half planes," in "Electromagnetic Scattering" (Uslenghi, P.L.E., ed.), Academic Press. Senior, T.B.A. (1979a), "Scattering by resistive strips," Radio Sci. 14, 911-924. Senior, T.B.A. (1979b), "Backscattering from resistive strips," IEEE Trans. Antennas and Propag. AP-27, 808-813. Senior, T.B.A. (1981), "Encased resistive sheet," unpublished notes. Stratton, J. A. (1941), Electromagnetic Theory, McGraw-Hill. Stevenson, A. F. (1953), "Solution of electromagnetic scattering problems as power series in the ratio (dimension of scatterer)/ wavelength," J. Appl. Phys. 24, p. 1134.

-105Tai, C. T. (1973), "Eigen function expansion of dyadic Green's function in electromagnetic theory," Math. Note No. 28, Weapons Laboratory, Kirtland Air Force Base, Albuquerque, N.M. Umashankar, K. R. and Butler, C. M. (1974), "A numerical solution procedure for small aperture integral equations," Interaction Note 212, Weapons Laboratory, Kirtland Air Force Base, Albuquerque, N.M. Van Bladel, J. (1964), Electromagnetic Fields, McGraw-Hill. Van Bladel, J. (1970), "Small-hole coupling of resonant cavities and waveguides," Proc. IEE 117, 1098-1104. Van Bladel, J. (1977), "Low-frequency asymptotic techniques," Chapter 1 in "Modern Topics in Electromagnetics and Antennas," Peter Peregrinus Ltd., England. Well, H. and Chu, C. M. (1980), "Scattering and absorption by wavelength size discs, in light scattering by irregularly shaped particles (Schuerman, D. W., ed.), Plenum Press.

UNIVERSITY OF MICHIGAN 3 901 5 03465 94851111111 3 9015 03465 9485