THE UNIVERSITY OF MICHIGAN 3673-1-F AFCRL-64-916 Ionospheric Disturbances-and Their; Effects on the Propagation of HF:.Electromagnetic Waves C-NMChu,, G. Hok and J.J.. LaRue Final Report Contract No. AF 19(604)-7257 Project 4649 Task 464901 September 1964 Prepared for Air Force Cambridge Research Laboratories Laurence G. Hanscom Field Bedford, Massachusetts 01731 Bedford, Massachusetts 01731

THE UNIVERSITY OF MICHIGAN 3673-1-F Requests for additional copies by Agencies of the Department of Defense, their contractors, and other Government agencies should be directed to: DEFENSE DOCUMENTATION CENTER (DDC) CAMERON STATION ALEXANDRIA, VIRGINIA 22314 Department of Defense contractors must be established for DDC services or have their "need-to-know" certified by the cognizant military agency of their project or contract. All other persons and organizations should apply to: U. S. DEPARTMENT OF COMMERCE OFFICE OF TECHNICAL SERVICES WASHINGTON 25, D. C. ii

THE UNIVERSITY OF MICHIGAN 3673-1-F ABSTRACT This report is concerned with the effect of ionospheric disturbances on HF radio waves whose propagation path is at a large distance from the source of the disturbance. The perturbation of the HF radio waves is assumed to be due to mechanical waves propagating in the atmosphere. The macroscopic equations governing the propagation of mechanical waves in the atmosphere, including the effects of viscosity and thermal conductivity, are developed. A simple model of the atmosphere is proposed and a technique, suitable for numerical analysis, is developed for determining the various modes propagating. The effect of the earth's magnetic field on the propagation of mechanical waves is evaluated and it is concluded that this effect is small. The modification of the disturbance due to the presence of charged particles is discussed in conjunction with a comparison of the microscopic and macroscopic formulation of the equations of motion. An idealized description of the source is given in the appendix and some solutions to some simple problems are discussed. iii

THE UNIVERSITY OF MICHIGAN 3673-1-F TABLE OF CONTENTS I INTRODUCTION 1 II PROPAGATION OF MECHANICAL DISTURBANCES 3 2.1 Equation of Motion 3 2. 2 The Unperturbed Atmosphere 10 2. 3 Propagation Constants 19 2. 4 Modes of Propagation 29 2.5 Mode Functions and Mode Expansions 45 2. 6 Formal Solutions of the Source Problem 55 III IONOSPHERIC DISTURBANCES 61 3. 1 Introduction 61 3.2 Effect of the Magnetic Field 61 3. 3 Extent of Travel of the Disturbance 64 3.4 HF Phase Delay 66 IV MODIFICATION OF ATMOSPHERIC DISTURBANCES DUE TO CHARGED PARTICLES 70 4. 1 Introduction 70 4. 2 Discussion of Fundamental Relations 72 V CONCLUSIONS 82 APPENDIX A: APPROXIMATE FORMULATION FOR LOSSES 83 APPENDIX B: IDEALIZED SOURCE DESCRIPTION 92 APPENDIX C: NUMERICAL ANALYSIS 99 APPENDIX D: SOME IDEALIZED SOLUTIONS 110 REFERENCES 119 -l I.. I I lVI

THE UNIVERSITY OF MICHIGAN 3673-1-F I INTRODUCTION This report summarizes the effect of ionospheric disturbances on HF radio waves propagating over long paths when the ionospheric reflection points of the HF waves are at a large distance from the source of the disturbance. In order to design detection systems for such ionospheric disturbances a realistic theoretical model capable of yielding some quantitative predictions regarding such disturbances is essential for providing system design and reliability criterion. The essential purpose of this report is to propose such a model and discuss the methods by which the numerical solution for problems involving such a model may be achieved. During the course of our research we have also investigated some idealizedproblems to emphasize the effect of the variation of the various parameters entering the problem. In order to maintain the continuity of the main context of this report, some pertinent results of such investigations are included as appendices. According to the difference in the nature of the physical principles (and their mathematical descriptions) governing the various stages of the evolution of the disturbance, we may subdivide the general problem into the following. a. The Source Problem The source of the disturbance will be assumed to be contained within a small volume moving along some path in the ionosphere. It will be further assumed that the disturbance is initiated by the release of mass and energy from this volume. The source problem is therefore that of describing the initial disturbance near the path of the source. During the initial stage of the disturbance ahighly non-equilibrium condition exists between the mass ejected by the source and the ambient atmosphere so that nonlinear effects must be taken into consideration. This initial disturbance,which lasts only a short time, is confined to the vicinity of the source. 1 1

THE UNIVERSITY OF MICHIGAN 3673-1-F b. The Propagation Problem The initial disturbance, after reaching near equilibrium with the ambient atmosphere, may be taken as the source of the long range disturbance. The dispersion and diverging of the initial disturbance into large scale weak disturbances far from the source can perhaps be adequately described by linearized equations of motion. c. Interaction with Electromagnetic Waves Problem This problem deals with the perturbation of the electron density in the ionosphere, and its effect on the forward propagation of HF waves. Due to the uncertainty of describing the energy and mass ejected from the source and the stratification of the atmosphere, adequate solutions for the initial nonlinear problem cannot be obtained easily. In Appendix B some discussions of the source are given. In the main report, however, we have chosen an idealized model, assuming that the source may be approximated by sections of horizontally expanding gases. The propagation of disturbances in a stratified atmosphere using linearized equations of motion is a well-defined problem. Restricted to the lower atmosphere with the source on the ground, such problems have been investigated by numerous authors (Weston, 1961; Pekeris, 1948) in connection with the propagation of pulses caused by meteors or explosions. In the present problem, however, we are interested in a disturbance in the ionosphere. Thus the investigation should be extended to include high altitudes with a source above ground, which complicates the problem considerably. At first glance, the propagating medium is ionized, so that the ordinary hydrodynamic equations may not be adequate. However, since the atmosphere is only weakly ionized (even in the peak of the F-layer) one may, as the first approximation, consider the medium as neutral. Some criterion as to the formulation of the problem if the effect of ionization is included is given in Section IV. The real complication comes from the fact that at high altitudes, the kinematic viscosity of the atmosphere increases almost exponentially and therefore must be included in the 2

THE UNIVERSITY OF MICHIGAN 3673-1-F formulation. An approximate formulation to incorporate the effect of viscosity and thermal conductivity in the equations of motion is given in Appendix A. Using this approximation the propagation problem, although somewhat complicated, can still be formulated in terms of the familiar mode theory. In Section II the basic problems of the propagation of disturbances in terms of mode theory including source excitation are discussed. The effect of mechanical disturbances on the electron density and the effect of HF forward propagation are discussed in Section III. In an earlier phase of this investigation some numerical calculations were carried out for a similar model by approximating the atmospheric parameters with smooth continuous functions. It was found that the first mode number could be obtained very easily. However, we were not able to calculate the second mode numbers due to the convergence problem in the numerical computations. It was later found that the continuous functions were not adequate near 100 Km where the derivative of temperature is discontinuous. In Appendix C the description of the numerical scheme used for the earlier calculation is discussed. Based on this experience the model proposed in this report has been designed specifically to avoid this difficulty. It is believed therefore that the model reported here is amenable to straightforward numeri cal computation, although actual computations were not done on this model. II PROPAGATION OF MECHANICAL DISTURBANCES 2. 1 The Equations of Motion The macroscopic description of the motion of a gas under the action of external forces may be adequately stated in the form of three conservation laws. The variables involved in these laws are: 3

THE UNIVERSITY OF MICHIGAN 3673-1-F p mass density p hydrostatic pressure T kinetic temperature P traceless pressure tensor Q heat flux v flow velocity e internal energy per unit mass h = e + enthalpy (heat energy) per unit mass P X external force per unit volume. All these quantities, in general, are functions of the space coordinate r and time t. In terms of these variables, the conservation laws may be stated as follows:,Mass conservation a p +V * (pv) = O (2.1) Momentum conservation (pv) +7V ~ (p + p v v + P)X (2.2) Energy conservation -( pv +pe)+V * p (v v+h)++v v' X (2. 3) In the above forms, the physical interpretation of the conservation laws are clear: the quantities involved in the time derivatives are the volume densities (of mass, momentum, or energy); the quantities involved in the gradient operation are the corresponding fluxes or flow rates; and the quantities which appear on the right hand side of the equations are the corresponding sources. I denotes the identify tensor. 4

THE UNIVERSITY OF MICHIGAN 3673-1-F Some alternate forms or combinations of the above laws are obtained by using the total derivative of a quantity such as mass, density, velocity, etc. in the moving fluid. The total derivative is usually denoted mathematically by the operator (- a +v.VV). Using this notation the alternate forms of the equations are: (i) the total derivative of the mass density (-a +v V)p=-pV. v (2.4) which is an alternate expression of (2. 1). (ii) the total derivative of velocity O X 1 1 (- +v'V)v - - Vp+ V' P (2.5) 5t P P P This can be obtained by taking the difference of (2. 2) and v times (2. 1). (iii) the total derivative of the internal energy+ p(t +v V)e =-V. -pV v-V v: P (2.6) This is obtained from (2. 3) by introducing the thermodynamic relation h = e + (2.7) p and equations (2. 1) and (2. 2). (iv) the total derivative of the entropy pT(- + v V)s =-V - -V v:P (2.8) where s is the entropy per unit volume. This is easily obtained by introducing the thermodynamic relation Tds = de _ dp (2.9) in (2. 6). Following the notation of Chapman and Cowling (1939), the scalar product of the two tensors is given by A: B = A.ijBij... 1J 1J _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _

THE UNIVERSITY OF MICHIGAN 3673-1-F Mathematically, these equations may be regarded as the first three moments of the Maxwell-Boltzmann equation for the microscopic description of a non-uniform gas (e. g. Hirschfelder, et al, 1954) and are, therefore, not complete. As is easily seen, the set of equations (2. 1) - (2. 3) contains five scalar equations with 14 unknowns. The set of 14 unknowns can be reduced if some local steady state properties of the gas are assumed to be known. For gases at moderate temperature, the perfect gas law is valid, i. e. p _R"' p T (2.10) p M where M is the molecular weight of the gas, and R-* = 8314 Joules/deg. Kelvin is the universal gas constant. If the heat capacity of the gas, C, is known as a function of the temperature, then de=C dT (2.11) v For the same range of variation, therefore, it is adequate to use a linear approximation and write (2. 11) as e=C T (2.12) v Equations (2. 10) and (2. 12) can be used to eliminate T through the familiar thermodynamic relation C= R*/M (2.13) v y-1 where y is the ratio of specific heats. The result is e 1 = P (2. 14) y-1 P Using this relation (2. 6) may be reduced to +Strictly speaking, the linear approximation yields (e-e) = C (T-T ) in the region 0 v of (e,T ). The linear approximation is valid in this region so thai (2. 12) is correct to within a constant. Since e appears in the equations as a differential, this constant is immaterial. i___ i 6

THE UNIVERSITY OF MICHIGAN 3673-1-F + V V) ( )+ -V p V =V (2. 15) In most problems involving weak disturbances one may assume that without the disturbances an equilibrium condition exists. If the equilibrium pressure, density and external force are denoted by p, p and X respectively, then (2. 5) yields VP0 = X (2.16) while (2. 15) yields V Q= 0. (2.17) Under the action of some perturbation force X, one may write p = p + p(2.18) P=Po 2 p=p +P' (2.19) o Using the above perturbation relations for p and p and equations (2. 1), (2. 5) and (2. 15), the following equations are obtained. at. + +v p+ V.(v (2. 20) (P +P at v+(P +')v Vv+VP+Vp +7 P =X (2.21) 0P P)V v +V/P +7 Pot0 a P +P -1)v 1) + 1 (po+P)V v+V V v P 0 (2. 22) Assuming that the perturbed quantities are much smaller than the corresponding equilibrium quantities, and neglecting second order quantities in the above equations, the following linearized equations for the perturbed quantities result. Vatp = 0 (2.23) 0t fo -~ -~'O,

THE UNIVERSITY OF MICHIGAN 3673-1-F av IPP=Xo V at +Vpv -=,x (2.24) Po at +VP -VPo+V' P =X_ aaJ Pp P +( )v- V (1 ) +'ypV * v+(y-1)V. Q = (2.25) The above set of equations cannot be solved without some knowledge of the relation of P and Q to the pressure, density and velocity. These relations may be expressed in terms of the viscosity and thermal conductivity, X, of the medium. For example in the quasi-static approximation, we may have MI p P, 2 P {v+ T - 3(V.* v) p (2. 27) | The direct introduction of these relations into (2. 23) - (2. 25) yields a set of equations which in principle can be used to solve for the perturbed quantities. These two equations are not used as such because: (1) The high order spatial derivitives introduced into (2. 26) makes the solution for the disturbance in an inhomogeneous medium extremely complicated. (2) As discussed in Appendix A, at a single freyuency in a homogeneous medium these erquations predict that the phase velocity of propagation approaches infinity as wv becomes large. This effect is not substantiated by experiment (Holway, 1963). In Appendix A a simplified, yet seemingly improved, model for the viscosity and thermal conductivity has been postulated. Essentially the following modifications were made to the loss terms: (1) The high order spatial derivatives were replaced by time derivatives obtained for a one-dimensional flow. (2) Extra factors were introduced to account for the lagging effect of the viscous flow. 8

THE UNIVERSITY OF MICHIGAN 3673-1-F The results are as follows. For a single frequency the Fourier transform of the perturbed quantities are governed by the following equations+ -ip+P * v+v _ _VP =0 (2.28) V Po -ip v + FV p -p =X(w) (2.29) p0 -ilp -+ V + G v = 0 (2. 30) where F and G are the factors which account for the losses. For non-viscous flow F = G = 1, while for viscous flow the argument given in Appendix A yields 1- 2 3 C c F (2. 31) wv 5 2 c 0 and.wv 7 1-i 22 c G = (2. 32) wv 5 1-i 2 2 c 0 Wo where c = \- is the local acoustic velocity in a lossless medium and v =/p| is the kinematic viscosity. The terms in the denominator of F and G account for the lagging effect of viscosity. In the solution of the problem involving disturbances we shall use F and G as given in (2. 31) and (2. 32) but the general mathematical formulation is valid for any forms of F and G that would adequately describe the medium. We have used the same symbol for the Fourier transform of a quantity and for the quantity itself since no confusion should result. 9

THE UNIVERSITY OF MICHIGAN 3673-1-F The discussion of the solutions of the set of equations (2. 28) - (2. 30) for the perturbation of the real atmosphere shall be given in subsequent sections of this chapter. At the present, however, it is proper to note the energy and power associated with the disturbance. From the interpretation of (2. 3) one may easily see the following: a) the mechanical energy associated with the disturbance is 1 2v (2.33) b) the mechanical power propagation is vp h - pv (2. 34) These equations indicate that in order to account for the effect of the variation of density on the energy and power propagation in an inhomogeneous medium, it is advisible to introduce the normalized quantities u=vp1/2 (2.35) ply (2.36) (y- 1)p/2 and R= (2.37) p0 These normalized quantities shall be used in the perturbation solutions of atmospheric disturbances. 2.2 The Unperturbed Atmosphere In discussing the propagation of mechanical disturbances in the real atmosphere, some properties of the ambient atmosphere must be known. Due to the fact that the atmosphere has, in general, a very complicated structure, only some average properties are tabulated in the literature (U.S. Standard Atmosphere, 1962). Since the atmosphere is always under the action of the gravitational field, 10

THE UNIVERSITY OF MICHIGAN 3673-1-F the unperturbed external force is given by X = -gp z (2.38) where the plane earth model is used and the ground plane is at z = 0. The vertical upward direction is taken as the positive z-axis, and g is the gravitational acceleration. The equilibrium conditions for the ambient atmosphere are given, therefore, by(2.2)and (2. 15) with v = 0. These are dPo 0dz Pg (2.39) dQo dz The perfect gas law may be assumed to hold, so that Po R* _ = -T. (2.41) p M o Due to the dissociation of molecules at high altitudes, the molecular weight M is not constant but a function of z. In general (2. 38) and (2. 40) permit the construction of a model of the atmosphere from a basic chosen set of values of M and T. If we denote Po A(z) - - T= (2.42) then the ambient pressure and density of the atmosphere are expressed by po(z) =p(0) exp AE dz (2. 43) and p(z)p() exp A dz] (2. 44) o A LoA 11

