Computer Experiments on Low Density Crossed-field Electron Beams

R H Levy and R W Hockney

Avco Everett Research Lab

June 1967

Research Report 273

ABSTRACT

The motions of thin and thick low density (q = ωp2c2 ≪1) crossed-field electron beams have been studied on the computer. Theoretical indications that such beams can be stable or unstable depending on the geometrical and electrical arrangements are supported, and the large amplitude behavior of a typical unstable configuration is found. In particular, thick, re -entrant beams can be quite stable against perturbations of the diocotron type.

I. INTRODUCTION

Crossed-field electron beams (i.e., electrons traversing crossed electric and magnetic fields} are characterized by the non-dimensional ratio

q = ωp 2 c 2 − nm/ε 0 B 2 (1.1)

It has been shown [1] that for thick, re-entrant beams having low values of q(<1/20), a substantial degree of stability is possible; that this stability can be achieved in practice is strongly implied by empirical considerations related to the success of the Vac-Ion Pump [2]. Two new applications for low density crossed-field electron beams, [3] [4] [5] have stimulated interest in this stability, and have given rise to new theoretical studies [1] [6] [7][8][9], experiments [10] [11] and to the computer experiments which are the subject of this paper.

For the calculations reported in this paper, we have chosen initial conditions susceptible of analysis, at least in the linear regime. The agreement between the computer calculation and the analysis in the linear regime gives us confidence that more complicated initial conditions can be chosen, and their development in time followed on the computer. In addition, our computer results in the non-linear regime are quite informative. The results of the computer calculations support the view that thick, re-entrant low density crossed-field beams are stable; in a case for which the linear theory predicts instability, growth is indeed observed, but the mode saturates without loss of electrons to the walls.

II. THE MODEL

The model used in our study is exhibited in Fig. 1.

Figure 1 - Schematic diagram of the model chosen for the computer calculation. The model is two-dimensional (no variations in the z-direction). It is also, as shown, periodic in the x-direction with maximum wavelength equal to the distance in the y-direction between the conducting electrode boundaries. The computer calculation is thus restricted to a square of side 2b. In the frame shown (center of the electron beam at rest), both electrodes are grounded, and the potential represents the mutual repulsion of the electrons. Note the shear in the electron drift velocity.

A beam of electrons of uniform number density n and width 2a is located symmetrically between two conducting plates separated by a distance 2b(> 2a). A uniform magnetic field of strength B points in the z-direction. The model is two-dimensional, no perturbations in the z-direction being considered. In view of this assumption, any potential applied between the conductors can be removed by transferring to a frame moving in the x-direction with a suitable velocity. Such a transformation can have no effect on the presence or absence of stability, or on any motion of the electrons in the y-direction. We can therefore consider both conducting plates to be at ground potential; the potential depression caused by the electrons in the unperturbed state is given by:

φ = − ne 2 ε 0 x 2a(b −y) 2ab − a 2 −y 2 2a(b+y) a ≤ y ≤ b −a ≤ y ≤ a −b ≤ y ≤ −a (2.1)

This potential depression represents the mutual repulsion of the electrons; the basic tendency of the electron beam is to spread out in the y-direction, with the electrons arriving eventually at the conducting boundaries.

In the steady state, the divergence of the electric field in the beam gives rise to a shear in the drift velocity. The magnitude of this shear is

ω 0 = d dy (E/B) = ne ε 0 B = ω p 2 c = qω c (2.2)

The shear is therefore independent of the electron mass, and much smaller than the gyro frequency. It is assumed that the quasi-static approximation is valid, that is, that E/B ≪ c.

We can form a characteristic gyro-radius by considering the electric field at the edge of the beam E = nea/ε0. The gyro-radius at the drift speed appropriate to this electric field is:

m(E/B) eB = mn ε 0 B 2 a = qa (2.3)

Thus for q ≪ 1 , and for perpendicular electron thermal velocities v1 on the order of E/B , the characteristic gyro-radius is much smaller than the beam width.

