THE UNIVERSITY OF MICHIGAN INDUSTRY PROGRAM OF THE COLLEGE OF ENGINEERING THE DETERMINATION OF THE ELECTRON DENSITY PERTURBATIONS RESULTING FROM THE MIXING OF TWO DIFFERENT PLASMAS Fo V, Schultz November, 1962 IP-593

TABLE OF CONTENTS Page STATEMENT OF THE PROBLEM.................................. 1 BOUNDARY CONDITIONS...o........... o............. 7 CALCULATION OF ELECTRON DENSITY PERTURBATIONS.............. 9 CONCLUSIONS............................................. 22 ACKNOWLEDGEMENTS..............o.......o...... o.. o o 23 APPENDIX A - DETERMINATION OF BOUNDARY CONDITIONS AT INTERFACE................................. 24 APPENDIX B - EVALUATION OF CONTOUR INTEGRAL................ 32 REFERENCES...o.... o.........................o............ o 53 ii

I STATEMENT OF THE PROBLEM The problem under consideration is that of determining the mixing of two plasmas, originally of different densities. It is assumed that, for t(time) < 0, the half-space z c 0 (Region 1) is filled with a plasma having electron and ion densities each equal to n10 particles per cubic meter, and the half-space z > 0 (Region 2) is filled with a plasma having electron and ion densities each equal to n20 particles per cubic meter. Here n10 and n20 are constants, independent of x, y, z, and t. For t 3 0 the two plasmas are allowed to mix and the problem is to determine how the perturbations in electron densities, n11 in Region 1 and n21 in Region 2, vary with t and z. The time interval considered is assumed to be short enough that the ion densities in the two regions remain essentially unchanged. Other assumptions are these: (a) The mean free path is assumed to be very large so that collisions can be neglected. (b) The initial electron temperatures (T ) are equal in the two regions. (c) The perturbations in electron densities, nll and n21, the z-directed electron streaming velocities, u1 and u2, and the induced electric field intensities, E1 and E2, are small enough that second-order, and higher, products of these 1

terms can be neglected. This is also assumed to be true of the derivatives of these quantities. (d) There are no externally impressed electric, magnetic, or gravitational fields. (e) The plasmas obey the perfect gas law. The electron gas pressure is, rather arbitrarily but necessarily, assumed to be a scalar given by p = nkT. (f) The various physical quantities vary only with z and t, there being no variations in the x and y directions. (g) Pressure changes occur adiabatically, so that, for the electrons, p=kap =ky kT n where 7 is the ratio of specific heats and k is a unit vector in the z-direction. Our working equations for the electrons, then, are these, written in two dimensions, z and t. MKS units are used. m n a- + kT n + n e E =0 (momentum equation) (1) an a at az a + a (n u) = 0 (continuity equation) (2) at az aE =~e a = (n - n) (Poisson equation) (3) az Eo o p =n k T (equation of state) (4) 2

There are five unknowns, p, u, n, T and E, in these four equations, plus the equation under (g), above. Here: m = mass of electron n = number density of electrons n = number density of ions o p = pressure of electron gas -19 e =1. 602 x 109 coulomb E = electric field intensity (z-directed) u = z-directed streaming velocity of electrons -12 o = 8.854 x 10, permittivity of free space -23 k = 1. 380 x 10 joule per degree K T = temperature in degrees, Kelvin. Now assume that, in Region 1, for t > 0, the electron density can be written nl= n10 + nll(z, t), where nl 1 n10 and n10 is a constant. Similarly for Region 2: n2 = n20 + n21 (z, t). 3

Also, in the same way, P1 P10 + P11' P2 P20 + P21' The unchanging ion densities are n10 and n20 in the two regions. Making use of the assumptions mentioned above, one obtains the following forms of the four basic equations for Region 1, with a similar set for Region 2: au, an11 m n10 at + kTo Lz +nl0eE =0, (la) anl au 11 1 aE +n -- 0- (3a) az e nll (3a) 0 11 =n10 k T11 + nll k T10 (4a) To obtain a single partial differential equation in nll alone, one may differentiate (la) with respect to z, differentiate (2a) with respect to t, and use this differentiated form of (2a), and (3a), in (la). There results: 2 2 a2nl ykT T an a 1 m - -— ~ - + 2 nll =0, (5) at2 m az2 P 11 4

where 2 2 e n10 2 = (6) pl mn 0 is the electron plasma angular frequency squared in Region 1. In order to solve (5) for nl1(z, t) one may take the Laplace transform, obtaining 2 an YkT 2 at m 2 pi 11 t=+0 (7) Here N (z, s) is the Laplace transform of n (z, t). We are assuming that both an 1 n (z, +0) and al are zero so (7) reduces to 11at tt=+O a2N -2(s+ 1) N1l =0, (8) where = jYkT /m (9) is the adiabatic acoustic velocity in the electron gas. The solution of (8) is z2 22 z 22 — s +w. -s +w N (z, S)=F (s)e + F2(s)e, (10) 11 1 2 where Fl(s) and F2(s) are undetermined functions of s, only. By a completely similar process the following equation is obtained for Region 2: 5

-z 2 z s 2 N1 (z, s)F ) + (s) e F (11) 21' 3 4 (11) In order to take the inverse Laplace transforms of (10) and (11) for obtaining nll(z, t) and n21(z, t), one must find the undetermined functions F (s),... F,(s). This is done by applying boundary conditions. 6