THE UNIVERSITY OF MICHIGAN 3673-1-F Values of po(z) and p (z) versus altitude are given in the United States Standard Atmosphere (1962). From (2. 42), (2. 43) and (2.44) one can easily see that the inhomogeneities of the atmosphere may be expressed in terms of A(z) and its derivative d (A(z)) and g. The following relations are often useful and can be derived from (2.43) and (2. 44) =-z g (2. 45) P0 VP 0 AA'+ g VP0 A ATtg (2.46) p A VPo g (2.47) - zg (2.47) p A 0 The other properties of the atmosphere involved in the equations describing the disturbances are: the ratio of specific heats y, the viscosity /u and the thermal conductivity X. The variation of y and the kinematic viscosity v=//Po have been calculated for the standard atmosphere. Since the thermal conductivity and kinematic viscosity are approximately linearly related, a separate calculation is not necessary for the thermal conductivity. The variation of A(-po/Po), the kinematic viscosity v, the ratio of the specific heats %y, and the gravitational acceleration g vs altitude are given respectively in Figs. 2-1 through 2-4 inclusive. From these figures, some genera features of the atmospheric properties can be observed. These are: (1) The effect of viscosity is practically negligible at low altitudes, perhaps < 100 KM. However, at high altitudes, especially in the ionospheric region, viscosity should play an important role in the propagation of mechanical disturbances except at very low frequencies. (2) Except at low altitudes, where the value A' shows some change of sign (corresponding to the temperature inversion in that region), the variation of 12

THE UNIVERSITY OF MICHIGAN 3673-1-F the temperature and therefore A' is quite smooth. Thus in the ionized region of the atmosphere, A may be approximated by a smooth curve. (3) The variation of g and y, in general, are very small, and as a first approximation, these quantities can be regarded as constant. For convenience in numerical work and to provide a guide for the range of values chosen for computation and analysis, the pertinent parameters A and v of the standard atmosphere have been approximated by simple piece-wise continuous curves. The approximate analytic expressions for these parameters are given in Table II-1. The values of the parameters computed from these expressions are also plotted in Figures 2-1 and 2-2 for comparison with the original data. It is apparent from these figures that a reasonable fit has been obtained. 13

THE UNIVERSITY OF MICHIGAN 3673-1 -F TABLE II-1 ANALYTIC REPRESENTATION OF ATMOSPHERIC PARAMETERS Kinematic Viscosity v (m sec ) Altitude Range Km -5.46 x 10-5 e14.7x0 z 0 - 135 -5 0. 662 e 135 - 190 -5 3 2.17 x 10 Z 4.24 x 10 e 190 - 400 -5 2.51 05 1.15x10 z 400- 700 2 -2 A( m sec ) Altitude Range Km 6.4 x 104 0 - 90 3.14z - 2.187 x 105 90 - 205.945z+ 2. 314 x 105 205 - 480.445z + 4. 665 x 105 480 - 700 14

THE UNIVERSITY OF MICHIGAN 3673 -1 -F 10 9 8 7:-I - 100 200 300 400 500 600 700 800 Altitude (Km) FIG. 2-1: A = p /p VERSUS ALTITUDE 15

THE UNIVERSITY OF MICHIGAN 36 73 -1 -F 9 10 8 10 - 7 10 6 5 -; - 10 4 3 0 - i 2 10 From Tables 10 10 oo -t -5 0 100 200 300 400 500 600 700 Altitude - Km FIG. 2-2. KINEMATIC VISCOSITY, v, VS ALTITUDE 16

THE UNIVERSITY OF MICHIGAN 3673 -1 -F From Table 1.8 1.6 ci 1.4 1.2 o 0.8. 0 o 0.6 0.4 100 t 300 400 500 600 700 Altitude - Km FIG. 2-3. RATIO OF SPECIFIC HEATS, I, VS ALTITUDE 17

THE UNIVERSITY OF MICHIGAN 3673-1 -F From Tables 10 9 cr 8 7 Ct D 6 0o = 4 0 o 3 2 0 100 200 300 400 500 600 700 Altitude - Km FIG. 2-4. ACCELERATION DUE TO GRAVITY, g, VS ALTITUDE

THE UNIVERSITY OF MICHIGAN 3673-1-F 2. 3 Propagation Constants The study of the effect of gravity, inhomogeneity and viscosity on the propagation of a disturbance in a real atmosphere is indeed a very complicated problem. The set of equations (2. 28), (2.29) and (2. 30) in general, can only be solved numerically. However some basic features of wave propagation in the atmosphere can be deduced by considering the single frequency propagation constant of a locally homogeneous medium. Hines (1960) has studied this type of propagation constant for the lossless case. In this section the propagation constant in a real atmosphere, including an approximate viscous term, shall be given. At a single angular frequency w the equations governing the propagation of small disturbances in the atmosphere may be obtained from (2. 75) - (2. 77) by introducing the relations: VPo A -z g p0 VPo = _ A +g p A The resulting equations are -iwp+poV 7 v- v P (2. A48) -iwp v + FV p~+ g= 0 (2.49) -iw-z vpg+ P GV v:=0 (2.50) where, as illustrated in Appendix A, 19

THE UNIVERSITY OF MICHIGAN 36 73 -1 -F F- 5 (2.51) I-13 2 c 0o and 7 wv 2 2 G 5 wvC(2.52) 2 2 c 0 are the terms containing the effect of viscosity and thermal conductivity, respectively. In accordance with the idea of normalization discussed in Section 2. 1, the following change of variables will be made 1/2 1/2 1/2; V O The set of equations (2.48) - (2. 50) then become A At -o -iwR+ V* V-z V 2A (2.53) -iw V+FVS- 2 z S+zgR 0= (2.54) -iw S +AGV- V-.V[g G(+g)=0 (2.55) The local plane wave propagation constant k is obtained by assuming that all the quantities are equal to an amplitude times a propagation factor exp(ik * r). Direct 20

THE UNIVERSITY OF MICHIGAN, 3673-1-F substitution into (2. 53) - (2. 55) yields the following set of homogeneous equations for the amplitude. A A'+g = 0 -iwR+k * V-z V 2A (2.56) 2A -iwcV+ik FS (2.57) -iwS+ i AG k. V-/. V -/G2 g)]=0. (2. 58) Here we have used the same symbols to denote the amplitude functions, since no confusion could arise. From the above equations it is evident that the vectors V, k, and z are coplanar. A A Thus, without loss of generality, the vectors V and k may be assumed to have x and z components only. In terms of these components, one may rewrite (2. 56) through (2. 58) as -iwR+ik V +ik V A+g V =0 (2.59) x x z z 2A z -iwV +ik FS = 0 (2.60) x x A'+G -icoV +ik FS- A S+gR=0 (2.61) z z 2A -iwS+iyAG(k V +k V )-V g G 2 A + = 0 (2.62) In order that the above set of homogeneous equations yield non-trivial solution for the amplitudes, the following dispersion relation must be satisfied. 21

THE UNIVERSITY OF MICHIGAN 3673-1 -F 42 G(A' +g) 2 2g vG(A?+gF]+gF[G(A?+gjkk2 4A + -AGF(k +k )-ik (F-i 4A x z z 2 x (2.63) If (2. 63) is satisfied, then the amplitudes, except for a coynmon proportionality constant, may be expressed as: S = -YAGk +i-(g G ) (2.64) z 2 R w=2 k + ik2 F(g-yG(A'+g)) +i A+g (2.65) V =w(w -yAGFk ) (2.66) z x V =wk F'vAGk +i ( G(A+g) (2. 67) x z 2 The dispersion relation (2.63) illustrates the local propagation characteristics of a stratified medium. Since both F and G are complex functions of t, the relation between the frequency and the propagation constant is, in general, quite involved. To illustrate the effect of the inhomogeneity, gravity and viscosity on the propagation constant, several special cases will be considered. For a homogeneous lossless medium with no gravitational force, F = G = 1, A' =g = 0, so that (2. 63) reduces to 4 2 2 wt -2o yAk = 0 (2.68) This yields 2 2 2 t t k =' = (2.69),yA 2 22

THE UNIVERSITY OF MICHIGAN 3673 -1 -F where c is the local speed of sound. This of course is the ordinary acoustic wave. For a lossless medium with a constant speed of sound, (2. 63) reduces to 4 2 "1-yg 2 2 2 Lo - 4 jA +/A k +g (y-l)k = 0. (2.70 This equation indicates that the medium exhibits anisotropic characteristics in the propagation of disturbances. For it we write k = k sin (2. 71 x where 0 is the angle between the direction of phase propagation and the vertical, the propagation constant k becomes W2( 2 2 2 a k 2(2.72)'yA(wi -w sin 0) g where for simplicity, we denote 2 2 22 =i gz -.2 (2. 73) a 4A 2 4c o0 and 2 2 2i ('y <( 2 (2. 74) g yA 2 a o This familiar anisotropic characteristic of plane wave propagation has been discussed in detail by Hines (1960). Briefly, depending on frequency,one may have the following modes of propagation. (i) For L> w 23

THE UNIVERSITY OF MICHIGAN 3673-1 -F k is real for any direction, i. e. plane wave propagation is possible in any direction. For w>> > the anisotropic character becomes much less prominent, and the disturbance behaves like an ordinary acoustic wave. Thus, such a mode is known as an acoustic mode. (ii) For w<wg g Plane wave propagation is possible in a limited direction given by ]sinoj >-. Such a mode is known as a gravity mode. (iii) In the frequency range w > w > wgo k is imaginary, so that unattenuated plane wave propagation is not possible. Such a frequency range is sometimes known as the forbidden region. In Fig. 2-5, oa and wo versus altitude are given for the values of A from the model atmosphere of Fig. 2-1 and Table II-1. The effect of the inhomogeneity, i. e. the variation of A, can also be roughly seen from (2. 63). Neglecting the gravitational effect, and letting F=G=1, we have i -t L- 4A /A = 0 so that 2 MA')2 2 4A k ='yA 2 (2. 75) Thus for w >( ), k is real and propagation is possible. A cutoff frequency 4A occurs at 24

gF Mr~nGA N rr 4 E'3673 -1 - 3.0.66 j.4 A.0 fZ4 0.6 0. ~2 iW?. 3AlitdeORACUSICAN 0 A0 GUAI,,FREUENIES pJE-R ~~~.O I~ —,,-v A - V TE TMS

THE UNIVERSITY OF MICHIGAN 3673-1 -F ( 2 1/2 Thus for w = 10.2-. (2. 76) The maximum temperature gradient occurs at an altitude of 120 Km and has a value of 200K/ Km, while the temperature at 120 Km is approximately 3500 K (U. S. Standard Atmosphere, 1962). For these values wto'10 radians/ sec and thus only very low frequency gravity waves are affected. The effect of the inhomogeneity when the gravitational force is considered is roughly reflected in the change of Lw and to. It is easy to see that if A is not constant, then the following modification of w and w should be made. a g 2 y/(A' + g) 2 A'2 | WI (4 A = w= (1+-) (2. 77) ~~a 4A a g and A' g 1 A + — )[g2yA g L'-' (2. 78) The effect of the inhomogeneity can be seen from the expression for t' and a'. In regions of negative temperature gradient both w' and t' are less than g a g t and wt respectively. In addition, the change in t is somewhat greater than int. Thus in regions of negative temperature gradient the forbidden band is both lowered and narrowed as compared with a homogeneous atmosphere. In regions of positive temperature gradient the bandwidth of the stop band is also narrowed but shifted to a slightly higher frequency. This effect is indicated by the dashed line in Fig. 2-5. 26

THE UNIVERSITY OF MICHIGAN 36 73 -1 -F The curves in this figure are based on the model atmosphere described by the values of A and A' given in Fig. 2-1 and Table II-1. As a result of the large discontinuity in A' at 100 Km the curves of to' and to' have not been drawn for a g altitudes from 0 to 100 Km. The effect of losses on the propagation constant, due to the involved functional relation of F, G andw is somewhat more involved. For a homogeneous medium which neglects the effect of gravity, one has 2 2 Lw -yAk FG=0 or v 5 iv 5 1 -- -.) ( -i- 23 22 cl /-c c k Li c/ HG Li-0 0 7. (2. 79) 2 2 2 c c o o0. Thus, the wave is always attenuated. For wv / c < < 1, one has 0 k=- + 2 i ( i + c2 ) (2. 80) o o The two terms in the imaginary component of k, corresponding to attenuation, are due to the effects of viscosity and thermal conductivity respectively. For 2 tov / c >> 1, one has k 1- -L; =. 62 ~ c 97 c O o As mentioned in Appendix A, this agrees roughly with experimental results. In principle, the effect of losses on the propagation constant when the gravitational 27

THE UNIVERSITY OF MICHIGAN 3673 -1 -F field is not neglected, can be calculated from the dispersion relation. However, due to the complicated functional form involved, it is not fruitful to discuss the effect in great detail. To show some effect of the viscosity, one may consider an ideal case of small viscosity and introduce the forms F = 1-i F, and G =1 - i a nG, where fF and AG are positive small quantities. For constant A, the propagation constant in the x(horizontal)direction (k = 0) can be obtained from (2. 63)as z 2(2 2 L2 F) 2 a k2 a a (2. 81) x2 y 2 2 22 ~/A -w - i(AF+ AG)(_ -I g 73-1 Due to the presence of the loss terms, k can never be zero or infinity. As the x frequency changes, however, k2 as a complex number changes both in magnitude and x argument. For instance, consider AF and AG to be of the same order of magnitude. Then the variation of the argument of k2 may be roughly stated as follows: x 2 2 2 For ow > 2, k as a complex number, is in the first quadrant, while for a x 2 2 o < t, it is in the second quadrant. For these ranges of w, both the real and g imaginary component of k are positive corresponding to an attenuated wave. For x the other range of frequency, the value of k calculated is physically unrealistic. This seems to indicate the broadening of the forbidden region as in the case of lossless atmosphere. The local propagation constant for the disturbance described in this section can be useful in performing ray tracing of the acoustic waves. The concept of the for bidden region for plane wave type of propagation, however, should not be strictly 28

THE UNIVERSITY OF MICHIGAN i 36 73 -1 -F adhered to, especially in a stratified medium with sources. In these cases, there may exist disturbances in the form of inhomogeneous plane waves, such as the surace waves well known in electromagnetic field theory. In general, if some source initiates a disturbance in the atmosphere, it may propagate in any direction. As the ave propagates in a stratified atmosphere, k is fixed and should remain fixed as he disturbance propagates. The value of kz however, should vary with altitude in a stratified atmosphere. At any altitude, the value of k can be estimated approximately from the dispersion relation appropriate to that altitude. For some alue of k, one may then find values of k to be locally imaginary even in the case of X z lossless medium. For such inhomogeneous waves, the average power propagation an be shown to be still in the x-direction and all amplitudes decay exponentially with altitude. Such a mechanism of propagation in the real atmosphere is possible and can best be considered directly from the solution of the differential equation..4 Modes of Propagation At a single frequency, the linearized differential equations of motion for a disturbance in the atmosphere may be written as (from equations 2.28, 2.29 and 2. 30) dp -iw+ pV *v+v -d 0 (2. 82) Po' v- z dz -iwpov+F V'+z pg =X(W) (2. 83) -iw -v P g+ yPo GV7 ~ v =0, (2.84) where F and G are the factors which account approximately for the effect of losses and X (w) is the Fourier component of the disturbing source. The solution of this 29

THE UNIVERSITY OF MICHIGAN 3673-1 -F set of equations must be obtained numerically,due to the variation of the parameters in the real atmosphere. Since the set of equations is linear, the complete solution of the system can be expressed in terms of a 6 -function source, i. e. the use of a Green's function. As is well known in mathematical physics, the Green's function may be expressed in terms of the source free solutions with appropriate boundary conditions. Most of the work done on disturbances in the atmosphere, therefore, has been concentrated on obtaining a solution of the homogeneous (source free)part of the equations. The general solutions of the homogeneous part of the above equations can, in turn, be obtained from the solution of a simpler problem, i. e. a problem with cylindrical symmetry. For such problems, cylinderical coordinates may be introduced, and thus each of the variables are functions of z, and r. By applying a Hankel transform to each of these equations, the above system may be reduced to a system of ordinary differential equations involving a radial propagation constant k. For a lossless atmosphere, such a reduction has been done consistently. Extension to a lossy medium, using the model proposed in this report, is somewhat lengthy but straightforward. A brief description of the steps of such reduction is given in the following paragraphs. The z-component and the divergence of (2. 83) are -iop v + F+ +Pg = 0 (2. 85) -,k0 Z See Section 2.6 for the solution of a more general case. 30

THE UNIVERSITY OF MICHIGAN 3673-1 -F dp -iwp v-iw v + FvV + dz v 2 (2. 86) 0- dz z a z dz a z Due to cylindrical symmetry, one may use the Hankel transforms of v, |p, p and V' v defined by z'M oD v' v | 0 =- rJ O(kr),dr P p The transformed variables then satisfy dp -iwD+po X+ dz W = (2.87) -iwp W+F dz + Dg = O (2. 88) -i~po dz dp0 d2P 2 dP dF aD -io -iw W+F 2 -Fk P+ +g z= 0 (2. 89) dz ~dz dz dz a -iwP -Wp g+' Po GX= 0. (2.90) Elimination of D and 2 in the above system yields i dP gpo 2 iwoF +Lz J 0 G =0 (2. 91) F w p FG Pw2 g G 0O0 2p dW P 0F2og 2 iPo dz VpG p is o b GP G (2.P92) It is to be noted that the parameter k used here should be identified with k x discussed in Section 2. 3. 31

THE UNIVERSITY OF MICHIGAN 3673 -1 -F These coupled equations can be used to eliminate one of the variables, yielding an ordinary second order differential equation. In accordance with the idea of normalization introduced in Section 2. 1, we may define the new variables U W 1/2 > P -1/2 so that from energy considerations, I 12 and f U2 should remain bounded everywhere in the atmosphere. The equationsfor ~ and U, then, are i wF (+ ( + yAFG) i+U) A9 A]= 1=0 (2. 93) Ldz 2A yAFG A yAG' iLdU(A+g U+ + - -kF = 0. (2. 94) 2A yAG'yAG These equations may be reduced to a self-adjoint form by a change of variables z exp 52 AFG (1 - F)dz (2. 95) 0 z U = U exp 2 (1 - F)dz, (2. 96) 1 j2,yAFG 0 and result in the following form: iW dt 1V 1+ U f 0 (2. 97) id U - VU1] + W2f 0 (2.98) 32