As stated in the Introduction, this model is susceptible of analysis in the linear regime. Thus, before discussing the computer calculations, we briefly describe the linear analysis, thereby indicating what is to be expected from the computer.

III. THEORETICAL CONSIDERATIONS

The theoretical treatment of the small amplitude motions of crossed-field electron beams that follows makes use of an Eulerian treatment. This treatment is justified when, as in this case, we are studying wavelengths much greater than the characteristic gyro radius.

In the limit of small q (small electron mass, low density, high magnetic field) the equations of motion of the electrons can be approximated by the equation

E ~ + v ~ x B ~ = − ∇ φ + v ~ x B ~ = 0 (3.1)

This is essentially the guiding center approximation. Since the three vectors in (3.1) are mutually perpendicular,

v ~ = − ∇ φ x B/B 2 (3.2)

Notice that the velocity vector of equation (3.2) is solenoidal. The divergence of equation (3.1) yields

n = ε0 B e [ curl v ~ ]z = ε0 e 2 φ (3.3)

The equation of continuity

Dn Dt = ε0 B e D Dt [ curl v ]z = 0 (3.4)

completes the system of equations. At the conducting boundaries

∂φ ∂x = vB = 0 (3.5)

Since v~ is solenoidal, it can be considered to represent the flow of an incompressible inviscid fluid. Similarly, Eq. (3.3) shows that the electron density is analogous to the vorticity. Eq. (3.4) is then the conservation of vorticity, and Eq. (3.5) the condition of no normal flow at the boundaries. Hence a complete analogy exists between the low density-crossed field electron beam, and the sheared flow of an incompressible inviscid fluid. The stability of such beams is then analogous to the Kelvin-Helmholtz instability of sheared fluids [12] [13] [14] . In the case of electron beams the instability is usually called slipping stream or diocotron. For perturbations propagating in the x-direction according to f(y) exp{i [kx - ωt]} (k real), the only waves of interest are two surface waves, associated with the two free edges of the beam. For these two waves one obtains the dispersion relation:

ω ω0 2 = ka − ch 2kb − ch 2ka 2sh2kb 2 sh 2 k(b−a) sh2kb 2 (3.6)

This dispersion relation is shown in Fig. 2, for b/a = 1, 2, 3 and 8.

The principal features of the dispersion relation are as follows: when the right hand side of Eq. (3. 6) is positive, the two waves propagate in opposite directions with real phase velocities of equal magnitude. When the righthand side of Eq. (3.6) is negative, the two waves have equal and opposite pure imaginary phase velocities, so that one wave represents pure exponential growth and the other pure exponential damping. This symmetry comes about as a result of the symmetric nature of the model in the chosen frame of reference. For b/a ≤ 2, stability is always present, while for b/a > 2, long waves are unstable, up to some critical wave number; beyond this wave number, stability is present.

Figure 2 - The dispersion relation (Eq. 3.6) appropriate to the configuration of Fig. 1. For b/a ≤ 2, stability is always present, but for b/a > 2 sufficiently long waves are always unstable. The curves for b/a = ∞ are indistinguishable from those shown for the case b/a = 8. When the wavelength is restricted to take on the values 2b, or sub-multiples of 2b, the abscissa is quantized. The first two permissible modes for the case b/a = 8 are shown; only the fundamental is unstable.

It was stated in the Introduction that a condition for stability of this type of crossed-field electron beam is that the beam should be re-entrant. The reason for this can be seen from Fig. 2, for if arbitrarily long waves are permitted, beams for which b/a > 2 will always be unstable. In practice, the re-entrant character of the beam is obtained by using cylindrical or toroidal arrangements. For our present model, the re-entrant feature is simulated by requiring the beam to have a finite periodicity in the x-direction. In particular, we require all potentials to be periodic in the x-direction with periods equal to 2b and sub-multiples of 2b. Thus the computer calculation is restricted to a square area of side 2b. Three squares of this kind are shown in Fig. 1. From the point of view of Eq. (3, 6), we require k to have only the values

