13630-9-T = RL-2262 AFCRL-72-0162 013630-9-T The Numerical Solution of Low Frequency Scattering Problems By Thomas B. A. Senior and David J. Ahlgren The University of Michigan Radiation Laboratory 2455 Hayward Street Ann Arbor, Michigan 48105 February 1972 Scientific Report No. 12 Contract F19628-68-C-0071 Project No. 5635 Task No. 563502 Work Unit No. 56350201 for public release; distribution unlimited. Approved Contract Monitor: John K. Schindler Microwave Physics Laboratory Prepared For Air Force Cambridge Research Laboratories Air Force Systems Command Laurence G. Hanscom Field Bedford, Massachusetts 01730

013630-9-T ABSTRACT The low frequency scattering of electromagnetic and acoustic waves by rotationally symmetric bodies is considered. By concentrating on certain quantities such as the normalised component of the induced electric and magnetic dipole moments, it is shown how the first one or two terms in the far zone scattered fields can be expressed in terms of quantities which are functions only of the geometry of the body. Each of these is the weighted integral of an elementary potential function which can be found by solving an integral equation. A computer program has been written to solve the appropriate equations by the moment method, and for calculating the dipole moments, the electrostatic capacity, and a further quantity related to the capacity. The program is described and related data are presented. ii

013630-9-T TABLE OF CONTENTS ABSTRACT ii ACKNOWLEDGEMENTS iv I INTRODUCTION 1 II PERFECTLY CONDUCTING BODIES 3 2.1 Formulation 3 2.2 Prooedure for P 9 2.3 Procedure for P33 and C 11 2.4 Procedure for M11 18 2.5 Procedure for M33 21 2.6 Disjoint Surfaces 22 m ACOUSTICALLY SOFT OR HARD BODIES 33 3.1 General Procedure 33 3.2 Soft Bodies 35 3.3 Hard Bodies 40 IV THE COMPUTATIONAL TASK 52 4.1 Integral Equations 52 4.2 The Kernels and Their Singularities 55 4. 3 The Body and Its Volume 62 V NUMERICAL SOLUTION OF THE INTEGRAL EQUATIONS 69 5.1 P / V Computation 69 5.2 P33/V0 Computation 73 5.3 ll /V Computation 75 5.4 Sampling Rate 78 VI CONCLUDING REMARKS 81 REFERENCES 86 APPENDIX: THE COMPUTER PROGRAM 87 iii

013630-9-T 1. INTRODUCTION When a plane electromagnetic wave is incident on a finite perfectly conducting body, or a plane acoustic wave incident on a finite acoustically soft or hard body, the scattered field in the far zone can be expanded in a power series in the wave number k if k is sufficiently small. The determination of the first few terms in these series requires the solution of certain elementary potential problems. We here consider the potential problems associated with the first (Rayleigh) term in the electromagnetic expansion and the first two terms in each of the acoustic expansions, and show how in the case of a singly connected body of revolution all of these terms can be deduced from the solutions of just five potential problems. If the body is not singly connected, only the axial component of the induced electric dipole moment is affected, and for a body consisting of two separate parts, an expression for the modified component is obtained. Each potential satisfies a simple integral equation. Computer programs are described for solving the equations by the moment method, and since most of the equations are of first order type, the computational procedures are rather similar to those of Mautz and Harrington (1970). The appropriate elements of the electric and magnetic polarisability tensors (Keller et al, 1972) are then computed, along with the electrostatic capacity and a quantity 7 related to this, and these are sufficient to specify the electromagnetic and acoustic scattering for any direction of plane wave incidence and any direction of scattering. For relatively simple geometries, the entire computation takes about 3 seconds on an IBM 360 computer. In our presentation we first examine (Section 2) the problem of a plane electromagnetic wave of arbitrary polarisation and incidence direction, and

013630-9-T isolate the potentials necessary for a complete description of the leading term in the far zone scattered field. This is followed (Section 3) by similar treatments of the acoustic problems, but here we seek the first two terms in the expansions. In Section 4 the integral equations satisfied by the potentials are cast into forms appropriate to digital solution, and the manner in which the body is specified is also described. Section 5 is concerned with various aspects of the computer program, a complete listing of which is given in the Appendix, and some of the numerical results obtained so far are presented in Section 6. This work was supported in part by the National Science Foundation under Grant GP-25321. 2

013630-9-T 2. PERFECTLY CONDUCTING BODIES 2.1 FORMULATION Let B be a finite, closed, perfectly conducting body of revolution about the z axis of a rectangular Cartesian coordinate system (x, y, z). In terms of the cylindrical polar coordinates (p, 0, z) where x2 + y, 0 arctan, x the surface will be described by the equation p ap(z) where p can be a multivalued function of z as, for example, in the case of a disk or a re-entrant shape, but is never infinite and is zero outside some interval in z. Let r be the radius vector to an arbitrary point in the domain 22 exterior to B and let n be a unit vector normal to the surface drawn into A. A linearly polarised electromagnetic wave is incident with electric and magnetic vectors A ^ ikk.r E1 ae A i A ik k.r Hi Ybe where k, a and b are mutually perpendicular unit vectors such that b X kA a; Y is the intrinsic admittance of the homogeneous isotropic medium (of permittivity E) exterior to B and a time factor e has been suppressed. 3