THE UNIVERSITY OF MICHIGAN 3673-1 -F where 2A 2 yAFG.99) 2 2 ow g(A'+g) + (2. 100) $ F A F yAGF and 1 k2 f = - F. (2.101) u yAG 2 Due to the introduction of losses, all the functions V, fu, f are complex and therefore do not vanish for any real value of z. Moreover, from the approximate expressions for F and G given in (A. 29) and (A. 30) as z->oo, the real part of the exponential in (2. 95) and (2. 96) is negative, so that the boundedness of ~ and U implies also the boundedness of ~1 and U1. A second order differential equation for (1 and U1 can then be easily obtained. The result is dz 2 d l I d + d V V - + f 0 (2. 102) dz f C dz dz f + f0 d 1 d d V V+ U + U + f 0. (2.103) dz f dz + U1 dz f =0. (2. 103) [fu daas U e] jdeU U f The above equations can also be reduced to the normal form, i. e. 33

THE UNIVERSITY OF MICHIGAN 3673 -1 -F d 2 d _ d1 1 d V V 3 dj dz f d dz2 dz f fp u 4 3 2 f2 (2.104) and 2d 2 2 2 U fUfU dV2 (_ f f __ 1 d V V 3 dz U 1 dz f d fF+ L f f2 0 2 f = (2.105) dz2 u- dZ f f f ) 4 f3 2 U U U Alternately, it is sometimes more convenient to obtain a differential equation in X? when solving some idealized problems such as was done by Pekeris (1948). This is accomplished by eliminating R and D and yields A-+X FG- = 0 (2. 106) AGI d k2 2 2 dz 2 AA [2 k2F Gg-(F-g-F AG- F(Ag) + l Wk2 2 (A' +) 0. (2.107) The use of these equations to obtain a formal solution of the equations of motion in a lossless atmosphere with a linear variation of A has been given by Pekeris (1948). To solve any set of these equations, boundary conditions are necessary. In order to satisfythese boundary conditions, it is generally found that non-trivial solutions of the homogeneous equation are possible only if the radial propagation 34

THE UNIVERSITY OF MICHIGAN 3673 -1 -F constant k, assumes a particular set of values. Solutions, corresponding to these values of k are then known as the different'modes' of propagation. The physical picture of the modes of propagation has been discussed in great detail by Budden (1962). In general, for each k, the corresponding variable has an r dependence for large r of the form 1 ikr J (kr)~ e Thus, the k's may be interpreted as radial propagation constants. For each value of k, the vertical propagation constant changes due to reflection at the boundary and refraction of the waves due to stratification of the medium. For certain values of k, the overall reflected and refracted waves corresponding to different vertical propagation constants, starting from the source, are in phase (or differ in phase by an integral multiple of 2r) in the far region so as to reinforce each other. For such values of k, then, the disturbances can be propagated from the source with small attenuation. For other values of k the waves are not in phase, and thus on the average, no net energy propagation is possible. In a lossy medium, where attenuation is inherent, the first mode, corresponding to the waves reflected (or refracted) only once, should have the least attenuation. The higher order modes, corresponding to the waves reflected and refracted several times, should have relatively larger attenuation. For the propagation of VLF radio waves in the atmosphere the mode calculations of Barron (1959) clearly demonstrate this fact. In our approximate formulation of the equations of motion, losses have been introduced. It is 35

THE UNIVERSITY OF MICHIGAN 3673-1-F therefore expected that only a very few number of modes will be necessary to describe the dominant solutions for disturbances far from the source. The boundary condition to be satisfied at the ground, z = 0 is given by W = 0 since the ground is rigid. The boundary conditions to be satisfied at higher altitudes is generally uncertain. In common practice, appropriate to the propagation of disturbances in the lower atmosphere, we may assume that the atmosphere is bounded at a plane at some arbitrary height z2. At z2 for the lossless case, one may postulate the following different conditions: (1) P = 0, corresponding to a free boundary. (2) W =, corresponding to a rigid boundary. (3) One may assume that beyond z = z2 the atmosphere is isothermal (A = constant) so that all the amplitude quantities should vary as 2 2 2 2 1 -W Li -Li a_ ~ ~2 e w YA 2 with the values of a, a) and y/A evaluated at z = z2 a g 2 Published results in the evaluation of the mode number k, limited to the lower atmosphere where the effect of viscosity can be neglected, generally fall into one of the following two categories. (i) The atmosphere is stratified into a very few layers, each idealized by assuming either a constant A (isothermal) or linear variation of A with altitude. The general solution for a lossless medium with constant A is given in Appendix D. 36

THE'UNIVERSITY OF MICHIGAN 3673 -1 -F For these layers, the general solution of the differential equations are known formally (given briefly in Appendix D). These solutions involve some function of k and z and arbitrary constants. By matching the boundary conditions at z = 0, and z = Z2' and by using the continuity of W and P at the interfaces between layers, these constants may be eliminated, resulting in an equation for determining k. Some of the results of such a calculation (limited to z21 0 Km) have been given by Pekeris (1948). (ii) The atmosphere between z = 0 and z = z2 is stratified into a multitude of layers each assumed to have constant properties. This stratification reduces equations such as (D. 26) and (D. 27) (in Appendix D) into difference equations relating the values of P and W of the i'th layer to those of the ( i+ l)'th layer. Thus, values of P and W at the ground can then be related to the values of P and W at z = z2 by a matrix whose elements are functions of the properties of each layer and k. Inserting appropriate boundary conditions at z = 0 and z = z2, the problem is then reduced to the condition that the determinant of some partition of the above matrix should vanish. The vanishing of the determinant then is the condition that determines the mode numbers k. Limiting the values of z2 to 108 Km and 220 Km and with the n 2 three types of boundary conditions mentioned above, Press and Harkrider (1962) have determined several values of k, or rather w/ k by computer techniques using 39 layers. Their results are quite valuable in explaining some basic features of the propagation of disturbances in the atmosphere but also point to the one important fact. The boundary condition at the upper boundary, z=z.................. 37

THE UNIVERSITY OF MICHIGAN 36 73 -1 -F may significantly affect the value of the mode numbers. Since the investigation in this report is concerned with disturbances in the F region, which probably extends to 300 Km, our choice of z2 should be at least 400 Km or possibly more. In this region, it is believed that the effect of viscosity should not be neglected. With the approximate formulation of viscosity included in the equations, it is felt that for the region z > z2 due to the rarefication of the gas, the propagation features and therefore the boundary condition at z = z2 should be inferred from the effect of viscosity. Unfortunately, due to the complicated variation of viscosity with altitude, no exact solutions are possible. In Appendix D, approximate solutions for a region with very high viscosity and nearly constant temperature are given. The results are' V i b, f Il......+2 dz (2.108) iT f0+1 d fVf 2 dz z z2 U1 iwf iU- f ~ iwf0 (2.109) (P V ilbf+2 dz f -Vz2 +One of the basic reasons for introducing the viscosity into our formulation is the hope of approximating the boundary condition at z = z2 more realistically. 38

THE UNIVERSITY OF MICHIGAN 3673-1-F These relations, therefore, are used as the boundary conditions at the upper boundary z = z2. Since the tabulated values of the atmospheric properties are available up to 700 Km, we may take z2 as 700 Km. The boundary condition at z= 0, because of the vanishing of the vertical component of velocity, is given by (2. 97) and (2. 98), as. =-V (2.110) and U1 = * (2.111) Due to the introduction of losses and the extension of the upper boundary to higher altitudes, it is our feeling that the use of a matrix formulation in the calculation of the mode numbers k may perhaps become involved and tedious. A method of numerical integration and reiteration, similar to those used successfully in calculating the mode numbers for electromagnetic wave propagation in the earthionosphere wave guide (Barron, 1959) seems therefore more convenient. A mathematical description of this scheme applied to the present problem is given below. The differential equations for ~1 and U1 given by (2.97) and (2. 98) are rewritten as d 1 dz 17 + V1i 1W 1 0 (2. 112) 39

THE UNIVERSITY OF MICHIGAN 3673 -1 -F dzd U1 -VU1 iW mi fue (2.113) From the expressions for V, f and f, we see that at a fixed frequency w, V is a function of z, f is a function of z, but f is a function of z and k = k/w. Thus, u 1 and U1 are functions of z and k. For fixed values of w and k, in principle,the solutions of (2.112) and (2.113) can be accomplished by numerical integration. But since these equations are homogeneous, considerable labor can be saved be solving for some ratios between the dependent variables, _ 4~~11 M(w, k, z) (2.114) _ 1 N(w, k, z) = (2.115) ~~~1~~ In terms of these variables (2. 112) and (2. 113) are reduced to M(w, k, z)+ Vc(w, z) = i- N(w, k, z)f (w, z) (2.116) and N(w, k, z)M(w, k, z)+dz N(w, k, z) -V(w, z)N(w, k, z)=iwf (w, k, z). dz u (2.117) If one considers (2. 112) and (2. 113) as generalized telegraphic equations, these ratios may be interpreted respectively as a reflection coefficient and an impedance, 40

THE UNIVERSITY OF MICHIGAN 3673-1 -F For fixed values of to and k therefore, we have only one nonlinear first order differential equation to solve. This equation is obtained by eliminating M from (2. 116) and (2. 117) to yield d 1 2 d N(w, k, z) = 2VN -i- f N + i f.(2.118) Suppose now that for a fixed to, we would like to find one value of k such that the solutions of the homogeneous equations (2. 112) and (2. 113) are possible in a region co > z > 0. From the approximate solutions of the equations for z> z2, we may obtain a boundary condition for N at z = z2 (eq. 2. 1 09), represented by a function N(w, k, z2) =N2(w k, z2) (2.119) Similarly, we may also choose a lower boundary, zl > 0 such that at z = zl N(w, k, z )=N(w k, z) (2.120) Both Nl(w, k, Zl) and N2 (w, k, z2) are known functions of k. By choosing a value of k = k, one may then integrate (2.118) numerically from z = z2 to z=z l The result of such integration yields a number N(w, ka, z ). If k is a characteristic value (or mode value) then the number N(w, ka z ) should be equal to Nl(w, k, zl) a' 1 If they are not equal then the numerical integration should be repeated for another choice of the initial value of k, kb, instead of k. A reiterative process to determine the correct choice of kb may be obtained by studying the effect of the Actually this represents two coupled equations for the real and imaginary parts of N. 41

THE UNIVERSITY OF MICHIGAN 3673 -1 -F initial choice of k on the values of N(w, k, z) calculated. If we differentiate a a (2. 118) with respect to k, and call aN(u, k, z) = J(w, k, z), the following differential equation for J is obtained. af d - U d J(w, k, z)= 2VJ-2iwof NJ+ iw ok (2.121) dz From the initial assumption of k, one may determine J(wt, ka z2) from (2.119). Thus by numerical integration of (2.121) we may obtain a value of J(wo, k Z). Thus, if for the new value of kb, N(t0, kb, Z) N1(to, ka, z1) then, by Newton's method we have kb a k (W k, z,) z N(W, k z Ca Nl(W. kP zl)-J(w, k, Z1S (2. 112) This reiterative process can be carried out until a sufficiently accurate value of k is obtained. To apply this procedure to the present problem, it seems appropriate to choose z2 = 700 Km, or a lesser value, perhaps 600 Km, so that N2(w, k, z2) is given by (2. 109), and, by differentiation, - aN J(o, k, z2) =-N (, k[, z2) 42

THE UNIVERSITY OF MICHIGAN 3673-1 -F For the lower boundary, one may choose z =0 so that by (2. 11), Nl(w, k, z )= 0. However, from our experience in a preliminary calculation on a similar model as described in Appendix C,it is felt that considerable improvement in accuracy could be achieved if we choose zl = 90 Km. This is because at z = 90 Km, the derivative of A changes very rapidly, and since both (2. 117) and (2. 121) involve the second derivative of A, the numerical integration may depend critically on the step size chosen. By choosing z1 =90 Km, this apparent'discontinuity' is avoided. Moreover, it is felt that the variation of A below 90 Km is slight in comparison with the large variation above 90 Km, and therefore we may assume A is approximately constant elow 90 Km. The general solution for the region between the ground and 90 Km, here viscosity is negligible, is given in Appendix D. From the results given in equation (D. 12) and (D. 13) we see that by choosing z1 = 90 Km., N( a k, )= i [wA 2] 1 t'coss X sinsz1] (2.123) | I 2 tO -k 11(W' k, Z1) cos s z2 W -W CO Z a g sin s z(2.123) and N1 (W, k, z1) -iw2k L cossz a sin s z1 (2.124) co sz v- /A A question may arise as to the identification of the order of the mode corresponding to any mode value k. This can be accomplished by appealing to the 43

THE UNIVERSITY OF MICHIGAN 3673 -1 -F phase integral method. For example, corresponding to any mode value ki = o ki, the propagation constant in the z-direction is given by k =V - f0 f, (2.125) where the imaginary part of the square root is chosen to be positive. At any height z, the phase constant is given by 3= Rek k /3z (, kn, z). Corresponding to each mode value k, there would be a height z due to the variation r of the parameters with altitude, for which / = O. This corresponds roughly to the height where a ray is reflected. By the phase integral method, the total phase shift for a ray starting from the ground and then reflected back to the ground is z 2 S /3 (o, k., z) dz. From the physical picture of mode theory, since U1 0 on the ground, one must have, for the first mode k1, z 2 5 /z(to, k z) dz r a, (2.126) nd for the jth mode, 2 S (to, k., z) dz' (2j+1)i. (2.127) 44

THE UNIVERSITY OF MICHIGAN 3673-1 -F Equation (2. 127) may also be used as a guide to the initial choice of k for the reiterative scheme. One may choose a set of values k, calculate the correponding values of zr, and obtain )1 by numerical integration of (2. 127). Then the value of k that yields a value of 0 closest to 7r should be used as the initial value of Re k in the reiteration scheme. Of course, this should be used only as a rough estimate, since the inclusion of losses renders the strict application of the phase integral method inadequate. 2.5 Mode Functions and Mode Expansions In the preceding Section 2.4, a physical explanation for the'modes' of propagation of a disturbance was discussed and a reiterative numerical method for computing the characteristic values ki for each mode was proposed. On physical grounds, we expect that due to the introduction of losses, only the first few modes, corresponding to lateral propagation with small attenuation, will be significant in the solution of the problem. In this section, the formal solution of the inhomogeneous problem due to a 6-function source with cylindrical symmetry will be formulated in terms of the mode functions. Assuming that the perturbing force has cylindrical symmetry, and denoting its Fourier transform by X(w, r, z) 9 X(t, r, z)e dt (2.128) -00 45

THE UNIVERSITY OF MI-CHIGAN 36 73 -1 -F the inhomogeneous equations of motion given by (2. 82), (2. 83) and (2. 84) can be reduced to the following inhomogeneous equations (corresponding to (2.102) and (2.103)): d 1 d dr 1 2 dzf 1 a1d'Z f,)+~+.G -k F =c (2.129) dz (f dz @ 1) bLhdz ( fp ) f'yAG 1 d 1 d d V +f dz (fl dz U1 )d+U1 V 2' (2. 130) dz ui 1L —(f 2 +fj 2' U u U where the excitation functions C1 and C2 are related to the perturbing force X (w, r, z). y a straightforward algebraic manipulation, it can be shown that if the Hankel transform of X (wo, r, z) and V ~ (wo, r, z) are given by 00 Sl(t, k, z)= rdr X (w, r, z)J (kr) (2.131) 1 0 z 0 nd O( S2(to, k, z) = rdr V ~X (L, r, z) J (kr) (2.132) 2 0 0 then 1 H(z) dS -, z 1/2 f2 +1/2 2 dz C1= d~. -V]0 2 +S ) ) 1 2 H dS1 (2.133) PO P0 and 46

