AFCRL-70- 0507 1363-6-T Scattering by a Periodic Surface by Tommy C. Tong and Thomas B. A. Senior The University of Michigan Radiation Laboratory 201 Catherine Street Ann Arbor, Michigan 48108 September 1970 Scientific Report No. 7 Contract F 19628-68-C-0071 Project 5635, Task 563502 Work Unit 56350201 Contract Monitor: John K. Schindler Microwave Physics Laboratory This document has been approved for public release and sale; its distribution is unlimited. Prepared For Air Force Cambridge Research Laboratories Laurence G. Hanscom Field Bedford, Massachusetts 01730

1363-6-T FOREWORD This report is part of a study whose ultimate objective is the development of adequate procedures for estimating the backscattering from a class of rough surfaces at oblique angles of incidence. In order to avoid a complete reliance on the physical optics approximation which is inaccurate at angles well away from specular and does not provide a reliable estimate of the de-polarizing effect of the surface, the approach that is adopted is to postulate a simple and deterministic base surface which is itself'rough'. By taking this base to be perfectly conducting, periodic and of infinite extent, the true surface field can be determined, as can the scattering at all angles of incidence. More general surface configurations such as might be appropriate to the sea are then treated by a modification or perturbation of the field on the periodic base. The first step in this process is to investigate the scattering of a plane wave by a two-dimensional, perfectly conducting periodic surface. This is the only topic considered in the present report. ii

1363-6-T ABSTRACT Numerical procedures are developed for the digital solution of the integral equations for the currents induced in a perfectly conducting, two dimensional periodic surface when a plane electromagnetic wave is incident. Data are obtained for both the surface and far fields for a variety of sinusoidal surfaces for plane waves of either polarization at oblique as well as normal incidence, and the results are compared with the predictions of physical optics. iii

1363-6-T TABLE OF CONTENTS Forewordii Abstract iii I Introduction 1 II Formulation of the Integral Equation for Numerical Solution 5 2.1 Formulation 5 2.2 Reduction of the Integral Equations for Periodic Surfaces 10 2.3 Numerical Procedures 18 III Application of the Numerical Technique and Presentation of Data 23 3.1 Discussion 23 3.2 The Scattered Field 24 3.3 The Current Induced on a Sinusoidal Surface 30 3.4 Energy Conservation 45 3.5 Convergence Check 52 IV Conclusions 57 References 59 Appendix A: Physical Opt i cs Approximation 62 Appendix B: Computer Programs 66 iv

1363-6-T I INTRODUCTION The problem of the scattering of an electromagnetic wave by a periodic surface is one which is easy to formulate but difficult to solve and it is only recently, with the advent of the high speed computer, that reliable results have been obtained. Prior to this time, a variety of approximate analytical treatments had appeared in the literature, most of them being based on the approach originated by Lord Rayleigh (1878). Assuming the surface to be infinite in extent and making use of the periodicity in (say) the coordinate x, Rayleigh expanded the scattered field in discrete spectrum of outgoing plane waves, which representation was assumed to hold right down to the surface. Application of the boundary condition leads to a single equation, valid for all x, from which to determine the complex amplitudes of the scattered waves, infinite in number. It is at this stage that approximations must be made, and the literature is replete with methods all of which are similar in character. Thus, for a perfectly conducting sinusoidal surface, Rayleigh (1878) obtained a solution by successive approximation based on the initial neglect of all attenuated waves; for the same surface, Tai ( 1948) proceeded via an orthogonal mode expansion, followed by matrix truncation, and others have pursued essentially the same path. An analogous treatment for a periodic interface between two homogeneous media was developed by Rayleigh (1907), and applied to a sinusoidal profile by Pavageau (1963) and to a triangular profile by Bousquet (1963). At best, all such solutions are valid only for corrugations whose height is much smaller than the free space wavelength, and in an attempt to overcome this restriction, Meecham (1956) used a variational method to find the angular distribution of scattered energy for a perfectly conducting grating. The scattered field was represented as a linear combination of known solutions of the 1

1363-6-T wave equation whose coefficients were obtained by a least squares fit to the boundary condition at the surface, and the procedure was then applied to a triangular (or saw tooth) profile. A somewhat different approach was taken by Eckart (1953) and Brekhovskikh (1952) who both used variations on the Kirchhoff approximation, and by Senior (1959) who employed physical optics. In either case the determination of the scattered field is reduced to quadratures, and Senior showed that for a plane wave at normal incidence on a sinusoidal grating, the physical optics integral can be evaluated exactly to give the complex amplitudes of the scattered waves. As indicated in Appendix A of this Report, the same is true for both polarizations, and for oblique incidence as well as normal. It should be emphasized, however, that the solution is still approximate by virtue of the postulated surface field distribution, and the failures of the physical optics estimate of the surface fields are examined in Chapter III of this Report. A method which is quite distinct from all of the above was developed by Sivov (1964) who used conformal transformation and a consideration of the static limit to analyze reflections from periodic surfaces with shallow and deep corrugations. The procedure is similar to that recently employed by Millar (1969, 1970) to investigate the inherent limitation of Rayleigh's method. As first noted by Lippmann (1953), it is not in general valid to assume that the expansion of the scattered field as a discrete spectrum of outgoing waves alone holds over the entire scattering surface, and this fact was later verified by Petit and Cadilhac (1964) in the case of a sinusoidal grating. In any general treatment of the grating problem it is therefore necessary to allow'ingoing waves in the immediate vicinity of the surface. Without exception, all of these analytical attempts to determine the fields scattered by periodic surfaces are subject to approximation, either implicit or 2

1363-6-T explicit, and it is only with the use of high speed computers that have permitted the direct digital solution of integral equations that reliable results have been obtained. Most of the initial work in this area was carried out by the French investigators, such as Petit, Cadilhac and Wirgin, and was motivated by the desire for more efficient optical diffraction. In his early papers, Petit (see, for example, 1963) followed Rayleigh's approach in expressing the scattered field as a discrete spectrum of outgoing waves alone, leading to a matrix equation for the determination of the spectral amplitudes. Since the matrix was then truncated and inverted numerically, it will be appreciated that the method is no more than a digitization of that originated by Rayleigh. However, in later papers (Petit, 1967), the Rayleigh assumption was circumvented by using an integral equation formulation. Series expansions were adopted for the incident and scattered fields and the integral equation converted to a matrix equation which was solved numerically. Specific results were obtained for plane wave incidence on gratings with triangular profiles, and the efficiencies computed. A rather different approach was taken by Pavageau (1967) who derived the integral equation directly in terms of the unknown surface current. The equation was cast in the form of a nonhomogeneous Fredholm equation of the second kind and solved by iteration. Methods which are very similar to that which we shall use have been employed by Neureuther and Zaki (1969) and by Green (1970). The former considered scattering by periodic structures, either dielectric or perfectly conducting. The integral equation was obtained from Green's theorem and the modified Green's function expressed in either of two ways depending on the parameters of the surface. The first (space harmonic) representation is analogous to that used by Petit; the second consists of an infinite series of 3

1363-6-T Hankel functions, and was computed using a Mellin transform and an asymptotic comparison scheme. Numerical data were obtained for sinusoidal profiles. Green (1970) also used the space harmonic representation, but improved its convergence by summation techniques. Data were presented for the surface fields and diffraction efficiencies of perfectly conducting gratings with triangular profiles. The present work also employs the numerical approach, and is concerned with the scattering of electromagnetic waves by infinite, perfectly conducting, two dimensional periodic surfaces of arbitrary but continuous profile. Plane wave incidence is assumed, with either E or H polarization, and both normal and oblique incidence are considered. A representation of the Green's function is employed that is similar to that used by Green (1970), and the convergence improved still further by subtracting the dc terms. This has the added advantage of making explicit the behavior of the Green's function in a neighborhood of its singularity. Digital programs have been written for computing the surface and scattered fields, and the program listings are included in Appendix B. Data are presented for a variety of sinusoidal profiles, and the results compared with the physical optics predictions. 4