013630-9-T For k small but kr large, the resulting scattered field E, H can be written as (Kleinman, 1965) ikr S e E.$ - - -4rr 4 7r r ikr S e H 4 - - 4r r k2 ( A(^1 r (A ) ^( rA^ ) + - (rAm) y ^ --- (2) where E and m are the electric and magnetic dipole moments respectively. As shown by Keller et al (1972), E. 11 Plla + (P33-P1 ) (a.)^z (3) mrn -Y M b + (M33 -M)(b.z)z where Pll P33' Mll and M33 are functions only of the geometry of the body. For a given body, P 1' P33 Mll and M33 are constants whose values are as follows: (i) Pll rB B x a(x ) dS (4)n (4) where 1s is such that V2 1 = 0 1 =x in V on B (5) 4

013630-9-T -2(r 1 O(r as r -4 oo ff (ii) p 3 33 z a ) dS an (z - 3 (6) B where $3 is such that 2 V2 ~3 + 0 3 = z+T in V on B; (7) y is a constant chosen to make dS = 0, B implying zero total induced charge on B, and ensuring that (8) -.2 3 = O(r ) B as r -oo (iii) n.x (x-,1) dS (9) where *1 is such that V2'l = in V on B a*l an ax -n On (10) 5

013630-9-T -2 2(r as r - oo M33 - f. (z-t3) dS B (iv) (11) where 3 3 is such that V2 43 = 0 a+3 a z an *4n 43 t O(r -2) ni3 = O(r ) in V on B (12) as r — oo. Although the values assumed by the potential function on B are quite distinct from thoseof a1, a 1 and 3 nevertheless, as shown by Karp an an ( 1956) and Payne (1956), 1 M = PI 33 2 11 (13) This obviates the need for solving the potential problem (iv) if the only purpose for finding 43 is to calculate M33. There is one other electromagnetic quantity of interest and this is the electrostatic capacity C of the body in isolation. If the body is raised to the potential unity, the surface charge density is a0 Ps An (14)an 6

013630-9-T where i0 is an exterior potential function satisfying the boundary condition 0= 1 onB. (15) The electrostatic capacity is then equal to the total charge induced on the surface and is C / -a dS (16) B Note, however, that if all portions of the surface are not in electrical contact with one another, charge can no longer flow freely over the entire surface, and additional (mutual) capacities can be defined. In particular, such electrical separation has a profound effect on the calculation of P33, and the modifications that result when the surface is disjoint are discussed in Section 2. The five quantities listed in eqs. ( 4 ), ( 6 ), (9 ), (11), and ( 16) can be computed by solving five separate potential problems of a rather standard nature, and the manner in which this is done is as follows. 2 Let V be some potential function satisfying V V = 0 outside and on B, and let V be the regular part of V. V is therefore an exterior potential and we can regard V-V Vs V (17) as an incident potential. Green's theorem applied to the function V in the region V then yields 1 () 1 a V(..)+ V(r) 47r? (2) n - -4 R 8n" V(- d) (18) 7

013630-9-T where R=r-r'1. If the boundary condition on the potential V is V(r) 0, r on B, (19) eq. (18) reduces to V(r) = vi(r) - 1 l 1 a(20) R'..... (20) B and since the integral exists for all r including points on B, we can allow r to lie on B and apply the boundary condition (19 ) to obtain i 1 1 a V (r) ~ 4T || R a 'V(r') dS', (21) B av which is an integral equation of the first kind for aV If, on the other hand, the boundary condition on the potential V is V(r) 0 r on B, (22) an 0 r eq. (18)reduces to V(r)V (r) 4 V(r ) dS' (23) B and because of the non-integrable singularity of a-r () at r = r eq. (23) is valid as it stands only if r is not on B. To obtain an integral equation for V, we allow r to approach a point on B in the direction of the 8

013630-9-T inward normal, in which case it can be shown that r V(r) a dS 27TV(r) + V(r') -{L dS B B where the bar across the integral signs denotes the Cauchy principal value. Hence V(r) a 2Vi(r)+ --- V(r t) a- ( ) dS (24) B for r on B, which is an integral equation of the second kind for V. 2.2 PROCEDURE FOR P 1 The solution of problem (i) and, hence, the computation of P is a straightforward application of the integral equation (21). If Vi x, VB = - then V x-1 = V1 (say) (25) with V 1x 0 on B. Since x X p cos 0 and the equation of the surface is independent of 0, the potential V must everywhere have the same 0 dependence as V, implying V (r) V1(p, z) cos 0. (26) avi This is true also of -, and we can therefore write an 9

013630-9-T a a a-, V ((r) -, V (p, z )cos0' an 1 - an, 1 T (s')cosO' where s' is arc length along a profile of the body. Moreover, 2 2 1/2 R (p-p) 2+ (z-z)2 + 2pp'(1 -cos ) ) with (27) (28) ' = 0 - 0' and since dS' a p' d ds' the integral equation (21) now takes the form 2 z s pcos0.... T1(s') p'd'ds' 1 cos 0) 4 j I T1(s ) R dds' 5 " 1 p' K T (s)ds cos 0 2" 1 1 i where the kernel is K 1 K (p, z; p 'I, z') 1 I cos dL R (29) 10

013630-9-T and the integration is along the profile of the body. The integral equation for T (s) is therefore 18 p K1Tl(s) ds 27rp 0 which can be solved to determine T (s). In terms of this quantity -x ' V(r) dS B (30) 2ir s 4 10 2 2 p cos 0T (s)dods 1 which reduces to Pll ' rI 2 p T (s)ds. 1 (31) 2.3 PROCEDURE FOR P33 AND C The solution of problem (ii), leading to the calculation of P33, involves two successive applications of the integral equation ( 21 ). In the first case we consider an incident potential i i V V2 (32) and seek the corresponding total potential V2 satisfying the boundary condition 11

013630-9-T i V = 0 on B. Since V is everywhere independent of 0, it follows that v 2 av 2 and are likewise independent. We can therefore write an e T 2 (s') and the integral equation (21) now becomes (33) 1 47 2r s 0 T2(s ) R p td0ds 2 R s 1 2r o p'K T (s')ds' 0 2 where the kernel is 7 K0 5 KO(PZ;p, Z') ~f ~0 R (34) The integral equation from which to determine T ( s) is therefore 2 s i p' K T2(s')ds' * 2rz (35) 12

013630-9-T The second of the two basic problems is that in which the incident potential is i i V=V = 1 3 (36) We again seek the total potential V3 satisfying the boundary condition V 3 0 on B, and writing av an' T81 (37) the integral equation (21) takes the form P 'K0T3(s) ds' 2i w0 (38) from which T ( s ) can be found. In problem (ii), however, Vi=z+ y and V X - 3 Thus, i i VUrV2 implying i + V3 + rV (39) V =V (40) 13

013630-9-T and if we write a n' T(s') a n'. (41) then T(s) = T2( s ) + T3(s) (42) where T (s) and T (s) are the solutions of the integral equations (35) and ( 38 ) respectively. The constant 7 is determined by the condition ( 8 ) for zero total induced charge on B, viz. * J n (v -v) dS X 0. (43) B n But JJi S an (z+ y)dS O0 aB B B B (44) as can be seen by application of the divergence theorem; moreover av dS = t2r 2 B s0 T(s) pd$ dS pT(s)ds * We are here assuming that the surface is not disjoint. 14

013630-9-T S s 27rJ pT2(s)ds+ 27r 7 p T3(s)ds 0 0 and hence, by virtue of eqs. (43) and (44), p T2(s) ds 0 -. (45) pT3(s)ds Since T ( s) and T3 ( s) can be found from the integral equations ( 35) and (38), the constant y given in eq. (45) now completes the specification of the surface field T (s ), and in terms of T ( s) P33 a J -n V(r) dS B 27r zpT(s) ds Hence s s P33 27r zp T2(s)ds + 2 7 Y zP T3(s)ds (46) A valuable by-product of the above analysis is the electrostatic capacity 15

013630-9-T C defined in eq. ( 16). This fact is apparent on recalling that the determination of C requires us to find the exterior potential 0 satisfying the boundary condition ( 15 ) on B, and this can be accomplished using the integral equation (21) with vi = 1 5 so that V X 1 -0 ( 0 on B). The problem is therefore identical to the second of the two basic ones considered above, and indeed j0 1-V3 implying a 0 av3 _ T (s) an an 3 where T3 ( s) is the solution of the integral equation (38). Hence, when the body is at unit potential, the surface charge density as a function of arc length is (47) (48) Ps = ET3(s) 5 3 and the capacity C is s C a 2x~ l 0 (49) pT3(s)ds. (50) We observe that the denominator of the expression (45) for 7y is simply 16

013630-9-T C /(2 7r E), which ensures that ' can never be infinite. S Some simplification of the preceding results is possible. Since V and S2 V S are both exterior potentials, Green's theorem can be applied to the region 3 22 exterior to B to yield the reciprocity relation 5 j s av3 fV 3 dS = V2 8n B (VanBladel, 1968). But avs a 2 A A an T2 (s) - n. z an 2 V s B 2 2 dS an (51) and s aV3 an - T (s), an 3 and from the boundary conditions on V and V 2 3 V2 -z 2 s V -1 3 on B. Substituting these into eq. ( 51), we have s zT3(s)dS = 2(s)-n. z dS B B which reduces to' pT2(s)ds a zpT3(s)ds (52) 17

013630-9-T With the aid of this result, the constant y of eq. (45) can be expressed in terms of the surface field T (s) alone: 3 S zpT3(s)d s O =' - (53) i.e. zpT3(s)ds 54) but whilst this reduces from four to three the number of separate surface field integrations involved in the calculation of P33 there is no way to avoid entirely the determination of the surface field T2(s). Indeed, the simplest expression for P33 is (55) P33 2r zpT2(s)ds - y (55) 0 2.4 PROCEDURE FOR M11 The solution of problem (iii) leads to the calculation of Ml and is a straightforward application of the integral equation (24). If V x and 18

013630-9-T then V ax- t1 V4 4 (say) (56) and aV4/av = 0 on B as a consequence of the boundary condition (10). 4/an Since x = p cos 0 and the equation of the surface is independent of 0, V 4: must have the same 0 dependence as V, namely, cos 0. In particular, on the surface V4(r) = V4(s)cos0 which enables us to write eq. (24) as (57) 1 V (s)cos 0 = 2p cos 0 + 2 -4 27r J V(s')cos '- a ()p'td0ds anlR But y sn 0') - z sn n cos a (x cos + ysin') - z sin a' (58) where ' = tan 1 a az' (59) so that a (i1 an' \RJ V1 1 n x - os a' (p ' -p cos ) - sin a (z' - z) 3 08 R (60) 19

013630-9-T with L ^ 0 - 0' as before. Hence 2 r 2 r os an i d x cos 0 cos ~ n ' i d 2cos 0-( p cos at Q + [(z' -z)sina'-p cosa' I } where 1 % (p, z;p', z') 3f d, (61) R 2 (62) Co s 2 =2(pz;p',z') d3 ' The integral equation from which to determine V (s) is now V4 (s) 2p+1- V4(s') (pcoofat2+ tz'-z) sina'-p'cos'aQ P tp'ds' 0 (63) and in terms of V4( s): M = ff ^ V4(r) dS B 27T s f cosaV4 (s) cos 0 p d0ds S0 J0 20

013630-9-T i.e. M 7r pV(s)cosads (64) 2.5 PROCEDURE FOR M33 Although it is not necessary to compute M33 directly because of the relation (13), the integral equation which the corresponding potential satisfies must be solved if the second term in the low frequency expansion for an acoustically hard body is to be evaluated. It is therefore appropriate to describe the determination of this potential function here. Once again we have a straightforward application of the integral equation (24). If V X z and V -1 3 then V z-3 a V5 (say) (65) with V5/ n 0 on B. Since V must be independent of e, eq. (24) 5 implies 27 S \ V5(s) 2z+ 2 f 5(') n ( a )P d0ds But 2 2 r an' n 0 0 foV / 21

013630-9-T x 2 cos a1 Q+ z1 -z)sina -pr cosa] Q ) where 1 is as defined in eq. (61) and 7T 20 a QO(Pz;p, Z). (66) The integral equation from which to determine V5 () is therefore 8 r V5(s) 2z+ +- V5 (s) co aO ' l+ (zt z))sinat -p'Cosa 1 20 p'ds' (67) and we note in passing that M33 J n z V (s)dS B 27T 8 a - J | sinaV (8)pd0ds i.e. M33 2 | pV (s)sinads. (68) 0 2.6 DISJOINT SURFACES So far it has been assumed that all portions of the surface are in electrical contact with one another, and if this requirement is not met, the analysis is no 22

013630-9-T longer valid. Thus, for example, an application of the above procedures to a body consisting of two separate spheres leads instead to the solution for the two spheres joined by an infinitesimal wire along the axis of symmetry, and though the presence of the wire (producing electrical contact) does not affect the values of M1l and P11 (and hence M33, by virtue of eq. 13), it does have a profound effect on P33. This is not unnatural since P33 is proportional to the longitudinal ( z) component of the induced electric dipole moment. The breakdown in our formulation when B has several distinct parts stems from the imposition of the zero induced charge criterion ( 8). If charge cannot flow freely between the n parts B1, B2...., B, eq. (8) must be replaced by the n equations JJ a dS ~ i= 1, 2...., n. (69) B i Since this obviously affects only the potential ( 3and leaves the procedure (and results) for P11, M11 and M33 unchanged, our efforts will be directed at P33 alone with the objective of finding an approach which is applicable when B consists of just two electrically isolated portions B1 and B2. So that we may use to the fullest extent the work that we have already done, it is desirable to have this new approach as similar as possible to that appropriate when the two portions are electrically connected. By analogy with problem (ii) of Section 2. 1, the task is to find an exterior potential 3 satisfying the equation V 3 O0 in the domain V exterior to B, together with the boundary conditions (3 = z + 71 on B (70) 3 z + 72 on B2 (71) 23

013630-9-T where the constants 71 and y2 are such that B1 If B2 - 3 dS X 0 an (72) a83 an dS 0 (73) The quantity P33 is then given by eq. ( 6 ) as before. Because the boundary conditions on B1 and B2 differ, it is no longer convenient to think in terms of incident and total potentials, with the difference representing the desired exterior potential. Let us therefore consider the basic potential problem in which ( ) is an exterior potential satisfying the boundary condition 1 on B1 - ( X(74) ~0 on B2 By application of Green's theorem to the domain ),, we have "(1) ) = 4 (:3 4 4' (r ')a R (1 1)(r t S an' -R 'ant ~3* - -4-I an dS --- ar n' )(r) dS' 4 n 47 R anB3 B1 B - 24

013630-9-T and the first integral is identically zero since B1 is itself a closed surface. If, now, r is allowed to approach B, application of the boundary condition (74) gives.f r on B 4 R an' 3 (,) dS' B B 2 -- which is an integral equation from which to determine / a. It can be simplified somewhat by observing that 1 3 and, hence, 3! are a/ n are independent of the azimuthal coordinate 0. When the 0 integration is performed, eq. (75) reduces to s 2 7, r on B p'KT3 (s )ds' (x 1 (76) 0, rr on B c.f. eq. (38), where T (s1) = a ( )(r') (77) and K0 is the kernel defined in eq. (34). It will be noted that the integration in (76) is over the entire profile of the body B - B + B2. Similarly, if 2 3 is an exterior potential function satisfying the boundary condition ro on B (78) 1 on B2 25

013630-9-T then PKT(2) )ds 1 (79) 2r, r on B2 where T(2)s) - n 3 (r' ). (80) Comparison of eqs. (76) and (79) with (38) shows that T ( s) + T3 (s) T (s) (81 3 3 3(81) where T3( s) is that surface field quantity which is appropriate when B and B2 are electrically connected. If T (s) has already been computed, it is 2 3 (1) (2) clearly necessary to compute only one of T (s) and T ()s). 3 3 Let us now return to the potential problem set forth in eqs. (70) through (73). As regards the boundary conditions (70) and (71), an exterior potential satisfying them is (1) 3 z -V2+ 3 + '1 3 (+ (82) where V2 is the total potential considered in Section 2.3. Hence 3 (1) (2) 's- '-TY ()s)- 1y2T (s) (83) a n (an) 23 ( ) and since T (s ) is given as the solution of the integral equation (35), it only 2 remains to specify the constants Y1 and 72' From the zero charge condition (72) and using the fact that 26

013630-9-T B1 or B2 a- dS = 0 an we have 71 pT3 1(s)ds + 2 PT3(2) (s)ds - pT (s) ds (1) (1) 3 2( (1) (1) (1) (84) where the symbol (1) below the integral signs shows that the integrations are carried out over the profile of the portion B1 alone. Similarly, from eq. (73), 1 I PT (1) (s)ds + y (2) (2) pT3(2(s)ds - - 3 (2 (2), pT2(s)ds (85) where the integrations are over the profile of B2 alone, and if we now define 2 C1 x 27r 11 (1) (1) (1) T13 (s)ds 3 C12 = 2 r I PT3(2)(s)ds (1) (86) p T3(2)(s)ds (2) 9 C2 = 2 r~ 21(84) a eqs. (84) and pT3(1)(s)ds, C22 = 27r 22 (2) (85) take on the more compact form YlC11 + 72C12 = -2re I pT2(s)ds (1) (87 a) 27

013630-9-T 1 C21 + Y2 C22 -2 r pT2(s) ds (2) (87 b) It will be observed that the quantities Cll, etc. all have the dimensions of capacity, and by virtue of eqs. ( 50) and (81), Cll + C12 + C21 + C22 = C (88) where C is the capacity when electrical contact is maintained. Rather than solve the eqs. (87) directly, it is more convenient to first eliminate the surface field quantity T2 (s) from the expressions. That this is possible can be shown by application of reciprocity to the exterior potential fuctions z - V, 3() and *. From the pair z - V2 and 3 we have 1) (-V2 )dS J (z 2) n dS B B Hence /f B1 T2(s)dS x IzT3(1)(s)dS B implying f (1) s pT2(s)ds J 8 pT2(s) ds = O 0 zpT (1)(s) ds (2) zpT (2(s)ds 3 (89) (90) Similarly, I (2) and we note that by addition of the last two equations we recover eq. (52). 28

013630-9-T Finally, from the function pair () and 32) 3 3 B1 T3(2)(s)dS (2 3 ffT3()(s)dS B2 implying (1) pT3 (s)ds d 012 = (2) pT3(1)(s) ds (91) i.e. C21 (92) as expected. Using eqs. (89 ) and (90), T ( s ) can be eliminated from the eqs. ( 87 ) and if we also eliminate T ) ( s ) using eq. (81), we obtain 3 ~S (1 -2) 0c11 +22E J pT3(s)ds - -2vre (1) rs (7Y1 72)C21+ '22 Tr J PT3(s)ds = 2rEf (2) 0 zp T3( )(s) ds zpT3(1)(s)ds s -2 r c zpT3 (s) ds 29

013630-9-T These can be solved to give s s s (1) 0 (93) (1) 0 s where s Ts3 A Js pT3(sd(s)ds pT (1) (J- pT3(s)ds (1)(s~dp S (s)ds (1) 0 (1) 0 (1) (95) We can now proceed to the calculation of P. If we write this quantity tv 33 as P33 to distinguish it from the 7P of eq. (55) for B1 and B in electrical contact, we have, from eqs. (8) and (83), P'' || T (s)+ Y1 T ()+ )(s)} dS (96) 30

013630-9-T But P33 2r f zP T2(s) + 7T3(s)) ds 0 (see eq. 46) where Y is given by eq. (53), and thus S s P33 P33 2r (y2 -) zPT3(s)ds+2ir (7Y1- 2) zpT3 (s)ds. 0 0 (97) Moreover, from eqs. (53), (94) and (95), after some manipulation, y2 =-(1 72) pT ()(s)ds 3 - i.n (98) pT3(s)ds which enables us to write eq. (97) as 33 33 P33 a P33 + 2 r(y1 -IY2) (z +7) pT ()(s) ds 3. (99) The factor (y -1 y2) is defined in eq. (93) and invoking yet again the expression (53) for 7 together with the identity (91), we have 71 -2 A" pT3(s)ds f pT3(s) ds + zp T3 (s)ds (1) 31

013630-9-T C 3 Z f ( >T (s)ds 0 giving,~-~ C p33 p3 e ( z+ y)pT3 (s)ds (100) where C and A are defined in eqs. (50) and (95) respectively. This is our final expression for P33. Compared to the situation when B and B are electrically connected, the only additional field quantity that (1) must now be found is T (s), which is given as the solution of the integral 3 equation (76); and since C/E and A are both positive, electrical separation decreases the longitudinal component of the induced electric dipole moment. 32

013630-9-T 3. ACOUSTICALLY SOFT OR HARD BODIES 3.1 GENERAL PROCEDURE Let B now be a finite, closed acoustically soft or hard body of revolution about the z axis of a Cartesian coordinate system (x, y, z). It is of no concern whether B is disjoint or not. A plane acoustic wave is incident and its velocity potential is written as ui ikk. r e (101) A where k is again a unit vector in the direction of propagation. If U is the scattered field that is produced, then U satisfies (V2 + k2) Us 0 ( a - ik U) - r ar in V (102) 0 as r —)oo (103) and the boundary condition Us -ui on B (104) if B is soft, or au8 Oau an An on B (105) if B is hard. Eqs. (104) and ( 105) are equivalent to * To avoid any possible confusion, we shall henceforth refer to U as a field. 33

013630-9-T U 0 on B (106) a U 0 on B (107) an respectively, where U U U + U is the total field. A general expression for U( ) at an arbitrary point in 2) is provided by the Helmholtz representation: i 1 ikB6 i1kR, a U(r) ) U (r)+ (r )-, R e B 4eRR n'r (108) where R = Ir-rl as before. For sufficiently small k, Ui, U and, hence, U can be expanded as power series in ik of the form 00_ Ui(r) = (ik)m U (r) (109) m=0 and when these are inserted into eq. (108), the coefficients of like powers of i k on both sides of the equation can be set equal to give - R al t U(rt())dS' (110) for m = 0, 1, 2,... By allowing r to lie on B, an integral equation is obtained from which U (r) can be found; and as is seen by substituting m 34

013630-9-T the power series for U (r) into eq. (102), 2 2 0 1 (111) V U m = U m-2 I m > 2, m showing that U (r) and U1(r) are potential functions, but U (r) is not unless U0(r) = 0. In the far zone (r -+ oa) the low frequency expansion of the scattered field deduced from eq. (108) is ikr U (r) e4U rr oD m (_-)m l+1 mO0 1=0 (ik)m"f)! r. rt)m(Bt) B X (n'. r) U -l(r) + U (r) dS (112) (Kleinman, 1965, with the correction of a sign error), provided U 1( r') is taken to be zero. Our objective is to calculate the first few terms in this series. 3.2 SOFT BODIES We now specialise the above results to the case of a soft rotationally symmetric body illuminated by the plane wave ( 101), and seek the first two terms in the low frequency expansion of the far zone scattered field. By 35

013630-9-T invoking the boundary condition U (r) = 0 m -- m = 0, 1, 2,.... we have on B, (113) ikr U (r) e (-1) 4 r r; U (r' ) dS' - ik fr. r BB n U(-r') B -an' U(r') dS'+ O(k2), 8ni 1 -J0 (114) showing that only the potential functions U0 ( r) and U1 ( r) are required. From eqs. ( 101 ) and ( 107 ) it follows that U0 (r) x 1 k. Ul(r) xk.r (115) and by inserting the boundary condition (113) into (110), the latter becomes U (r) x U (r) - 1 m - m- 4r m y1 (m-l)! B m-l-1 R B an' U(rt')dS' ant which, for r on B, reduces to U (r) - m - 47r m two 1 (re'l)' ffR On' U1(r') dS' B (116) 36

013630-9-T When m 0, eq. (116) gives r41r R B a an U (rf) d8t (117) This is identical to the integral equation satisfied by the potential V (r) of Section 2.3, and hence U0(r) = V3(r) a U0(r) = T3(s) an - (118) We note that IJfa U0(r)dS B s 27r PT3(s)ds = Jo C E (119) (see eq. 50), where C is the electrostatic capacity. From eq. (116) with m =1, A 1 k. r - - 4vR Uo(r )dS' 4; n, U (r) dS' B B (120) and using eq. ( 119), the left hand side can be written as A A) A A C AO C (k. x)pcos (kysin + (k )psin. z)z - C Since the surface of the body is independent of 0, it follows that U (r) must have the form have the form~ 37

013630-9-T U1(r) ( k. x)o 0+ (k. y) sin U] where the individual U1 (r), j 1, (1)(r) L ^ +(^. (2) - c u (3).z)U1 4r)'4c 1r) 1, 47rc I (121) 2, 3, satisfy 1p P 4r 4R B Zx 4& R B a (1) d an'U (r)dS (122) a U(2) S ' U1 (r)dS' (123) 1 1 4 r f-1 RI a an U()(r') dS'. 1 (124) B Comparison of eqs. (124) and (117) shows U1 (r) U(r) = V3(r) U1 _ _ 3 a u(3) an 1 9 = T (r) 3 (125) Similarly, U(1 (r) cos p is identical to the potential V (r ) of Section 2.2, implying Ul(1)(r) - V (p, z) 1 1 (126) so that a (1)(r) = T(s) an 1 1 9 (127) 38

013630-9-T (2) and U( (r) is identical to the potential V (r) of Section 2.3, so that 1 2 an U1 (r) T2(s). (128) It is now a trivial matter to evaluate the right hand side of eq. (114). The first integral is clearly C/e, and the second can be written as j r. r') T3(s )-.cos +(k.y)sin3 Tl(s') B A A C + (k.z) T2 (s) ' Ts3( dS' 27 r 1 Jo J 3 pds' But pIT(st)dst z'p'T3(s) ds' (52) 0 ^0 and hence the second integral on the right hand side of (114) is 2ir T 4 (r - y k) z p'T 3( s ds' 4e - ( r - Jk) where y is as defined in eqs. (53) and (54). 39

013630-9-T The low frequency expansion of the far field is therefore 4rr (4r e- - (r - ). + O(k ), (129) showing that a knowledge of C and y alone is sufficient to specify the first two terms. As demonstrated by Van Bladel (1968), a similar result obtains even for a body which is not rotationally symmetric. 3.3 HARD BODIES The final case to be considered is that in which B is a hard rotationally symmetric body. The boundary condition on U (r), m = 0, 1, 2,.... m is then U (r) 0 on B (130) an mand when this is inserted into eq. (112), the low frequency expansion of the far zone scattered field becomes ikr e Ie 3 n 2A U (r) (1i) ik. r)U0(r')dS'+k2r. r U (r) B B -U1 (n r) dS ik r r)2 U0(r') B -(r. r U(r ) + U (r ) (n'.) dS' + O(k4). 1 2 -1 (131) 40

013630-9-T As we shall see later, the first term O(k) is identically zero, and we therefore need U0 ( r), U1 (r) and U2 (r) to compute two non-zero terms in the expansion. From eq. (110) and the boundary condition (130), an expression for U (r) at an arbitrary point r in V is m 4.U ( r)Ui (r) + -m+l) R Ur dS t=0 B (132) and in particular, when m 0, U0(r) U0 (r) + Uj (r ) an- dS'. (133) B s i Clearly U (r) a U (r) - U0 (r) is an exterior potential function a a i and -- U (r) 0 on B since U (r) = 0. In addition, an n 0 5 -1 U0( r ) vanishes more rapidly than r as r -4 oo since there is no term 0 O( k ) present in the expansion (131), and hence U (r) - 0 (134) implying U0(r) = U(r) 1. (135) From eq. (132) with m - 1, we have u(r) u(r) + ffU (r) an' RdS (136) B 41

013630-9-T which can be converted into an integral equation for U1(r) by allowing r to approach B. Because of the non-integrable singularity of the kernel for r on B, it is necessary to apply a limiting process, and if a bar across an integral sign is again used to denote the Cauchy principal value, we obtain U (r) 2 U (r) + 2 f U(r) an (), dS' (137) B i for r on B, where U (r) is givenbyeq. (115). Interms of the 1 cylindrical polar coordinates (p, f, z), i A (A U(r) {(kx) co + (k.y)sin + (k.)z and since the surface of the body is independent of O, it follows that U (r) can be split up into three parts each of which has the 0 dependence of that i part of U (r) giving rise to it. In particular, on the surface, 1 -Ul(r) (= (.)cos + + (. )sin V4(s) + (.z.)V5(s) (138) where V4(s) and V ( s ) are the potentials introduced in Sections 2.4 and 2.5 respectively and satisfying the integral equations (63) and (67). For the remaining function U2(r) an expression at an arbitrary point r in V is given by eq. (132) with m = 2 and is U r) 2 (r) + 4 f U2(r.E) n' )dS (139) B 42

013630-9-T where (see eqs. 101 and 109) U (r) 1 i (. r)2 (140) 2 2 An integral equation for U2 (r) can be obtained by allowing r to approach B, but it proves unnecessary to determine U (r) explicitly if the only pur32 pose is to calculate the term O(k ) in eq. (131). To see this we first note that since U (r) - 0, the eqs. (111) imply 0 V2 U 2 showing that U2 is a potential function. Moreover, from eq. (139), 2. U 5 U - U2 is an exterior potential, being of double-layer type, and 2 2 2 since a A fAA U2 (r) (k. r) (k. n), (141) a n 2 the boundary condition on U (r) is 2 - a s A AA a U (r) - (. r) (k.n) 142) a n 2 for r on B. U2 clearly depends on the direction of incidence as well as 2 that of the normal to the surface, and in principle nine separate but elementary potential problems must be solved to find U. In terms of these potentials, 2 3 3 U2 (r) k G(r) (143) i=1 jul where, for convenience, we have put 43

013630-9-T x=x1, ywx2, zx3 X and the potential functions Gi () are such that a AA A 8 G (r) = -(n xi)(r. xj) 8n ij ' - for r on B. In like manner we can write 3 U (r) ' k.Fi(r) io1 where the functions F (r) are such that i i, j 1, 2, 3 (144) (145) aA A an (r) -n.xi (146) on B, and comparison of (145) with (115) and (138) shows that on the surface F (r) (V4(s) -p os0 F2(r) (V4(s)-p) sin (147) F3(r) - V5(s)-z. Following Van Bladel (1968) we now apply reciprocity to the exterior potentials F, (r) and Gi (r) i,, j= 1, 2, 3, in the region V to get Fr) an ijSB BF')i B B 44

013630-9-T which reduces to I Gj (r')( '. )-dS' = I( r n ')(.xi)(r'..)dS' B B when the boundary conditions (144) and (146) are employed. Hence (148) ff 2(r') (n'. xI)dS'= B I J - I J kikj f Gi (r)n. x)dS B kij J (r)(n ' )( r'.x )dS' B B FL()(k. n )(k. r ') dS' B implying Jf U(r)(' r )dS' B - U (r')('.r)dS' B I F r )( r x,)(k. n )(. ) dS' B and Bf I B F(rt)(r. x) (k. n, ) (. r_') dS' ff (r^')2(,. )dS B (149) 1 2 45

013630-9-T This integral is the only form in which U2 ( r ) enters the far field expansion through terms 0 (k ), and since the F (r') are known by virtue of the eqs. ( 147 ), the integral can be computed without the explicit determination of U 2 ( r ) itself. We are now in a position to evaluate the individual terms shown in eq. ( 131). Since UO (r) * 1, we have J (n'. r) U0(r') dS' = 0, (150) B 2 verifying that the leading term in the far field expansion is 0 (k ), and (r. r) U( )(n. r) dS' r. ^n(r r') dS' B B r(. | V(.r')d T' V0 V0 (151) where V0 is the volume of the body. Also, from eq. (138), U1 (r') (n'. r^)dS' r. J {cos a' cos 0' x + cos a' sin 0y-sina'm z B1 X (k. x)cos0'+k.ysin V4(s') + (k. z)V( s' )p'd0'ds' r. {( k. x) + y(k. y ) P' V4 (s) cosads 46

013630-9-T S -(r.z)(k.z ) 2 p'V5(s') sin a ds' 0 where a' is the angle defined in eq. (59). Hence, from eqs. (64) and (68), J U ( (n. ) dS' k M1-(k. z)(M -M 33) B (152) where Ml and M33 are the elements of the magnetic polarisability tensor 11 33 discussed in Sections 2.4 and 2.5 respectively. As we have previously noted, for a body of revolution M33 is related to Pll (see eq. 13). When the results of eqs. ( 150) through ( 152) are substituted into eq. (131), the low frequency expansion of the far zone scattered field is found to be ikr ( 2 A: A. U8() 4 Ir k [ M ^)(M) -(M M V + (k3)) (153) 3 where the actual term involving k is ik3 [1 (r r U 0 r r) U1(r')+ 2( (r dS. B (154) Unfortunately, the evaluation of this is rather a messy task. Since U (r) = 1 (see eq. 135), Jr2 r.' U()(n. r)dS. n'( r. dS' B B 47

013630-9-T l A 2'r. V 0 V r A Qr. dTI r') }1T ( 155 ) f ff V0 To simplify the treatment of the next two integrals in ( 154 ), write so that (see eq. 145 ) ( 156 ) U (r) I k A F I (r) ( 157 ) Using eq. ( 149 ) we then have f A. 1 (r' B a A )+ u2. IB A A ( n,r) dS' AA A A AA tA~ (n'.k) - (k.x1)(r. r)(n,r)J d St ( 158) B l A A IA -1k. r )(n' r) -( 2. But ff (~. r t)(f n )dSt (k.r) kfffr? B V0 48

013630-9-T and | (r.r')(n'.k)(k.r)dS' r +(k. r)k r J 'dT' B V 0 as may be shown by analyses similar to that performed above. Hence (k. r')(n. r) - (r.')(n k)'. k (k. r')dS' = -r. r d T' B V 0 (159) which cancels the contribution (155) of the first term in the integrand of (154). The complete integral (154) is therefore 3 A A A ik3 Fi (r)( r. x)(k. r')(.) ) (k.x)(r.r')(n'. r dS' (160)! B and to simplify this we now invoke the rotational symmetry of the body. From eqs. (147) and (156) we have F (r') x V (st)cos, F2(r') V(s')sin0', F (r ) V (s') (161) When these are inserted into (160) and the azimuthal integration performed, the contribution of the first term in the integrand is r (. 3. r - (k.z)(r.z) p (z cos a - p sina) V4(s) ds {^ ')f4 49

013630-9-T A A ^ A^ 2 + 2 (r.z) 1-(k.z)2 - 2, (r. z)(k. z) I p V5(s) cos ads p z V5(s) sin ads, The contribution of the second term in the integrand of (160) differs only in A% A having r and k interchanged, and when the two are subtracted, the final 3 expression for the term in k in the far field expression ( 153 ) is ik T (k - r). z. r - (k. z) (r. z) P(zcosa-psina) V4(s)ds 5 -f +(i.z)(,.Z r w0 p V (s) cos a ds A A A(k. z -2 (k. z)(r.z) s I pz V5(s) sina ds 0 (162) Although this is only the second non-zero term in the low frequency expansion, it is much more complicated than the second term in the expansion for a soft body. The surface field quantities involved are the same as those associated with M and M3, but there is now no simple relationship analogous to (13) which enables us to dispense with V5 (s). If the direction of incidence or 5 50

0136309-T A A A A observation is parallel to the axis of symmetry, i. e. k = + z or r = + z, the integral containing V4(s) disappears, but there is no comparable situation where the integrals containing V5 ( s) are absent except for the special case of forward scatter, k, when e entire expression (162) vanishes. of forward scatter, r ak, when the entire expression (162 ) vanishes. 51

013630-9-T 4. THE COMPUTATIONAL TASK When this study was first undertaken the main objective was to develop an effective program for computing the quantities P l, P33 and Mll specifying the low frequency scattering behavior of perfectly conducting rotationally symmetric bodies. The realisation that the calculation of P33 produces as a by product the electrostatic capacity led us to add this to the list of quantities considered, but it was only later that the question of acoustic scattering came up. Since the first two terms in the low frequency expansion for a soft body are expressible in terms of C / e and y, and 7 is implicit in the P3 33 computation, it was only natural to add this to our list, and for a hard body the first term involves no additional work. But the second term, (162 ), is another matter. In particular, it requires the explicit calculation of the surface field V5 ( s) that had hitherto been avoided by virtue of the relation (13), and even if this were done, the nature of the k term is almost such as to preclude any physical understanding of the data. For these reasons it was decided not to implement the computation of V5 ( s) and, hence, to ignore the second term (162) in the hard body expansion. The quantities which we are now left with are all ones which are needed for the electromagnetic problem. 4.1 INTEGRAL EQUATIONS It is convenient to begin by listing the integral equations which have to be solved and the quantities to be computed from their solutions. Assuming that the profile p = p ( z) of the finite, closed, rotationally symmetric body has been specified in some manner and its volume V0 computed as a preliminary step, then: 52

013630-9-T (i) solve I p' K1T1 (s')ds' = 2xp 0 (163) where the kernel K is defined in eq. (29); compute 1 p11 VO 5 7T Vol v 0I 2 p T (s)ds (164) (ii) solve s r. p'KOT2(s')ds' X 2wrz (165) p'KoT3 (s')ds' X 2ir (166) where the kernel K is defined in eq. (34); retain the option to print out T (s ); compute 3 C E = 2 s pT3 (s)ds, (167) (168) s 7 = - - 2 i zPT3(s)ds C 3 53

013630-9-T P P33 Vo 01 s fzpT2(s)ds - V V0 (169) (iii) if and only if B consists of two separate closed solve s p2K p' K0 T3( (s') ds' x 0 0 Iparts B and B2, r on B; (170) r on B 2 compute S 2 6 P33 V0 2r Vo (1) (z+ y) pT3 (s)ds 3 (171) 2 T(1) ds - 2 p T (s)ds 27r 3 c p T3(s) ds (1) (1) where the symbol (1) below the integral sign is carried out over the profile of B1 alone means that the integration (iv) solve s f V4(s')(pcos ' f2 + (z- z) sin a' -p cosa' 1 ' ds' x r v (s) - 2p (172) 54

013630-9-T where S 1, 22 and a ' are defined In eqs. (61), ( 62) and ( 59 ) respectively and the bar across the integral sign denotes the Cauchy principal value; compute pV(s)cosa ds. (173) 0 0 0 We therefore have four (five) integral equations to be solved, three (four) being of the first kind and one of ethe second, and five (six) derived quantities to be computed from their solutions: the numbers in parentheses refer to the unusual situation where B is disjoint. Before attempting this task, there are certain features of the equations to be examined. 4.2 THE KERNELS AND THEIR SINGULARITIES The kernels K0 and K1 of the integral equations (163), (165), (166), (170) can be expressed in terms of complete elliptic integrals of the first and second kinds. From the definition of R given in eq. (28), we have R = (p+P)2 + (z-z')2 (1 msin2 ) 2 (174) where 4p0' m = 'E (175) (p+p')2 + (z-zt)2 and 2 (Tr-. (176) 55

013630-9-T Hence 1 ( ( )1/2 2 -1/2177 (l-m sin e O) (177) and when this is substituted into the definition (34) for K, we immediately obtain (m 1/2 K ( 1P K(m) (178) where r /2 K(m) = (1-msin2 e1/2 de (179) 0 is the complete elliptic integral of the first kind (see, for example, Abramowitz and Stegun, 1964, p. 590). By a trivial manipulation, we also have cos 2 1 R - - 1 - -- (180) R m R 2pp' implying cos ~ -1/22 12 1/2. (mpp ) 1/2 (1-_m -msin2 -)l (181) and hence, from the definition (29) of K1, 1 (mpp )1/2 (( 2- ) K(m)- E(m) (182) wherer /2 E(m) J (1 -m sin2 e) 1/2 de (183) 56

013630-9-T is the complete elliptic integral of the second kind (loc. cit. ). The above representations of K0 and K are exact. Since p, z, p, z' are all real with p, p' > 0, it can be verified that 0 < m < 1. Over this range E (m) is a finite slowly-varying function, having the values 7r /2 for m = 0 and unity for m = 1. A finite polynomial approximation sufficient for -8 computing E (m) with an error of less than 2 x 10 is given in Section 17. 3. 36 of the above reference. Through the first three terms the precise expansion is ( Jahnke and Emde, 1945): E(m) 1- 4 ml + m (m, m2 i2 (184) 4 121 O ml with ml1 1-m, (185) 2 2 i.e. (p-p)2+(z )2 (186) (p +') + (z- z ') and In. (187) We observe that m1 a 0 if and only if p I' p, z ' z, that is, when the integration and observation points coincide. For an integration point in the immediate vicinity of the observation point, ( - 2 m1 Cl 2 ) (188) where s is to a first order the arc length between the points. The elliptic integral K(m) also has the value vr/2 for m = 0 but becomes logarithmically infinite as m -- 1. A finite polynomial approximation sufficient 57

013630-9-T to compute K( m) with an error of less than 2 x 10 is given in Section 17. 3. 34 of Abramowitz and Stegun (1964), and a precise expansion through the first three terms is (Jahnke and Emde, 1945 ): K(m) r+ m1 r - 1 + O(m, m2 ). (189) Because of the infinity of K(m) as m- 1 (m-+ 0), K0 and K1 are also infinite in this limit, but their behavior in the vicinity of the singularity is easy to determine. Using (184) and (189) we have 0 (p - + -o (r m ) (190) and K-1 1 2 + O (ml (191) K1 =' itp,)1/ showing that the singularity at p p z, z ' z is an integrable one in each case. The contributions of the singular (or self ) cells to the integrals in eqs. (163), (165), (166), and (170) are therefore finite and can be analytically approximated as follows. Consider for example the integral equation (165). If the self cell in the sampling procedure is centered on s = s (where p =p ) and is of arc n n length As, then s + 1 n 2 J p'K0T2(s')ds' p' K T (s)ds' self s - - s n 2 58

013630-9-T n T 2(s) B +- AB n 2 1 s -- As n 2 K ds' 0 - T2 (s) 2- n 1 T — 2 A - I A s 2 8p ds In Is'( ds' and hence 16 p n pt' K T2(s')ds' I T 2(s) (n +1 A self (192) It is desirable to retain the first correction, unity, to the logarithmic term to ensure the necessary accuracy when the sampling is relatively coarse and / or p is small. For the integral equations ( 166 ) and (170) the results differ from the above only in having T (s ) and T ( s ) respectively in place 3 n 3 n of T2 (s ); and for the integral equation (163): 2 n self p' K1 T1 (s')ds' - p T1 (s ) 1 1 ni n s +1- As n 2 n 2 K dst I 2l 8p n -9 As 592

013630-9-T giving 16Pn 'K1iT (s)ds' - T1 () n -1 ( -i AS self (193) For the integral equation (172) the computation of the kernel is a more complicated task due partly to the presence of the functions Q and Q. However, these also can be expressed in terms of complete elliptic integrals, and the resulting method of computation is much less time consuming than a direct numerical evaluation of the integral expressions for Q 1 and 2 The definition of 1 is given in eq. (61), and using eqs. (177 ) and (181), the integrand can be written as / 1 -3/3/2 COS = 1. m 2 (1 - -- )(1-m sin20). 3 4m pp1 2 R 2 -(1-m sin 0 -1/2 from which we have \3/2 r 3/2 Q -1 Tm1 -) ( 1 - m sin20) 1 2m pp 2 v / L9) -3/2 dO-K(m). (194) To evaluate the remaining integral, differentiate the expression (179) for K( m) with respect to m to get 1 K'(m) - 2 T/2 0 2 sin2 0 23/2 (1 -m sin 0) 60

01363 0-9-T 1 2m 2m ir /2 0 x /2 Jo 1 2 3/2 ((1.-msin20) _____ dO0 (1.-msin 20)1I) (i~ sn20 )- 3/2 d — 1 P M (1-.sin ) dO 2m K'm Hence 4w 2 (1 -Msin 2e). 3/2 dO e= K(m) +2mK'(m) ( 195 ) and when this is substituted into eq. ( 194), the result is ( 3/2 m - -m 1pp / 1 2 4 ( 196 ) The procedure for 02 is similar. From eqs. (177) and (181), 2 2 Cos '= R 3 1 2 2m ) 3/2 (-2 2 - 3/2 l~l-m)(1-rn sin 0) 2 " 1/2 2 / - (2 - m) (1 -msin 2O) + (1 -msinO20)) implying 1 2 = 2 m ) 3/2 2 4r/2 Jo -3/2 (1 -msin 20) dO0 -(2-MK~m) E(m)j. 61

013630-9-T and when the expression (195) for the integral is substituted into this, we have 3/2 2 2 2 2m(1- K'(m)- - 4 ) K(m)+E(m). m I(197) The finite polynomial approximations to E (m) and K (m) were mentioned earlier, and in particular, for the latter, 4 4 1 K(m)=(a + a m +....+ am1 ) + (b + b m +....+b m ) 0 11 4 1 0 1 1 4 1 m 5 5 1 + O(m, m1 n m ) (198) 1m where values for the coefficients ai and b., i s 0,...., 4 are given in Section 17.3.34 of Abramowitz and Stegun (1964). Since d/dm X -d/dm1, it follows that ^0 2 K (m)1 +(b- al )+(b2- 2a2)m +(b3-3a3)ml +(b -4a )m 2 3 1 4 4 1 -(b + 2b2m + 3b3m + 4b4ml ) In + O(mi m14n -) (199) which can be used to compute K' (m). We note the pole-like behavior of K'(m) when m 1 1 - m X 0, and this is reflected in the non-integrable singularity of the kernel of eq. (172) at p ' p, z' z. 4.3 THE BODY AND ITS VOLUME One of the many factors motivating the present study was the need to 62

01 3630-9-T compute the low frequency scattering behavior of missile-like targets. These are generally rotationally symmetric bodies (or can be approximated as such to an accuracy which is adequate at low frequencies), and are often made up of several distinct parts, e. g. a cone mated to a cylinder which is terminated in a spherical cap. Although the complete profile of such a body is certainly not an analytic curve, each individual segment has a relatively simple equation whose form can be used to advantage in the numerical process. It is therefore assumed that the profile is a finite piecewise smooth curve composed of straight line and circular arc segments. For definiteness, the number of segments is limited to 15 or less. At the end points of the profile where it intersects the z axis of rotation of the body, p = 0 (of course), and the nature of the program is such that segments which are perpendicular to the z axis can be handled, as can a 'disjoint' body having two separate parts provided each portion of the complete profile terminates on the axis. Every segment contributes to the total volume V which can be found by adding the individual contributions 6 V0. In certain cases, a volume contribution can be negative and subtract from the volume attributable to the other segments. Where this occurs, it must be noted as part of the input specification for the segment in question. In the following we list the input specifications of circular arc (Types 1 and 2) and linear (Type 3) segments, and give expressions for the corresponding volume contributions (assumed positive). The segments must be described sequentially starting at the intersection of the left hand segment with the axis, and the ordered sequence of segments defines the profile of the body. In some cases it may be desirable to regard a single linear or curved portion of the profile as two or more segments to permit a non-uniform 63

013630-9-T spacing of the sampling points over the whole. Tvpe 1 Segment (Circular Arc Concave Down) - 1 171 -- -- - Z-2 --- - - -lb --- - --- - -- - z - ---- - - - -.. - r Specification: z ' 2 1' 2 (ze, pr ) p = p(z ) P1 1 0 ( degrees ), volume sense P2 = (2) 0 < 0 < 180 XX\ \ x \ z( zi Sp ' z (z\,A If 0 is the angle subtended by the arc at the center of curvature, then the radius a is (1 (.p 2 + (2 = 1/2 a=... -P2 -Pi ) + (z2 - z1)2 2 sin -2 (200) Since we permit the specification of re-entrant circular arc segments we do not require z2 > z1. In order to obtain correct results for both standard ( z2 > zI ) and re-entrant segments define the quantity d -zl - z2 1 d X 1 ~2 - 1 (201) 64

013630-9-T Then, the coordinates (z0o p ) of the center curvature are O 2 1+ 2 -d (pl -p2) cot ) 2 ) (202) + P2 d (z -z ) cotj 2 The volume of rotation is given by z2 6Vo= T p (z)dz z1 9 and since the equation of the circular arc segment is )2 + ( - z )2 2 (p-p0Z) +(z- ) a the incremental volume 6 V is 0 0 3 6Vo= 7 (z2~-z1) (o+a-~' u22+UlU2+U)~ 2 2 0 +0 P ( p2-p0)-iul(pl-pO0) + d -a2 U2 z2 - z0 (203) where u1 ' Z1 -z0 65

013630-9-T Type 2 Segment (Circular Arc. Concave U ) - -- - - - -. - - jr- 9 Specification: same as for Type 1 (Zo, Po) )z \\ \ (z2' P2 Eq. (200) ( ZO p ) of the gives the radius a of the type 2 segment, but the coordinates center of curvature are now {Z +Z2 +d(pP- P2 ) cot+) z0 2 2 PO 2 1 2 2 - 1 (204) The incremental volume of the type 2 segment is 6V (z2-zl)P + a - (u22+ulu +u2) - U 2(P2 -Po) - ul(pl-pO)+ da20 (205) Note that only simple sign changes distinguish (202) from (204) and (203) from (205). Relationships that hold for both type 1 and type 2 segments may be derived by using a constant g defined as follows: 66

013630-9-T -1 type 1 segment (206) 1 type 2 segment The center of curvature (z0, p ) and the incremental volume 6 V0 for circular arc segments of types 1 and 2 are then z0 ' i Zl+z2+ dd(pl-2)cot2j 207) P0 1 +p2 + 2 d(z2- l)cot2 6VO = 7 (z2-l) p+a - (u1 + +uu2+ u2j -P0 2(P2-P0) (208) -Ul(Pl-p) +da (208) where, as before, U2 = Z2- Z U1 1 z 0- z0 Type 3 Segment ( Linear) Specification: (z2' P2) 2' 22 (z1 p1 ) P1 P(1), P2 P(z2) __ 67

013630-9-T The equation of the segment is clearly P2 -Pi p p1 t+ (z-zo ) Pl z2 - z1 and the volume contribution is (209) 6Vo z1 p2 r p dz i.e. 7 2 2 6Vo -(z -z ) (p +p2p ) 2 1 2 P l+ pi (210) which is positive or negative according as z2 > z, z2 < z1, respectively. 68

013630-9-T 5. NUMERICAL SOLUTION OF THE INTEGRAL EQUATIONS The numerical procedures involved in finding 33 /V, 11 /V 0 and where appropriate, 33 / V are quite similar to those required for 11 /V0, and it is therefore sufficient to give full details only for 11 /V. 5.1 P11/V0 COMPUTATION The primary task is the solution of the integral equation ( 163 ) for the function T1 ( s ) and this entails the determination of a sequence of values T ), i = 1, 2,...., N, approximating T1 ( s ) at the sampling points s = s. on the profile p x p( z ). For this purpose the profile is divided into N cells C. of arc length As i and midpoints s. corresponding to the coordinates (Pi zi ) Within each cell we also define the points si and si+ where s. s -i a Asi 1 0 i (211) s. = - a As. 1+ i 0 1 with the restriction 0 < o <1 By assuming that T ( s) has the constant value T1 over the i th cell, 69

013630-9-T the integral on the left hand side of (163) can be evaluated as a linear com(i) bination of the T1 whose coefficients depend on the position (p, z) of the field point, leading to a linear system of N equations in N unknowns, viz. (212) i l 2,...., Ni Hence, the system to be solved is At b 1 (213) where t is a column vector with elements 1 tli = T1 1 i x 1,l 2,...., N, (214) A is a square matrix with elements X 0 p'Klds' ai. 1J ir-j 1, 2,.. *.., N, (215) Cj and b is a row vector wth elements and b is a row vector with elements b = 2rpj, j,1, 2,...., N. (216) Increasing the complexity of the quadrature technique used to evaluate the integrals I will generally improve the accuracy but will almost C. J 70

013630-9-T certainly increase the computational cost. What is therefore desired is the least expensive procedure capable of giving the required accuracy. The two simplest approaches are to integrate first and second order approximations to give (i 9 j): ai p K (i, j) asj (217) a w [j Kl(ii p-)+ pj+K (i, j+ +Wo(o) PKl(i, j) As (218) respectively, where the subscripts j- and j + correspond to the points s and sj+ of eqs. (211), and K (i, j) isthekerneldefined ineq. (29) J- J' 1 1 and evaluated at the points (pi, zi ) (Pj, zj). By requiring a o < 2 we ensure that the sampling points s. and s+ do not coincide with the end points of the cell C., and thereby avoid any difficulty in the computations of Q 1 and n2 (see eqs. 196 and 197). When -- (219) a0 2 5 eq. (218) reduces to the three-point Gaussian formula for which 4 5 w 9 w1 1 (220) l9 18 With this choice of w0 and w, the advantages of eq. (218) vis-a-vis eq. (217) were now determined by computing P11/V for a sphere using various values of N. Fig. 1 shows percent accuracy and C. P. U. time versus N for each integration scheme. It is apparent that for a given expenditure of C. P. U. time the Gaussian three-point technique is much more accurate than 71

0 5.0 4.0 4.0 percent 3.0 error 2.0 3.0 (sec.) 2.0 1.0 0 C 0 3 40 Fig. 1:pgroenterror aindaCPrule calpsoulation for a BPLre T I~ rv-sorwe O0~1 tion Gzd d'S~tS5 thrW P 72

013630-9-T the trapezoidal method, though the accuracies of both are severely degraded if N is too small (N < 5). Since the Gaussian scheme with N = 10 produces an accuracy of better than 99. 8 percent for a sphere, there is no point in going to a more complicated procedure, and the computer program was therefore written using three-point Gaussian quadrature to determine the matrix elements aij In summary, the integral equation- (163) is solved by conversion to the matrix system (213) in which. [ (P K(i. j-)+ pK1 (i, ji )+ PjK(i.j)] As i,j X 1, 2,...., N; iA j (221) I 16pi a a [In -(\ IA] s i i = 1, 2,...., N. (i) T() P/ V is Having determined the sampled values T T (s ), 11/V is 1 1 i 0 computed from eq. (164) by integration over each segment of the profile using a second order integration procedure (subroutine INTEG, described in the Appendix). 5.2 P33/ V COMPUTATION The point sampling method of solution of the integral equations (165) and (166) requires us to find the sequences T2 = T2(s ) and T3 = T3( s ) 2 2 1 i 3 3 1 (i) i = 1, 2,...., N, from these equations. To determine the T, choose, and w in accordance wi (219) and (220) and thence solve the O 0' w0 and w1 in accordance with (219) and (220) and thence solve the

013630-9-T matrix system At2 = b where 21 t = T (i) 2i 2 bi 2r i i i= 1, 2,...,N (222) and aii 18 (p K0(i,j-) + pj+Ko(i j+)j+ p K0(i j) As. -J i, j = 1, 2,...., N; i j j, (223) a= 9n 16pi + 1 As ii Asi i = ls 2p o 9 o * * No The T are similarly determined by solving the matrix system At = b 3 3 where the elements a are again given by (223), but t TX () 31 3 i = 1, 2,, N (224) b X 2r. i The quantities C/e, and 33 /V defined in eqs. (167), (168) and (169) respectively are computed using the same second order integration procedure p employed in calculating 11/V0 If the body profile consists of two discrete parts, it is also necessary to solve the Integral equation (170). The corresponding matrix system is 74

013630-9-T almost identical to that in (224), and from the sampled values T3 (s ) 3 i andT (si), 6 P33/V0 (seeeq. 171) is computed and, hence, 33/V. 5.3 M11/V0 COMPUTATION The basic approach is similar to the above in spite of the more complicated integral equation (172) that must now be solved. The matrix equation for the sampled values V (s ) V( is Av * b where 41 i 4 4 2.... N (225)V () v4i 4 i = 1, 2,...., N (225) bi 2rp. i i and ai - [ f(ij-) + p+f(i j +)+ 4 f(ij) As. i j -9 ' I,,}4 i,j = 1, 2,....,N; i A j, (226) 1 i 2 i a 7 - ii p' f(s') ds' in which f(i, j) P cosa.j 2 (i, j) + (zj - z.)sin.j - p cos a j..., N; i j. Ql (i, j ) (227) i, j = 1, 2,. 75

013630-9-T We observe that the computation of each diagonal element of A requires the numerical evaluation of a Cauchy principal value (denoted by the bar across the integral sign in the above expression for aii ). As an approximation to this principal value, we remove from the cell C. a slice defined by the interval I 1 (s - 3 A si, si + - B As ) where 3, 0 < < 1, is the fractional exclusion; /3 1 implies no exclusions, i. e. that the principal value is not taken. We now have S - S. s + As i 2 si +2 i a Ch - J P p f(s')ds' - p'f(s')ds' 1 + s - -As + A s i 2 i i 2 i (228) and these integrals are also computed using three-point Gaussian quadrature. Defining 1 si2 si 4 i(1+3) As 1 il si2 2 o(10 3) Asi 1 si3 si2 + 2 0(1-. ) As1 (229) si5 - S+ 4 (1+1) As si4 si5 2 0( 13) As. 1 -i6 si5 + a0(1-1) A si 76

013630-9-T we obtain 1 ii ii I ii 2LI Pil f f(i, i)+Pi3 f(i,i )+ Pi 4f(i i4) +p. 6f(ii6)j + (Pi2f(lIi2)+Pi5f('57}9j i (230) Equations (225) through (230) completely describe a system of N linear equations in N unknowns V4i, i = 1, 2,...., N. Their solution and subsequent integration of the V4i according to eq. (173) yield M1/V. Experiments were performed to find an appropriate value for the fractional exclusion. As an example, for a sphere(M l/V 1.5) with N X 20, the data in Table 1 were computed. If we exclude the fortuitous (?) Table 1 1.0 0.1 0.01 0. 001 11/V, 1.480 1.516 1. 503 1.501 percent error -1.33 1.07 0.18 0.08 error zero occurring for fB somewhere in the range 0. 1 < 3 < 1, these data indicate that the choice fS z 0. 001 is sufficient to keep the error less than 0.1 percent. 77

013630-9-T 5.4 SAMPLING RATE Increasing the number N of points at which the surface is sampled will generally increase the accuracy of computation, but since the number 2 of matrix elements increases as N and the cost of a linear system 3 solution increases roughly as N, this improvement is obtained at the expense of an increase in computation cost. Unfortunately, there is no rule for specifying the minimum value of N sufficient for a given accuracy, and the information which follows is based only on our experience in using the program. The results in Fig. 1 and Table 1 show that for a sphere Il/V and 1/0 11 /V0 are accurately determined with N as small as 20, and this is also true of P33/ V0. On the other hand, if the body has a discontinuity in dp/dz lying off the axis, it appears necessary to increase N to 50 or more to maintain the same accuracy (error < 0.5 percent) in the Pll /V "p and 33/V0 computations. This is illustrated by the results in Table 2 for Table 2 N N1 N2 ll/V0 11 7 4 2.752 20 15 5 2.801 40 30 10 2.872 70 50 20 2.888 Note: N1 is the number of sampling points on the generator of the cone (linear segment) and N2 is the number on the (half) base (circular arc segment). 78

013630-9-T a rounded cone with half angle 150. The small but not negligible (0. 58 percent) change in Pll/V0 as N increases from 40 to 70 suggests that such large values of N may be essential for bodies such as this for which T1 ( s ), T2( s ) and T3( s ) have infinities at one or more points on the profile. For the same rounded cone, the results for 11 / V are given in Table 3. Since an increase in N from 17 to 35 produces only an insignifiTable 3 N N1 N2 M11/V0 17 10 7 1. 680 35 25 10 1.678 cant change in M11 / V0 the choice N = 20 is now adequate. Observe that the surface field V (s) associated with 11 /V0 does not become infinite at a discontinuity in dp /d z, and this is undoubtedly the reason why in many cases a small value of N now produces the same accuracy as does a much larger value in the P11/V and P33/V0 computations. No attempt has been made to exploit this finding in the general program. When treating bodies composed of several segments, a strategy which has proved successful is to divide all segments into cells of approximately equal length. This serves to fix the allocation of any given number N of sampling points among the various segments. Tests so far performed have not conclusively shown the advantages of dividing a single segment into two or more smaller segments so as to effect a non-uniform sampling. It is, however, believed that such a sub-division may, for a given N, improve 79

013630-9-T accuracy in the 11 /V0 and 33/ V0 computations for bodies like the rounded cone having infinities in the surface field quantities. 80

013630-9-T 6. CONCLUDING REMARKS We have here considered the low frequency scattering of electromagnetic and acoustic waves by axially symmetric bodies, By concentrating on certain quantities such as the normalised components of the induced electric and magnetic dipole moments, we have shown how it is possible to arrive at rather elegant expressions for the far zone scattered field in terms of quantities which are functions only of the geometry of the body. Each such quantity is expressible as a weighted integral of an elementary potential function which can be found by solving an integral equation. A computer program has been written to solve these equations by the moment method and to calculate the dipole moments, the electrostatic capacity, and a further parameter y related to the capacity. Any body can be treated whose profile is made up of straight line and circular arc segments and it is even possible to have two distinct bodies with or without an electrical connection between them. Although no serious attempt has been made to optimise the program, only a few seconds are required to compute all of the above quantities to an accuracy of better than one half percent. We have already used the program to compute the scattering from a variety of shapes, and it may be helpful to list some of the results obtained so far. Data for a rounded cone consisting of the intersection of a cone of half angle e with a sphere centered on the apex are given in Table 4. A/w is the length-to-width ratio of the body. For 0 < 900, the values of 11/V0 and 33 /V are quite similar to those previously computed by Senior (1971) using a mode matching method, but since 11 /V showed significant discrepancies, this quantity was determined for a variety of 0. Detailed checking has confirmed that the present data are accurate 81

013630-9-T Table 4: Rounded Cones V0 P33/V0 (deg.) I/w 3 7.5 15 30 45 60 80 90 93 99.2 108 120 140 151.7 180 9.554 3.837 1.932 1. 000 0.7071 0.5774 0.5321 0.5000 0.5262 0.5799 0.6545 0.7500 0.8830 0.9402 1. 0 0.002870 0.01792 0.07137 0.2806 0. 6134 1. 047 1.731 2. 094 2.204 2.429 2.742 3.142 3.699 3.938 4.189 2. 865 3.664 8.147 3.494 4.520 1.931 11 /V 1.884 1.813 1.678 1.484 1.366 1.312 1.334 1.373 1.386 1. 416 1.458 1.507 1. 547 1. 540 1.5 C(e ) 5. 406 6. 386 7.428 7.303 7.123 6. 889 6. 586 6.441 6.283 4.428 4.368 4.261 4.071 3.789 3.370 3. 187 3.0 2. 184 2.242 2.372 2. 553 2.769 3. 006 3.042 3.0 82

013630-9-T to three significant figures. When 0 = 900 the body is a hemisphere for which precise values of P11/V0 and P33/V0 are available: 11 /V 4.430...., 33/V 2.189.... (SchifferandSzego, 1949, p. 152). The corresponding values in Table 4 are within 0. 2 percent of these. For 0 > 900 the cone is a re-entrant one, i.e. a sphere with a conical region removed, and when 0 = 1800 the body is a sphere for which exact data are also known. Results for ogives and symmetrical lenses whose arcs subtend an angle 0 at their centers of curvature are shown in Table 5. The transitional shape is a sphere for which 0 X 180. To illustrate the computations when two bodies are present, Table 6 gives data for two identical spheres separated by a distance c d where d is the sphere diameter. When the two spheres are touching (E X 0) it is known that P33/Pll 8/3 (Schiffer and Szego, 1949, p. 154); the ratio deduced from Table 6 is 2. 678, which is within 0.4 percent of the exact value. As E increases, P11/V0, M11/V and P33/V 0 3/V0 rapidly approach the values appropriate to a single sphere in isolation. 33 /V0, on the other hand, is proportional to the axial component of the induced electric dipole moment for two spheres which are electrically connected by an infinitesimal wire, and with increasing e this increases indefinitely, as expected (Kleinman and Senior, 1972). The same is true of C /(. The parameter y has also been included in Table 6, and since its exact value can be shown to be -(1 + e/2 ), the accuracy of computation can be judged. 83

013630-9-T Table 5: Ogivee and Lenses Shlape 9(deg. ) t/w ogi. 22.8 36 ye 56 88 132 150 iere 180 64.4 43. 6 is 28 17. 2 10. 02 6. 314 4.011l 2. 475 1. 540 1. 303 1 0. 6297 0.4000 0. 2493 0. 1512 V 0 0. 004146 0. 01058 0. 02650 0.07148 0. 1966 0. 2847 0. 5236 1. 212 2. 586 6.448 17. 30 1/0 2. 089 2. 100 2. 189 2. 363 2. 647 2. 775 3 3. 779 5. 182 7. 649 11. 88 17.49 p33/V0 49.88 24. 15 12. 57 6. 778 4. 136 3. 595 3 2. 161 1. 674 1. 390 1. 225 1. 144 1.943 6. 128 1. 910 5. 696 1. 843 5.451 1. 739 5. 424 1. 611 5. 696 1.564 5. 880 1. 5 6. 283 1. 368 6. 570 1. 252 7. 375 8. 758 1. 098 10. 759 1. 061 12. 982 -M i/ i/0 c / (c/7) sp1h len 11.4 0.09981 39.55 84

013630-9-T Table 6: Two Spheres eI/Vo P33/Vo 33/VO I1/V0 c/(E ) 7 0 2.702 7.237 7.237 1.605 6.153 -1.0004 0.005 2.706 7.280 5.021 1.586 6.151 -1.0029 0.01 2.709 7.319 4.800 6.150 -1. 0054 0.02 2.715 7.402 4.557 1.607 6.148 -1.0100 0.05 2.732 7.655 4.210 1.592 6.141 -1.0250 0.1 2.759 8.026 3.922 1.579 6.130 -1.0503 0.5 2.891 12.02 3.299 1.528 6.015 -1.2502 1.0 2.950 18.19 3.142 1.511 5.822 -1.5001 5.0 2.994 120.6 3.045 1.500 4.384 -3.5001 10.0 3.002 383.1 3.046 1.501 3.470 -6.0000 85

013630-9-T REFERENCES Abramowitz, M. and I.A. Stegun (Eds. ) (1964), Handbook of Mathematical Functions, NBS Appl. Math. Series No. 55. Forsythe, G. and C.B. Moler (1967), Computer Solution of Linear Algebraic Systems, Prentice Hall, Englewood Cliffs. Jahnke, E. and F. Emde (1945), Tables of Functions, Dover, New York. Karp, S. N. (1956), Relation of Electric and Magnetic Dipole Moments (abstract only), Quart. Status Report No. 1, Contract AF 19(604)1717, Institute of Mathematical Sciences, New York University. Keller, J. B., R. E. Kleinman and T. B. A. Senior (1972), Dipole Moments in Rayleigh Scattering, J. Inst. Math. and Applics. (to be published). Kleinman, R.E. (1965), The Rayleigh Region, Proc. IEEE 53, No. 8, 848-856. Kleinman, R. E. and T.B.A. Senior (1972), Rayleigh Scattering Cross Sections (to be published). Mautz, J.R. and R. F. Harrington (1970), Computation of Rotationally Symmetric Laplacian Potentials, Proc. IEE (London) 17, 850-852. Payne, L. E. (1956), New Isoperimetric Inequalities for Eigenvalues and Other Physical Quantities, Comm. Pure Appl. Math. 9, 531-542. Schiffer, M. and G. Szego (1949), Virtual Mass and Polarization, Trans. Amer. Math. Soc. 67, 130-205. Senior, T. B. A. (1971), Low Frequency Scattering by a Finite Cone, Appl. Sci. Res. 23, 459-474. Van Bladel, J. (1968), On Low Frequency Scattering by Hard and Soft Bodies, J. Acoust. Soc. Amer. 44, 1069-1073. 86

013630-9-T APPENDIX: THE COMPUTER PROGRAM The program computes ll/Vp, C/e, 7, 33/Vo 11 / V0 and, where appropriate, 33 /V0, and consists of a main program and six subroutines. A. 1 DATA SET A data set is made up of one control card and a number of segment specification cards, one for each segment (or sub-segment) of the profile. The segment specifications conform to the convention stated in Section 4. 3. Control Card Columns Description 1 The number (1 or 2) of bodies. 3- 4 Two digit integer (right justified): the number of segments on the first body (the body to the left). When there is only one body, use these columns. 6- 7 Same as columns 3 - 4, but for body to the right. 9 A printing key: 1: print T from 33 / V com3 o putation. 0 or blank: do not print T. 3 87

013630-9-T 11 A computation key ( 0, blank or 1 ) 1: suppresses computation of 11 / V C/ T P33/V0 13 A computation key (0, blank or 1) 1: suppresses computation of 11 / V 21-30 A real number: the fractional exclusion /3. If these columns are blank, 13 defaults to 0.001. Segment Specification Card Columns Description 1 - 2 Two digit integer ( right justified): the number of sampling points or cells on the segment. 4 Segment type key: 1: circular arc, concave down 2: circular arc, concave up 3: linear. 6 Volume sense: + or blank: additive volume -: subtractive volume. 11-20, 21-30 Two real numbers: respectively, the end coordinates z 1 and z2 of the segment. 88

013630-9-T 31-40, 41-50 Two real numbers: respectively, the end coordinates p1 and p2 of the segment. 51-60 A real number: for circular arcs, the included angle in degrees. There are the following restrictions: (i) the total number of segments must not exceed 15, and (ii) the total number of cells over all segments must not exceed 80. The profile is specified in the direction of increasing profile-length, beginning at its left-hand intersection with the z-axis and ending at its right-hand intersection with the z-axis. Re-entrant segments are permitted, allowing z > z2 i. z2 89

013630-9-T A. 2 MAIN PROGRAM The main program reads and prints data and supervises all computations. A rough flow chart showing the interaction of the subroutines is given below For each segment 1) read speclfication card 2) 90

013 630O-9-T kF AI. AP Il( 80,t80) pAP 33(80,8CHO) -AMII(R80,80),X R t(80),K(UI IN\TEGE)P, NUMPTS (15 ),tPLUS/ I+'/,kL/' /I NOX(2 C,f1mm(1N P.HO(80,9 ) Z(80,99),gARC,(80),gC(80,9),S(P0,9)/SU)L/'IPS(8')) DATA MINI,TwnPI,PI,WO,Wl/'-',F.?83185,3.141593,.4444'444,.?777-117/ 37 kFAD(5,34,FND=999 )NBOD,NS1-,NS2, IPklINT,KFYP~iI,KFYMIIt,Fk~ 34 FF)PMAT( I1,1X,?(I2,1X)93(I1,1X),6X,F1?.7) WR I TE ( 6,t4)N80F0vNS I 4 FOR MA T('1 R*F(-FIGINNI NG0 F DATA S ET:'/0 v5X v F~(]Dl ES',X,' -,1 2/ I F ( NP0D.L-F. 0 ORk. NBOI) GT. 2) GO 10 9904 N\J S F (; = NS SI + N 5,'? I F ( NR.)FQ. 2) W RI T F( 6, I) NS I FnRMAI( II,5X,'8fRODY h? I,5-X,'=I',1 2 I F (NKI).FQ. 1.) GO) TO 1 0 IF(NSI1.L-F. 0 n0\. ~s2 L.LF. 0) GO TO (90( 10 TF(NSFGS,'.[.F. 0 f0R. NSE-GS.GT. 15) (;n TO 4990 VIRTF( F,935?) IPRINT,)<FYP11,KFYMI1I 52 FOR M A T( 5 X,'P RI NI K F Y'I93 X v',I2/ c',5X v ICOMP K:Y P''z,? I,5 X'C (M P KEFY M''S=I'I,12) IF( FP.L-F. 0. nOR. FR.(;T. I * FR=.001 IF(KFYM11l FQ 0) WRITF(6,5)FP\ 5 FOPMAT(' ',v5X,'FXCl-(JSI0Nr =',F7.4) IF(KFYP1.1.NF. 0 SAND. KEYM1I..NF. 0) G0 1(- 990) N =0 NC I = 0 V 0=0 *0 DOl 11 1=1,NSEGS RFAD( 5,1?) NIJMPTS( I),ITYP,ISIGNI,ZFP.,RHnFP,THFTA 1.2 FF)kM AiT( T?12,IXIT1.tI IX,9AIt4XI5F1l0.7 ) TIF(NUMPITS(I ) LF. 0.OR I TY P FL. 0.OR I TY P.61. 3) G To 990 I F(I SI(;N.F('. K~L) I SI CN= PL-US WlRITTF ( 6 T3 )I,-pNUMPT S(I),ITYPrI SI GN Y7 F PPHO-F P 1 3 FnPMA1'O SEGMENT #1',2,':/ 5X (-,FL I.-,7X,'='2/ ~,Xv 6 'TYPF K(FY I,4X '=,I?/',5X p VOLLMF SFNSE'= I',A4I',X, 8 7 7- C 0 RDI N AT E FN D P 0 1 NTS =S,12*7',',l2*1, )/'5X, 5 'PHn-CnflRDINATF ENID PnOINT S=,F?7,'1?1)) IF(ITYP.NF. 3) WRITF(6,1l4) THFTA 14 F(1RMAT(I 1,5X,'THEIA (DFG) =',FI0.5) IF(KFYM11 F() 1) FR=1.0 N1 v TH11F T A= P I *T H FT A 180. I(I.F. NS1 ) NCI =NC I +N\IiMPTSQ~( I) I F(M.61. 80) GO TO 99 0 (IAl L 1DAJlA( ITYPNM,/ZFPRHOFPtTHETAFRgVO-LINC,) IF( ITYP.NF. 3.AND. ISIGN.FO. MIN) VOLINC=-VOLIIg~C I11 v V= v o+ vOLt-I ~c W I TEF(F,v5-2)VO F DRM A T ( IO0C IIM P TEFO R F SIJ)L T S:'/' X,5xI'V ()L HME'F 6X,''F 1. 91

013630-9-T OF) 2 Ni=l, m I NnX (1 I=N~ A N =A RC (N) TN=RHO( N,8) nn 3 L=N)?M IF(L FEO. N) GO1 TO 82 A L =A RC ( L ) TL=RHn(iL,8R) IF(KFYMll-l)l10,n109,1.09 1 10 AMIlI(NL-)=0.0 A\M.1 1-(,N) =0 * I. 09 IlF (KFYP 11-1)11 1.,112,1 1.2 Ill APII(N,L)=0.0 A PI11(L,v N) = 0 * A P 33 ( N? L ) = 0. 0 A P33 (LvN )=0,0 112 INDX(2)=L 92

013630-9-T nr! 103 J=l,3,2 J P 6=1 + 6 )r) 104 LL=1.? T =3-LL I 1=INOX( LL) J?=INN)X( I ) TI2=RHf( I2,JP 6) CALL SFTU P(KFYP11,KFYMl, 1, I? JP6,API 1,API33,AM l 11, 1 ) I F (KFYP11-1)105,106,106 105 AP1 I ( I 1, I 2) =APll ( I ], I 2 )+AP I 1 ),T I 2 AP33( I 1, I2 )=AP33( I 1, I 2 )+AP I 33'T I?2 106 I F (KEYMl-1) 107,104 104 107 Ar 11 ( I, I2 2) =AM11 ( I 1, I 2 )-AMI 11 ^T I2 I 04 CON! TI NIlF 103 C(INT I IlJF CALL SFTIJP(KFYP11,KFYMll,N,L,T,API 11,API33,AMI I1, I ) F ( KFYP 1 1-1 ) 1 08, 09,?09 1 OR U=-^JO'AP I 1. 1 AP1 1 (., L ) =AL (Wl AP11 (N,L ) +1 TI ) AP I. ( L,N) =AN ( W 1 AP 1 1 ( L,N)+t' TN) JI=W^0,:API 33 AP33(N,L )=AL (W 1 AP33(N,l )+UtTL) AP3 3(, N ) =AN' (W 1AP33 ( L, N ) +lJTN)?09 IF(KFYM11-1 )210,3,3 210 A 411 (N, L )=AL* (W1 AM11 ( N, L )-AMI 1 W TL) CALL SFTlJP(,0,L,N 8,AP I 11 API 33,AMI 11,0) AM11 ( L, N ) =AN,* (W 1:AM11 ( L, N )-AM I 1 -"WO'TN ) GfO TO 3 2? IF(KFYP11-1 )R3,84,84 R3 Ij= AL On( 16.lTN/AN) AP11 (N,N)= ( J-1.0) AN AP33( N,N) = (U +.0 ) -:AN a4 I F (KEYM 1-1 )85,3,3 2 5 TF(FR.FQ. 1.0) c0 TO 3 nIn 6 I=1, 6 a6 CALl SFTUtP( 1,0,,N, I,U,tJ,ST3( I ), 1 ) I1=. 5'( 1.O-FR) ':AN AMII( N,N)=PI-tJ (WO ( RH ( (N,2 )ST3(2 )+RHO(N,5) S: T3(b) ) 3 +W1*(RHF(N,1 )*ST3( 1 )+RH (N,3 ) ST3( 3 ) +RHO(N,4) S ) 3(4)+ 4 RHn(N, 6) 'ST3( 6) ) ) 3 CONTINUE 2 CI -N TI Nl I F 93

01363 0-9-T 20 R(I)=TWOPlI-,RHO(1I,8) IF( KFYPL 1-1)21,24,24 2 1 CA LL D Fcnmp(AP1I1vI,M) CALL SnLVE(APL],qXKM) nO 22 1=1,9M CALL INTFG( X NSEG7SNUMPTS vP11 P1 1=P1, —I )Z/VO 24 I F (K FY MlI-1)?,2 25 CALL. DFFcnmp(Amil1, m) C A LL SlL V F (AMII vXv PvM) )() 26 I1,11M CAI.L I NTFG ( X IN1S FG;SN1JM PTISM 1l) M I11 = P I,- Ml I /vo0 2 8 I F ( KEFY P I 1i ) 32,9 45,v4 5 3i2 nn 29 i=i1,m 29 Ki( I ) = TV, IP I1*7 ( I, 8 CAL.L DFcnmp (AP33,M) (CALL SF)LVF (AP33,qXqR, M) DO 134 J=1,M X ( I =Z (I, 8 ) RHO)( 8 ) X(I 134 PII)=TWOPI CALL INTEG(XNSFGSNUMPISP33) CALL SOLVF( AP33,130,8 M) Of] '35 I=1,M X( I)=RHn( i,,8)*13(lI) 35 KI=(,)XI C A LL I NTFCG( X vN,NSE(GS \1IJM PITSy C,CAP) CAP= TWOP I CAP CA4L L I NTF(; ( K, N1SEFGS 1\IJM PT Sv GAM) C A VI = -1wniplI * GA M / C A P) P'33~= TW40P I *P 33 - C1\P G A M* (A M)/V() IF(NRnn - 1) 39,t39,54 5 4 1)0C 36 = I =1,M I F( I - N CI) 3 05, 3 05,- 3 06 305 K~( I)TWOP I ((I) T10 36f 306 ( I )=0.O 36 (JiNT I NI JE 94

013630-9-T CAII '~f \/ F (AP33,vTI'?,RVM) 0(1C 3 07 1 =I,9M 30O 1 X ( I )= H H,8)T11 (I) CA LL INKITF(, XN)SF(;SNIJMPTS 9iTL CALI-L IN11TEC,( X,-NS1,NUMPTS y TN) DO '308 I lM 308 R X ( I ) =Z7 (I 8,8)- X( C A LL I NTFG(, X,9NS- F G SNU()M PT S, I) DE FL TA4P=-( VIWOP I /V ) t, ( U+CA M TL )2 / (J I- TW OP I ~TL I L/CAi\P t I= P3 3 + DF LITA P 40 FnRIMAT(l ',5X, 'CAPACITANCE =',F1O.5/1 ',5X, G'(MMA',7X, '~',FlO.5/ 8 ''5X, *P11/V1 7X, = ' FIO.5/1 ',5X, 'P33/V',7Xt, '' 'Ft. 5) ILF (NJR n.F (. 2 ) WRI TE( 6, 3 09 )DE L TA PU. 309 F(IPM A T(I 1,5X,'DFLT P33/V',2X,'=',FlO.5/',5 Xv LI SJ NT P 3 3/Vz' 6 FIO.. 5) 45 I F (K FYNMIli1 -I) 4?,9337,v3 37 42? PRI TF ( (-,4'3 )MI11 43~ F0PMAT(' ',5X,'M11./V',7X,'1=',F10.5) ~i337 I F (I P RI NT FEQ. 1.A ND KEY P I.1EO. 0) W RI TE('6,44) 44 F(RZA(1,8) 5,RH9O( 18),1(R)HI=,)'t2gT'( G Ii TO '37 99 0 WPJR ITE ( 699 1 9 91 FOnR MA T (I 04: --- ERR OR I N D)A TA' 999 CALL SYSTEM E ND 95

013630-9-T A. 3 SUBROUTINE DATA (IN, MX, MY, ZEP, RHOEP, THETA, B, VOL) This subroutine is called once for each segment of the profile. From the input specification for the segment, DATA computes the (z, p) coordinates of the necessary sampling points on the profile, the quantities cos a and sin a at these points and the incremental volume of the segment. Arguments: IN Type key for segment. MX Total number of cells in segments to the left. MY MX + (number of cells in this segment). ZEP z-coordinate end points of segment: ZEP ( 1 ) z ZEP(2) X z. RHOEP p- coordinate end points of segment: RHOEP(1) p1, RHOEP (2) p2 THETA Angle (in radians) subtended by a circular arc at its center. B Fractional exclusion, 3. VOL Incremental volume of segment. Comments: Stored in COMMON are the arrays RHO (80, 9), Z(80, 9 ), ARC (80), C (80, 9) and S (80, 9) which contain the numbers computed by DATA. For the Ith cell, the subscripts (I, J) correspond to the points s. of (229) when 1 < J < 6. For J x 7, 8, 9, the subscripts (I, J) refer to the points s, s and s + respectively of (211). 96

013630-9-T SUBROUT IN E DATA( I NtMX, MYvZEP,,RH OEP,TH ETA IB, VOL) DIMENSION ZEP(2),RHOEP(2) COMMON RHO(8099),Z(80,9),ARC(80),CU80,9),S(8099) DATA STEP/.3872988/ MXPl =MX+l FN=FLnAT( MY-MX) TF(B.NF. 1) SliBSTP=.5*(l.0-B)*STFP IF (IN-2) 1,2,3 cc=-100 (;( TO 1 0 CC = 100 0 T12=S I N(T HFT A /2.O) A= /F P (2 )-7FP ( I1) RAD=0 5*SORT((RHOFP(l)-RHOEP(2))**2+A*A)/ST2 DD)=A/ABS (A) T=GCC"DD%`COS (THE-TA!?.) 1ST2 Z CN T=0. 5-( Z EP (I) +Z EP (2 )+T* (RHOFP (1I)-RHOFP (2) RHOCNT=0.5,-(RHO)EP(1)+RHO(EP(2)+T*A) 02= E-P (2 )-ZCNT U1l=ZEP (1)-ZCNT VOL= 3. 141 59 3*ABS (A"'(R HOCNT4*2+R AD*R AD- ( t2,*2 +UJ1 *02+ tj I2)/3. 0) 3 -CCx,*RHOCNT4l,(Li2*(RHOEP(2)-RHOCNT)-t)l*(RHOEP(1)-RHU-CNT) +kRADf)-R A D 3 )':-DD*IHFTA)) KF TA=CC4DD —THF-TA/FN TH FIL=1A TA N2 (RHO )E P ( 1R HOC NT Z E P(1Z C NT) U.=ABS (KFTA*RAD) K3=S TFP*8ETA DO 902 I=MXP1,MY PH I =THF Tl ( I-MX-. 5) *RFTA IF(B.FO. 1.0) GO0 TO 1905 DOl 1 9 02 J = 1f2 ANG=PHI+.5*LJ-1,5)*BETA*(l.O+B) DO 1903 L=1,3 PSI=ANG+(L-2 )*SUBSTP*BFT4 M=L+3'(J-1) C( I M)=-CC*SIN(PSI) s(1, M) =cc',cOs (PSI) 7(1I,M )=ZCNT+RAD*CC*S ( I,M) [903 R HO((IvM) =RHOC NT.-CC*R AD*C(IM) L90o? CON TI NU [905 DOM 903 J=7,t9 ANG=PHil+( J-8 )*fi3 C( I,J)=-CC*SIN(ANG) S( I,J ) =CC*COS ( ANG) 7 (1 J J) =7CNT4-RAD*CC'*S ( I, J) (03 RHU(I J)=RI-iOCNT-CC*RAD*C(IvJ) 0O2 ARC( I)= P. F TURN 97

013630-9-T 3 DX=( /LP(2 )-ZEP( 1) )/EN DY= (RHnEP (2 )-RHOEP (1)) /EN U =SORT(DX*DX+DY*DY) SI =DY/0 CLI = DX/ U DC) 917 I=MXP1,MY PHI=FLOAT( I-MX) —.5 IF( B.EO, 1.0) GD0 TO 1800 DO 1802 J=192 ANG=PHI+,5*(J1-1.5)*(1.0+B) C)n 1803 L=l,3 M=L+3', ( J-1) PS I=ANIG+ (L-2 )*4SUBSTP Z (1,M)=ZEP( 1)+PSI*DX RHC( IIM)=RHOFP( 1)+PSI*'DY S(I,M)=SI 1803 C(IVM)=CI 1802 CONTINUE 1800 Dn) 913 J=7,9 ANG=PHI+( J-8 )*STEP Z (IJ )=Z EP (1) +ANG*DX fkH0(I,J)=RHOEP(1)+ANG*DY C( I,J)=CI 913 S(I'VJ)=SI 917 ARC(I)=U VO)L=1.047198*(ZEP(2)-ZEP(l))*(RHOEP(1)**2+RHOEP(1H*RHOEP(2)~ 8 RHI1EP(2)**2) REFTUR N FND 98

013630-9-T A. 4 SUBROUTINE INTEG (V, NSEG, NUMPTS, SUM) INTEG numerically integrates quadratic interpolating polynomials approximating the data on each segment of the profile. When the profile is composed of several segments, no interpolation is performed across segment boundaries. Hence, the integration is accurate even for disconnected segments, e. g. the circular arcs of two spheres. Arguments: V Real vector of function values, ordered as the cells. NSEG Total number of segments in the profile. NUMPTS Integer array containing in NUMPTS (I) the number of cells on the Ith segment: I = 1, NSEG. SUM Integral of V across the profile. Comments: Stored in COMMON are the arc lengths ARC (I), I = 1,...., N required to compute the integral. 99

013630-9-T SUFRRMJTINE INTEG-( VNSFG',.NUMPTSSUJM) cOMMONi RHO(80,9 ),Z(80,9 ),ARC(80),C(80,9),S(80,Q) DI MFN)SiON) V (80),NUMPTS (15) SUJM=0. 0 J ACC,= 1 D~O 3000 I=1,NSFG' T=ARC ( JACC ) i =N~tJMP TS( I ) SiN=IS~JACC-i,65:( (JC )V JCl)+(I I F(L/-2 t\NF. (L+l )12) GO TO 3001 SIM SJ+ (0 6 0 7* (N 1)-, 8 3 3 3-(N 2 +. 1 6 6 -V 3001- I-Ml=N-1 JL-O=JACC+ 1 00 -300? j=j.LnOLMl,? 300 2 SIM ()+.3 3 3 "- J 4 "V(J)+ J 3000( JiACC=JA /CC+LP F TI RN FN[) 100

013630-9-T A. 5 SUBROUTNES DECOMP (A. N) AND SOLVE (A, X, B. N) Used together, DECOMP and SOLVE solve the linear system AX = B. DECOMP performs a L-U decomposition of the N x N matrix A and SOLVE performs back-substitution. These routines are adapted from Forsythe and Moler (1967, pp. 68 -69). SIR ROUTINE ECnOMP(tJ_,N) DIMENSION t)L(808,O) CnMMON /SOL/IPS(8O) DO 5 1=1,N 5 IPS( I )=I NM1=N-1 DO 16 K=1,NM1 RIF=O.O DO 11 I=K,N IP=IPS( I ) IF(ABS(I)L(IPK)).LE. BIF) GO TO 11 RIF=ARS(UL( IP K ) ) IDXPIV=I 11 CON TI NUE IF(IDXPIV.EO. K) GO TO 15 J=IPS(K) IPS(K)=IPS( IDXPIV) IPS( IDXPIV)=,) 15 KP=IPS(K) PIVOT=tL (KP,K) KP1=K+1 DO 16 I=KP1,N IP=IPS( I) EM=-JI ( IP,K)/PIVOT tl( IP,K)=EM DO 16 J=KP1,N UtJ ( IP,J )=UL( I PJ)+EM*tUL(KPJ) 16 CON TINUE R E TR N END 101

013630-9-T SUBROUT INESOLVE MIL,X ptN) DIMENSION UL(RO8O,8) B(80),X(80) COMMON /SOL/IPS(80) NPI=N+l I P= IPS (1) X (1) =R (IP) DO 2 I=2,N I P= IPS (I) IM1=I-1 SIJM=0. DO 1 j=1iTMI SUJM=SUJM+UL (I IP,j)*X(,J) 2 X(I)=B(IP)-SU)M I P= IPS (N) X(N)=X(N)/U)L(IPpN) DO 4 IBACK=2,N I=NPl-IBACK I P= IPS (I) I P1=1+1 S(UM= 0.0 00 3 J=IPIN 34 SJM= SUM+UJL (I PJ )*X(J) 4 X( I )=(X ( I)-SUJM)/L (I IP I) R ETURN END 102

013630-9-T A. 6 SUBROUTINE ELLI (M1, K, E, KPR, KEY) This computes the elliptic integrals K ( m) and E (in) and the derivative K' (m) from their power series approximations (see Section 4.2). Arguments: M 1 Real, the quantity ( 1 - m). K Real, K (m). E Real, E (m). KPR Real, K' (m). KEY Integer: 0 Compute K, E and KPR; 1 Compute K, E but omit KPR. SIJRRnFtRlTINF ELLI (MI,K,E,KPRKEY) REAL M1,KKPR T=-ALOG(MI) K=I.386294+.5*T+Ml*(9.o666344E-2+.l249859*T+Ml*(3.590092E-2 5 +6.880249E-2*T+Ml*(3.742564E-2+3.328355E-24I+Ml*(1.451196F-2 5 +4.41787E-3*T)))) F=1,0+Ml'"(,443?514+.2499837*T+Ml*(6.260601E.-2+9,20018E-2*T~MI4 R 4.757384E-2+4.069698F-2*T+Ml*l1736506E-2+5.264496E-a*T)))) TF(KFY.EQ. 1) RETURN KPR=.5/MI +?.832?5F-2 -.1249859*T + M].*(-2.999362E-3-.137605S'T' A, +Ml*(-7.R9933-,E-? - 9.98506F,-2*T + Ml*(-5.362998E-2 - 5 1.7-714RF-? ~ T )) RF TURN FND 103

013630-9-T A. 7 SUBROUTINE SETUP (KEYP1l, KEYM11, I, J, L, APIll, API33, AMI11, IJ) This is essential in computing the linear systems. Specifically SETUP, after calling ELLI, computes the quantities API 11 (K1 of eq. 182), API33 (K0 of eq. 178), Q1 (eq. 196) and n2 (eq. 197). The quantities Q1 and 2 are used to compute f(i, j) (AMI 1) of eq. (227). Arguments: KEYP11 0 when computing API11 and API33, else 1. KEYM 11 0 when computing AMI 11, else 1. I Subscript of observer (unprimed) cell. J Subscript of remote (primed) cell. L Index of the point within remote cell for which the kernels are to be computed (see DATA, Comments). API11, API33 Described above AMI11 IJ 0: use last value of M 1 in kernel computations; 1: compute new M1. 104

013630-9-T SUBPROUTINE SETIJP(KEYPI11KEYM11,I,J,,LAPIIIAP133,YAMvI 1v,I)) CnmMON RHn(80,9),Z(8O,9),ARC(80),C(8O,9),S(80,9) PEAL MMlKKPR 7D=Z (JL )-Z (I,8) R=RHnO( i 8 ) RP=RHO( J,L) I F (IJ.FQ, 0) 0 TO 11 5 A1=RRP+RRP A 2= R *R +R P *R P +0*z0 MI= (42-Al )1/(42+41) M=1.-Ml CALL FLLTI(M1,KvEK PRvKEYM11) AO=M/RRP Al=SOR T(A40) A?=M+M A 3 =2.- M IF (KEYPII-1) 113, 114, 114 11~3 API1l1=Al*(A3*K-E-E) /M AP 13-3=41I *K 1 14 A I= A1I* A0 4 3= * 5*43 115 AO=C(JL) nmi=-Al*(.25*K-A3*KPR) nM2=Al*(EF-A3*((A3+M)*K - A2*A3*KPR))/(M*M) AMI11=R*AO*OM2+(ZD*S(J,L)-RP*AO)*nmI 1 hA4 RETURN END 105

_UNC ASSrFIE n. DOCUMENT CONTIOL DATA H & D /-o' l ti,o It,. i (Iit. lt O,f t/t I.', (it>1(-'!I e ietfrtt P t ( wI tl $ It, I e(! i It i tt(l ori r lIt ' tf t eI! tIft. rd wht f Irl n ivt o v ll r If, rfI t r( If ~ei1t (Il) ' Ik,..,N 1>NK; tt 1 >VA~ I I G ArI Ir V I /IC.~ ) 1 t- HI -PO FT Si:CU1I1 V CLASI I CA T ION 'Phc1 Ull iv ' sity of Mij'i/;', llali 'tion Laboratory UNCLASSIFIED 22: 1 S('pace ( c Rsea(lc lll 1., North Campus 2b. CGOUI An1 Arr:\l)M. lic'lli: 11 t8105 ----- 4: t ' t0 I II L t' THE NUMERICAL SOLUTION OF LOW FREQUENCY SCATTERING PROBLEMS OESCRIPTIVE NOTES (Ty)pe a ( rtport nnd Inclusive datlo) Scientific Interim. AU THOq"IS) ('Frst name, mlddJl Initial, I lat nlome) Thomas B. A. Senior David J. Ahlgren RLPOR T DATe 7li. TO"l AL NO. OF r AGE S 7t. NO. Or lncrFS February 1972 105 12 CONTRACT OH GRANT NO. Ou. OHIGINATOfR'5 RLIJORT NUMAr.RIS) F 19628-68-C-0071 013630-9-T. PROJECT NO. Project, Task, Work Unit Nos. Scientific Report No. 12 5635-02-01 ),h. OTHER REPORT NOIS) (Any ot/icrt Imhcrli, that nuiy he. uo^.irncv. DoD Element 61102F AFCRL-72-0162 AFCRL-72-0162 DoD Subelement 681305___________________ D STRI UT ON STATEMENT A - Approved for public release; distribution unlimited. SUPPLEMENTARY NOTES 12. SPONSORING MILITARY ACTIVITY Air Force Cambridge Research Laboratories (LZ TECH, OTHER L.G. Hanscom Field.Bedford Massachusetts 01730 * ABSTRACT The low frequency scattering of electromagnetic and acoustic waves by rotationally;ymmetric bodies is considered. By concentrating on certain quantities such as the lormalised component of the induced electric and magnetic dipole moments, it is shown Low the first one or two terms in the far zone scattered fields can be expressed in terms )f quantities which are functions only of the geometry of the body. Each of these is the veighted integral of an elementary potential function which can be found by solving an ntegral equation. A computer program has been written to solve the appropriate equations )y the moment method, and for calculating the dipole moments, the electrostatic capacity, md a further quantity related to the capacity. The program is described and related data ire presented. \ D F' 1 473.. -—.. —. -I....... FORM 1 A73 1 NOV 5l US.....UNCLASSIFIE D %'t' lrmtl ('I Issi -.,toi, ll

TITNC T XASRTIR....:ntrt lasfcto L I t A LI KJ H C; C FV WONrL)DS. A=_L I......4..... ~ w -T_,, -f...7! Low Frlequency Scattering Electromagnetic Acoustic 1Rotational Symmetry Computer Program Numerical Data I I 1,I ll -U NCLASSUIFIEDn S,", itrIt V ('1,t', I I... l i,,;