THE UNIVERSITY OF MICHIGAN 3673-1 -F:~ ____ dS SH (z) C L =.+ -d (S + i 1/2 (2.134) C2it Ldz'JL f1/2 2dz j p LfU P c where, for simplicity, we denote z H(z) = exp 2yAFG (-F)dz. o From the definition of S2 and S1, it is evident that C1 may be considered as the Hankel transform of X H(z) dX z H(z) C(w., r, z) ) -V)( + V ) (2.135) dz P1/2f Li2 1/2 dz 0 0 while C2 cannot be interpreted as such due to the implicit dependence of f on k. For convenience in interpreting results, therefore, we shall concentrate on the solution of (2. 129). The homogeneous part of (2. 129) has been discussed in Section 2.4. Corresponding to each characteristic value k., one may find a solution P. (the 1 1 corresponding characteristic function) satisfying 2 (d V V 1 -2 __ -+ P.-k F)= O (2.136) | dz Lf~ dzd' i dz fp fp'yAG i The value of P.(w, k., z), corresponding to a characteristic value k,, can be ob1 1 1 ained by a simple numerical integration. In the reiterative procedure for the alculation of the characteristic value k., we have calculated N(w, ki, z) so that by 47

THE UNIVERSITY OF MICHIGAN 3673-1 -F (2.116), the value of M(, k' z) i -N(W, k, z)f (, z) -V(w, z) can be easily evaluated. But by definition d P.(w, k., z) =P.(w, k., z)M(w, k., z), dz 1 1 1 1 so that the functions P., corresponding to each k., can be calculated. For 1 1 z x Z1 -90 Kmin, the functions P. have been obtained analytically in Appendix D. From (D. 11) we have for z z1 = 90Km, 2 2 Ud -W a s cos sz - sin sz P.(w, k, z) = A (2.137) -k) yA where 2 2 a 2 2 2 S a (- W )k2 (2.138) YA g Due to the self-adjoint form of (2. 129) it is easy to show that 00 1 - d 1 d -P(w,.i Z) Pi(w>, ki, z) f- P.(w,k z) Pj(w, Z) ff~ dz if I i' dzj 2 2 oo (k. ) F PiPj dz (2.139) 1 j 1 j 0 Due to the homogeneous boundary conditions satisfied by P. and P. at z = co and z = 0, we have formally 48

THE UNIVERSITY OF MICHIGAN 3673-1-F 0 For i = j, we shall denote FP. dz = A (2.141) 0 To calculate A. we rewrite (2.139) as 00 00 (k2 2 FP dz P.M. -M] O (2.142) iki f j where M is defined in Section 2.4 as dP.(w, k, z) Mi(w, k. z) = 1 Pi(W, k., z) 1 1' dz 1 1 Due to the introduction of losses, P. and P. decay exponentially as z o0, and 1 thus 0 FP P.odz) P.(w, k.i 0)Pj(Wo, kk. O) k-2 M2 j o k. -k 1 i (2.143) Assume first that k is not a characteristic value, and let it approach the characteristic value k.. The limiting value of both sides of (2. 143) yields 49

THE UNIVERSITY OF MICHIGAN 3673-1 -F P (w, k., 0) a= -1 ( k,, 0) 2k (2.144) i (w 0) 2k k 0k The value k M(w, k, O) ak can be calculated analytically for the simplified model from the value of J(o, k., zl)determined in Section 2.4. Briefly,since in the region between z = 0 and z =90 Km all parameters are assumed to be constant, the following relations may be deduced: (i) From (2. 116), one has aN(w, k, z) =2V(w)N(w, k, z) - f (w)N (w,k,z)+iwf (w, k) (2.146) Equation (2.146) can be integrated analytically to yield N(w, k, 0) = N(w, k, z f s tan s z (2.147) where s is defined in (2. 138). (iii) Differentiation of (2. 147) with respect to k yields 50

THE UNIVERSITY OF MICHIGAN 3673-1 -F dN(ow, k, 0) Z + iwk tan s-1 sec s a k k 1 s 1 (2. 148) Combining (2.144) through (2.148) we have 2 _ P(w, k, 0) tan s z1 = - E(w, k. v z)+i ki -z sec sz) 0 W 2 k. i 1 (2. 149) Thus, A0 can also be calculated with ease. The formal solution of 1 due to a vertical section of a radially expanding source, at z = z' of length dz, which begins to expand at t = t', can now be easily obtained. Taking the z-axis as the axis of the cylindrical expansion, we have from Appendix B 1/2 V X (t, r, z) 2r So 6(x, 0)6(y, 0)6(z-z')T(t-t'(z')) The Fourier Hankel transform of the above yields 1/2 S (w, k, z) = 2r S Q P ( iwt((zT()) 2' (Q0' e and S1(W, k, z) = 0, here T(w) is the Fourier transform of T(t). Introducing these source functions into 2. 129) we obtain formally 51

THE UNIVERSITY OF MICHIGAN 3673-1 -F T -2 = 1 i 1/2 t(z' ) L 1 k F 1 2 H(z') Q 0(z') e (z -z') dz (2.150) where L denotes the linear differential operator L d 1 d d V (2.151) Ldjj+ ~ y-I (2.151) 1`dz fi dz 1z'f yAAGThe solution of (2. 150) is equivalent to finding the Green's function satisfying - 2 G(w, k, z -z') - k FG(z -z') = 6 (z -z') (2.152) If we expand the Green's function in the form G(w, k, z-z') = A P.,. i 1 1 then it is easily seen that, formally G(w, k' z-z') Pi(z')p 2(z) i A 2_ 2) 1 1 The Hankel transform of the above is G(w, r, z-z')Z = A 2 k2 i i 0 k.-k y using a contour integral with a proper deformation of the contour, it can be shown Sommerfeld, 1949) that 2 _1In i kr o2 2 = 2 H (k r),-i e i 0k -k o 1 l~1 52

THE UNIVERSITY OF MICHIGAN 3673 -1 -F This yields 2 2 P-1 P.(z')P.(z) ik.r G(wo, r, z-z')-,J-i Wrr 1 1 i A1i The Fourier component of the excess pressure at any altitude z, and distance r from a section of the path of the source of length dz' at z' on the z-axis due to a vertically moving'source is therefore given by 1/2 H (Z)1/2 (wo r, z-z')= 2rS T(w)e dz x p QOP H(z) / 2 -\ Pi(z')P.(z) ik r 2l /. e 1 (2. 153) The interpretation of (2. 153) is quite obvious, and may be explained as follows: (i) T(w) is the Fourier spectrum of the source. it t' (z') (ii) e is the modifying factor for the spectrum due to the motion of the source. (iii) 2ir S [Q p 1/2 is the strength of the source. 53

THE UNIVERSITY OF MICHIGAN 3673 -1 -F iir ik r(iv) e is the radial propagation faptor. P.(z')2 th (v) 1 is the excitation factor for the i mode, ( 2 A. 1 and H(z') Pi(z) (vi) H(z) P(z')- may be interpreted as the altitude factor, which is unity if H(z) P.(z') the source and the point of observation are at the same height. The time variation of the excess pressure, i. e. the inverse Fourier transform of (2. 153) is quite an involved problem. In general this is to be carried out numerically. Due to the introduction of losses, it is felt that corresponding to each frequency band, probablyonly the one value of k. that has the least imaginary component will be necessary. Thus the inverse integral may be performed numerically. Associate each frequency band with one value of ki, for which the values of Pi and A. are evaluated. From these values the time function jp(t, r, z-z') can be obtained as a summation over the significant range of frequencies generated by the source. A rough indication of the duration of the disturbance observed at z due to a source at z' then depends on the difference of the group travel times (distance divided by group velocity) between the band with the highest group velocity and that ith the lowest group velocity. 54

THE UNIVERSITY OF MICHIGAN 3673-1 -F The disturbance observed at a distance away from the vertically ascending source depends on the summation of the effect over all the segments of the path of the source. Here again a numerical integration has to be carried out over z'. The significant contributions would be those at altitudes where P.(z') and the source term have the greatest values. 2.6 Formal Solutions of the Source Problem A formal solution for the disturbance initiated by a source moving along a known path can also be expressed in terms of the mode functions discussed in the previous section. Consider the path of the source as shown in Figure 2-6. In cylindrical coordinates, the path may be expressed as functions of z by t =t (z) 0 0 r =r (z) 0 0 0 0 (z). O O If the idealized model of the source discussed in Appendix B is employed, the'equivalent source' at any height z is given by: X =0 (2.154) z X 27rs1Q p 1 /2 T(t-t (z) (r, z) )((2.155) This source function can be formally expanded into a Fourier Bessel series. The result is 55

THE UNIVERSITY OF MICHIGAN 3673 -1 -F Z =~expansion of gas (r, 0 ) o o z = 0 flat earth FIG.2-6: GEOMETRY FOR A STRATIFIED ATMOSPHERE 56

THE UNIVERSITY OF MICHIGAN 3673 -1 -F n =cc Oo oo -V *- = > e 1 e dw 2 kS (w, k, z)J (kr)dk 2r jn n n =-o O 0 (2.156) where 1/2 i o t 0(z) -inO 0(z) S (w, k, z) = T(w, z) 27rS [ Q P] e e and T(w, z)= T(t, z)e dt (2.157) -00 In the frequency domain, the equations of the disturbance take the following form. dp -iwP+P0 PV7 * v+vz = 0 (2. 158) -iwp v+ FV p + Azpg =X(w, r, 0, z) (2.159) -iwo'-v p g+ vy P GV v = 0. (2.160) Since the assumption of cylindrical symmetry is no longer valid, we define the Fourier-Bessel transforms of v, V v,,' by the following: v (w, r, 0, z) W (w, k, z) z n V v(w, r, 0, z) n=oo Go X(L, k z) ekdk D k n-oo 0 n rw, 0, z) P (to, k, z) (2.161)...... 57

THE UNIVERSITY OF MICHIGAN 3673-1 -F Taking the Fourier-Hankel transform of (2. 158) through (2.160) and the z-component and divergence of (2.159), the following scalar equations for the Fourier-Bessel components are obtained: dp -iwD+ +P ~ + d W =0 (2.162) n o n dz n dP -iwPo W + F n +D g=0 (2.163) o n dz n dp d2P dP dD 0 n 2 n dF n -iw pon-iw d W +F -Fk P + =S (, k, z) (2.164) i~ A dz n 2 n dz dz dz n dz -iwP -o p g+yP GX = 0 (2. 165) n no o n These equations, except for the indices n and the source term S are identical to the set of equations (2. 87) to (2.90); thus a similar reduction of the equations is possible. Parallel to the reduction of the equations given in Sections 2.4 and 2.5, one may define, for each n, z n=P Pn p /2 exp AFG (1 -F) dz (2.166) n n o 0 2,yAFG z U =W /2 exp 5 (1l-F)dz (2.167) nn 0 o 2,v AFG and obtain the following set of differential equations 58

THE UNIVERSITY OF MICHIGAN 3673-1 -F d 1 d V2d V V 1 2 H(z) and d z'f dz_ V + Uf'+ UV+H(z) s /(2.169) u'o The mode functions developed in the previous sections can therefore be used without modification. By use of the Green's function, P.(z')P.(z) G(w, k, z-z') =(2.170) i k2 1-k 1 1 the far region disturbance caused by the source can be obtained by multiplying (2. 170) by the strength of the source and integrating over the whole length of the path. From the path of the source we have, in a distance dz, the length of the path given by dr dO dL = dz r2) +r o2 (2.171) o dz o dz Introducing the notation dL L(z) = (2. 172) dz the expression for the Fourier Bessel component of the excess pressure is then formally given by +This holds except for a source moving in a horizontal direction. In that case, z is fixed and the integration should be carried out along dL. 59

THE UNIVERSITY OF MICHIGAN 36 73-1 -F [P(Z) 1/2ZT) P (w, k, z) 2 irS l dz' ()1 T(w, z' [QoPo length of path iwt (z) -inO (z')\ JP.(z')P.(z) L(z')e 0 e ~ (2.173) i A.(k2 -k ) 1 1 A partial inverse transform of the above can be carried out, yielding (w r, 0, z) 2?rS o( 1/2 dz H(z') 1/2 (r0z2R LS(z]/ length p(Z'1/2 T(z' e x of path P.(z')Pi(z) ik.R -iL(z') Z 1 1 e l 2 (2.174) i i i Ai where 2 2 R = r +r (z) -2rr (z')cos(0 -0 ) Equation (2.174) of course can be obtained directly from (2.153) by integrating over the contribution of all sources. 60

THE UNIVERSITY OF MICHIGAN 3673-1-F III IONOSPHERIC DISTURBANCES 3. 1 Introduction The neutral disturbances predicted by the linearized theory of wave propagation is usually of small intensity with the fraction of excess pressure, p/ P0 (probably on the order of a percent or less) and covers a large region of space. Such disturbances cause the electron density of the ionosphere to change slightly. This slight change of electron density modifies the phase delay of radio waves passing through it, thus causing a detectable signal. Our proposed model is primarily oriented toward the calculation of such effects, but since no concrete numerical results are available now, only a general discussion will be given in this chapter. 3.2 Effect of the Magnetic Field We shall suppose that a weakly ionized gas which is neutral as a whole is bound together not only by intermolecular attraction, but also by the attraction forces between ions and electrons. Thus if this ionized gas as a whole is moving in a direction perpendicular to the magnetic field, the magnetic forces tend to separate the ions from the electrons. The separation of these charges results in a force against the action of magnetic field thus minimizing its effect. To deduce a mathematical relation for such a mechanism, we shall assume that in a disturbed medium, the change of electron density is primarily due to collisions with neutral 61

THE UNIVERSITY OF MICHIGAN 36 73 -1 -F particles. The average transfer of momentum per unit time to the electrons from the neutral particle motion is therefore v(v-v ) -e where v is the collision frequency v the velocity of neutral particle motion and v the velocity of electron motion. The average equation of motion of the electrons is therefore governed by the linearized relation av - -e v -v wt mh -eree where e is the charge (numerical) of an electron E is the electric field, and B is the earth's magnetic field. To simplify the analysis without loss of the essential feature of this model, let us assume that due to its large mass an ion is less sensitive to the electromagnetic 62

THE UNIVERSITY OF MICHIGAN 3673-1 -F forces and is therefore moving with the same velocity as the neutrals. The'binding' force between the charged particles can be estimated from E. If the medium is neutral as a whole, then eE =V [eN (v -v] (3.2) at * o -e where N is the number density of electrons (or ions) in the unperturbed medium. This yields E = eN (v v), (3.3) which means physically that as the ions and electrons move apart, the bunching of ions and electrons creates an electric force that tends to pull them together. Combining (3. 1) and (3. 3) and taking the Fourier transform yields 2 eN e -t v + (v -v)-iwe v xB=-iwv (v-v ). (3.4) -e e m -e m -e -e o e e Introducing the notation e2N eWP N (plasma frequency)'b =. B (cyclotron frequency) b = (direction of magnetic field). the motion of the electrons caused by the neutrals is governed by 63

THE UNIVERSITY OF MICHIGAN 3673-1 -F 2 2 A 2 (-to + w - iw )v ) -iwob v xb = ( -iov )v. (3.5) p -e b-e p or 2 2 JL-1-i v % v A.i V -i-)v x b v (3.6) 2 o -e t -e 2 w tw Solving the above, one has e 2 2 2 -i -1 v + i -i v (vx) ( b) (b v) 2 |U2 CiV -2 Le- iW -1) — 2 0 (3.7) For a typical case near the F-region, from which we see that for to up to a few mc/ sec., 2 2 -it v v (3.v 8) -e 2 2 2 ~~-i- --- -- 3. 3 Extent of Travel of the Disturbance(37) Fo r a typical cedisturbances it is hoped that by cal-region, culmotion, but intrfew mode numbers andof the electric fore tespondings to'bind' the gaions for a typicalhere64

THE UNIVERSITY OF MICHIGAN 3673-1-F atmosphere, one may be able to calculate the dominant contributions to the variation of pressure in the ionospheric region. In general, from the frequency spectrum of the source and the group velocity associated with a dominant mode, one may estimate the lateral extent of the disturbance. On the other hand, the source spectrum and the mode function P. would yield information on the vertical extent of the disturbance. Crudely, the simplest manner by which the spectrum of the source may be estimated is by the piston model, assuming the mass and energy flow from the source expands laterally into the surrounding atmosphere until the energy released by the source is expended in the expansion. This yields the initial extent of the disturbance Q r = o (3.9) o 0 o and the initial time of expansion T =- X. (3.10) n the other hand, Holway (1960) used a blast wave analysis of a cylindrical expansion and assumed the initial expansion was limited to the region where the excess pressure is three times that of the ambient atmosphere and thus obtained the relation r =.392 a (3.11) o ~6P 65

THE UNIVERSITY OF MICHIGAN 3673-1 -F and PIQ o T =.160 (3.12) (+l 1)p Except for a numerical constant, these two sets of formulae have the same functional dependence and order of magnitude. The frequency spectrum of the source in terms of T is therefore proportional to sin2 2 This indicates that the dominant frequency band is limited to 1 In terms of 2z mode analysis, this source excites various modes of propagation. If losses are neglected, the source excites both the gravity mode and the acoustic mode. These modes propagate with different velocities, so that as the disturbance propagates far away from the source, the extent of the disturbancelbecomes larger. The extent of the disturbance also increases,due to the distribution of sources and the effect of viscosity. However, as a rough estimate in the ionospheric region, the extent of the disturbance is perhaps less than 100 Km. 3.4 HF Phase Delay From the analysis of Section 3.2, it is seen that due to the attraction between positive and negative charges, the motion of the charged particles follows quite closely the motion of the neutrals in the ionospheric region. Thus, the 66

THE UNIVERSITY OF MICHIGAN 3673-1-F approximate electron density change in the ionosphere may be written in the form An.'ne, = (3.13) n p`yp e 0o 0 Thus the calculation of the excess pressure at any point in space, to a first approximation, is sufficient to describe the perturbation of the electron density of the ionosphere. Due to the change of the electron density, the characteristics of radio waves propagating via ionospheric reflection are necessarily modified. Since, in the traveling disturbances, the fractional change is small while the range of the disturbance is large, the effect of the disturbance on HF propagations can be estimated approximately from ray theory. The fractional change of electron density causes a change of the equivalent dielectric constant of the ionosphere. Since E 3180 ne lo 2 (3.14) the change of the dielectric constant is 3180 A= 2 8 Ann (3.15) 2 e From the ray theory of propagation, the phase delay for a ray passing through a region of length dL with relative dielectric constant e is d =W FE- dL (3.16) C where C= 3 xl 0 m/sec is the velocity of light. Due to the change of electron density 67

THE UNIVERSITY OF MICHIGAN 36 73 -1 -F d(Az) o l e) dL (3.17) w dL 1 3180 (3.18) C 2 2 e Thus, we have (1 (l-E) An w 1 (l-e) " d(A)) 1 ( -) -d- dL. (3.19) C2 F n h d C 2 E yp If the variation of /Po along the total path of the ray is known, we then have A S-2C (1- 1 dL. (3.20) ray path This change of phase shift may be reflected as an instantaneous frequency shift given by f d (1- ) 1 (3.21)| Af = + - -P-) d L (3.21) 2C dt E v Po ray path or A 2 51=+)f (I a ) dL (3.22) ea at pt ray path provided that we assume that the ray path is essentially the same even with the perturbation. If one assumes that along the part of the path of the ray that is influenced by the disturbance, that e, ry, 5j/p are nearly independent of distance, then we see that the doppler shift is proportional to the time derivative of the pressure 68

THE UNIVERSITY OF MICHIGAN 3673-1 -F variation (Barry, 1962). For a more detailed description of the variation, the information of the pressure variation along a ray or bunch of rays as a function of time must be known. Then the doppler frequency for the composite signal can be calculated. Quite apart from the equation (3.22) one may also estimate approximately the time interval during which the perturbation can be detected. Roughly, if the extent of the disturbance is R moving with the effective velocity v and making an angle ac with the ray, then the time that the ray intersects the disturbance is calculated by simple geometry to be approximately R (tan c + coto ) 2v 2 2 More detailed calculations, including both the change of ray path and the variation of pressure as a function of space and time would have to be obtained for individual cases. In general the ray path of the unperturbed ionosphere has to be determined, and the effect analyzed. 69

THE UNIVERSITY OF MICHIGAN 3673-1-F IV. MODIFICATION OF ATMOSPHERIC DISTURBANCES DUE TO CHARGED PARTICLES 4. 1 Introduction When a disturbance takes place in or propagates through a region of the atmosphere which is partly ionized, electric and magnetic force fields and electric and magnetic energy storage deserve consideration in the theoretical study of the behavior of the atmosphere exposed to the disturbance. The objective of the present chapter is a critical evaluation of the modifications of this behavior caused by the presence of the charged particles in the atmosphere. The theory of wave propagation in an ionized gas, just as in a neutral gas, may be approached either from a macroscopic or from a microscopic point of view. In the former case the gas is treated as a continuous medium, forwhichmacroscopic characteristics, such as mass and charge densities, thermal conductivity, viscosity, etc, are postulated from experimental results. In the latter case, the statistical mechanics of the discrete particles is considered, based on the more fundamental postulates of particle characteristics, such as mass, charge, number density, and collision parameters. A discussion of the merits of both approaches will be included below. The decomposition of the atmosphere into several classes of particles, such as molecules, atoms, ions and electrons, raises the question of whether or not the gas behaves like a single homogeneous medium, the motion of which integrates the 70

THE UNIVERSITY OF MICHIGAN 3673 -1 -F effects of all the pressures, force fields, etc, acting on the various components of the gas. At the opposite extreme, a model of the atmosphere would be a superposition of a number of widely different media, each with its own approximately independent motion. In a whole spectrum of intermediate alternatives, the atmosphere may be pictured as a sum of distinct but more or less strongly interacting media. Because of the wide variation of density and composition of the atmospheric strata involved in the wave propagation considered in this report it does not seem possible to pinpoint one of these particular models of the atmosphere as the single appropriate one on which to base the analysis of the far-field propagation of ionospheric disturbances. The quantitative diversity of particle interactions also effects the behavior of each component by itself, placing it anywhere between the extremes of adiabatic vibrations, pure diffusion, and collisionless expansion, depending on density and collision parameters. Obviously only approximately adiabatic vibrations offer low attenuation and reasonably far reaching propagation; all other disturbances will be dissipated within a fairly short distanice of the source. In the neutral gas the *particle interaction and medium like behavior is entirely due to collisions; adiabatic conditions at a given frequency require a certain minimum density. However, since important components in the spectrum of the disturbance are of extremely low frequencies and long wavelengths, the maintenance of approximate local and temporal equilibrium throughout a complete cycle is possible even at fairly high altitudes. In 71]

THE UNIVERSITY OF MICHIGAN 36 73 -1 -F the neutral gas the particle interaction and medium like behavior is entirely due to collisions; adiabatic conditions at a given frequency requires a certain minimum density. However, since important components in the spectrum of the disturbance are of extremely low frequencies and long wavelengths, the maintenance of approximate local and temporal equilibrium throughout a complete cycle is possible even at fairly high altitudes. In a gas of charged particles, on the other hand, the interactio is due to the electromagnetic field, the longer range of which makes a medium-like behavior possible at much lower densities. In a thin, fully ionized plasma an adiabatic compressional electron wave is possible for a limited frequency range above the plasma frequency. This frequency band is, in general, too high to be of any interest for ionospheric disturbances. A compressional positive-ion wave has no cut-off frequency and will be affected by Landau damping - in other words, no adiabatic wave of this type exists. 4.2 Discussion of Fundamental Relations A macroscopic theory of the dynamics of a thin, completely ionized gas differs from that of a neutral gas by the addition of electromagnetic accelerating terms in the force equation and by the introduction of Maxwell's equations relating the fields to the charge and current densities in the gas. Under adiabatic conditions the results are reasonably simple and have been thoroughly investigated by many authors. In the absence of a terrestrial magnetic field and neglecting discrete collisions, linear propagating disturbances may be separated into two classes: 72

THE UNIVERSITY OF MICHIGAN 36 73 -1 -F transverse waves, which are similar to electromagnetic waves in a vacuum except; for a modified dielectric constant, and compressional or acoustic waves. The latter are those of interest here. In these modes of vibration the behavior of the plasma may be described as that of two separate vibrating fluids, the ion gas and the electron gas, respectively, mutually coupled by the electric field. The differential equations may easily be brought to the same form as the conventional equations for two transmissio lines with mutual coupling. au l n i (4. 1)| p =-jwL *p (4.1) ax'y' n n n'/n Pno apn 1 2 1 =-_jp.u - p (u -u ) =-jwC u - 1 (u -u) (4.42) ax pno n jw no n n p nn JwL n p where u, p and p are velocity, pressure and mass density, respectively. The su3) nx p p m wh(see psketch A of Fig. 4 -1). They, pre issure a dimensional inconsistency in the symbols, n Y Pmp 3 p p7 (see sketch A of Fig. 4 -1). There is a dimensional inconsistency in the symbols, since jwL n, jwL are reactances per unit length, while (jwL ) is susceptance per 73

THE UNIVERSITY OF -MICHIGAN 3673-1-F unit length. The choice of velocity as equivalent voltage and pressure as equivalent current rather than vice versa avoids the use of an ideal transformer in a capacitive coupling link between the various branches of the networks. The equivalence provides a simple description of the different elements of energy storage involved in the propagation. In Fig. 4-1, B and C illustrate the two types of waves possible. In the electronic waves, the two lines vibrate approximately in phase opposition, and the resonance between L and the capacitive elements cause m the cut-off frequency o =(2 + ) 1/2 In the ion oscillations, on the other hand, c n p the two systems are nearly in phase; consequently there is very little energy storage in L, which has a negligible influence on the propagation, and no low-frequency cut-off occurs. If the gas is not fully ionized, the neutral gas particles form a third component in the system which depends entirely on interparticle collisions for exchange of energy and momentum. If the neutral-particle density is high enough for adiabatic ehavior and if the masses of the neutrals and the positive ions are approximately equal, the neutral gas component and the positive-ion component are necessarily very strongly coupled, and will nearly act as a single medium. Statistically speaking, every heavy particle will have roughly the same history of collisions, gain and oss of momentum and energy by collisions, regardless of whether it is neutral or ositive. Tentatively, the behavior of the gas may be described by the same equaions (4. 1 - 4. 4) if the subscript now refers to the total population of heavy particles 74

THE UNIVERSITY OF MICHIGAN - 3673 -1 -F L L L n n n vn n_ A Lm L P P P L L L n n n p- C L L L p p P I I FIG. 4-1 75

THE UNIVERSITY OF MICHIGAN 3673 -1 -F and w2 is obtained by spreading the same total charge uniformly over all the heavy p particles. Clearly C and 1/C increase with the total number of heavy particles while L remains constant; the result is reduced coupling to the electron component. m Despite the large difference in mass, the electron behavior will not be completely unaffected by collisions with neutrals; a certain amount of scattering and attenuation results. If the behavior of the neutral component is not adiabatic, the more complete macroscopic theory presented in previous chapters should be introduced, with appropriate additional coupling terms representing exchange of energy and momentum between components. For several reasons it does not seem fruitful to pursue this course of investigation; these will be mentioned below in the discussion of the ionized atmosphere from the microscopic point of view. A macroscopic theory of the behavior of an ionized gas in the presence of a magnetic field has been the subject of a large amount of literature, conventionally classified under the name of magnetohydrodynamics. The gas is characterized by its fluid-mechanical constants and its electric conductivity; its wave equation is derived from the hydrodynamic and electrodynamic equations. The acoustic wave of the neutral gas breaks up into two waves, one faster and one slower than the velocity of sound; in addition, a transverse wave occurs, which has the Alfven velocity vA'=1 __ (4.5) A Pp 76

THE UNIVERSITY OF MICHIGAN 36 73 -1 -F where p is the permittivity and p the mass density of the gas. For the small geomagnetic field the magnetohydrodynamic modification of the acoustic propagation of disturbances in the atmosphere appears to be insignificant. Statistical mechanics furnishes the alternative approach, the microscopic approach, to the theory of propagation of disturbances in a gas. At the price of considerably increased complexity and statistical rather than completely deterministic conclusions it. is now possible to express the results in terms of the more fundamental properties, such as mass, number density, charge and velocity distribution of the molecules, atoms, ions and electrons of which the gas is composed. The unperturbed gas is assumed to be in equilibrium, so that the Maxwell-Boltzmann distribution can be assumed. The continuity equation in phase space, i. e. the Boltzmann transport equation or the Boltzmann-Vlasov equation, combined with Maxwell's equations for the electromagnetic field constitute the basis for calculation of the propagation of disturbances. In a gas where discrete collisions are the dominating means for exchange of energy and momentum, the random-flight theory and the Fokker-Planck equation offer an alternative microscopic point of view instead of the one obtained from the Boltzmann transport equation. For the range of parameters with which the present study is concerned, however, the Fokker-Planck approach is not felt to be sufficientl appropriate to justify any further discussion. The thin completely ionized plasma, where discrete particle collisions can 77

THE UNIVERSITY OF MICHIGAN 36 73 -1 -F be neglected, is the simplest, most completely investigated problem in this area. If the Boltzmann equation is multiplied successively by powers of the velocity from zero and up and integrated over velocity space, a set of equations is obtained corresponding to the macroscopic conservation relations, but it does not terminte. It relates by iteration all the infinite set of moments of the velocity distribution with each other. By instead performing a Laplace transformation in time and a Fourier transformation in space, Landau (1946) obtained a dispersion relation which permitted a general analysis of all the possible solutions. He showed that approximately adiabatic solutions must have wavelengths large compared to the Debye wavelength. This condition is satisfied in a range of frequencies immediately above the plasma cut-off frequency, where the phase velocity of the perturbation falls outside the main range of thermal particle velocities. Approaching the limit T = 0 this frequenc range shrinks to a single frequency, the plasma frequency, the phase velocity becomes arbitrary (0/ 0) and the group velocity zero, in agreement with the classical theory for a tcold' plasma. Landau's conclusions have been studied and extended by a number of other authors and are, with minor modifications, universally accepted. If classical statistical mechanics is postulated to give a satisfactory model of a thin plasma, the theory of Landau and his successors gives the most fundamental and rigorous results as far as the propagation of small perturbations is concerned. 78

THE UNIVERSITY OF MICHIGAN 3673-1-F e may then use this theory for a comparison with and a criticism of the macroscopic theory discussed previously. The conservation equations in the macroscopic theory correspond one by one to the moment equations obtained by integration of the microscopic transport equation over velocity space. The microscopic theory implies that also the macroscopic theory should give an infinite set of equations between an infinite set of moments. Is it justifiable to cut the number down to a finite set of equations, and if so, where should the cut be made? The most drastic cut is obtained by assuming adiabatic variations only. The Landau theory concludes that approximately adiabatic, i.e. reversible, unattenuated variations can occur only in a relatively narrow frequency range above a plasma cutoff frequency, if such a frequency exists. This suggests that the electron waves obtained from the macroscopic adiabatic theory of a completely ionized plasma discussed above exist in a frequency range above this low-frequency cutoff, while the eavy ion gas is unlikely to show any appreciable medium-like propagation of comressional waves, since it has zero cut-off frequency. Within the adiabatic frequency range the same solution for the electron waves is obtained from both macroscopic and microscopic theory; the two theories thus are supporting each other in this ange. Any attempt to extend the scope of the macroscopic theory by including a larger number of equations and moments appears to be on shaky ground, as far as he thin plasma is concerned. The number chosen appears somewhat arbitrary and 79

THE UNIVERSITY OF MICHIGAN 3673-1-F difficult to justify from convergence considerations. The microscopic theory views the disturbances in a different light: An infinite number of solutions must exist in order to make it possible to match initial and boundary distributions of any number of degrees of freedom. Most of these solutions are of a local or transient character, and the problem is to find those solutions which reach far andlast a long time because of their wave-like character and low attenuation, if such solutions exist. Arbitrary termination of the infinite set of moment equations does not necessarily pick out physically significant solutions from the infinite set. At higher density and in the presence of neutral components in the gas, the collision integral term in the transport equation cannot be neglected, and the microscopic theory becomes very complicated to handle. However, it has been shown (Bhatnager, et al, 1954) that a similar general disperion relation study to that of Landau can be carried out if a somewhat simplified model of the collision integral is used. Thus, any arbitrary inclusion of higher moments is avoided but the adoption of a simplified collision-integral model admittedly involves a certain arbitrariness. The design of the model so as to be consistent with a finite number of conservation conditions is somewhat analogous to the moment truncation. Since it does not as severely restrict the degrees of freedom of the system, however, it is a less drastic and intuitively more appealing procedure. The results indicate that waves with low attenuation occur when their wavelength is large compared with both the mean free path and the Debye length. If the wavelength is of the order of the 80

THE UNIVERSITY OF MICHIGAN 36 73 -1 -F Debye length, the same result is obtained if the collision frequency is large compared to the plasma frequency. At smaller collision frequencies, very large'phaseslip' damping of the same nature as Landau damping is indicated. The latter condition prevails, of course, even more at shorter wavelengths. Neither the electric field nor the collisions are then capable of maintaining phase coherence between the particle vibrations under these last-mentioned conditions. So far we have here discussed only homogeneous media. If the stratified character of the upper atmosphere is taken into account, the possible modes of wave propagation multiply considerably. Variable composition, density, and temperature modify the characteristics of the waves described above, and'waveguide modes' surface waves, etc., associated with different strata or boundaries between strata are added to the picture. A detailed quantitative investigation of the contributions of ions and electrons in this area would be a very large project in itself. However, some general inferences may be reached by comparison with the results for uniform media. It is safe to extrapolate the fact that pure electron waves are far below their cut-off frequency in the frequency range of ionospheric disturbances to all kinds of such waves. As far as the ions are concerned, on the other hand, it is apparent that they will play a similar part as in the homogeneous media; they will in general be tightly coupled to the neutral gas, whenever the gas has propagating properties, and merely contribute a small amount of electrical energy storage, which will very slightly modify the propagation. 81

THE UNIVERSITY OF MICHIGAN 3673-1 -F V CONCLUSIONS A realistic linearized model for the calculation of long range ionosphere disturbances due to a source moving through the atmosphere has been formulated. This model incorporates the stratification of atmosphere, the effect of the gravitational force, viscosity and thermal conductivity. The form of loss terms including the effect of viscosity and thermal conductivity is believed to be somewhat better than the Navier-Stokes form, in that it reduces to the Navier-Stokes result for small values of viscosity, but for large viscosity, leads to results closer to those predicted by the kinetic theory of rarefied gases. The formulation is based on the mode theory of propagation, thereby reducing the mathematical problem to the calculation of a few mode numbers and mode functions. The disturbances due to any source may, in principle, be evaluated. Based on the experience gained in a previous similar calculation, a procedure (partly analytical and partly numerical) has been formulated to calculate the mode numbers and mode functions. The effect of a neutral disturbance on the ionospheric perturbation and HF forward propagation has been discussed. An explanation clearing up the long unanswered question of the effect of the magnetic field on the ionospheric disturbance due to a neutral source has been given. The effect of space charge forces on the propagation of disturbances has been discussed in general and it has been concluded that in the low frequency range of interest in the present problem, the space charge forces can be neglected. An idealized model for the source of excitation, incorporating the nonlinear effect in the initial stage of disturbance has been suggested. Although no actual calculations were made using the present model, it is flt that the computation is feasible and the results would be useful in the prediction of traveling disturbances. 82

THE UNIVERSITY OF MICHIGAN 3673-1 -F APPENDIX A: APPROXIMATE FORMULATION FOR LOSSES The hydrodynamic equations governing the motion of a fluid are given by at p +7V (pv) = O (A.1) P (O- v+-.7v)+V P+V' _ (A.2) - at-1 +- / /-1 TPV' V+ V7 Q+VV: P 0 (A.3) Since these equations can be interpreted as the first three moments of the Boltzmann equation, numerous formulas have been derived relating _ and Q to the macroscopic quantities, p, p and v by using different approximate solutions of the Boltzmann equation and different models of collision. A particular class of the formula applicable to rarefied gases that has been verified qualitatively by experiments was given by Holway (1963) and Kahn (1963). Since the upper atmosphere consists of a very rarefied gas several hundred kilometers above the ground, it seems important that some features of their result should be incorporated into the equations of motion for the atmosphere. On the other hand, their results are derived from the kinetic theory of homogeneous gases; thus, the direct application of their results to an inhomogeneous atmosphere is not easily justified. Moreover, the propagation of disturbances in the atmosphere is, in itself, quite a complex problem, and the introduction of their results into the equations of motion would complicate the problem considerably. 83

THE UNIVERSITY OF MICHIGAN 36 73 -1 -F It is, therefore, the purpose of this appendix to derive a simple formulation based on the macroscopic description of a gas and yet retain some of the features of their analysis. Our simple approach is based on the following intuitive assumptions. a) An inhomogeneous medium may be approximated by a continuum of locally homogeneous media with varying local properties. b) The macroscopic equations of motion and the static approximations for viscosity and thermal conductivity, i.e.: -/ (VVT +V v-7. _vI) (A.4) and Qe = - V (e-1) p)' (A. 5) perhaps with some modification, may be adequate. c) Since all the losses are associated with the motion of the fluid, all loss terms may be expressed in some form related to the time derivatives of the macroscopic variables, p, v, and p. Thus, for perturbation solutions at some frequency the effect of losses may be expressed in terms of some function of w times the amplitude of the perturbation. These functions are to be evaluated from the one dimensional acoustic wave propagation in a viscous medium. The usual small signal perturbation equations for a sound wave of angular frequency w are easily obtained from (A. 1) through (A. 3). We shall write them This simplifies the mathematical formulation for the wave propagation in stratified media considerably. 84

THE UNIVERSITY OF MICHIGAN 3673 -1 -F in the forms -ofp+ P o * v 0 (A.6) -iw p v +tvp+ f =0 (A. 7) -iW''+YPV.v+g =0 (A.8) -wp + y P7o7' v + g= =0 where, in view of (A. 3) and (A. 4), and the one-dimensional approximation ft =V -P VV V (A. 9) g=('-1)V * =_ - C p V p V P (A. I0) -p0 0 Some remarks on p and X are to be noted before the derivation of the approximate equation. (i) It is well known that p/ p is of the order of the mean free time of a homogeneous gas; therefore, the parameter determining the viscous loss should be proportional to w p / p (ii) The thermal conductivity may be expressed approximately in terms of the viscosity. As a first approximation, we may have 3 = (9 Py -*5) (A. 11) Thus, we may write Cfu 3 Po=J 4v (A.12) and 85