1363-6-T II FORMULATION OF THE INTEGRAL EQUATION FOR NUMERICAL SOLUTION We consider here an infinite, perfectly conducting periodic surface illuminated by a plane electromagnetic wave. Since the surface is assumed two dimensional in the sense of being independent of a Cartesian coordinate z, the entire problem is two dimensional, and the most general solution can be deduced from the particular solutions appropriate to incident plane waves having either E or Hin the z direction, i.e. parallel to the corrugations. In either case, the problem is essentially scalar. 2. 1 Formulation It is convenient to develop first the integral equation in the somewhat simpler case of a scattering surface of finite extent. Let S (see Fig. 2-1) be this surface, and surround it by another closed surface SR. Let i (r) and G(r r') be two scalar functions which are continuous, together with their first and second derivatives, on S and SR and throughout the volume V enclosed by them. Assume, moreover, that 0(r) satisfies the homogeneous wave equation (V 2+k r)= 0 (2.1) inside and on the boundaries of V, whereas ( V2 + k2) G(rr') = - 6 (r-r'). (2.2) Applying Green's theorem to the volume V, we obtain 5

1363-6-T S Ftictnofre V FIG. 2-1: Geometry for the application of Green's theorem. 6

1363-6-T G(rlrt) n' )r' )-(r) an, G(tr') dS' S+SR G(rr) (V k (r')-() ) = I G(rjr') (V2k2)O(r')-(r')(V2k2)G(r') dV' V = (r) (2.3) by virtue of the properties of the delta function. We now identify b (r) with a total field tot tot (r) = b i(r) + s (r) where i'(r) is an incident plane wave originated by a source at infinity (and therefore outside V), and Os(r) is the field scattered by the surface S. Since Oi S(r) must satisfy a radiation condition at infinity, its contribution to the integral over SR decreases to zero as the surface SR recedes to infinity, whereas,i(r) contributes itself. Equation (2.3) then becomes tot i Ca a(r Gr i/~ (r)= pi (r) + ( GirInr') o n (r)-(anr Gr) dS' (2. and in spite of the assumptions of an incident plane wave and a surface of finite extent, Eq. (2.4) is also valid for an arbitrary incident field and for a surface S extending to infinity. In the latter case, however, the proof is by no means trivial (Jones, 1952). The particular situation of concern to us is that in which S is doubly infinite and divides space into two regions. It is then sufficient to integrate over only the upper ('illuminated' ) side of the surface. We also assume that S is independent of the coordinate z, and this allows us to distinguish two particular cases according 7

1363-6-T to the polarization of the incident field. Mks units are employed and a time factor e suppressed throughout. 2.1.1 E (or horizontal) Polarization i lA If E = E z, the electric vector in the scattered field will also be confined - z to the z direction, and we can make the identification p(r) E t(r). (2.5) Since the normal derivative of Et is related to the surface current density K z z by the equation a tot 3 Et (r') -jw K K ), (2.6):n z - z - where u is the permeability of free space, Eq. (2.4) can be written as E (r) = E(r) (r))+ Er') SG( Z z z - antS (2.7) At the perfectly conducting surface S the boundary condition is E (r) - 0 (2.8) Z - and hence, on allowing r to approach the surface, Eq. (2.7) gives rise to the integral equation 8

1363-6-T fG(r.r')K (r')dS' = (r). (2.9) - - Wmj zS For a two-dimensional problem the free space Green's function is G(rlr') - H'(k r-r'l ) (2.10) (2) where H is the zero order Hankel function of the second kind. The final integral equation is therefore f (2) 4 i S (-r)H' )di I E 1 (r) (2.11) Z- 0 Wm(A Zfrom which K has to be determined. z 2.1.2 H (or vertical) Polarization i lA If H = H z, the total magnetic field is likewise in the z direction, and on z making the identification (r) = H (r), (2.12) _z - Eq. (2.4) becomes t(r)=H (r) + G(rr) a Htot(r)- Htot( rr) aG(rr dS Z - - n z z an S (2.13) The boundary condition at the perfectly conducting surface S is 9

1363-6-T a Htot (r) =O (2.14) ant - and on allowing r to approach the surface, Maue (1949) has shown that (2.13) reduces to 1 tot tot 1 H (r) = H (r) - Ht (r) G(r |r) dS. (2.15) 2 z z- - -ant - S The quantity H is the induced current density, and since this current flows z tangential to the surface in a plane perpendicular to the z direction, we can write tot A K (r) = H (r) t (2.16) t- z - where t is a unit vector tangential to S. Inserting finally the expression (2.10) for the free space Green's function, we have -a }. J K(r) H ( ) r-r ) da = 4j H (r) K (r), (2.17) -S i.e. I K(r)H(2) r-r' )cos( r-r')di 4j (r) - Kt(r). (2.18) t- 1 z — - - 2 S This is an integral equation from which to determine K. 2.2 Reduction of the Integral Equations for Periodic Surfaces We now make use of the fact that the incident field is a plane wave and that the surface y = f(x) is periodic with period d, i.e. 10

1363-6-T f(x + md) =f(x), m + 1, + 2, 3,.... (2. 19) As before, it is convenient to consider separately the cases of E and H polarization. 2.2. 1 E Polarization Let us assume that Ei -jk(x sinO -y cos) (2. 20) E (r) =e (2.20) z - a where 0 is the angle of incidence with respect to the normal to the mean surface. The integral equation is as shown in Eq. (2. 11), but since the integration extends from - oo to oo, this is not appropriate for a numerical solution. However, by invoking the periodic property of the surface, the integral can be reduced to one over a single period alone at the expense of a more complicated form for the kernel. From Floquet's theorem, we have that K (r + md) =K (r) ejkmdsine, m=,3.. (2.21) z- z - Moreover, d'= 1+ ft(xt) dx' where f'(x) is the derivative of f with respect to x. This allows us to express the integral as one along the x' axis, and since -rl = (x-x'-md)2+ (y_-y)2 11

1363-6-T we obtain f Ho (k x-x'-md)2+(y-y')2' -)e jkmdsinO l+ {'' 0 m`-o 4 -jk(xsin-ycos)2.22) valid for 0 < x < d with y f(x), y - f(x'). In arriving at (2.22) we have implicitly assumed that interchanging the order of integration and summation is valid (a fact which is by no means obvious), but even so Eq. (2.21) is not a very promising equation for numerical purposes because of the extremely poor convergence of the infinite series constituting the kernel. To rectify this situation, consider P m H(2) (k (x-x-md)2+(y-y)2' )ekmdsin. (2.23) mn =-oo Using the Fourier integral representation of the Hankel function, the Poisson summation formula (see, for example, Morse and Feshbach, 1953) applied to (2.23) gives 2m y-y' 2 x) 1 j -j( - +k sinO)(x-xt) -j y-I X P d X e e (2.24) m=-oo m where X =k (2mr + k sinO)2 (2.25) m d and the chosen branch of the square root is that for which ImX <0 m12

1363-6-T in order to satisfy the radiation condition. The expression (2.23) for Plhas been used by Yen (1962) and Green (1970), and is equivalent to that obtained by Neureuther and Zaki (1969) using a space harmonic expansion. The features that should be noted are that if yj y' the successive terms decrease exponentiallywithincreasing Iml for all m +kdsin > kd T27i 2 27 and that in addition there is an algebraic decrease (proportional to m-1) provided by the factor X in the denominator. m We can produce still a further improvement in the convergence properties of the series by separating out the zero-frequency (k=0) terms, and this has the added advantage of making explicit the behavior of P in the vicinity of the singular point x'=x, y'=y. Putting k = 0 in (2. 24), we have c2m 2m 2m(y7y, ) P2 d 1 -XI) -- - m(y-y 21 e e (2.25) 1 k= d m-oo jiml 27r and since d1o 27- (x-x') 27T I y-y', _? d' T - o log 1-e i, — 27rj m=1 -1 j -e (x-x) -r- I-Y' 2m =^ 2 -log 1-e 27rj (see, for example, Collin, 1960), provided y' y, it follows that 13

1363-6-T Pl = (m=O term)- 2 log[ {-e 1-e _d (2.27) Hence, by subtracting the zero frequency value of each term of (2.24) and then adding (2. 27), the integral equation (2.21) takes the form dr~~ ~jkx'sine 2d jkycosO f G(x,y; x',y')K') e (x') = e 0 <x<d, (2.28) where G(x,; xt,' = k 1 e-jk ly-y' Icos0 G('' T' kcose -d! g 1e -2ed [y-y'l j(x-x' | |y ty' -j(x-x' jdoo - Y-YT I XF 2 -2 log 1 -e Jd m 1 X + 2m7r oD.2m7 -)' ( -jl Ix yY-Y'1 m d _ ~_ _a (x-xx) j yy m d + d e e m d (2.29) 14

1363-6-T with 2 2mr 2 X+ k2 - (- +k sin e) m d X =k2 2m -ksine)2 (2.30) m d We observe that the Green's function has now been expressed as the sum of a term involving explicitly the amplitude and the angle of incidence, a logarithmic function representing the true behavior of G, in a neighborhood of its singularity, and two series which are themselves differences of two convergent series. 22. 2 H Polarization The incident field in this case is taken to be i -jk(x sin0-y cosO) (2.31) H (r) - e (2.31) Z - and the reduction of the integral equation (2. 17) to a form suitable for numerical solution proceeds in much the same way as for E polarization. Since a 1 -f'(x) a + a Eq. (2. 17) can be written as P2Kt(x')dx'- 4jHz(r)- Kt(r)} (3.32) 0 ~t (iz ~Kt~_15

1363-6-T where 00 -d - dLf X x)(d +ksin0)-X sgn(y-y') i X m=-oo 2j 1 f(x 2m'r -j( —-+ksinO)(x-x') -j |y-Y'l X d m e e (2.33) in which sgn(T) is the signum function and Xm is given in (2.25). The zero frequency (k=O) limit of P2 is 2|k=p d 2 f'lx')+j sgn(y-y') e e m=-oo 00oo-j 2mir O.2m2r 2f i d-77 -j-(x I' }-x' - y-y'I =d f(xI)+j sgn(y-y') + fe(xt)+jsgn(y-y9 e m=l 2f 2-" m{(x-x'~-jly-yi' + -f'(x')+j sgn(y-y') ex')+j m=l and since oo jmT _1 j T 2 ei = - 2- 2 cot m=l 16

1363-6-T we have P2 k= f' (x)+ { f'(x')+j sgn(y-y')} cot{ (x-x')+j y- 17 jk~ru + -f?(x)+j sgn(y-y') cot d (x-x')+jd ly-y' (2.34) Hence, by subtracting the zero-frequency value of each term of (2. 33) and then adding (2.34), the integral equation (2.32) takes the form G2(x, y;x', y') Kt(x)dx' =jd {2H ()Kx)-K ), O<x <d, (2.35) where G2(x, y;x', y')= -j sgn(y-y')-j {f'(x')tan0-sgn(y-y') e-ksin(x-x')-jkcos Iy-y'1 - 2m xx,) -j I Y-y' IX I f' x (-dW+ksinO) -X s gn&y-} -e 2mwy-y { jf'') sgnM yy)}] -jjksine(x-x') e d x') m -je e [ L X+ m=l m f'(x')(2 d7+ksin) -X+ sgn(y-y')-e d jf(x')-sgn(y-y') d nm _je-jksin0(x-x') e d m=l e 2mnr i f.(x')(-2m7r^+ksin) -Xmsgn(y-y')t d Y Yj(x )+sgn(Y-y 9 + f't(x')+jsgn(y-y')| cot f- d (x-x')+j y-y'l I 2 { -f'(x')+j sgn(y-y')} cot { d(x-x')+j d l y-y' | (2. 36) 17

1363-6-T In spite of the obvious complexity of the Green's function G2, its form is directly analogous to G1 in consisting of certain explicit terms, a pair of cotangent functions representing the true behavior of G2in a neighborhood of its singularity, and two series which are themselves differences of two convergent series. Once again, therefore, the singularity of the Green's function has been separated out. We further note that the integral equations (2. 28) and (2.35) have been derived without approximation. In consequence, the formulation so far is exact for any two dimensional perfectly conducting periodic surface which is smooth in the sense of having a continuous first derivative. 2.3 Numerical Procedures Methods for the numerical solution of integral equations have been extensively discussed in the literature (see, for example, Harrington, 1968). The general procedure consists of reducing the equation to a finite set of algebraic equations, i.e. to a matrix equation, and can be illustrated by considering d G(x,x')K(x')dx' = F(x), O<x<d. (2.37) -0 We assume that the unknown function K(x') can be expanded in terms of linearly independent base functions n (x') such that N K(x') an (') (2.38) n=1 where the a are the associated constants. Substitution of (2.38) into (2. 37) n gives N d anI G(x,x'))n(x')dx' = F(x) (2.39) n n n=1 18

1363-6-T and the solution of the integral equation has now been reduced to the determination of the constants, a, n=l, 2,... N. n There are several possible ways of finding the a, e.g. least squares fit, Galerkin's method, and the collocation method, and it is the last of these that we shall employ. The collocation method converts (2.39) into a system of N linear equations by forcing the two sides of (2. 39) to be equal at N sampling points in the interval (0, d). This is simply a point matching procedure and results in the matrix equation a r G(x,x')0 (')dx' = F(x ), 0<x <d, m=l,2,...N. (2.40) /, nj m n m - mn-1 0 There now remains the problem of choosing the base functions n (X'), and here again there are several possible choices, e. g. rectangular, quadratic and sinusoidal. By appropriate choice, we can economize in the number of sampling points required for an'accurate' approximation to the solution K(x), and experience has shown that a rectangular function is not in general a good choice, whereas sinusoidal interpolation often works rather well. The particular form of sinusoidal interpolation that we have adopted is predicated on the use of sampling points which are uniformly distributed in 0 < x' < d. The range of integration is therefore broken up into N increments of length A= d/N (see Fig. 202). Furthermore, let x' be the midpoint of the 1n A A n'th cell, i.e. x = (n- ) i, and let Ax' denote the interval x' < x' <x'+ n 2 n n 2- -n We assume that A +B sin k(x'-x )+C cosk(x'-x ), if x' Ax' K (x) n n n n n n (2.41) 0 otherwise 19

1363-6-T yI,A-n n'th cell _____ _^^. — ___^I I x~ I I FIG. 2-2: Illustration of the interpolation procedure. 20

1363-6-T and specify the constants An, B and C by continuing the form appropriate t: n n x'EAx' out to the centers of the adjacent cells and imposing continuity. Defining n K K(x ) and K K(x ) we have n n n 1 n+ 1 K =A + C, n n n Kn+l =A +B sin k A+C cos k A n+l n n n Kn = A - B sin k A+ C cos k A, n-i n n n from which we obtain -K +2K cos k A- K A n+1 n n-1 A = n 2(cos k A- 1) Kn+ -Kn-I B = n+: n-i (2.42) n 2sin k A K + -2K + Kn_ nei n n-i n 2(cos kA- 1) Substitution of (2.41) and (2.42) into (2.40) now gives N' K { cos k A-cos k(x'-x ) >,J mG(XmX) [Kn cos kA- 1 n=l Ax n -sinkA+(cos k-l)sink(x'-x )+sin k Acos k(x'-x ) n n + n+1 22sin kA(cosk A c- 1) -sink A-(cosk A-l)sink(x'-x )+sink Acosk(x'-x ) n n n-l 2 sin k A (cos k A- 1) =F(x ), (2.43) m where x =(m — ) A, m = 1,2,... N. We note in passing that when n = 1 or N m 2 21

1363-6-T the periodicity of the problem must be used to determine the constants K_1 or KN+1 required in (2.43) The above procedure is immediately applicable to the integral equations for E and H polarizations on inserting the appropriate values for the Green's function and the forcing function, but a brief comment is desirable concerning the treatment of the singular cell. The Green's functions of concern to us are singular when x = x', the singularity being logarithmic for E polarization, and a first order pole for H polarization, and it is therefore necessary to modify the numerical scheme when x = x In line with the usual practice, we divide the m n singular cell into three portions: A E E ~ A (x —,x x - )x +2), (x +2, xn+- ) n 2 n 2 n 2n 2 n 2' n 2 with 0< c < A. The first and third segments are handled by the standard numerical technique, whereas the central portion is treated analytically by means of a limiting process (Andreasen, 1964). 22

1363-6-T III APPLICATION OF THE NUMERICAL TECHNIQUE AND PRESENTATION OF DATA 3. 1 Discussion Although the only data that we shall present are for sinusoidally corrugated surfaces, it is well to remember that the formulation given in Chapter II is quite general and appropriate to any perfectly conducting,two-dimensional periodic surface. It is not even necessary that the surface have continuous first derivative, and providing care is taken of any singularities at sharp edges, the procedure is applicable to surfaces whose first derivative is only piecewise continuous as in the case of a saw-tooth profile. The particular numerical scheme that has been adopted is especially suitable for surfaces having 2ra/d > 1, where a is the peak deviation of the surface from its mean, and d is the period. In the infinite series representation of the Green's function (see Eqs. (2.29) and (2. 36) ), the first few terms (for which i" + sin e|< 1 ) d - are oscillatory, but for larger m the terms decrease exponentially at a rate which is ultimately determined by the quantity 27ra/d. Whereas the earlier terms correspond to propagating modes in the scattered field, the later ones correspond to evanescent modes, and it is only necessary to include a few of the latter in order to achieve adequate accuracy. It is therefore apparent that for surfaces of relatively small period ( d X) the infinite series can be approximated by the first few terms alone. The unknown function in the integral equation is the current induced in the surface, and this is the quantity that is computed initially. Knowing the surface field it is, of course, a trivial matter to determine the scattered field, and the means of obtaining this in the form that is most convenient for our purposes is described in Section 3.2. 23

1363-6-T As noted earlier, all of the numerical results given in this report are for sinusoidal surfaces whose profile is taken to be 2wx y = a cos -d (3.1) Several different combinations of d/X, a/X and incidence angle are considered, and the distributions of surface and scattered fields computed. Comparisons with the predictions of physical optics are also given. 3.2 The Scattered Field When a plane wave is incident on a periodic surface, the scattered field can be represented as an angular spectrum of plane waves, which spectrum is discrete by virtue of the periodic nature of the boundary condition at the surface. Each of the infinity of waves making up the spectrum has associated with it a diffraction angle which may be real or complex and is determined by the grating law. Whereas the amplitude of the wave is a function of the profile size and shape, and the directions of incidence and diffraction, the diffraction angle depends only on the value of d/X and the direction of incidence. A finite number of the diffracted waves represent propagating modes and these are the important ones far from the boundary. The remaining modes are evanescent and though these do not serve to carry energy away from the surface, they do play a vital role in affecting the amplitudes of the propagating modes. The number of modes that propagate can be determined from the expressions for X given in Eq. (2.30): if X- is real, the corresponding mode propagates m M without attenuation, whereas if Xi is pure imaginary, the mode is evanescent. m To find the (complex) amplitudes of the diffracted waves we proceed as follows. Let y = f(x) be the profile of the surface and, for convenience, assume the incident plane wave to be E polarized. We then have 24

1363-6-T i -jk(x sinO-y cosO) z (see Eq. (2. 20) ), and by invoking the periodicity of the surface, the scattered field can be written as oD x + XmY) Es = A e m (3.2) z Im m=-o for y > max. f(x), where 2m7r 3 = + k sin O, m d and 22 x k =k-1f3 m m The field arising from the currents induced in the surface is given by d Es(x, y)= Li Pl 1+ (If'(xI)2 K (x') dx' (3.3) where P1 is as shown in (2. 24). In particular, this is valid for y > max.f(x), and hence, by combining (2. 22), (2. 24), (3. 2) and (3.3) we have 2 -j(3x+Xy) fd e 1m +m - 1+(x-xl)+Xm y') f 2K d A e e +'(x' K (x.)dx. m m -oo from which we obtain _ Fe j(1m x'+XmY') v' 1+) f 2 m rf'') K 2x m J Kx')dx'. (3.4) m 0 If the incident wave is H polarized, the procedure is directly analogous. For the incident magnetic field shown in (2.31), we write 25

1363-6-T oo -jx3 x+X y) Hs= Am m (3.5) z m m=-oo with Ad jOfI+X, _ Al — I e X- f(x') K(x')dx (3.6) 2dX | eJ(mX'+X mY') -l{) Kt(x')dx'. (36) m 2dX j | m t Having determined the current distributions K (x') and Kt(x'), it is therefore a trivial matter to compute the amplitudes of all the diffracted modes, both propagating and evanescent, and the amplitudes of the propagating waves for several different values of d/X, a/X and 0 are given in Tables III-1 through III-4. A widely used technique for estimating the field which is scattered by an object is the physical optics method in which the surface field is approximated by its geometrical optics value. Although the estimates are deficient as regards the polarization characteristics of the scattering, the method often provides an adequate approximation of the scattering if all dimensions (including radii of curvature) of the scattering surface are large compared with the wavelength, and because of its great convenience, the method is commonly employed in rough surface scattering analyses even when all radii are not large. A matter that is then of some debate is whether shadowing should be taken into account, and, if so, how. According to the physical optics method the current induced in the surface is K = 2 ^H (3.7) in the illuminated region, and K = 0 (3.8) A in the geometrical shadow, where n is the outward normal to the surface. 26

TABLE I-l: Amplitudes of diffracted waves: d= 0.2X, a = 0. 1X. E WAVE H WAVE 0 = 0~ 0=30~ 0=60~ 0=00 0=30~ 0=600 Mag. Phase Mag. Phase Mag. Phase Mag. Phase Ma. Phase Mag. Phase IAO1 0 (O) 1 | (0) IIAo l 00() IAl 0(OI ( A 0 (~) i umerical - oNumerical 1.0050 50.81 0.9940 44.33 1.0000 25.90 0.9830 -12.45 0.9910 -23.82 1.0640 -67.24 Solution _. _.. Physical optics 0.6425 0 0.7251 0 0.9037 0 0.6425 0 0.7251 0 0.9037 0 w/o Shadow ____ __ _ __ PhysicalOptics 0.6425 0 0.1253 15.96 0.1770 19.41 0.6425 0 0.1550 -65.17 0.5440 13.69 w/ Shadow...,,_ TABLE III-2: Amplitudes of diffracted waves: d = 1. 9X, a = 0. 25 X, 0 = 0~. E WAVE H WAVE _ Mag. Phase Ma Phase Mag. Phase Mag. Phase 1A01 o1 v =(~) | I 0(~), =-1________4A...........l...... = 11,, ~L) Nume ri cal Numerical 0.4920 19.73 0.6630 -82.95 0.8820 -77.78 0.3405 52.93 Solution Physical 0.3042 0 0.4389 90 0.3042 0 0.4389 90 Optics _____

TABLE 111-3: Amplitudes of diffracted waves: d = 0.4X, a = 0. 2X. E WAVE H WAVE 00~____ _ = 60_ 6 =0~e 0 = 60 Mag. Phase Mag. Phase Mag. Phase Mag. Phase |A j?(~) |A ( (o) (At of(0) |Atj| o 0'() Numerical 0.9998 -80.81 1.0230 49.89 0.9350 72.81 1.0140 34.72 Solution PhysicalOptics 0.0549 0 0.6450 0 0.0549 0 0.6450 0 w/o Shadow ____I ______________ PhysicalOptics 0.0549 0 0.4750 60.37 0.0549 0 0.1290 48.68 w/ Shadow __________ —_-__________^__ _____________ I.-' 00 E WAVE H WAVE H Mag. Phase Mag. Phase 1A01 0 (0) IAl| _(_) Numerica8.2098-.5 Solution 10000 8_120.96 Physical 0.9647 0 0. 9647 0 Optics TABLE 111-4: Amplitudes of diffracted waves: d=0. 2 X, a=0. 03 X, 6 =0

1363-6-T In many instances, however, some of the failures in the resulting estimates of the scattering can be traced to the discontinuity in the surface field created at the shadow boundaries, and improved estimates can be had by using (3.7) over the entire surface, shadowed as well as lit. In the present case of a periodic surface, the amplitudes of the diffracted modes predicted by either version of the method can easily be obtained by inserting the postulated expressions for the current into Eqs. (3. 4) and (3. 6). A computer program has been written to carry out the calculation, and the appropriate results are included in Tables III-1 through III-4. It will be noted that the results produced by either version of the method show very little agreement with the exact data. In general, the physical optics values, with or without shadowing included, tend to be too small, and we further note that if shadowing is excluded (or is not present), the physical optics esimates are the same for both polarizations (see Appendix A). It must be admitted, however, that for all of the surfaces considered, the minimum radii of curvature are too small to provide reasonable hope for the physical optics method to be accurate. With a sinusoidal surface, the minimum radius of curvature, Pmin, occurs at the peaks and troughs, and the radius increases to become infinite where the surface crosses its mean. For the four different surfaces embraced by Tables III-1 through III-4, the values of kp. are as follows: d = 0.2X, a= 0.1X, kp.=0.064, mmin = 1.9X, = 0. 25X, = 2.30, = 0.4 X, = 0.2X, = 0. 127, = 0.2 X, = 0.03 X, = 0.212 To judge from this listing, one might expect the physical optics method to be most successful for the second surface (Table 111-2), but in practice it works best for the fourth surface, which has a much smaller radius, but a small 29

1363-6-T enough value of a/X to be nearly plane. Before leaving this discussion of physical optics, we note that if shadowing is ignored or, because of the incidence angle (0< tan -1 ), is not present, a'- )snt p27ra simple analytical expression for the amplitudes of the diffracted modes can be obtained by an extension of the treatment in Senior (1959). An outline of the procedure is given in Appendix A, and the results obtained from these formulas are in agreement with the computed values in the last line of each of Tables III-1 through I-4. 3. 3 The Current Induced on a Sinusoidal Surface We here examine the nature of the surface field distribution for five different sinusoidal surfaces, considering first the results for E polarization and then those for H polarization. In each case, the modulus and phase (in degrees) of a normalized current are plotted as functions of x/X over a single period of the surface running from peak to peak. The horizontal scale is therefore the horizontal distance in wavelengths and not the distance along the actual surface. The particular current normalization chosen is such that for E polarization the quantity plotted is k K, the physical optics value of which is ^W2 2 (cos 0 - f'(x) sin0)e-jk(x sinO-f(x) cos 0) whereas for H polarization the quantity plotted is simply Kt whose physical optics approximation is 2-jk(x sin0 -f(x) cos 0) In each of the following figures, the exact computed values are shown as circled points, and are joined by a broken line only to guide the eye; the physical optics approximation is shown as a solid line. 30

1363-6-T 3.3. 1 E Polarization For a surface having d = 0.2X and a = 0. 1X, the results for the three incidence angles 0 = 0 (normal incidence), 300 and 600 are given in Figs. 3-1 through 3-3, respectively. We observe that for this relatively small period most of the current is concentrated in the vicinity of the surface peaks, with the current being almost zero in the troughs. There is, indeed, almost an exponential decrease in the current modulus away from the peaks, and the main effect of increasing the incidence angle is to scale the curves, leaving the general shape unchanged. The phase is somewhat more sensitive to 0, and whereas the curve is almost flat for 0 = 0, the structure increases noticeably with 0. Since kp = 0.064 for this surface, it is not surprising to find that the min physical optics approximation bears no resemblance to the exact data. This is particularly true of the modulus; the exact curve is a great deal more simple than in the physical optics one. The physical optics phase is also poor for 0 = 0, but agrees better for the larger 0, at least in an average sense. The effect of decreasing the amplitude of the surface while keeping the period constant is illustrated in Figs. 3-4 and 3-5 for the case of normal incidence. It is observed that as the height decreases, so does the current concentration near the peaks, but even for a as small as 0.01X (Fig. 3-5) there is still amost a 2: 1 variation between the peak and trough values. The phase, on the other hand, is much more nearly constant, and is more akin to that for a flat surface than a sinusoidal one. Results for d = 0.4X and a = 0. 2X with 0 = 0 and 600~ are shown in Figs. 3-6 and 3-7. By comparing Fig. 3-6 with 3-1 it is seen that doubling both d and a has little effect on the modulus of the current, but has a marked effect on the phase. This is true also of the curves for oblique incidence (cf. Figs. 3-7 and 3-3). 31

Phase (degrees) Normalized Current Modulus o(S o o o or I0 W l —, l — I _ l l-l-l-lI at01 I I I / /. M s f f E p w, ot. X I /I ^ - I ^ M- I FIG. 3-1: Modified surface field for E polarization with d = 0.2X, a = 0. I. 6 = 00: — o —, exact; ---—, physical optics.

Phase (degrees) Normalized Current Modulus! t3 0 tI3 C3 4 U 0 0 0 0 K ~\ W 1 ~ >0e 1 - - o / x / 0 I FIG. 3-2: Modified surface field for E polarization with d = 0. 2X, a= 0. 1X and 0 = 300: — o —, exact; /-,physical optics. -, physical optics.

Phase (degrees) Normalized Current Modulus I II o00 0 o0cs / 0' P1, physical optics.

Phase (degrees) Normalized Current Modulus 0 I I 0I - I~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~I I / I / Q FIG.304:Modiied urfae fild or Epolaizaton wth d a0 and e. 0, physical optI c ^~~~~~~, C I~~~~~~~~~~~~~~ o h~~~~~~~ x <> ~ Q ~~~~~~~~~~~~~~~~I ^ N \^ ~~~~~~~~ FIG. 3-4: Modified surface field for E polarization with d = 0. 2X, a = 0. 03X and 0 = 00: o —o, exact; -—, physical optics.

Phase (degrees) Normalized Current Modulus | I:~ I 0 0 0 00 IS3 h-~. - - -'; / oo.^ \ ~~ ~~~~~~~~~~~~~~~~~~~~~~~~/'a FIG.3-5: Modifed urface eld for E polarization withd=0.2X, a = 0O1X and e00: - exact;.physical optics. w~^ o \ \. tO M FIG 3-:Mdfe ufc il o oaiainwt.X.1 n ~ o-,ea --- phscl pis

Phase (degrees) Normalized Current Modulus 0 It~~~ ~ I'' vt2;, o o O O ~- 1 Co r> IO~~ / F. 3-6 Modified surface field for E polarization with d. 4X a. 2 and 0 0: — o —, exac I I / i \/ ~ - ^xI \ i FIG. 3-6: Modified surface field for F polarization with d = 0.4X., a = 0.-2o and 6- 00: -o-, exact; -, physical optics.

Phase (degrees) Normalized Current Modulus c o o O O 10~~~~~~~~~~~'I.lo o p..0 I. 0 0 0- -^, 0^ / ^><~ j/A 0 FIG. 3-7: Modified surface field for E polarization with d =0., a =0. 2A and. 600: — o —, exact; /, physical optics. CO^" C)! \ \ FIG. 3-7: Modified surface field for E polarization with d = 0.4X, a = 0. 2X and = 600: — o —, exact; ---—, physical optics.

1363-6-T All of the surfaces so far considered have had rather small values of d/X and a/X, and have not been such that one could expect physical optics to provide a reasonable approximation. For a surface with a much larger period, d = 1. 9 X with a = 0. 25 X, at normal incidence, the appropriate curves are shown in Fig. 3-8. Physical optics now constitutes a more accurate approximation, and though there are still noticeable discrepancies, qualitatively as well as quantitatively, the physical optics curves do approximate the exact data in a mean sense. The exact current modulus still has a maximum at the surface peaks, but the minimum now occurs where the surface slope changes sign, i. e. where the surface crosses the x axis. The phase remains fairly constant over the lower part of the concave portion of the surface, but changes rapidly from positive to negative and vice versa where the modulus has its minimum. Since the particular case considered in Fig. 3-8 has also been treated by Neureuther and Zaki (1969) using their numerical program, it is appropriate to compare the data. This is done in Fig. 3-9. The agreement is very good in spite of the rather different numerical procedures involved in the two methods. 3.3.2 H Polarization For each of the situations considered above, we have also computed the surface fields when the incident wave is H polarized, and we now present these data in the same sequence as we did for E polarization. The results for incidence at 0 = 0, 300 and 600 on a surface having d = 0. 2 X and a = 0.1 X are shown in Figs. 3-10 through 3-12. Eight sampling points were used in the computations. We observe that the curves for the modulus are roughly sinusoidal in form and are relatively insensitive to changes in 0, with at most a slight reduction in amplitude as 0 increases. In contrast to the behavior for E polarization, however, the maximum now occurs in the surface trough, with a minimum at the peak. There is again no agreement between the exact data and physical optics. 39

Phase (degrees) Normalized Current Modulus II 000 F0 LIGx. 3- 8: ---—, physical ol0~ ~ ~ ~~~a s. ~0\ -r\ r \ ftwIft FIG. 3-8: Modified surface field for E polarization with d = 1. 9X, a = 0.25X and 0 = 0~: — o —, exact;, physical optics.

Modulus o / — I )' <\ //'\ // S \\ // \\ // \\ //\// 40 1 Q) I Q) <D \ ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~/ ~~~ N~ %.475.950 1.4751 1.900 I 0 Cd I / ri\ / H -40 -- z~ I -40S ~ -1~~- 1- 7'Phase -80 — 2 * Tong and Senior X o Neureuther and Zaki (1969) FIG. 3-9: Modified surface field for E polarization with d = 1. 9X, a = 0. 25X and 0 = 00.

Phase (degrees) Normalized Current Modulus Lpi' F.. i.... I.../ 0 \ \tC / O~\/ v -/ J / FIG. 3-10: Modified surface field for H polarization with d = 0. 2X, a = 0. 1X and 0 = 00: — o —, exact; ---, physical optics.

Phase (degrees) Normalized Current Modulus 0 0 0 0~ t t / \ / \ / \ I x. I I I I I If\ / / FM. 3-11: Modified surface field for H polarization with d = O. 2A, a = O. 1X and 8 = 30~: — o —, exact; /-, physical optics./ ^ I > ^~ ^ ^ I// \/ / \ / FIG. 3-11: Modified surface field for H polarization with d = 0. 2, a = 0. IX and 6 = 300: — o, exact; ---, physical optics.

Phase (degrees) Normalized Current Modulus I I I I 0! ~" —- I I p /J 0 / \ / \ / / \ / x I \\ I / \ I I I I I!-' FIG. 3-12: Modified surface field for H polarization with d = 0. 2X, a = 0. IX and e = 600; — o —, exact; -—, physical optics.,physical optics.

1363-6-T The effect of changing the surface amplitude while keeping the period constant is illustrated in Figs. 3-13 and 3-14 for the case of normal incidence. As the surface height decreases, the numerical results approach the physical optics values, a fact which is not surprising since physical optics is most effective for a planar surface. As in the case of E polarization, the phase curve is relatively flat even for a = 0.03X, and does not show the variation displayed by the physical optics plane. Figures 3-15 and 3-16 show the results for d = 0.4X and a = 0. 2X with 0 = 0 and 60~, respectively. By comparing Fig. 3-15 with 3-10 it is seen that doubling both d and a has a somewhat larger effect on the current modulus curve than was the case for E polarization. In particular, there is a notable change of shape, with the minimum no longer occurring at the surface peak. Similar changes occur for 0 = 60~ (cf. Figs. 3-16 and 3-10). Once again, there is no agreement between the exact data and physical optics, either in modulus or phase. The current distribution for a surface with larger period ( d = 1. 9 X and a = 0.25X) at normal incidence is presented in Fig. 3-17. The modulus has a marked oscillation, quite distinct from that found for E polarization, with three maxima per period. One of these is at the surface trough, and the current is again a minimum at the surface peak. The agreement with physical optics is poorer than for E polarization, with neither modulus nor phase being approximated to any real extent. The exact data are, however, in agreement with those obtained by Neureuther and Zaki (1960) for this case, as evident from Fig. 3-18. 3.4 Energy Conservation In solving problems of electromagnetic scattering from perfectly conducting surfaces, a common procedure for checking the accuracy of the solution 45

Phase (degrees) Normalized Current Modulus 0cn ocn0 0 o t0 I I ^~~~~~~~~ C b r ^^ 0^~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~I y^ A co FIG. / -13: <odified vurface JLeld for ii polarization with d = 0. 2X, a E 0.03~and 8 = 00: — o —, exactC --- I I I. I.( FIG. 3-13: Modified surface field for H polarization with d = 0. 2X, a = 0. 03X and 60=00: — o —, exact; —,physical optics.

Phase (degrees) Normalized Current Modulus I CD 0 I.... 0 CD 0 0' ^I I o 0 I c ~ x. I FIG. 3-14: Modified surface field for H polarization with d = 0. 2, a = 0. O1X and 0 = 0~: — o- exact;, physical optics.

Phase (degrees) Normalized Current Modulus \ /J / I I I \ I0 s3 // tP ^ i3 o) y1-' 0/ 0I co. i^"^ J j I \ ^x. \ IIJI / / N! \ FIG 3-15: Modified surface field for H polarization with d FI. 3-15: Modified surface field for cs. polarization with d = 0.--, physical optics.

Phase (degrees) Normalized Current Modulus 00 rf^ ~4^ 00! I 00a~~, s, I I i!0010 I^ / FIG. 3-16: Modified surface field for H polarization with d = 0. 4X, a = 0. 2X and 0 = 600: — o —, exact; physical optics. a,~~ FI. -6:Mdiie urae iedfo i oarzain ih =04X =0 21 n 8=60: --- xJ t ~ —---, physical optics.~~~~~~~~~~~

Phase (degrees) Normalized Current Modulus 0~ / MO 0 I II,woo \ — I. — physical optics. c)0) FIG. 3-17: Modified surface field for H polarization with d = 1. 9X, a = 0. 25X and 0 = 00: — o —, exact; physical optics.

80 -4, Modulus / \ ~/ 30 - o. L L5/ -, x / / X 80 /- \ \ FIG. 3-18: Modifd s e f475u 950 1. 4250 0 4I /4 \ 0,vl \/ \ Cd ~~~-40 - J ~~~~/ \\ c.k ~. ~ 0 I I I \ p. 0 0 0 P47.950"1.ase - 1-9 -80 * Tong and Senior X o Neureuther and Zaki (1969) FIG. 3-18: Modified surface field for H polarization with d = 1. 9X, a = 0. 25X and 6 = 0~ FIG. -18:Modiied urfac fied fo H plariationwithd =. 9X a =0. 5 n = 00.

1363-6-T is to determine the degree to which the law of conservation of energy is fulfilled. This is essentially the check employed by Petit and Cadilhac (1964), Neureuther and Zaki (1969) and Green (1970), as well as by many others. However, as pointed out by Amitay and Galindo (1969), energy conservation does not provide a measure of accuracy of a solution found by the Ritz or other related methods; and hence, in order to use energy conservation to check other than computational round-off errors, care must be exercised in choosing the method of solving the integral equation numerically. The method that has been adopted in this Report is a point matching and sinusoidal interpolation scheme, and since the solution does not automatically satisfy energy conservation, we are able to use this as a check. The relation that must be fulfilled is N I |Am Re (Xm) = k cos 0 (3. 9) m=o where N is the number of propagating modes. The extent to which this is satisfied by the numerical results that we have obtained is shown in Table 111-5. For comparison purposes, the percentage errors implied by the physical optics method, with and without shadowing, are also included. Whereas all of the exact data satisfies the energy conservation law to within (about) one percent or better, the only case where the physical optics solutions does so is when d = 0. 2X and a = 0.01X. Such a surface is very close to being planar. 3. 5 Convergence Check In solving an integral equation numerically it is always desirable to carry out a convergence test, a) to determine the number of sampling points necessary to achieve the required accuracy, and 52

TABLE 111-5: Conservation of Energy E WAVE H WAVE Physica Optics Physical Optics Numerical w/o Shadow w/Shadow Numerical w/o Shadow w/Shadow co I tI to e; Io; to 0; 1 d a of of of of of of k o(~) Energy Error Energy Error Energy Error Energy Error Energy Error Ener Eor 0.2 0.1 0 0.9970 -0.70 0.4128 -58.72 0.4128 -58.72 0.9848 -1.52 0.4728 -58.72 0.4728 -58.72 0.2 0.1 30 0.8564 -1.18 0.4553 -47.50 0.5075 -41.45 0.8542 -1.43 0.4554 -47.30 0.0208 -98.00 0.2 0.1 60 0.5013 0.26 0.4084 -18.32 0.6172 23.44 0.5075 1.50 0.4083 -18.34 0.1479 -70.42 0.2 0.01 0 0.9990 -0.10 0.9921 - 0.79 0.9921 - 0.79 0.9980 -0.20 0.9921 - 0.79 0.9921 - 0.79 0.2 0.03 0 1.0042 0.42 0.9308 - 6.92 0.9308 - 6.92 0.9987 -0.13 0.9308 - 6.92 0.9308 - 6.92 0.4 0.2 0 0.9976 -0.76 0.0030 -99.70 0.0030 -99.70 0.9823 -1.77 0.0030 -99.70 0.0030 -99.70 0.4 0.2 60 0.5068 1.36 0.2061 -58.78 0.1147 -76.47 0.5042 0.84 0.2061 -58.78 0.0083 -98.50 1.9 0.25 0 0.9996 0.04 0.4202 -57.98 0.4202 -57.98 0.9974 -0.26 0.4202 -57.98 0.4202 -57.98 _~~..J

1363-6-T b) to test whether or not the numerical solution approaches a stable value as the number of sampling points is increased. Although it has not been proved mathematically that the accuracy can be improved by increasing the number of samples, it would appear reasonable to believe so on a physical basis. Figures 3-19 and 3-20 show the results of a convergence test applied to the case d=0.2X, a=0. 1X and 0=0 for E and H polarizations respectively. It can be seen that the solution does appear stable, and that surprisingly good results are obtained even with as few as 4 sampling points. Two other aspects of the computer program that should be mentioned are the computational time consumed and the number of terms included in the summation of the infinite series for the modified Green's functions. For E polarization, a typical figure for the computational time on an IBM 360 machine using 12 sampling points is about 60 seconds. This time includes the computation of the amplitudes of the diffracted modes, the energy conservation check and the calculation of the physical optics results. A typical figure for H polarization using 8 sampling points is 28 seconds. The time increases roughly as the square of the number of points, and most of it is eaten up in the computation of the matrix elements. Except for the large period surface (d=l1. 9X, a=0. 25X), where the number of terms retained in the Green's function series was 6, only the first three terms were required. Tests that were performed indicated that three was sufficient. 54

'sulod 2ulIduiBs Zg put 01'c 2uisn poandudoo 00 = 0 pUB YXT *0 = V'XZ'0 = P t^M uojvziielod a joj plaFij;aoJns paTJTpo~A:6T-'013I X/x 1._____________^____________, 0 OI / \V / Z _ 0 0 (D * 1~~5_1~~ ~~ /v \ ~~/x V I - I/ \ II _ I\ 1 I 1,II fod ]IHT]rluOd TTd -- E 0 / s~u~od ufnduxBs a - V S L-9-9~ET

1363-6-T 3 b —-....-~~. — 2 ^ ~. / ) A// A - 4 sampling points 2 X - 6 sampling points m^^ eL - 8 sampling points ~ 1 _ 0 l x.1 <.2 x/\.1.2 ~o vx/X bn i -5 -5z.4 to — 10 -10 FIG. 3-20: Modified surface field for H polarization with d = 0. 2, a = 0. 1X and 0 = 0~ computed using 4, 6 and 8 sampling points. 56

1363-6-T IV CONCLUSIONS We have here considered the problem of a plane electromagnetic wave incident on a perfectly conducting, two dimensional periodic surface, and have developed numerical procedures for the direct digital solution of the integral equations for the surface fields. By using special summation techniques followed by the subtraction of the dc term to improve the convergence of the series for the modified Green's function, a relatively efficient procedure has been arrived at, and this has been programmed for a computer. Data have been obtained for the surface fields on several different sinusoidal surfaces when illuminated by E and H polarized plane waves at oblique as well as normal incidence. It is found that the polarization has a marked effect on the field. Since most of the surfaces considered had relatively small values of d/X and a/X, there is little agreement with the physical optics estimate of the surface field. However, this is also true for the one surface of much larger period that was examined. Knowing the surface field, the amplitudes of the diffracted waves in the discrete angular spectrum representation of the scattered field can be computed, and this has been done using the exact surface fields as well as the physical optics estimates with shadowing either included or ignored. Here again the physical optics predictions are deficient, and whereas the results derived from the numerical program satisfy the conservation of energy law, the physical optics values do not. In the continuation of this study, the numerical program will be applied to a wider variety of periodic surfaces in order to build up an understanding of the surface field behavior, including its dependence on profile size and shape, and the direction and polarization of the incident field. It is our hope that this will 57

1363-6-T enable us to develop synthesis procedures for the efficient computation of the fields scattered by periodic surfaces of large amplitude and period, leading ultimately to better prediction techniques for rough surfaces in general. 58

1363-6-T REFERENCES Adachi, S. (1965), "The Nose-on Echo Area of Axially Symmetric Thin Bodies Having Sharp Apices," Proc.IEEE, 53, 1067-1068. Amitay, N. and V. Galindo (1969), "On Energy Conservation and the Method of Moments in Scattering Problems," IEEE Trans., AP-17, 747-751. Andreasen, M.G. (1964), "Scattering from Parallel Metallic Cylinders with Arbitrary Cross Sections," IRE Trans., AP-12, 746-754. Bousquet, P. (1963), "Diffraction des ondes electromagnetiques par un reseau a profil triangulaire," C.R.Acad. Sci., Paris, 257, 80-83. Brekhovskikh, L. M. (1952), J.Exp.Theor.Phys.,USSR, 23, 305. Collin, R.E. (1960), Field Theory of Guided Waves, McGraw-Hill Book Co., Inc., New York, NY. Eckart, C. (1953), "The Scattering of Sound from the Sea Surface, " J.Acoust. Soc.Am., 25, 566-570 Green, R. B. (1970), "Diffraction Efficiencies for Infinite Perfectly Conducting Gratings of Arbitrary Profile," IEEE Trans., MTT-18, 313-318. Harrington,R. F. (1968), Field Computation by Moment Methods, The Macmillan Co., New York, NY. Jones, D. S. (1952), "Removal of an Inconsistency in the Theory of Diffraction," Proc.Camb. Phil. Soc., 48, 733-741. Lippmann, B.A. (1953), "Note on the Theory of Gratings," J. Opt. Soc.Amer., 43, 408. Maue, A.W. (1949), "Zur Formulierung eines Allgemeinen Beugungsproblems durch eine Integralgleichung, " Z. Physik, 126, 601-618. Meecham, W. C. (1956), "Variational Method for the Calculation of the Distribution of Energy Reflected from a Periodic Surface, " J. Appl. Phys., 27, 361-367. 59

1363-6-T Millar, R. F. (1969), "On the Rayleigh Assumption in Scattering by a Periodic Surface," Proc. Camb. Phil. Soc., 65, 773-791. Millar, R. F. (1970), "On the Rayleigh Assumption in Scattering by a Periodic Surface-II, " Proc. Camb. Phil. Soc. (to be published). Morse, P.M. and H. Feshbach (1953), Methods of Theoretical Physics, Part II, McGraw-Hill Book Co., Inc., New York, NY. Neureuther, A. and K. Zaki (1968), "Numerical Methods for the Analysis of Scattering from Nonplanar Periodic Structures," Alta Frequenza, Special Issue on the URSI Symposium on Electromagnetic Waves, Stresa, Italy. Pavageau, J. (1967), "Equation integral pour la diffraction electromagnetique par des conducteurs parfaits dans les problemes a deaux dimensions; application aux reseaux," C.R.Acad. Sci., Paris, 264, 424-427. Pavageau, J. (1968), "Sur la methode des spectres d'ondes planes dans les problkmes de diffraction," C.R.Acad. Sci., Paris, 266B, 1352-1355. Petit, R. (1963), "Diffraction d'une onde plane monochromatique par un reseau periodique infiniment conducteur," Comptes Rendes, 257, 2018-2021. Petit, R. and M. Cadilhac (1964), "Etude theorique de la diffraction par un reseau," C.R.Acad. Sci., Paris, 259, 2077-2080. Petit, R. and M. Cadilhac (1967), "Sur la diffraction d'une onde plane par un reseau, infiniment conducteur, " Comptes Rendes, 264, 1441-1444. Rayleigh, Lord (1878), The Theory of Sound, Vol.II, Dover Publications, New York, NY (1945). Rayleigh, Lord (1907), "On the Dynamical Theory of Gratings, " Proc. Roy. Soc., 79A, 399-416. Senior, T. B. A. (1959), "The Scattering of Electromagnetic Waves by a Corrugated Sheet," Can. J. Phys., 37, 787-797 and 1572. Sivov, A. N. (1964), "Reflection of Electromagnetic Waves from Corrugated Surfaces with Small Periods, " Radio Eng. and Elect. Phys., 10, 1509-1515. 60

1363-6-T Tai, C-T (1948), "Reflection and Refraction of a Plane Electromagnetic Wave at a Periodic Surface, " Harvard University Cruft Laboratory Technical Report No. 28, Cambridge, Mass. Wirgin, A. (1966), "Diffraction par un reseau de profil rectangulaire illumine par une onde plane associee a un vecteur champ electrique polarise' parallelement aux sillons," C.R.Acad. Sci., Paris, 262, 385-388. Yen, J. L. (1962), "Multiple Scattering and Wave Propagation in Periodic Structures," IRE Trans., AP-10, 769-775. Zaki, K. (1969), "Numerical Methods for the Analysis of Scattering from Nonplanar Periodic Structures," Ph.D. Dissertation, University of California, Berkeley. 61

1363-6-T APPENDIX A: PHYSICAL OPTICS APPROXIMATION Although the main purpose of this report is the development of a numerical technique for the determination of exact values for the surface and, hence, scattered fields, we have found it desirable to compare the results obtained with those provided by the physical optics approximation. This approximation is a physically-based one which postulates an explicit form for the surface field arrived at by assuming each element of the surface to bear that current which it would were it part of the local tangent plane. The calculation of the scattered field is then reduced to quadratures. In many instances, however, and a periodic sheet is one, an analytical evaluation of the physical optics integral is a difficult procedure, particularly in such cases where part of the surface is shadowed, and some of the shortcomings of physical optics estimates in general certainly arise from sloppiness in the evaluation of the integral. In still other cases, the physical optics result proves to be more accurate if shadowing is ignored* (see, for example, Adachi, 1965) and it is of interest to observe that for a sinusoidal sheet it is then possible to produce an exact evaluation of the integral. Needless to say, however, the resulting scattered field is still subject to the unknown errors inherent in the use of the physical optics approximation, and to the neglect (if present) of all shadowing effects. The procedure is directly analogous to that given by Senior (1959) for the particular case of an H polarized plane wave at normal incidence on a sinusoidal sheet and consists of three steps: This modified method is sometimes called'extended physical optics' 62

1363-6-T (i) writing down the physical optics integral for the scattered field without any shadowing included; (ii) asymptotically evaluating this expression appropriate to an observation point at large distances from the sheet; and (iii) matching this expression to a discrete angular spectrum of waves to obtain their amplitudes. Since these amplitudes are independent of the field point, we have, in effect, produced an exact evaluation of the integral valid certainly in the half space above the sheet. Let us take the equation of the perfectly conducting sheet to be y = a cos Kz (A. 1) where a is the amplitude of the corrugations and 2ir/K = d is the period. If the incident field is E-polarized, we write i ^ -jk (x sin -y cos0) (A - ze (A. 2) (cf. Eq. 2.20), implying i A A -jk (x sino -y cos 0) H -Y (cosx + sin 0 y)e- - (A.3) where Y is the intrinsic admittance of free space By virtue of the periodicity of the sheet and, hence, of the problem as a whole, the scattered field can be expressed as a discrete spectrum of waves which waves are certainly outgoing as regards y> a. Thus, we have o -jk(x sinO + y cos ) E- z A e m (A. 4) m m=-oo where 63

1363-6-T ksin = mK+ksin, In m k cos = k2 -(mK+k sin0), (A. 5) m and that branch of the square root is chosen having imaginary part non-positive. Application of the physical optics approximation, followed by steps (i) through (iii) above now gives.m k Ka A = _ cose-mKa sin J ( m) (A. 6) 2m m yk-(mK+ksine) m wrhere m = a cos +{k2-(mK+ksin 0 ).(A. 7) In particular, A =-J (2a k cos 0), O o and as a -- 0, A -1, Am- 0, m 0, in agreement with the known solution for a flat sheet. If the incident field is H-polarized, we take Hi A -jk(x sine -y cos ) (A.) H z e (A. 8) implying i 1 A A -jk(x sin -y cos 0) E = (cos 0 x + sin ) ejk( sin -y cos (A. 9) and expand the scattered magnetic field as a o -jk(x sin 0 +y cos 0 ) H= z A1 e (A. 10) m=-oo 64

1363-6-T where 0 is as before. On evaluation of the physical optics integral, we m now obtain A' =A. (A. 11) m m The equivalence of the results for the two polarizations is consistent with the known fact that the physical optics approximation for a perfect conductor is inherently polarization insensitive. Results computed from Eqs. (A. 6) and (A. 11) agree with those given in the appropriate columns of Tables ID-1 through III-4 and obtained from a numerical evaluation of the integrals in Eqs. (3.4) and (3.6). 65

1363-6-T APPENDIX B: COMPUTER PROGRAMS Two separate computer programs have been written, one for E polarization and the other for H polarization. Each program consists of the following parts. 1. Main Program The functions which this performs are a) reading in the input data consisting of the amplitude and period of the surface (measured in wavelengths), the angle of incidence, the number of sampling points and the number of terms employed in the series summation; b) calling the subroutine which computes the matrix coefficients, including the analytical treatment of the singular term; c) computing the physical optics current; d) calling the subroutine which inverts the matrix and hence obtains the induced surface field; and e) printing out the appropriate numerical results. 2. Subroutines a) Subroutine CAVE This determines the appropriate constants in the sinusoidal interpolation formula and calls the integration subroutine CSIMEQ to integrate the kernel function over the sampling cell. 66

