THE UNIVERSITY OF MICHIGAN AFCRL-746 4134-2-F STUDIES IN NON-LINEAR MODELING - III ON THE INTERACTION OF ELECTROMAGNETIC FIELDS WITH PLASMAS by Kun-Mu Chen R.E. Kleinman J. Meixner O.G. Ruehr D. Sengupta 31 July 1961 Final Report The work described in this report was partially supported by the ADVANCED RESEARCH PROJECTS AGENCY ARPA Order Nr. 147-60, Project Code Nr. 7600 ARPA Order Nr. 147-60, Project Code Nr. 7600 Contract AF 19(604)-7428 prepared for Electronics Research Directorate Air Force Cambridge Research Laboratories Office of Aerospace Research Laurence G. Hanscom Field Bedford, Massachusetts - -

--- THE UNIVERSITY OF MICHIGAN 4134-2-F Requests for additional copies by Agencies of the Department of Defense, their contractors, and other Government agencies should be directed to: ARMED SERVICES TECHNICAL INFORMATION AGENCY ARLINGTON HALL STATION ARLINGTON 12, VIRGINIA Department of Defense contractors must be established for ASTIA 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 4134-2-F TABLE OF CONTENTS ABSTRACT I. INTRODUCTION II. PLASMA SHEATH FORMED BY A STATIONARY PLASMA ON AN INFINITE PLATE iv 1 6 III. PLASMA SHEATH FORMED BY A MOVING PLASMA ON AN INFINITE PLATE IV. THE INTERACTION OF A HIGH INTENSITY ELECTROMAGNETIC FIELD WITH A LOW DENSITY PLASMA V. THE INTERACTION OF A HGIH INTENSITY ELECTROMAGNETIC FIELD WITH A WEAKLY IONIZED GAS VI. THE ELECTRICAL CONDUCTIVITY OF A PARTIALLY IONIZED GAS VII. NON-LINEAR ELECTRICAL CONDUCTIVITY OF A FULLY IONIZED GAS APPENDIX A: NON-LINEAR MODELING OF BESSEL FUNCTIONS APPENDIX B: LOCAL ALTERATION OF ATMOSPHERIC DENSITY WITH ELECTROMAGNETIC ENERGY 18 31 51 74 90 107 125 iii

THE UNIVERSITY OF MICHIGAN 14134-2-Fl ABSTRACT The results of an investigation of the interaction of electromagnetic fields with plasmas, including the development and use of appropriate non-linear modeling techniques, are presented. Expressions for density distributions in the sheath and potential induced on a conducting plate in the presence of a plasma (moving or stationary) are found. The potential at one temperature in the sheath is non-linearly modeled to give the potential at another temperature. Correction terms are found for the standard treatments of the interaction of electromagnetic fields with collisionless plasmas and weakly ionized gases. These corrections are functions of field intensity and give significant results for high intensity fields. The current induced in the plasma with a low intensity field is non-linearly modeled to the high intensity case. Expressions for the conductivity (both a. c. and d. c.) are found for various plasmas and incident field intensities. iv

THE UNIVERSITY OF MICHIGAN - 4134-2-F I INTRODUCTION This is the final report on contract number AF 19(604)-7428, ARPA Order No. 147-60, concerning an investigation of the interaction of electromagnetic fields with plasmas. An understanding of this phenomenon is essential in treating many of the vital problems associated with the coupling of high intensity electromagnetic energy and plasmas where it is desired to optimize the attenuation of electromagnetic waves by a plasma or to minimize this attenuation without resorting to brute force techniques. This will allow us to optimize systems for telemetering through plumes and re-entry plasmas and also to optimize radar systems for terminal guidance of maneuverable I. C. B. M. 's. In attempting to study these physical problems in the classical tradition by combining knowledge gained from both theoretical and experimental investigations to construct a physical picture, one is confronted with two major obstacles. On the one hand the theoretical problem, the coupling of the Boltzmann equation with the Maxwell equations, is extremely difficult in general and despite much activity in this field, results even in the simplest idealized situations are quite rare. On the other hand useful experimental results are equally, if not more, difficult to obtain. Full scale laboratory experiments of physical importance are, because of the large field strengths, if not impossible, very costly. Furthermore linear modeling, the usual technique for scaling physical phenomena to dimensions feasible for r I. 1

i THE UNIVERSITY OF MICHIGAN 4134-2-F laboratory experiments, is inapplicable since the basic particles of a plasma cannot be modeled correctly. If reasonable laboratory experiments are to yield significant full scale information for high intensity electromagnetic fields this must necessarily be accomplished by application of suitable non-linear modeling techniques. The concept of non-linear modeling as developed in [13 through t(3 provides the basis for theoretical investigations which will enable data obtained in low power laboratory experiments to be transformed to give results in cases of physical significance. In order to develop the non-linear modeling techniques necessary in this problem, work in the Radiation Laboratory under this contract has proceeded by considering a number of special problems involving the interaction of plasmas with electromagnetic fields. These problems were chosen for a variety of reasons. In the first place they are of considerable interest in themselves. Secondly they serve, in a sense, as canonical problems on which the non-linear modeling techniques can be developed and demonstrated. Thirdly the theoretical solutions can be found and used to guide and check the modeling attempts. These solutions themselves represent new results in the field and are presented here for the first time. In chapter II, the plasma sheath formed by a stationary plasma on the surface of an infinite plate is considered. With the assumption of a low density plasma and a Maxwellian distribution of both ions and electrons up to the plate, exact expressions for the density distributions of electrons and ions in the sheath - - 2

THE UNIVERSITY OF MICHIGAN 4134-2-F and the potentials induced in the plate and the sheath are found. The potential in the sheath for one temperature is also non-linearly modeled to yield this data at other temperatures. Chapter III treats the same problem when the plasma, instead of being stationary, moves toward the plate with constant mainstream velocity. Although the sheath formed is in many respects similar to that formed in the stationary case, the analysis is more complicated, hence approximate rather than exact expressions for the potential and density distributions are found. These expressions include a first order correction to results available in the literature. In chapter IV, the interaction of a high intensity electromagnetic wave with a collisionless plasma is investigated. The coupled Boltzmann and Maxwell equations are solved without small signal approximations. The zeroth-order analysis yields conventional results; however, the first-order analysis produces new results for the velocity distribution function and quantities derived from it, i.e. conductivity, permittivity, current, and energy density of the plasma. These are evaluated as functions of the intensity of the EM wave. An attempt to use nonlinear modeling in this problem is unsuccessful and the difficulty encountered is explained. The interaction of a high intensity EM wave with a weakly ionized gas is studied in chapter V. A rigorous mathematical treatment is presented with a simple collision model. The zeroth-order analysis gives the plasma parameters 3

THE UNIVERSITY OF MICHIGAN 4134-2-F as functions of the intensity of the incident wave. The results obtained can be reduced to the well-known results when the intensity of the EM wave is assumed to be very small. A first-order analysis is also presented. Non-linear modeling of the current induced in the plasma by the electric field is found to be possible in this case and the modeling function is presented explicitly. In chapter VI the electrical conductivity of a low density and partially ionized gas where both electron-neutral particle and Coulomb type collisions play important roles is discussed. Assuming that the ionized gas is perturbed by a weak electric field, the velocity distribution function for the electrons is obtained by solving the Boltzmann equation; the collision between neutral and charged particles is accounted for by the hard sphere model for the particles while the collision between the charged particles is taken care of by the Fokker-Planck equation. Explicit expressions for a.c. and d.c. conductivity are given for various cases. To the extent that the assumptions made for the collision models are valid, the expressions for conductivity are quite general and can be used for any degree of ionization. The conductivity of a fully ionized gas for arbitrary electric field intensity is investigated in chapter VII. In the d. c. case an instability phenomenon, the runaway effect, limits the intensity of the nelectric field that may be applied and under this restriction an expression for the conductivity is obtained. In the a. c. case, above a certain critical frequency, the problem of finding an expression l id IIII J - _ I I I I I I _ 4

-- THE UNIVERSITY OF MICHIGAN 4134-2-F for the conductivity becomes one of solving a non-linear first-order differential equation. Approximate expressions are presented in this case. Not included here is the theoretical treatment of propagation of electromagnetic waves in a plasma column surrounded by an annular isotropic medium in a cylindrical wave guide 016]. With the source of electromagnetic waves taken to be a thin ring source, the formal solution has been obtained for arbitrary angular variation of the current intensity and the case of constant intensity is treated at length as a special case. In this case the excited fields will also be axially symmetric and this problem is solved where the anisotropy of the medium is taken into account by considering both the permittivity and permeability to be tensors. The most general problem treated in [16] is too complex to be amenable to non-linear modeling techniques at present. However, suitable techniques have been developed to handle some limiting cases where the physical situation is described by a Bessel function. The details of this analysis appear in Appendix A. A problem of considerable interest involving the interaction of electromagnetic waves with a plasma to alter the density of the atmosphere is discussed in Appendix B. Calculations of the power needed, both on the ground and at various altitudes in the atmosphere, to effect particular temperature (hence density) variations are presented. Appreciation is extended to R. W. Larson for these calculations. ] - 5 --

I I I THE UNIVERSITY OF MICHIGAN 4134-2-F II PLASMA SHEATH FORMED BY A STATIONARY PLASMA ON AN INFINITE PLATE In this chapter the plasma sheath formed by a stationary plasma on the surface of an infinite plate is considered. This is a problem frequently met in plasma physics and has been studied by several authors (e.g. (4] ), however, the solutions are mostly approximate. It is the purpose of this study to treat the problem more rigorously and an exact solution is attempted. The problem is then examined from the non-linear modeling viewpoint. Assume a stationary, uniform plasma fills the half space with an infinite plate located in the x - y plane as the boundary. If the potential of the plate is allowed to float, the plate gets charged and the electrical neutrality of the plasma is not preserved in the vicinity of the plate. Due to the fact that the root mean square velocity of the electrons is much higher than that of the positive ions (assumed to be singly charged), the plate must be charged negatively in order to adjust the amount of the electrons and the ions which hit the plate per unit time at equilibrium. This negative potential of the plate results in a low electron density and a high positive ion density in the vicinity of the plate. The deviation from the uniform densities of the electrons and the ions, or the deviation from the electrical neutrality, dies out in the direction away from the plate. The thickness of the layer where the electrical neutrality 6 6 U

I THE UNIVERSITY OF MICHIGAN 4134-2-F deviates, or the plasma sheath, and the density distributions of the electrons and the ions in the plasma sheath are found. Distribution Functions The velocity distribution function of the electrons in a plasma is m c2 3/2 - e m / 2 KT f =ne 2iKT ) e (1) e e 27r KT the velocity distribution function of the positive ions in the plasma can be written as mic m. 3/KT 1 f i= n 2 ) e (2) where n and ni are the number densities and can be expressed as functions e 1 of the space. In general, we can write f f + f (3) e eo e f. = f. + f.' (4) 1 10 1 where m c2 3/2 - e /m 2 KT e e eo no 27r KT e (5) e m.c2 m=. 3/2 KT. f. - n e2 KT ) (6) 10o0 27T KT. 1 I - -- -- 7

THE UNIVERSITY OF MICHIGAN 4134-2-F f ' and f.' are the perturbed terms and n is the unperturbed density of the e i 0 plasma. If the velocity distribution functions for the electrons and the ions are assumed to be Maxwellian at any point in the physical space, the Boltzmann equations lead to af e + c. Vf =0 (7) at e m c e e af. + c * Vf. + (8) at 1 m c 1 1 Since the plasma is stationary, the time derivative terms should drop out. Due to the one-dimensional geometry, the variation is only in the z-direction. Therefore af e c Vf c e e z az '_ af af E Vf e E E f z = E ce ac ac z z 3f 2 m E e 3c _ e E c f = 2 3c KT z e ac z e With these relations, Equation (7) becomes af e eE + e f =0 (9) az KT fe ( e - - 8

THE UNIVERSITY OF MICHIGAN 4134-2-F Equation (9) can be solved subject to the boundary condition f -- f as e eo x -* co, and the relation E = - 0_ = - A. The final solution for f is z ' e e0 3/2 mc2 e0 KT m 2 KT e e_ e e f =f e n e e (10) e eo o 27r KT e Where 0 is the potential in the space and is normalized to zero at infinity. Thus the velocity distribution function of the electrons is determined as a function of space and velocity. Following similar reasoning, Equation (8) can be reduced to af. 1 eE az KT. i 1 and the final solution for f. is 1 m 2 - ~3/ 2 -i KT. m. / 2KT. KT. f. = f. e e (12) I 10 o 27r KT. Equation (12) gives the velocity distribution function of the positive ions as a function of space and velocity. It follows from Equations (10) and (12) that the number densities of the electrons and the positive ions can be expressed as follows: e0 KT n n e e (13) n -- n e (13) e o 9

THE UNIVERSITY OF MICHIGAN 4134-2-F _ e KT. n. = n e 1(14) 1 o Potential Distribution in the Plasma Sheath The number densities of the electrons and the positive ions in the plasma sheath are thus determined as functions of the potential in the space. It is now necessary to evaluate the potential distribution in the plasma sheath. The potential in the plasma sheath is the solution of a Poisson's equation as follows: 02 - e (n. -n) ~ 1 e o _ e e n e KT KT o 1 e = - (e e ) (15) 0 o Assume T = Ti= T for simplicity, and owing to the one dimensional variation, e i Equation (15) is reduced to the differential equation 2 2n e = -h sinh ( ).(16) dz o The boundary conditions for 0 are 0 = 0 at z = O (on the plate) 0=0, and = 0 as z - o aaa Equation (16) can be solved as follows: I 10

THE UNIVERSITY OF 4134-2-F MICHIGAN Rewrite Equation (16) as 2 = a sinh (B0) dz (17) with 2n e 0 o E 0 e KT Let = p() dz 2 dz2 dpd Pd ' then Equation (17) becomes p d -= a sinh (30) 2 p- - cosh 30 + c or d + 1 dz 2 cosh 0 + 2c - oosh 10 + 2c (31 (18) Subject to the boundary condition of -~~ 0 as z - oo (or as 0 - 0) dz ai c is found to be - - 1 P Equation (18) then becomes dz + -J 13 cosh 10 - 1 -+ 2V sinh 2 2 (19) -- -- ~~ 11

I THE UNIVERSITY OF MICHIGAN --- 4134-2-F 0 is always negative and is always positive in the plasma sheath. Theredx fore, the negative sign in Equation (19) is taken. The integration of both sides of Equation (19) leads to 0; dt 2_Vd sinh - 2 z 00 or or - logtanh - log tanh 2 0 4 4 This can be simplified to tanh -- (tanh ) e a 4 4 or 0 = - tanh1 e tanh (20) 2n e o e The substitutions of c -- and 3 -- in Equation (20) give a final solution e KT 0 for 0 as follows: 2n e2 O e KT Z 4KT 1 V 4KT tanh e tanh 7 (21) n 4KT The solution of 0 can be checked easily because 0 -- 0 as z -+ o, and = 0 as z = 0. o I 12

THE UNE IVERSITY OF MICHIGAN 4134-2-F Potential of the Plate The potential distribution in the plasma sheath is found as a function of some parameters and the potential of the plate, 0. 0 can be determined from the condition that the same amount of electrons as positive ions will hit a unit area of the plate per unit time at equilibrium. Owing to the fact that the rms velocity of the electrons is much higher than that of the positive ions, more electrons than ions tend to hit the plate per unit time except the plate is charged negatively so that only high energy electrons can reach the plate. Assume the potential of the plate as 0, which is negative. A critical velocity of the electron is defined as 1 2 meCe = ce 2 e eo 0 or 2 Ie0oI C = (22) eo m e c means that only those electrons whose velocities are higher than c and eo eo point toward the plate can overcome the potential barrier on the plate and reach the plate. Similarly a critical velcoity for the positive ion is defined as c,.= \ e (23) lli 13