THE UNIVERSITY OF MICHIGAN 3673-1 -F -- (9 - 5) 2 - V (A.13) Now, as an approximation, we may replace the space derivatives of p and p in the expressions for f and g. by the functions obtained from the solution of the one dimensional non-viscous flow problem. From (A. 6), (A. 7) and (A. 8), one notes that if f and g are neglected, P o. (A. 14) P o. iwr V' v p (A. 15) -/p0 and V -— ip oV ~ (v (A. 16) Introducing these relations into (A. 12) and (A. 13), one obtains an approximate formulation for loss terms as f -i 4 Vp (A. 17) -p Co and \ -iWy ( )2 y)Adv (A.18) 4Cy 2 7 o here C = is the local velocity of sound in the lossless case. o P86

THE UNIVERSITY OF MICHIGAN 3673-1 -F To check the validity of this approximation, (A. 18) and (A. 19) are introduced into the equations of motion (A. 7) - (A. 8) to yield -iWop v+ VP i2 (A. 19) | — i(+^YP v -'PO V v 0~ (A. 20) ik- r For plane wave propagation according to e -, these equations yield kCr1 _I1/2 W C- 2 3 2 (A. 21) C C O ) For w v < < 1, one finds that k C3 C P 2 +..) o This checks with the conventional attenuation constant of sound waves. For v > > 1, one finds that kC - O0 (A. 23) ince the phase velocity of propagation is C = a/k, (A. 23) indicates that C/ Co- co (A. 24) ~s wv - co. This is the usual conclusion that is obtained from the Navier-Stokes equation of motion. In general, however, this conclusion is not verified by experiment (Holway, 1963). 87

THE UNIVERSITY OF MICHIGAN 3673-1 -F Thus in the investigation of disturbances in the ionosphere where the atmosphere is quite rarefied, the effect of viscosity must be examined more closely. This indicates the necessity of introducing a more adequate model to account for the viscous forces than is given by (A. 4) and (A. 5). To propose an improved model, we recognize that (A. 4) and (A. 5) are only approximate in the quasi-static case, i. e. assuming very slow variations. Physically, the viscous effect takes place through molecular collision, so that for wv not very small, the viscous effect may'lag' somewhat behind the velocity or pressure changes. To account approximately for this lagging effect we resort to the more complete description obtained from the higher order moments of the Boltzmann equation. For example, by Grad's 13 moment approximation (Grad, 1959) it can be shown that the pressure tensor P and heat flux Q are interelated by the following: + +P +' -1 Q+ 7T 2 q - p at - V 3 9 + P vv+ vv -3 V' v I 0 (A.25) +' - +' P++ Fp- V - y-1 P at y —1 P0 4e(y-1) p[V P0 o (A. 26) hese equations show clearly that (A. 12) and (A. 13) are only approximate, and for w- large, the time derivative term should be considered. A more precise po approach to account for the effect of losses is, of course, to solve (A. 25) and + These equations were originally derived for mono-atomic gases. We have made the proper corrections to the coefficientsby introducing the factor y. 88

THE UNIVERSITY OF MICHIGAN 3673-1-F (A.26) together with the equations of motion (A. 1) through (A. 3). For a onedimensional flow in a homogeneous medium, this approach has been used by Al (1961). However, for an inhomogeneous medium, we find that this precise approach is somewhat inconvenient. Therefore, as a simplification, we again consider a one-dimensional problem. From the Fourier transform of the linearized form of (A. 1) through (A. 3), (A. 25) and (A. 26), the following relations between and P are obtained 4 v -i 3 C 0 V. p.wv vP (A. 27) and V 9y - 5)) P -i ( 4 ) i- i VV V (A.28) ( -1) C Using (A. 27) and (A. 281 the expressions for g and f are v 4 2 3 f v 5 VP (A.29) 1-i ~oyW 5 P 2 3 0 and Since these expressions are only approximate we shall take v=5/3 in (A. 27) and (A.28) 89

THE UNIVERSITY OF MICHIGAN 3673-1-F -i 2V.wv% 5 yp_ v (A.30) 1-i2 2 0 Since these forms cannot be directly interpreted as time derivatives, they can nly be used in the Fourier transform of the perturbation equations. Introducing (A. 29) and (A. 30) into (A. 6) through (A. 8) we obtain |iwP'p7V 0 (A.6) -iw p +po ~Jv - -iwpL v + FV P=O (A.31) and -iwjp+ y/p GV * V = o (A.32) here the factors used to account for viscosity are W V l-i —-2 3 2 C F 1ix0 (A.33) xxjv 5 23 0 1 -i- C2 090

THE UNIVERSITY OF MICHIGAN 3673-1-F These equations yield kC (A.35) For wv <<1, (A. 35) yields kC ( ~ ) +i 2 2C2 0 which is the same as (A. 22) if we assume V= 5/3. On the other hand, for wv>> 1, (A. 35) yields the limiting value C k -0.62 This is comparable to the experimental value of.42 given by Holway (1963). Although this model accounts for the effect of viscosity in a more realistic manner than the Navier-Stokes equation in that the phase velocity remains finite in the high wv limit, it is recognized that it has its limitations such as zero damping in the high wv limit. It is perhaps too much to expect more of such a simple model. It is felt however, that this model should adequately account for the effect of viscosity for the purpose of this report. APPENDIX B IDEALIZED SOURCE DESCRIPTION The initiation of a disturbance by a source of mass and energy moving through the earth's atmosphere is a very involved non-linear problem. We shall assume that the primary source of the disturbance is the release of a gas which is at a very high temperature and pressure with respect to the ambient atmosphere 91

THE UNIVERSITY OF MICHIGAN 3673-1-F and which expands rapidly, thus producing a shock wave. The general solution for a gas expanding into an inhomogeneous medium, including the effect of viscosity, has never been attempted. A rough estimate of the parameters describing the initial stage of the disturbance obtained by similarity solutions of a cylindrical blast wave are available in the literature (Holway, 1960). Assuming that the pressure of the expanding gas is much larger than the ambient, the growth of the radius of the initial disturbance depends on the total energy Q deposited by the source per unit path length and the ambient density p according to r _=S~z-1)(S0l) /(') / t1/2 (B.1) Of course, as the gas expands, the pressure begins to drop, thus violating the assumption of a large pressure ratio and thus rendering (B. 1) inadequate. An order of magnitude analysis using (B. 1) has been carried out by Holway (1963) for two typical sources. Assuming that the final radius of the source gas is determined by the equilibrium of pressure of the gas and ambient atmosphere, the final radius can be estimated from a knowledge of the radius of the source and the energy released by the source. For the purpose of this report we shall assume a final radius of the cylinder of gas as less than 1OKm while the initital time of expansion is a few seconds. Since we are interested in long range disturbances, of the order of a few undred kilometers from the source, the initial stage of expansion may be considered 92

THE UNIVERSITY OF MICHIGAN 3673 -1 -F to be confined to a cylinder of very small radius. From (B. 1) one finds that the radial velocity of expansion of the cylinder is r 1dt 2 2(y -1)(y+1 /2 Q 1/2 1 v (B.2) r dt 2 7 r During the initial stage of the expansion, therefore, yrr 2 II N ] (Q =S() (B.3) = 2[2 (z1)(/+l)] l/2 Q 1/2 0 0 is a constant. The fact that the source may be enclosed in a volume that is small compared with the distance at which the disturbance is observed, and that the product of v r is a constant, leads to our idealized formulation of the source. For a source that ascends vertically, one may use cylindrical coordinates with the path of the source as the axis. The position of the source can then be considered the composite of small vertical sections. Consider a section of gas at z' which starts its initial expansion at t'(z'), with the radial component of the velocity satisfying (B.2). Confine this initial expansion to a small cylinder of cross-section area da. Then since da is small, to an observer far away from the source this velocity may be considered as caused by a distributed force per unit volume X satisfying t' +T' p (z') dt = v, (B.4) 93

THE UNIVERSITY OF MICHIGAN 3673 -1 -F where I' is the time for which the blast wave analysis is valid. In view of (B. 2) one has x z (B. 5) Moreover, from the fact that 27rrv ='v da, (B.6) one has tt + T It 1/2 x 2 ~2r rv Q (zt) r o 1 V (- ) dt= V' vS L ) =Sa (B. 7) P, da P (Z da Since to a distant observer da is small, the source may be considered as a line source, and therefore may be idealized as a 6 -function. Therefore t' +'" r 1/2 X FQ (z') ( ( —) dt = 2 ir S j 6 (x, 0) (y, 0). (B. 8) T is generally a function of the height, the energy and the criterion by which one may assume that blast wave analysis is valid. Assuming that the blast ave analysis is valid when the pressure behind the shock is 3 times the ambient pressure, Holway (1960) has shown that QP T' (z') =.160 ( -lo (B. 9) (y+l)p0 From (B. 8) one may see that as a function of space and time 94

THE UNIVERSITY OF MICHIGAN 3673-1 -F V - = =2rS XP0 =2L7rS 0 (z') 6 (x, 0) 6(y, 0) T(t), (B. 10) where the time function T(t) describing the build-up of the initial shock must satisfy t'+ T' \,. T(t)dt = 1 (B. 11) As a crude estimate, since I' is usually small compared to the time required for the disturbance to propagate from the source to a distant observer,we may approximately take when T(t)=,(z,)'+t' >tI t' (z') (B. 12) = 0 otherwise. For the limiting case of T' —- 0, of course, T(t) approaches a 6 -function, given by T(t) = 6 (t - t'(z') ) (B. 13) Thus, for a vertically moving source with a path given by t(z), to an observer at a large distance from the source the disturbance may be assumed to be caused by a ictitious volume force generated from the path. This fictitious force satisfies (B. 5) nd (B. 10) everywhere along the path. The above argument can easily be applied approximately to a source that moves along an arbitrary path. If one denotes the path of the source by x =x(t), =y (t) and z = z(t), then at any time t = t', one has approximately 95

THE UNIVERSITY OF MICHIGAN 3673-1-F o)1/2 Q (zt) 7* = 2r S 6(x-x) 6 (y-yo) T(t) (B.14) Lp (zt _') 0 and, on the average X = 0. (B.15) z Both of these relations are approximate in the sense that we have assumed that the initial expansion has cylindrical symmetry. Thus,in replacing the initial disturbance by a volume force, the average volume force in any direction is zero, which implies X = 0. The volume force is therefore described by its divergence relation which is (B.11). A similar conclusion can also be obtained if we assume that the initial expansion takes the form of a series of spherical blast waves originating at each point along the path of the source. It is to be noted that the use of an'equivalent source' described by a scalar quantity, i. e. ~ (X / po), is equivalent to considering the initial disturbance as a ressure source. To demonstrate this, let us consider the linearized equations of motion for a non-viscous homogeneous fluid. These equations are: 0o at + P x'+ /po. v=0 at0 Eliminating v yields 96

THE UNIVERSITY OF MICHIGAN 3673 -1 -F YP 2 a2 P 2 P Po at o Thus the source we specified is, in reality, a source of pressure. Since it is well known that the solution of the inhomogeneous wave equation is equivalent to that of a homogeneous equation with initial values specified on the surface enclosing the source, the initial stage of the disturbance may also be described by a set of boundary conditions in terms of the excess pressure on a closed surface at every instant. The solution of the linearized equations with a given set of initial values or pressure on an imaginary surface is an approach which has met with considerable success (Weston, 1961) in obtaining a description of the acoustic disturbance due to a nuclear explosions. The above discussion suggests an alternative possibility for obtaining a more detailed description of the source. For a more realistic description of the source, one may carry out the following near field analysis: (i) A blast wave analysis for a very short period of time, assuming spherical or cylindrical symmetry. This yields a velocity distribution at every point in space at time t. (ii) The results of the blast wave analysis are then used as the initial values for the geometrical acoustic analysis described in detail by Lighthill (1956). This approach, using the solution of the Burger's equation along various rays determined by the linearized equation incorporates both the non-linear and viscous effect and 97

THE UNIVERSITY OF MICHIGAN 3673-1-F hence is more realistic than a direct application of the linearized equations. The result of this analysis could yield the values of p (as functions of time) on an imagineary large vertical cylinder surrounding the path. The size of this' cylinder is determined by the condition that the value of p on the cylinder is a small fraction (say, 0. 1) of the ambient pressure. The values of pressure on this cylinder can then be used as the boundary conditions for the linearized equations. Calculations of this sort, which depend on detailed information regarding the source energy and path are necessarily complicated and have not been attempted in this report. 98

THE UNIVERSITY OF MICHIGAN 3673 -1 -F APPENDIX C NUMERICAL ANALYSIS C. 1 Description of the Problem An iterative procedure for obtaining the characteristic values of the boundary value problem describing the propagation of a mechanical disturbance in the atmosphere was described in a paper by Chu and LaRue (1963). The problem outlined in the above paper differs from that described in the present report in the following manner: 1. The effect of thermal conductivity was neglected. 2. The variation of the gravitational acceleration g and the ratio of the specific heats y was considered. 3. The approximation for the effect of viscosity was that described in the first part of Appendix A. i. e, the lagging effect due to collisions was not considered. 4. The variation of A (temperature) was approximated by quadratic rather than linear functions of altitude. 5. The numerical integration necessary to determine the eigenvalues was carried out from 700 Km to the ground. The analytic functions used for A, v, y and g are given in Table C-I and these functions are compared with values from the U. S.Standard Atmosphere(1962) and the Astrophysics Research Corporation(1963), in Figs. C-1 through C-4 respectively. The linearized differential equations of motion describing the motion of the atmosphere at a large distance from the source were expressed as dz +Vp =fl Z (C. 1) The notation used in this appendix is the same as that used in the reference given and, in general, differs slightly from that used in the rest of this report. 99

THE UNIVERSITY OF MICHIGAN 3673-1-F Table C-I FUNCTIONAL REPRESENTATION OF ATMOSPHERIC PARAMETERS Function Altitude Range (km) R * T A = M 6x104 0-100 3.43Z - 2.83 x 105 100 -135 -1.28 x 10-5 Z2 +21Z - 1.685 x 106 135 - 190 -1. 18 x 10-6 2 + 1. 82Z + 8.7 x104 190 - 700 Kinematic Viscosity, v Altitude Range (km) -5 14. 7 x 10 Z 1.46 x10 e 0 - 135 -5 6.77 x 10 Z 0. 662 e 135 - 190 -5 4.24 103 2.17 x 10 Z 4. 24 x 10 e 190 - 400 5 1.15 x 105 2.51 x10 e 400 - 700 Ratio of Specific Heats, y Altitude Range (km) 1.4 0 - 135 4.78x 10-7 Z + 1.3355 135 - 700 Acceleration due to Gravity, g Altitude Range (km) -2.63x10-6 Z+9.81 0- 700 All Functions are in Mks Units. 100

THE UNIVERSITY OF MICHIGAN 3673 -1 -F 80 70 60 50 Cq S 40 X / - From Tables ~ 30.. Approximation 20 10 0 I 0 100 200 300 400 500 600 700 Altitude -km FIG. C-1: A VS ALTITUDE 101

THE UNIVERSITY OF MICHIGAN 3673-1 -F 9 10 108 0 10 -2 5 107 0 10 Altitu d -A 101 -3 0 0 100 200 400 500 600 Tables Altitude -km Approximation 10~~~~~102

THE UNIVERSITY OF MICHIGAN 3673-1 -F 1.8 1.6 t.4; 1.2 1.0 - From Tables 0.8 m I - ~ Approximation o 0.6 ~0.4 0.2 0 100 200 300 400 500 600 700 Altitude — Km FIG. C-3: RATIO OF SPECIFIC HEATHS, v, vs. ALTITUDE 103

THE UNIVERSITY OF MICHIGAN 3673-1 -F 10 9 C\1 7 6 5 4 4 - From Tables a-)..- Approximation A3 a)., C,) o 0 100 200 300 400 500 600 700 Altitude-Km FIG. C-4: ACCELERATION DUE TO GRAVITY, g, vs. ALTITUDE 104