ka = a/b m π, m = 1, 2, . . . (3.7)

For the case b = 8a, these wave numbers are marked on Fig. 2, for m = 1 and 2; it will be seen that only the fundamental mode m = 1 is long enough to be unstable, with a growth rate equal to 0.2 ω0.

The inertial terms ignored in Eq. (3.1) are characteristically smaller than the electric and magnetic terms by the factor q . The inclusion of these terms has been shown (in a cold approximation) to yield a short-wave instability with a growth rate ∼ ω0e-2/q. In practical cases, this growth is negligible for about q < 1/20.

The linear theoretical analysis gives no information on the behavior of the system after saturation. However, the aerodynamic analogy suggests the formation of a vortex street. The formation of such a street is shown by the computer calculations (Section 5).

IV> COMPUTATION TECHNIQUE

In the computer calculation, the continuous distribution of charge implied by Eq. (2.1) is replaced by 1,000 discrete rods of charge whose motion is calculated step-wise in time subject to a self-consistent electrostatic field. This is essentially a Lagrangian technique. Such a technique is, of course, perfectly sound; however, unlike many neutral plasma and high density (q ~ 1) electron beam calculations, it is not positively required. This point was made in justification of the Eulerian approach to the analysis of the previous section.

For our calculation, the Lagrangian approach was chosen because of the existence of a substantial program in this form [15], developed by one of us (R.W.H.) for use in the study of the anomalous diffusion of plasmas. This program required only slight modification before it could be used for the present studies.

The positions and velocities of the 1,000 charged rods are advanced stepwise in time as follows: (1) from the position of each charge, one derives a charge distribution; (2) from the charge distribution, using Poisson's equation, one calculates the electrostatic potential; (3) from the potential, one derives the electric field acting on each charge; (4) one then integrates Newton's law of motion for short time step DT, using the local electric field and a constant external magnetic field. This calculation gives a new position and velocity for each rod, and the cycle repeats at step (1). As an initial condition the 1,000 rods are distributed randomly within a band in the interior. The initial velocities are selected to give a Maxwellian distribution with the temperature as an input quantity. The errors inherent in such a model are discussed in Ref. 15, and the special rapid technique used for the solution of Poisson's equation is described in Ref. 16.

An interesting feature of the present calculation arises in connection with step (4), the integration of Newton's law of motion. The full equation of motion

m d2x~ dt2 = − e E ~ + dx~ dt x B ~ (4.1)

is represented by the central difference approximation

m [ x~ (t+DT) − 2x~ (t) + x~ (t − DT)] DT2 = − e E ~ (t) + {x~ (t+DT)−x~ (t−DT)} 2DT x B~ (4.2)

where x~ is the position vector of a gyrating rod. This scheme was recommended by Buneman [17], and has the important property of being reversible in time; it can also be used for time steps greater than a gyro-period, provided that the net displacement in each time step is small compared with the scale size for changes in the electrostatic potential. To see this we note that the solution in constant E and B fields to both the differential Eq. (4.1) and the difference Eq. (4.2) can be expressed as a gyration with Larmor radius R and gyro-frequency ωc about a guiding center drifting with a velocity v~D . That is to say, taking axes such that Ex = 0 ,

x(t) = vDt + R cos(ω ct + φ) + constant y(t) = R sin(ω ct + φ) + constant (4.3)

where φ is a constant phase determined by the initial conditions.

It is easy to relate the parameters of the motion derived from the difference equation (distinguished by a star) to the parameters of the exact motion derived from the differential equation. For simplicity we take the case that the electron has initially zero velocity (the cold case of Section V) and obtain the relations:

VD* = VD = E/B R* = R 1+ ωcDT 2 2 1/2 = VDDT 2 1 + ωcDT 2 2 1/2 ωc* = 2/DT arctan(ωcDT/2) (4.4)