II BOUNDARY CONDITIONS Boundary conditions on n (z, t) and n21(z, t) must be satisfied at z = 0 and at z =+ co. These last conditions are: lim n21(z, t) =0, (12) z -- oo lim n l(z,t) =0. (13) z ---- - 0o The boundary conditions at the interface (z = 0) are obtained by using the conservation equations for mass, momentums and energy, together with the following assumptions: (a) Viscous effects disappear because of the uniform unidirectional drift motion of the particles. (b) Second, and higher, order terms are neglected. (c) Heat conduction is negligible, because of the very low density of the gas. The resultant boundary conditions at z = 0, then, can be shown to be these: T=T2 (14) u =u2 (15) n1 n2 (16) 7

From (16) one obtains (for z = 0) n21 = nll+ (n0-n20 (17) which shows that our results will be valid only for values of (n10-n20) of the same order of magnitude as nll and n21.

III CALCULATION OF ELECTRON DENSITY PERTURBATIONS We are now in a position to apply the boundary conditions listed above and thus to determine F (s),... F4(s) of (10) and (11). It is necessary to use also the conservation of momentum equation, and (12) and (13), in the process. When this is done the following relations are obtained: F (s) =0 (18) (n -n S2+ 2 Fn20 10 s 2 (p2 F 2(s) = -, (19) n20 22 2 2 51s + +S + L n10 i p2 F3(s) = F2(s) + (n - n) (20) 3, 2 10 20s F4(s) =0. (21) We may now determine n l(z,t) and n21(z,t) by using the results of the last paragraph and then taking the inverse Laplace transforms of (10) and (11). We obtain T+ i3 1 1 i ra st nl(z,t) = l2im i- N (z, s) e ds, (22) 11and a similar euation for n 2(zt). and a similar equation for n 21(z,t). 9

By using the expression for F2(s) given in (19), we obtain z S 221 y+ i 2 2 c PC st 20-n10 ____ (z t)- hlirm /- - - -- - - - - ds. 12i --- oo 0f [ 2 2 + j 2n 20 s1 ---- +_ 4S +~p2- -- i n 0 pi s p2 (23) An eixaminaation of the integrand of this integral shows that there is a simple pole at the origin, branch points at + i w and + i w2, and only these. Branch p1 p2' cuts can be chosen as shown in Figure 1. It turns out that the Riemann surface of the integrand h; as four sheets. On the branch cuts between i op and i p2' and p p 21 between - iC) p and — i U 2- the first and second sheets are connected and the third and fourthrv sheets are connected. Between it 2 and - is 2 the first and third sheets p2 p2 are connlected and the second and fourth sheets are connected. One nay, as is scustomary, alter the path of integration of (23), as is indicated in Figure 1, in order to facilitate the integration process. In choosing the alte:rnatte path -we must make sure that, in traversing this path, we remain on the samte Ri:xn.emann surface as is used for Path 1 from TY - i 3 to 7+ i f, and either that no poles arel. ecilosed between the two paths, or that the values of the residues at the poles are properly taken into account. In order to make a proper choice of alternate path, let us rewrite (23) as follkows 10

D' 7+i 1 11 / (^"pi I \. / 16t1.13 2 / 16 112 1 3 8 Path l p2 7 - iw A' X Branch Points O Simple Pole.- Branch Cut - - - Paths of Integration FIGURE 1: S-PLANE FOR INTEGRAND OF EQUATION (23) 11

2 s(t+ Z 1+ s-2 +' 2 e ds n1(z, t) 20-n10 l im / p2 tl IzI/c (24) 11lim 2ni(z, t)=0. (26) Re s- -o -s -- 2s +u 2 s + W 2 T-i~ Is2 p 2 i_ remembering that z < 0. If we assume that (t + ). 0, which means that c t - Iz|/c (25) for the real part of s very large, we see that lim nll(z, t) - O. (26) Re s -oo Consequently, since no singular points of the integrand are located on Paths 1 or 3, or within the region bounded by these two paths, we see that lim ( / /- ) = lim / =0. (27) oO Ph ath 1 Path 3 1 Therefore, n (z, t)= 0 for t - Izl/co (28) For the case when (t + - ) > 0, or c t;lzl/c, (29) we have devised Path 2 of Figure 1, since in this case the integrand becomes infinite on Path 3, as Re s-co. Path 2 has been so chosen that one remains on one sheet by crossing no branch cuts. Also, no poles or other singularities are 12

enclosed between Path 1 and Path 2. Therefore, the value of the integral as determined along Path 2 is equal to the value determined by following Path 1. In order to carry out the integration it was found necessary to make some simplifying approximations in the integrand, which limit the validity of the results to values of z less than about 1. 5 Debye lengths. Then for |z 1<1.5 Debye lengths, t > zI/c the following approximate result was obtained: W (cosh z ) W n- (zt)=(n -n) [ p2 c P1 + 1 z Si(w t) 11 (n 20-n) Ln20 7 c p2 nO1 pi p2 1 z W 2 z sin w2t + -cos w t 2w 2 c p2 2 7r t c p2 p1 p 16 n2J 1 c12 1 -+ ( —n- 2u'Plp2 j - (cos W tt+cos wt) 2 %2 Pi P" 7p2 (n20-n10 )(plwp2) + 21 (sin wlt - sin wp2t) p2P pp2 Wp2- pl t (Eqn. continued on next page) 13

4n A Bt C1 p - ( - + — + ) Cos 7 t0 7r 2 (n 0-l)( +) t A B C D E +_2) cos w t+ ( +- ) sinpl t t5 t3 t p2 t4 t2 p1 D E + ( - + 2 ) sin p2t. (30) t4 t2 P Here Al, B1,.. E2 are functions of wp and t p2 It is obvious that the first four terms in the expression for n (z, t), given by (30) go to zero as (n20-n10 ) — 0, as they must. The remaining terms in this expression also go to zero as (n0-n" 0)-> O, since they arise from evaluating integrals between the limits of wp2 and op 1 and these limits are equal for n20= n10 The expression for nll(z,t) given by (30) was obtained by assuming that n10 n20. In order to complete the analysis it is necessary to either compute the perturbation n21(z, t) in Region 2 for n10 > n20 or to recalculate nll for n10L n20. The latter procedure was followed and the following result obtained (for n10 n20): 14

% cosh ( ~p p ) n (z,t) (n20-n 1) p2 c pl + Si (w 11 2010 n20 2 C - sto w n10 pi p2 1 z sin(w t)+ 1 z cos ( t) 2n 2 t2 c p1 27rt c pl pi 16 n o t + ( 20 p1p p2+ (cos w t + cos wp2t 7T W2 (n )(W +w2 Pi + w2 (n20 n10)( + to2 t (w - ) p2 1 2 "pi1 4n A B C1 -20 pl z_ ( + 1+ Cosw t 2 c Vt5 t3 t P2 r wp2 (n20- n10)(wpl + wp2) At B2 C? DI E + ( -- + - + ) cos t t+( —+ -- sin w t t t t p1 t4 t p +( + 2) sin w1t (31).4 2 P tB Here A, B',.... E' are functions of to and o. 15

It is now in order to investigate the physical significance of the two different expressions for nll(z, t): equation (30) for n10 > n20 and equation (31) for n10 < n2 This is done only for the case of n1o > n20. It will be remembered that n (z,t) = 0 for Izl/c > t. This means, of course, that the electron density perturbation n11 has a front which is propagating in the direction of decreasing z with a velocity c, the adiabatic acoustic velocity in the electron gas. This is characteristic of weak discontinuities, as discussed in Section 93 of j1 ]. We are, of course, dealing with a weak discontinuity since the linearization of the original partial differential equation, (1), is valid only for ill/nlO (n10-n20)/n10 L O. 1. Another restriction on the solutions is that z must not exceed about 1. 5 Debye lengths. The following are the values of A,..... E2 which appear in (30): A =- 24a /p 1 pi 2 B =2a(5- -2 )+4b pl C1 =0 24a + 6b 1 o pl 16

/2 a+b\ 2 E w = p p2 A2 24a/ 1 2 Pi 2 B- l1 - a + +b( p 2 + 2c L pl p pi C2Here a -90.1 D = ~ (24 a 2+ 6b) pi pi E2 -2 a( +b( p c( Pi Pi pl pl W Here a = -90.1 b =169.5 c =-79.3 a + b + c = 0 (exactly), these being numerics arising in a curve-fitting procedure used as an approximation in evaluating the integral of (23). If we now consider the following numerical values: 12 n10 =10 particles per cubic meter, n20 = 9 x 1012particles per cubic meter, T = 1000 degrees, Kelvin, 17

we find that -3 h = 2.18 x 10 meter, (Debye length) c = 1. 590 x 10 meters per second, (adiabatic acoustic velocity of the electron gas) f = 8. 98 megacycles per second, (plasma frequency of Region 1) f = 8. 52 megacycles per second, (plasma frequency of Region 2) 2 7 w0p = 5. 64 x 10 radians per second = 2r fl, pi 7 wt 2 = 5.35 x 10 radians per second = 27r f p2 2 p p = 0. 9486, p2 p1 (W/p ~ 0.900. If these values are used in (30), there results, after collecting terms: n1(z, t) 7 -11(Z- -- 0. 512 cosh 355z - 113 z Si(5. 35 x 10 t) 10- n20 +5. 03 x 1013 z +5.52 x0-21 1 2.62 -25 z sin(5. 35 x 7t) t2 t t 10 - +5.52 xt10 - +2.62 x 10 + - 0. 949 x 10-6 + 802 10-14 1 003x 1-19 z L t t2 t3 + 9. 74 x 1033 z cos(5. 35 x 107t)+8 65 x 10-13 -5. 52 x 10 2 1 t5 t 2 -25 14 1 -19 z - 2. 91 x 10 sin( 5.64 x lO t)i 0. 802 x 114 + 8.60 x 109 t4- t2 t3 - 9. 74 x 10~33 z- 1 cos(5. 64 x 107 t). (32) t5 t018 18

If one now attempts to use this expression to determine numerically the values of nll on the wavefront progressing in the direction of decreasing z, it is found that, because of the very small values of t involved, it is very difficult to obtain accurate results. Consequently the integral of (23) was re-evaluated under the assumption that t was small enough to allow us to replace sin rt by rt, and cos rt by unity. The resulting required integration is straightforward (after making other simplifying assumptions as before) and the final result is (for n0 > n20 and very small t): z 2 (Z, _ O W 2(cosh10 pl p Z t 36_ (n - c 6 2 (10- n20) n20 + (pl pi pp p2 pl pl p p2 ~0~(19 3r (n 10-n20) + 2 i (p +p2 p2 p1 ~4 n1pl p 1 Z 2 c (nlO -n20) (W1 +Wp2) p2 c 6 5 4 - 6~ )+ (1- (1 4 pl p Wpl 4 3 2 2 ~p2 b p2 c p2

If we now use (33) to calculate n11/n10 -n20), as a function of z for t = Iz |/c, we get information concerning the amplitude of the wavefront of the electron density perturbation as it moves in the direction of decreasing z. The results are as follows: TABLE I z (mm) n 1/(n -n20) 0. 0 -0.510 -0. 1 -0. 510 -0. 2 -0. 560 -0. 3 -0. 643 -0.4 -0. 759 -0.5 -0.908 -0.6 -1.091 -0 7 -1. 306 -0.8 -1. 555 -0. 9 -1. 837 -1.0 -2.153 Another result of physical interest is the variation of n 1/(n10-n20) as a function of z for values of t large enough that the oscillatory terms of the electron density perturbation have essentially damped out, but t small enough that the ions have not yet started to move. This result is easily calculated by using only the first two terms of (32) and one obtains the following values: 20

TABLE II z (mm) nll/(n10-n20) 0.0 -0.512 -0. 1 -0.494 -0.2 -0.477 -0.3 -0.461 -0.4 -0.445 -0.5 -0.431 -0.6 -0.418 -0.7 -0.404 -0.8 -0.391 =0.9 -0. 379 -1.0 -0. 368 -1. 2 -0. 347 -1.4 -0. 328 -1. 6 -0.313 1. 8 -0. 301 -2.0 -0. 293 It also is of interest to obtain an indication of the rapidity with which the electron density oscillations, as shown by (32), damp out. The t term will damp out the most slowly. Our results are valid only out to | zl 2 millimeters and it takes the perturbation about 10 second to reach this distance. Consequently, -7 -1 at the end of 10 second this t term will damp out to about 10 percent of its original value. Since the frequency of oscillation is about 8. 52 megacycles, only a cycle, or so, of this oscillation would occur. The Si(5. 35 x10 t) term is also a damped oscillating term, and a plot of it shows that several cycles of this exist. It is considered unnecessary to make a detailed analysis of the behavior of n21 since the equation for it is very similar to that for n 1, so very similar phenomena will occur. 21

IV CONCLUSIONS From (32), and Tables I and II, one can conclude that as t increases from t = 0 the following phenomena occur: 1. Immediately the electron density at z = 0 assumes a value approximately midway between n10 and n20. 2. A rarefaction wave of electron density moves in the direction of decreasing z, with a velocity c and with an increasing amplitude. 3. At each point of the medium, after the passage of the rarefaction wave, the electron density oscillates at two different frequencies, one corresponding to w and the other to w These are damped pi p2 oscillations consisting of only a few cycles. 4. After the oscillations die out the electron density varies smoothly from a value of about 0.5 (n10-n20) at z = 0 to a value of n10 for large negative values of z. 22

V A CKNOWLEDGEMENTS The writer wishes to acknowledge the helpfulness of discussions with several members of the staff of the Radiation Laboratory, particularly Dr. R. J. Leite and Dr. R. K. Ritt, and with Dr. R. A. Holmes of Purdue University. 23

APPENDIX A DETERMINATION OF BOUNDARY CONDITIONS AT INTERFACE The boundary conditions at the interface, z = 0, are determined by first using the applicable form of the momentum equation, 3u au q At = -'-k- kn 3-a (nV.Vk) + -q E (A-l) at uk Tt ^ n1. c (nik m i 3' k axk k ak " where V and V are components of the peculiar, or thermal, velocity of a par1 k ticle, and q is its charge, plus the continuity equation, - = - (pu) (A-2) t -k aXk where p = mn, the mass density of the fluid. To obtain useful results from (A-l) and (A-2), one uses these equations in connection with the equation 3u. 3p (pu) = p -t + u. (A-3) ti 1 T' By substituting (A-l) and (A-2) for the right-hand terms of (A-3) we find that 3__ (pu.) = -p EUk' u P 3 ---- it 1 k a uk n i (nV i Vk +qp E. u... (A-4' n o k aX In our case of uniform unidirectional drift motion of the particles, viscous 24

effects disappear, with the result that Vik = k Vvk (A-5) I k ik i k' where 6ik is the Kronecker delta. Then - (nVV ) (nV2) == _ - P (A-6) k axk 1 k ax. 1 x m a m ik a x By using this result in (A-4) we obtain a' a (pui): ( pu.) Z~Puu. uk 89]+nE. q(A-7) a (PUi) 1 = P axk i axk ik 1 a Equation (A-7) can be made more compact and its physical significance made more evident by using this result: au i ( puk) ^ — (P U. P= u + U a (p uk) = k x u a xk Then (A-7) becomes (pu.) = - (Pu.+ i p) + qn E. (A-8) To get (A-8) into its final useful form, one expresses the electric volume force in terms of the Maxwell stress tensor, by pages 95-97 of [2]: a Tik qn Ei= Z axk (A-9) k k 25

In our case, the Maxwell stress tensor is 1 21 T E E 1. E2 (A-10) ik o i k 2 ik since we are dealing with a plasma in otherwise empty space, with a constant electric permittivity c By using (A-9) and (A-10) in (A-8) we obtain a pu.) = 7, FE (A-li) a1t k i E Lx [0 P k + ik i k ik 2 (A-l1) k Xk By pages 12-15 of 1, we can write T 1ik P U., ^-i~ kk (pu.) = - a at 1 i Xk where Tik is the momentum flux density tensor, given in our case by i P u k s- i o E Ek + 6 E2 (A-12) ik 1kPPu k oik ik 2 Now 7Ik n kis the flux of the i component of momentum through a unit surface area having the unit vector n along its outward normal. The flux of each component of momentum must be continuous at an interface between a Region 1 and a Region 2, so E (c ikl nk = (Tk nk' i = 1,2,3. (A-13) k k 26

At the interface z = 0, we have k = 3 in (A-13) with the following three resulting scalar equations for i = 1, 2, 3, respectively. p u u - c E E = p u u -e E E (A-14) 1 x lz o ElxElz P22x 2z o 2x 2z p u u - cE E = p2u u - EyE E(A-15) 1 lyulz o ly 12 z 2y z A-15) Pl+ Plz z E1 2 + 2u2z - E2 + E2 (A-16) Pi l 1Z 0 z 2 1 2 2 2z 2 o 2z 2 2' Because our regions are unbounded and homogeneous in the x and y directions, the x and y components of u and E are zero in both Region 1 and Region 2. Also, we know that E1 = E2z Therefore, only (A-16) gives a useful result: p+p u2 =p2+pu 2 P1 +1 Ulz 2 2 2z or p + 1) + P1+P11) = (P20 +P2 )+(P2 +P2 )u 2 P10 p11 (P 112 1z 0 1 20 21 2z Dropping terms of second, and higher, order gives PP1 0 +p1 P = P2 = P20 +p21 (A-17) A second boundary condition at z = 0 is given by the fact that the normal flux of electrons must be continuous at the surface z = 0. That is, u1 n1 = 2 n2. (A-18) 27

W\e obtain a third boundary condition here by recognizing that the normal conmponent of energy flux must be continuous at this plane interface between the two regions. By sections 51 and 53 of [3] we have for the conservation of energy, in our case, _ 1 1 3( -1 pu2 + Pe + E2 + H2) at 2 2 o 2 o -V -p( 2 u2 +w) -u - T' - wVT+ EXH - -V.. (A-19' Here: e = internal energy per unit mass, w = enthalpy per unit mass = + p/ p, a' = viscosity stress tensor, 6 = thermal conductivity. As mentioned previously, viscous effects do not enter in our case, so the term u * o' in (A-19) is dropped. Also, we are assuming very long mean free paths so collisions can be neglected, allowing us to drop the heat conduction term XVT. Because we are dealing with longitudinal oscillations, there is no magnetic field and the Poynting vector term disappears. By Section 85 of [1 we have for the enthalpy per unit mass of a perfect gas: = pV = e + - (A-20Y 7-1 y-1 p P 28

whence 1 P (A-21) y-1 p Now if these expressions are used in (A-19), together with the results of the last paragraph above, we obtain 1 P 2 + P2 mnu1 (- u2 + e-1 n ) = mn u ( 2 -1 mn2- (A-22) 112 1 y-1 mn 2 2 2 2 y-1 mn2 By using the facts that P= P2 and u n = u2n2 given by (A-17) and (A-18), respectively, (A-22) above becomes 1 u2 + Y 1 = 1 u2+ (A-23) Ul - (A-23) 2 1 y-l mnn1 2 2'y-l mn2 1 2 and 1 2 in(A-23) weobtain By dropping the second-order terms, u2and u2in (A-23), we obtain 2 1 2 2' n = n2 (for z = 0) (A-24) By using (A-24) in (A-17) and (A-18), we get T1 = T2, (z = ) (A-25) 29

u1 = u2, (z=0). (A-26) From (A-24) we have n 10 + n 2 n + 2 n10o +nll n20+ 1 n21 = nll+(n10 -n20) (Z0) (A-27) This equation tells us that our results will be valid only for values of (n10 - n20) of the same order of magnitude as n1 and n21. Taking the Laplace transform of (A-27) gives nl0 - n20 N (0, s) = Nll (0,s) + 20 (A-28) 21 11 s A second relation involving N 1(0, s) and N21(0, s) can be obtained by taking the Laplace transform of (la), obtaining 8N (z, s) m n10 sU (z,s)-ul(z,0)J + n kT ) = 01 (A-29) where U(z, s) and f(z,s) are the transforms of u (z,t) and El(z,t), respectively. One of our initial conditions is that ul(z 0) = u 2(z,) = 0 Consequently (A-29) can be rewritten ~e y~Tk T 3N1 (z s) U (z,s) = -- (zs) - (A-30) 1 ms 1 mn 0s a( 30

with a similar equation for U2(z, s). Because of (A-25) and (A-26), we get from (A-30): a N21 ( 0, s) n20 aNll(11s ) E) 8~z ~ ~ n^ 3z It-~ (A-3 1) n10 where we have made use of the fact that El(0,t) = E(0,t). We can now use (12), (13), (A-28), and (A-31) to determine F (s),... F (s) in (10) and(11). From (13), lim N (z,s) = 0. (A-32) z -- - oo Now F ( s) and F (s) are not functions of z, so (A-32), in conjunction with (10), shows that it is necessary that F (s) = 0. Similarly, it can be shown that F (s) = 0. Lastly, wecandetermine F2(s) and F3(s) by 4 2 3 using (A-28) and (A-31). This is routine algebra and the results are: (n20 - n ) s2 + p2 F (s) - S20 —10 v p2, (A-33) 2 i20 2 2 2 2 s20- s + + +W nn P1 -n F3(s) = F2(s) + 1 0 (A-34) 3 231 s 31

APPENDIX B EVALUATION OF CONTOUR INTEGRAL In order to evaluate the integral in (23) we use Path 2 of Figure 1 and rewrite the integral thus n20 - n10 all(z, t)= -- - 2_i (011012)i n (zt) i (021 +22) c e ePllP1 2 eit i) / P21P22e e reit d(reio) 0im. 4 + \ Be Pao 2 (rei0) nio 1 e2 Pll 12)+ IP21P22 e2 Path 2 where (B-l) in11 s- iw P1le i021 s-idp2 =21e s+ipl =P12e i022 s + iw2 =P22e s = re, and -7/2 0 mn 3 (3/2) r It is to be noted that y > 0, but y can be a small quantity. Figure B-1 illustrates the geometry. The value of the contour integral evaluated on Path 2 is, of course, equal to the sum of the integrals evaluated along the various separate paths which make up Path 2, as shown in Figure 1. We now evaluate our integral along the various separate paths. 32

S 1 wpl P2 12 -i -) -l~pli1 FIGURE B-1i POLAR REPRESENTATION OF s + i 33 33

First, we note that our integrand is analytic along the negative axis of reals, so the integrals along paths B-O and O-C are equal in magnitude and opposite in sign. Also, then, the integrals from 1 to 2, and from A' to D' can be evaluated along continuous paths between these points. Next, let us consider circular arc A-A'. On this path r = r and r — oo, O O p - oo, 0- -7r/2, and 0-m -Tr/2 as j —- oo. Then.Z2 211+ 12) -i/2 p 21p22e e e2 e id2 lim IA'- lim /: A-Alo i( 3 — *o [-~ODn +,3-po 1-oo / rn_2_ (011+012)2 00-r/2 J L 11 12 O n02122e J o=0o (B-2) i ty+itB lim I, =lim e d =O OaD A-A' (n20/n10)+ 1 -- ao J 0 --— /2 0=0 Similarly, lim ID-D 0 (B-3) 34

First, we note that our integrand is analytic along the negative axis of reals, so the integrals along paths B-O and O-C are equal in magnitude and opposite in sign. Also, then, the integrals from 1 to 2, and from A' to D' can be evaluated along continuous paths between these points. Next, let us consider circular arc A-A'. On this path r = r and r — oo, O O p -— oo, 0 — -7r/2, and -mn -Tr/2as 3-mco. Then -7r/2 1 21P22ee e id lim I A-A' lim i.3 — *.oD i - 0.........:3 —00A^ oo /f F _2_ (011o012) - 21 2 - /2 n10 P11P12 o P21P22e no 0 (B-2) —,_/ 2 ~i'~~-/ t7+it3 lim IAA- lim f e d0 = O-DAo (n20/nO)+l 1 -0 0 —7-r/2 0=0 Similarly, lim ID'-D =0 (B-3) ] —-3-oo 34

On the circular arc, A' - D', we use the form of the integral given by (24), remembering that (t + c ) > 0. In the second and third quadrants (Re s) <0, c and as 3 — oo, (Re s) — - co, so lim IA'-D =0 (B-4) /3-> a) We next consider the integrals over the paths 1-2 and 9-10. Here r= r, a constant during the integration, and Pll = Wpl lim P21 p2!im r —,-0 P n 22 p2 P12'pl On the path 1-2: 011 = - lim 21 2 0o ~ a 022 2 1 0 = 17r 012 2 For the path 9-10, 011 - 2 r 1 0 = 2 7r ~21 2 lim 2 r -r 0 j - - Ir O 022- 2 7r iT 012= 2 i35 35

Then Z tr ei 3 c Wple o 2 i ip2e I12 = lim (B-5) r — 3O r 1 20 - + ~2 e ~0 pi p2 z c pi iwp 2e -/ Er- roteil0.I -lim e do ) 1-2n r -- -O 20+w +w 0=n. p1 p2 E —0 2 — r 1i 2e — 2 n n20 - lo i + %2 1 10 pi p2 For the path 9-10, zw rotei0 c pi 9-10 = limr ^0 pl "iOpl e p2 9-1 r — "0 n.20] 0= I -w+E 4W + E 02 Ln 1n10 p 2 (B-7) z topl p2 2*^ ~ ~~E rot ei' *Io3~~ ~ —-- lim e d, 9-10 n20 r - 0 f to + w 0 0 o = - - n10 p1 p2 e —0 2 36

zw c pl i 7r p2 e I91 = - (B-8) 9-10'n20 nlo pl p2 27ri w cosh - w p2 c pi (B-9) I,1 1-2 + -10 (B-9) 20 -- 1 _ + W n10 pl p2 fort > zl/c. 7r 3ad For the circular path 3-4 we let p 2- 0 and 2+ e 0 22 2 7r - e, and then let -- 0. On this path as p2- 0: P11 =W1 +Wp2 011 2 21 = 2 p221 2 = P12 p1 p2 012 =2 s iw + P e'022 p2 + P22 i022 ds =i P2e dO22 22 Then lim I = lim P ~ P i22~ 3 7r+ 0 2 2 2 r(-i p2 P22e )t 22. ~-e E 2;e e eip e r'2, O-p2P22.P22e f'r-'-'_ i "3' ) d322 (i +e 20 2 2 (~ 1+ 2 (-iwp2+P i~22 n~ i2 32 l -p2 22e (B-10)..... —------- 37

lim I3_40. (B-11) P22'0 Similarly, lim I7-8 = 0, (B-12) p22 0 lim I11 12 = (B-13) lim I1516= (B-14) P21-* 0 On the path 5-5', we let p 0 and e 012-. Here, as p - 12 2 12 0n+ 12 l'211 Pi 011 2 _ 3 P21 p2 + p1 021 - 7r P22'pl P2 022 2 i,12 s=-io ++pl2e i0 12 ds =ip12e d 12 38

Then we obtain lim I = lim 12 ~ 12 -— 0 (- +012) i012 * I i. d0 2 2 / (-i +p +P12e )t --' (-i P12 c w p12le2 12 e p12e 2 (B-i5) lim I5_5, = 0. (B-16) Similarly lim 2 =2, (B-17) lim I =0. (B-18) 1I3-14 1114 P 0 We now consider the following four integrals as a group, since they combine well. 0 z 2 2 f 2 c p1- irt / -2r e e e dr 16-1 20 2 2 2 2(Br= r2 n U-1p r + -p2 - r Wp2 "10 ^ p1P2]39 39

z 2 2 _ - Iw - r p2 2 2 e c Pl -irt d U | 2 2 c pl - -irt r2-3 wp 2 (B-20) O/ [2"2op 2 r2 2 2 r tow r + w r r 0 "^Kl^n P' P2 Lo w I to2 -r e dr p1 p2 2 J p2 1~8-9 lp t r2e eirt dr I ir20 2 2 (2 2 ) r = r - w - r w -r 4 pi 2 p wr2 r Z j 2 2i "P 2 2 -c -p irt p p- r e P e dr I10 ( -24) L 20 2 2 2 2 r =0 r t- — r +tI- r r = 0 = 0 p1 i p2 n1040

We can simplify this integral by considering the first term in the denominator n0 20 2 20 2 - w -r = - - r n10 pl 10 p2 n10 Following (17), we noted that our results will be valid only for (n10 -n20 ) of the same order of magnitude as nll and n21, and we have restricted our work to the case of n!!' n10 and nP n2. Therefore we must restrict n20/n10 as follows: n20/n10 0.9, n / ^ 0. 95, n20/n10 ~0.95, and we can say that A r -2 A 2 2 2 n/n. W - r = c -r p2-n10 n p2 20/ 10 | p2 n p210 Then we can write 2 P2 (sinh - 1r (sin rt) dr I 2i. (B-25) r=O This integral still is untractable so we must look for further simplifications. Now 3 5 sinh x =x+ + +..,. (B-26) 3! 5: We can use only the first term with fair accuracy if 3 x410x 41

or if x,_ 0.774. (B-27) Sinh - W - r has its maximum value, sinh -w, when r = 0. Therefore, to C I pl pi represent sinh z J by only the first term in the corresponding series we must c ppl satisfy the relation z z- W 0.774. c pl By using the numerical values for w1 and c, which are given just before equation p1 (32), we find that our approximation for sinh z W1 will be valid for c P1 z 3 millimeters, (B-28) which is of the order of a Debye length. We then assume that z | 2 2 2 z 2 2 sinh 1- r 2p1- -r (B-29) Our last simplifying assumption is this, 2) 2 z u 1 -r- (1 2 ) (B-30) pl In our integral of (B-25) the upper limit on r is wp2 = 0. 948 w. For this value of p2 pl r we have: 22 2 2 1- (r2/ )= - ( 2/w ) 0. 316, 2 (1 - ) =0o550. 2 2 pl 42

Our approximation (B-30) is thus seen to be only fair at this limit, however, it is much better for smaller values of r. By using the above-discussed approximation in (B-25) we get p2t W2t P I 22i p d(rt) - (rt)sin(rt)d(rt), 12 C PI d(rt)-i t rt=0 P t0 rt=O (B-31) cY Z 1 I z 2 p1 (B-32) W t Here S i( t) / sin x dx (B-33) 0 In order to obtain some idea of the validity of (B-32), we have numerically evaluated I2/2i, as given by (B-25), for the following values of the various parameters: meters: ~z = 1. 0 millimeter, 5 c = 2 x 10 meters per second 7 o = 5 x 10 radians per second, pl 7 = 4.74 x 10 radians per second, p2 and for w 2t taking the values 0. 5, 1.0 1.5,.... 20. 0 radians. Then I2/2i was calculated by using (B-32). The results of the two different calculations are shown 43

in Figure B-2. The numerical integration of I2 was carried out by using Simpson's rule, with the number of subdivisions increased at w 2t = 5, 10, and 15. This p2 accounts for the jumps on the curve of yl at these points. Our last integral to evaluate is 3 45+ I6-7+ I12-13+ 114-15 (B-34) where p i r 2 c P-r e-irt dr ----, (B-35) 4-5 n/ = J r 20 12 2 2 2 r -w=u r- w -r +ir -w0 =p2 n10 p1 - r 2 w, w2z -r p2 o p2 P2. 2 c p -irt dr ir-w2 e e I_7 = / --- p (B-36) pi 2 e c pi irt 213( i rwP2 e2 ee dr I12- = / ---------------— _____ ____________ _____(B -37) n20 2 2 2 2 r =u r CO 10 w -r +i r - op2 p2 rLn 2] ir

0.5 - -- - -- 0.4 0.3 C 0.3 - 0.2 I 0.1 0 5 10 15 20 w t (radians) P2 P (sinh. - r )(sin rt) d(rt) FIGURE B-2: PLOTS OF y / pi( 1 } (rt) rt=0 z 1 z ~p2 z AND wY2: Si( t)- - sin(wp2t)+ cos(w 2t), AND y pi p2 1 p 22 c p2 2w t c p2 pl -3 5 7 7 FOR z=10, c=2 x10, w =5x10, w =4.74 x10 1 45 45

z 2 2 22- ~ 2 pl-r irt p2 i r-p2 e e dr p2I~~~~~~ =/ ------— (B-38) I14-15 = (B-38) / r 2 20 2'2 r2 2 r - n l-r +i jr - pr=p l 1 p2 After considerable algebraic manipulation, we find that pi 2 2 2 2 z2 2 4 p (r -W )(w -r )(cosh - r )(cos rt) dr 4i20 p p c pl T= p2 p1 - 3 0n 0 j n/ 0 n20 2 2 gp r nln) r-t~p2 % 2r2 z 2 /pi n (r w 2)(sinh - w - r )(sin rt) dr 4i n10 /p2 c P + (B-39) n10-n20 J r n+n 2 (B-3 r top2 r K1 0n2) 20 2 2 P2 _r n ~r-w p2 In order to handle these integrals we make use of the fact that our range of integration is very small, since o 2 =0.948 o. p2 pi First, we consider the denominator of the two integrals. 1/r is nearly constant throughout the range of integration, so we replace it by its value at the midpoint of the integration range: 1 2 1- ^ -- = ---' (B-40) r (W o- 2) /2 p +Wp2 46

Next we consider the other factor in the denominator, at the end points of the range of integration nlO+ n20 2 2 20 2 r= ( )r - - o p2 nlo p2 n10 p2 n0 + n20 2 2 n10 2 r=W,: ( )r — P~1( nlo ) p2 n20 p2 2 Since n20 = 0.9 n10, the above factor is very close to w 2 throughout the range of 2O p2 2 integration, and will be taken equal to w 2 p2 The first of the integrals in (B-39) can then be written: pl T4..i n20 2 1 2 2 2 z 2 2:=n0 —20 l+p2 / (r -w 2)(W -r )(cosh W-r (cos os rt)dr. 3 "n0-n20 WJpl+JP2 L2 P2 pi c pi p2 r = p2 (B-41 In order to handle this integral one must resort to further simplifications. The cosh term will be equal to unity at the upper limit of integration and equal to z 2 2 (cosh - jiw - W2 ) at the lower limit. At this lower limit, the cosh term will be c p1 p2 equal to 1. 1 if c Wl 2 = 0.445 (B-42) By using the previously listed numerical values of c, W 1, and w 2, one finds that the relation of (B-42) will be satisfied if z = 5. 56 millimeters. Therefore, one is justified in replacing cosh z 1 - r2 by unity, if z /5. 56 millimeters. 47

Last, we consider the factor 2 (r - 2 )(W - r in the integrand of * p2 p1 (B-41). To handle this radical we do the following: 2 2 r-w r-w r+2w ( 2w r-w p2 J p2 p2 | p2 p2 2 22w -r -- -r +r + 2 w -r p p 1 P1 pi p so 22 2 2 (r - )(w r ) 2pw (r-w )( -r) (B-43) \ p2 pi pl p2 p2 pi Now the function J (r-w 2)( pl-r) is zero at each limit of integration, has a maximum at r = (w 1+wp2)/2, and is symmetrical about this maximum. We replace the radical by a parabola which passes through these three points. The parabolic function is 2 y r - p + ( 2 P2 (B-44) -w 2 2 To check the validity of this substitution, a normalized numerical check was carried out to compare y1, = (r-l)(1. 0548 - r) with 2 = -36. 5 (r-1. 0274) + 0. 0274 48

as r varies from 1.00 to 1. 0274. The results are shown below: r Y1 Y2 1.000 0.0000 0. 0000 1.005 0.0158 0.0091 1.010 0.0212 0.0163 1.015 0.0244 0.0218 1.020 0.0264 0.0254 1. 025 0.0273 0.0272 1. 0274 0. 0274 0. 0274 These values of y1 and Y2 check quite well. One feels fairly well justified, then, in replacing the radical by (B-44). We finally arrive at pi 8in20 2 2 - 1+ 2 )2+(L 2 n2 )(r- -t 2 co22(nlo-n20)(W pl -P2 P2 p2 cos rt dr. (B-45) The integration is straightforward now, and we obtain 32i n w w r._____ 1'_= 22 p-[p2 2 p1p)(cos pt+ cos w t) + 3 (sinwlt-sinw t 3 p2 pl WP2 (n 10 )(W W p4pl p (B-46) It is to be noted that I3 = 0 if 2 =W 49

By using the simplifications discussed between (B-39) and (B-41), we can write the second integral of (B-39) thus: Wpl u 4i n10 z2 f2 2 2 2 4Ill i n1- 2- - - (r- u )(sinh - - r )(sin rt) dr. 3 n10-n20 W +W 2 p2 c p p p p2 r = p2 (B-47) By using the same reasoning as was used in obtaining (B-29), we find that, in the integrand of (B-47), we can use this approximation, z 2 2 z 2 2 sinh C w - r - r if z is not greater than one centimeter. With this simplification the integral still is intractable because of the factor 2 2 2 - r =) 1 - (r/w ).The series approximation of (B-30) is not satisi p l p1 factory here, because, throughout the range of integration in this case, r/Wpl is close to unity. Because of this, the radical was approximated by the following polynomial: - (r/W)2 =a(r/wp)2 + b(r/w )+ c, (B-48) where the constants of the polynomial were determined by matching the two functions near the end points and center point of the range of (r/w ). Following are the values of a, b, and c for our particular range of values of r/pl: 50

a =-90.13 b =169.46 c =-79.32 The goodness of the fit of the parabola for these values of a, b and c is shown by the following table: r/wpl - (rp a(r/l) 2+b(r/W )+ 0. 95 0.3122 0. 3150 0. 96 0.2800 0.2881 0.97 0.2431 0.2431 0.98 0.1990 0.1801 0.99 0.1411 0.0991 1.00 0.0000 0.0000 With these simplifications, the integral of (B-47) becomes: 8 i n0W 1 0 IT -= _ 10 lpi __ z 3 W2 p2( nl0 -n20)(pl wP2) pl (r2-w2) a(r/W 2 + b(r/ 1) + c (sin rt) dr r =w p2 (B-49) 51

This integration is straightforward, giving ^ r 8 0p za 1 4 4 4,3 I3' -T2. - pl z _a___ ( 2_ t,_2cos p2t-w lcoswplt)- (2sin wp2t 3 12 2 24 -p2 (nlO-n20)(wpl+wp2) tw - wpi sin wpl - 1i (w2 cos w2t- WOpl cosW pt) + 4 (wp2sinwp2t - sin w plt) + (Cos W t- Cosw t) +w [ (W3 Cosw t-w3 cosw t) 3 2i5 p2 pi p2t p2 pi pi +c- p p) tp2 2 PiPt P 32. 2 6 6 - 2 sin w 2t- i w 1 -c t- (w cos wp2t-) cos wplt)+(sin tsinlt) t ttt5 2 ~2 t3 (cos p2t-cos it) -b ( p2Cos pt- plcoslt)3 p2 P1 2 p2 P1 p1 12~ ~2 p1 - 1 ( si n t t-sin l t -c 2os t - cos (B-0) ( p1 jPiS pi2 Lt p2 pt jj(B-50) t2

REFERENCES 1. Landau, L.D., and E. M. Lifshitz, Fluid Mechanics, Addison-Wesley, 1959. 2. Panofsky, W. K. H., and M. Phillips, Classical Electricity and Magnetism, Addison-Wesley, 1955. 3. Landau, L.D. and E. M. Lifshitz, Electrodynamics of Continuous Media, Addison-Wesley, 1960. 53