THE UNIVERSITY OF MICHIGAN 3673-1-F and dT dz V Z= f2 0(C.2) where the functions ~ and zcan be compared with (1 and U1 in (2. 97) and (2. 98), respectively. The functions fl, f2, V and F are 1 2 f - (Fg-w ) (C. 3) r i2] f2=k(1- A ) v (C. 4) %(: -1) Al' A V 2A 2y(y-1) - i 2A(C.5) A and F AI'71 + g(,y- (C. 6) A (P-1l)Py yA From (C. 1) and (C. 2) the differential equation satisfied by ~ can be easily obtained _s_ d'] V V2 1 dz Lbf dz j Ldz f f (C. 7 The boundary conditions which J must satisfy are z =O,= -v(O)) (C.8) z =b(700 km) =' = + u(b)] (C.9) 105

THE UNIVERSITY OF MICHIGAN 36 73 -1 -F where 2 A k iwv + F = u + IV yA (C. 10) 2 i~ ) A An iterative procedure similar to that outlined in Section 2. 4 was used for calculating the characteristic values k and for convenience will be repeated briefly here. We begin by defining the quantities N(k, z, w) = (k, z, w) (C.12) and aN J a= (C.13) ak The differential equations satisfied by N and J can be obtained from the definitions (C. 1) and (C. 2) and are N' = f2 - f N2 + 2VN (C. 14) J' =af2 -2f J+2VJ. (C.15) ak t is to be noted from the boundary conditions that 106