R* is the radial distance of the vertices of the polygon which the electron follows in this numerical approximation. It differs as indicated from the true gyro-radius R. ω*c is the frequency with which this polygon is traversed. This frequency differs as indicated from the true gyro-frequency ωc. We note that the difference equation reproduces exactly the guiding center drift. The only errors introduced by differencing are in the gyro-frequency and Larmor radius.

Figure 3 - Growth in amplitude of the unstable fundamental mode for the thin beam (b/a = 8) case. The initial amplitude results from the random initial position of the electrons in the beam, and the growth rate is derived from Fig. 2. The time step used for this calculation was about 1/5 of a gyro period, or 1/100 of the shear period τ0 = 2 π/ω0. The indicated exponential growth rate is approximately 0.2ω0.

For the study of the linear growth regime in Fig. 3, covering five growth times, we have been able to use five time steps in a cyclotron period. This leads to R* = 1.18R and ω*c = 0.89 ωc.

The angle turned through between samplings of the electric field is ω*cDT=66°. With these values one has a reasonable representation of the full dynamic motion of the electrons.

Figure 4 - Frames from the movie exhibiting the growth and development of the thin unstable electron beam. The time step used in this calculation was approximately one gyro period, or 1/20 the shear period τ0 Times on this figure are in units of 0.98 τ0. The relatively long time step used in the calculation should not affect the long wave behavior exhibited in this figure.
Figure 5 - Energy spectrum associated with the computer run exhibited in Fig. 4. The unstable fundamental mode grows and then saturates. The total energy remains constant to within 3-1/2%. The energy residing in the higher modes starts out at about the same level as the fundamental mode ( about 10-4) and does not increase significantly during the run.

In the study of the long term non-linear behavior in Fig. 4 and Fig 5, covering 25 growth times, it was necessary to use a longer step which was approximately equal to the gyro-period. In this case

R* = 3.3R = 0.52VD DT ω*c = 0.4ωc ω*c DT = 144°

The Larmor radius and gyro period are now badly distorted, hence such a large step would not be acceptable in the study of effects which depended on electron inertia. However, since we are studying the diocotron instability which is independent of the electron mass, this is unimportant. It is only necessary that the correct guiding center drift motion be obtained. This will be the case if the distance moved by the guiding center in a time step and the Larmor radius are small compared with the wavelength, λ, of the diocotron wave. The relations (4.4) show that for large time steps the first condition guarantees that the second is satisfied. In our case we have taken q=1/20, which leads to

VDDT λ = 0.2 , R* λ = 0.01

where λ=2b is the wavelength of principal interest. The description of very short waves is, however, likely to be in error, The shortest wavelength available to the program is ~2b/24 . The sampling of the field every 144#176; ensures that the field is evenly sampled on either side of the line of motion of the guiding center, Thus we expect the guiding center motion to be correctly represented even though the time step is too large to represent the inertial aspects of the electron motion correctly.

The main output from the computer model is produced on a Stromberg-Carlson 4020 computer recorder. This may be a complete record of the rod positions and potential distribution or alternatively simply a record of the rod positions. The latter form of output is much less costly in computer time and is shown in Fig. 4. Either of these displays can be brought out on movie film for projection purposes. A time step, with the record of rod positions only, takes 2.1 secs on an IBM 7090 computer; the total computer time for a typical run of 250 steps with movie output was 17 minutes

V. RESULTS

Linear Growth Regime

Two cases have been run, namely b/a = 2 (the thick stable beam) and b/a = 8 (the thin unstable beam). For both cases q was taken as 1/20, which is the same value used by Byers [18] in a similar experiment. This is sufficiently small to ensure that we are studying the long wavelength phenomena, since the short wave growth rate is estimated to be ω0e-2/q≈4×10-18ω0 - leading to a growth time roughly 1018 times longer than the calculation time. On the other hand q is sufficiently large for us to take approximately five steps of the calculation in a gyro period of the electron. Five steps per gyro period has been used successfully previously by Hockney [15] and enables the full orbit of the gyrating electron to be followed. These results thus include the inertia of the electrons and we do not make the approximation of guiding center motion.