I I THE UNIVERSITY OF MICHIGAN 4134-2-F Thus, the ion having a velocity lower than ci and points away from the plate may be attracted back to the plate. The application of the boundary condition at equilibrium yields an equation as follows: C eo Cio oo - c d ccd \c d dc \ d c +\ c d dc [f j z xo L d zO 00 -00 -00 00 0 o00 The integrals can be carried out to be I eol le0 l n e = n 2 (2 - e KT o V 27rm o V 27rm. e i or |e0 |M F 1 -- = log + (24) KT 2 It is learned that the potential of the plate, 0, is a function of the temperature o and the mass ratio of ion and electron. Density Distributions of the Electrons and the Ions in the Plasma Sheath Up to this point the potential distribution in the plasma sheath is completely determined. The density distributions of the electrons and the ions in the plasma sheath are obtained by substituting the potential 0 of Equation (21) I - -- - 14

THE UNIVERSITY OF MICHIGAN 4134-2-F in Equations (13) and (14). The final solutions are ne = n exp 4 tanh (e I tanh 4KT ) (25) 2n e2 0 - eKT e0 e KT n. = n exp - 4 tanh (e o tanh (26) It is noted that 0 is negative. From Equations (25) and (26) it is observed that both n and ni approach n quite rapidly as z increases. The thickness of the plasma sheath, H, is obtained by finding a value of z at which n is 0. 95 n eo~4KT~e (arbitrary), and after that point n and n. are very close to n. c? 10 H= cKT log tanh 4KT 2H = lo ^ tanh(0. 0128) ( 2n e 2 0 to be as follows Veo s an exam3.08, a pasa 0. 266 volt H =1 0. 605 cm n = 10 1/ considered. In this case, the numer ical results are found. to be as follows Ie0O -- = 3.08, 0 =-0.266volt, H= 0.605 cm. r~~~~~~~i 0~~~~~~~~~ 15

THE UNIVERSITY OF MICHIGAN 4134-2-F Non- Linear Modeling the Potential Distribution Equation (16) could also serve as the starting point in a non-linear modeling attempt, where the potential distribution for two different temperatures is sought, with but one experiment. This equation is a representative of a general class of equations treated in some detail in W. Specifically, simplifying the equation somewhat as in (17) we have the two equations 2 a - a sinh (130) = (28) dz 2 d2 - c sinh (i32) = 0 (29) dz where the only difference is a change in the constant f3. We wish to find 0 as a function of W. With W3 we can immediately write, without employing any boundary condition do d= (30) Sa sinh/31 0d0 a dinh /32~d~ or, do = (31) o1 / cosh 31O+c i2 cosh32V+ c where cl and c2 are constants of integration. In general the integrals in (31) will be elliptic but in the special case when c l c2= -1, (31) becomes * 16

I THE UNIVERSITY OF MICHIGAN - 4134-2-F J ds 2 ds sinh sinh 2 2 (32) (33) or 3 1 log tanh -= - log tanh 2 + c pi 4 f2 4 where c is a constant of integration. To evaluate c we tion of the form I- = i/ where 0 = 0 hence 1 0o c = -, log tanh - log tanh - P/ 4 P2 apply an initial condi 12%o 4 (34) Substituting (34) in (33) yields 110 tanh 1 4 1 log - tanh 4 tanh 2 1 log tanh 4 (35) or explicitly 2 2 tanh 4 Thus in this case it is possible to find a non-linear modeling function and in fact this does correspond to the physically significant case as can be verified by eliminating z from (20) and a similar equation for different values of 3. 17

THE UNIVERSITY OF MICHIGAN 4134-2-F III PLASMA SHEATH FORMED BY A MOVING PLASMA ON AN INFINITE PLATE In this chapter the plasma sheath formed by a moving plasma on an infinite plate is investigated. The significance of this study is to investigate the behavior of a moving plasma on a plane boundary or the behavior of a stationary plasma on the surface of a moving plate. This analysis may be useful in the investigation of the antenna on a space vehicle which is wrapped by a plasma sheath. The analysis is approximate, however some significant correction terms which are ignored in the conventional studies are evaluated. Assume a moving plasma, having a mean stream velocity V, moves toward an infinite plate located in the x-y plane. If the potential of the plate is allowed to float, the plate gets charged negatively and a plasma sheath forms on the plate. In this analysis, the distribution functions and the densities of the electrons and the positive ions and the potential distribution in the plasma sheath are obtained. Distribution Functions A moving plasma moves toward an infinite plate which is located in the x-y plane. The unperturbed velocity distribution functions for the electrons and the positive ions can be expressed as I - --- -- 18

THE UNIVERSITY OF MICHIGAN 4134-2-F m e 2 2KT /- [^ ] e 1 (2) f.t = n (/. 2KT (2) with V = -Vz. As the moving plasma hits the plate, the distribution functions are disturbed in the vicinity of the boundary. In general, the distribution functions can be written as f = f + f e eo e n \H~e (2<)~(3) f. = f. + f.' 1 10 i f ' and f.' represent the perturbed terms. f and f. can be determined from e 1 e 1 the Boltzmann transport equation. In the case of a moving plasma, f and f. e 1 are functions of space and velocity. Moreover, the dependence on the space and the velocity coordinates are found inseparable. Therefore, an approximate approach is devised to determine f and f.. e 1 It is necessary to make some reasonable assumptions here. (1) As the moving plasma hits the plate, there occurs only the diffuse reflection of the surface. The specular reflection is ignored. (2) Equal quantities of the electrons and the positive ions hit the unit area of the plate per unit time at equilibrium. " ' '..-... --- " "19

THE UNIVERSITY OF MICHIGAN 4134-2-F (3) The neutralization taking place on the plate is assumed to be complete. This implies that those positive ions and the electrons reflected from the plate are assumed to be neutralized. (4) The plate is charged negatively. With these assumptions the pile up of the reflected particles in the front of the plate can be ignored because they are neutralized particles. Use the Boltzmann transport equation, af Vf+ - vf y) at m C St (4) ol.. and assume the distribution functions as m / m \3/2 - e (' V)2 m 2KT f n (z) e (5) e e fe KT \ e/ m..2 1 1 m.(2 T \3/2 - - V f. = n. (z ) e (6) i 1 7r KT Note that Equations (5) and (6) are approximate expressions, because the space and the velocity coordinates are not separable in this case. The substitution of Equations (5) and (6) in (4) for the steady state give af e eE V + +- (1 + )f =0 (7) az KT c e e z af. 1i eE V (+ )f.= 0 (8) az KT. (1+ c 1(8 1 z _ — - ---- -20

THE UNIVERSITY OF MICHIGAN 4134-2-F The integration of Equation (7) in the velocity space yields an equation for n if f is approximated with f in the third term of (7). That is e 3 eE 3 eE V 3 ed c + -- \f dc -- c az KT e KT c' eo e e z or an eEn e + e eE V f d3 "a + KT-f d c (9) az KT KT c eo e e z If the r. h. s. of (9) turns out to be small, the approximation is valid. It is found to be V' 3f 3 -V'2 t2 f d c = - n 2V e dt n (10) f c eo o e o o 0 where V' = e V 2KT e In the usual case, ye turns out to be very small. Therefore, the approximation made in Equation (9) is valid. Now (9) can be solved subject to the boundary condition of n -> n as z -oo, and the relation E = as follows: e o az KT e eE n = n e KT (11) e o oe KT e Equation (11) gives the electron density as a function of space because 0 is a function of space. Although (11) is an approximate solution it carries an extra 21 -- --- ----

THE UNIVERSITY OF MICHIGAN 4134-2-F correction term which is usually ignored. For the case of the positive ions, the approximation made in solving (7) is not valid because the effect of the third term of (8) is enormous. Equation (8) is solved approximately as follows: f i 3 * eE IV 3 _ eE (1+ | d c(1 + ) fi d c J zKT. c 10 or an! i. eE eE V f.3 d1 az KT. n + io (12) 1 1 Z It can be shown that 2 set2dt - Y n (13) VTT V 3 -V"2 t2 V f. d- - n 2V" e e dt -7in (13) 10 o L 0 0 where V" = V Yi turns out to be near unity in the usual case. Subject to the boundary condition of n. -4 n as z -> oo and the relation E = - the perturbed term of the 1 o z positive ion, n.', can be determined as n.' = - (1 - ) n (14) iKT. 1 1 KT, (~1 y n 22 1

I I I THE UNIVERSITY OF MICHIGAN 4134-2-F Therefore, n. = n. + n' = n n( (15) i 10 i o o i KT. 1 Equation (15) gives the density of the positive ions as a function of space. The factor, (1 - -yi) is very close to zero for V" greater than 1. Therefore, the second term of (15) serves as a correction term. To get an idea about the magnitude of ye and.,i these are calculated for the case of T = 1000 K, and 'e m = (16 x 1825) m. The results are: Y = 3.2 x 10 and. _i 1 for i e e e -5 V = 7Km/sec, ye = 6 x 10 and yi =0. 97 for V = 1 Km/sec. This implies that when V is around or higher than the rms velocity of the positive ions, but much lower than that of the electrons, the assumption of n = n e KT and e o ni = n is quite accurate. This is equivalent to stating that at this mean stream velocity of the moving plasma the potential of the boundary gives a small effect on the distribution of the positive ions but the distribution of the electrons is entirely governed by the potential of the boundary. For the mean stream velocity of the moving plasma lower than the rms velocity of the positive ions the correcx 2 tion terms should be taken into account. e dt is tabulated to facilitate computations. (See Table I) Potential Distribution in the Plasma Sheath The density distributions of the ions and the electrons have been determined as a function of the potential or of the space in the preceding section. The potential distribution in the plasma sheath can be found by solving a Poisson's I 23

I I THE UNIVERSITY OF MICHIGAN 4134-2-F TABLE I x t2P xe dt = AxlOP x A P x A P x A P 0 0 1.70 6.7035 75983 0 4.8 1.0811 27786 9.05 5.0041 6992 -2 1.75 7.6850 89817 0 4.9 2.7909 92389 9.10 1.0033 43383 -1.80 8.8543 99688 0 5.0 7.3541 89348 9.15 1.5113 26386 -1 1.85 1.0254 42272 1 5.1 1.9778 61774 10.20 2.0269 89793 -1 1.90 1.1939 08605 1 5.2 5.4291 66458 10.25 2.5530 74677 -1 1.95 1.3976 42572 1 5.3 1.5210 31469 11.30 3.0924 83086 -1 2.0 1.6452 62808 1 5.4 4.3491 39612 11.35 3.6483 25877 -1 2.1 2.3190 52389 1 5.5 1.2691 74669 12.40 4.2239 76023 -1 2.2 3.3452 50713 1 5.6 3.7799 45897 12.45 4.8231 28889 -1 2.3 4.9398 02230 1 5.7 1.1489 19440 13.50 5.4498 71235 -1 2.4 7.4676 21621 1 5.8 3.5639 15749 13.55 6.1087 61080 -1 2.5 1.1556 02507 2 5.9 1.1282 19327 14.60 6.8049 20923 -1 2.6 1.8302 26303 2 6.0 3.6448 74848 14.65 7.5441 47374 -1 2.7 2.9659 42341 2 6.1 1.2016 82417 15.70 8.3330 40927 -1 2.8 4.9165 86436 2 6.2 4.0430 59123 15.75 9.1791 60408 -1 2.9 8.3346 99927 2 6.3 1.3881 59382 16.80 1.0091 20769 0 3.0 1.4445 45766 3 6.4 4.8637 72076 16.85 1.1079 24967 0 3.1 2.5591 06616 3 6.5 1.7390 32700 17.90 1.2154 98595 0 3.2 4.6331 25068 3 6.6 6.3451 09838 17.95 1.3332 07308 0 3.3 8.5706 33926 3 6.7 2.3624 56410 18 1.00 1.4626 51863 0 3.4 1.6197 22651 4 6.8 8.9759 06314 18 1.05 1.6057 16168 0 3.5 3.1268 14157 4 6.9 3.4800 01325 19 1.10 1.7646 26158 0 3.6 6.1652 40644 4 7.0 1.3767 78248 20 1.15 1.9420 22025 0 3.7 1.2414 92189 5 7.1 5.5581 34208 20 1.20 2.1410 47271 0 3.8 2.5529 88287 5 7.2 2.2896 64923 21 1.25 2.3654 58828 0 3.9 5.3608 52644 5 7.3 9.6247 57736 21 1.30 2.6197 63680 0 4.0 1.1494 02557 6 7.4 4.1283 84333 22 1.35 2.9093 88950 0 4.1 2.5161 68785 6 7.5 1.8069 26236 23 1.40 3.2408 94369 0 4.2 5.6236 11905 6 7.6 8.0699 03889 23 1.45 3.6222 38609 0 4.3 1.2831 56637 7 7.7 3.6775 76227 24 1.50 4.0631 14269 0 4.4 2.9889 27772 7 7.8 1.7100 87084 25 1.55 4.5753 70719 0 4.5 7.1073 35497 7 7.9 8.1140 17191 25 1.60 5.1735 49665 0 4.6 1.7252 00132 8 8.0 3.9283 73639 26 1.65 5.8755 65931 0 4.7 4.2746 46980 8 L __ __ 24

-- THE UNIVERSITY OF MICHIGAN - 4134-2-F equation as follows: V20 _ e (n.-n) ~ i e 0 e eno KT e e 0 e0 'e KT e 1+ - i) KT ) (16) For simplicity, let T = Ti = T, (1 - yi- ye) = y. e i 1 e Equation (16) becomes 2 O _ az2 eno KT -'le o \ -l+ KT 'KT (17) The boundary conditions for 0 are 0 = at z = 0 (on the plate) 0 0 as z- oo. Equation (17) can be solved as follows: 2 o - a(e -1+ S0) az (18) with en 0 a -- 0 e -KT Let az. 2 az2 dP dO f _ ~ ~ ~ ~ ~ 25 --

i THE UNIVERSITY OF MICHIGAN 4134-2-F then Equation (18) becomes P d a(e - 1+ S0) 1 2 1 30 + 1 2) 2 2 1 or = \/2a( 1 e-z0 - 0 + 2C i (19) dz o 2 1 The arbitrary constant C1 can be determined subject to the boundary condition of - O when z-> oo (or as 0 -- 0) as az C1 - - Therefore, (19) yields d d dz = 2-(20) \/a t(e 3 1)- 0+ 2 so and 1 dt F Ax\2 V (e t- 1) -t +1 t2 C2 0VP 2 The abritrary constant C can be determined by the condition of 0 = 0 at z =0. Finally, an implicit solution of 0 as a function of z can be obtained as 0 Vz = a V- (e\ 1) -dt (21) (e -l)- t+ St 0 2P I 26

I THE UNIVERSITY OF MICHIGAN 4134-2-F In the usual cases the integral can be evaluated approximately by expanding the denominator into a power series. The result is z Vne 1 log e n e2(1+ y) KT 1 + + ieo2 + 1 | e0- 2 6KT 24 (KT \ 2 + 6 (1 Gi ~e( 7 01 2 __T 2 ) (22) For the case of K < 1, a zeroth order approximate expression of (22) becomes KT n e2( +y) O z e KT 0 0 e 0 (23) Potential of the Plate The potential distribution in the plasma sheath is found as a function of some parameters and the potential of the plate, 0. 0 can be determined from the condition that equal quantities of the electrons and the positive ions hit a unit area of the plate per unit time at equilibrium. Owing to the fact that the rms velocity of the electrons is much higher than that of the positive ions, more electrons than ions may hit the plate per unit time except the plate is charged negatively so that only very energetic electrons can reach the plate. Assume the potential of the plate as 0, which is negative, and critical velocities for the electrons and the positive ions are defined in the same manner I - - 27

THE UNIVERSITY OF MICHIGAN 4134-2-F as in Chapter II. C eo 2 leoI m e 2le0o C io m. 1 The condition at equilibrium yields an equation as follows: -C r e -OD do OD c dc \ z Z \ -o 00 x -o y o Z =00 10 - cdc + 0 OD d'z \ x \ y 1fi z = 0o -- OD (24) Under the conditions 2KT ~, m e 2>> (2KT) k, mi (25) Equation (24) becomes after the integration as follows: 1 2 V"V 2KT m e eoI KT e 1 i - 2 2 - erf Kin. 0K 2KT KT (26) The potential of the plate, 0, can be determined numerically from (26). How W 28

THE UNIVERSITY OF MICHIGAN 4134-2-F ever, when the conditions in Equation (25) are valid, (26) may be approximated as follows eol 1KT KT0- log V m (27) KT T7- V e Density Distributions of the Electrons and the Ions in the Plasma Sheath Up to this point there is enough information for the determination of the density disbributions of the electrons and the ions in the plasma sheath. The procedures are as follows: (1) For a specific plasma with m, mi, no, V and T given, the potential of the plate is determined from Equation (26). (2) The potential distribution in the plasma sheath can be calculated from Equation (22). (3) After the potential in the plasma sheath is found the densities of the electrons and the ions are obtained from Equations (11) and (15). As an example, a plasma with the following parameters is considered. m. = (16 x 1825) m (oxygen ion) T = 10000K V = 7 Km/sec n = 106 1/cc O In this case KT n - n e n.= n. e o 1 o 29

- THE UNIVERSITY OF MICHIGAN 4134-2-F eo T - - 2.61 KT A distance, H, at which n = 0. 95 n can be calculated as follows: e o If n = 0.95n, e o = - 0.0513. KT So from Equation (22), H = 0. 925 cm. W 30 - --

THE UNIVERSITY OF MICHIGAN 4134-2-F IV THE INTERACTION OF A HIGH INTENSITY ELECTROMAGNETIC FIELD WITH A LOW DENSITY PLASMA The purpose of this study is to explore the basic properties of a plasma when it interacts with a high intensity EM wave. There are many publications on this subject but they are concerned mainly with the small signal case. The conventional result valid for the small signal case is no longer accurate for the large signal case. In this study it is found that some basic parameters of the plasma vary as the functions of the field intensity. The approach used in this chapter is to find the velocity distribution function by solving some basic equations exactly. After the velocity distribution functions of a plasma are obtained many basic properties of the plasma are readily found. The conventional approach adopted in many papers 5 through 8 is to assume the velocity distribution function as the sum of an isotropic part and a nonisotropic part varying with the frequency of the incident wave. In the case of a high intensity incident wave this assumption is not valid and the velocity distribution function should be found directly from a Boltzmann equation without any approximation made before solving the equation. The zeroth-order velocity distribution function obtained exactly from a simplified Boltzmann equation produces some conventional properties of the plasma. The first-order velocity distribution function obtained from a Vlasov's I. 31

THE UNIVERSITY OF MICHIGAN 4134-2-F equation can produce some significant results. Some parameters of a plasma are evaluated as functions of the field intensity. Formulation of the Problem The interaction of an EM wave with a plasma can be described with two sets of equations, namely the Boltzmann equation and the Maxwell equations. These two sets of equations are mutually coupled and results in the non-linear character of the partial differential equation which is to be solved. Assuming an incident EM wave of high intensity interacts with an infinite plasma, the Boltzmann equation is af e (% e ) af t+ v Vf + - (E + v x B). V f + ( + v x b) Vf- ( ) at m v m v at coll (1) where f is the velocity distribution function, v is the velocity of the charged particle of the plasma, E and B are the electric and magnetic fields of the incident EM wave, and e and b are the internal fields induced in the plasma. E and B can be represented as follows: E = Ecoswt x - E (2) B = Bcoswt y = - coswt 9 0 g and b are the solutions of the two wave equations which are derived from the Maxwell equations. -- 32

THE UNIVERSITY OF MICHIGAN 4134-2-F 3-a 1 3 i~ 77 e AfAd3v c-r E 2 rU j + Vp = p LA f d v+ A V o at2 o At e o at A A A E eA v at o A o A (3) 9 12h i _ C ' f A 'R V - bo b - MO V XJ = -1 Vx, A v^dv (4) at A To make the problem more specific the following assumptions are made. (1) The plasma consists of electrons, positive ions (singly charged) and neutral particles. (2) For the low density plasma the collision term is ignored. This is the case of the ionosphere. The collision term is important in the case of high density plasma which will be investigated in later chapters. (3) The incident EM wave is assumed to be of high intensity but still low enough that a non-relativistic analysis is valid. (4) The frequency of the EM wave is higher than the electron plasma frequency of the plasma. (5) The plasma is of infinite extent and homogenous in its unperturbed state. (6) In the analysis the spatial variation is neglected. This implies that the velocity function is to be determined as a function of the velocity and the time only. (7) Although the intensity of the EM wave is probably limited by the v 33 - -

I - THE UNIVERSITY OF MICHIGAN 4134-2-F breakdown condition of the plasma, this phenomenon is not considered in the present study. Based on these assumptions, two groups of basic equations for the electrons and the positive ions of the plasma can be formulated as follows: (A) For the electrons: af _- - e(E+v x B) V f - -(+v x b) V f = 0 (5) at m v e m v e e e 2-. ~V a a ' Ye f v +((f (f- f )d3v 0 0,2 t o t e ~ \ E e at o (6) V2- = 2 eV x - V(f-f)d v (7) 'o o.,2 -o e at (B) For the positive ions: af. - + e (E+vxB) 'V fb) f. + (8) at m. vi m. v 1 1 with Equations (6) and (7). These equations are solved in the following sections. After f and f. are evaluated many basic properties of the plasma, e. g., the e 1 conductivity, the permittivity, and the energy density etc., can be found easily. Zeroth-Order Velocity Distribution Function In the zeroth-order analysis, the incident magnetic field and the internal fields are neglected. The assumption is made on the ground that the effect I -- --- ~ ~ ~ ~ ~~~~ ---~ 34

THE UNIVERSITY OF MICHIGAN 4134-2-F due to the incident electric field is much greater than that caused by the other fields. This assumption leads to a simplified Boltzmann equation and a zerothorder velocity distribution function. A refined first-order solution is attempted in a later section. The simplified Boltzmann equation for the electron is (0) (0) Of af e e e at m av e - -E coswt aV = ~ <~> e x and for the positive ions, (0) (0) af. af(0 i e 1 + E coswt 0 at m. av i x Equation (9) is solved exactly as follows: (0) (0) af() af(o) e e - ycoswt = 0 (9a) at av eE x With 7 = The value of v is between -oo and co and a steady state solurn x e e (0) tion is sought. One more condition is provided by assuming that f) is a e Maxwell-Boltzmann distribution in the absence of the external field. Define (0) ikv F(, t) = f) e dv e x and take the Fourier Transform of (9a). This yields the following equation. 35

i THE UNIVERSITY OF MICHIGAN 4134-2-F 3F - + iyk coswt F = 0 at Equation (10) can be solved and the solution for F is -i-yk- sinwt F(k,t) = F (k) e O (10) (11) F (k) is a constant with respect to t and is to be determined from a boundary condition. Since in the absence of the external field (y is zero) f( is MB e distribution, F (k) can be found as follows: O 5 / m / 2KT v e e F (k) = no e o o r KTe ikv x e dv x m e = n 1 I e e m (v + v ) e y z 2KT e KT e 2 - k 2m e e (12) f )is obtained by inverting F(k, t) as follows: e Go\ 1n s -ikv f() 1 F(k t) e dk e 29 -M -n~9 9. m (v + v") KT e v z e k2 y k oo k - -iy sinot - ikv 27T 2KTKTK 2m ( e n1 n(me e no2 0 2; KT e e e 2 KT e(v + sinwt)2 v+ v2 v (^ VT ^^ l^ uy z o \2rKT e This is the exact solution of (9a). This is the exact solution of (9a). K dk (13) - - - --- -- 36

I I THE UNIVERSITY OF MICHIGAN 4134-2-F Similarly, for the positive ions the zeroth-order velocity distribution can be found as m. where y' stands for eE 12:1 mi These velocity distribution functions are obtained from a simplified Boltzmann equation without any approximation made in solving the equation. These functions can be used to calculate some approximate properties of a plasma. In effect, these functions produce the conventional results which have been derived by other methods. Applications of Zeroth-Order Velocity Distribution Function In this section the zeroth-order values of the current, the conductivity, the permittivity and the energy density of the plasma are obtained by using the zeroth-order velocity distribution functions found in the preceeding section. Zeroth-Order Current: By definition, the zeroth-order current is J = -e vf dv + e vf. d v 5-~0~vf(0'd3veS f!0d3v e A 1 After carrying out the integration, it is found that I 37

THE UNIVERSITY OF MICHIGAN 4134-2-F - 2 2 e n e n (~)=J x= i-~ + (~ Esinwt (15) x Lom urn. This agrees with the known result. Zeroth-Order Conductivity: The conductivity is defined as the ratio between the current and the applied electric field. That is?r 2 2 n en en =-a L-i -- + -o (16) Zeroth-Order Permittivity: If the permittivity is defined as = (1 + ), the zeroth-order 0 jWo 0 permittivity of the plasma in the absence of a constant magnetic field is found to be - 2 2 W W e 4eA 1 -2 (17) 2 2 en en 0 0 Where -- w --. This result is a well-known one. pe m e pi m.e e o 1 o Zeroth-Order Energy Density: With the energy density of the plasma defined as 1 2 (0) dv3 I 2 f(0) 3 u 2 [ i fA 0d v + 2 V i-my f dvv (18) u= m v2f()dv I 38

THE UNIVERSITY OF MICHIGAN 4134-2-F The zeroth-order result is 3KT e 22 __~e me / e2E2>. 2 l u= no t + 2 2')sln wt -2 2 2 2/ J 3KT. mi e2E2 2 + n[ + sin * (19) 1 o2 \m2 / 2..J The time average value of u is 3KT m / 2 u = % ~ + n + 22 ) e 3KT. m2 / 22 2 + n +?) (20) o -.2 4 2m 2 First-Order Velocity Distribution Function The zeroth-order velocity distribution functions of the plasma have been found by neglecting the incident magnetic field and the internal fields. These approximations are justified in the small signal case. However, for a high intensity incident EM wave these approximations are not justified, because the velocity of the charged particles of the plasma induced by a strong EM wave could be very high. In order to find some basic properties of the plasma with a higher accuracy and valid for the strong signal case, a first order velocity distribution function is attempted in this section. The basic equations to be considered are the Vlasov's equation and the Maxwell equations. I - 39

THE UNIVERSITY OF MICHIGAN 4134-2-F af(1) e - + -- +vx B) V f(l)+ e ( +v x b) - f = O at m v m v (21) Equation (21) is the Vlasov's equation without the space variation term. The internal fields are found as follows: By neglecting all the space variation terms, the internal fields can be written down from Equations (6) and (7) as b= 0 (22) and a e2 —\v(f. - f ) d v at2 E 1 e ~~a~~~ ~it With the assumption that t varies as e, e can be found approximately as - * -e a.,(0) d3 vf d v 2 at e w o 2 W - Ecos wtx (23) 2 W The substition of Equations (22), (23), and (2) in (21) give two basic equations for the electrons and the positive ions as follows: af Qaf ( v1) N af() e e z e x e -a - y coswt av + ycoswt - -cost = 0 (24) a~t ~ av ~xc avx 0 Z X 0 X 0 Z 40

THE UNIVERSITY OF MICHIGAN 4134-2-F and (1) (1) (1) (1) af( af( v af() v af() 1 1 1__ at + -y'coswt - ost + y'coswt = (25) at 8v c 8v c av X 0 X 0 Z 2 L where 3 = (1 + ), c = velocity of light, 7 and y' are defined as before. 2 o Li Equations (24) and (25) can be solved exactly and the method used in solving (24) is briefly outlined as follows: Write (24) as (af) af af af() 1 (l +A e e(ye) e s + a A -e+B v a -v a- = 0 (26) x x Z where A= -^y, B= - c 0 Define a new variable 7'. '7 = A sinwt and a new constant 1 C = B/A — Equation (26) then becomes Equation (26) then becomes Of(1) Of(1) Of(I) af af at " e e e d? ( z + cv -v 0. (27) d' + +vz av x av x z The initial condition for f) is assumed to be e I -- -- 41

THE UNIVERSITY OF MICHIGAN 4134-2-F (v2+ + 2 v2 ) (1) (2KT x y z (1)(vm e fe (v, O) = n / 2KT (V+ v (28) Equation (27) implies, [9), dv dv x z 1 + cv cv z x Using a well-known technique of integrating first order partial differential equations, it can be shown that (27) has a general solution as (1) e 12 fe = g(ul, u2) (29) u= c(v2 + v)+v (30) 1 2 x z z 1 cv \ u = T +-sin 1 \1+2cu (31) U2 cV + 2cu1 where g is an arbitrary function of u1 and u2. To construct an appropriate (1) solution for f, one can make use of the condition that as t = 0 (or 7 = 0), e 2 2 g(ul, u2) should reduce to (28). By finding an expression for v + v in terms 1'j.~~~ 2^ ~x z of u1 and u2 at 7 =0, fe) can be constructed. ~1 ~2 e v + v = u - 2 (I1+ 2cu1 sincu2- 1), at T = 0 x z c 2 c Using (30) and (31), it follows Using (30) and (31), it follows - -- -- 42

----- THE UNIVERSITY 4134-2-1 OF MICHIGAN -? 2 2 cu!- 2 c ( T+2 cu1 2 2 (1- cosc7)(l+cv ) sincu2-1) = +v + 2 x z 0c/2 - cv sin c x Hence f(1) _ e / m 3/2 e o2irK KT) 3Xp m e 2KT e 2 2 2, v +v +v + x y z c2/2 (32) A further simplification is made as follows: 2 2 2 2 v +v+v + (1- cosc)(l + cv )- cv sincT x y z 2 z x c v+ v+ v+ 23c v sin x y z ox ( 7 sinwt) + WC o 22 z 223 c (1 - - ) 00 /3co - cos ( -sin wt) c o _ = v + 3sc sin( t) + v2+ v - x o Wc e y z o 3Co(1 - cos ( 2L 0c o sint) 2 sin wt) ) Therefore, the final solution for f() is e (1)me 2 f = n exp + c sn sin t) e o 2r 2KT 2 - x o eC \ e/ L e o0 + v2 y + v - 3c (i-cos ( - sinwt) z o WC L-~~~~o )jfj (33) with 2 w = (1+ - ), 2 u eE 7= - m e c = velocity of light. o By using the same technique, f can be found-as 1 43 - -- --

-THE UNIVERSITY OF MICHIGAN 4314-2-F f(1) mn. fo 2= n 0 KTi i V +exp - i cV1- Cos( sin wt) + v2 2KT. x o0 (WC y + v z- c (1 - cos (C' sin wt) ) L%- O CJ (34) eE with '= - m. Equations (33) and (34) are the exact solutions of (24) and (25). These solutions are mathematically very neat and physically very plausible. Some interesting results are found with (33) and (34) and are presented in the next section. Applications of First-Order Velocity Distribution Function In this section the first order values of the current, the conductivity, the permittivity and the energy density of the plasma are obtained. First-Order Current: With the current defined by: 5(1) -e,(1) d3 f(1) 3 - =r-e vf d- v+ef Vf. d v e I1 the result is found to be — (1). J Jx z + J x z (35) 44

I I I THE UNIVERSITY OF MICHIGAN 4134-2-F 2 2 W eE eE J = en (1+ — ) c sin( sinwt) + en (1+- e) c sin( sinwt) (36) x o 2 o [m c o 2 o wm.c w e o W 10 2 2 J =en (1+ -- e) cos( sinwt)- en (1+ - )c cos(e — sinwt). (37) z 2 o o m co 2 o w2 m.c D e o co 10 eE eE The factor e is always smaller than unity, because - is the velocity wm c mrn e o e of the electron induced by an incident electric field E and this cannot exceed the velocity of light, thus (36) and (37) can be expanded into Fourier series as follows: 2 ) 3 5 eE eE 1 1 eE J en (1+ ) ) ( ) mx o 2 m 2 cm 320 4 om o e 8c e c e o o eE 1 eE 1 1 eE + -— (- — (-) + --- sint om 2 w m. 320 4 1m. i 8c i c i 0o 2 2 3 5 + en (1+ ( eE 1 1 (eE) 1 1__ eE o 2 24 2 Wm 640 4 wm o c e c e O O 3 5 1 1 eE 1 1 eE 5 + - -2 ( --- ) - - - ( - ) sin3cot 24 2 corm. 640 4 cm. C 1 C 1 0 0 W leE 51 _ + en (+ - (eE sin5 t (38) o 2 [1920 4 wm 920 4 Dm 1920 4 rm. W c e c o o 45

- THE UNIVERSITY OF MICHIGAN 4134-2-F 2 co 2 4 j -en(l+ i 1 ( I1 eE 2 1 1 eE 4 Jz = - ( + 2 ) 4 c (Wm ) - 64 3 (wm tj _ o e c e o 4 c om. 64 3 com. O 1 c i O 2 ~2o c 2 4 I _ 1 1 eE + en (1+ ) - (e ) eE o0 2 4 c uWm 48 3 Wm W L- o e c e O ~2~~~~o 1 1 ( eE ) 1 1 E 4 4c (~ ) cos2cot ~ 4 c Wm. 48 3 (cm. o 1 c 1 0 +enp(1~+-e( I1 1 eE cos41 3 o 1 o( W2 )[192 3 ( m 192 -3 (Wm) 9 2 * (3 ) 0 0 Some significant points are summarized as follows: (1) In the absence of a constant magentic field, a strong EM wave produces a drift current in the plasma in the direction of the propagation of the wave. This phenomenon can be visualized from the physical point of view. (2) The current produced in the direction of the E field of the incident EM wave has odd harmonics and that in the direction of propagation of the incident wave has even harmonics. (3) The current induced by a strong EM wave is not a linear function of the field intensity of the incident EM wave. -- -- 46

THE UNIVERSITY OF MICHIGAN 4134-2-F First-Order Conductivity: If the conductivity of a plasma is defined as the ratio between the fundamental component of the current to the incident electric field, the conductivity is found to be 2 o- =-jen (1 + -) + ) ) E xx -joen Wm wme 2 m Wm W LQ e 8c e i 0 + 1 e e4 (40) + — (..) - (e) E4 (40) ~-. 4 Dm 't Ym om. 320c L e wi 0 The conductivity in the z-direction cannot be defined since I does not have a z fundamental component. First-Order Permittivity: According to the definition of the permittivity of the plasma, e - e (1 + -xx xx o jwe the result is found to be 2 2 2 L L2 2 xx LoD){ +4 - e e [ D 2 D 2 ~e 2tow w 2 e 2 1- (1 2 2 2 m 2 m 0 8c w e U 1 r - 2 2 1 _ _ _ _ _ _ - 2 _ +30 [ p(e )4+ (i )n E4 } (41) _320c L i8 e t 0 2 This is a function of the intensity of the incident EM wave. ---- -- --- _.,, — 47

THE UNIVERSITY OF MICHIGAN 4134-2-F First-Order Energy Density: With the energy density of the plasma defined as 1 f 1 2 f(1)d3v 2 e e 2 1 i the result is 2 2 u= n [KT +m (1+ e) c (1- cos ( sinwt) 0o L2 e e 2 o om c L-U cje o c 2 2 + n [ KT+ mi(l+ e) c2(1- cos (eE sint) (42) o 2 1 1 2 o Lim C L ~ ~~o wmcU i o J The time average value is o 2 e 4 e 2 tm 64 e 2 2 wm w e w c e 2 2 3 to 2 eE 2 2 2 4 + TKT.+m.(l+-) - -m.(l+ ) 1 ( (4) 2 1 4 i 2 tm. 64 i 2 2 w m. U i u c i o Non-Linear Modeling the Velocity Distribution Function Equation (9a) can serve as the starting point of a non-linear modeling attempt in the following way. For different values of 7y, say y1 and y2, we write af af at - 71costt - = (44) x1 v x 48

THE UNIVERSITY OF MICHIGAN 4134-2-F and af2 af2 2 2 at - '2C wt av = 0 (45) x Now we ask if it is possible to find f2 as a function of f alone. Assuming that such a relationship does exist we try to determine it as follows. We assume f = f(f) hence af2 df2 af at dfl at and af2 df2 av df1 af av x and (45) becomes df1 [a - 2 av or making use of (44) df2 af df cos t a- (Y1-Y)= 0 1 x (46) (47) af1 Since a- is not identically zero and y is chosen 7Y2 it follows that x df df1 or f = constant. This result is physically impossible thus we conclude that 2 f2 = f2(fl) is invalid. We might then ask, more modestly, if we can find f2 as 49

THE UNIVERSITY OF MICHIGAN 4134-2-F a function of fl and t, however for these modeling attempts it seems more appropriate not to work directly with the velocity distribution function but an integral of this function as will become evident in the non-linear modeling discussion of the next example, the high intensity field acting on a weakly ionized gas. 50

---- THE UNIVERSITY OF MICHIGAN 4134-2-F V THE INTERACTION OF A HIGH INTENSITY ELECTROMAGNETIC FIELD WITH A WEAKLY IONIZED GAS The previous chapter treated the interaction of a high intensity EM field with a very low density plasma where the collision effect can be neglected. In this chapter the interaction of a high intensity EM field with a weakly ionized gas, or a plasma of low ionization, is considered. The collision effect is important in this case. A partially ionized gas is assumed to be composed of electrons, the positive ions (singly charged) and neutral particles. The possible collisions are e-n, i-n, e-i, e-e, i-i, and n-n. If the degree of ionization of the gas is low, only the e-n and the i-n collisions are important in the analysis. Fortunately, the mathematical models for these two types of collision are simple, therefore, a rather rigorous analysis is possible. The Boltzmann equation is solved exactly with a simple mathematical model for the collision. Some reasonable approximations are made to derive some useful parameters only after the exact solution of the Boltzmann equation is obtained. The zeroth-order and the first-order solutions for the velocity distribution function are found. The basic parameters of a plasma are evaluated as the functions of the intensity of the incident EM field. 51

THE UNIVERSITY OF MICHIGAN 4134-2-F Formulation of the Problem The Boltzmann equation is [f e - e - Of + - f+ (+ f+ + f = () (1) at m v m v at coll. where E and B are external fields and can be expressed as E = E coswt x A. E A (2) B = Bcos wt y= -cos t C 0 The internal field g and b are found in the previous chapter as 2 w - E cos wt x < a w a(3) T bs- 0 The spatial variation is neglected again in this part of the study. The collision term is approximated as follows: For the electron case, among all types of collision the e-n collision is predominant when the degree of ionization is low. To take into account this collision the simplest mathematical model is to assume af (- )cl -v(f-f ) (4) at coll. o v is the collision frequency of the electrons and the neutral particles. v is known to be nearly constant for the low velocity electrons and nearly proportional to the velocity for the high velocity electrons. A more accurate behavior of v was i -- -- 52

THE UNIVERSITY OF MICHIGAN 4134-2-F found in experiments but is too complicated for the theoretical analysis. In order to improve the accuracy of the collision model expressed in (4), f must be modified. f is the equilibrium state of f and can be assumed to be a Maxwellian distribution in the small signal case. However, if the incident EM field is strong, f may be quite different from Maxwellian. It has been found that the isotropic part of f which is constant with respect to time varies as a function of the intensity of the incident EM field (5]. After the substitution of equation (4) in (1), the comparison of both sides of (1) suggests that f in (4) should be the isotropic part of f because the left hand side of (1) is dependent on time when no constant external field is present. By modifying f and specifying v carefully, the collision model expressed in (4) is acceptable for the present analysis. The zeroth-order analysis is made by neglecting the incident magnetic field and the internal fields. These neglected fields are taken into account in the first-order analysis. Zeroth-Order Velocity Distribution Function For the zeroth-order analysis, the effects due to the incident magnetic field and the internal fields are neglected. The simplified Boltzmann equations to be solved are as follows: For the electrons: af(1) (o) e - y cos t = (f0)- f (5) at av e eo x -.......1 -53

THE UNIVERSITY OF MICHIGAN 4134-2-F where 7 --, v is the collision frequency of the electrons and the neutral e (0) particles, and f is the equilibrium state of f. f will be specified careeo e eo fully later. For the positive ions: (0) (0) 1 + y'cos t a1 = (f()fio (6) at av i 10 x eE where 7' --, v' is the collision frequency of the positive ions and the m. IL~~~~~~~~ ~~(0) neutral particles, and f. is the equilibrium state of f.. f. can be assigned 10 1 10 to be Maxwellian due to the large mass of the positive ions. The analysis is carried out for the electron case only. Equation (5) implies the following relations: dv df) x e dt = x (7) -ycos wt _v(0)_ f e eo From (7) two equations are obtained as follows: v + sinwt = C = constant (8) x U df(0) + vf 0 f (9) dt e eo Equation (9) leads to 5dt Svdt t dt e f f) vf e dt. e L ej t-t- J eo -00. - -- 54

THE UNIVERSITY OF MICHIGAN 4134-2-F If the collision frequency v is assumed to be independent of time, t, f can be e written as vt (0) Z/vt e f = f e dt (10) e eo -oo In evaluating the integral, f is considered to be a function of t. Originally eo f is a function of v, v, and v. However, with (8) f in this integral eo x y z eo should be rewritten as f (v, v, ) f (C- 1 sinwt, v, v) eo x y z eo w y z or (10) is rewritten as t ef = vf (C- Y sinws,v,v ) eV ds. (la) e eo w y z In order to put the integral on the right of (1Oa) more explicitly, a new variable u = s-t is introduced. This changes (lOa) to 0 e f -) vf (C- Y sinw(u+t), v, v) e du (11) e eo w y z The constant, C, is eliminated by substituting (8) in (11). Thus f() \ f v +2-sinwt - sinw(u+t), v v e du e eo w ' y z Changing the variable u into -u, the final solution for f) is obtained as e. --- 55

THE UNIVERSITY OF MICHIGAN 4134-2-F 00 f() = vf rv + - (sinwucoswt + (1 - coswu) sinwt), v,v e du e eo LX w y z L o (12) Equation (12) is the exact solution of (5) subject to the assumption that v is independent of time. The argument of f appearing in (12) means that the v eo x term in the original form of f is to be replaced by eo v + [sinwucoswt+ (1 - coswu) sinwt] x oe Equation (12) is a definite integral and can be evaluated immediately once feo and v are specified. In order to have an accurate form of feo which is also appropriate for eo the formulation (5), the following cases are considered. As the first choice, the result obtained by Margenau [53 will be used. That is 2 m /2dc2 f = A exp - \ --- - (13) eo ^ KT+My2X2 /6(c+)) 2 2 where M is the mass of neutral particles, X is the mean free path of the electron and the neutral particles. A is a constant defined in such a way that f d3v = n i eo o If (13) is used, f) can be expressed explicitly as e I 56

I THE UNIVERSITY OF MICHIGAN - 4134-2-F 2 22 (v + H(u,t)) +v +v x I y z OD,.~~~~~~~ 0?~~mm /2dc2 (0) _ u e f = ve-ve du A exp - (14) e J _ J KT+M" X /6(c + X ) <e =5e dUexKT + My2X2 / 6(c2+w2x2) j 0 ~~~0 where H(u,t) -= [sinwucoswt - (1- coswu) sinwt]. In the evaluation of (14), v can be assumed to be a constant or v/X depending on the circumstance. As some special cases the following are considered: M(y 2 2 (1) For the case of KT > v2 W >, or when the frequency of the EM field w is higher than the collision frequency and the kinetic energy added by the incident EM field is smaller than the thermal energy of the electrons, an approximate (0) expression for f) can be found as e (0) -vu e 2 2 f() = 5ve du A exp- 2KT [(v +H(u,t)) ++v ]} (15) withQ x ye J with with My2 T =T+ e (16) 6K w2+ 3 (KT6 + L X2m 6w e 2 (2) For the case of - >> KT, or the very strong EM field case, f() can be 2 e w expressed as = vedu A exp 2 v+H(,t)) +v f(0) v e du A exp - eH 2 e 2M22 x y z + 2w2 2 L(v +H(ut)) 2+ +v2 L x y zY (17) 57

THE UNIVERSITY OF MICHIGAN 4134-2-F corresponding to the limiting case of f e eo r 3m (v4 + 2 2v2) ] f =A exp M- v eo 2M-2 __ 2 2 (3) For the case of v > w, or as the frequency of the EM field is lower than the collision frequency, the result obtained by Chapman and Cowling [10] seems more appropriate. In this case v 3m cdc fe — exP(18) (0 y \du A exp / 2 22 (19) -v0e 0 and 2 2 2 1/2 O (4) For the case of v > w and My - ~ KT, the limiting case of (18) or the 2 c result obtained by Druyvesteyn is used. That is 3m v4 feo = A exp M- — 2 j (20) - 4My 2 20 and (0) -vu e 3m 2 f v e du A exp v (21) e duAexp (v x+H(ut))+v y+v (21) ~0 I 58

THE UNIVERSITY OF MICHIGAN 4134-2-F As a matter of interest and with the purpose of checking the theory, a well known result of f for the small signal case is reproduced from (12). e When the incident EM field is small, it is accurate to assign f as eo 2 m V e / me 3/2 _ e -f n Te 2KT eo o27r KT and it is also reasonable to put f [v + V' (sinwucoswt + ( - coswu)sinwt), v vI eo x w ' x z 2 m v tin me — e e - e n o e e K1 - Y v (sinwucost + ( - cosu)sint The substitution of this expression in (12) gives,_. r m Inn -(0) fp e- - Z Cos t- e sin e o 1 KT x 2 2 KT x 2 2 UO +v tU +1 yv x af af eo eo f x 2+7 - ( -v cos t+^o 2 2 2 v sinwt (22) eo x 2 2 2 X 3v x 2 222 3v / v +w X v +Co X This result was derived by Margenau [5] Application of Zeroth-Order Velocity Distribution Function In this section the zeroth-order values of the current, the conductivity and the perittivity are obtained by usin the reeding zeroth-order velocity distribution function found in the preceeding section. 59

THE UNIVERSITY OF MICHIGAN - 4134-2-F 1 Zeroth-Order Current: With the zeroth-order current J(0)). J~~~\ defined as " (0) 3 vf dv e it can be found that ^ -110) A J =Ix x and J = -e v dv dv dv A expx \ \ l x x y z Vxdxd d~ Ap-W x~~~0 (v +H(u,t)) +v2+v2 x y z 2 m /2dc e KT+My X /6(c 0 du 2+w2x2) 00 so -Zt ve (23) To evaluate this complex integral approximately, a new variable v' = v + H(u, t) x x is introduced and v is held constant at this step. Equation (23) then becomes 2 e 5 5dv2 22/c2 22) ]H(utve J = 47r e \du \dv v A exp - 2 2 22 H(u O e )2XL _ KT+M7 2X/6(c2+ X2 ) '0 ^0 t) Think of v as v/X, and find a value of v which corresponds to a point near the peak of the integrand of the above integral. This value of v is a root of the following equation. - -- 60

--- THE UNIVERSITY OF MICHIGAN 4134-2-F 2 a 3 m /2de2 exp - - 0 v L J KT+M X /6(c +w X 0 The root of this equation is found to be I 1 \ 22 2 2 2 2m 2 /2KT *V- — 2 2e V 8m My X + (2m w AX+ 6 KT) - 2m WX +6KT 2 Vmme e e ~~~~~~~~~~~~~e ~~~(24) J is approximated as x 2 m 2/ C 2 2dc2 J 47re v KdvAT+2 e H(u,t)- e du x \ \dKT+ My2h2 /6(c2 +w2 2 o 0 = en -iA2 coswt + en y sin t (25) i + (2/X) wi +t (ic/X) It is noted that (25) is the approximate value of the current. If the harmonic components of the current are needed, (23) must be evaluated more accurately. Zeroth-Order Conductivity: The conductivity is defined as the ratio between the current and the applied electric field. It is found to be o = - - io. (26) r 1 61

THE UNIVERSITY OF MICHIGAN 4134-2-F n e 2m X 2m X 2m X 2(2 n e 2 m e X 2 me x 2 m e j (27) r r2 2+ 2 2 e1e (+ ) + 3 K 3KT 2m XZ2 2m 2 2m X 2 n e (28) Oj= m 2 2 2 < e wT 2 3KT +w 3KT 2m X 2 2 2m X 2 2m X2 e e e Zeroth-Order Permittivity: If the equivalent permittivity is defined as e = e (1- - ), the zeroth0 WE 0 order permittivity of the plasma in the absence of a constant magnetic field is found to be 2 E = E ____pe (29) E m 1 2 0 | /2 3KT\2 2 3KT - 2 2 -- 2 2 2 2m X2m X 2 2m X e e e 2 2 n e where w = - pe me e o It is interesting to show that a well known result for the conductivity 2 in the small signal case can be derived from (23). When KT >>, (23) 2 can be written as 62

T THE UNIVERSITY OF MICHIGAN 4134-2-F 2 ODI m V e J = -e v dv dv dv A e x \\ x x y z x^ [ - m v H(ut) KT x \ ve du -0 41re me \ 4 - 3 KT \ dv 3 K 2 mv 2KT A e O H(u, t) v e du Letting u = v, J becomes u 2KT o x 8 Jx- 3 V 0 ), noe 2 4 -U /YU -u i/u e 2 2 +/, 8 + 3 00 r 2 \ -u 7ynoe \ e 4 du cos wt 2 wu 2 2 du sinwt 2 +2 t0+ v (30) and the conductivity can be written as 2 n e O Y = o r m e 2 n e o i m e 8 3 V 2 2 -u vu 2 2 W +V du du 8 3 V 00 2 2 -U WU e 2 2 w +v 0 O These are the results obtained by Margenau C83. 63

THE UNIVERSITY OF MICHIGAN 4134-2-F First-Order Velocity Distribution Function The Boltzmann equations to be solved for the first-order analysis are as follows: (A) For the electrons, af(l) af(1) (1) V f(1) o~f 3~f ~ v af v f e e z e x e -at - y37cost + cost v - Coswt at av c av c 3v x 0 x o x - v(f-f ) (31) 2 e. eo w where 3 = (1+ ~2), c = velocity of light, y and v are defined as before. O (B) For the positive ions, (1) I(1) (1) (1) af af) v af) v af) 1 1 Z 1 X 1 + 3y' cos wt - -' cos t -z - - + 7'coswt - at av c av c av x 0 x 0 Z - -v(f( - f ) (32) 1 10 where y' and v' are defined as before. Equation (31) implies dx(1) dx - dz e dt =(33) v - v z )x - v( - f -( - - ) Co cos wt e cosLf eo c c 0 0 Three equations are obtained from (33) as follows: v v + (Q -z ) Y sinwt = c (34) x c w 1 0 I 64 -- --

THE UNIVERSITY OF MICHIGAN 4134-2-F V v + Y sint = c (35) z c W 2 o w O (1) e (1) -e + vf vf (36) dt e eo or C2 c- (- ( ) sinwt 0 v 2 sin (37) 1+ 1__ sin2wt 2 2 0 C t Fcw1 _ L Ys intt) Y sinwt 2 c c w w 0 0 v (38) +1 y _ sin wt 2 2 C W 0 Following the same procedure as in the zeroth-order analysis, an expression for f(1) can be found to be e G0 f(1) r v e-vu f = ve duf (vi, V, v') (39) e eo x y z f (v', v', v') means that v, v, v terms in the original form of f are to be eo x y z x y z eo replaced by v', v', V v' v', v' are expressed in terms of v,v, v, t, u, xd s z x y zr r a z and some other parameters as follows: W 65

THE UNIVERSITY OF MICHIGAN 4134-2-F vt = 1 (sinwu sin2wt-cosuu(1-cos2wt)) V 1+ i sin (t- u) o 2 2 C kC v W zV ~ 2 C- - - in u cos + (1 - cos nu) sin4i c 22 0 2 o v' = v (41) Y Y -z -z 1+ 22 2 { ] c 00 1- -L _ - siwu sin2wt c 22 o D -: J3 1 (1- coswu)(l- cos2ut) (42) o CO After the substitution of (40) - (42) in (39), f) can be theoretically e evaluated as a function of velocity, time, and the intensity of the EM field. The actual evaluation of (39) may be quite impossible without some approximations. The applications of the first-order velocity distribution function are straightforward but are omitted here to avoid the lengthy mathematical formulas. I. 66

--- THE UNIVERSITY OF MICHIGAN I 4134-2-F Non-Linear Modeling the Induced Current We may use this example in a non-linear modeling attempt in the following way. In one system we measure the effect of applying a weak EM field to a plasma, while in another system we would like information when a strong field is applied. In this example, our purpose is to predict the result of the strong field experiment from the data obtained in the weak field experiment. This problem can be solved if we can model the basic differential equations which govern the two systems. That is to say we aim to obtain a quantity in system II as a function of the quantity in system I theoretically. The quantities of interest are the induced currents in two plasma systems to which two different electric fields are applied. We assume that the two systems have identical plasma of weakly ionized gas type. If an electric field of E cos wt is applied to the first system, there will be an induced current i1. The question is what will be the current i2 in the second system if a strong field, JE coswt, is applied to it. e is a constant much bigger than unity. If we can succeed in solving for i2 as a function of i theoretically, we can predict i from i which can be obtained from a much easier 2 1 experiment. For the two systems of plasma assumed above, we can formulate two Boltzmann equations to describe them as follows: 67

THE UNIVERSITY OF MICHIGAN - 4134-2-F 3fl eE af1 nl(t) - - cost av = -v (f - f (43) m av 11 n 10 x o f2 feE n 2 (t) - t - — lcos wt -- (f(44) at -tm av 22 n 20 X 0 f and f represent the velocity distribution functions of the electrons in system 1 2 I and II. To system I, an electric field, E cos wt, is applied in the x-direction and to system II, a strong electric field of I E cos wt is applied in the same direction. I is a constant much larger than unity. vl and v2 are the collision frequency of the electrons with the neutral particles in the two systems. v1 and v2 are the functions of the gas temperature and the intensity of the applied field. The collision model adopted in this analysis is the elastic collision type with conservation of particle during the collision. The conservation of particle is a necessary condition for the formulation of (43) and (44) in which the spatial variation terms are dropped. nl(t) and n2(t) are assumed to be the fluctuating densities and n is the unperturbed density for both systems. It can be proved later that nl(t) and n2(t) are independent of t from (43) and (44). fl and f2 are functions of velocity and time and the attempt to solve f2 as a function of f is difficult. However, there is a trick to sidestep this difficulty by eliminating one of the independent variables before solving (43) and (44) in the velocity space. 68 -- -- -- - II II I ~ I I I' I

THE UNIVERSITY OF MICHIGAN 4134-2-F Thus, from (43) at fl d atSldv eE af cos t \ m I av x 3v - d v = 5(fl n1(t) n o f10) dv or, since af3 d3v = 0 -avx I d3 f dv =- n, and f d3v = n 1 1 a at nl (45 The same operation on equation (44) gives a at 2 =0 (46 If the unperturbed densities in the two systems are the same and equal to no then we can put n= n2 n. The operation, Sv d3v, on equation (43) yields >Jx -) I) a 5 3 at f1vdvx eE 1 - cos wt m \ v x v d3v -V (f x 1 1 n1 (t) n o 3 f 0)v dv 10 x or a eE (nu) + os wt n = -Vn u at 11 m 1 111 (47) where nl = flv d v 1 1 Jx - - 69

I THE UNIVERSITY OF MICHIGAN - 4134-2-F and f10 is the equilibrium state of fl and is assumed to be isotropic in the velocity space. The same operation on (44) gives - (n u )+ eE - cswtn -v n u4 at 2 2 m co tn2 = -nu2 (48) The combination of equations (43)-(48) gives two differential equations as follows: d eE dt 1 + vu =- - cos wt (49) dtl 11 m +vu - cos Wt (50) dt 2 22 m Now, we have transformed (43) and (44) to (49) and (50). The quantities u2 and u2 are to be found instead of fl and f2. The problem is simplified because u and u2 are functions of t only. The next step is to try to solve u2 as a function of ul. An equation relating u1 and u2 is immediately found as follows: dt u vu Ul + vU) (51) dt 2 22 dt 1 1 An expression of u2 as a function of u1 can be produced from (51) in the following way: -v2t d V2t u= e u + u) e dt =1 dt2 = -2(ue V2t uV2t = e (ule )dt+ (v1-2 ule dt I. 70 --

THE UNIVERSITY OF MICHIGAN 4134-2-F That is -v2t t2t u2 ieu + (vl-V2)e 2 ue 2 dt (52) Equation (52) shows that the mean velocity of the electrons does vary non-linearly with the applied electric field. There are the following ways to express u2 as a function of ul more explicitly. (1) If we can determine u from an experiment, u2 will be known immediately after u1 is substituted in the integral appearing in (52). (2) Solve (52) by successive approximations. That is, consider the second term on the right of (52) as a correction term and substitute u with 1 j U2 in the first approximation. (3) If we can determine theoretically an approximate solution of ul (this is possible in many cases) then u2 can be obtained after the approximate solution of ul is substituted in the integral in (52). In this case u2 will be an approximate solution. It is perhaps interesting to show the result for this particular example. From (50), a solution (exact in this particular example) of u1 is found to be eE 1 u1= m 2 2 (v coswt + sinwt) (53) 1 m 2 2 1 v +w The substitution of (53) in (52) gives - - 71

I THE UNIVERSITY OF MICHIGAN 4134-2-F U2 1= + L2 1 m (v- 1) (21+ )(2 2 2) (v +w )(v + w) 1 2 (v v2- w )coswt+ (v +v2 )wsinwt (54) If the electric currents are defined as i1 = -enlu, i = -en u 2 22 and n = n = n 1 2 o, then the relation between two currents can be expressed as n e2 i2= Lil- oL =2 1 m (v2- V) 2 (2 2 12 2 L v1v2- w2)cos wt+ ( + v2) w sinwt (v +w )(v +t )L '. IL(55) Thus i2 is determined as a function of i and t. This time dependence will not appear in the actual experiment. The quantity we can measure in an experiment is the magnitude of the current. For this quantity the relation is n e2E 9\ + i- Jn n e2E + ( - i-~ ~m 1 ms 2212 (v2- v1)(v1Z2- + ) (V1 + 2 )(v2+ o ) 1 2 2 1/2 (v Vl )(Vl+ v2)wU (v1+ 2)(v2+ W ) 1 2 (56) This relation can be checked experimentally. 72

THE UNIVERSITY OF 4134-2-F MICHIGAN - Note that 1 l3KTle v - V 1 x\ m 1 3KT2e 2 X m 2 My T =T+ le T e -f 2e 6K 2 + -2 (KT+ M2) X m 62 t2 Mv2 6K 2 + 3X2 m L2 2 (KT + le 6w J where X = mean free path of the electrons and the neutral particles, T = temperature of gas, eE 7 v m- - M mass of neutral particles. M -= mass of neutral particles. 73

THE UNIVERSITY OF MICHIGAN 4134-2-F VI THE ELECTRICAL CONDUCTIVITY OF A PARTIALLY IONIZED GAS In this chapter an expression for the electrical conductivity of a partially ionized gas is derived where both electron-neutral particle and Coulomb type collisions between particles play important roles. As mentioned previously, considerable work has been done along this line. Margenau C5], 0C6 and his group published a series of papers dealing mostly with low intensity fields and weakly ionized gas. Spitzer iD and his co-workers on the other hand dealt with the small signal static conductivity of a fully ionized gas. In a weakly ionized gas Coulomb collisions between charged particles are neglected. Boltzmann's equation is solved in these cases by considering only collisions between the electrons and neutral particles. In the fully ionized case, however, it is the Coulomb collision which determines the electron velocity distribution function. This chapter deals with an intermediate case of the interaction between a low-intensity electromagnetic field and a partially ionized gas where neither the electron-neutral particle nor the Coulomb collisions can be neglected. The term conductivity, as used here, is defined in the usual manner as the ratio between the current produced in the ionized gas to the amplitude of the incident electric field. The effects of inelastic collisions and a steady magnetic field are neglected in the present analysis. 74

THE UNIVERSITY OF MICHIGAN -- 4134-2-F The contribution to the collision term in Boltzmann's equation due to collisions between neutral particles and electrons is accounted for by the standard analysis ClO, 2i. The contribution due to Coulomb collisions, however, is rather complicated to analyze. Here, we shall use the Coulomb collision model derived by Dreicer A3i from the Fokker-Planck equation. After assuming that the interacting electromagnetic field is of low intensity, an expression for the electron velocity distribution function is derived by solving Boltzmann's equation with the above two models for the collision terms. Expressions for the electrical conductivity in various cases are then derived. To the extent that the assumptions for the collision models are valid, the expressions derived here for the conductivity are quite general in nature. The expressions reduce to the wellknown relations for both the limiting cases of fully and weakly ionized gases. Basic Formulation of the Problem In this section we shall formulate the problem in general terms. The basic parameter which has to be determined before one can obtain any information about the ionized gas is the electron velocity distribution function. Let us define the electron velocity distribution function F(v, t) such that Fd v gives the number of electrons whose velocities lie in the element of volume d v situated around the point v in the velocity space. It is assumed that the macroscopic properties of the gas do not vary from point to point. Then the distribution function F satisfies 75

- THE UNIVERSITY OF MICHIGAN 4134-2-F Boltzmann's equation: - + a(t) F (7 (1) at v t where a (t) is the force per unit mass on the electrons, V is the gradient operator in velocity space, and (t ) is the rate of change of the distribution function due to various ' coll. types of collisions. Here we assume, t oll (atw c at (2) 1coll. -cn /cc where the subscript en means collision between electron and neutral particles, and cc means Coulomb collisions between the charged particles. Explicit expressions for the collision terms depend on the type of model one assumed for the gas. In general there should be considered two other equations similar to (1) in order to account for the velocity distributions of the heavy neutral particles and positive ions. However, for simplicity of analysis it is assumed here that the heavy particles are stationary relative to the motions of the electrons. This is a reasonable assumption compatible with the physical cases where the mass of a heavy particle is about 1830 times heavier than that of an electron. In the present problem the acceleration of the electrons is assumed to be of the form -- -- - -~ 76

THE UNIVERSITY OF MICHIGAN 4134-2-F eE coswt a(t) - - m (3) where E is the amplitude of the incident electric field (polarized along the z-axis), w is the angular frequency of the incident field, -e is the charge of an electron, m is the mass of an electron. The effect of the alternating magnetic field is neglected. After the distribution function F is determined from (1), the current produced in the gas may be obtained from the relation, 3 J = -n ev = -e Fd vvcosO (4) o z where, J is the z-component of the current, assumed positive along the positive direction of the z-axis, Vv is the mean velocity of the electrons along the positive z-direction, z 0 is the angle between the z-axis and the velocity vector i, n is the electron density. The integration in (4) extends over the entire velocity space. entire velocity space. 77

THE UNIVERSITY OF MICHIGAN 4134-2-F Discussion of the Collision Terms As mentioned in the introduction we shall consider only two types of elastic collisions-electron-neutral particles and electron-heavy positive ion. If the gas is fully ionozed and of high density, then collisions between charged particles of the same kind (for example electron-electron, ion-ion) should also be considered. The electron-neutral particle interaction is a short-range phenomenon and explicit expressions for this are taken from the standard analysis C0l, A.2J. Since Coulomb force is a long-range phenomenon, the cumulative effects of small deflections suffered by the electrons at large distances become very important [1 and they cannot be accounted for in the same way as is done in the first case. Random two-body encounters associated with distances smaller than the Debye length X are assumed to be the sole mechanism for Coulomb interaction considered here. This type of interaction is in general accounted for by the Fokker-Planck equation. The detailed derivation of the collision term from the Fokker-Planck equation is given by Dreicer 013J. Here, we take our collision term after applying the relevant approximations to the general Coulomb collision term give by Dreicer. Evaluation of the Electron Distribution Function In the analysis we shall make use of the Lorentz approximation which states that the collisions between various particles produce a spherically l 78

THE UNIVERSITY OF MICHIGAN 4134-2-F symmetric velocity distribution-small deviations from spherical symmetry are then explained accurately enough by the second coefficient F (v) in the spherical harmonic expansion of F: oo F(cos 0, v) = (v) F (cos) PF ) F ) (cos5) n where, the polar axis is chosen to be the z-axis, P' s are the Legendre polynomials. n 1 0 For the small signal case it is assumed that F << F This condition physically means that the average velocity of the electron gas is small compared to the rootmean square electron speed. This is reasonable if the perturbing field is small. By using the orthogonality property of spherical harmonics the following relations are obtained after the integrations over the angles are carried out. electron concentration n = \F4rv2dv, (6) the electric, current J =- 4 \ F v dv. (7) Z 3 After expanding the term.(t) V F in equation (1) into its spherical coordinate v components in the velocity space and substituting (5) into (1) and equating the terms independent of cos 0 and those dependent on cos 0 the following relations are obtained. 0I I I 79

THE UNIVERSITY OF MICHIGAN 4134-2-F 0 F 0 0 21 (aF>F(aF (8) aF _ ycost a 21 ( )+ ( ) (8)v at 2 av av a t \ at o3v v 'n cc F cont s t lctnntal atcl cll n a ) 2M a (V F v (v)) + 2T a (vv (v) a ) (1 0) at n = /\2at / (9)\ at ' v 'en //cc where eE the mass of the neutral particles, The contributions due to electron-neutral particle collisions are [10, 2,2] i iin n n a 3ion of v; KT a,2 Fhydr at) 2 v(V F - (v))+ -(vv^v)-) (10) 2 Mav e 2av e av =-Fv (v) (11) where at v (v) = 2irNv S (1-cosp)sinPc(I (v) dI, (12) N = the number density of neutral particles, uj (, v) = differential cross section for elastic scattering through the e angle f, M = the mass of the neutral particles, T = temperature of the electron gas K =- Boltzmann constant. The factor v (v) may be identified with the collision frequency of electron-neutral particle collisions. In general ve (v) is a function of v; but for a hydrogen e I. 80

THE UNIVERSITY OF MICHIGAN 4134-2-F plasma v (v) is independent of v for v above several electron volts. The Coulomb collision terms, subject to the Lorentz approximation and the condition << 1, are given by (13], (_F1 n aF 0 ei (13) at 3 c v (at (14) at c where e ee. 2 rei = 47r m ln(X/p), (15) p = the average impact parameter (distance of closest approach between the two colliding particles), e. = charge on a heavy ion = Ze, Z being the degree of ionization, X = the Debye shielding distance, 47lr. 1/9T 109 Coulomb-volt -meter 0 After substituting the relations (10), (11), (13) and (14) into (8) and (9), the following two equations are obtained. taF _ -ycos wt a 21. 1 m a 3 0 (vF)= -(v F (v)) \at 2 av 2 M av e 3v v KT a aF + - (v v (v ) (15) M2 av e av Mv - --- 81

THE UNIVERSITY OF MICHIGAN 4134-2-F aF~~1 0 ycswtit~ CI nV F1 Equations (15) and (16) together determine the distribution function F. it can be shown that in the solution F of () n i terms of F is given by theion F Equations (15) and (16) together determine the distribution function F. In order to solve the two simultaneous equations (15) and (16) we assume that the isotropic part of the distribution function is independent of time. After it can be shown that the solution F1 of (16) in terms of F~ is given by the following n - (v)+ el t 3 1 1 v 3 I F F (at t = 0) e w + W] sinwt av/) [ ~ e () no sine2 v (v) + n /v ~ e o ei +[v -(v)+ 2 cos wt v v + n(v)+ o e 2 3V o el v (v)+ n e /v (v)+ n e o ei v 1 w- ~ Le- T --- n ~i]2 e j (17) 2 + k(v) + 0 e -e 3 v l __ 82 - -

THE UNIVERSITY OF MICHIGAN --- 4134-2-F I It can be seen from (17) that the combined collision effect (electron-neutral particle and Coulomb) relaxes the perturbed distribution function F to the following steady state value 0 o~ r v (v)+ n Ci/v 3 __ __________ e e i F = ( )[ 2 n r i 2 sinwt + v co Vvn 2 n / w+ (V)+ o 3 VV V 3 ut I n% tLu) After substituting (18) into (15) and taking the time average, the following is obtained, 2 aF 2 aF _ 2 1v 2 L- ~ 6v _ 1 m a 2 M av v v (v)+ n F./v3 e 0 el e ____ o i +[v (v)+ n r./V 32 -e o ei 3F (v)+ KT a (2 aF Fve(v + 2 v v(v) - Mv a (19) From (19) F can be written in the following form, F = AexpO mvdv (20) y2 6v (v) e v (v)+ n r/v3 e r no e 2+ KT 2 + ()+ orei] V where A is a constant. If Coulomb collisions are neglected the distribution function (20) reduces to the one given by Margenau [5, [6) for the constant mean free path case. If the thermal energy of the electrons is large and the electric 83 --

-- THE UNIVERSITY OF MICHIGAN 4134-2-F field E is sufficiently small, then the term KT in the denominator of the F0 integrand in (20) is predominant and F reduces to the Maxwellian distribution 0= 2 (2 -mv2/KT F n e (21) 0 \21rKT/ where the constant A has been determined by using the relation (6). Equation (20) indicates that, strictly speaking, F cannot be assumed to be Maxwellian and independent of the field intensity E. As we shall see later this fact makes the electrical conductivity of the ionized gas a non-linear function of the field intensity. We are now in a position to calculate the conductivity of the gas. This is done in the next section. Expressions for Electrical Conductivity After using (7) and (18) the following relation is obtained for the electric current, 47re E 3 = * W ~B 26+ v3-(v +n1i2 sinwt 3m L.2v6+ v v(v) + n ]2 e o ei Defining complex conductivity as = - i t can be shown that + Lve(v)v- +nr.] coswt v6 F0 )v (22) 6 v (v)v3 + +no Defining complex conductivity as or = a - icra., it can be shown that, r 1 - -~ ~ 84

i THE UNIVERSITY OF MICHIGAN 4134-2-F 4e2 3m 47r e 2m r 3m 00 I 00 s wv3 6 (aF 2 6+ v3v (v)nr-2 v dv vv(v)v+n ei. 2 e o e w2v6+ v(V)3+n r.2 v )dV e o eI (23) (24) For the d. c. case w = 0, a. = 0, a = aD. If Coulomb collisions are neglected, and F is assumed to be given by (21), then (23) and (24) can be written in the 2 following forms after introducing the dimensionless parameter x = 2K 2KT 42 00 2 2 47re n 2 o 2 w -x 4 i 3m r3/2 2 2 x dx (71) w +v (x) e (25) 2 0 2 r 3m, )3/2 (7Tr 00 v (x) e 2 2 + v2 (x) e 2 -x 4d e xdx (26) In general v (x) is a complicated function of x and the integrals in (25) and (26) e are not always amenable to integration in closed forms. For constant v, ca. e 1 and ca as given above reduce to the familiar forms. r Conductivity of a Fully Ionized Gas In this section we shall derive expressions for the a. c. conductivity of a fully ionized gas. In order to simplify the analysis it is assumed that F is 85

--- THE UNIVERSITY OF MICHIGAN 4134-2-F given by (21). Since the gas is fully ionized v (v) = 0 and hence (23) and (24) e may be written as follows: eo 2 47re n ( KT\3/2 nei 7 -x 2 o 2 ___ ____ x cr 3m 3/2 Km 22KT)36+2 2 x e dx (27) (r t m ) o ei 0 2 2 4iren3/2 K 10 -x 1i 3m t 3/2 m 2/2KT\3 6_ 2.,2 (x2 1,7r _/2 %\7 / | W \ A *I 11 I. </~~ m / o el For the d. c. case, w = 0 Consequently, a. = 0 and the d. c. conductivity is 1 given by: 4-re 2 (2KT3/2 DD.C 3mi ( 3/2 im ei (ir) 00 7 -x2 xe dx. S X0 - (29) Bu 7 -x2 But x e dx = 3 and therefore, 0o - i - 2 3/2, 47re2 2 /2KT / aD.C. m 3/2 m (30) ei (7r) Expression (30) can be transformed into the familiar form for resistivity (I/or ) D. C. of a fully ionized gas a given by Spitzer [14. We shall now calculate o as given by (27). For this we need to evaluate r the integral -GO - 2 OD 7 -x2 I \ x e d 6 3 dx x +a (31) 86

-THE UNIVERSITY OF MICHIGAN 4134-2-F where 2 3 o ei a = (32) 2 /2KT)3 (3 \ m 2 After a change of variable t x, the integral (31) can be written as follows: oD 1 C t3 t I= 2 3 dt. (33) t+ a 3 3 After breaking the term l/t + a into partial fractions and making use of the 0D tn-t known result t + t- dt = Kn(tl), where K is the modified Bessel function of the second kind and order n, the following is obtained, I 1 1 7r H(1) -7ri/6 (2) aei/6) where H ) and H ) are the usual notations for Hankel functions. Using the 3 3 result C15ithat H3) (ae ) is the complex conjugate to H (aei/6 (34) may r l34 3a2 2( I= - K(a) + 2 pcos0+ psin0, (35) 6a 12a 4Via where H ()(ae /) = pcos0- ipsin. 3 87 -- --

-- THE UNIVERSITY OF MICHIGAN 4134-2-F Thus the following is obtained for ar: 4re 2n 3 /2KT32 a 3(a)_ 1 a r ~J32 2KT) n3/ a + + s0in 7r 0 (36) r 3m J/2 M n 6 Or ) o el L- -i The integral in (28) can be similarly evaluated and the following expression is obtained for o.: 1 2 4ire n 0 2w 2KT i m 3/2 m (r) )/2 9 9/2( 6 a 2 2 np. o ei 1ri/3 -i/ Kg/(-ae )+K/2(-ae T i/3) 12 F (-ae )- K/2 (-ae -i/ 3] (37) If F is assumed to be non-Maxwellian, in general, one has to take recourse to numerical methods in order to evaluate a and a.. r 1 Discussion In the above the electrical conductivity of a partially ionized gas has been discussed in detail. The general expressions derived for conductivity may be applied to specific cases keeping in mind the assumptions made in the analysis. For non-Maxwellian distribution for the unperturbed distribution function and collision frequency depending on velocity, the conductivity has to be calculated by numerical methods. Depending on the type of gas and the degree of ionization, - - - - W 88

I - THE UNIVERSITY OF MICHIGAN 4134-2-F one can, however, make some reasonable physical assumptions which will greatly simplify the calculation. It is evident from the above analysis that even in the small signal case and for weakly ionized gas, strictly speaking the electrical conductivity turns out to be a non-linear function of the field intensity. From this, at least qualitatively, one would expect that in the case of high intensity field the electrical conductivity, if it can be defined in the usual sense, will heavily depend on the field amplitude. This conclusion should be true for weakly as well as strongly ionized gas. The case of the interaction between high intensity EM field and strongly ionized gas will be investigated in the next chapter. 89

I I THE UNIVERSITY OF MICHIGAN 4134-2-F VII NON-LINEAR ELECTRICAL CONDUCTIVITY OF A FULLY IONIZED GAS The electrical conductivity of a fully ionized gas was investigated by several authors [113, C13 for a small d. c. electric field case, in which case the analysis is linearized and small perturbation approximation is justified. However, when the impressed electric field is not small the analysis can be much more complicated and the conductivity may have non-linear behavior. The purpose of this chapter is to study the general behavior of the electrical conductivity of a fully ionized gas based on the Fokker-Planck equation. A method which may be useful in the study of electrical conductivity is developed. An electric field of arbitrary intensity is assumed to be applied uniformly throughout a fully ionized gas which is of infinite extent and spatially homogeneous. To find the electrical conductivity the mean velocity of electrons induced in the plasma is first obtained and from that the induced electrical current and the electrical conductivity of the plasma are determined. In the d. c. case, an instability phenomenon which is called the runaway effect can be observed. This instability automatically restricts the intensity of the impressed electric field to be lower than a critical value. Under this restriction, the d. c. electrical conductivity is obtained as a function of the intensity of the electric field and other parameters. In the a. c. case, the runaway effect loses its significance if the impressed frequency is higher than a critical value. For high intensity micro I m 90

THE UNIVERSITY OF MICHIGAN -- 4134-2-F waves, the electrical conductivity can be very non-linear. In this chapter a simple way of obtaining an approximate electrical conductivity is also presented. The Basic Equations An electric field, E cos wt, is assumed to be applied along the x-direction to a fully ionized gas which is of infinite extent and spatially homogeneous. The basic equation which describes the system is a Fokker-Planck equation as follows: af f --- fe eE fe a 1T 2 - -- coswt = -e-(f <AV.>) +- - (f <Av.v.>) at m av.-av. e 1 2 - av av. e 1 j e x 1 i,j i j (1) f is the velocity distribution function of the electrons of the plasma. The spatial variation term is neglected so that f is to be determined as a function of velocity, v, and time, t, only. The right hand side of (1) represents the collision effects of electrons with positive ions and electrons due to the friction and the dispersion in the velocity space. These two terms can be represented I by - v e (fe<Av>) = - - v (fe a ) (2) 2 2 a2G (f (3) E a2^ (f <<v.Av.>) = a (f e 2,j.i j. 1v.v. e jv. v. iJ i J i,J i j i j and "~f (t- v ', f t) 3" ()' H = H + H (= r -_ i d + 2 \ d V (4) e ep ee - v - -~- -- - -- -- --- - -- - -- -- 91 ~_^ ___ _^^ ^__^.__ ___ _^_^~~

i THE UNIVERSITY OF MICHIGAN 4134-2-F G e =G+G = rev- v' fp('d v f 't)d3 (5) e ep ee e p e 47r 4 h ' e 2 p m o e h = Debye shielding length p = average impact parameter for a 900 Coulomb deflection. The positive ions of the plasma are assumed to have a Maxwellian distribution and undisturbed by the impressed electric field. The velocity distribution function of the positive ions is 3 22 a -a v f n P — e (6) p 0o 3/2 m with a =, n = unperturbed density of the plasma. p 2KT' o With the information expressed in (1) to (6), we are, in principle, able to solve for f. However, it is hopelessly complicated, especially if (1) is e considered without making any approximation. Fortunately, for the purpose of studying the electrical conductivity we can find it much easier to solve the problem by integrating (1) in the velocity space and solving some moment equations which are thus obtained. First of all, the operation Ed v on (1) gives I 92

THE UNIVERSITY OF MICHIGAN 4134-2-F an (faf\ eJ \ V e) d v = ~ (7) at - at /coll. This is valid as long as the conservation of particles during the collision is asserted. Equation (7) implies that n is independent of time. This is actually e a necessary condition if the spatial variation of f is assumed to be zero. We e can, therefore, let n = n (8) e o Secondly, the operation v d3v on (1) gives x a eE af> 3 - n u+ coswt n =\ - vdv (9) at e m- e at C x e coll. where neu f v d3v e e x In carrying out the right hand side of (9) the following facts are used. (1) The collisions between like particles do not alter the total momentum of the parent gas. (2) The dispersion in velocity space leaves the momentum of a gas element unchanged. Based on these reasons and using (2) to (5) the right hand side of (9) is simplified to the following form. 93

THE UNIVERSITY OF MICHIGAN 4134-2-F d 3 3 ad v = a (f - ) dv v = f v d - re e av. r IV t7 3 a $ ~f ( t) 3 T h e c r u c i a lf evt o f jn g \ i t ev a i r d v s I e dv Iva In deriving (1) is to be solvsumed. From it the mellian andveloity of electrons, u, c be an beven function of vand v The subSimplificatution of Equation (119), with (8), yields f( e' r, f ("t) t) M(at)t= 1 e d v' p e Simplification of Equation (1 1) f ('~" t) M ('V, t) b- -,[ d3v' L - -- 94

I I I i I THE UNIVERSITY OF MICHIGAN 4134-2-F can be evaluated exactly if f is assumed to be Maxwellian as expressed in (6) That is 2 22 3 -a2V2 0 " 3/2 \d v' This integral is analogous to the potential integral in electrostatics. Analogously the integral gives the potential at the point v in the velocity space produced by a spherically symmetrical charge density expressed as 22..3 3/2 -a V 47rn (a 3/r 32) e o p This potential is well known in electrostatics and the exact answer is 3 2 2 2 oo 2 2 - a -a vI -a vt M(v) = 47rn -7 I i e P v'2dv' + e P v'dvI 0 (0)3/2 = n - erf(a v) (12) ov p The right hand side of (11) then becomes r.h.s. = r d3vfe a erf (a v) v(vex r r -V~~~~~~~~2 2 ~3 1 2 ap -av - \ d v f erf (a v)- e (13) e xe 3 p vL P 2 VV 95

I THE UNIVERSITY OF MICHIGAN 4134-2-F Therefore, (11) finally becomes 22 au eE co3 1 2 aj - P - + m — S \t =- r d vvef e erf(ca v)- 2 e 3t m exe 3 p V7 2 e L-v v ~ (14) Approximate Electrical Conductivity For the purpose of showing the usefulness of (14), the approximate electrical conductivity will be found in this section. We can interpret the r. h. s. of (14) as ca -a2 v 2 vv f erf(ca v) - 2- e -— n duR(v) xe x 3 P 3v2 m Lv e (15) where 1 3 du - f v d v = density of the mean velocity of electrons in the n e x o velocity space, 1 r 3 u - \ f v d v - total mean velocity of electrons, n e x ~v 0V)= 2 03efo )-A 0 \e 2 _12 P ["m r 2 a-av R(v) = e erf(v) - 2 f e p (1f 2 p 2 i) e Lv v J - electrical resistivity in the velocity space. For practical purposes, in many cases, we can write (15) approximately as 96

THE UNIVERSITY OF MICHIGAN 4134-2-F 2 \ ~~~2 2 e2. e 3e - -- d vv f R(v) - - R(V) d vv f =- -R(v)n u (17) m xe m xe m o ^ev e e where v is an appropriate value of v which is to be substituted in R(v) so that (17) is valid. For example, v can be assigned as the rms velocity of electrons. If the r. h. s. of (14) can be approximated in the way expressed in (17), (14) becomes -t + coswt - - n R(v u (18) at m o e e Equation (18) gives a steady state solution immediately as eE 1 e2 u= -m L o i2 — n R(, ) coswt + wsinwt (19) n R(~) + ~ m o L- e _ From this the induced current, J = -en u, and the complex electrical conductivity, o = J/E = - r.- are easily determined. The approximate electrical conductivity is then 2 e 2 n R(v) n e m o -r = L-e oR(- )+w (20) m o L e J 2 n e 0 W ori = F2 (21) e e _,J 2 p 97

I THE UNIVERSITY OF MICHIGAN 4134-2-F and 2_2 R(v ) = - erf(a - 2 e 2 -3 p VP -2 e v v 2KT a \= m e p 2KT In the next two sections, the electrical conductivity for the d. c. and a. c. cases are analysed more exactly from (14). The non-linear behavior will then appear. D. C. Electrical Conductivity of a Fully Ionized Gas For the d. c. case, or w = 0, (14) can be reduced to rr7 P r2 22 eme i 2 a -a E = - -e \ dvv e rf ( v)- e fpv (22) e xe p 2 J If the electrons acquire a mean velocity, u, after a d. c. electric field, E, is applied to the plasma, f can be assumed to have a form as e3 y2 e 2 2,2 (23)2 f = n e (23) e (o )3/2 \m with a \2K e 2KT Upon the substitution of (23) in (22), the mean velocity, u, can be determined. In evaluating the integral in (22), it is learned that I w 98

THE UNIVERSITY OF MICHIGAN 4134-2-F m a v aplae = \| >> 1 p p e m e Hence, (22) can be approximated as 3 00 oo -a2 rv u)2+v2+ v2 e e pe eex y zJ E E = - n Ldv udv \ -dv v ] (24) e 0 3/2 z x 2+v v+ v2 /2 E ton ) so tt te e elot u i o-lerl y zd with a very high degree of accuracy. The integral can be evaluated exactly (details are included at the end of this chapter) and (22) finally becomes 22 I m a -a u E -- n erf(a U)- e 24) Equation (24) shows that the mean velocity, u, is non-linearly dependent on E and other parameters. u can be solved as a function of E from (24) at least numerically. After u is determined as a function of E, the current and the electrical conductivity are obtained immediately The behavior of the induced current, J, as a function of the electric field, E, is shown graphically in Fig. 1. We observe a very important phenomenon at this point. That is when E is increased higher than a critical value, E, the induced current increases monotonically and shows instability. This effect is called the runaway effect. If can be found numerically that E corresponds approximately to a u = 1 and its value is e __ 99

I THE 3.0 — 2.5. UNIVERSITY OF 4134-2-F MICHIGAN to a lV 2. 1. 1. 0. 0.1 0.2 0.3 0.4 eme o e e o e FIGURE 1: THE INDUCED CURRENT IN A FULLY IONIZED GAS AS A FUNCTION OF D.C. ELECTRIC FIELD 0. 5 100

THE UNIVERSITY OF MICHIGAN 4134-2-F 't m m 3, n E =0.43n)( ) = 0.86 T (25) C e 2K K p T 0 Because of this phenomenon the analysis of the d. c. electrical conductivity of a fully ionized gas should be restricted to the case of E < E, or u < \2 (26) c m e if sensible results are expected. Under the condition expressed in (26), (24) can be written as em 2 3 2 2 23 1 4 5 E =n OF o L3 U - 7u 7 e If the d. c. electrical conductivity is defined as J u u ^.Cr = E"'-en - 'D. C. E -en E OD.C. can be expressed as follows: 3V e2 (2KT 3/2 277 ( e 2 /2KT 2 + 3584 -2 m non [(m) E (27) T eo o e ad o This expression shows that ogD C is dependent on E, n and other parameters. When E is very small cD C reduces to 101

THE UNIVERSITY OF MICHIGAN 4134-2-F DC _ e2 (2KT\3/2 D.C. 4 me \me e e e which is about 1/2 that obtained by Spitzer [1. A. C. Electrical Conductivity of a Fully Ionized Gas If f has a form as expressed in (23), (14) can be written as e u eE 12 e -au a+ - cos wt = - enI \- erf(a u)- e (28) stable solution for u is to make the linear term always larger than the non-linear term in (23). That is!2 rQ 2 22 -au p 1, 2 e e a- > n l erf(a u) - \ e et eo 2 e V ur ~au Lu Since at is roughly equal to wu and the right hand side of the above inequality has a maximum at u, the following relation can be obtained. That is a e r 2 e e22 a -aU ] Wu > Fn - erf(a u) - e - e e o e e 7 u 3u and a critical value of o is obtained if u is replaced with as follows: a e - - ~ ~ ~ ~ 102

I THE UNIVERSITY OF MICHIGAN - 4134-2-F 0 2 ( 3 /2 el > w = 0.43n ['3 = 0.43n e h471e ~c 0 'e 3 02KT 4 3 2 /\pe em / 0 e (29) Therefore, we can state that in the analysis of the a. c. electrical conductivity of a fully ionized gas the electric field E can be of arbitrary value if w > w. In the actual case, the determination of u from (28) is quite involved and only special cases are discussed here. 22 When E is small and a u << 1, (26) can be written as e _u + 4 3n eE u -- + (4- 'na u ==- -- acos wt (30) at \3VF eoe m e The steady state solution for u is eE 11 = - 1 r r 3 ] n a coso t + wsin wt oe J (31) m /41 2 3VT e e (4 ' n a3) + c t3V- e 3 e o e/ J With the definitions J = -n eu, and o = = - - i o E r conductivity can be found as 2 4, 3 n e, naoe o 3 Vt- e o e r m 7 4 n3\2 2 e +- *, n a + \3Vt- e o e., the a.c. electrical (Yi (32) 2 n e o W a=i m n4 rna3 ~2 \3VT eoe (33) J -- 103

I THE UNIVERSITY OF MICHIGAN 4134-2-F When E is large and a u 2 1 but w > w, (28) can be written as e c au eno at 2 U eE oswt. m e (34) u can be obtained if this non-linear differential equation can be solved. Evaluation of a Definite Integral A definite integral appeared in deriving (24), the evaluation of which was omitted to avoid obscuring detail, but is presented here for completeness. The following integral is to be evaluated as a function of the parameters a and c: I = x dx dy \ dz / -0o G-O - aO -a 2x - c)2+y +z e [2 2+ 2 3/2 [x ++zJ (35) Let: z' = ax x' = az s = ac (36) y' = ay 00 ~ jJOD 4JO e- z' s)220,2 ] -00 -00 +00x I I- z'dz' fy de y ' 3 z+ (y2+ 3/2 z')2+ (Y ) (x') Note that aI is a function of s only. Define Il(s) = aI Introduce spherical 1eieI~)= cI nrdc peia coordinates: --- W 104

m THE UNIVERSITY OF MICHIGAN --- 4134-2-F x' = rsinecos0 y' = rsinOsin0 z' - rcos (37) One obtains: 2 7T I (s) - do 0 7r sin Ocos OdO 0 00 \dr 0 - r- scos0)2+ s sin2] e (38) Note that the integrand is independent of 0. Let t = r- s cos 0. It follows that: rT 2 2 I (s) = 2r es2sin 0 0I Integrate by parts: s cos 0 2 t2 sin 0 cos 0 d 0 + e [~~^ dt] (39) = s cos 0 2 + et dt 4) 2 2 du * - s sinos 0 du - -ssin0 e dO 2 2 -s sin 0 dv - e sin cos 0 dO 2 2 -s sin 0 -e v = 2 2s I1(s) = 27r -2f2 2 2s 2 2 -i -S2 sinO eS d (40) 105

"" THE UNIVERSITY OF MICHIGAN 4134-2-F 2~s)- ~ Fl-fs2S1 I1(s) = 2 erf s - 2s eys (41) In terms of the original parameters we write I: r2 2 -I(a, c) = 3 2 r erf(cac) - 2ac e a(42) a c 106

THE UNIVERSITY OF MICHIGAN - 4134-2-F APPENDIX A NON-LINEAR MODELING OF BESSEL FUNCTIONS This appendix treats the non-linear modeling of two systems described by Bessel functions. The domain of regularity of the modeling function or transformation between the systems is investigated and an explicit series representation of the modeling function is obtained. The corresponding results 123 for trigonometric functions are included (and expanded), not only for comparison, but also because they suggest the structure of the coefficients in the modeling function series. Statement of the Problem Suppose a (physical) system is described by the following differential equation with initial conditions: 2 e+ + + K2y = (1) dt2 t dt 1 at t - 0: y(0) =1 t (0) = 0 (2) dt Another system is described in terms of the independent variable s: d2x 1 dx 2 ds + - + Kx x=0 (3) a2 s ds K2 ds at s = 0: x(0) = 1, (4) ds 107

I THE UNIVERSITY OF MICHIGAN ---- 4134-2-F Let the independent variables of the two systems be related by a change of scale: s = at (5) The initial conditions on x remain unchanged, (4), and the differential equation for x as a function of t becomes: dx 1 dx 2 2 2 t + a Kx = 0 (6) dt2 t dt 2 A non-linear modeling function is a transformation, y Ix0,giving the solution y(t) in terms of x(t). We assume that y (Cx is not explicitly a function of t. Thus dy d dx dt dx Yx ' dt (7) Essentially, finding y [xD means eliminating t from the two systems. Suppose the solution of (1) and (2) is: y = g(Klt) (8) It is not difficult to show that x(t) is given by: x = g(K2s) = g(K2at) (9) We can obtain y (xJ formally or symbolically as: K1 ' LK2a -1 (10) y[xj=g - g-l '(x] = g[Kg(x] (10) where K K - (11) 2a I 108 --

- THE UNIVERSITY OF MICHIGAN 4134-2-F The characteristics of the transformation y [x] thus depend on the parameter K. This can be seen also by making the change / = aK2t in (1) and (6). Restrictions of K which make y (x) single-valued are called similitude conditions C1] We will attempt to investigate these conditions by expressing y 0x in a power series whose coefficients depend on K. Mappings Before proceeding to the development of the modeling function, y [x), we digress to a discussion of the systems (1)-(4) and the corresponding situation for g = cos t in which we suppose that the solutions are known. Strictly speaking, if non-linear modeling is to be fruitful, we cannot solve the original equations analytically, at least not in a practical representation. We hope to find the transformation y (x] between two "unsolvable" systems and by performing one, x(t), experimentally, then to deduce a solution, y(t), of the other through y (Cx. However, one must distinguish between determination of modeling transformations on the one hand and the theory or study of such determinations on the other. The latter, of course, properly includes some questions of technique as well as of existence and limitations. By considering modeling of "known" equations we hope to develop general principles and methods more readily. In later -- - 109 -----

THE UNIVERSITY OF MICHIGAN 4134-2-F sections we determine modeling functions for systems corresponding to Bessel and trigonometric functions without solving the original differential equations. In this section, by assuming known properties of the solutions, i.e. y(t) and x(t), we obtain information which allows us to interpret and extend these later results. Consider then, the solutions of the systems (1)-(4). (We assume the simplification of parameters introduced previously.) y(-) J (K7) (12) x(7) J (7-) (13) The modeling function is obtained symbolically as y [x = J [KJ1 ()] (14) It is clear that y (x] is not necessarily uniquely defined for all values of x and K since J1 is a multi-valued function. Even so, it is conceivable that y IxJ 0 is single-valued for some choices of K. Consider systems defined by: Yl(T) = cos(KT) (15) x1(7) = cos(7') (16) Y X co- [ cos l(x) (17) 1 When K is an integer, n, (and for no other choice) yl Cx} is single-valued 1 2.23 K = n =4> y [x] = T (x), Tschebyscheff polynomial (18) 1 n m i _ _ ~ ~ ~ ~ 110

THE UNIVERSITY OF MICHIGAN 4134-2-F We see that not only is y1 xD single-valued, it is actually an entire function of x considered as a complex variable. For example, the case n = 2 is: 2 Y1 x = T2(x) = 2x- 1 (19) To obtain all values of x and y, however, we must consider that 7 in equations (3) and (4) also takes on complex values. For suppose 7 is real. Then Ixl does not exceed unity and we obtain only a portion of the function defined by (19) 2 Y1 X = T2(x) = 2x - 1 Ixl < 1 (20) Suppose T' varies over a strip: t = o + ir7 (21) 0. < r0 < 7 < 0o (22) 0 < oa < r, -oo < r7 < (23) Now we obtain all complex values of x just once, i.e. the function cos7 maps the strip onto an entire x plane. Hence cos x can be defined so that it is one-valued onto the strip. For this restriction we would expect yl (x] to be unique for all K (and not just integer cases). What is the behavior of this function? To investigate singularities we compute the derivative of y1 [x]. Here, of course, we are again freely using knowledge of the solutions of the original systems. dy (7) dy, _ d7~ -Ksin K7' Yl -T - (24) yj dx dx(r) -sin 7 d1' id 111

THE UNIVERSITY OF MICHIGAN 4134-2-F Singularities will occur only when sin 7' is zero but not necessarily even there. If sinK7' is also zero the singularity is removable. This happens when K is the ratio of two zeros of the sine function, in this case, when K is an integer. The question arises: suppose we consider a region analogous to the strip for the Bessel function, i.e. a region of the complex T plane which maps under x = J (T) onto the x plane one-to-one. We have here again a proper O (i.e. single-valued) definition for J (x). Thus the function y Cx) can be O properly defined for all complex x (and for any K) as long as x results from taking 7" within the fundamental region. This gives rise, in general, to many different possibilities (one for each fundamental region) for y x], each restricting the range of 7. For example the fundamental region containing 7 = 0 and small positive values of 7' extends along the real axis up to 1= 1 such that Jl(7)= 0. If we represent y x]J for this 7T-region, say in a power series, we would expect a singularity (divergence) when x = J (71) unless K were a ratio 0 of two zeroes of the function J (71) since dy ~ -KJ 1(K7') (25) dx J( ) (25) X. We expect these numbers K. --,where J (X.) = 0, to be significant in the power series representation for y just as in the trigonometric example. power series representation for y [x] just as in the trigonometric example. 112

THE UNIVERSITY OF MICHIGAN 4134-2-F Taking these K values is no guarantee that various branches of y [x] will have the same representation, however. The periodicity of the trigonometric function is apparently what allows this for yl Cx. Derivation of Differential Equations for the Modeling Function No loss of generality (see (11)) is suffered by lumping all parameters into one constant K in terms of which the similitude conditions are described. Consider, then, y(t) and x(t) defined by: y + J + K y = 0 (26) t x + + x 0 (27) t at t = 0: y(0) = x(0) = 1 (28) y(0) = k(0) - 0 It is desired to find the modeling function y [x] in terms of K (i. e. by eliminating t) without solving the differential equations. From (26) and (27) we determine a differential equation for y x]. From (28) we can obtain sufficient boundary conditions to determine y Cx] uniquely. Let primes denote derivatives of y with respect to x. Since y(t) = y[x(t)] we have: y = yt' (29) y - y"k2 + y'x (30) - __ — ~ 113

I -. THE UNIVERSITY OF 4134-2-F MICHIGAN From (1), (2), (4), and (5):.2y - xy'+ K2Y-~ (31).2, - +2 = 0 (31) Since the modeling functions are non-linear in general (y" # 0) we see that k is a function, say f, of x only:.2 K 2 x2 = y' - = f(x) (32) y" f(x)y" - xy' + K2y 0 (33) From (27) and (32) we determine an equation for f(x): 2x f' = = 2x = -2x- - (34) x t f + 2x =- (35) t Differentiate again with respect to x: f" +2 = - f(f"+ 2) = 2 2 t 2x xt (36) (37) 2(k)2 _ 2 t - f' t Substitute from (35): 1 2 f' f(f" + 2) = - (f'+ 2x) + - (f'+ 2x) f(f"2) = (f' 2x)(f' f(f"~2) = (f'~2x)(f'~x) (38) Finally: (39) 114

THE UNIVERSITY OF MICHIGAN 4134-2-F From (39) we can obtain f as a function of x. Then (33) determines the modeling function y fx. From the conditions (28) on y and x as functions of t we find at x = 1 (i.e. at t = 0): y [1] 1 (40) f(l) 0 (41) f'(1) = 2x (42) from which (see below) f'() -1 (43) In (42) we have encountered a removable singularity in the term -. A similar situation occurs in determining y' at x = 1: y' = (44) From (1) and (2) we find that, as t -- 0: x + _- + x-> x + + x x-> 0 (45) t x 1 '.. x - 2 - - (46) 2 2 Similarly 2 2 y + + Ky ->. 2y + K y - 0 (47) t K2 K2 y -~ -- (48) 2 2 Hence from (44), (46), and (48): y'(1) = K2 (49) 115

THE UNIVERSITY OF MICHIGAN - 4134-2-F Summarizing from (33), (39), (40), (41), (43), and (49) we have the following system for the modeling function y [x3: fy" - xy' + K2y = 0 (50) f(f" + 2) = (f' + 2x)(f' + x) (51) y(l) = 1 y'(l) = K2 (52) f(l) 0 f'(1) =-1 (53) Series Solutions for the Modeling Functions A brief description of the situation corresponding to g = cos t will be helpful. These results have been obtained elsewhere, 13, [2] but the emphasis here is on the power series representation. Let: y y() = 1, y( 0 ) = 0 (54) + x = 0, x(O) = 1, (0) 0 (55) y = y ix) again implies (31), fy"- xy'+ K2y. (56) Since f = x we obtain, using (55), a differential equation for f as a function of x: f'= = - 2 x = 2 -2x (57) f = c- x2 (58) 116

THE UNIVERSITY OF MICHIGAN 4134-2-F when t = 0, x = 1, and f= 0. Hence f(l) = 0. f = 1 - x (59) Boundary conditions on y [xJ near x = 1 from (54) are: yl] = 1, y10 = - = - K2 (60) x e e A power series solution is now written down for the equation for y Cx] using conditions (60). (1-x)y" - xy' + K2y = 0 (61) (62) y x = + s (K2)(x- 1) (62) S = K (63) K2 (K2- 1) 2 3 -32: (64) K2 (K2- 1)(K2- 4)( S3 = 3 5- 3! (65) m-= (K2_ i2) Sm 7 m 3 5....(2m-l) (66) i=0 The zeroes of s are just the integers with absolute value less than m. Thus m they are cumulative and trivially approach these values as m increases. Consider the zeroes, ri., of the derivative of g, 0 W 117

THE UNIVERSITY OF MICHIGAN 4134-2-F g(t) = cost (67) g(t) = -sint (68) sinr-i = 0 (69) = Ti i = 0O, +1 2,... (70) r)= 7T (71) i (72) r1 We find that the coefficients of the power series expansion of the modeling function are polynomials whose zeroes approach the ratio of the zeroes of the derivative of the function being modeled to the first such zero. If we take K = i the series (62) becomes finite. It is the series for T.(x) (Tschebyscheff polynomial) expanded near x = 1. A series solution for y x)] corresponding to the Bessel function is now obtained near x = 1 by first obtaining a series for f(x) from (51) and (53), then substituting this into (50), and using (52). Let: f(x) = a(x1 (x73) n=o Inserting this series into the equation ff" - (f) = 3xf' - 2f + 2x, we obtain, near x = 1:. 118

I THE UNIVERSITY OF 4134-2-F MICHIGAN - I a (x- )n n a m(m- )(x- 1 )m-2 m m=o n n=o n=o a (x- l) n n a (x- l)n n m=o o3 n=o m a (x- l)m-1 m n a (x- 1) n Iaxu - 2=o n=o a (x- )n+ 2(x-1)2 + 4(x-1) + 2 n (74) After computing the Cauchy products and changing indices we have: n=o =o a. a _ -(j-l)-j(n-j) (x-l) -2 -- (3n-8)a (x-1)n 2 n-2 HE,,,_,,, + 3 (n- 1)a nn=i n-2 (x- 1) 1 + 2(x-1)2+ 4(x-1)+2 Equating like powers of (x- 1) yields: a2 0 0 ' ala = 3' 0' a0 2 -a1 = -2a0+ 3a1 + 2 (75) (76) (77) (78) -- 119

- THE UNIVERSITY OF MICHIGAN 4134-2-F -2aa+ Oa2a + 6a3a = a + 6a2+ 4 (79) -3aa -2a2 + 3a3a + 12a a0 = 4a +9a + 2 (80) aa3 2 3 1 40 2 3 Finally for n > 5: jaja j(2j-n-1) = (3n-8)a 2+3(n-l)a (81) 3= The conditions on f(x) obtained in (53) imply a0 0 (82) a = -1 (83) Clearly (82) and (83) satisfy (76), (77), and (78). Equations (79) and (80) determine a2 and a3 uniquely: a2 = -3/4 (84) a3 -1/72 (85) Noticing that (81) is a linear equation in the highest subscript of a for n > 5, we see that the coefficients are uniquely determined in a recursive way. Since we need them later we list the values of a4 and a5: a 1 2 82 a = ---- = (86) 4 288 (4)2 -13 4 - 52 a= 2 2 (87) 5 3.(5!)" 3 (5)2 120 --

- THE UNIVERSITY 4134-2-F OF MICHIGAN For this application (i.e. starting with a = 0 and a = -1) (81) can be simplified: 2 a. a j(2j-n-2) -= (3n-5) a + n a j n-(j-1) n-1 n n >4 (88) Assuming a series solution for f(x) near x = 1 is known, we can solve fy+ Ky = near x 1 as follows. Let fy" - xyt + K y = 0 near x = 1 as follows. Let yx] -= b (x-l)m m=o (89) From the equation for y x]: n=o a (x-l)n n ~ m(m-l)b (x-l)m-2 m m=o m=o m=o k mb (x- l)mm mb (x-l1) m+ K m ~ m m-2 (90) i(i- 1)b..(x-l)m 2 =0. 1 m-i - > (m-2)b ]~~~~~~~ -> (m- 1)b (x- 1)m-2+ K2 m -1 b (x-1) )m-2 0 m-= (91) Equating like powers of (x- 1) we obtain: 0 - bao =0 O bal + O' blaO- b- 0 (92) (93) 121

THE UNIVERSITY OF MICHIGAN 4134-2-F Equations (92) and (93) impose no restrictions. For m 2, i(i -)b a - (m-1)b + (m-K2-2)b (94) i rnm-i m-1 m-2 We can determine the coefficients b uniquely in terms of the a. A m m few terms are calculated for reference: (We assume a0 is zero for all applica tions ) 2 b = K b (95) (This agrees with (52).) 2 (K - l)b1 b2 -- (96) 2 2(1-a ) (K - 2 + 2a )b 2 (97) 3 3(1 -2al) (K2- 3 + 6a b 3 + 2a b b =2 3 32 (98) 4 4(1- 3a) Taking the values of a previously found, putting them into (96) - (98), and continuing the process we find an expansion for the modeling function, y x, near x= 1. m=o m=o (mi) - -- -- - 122

I THE UNIVERSITY OF MICHIGAN 4134-2-F R1 R2 R3 R4 R5 R6 K2 -K 2 2 = K2(K- 1) - K2(K2- 1)(K - 7/2) - K2(K2- 1)(K4- llK2+ 26) - K (K2- 1)(K6 24K4+ K24 997 3 3 22 2_2 8 6 1287 4_ 7263 K2+ 13, 001 K(K 1)(K 44K+ K 2 2 - ~2 2 2 (100) The zeroes of these polynomials possess some very interesting properties. Up to m = 6 they are all real. Their squares are listed in Table II. TABLE II m Squares of Zeroes of R Ordered in Magnitude 1 0 2 3 4 5 6 0 0 0 0 0 1 1 3. 5000 1 3.4385 1 3.4103 1 3.3944 7.5615 7.3735 7.2799 13.2162 12.8431 20.4826 (>3 0 1 3.3523 7.0493 12.091 18.4772 - ~ ~ ~ ~ ~ ~ ~ — ~ 123

THE UNIVERSITY OF MICHIGAN 4134-2-F X2 The last line is a listing of the quantities, where the X. are defined as the non-negative zeroes of the Bessel function J1. The natural and X. obvious conjecture is that the zeroes of the polynomials R approach -. This conjecture has not, as yet, been proved. However, similar results are obtained for other cases and a result of this nature is believed valid for a rather large class of modeling functions of the type: y = gKg (x)] (101) 124

THE UNIVERSITY OF MICHIGAN 4134-2-F APPENDIX B LOCAL ALTERATION OF ATMOSPHERIC DENSITY WITH ELECTROMAGNETIC ENERGY Summary It is assumed that by some means a spherical domain in a homogeneous atmosphere can be uniformly heated so that the density reduces by a factor 1/7y at the center of the sphere. First the steady state with this density ratio is considered and the power needed to maintain the steady state is evaluated. Second there is assumed uniform atmosphere for t < 0, the power supply starting at t = 0. Then the time is evaluated after which a density ratio of l/y at the center of the sphere is reached. The density and temperature distribution are given and a simple formula is derived for the power supply needed to achieve a density ratio l/y in a spherical domain of radius a and in a time interval t. For 7 = 10 and a = 100cm, the total energy to be delivered until this density ratio is obtained is of the order of 3.10 Joules. Problem I In an infinitely extended gas a spherical volume of radius a is constantly heated with a heat production P per unit volume and unit time. Determine the steady state and the density distribution. Make the center of the sphere the origin of a coordinate system. Since we are concerned with the steady state the continuity equation is trivially satisfied. J 125

I i THE UNIVERSITY OF MICHIGAN 4134-2-F The equation of motion gives a constant pressure throughout the whole space. The energy law gives - V. (VT) = P, (1) where a is the heat conductivity. Finally we have the equation of state which is pT = p T. (2) 00 00 The suffix o refers to infinite distance, p is the density. The heat conductivity is practically independent of the density, however 1/2 varying with temperature as T /. We neglect this variation, which will not change the qualitative picture. We have then in polar coordinates dT 2 dT 2 + r dr -P (r < a) (3) =0 (r >a) We have the boundary conditions T and - continuous at r a, T T at r = oo. (4) dr co The unique solution is P a3 T T + - (r > a) (5) oo 3S r T T + (3 a2 r2 (r < a). (6) oo 6is The density distribution is then given by. 126

- THE UNIVERSITY OF MIC] 4134-2-F 3 P =P [1+ P a 3 p p I + T- (3 a2- r HIGAN (r > a) (r < a) (7) (8) The problem contains only one essential parameter which is Pa2 a 66-T 00 (9) Put r = as, then P_ P =+ 21+ s 1 (s > ) (s < 1). (10) p0 1 + a(3- s2 (11) Figure 2 gives a rough sketch of p/p for a = 1 and a = 3. 1 (I — 0.8 0.6. P/P00 0.4. 0.2 0 a = 1 _.^ ^ ^ ^ a= 3 0 1.0 2.0 3.0 4.0 5. 0 S = r/a FIGURE 2 -- 127

----- THE UNIVERSITY OF MICHIGAN 4134-2-F As is seen, one has p(0)/p = 0.1 if a = 3. Then p/p = 0.5 for s = 6. The value a = 3 is equivalent to 2 Pa = 18oT.(12) 00 The total production of heat per unit time is 47T 3 a P = 247r auT.(13) 3 oo The heat needed to make p(0)/p = 0.1 is therefore proportional to a, that is, 00 to the radius of the sphere inside which an appreciable drop in density is desired. It is remarkable that a, C, T are the only parameters entering the necessary power. A rough estimate is given with o " 2. 10 erg/cm sec K, T = 300 K, and a= 100cm. 47 a3 P 5.O1 ergs1 = 5000 J s1 (14) 3 Of course, one has also to compensate for the loss by black body radiation, which is certainly not negligible for T(0) 3000 K. Hence the above result gives a lower limit of the power needed to maintain the steady state. Problem II In an infinitely extended gas a spherical volume of radius a is constantly heated with a heat production P per unit volume and time. Determine the temperature, the density and the velocity field in the initial stage after the power 128

""" THE UNIVERSITY OF MICHIGAN 4134-2-F supply has started, if the gas was at equilibrium before. Assumptions: 1. viscosity can be neglected 2. heat conduction can be neglected 3. the pressure is practically constant While the first assumption is certainly a good approximation, since there are no boundaries and hence no boundary layers, the second assumption is certainly not generally correct. It is, however, a good approximation in the initial stage and the duration of the initial stage can be well defined. The temperature conductivity k =- is approximately k = 2 cm /sec; this says that heat spreads over a Ppc distance of 1 cm in about 1 second, or over a distance of 100 cm in about 10, 000 seconds. In our problem the initial stage is determined by the condition that. heat conduction takes place only over a distance which is small as compared to a. This is equivalent to t a 2/k. (15) As long as t satisfies this condition the neglect of heat conductivity should not alter the results essentially. The third assumption is usually made in the theory of flames and its justification is given there. It can certainly be applied in our case as well. As a consequence we can disregard the equation of motion which would only give the small accelerations in which we are not interested. I —1 129.

THE UNIVERSITY OF MICHIGAN 4134-2-F First we have the equation of state which is written pT = p T (16) 00O I p and T are the initial values of the density and of the temperature respectively. Obviously we have p = p and T =T, where p and T refer to infinite o oo o Too 0ro oo distance. The continuity equation is dP + p V v dt (17) with d/dt being the substantial derivative. d a. 7 d = - + v dt At (18) The temperature satisfies the equation dT pc d P(r < a), ~0 p dt and c is the specific heat at constant pressure. From (16), (17), and (18) we obtain Pd dt- - - Tp v (r < a), dt P (r > a). (19) I 0 = -Tp 7- v (r > a) (20) or V 'v - P (r < a), V v = O (r > a) 0 o p Because of spherical symmetry we may assume that v has only a radial component which depends only on r. Then (20) can be easily integrated to give 130 --

THE UNIVERSITY OF MICHIGAN 4134-2-F 3 r a v = <r. v (r <a), v = (r > a) (21) Tr, which has the dimension of a time, is defined by 3 = c Tp. (22) P p oo In (21) we have assumed that v is continuous at r = a. This can be concluded from T and hence p being continuous as well as the mass flow p v. With the velocity field as given in (21) we integrate the continuity equation (17), which we write now in the form + V. (p) = 0 or apf + 1 (r pv) = 0. (23) at at 2 ar r We have to distinguish three domains of the variable r. The first domain is 04 r. a. Upon insertion of (21) into (23) one finds p = p e 3t (0. r a). (24) The second domain extends from a to ~(t) where V(t) is the distance from the center which has been reached by the gas originally at r = a in the process of expansion. Beyond r = E(t), the gas is at temperature T, since heat conduction is neglected and therefore the density is also unchanged. We have now to determine (t). In the interval a < r < I, we have from (21) and (23) - 131

THE UNIVERSITY OF 4134-2-F MICHIGAN 3 ap + a _ = o at -r2 ar;r (25) Hence p is a function of t we obtain pr = -3t/U p(r) = p e o - 3- which assumes the value (24) at r = a. 3a Thus + r /a - 1 (a, r 4 (t)). (26) Now we are able to determine V(t). It is given by p() =- p or = a3 [ + t (27) 3 Figure 3 gives the density as a function of r/a for -/t = 3/4, p(0, L )=0.1 Po; 1 1.5a, and for t/r = 1/3, p(0, -U) 0.36p, 3 0 t = 0 1..... _ 0.8 — 0. 6 e 1.25a. P 0 T o T 0.4 0.2 0.1 0 0 0.5 /, 1.0 1.5 FIGURE 3 -- - - - 132

THE UNIVERSITY OF MICHIGAN - 4134-2-F The whole picture is independent of the special values of P and a. Their values enter only in the scales. Of course the picture will change if the condition (15) is violated. It is convenient to express P in terms of the P in (9) which is the power per unit volume needed to maintain the steady state with a density ratio (1 + 3a) at the center. (see (11)). This P will be denoted by P t We have then P Pst (28) with some numerical constant 3. It is clear from the outset that a density ratio (1+ 3ac) will not be ultimately reached if 3 < 1. It is also clear that 3 1 will bring us ultimately down to this density ratio; but the neglect of heat conduction will not permit us to apply our results down to this density ratio. If, however, 3 >> 1, let us say B = 10, then we reach this density ratio (1 + 3a) during a time interval 0 < t < t which complies with the inequality (15). The ultimate density for t -> oo will, of course, be lower than (1+3) a) To prove this, we insert (28) with the value (9) of Pt into (22) to obtain st 2 1 a2 / = -- -.(29) 2a: k Then we have from the definition of t and from (24) O 3to/ 3a2 e =1+3a = 1 + (30) 20k C 133

THE UNIVERSITY OF MICHIGAN 4134-2-F or kt,_2 o- 1 kC k(1+ 3a (31) 2 3 2 2f3kU a a The right hand number is always less than (213) so for 3 = 10, (15) is very well satisfied. So far we have only proved that for 3 ~> 1 the neglect of heat conductivity is permitted up to t = t. Now we give an expression for the time t needed to 0 0 -1 -1 achieve a density ratio (1 +3a) =. From (29) and (30) we obtain 2 t = anT a (32) o 23 (y - 1)k Furthermore we evaluate the total power delivered up to t. It is, from (22) and (24) _ 3 4 3 4r 3 3 a Pt = - a 3p c T = a- c p T (33) 3 o 3 opot 3 poo 0 3-3 -3 7 or, with T = 300 K, p 10 gmcm, c ~0.24cal/gmgrad = 10 erg/gm~K 47rT 3 47r 3 6 -3 47ra3Pt = 4- a a-y. 3.10 ergcm. (34) 3 0 3 We obtain in particular for a = 100 cm and y = 10, 47a a Pt = 3.10 Joule (35) 3 o It is, of course, evident that/3 should not enter the equation (34). This quantity /3 served only to indicate the power needed in this problem to achieve a certain - - -- 134

I THE UNIVERSITY OF MICHIGAN 4134-2-F density ratio in relation to the power needed to maintain the steady state at this density ratio. As is seen, by comparison with (14), in this numerical example Pt is of the order of 600. 7 Using the equation of state of air and assuming a specific heat = - R/M 2 (M = average molecular weight of air), we can simplify (33) into a t = 4 a3 - 7/2 p L, (36) 3 o 3 which reveals clearly the dependence on pressure. Our assumption of a constant pressure over the whole space is, of course not exactly correct. One important consequence of this assumption is the steady velocity field which we have obtained in (21) while we would expect that a disturbance travels at a finite speed and does not affect the initial state at some distance until after some time. However, if the process of heating is relatively slow so that t is large compared to the time during which a disturbance travels over a o distance a, we can expect our results to be valid though (21) would no longer be correct at least outside the zone reached by the disturbance. But it is certainly a good approximation in the domain of interest, that is inside the sphere of radius V(t) which is less than 1.5 a in the numerical examples. It is possible to give a justification by integrating the equation of motion p d = - gradp (37) for the pressure by using (21). We restrict to r < a and find I. - - — ~ ~ 135

- THE UNIVERSITY OF MICHIGAN --- 4134-2- F 3t' r 2 p(r,t)= p(O,t) - p e-3t r (38) 27 p(0, t) does, according to our assumptions, not depend on t. Inserting now the equation of state p= R/M oTo (39) we obtain 2 P(r M r - p(r) M r. e-3t/t (0<r<a) (40) Po RTo 2~ Now 1.4RT /M is the square of the sound velocity v. The correction term in (40) is therefore 0. 7(r/Zvs) e -3 which is in fact small if the time sound needs to travel over a distance a is small compared to ZT. The preceding analysis is applied to find power required to heat a 1 meter sphere in the ionosphere between 60 and 100 Km so that the density at the center is 1/4 the ambient value. The results are shown in Fig. 4. The power required on the ground was computed assuming the heating is accomplished by high frequency power absorption and that the power given by (9) is increased by two factors. The first factor accounts for the small fraction of the power which will be absorbed in traveling a distance a (taken as 1 meter). For frequencies greater than 3 Mcs (greater than expected plasma and collision frequency), this factor is given approximately by: - -- -- 136 - --

THE 14 4 -I l1 UNIVERSITY OF MICHIGAN 4134-2-F P = absorbed power = 8acrauT a 00 I Co 4a) 4 —j a) 0) 03) "0 0( PI T, Cd 1.0 -0 T6 Ei 0 I q a = density factor = 1 a = sphere radius = 1 meter a = heat conductivity = 1. 51 x 10 sec-rg sec-cm-OE T = ambient temperature = 1660K h22 F = pattern factor = 422 w = microwave operating frequency = 30 Me D = antenna diameter = 10 meters F P =P - c = velocity of light G a FA A F = absorption factor = 2 awo v 2 cw w = plasma frequency =.28 Mc collision frequency = 3 M v = collision frequency = 3.3 Mc 1013 0 PLq I I 66 7b 8b 9b h, altitude (kilometers) FIGURE 4: POWER REQUIRED TO EFFECT A 750/ REDUCTION IN DENSITY AT THE CENTER OF A 1 METER ATMOSPHERIC SPHERE AS A FUNCTION OF ALTITUDE 100 137

i I THE UNIVERSITY OF MICHIGAN 4134-2-F 2 W V F a 2ka, b -a-P2 2 cw where k is the absorption coefficient, w the plasma frequency, v the collision P frequency and w the microwave frequency. The second factor accounts for the fact that at these high altitudes, even very large antennas with relatively narrow beams will illuminate a much larger area of the ionosphere than the desired ir square meters, and is given by h2X2 F - P ~4D2 where D is the antenna diameter, h the altitude, and X the wavelength. Using (9) with a = 1, values of ionospheric temperature and heat conductivity given in the 1959 ARDC model ionosphere 18], and the above factors, the ground based power requirements can be quickly estimated, and are given in Fig. 4. It is found that this power is relatively independent of frequency in the frequency range indicated. This power requirement is seen to be large even though the value given by (9) is not. w 138

THE UNIVERSITY OF MICHIGAN 4134-2-F REFERENCES 1. Ritt, R.K., IRE Trans. AP-4, 3, (July 1958). 2. Belyea, J.E., et al, "Studies in Radar Cross Sections XXXVIII - NonLinear Modeling of Maxwell's Equations", University of Michigan Radiation Laboratory Report No. 2871-4-T (December 1959). 3. Belyea, J.E., et al, "Studies in Radar Cross Sections XLV - Studies in Non-Linear Modeling II", University of Michigan Radiation Laboratory Report No. 2871-6-F (December 1960). 4. Tonks, L. and Langmuir, I., Phys. Rev. 34, 876 (1929). 5. Margenau, H., Phys. Rev. 69, 508, (1946). 6. Margenau, H., Phys. Rev. 73, 309, (1948). 7. Hartman, L.M., Phys. Rev. 73, 316, (1948). 8. Margenau, H., Phys. Rev. 109, 6, (1958). 9. Webster, A.G. Partial Differential Equations of Mathematical Physics, (New York: Dover Publishing Co.) p. 59. 10. Chapman, S., and Cowling, T. G., The Mathematical Theory of NonUniform Gases, (Cambridge: Cambridge University Press, 1939) p. 348 11. Cohen, R.S., et al, Phys. Rev. 80, (2), 230-238, (1950). 12. Holstein, T., Phys. Rev. 70, 367, (1946). 13. Dreicer,H., Phys. Rev. 17, (2), 343-354, (1960). 14. Spitzer, L. Jr., Physics of Fully Ionized Gases, (New York: Interscience Publishers, Inc., 1956) p. 83. 15. Jahnke and Emde, Tables of Functions, (New York: Dover Publishing Co.) p. 134. El 139

-- THE UNIVERSITY OF MICHIGAN 4134-2-F 16. 17. 18. Samaddar, S. N., "Wave Propagation in an Anisotropic Column with Ring Source Excitation", (To be published). Rosenbluth, N.M., et al, Phys. Rev. 107, 1, (1957). Minzner, R.A., et al, "The ARDC Model Ionosphere", Air Force Survey in Geophysics No. 115, AFCRC-TR-59-267, (August, 1959). - -- 140

UNIVERSITY OF MICHIGAN 3 9015 02829 6138 III II 3 9015 02829 6138