THE UNIVERSITY OF MICHIGAN 3673-1-F N(k, bi, at) (ac-u + V) iwoG (Fg-w2) itoG Laa J(k, b, ) = 2 ak Fg -o are known functions of k. Moreover, N(k, 0, o) = 0 is independent of k. Assuming initially any fixed value of k and using the boundary values, (C. 14) and (C. 15) can be integrated from z =b to z = 0. Denote these values of N and J by N(k, 0, w) and J(k, 0, o). Now if k is an eigenvalue, say k, then n N(k, 0, to) = 0. n But by definition of J one has N(k, 0, Lo) J(k, 0, tw)(k-kn) (C.16) Hence, as a first approximation kk- N(k 0 N w) (C. 17) n J(k, 0, 0o) The process can be reiterated until it converges. C. 2 Computer Program The iterative procedure outlined above was programmed for numerical analytis on an IBM 7090 digital computer. Considerable difficulty was experiencedwith convergenceproblems inthe originalprogram, particularly at points of discontinuity of A' (proportional to T', the temperature gradient). Accordingly an error criteria was established and inedincorporated into the program. This criteria consisted essentially 107

THE UNIVERSITY OF MICHIGAN 3673 -1 -F in comparing the absolute value of the real and imaginary parts of N and J at a new point with those obtained at a previous point. If the differences in these values were too large,the integration'step size' was automatically halved and a new point estimated. This procedure was continued until a satisfactory estimate for the new point was established. After the modification of the computer program as described above, convergence was obtained with little difficulty. The procedure used in the numerical calculations was as follows: A value of the angular frequency w was chosen such that w <fog = gF. This value (10 ) was chosen in order to prevent any difficulty in the solution due to the factor fl in (C. 7). An estimate for the initial value of k was obtained by the phase integral technique described in Section 2.4. For = 10 corresponding values of Re k and the phase integral were as follows: Re k Phase Integral -6 10.43 2x10-6 6.82 These approximate values of Re(k) serve as an initial value of k and also indicate the separation of the real part of the characteristic values, i. e., in this case one would expect two characteristic values to lie between 10 and 2 x 10 or rather within a similar interval since the values obtained for the real part of k are approximate. 108

THE UNIVERSITY OF MICHIGAN 3673-1 -F Thus, for w = 10 and an initial value of k = 10-6 + IO, the interaction process converged after 11 iterations to a value of k =6.86 x 10 + i2.81 x 10 In view of the values of Re(k) estimated above several new values of Re(k) were selected such as 5 xlO, 6x10, 7x107, and 8 x 10 again assuming ki = 0. In each case the numerical process converged to the same value as indicated above. Examination of the values of Re(k) estimated at the end of each integration showed that these varied over a range of 10-5 to 10-8 and thus it was probably necessary to assume a reasonable initial value of I(k) in order to separate the various modes of propagation. Thus, initial values of Re(k) were assumed as above but with i(k) = 3 x 10. This reduced the number of iterations required for convergence but only the single eigenvalue given above was obtained. In addition, although the fluctuation in the estimated values of the real part and imaginary part of k were not as great as before, they still exceeded the difference in the initial values of Re(k) assumed. Because of this, it was felt that the difficulty was due to the large discontinuity in the temperature gradient at 90 Km. This was at least partially substantiated by the fact that the integration step size was reduced from 10 Km (the initial value) to the order of 1/2 Km at an altitude of 100 Km. Therefore, it was ecided to modify the model of the disturbance which initiated the investigation discussed in the body of this report. 109

THE UNIVERSITY OF MICHIGAN 3673-1-F APPENDIX D: SOME IDEALIZED SOLUTIONS Some special solutions of the homogeneous equation describing the disturbance that may be useful in fixing the boundary conditions for the problem investigated in the text, or to explain physically the nature of the propagation of disturbances, shall be given in this appendix. The simplest problem that can be solved analytically is that of an isothermal region, i. e. A is constant. For such a medium, if losses are neglected, the equations for the normalized pressure and z-directed velocity are given by (from 2.112 and 2.113) A +\ + U[W2_W] 0 (D.1) dz yA + idU a- g + -k2 iw[dzU yAjU] +L kj 0. (D.2) Eliminating either U or, the following equation is obtained, d2 ~U [W2- 2 {(w2-w2) - A where 2 -1 g (D.4) g yA and a2 g_ y (D.5) a 4A 2~~~~~~1

THE UNIVERSITY OF MICHIGAN 3673-1-F have been defined in Section 2. 3 of this report. From (D.3) the mathematical significance of wa and wg become quite clear. The condition w=wa or W=wg defines the regions where the nature of the solution changes, i.e. from oscillatory to exponential in this case. Formally, the general solution of (D.3) may be written as U=U+a gZ k2) +u a g k2) z (D.6) U=U+ eXP~i k exp TA 02 yA z 2 where U+ and U are two arbitrary constants, and to avoid ambiguity, one may choose the branch of square root, when complex, to have a positive imaginary component. Equation (D.6) indicates that in the'forbidden region', i. e. wo > >w (D.7) a g U must be non-oscillatory, while outside this range, U can be oscillatory for some real k. Alternately, the same equation indicates that if k2 1 k2 > 1 (D.8) 2 yA then, U is exponential for all frequency ranges. For a region of the atmosphere between the ground z=O to some height z=zl, if we assume A is nearly constant and that the viscosity is negligible, then we may rite U=U sinsz (D.9) here W2_W2 | a -(w2-02)k2 (D.1O)