For the long wave diocotron instability the growth time is of the order of τ=2π/ω0, (see Fig. 2). In the case of the thick beam there was no measurable growth of any instabilities over a period of five growth times, and the beam remained well confined to its initial width.

In the case of the thin beam we find, as predicted theoretically, that the m = 1 wave is unstable. Fig. 3 shows the excellent agreement that was obtained between the theoretical growth rate and that observed in the computer model. This agreement gives us confidence in believing the new results which are found with computer model and which are not predicted by linear theory. These relate to the large amplitude, non-linear behavior of the system. We observe in Fig. 3 that the amplitude of the m = 1 wave saturates after about 3τ0 At this time the beam has twisted itself into a vortex but the amplitude of the vortex is not sufficient to drive any of the electrons to the wall. This result is in agreement with that of Byers [18] who found that the Kelvin-Helmholtz instability gave rise to confined stable vortices in a two-fluid Eulerian model of a plasma.

The above runs were performed for two initial conditions. In the first, or cold case, the electrons, which are given random positions within the beam, have zero initial velocity. The kinetic energy which they quickly acquire is due entirely to the electrostatic field of the beam itself. In the second case a small initial velocity, selected from a Maxwellian distribution, is given to the electrons. The characteristic thermal velocity , Vth=(kT/m)1/2 was taken as approximately one half of the maximum shear velocity of the beam. We observed no difference in the macroscopic long wave behavior of either the thick or thin beams. This result is in agreement with the fact that the long wave instability is predicted from guiding center theory which is not influenced by the electron thermal velocity.

Long Term Non-linear Behavior

For the long term behavior we have used a time step approximately equal to the gyroperiod, and the results have been displayed in a short 16mm movie. Still frames from this movie are reproduced in Fig. 4 in which each electron-rod is represented by a dot. As expected from the analogous situation in hydrodynamics, the beam breaks up into a vortex street. The vortices are completely formed at t= 5.0τ0, and persist without change for at least the next twenty shear periods.

As a check on this long term behavior we have also computed the variation of electric-field energy during this time, and its division between the different wave numbers. This is shown in Fig. 5. We find that the total field energy drops by 3-1/2% during the 500 time steps represented in Fig. 4 and Fig. 5. This represents the cumulative effects of all the finite difference errors arising in the computer model. These are the finite size of the space-mesh used in the solution of Poisson's equation, and the finite size of the time step. The energy in the m = 1 unstable wave is seen to increase by 3 orders of magnitude during the first three growth times and then to saturate with about one eighth of the total field-energy. The energy which goes into the m = 1 wave is taken from the electrostatic energy of the initial beam which is the energy in the m = 0 wave. An interesting feature of the computer results is that this interchange of energy between the m = 0 and m = 1 waves persists on a smaller scale for the duration of the run. One can see a periodic exchange of energy with a period approximately equal to 2τ0 . The amplitude of this oscillation appears to decrease by a factor two during the run. The effect of the oscillation on the form of the vortex is minor and can be seen by comparing the vortex in Fig. 4 and t = 5τ 0 when the m = 1 wave is a maximum and with that at t=10τ 0, when the wave is a minimum.

VI. CONCLUSIONS

A Lagrangian computer model consisting of 1,000 electron rods has been used to investigate the long wave diocotron instability in a crossed-field electron beam. A thick and thin beam were found to be stable and unstable respectively, as predicted by theory. Excellent agreement was obtained between the theoretical growth rate and that observed in the model. The computer model was able to follow the non-linear behavior of the system. The unstable beam broke up into a vortex street which remained for a long time, The size of the vortex was not enough to transport electrons to the wall of the system and the beam remained confined to the interior. The results of this paper have been made into a 4-minute, l6mm movie.