1363-6-T b) Subroutine FUN The constants in the sinusoidal interpolation formula are here computed. c) Subroutine SUBC This computes the kernel function of the integral equation. d) Subroutine CIMPS The numerical integration is performed using the Newton-Cotes method. e) Subroutine CSIMEQ This performs the matrix inversion using a Gauss reduction method. f) Subroutine ENERGY The amplitudes and phases of the diffracted waves are computed and the scattered energy found. 67

-~ —----.... —y-' _ _ - --- -- _ _^ / - - -.. T'='... -... 7 0.... aW_0'-7 0.ST GNON S44,S T= 1ROS TONN, $RIIN *WJATFnR.tnMPI __ __ __F C THTS 1S THE MAIN PRnGRAM FOR THF STUDY OF EM SCATTERING FROM A PERFECTLY C CONNnJCTING,PERIOnTC SIRFACE FOR THE E WAVE DnMENSION A(?n,?n),R(?O),W(?n),Z(20),ZMAG(20), PHASE(20),IZW(20) 1,n T;T(?n),n(2n),PHAS(20),ZM A(0),nnST(20). COMMnN /DATA&/X,N,iAMP,MM,THFTA /DATA2/ XI,XK,PI I,CS,SS,CSN,DCS 1,nDSS,CTH,STH,C.SS,SINC,PIDPSO,NFA, TKAME/DATA3/COY, FUNX/DATA4/A ___?,W CnMPLEX A,,W,TE,CEXPCOY, DAG, FUNX,TG, SUM1, TPO,CAB,STORCONJG _ __ EXTERNAL CIGUN _ RFAD 1, XDAMP,THFTA,N,MM_ TF (N,GT. 1nn) GO Tl 9 PRINT 1, XDAMP,THETA, NMM_ __ ____ 1 FORMAT( 3Fin.5, 2 T 5 -_ C TO CALCIULATE OUT SOMF CONSTANTS FOR COMPUTATION ACIJY=n.02 __ I TN T=__ PI=3. 14159265 PII=2.*PI _______ XN=N c..g_ nFL=XD/XN __ rS=cs (PTI*nEL)' SS=S N( ( P I I *OLnFL ) _. ---- rC'N=CS-. - S TH=SIN ( THETA ) PT __ r C.TH=CnS(THETA)*PII - CSS=?. *CSN*SS P S=P I **2 FL TA=DEL/ T10.______ PIn=PI T/XnI - ---- - - ----------- TY= (2.,. )* nFLTA/CTH TG.=(n.,1. )*STH R7= XD/( 3.*PI**?)*25.___ PnCS=COS(?.*PTI*nFL) nSS=SIN(2.*PT I*rFL)_ AME=AMP*P l

C T TR TFAT THE SERIES SUMMATTnN IN THE SINGULARITY TERM R!IMl=(O.,n. 1 M 7 = I n MM nn 21 I=i,MZ X TI - Al=PSO A?=(-PTTIXI/XOn+STH) **? IF (Al.LE. A? ) rn Tn o01 StIMl=SIJMl4+ I./SORT(Al-A?)+Xn/XI/PTII*(0.,l.) Gn Tn 2n3?01 SIJM1 =RIJM1-(1./SORT(A?-A1 )-Xn/XI /PIT)*(O.,..) 2n0 A4=( P II*XT /Xr+STH) **? TF ( Al.LE. A ) GOn TO 2n? SIJMl=SlIMI+l. /SORT( A-A4)++Xn/X /PI*(0.,1.) Gn Tn?1?n? I IMI = ijMl - (. /SORT ( A4- 1 )-Yn/x /P T l )* i ( o. ti. )?1 CN TIN IE )nAG=TY+XD/P T*(n.,1.)*(4.*nFLTA*ALnr,(OELTA)-DFLTA*?.) +SI-Ml*2.* i nFLTA C Tn CLEAR THE MATRTCFFS nn T=1,N _ __' nn 7 K=1.,N A( T,K =( n.,. n r r C'nNTINIIE C F Tn FTLL UP THE MATRICES nn 3 I=lN X=T. =....__ ____. ______ XI I=XIDEL XT= ( XTI-. 5 )*nEL F 7( T)=36S.*XT XIX=XI-HDEL CA-1_ CITMPS( XTX,XI I, ACIY, INT,,lIJN, STnR, INDEX) TF (I.FO. i) Gf T n ____ _ _ _____ _ ____ __ nIT(:T )=nnSTT-I) +ARS( STnR) /?. nnST( ) =nnST( T-1. )+CARS(STOR) nn TO 706?n5 nns T( T )=CARS (STnR) OTST( )=DOST(T )/?. 2n6 SINC=CnS(PID*XIT _ R ( T ) =RZ*CEXP ( ( n -1,-.)*-STNC*CTH-*AMP ) nn 3 K=I,N XK=K XK= ( XK-0. 5 )*nFL

C, Tn CAL,C.LA TE PHYSTCS OrPTTCS CI.RRENT CAR=(.,-1..1) XK*STH-AMPCOnS PI DXK)*CTHN) PO(K)=-lnon,/6./pT/SORT(1.+(AMF*SIN(PID*XK) )$*2-)(AME$SIN(PID*XK)* ISTH-C TH) /P I T*CFXP(CAR) IF (K.EO. I) cGn TO 4 R1 =XK-HnEL R?=XK+HDFL CALL CAVE ( R1,R?, CIJY, NT, 1.2,3 ) rn Tn 3 4 1R=XK-HnEL R?=XK-nELTA CALL CAVE (R,R?, AC1Y, INT,1,?, 3) R1=XK+OELTA R2=XK+HDEL CALLt CAVE(R R?,ACIlY,INT, l,2,3) TF= TG*XI ____ _, A(I,I)=A(I,I)+SORT(1.+(AME*SIN(PID*XI))**2)*CEXP(TE) *DAG c.3 CONNTINUE _ C Tn INVERT THF MATRIX TO GFT SOLJT ION -- -- --- - ---- - CALL CSIMEO(N,N,A,W,R, IZW) -- nn-?51- i —1-..... W(T)=CnNT,.G(W(I )?51 CONTINlUE CALL ENERr.Y C Tn PRINT nIOT THE NFCFSSARY RESUILTS PRINT 3? 3?7 FnRMAT (lnX,?4 I,RX-,5HZ ELFCT. DECREE,RX,6HI REAL,12X, IHI IMAGIN IARY) PRINT 33,(T,Z(I),W(I),I=1,N 33 FnRMAT ( OX, 3, lX,F7. 2,3X,F?0.5,2X,F17.5) nn 35 TI=,N 7MAGC(I)=CAfS(W(I ) ( ) _ __ __ _ _____ 7 M Ar, ( T =C A R S ( W ( I ) ) 7MA(I ) =CARS (P ( ) PHAS;F( I =( Rn./PT )'ATAN(AIMAC(W( I ) )/RFAL(W( I ) )) PHAS(I)=(lRn./PI)*ATAN(AIMAr,(Pr(I))/REAL(PO(I))) 35 CONNTINUE PRINT 36..,

3, FnRMAT (lnX,?H T.RX,-1H7 FLFC.DEF RFFRX,l lHI MA NITIJTDE,RX,7HT PHAS IF:) PR ThT 37, ( 1,7 ( t ),7MAC,( T.PHASE( T ), 1=1,N) 37 FnRMAT( lnx, n, 1 nX,F7.?,,qX,F 0. 52XF.2) PR INT 43 43 FORMAT (?7HRRIILTC FnR PHYSICAL nPTICS) PRTNT 4? 4 FnRMAT (n1X,?H T,PX,l HSIIRFACF I STANCF, RX,6HI REAL,12X,11HI IMAGI 1 NAR V ) DRTNT 44,(T,nT T (T1),Pn( I ) I =1,N) 44 FnRMAT( lnXT, X,FR.,3X, F?n.5,X, Fl7.5) PRIN\T 4f 4h FNRMATfTlOX,?N TIRX,16HSIJRFACF DISTANCE,-8X, 1HI MAGNITIJDF,arX,7HI PH 1 aF) PRINT 47,(i,nTST(T ),Z A( I ).PHAS( I ),1=1N) 47 FnRMAT(InX,ITqnX,FR.q,~X,F2), S2XFll.2) ~ I' Gn Tn 19 9 CTnP FNO

U11RRnlI TINE CAVF( AA, RR,CC,JN,NA,NR,NC) nT MFNS TnN A(?n,?),W (?0) r.nMMnN /DATA1 I/Xr,N,AMP, MM,TRFTA / ATA/ X I,XK,PT,P I,CS,SS,CSN,ODCS 1,nSSCTH,STHrSS,STNCPID,PSO,NFA, I,K,AME/DATA3/COY,FIJNX/DATA4/A 7,w COMPLFX A,cnY,FI'NX, TAl, TA2, TA3,W ___ ___ ____ __ _ FXTFRNAL FUN Q1 =AA AC II=CC TNT=JN N F A=NR CAIL (,TMPS(R 1, R?, Arl!Y, I NT, FtIJn-TA- I. INDEX-) NlF A=NR CALL CTMPS(R,R7?,Art!Y,INT,FIJ)N,TA2, INDEX) NF A=NC __ CALL CTMPS(RA,R?,ACr. IIY, NT,FUN!,TA3,INDEX) A( T,K)=A(I,K)+TA1 _ -'T'IF (K.EO. i) cn Tn - I TF (K.EO. N) C,n TO? A( I,K+I )=A (I,K+1 )+TA? - A(T,K-l)=A(T,K-l)+TA3 _ GOn TO 3 A( T,?)=A( I,?)+TA? _ A(T,N)=A( I,N)+TA3 Gn Tn 3? A(T,N-1)=( (T,N-1 )+TA3 A(T,1 )=A-(,1)+TA? 3 RFTI RN FNO

rnMPl FX FIINCTT nN FIIN(X) nTMENSTnN A(2?n,?O),W(20) rnMMN/ /DATA/ /XD,N1,AMPMg,TFTA /DATA2/ XI,XK,PIPIT,CSSSCSN,DCS 1,nS,CTH,STH,.CfSS, ST C, PI n,PSO, NFA, T,K,AMF/DATA3/COY, F(NX/nATA4/A?,w CnMPLFX A,W,FIINX, CnY XW=PT T*(X-XK) TF (NlFA. FO. 1 ) rn Tn i TF (NF.F:O. 2) GC Tn. TF (NFA.En. 3} Gn Tn 3 1 FN= (CS-C nS(XW)) /CSN nO TO t1 2 FN1 = (-SS+ST N( Xw)cJ C +SCns ( x r ) / (CSN-2?."SS),n" Tn 11 - FN1=(-SS-STN(Xwt —CfSN+SS*Cfns(XW)) /(CSN*2.*SS) 11 V=X rACl. stICRr, (V) I IN= FN1 *Fl INX._ RFTIIRN.. FNO i H

S, iRRn TINE S iR, (V) nTMENcTNN A(?n,?0),W(?n) CrnMMnN /DnATA1/Xrn,N,AMP,MM,THETA /DATA2/ Xl,XK,PIPIT,CSSS,CSNDCS I nfSSCTH,STH,CCC,SINC, PI DPSONFA,I,K,AME/DATA3/COY,FI JNX/DATA4/A 2,W r.nMPLEX A,FIJNX,CnY,XXl,XX2,C1,C2,F?,F3,YV2,TETlYYZ IYZ2,ZZl ZZ2 l,CEXP,CLnG, XIl, F F,W YY=AR ( S I NC-C nS ( P ID V) )*AMP v=Cny*YY y?=( n., 1. ) *YY XI-=(n.,1. )*(V-XT ) nI=-PTD (((!.,n.)*YY+(n.,1.)*(XT-V)) c?=-Pin ((.,. )*YY-(O.,1. )(xI-v) )) T=( (.,. )-CFXP(C1 ) )*( (1.. -CFXPC2) ) F2?=CFXP(Y)/CTH+(n.,1. )*X/PI T`CLn( Tl) F3= (.,n. ) nn I M=I,MM XM=M XX=-PT F)XM y TH= ( XX+STH ) **2 TF (PcO.LE. XSTH) CnA Tn? 771 =SORT( PSO-XCTH) (1.,0. ) r.n Tn 3 771=SORT(XSTH-PSO)(0.,1.) ) 3 XSCH=(-XX+STH)**7? TF ( PO.lF. XS'H) cr TO r,. 6 77 =SORT(PSO-XSSH)*( 1.,. ) (nn Tn 5 4 772=SORT(XSSH-PS.)*((n.,1.) 5 xX1=XX*XIJ XX?= (0.,-1. )*XX YY1 =YY*XX Y71=YY?*ZZ1 Y7Z=YY7*ZZ? TF=EXP(YY )/XX? FFR=CEXP(-XXI )*(CEXP(YZ1)/ZZ1-TF)+CEVP(XXI)*(CEXP(YZ2)/ZZ2-TE) F3=F3+FF3 _ 1 CnNTTNIJE 6 H=1.+(AME*STN(PIO*V))**? FIINX=( F2+F3)ORT (H ) *CFXP( (n.,-. ) STH*V) R F TI iRN FNn

CIRP nlI TT NF C r MPS( A,, P F. MAY I Y, NTr,, Y P I NDEX) C. HIH Srll RnlilTINF TEn nn THF TNTGRPATInN C.n Pl FX FAFR, X TP,. X,I, FI4JFX,XTS,, X INTC, F F R = X TNTr (A) + X T T (R) YH=R-A YITR=XH.5 X,j=XTR*FAFR Y MFW = A+X I R YHA=XIR/ 3 TNnFX=n 10 ( FNPWX=XTNTG,(XNFW) Onl IF ( INDFX.GT. } rn Tn lnn0 100n ThinFX=l. yT=XHA*(FAFR+FNFWX*4. 1004 Xl=(X,I+XI'*3. 1.25 TNDFX=TNnFX+1 TF ( INnX.,T. MAX! ) GO TO 1009 1010 XH=XH*.5 XhIFW=A4^XH*. 5 _ c=(r).,0.) nn5 F( XNFl.LT. P ) GO TO 1 00r 1nn7 XTP=(XJ+XH2.*)/3.`c c 0T F ( CARS (XIP-XT.LF. CARS(F*XIP) GO TO 1009 0nnP yI=XIP Gn Tn Inn4 n103 S=S+FNEWX XNEW=XNEW+XH rn Tn 0nn5 10nn RFT)IJRN FNn

S1IRR nTl TINEr CrTMFO ( TnI M,N, A, X, JPRM). THTS StIRRFIJTTNF TS Tn INVERT THF MATRIX CnMPl_EX A, X,R,nFT,9I IM,C nTMEN cnN A(TniM,TnTM m) X ( I)),R( M),JPRM(IDIM) nFT=(!.,n. ) nn 1! T=l,N 13 IPRM(T)=T - C F IN THE ELEMENT nF MAX MIIM ARSfnLUTTE VALtJE nn I K=1,N r.=A(K,K) TT=K I.,I=K nn 3 J=K,N 0nn? I=K.N TF (CARSI (C,)-CARS ( A ( T,)) ) )3,, 2 3 C=A(T,J) TT=T Ir = A J T. T T = t,t l=J 2? CnTTNIE ) OnFT=nET*C - TF (CARS1 (nFT))?n,?n,30?n PRTNT inn inn FnRMAT (9H SINCGuLAR ) r.ALL EXIT 3n R()TT=R (II)/C C DFVTiF EACH FLEMENT nF THF TI=TH RnW RY C K Pn=K+1 nn 4 J=K,N 4 ( TT.T JI=Af TT,t)/ A( T!, JJ l )= ( I,n. )

TF (IT P. FO. K ) n.n Tn hn C CWITCI-T THF K=TH RnW Amn THF TI=TH RnW nn S.I=K,N C =h (K,,1 A(K, tI =A(, S /X(IT,J)=C C=R( K ) (K )=R( I T ( T )=C An TF (j1.FO,. K ) C- Tn 7n C STnRF THF InC A T TN OF THF.CLUtMN OF THF MAX PIVOT TT=JPRM(JI).PRM (.I) I =1 PRM ( K ).JPRM( K)=TT r. SWT TCH THF K=TH ANM, THF J=TH C(nLUIMNS n ^ T= 1,N C=A( T,K) A( I,K ) =A( I,1, 1 6 A(T,JJ)=r. C, GrF SFIiR=CnLlJMN nF ZFRnS IN THE K=TH CnLtJMN 70 TF ( K.FO. N) Gn Tn 1 nn 7 I=Kpn,N 0n R,j=KPO,N -i ( A(,J)=A( I,Jv -A( T,K)*AK,JK) R( T I)=R( I )-A(,K )R (K ) 7 A(T,K)=(n.,n. 1 CnNTT NIE C nRTATN THF X=S RY PACK SIRSTITI TITnN nn 9 T=I,N K=N-T+I K Pn=K+411M=(., n.) TF (KPn -N) 24,?4,?6 24 nn 1n J=KPnN S.tM= S IM+X ( ) A ( K,, J 1n rnNTTIIUF pS Y( K ) =R(K )-SIJM 9 rnNTTNIlE C nRnFR THE X=S nn t11 T=,N 11 R(T)=X(I) nn 1? I=,N TTI=,PRM(I) 17 X(I )=R(I) RF TIRN ENn RFAL FtUNCTTIN CARS1(Z)

.nMPLEX 7 TF ( RFAL (Z. Fn.n. AK n. TAMA.(Z ). FO.. ) GO TO 10 CARl 1=CARS( 7 Gn Tn?0 In rARl-=-n. 2n RF TIRN ENO ~~~~~~~-.~~~~~~~~~~~~~~~~~~~~~~3 y~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 00 0)

SiJRRnlU TINE ENER'GY C THT, CIJRRNOITTNE TS TO CHECK THF ENERGY CONSERVATION nTMFNSTnN A(?n,n),W(?o) nCMMNN /nATA1/Xn,i,AMP,MM,TFHTA /nATA2/ XIXKPI,PI ICSSSCSN,DCS 1,rS S,CTH, STH,C. SS,SI C,PI D,PSn,NFA, TK, AME/DATA3/C 1Y,Ft'JNX/IATA4/A,W rnMPl FX A,RR,CC, RBW1,CRW?,CFXD,Y1, v2,,Y3,Y4,CMPLX,CSORT,YY1,YY2, 1Q1 r1, S,W,C nYF, FI)XCl, r nCNJr, XN=N nFl =XO/yN nn 3 i1 N=!,I M=MN-1 XM=M YV =( PSO-(PTnIXM+PI T N (THFTA ) *? )( 1., 0. ) YY 2 =(P O-(-PT nXM+PT TSSIN(THTA ) *2)( 1.,0. R1 =CSORT (YY1 SRP=CSOR T ( Y?) Pnw!=(n.,. ) CR,?=(.,na. nn 3?? K=I,N -3 XIX=K Co YK= (XK-. 5 ~X *nFL Y1=(O.,1. )*(PTnr!xM+PI SIN(THETA))*XK Y?=(n.,1.)*(-PTDlOXM+PTT*SIN(THETA ))}XK - FP= COns( P I n*x) Y3=^MP*FR*SRI*(,l.,1...........)Y4=AMP*FR*SR2?(n., 1. PPR=l.+( hMP*PT n*I N (P In*XK) )**2 PRWI=CEXP(Yl )'CEXP(Y3)*SORT(PPR)*DEL/CSORT(vYY1)"W(K)+BRW1 CRW7=CFXP(Y2)*CFXP(Y4)*SORT( PPR) DEL/CSORT(YY2)*W(K)+CBW2 322 COlN TTNIJE RR=CnNJG(RRW ) CC=CnNJC (CRW?) AS=(PTI*2*3. /5./X D)**/P7 I WI FRR*RRWl ISR 1 AS W?=CCr*CRW2*SR,7*a PR TNT 333,W1,W? 3 3 FnRMAT(?F2n.7) TF (FAL (YY?).LT.. ) Crn Tn 444 31. 1 CONTINlUE 444 RFTIJRN FNn $nATA n.? n.0I 0.523h 2nO 3 S rTnP

$STr.NnN S44h T=1RnS' TnNG SRIIN *WATFnR nCOMPILE C THIS IS THE MATN PRnGRAM FOR THF STUDY OF EM SCATTERING FROM A PERFECTLY. CnNnICTTNG,PERInnTrC R,1RFACF FOR THE 14 WAVE nTMENSTnN A(?,?n),P(?0),W(?2),Z(?0),t7MAG(201,PHASF(20), IZW(20) I,nAG(?0),P" (n P (?0),R (? )RT(0), 20 MAG ( 20) tPASE(20) CNMMON /nATAl/XY,N,AMD,MM,THFTA /DATA2/ XI,XK,PI,PII,CS,SS,CSN,DCS 1,nSS,CTH,STH,CS 9,SITNr.,PInD,PSO,NFA,I,K,AME/DATA3/COY,FlJNX/DATA4/A rnmMPLFX A,R,W,TF,CEXP,CnY,DnAG, FJNX,TG, SIJMI,,TY,CMPLX,GC, CnNJG I,CAR,PI,P? 19 REAr I., XD,AMP,TH;TA.M,MM TF (N.rT. inn) Gn Tn 9 PRTNT 1., Xn,AMP,THFTA,N,MM 1 FnRMAT(3F0n.5,?T 5 C. Tn CrALC.ll ATE OiT SOMF C.ONSTANTS FOR CnMPUlTATInN T N T=? A.C IIY=n,. 05 PT=3.141 59? 5 6 PT T=?.*PT co YN=N nFL-=xn/N C-,=rcns( PI* lnFt) SS=SIN( PT I*nFL ) CS = C S- 1. CTH=SIN(THFTA)IPTl I C TH=Cn ( THETA )*PT _ _ _ CSs=?.*cSN*Ss P' 0=PT IT **7 nFL TA =FL / n. CnY=(0.,1. )*CTH HnfFL=nEL/2?. PT n=PI I/XD Tr,=(nI.,1. )*ST ncs=Ccns (2.*PTT*nFL) nsS=S TN(?.*P TT*FL) AMF=AMP^PID C Tn CLEAR THE MATRI CES nn? I=1,N __ __ nn? K=I,N.A T,K)=(n0.,. 2 CrnNTINIUE

C THTS PAP T F FTHF PPFlG, AM IS TO CA!C IIIATF TH-F PHYSICAL UlTICS CIIRRFNT nn 4 =+,N XK =K y< (y XK-. 5) ^rFL CAt =(.,-I. )'-(XK ST C f- t.AMP'.cn,(PIDn-YK)I,CTH) p I ( K ):-. / S R T ( 1. +( AMF ( P i n K ) )?) *CFXP (CAB) P7 ( )= PI (K ) ^ AtMF ST ( ( P D XK Qi (K) =SqNT(RFAL (Pl (K ) c2?PFA_ (P2(K ))^*2) RT(K)=CSOPRT(ATMAC,(P1(K))'-:?+ATMAG(P?(K))4*?) AMACA ( K ) =SOPT ( Rl ( K ) * 2+R ( K ) *72 ) PA F( K )-=tO./PTATAN(AIMAC (Pl(K) )/R AL(P1(K) ) 4 rnNTTNIIF PRITMT 5,(K,P1 (K),P?(K),AMAC,(K),PASF(K),K=1,N) F FnR MAT( (HK=,T 5, F1 5.5,?F 5.,P 5., F 15.3)) C TN FTLL tUP THE MATRTCFS _ _ nn 3 T=I,N XT=T rC=(O.,-1. )*STH*XT A( T,T ) =C, XP (GC)XD (.,.) 7 T )=3-.*XT_ _ ______ - o S T TNC =cn(P( P T n*x I - R ( ) =CFXP( (.,-1. ) I NCCTH*AMP) (., 2.)*XD nn 3 K=19,N XK=K YK= (XK-O. 5) ~*OFL R1 =XK-HfFL R?:=XK-DFLTA CALL CAVE(P,R7,AC.IIY, TNT,1,?,3) Ri =XK+OFLTA R?=XK+HnEi CA L C.AVF ( R,R, A C. I, NT,1,,3 ) 3 C.NtTINlUE r Tn TNVFRT THF MATRIX Tn GFT SOLIITINN CALL CrSIMEO(l,N, A,W,R, 17W) 0n 34 T-1,N w( = )=CNN, (W(TI) 34 rCnNTTNIJE

CA. I. FNFRGYC DR TNT 3? 3? FORMAT (Inn,?, p,RX,15H7 FLFCT. DFCRFF,RX,6 HI RFAL,1 X,1 HI IMACGIN I ARY) PRINT,( T, ( T ).W(T ), T1,N) 33 FnRMAT (nX, T3,1OX,F7.?,X,P?0.5,?X,F17.5) nn 3S I=1,N 7MAr( T )=C.ARC(W( T I PHAF( T )=( 1n. /P )*ATAN(A IMAr,(W( /RFAL(W( ) ) 35 c.nNTTNIIF PRTN T 3' 3A FnRMAT (lOnX H I,tX,I1'HZ FL C. OFGRFE, X, 11HT MAGNI TJDEF, 8X, 7H PHAS IF) PRINT 379 ( T HZ ( T 7MG( I T PHAF( 1 ), I=1, N).- -- 37 FnRMA T( I nX, T3, I nX,F7. 2, 3X, F0. 5,X, F1 1.2) nn Tn 19 9 S Tn FNO co tP\~~~~~~~~~~~~~~~~~~~~~~3 Q ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ C0 o*I 10 CT)~~~~~~~~~~~~~~~~~~~~

S IRQ ni TT NF CAV\/F ( AA, PR,CC r,, J. NA, NR, NC) nrTMFNSInN A(2n,)},W(?n) cnMMnN /nATA1/Xfn,N,AMPMMTHFTA /nATA2/ XI,XK,PI,PII,CS,SS,CSN,DCS 1,nFS,CTH,STH,C.<,SINC,P T n, PS,NFA,I,K,AME/DATA3/CO)Y, FIJNX/DATA4/A?,w rnMPLFX A,CnY,FIiJhX,TA1,TA?2,TA3,W FXTFRNAL FIIN sl=AA R?=RR Mh, IY=f C. TNT=JN NFA=Ng CALL CTMPS (R1, RACY,INT, FIN,,TA2,INDEX) NFA=NC. co CALL CIMPS(R1,P7,ACuIY,TIT,,FIINTA3,,INDFX) >A /( (I,K)=A( I,K)+TA1 TF (K.EO. 1) rCn Tn i _ IF (K.FO. N).n Tn? A(T,K+I)=A(T,K+I)+TA? A( I,K-1 )=(I,K-1 )+TA nn Tn 3 ( (I,?2)=A(T,?)+TA? A(T,N)=A(I,N)+TA3 r.n Tn 3 A(TiN-l)=A(IN-1 )+TA3 A(I,1)=A(I,1)+TA? 3 RFTtIRN FNn;

C.MPIFX FIINr(TTON FllJ(X) nTMFNiN A(l? n,7n),W(7n) CnMMnN /nATA!/yr),NAMPMM,T-FTA /nATA2/ XIXKPI,PTI,CSSS,CSN,)CS I.nS,CTH, CTH,CR, S INC,P I), PS.O,NFA,I,K,AMF/nATA3/COY, FINX/nATA4/A r.OMPLFX A,FIINX,CnY,W XW=PT T ( X-XK) TF (NFA.FO. ) rfn Tn 1 TF (NF.FO. 2) rn TOn TF (NFA.FO. 3) rn Tn 3 1 FMl = (CS-cnS(XW))/CSN r.n Tn 11 2 FNI=(-S+STN(XW )*r.sN+SSsn-Cfl(W) )/(CSN*?.*SS) On Tn ii 3 FNhll =(-SS. -STN(XW)*C sN+SS*CNs(yW )/(CSN*?.*SS) It V=X CALL SlIPC(V) FI IN=FNl F INX cI oo RFTt IRN_ _ F4 FNn H

SII RR nlITI NF I JRC (V) nTMFNI'nN A(?n,,n),w(?0) C("MMnN /f)ATA I /X,N,AMP,MM,TFH E TA /nATA2/ XI,XK,PI,PII,CSSS,CSN,DCS 1,n SS,CTH,STH,CS, SINC.,PI n,PSO,NFA, I,K, AME/DATA3/COY,FU NX/DATA4/A 7, J CMPi EX A,FIJNX.CnY, XXl, F?,P3,YY?,TE,,YZ1,YZ2,CEXP,XIJ,FF3 - -- 1S;N, AR P, ARN,FT3,CST rcn.S, ZZ, ZZ2,W yN=N TF (XK-XN/?. 17,17,19 17 AM1=AMP-CrnOS(PT1-V) nn Tn n _____ ____ _ 19 /M1 =AMP, nOS ( PT D*V) 2n AM?=AMPCO,,S(PTn*XT TF (AM1-AM2) 14, 13, I 13 S.TnP 14 SCN=(O.,-1. ) rn Tn R ___ 16 SrN=(n.,.) P OFK=-AME*STIN.(PT n*v ) YY=ARR( SINC-C nc ( PlBVi )vAMP Y=C.Y*YY YY?=(n.,1.ltYY F?=rEXP ( Y /CTH ( (.,1. )I FN-STH+CT*Sr.N) -SGN. F3=(n.,n., nn i M=.,M'Mo XM=M XX -P Tn*xM XS TH= (XX+STH ) T*7 TF (PrO.1F. XSTH) r0 TO 2 771=SORT(PSO-XSTH)*(1.,0.) r,n Tn 3? 771,=SORT(XSTH-PSO)*(n..1. 1.) 3 XSSH= (-XX+STH)*? T f PSO..F. XSCH) GO n Tn 4 777 =SnRT( PSn-XS *( 1..,nv rn Tn 5

4 772=SORT(XSSH-PRO)^( r., t.) 5 XX1=XX*XIJ YY1 =YY*XX Y71 =YY*ZZ1 v7?=YY?Z*Z2 TF=FXP(YYI) FF3 =CFX( -X ( CF X P ( Y ) / Z 1*( ( 0.,. ) F N ( X X + S T H ) + SGN Z I )-TF ((rJN-nFN) )+(CFXP(YZ?)/Z7?*((n.,. )*DFN*(-XX+STH)+SCN* 2 77? )-TF*(nF+tN+rGM)},CFXP(yX ) F3=F-3+FF3 1 CnNTTNtF AQP=PT/xnr,(v-!T )t~!. qn. )+(O.. }*PT/Xn*YY APR=PT/Xn*(XT-V)*(1.,^.)+(n.,1.)*PI/Xn*yy ApN= p!/xrn (XT- v)(.,n. )+ (.,l.),Pl/x4YY FT3=(n.,O.5)*((CCnS(ARP)/CSIN(ARP) *(-nFN+SGN)+CCOS(ARN)/CSIN(ARN)*( 1 nhFN+SN) ) Fl:NX=(F?+Ft+FT3)CF:XP( (n.,-1. )*STH*V) RFTIIRN PNP CJ) I H

,iRRnITITNF F FPP GY COnMM nN /nAT / JXn,N,AMP,MM,THFTA /nATA?/ XIXK,PT,PI,CS, SS CSN,DCS 1,nSS,CTH,STH,rS, cTNr,,PT n,PSO,NFA,I,K,AMF/DATA3/CnY,FIJNX/DATA4/A?,W nTMFNTInN A(?n,9n),w(?( - C.nMPl FX A,RR, CR,RRW,CRW?,CF XP,Y1,Y?,Y3,Y4,CMPLXCSORT,YY1,YY2, lt R, SR?,W,CY,FIJIlX, CnNJ rTHTS S1iRRlIJTINE TS Tn CHFCK THF ENERGY CnNSFRVATION YN=N nFL=Xn/XN nn 311 MN=t,15 M=MN - YM=M I1t=PT I*XM+DT T.STNfTHFTH tI?=-PTDXM+PTT *!IN(THFTA) YY1=(PSO-I1*4r)*(1.,.) - _-___ -_____- ___ —- -- YYI=( P. *~? ) * ( 1., ~ i1, ) YY?= ( PSO-1i? )$ ( 7.,n. SP1=CSORT(YY1 CR?=C SORT (YY7) 00 RRW1=(n.,n.) - 1nn f32?? K=I,N __ YK=n ~< K - _,,..,, _.. —.. - -...- ----------------- --- —. —........ -.-.. XK=K c XK=(XK-. 5) rOFl YV=(n.,1.)*(PTnOXM+PTI*SIN(THFTA))*XK H Y=(fn.,!.)*(-PTinXM+PTT*SIN(THFTA))*XK FP=AMP cns( Pi n XK) F\/=- AMP*S TN(PI nXK_ Y3=FR*$CRP(*n.,1.) Y4=F R ( O.,. 1. RRW=C FXP(Y! ) Ir, PX P{(Y,) K nFL/SRlR( s 1- - I I PDT*FV) -W (K )/XD/2'.+BBW1 CPW?=CFPXPXp(Y)CFXP(Y4)*nFL/SR?^(SR?-tJ2*P IFV) W(K)/xn/?.+CRW2 3? rnNTTltf\lF RR=CnNJr,(RRW1 ) rr =nNJ, (CRW7) W1I =R*RRW SR1 /PT T ^?=CC rCRW 2?*SR- P /P T PR TNT 333,W1,W? 3'3 FnRMAT(FP2n.7) TF (RFA (YY7).IT. n.) rn Tn 444 311 rnNTTNIIF 444 RF T IRN FNn

St I RRflTIN CTMP ( A,, FMAYI. IN TC,,XIP, INDEX) C THIST- S1IRRn ITTNF TIS FnP TNTF.rRATTON CFMPLFX FAFR,X TP,X.J, FNFWX,XT,., XI NTG.. FAFR=XTNTG( ^A)+XINTG(R) XH=R-A X TR= XH*. Y,= X IR*FAFR Y IFW = A + XI R XHA=XIR/.f. TNrtFX=0 1jnni FNFWX=XINTr,(XNFWI 1 nl TIF (INnFX rT. n0 ) GO Tn Inn00 n100? TnDFX= T=X-HA*(FAF R + F N W X*4.) 1nn4 Xj=(X.l+xT*3. *.-?5 T n FX = T N n FX + 1" TF (ITNDEX GT. MAXI ) rn TO 1009 co lO l n H=XHI*.5 XnFW=A+XH*n.5 1005 T F f XNFW. L T. R ) Gn TnO 100 1007 XTP=( XJ+XH*?. S) /3. TF ( CARS (XTP-XT).LF. CARS(E*XIP)) OC TO l009 i oR YT = X T P C,n TO nn04 10on03 =S+FNFWX XM FW=XNFW+XH c,n TnO 1005 100nn9 R F T IRN FNn

FNn SuI Pi TINE T MFn ( T n I M,,A,X,R,JPRM) C TH I Cl RR'nilTINFP Tr, FnR MATR X T i\IRST nN (f,rnmPt FY,..YR-,'nFT, S tM,.... nITENTFON (TintM,TnM), (TinT M) R( Ii, M,JPRM( InM) I)F T= (1.,. ) nn 13 T=1,N 13 IPRM( )=I r F Nln THF EL EMFNT nF MAX MIMIM ARSnLTTE VALUIE nn 1 K=l,N r=A (K,KI) TT=K,IJ=K nn?,=K,N nnr? T=K,N TF (CARSI(C)-CARs1 (A ( I,. J) ) ) 3.2,2 q r.-A(TJ) llTT,JJ=J 2 CDNTTTNtJE nF Tr)E FT*C, TF (CARSl ()nET))? n,?n,n _?n PRINT Inn 1nn FnRMAT (9H STNGUILAR ) oo CALL FXTT o') ~ n R(TT)=R(II)/C C nFVTF EFACH FLEMENT nF THF T =TH ROW BY C c K Pn=K+1 l nn 4 J=K,N 4 fTT,, ) =A (T, J )/C A( TIJJ)=(1.,.) TF (TT.EO. K ) Gn TO 60 C SWTTCH THE K=TH RnW Ain THE TT=TH ROW )O 5,=K,N r.=A (K,t A ( K,, J ) =A ( T T, ) 5 A T T j)=C C=R K) R(K)=R(TI) R I T)T=C An TF (JJ.FO. K ) Cn TO 7n C STnRF THE LnCATTnN OF THF COLIUMN nF THE MAX PIVOT IT=JPRM(JJ) JPRM(J)=I = PRM(K) IPRM( K )=I T

rC:, TCH THF K=TH ANn THF I=TH (,OLI MNS nf h T=1,N C=A( T,K) A( I,K)=A( T,.1.1) 6 A( I,TIJI=C C CF T SlIR=CnLIIMN nF Z FRllR TN THE K=TH CnLIJMN 7n IF ( K.EO. N) Gn Tn 1 nn 7 I=KPn,N r)n R =KPn,N RP ( I,J)=(f I,.1-A( I,K )hA(K,J) R( T )=R( T )-A( T,KI R ) (K 7 A( T,K)=(n.,n0 1 C.nh TINlJE C nRTATN THE X=S RY PACK SIJRSTITITTiONN... nn T=1,N K=N-I+ I Kpn=K+i st M=(n.,,. T'F (KPO -N)?4,2426 _ 24 nn 10 J=KPI1,N -.... __ __ —_ —S!JM=SUM+X ( )A (K, J) in rnNTiNIj)E? X(K )=R(K )-SIM. 9 C nNTINtE, C ORnFR THE X=S 3_ nn 0 1 =1,N... _. ___ 11 ( I)=X(I) nn 1? I=1,N TI=JPRM( I) 12 X(II)=B(I) R TIRN ____ RFAL FUINCTInN CARSI(Z) CnMPL X Z TF (RFAL(Z).FrO.n..AKhn.AIMA r,(Z).E O.O. ) GO TO 10 CARSl=CARS( Z Gn Tn?n 10 C A RS 1 =0. -........................... In rCAsl=n. 2n RFTtIRN FN) sn TA n.? 0.01 0.5?R3 200 3 7 TnP 763 LINES PRINTED

UNCLASSIFIED Security Classification DOCUMENT CONTROL DATA - R & D (Security classification of title, body of abstract tnUd ilndexi.nl iannotatlon must be entered when the overall report is claaalled) 1. ORIGINATING ACTIVITY (Corporate author) 2a. REPORT SECURITY CLASSIFICATION The University of Michigan Radiation Laboratory, Dept. of UNCLASSIFIED Electrical Engineering, 201 Catherine Street, 2b. GROUP Ann Arbor, Michigan 48108 —.... 3. REPORT TITLE SCATTERING BY A PERIODIC SURFACE 4. DESCRIPTIVE NOTES (Type of report and Inclusive dates) Scientific Interim 5. AU THOR(S) (First name, middle initial, last name) Tommy C. Tong Thomas B. A. Senior 6. REPORT DATE 7;I. TOTAL NO. OF PAGES 7b. NO. OF REFS September 1970 96 30 8a. CONTRACT OR GRANT NO. 9.. ORIGINATOR'S REPORT NUMBER(S) F19268-68-C-007 1 b. PROJECT NO., Task No. and Work Unit No. 1363-6-T 5635-02-01 Scientific Report No. 7 c. DoD Element: 61102F ()h. OTHER REPORT NO(S) (Any other numlbcrs thao may be assionedj I(is report) d. DoD Subelement 681305 AFCRL-70-0507 10. DISTRIBUTION STATEMENT This document has been approved for public release and sale; its distribution is unlimited. II. SUPPLEMENTARY NOTES —- | ----— 12. SPONSORING MILITARY ACTIVITY TECH, OTHER Air Force Cambridge Research Laboratories(CR]) Laurence G. Hanscom Field Bedford, Massachusetts, 01730 13. ABSTRACT" -- Numerical procedures are developed for the digital solution of the integral equations for the currents induced in a perfectly conducting, two dimensional periodic surface when a plane electromagnetic wave is incident. Data are obtained for both the surface and far fields for a variety of sinusoidal surfaces for plane waves of either polarization at oblique as well as normal incidence, and the results are compared with the predictions of physical optics. DD FR.M.1473U-NCASF *_ _; _ — UNCLASSIFIED S,.',riiv' (."i s.ilic,.-....

1363-6-T UNCLASSIFIED Security ClasHification " 4-. | LINK A LINK 0 LINK C KEY WORDO ROLE WT ROLE WT ROLE WT Periodic Surfaces Diffraction Gratings Surface Currents Diffraction Modes Numerical Techniques Exact Data Physical Optics UNC LkASSIFIED SzcrilrIy (Iassifi:al Io1n