THE UNIVERSITY OF MICHIGAN 3673-1-F This yields i a g:!U sc Cos sz asins[ k (D.11) so that U __ _ _ _ __ _ _ _ 2 io[7A -k-2, (D.12) Z=Z1 ECS os sz- a 9sin sz and a- g -s +ssin szl+ A cossz1 ~L 1~~1 (D.13) 1 sos szl- sin. sz These relations are used as boundary conditions in the calculation of characteristic values. For a lossless medium with the parameter A varying linearly with altitude, an analytic solution for disturbances has been given by Pekeris (1948). Briefly, if A is represented by A = A'(z-za), (D. 14) then, by a change of variables = z-za, (D. 15) equations (2.106) and (2.107) may be written as w3 Y k2gW+'{Aek2-L2] = 0O (D.16)

THE UNIVERSITY OF MICHIGAN 3673 -1 -F 2 d;( 2 2 4_ 2 d2 2 TA'T + [A'k2g-2 yg +W [ k g = 0 (D. 17) The equation for, can then easily be shown to be dA x d + ) d 2 yA'k 2g+k 2g 2(-1) 2 0 (D 18) TAI' 2 + AA I -g) d~ + ~ Li 2'A = 0 (D. 18) According to Pekeris (1948) the solution of this equation can be represented in terms of confluent hypergeometric functions. In our notation, if we introduce y = 2k~ (D. 19) c=(-A ) 2k (D. 20) A4 2k (1- ) 4_Ak2g+g2(Ty-l)k2 a-4k 2. (D. 21) 4k 2 yA'(2k)w the above equation is reduced to c 2 y dy 4 Y dy If we let,I =e-j/2Y F, (D. 23) F satisfies the confluent hypergeometric equation 2F dF d + (c-y) dF- aF = 0. (D. 24) y 2dy Therefore the general solution of (D. 18) can be written as.=y1. e-1/2y o(a,c,y)+X; e1/2 Yy y),l (a+1 D. 25) where ~(a, c, y) is the standard confluent hypergeometric function.

THE UNIVERSITY OF MICHIGAN 3673-1-F In terms of'-, we then have 2 2 D. 26) W oyAl~ d'- +( (D 4 2 2 d + 2 A') w-k g L and Pog w2 yAl ~ 2t2A,2 P t2A- id 2 _2_ k2 1 22 Pog 4 2d2 (2 A (D. 27) O -k g i -k g to For the purpose of numerical evaluation, since the confluent hypergeometric functions are not extensively tabulated, one may obtain some approximate solutions directly from the differential equations. For a constant temperature gradient, (2. 99) through (2. 101) yield 1 F_1 (2-) g1 V=+ 2'+' yA' g g- -— + f = A - A g i - 1 -2; = yA'Equation (2. 102) then indicates that for a fixed o, the differential equation has a singular point at 2 + A' 1 g2(y- + A' 7 g Evidently, at ~ = 1 2 ()2 g asdefinedby(2.78). Thebehaviorof i (or; 1 since losses have been _________________ _ -114 __

THE UNIVERSITY OF MICHIGAN 3673 -1 -F neglected) near this singular point can be obtained by writing (2. 104) as 2 a a-2 a21 d + + [+a +. 2 v1 (2 -2 dz 1 where, it can easily be verified that a_2 =- 3 4 -t + /a' -g A a) By neglecting the other g terms that are small near ~ = 1' the solution of ~ is given by a Bessel function where Z2 is a Bessel function of order 2. To study the effect of the variation of viscosity on the propagation of disturbances, one may consider the idealizedproblem of a medium for which A and g are approximately constant but v varies exponentially. This model is perhaps quite realistic for the model of the upper atmosphere, say above 400 Km. Thus let v =v es(Z 2) (D. 28) The effect of losses may then be expressed in the forms s 1 -i3 a e F= -5 (D. 29) 1 -i sae

THE UNIVERSITY OF MICHIGAN 3673-1 -F 1 -i-ae l G 5 s (D. 30) 1 ie where z- z2 (D.31) and 2 c- V (D. 32) YA o Due to the involved variation of F and G with altitude, only an approximate solution is possible. Since as -- co, e becomes negligibly small, we may obtain an approximate solution of a differential equation such as (2.105) by expanding the coefficients of in terms of e-S and neglecting all the powers of e e beyond the first. Formally, then, we can write d2 U d2 (U )+(.U. ) bl+ib 2 e =0 (D. 33) where

THE UNIVERSITY OF MICHIGAN 3673-1 -F 2 5 2 g(A+g) + 5 5 k 2 A+g5 g+ 19 L A 7 2A'7YA W2 5 2A 9 Aj (D. 34) 240 2 4 (W2 -g+ ( 2g +9 A 5 22 9 L 2147 2A 15 A 2A 9 YAA- E + 2 189 ^yA c g(A J 2 2 g(A' f+g)+5g2 A 7yA (w2 gA?+g)~ )( 4 )+ +4 2 A 7-yA 49 A 2 25 147 -A 15 -2i 5S k2 ~ 2 L g74 2 4 2 g(A+g) 5 k 9 s 47 -A 5 A _y 25 a, 2 7 WA 2 5 a 2 g(A' +g) 5.2 A 7 yA -2 ~_A'+g +5 0g] 22 _ (D. 35) a L 2A 9 -yA 189 (A The solution of (D. 33) can be expressed in terms of Bessel functions of order s and argument Due to the boundedness of U1 as — >oo, we choose to use the Bessel function, and write

THE UNIVERSITY OF MICHIGAN 36 73 -1 -F U^ = J2iV ei2 (D. 36) s For large I, one may replace the Bessel function approximately by its first term and write i f1(z-z2) U1 =A exp e This yields -iA 1 d V f 1 (z-2 (D. 32) 1 dz so that U1 iw f (D. 38) Also, from (2. 97) the following ratio is easily obtained: U 1 i-V fi b f f | id e l faio b f +a — f- f zz2 These relations are used in the upper boundary condition in evaluating the mode numbers in Section 2. 4.1 ____ ___ ____ ___ ____ ___ ____ __ 18

THE UNIVERSITY OF MICHIGAN 3673-1-F REFERENCES Al, D.K. (1961), "Small Perturbations in the Unsteady Flow of a Rarefied Gas Based on Grad's 13-Moment Approximation," Rarefied Gas Dynamics, Academic Press, p.345. Astrophysics Research Corporation(1963),"Memorandum on Properties of Acoustic Waves in the Atmosphere". Barron, D.W. (1959),"The'Waveguide Mode' Theory of Radio Wave Propagation When the Ionosphere is not Sharply Bounded," Philo.Mag., 7, No. 45, pp. 1068-1081. Barry, G. (1962), "Acoustic Waves in the Ionosphere," Stanford Electronics Laboratory Report TR-50. SECRET Bhatnagar, P.L., E.P.Gross and M.Krook(1954),"A Model for Collision Processes in Gas es -I:Small Amplitude Processes in Charged and Neutral One-component Systems," Phys.Rev., 94, No.3, pp.511-525. Budden, K.G. (1962), The Waveguide Mode Theory of Wave Propagation, PrenticeHall. Chapman, S. and T. G. Cowling (1939), The Mathematical Theory of Non-uniform Gases, Cambridge University Press. Chu, C-M and J.J.LaRue (1963), "Propagation of Mechanical Disturbances in a Stratified Media(U),". Proc. ARPA-OHD Meeting, 19-20 June 1963, Stanford University Research Institute Report SRI-3-1384. SECRET Grad, H. (1959), "On the Kinetic Theory of Rarefied Gases, " Comm.Pure and Appl. Math., 2, No. 4, pp. 331-407. Hines, C.O. (1960), "Internal Atmospheric Gravity Waves at Ionospheric Heights, " Can. J. Phys., 38, p. 1441. Hirschfelder, J.O., C.F. Curtis and R.B.Bird (1954), Molecular Theory of Gases and Liquids, Wiley and Sons. Holway, L.H. (1963), "Kinetic Theory of the Neutral Disturbance," OHD Meeting of Theoretical Working Group, 20-21 Nov. 1963 at SRI, p. 61. SECRET Holway, L.H. (1960), "Self Similar Flow of Blast Waves, " Raytheon Equipment Division Report ESRN-10. Kahn, D. (1963), "Kinetic Theory of Sound Propagation in Rarefied Gases, " OHD Meeting of Theoretical Working Group, 20-21 Nov. 1963 at SRI,. SECRET Landau, L.(1946), "On the Vibrations of the Electron Plasma, "J.Phys.USSR, 10 p.25. 119

THE UNIVERSITY OF MICHIGAN 3673-1-F Lighthill, M.J. (1956), "Viscosity Effects in Sound Waves of Finite Amplitude, " Surveys in Mechanics (Ed. C.K.Batchelor), Cambridge University Press, p. 250 Pekeris, C.L. (1948), "The Propagation of a Pulse in the Atmosphere-It, "Phys. Rev., 73, pp.145-154. Press, F. and D. Harkrider(1962), "Propagation of Acoustic-Gravity Waves in the Atmosphere, " J. Geophys.Rev, 67, No. 10, pp. 3889-3908. Sommerfeld, A. (1949), Partial Differential Equations in Physics, Academic Press. U.S.Weather Bureau(1962), U. S. Standard Atmosphere, 1962, U. S. Government Printing Office, Washington, D. C. Weston, V.H. (1961), "The Pressure Pulse Produced by a Large Explosion in the Atmosphere, " Can. J. Phys., 39, pp. 993-1009. 120

TTnrlss ifieP. Security Classification DOCUMENT CONTROL DATA- R&D (Security classification of title, body of abstract and indexing annotation must be entered when the overall report is classified) 1. ORIGINATING ACTIVITY (Corporate author) 2a. REPORT SECURITY CLASSIFICATION Unclassified The University of Michigan GROUP Radiation Laboratory 3. REPORT TITLE IONOSPHERIC DISTURBANCES AND THEIR EFFECTS ON THE PROPAGATION OF HF ELECTROMAGNETIC WAVES: FINAL REPORT 4. DESCRIPTIVE NOTES (Type of report and inclusive dates) Final Report June 1960 - June 1964 S. AUTHOR(S) (Last name, first name, initial) Chu, Chaio-Min; Hok, Gunnar; LaRue, John J. 6. REPORT DATE 7a. TOTAL NO. OF PAGES 7b. NO. OF REFS September 1964 120 20 8a. CONTRACT OR GRANT NO. 9a. ORIGINATOR'S REPORT NUMBER(S) AF 19(604)-7257 3673-1-F b. PROJECT NO. 4649 c. TASK s9b. OTHER Pt pORT NOS) (Any other numbers that may be 464901 assigned aths report) d. AFCRL-64- 916 to. JVAILABILITY,/ LMITATION NOTICES Agencies of the Department of Defense and other Government agencies may obtain copies of this report from DDC. Other persons and organizations should apply to U. S. Department of Commerce Office of Technical Services. 11. SUPPLEMENTARY NOTES 12. SPONSORING MILITARY ACTIVITY Air Force Cambridge Research Labs. Laurence G. Hanscom Field Bedford, Massachusetts 13. ABSTRACT This report is concerned with the effect of ionospheric disturbances on HF radio waves whose propagation path is at a large distance from the source of the disturbance. The perturbation of the HF radio waves is assumed to be due to mechanical waves propagating in the atmosphere. The macroscopic equations governing the propagation of mechanical waves in the atmosphere, including the effects of viscosity and thermal conductivity, are developed. A simple model of the atmosphere is proposed and a technique, suitable for numerical analysis, is developed for determining the various modes propagating. The effect of the earth's magnetic field on the propagation of mechanical waves is evaluated and it is concluded that this effect is small. The modification of the disturbance due to the presence of charged particles is discussed in conjunction with a comparison of the microscopic and macroscopic formulation of the equations of motion. An idealized description of the source is given in the appendix and some solutions to some simple problems are discussed. 1DD JAN 6473 Security Classification

UNIVERSITY OF MICHIGAN I IrIi II iII rlJ 111 11111 lll ll 3 9015 02828 5651 Unclassified Security Classification 14, LINK A LINK B LINK C KEY'WORDS ROLE WT ROLE WT ROLE WT Ionospheric 5,6 1 Propagation 8 1 Electromagnetic 10 1 Disturbances 5, 6 1 5, 6 1 HF 10 1 Viscosity 5, 6 2 Perturbation 5, 6 1 5,6 1 Detection 4 2 INSTRUCTIONS 1. ORIGINATING ACTIVITY: Enter the name and address 10. AVAILABILlTY/LIMITATION NOIICES: Enter any limiof the contractor, subcontractor, grantee, Department of tations on further dissemination of the report, other than those Defense activity or other organization (corporate author) imposed by security classification, using standard statements issuing the report. such as: 2a. REPORT SECURITY CLASSIFICATION: Enter the over- (1) "Qualified requesters may obtain copits of this all security classification of the report. Indicate whether report from DDC." "Restricted Data" is included. Marking is to be in accord- (2) "Foreign announcement and dissemination of this ance with appropriate security regulations. report by DDC is not authorized." 2b. GROUP: Automatic downgrading is specified in DoD (3) "U. S. Government agencies may obtain co ies of Directive 5200. 10 and Armed Forces Industrial Manual. this report directly from DDC. Other qualified DDC Enter the group number. Also, when applicable, show that users shall request through optional markings have been used for Group 3 and Group 4 as authorized. (4) "U. S. military agencies may obtain copies of this 3. REPORT TITLE: Enter the complete report title in all report directly from DDC. Other qualified users capital letters. Titles in all cases should be unclassified. shall request through If a meaningful title cannot be selected without classification, show title classification in all capitals in parenthesis " immediately following the title. (5) "All distribution of this report is controlled. Quali4. DESCRIPTIVE NOTES: If appropriate, enter the type of fied DDC users shall request through report, e.g., interim, progress, summary, annual, or final.., Give the inclusive dates when a specific reporting period is covered. inclusive dates when a specific reporting period is If the report has been furnished to the Office of Technical ~~~~~~~~~~~covered. ~Services, Department of Commerce, for sale to the public, indi5. AUTHOR(S): Enter the name(s) of author(s) as shown on cate this fact and enter the price, if known. or in the report. Enter last name, first name, middle initial. If military, show rank and branch of service. The name of 11. SUPPLEMENTARY NOTES: Use for additional eplanathe principal author is an absolute minimum requirement. tory notes. 6. REPORT DATE: Enter the date of the report as day, 12. SPONSORING MILITARY ACTIVITY: Enter the name of 6. REPORT DATE: Enter the date of the report as day, the departmental project office or laboratory sponsoring (paymonth, year, or month, year. If more than one date appears the departmental project office or laboratory sponsoring (pa on the report, ormonth'year.Ifmrenuse date of publication. ing for) the research and development. Include address. 7a. TOTAL NUMBER OF PAGES: The total page count 13. ABSTRACT: Enter an abstract giving a brief and factual 7a. TOTAL NUMBER OF PAGES: The total page count summary of the document indicative of the report, even should follow normal pagination procedures, i.e., enter the sugmary of the document indicative of the report, even number of pages containing information. though it may also appear elsewhere in the body of the technica report. If additional space is required, a continuation 7b. NUMBER OF REFERENCES: Enter the total number of sheet shall be attached. references cited in the report. It is highly desirable that the abstract of classified re8a. CONTRACT OR GRANT NUMBER: If appropriate, enter ports be unclassified. Each paragraph of the abstract shall the applicable number of the contract or grant under which end with an indication of the military security classification the report was written. of the information in the paragraph, represented as (TS), (S), 8b, 8c, & 8d. PROJECT NUMBER: Enter the appropriate (C), or (U). military department identification, such as project number, There is no limitation on the length-of the abstract. Howsubproject number, system numbers, task number, etc. ever, the suggested length is from 150 to 225 words. 9a. ORIGINATOR'S REPORT NUMBER(S): Enter the offi- 14. KEY WORDS: Keyv words are technically meaningful terms cial report number by which the document will be identified or short phrases that characterize a report and may be used as and controlled by the originating activity. This number must index entries for cataloging the report. Key words must be be unique to this report. selected so that no security classification is required. Identi9b. OTHER REPORT NUMBER(S): If the report has been fiers, such as equipment model designation, trade name, miliassigned any other report numbers (either by the originator ary project code name, geographic location, may be used as assigte anysothe rport enerts (eithbers key words but will be followed by an indication of technical or by the spnsor), also enter this number(s). context. The assignment of links, rules, and weights is optional. Unclassified Security Classification