ACKNOWLEDGMENTS

It is a pleasure to acknowledge numerous helpful discussions held with Professor O. Buneman, and also some useful conversions with Dr. J. Byers and Professor C. K. Birdsall and Dr. H. E. Petschek.

This research was supported by the following contracts: George C. Marshall Space Flight Center, National Aeronautics and Space Administration, Contract NAS8-20310; Air Force Office of Scientific Research, Office of Aerospace Research, United States Air Force, Washington, D.C., Contract AF49(638)-1553; National Science Foundation Grant GK 625 and Tri-Service Contract NONr-225(83). The development of the program used in the calculations was supported by the Tri-Service contract.

REFERENCES

1. O Buneman, R H Levy, L M Linson, Stability of Crossed-Field Electron Beams, J. Appl. Phys. 37, 3203-3222 (1966).

2. J C Helmer, R L Jepsen, Electrical Characteristics of a Penning Discharge, Proc. I.R.E. 49, 1920 (1961).

3. R H Levy, G S Janes, Plasma radiation shielding, AIAA J. , 1835 (1964).

4. G S Janes, R H Levy, H E Petschek, Production of BeV Potential Wells, Phys. Rev. Letters 15, 138 (1965).

5. G S Janes, R H Levy, H A Bethe, B T Feld, New Type of Accelerator for Heavy Ions, Phys Rev. 145, 925 (1966).

6. R H Levy, Diocotron Instability in a Cylindrical Geometry, Phys Fluids 8, 1288 (1965)

7. R H Levy, Effect of Coherent Radiation on the Stability of a Crossed-Field Electron Beam, J. Appl. Phys. 37, 119 (1966).

8. R H Levy, J D Callen, Diocotron Instability in a Quasi-Toroidal Geometry, Phys. Fluids 8, 2298 (1965).

9. J D Daugherty , R H Levy, Equilibrium of Electron Clouds in Toroidal Magnetic Fields, Phys. Fluids 10, 155 (1967).

10. G S Janes, Experiments on Magnetically Produced and Confined Electron Clouds, Phys. Rev. Letters 15, 135 (1965).

11. G S Janes, Experiments Relating to the Confinement of Electron Plasmas in Toroidal (HIPAC) Type Geometries, presented at the 8th Annual Meeting of the American Physical Society, Division of Plasma Physics, November 1966.

12. S Chandrasekhar, Hydrodynamic and Hydromagnetic Stability, Clarendon Press, Oxford, 1961), Chapter XIII.

13. H Helmholtz, Ueber discontinuirliche Flussigkeitsbewegungen, Wissenschaftliche Abhandlungen (J. A. Barth, Leipzig, 1882), pp. 146-57; translation by Guthrie in Phil. Mag., Ser. 4, 36, 337-46 (1868).

14. Lord Kelvin, Hydrokinetic Solutions and Observations, "On the Motion of Free Solids Through a Liquid," 69-75, "Influence of Wind and Capillarity on Waves in Water Supposed Frictionless, 11 76-85, Mathematical and Physical Papers, iv. Hydrodynamics and General Dynamics (Cambridge University Press, Cambridge, 1910).

15. R W Hockney, Computer Experiment of Anomalous Diffusion, Phys. Fluids 9, 1826 (1966).

16. R W Hockney, A Fast Direct Solution of Poisson's Equation Using Fourier Analysis, JACM, 12, 95 (1965).

17. O Buneman, Technical SUIPR Report No. 110, Institute for Plasma Research, Stanford University, October 1966 (to be published in J. Computational Phys. ) .

18. J A Byers, onfined and Nonconfirmed Interchange Instabilities Obtained from Nonlinear Computer Models, Phys. Fluids 2, 1038 (1966).

More Computer Animation Papers 1964-1976