Nonlinear Phenomena in a Pure Electron P

Nonlinear Phenomena in a Pure Electron Plasma Studied with a 2-D Fluid Code - CaltechTHESIS
CaltechTHESIS
A Caltech Library Service
About
Browse
Deposit an Item
Instructions for Students
Nonlinear Phenomena in a Pure Electron Plasma Studied with a 2-D Fluid Code
Citation
Bachman, David Alan
(1998)
Nonlinear Phenomena in a Pure Electron Plasma Studied with a 2-D Fluid Code.
Dissertation (Ph.D.), California Institute of Technology.
doi:10.7907/72z8-gm50.
Abstract
This thesis presents a new computational tool for studying cylindrical pure electron plasmas. Previous research used linear methods to describe the evolution of small plasma perturbations. This new tool numerically solves the nonlinear two- dimensional (2-D) fluid equations in cylindrical coordinates, allowing the exploration of many phenomena that have been observed experimentally in these plasmas. Most experimentally observed phenomena are large in amplitude and follow nonlinear dynamics. Clearly, codes based on the linearized equations can not reproduce these nonlinear phenomena.

The plasma was first studied for small amplitude perturbations using the nonlinear fluid code, producing results that agree with linearized calculations. The perturbed electric field decays in time, due to shear in the flow of density perturbations at different radii, which results in the total contribution to the perturbed electric field phase mixing away.

At higher amplitudes, the decay envelope becomes modulated, which is a result that has been observed experimentally and is also reproduced by this fluid code. The modulation is caused by nonlinear trapping of fluid elements within the wave, which is illustrated in images of the perturbed density.

For two applied pulses separated in time, a third echo response is observed after the responses to the two applied pulses have decayed away. Echoes have been observed experimentally in neutral plasmas, but have not yet been observed experimentally in non-neutral plasmas.

At very high amplitudes, a nonlinear decay instability occurs. A high amplitude wave decays into a wave with lower azimuthal symmetry number due to an interaction occurring at the beat frequency between the two waves. The beat-wave decay instability has been observed experimentally, and is also observed using the 2-D nonlinear fluid code.

The 2-D cylindrical nonlinear fluid code is capable of reproducing a wide range of non-neutral plasma phenomena, and is an important new tool for future research on non-neutral plasmas. This research is also relevant to ordinary fluids, since the equations describing a non-neutral plasma are analogous to the 2-D Euler equations for an inviscid fluid.
Item Type:
Thesis (Dissertation (Ph.D.))
Subject Keywords:
(Applied Physics)
Degree Grantor:
California Institute of Technology
Division:
Engineering and Applied Science
Major Option:
Applied Physics
Thesis Availability:
Public (worldwide access)
Research Advisor(s):
Gould, Roy Walter
Thesis Committee:
Gould, Roy Walter (chair)
Corngold, Noel Robert
Bellan, Paul Murray
Defense Date:
1 September 1997
Record Number:
CaltechETD:etd-01232008-081243
Persistent URL:
DOI:
10.7907/72z8-gm50
Default Usage Policy:
No commercial reproduction, distribution, display or performance rights in this work are provided.
ID Code:
299
Collection:
CaltechTHESIS
Deposited By:
Imported from ETD-db
Deposited On:
15 Feb 2008
Last Modified:
17 Jul 2023 23:30
Thesis Files
Preview
PDF (Bachman_da_1998.pdf)
- Final Version
See Usage Policy.
21MB
Repository Staff Only:
item control page
CaltechTHESIS is powered by
EPrints 3.3
which is developed by the
School of Electronics and Computer Science
at the University of Southampton.
More information and software credits
Nonlinear Phenomena in a Pure Electron Plasma

Studied with a 2-D Fluid Code

Thesis by

David A. Bachman

In Partial Fulfillment of the Requirements
for the Degree of
Doctor of Philosophy

California Institute of Technology

Pasadena, California

1998
(Submitted September, 1997)

David A. Bachman

ill

Acknowledgements

It has been a pleasure and a privilege to work towards my Ph.D. under the ex-
pert guidance of Professor Roy Gould. He provided clear explanations of concepts I
struggled over, and more resources than any graduate student deserves. Along with
technical knowledge, I have also gained much knowledge outside of plasma physics
from our many conversations. I wish to thank Professor Gould for all the effort he
has invested in me.

I also want to thank Professor Noel Corngold and Professor Paul Bellan. They
were both involved in teaching me many concepts in plasma physics, and greatly
assisted in the completion of my degree.

My experiences at Caltech have been interesting and positive. I wish to thank
the outstanding faculty and staff for providing me with this positive learning environ-
ment. The Caltech staff have always been friendly and more than willing to go out
of their way to assist me. Special thanks to Connie Rodriguez for helping me with
administrative matters and providing pleasant conversation.

I would also like to thank the Office of Naval Research and the contract monitor
Dr. Charles Roberson for funding this effort in non-neutral plasmas. Without Dr.
Roberson’s strong support for the program, non-neutral plasma physics would not be
where it is today.

Fellow graduate students Jimmy Yee, Freddy Hansen and Steve Sanders provided
very useful feedback on the rough draft of this manuscript. This thesis was signifi-
cantly improved by their involvement.

I want to thank my parents for their constant support of my educational goals.
Without their encouragement, I could not have progressed this far in my education.
I also received emotional support from my girlfriend, Yun Huang. She provided the

final encouragement I needed to complete this thesis.

Abstract

This thesis presents a new computational tool for studying cylindrical pure elec-
tron plasmas. Previous research used linear methods to describe the evolution of
small plasma perturbations. This new tool numerically solves the nonlinear two-
dimensional (2-D) fluid equations in cylindrical coordinates, allowing the exploration
of many phenomena that have been observed experimentally in these plasmas. Most
experimentally observed phenomena are large in amplitude and follow nonlinear dy-
namics. Clearly, codes based on the linearized equations can not reproduce these
nonlinear phenomena.

The plasma was first studied for small amplitude perturbations using the nonlinear
fluid code, producing results that agree with linearized calculations. The perturbed
electric field decays in time, due to shear in the flow of density perturbations at
different radii, which results in the total contribution to the perturbed electric field
phase mixing away.

At higher amplitudes, the decay envelope becomes modulated, which is a result
that has been observed experimentally and is also reproduced by this fluid code. The
modulation is caused by nonlinear trapping of fluid elements within the wave, which
is illustrated in images of the perturbed density.

For two applied pulses separated in time, a third echo response is observed after
the responses to the two applied pulses have decayed away. Echoes have been observed
experimentally in neutral plasmas, but have not yet been observed experimentally in
non-neutral plasmas.

At very high amplitudes, a nonlinear decay instability occurs. A high amplitude
wave decays into a wave with lower azimuthal symmetry number due to an interaction
occurring at the beat frequency between the two waves. The beat-wave decay insta-
bility has been observed experimentally, and is also observed using the 2-D nonlinear

fluid code.

The 2-D cylindrical nonlinear fluid code is capable of reproducing a wide range

of non-neutral plasma phenomena, and is an important new tool for future research
on non-neutral plasmas. This research is also relevant to ordinary fluids, since the
equations describing a non-neutral plasma are analogous to the 2-D Euler equations

for an inviscid fluid.

Contents

Acknowledgements

Abstract

1 Introduction

1.1
1.2
1.3

1.4
1.5

2 The
2.1
2.2

2.3
2.4
2.5
2.6
2.7

Overview... 6
Background ... 2... ... 0. ee
Pure Electron Plasma Geometry. ............--..000.4
1.3.1 Plasma Confinement ...............-.22.+200-
1.3.2 Trap Operation ..........0. 0020020 eee eee
1.3.3 Parameters in the Caltech Pure Electron Plasma Experiment .
Previous Work at Caltech ..........0...02.2000 000004
Thesis Outline. 2... 2.

Fluid Model

Introduction... 2... 2 ee
Model Assumptions and Approximations ................
2.2.1 Approximations to Maxwell’s Equations ............
2.2.2 Approximation of a Two-dimensional Plasma .........
2.2.3. Approximation of Zero Pressure Gradient. ...........
2.2.4 Drift Approximation ..................20-04
Steady-State Angular Velocity .................0.024.4
Nonlinear Model .. 2... 0... 2
Normalized Variables... ........0.0 0000008 ee eee
Linearized Model ... 2... 2 2

Summary .. 0... ee

ili

oN Oo F A BR YO HS &

vil

3 Computational Methods

3.1 Introduction. ........2..020200.000..000000000000.
3.2 Nonlinear Fluid Code... ........2.0.0....000 000000.
3.2.1 Poisson’s Equation ...........0.0.0.0.000 000004
3.2.2 Density Derivative Routine... ........2.0.2.020000.
3.2.3 Leapfrog Time Advance ................0.0000.
3.2.4 Diagnostics and Outputs... .......0...0.2.00000.
3.2.6 Accumulated Errors ................20000004
3.3 Linearized Fluid Code ..........0...0.......00000.
3.4 Conclusion... 2... 2 ee

4 Small Amplitude Disturbances
4.1 Introduction. ............000.0 0000000000 eee
4.2 Collisionless Damping ..............2.....000.000.4
4.2.1 Fast Decay Profle ...................00048.
4.2.2 Slow Decay Profle ...........0...000.0.000.0.
4.2.3 No Decay Profle ...............0....0000.
4.2.4 Other Profiles. ...........20.0.0.000000000040.
4.2.5 Briggs, Daugherty, and Levy Decay Rate... ......2..
4.3 Density Perturbations ................0..00.2..040.
4.3.1 Fast Decay Profle ............0..0..00.0.2.00..
4.3.2 Slow Decay Profile .....................0..
4.3.3 No Decay Profile .............00. ree
4.4 Linear Model ..........20.0...0..00 000000000 eae
4.5 Conclusion... 2... 0.000000 000 0c ee
5 Fluid Trapping
0.1 Introduction... .... 0.020.000.0000 00000000 ee
5.2 Fluid Code Results... 2... ..0.0.0.0000002000000000004
0.2.1 Perturbed Electric Field ..........0.0.0...0.002.
5.2.2 Bounce Frequency Dependence .................

23
23
24
29
29
dl
32
36
37
39

Al
41
42
42
49
52
59
57
58
58
60
63
63
69

5.2.3
5.2.4

5.3. Conclusion

Density Perturbations ..........

BGK-Like Modes .............

6 Echoes in a Pure Electron Plasma

6.1 Introduction
6.2. Echo Concept
6.3 Echo Results

6.3.1
6.3.2
6.3.3
6.3.4
6.3.5

6.4 Conclusion

Demonstration of the Echo. .......
Control Sets. ............-..
Amplitude of the Echo ..........
Higher Amplitude Effects ........
Slow Decay Density Profile. .......

7 Beat-Wave Decay Instability

7.1 Introduction

7.2.1
7.2.2
7.2.3
7.2.4

7.3. Conclusion

Choice of Density Profile.........

Demonstration of the Decay Instability

Growth Rate Versus Amplitude .... .

Comparison to Experimental Results . .

8 Conclusion and Future Work

8.1 Summary of Results

8.2. Further Work

8.2.1

8.2.2 Further Research on Non-neutral Plasmas

Bibliography

Improvements to the Nonlinear Fluid Code

92
92
93
99
99
102
107
113
121
122

124
124
126
126
128
128
132
134

136
136
137
137
139

141

List of Figures

1.1
1.2

3.1

4.1

4.2

4.3

4.4

4.5

4.6

4.7

Schematic diagram of a cylindrical Penning trap. ...........
Plot of the envelope of the perturbed electric field at the wall ver-
sus time on a semilog scale for the experimental data collected by N.

Sateesh Pillai... 2... ee
Flow chart describing the logic of the nonlinear fluid code. ......

Plot of the normalized unperturbed density and normalized angular
velocity versus radius for the Fast Decay Profile no(r) = no(0){1 —
3(£)? + 3(£)* — (£)°} with the edge of the plasma at a = 0.8.

Plot of the perturbed electric field at the wall versus time for the Fast
Decay Profile no(r) = no(0){1 — 3(£)? + 3(£)* — (£)°} with a = 0.8. .

Semilog plot of the amplitude of the Fourier components of the electric

field at the wall versus time for the Fast Decay Profile no(r) = no(0){1—

3(£)? + 3(£)* — (£)°} witha=08........0...20.000.2.
Fourier components of the electric field at the wall versus time for the
Fast Decay Profile no(r) = no(0){1—3(£)?+3(£)*—(£)°} with a = 0.8
and a pure m = 2 wall excitation. .. 2... 2.2... ...2.2204.
Plot of the normalized unperturbed density and normalized angular
velocity versus readius for the Slow Decay Profile no(r) = no(0){1 —
6(£)* + 8(£)® — 3(£)8} with the edge of the plasma at a=0.8.... .
Fourier components of the electric field at the wall versus time for the
Slow Decay Profile no(r) = no(0){1 — 6(£)* + 8(£)® — 3(£)8} with
a = 0.8 and a pure m = 2 angular wall excitation... .........
Plot of the normalized unperturbed density and normalized angular
velocity versus radius for the No Decay Profile no(r) = mo(0){1 —

10(£)® + 15(£)® — 6(£)°} with the edge of the plasma at a= 0.8. . .

4.8

4.9

4.10

4.11

4.12

4.13

4.14

4.15

4.16

5.1

5.2

5.3

Plot of the perturbed electric field at the wall versus time on a semilog
scale for the No Decay Profile no(r) = no(0){1 — 10(£)® + 15(£)8 —
6(£)}, witha=08. 2... ee
Plot of the normalized density versus radius for the normalized Fermi
function with 6 = 2,4,6,8,10,12,14,and16...............
Density perturbations at six successive times for the Fast Decay Profile
no(r) = no(0){1 — 3(£)? + 3(£)* — (£)°} witha =08. ........
Density perturbations at six successive times for the Fast Decay Profile
no(r) = no(0){1 — 3(£)? + 3(£)* — (£)°} with a = 0.8, assuming the
plasma convected in the theta direction at the unperturbed angular
velocity, ignoring the self consistent perturbed electric field... ....
Density perturbations at six successive times for the Slow Decay Profile
no(r) = no(0){1 — 6(£)* + 8(£)® — 3(£)°} witha =08.........
Density perturbations at six successive times for the No Decay Profile
no(r) = no(0){1 — 10(£)° + 15(£)8 — 6(£)"9} witha =0.8. ......
Magnitude of density perturbations versus radius at six successive
times for the Fast Decay Profile... .............-..200.
Magnitude of continuum density perturbations versus radius for the
Fast Decay Profile obtained after the quasimode has decayed away.
Resulting density perturbations for six different times for the Fast De-

cay Profile using the Quasimode and Continuum mode linear model. .

Plot of the envelope of the perturbed electric field at the wall ver-
sus time on a semilog scale for the experimental data collected by N.
Sateesh Pillai. 2... ee
Plot of the perturbed electric field at the wall for the Slow Decay
Profile, using an applied three cycle pulse amplitude of 4.28 x 107°¢p.
Plot of the perturbed electric field at the wall demonstrating bounce

oscillations. 2... 20.0.0.

5.4

5.9

5.6

o.¢

5.8

5.9
5.10

5.11

6.1

6.2

xi
Plot of the envelope of the perturbed electric field versus time on a
semilog scale for various amplitude pulses applied to the Slow Decay
Profile... 2... 77
Plot of the angular bounce frequency versus the perturbed electric field
at the wall for 12 high amplitude traces. ........0....0.. 80
Plot of the envelope of the perturbed electric field at the wall for a
typical computation displaying bounce oscillations. .......... 81
Images of the total density at six different times. An external Gaussian
pulse with a maximum amplitude of 0.005¢,, and a FWHM of 1.06t¢,
is applied to the Slow Decay Profile. The resulting perturbed electric
field at the wall exhibits strong bounce oscillations. .......... 83
Images of the perturbed density at six different times for an external
Gaussian pulse with a maximum amplitude of 0.005¢), and a FWHM
of 1.06¢., applied to the Slow Decay Profile. The resulting perturbed
electric field at the wall exhibits strong bounce oscillations... ... . 84
Magnification of the trapped fluid elements causing bounce oscillations. 86
Plot of the envelope of the perturbed electric field at the wall resulting
from a slow Gaussian pulse applied to the Slow Decay Profile and
exciting a BGK-like mode. ..........0..000..0..0000. 89
(a) Perturbed density resulting from slowly applying a Gaussian pulse
to the Slow Decay Profile, which excites a BGK-like mode. Image (b)

is a magnification of the outer plasma edge, stretched out flat. .... 90

Plot of the perturbed electric field at the wall versus time for a calcu-
lation which demonstrates the echo.................0.. 100
Plot of the perturbed electric field at the wall versus time for a one
cycle m, = 2 pulse at time t = 0, and a one cycle mz = 6 pulse at time
T = 25t., applied to the Fast Decay Profile. An echo occurs at time
t = 37.5t¢r with m3 = 4 azimuthal symmetry............... 101

6.3

6.4

6.5

6.6

6.7

6.8

6.9

6.10

6.11

6.12

6.13

6.14

6.15

6.16

xii

Perturbed electric field at the wall versus time for a one cycle pulse
with m = 2 azimuthal symmetry applied to the Fast Decay Profile.

Perturbed electric field at the wall versus time for a one cycle pulse
with m = 4 azimuthal symmetry applied to the Fast Decay Profile.

Perturbed electric field at the wall versus time for a one cycle pulse
with m = 6 azimuthal symmetry applied to the Fast Decay Profile.

Perturbed electric field at wall for a one cycle m, = 4 pulse at time
t = 0, and a one cycle m2 = 2 pulse at time 7 = 25t,, applied to the
Fast Decay Profile. Both pulses decay away, with no echo produced. .
Perturbed electric field at wall for four echo computations varying in
applied amplitude. ........0.00020.000000000. 000004
Plot of the maximum amplitude of the echo versus applied m, = 2
amplitude of the first pulse on a log-log scale. . 2... 2...
Plot of the maximum amplitude of the echo versus applied m,; = 4
amplitude of the second pulse on a log-log scale... 2... 2. 2...
Perturbed electric field at wall for three echo computations, varying
the interpulse spacing, T...........0.0.0.. 00004 eee ee
Plot of the amplitude of the echo versus interpulse spacing of the ap-
plied pulses on a log-log scale... 2... ee
Response of the m; = 2 pulse, mz = 4 pulse, and m3 = 2 echo response
as a function of pulse 1 amplitude and pulse 2 amplitude... .....
Perturbed electric field at the wall versus time for two high amplitude
applied pulses... 2.2...
Semilog plot of a computation demonstrating higher order echoes. . .
Plot of the perturbed electric field at the wall versus time for an m = 2
applied pulse at time t = 0, and an m2 = 6 applied pulse at time
T = 25t-r, demonstrating an echo of anecho...............
Plot of the perturbed electric field at the wall versus time for a com-

putation demonstrating the echo while using the Slow Decay Profile. .

103

104

105

107

109

110

111

112

114

116

122

7.1

7.2

7.3

7.4

7.5

7.6

7.7

xill
Experimental results of Mitchell and Driscoll for the beat-wave decay
instability in a pure electron plasma. .........-.....2056.
Plot of the normalized density and normalized angular velocity for the
Beat Profile... 2.0.0.0. 0202 ee
Plot of the envelope of the perturbed electric field at the wall versus
time on a semilog scale for the m = 1,2 and 3 Fourier components
due to a large applied pulse with m = 3 azimuthal symmetry. The
beat-wave decay instability causes the decay of the m = 3 mode and
growth of them =2 mode.............. 0000005 2 eee
Images of the total density at times 2.47¢.,, 6.92te, 11.4t¢,, 15.8te,,
20.3t,,, and 24.7t,,, after excitation from a large pulse with m = 3
azimuthal symmetry. Due to beat-wave resonance damping, the m = 3
mode decays into anm=2mode........----+++2 00005
Plot of the envelope of the perturbed electric field at the wall for the
m = 3 and m = 2 mode versus time on a semilog scale. ........
Estimated growth rates 7. for daughter mode m = 2 normalized to the

central rotation time versus normalized amplitude A3 of parent mode

Plot of the experimental and computational growth rates for the m = 2

daughter mode versus parent mode m = 3 amplitude. .........

126

127

129

130

131

134

List of Tables

1.1

2.1

4.1

4.2

6.1

6.2

7.1

Typical physical parameters for the Caltech Pure Electron Plasma Ex-

periment. 2...

Description of the normalized variables used in the numerical calcula-

Decay Rates for Various Profiles. ................0004
Comparison with Briggs decay rate for various profiles. For density
profiles which deviate strongly from a step function, the Briggs decay

formula is not very accurate... 2... 20.0... ee ee

Echo amplitude versus applied pulse 1 and pulse 2 amplitude. Pulse
1 is a one cycle pulse with m, = 2 symmetry and is applied at time
t = 0. Pulse 2 is a one cycle pulse that has mz = 4 symmetry and is
applied at time 7 = 25t,,. The echo has m3 = 2 symmetry and occurs
at time t = 50t,,. The echo amplitude is given as a percentage of the
unperturbed electric field at the wall... 2... 2... 2... ee,
Echo amplitude as a function of interpulse spacing. A one cycle pulse
with m, = 2 symmetry and amplitude 1.0 x 107d, is applied at time
t = 0. A second one cycle pulse with m2 = 4 symmetry and amplitude
2.0 x 1073, is applied at time t = 7. The echo amplitude is measured

as a percentage of the unperturbed electric field at the wall. .....

Growth rate of the m = 2 mode as a function of m = 3 mode initial

amplitude... 2... ee

108

Chapter 1 Introduction

1.1 Overview

This thesis presents results obtained using a new computational tool for studying a
cylindrical pure electron plasma. Previous research has studied these plasmas using
linear methods to describe the evolution of small perturbations of the steady-state
density and electric potential. This newly developed tool numerically solves the full
nonlinear two-dimensional (2-D) fluid equations in a cylindrical geometry. The non-
linear code allows the exploration of many phenomena that have been observed ex-
perimentally in these plasmas. Most of the experimentally observed phenomena are
of large amplitude and therefore follow nonlinear dynamics. The previous codes based
on the linearized equations can not reproduce these high amplitude phenomena.

To test the nonlinear fluid code, the plasma response was studied for small ampli-
tude perturbations and compared with results from the linearized theory. For small
amplitude perturbations, the results from the nonlinear fluid code were found to agree
with the results from linearized calculations. The perturbed electric field at the wall
was found to decay exponentially in time. This decay is due to a process where density
perturbations at different radial layers undergo shear, so that the total contribution
to the electric field at the wall phase mixes away.

At higher amplitudes, a modulation of the decay envelope occurs. This modulation
has been observed experimentally and is qualitatively reproduced by the fluid code.
The modulation is caused by the nonlinear trapping of fluid elements within the
excited wave. The trapped fluid elements can be observed in images of the perturbed
density that are produced by the nonlinear fluid code.

For two applied pulses that are separated in time, a third echo response is observed
which occurs after the responses to the two applied pulses have decayed away. Echoes

have been observed experimentally in neutral plasmas, but have not yet been observed

experimentally in a pure electron plasma.

At very high amplitudes, a nonlinear decay instability can occur. A high amplitude
wave can decay into a wave with a lower azimuthal symmetry mode number due
to an interaction at the beat frequency between these two waves. The beat-wave
decay instability has been observed experimentally, and is also found to occur in
computations using the 2-D nonlinear fluid code.

The 2-D cylindrical nonlinear fluid code is capable of reproducing a wide range of
non-neutral plasma phenomena. It is expected that this nonlinear fluid code will be
an important new tool for future research on non-neutral plasmas. This research can
also be applied to other single species plasmas such as ion, positron, and antiproton
plasmas, and also to ordinary fluids, since the 2-D equations describing a non-neutral

plasma have the same form as the 2-D Euler equations for an inviscid fluid.

1.2 Background

Plasma is common throughout the universe, comprising more than 99% of the lu-
minescent matter within it. On Earth, plasma may be found in neon signs or in
laboratory experiments attempting to harness fusion energy from a confined plasma.
These plasmas are neutral in total electrical charge. Non-neutral plasmas can also
exist, usually comprised of a single species of charged particles; e.g., only electrons,
ions, positrons, or antiprotons.

Non-neutral plasmas have numerous applications. Non-neutral plasmas are similar
to the beams in microwave devices. Non-neutral plasma traps allow antiprotons
produced in a collider to be transported to experiments located a large distance away
[1][2]. Such a trap may eventually hold antimatter fuel for an extended space flight.
The National Institute of Standards and Technology is pursuing research directed at
using a non-neutral plasma to produce a clock that is more precise than the atomic
clocks of today [3]. Another research effort is attempting to obtain fusion in a non-
neutral plasma [4]. These applications, as well as unforeseen future applications, will

benefit from a better understanding of the physics of non-neutral plasmas.

The history of pure electron plasmas begins with work on magnetized electron
beams. Brillouin [5] first noticed that an electron beam could be focused with a
magnetic field, where the repulsive force from the space charge of the electrons could
be balanced by a Lorentz force due to the electrons passing through a magnetic
field. Further research by Webster [6] found a magnetically focused hollow beam to
be unstable. The hollow beam instability was attributed to shear in the azimuthal
electron flow and named the slipping stream instability [7], which later became known
as the diocotron instability [8].

Research performed by Levy [9] found that an analogy existed between the ap-
proximate equations describing a magnetized electron beam and the equations for 2-D
incompressible inviscid fluid flow. Further research by Briggs, Daugherty, and Levy
[10] predicted collisionless damping of wave perturbations for plasmas with density
profiles that. are monotonically decreasing in radius. This damping was found to be
analogous to Landau damping in neutral plasmas [11].

Malmberg, deGrassie [12], and Driscoll developed the cylindrical Penning trap to
confine a non-neutral plasma and study the basic physics of these plasmas. Many
results have come from their research, including experimental results for 2-D vor-
tex dynamics [13]. They have shown that the 2-D drift-Poisson equations suffice to
describe a variety of phenomena in non-neutral plasmas, and the 2-D drift-Poisson
equations will be used in this thesis. Theoretical work using the linearized 2-D drift-
Poisson equations has been performed by Pillai [14], Corngold [15], and Spencer and
Rasband [16]. Recent research on non-neutral plasmas has been summarized in two
AIP Conference Proceedings [17][18].

The pure electron plasma laboratory at the California Institute of Technology is
actively studying the physics of non-neutral plasmas with both experimental and com-
putational efforts. Those experimental results relevant to this thesis will be summa-
rized below, after explaining more clearly the geometrical configuration of a confined
pure electron plasma in a cylindrical Penning trap. The results from the experiment

motivated the computational research described in this thesis.

1.3. Pure Electron Plasma Geometry

1.3.1 Plasma Confinement

Pure electron plasmas are typically confined in a cylindrical Penning trap [12]. A
schematic diagram of this trap is shown in Figure 1.1(a). The plasma is confined
axially by electrodes biased to a sufficiently negative electric potential. The plasma is
confined radially by a uniform externally applied axial magnetic field. In this configu-
ration, there is a radial electric field due to the space charge of the electrons. Because
the electrons are in an externally applied uniform magnetic field, they undergo cy-
clotron motion.

This thesis will focus on the low frequency phenomena related to rotation of the
plasma. The plasma rotation frequency is much smaller than the cyclotron frequency
for these plasmas. The cyclotron orbits of the electrons are not important, and the
plasma will be discussed in the guiding center drift approximation.

The radial electric field causes an azimuthal & x B drift in the motion. The
azimuthal drift occurs at every radius, and thus the plasma rotates. The rotational
angular velocity is not necessarily the same at every radius, so the flow of electrons
may be sheared. With the flow of electrons in the azimuthal direction, the electrons
are radially confined by the external magnetic field. The radial space charge repulsion
force is balanced by the inward Lorentz v’ x B force due to the plasma rotating in

this uniform magnetic field.

1.3.2 Trap Operation

The procedure for operating the cylindrical Penning trap is as follows. A uniform
external magnetic field in the axial direction surrounds the trap. Electrons are emit-
ted from a hot filament and injected into the trap. The dump electrode is biased
negatively to prevent electrons escaping from that end of the trap. After the trap is
filled with electrons, the inject electrode is biased negatively and the electrons can

not escape. During the confinement time, experiments are performed on the plasma.

(a) magnetic field

trap octupole octupole trap

collector

f sma,
filament ; | I

r.f. switch :
_ low noise
inject 5 amplifier dump

0.3-0.6 MHz data acquisition
source

vyTX Vx
Vv jv tH He
V v

low noise amplifier

(b) (c)

Figure 1.1: (a) Schematic diagram of a cylindrical Penning trap. (b) Diagram of the
first sectored electrode structure in a configuration to produce a perturbation with
m = 2 symmetry. (c) Diagram of the second sectored electrode in a configuration to

receive the response of the plasma at the wall.

Parameter Value
Magnetic Field 50 Gauss
Central Electron Density ~1x 10° cm?
Electron Temperature 2eV

Trap Radius 2.5 cm

Trap Length 40 cm
Cyclotron Frequency 140 MHz

Central Plasma Frequency | ~ 9 MHz
Wall Rotation Frequency ~ 0.15 MHz
m = 2 Diocotron Frequency | ~ 0.5 MHz
Axial Bounce Frequency 1 MHz
Larmor Radius 1mm
Debye Length lcm

Table 1.1: Typical physical parameters for the Caltech Pure Electron Plasma Exper-
iment. ,

Finally, the dump electrode is grounded, and the electrons are emptied from the trap.
The experiment is repeated periodically in an inject, confine, and dump cycle.

The plasma is surrounded by a conducting cylindrical wall, which in the Caltech
device is divided into two regions containing eight sectors each. A diagram of the
sectored electrodes is shown in Figure 1.1(b) and (c). While the plasma is confined,
voltages can be applied to the first set of wall sectors to excite perturbations in the
plasma. Figure 1.1(b) demonstrates the sector configuration to excite a perturbation
in the plasma with m = 2 azimuthal symmetry. As the plasma rotates, perturbations
in the plasma rotating past the second set of wall sectors cause an oscillating perturbed
electric field at these sectors. This signal is then amplified and digitized, as shown in
Figure 1.1(c). Sensing the oscillating perturbed electric field at the wall is the primary

diagnostic measurement for the Pure Electron Plasma Experiment at Caltech.

1.3.3. Parameters in the Caltech Pure Electron Plasma Ex-

periment

Physical parameters for the Caltech Pure Electron Plasma Experiment are summa-
rized in table 1.1. This table is based on quantities measured in the experiment

[19].

This thesis will focus on the low frequency diocotron modes that are near the
average plasma rotation frequency. The m = 2 diocotron frequency of 500 kHz
corresponds to a resonant layer of plasma with a rotation frequency of 250 kHz.

This diocotron frequency is less than the plasma frequency, and much less than the
cyclotron frequency, since the cyclotron frequency is much greater than the plasma
frequency for this pure electron plasma. The cyclotron frequency being greater than
the diocotron frequency justifies the approximation that inertial effects may be ig-
nored, and the average electron velocity is given by E x B drift. The plasma model
will be discussed in more detail in Chapter 2.

Additionally, the axial bounce frequency is greater than the diocotron frequency,
and for many phenomena the plasma may be treated as two-dimensional.

This is the basic geometry for a cylindrical pure electron plasma. The next section
will describe some previous research performed on the Caltech Pure Electron Plasma

Experiment, which motivated much of the effort in this thesis.

1.4 Previous Work at Caltech

The Caltech Pure Electron Plasma Experiment described above was used to study
the response of the plasma to an externally applied sinusoidal burst. The experiment
performed by N. Sateesh Pillai [14] found that for a small amplitude applied pulse,
the oscillating electric field of the plasma at the wall decayed exponentially, experi-
mentally verifying the exponential decay results of Briggs, Daugherty, and Levy [10].
For larger pulses, a modulation of the decay envelope was observed. ‘These results are
summarized in figure 1.2 which shows the envelope of the perturbed electric field at
the wall versus time on a semilog scale.

Linearized fluid theory, reviewed in Chapter 2, was used to explain the exponential
damping of the electric field at the wall for small applied pulses [14]. The damping
is due to a collisionless phase mixing process [10]. The decay rates for power law
density profiles has been studied analytically by Corngold [15]. The modulation of

the decay envelope for larger pulses was attributed to trapping of a small fraction of

Amplitude (dB)

Sore ee, Ci Sc, Cer cecerrre | Cron ser nenen | Onn, © Senne nennet iby Seer eeteeeeeeeereas, % POtiSC Sitti Peete’. © See renreeres, +, Shine

Time in pS 300

Figure 1.2: Plot of the envelope of the perturbed electric field at the wall versus
time on a semilog scale for the experimental data collected by N. Sateesh Pillai. The
applied m = 2 excitation voltage ranges from 10 mV for the lowest amplitude pulse,
to 1000 mV for the highest amplitude trace. The vertical scale is calibrated in dB.

plasma by the excited wave.

The experiment did not permit direct measurement of plasma density perturba-
tions. Therefore the hypothesis of trapped plasma at higher amplitudes could not
be verified by direct observation. This motivated my use of computational meth-
ods to better understand these phenomena. Computation allows access to quantities
which can not be diagnosed in the Caltech experiment, and provides verification of

the nonlinear fluid model used in the computation to describe the plasma.

1.5 Thesis Outline

This thesis is based on the results of a nonlinear 2-D cylindrical fluid code. The
fluid code is used to study the response of a pure electron plasma to one or more
externally applied pulses. First the plasma model will be discussed in Chapter 2,
followed by a discussion of the numerical techniques applied to the model equations
in Chapter 3. Chapter 4 and Chapter 5 will describe the results from the fluid code for
initial plasma conditions that have been studied in previous experiments at Caltech.
Chapter 4 examines the plasma response for a small amplitude applied pulse and
Chapter 5 examines the response for larger pulses. Chapter 6 will discuss two applied
pulses interacting to produce an echo, a result seen in the fluid code that has not
yet been observed experimentally. Chapter 7 will discuss preliminary research on
a decay instability that occurs at very high amplitudes. This decay instability has
been observed experimentally and also occurs in computations using the nonlinear
fluid code. Chapter 8 will present a summary of results and discuss further research
opportunities.

Chapter 2 will discuss the fluid model used to describe the pure electron plasma.
Under certain approximations, a pure electron plasma can be modeled by the equa-
tions for an incompressible and inviscid fluid. The dynamics of the plasma are gov-
erned by Poisson’s equation, the continuity equation, and the equation for Ex B
drift. These equations are then linearized around steady-state quantities to obtain a
more tractable equation for the perturbed quantities, valid when the perturbations are
small. A set of normalized variables will be introduced, that will be used throughout
the remainder of the thesis.

Chapter 3 will discuss the application of numerical techniques to the fluid model.
A nonlinear 2-D cylindrical fluid code will be presented. Numerical solution of Pois-
son’s equation and the time derivative of the density will be discussed, followed by
the method for advancing the density in time. The diagnostics in the code will be
mentioned, and errors in the fluid code will be discussed. The nonlinear 2-D cylindri-

cal fluid code is the primary tool used to obtain the results presented in this thesis. A

10
model based on the linearized equations presented in Chapter 2 will then be discussed.

Chapters 4, 5, 6, and 7 contain results from the nonlinear fluid code. Chapter
4 describes the response of the plasma to a small amplitude applied pulse. The
response in the electric field at the wall is found to decay exponentially in time. This
is consistent with experimental results for small applied pulses. The decay is due
to a collisionless phase mixing process where density perturbations at different radii
convect at different angular velocities. This mixing process is clearly illustrated by
successive images of the perturbed density. An approximate analytical solution to the
linearized fluid equations is then presented, further illustrating the plasma dynamics.

Chapter 5 discusses higher amplitude pulses which cause a modulation of the decay
envelope due to trapping of fluid elements within the excited wave. The experimental
results showing an exponential decay at lower amplitudes and the modulation of the
decay envelope at higher amplitudes are qualitatively reproduced using the nonlinear
fluid code. Images of the perturbed density show regions of plasma trapped by the
wave, causing the modulation of the decay envelope at higher amplitudes. Thus
the explanation of modulation of the decay as being due to trapping is confirmed.
Applying the external pulse slowly excites a wave that has a greatly reduced decay
rate with no apparent modulation in the decay. This suggests a similarity to BGK
modes in neutral plasmas [20].

Chapter 6 describes the phenomena of echoes which is predicted to occur in a
pure electron plasma [21]. An echo can occur when a single pulse that decays in time
is followed by a second pulse of higher azimuthal symmetry that also decays in time.
These two perturbations then undergo a nonlinear interaction to generate a third
echo response at a later time. An approximate theoretical description of echoes is
first presented, followed by the results from the fluid code. The time of the echo and
its amplitude dependence are studied, as well as higher order effects for larger pulses.
These effects have not yet been observed experimentally in a non-neutral plasma.

Chapter 7 discusses the effects of very large pulses on the plasma. Such large
perturbations are found to be subject to a “decay instability.” A large perturbation

with m = 3 azimuthal symmetry is found to decay into a perturbation with m = 2

11
azimuthal symmetry due to an interaction with the plasma at the beat frequency of
the m = 3 and m = 2 waves. The m = 2 growth rate is systematically studied as
a function of m = 3 mode amplitude. The beat-wave resonance damping results are
compared to the experimental results of Mitchell and Driscoll [22].
Chapter 8 presents a summary of the results of this thesis. The summary is

followed by a discussion of questions which remain for future research to answer.

Chapter 2. The Fluid Model

2.1 Introduction

The geometry of the experimental setup for creating and trapping a pure electron
plasma was previously shown in Chapter 1. This chapter presents a mathematical
model for a pure electron plasma which will be used in subsequent chapters.

A pure electron plasma in the discussed geometry and under certain approxi-
mations can be described as a cylindrical column of rotating fluid surrounded by
a conducting cylinder. A strong axial magnetic field confines the plasma radially,
while electrostatic trapping end electrodes confines the plasma axially, as previously
discussed in Chapter 1.

The first section will discuss the approximations made in the fluid equation model.
The plasma may be taken as electrostatic, two-dimensional, and with zero pressure
gradient.

After the discussion of approximations, the fluid model will be discussed. Under
the approximations, the plasma may be treated as an incompressible, inviscid 2-D
fluid. The 2-D nonlinear fluid equations which describe a non-neutral plasma will be
presented. These equations will be used to model a pure electron plasma, and are
analogous to the Euler inviscid fluid equations. The presented set of equations are
the basis for the Nonlinear Fluid Code presented in Chapter 3.

A normalized set of variables will be introduced. These normalized variables will
absorb such quantities as the radius of the wall, magnitude of the magnetic field,
charge of the electron and the permittivity of free space. The use of normalized
variables simplifies the equations, and these new variables will be used throughout
the thesis.

The nonlinear 2-D fluid equations can be linearized for small perturbations around

steady-state quantities. Such a linearization reduces the 2-D fluid equations to a

13
single equation for the perturbed electric potential. This equation is valid for small
amplitude perturbations to the plasma, and forms the basis of the Linearized Fluid

Code presented in the next chapter.

2.2 Model Assumptions and Approximations

2.2.1 Approximations to Maxwell’s Equations

The fields due to the plasma are taken to be entirely electrostatic. Currents produced
by the motion of electrons are assumed to be small, and the induced magnetic fields
due to these currents are taken as negligible compared to the strong external magnetic
field. The magnetic force between moving electrons is taken to be small compared to
the electrostatic force. Electric fields induced inductively by time dependent currents
in the plasma are negligible compared to the electrostatic space charge fields of the
plasma. Additionally, changes in the electric field are assumed to propagate much
faster than any motion of the plasma.

The surrounding cylinder is taken to be perfectly conducting, so that no dissipa-
tion from resistance in the wall will be considered. This corresponds to a boundary
condition for the electric potential at the wall. The wall will be held at zero electric
potential, except during application of an external pulse of voltage along wall sectors
to induce plasma perturbations.

Under these approximations, Maxwell’s Equations reduce to

=> oO e

V:-E = -—n
E0

Sl oO

V:-B = 0

= —=

Vxh = 0

=>lUlUlC CO

VxB = 0

for the electric field E, and the magnetic field B, with n representing the electron

14
density. The quantity e is the magnitude of the electron charge, and €o is the permit-
tivity of free space. The electric field may be written as the gradient of an electric

potential

—>
E=-Vé
so that the electric potential is given by Poisson’s equation

Vo=—
o a
with the boundary condition at the wall discussed above. The magnetic field is a

—>
constant B in the Z-direction.

2.2.2 Approximation of a Two-dimensional Plasma

A large number of experimental results are essentially 2-D in nature, and this thesis
will focus only on 2-D phenomena. A possible reason for the 2-D behavior of these
plasmas can be argued from the axial motion of the electrons along the magnetic
field lines. The electrons rapidly bounce back and forth between the electrostatic
trapping electrodes at each end of the cylinder, with a bounce frequency higher than
the plasma rotation frequency, effectively averaging over any 2-direction asymmetries.

For the rest of this thesis, only two-dimensional effects will be considered. The
z-direction will be ignored, treating the plasma column as infinite in length, and only

components in r and @ will be considered.

2.2.3 Approximation of Zero Pressure Gradient

The electron density is assumed to have a monotonically decreasing radial distribution
with the density highest at the center and approaching zero near the wall. Related to
the radial density gradient in a warm plasma is a radial pressure gradient which will
cause a diamagnetic drift. Previous research [13] has found that the diamagnetic drifts

and thermal effects appear to be of minimal importance in describing phenomena

15
observed in non-neutral plasmas. Accordingly, all thermal effects will be ignored in

this thesis, and the plasma will be taken to have a zero pressure gradient.

2.2.4 Drift Approximation

This section discusses the drift approximation. Neglecting viscosity and thermal

effects, the equation for the fluid velocity is given by

<= —ne(E +0 x B) (2.1)

where me, is the mass of the electron, W'(r, 0, t) is the fluid velocity, and 2 denotes
the convective derivative. This equation may be rewritten as
Dw e

>» ~
=-— FE -w.v x2
Dt Me

where the cyclotron frequency w, = ek and the magnetic field B = Be. For low

density electron plasmas in a strong magnetic field, the electric field E from the
space charge of the plasma, and the v x B term are large compared to the inertial

Dw
Dt

term. The low frequency modes studied in this thesis are of the same order as

the fluid angular rotation frequency wo, which is assumed to be much less than the

>, term is of the order wov, which is therefore

cyclotron frequency. The inertial
small compared to w.v in the magnetic term, and can be neglected. By neglecting

the inertial term in equation 2.1, the fluid velocity is given by

+ ExB
3/2
which is the standard expression for the velocity due to E x B arift.
Applying equation 2.2 to a cylindrical pure electron plasma, the strong radial
electric field causes a drift velocity in the 6-direction, with central rotational angular

frequency wo(0) = =) where w? = ne is the plasma frequency (A discussion of the

steady-state angular velocity wo(r) will be left to the next section). Therefore the

inertial term in equation 2.1 is of order “ compared to the magnetic term. If the
plasma frequency is much smaller than the cyclotron frequency; i.e., “3 <1, then at
every radius the rotational angular frequency wo is small compared to the cyclotron
frequency, satisfying the above assumptions. In fact, the criterion “ < 1 is well
satisfied for a low density electron plasma in a strong magnetic field, as seen in table

1.1 for the Caltech Pure Electron Plasma Experiment.

2.3 Steady-State Angular Velocity

This section describes the steady-state angular velocity for an unperturbed pure elec-
tron plasma. The angular velocity of these plasmas is a very important quantity
responsible for many plasma phenomena. For a nonuniform density, the angular ve-
locity is sheared. The shear in the angular velocity is responsible for the collisionless
damping described in Chapter 4.

For a confined, unperturbed cylindrical pure electron plasma, the electron density

n(r) is only a function of the radius. The radial electric field is given by

E,(r)=-22 [ "rin (r!)dr!, (2.3)

eo r

and the plasma motion is dominated by E x B drift. Inserting equation 2.3 into
2.2, the azimuthal velocity is given by
el [",

ve(r) = =rB Jp r'n(r’)dr’

where B is the magnitude of the uniform magnetic field that is in the 2-direction.
The angular rotation frequency wo(r) = ve(r)/r can be written as
e 1 f["

wo(r) = a rB F r'n(r’)dr'. (2.4)

For a plasma with uniform electron density that is independent of radius, equation

2.4 can be integrated to give an angular rotation frequency which is independent of

17
radius. Thus, in this special case, the plasma rotates as a rigid body with angular

frequency

len
Wo = -—=.
7 26,B

On the other hand, for electron density profiles that do vary with radius, the
angular rotation frequency becomes dependent on the radius and can be calculated
by equation 2.4. This thesis will focus on non-uniform density profiles that decrease
monotonically with radius. For these density profiles, the plasma flow is sheared. The

angular velocity at the origin is given by

2.4 Nonlinear Model

The plasma is described in two dimensions by the electron density n(r, 0, t), the elec-
trostatic potential ¢(r,@,t), and the fluid velocity w’(r,0,t). The 2 dependence is
ignored. Under the assumptions of the previous section, the instantaneous relation-

ship between the density and electric potential is given by Poisson’s equation

V9(r; 0,t) = —n(r, 6,2) (2.5)

where e is the magnitude of the electron charge, and € is the permittivity of free
space.
Under the approximations discussed in the previous section, the fluid velocity is

given by

9 xB
By

The final equation for this system is the continuity equation for 2-D inviscid fluids,

(2.6)

written as

Variable Name | Normalized Variable Symbol | Dimensionless Variable
time t sot
radius r Ta
electron density n m(0)
electric potential o) EOWA
electric field E ono (OVFa E

Table 2.1: Description of the normalized variables used in the numerical calculations.

= +V-(nv’) =0.
Because the drift form of the velocity is divergence free, the fluid is incompressible,

and the continuity equation reduces to

— +0 -Vn=0. (2.7)

Equations 2.5, 2.6, and 2.7 describe the full nonlinear time evolution of the plasma
under the assumptions of the previous section. These equations are isomorphic to
the Euler equations for incompressible, inviscid fluids [9][10], where the vorticity is
analogous to the electron density, and the electric potential is analogous to the stream

function.

2.5 Normalized Variables

The three equations governing the full nonlinear dynamics of the plasma are Poisson’s
equation 2.5, the incompressible fluid continuity equation 2.7, and the drift equation
2.6. To simplify the calculations and to make the results of the fluid code easily
scalable, a normalized set of variables will be introduced which eliminate the physical
constants from the equations. The definition of these new normalized variables are
summarized in table 2.1. The quantity r,, is the wall radius, B is the magnitude of

the magnetic field, and no(0) is the central electron density.

Using these dimensionless variables, the fluid equations 2.5, 2.6 and 2.7 reduce to

Vo = n (2.8)
On =
— = -v-Vn (2.9)
ee xe
vo = -VOXxzZ. (2.10)

2.6 Linearized Model

It is useful to analyze the case in which the plasma deviates slightly from a cylindrical
equilibrium [10][13][14][15][16]. Separating the electron density, electric potential, and

fluid velocity into steady-state quantities and perturbed quantities gives

n(r,0,t) = no(r) + ni(r, 6, t)

dlr, 0, t) = bo(r) + dy(r, 6, t)

U(r, 0,t) = U(r) + Wr, 6, t)

where the subscript 0 indicates large steady-state zero order quantities and the sub-
script 1 represents perturbed quantities that are first order in amplitude. In a later
chapter describing echoes, second order quantities will be considered to produce a
source term. In this discussion, quantities higher than first order are neglected.
Poisson’s equation is linear, so it can be separated into unperturbed and perturbed

components
V 7
Po = No

and

V7¢, = 71. (2.11)

For the drift equation 2.10, the magnetic field is in the Z direction and V¢@p is in
the 7 direction, so the only unperturbed velocity component is in the 0 direction, and

is given by

u(r lo = Se,

The steady-state angular velocity is given by

1 10¢o
wo(r) = —vo9 = Or”
The continuity equation is simply
Ono
—~ =9
ot

as this is a steady-state quantity, and for the perturbed quantities

On,

apt (Wi + G)- V (no + m1) = 0.

Neglecting the higher order term v7 - Vn, the perturbed continuity equation

becomes

a Un + Vm =0.

As only the t) component of % and the 7 component of Vino is nonzero, the

perturbed continuity equation reduces to

On, Voe On, dno

a or 00 ar

where v4, = — 126 The term vog/r is equal to the steady-state angular velocity

21
wo(r). Making this substitution, the continuity equation for the perturbed quantities
becomes
On, On, dno

ae tl" ag + MG =O (2.12)

The next step is to assume an angular dependence of the perturbation. Complex
exponentials will be used to represent this angular dependence. Since the perturbed
density and perturbed potential are real quantities, it will be necessary eventually
to either add the complex conjugate to the results, or to take the real part of the
quantity to obtain real values. Taking an angular dependence of e’”® for the perturbed
quantities, ny(r,0,t) = nam(r, te”, 1(7,0,t) = o,(r, te" and [v;(r,6,t)|, =

[v1.m(r, t)]-e"? the continuity equation 2.12 for the perturbed quantities is

Oni m(r,t)
Ot

dno(r)
dr

+ imwo(r)nim(7, t) + [vim(7, t)]r = 0. (2.13)

A Laplace transform is applied to the equation 2.13. The Laplace transform
variable s = —iw, where w is complex and includes a small positive imaginary part.
The variable s lies in the right hand plane, so w is in the upper half plane. The
plasma is initially unperturbed, leading to the initial condition n1,,(r,0) = 0. The
transformed equation may be solved for the transformed perturbed density with the
solution

m dro(r) Prm(T, ¥): (2.14)

Mm(r,W) = r(mwo(r)—w) dr

Using this form of the perturbed density, and substituting back into the perturbed

form for Poisson’s equation 2.11 gives the result

do; (r, w) 1 dd, (r, w) m m dno(r) 7 _
dr? dr 7 Le + r(muwo(r) —w) dr JOrm(rw) = 0. (2.1)

This is an ordinary linear differential equation that can be solved for the perturbed

potential Drm (7; w). The boundary condition at the origin is given by 1 m(0, w) =0,

22
and the boundary condition at the wall is given by the Laplace transform of the wall
potential. Once the solution for Pi ml? w) is found, the perturbed density is given by
Poisson’s equation or by equation 2.14. The time dependence of the perturbed density
and perturbed potential is given by the inverse Laplace transform of 71 ,.,(r,w) and
bi m(T5 w). This is the equation used in the linearized computation to obtain results

for small amplitude applied pulses.

2.7 Summary

A pure electron plasma under certain approximations can be treated as an incom-
pressible, inviscid fluid. Under these approximations, the equations of motion are
given by Poisson’s equation, the continuity equation, and the equation for Ex B
drift. These equations can be linearized to obtain a single equation for the perturbed
potential. A solution to this linearized equation should be valid for small ampli-
tude applied pulses, where the resulting perturbed potential is much smaller than the
unperturbed electric potential of the plasma.

These are the model equations used throughout this thesis. The next chapter will

discuss the methods for applying numerical techniques to these equations.

Chapter 3 Computational Methods

3.1 Introduction

The fluid model equations discussed in the last chapter are difficult to solve analyti-
cally for anything but the most trivial choices of plasma density and wall potential.
However, Corngold has made substantial progress for certain power law density pro-
files [15]. In order to allow exploration of a variety of different nontrivial phenomena,
computational methods were used to solve the 2-D fluid equations. Two different
types of codes were developed to solve the equations, each applicable to a particular
regime of plasma behavior.

The first code solves the full 2-D nonlinear fluid equations in cylindrical geometry.
This code will be referred to as the Nonlinear Fluid Code throughout the remainder of
this thesis. The Nonlinear Fluid Code presented here is a new tool for studying pure
electron plasmas. The primary results of this thesis are based on results obtained
from the Nonlinear Fluid Code which was written in Microsoft PowerFORTRAN.

The second code numerically integrates the ordinary differential equation resulting
from the linearized fluid equations (equation 2.15) for infinitesimal perturbations
in a single azimuthal symmetry eigenmode number and as a function of frequency.
Throughout this thesis, the program which solves equation 2.15 will be referred to as
the Linearized Fluid Code. This code was written in Microsoft Quickbasic. Similar
codes to solve equation 2.15 have been used successfully to study pure electron plasmas
[14] [16].

Results from these two codes were compared for small amplitude perturbations.
The results from the Nonlinear Fluid Code agree with the results from the Linearized
Fluid Code. This provides evidence that the computational methods are valid, since
each code uses a quite different method to solve the equations.

This chapter will first present a discussion of the various components of the Non-

linear Fluid Code. The density and potential are discretized on a cylindrical mesh,
and numerical techniques are used to solve Poisson’s equation and the nonlinear equa-
tion for the time derivative of the density. The density is then numerically integrated
forward in time. The diagnostics from the nonlinear code are discussed, including the
effects of accumulated errors.

The next section will describe the Linearized Fluid Code. This code numerically
integrates equation 2.15 for 500 points in radius and 512 discrete frequencies, leading

to a solution for the perturbed potential as a function of time.

3.2 Nonlinear Fluid Code

This section will describe the two-dimensional cylindrical fluid code used to obtain
most of the results in this thesis. Originally this code was written in Cartesian
coordinates by R. W. Gould and produced some useful results. However, it was
found that the Cartesian code generated some numerical artifacts, caused by the
Cartesian grid failing to properly represent the cylindrical symmetry of the problem.
By rewriting the code in cylindrical geometry, these numerical artifacts were avoided.
A grid resolution of 512 points in radius and 128 points in angle was used for most
computations. For several computations, a grid resolution of 1024 points in radius
and 256 points in angle was used, but this increases the necessary computation time
by a factor of 8. This section will focus on the methods used in the cylindrical code.
The specific numerical methods for solving each part of the system will be discussed,
with the final section describing the diagnostic outputs of the code and numerical
errors in the code.

The nonlinear fluid code has three main subroutines. These three subroutines
solve Poisson’s equation, calculate the time derivative of the density, and advance
the density in time. Given a density input, the Poisson solver computes the electric
potential. Given the density and electric potential, the density derivative routine
computes the time derivative of the density. Given the density, time derivative of the

density, and a time step, the time advance routine advances the density in time. For

25
computational stability, the density is advanced using a half-step, full-step leapfrog
time integration method. The program logic is summarized in the flow chart diagram
of figure 3.1. The details of the individual subroutines will be described in the next

sections.

3.2.1 Poisson’s Equation

The first equation to solve is Poisson’s equation (equation 2.8). The following proce-
dure is partially based on the method presented by Birdsall and Langdon [23]. Given
a specific electron density distribution n(r,0), this allows the calculation of the elec-
tric potential ¢(r,6) throughout the plasma, subject to the boundary conditions at
the wall. In cylindrical coordinates, Poisson’s equation is given by

2 2
Oo 106 18H _

Or? rr | 2 Aer

The electric potential is decomposed into two parts, the electric potential from the
plasma, ¢p1asma, and the electric potential from the applied sectors in vacuum without
the plasma, @appieq- The total potential is then just the sum of these two potentials,
® = Pplasma + Papptied- This is done so that the details of the plasma potential can be
observed without being hidden by the applied pulse potential. In this manner, the
buildup of the plasma potential can be observed during the application of the applied
pulse. The plasma potential ¢,ja5ma May be solved with the following procedure by
setting the wall boundary conditions to zero, and the applied potential ¢pji:eq May
also be solved with this method by using a density of zero, and only considering the
boundary conditions at the wall. The applied potential is non-zero for only a short
period of time to initiate a perturbation in the plasma.

To solve Poisson’s equation, the first step is to represent the angular part of the
variables as a Fourier series (7,0) = S>¢,,(r)e” and n(r,@) = > mm(r)e"™, with
each azimuthal symmetry eigenmode number represented by the variable m. The

Poisson equation then becomes a set of equations, one for each m value given by

Poisson Solver

(1,9)

Dertvative Routine

on'(r,6

Time Advance

K—

bole

Tramp (1,8)

Poisson Solver

Gren (1,8)

Derrvative Routine

Nera '(£,0)

Time Advance

K— At

n*“"(r,6)

Figure 3.1: Flow chart describing the logic of the nonlinear fluid code. 'The superscript
represents time, with the time step given by At.

temporary.

The subscript temp stands for

Pom, 1d — m?
dr2 r dr r2

Pm = Nm

where ¢,,,(r) and n,(r) are the Fourier series coefficients for the electric potential and
electron density, respectively.

The next step is to discretize the system in the radial direction, setting r; = zAr
where Ar is the radial grid spacing and the subscript 7 represents the radial integer
index for the discretized system. Applying second order finite difference methods the

equation becomes

it+1m ~~ 2 i,m, + i—1,m it+1,m ~~ Yi-1,m m
bisa 2m t Germ , Gorin cum My
(Ar)? 2r;(Ar) re"

which can be rewritten as

1 1 9) m2 1 1

The index 7 has a range of 1 to N;, where N; is the total number of radial grid
points. The origin at index 7 = 0 is a boundary condition to be discussed below.
This system is a set of linear algebraic equations, one for each index 7, that can be
solved by linear algebra matrix techniques. Specifically, the system can be rewritten
in vector form as Ag = 7” where A is a tridiagonal matrix and o and 7’ are vectors
with index 7 ranging 1 to N; —1. The boundary conditions at the origin and at 1 = N;
have been included in 1m and ny, m-. A tridiagonal solve routine, IMSL library call
DLSLTQ, was used to solve for each Pim for the given set of nim. The set of nim is
found by performing a discrete Fourier transform in angle on the density n,;, where
a is the index in radius and 7 is the index in angle.

For the plasma potential ¢yjasma, the boundary condition at the wall is zero for
the electric potential due to the plasma since any potentials at the wall are contained
in the applied potential. The boundary condition is satisfied by setting the Fourier

coefficients dy, » = 0 for all m values. For the applied potential @,,pi:¢q, the boundary

28
condition at the wall for the electric potential from applied voltages at the wall in
vacuum are determined by Fourier analyzing the applied wall potential in angle @ to
obtain the Fourier coefficients ¢y, ,,. Using this wall boundary condition, the vacuum
potential is found by using a value of zero for the density, ie., Nim = 0 for all ¢ and
mM.

The boundary condition at the origin is simply given by ¢,, = 0 for all |m| > 1,
since a well behaved solution to Poisson’s equation can not contain a finite angular
component at the origin. However, the m = 0 case must be treated separately, using
a different boundary condition at the origin.

The origin boundary condition for the m = 0 component can be understood by a
simple application of Gauss’ Law. For the plasma potential ¢yjagma, the origin m = 0
density component noo represents the density in an area of a circle with a radius of
sAr. Therefore, if we draw a Gaussian surface of a circle around this disk with a
radius of }Ar, the enclosed charge is given by —{7o0(Ar)?. This leads to a radial
electric field at the radius r = $Ar, of E,($Ar) = —}nooAr. This electric field

at r = 3Ar gives rise to a boundary condition relating to the radial derivative of

the electric potential at r = sAr, which can be computed using the points at ¢ 9

and $19. Using the finite difference form of the equation E,($Ar) = — oe r=iAr? We
obtain the result — £1,0- 20.0 = —4nooAr, or do — b19 = —4N00(Ar)?. This is the

boundary condition used at the origin for the m = 0 component of the system for the
plasma potential. For the applied potential ¢gppiieq, the m = 0 boundary condition
at the origin is zero.

After solving for the set of ¢;,,, for each value of m, iterative refinement is per-
formed to remove errors accumulated in the tridiagonal solve routine. Iterative re-
finement takes the set of @;,,, for each m and multiplies it by the tridiagonal matrix.
This should give back the set of nim. However, due to numerical accumulation of
errors, this new set of 7; differs from the original by a small difference. This small
difference An, is recorded and reapplied to the tridiagonal solve routine, resulting
in a set of small corrections to the potential, A¢,,,,. Subtracting off this error from

the original potential set of ¢;,,, gives the refined potential solution that should have

29
near machine accuracy precision. A further discussion of iterative refinement may be
found in [24].

This process for finding the potential solution set of ¢;,, is repeated for each
value of m. After all @;,, for all values of m have been computed, the Fourier series is
summed, obtaining the real space electric potential solution ¢; ; = >> Lime” where
zi is the index over radius and 7 is the angle index which ranges from 1 to N;. The
electric potential at the origin is represented by @ o-

This is the method used to solve Poisson’s equation for the electric potential. It
is second order accurate with respect to Ar. The electric potential found from this

solver is then passed to the density derivative routine, discussed in the next section.

3.2.2 Density Derivative Routine

The equation for the time derivative of the density is a combination of the continuity
equation 2.9 and drift equation 2.10. Inserting equation 2.10 for the drift velocity v”

into the continuity equation 2.9 gives

ot = (Vx 2)-Tn
which can be rewritten as
on = (Vax V9) -2 (3.1)

Writing this equation in cylindrical coordinates gives the result

On _ 19nd _ Onde 2)

ot r\drd00 000r™
One might consider applying standard second order finite difference techniques
to this equation to obtain a value for the time derivative of the density. However,
this approach will lead to computational instability first demonstrated by Phillips
[25][26]. Akio Arakawa recognized the cause of this instability, and devised a method

to compute the spatial derivatives while maintaining computational stability [27].

30
Specifically, Arakawa took four second order finite difference representations for the
spatial derivatives that used neighboring grid points, and found that a linear combina-
tion of these would produce a method that conserved mean density, total energy, and
mean squared density, while maintaining computational stability. Using this method,

the finite difference representation for the density time derivative is

aE = Tar ApAg (bid + Pittj—1 _ Pijt _ Pigs jt1)(Mi41,3 — i,j)

Pi-15-1 + Pi i,j—1 — Ogu — Pi j41) (M4, — Ni-1,3)

1j + Pisrje. ~ Pi-1j ~~ bi-1,j41) (Mig ~ ni.)

+(@

(9:4
+ (Giga j—1 + bist, — Pi-13-1 — Gi-1,5) (Mag — M451)

(9:4

(6,5

15 ~ Pigta)(Mt1,941 — Nay)
-1 7

bi-1,5)(M4,3 — M-1,j-1)
+( bi 541 — O15) (Mi-1,541 — N45)
(ni,

+($i41,5 — Pi 5- -1) — Ni41,5-1)|

where 7 is the index over radius, 7; = i7Ar, and 7 is the index over angle. The radial
grid spacing is represented by Ar, and the angular grid spacing by AQ.

The boundary condition at the wall is straightforward. For the electric potential,
the wall potential is used which is either zero when no external pulse is present, or
equal to the applied voltage at the wall during an applied pulse. For the electron
density, it is assumed that any electrons that hit the wall will be lost, and the density
is therefore zero at the wall. At the origin, the equation must be modified since the
area must be calculated differently, and there is no angular component at the origin.

The derivative for the density at the origin is given by

N;
Ono,o ;

Dt 3n(Arp S I (Puja a 1,3-1) (Noo +715)

j=l
an — ;)(No,0 + 1,541)

+(1,; — %0,0)(,3-1 + N0,0)]

where noo is the density at the origin and ¢o is the electric potential at the origin.
The origin does not have an angular component and therefore does not have an angle
index j.

The system has cylindrical symmetry, and therefore is periodic in angle. To handle
the boundary conditions in angle, the indices are “wrapped” around. For example,
when using the value ¢; y,,1, the value for @;, is used, and when ¢;, 9 is encountered,
the value for ¢; y, is used.

This routine will give a second order accurate time derivative of the density at
each grid point, given the density and electric potential at each point. With this
derivative, the density at each grid point can be advanced in time by the procedure

described in the next section.

3.2.3. Leapfrog Time Advance

With the above two procedures, the electric potential can be calculated given the
electron density distribution, and the density time derivative can be computed using
the electric potential and the electron density. The leapfrog method then allows the
density to be advanced by a time step At.

The leapfrog method begins by calculating the electric potential and density time
derivative from the given density distribution. A new temporary density is calculated

at half a time step. This can be written as

nitah — nt 1 Omi
tempi,g ~ ' "4,7 9 Ot

with the superscripts representing the time variable.

This temporary density is inserted into the Poisson solver and derivative algo-
rithm, resulting in a temporary potential ¢jemp,, aud derivative OnterPid at a time of
t+ 5A. This derivative is then used to advance the original density a full time step,

written as

nitadt
t+At ot tempi, j

‘J at

where the superscripts ¢,¢+ At, t+ At represent the original time, time advanced
half a time step, and time advanced a full time step, respectively.

This procedure is then repeated to advance the density at each grid point in time.
The routine is second order accurate with respect to the time step At. For the
procedure to be numerically stable, the Courant condition requires that the time step
be sufficiently small [24]. The maximum stable time step Atmax must satisfy Atmax S
ar and Atmax S ae for all r. For the typical resolution of 1024 points in radius and
256 points in angle, Atmax was experimentally determined to be around 1 x 107%t¢,

_ 4neoB

where t,, = eno(0)

of At somewhat below Atma, were used for the time step. For a resolution of 1024

is the time for one central rotation. For most computations, values

points in radius and 256 points in angle, a value of At = 5 x 10~‘t,, was used.

This now completes the system necessary to advance the density as a function
of time and observe the time evolution of the system. The next section discusses
diagnostics and output data in the code so that the time evolution of the electric

potential and electron density can be conveniently displayed and analyzed.

3.2.4 Diagnostics and Outputs

There are several computer data files produced by the code at periodic intervals during
a computation, with this interval denoted as a frame. The time step is too small for
data to be stored at every time step. Therefore successive frames are separated by
hundreds of time steps, with the perturbation rotating approximately ; of a revolution
between frames. These files contain data describing the density, perturbed density,
perturbed potential, Fourier analyzed density and potential, electric field at the wall,
and conserved quantities. These output files will be described below.

The plasma density is stored internally as a double precision real value for each

grid point throughout the computation. When this file is stored to disk, each grid

point is scaled and rounded into a one byte integer. Thus the density is stored with 8
bit resolution, and the value for the density varies between zero and 255 at each grid
point. This scaling is necessary due to the limited disk capacity, where even at 8 bit
resolution, for a typical run of 200 frames at 1024 points in radius and 256 points in
angle, the total storage requirement is 50 megabytes. Storing the full double precision
values at each grid point would require 400 megabytes for this single file (significant
in the year 1997). This byte resolution storage method is also used for the perturbed
density and perturbed potential.

Often the perturbations are small, and thus no changes are perceptible when
observing the total density. For this reason, the perturbed density is also stored. The
perturbed density is calculated by subtracting off the unperturbed density at each
radius from the total density. This leaves the perturbed density which varies from
the background average density. The perturbed density is then scaled so that the
perturbation has a value ranging from -127 to +127. The initial scaling factor for
the perturbed density is not known ahead of time. If too large a scale factor is used,
the perturbed density will exceed the limit of 127, and thus the stored file will be
inaccurate. If too small a scale factor is used, the full 8 bit dynamic range will not be
exploited. Therefore, a variable scale factor is used. The scale factor is normalized so
the highest absolute perturbation has a value of 127 at each instant of time during
the application of the applied pulse. After the applied pulse is turned off, the scale
factor remains fixed for the remainder of the computation. The value of the maximum
absolute perturbation is stored in a separate file, so the scale factor is always known.
As the maximum perturbations should not grow after the applied pulse is turned off,
this method maximizes the use of the 8 bit dynamic range without saturation. The
perturbed density is then stored in byte resolution.

Usually, variations can not be seen in the total potential either. Therefore the
same method used for the perturbed density is also used for the perturbed potential,
where the average potential is subtracted off, and the difference scaled to byte reso-
lution. Again the scale factor varies during the applied pulse, and then remains fixed

throughout the rest of the computation.

The exact rotation frequency of the perturbation is not known before a calculation
is performed. The angular frequency profile wo(r) can be computed, but the resonant
frequency of the wave remains unknown. However, a reasonable guess can be made
which will be used to calculate the time between frames and the frequency of the ini-
tial applied pulse. After the calculation is finished, the actual perturbation rotation
frequency can be determined and used to calculate an accurate time between frames
and used for the applied pulse frequency in subsequent calculations. The perturba-
tion rotation frequency is less than the central rotation frequency for monotonically
decreasing profiles. Thus the time between frames may correspond to ~ 0.3t,, for
a particular density profile, and a total run of 200 frames would correspond to an
elapsed time of 60 central rotation periods for this density profile.

After the density, perturbed density, and perturbed potential files have been com-
puted and stored, a conversion program is used to produce 2-D contour plots of these
variables. The converter program takes the cylindrical input data, and converts it
into rectangular coordinates with 480 by 480 resolution using bilinear interpolation.
This resolution is used since the typical resolution for a VGA screen is 640 by 480.
This file can then be read by an animation program which will display the files frame
by frame in rapid succession, showing smooth flow of plasma in time. The animation
program divides the data by 16, giving 16 contours. In the total density picture,
the white regions represent high density, and the black regions low density. For the
perturbed pictures, the white regions represent areas of high positive perturbation,
the middle grey transition represents a contour of zero perturbation, and the black
regions represent an area of high negative perturbation.

The sin(m@) and cos(m#@) components of the density and potential at each radius
are stored for the values m = 0,1,2,3,4 at each frame. The components are stored
in full 8 byte double precision. These files allow analysis of the Fourier component
variation in both amplitude and phase as a function of radius and time.

Additionally, the electric field at the wall is written to a file in full 8 byte double
precision. The electric field is computed by taking the negative of the derivative of

the potential at the wall. A second order finite difference is used to approximate the

35
derivative, using the neighboring grid points to the wall. The electric field is stored
for each grid point in angle. The time between storing the electric field is TH the time
between storing frames. Thus in the time between each density picture, the electric
field has been stored 10 times.

Finally, for every frame a file is written which contains several calculated conserved
quantities. The finite difference method used in the code should ideally conserve mean
density, mean squared density, and total energy. These values are computed over the
area of the plasma and written to a text file. It is found that for a small amplitude
perturbation, the mean density and mean squared density are found to vary by less
than 1 part in 10°. For a very large amplitude perturbation, the mean density may
vary by 1 part in 10’, and the mean squared density may vary by 1 part in 10°. At
very large amplitudes, some plasma may hit the wall and be lost.

There is one more feature of the code worth mentioning here. It is often useful to
observe the plasma dynamics in the rotating frame that is moving with the plasma
perturbation. For the system of equations we are studying, moving into the rotat-
ing frame is equivalent to including a uniform, time independent, positively charged
background density in the calculation. This can be seen by using the result that
the rigid rotor solution to the fluid equations has a uniform density distribution as
discussed in Chapter 1. The rotation frequency of the rigid rotor is proportional to
the value of the uniform density. This value for the background positively charged
density is placed into the Poisson solver, where the density used in the solver is the
electron density minus the uniform positively charged background density. This will
then remove the rotation of the wave, and allow the plasma dynamics to be observed
in the rotating frame.

The value of the uniform positive background density is chosen to remove the
rotation of the perturbed potential. This wave angular frequency value may not be
known ahead of time, as discussed earlier. Therefore, a feedback loop is used to track
the perturbed potential position in angle, and adjust the value of the uniform positive
background density to keep the perturbation from rotating. This adjustment is made

at each frame, and the new perturbation frequency is written to a file. This tracking

method is also useful for determining if the wave frequency changes with time.

3.2.5 Accumulated Errors

There are several sources of errors in the Nonlinear Fluid Code. First, all spatial
derivatives are second order accurate in Ar and A@. In time, higher order error
terms accumulate, producing fine scale numerical structure. The error becomes more
significant when spatial features of the plasma approach the dimensions of the grid
resolution. The grid can not accurately represent features that are smaller than the
grid spacing.

Additionally, errors exist in the time integration method. The leapfrog time step
method discussed above is second order accurate in At. These errors also accumulate
in time, and limit the accuracy of the code for later times.

The nonlinear fluid code also suffers from aliasing errors [27]. These errors manifest
themselves as phase inaccuracies and distortions in the spectral distribution of energy.
However, the total energy and average scale of motion do not suffer from these aliasing
errors.

Due to these accumulated errors, the results from the code become less accurate
for long integration times. These errors may be observed in images of the perturbed
density as small points of perturbed density with dimensions on the order of the
grid resolution. As the computation proceeds, this high spatial frequency component
becomes more evident, and contours of constant perturbed density become “fuzzy”
due to accumulated numerical noise.

Therefore, the computation is valid for a limited range in time. ‘Typically, a
computation is stopped at 200 frames, corresponding to 60 — 100 central rotation
periods. At 200 frames, the effects from numerical noise become noticeable.

It is possible to reduce the accumulated errors by increasing the grid resolution
and decreasing the time step. Increasing the radial and angular resolution by a factor
of 2 requires decreasing the time step by a factor of 2, and the required computation

time increases by a factor of 8. Since increasing the grid resolution substantially

increases computation time, it is not always practical to use higher resolutions.

3.3 Linearized Fluid Code

The final linearized equation from the previous chapter, equation 2.15, is a differential
equation for the perturbed electric potential bi m(T, w) as a function of the radius r
and angular frequency w, with an assumed angular dependence of e””® where m is the
azimuthal symmetry eigenmode number. The plasma density is then related to the
potential through the Poisson equation. To obtain the time evolution of the system, it
is necessary to solve bi m(T, w) for various frequencies, and perform an inverse Laplace
transform. In the linear equations, different m modes are clearly independent of each
other. Therefore, for simplicity a single m mode is studied.

A range of w is chosen, and this frequency range is then divided into 512 discrete
frequencies. This range is selected to include the azimuthal symmetry eigenmode
number m times the fluid rotational velocity wo(r) at each radii and neighboring fre-
quencies at which the plasma has a significant response. For each discrete frequency,
the equation is solved for 500 points in radius using a fourth order Runge-Kutta
method. There is a pole in the equation when w = mwo(r), and numerical difficulties
with this pole are avoided by including a small positive imaginary component to the
frequency.

Boundary conditions for the transformed potential are specified at the origin and
at the wall for the integration. The density perturbation is expected to be very small
at the origin, and therefore the Laplace equation vacuum solution bamlt w) = Ar™
is used for the first few integration points, where A is an unknown constant which is
initially assumed to be one. Equation 2.15 is then integrated outward to the wall via
the Runge-Kutta method for a discrete set of 512 angular frequencies. The boundary
condition at the wall is set by assuming a potential perturbation which is described
by the discrete transform of a delta function in time at t = 0. The delta function
solution can be treated as a Green’s function in order to construct a solution for an

applied pulse with a finite time extent. A delta function impulse in time corresponds

to a flat frequency spectrum with all frequency modes of the applied potential at the
wall equal to one. When the Runge-Kutta scheme reaches the wall, the potential in
general will not be equal to one. To remedy this, the solution is then renormalized
by dividing the potential at each radius by the final wall potential. This procedure is
valid because the equations are linear. This then gives the final solution for Dt nl"; w)
with the correct boundary conditions at the origin and at the wall. This solution
procedure is applied to all 512 discrete frequencies.

The discrete numerical solution 1 m(T;W) is the approximate delta function im-
pulse response of the system to a potential applied at the wall. This solution is
“windowed” by multiplying the solution by a function that approaches zero near the
boundary frequencies to remove aliasing errors due to the finite extent of the dis-
crete frequency range. A Gaussian in frequency which rolls off quickly near the edge
frequencies is used to window the solution Pt ml? w). This Gaussian function is mul-
tiplied by the solution bi m(T; w) for all radii at each specific angular frequency. This
is equivalent to using a burst with a short Gaussian envelope to excite the perturba-
tion, instead of a delta function. Again this is valid due to the fact that the system
is linear.

To find the time dependence, an inverse Laplace transform of Pt m(T; W) is re-
quired. However, the function to be transformed has only been computed at 512
discrete frequencies. We would like to perform this inverse transform using discrete
fast Fourier transform (DFT) methods, except that the frequencies used in evaluating
the function are not real, but are complex and contain a small imaginary part. If a
DFT, which assumes data at real frequencies were to be used, an artificial exponential
decay would be introduced. However, as the argument below shows, this is readily
corrected.

To correct for this artificial decay, the solution ¢, ,,(r,t) is multiplied by e”‘* where
w; is the constant imaginary part of the angular frequency. This correction may be

understood from the definition of the inverse Laplace transform

fit) =5- [Flue aw

where w is in the upper half plane, i.e., w = w, + iw;. This can be written

1 oO . .
f(tj= a / F (Wp + iu, eer duo,

Since the imaginary part of the frequency w; is a constant, the exponential can be

taken outside of the integral giving

1 oo ,
f(t) = ae I. f (wp + iwi )e"dw,.

This is a routine to perform the inverse transform using real values for w even
though the function f(w) has been evaluated by including an imaginary part for w.
Therefore DFT routines can be used, provided that a correction is included after the
inverse transform has been completed.

Now that the time dependent solution of the potential has been found, the time
dependent electric field at the wall is equal to the negative of the derivative of the
electric potential at the wall. Additionally, the time dependence of the perturbed
density is given by the Laplacian of the time dependent potential through Poisson’s
equation 2.8, and is computed in the code with the electric potential. Images of the

perturbed density are generated by including the angular dependence, and dividing

the density into 16 equally spaced contours.

3.4 Conclusion

Both the Nonlinear Fluid Code and the Linearized Fluid Code are useful for studying
pure electron plasmas. The Nonlinear Fluid Code is useful for studying higher ampli-
tude effects where nonlinear terms are important. The required computational time

varies, depending on the complexity of the case running, but a typical calculation

AO
requires one week of computation on a 133 MHz Pentium PC. The Linearized Fluid
Code is convenient for studying a variety of small amplitude effects quickly, since the
computation requires only 20 minutes on a 66 MHz 486 PC.

For small amplitudes, the results from the Nonlinear Fluid Code are found to
agree with the results of the Linearized Fluid Code. This is a useful diagnostic, as
each code uses a different method for solving the equations. Obtaining the same
results with each method for low amplitude applied pulses provides evidence that the
computational methods are valid.

These two computational methods are used to study a variety of aspects in pure
electron plasmas. The results from these computations are the subject of the next

several chapters.

Chapter 4 Small Amplitude

Disturbances

4.1 Introduction

Before attempting to solve large amplitude nonlinear problems, the Nonlinear Fluid
Code was first used to study small amplitude disturbances in order to verify that the
well understood results from linear calculations could be reproduced using the Non-
linear Fluid Code. For small amplitude perturbations, the results from the Nonlinear
Fluid Code were found to agree with the results from the Linearized Fluid Code.

The linearized differential equation approach has been used many times in previous
research on non-neutral plasmas. In particular, computational work by Pillai [19]
using a Linearized Fluid Code found that, for a small amplitude applied pulse, the
response of the plasma at the wall decays exponentially for early times, and transitions
to an algebraic decay at later times. He also found that by choosing a plausible density
profile, he could match the experimentally observed decay rates. However, there is
no reason to believe that the density profile he chose is the only profile that has a
decay rate which matches the experimental decay rate.

The first section will describe perturbations in the electric field at the wall resulting
from application of an external pulse. The perturbed electric field at the wall is found
to decay exponentially in time. Several density profiles are considered, and the decay
rate is found to be strongly dependent on the radial density profile. The decay rate is
studied as a function of the density profile, and compared to the theoretical prediction
of Briggs, Daugherty, and Levy [10].

The second section will describe the density perturbations caused by the external
pulse. The perturbed density images illustrate the sheared flow of the plasma, and

explain the collisionless damping as a phase mixing process caused by the different

42
radial layers rotating at different rates. Images of the perturbed density will be
compared for plasmas with different radial density profiles and different decay rates.
Finally an approximate analytic solution to the linearized fluid equations will be
presented and compared to the results from the computation. The approximate solu-
tion models the perturbed potential as a decaying exponential function, and models
the perturbed density as a sum of two exponentials. Using parameters obtained from

the computation, the images of the perturbed density can be reproduced.

4.2 Collisionless Damping

4.2.1 Fast Decay Profile

The Nonlinear Fluid Code is used to study the response of a pure electron plasma to
an externally applied pulse. The response of the plasma is observed in the perturbed
electric field at the wall. After application of the external pulse, the perturbed electric
field at the wall is found to decay exponentially in time. The decay rate is a strong
function of the unperturbed density profile. These results are consistent with those
found by N. Sateesh Pillai [19].

The Nonlinear Fluid Code is used to study the unperturbed normalized density
profile no(r) = no(0){1 — 3r? + 3r4 — r6} for r < 1. This normalized density profile
has a value of one at the origin, with the first derivative vanishing at the origin.
Furthermore, the density is monotonically decreasing to zero at the plasma edge at
r = 1, with two derivatives vanishing at the plasma edge. The density profile may
also incorporate a parameter @ which is the radius of the edge of the plasma, beyond
which there is vacuum until the wall. For values of a < 1, the normalized density
profile is modified to become no(r) = no(0){1 — 3(£)? + 3(£)* — (£)®} for r < a, and
no(r) = 0 for r >a.

Each density profile has a corresponding radial electric field due to the electrostatic
space charge of the plasma. The radial electric field causes an E x B drift velocity

in the plasma perpendicular to the electric field, creating a rotation of the plasma.

43
This drift velocity at each radius corresponds to an angular velocity at that radius.
For a plasma with uniform density, the angular velocity is constant in radius, and the
plasma rotates as a rigid rotor. For a monotonically decreasing profile as presented
above, the angular velocity is different at different radii, and the flow is sheared.
For the normalized density profile presented above, the normalized angular velocity
profile is given by wo(r) = wo(0){1 — 3r? + r4 — 4r®} where wo(0) = 2) = pee is
the central angular rotation frequency. For values of a < 1, the normalized angular
velocity profile is modified to become wo(r) = wo(0){1 — $(£)? + (£)* — Z(£)°} for
r a. For the computation in this section, the
value a = 0.8. is chosen. A plot of the normalized density profile and normalized

angular velocity is shown in figure 4.1.

Unperturbed Density and Angular Velocity vs. Radius

angular velocity

resonant radius

Normalized Profiles

0.2 o.4 0.6 0.8 1

Normalized Radius

Figure 4.1: Plot of the normalized unperturbed density and normalized angular ve-
locity versus radius for the Fast Decay Profile no(r) = no(0){1—3(£)? +3(£)* — (£)°%}
with the edge of the plasma at a = 0.8. The m = 2 resonant radius at r, = 0.72r,, is

indicated by the dashed vertical line.

A perturbation in the plasma will propagate with a certain phase velocity. For a

monotonically decreasing density profile, at a single radius, the E x B drift velocity
of the plasma may be equal to the phase velocity of the wave. This radius is denoted
as the resonant radius r,. For the density profile in figure 4.1, the m = 2 resonant
radius is indicated by a dashed vertical line.

The computation begins with a three cycle sinusoidal pulse of voltage applied to
sectors along the wall to excite a perturbation in the plasma with m = 2 azimuthal
symmetry. The configuration of the sectors is modeled after the experiment, which
was shown in figure 1.1(b) of Chapter 1. The pulse is applied at the m = 2 resonant
angular frequency, which is found to be w, = 0.61w9(0) for this profile. Because
the resonant angular frequency is not known before the computation is performed,
each computation is performed with an initial guess of the resonant frequency. From
the response of the initial computation, a resonant frequency is calculated and used
for the applied pulse frequency in subsequent calculations. The resonant angular
frequency of w, = 0.61w9(0) corresponds to a resonant radius at r, = 0.72r,, for this
profile, where r,, is the wall radius. At this radius, the plasma has an unperturbed
density of 0.0060n9(0). The amplitude of the applied pulse is 5 x 107°¢) where gy is
the electric potential at the origin. For this profile ¢y = —0.101185£n9(0)ri,.

Figure 4.2 shows a plot of the perturbed electric field at the wall as a function
of time. Time is presented in units of the central rotation period t,, = OR The
perturbed field oscillates as the plasma density perturbations rotate past the wall
sectors. In this plot, you can see the perturbation increasing in amplitude, approxi-
mately linearly, for the first three cycles while the applied pulse is present. After the
applied pulse has terminated, a decay of the perturbation, which is approximately
exponential, is observed. The perturbed electric field reaches a maximum of 0.0086%
of the unperturbed electric field at the wall for this three cycle applied pulse. For this
profile, the unperturbed electric field at the wall has the value —0.08279(0) rw.

Figure 4.3 is a plot of the envelope of the perturbed electric field at the wall
versus time on a semilog scale, separated into angular Fourier components. This plot
more Clearly illuminates the exponential decay of the perturbation, and also shows a

transition to an algebraic decay later in time [10][19][28], with interference between

0.01 T T T T T T T T T T T T T T
w= 0.610,(0)

co 606f ¥ = 0.0500,(0) |
o | wrly = 1241
Nene”
us) |
Lv r 7
om
2 7
Sm
—_
3 0
fx) L J
ar
2 — +
= L |
cP)
PA nm 4

=- 1 L 1 L | L L 1 L | 1 1 Ll L

0.01, 25 50 75

Time (Central Rotation Periods)

Figure 4.2: Plot of the perturbed electric field at the wall versus time for the Fast
Decay Profile no(r) = no(0){1 — 3(£)? + 3(£)* — (£)°} with a = 0.8. A three cycle
burst of voltage is applied at the wall with an amplitude of 5 x 10~°¢). The response
in the perturbed electric field at the wall reaches a maximum amplitude of 0.0086%
of the unperturbed field at the wall, and then decays exponentially in time.

a= 0.610,(0)
~y = 0,050«,(0)-
fy = 12.1

Log of Envelope of Perturbed Electric Field
(each tick is a factor of 5)

25 50 75

Time (Central Rotation Periods)

Figure 4.3: Semilog plot of the amplitude of the Fourier components of the electric
field at the wall versus time for the Fast Decay Profile no(r) = no(0){1 — 3(£)? +
3(£)* — (£)°} with a = 0.8. The exponential decay of the m = 2 perturbation i is
observed, ‘followed by a transition to an algebraic decay. The perturbation is excited
by a three cycle burst of voltage on wall sectors with amplitude 5 x 10-°¢). The
perturbed m = 2 electric field at the wall reaches a maximum amplitude of 0.0086%
of the unperturbed electric field at the wall. The m = 4 component is a harmonic of
the m = 2 pulse that is generated within the plasma, while the other higher m values
are excited by spatial harmonics of the sectored electrode. Each tick on the vertical
axis represents a factor of 5 in amplitude.

AT
the two contributions during the transition. The initial decay rate for the m = 2
component is found to be 7 = 0.050w9(0) leading to a value of “* = 12.1 for this
density profile.

In addition to the m = 2 mode, higher Fourier components are also observed
in figure 4.3. The components m = 2, m = 6, m = 10, and m = 14 are excited
directly by the sectored electrode structure (shown in figure 1.1(b) of Chapter 1),
which has angular Fourier coefficients of aul? ave 2/2 and ava respectively for this
m = 2 sector configuration. In contrast, the m = 4 Fourier component is not excited
directly from the sectors. Instead, this component is a harmonic of the m = 2 mode
that is generated within the plasma. These higher Fourier components are not of too
much concern for this chapter, but they will become more important in Chapter 6
where plasma wave echoes are discussed.

To demonstrate the effect of a sectored electrode structure, a second computation
was performed using the same profile, but with a pure angular m = 2 electric potential
perturbation rather than a sectored electric potential configuration applied to the wall.
Additionally, a Gaussian pulse shape in time is applied in a frame rotating with the |
wave, rather than using a three cycle sinusoidal pulse. The Gaussian pulse reaches a
maximum amplitude of 1.18 x 10~*¢, and has a full width at half maximum (FWHM)
of 1.06¢., where t,, is the time for one central rotation. The perturbed electric field
at the wall reaches a maximum value of 0.015% of the unperturbed electric field at
the wall.

Figure 4.4 is a semilog plot of the Fourier components of the electric field at the
wall versus time for this applied Gaussian pulse. As is seen from the figure, the initial
exponential decay rate is unchanged for the m = 2 component of the electric field at
the wall. However, since the wall perturbation has a pure m = 2 angular symmetry,
the m = 6, m = 10, and m = 14 components are no longer directly excited from the
applied pulse at the wall. The m = 4 component remains, however, as it is a harmonic
that is generated by the plasma. Higher order harmonics are also generated, but are
too small in amplitude to be visible in figure 4.4.

As we will see in the next section, the decay rate depends strongly on the unper-

a oo \ pte fanindansoustisnuafenmntannniann ac toeein’ = 0.0500,(0).
ON tae = 121

(each tick is a factor of 5)

Log of Envelope of Perturbed Electric Field

0 — — 35 — — 50 — — 75
Time (Central Rotation Periods)

Figure 4.4: Semilog plot of the amplitude of the Fourier components of the electric
field at the wall versus time for the Fast Decay Profile no(r) = no(0){1—3(£)?+3(£)*—
(£)°} with a = 0.8. Each tick on the vertical axis represents a factor of 5 in amplitude.
This calculation uses a pure angular m = 2 excitation at the wall. The applied pulse
has a Gaussian shape, reaching a maximum amplitude of 1.18 x 10~4¢,) with a FWHM
of 1.06¢,,. The m = 2 perturbed electric field reaches a maximum amplitude of 0.015%
of the unperturbed electric field at the wall. The m = 4 component is a harmonic
that is generated within the plasma.

49
turbed density profile. The perturbations of the electric field at the wall decay rather
quickly for the density profile discussed in this section. Accordingly, I will refer to
this density profile, no(r) = no(0){1 — 3(£)? + 3(£)* — (£)°} with a = 0.8 as the Fast
Decay Profile. Since this profile exhibits strong damping, it is useful for studying the
mechanism that causes the decay of perturbations in the electric field at the wall.

The Fast Decay Profile will be the primary density profile used in Chapter 6.

4.2.2 Slow Decay Profile

The previous section looked at the decay rate of the electric field at the wall for the
Fast Decay Profile. This section will examine a density profile that has a much lower
decay rate. The computation in this section uses the same Gaussian pulse shape
excitation in a pure m = 2 electric potential perturbation at the wall, reaching a
maximum amplitude of 1.18 x 10~*¢, and a FWHM of 1.06t¢,.

This section discusses the normalized density profile no(r) = no(0){1 — 6(£)* +
8(£)° — 3(£)§} with a = 0.8. The corresponding normalized angular velocity is
wo(r) = wo(0){1 — 2(£)* + 2(£)® — 2(£)8}. This density profile is flatter near the
origin, with three derivatives vanishing at the origin, and smooth near the edge with
two derivatives vanishing at the plasma edge. A plot of this normalized density
profile and normalized angular velocity profile is shown in figure 4.5. This profile has
a central electric potential of ¢) = —0.142696£no(0)rz,, and an unperturbed electric
field at the wall of —0.1282n0(0)rw.

For this profile, the m = 2 resonant angular frequency is found to be w, =
0.843w (0), which corresponds to a resonant radius of r, = 0.779ry. The plasma
density at this radius is 0.00053n9(0). For the given applied pulse, the electric field
at the wall reaches a maximum perturbation of 0.027% of the unperturbed electric
field at the wall. The initial decay rate of the perturbed electric field at the wall is
found to be y = 0.0089w9(0) which leads to a value of “s = 95. A semilog plot of

the Fourier components of the electric field at the wall versus time for this profile is

shown in figure 4.6.

Unperturbed Normalized Density and Angular Velocity vs. Radius

angular velocity

resonant radius
0.42

Normalized Profiles
ip

Oo
bo
iti ee

ao
HB

O.2 0.4 0.6

Normalized Radius

Figure 4.5: Plot of the normalized unperturbed density and normalized angular veloc-
ity versus readius for the Slow Decay Profile no(r) = no(0){1—6(£)* +8(£)° —3(£)8}
with the edge of the plasma at a = 0.8. The m = 2 resonant radius at r, = 0.779r.,
is indicated by the dashed vertical line.

m2 5 v= 8.008960)

Log of Envelope of Perturbed Electric Field
(each tick is a factor of 5)

| rn
Time (Central Rotation Periods)

Figure 4.6: Fourier components of the electric field at the wall versus time for the Slow
Decay Profile no(r) = no(0){1—6(£)*+8(£)° —3(£)8} with a = 0.8. This calculation
uses a pure m = 2 excitation at the wall. The applied pulse has a Gaussian shape,
reaching a maximum amplitude of 1.18 x 10~‘¢y) with a FWHM of 1.06t,,. The
perturbed m = 2 electric field at the wall reaches a maximum amplitude of 0.027%
of the unperturbed electric field at the wall. Each tick on the vertical axis represents

a factor of 5 in amplitude.

As seen from figure 4.6, this profile has wall electric field perturbations that decay
much more slowly than those of the Fast Decay Profile. Therefore, I have named the
density profile no(r) = no(0){1 — 6(£)* + 8(£)*® — 3(£)8} with a = 0.8 the Slow Decay
Profile. This profile was chosen because, in the low amplitude limit, the wall electric
field perturbation decay rate of this profile is similar to the decay rate found in the
experiment of N Sateesh Pillai [19]. Since the linear decay rate of the experiment
matches the linear decay rate for this profile, the Slow Decay Profile is expected to
resemble the experimental density profile. The Slow Decay Profile will be used for
the computations of Chapter 5 where higher amplitude applied pulses cause trapping
oscillations.

An m = 4 perturbation is also visible in this computation. As mentioned previ-
ously, the m = 4 perturbation is a harmonic that is driven by the m = 2 perturbation
in the plasma. The m = 4 mode has its own resonant frequency and decay rate that
differs from the m = 2 mode. Since the m = 4 mode is being driven by the m = 2
perturbation at the m = 2 resonant frequency rather than at the m = 4 resonant fre-
quency, there is an interference effect present that causes a modulation of the m = 4

perturbation.

4.2.3. No Decay Profile

Another density profile with a still different behavior is given by no(r) = no(0){1—
10(£)° + 15(£)8 — 6(£)°} with a = 0.8. The corresponding angular velocity profile is
wo(r) = wo(0){1 — $(£)® + 3(£)8 — (£)!°}. This density and angular velocity profile
is shown in figure 4.7. The central potential is ¢ = —0.165036£no(0)rz,, and the
unperturbed electric field at the wall is —0.16£7n9(0)rw.

Again a Gaussian pulse shape excitation is applied with pure m = 2 symmetry,
reaching a maximum amplitude of 1.18 x 10~*@) and having a FWHM of 1.06t,,. For
this profile, the m = 2 resonant angular frequency is found to be w, = 0.98w (0),
which corresponds to a resonant radius of r, = 0.81r,,.. The resonant radius lies

outside of the plasma edge, at r = 0.80r,, since a = 0.8. For the given applied

Unperturbed Normalized Density and Angular Velocity vs. Radius

density angular velocity

oO
an

resonant radius |
r=0.81

Normalized Profiles
o.Ucs
* ft
ants)

Oo
ba

ne or

O.2 o.4 0.6 QO.

Normalized Radius

Figure 4.7: Plot of the normalized unperturbed density and normalized angular veloc-
ity versus radius for the No Decay Profile no(r) = no(0){1— 10(£)® + 15(£)® —6(£)1°}
with the edge of the plasma at a = 0.8. The m = 2 resonant radius at r, = 0.817, is
indicated by the dashed vertical line, and lies outside of the plasma for this profile.

54
pulse, the electric field at the wall reaches a maximum perturbation of 0.038% of the
unperturbed electric field at the wall. For this profile, the electric field perturbations
at the wall are found to not decay at all. This profile is therefore named the No
Decay Profile. A semilog plot of the envelope of the perturbed electric field at the

wall versus time for this computation is shown in figure 4.8.

Zz eceesegadanesseses See ceeceeSeesecechccaeaneecgeeeteees laces cesetlseisiee iuessiisisissiissies i veces hee er eeeeeshs cesses tossessssel
OR pm Bef be eee ceed eceee der sntneeedeeeeseeeeibessssteecbeeeeeeeeelennnraeea be eeeeet ec ceecee less eeeereebesreeee denen
eo a
ia *s W,= 0.98a,(0)
2 2 a a bcestensduseussetteleenees ¥ eS
OP ccc hicsccsccedecsecressfevseseesstesseattes besesevsslevetsesees : bevbb bees cttetenteslevtrttttetlettveesedereeeseeed
Lo) :
a = jcseesedsscetsonaerisesseeedeeesicsibecesite beste besssresebesrsereediesersieedereeeses beeeseesdiesiesseefeeseresssberessseeedecrereeen
— YY *
2 § :
o oe EE A ee ee oe ee oes ; cee ee ee ea dh eteeseeeedeeeseesteeheessesesdeeeeeee ced
Se :
o :
6D a : ne Soe See
g :

0 15 30 45

Time (Central Rotation Periods)

Figure 4.8: Plot of the perturbed electric field at the wall versus time on a semilog
scale for the No Decay Profile no(r) = no(0){1 — 10(£)® + 15(£)® — 6(£)!°} with
a = 0.8. This calculation uses a pure m = 2 spatial excitation at the wall. The
applied pulse has a Gaussian shape, reaching a maximum amplitude of 1.18 x 10~*¢,
with a FWHM of 1.06¢,,. The perturbed m = 2 electric field at the wall reaches a
maximum amplitude of 0.038% of the unperturbed electric field at the wall. Each
tick on the vertical axis represents a factor of 5 in amplitude. The m = 2 resonant
radius lies outside the plasma for this profile, and thus the perturbed electric field at

the wall does not decay in time.

+535)
4.2.4 Other Profiles

To study the effect density profiles with differing degrees of plasma edge “sharpness”
has on the decay rate in a more systematic manner, a standard profile is chosen with
one free parameter that varies the profile. This density profile is given by the Fermi

function

nr) = no(0) aay

Renormalizing the equation so that n(0) = no(0) and n(1) = 0, the equation

becomes

e388 _ e(4r?—-1)
(e389 — 628 + e® — 1)(1 + e8(4r?-1))

no(r) = no(0) (4.1)

This profile smoothly approaches the origin and the wall, with the maximum

change in density near r = 5r,. Computations were performed with this profile
for the values of G = 2,4,6,8,10,12,14,16. The profiles are shown in figure 4.9.
All computations use the same applied Gaussian pulse as described in the previous
sections.

For small values of (3, the profiles look similar to the polynomial profiles discussed
earlier. As @ increases, the profile edge becomes sharper, and the profile resembles a
step function, as can be seen in figure 4.9. The sharper density profiles are found to
have a higher resonant frequency and a lower decay rate compared to density profiles
that change more slowly. The corresponding resonant frequencies and decay rates are
summarized in table 4.1. For G = 2, the value for a is similar to that for the Fast
Decay Profile, and for 6 = 6, the value for os is similar to the value for the Slow

Decay Profile. The Fermi function profile with @ = 8 is used in Chapter 7, and is
named the Beat Profile.

Unperturbed Density vs. Radius

o o a

Normalized Density

Normalized Radius

Figure 4.9: Plot of the normalized density versus radius for the normalized Fermi
function with 6 = 2,4,6,8,10,12,14, and 16.

Wr Ys y er

0.870 | 0.825 0.0741 i1.7
0.916 | 0.746 0.0314 29.2
0.977 | 0.716 0.0105 93.0
1.012 | 0.703 0.00302 335
1.030 | 0.697 | 7.66 x 10-* | 1340
1.040 | 0.693 | 1.76 x 10-* | 5910
14 | 1.046 | 0.691 | 3.68 x 10-° | 28400
16 | 1.050 | 0.690 | 6.40 x 10~® | 164000

looedl ileal
bo] | PO] @] ALM] @

Table 4.1: Decay Rates for Various Profiles.

57
4.2.5 Briggs, Daugherty, and Levy Decay Rate

Briggs, Daugherty, and Levy [10] derived an approximate analytic expression for the
perturbed electric potential decay rate of a pure electron plasma subjected to a small
excitation with angular symmetry m. The analysis made use of results derived from
a step function density profile, and assumes that the derivative of the density is small

except near the radius r = b. Their formula is

* wp — mw(b)]?

b |wo(rs)|

(Ts)
o(0)

mu bn4(rs)
m. no(0)

where y is the decay rate, no(r) is the unperturbed density profile, ng(r) is the
derivative of the density profile, wo(r) is the unperturbed angular velocity profile,
wo(r) is its derivative, d(r) is the perturbed electric potential, m is the azimuthal
symmetry eigenvalue, w, is the resonant angular frequency, and r,, the corresponding
resonant radius. Unfortunately, the quantities appearing in their expression can only
be obtained from a full solution to the problem. In particular, the resonant frequency
w, and the corresponding resonant radius r, are not known for a specific profile before
a numerical computation is performed. Furthermore, the parameter b is not clearly
defined, but will be taken as the radius where the absolute value of the derivative of
the density profile is maximum. The value of ¢(r) varies with time. It is expected,

though, that ¢(r,) will have the same time dependence as ¢(b), so the ratio OF

will be independent of time. The value for ¢(r) is taken at time t ~ 7.5t,,, which is
within the exponential decay range for all of these density profiles. The results from
the formula for the various Fermi function profiles are summarized in table 4.2.

In this table, y represents the decay rate of the perturbed potential. The column
% discrepancy represents the percentage error from the measured values in table 4.1.
It is found that in the best case, the formula can be within 7% of the value calculated
from the Nonlinear Fluid Code. When the density profile is not similar to a step
function, the formula is less accurate. The percent discrepancy in the computation
with @ = 16 may be caused by inaccuracies in the fluid code, where the finite grid

resolution does not accurately represent the sharp transition in density at r = b.

B b y “oe % discrepancy in “*
2 | 0.550 0.199 4.37 168%
4 | 0.515 0.0462 19.8 47.5%
6 | 0.507 0.0134 72.9 27.6%
8 | 0.504} 0.00352 288 16.3%
10 | 0.503 | 8.39x10-4 | 1230 8.94%
12 | 0.502 | 1.87x10-* | 5560 6.29%
14 | 0.501 | 3.96x10-° | 26400 7.58%
16 | 0.501 | 8.12x10-® | 129000 27.1%

Table 4.2: Comparison with Briggs decay rate for various profiles. For density profiles
which deviate strongly from a step function, the Briggs decay formula is not very
accurate.

4.3 Density Perturbations

The previous section discussed the damping in time of the perturbed electric field
at the wall for a variety of different unperturbed density profiles. This section will
attempt to relate the mechanism for this damping to the time dependence of the
density perturbations in the plasma.

The density perturbations are small and can not be observed in images of the total
density. Therefore, the average density at each radius is subtracted, leaving only the
perturbations in density. In the following images, the white areas represent regions
of positive density perturbation, while the black regions represent regions of negative
density perturbation. The various shades of grey represent a linear scale of density

perturbations.

4.3.1 Fast Decay Profile

As mentioned in the previous section, the density profile no(r) = no(0){1 — 3(£)? +
3(£)* — (£)°} with a = 0.8 exhibits a strong damping of the perturbed electric field
at the wall. A Gaussian pulse in time is applied to the wall with m = 2 azimuthal
symmetry. Figure 4.10 illustrates the resulting perturbed density for six successive
times.

As is seen from the figure, the initial density perturbation due to the applied

Figure 4.10: Density perturbations at six successive times for the Fast Decay Profile
no(r) = no(0){1—3(£)?+3(2)*—(£)°} with a = 0.8. The white represents a high pos-
itive perturbation, while black represents an equal magnitude negative perturbation.
Each successive shade represents an equal change in perturbed density, with lines
between shades representing constant perturbed density contours. Between frames,
the plasma rotates clockwise slightly more than half a revolution. As the perturba-
tion phase mixes, the contribution to the perturbed electric field at the wall decays
exponentially.

pulse is convected in the theta direction at different rates for different radii. This
shearing causes the density perturbation at each radius to become out of phase.
As time progresses, alternating bands of positive and negative density perturbation
become thinner and more evenly distributed, and thus the collective contribution to
the perturbed electric field at the wall rapidly decays away in time. The density
perturbations persist, even though the perturbed electric field at the wall decays
because of this phase mixing.

It should be noted that there are density perturbations which remain in phase,
visible near the central region of the plasma. Even though it is difficult to see,
there is a coherent mode extending throughout the plasma. The outer region that is
mixing has a perturbed magnitude that varies with the same phase as this coherent
mode. The positive and negative density perturbations are slightly stronger in the
areas in phase with the wave, even though the perturbation is being mixed. As time
progresses, the coherent central region does slowly mix.

For comparison to the results from the computation, figure 4.11 is a plot of the
density perturbations for different times if the layers merely convected in the theta
direction at their unperturbed angular velocity, wo(r) = wo(0){1 — 3(£)? + (£)4 —
+(£)°}. This ignores the requirement for a self consistent perturbed electric field,
and thus this figure does not display a collective coherent effect. Comparing the two
figures, the coherent portion of figure 4.10 should be clear. Further discussion of the
coherent portion versus the shearing portion of the density perturbation will be left

to a later section in this chapter.

4.3.2 Slow Decay Profile

Figure 4.12 illustrates the density perturbations due to an m = 2 Gaussian pulse
of voltage on the wall for six successive times for the Slow Decay Profile, no(r) =
no(0){1 — 6(£)* + 8(£)® — 3(£)8} with a = 0.8.

In comparing figure 4.12 with figure 4.10 for the Fast Decay Profile, it is clear

that the collective coherent portion of the density perturbation is much stronger for

Figure 4.11: Density perturbations at six successive times for the Fast Decay Profile
no(r) = no(0){1 — 3(£)?+3(£)* — (£)°} with a = 0.8, assuming the plasma convected
in the theta direction at the unperturbed angular velocity, ignoring the self consistent
perturbed electric field.

Figure 4.12: Density perturbations at six successive times for the Slow Decay Profile
no(r) = no(0){1 — 6(£)* + 8(£)® — 3(£)°} with a = 0.8. The white represents
a high positive perturbation, while black represents an equal magnitude negative
perturbation. Each successive shade represents an equal change in perturbed density,
with lines between shades representing constant perturbed density contours.

63
the Slow Decay Profile than for the Fast Decay Profile. As time progresses, the
coherent density perturbation slowly decays, and the shearing component becomes

more evident, slowly progressing towards the center of the plasma.

4.3.3. No Decay Profile

As a final comparison, figure 4.13 illustrates the time evolution of density perturba-
tions for the No Decay Profile, no(r) = no(0){1 — 10(£)® + 15(£)® — 6(£)!°} with
a = 0.8 after application of an m = 2 Gaussian pulse of voltage applied to the wall.
For this profile, the coherent collective density perturbation does not decay in time.
However, the shearing component still has an observable effect in the plasma. The
shearing mode does break the density perturbations into filaments, but the collective
mode does not allow these filaments to wrap around and mix as with the previous

two density profiles.

4.4 Linear Model

This section presents a simple approximation based on two exponentials which helps
understand the behavior of the perturbed density as a function of time. In section 4.2,
the perturbed electric field at the wall was found to decay exponentially in time after
excitation from an externally applied pulse. This exponential decay after the applied
pulse can be modeled in the transformed space as a simple pole in the complex plane
[19]. This pole is located at

w* = Wp — 1y

where w* is the complex pole, w, is the real resonant frequency, and ¥ is the decay

rate. Rewriting the linearized fluid equations 2.14 and 2.15 from Chapter 2 we have

m dno(r) ~ (4.2)

Mym(1,w) = r(mwo(r)—w) dr "*

Figure 4.13: Density perturbations at six successive times for the No Decay Profile
no(r) = no(0){1 — 10(£)® + 15(£)® — 6(£)!°} with a = 0.8. The white represents
a high positive perturbation, while black represents an equal magnitude negative
perturbation. Each successive shade represents an equal change in perturbed density,
with lines between shades representing constant perturbed density contours.

29, 1b im _ me m dno(r) 5 9
dr? r dr r2 r(mwo(r)—w) dr 77 h™

which can be solved for the transformed perturbed density 71.(r,w) and potential
bi m(t, w). Using the simple pole approximation for the transformed perturbed po-
tential, the approximate solution takes the form

Am(r,w)

W— wt”

Pim(,w) = (4.3)

This approximation is useful after the applied pulse has been turned off, and
during the time the decay is exponential. The function A,,(r,w) is the solution to
the equation with the simple pole at w* removed. It is assumed that the pole is the
dominant feature, and the features of the function A,,(r,w) are small compared to
this pole.

Inserting this form for the perturbed potential 4.3 into the equation for the per-
turbed density 4.2 we get

m dno(r) Am(r, w)
r(mwo(r)-—w) dr w—u*

Nim(T, Ww) =

or absorbing some of the constants and functions of r into a new function C(r,w) we

obtain
1 1

w — w*w — muo(r)

Ni m(r, Ww) = C(r,w)

Performing an inverse Laplace transform on the perturbed density results in
Nim(7, t) = D,(r)e™"* + D.(r)e—mort
which can also be written as

Ni m(r,t) = D(r)eNer" + D (rye reo", (4.4)

This form for the perturbed density is comprised of two components. The first
term is a collective coherent mode with resonant frequency w, which decays expo-
nentially. This collective coherent mode which decays will be taken as the quasimode
[15]{19]. The radial dependence of the quasimode has been lumped into the coefficient
D,(r). The second term does not decay in time. It is also not coherent, but rather has
a frequency that is a function of the radius, and therefore a resulting shear. Since the
frequency is a continuous function of the radius, this term describes the continuum
modes. The radial dependence of the continuum modes has been incorporated into
the coefficient D,(r).

This model contains several features. If the quasimode term is larger than the
continuum term, the quasimode will dominate the perturbed density early in time.
However, as the quasimode decays, the continuum modes will become more apparent
in the perturbed density. At some time and at some radius, the two terms should
be comparable, and an interference pattern should be visible. Looking back at figure
4.10 and figure 4.12, there is a visible interference between these two modes occurring
near the radius that separates region that is sheared and region that is coherent.

To further illustrate the interference between the two modes, figure 4.14 plots the
magnitude of the m = 2 density perturbations versus radius for the computation
using the Fast Decay Profile. The frames (a)-(f) represent the same times used in
figure 4.10.

From this plot, the interference pattern is clearly visible. Additionally, the decay
of the quasimode is also visible as the maximum amplitude decays from frame (a) to
frame (f). After the quasimode has decayed, the remaining contribution is from the
continuum modes, which do not decay in time. In frame (f), the quasimode has nearly
decayed away. Figure 4.15 is the same plot of the perturbed density magnitude versus
radius, but for a time much later than frame (f), where the quasimode has completely
decayed away.

Figure 4.15 shows the form of the continuum mode as a function of radius, and
figure 4.14(a) shows the sum of the quasimode and continuum mode early in time. It

is therefore possible to use the function of figure 4.15 for the continuum term D,(r)

(a) (b)

Amplitude
Amplitude
a —

01 0.2 0.3 0.4 0.5 0.6 0.7 oT 0.2 0.3 0.4 0.5
Radius Radius
() (d)
/\ ;
a \ aw \
aa \ ae \
0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.1 0.2 0.3 0.4 0.5 0.6 0.7
Radius Radius
(e) (f)
oO vo
5 F

fl i

/\f \ J
YN \ \ ava \
03.04 03 : : :

0.6 0.7 0.1 0.2 0.3 0.4
Radius

0.1 02 4
Radius

Figure 4.14: Magnitude of density perturbations versus radius at six successive times
for the Fast Decay Profile. As time progresses, an interference can be seen between
the quasimode and the continuum modes.

Amplitude

0.1 0.2 03 05 0.6 0.7

0.4
Radius

Figure 4.15: Magnitude of continuum density perturbations versus radius for the Fast
Decay Profile obtained after the quasimode has decayed away.

69
in the linear model. Additionally, the function of figure 4.15 subtracted from the
function of figure 4.14(a) can be used for the quasimode term D,(r). Figure 4.16
shows the resulting density perturbations at six successive times from using these two
terms in the linear model. A good qualitative agreement is seen by comparing figure
4.16 to the results from the full computation in figure 4.10.

This approximate two-exponential solution for the perturbed density is successful
at reproducing images of the perturbed density that agree with the results from the
code. However, the decay rate of the quasimode term depends on the phase mixing
caused by the continuum modes, and can not be determined for an arbitrary density

profile without performing a numerical computation.

4.5 Conclusion

This chapter examined the response of a pure electron plasma to a small amplitude
pulse applied at the wall. The response of the perturbed electric field at the wall is
found to decay exponentially. The decay rate is a strong function of the unperturbed
density profile, and roughly matches the approximate decay rate expression derived
by Briggs, Daugherty, and Levy. However, their expression for the decay rate can not
predict the decay rate for a given profile since it requires values that are not known
before a computation is performed. Their expression does not agree well with the
computation for density profiles which strongly deviate from a step function. For a
certain class of power law profiles, Corngold has found an analytical expression for
the decay rate [15].

Images of the perturbed density illustrate that the density perturbation undergoes
shear. The shear in the plasma flow causes the density perturbations to become out
of phase at different radii. This phase mixing causes the perturbed electric field at
the wall to decay exponentially.

Application of a linear model suggests that the perturbed density can be separated
into a coherent quasimode term, and a shearing continuum mode term. However, the

parameters for this model must still be obtained from a complete linear calculation.

Figure 4.16: Resulting density perturbations for six different times for the Fast Decay
Profile using the Quasimode and Continuum mode linear model.

71
Using the parameters from the linear calculation in this model enables the images of
the perturbed density to be reproduced.

The Nonlinear Fluid Code successfully reproduced the results from the Linear
Fluid Code for small amplitude disturbances. This provides verification of the Non-
linear Fluid Code method. This chapter studied several different density profiles.
However, a systematic study of the effect of profiles on the plasma response has
not been completed. The next chapter will examine the plasma response to larger

amplitude pulses, where nonlinear effects become important.

Chapter 5 Fluid Trapping

5.1 Introduction

Chapter 4 examined the effect of a low amplitude externally applied pulse on a pure
electron plasma. The response in the perturbed electric field at the wall due to this
applied pulse was found to decay exponentially. This decay was explained using
linear theory, which described a phase mixing process for the perturbed density. ‘This
chapter will examine the effects of single, larger amplitude pulses, where the linear
theory is no longer valid and nonlinear self-interactions become important.

Briggs, Daugherty, and Levy [10] predicted that with higher amplitude pertur-
bations, fluid elements might be trapped in the large amplitude excited wave. Such
trapping would cause a quasi-periodic modulation of the decaying amplitude of the
perturbed electric field at the wall as the trapped fluid elements bounce within the
potential well of the wave. According to Briggs, Daugherty, and Levy, the quasi-
periodic modulation should have a bounce frequency proportional to the square root
of the amplitude of the perturbed potential.

Pillai’s [19] experimental data, showing the envelope of the oscillating electric
field at the wall versus time on a semilog scale, is reproduced in figure 5.1. He
found that for small amplitude applied pulses, the perturbed electric field at the wall
decayed approximately exponentially. For higher amplitude pulses, a modulation of
the decay envelope occurred with a bounce frequency dependence proportional to the
square root of the amplitude, agreeing with the bounce frequency dependence on wave
amplitude predicted by Briggs, Daugherty, and Levy [10]. Because of this agreement
with the Briggs prediction, Pillai attributed the envelope modulation to trapping of
fluid elements. Since there was no diagnostic in the experiment to measure the density
distribution directly, these trapped fluid elements could not be directly observed.

In order to show that trapped fluid elements are responsible for the modulation

he
ate

Amplitude (dB)

core Gees, ae ne, Senne rere Cleeneenrnenh | Swann © Satire rrrrer ens Sener Seer erereeneean, © MPTrrerrrrrttrertirrt. “tetrceisnrre, +, “Tltrn

Time in pS 300

Figure 5.1: Plot of the envelope of the perturbed electric field at the wall versus
time on a semilog scale for the experimental data collected by N. Sateesh Pillai. The
applied m = 2 excitation voltage ranges from 10 mV for the lowest amplitude pulse,
to 1000 mV for the highest amplitude trace. The vertical scale is calibrated in dB.

of the decay envelope, the nonlinear 2-D cylindrical fluid code is used to study the
response of a pure electron plasma to higher amplitude pulses. The results from
the fluid code qualitatively agree with the experimental results. For low amplitude
perturbations, the envelope of the oscillating electric field at the wall decays expo-
nentially. For larger perturbations, a bounce oscillation is observed in the decaying
envelope. The bounce modulation is caused by fluid elements that are trapped in the
wave. The trapped fluid elements are directly observed in images of the perturbed
density. For a perturbation that is slowly applied to the wall, a non-decaying wave
without bounce oscillations is excited, which is similar to a BGK mode [20]. A BGK

mode is a large amplitude mode that contains regions of plasma trapped within the

excited wave, and is stationary in a frame rotating with the wave.

5.2 Fluid Code Results

This section describes results from the Nonlinear Fluid Code for pulses applied to
the Slow Decay Profile. The results are observed in the perturbed electric field at
the wall, and in images of the perturbed density. The Slow Decay Profile, no(r) =
no(0){1 — 6(£)* + 8(£)® — 3(£)8}, with a = 0.8, previously discussed in Chapter 4,
was chosen because its decay rate closely matches the experiment above for small
amplitude applied pulses. This density profile may resemble the density profile in the
experiment. This profile has a central electric potential of d) = —0.142696 2 n9(0)ri,,

and an unperturbed electric field at the wall of —0.12851n9(0)rw.

5.2.1 Perturbed Electric Field

This section describes the perturbed electric field at the wall resulting from pulses
applied to the slow decay profile. The electric field at the wall is systematically
studied as a function of applied pulse amplitude.

The 2-D cylindrical fluid code is used to study the response of a pure electron
plasma to an externally applied pulse. A low amplitude sinusoidal pulse was applied
for three cycles with an amplitude of 4.28 x 107°¢9. The resulting perturbed electric
field at the wall reaches a maximum amplitude of 0.00421% of the unperturbed electric
field at the wall and decays exponentially in time. The perturbed electric field at the
wall versus time is plotted in figure 5.2.

As seen in the figure, the response from this applied pulse decays exponentially.
The pulse amplitude is in the linear regime, so that the results of Chapter 4 are
valid for an applied pulse with this amplitude. This computation is similar to a
computation in Chapter 4, where a semilog plot of the envelope of the perturbed
electric field versus time for this computation was shown in figure 4.6.

For comparison, figure 5.3 is a plot of the perturbed electric field at the wall versus

time for an applied pulse amplitude that is 200 times higher, with a three cycle pulse

0.005 T T T T T T T T T T T T T T
_~ | = 0089040) |
= 0.0089w,{0

x y=0. 0
a L arly = 95

o : ,
ca iI

8 0

mer}

' L
ym | 1
oO

-0.005 a ee ee ee ee
0 20 40 60

Time (Central Rotation Periods)

Figure 5.2: Plot of the perturbed electric field at the wall for the Slow Decay Profile,
using an applied three cycle pulse amplitude of 4.28 x 10~-°¢). The perturbed electric
field at the wall reaches a maximum amplitude of 0.00421% of the unperturbed electric
field at the wall and decays exponentially in time.

amplitude of 8.55 x 1073. For this applied pulse, the perturbed electric field reaches
a maximum amplitude of 0.851% of the unperturbed electric field at the wall.

The plasma response at the wall is quite different from the previous lower am-
plitude example. For this higher amplitude applied pulse, the decay rate has been
reduced. Also, a bounce modulation of the decay envelope is evident.

A systematic study of the dependence on amplitude was made for this profile,
varying the amplitude of the initial three cycle pulse. The highest amplitude case has
an applied pulse amplitude of 0.0214¢) causing a perturbation in the electric field
at the wall of 2.07%. The lowest has an applied pulse amplitude of 1.07 x 107d,

1% T T T T T T T T T T T T T T

Perturbed Electric Field (% )

i i 1 1 I 1 1 1 L

20 40 60
Time (Central Rotation Periods)

-1% L 1 \ \

Figure 5.3: Plot of the perturbed electric field at the wall demonstrating bounce
oscillations. A three cycle pulse with an amplitude of 8.55 x 10~*¢, is applied to the
Slow Decay Profile. The perturbed electric field reaches a maximum amplitude of
0.851% of the unperturbed electric field at the wall.

leading to a maximum response of 0.0105% of the unperturbed electric field at the
wall. The logarithm of the envelope of the perturbed electric field of the wall versus
time is plotted in figure 5.4.

As seen in the figure, the lower amplitude traces decay exponentially. At higher
amplitudes, the perturbed electric field at the wall has a lower decay rate. Addi-
tionally, a modulation of the decay envelope occurs. The time scale for this plot
was chosen for easy comparison with the experimental data in figure 5.1. The cal-
culation is terminated after 60 central rotations because of computational time and

accumulation of errors.

Log of Perturbed Electric Field
(Each tick is a factor of 2)

Time (Central Rotation Periods)

Figure 5.4: Plot of the envelope of the perturbed electric field versus time on a
semilog scale for various amplitude pulses applied to the Slow Decay Profile. The
highest amplitude case has an applied three cycle pulse amplitude of 0.0214¢) causing
a perturbation in the electric field at the wall of 2.07%. The lowest has an applied
pulse amplitude of 1.07 x 10~*¢, leading to a maximum response of 0.0105% of the
unperturbed electric field at the wall. The lower amplitude traces show an exponential
decay, while the higher amplitude traces show a modulation of the decay envelope.
The computations were terminated after 60 central rotations due to computational
time and accumulation of errors.

In comparing the computational data of figure 5.4 to the experimental data of
figure 5.1, these two figures exhibit remarkable similarities. Both traces have the
same low amplitude decay rate. This is expected as the slow decay density profile
was chosen to have the same low amplitude decay rate as the experiment. As the
amplitude is increased, both the experimental data and the computation show a
decrease in the decay rate and a modulation of the decay envelope.

The modulation frequency does, however, differ between the results from the com-
putation and the results from the experiment. Additionally, the decay rates differ at
higher amplitudes. The difference is probably due to imperfect knowledge of the ex-
perimental density profile. In the experiment, no direct measurement could be made
of the radial dependence of the density. The Slow Decay Profile was chosen to match
the experimental decay rate at low amplitudes, but this profile is not unique. Addi-
tionally, other effects such as viscosity which are not included in the fluid code may
also be responsible for the discrepancy.

Since there is a good qualitative agreement between the experiment and computa-
tion, it is expected that the fluid code contains the physics responsible for the bounce

oscillation found in the experimental data.

5.2.2 Bounce Frequency Dependence
According to Briggs, Daugherty, and Levy [10], the quasi-periodic modulation should

have an angular frequency given by the formula

2 291(Ts) dwo(r) (5.1)

T=Ts

where wy, is the bounce oscillation angular frequency, m is the azimuthal symmetry
eigenvalue, r, is the resonant radius, ¢,(r) is the perturbed electric potential, and
wo(r) is the unperturbed angular velocity profile. This formula predicts that the
bounce frequency should be proportional to the square root of the amplitude of the
perturbed potential. As with their linear decay formula, the other factors are difficult

to determine, but they are assumed not to change when the potential is varied.

Some features of this bounce formula may be easily understood by considering
a one dimensional model of a particle in a potential 6 = ¢)cos(m@). The electric
field E = —mdysin(m@), and the force F = gE, where q is the electric charge of
the particle. Using Newton’s second law, linearizing around @ = 0, and solving the
differential equation results in periodic motion with a frequency w, ~ m/ bo .

In figure 5.4, the 12 higher amplitude plots complete at least one half period of
bounce oscillation. From these plots, a bounce angular frequency can be determined,
which can then be compared to the square root relationship above.

Some care must be taken when calculating the bounce frequency. The modulation
of the decay envelope is imposed on a decaying exponential. To obtain the maxima
and minima of the oscillation, a straight line was fit to the logarithm of the perturbed
electric field at the wall for each trace to obtain a decay rate. Then this exponential
decay was subtracted from the data, leaving the modulating envelope. The time
between the first minimum and the next maximum, corresponding to half a bounce
period, was used to calculate the bounce frequency.

Additionally, the perturbed potential used for each trace is calculated in a specific
manner. Equation 5.1 requires the electric potential at the resonant radius. The
perturbed potential is assumed to be proportional to the perturbed electric field at
the wall. During a bounce oscillation, the perturbed electric field at the wall decays
exponentially. Therefore, the perturbed electric field used for comparison to equation
5.1 was averaged over the time of the half bounce period between the first minimum
and the next maximum.

Using the calculated bounce angular frequency and the averaged perturbed poten-
tial during this bounce oscillation, a plot was made on a log-log scale of the bounce
frequency versus perturbed amplitude for the 12 high amplitude traces of figure 5.4.
This plot is shown in figure 5.5.

On this log-log scale, the points fall along a straight line. The dashed line has a
slope of 0.5. A least squares fit through the data points gives a slope of 0.493 + 0.008.
This indicates that the bounce frequency is proportional to the square root of the

perturbed electric field at the wall.

Bounce Angular Frequency
(uruts of central rotation angular frequency)
as

0.1 1 10

Perturbed Electric Field at Wall (%)

Figure 5.5: Plot of the angular bounce frequency versus the perturbed electric field at
the wall for 12 high amplitude traces. The dashed line has a slope of 0.5, indicating
that the bounce frequency is proportional to the square root of the perturbed electric
field at the wall.

Agreement of the computation with the square root relationship suggests that the
modulation of the decay envelope is due to the trapping of fluid elements. The next

section will describe a search for these trapped fluid elements directly by observing

images of the perturbed density.

5.2.3. Density Perturbations

This section describes images of the density and perturbed density resulting from
application of an external pulse to the Slow Decay Profile. The pulse has an amplitude
sufficiently large to cause a modulation of the decay envelope.

A slightly different external pulse is applied to the Slow Decay Profile. A Gaussian
pulse shape is applied in the frame moving with the wave, having an m = 2 resonant
frequency of 0.851w 9(0). The applied pulse has a maximum amplitude of 0.005¢p,
and a FWHM of 1.06¢,,. The resulting perturbed electric field at the wall reaches a

81
maximum perturbed amplitude of 1.18% of the unperturbed electric field at the wall.
A plot of the envelope of the perturbed electric field at the wall for this computation

is shown in figure 5.6.

1.5% t qT T T T T T T T T T T T T

Perturbed Electric Field ( % )

20 40 60

“1.5963 4
Time (Central Rotation Periods)

Figure 5.6: Plot of the envelope of the perturbed electric field at the wall for a
typical computation displaying bounce oscillations. The applied pulse has a Gaussian
pulse shape, reaching a maximum amplitude of 0.005¢), and a FWHM of 1.06¢,,.
The letters (a)-(f) correspond to times 2.65t., 7.36te, 12.1tgp, 18.3ter, 24.4tep, and
31.5te, respectively. The perturbed field reaches a maximum amplitude of 1.18%
of the unperturbed electric field at the wall. The bounce oscillation has an angular
frequency of wy = 0.0407w9(0).

As seen from the figure, the perturbed electric field at the wall exhibits a strong
bounce oscillation. The bounce oscillation modulates the decay envelope with a max-
imum modulation of 52.9% at an angular frequency of w, = 0.0407w (0). The next

section will look at images of the density to see if this modulation can be observed.

82
Total Density

This section describes images of the density resulting from the application of the
pulse described in the above section. In the described computation, the perturbation
is large enough for nonlinear effects to occur. The perturbed electric field at the
wall is greater than 1%, and nonlinear bounce oscillations are observed. It might be
expected that an observable effect could be seen in the total density of the plasma.
Figure 5.7 shows images of the total density of the plasma at six different times in
the reference frame that is rotating with the wave.

Images (a) - (f) in figure 5.7 correspond to the total density at times 2.65t,,,
7.36ter, 12.1tep, 18.3t¢-, 24.4¢.,, and 31.5t,,, respectively. These times correspond to
covering the duration of the first bounce cycle, labeled (a)-(f) in figure 5.6.

As seen in the figure, very little difference can be seen between the different frames
for the total density picture. Even the initial strong perturbation is barely visible,
giving the total density a very slight elliptical shape. This indicates that even though
nonlinear bounce oscillations are observed, the amplitude of the perturbation is not
so high as to significantly distort the shape of the plasma from the unperturbed
configuration. In other words, nonlinear effects may occur at amplitudes that hardly

distort the plasma.

Perturbed Density

Since the density perturbation is too small to be seen in the total density, images
of the perturbed density are plotted in figure 5.8. These images are generated by
subtracting the unperturbed density, revealing the difference of the density from the
unperturbed density. The times for each frames are the same times used in figure 5.7.

The perturbed density in figure 5.8 shows structure and variation that could not
be seen in images of the total density. Some of the structure is similar to that found
in figure 4.12 from Chapter 4, which shows images of the perturbed density (without
the gap to the wall) at different times after excitation from a low amplitude applied

pulse. This similarity in structure is due to the same decay mechanism that exists

Figure 5.7: Images of the total density at six different times. An external Gaussian
pulse with a maximum amplitude of 0.005¢), and a FWHM of 1.06¢,, is applied
to the Slow Decay Profile. The resulting perturbed electric field at the wall exhibits
strong bounce oscillations. Figures (a)-(f) correspond to times 2.65t¢,, 7.36ter, 12.1ter,
18.3t¢-, 24.4t,, and 31.5t,,, respectively. The white region represents an area of high
density, and the black area represents a region of low density. Very little change is
visible in the total density pictures.

Figure 5.8: Images of the perturbed density at six different times for an external
Gaussian pulse with a maximum amplitude of 0.005¢), and a FWHM of 1.06t,, ap-
plied to the Slow Decay Profile. The resulting perturbed electric field at the wall ex-
hibits strong bounce oscillations. Figures (a)-(f) correspond to times 2.65t,,, 7.36te,,
12.1t, 18.3te, 24.4¢-, and 31.5t.,, respectively. The white area represents an area
of positive density perturbation, where the black region indicates an area of negative
perturbation. Near the edge of the plasma, two lobes of plasma form. These lobes are
regions where fluid elements are trapped in the wave. The oscillation of the trapped
fluid elements in the potential wall causes the bounce modulation of the perturbed
electric field at the wall.

in both computations, which causes the perturbed electric field at the wall to decay
exponentially. However, the mixing seen in figure 5.8 is less than in figure 4.12. This
decrease in the mixing at higher amplitudes accounts for the corresponding lower
decay rate at higher amplitudes. Essentially, the fields from the higher amplitude
perturbation, caused by the higher amplitude applied pulse, are better able to hold
the density perturbations in phase, causing the fluid elements to follow nonlinear
orbits which do not mix as much at different radii. This leads to the lower decay rate
in the perturbed electric field at the wall for higher amplitude applied pulses.

Another significant difference is visible in the images of the perturbed density.
Near the outer region of the plasma, a tail of plasma is left behind. As time progresses,
this tail stretches past the region of negative perturbation, and joins the other lobe
of positive perturbation.

To see this effect still more clearly, an annulus containing the outer edge of the
plasma was magnified and stretched flat for each of the images in figure 5.8. These
strips of perturbed density are shown in figure 5.9.

These strips displaying the perturbed density are a magnification of the outer
plasma edge. The vertical axis represents the radial direction of the previous figure,
while the horizontal axis represents the theta direction.

In figure 5.9, one can actually see the trapped fluid elements oscillating in a
potential well. The oscillation of the trapped fluid elements is the cause of the bounce
oscillation seen in the perturbed electric field at the wall. The perturbed electric field
at the wall is found to be minimum during the time the trapped fluid elements are
traversing the region of negative density perturbation. The perturbed electric field
at the wall is maximum when the trapped fluid elements are bunched up near the
region of positive density perturbation. Some fluid elements near the outer edge are
observed not being trapped, but instead are passing fluid elements. The figure is
displaying a physical space image of Kelvin’s “cat’s eyes” [29] for trapped orbit fluid
flow.

A final observation is that the density perturbation throughout the entire plasma

oscillates with the trapped fluid elements. This feature is visible by comparing image

Figure 5.9: Magnification of the trapped fluid elements causing bounce oscillations.
The images shown in this figure are a magnification of the edge of the plasma in
figure 5.8, stretched flat. The vertical axis represents radius, while the horizontal axis
represents theta. As time progresses from frame (a) to frame (f), the nonlinear orbits
of trapped fluid elements are visible. The oscillation of the fluid elements within the
potential well cause the modulation in the perturbed electric field at the wall. The
image is suggestive of Kelvin’s “cat’s eyes” for trapped orbit fluid flow.

87
(c) and (e) of the figure 5.8. The coherent positive density perturbation region of
image (e) is brighter than the coherent positive density perturbation region of im-
age (c). Thus a type of feedback mechanism exists. As the trapped fluid elements
approach the region of positive density perturbation, the field due to the trapped
fluid elements add to the field due to the coherent density perturbation. The in-
crease in the perturbed electric field increases the perturbed density of the coherent
density perturbation, again amplifying the perturbed field. Thus the trapped fluid
elements do not oscillate independently within the potential well, but are coupled to

the amplitude of the well by the corresponding fields.

5.2.4 BGK-Like Modes

This section describes a pulse that is slowly applied to the plasma. The plasma
response has features resembling a BGK mode [20]. In the wave frame, a BGK mode
is a steady nonlinear wave that is supported by a small class of particles trapped in
the wave. A BGK mode would create an oscillating perturbed electric field at the
wall that does not decay in time. In the wave frame, the perturbed density would
not change in time, and the density is a function of the potential.

A high amplitude, short duration applied pulse is found to cause a decrease in the
decay rate of the perturbed electric field at the wall and a bounce oscillation in the
decay envelope. The bounce oscillation of the perturbed electric field is caused by
discrete trapped fluid elements bouncing back and forth in the potential well of the
wave. If the potential well was instead filled with fluid elements distributed evenly
throughout the well, there would not be a bounce oscillation of the envelope of the
perturbed electric field at the wall.

The discreteness of the trapped fluid elements is believed to be caused by the short
time duration of the applied pulse. Such a short, high amplitude pulse perturbs the
plasma in a manner that does not allow enough time to fill the potential well with
fluid. Exciting the plasma with a smaller amplitude pulse for a longer duration in

time may allow the plasma to distribute itself evenly throughout the potential well,

88
and eliminate the bounce oscillations, resulting in a perturbation resembling a BGK
mode.

A long duration low amplitude external pulse is applied to the Slow Decay Profile.
The external pulse has a Gaussian pulse shape and is applied in the frame moving
with the wave with m = 2 azimuthal symmetry. The pulse has a maximum ampli-
tude of 2.755 x 10~*¢@ 9, and a full width at half maximum of 31.8¢,,. Note that the
applied pulse has a duration that is 30 times as long as the pulse in the previous
section. To obtain the same amplitude for the perturbed electric field at the wall,
the applied pulse amplitude is greater than 30 the amplitude of the pulse from the
previous section. This is necessary because of the linear decay mechanism discussed
in Chapter 4. Over the longer duration of this applied pulse, the decay mechanism
prevents the perturbed electric field at the wall from reaching the same amplitude.
Therefore the applied pulse amplitude is increased to compensate for the decay. The
resulting maximum perturbed electric field at the wall is the same in both computa-
tions. However, in this computation, the energy is focused into a smaller frequency
spectrum, centered around the m = 2 resonant frequency. A plot of the envelope of
the resulting perturbed electric field at the wall is shown in figure 5.10.

The perturbed electric field at the wall reaches a maximum amplitude of 1.18%
of the unperturbed electric field at the wall. By using a slowly applied pulse, the
bounce oscillations have been nearly eliminated. The decay of the perturbed electric
field at the wall has also been greatly reduced compared to the decay rate for a small
amplitude applied pulse that was shown in figure 5.2.

Figure 5.11(a) is a plot of the perturbed density at time 147.2¢,,, and Figure
5.11(b) is a magnification of the outer plasma edge. The region of trapped fluid
is visible near the plasma edge, and more evenly distributed for this slowly applied
pulse, forming very distinctive “cat’s eyes.” With the plasma in the potential well
evenly distributed, there is no bounce oscillation in the perturbed electric field at the
wall. Additionally, the flow within the potential well appears to be sheared, causing
the perturbed density within the well to form a spiral structure.

For an ideal BGK mode, the perturbed density in figure 5.11 would be independent

1.5% qT T v T T T T T T T T T T T

Perturbed Electric Field ( % )

200 300

-1.5% po
0 100
Time (Central Rotation Periods)

Figure 5.10: Plot of the envelope of the perturbed electric field at the wall resulting
from a slow Gaussian pulse applied to the Slow Decay Profile and exciting a BGK-
like mode. The applied pulse has a maximum amplitude of 2.755 x 10~*@ , and a
FWHM of 31.8¢,,. The perturbed field reaches a maximum amplitude of 1.18% of the
unperturbed electric field at the wall. By using a slowly applied pulse, the bounce
oscillation has been nearly eliminated.

of time. For this computation, the perturbed density is steady, with the exception
of slight motion of trapped fluid elements within the well forming a spiral structure.
It should also be noted that this computation is pushing the limits of the fluid code
by computing to a time of nearly 300¢,,, and may have more accumulated error than

other computations in this thesis.

Figure 5.11: (a) Perturbed density resulting from a slowly applied Gaussian pulse
exciting a BGK-like mode perturbation on the Slow Decay Profile. The applied pulse
has a maximum amplitude of 2.755 x 10~“¢,, and a FWHM of 31.8¢,,. Near the edge
of the plasma, two lobes of plasma form. These lobes are regions where fluid elements
are trapped in the wave. Image (b) is a magnification of the outer plasma edge where
trapping occurs, stretched out flat. The “cat’s eyes” are clearly visible in (b). The
flow within the potential well appears sheared, causing the perturbed density within
the well to form a spiral structure. The density within the potential well is evenly
distributed, thus there is no bounce oscillation in the perturbed electric field at the
wall.

91
5.3 Conclusion

In Chapter 4, low amplitude applied pulses were shown to result in an exponentially
decaying perturbed electric field at the wall. This chapter showed that for higher
amplitude pulses, a modulation of the envelope of the perturbed electric field is ob-
served numerically, consistent with the experiment [19] and the Briggs prediction [10].
The modulation frequency is proportional to the square root of the amplitude of the
perturbed electric field at the wall. The modulation of the decay envelope is caused
by trapped fluid elements which are directly observable in successive images of the
perturbed density. The oscillation of the fluid elements within the well corresponds
to the bounce oscillations in the perturbed electric field. The bounce oscillations can
be minimized by producing the same magnitude density perturbation with a lower
amplitude external pulse applied over a longer period of time, exciting a perturbation
within the plasma that resembles a BGK mode.

Even though nonlinear effects are occurring at the applied amplitudes discussed
in this chapter, the total perturbations are not large, and can not be easily seen in
images of the total density. Chapter 7 will look at nonlinear effects that can occur
with very large applied amplitudes, where the perturbations can be seen in images of

the total density.

Chapter 6 Echoes in a Pure Electron

Plasma

6.1 Introduction

Chapter 4 examined the response of the plasma to low amplitude applied pulses. The
perturbed electric field at the wall was found to decay due to phase mixing of various
radial layers of the perturbed density. This phase evolution of the perturbed density
can be reversed by a nonlinear interaction with a second applied pulse. Application
of a second pulse results in a subsequent reappearance of a perturbed electric field
at the wall, many damping periods after the application of the second pulse. This
reappearance of the perturbed electric field at the wall is referred to as an echo.

Echoes occur in systems where perturbations decay by a non-dissipative process,
and many different echo phenomena have been observed [30][31][32][33]. In neutral
plasmas, the plasma wave echo has been studied theoretically [34], and observed ex-
perimentally [35]. Since the perturbed electric field at the wall decays by a collisionless
process in a non-neutral plasma, it is reasonable to expect that this system might also
exhibit echoes. At this time, echoes have not been experimentally observed in non-
neutral plasmas. Echoes in non-neutral plasmas have been analyzed by an extension
of the linear analysis of Chapters 2, 3, and 4 [21], but those results are restricted to
small amplitudes. A similar quasi-linear argument is presented in section 5.2. The
Nonlinear Fluid Code is not restricted to small amplitudes and can be used to study
echo phenomena.

This chapter will examine the effects of applying two pulses which are separated
in time to a pure electron plasma. The first section presents a rough argument for
the interaction of the two pulses. The argument describes the angular dependence of

the echo, the time when the echo occurs, and the dependence of the echo amplitude

93
on the two applied pulse amplitudes and as a function of the time separation between
the two applied pulses.

The next section describes echo calculations using the Nonlinear Fluid Code. Two
external pulses, separated in time, are applied to a pure electron plasma. An echo
response is found to occur long after the responses from the two applied pulses have
decayed away. For small applied pulses, the echo is found to have a time and am-
plitude dependence that agrees with the rough argument for the echo. For higher
amplitude pulses, the echo amplitude saturates, and higher order echoes caused by
the interaction with harmonics of the applied pulses are observed. Echoes of echoes

are also found to occur.

6.2 Echo Concept

This section presents a simplified analysis for the interaction of two applied pulses
generating an echo. The echo occurs from a nonlinear interaction between the two
applied pulses, where the second pulse partially unmixes the first pulse. The two
applied pulses are treated in the linear approximation discussed in Chapter 2, and
the echo is generated by a nonlinear interaction between pulses that is second order
in amplitude. This quasi-linear description will lead to the following approximate
results: symmetry requirements for an echo to occur, the angular symmetry of the
echo, the time the echo occurs, and the amplitude dependence of the echo. These
results are analogous to the results for the plasma wave echo [34]. Since this is a quasi-
linear argument, the results are limited to small applied pulses. The results from the
Nonlinear Fluid Code, presented in the next section, contain the full nonlinear effects
and are not limited to small pulses.

The following argument is similar to the treatment presented by Gould [21]. Using
the fluid model presented in Chapter 2, we will take the density n(r,6,t) and elec-
tric potential ¢(r,0,t) and divide them into steady-state quantities and perturbed

quantities

n(r,6,t) = no(r) + ni(r, 6, t) + no(r, 0, t) + ng(r,6,t) +--- (6.1)

ar, 0, t) = do(T) + dy(r, 9, t) + b2(r, 0, t) + 3(7, 8, t) rete (6.2)

The quantities no(r) and ¢o(r) are steady-state quantities with no angular depen-
dence. Subscript 1 corresponds to the perturbation due to the applied pulse 1 at time
t = 0, and is first order in amplitude. Subscript 2 corresponds to the perturbation
due to the applied pulse 2 at time ¢ = 7, and is first order in amplitude. Subscript 3
represents the echo perturbation from the nonlinear interaction between pulse 1 and
pulse 2, and is second order in amplitude. Higher order terms exist, which could be
included to form a theory on higher order echoes. In this argument, these higher
order terms will be neglected. The response for the first and second pulse will be con-
sidered to first order in amplitude, and the echo will be considered to second order
in amplitude.

Equation 6.1 and equation 6.2 for the density and potential are inserted into the
equation 3.2 for the time derivative of the density in cylindrical coordinates. Keeping

terms to first order for pulse 1 and pulse 2, and terms to second order for the echo

results in
where

_ 1,On, O¢2 _ Oni Ob, , On2IG, _ Anz Og,
S3(",8,t) = "(aap — a9 ar + Or 00. OO Or)” (6.6)

Equation 6.3 and 6.4 are the same linearizations discussed in Chapter 2 and de-
scribe two independent small amplitude pulses that will undergo mixing and decay
exponentially in time. Equation 6.5 describes a third echo response that will undergo
mixing and decay exponentially in time, with an additional source term S3 from the
nonlinear interaction of pulse 1 and pulse 2.

If the time between the two applied pulses 7 is much greater than the decay time
of the response to these pulses, it can be shown that the first term of equation 6.6 will
dominate the other terms. This is because ¢,, the perturbed potential due to pulse 1
decays exponentially in time, but n,, the perturbed density due to pulse 1 has a term
that persists and does not decay in time. Additionally, the radial derivative of the
perturbed density of pulse 1, ora. has a term that is proportional to t, and so it will
dominate for a large time separation between applied pulses. Other terms exist, but
these smaller terms will not be considered here. For + large compared to the decay
time, the dominant component of the source term may be written as

1 Oni 0

Thus, source term is proportional to the radial density derivative of pulse 1, and
the angular potential derivative of pulse 2. The real form for the perturbations rather
than the complex exponential notation will be used since the source term involves the
multiplication of two perturbed quantities.

We first give an approximate expression for the density perturbation caused by
the first pulse. Pulse 1 is applied to the wall at time t = 0, and is approximated
by V; cos(m,6)6(t), where V, is the amplitude of pulse 1, and 6(t) is the Dirac delta
function. Treating the applied potentials as a delta function is a convenient approx-
imation which is not too bad if the pulse separation is much greater than the decay
times. At time t = 0, the impulse at the wall causes a density perturbation which

begins to mix. The time dependent form of the density perturbation is given by

ny(r,6,t) = ViD-.(r) cos[m1(9 — wo(r)t)] +--- (6.8)

96
where only the term that does not decay in time is considered. The variable D,(r)
describes the radial dependence of the perturbed density, and its exact form is not
important for this discussion. This phase mixing process for the perturbed density
was described in Chapter 4, and is due to the radial dependence of the angular
velocity. Pulse 2 is applied to the wall at time t = 7, and will be approximated by
V2, cos(m26)6(t — 7), where V2 is the amplitude of pulse 2. The resulting perturbed

potential is given by

do(r, 0, t) = Va Ae(r) cos(m29)6(t — 7). (6.9)

The radial dependence of the perturbed potential is included in Ao(r). The exact
form of A2(r) is also not important for this discussion. The interaction between
pulse 1 and pulse 2 is assumed to occur at the time t = 7 for simplicity. After this
interaction, the plasma evolves in time by the phase mixing process. The source term
is proportional to the derivatives of n, and $5. The radial derivative of n; is given by

On,

an Vir D(r)m4

sin[(m,(@ — wo(r)t)] +--- (6.10)

where only the most important term is shown. The angular derivative of ¢, is given

a = V2Ao(r)mz sin(m29)6(t — 7). (6.11)
Inserting equation 6.10 and equation 6.11 into equation 6.7, the source term $3 is

approximated by

Owo(r) 1
Or 2
— cos[(m, + m2)O — mwo(r)t] }6(t — 7).

S3(r,0,7T) = ~ViVar De(r) Aa(r)mime cos[(m, — m2)8 — mywo(r)E]

The source term $3 exists only at time t = 7. This term creates a density pertur-

bation n3 which evolves in time due to phase mixing, and is given by

Owo(r)
Or

5 {cosl(m, — mM )0 — Mmywo(r)T — (m1 — M2)wo(r)(t — 7)|

n3(r,0,t) = “Vi Var De(r)Aa(r)mmyrtg

— cos[(m + m2)9 — mywo(r)T — (mM + M2)wo(r)(t — T)]f

for time t > 7. This is the resulting density perturbation which generates the echo,
due to the nonlinear interaction of pulse 1 and pulse 2. After the interaction at time
t =7T, the plasma phase mixes due to the radial dependence of the angular velocity.
The next step is to examine the time dependent phase for the two terms in n3. It
is expected that the echo will reach its maximum amplitude when the phase is inde-
pendent of radius and the perturbation is aligned to give the maximum contribution
to the electric field at the wall. When the perturbation is aligned, the second applied
pulse has partially unmixed the first pulse. The time dependent phase is independent

of radius when
[—myT — mt + met + MT — MeT|wWo(r) = 0
for the cos[(m — m2)6] terms and
[—myT — mt — met + mT + MgT|wWo(r) = 0
for the cos[(m + mz2)6] terms. Solving the equation for the m, — mz terms for ¢ leads

mg
t = —————-T.
mz — ™y4

In this case, the echo occurs after the second pulse as long as mz > mj, and has
azimuthal symmetry m2 — m,. At this time, the density perturbations at different

radii are aligned, and the echo has a maximum contribution to the perturbed electric

98
field at the wall. If m,; > me, then the echo would occur before both pulses. This
violates causality, and therefore an echo does not occur when m, > mz.

Solving the equation for the m, + mz terms for t leads to

mag

= —T.
Mz +m

This is a non-physical result, since this equation states that the echo would occur
at a time t < 7, since mg +m, > mg. The echo obviously can not occur before the
second pulse. We can conclude that the m2 +m, terms do not contribute to form an
echo.

An echo can only occur when mz > m,. Since only the m; — mz terms contribute,

the density perturbation term that generates the echo is given by

M Mz Owo(r)
2r Or
+muo(r)T — (mz — m1)wo(r)t}.

n3(r,0,t) = VyVerD,(r)Ao(r)

cos[(™mz — ™1)@ - (6.12)

The density perturbation evolves in time from phase mixing, with different radial
layers rotating at different rates. At time t = —™—7, the density perturbations at
different radial layers will align, producing a coherent electric field in the plasma that
can be measured at the wall. In equation 6.12, the perturbed density is proportional
to Vi,V2, and r. Therefore, the echo is proportional to the amplitude of pulse 1, the
amplitude of pulse 2, and the interpulse spacing T.

To summarize, there are six important conclusions from this argument describing
the echo.

eAn echo can only occur if the azimuthal symmetry number of the second pulse
is greater than the azimuthal symmetry number of the first pulse, i.e., mz > m1.

eThe echo will have azimuthal symmetry m3 = m2 — ™}.

eThe echo will reach a maximum amplitude at the time ¢t = mom"

eThe echo amplitude is proportional to amplitude of the first applied pulse, V.

eThe echo amplitude is proportional to amplitude of the second applied pulse, V2.

99
eThe echo amplitude is proportional to the time between the two pulses, T.
The next section will look at results from the Nonlinear Fluid Code for two applied
pulses. These results will be compared to the results of the argument presented in

this section.

6.3 Echo Results

The previous argument suggests that an echo should occur if two pulses are applied
to a pure electron plasma, separated in time, with the azimuthal symmetry number of
the second pulse greater than for the azimuthal symmetry number of the first pulse.
The Nonlinear Fluid Code is used to study the echo phenomena, and confirm the
results of the argument presented in the previous section.

For these computations, the Fast Decay Profile, no(r) = no(0){1—3(£)? +3(£)* —
(£)°}, with a = 0.8, was chosen. The reason for this choice is that the perturbed
electric field at the wall due to an applied perturbation will quickly decay away in
time. Therefore, the response to the first pulse and the response to the second pulse
will have essentially vanished at the time of the echo, allowing the echo to be clearly

observed.

6.3.1 Demonstration of the Echo

This section presents the results of the Nonlinear Fluid Code for a computation that
satisfies the requirements to produce an echo. Two pulses, separated in time, are
applied to the plasma, with the azimuthal symmetry number of the second pulse
greater than the azimuthal symmetry number of the first pulse.

A one cycle sinusoidal pulse with a frequency of w, = 0.62w 9(0) is applied at time
t = 0 with m, = 2 azimuthal symmetry and an amplitude of 1.0 x 10734). A second
one cycle sinusoidal pulse at twice this frequency is applied at time T = 25¢,, with
mg = 4 azimuthal symmetry and an amplitude of 2.0 x 10~%¢,. An echo is expected

to occur at time t = mT = 755 (25ter) = 50t., and have mg = Myg—m, = 4-2 = 2

100

azimuthal symmetry. The resulting perturbed electric field at the wall is shown in

figure 6.1.
0.1 v t T T T T T T T T T T T T
—_ - +
\wenee” r 4
ue)
=|
"o
fy
3 |
oO 0
Le |
wor
oO
= L
o |
oO
ae
-0.1 1 1 1 1 I L 1 L 1 l 1 L L L

0 25 50 75
Time (Central Rotation Periods)

Figure 6.1: Plot of the perturbed electric field at the wall versus time for a calculation
which demonstrates the echo. A one cycle pulse with m, = 2 azimuthal symmetry
and amplitude 1.0 x 10734, at time t = 0 is applied to the Fast Decay Profile. A
second one cycle pulse with m2 = 4 azimuthal symmetry and amplitude 2.0 x 1073¢,
is applied at time 7 = 25t,,. An echo then occurs at time t = 50t,, with m3 = 2
azimuthal symmetry, obtaining a maximum amplitude of 0.0192% of the unperturbed
electric field at the wall.

As seen in figure 6.1, the first m; = 2 pulse at t = 0 excites a response which decays
exponentially. The second mz = 4 pulse at time T = 25t,, also excites a response
which decays exponentially. A nonlinear interaction between these two applied pulses
causes an m3 = 2 echo to occur at time t = 50t,,. The echo reaches a maximum

amplitude of 0.0192% of the unperturbed electric field at the wall.

101
As a second example, a one cycle sinusoidal pulse with a frequency of w, =
0.62w (0) is applied at time t = 0 with m; = 2 azimuthal symmetry and an amplitude
of 1.0 x 107345. A second one cycle sinusoidal pulse at three times this frequency
is applied at time 7 = 25t,, with mg = 6 azimuthal symmetry and an amplitude
of 3.0 x 10~%¢,. For these initial conditions, an echo is expected to occur at time
t = 37.5ter with azimuthal symmetry m3 = 4. The resulting perturbed electric field

at the wall is shown in figure 6.2.

0.1 T T T T T T T T T T T T T T
—_ nu |
\weee” rc 4
aS)
[oy
4 0 Ww aT PaavAVAVAN:
aa) i 1
ne
oD
= L J
o L |
- 1 1 1 L L 1 1 L 1 l 1 L 1 1
O15 25 50 75

Time (Central Rotation Periods)

Figure 6.2: Plot of the perturbed electric field at the wall versus time for a one cycle
m, = 2 pulse with amplitude 1.0 x 10~%¢ at time t = 0, and a one cycle mz = 6
pulse with amplitude 3.0 x 10-°¢, at time 7 = 25t., applied to the Fast Decay
Profile. An echo occurs at time t = 37.5t,, with m3 = 4 azimuthal symmetry and
has a maximum amplitude of 0.00347% of the unperturbed electric field at the wall.
Near time t = 75t,,, a higher order echo occurs.

As seen in figure 6.2, an echo occurs at time t = 37.5t,,. The echo has the expected

102
m3 = 4 azimuthal symmetry, and reaches a maximum amplitude of 0.00347% of the
unperturbed electric field at the wall.

There is also another echo occurring at a later time t = 75t,,. This echo is actually
a higher order echo of an echo, where the m3 = 4 echo is interacting with the m, = 2
pulse at ¢ = 0 to produce an m = 2 echo of an echo at time t = 75t,,. Higher order
effects will be discussed in more detail in a later section of this chapter.

The Nonlinear Fluid Code successfully demonstrates the phenomena of the echo.
The echo has the correct azimuthal symmetry, m3 = m2 — m1, and occurs at time
t= maT) in agreement with the argument of section 5.2. The next section provides
a set of control cases, verifying that two pulses are required to produce and echo, with

the mode number of the second pulse higher than the mode number of the first pulse.

6.3.2 Control Sets

This section provides the results from several control cases. These sets verify that
the observed echo is not a numerical artifact, and that the echo only occurs when
the azimuthal mode number of pulse 2 is higher than the azimuthal mode number of

pulse 1.

Single Pulses

This section verifies that the response to a single applied pulse decays exponentially.
For the amplitudes used in this section, no nonlinear effects such as bounce oscillations
are observed. This is verified for the m = 2, m = 4, and m = 6 pulses that were used
in the previous section.

A one cycle sinusoidal pulse with a frequency of w, = 0.62w9(0) is applied at time
t = 0 with m = 2 azimuthal symmetry. The pulse is applied at the m = 2 resonant
frequency for this profile, with an amplitude of 1.0 x 10~%¢, and a pulse length of
1.613¢,,. The resulting perturbed electric field at the wall is plotted in figure 6.3.

This result is similar to the result of the previous chapter for this profile, shown

in figure 4.2, but at a higher applied amplitude. The maximum perturbed electric

103

0.1 ——— rr
~ m= w,= 0.610,(0)
xe y= 0.050a,(0) 7
fa} _

a w,/y = 12.1

uc)

—_

Ay
- ! 1 Ll 1 1 ! 1 ni 1 1
0.1 25 50 75

Time (Central Rotation Periods)

Figure 6.3: Perturbed electric field at the wall versus time for a one cycle pulse with
m = 2 azimuthal symmetry applied to the Fast Decay Profile. The pulse has an
applied amplitude of 1.0 x 10~%¢,, and the perturbed electric field at the wall reaches
a maximum amplitude of 0.0916% of the unperturbed electric field at the wall. The
plasma response decays exponentially.

field at the wall is found to be 0.0916% of the unperturbed electric field at the wall.
The response to this applied pulse quickly decays in time.

A similar computation was performed using an applied pulse with m = 4 sym-
metry. The m = 4 pulse is applied at time t = O for one cycle at a frequency of
w = 1.24w9(0). The amplitude of this pulse is 2.0 x 107345, with a pulse length of
0.806¢,,. The resulting perturbed electric field at the wall for the m = 4 pulse reaches
a maximum amplitude of 0.0188%, and is shown in figure 6.4. Note that the scale for
this figure has been increased by a factor of 5 compared to the figure for the m = 2
pulse.

For the m = 6 pulse, the applied frequency is w = 1.86w (0), leading to a one
cycle pulse length of 0.538¢,,. The amplitude of the m = 6 pulse is 3.0 x 10-°¢p. The

resulting perturbed electric field at the wall for the m = 6 pulse reaches a maximum

104

0.02 v T T T T T T T T T T T
ce mm an 20100)
x Y HULU wy
S. w,ly = 6.5
MS}
ac}
fy
a) 0
ea
er]
0)
Ay
. | \ ! \ ! | ! \ A \
0.02, 25 50 75

Time (Central Rotation Periods)

Figure 6.4: Perturbed electric field at the wall versus time for a one cycle pulse with
m = 4 azimuthal symmetry applied to the Fast Decay Profile. The pulse has an
applied amplitude of 2.0 x 107¢9, and the perturbed electric field at the wall reaches
a maximum amplitude of 0.0188% of the unperturbed electric field at the wall. The
plasma response decays exponentially.

amplitude of 0.00509% and is shown in figure 6.5. In this figure, the scale has been
increased by a factor of 20 compared to the figure for the m = 2 pulse.

The frequencies used for the m = 4 and m = 6 pulse are not the exact resonant
frequencies, but rather twice and three times the resonant frequency for the m = 2
pulse. However, since only a one cycle pulse is used, the frequency spectrum is quite
broad, so it is unnecessary to apply the exact resonant frequency. Additionally, the
applied amplitude for the m = 4 and m = 6 pulse has been increased by a factor of two
and three, respectively, to compensate for the shorter pulse length, as compared to the
m = 2 pulse. Thus the total energy applied to the plasma in all three computations
should be equal.

As can be seen from the figures, the single pulses behave as expected. The per-

turbed electric field at the wall rapidly decays away for each pulse, consistent with

105

0.005 —~—>—+—+ —} + — + + +> + 1+ —

~ m=6 w= 1.90,(0)
=0.2640,(0) |

x " 3

a wry = 7.2

3 ,

2 0

eal

7 Ll 1 ' A 1 l 1 A 1 \

0.0055 25 50 75

Time (Central Rotation Periods)

Figure 6.5: Perturbed electric field at the wall versus time for a one cycle pulse with
m = 6 azimuthal symmetry applied to the Fast Decay Profile. The pulse has an
applied amplitude of 3.0 x 10~°¢9, and the perturbed electric field at the wall reaches
a maximum amplitude of 0.00509% of the unperturbed electric field at the wall. The
plasma response decays exponentially.

the results of Chapter 4. No nonlinear phenomena such as bounce oscillations are
observed at these amplitudes.

It should be noted that the resulting perturbed electric field at the wall is lower
for the m = 4 pulse and m = 6 pulse, even though the total energy for each pulse is
the same. The reason for this discrepancy is that for the applied pulse, the perturbed
potential at the resonant radius r, is proportional to (== )™ where r,, is the wall radius.
Thus for higher m values, the field falls off more rapidly as the plasma is further away
from the wall sectors than for lower m values. Additionally, the resulting density
perturbations at the resonant radius create a perturbed electric field at the wall that
is also proportional to (7+). Therefore, the perturbed electric field at the wall will
be smaller for higher m valued density perturbations. The net effect is that an applied

pulse of voltage at the wall will induce a response in the perturbed electric field at

106

the wall which will have an amplitude response that is proportional to (£*)?”. Thus

r.
rw

higher m valued applied pulses will have a lesser contribution in the perturbed electric

field at the wall.

Two Pulses with m,> mz

An echo can only occur if the azimuthal mode number of the second pulse is higher
than the azimuthal mode number of the first pulse. This section searches for a counter
example by applying an m, = 4 pulse followed by an mz = 2 pulse. This computation
uses the same pulses that were used in the echo computation above, but with the order
of the two pulses reversed.

A one cycle pulse is applied at time t = 0 with m, = 4 symmetry and an amplitude
of 2.0 x 10-3¢9. A second one cycle pulse is applied at time tT = 25t., with m2 = 2
symmetry and an amplitude of 1.0 x 107¢ 9. The resulting perturbed electric field at
the wall is shown in figure 6.6.

As seen in figure 6.6, a nonlinear interaction between these two pulses is not
observed. The first pulse is applied, and it decays away. The second pulse is applied
later, and it also decays away. No echo is generated. This result is consistent with
the argument from section 5.2. The first pulse has an azimuthal symmetry eigenvalue
greater than the azimuthal symmetry eigenvalue of the second pulse. In the argument
presented in the previous section, no echo can occur if m, > mz.

To summarize the results of this section, these computations using the Nonlinear
Fluid Code verify the existence of t:1e echo in non-neutral plasmas. The echo is found
to only occur when mz > my. Adilitionally, the echo is found to have azimuthal
symmetry of mg = m2 — mj, and occur at a time t = mom T" These results are
consistent with the argument presented in section 5.2. The next section will explore

the dependence of the echo amplitude on the amplitudes of the applied pulses and

the interpulse separation.

107

0.1 ' T T T T T T T T T T T T T
~~ s
Nae 4
fy
ond
8 oH
ea I 1
mer)
oO
= L J
% L ]
oO
ae

-0.1 i 1 L 1 | L 1 1 1 l L L L L

25 50 75
Time (Central Rotation Periods)

Figure 6.6: Perturbed electric field at wall for a one cycle m; = 4 pulse at time t = 0,
and a one cycle mz = 2 pulse at time 7 = 25t,, applied to the Fast Decay Profile.
Both pulses decay away, with no echo produced.

6.3.3. Amplitude of the Echo
Dependence on First and Second Pulse Amplitudes

The previous section verified the existence of the echo, as well as the azimuthal
symmetry of the echo and the time the echo occurs. This section will explore the
amplitude of the echo in relation to the amplitude of the two applied pulses. For the
following computations, a one cycle pulse with m, = 2 azimuthal symmetry is applied
at time ¢ = 0, and a one cycle pulse with mz = 4 azimuthal symmetry is applied at
time T = 25t,,. The amplitude of the two pulses is varied, and the resulting change

in the echo amplitude is observed. The results are summarized in table 6.1.

108

Pulse 1 Amplitude | Pulse 2 Amplitude | Echo Amplitude
1.0 x 10-745 2.0 x 10°-%h 0.0192%
1.0 x 10~°d, 1.0 x 107%¢5 0.00972%
5.0 x 107*¢5 2.0 x 107° 0.00959%
5.0 x 10-4, 1.0 x 107°¢5 0.00485%

Table 6.1: Echo amplitude versus applied pulse 1 and pulse 2 amplitude. Pulse 1 is a
one cycle pulse with m, = 2 symmetry and is applied at time t = 0. Pulse 2 is a one
cycle pulse that has mz = 4 symmetry and is applied at time T = 25t,,. The echo
has m3 = 2 symmetry and occurs at time t = 50t,,. The echo amplitude is given as
a percentage of the unperturbed electric field at the wall.

With the amplitude of the second pulse halved, the echo amplitude is halved.
For the amplitude of the first pulse halved, the echo amplitude is halved. With the
amplitude of both applied pulses halved, the echo is reduced by a factor of 4. This
is clearly illustrated by figure 6.7, which plots the perturbed electric field at the wall
for these computations. The echo is proportional to the amplitude of the first pulse
multiplied by the amplitude of the second pulse.

The amplitude of the echo is systematically studied as a function of the first
m, = 2 pulse amplitude, occurring at time t = 0. The amplitude of the m2 = 4 pulse
at time T = 25t,, is held fixed at 2.0 x 107%¢9. The echo amplitude versus applied
pulse 1 amplitude is plotted on a log-log scale in figure 6.8.

The points of figure 6.8 lie along a straight line of slope one. This indicates that
the echo amplitude is proportional to the amplitude of the first pulse.

The amplitude of the echo is also systematically studied as a function of the second
mz = 4 pulse amplitude, occurring at time T = 25t,,. The amplitude of the m, = 2
pulse at time ¢ = 0 is held fixed at 1.0 x 10~°¢9. The echo amplitude versus applied
pulse 2 amplitude is plotted on a log-log scale in figure 6.9.

The points of figure 6.9 lie along a straight line of slope one. This indicates that
the echo amplitude is proportional to the amplitude of the second pulse.

The amplitude of the echo is proportional to the amplitude of the first applied
pulse, multiplied by the amplitude of the second applied pulse. This amplitude de-
pendence of the echo agrees with the results predicted by the argument of section

5.2. The next section will study the amplitude dependence of the echo as a function

109

cms

co)

jan

‘SS Second

pulse

§ halved

RP L

owt

Poy

ae)

eal

ae]

8 Fst |

5 pulse

fy halved
Both
pulses Mn ni\\f\Wornre >
halved

: f 1 ‘ f ! : : : 1 ! f \ 1 \
0 25 50 75

Time (Central Rotation Periods)

Figure 6.7: Perturbed electric field at wall for four echo computations varying in
applied amplitude. The second trace has the amplitude of the second pulse halved.
The third trace has the amplitude of the first pulse halved. In the fourth trace, the
amplitude of both applied pulses are halved. The echo amplitude is found to be
proportional to the amplitude of the first pulse, and proportional to the amplitude of
the second pulse.

110

0.1

Maximum Echo Amplitude
(% unperturbed elecric field at wall)
ss,

ed
0.001
0.00001 0.0001 0.001

Applied m=2 Amplitude
(units of the central potential)

Figure 6.8: Plot of the maximum amplitude of the echo versus applied m, = 2
amplitude of the first pulse on a log-log scale. The m2 = 4 amplitude is fixed at
2.0 x 10-8. The first pulse is excited with a one cycle burst at time t = 0, and
the second pulse is excited with a one cycle burst at time 7 = 25t,,. The solid line
is a best fit with slope 1, indicating that the echo amplitude is proportional to the
amplitude of the first pulse.

of the time between the two applied pulses.

Dependence on Interpulse Spacing in Time

The previous section found the amplitude of the echo to be proportional to the am-
plitude of the first applied pulse, and proportional to the amplitude of the second
applied pulse. This section will examine the amplitude of the echo with respect to
the time between the two applied pulses, 7. For the following computations, a one
cycle pulse with m, = 2 azimuthal symmetry is applied at time ¢t = 0 with an ampli-
tude of 1.0 x 1073¢). A second one cycle pulse with mz = 4 symmetry is applied at

a later time 7, with an amplitude of 2.0 x 10~°¢9. The time 7 for the second pulse is

111

0.1
BS
23
5 I
2 6
5S 3 oo 4
] es ;
Es Z
E¢
qo Li
Ss 5
0.001
0.0001 0.001 0.01
Applied m=4 Amplitude

(units of the central potential)

Figure 6.9: Plot of the maximum amplitude of the echo versus applied mz = 4
amplitude of the second pulse on a log-log scale. The m; = 2 amplitude is fixed at
1.0 x 10-°¢9. The first pulse is excited with a one cycle burst at time t = 0, and
the second pulse is excited with a one cycle burst at time T = 25t,,. The solid line
is a best fit with slope 1, indicating that the echo amplitude is proportional to the
amplitude of the second pulse.

varied, and the resulting change in the echo amplitude is observed. The results are
summarized in table 6.2.

As the time between pulses is increased, the amplitude of the echo increases. The
perturbed electric field at the wall versus time for three cases is plotted in figure 6.10.
The echo amplitude is found to be proportional to the temporal interpulse spacing,
T. This result is also consistent with the echo argument presented above. The echo
amplitude results of table 6.2 are plotted on a log-log scale in figure 6.11.

The straight line has slope one. For the first several points, the echo amplitude
is proportional to the interpulse spacing. For the longer times between pulses, some

saturation of the echo amplitude occurs. Saturation effects will be discussed in the

112

Perturbed Electric Field (each tick represents 0.02%)

I»

0 40 30. ~~120
Time (Central Rotation Periods)

Figure 6.10: Perturbed electric field at wall for three echo computations, varying the
interpulse spacing, t. The first trace has half the interpulse spacing of the second.
The third trace has twice the interpulse spacing of the second. The echo amplitude
is found to be proportional to the interpulse spacing.

113

Interpulse Spacing 7 | Echo Amplitude
12.5te, 0.00932%
17.7tep 0.0137%
25ter 0.0192%
39.4 bop 0.0261%
50ter 0.0309%

Table 6.2: Echo amplitude as a function of interpulse spacing. A one cycle pulse with
m, = 2 symmetry and amplitude 1.0 x 107°¢, is applied at time t = 0. A second
one cycle pulse with mz = 4 symmetry and amplitude 2.0 x 10~%¢, is applied at time
t = T. The echo amplitude is measured as a percentage of the unperturbed electric
field at the wall.

next section.

To summarize this section, for small amplitude pulses, the amplitude of the echo
is proportional to the amplitude of the first pulse, multiplied by the amplitude of
the second pulse. The echo amplitude is also proportional to the time between the
two applied pulses. These results are consistent with the echo argument presented in
section 5.2. The next section will describe the effect on the echo when larger amplitude

pulses are applied to the plasma, and higher order terms become important.

6.3.4 Higher Amplitude Effects

Saturation at Higher Amplitudes

The previous section examined the effect of two small pulses interacting. In this low
amplitude limit, the echo was found to be proportional to the amplitude of the first
pulse multiplied by the amplitude of the second pulse, and also proportional to the
interpulse spacing. This section will look at the effect larger amplitude pulses have
on the echo. The effect from larger amplitude pulses is not contained in the echo
argument presented in section 5.2 because higher order terms were neglected. The
Nonlinear Fluid Code contains all the higher order terms, and therefore is not limited
to studying low amplitude effects. For large applied pulses, the amplitude of the echo
can saturate.

The amplitude of the echo was systematically studied as a function of pulse 1

114

0.1
3 8
Bs
as
5 a
o o
ay
=e vs
E¢ ZL
pa &
Ss 5
0.01-—]
10 ; 100
Time Between Pulses
(central rotation periods)

Figure 6.11: Plot of the amplitude of the echo versus interpulse spacing of the applied
pulses on a log-log scale. The one cycle m; = 2 pulse has a fixed amplitude of
1.0x 10-3¢,, and the one cycle mz = 4 pulse has a fixed amplitude of 2.0x 1073¢9. The
straight line has slope one. For the first few points, the echo amplitude is proportional
to the interpulse spacing. For longer spacings, the echo amplitude saturates.

amplitude and pulse 2 amplitude. Pulse 1 is excited by a one cycle sinusoidal burst
with m, = 2 azimuthal symmetry at time t = 0. Pulse 2 is excited by a one cycle
sinusoidal burst with mz = 4 symmetry at time 7 = 25¢,,. An echo with azimuthal
symmetry m3 = 2 occurs at time t = 50t,,.. The amplitude of pulse 1 and the
amplitude of pulse 2 were varied together.

The lowest amplitude case has a one cycle m; = 2 applied pulse amplitude of
1.25 x 10-*¢5, and a one cycle mz = 4 applied pulse amplitude of 2.5 x 10~*¢.
The echo reaches a maximum perturbed amplitude of 0.000305% of the unperturbed
electric field at the wall. For each successive case, the applied pulse amplitude of

the m,; = 2 pulse and the applied pulse amplitude of the mz = 4 pulse was varied

115
together, increased by a factor of /2. The applied amplitude of the mz = 4 pulse
always has twice the applied amplitude and half of the pulse duration of the m, = 2
pulse.

The response of the electric field at the wall for the first pulse is proportional to
the applied amplitude of the first pulse. The response of the electric field at the wall
for the second pulse is proportional to the applied amplitude of the second pulse. For
the lower amplitude pulses, the echo amplitude is proportional to the amplitude of
the first pulse multiplied by the amplitude of the second pulse. For higher amplitudes,
the echo saturates. These results are summarized in figure 6.12. The response of the
first pulse, response of the second pulse, and the echo response are plotted versus
applied pulse 1 amplitude (also applied pulse 2 amplitude divided by 2) on a log-log
scale in figure 6.12.

On this plot it is possible to see the heavy saturation occurring at the higher
amplitudes. The slope is one for the response of the m,; = 2 pulse and for the
mz = 4 pulse, indicating that the response amplitude is still in the linear range for
each individual pulse. For the lower amplitude points, the slope for the echo response
is two, illustrating that the echo is proportional to the amplitude of the first pulse
multiplied by the amplitude of the second pulse. For the larger applied amplitudes,
higher order effects are dominating, limiting the response amplitude of the echo. The

next section will examine higher order echoes that can occur at larger amplitudes.

Higher Order Echoes

The approximate argument in section 5.2 describes the echo phenomena for small
applied pulses. At larger amplitudes, higher order effects become important. Specif-
ically, harmonics of the applied pulses can cause higher order echoes.

If we describe the harmonics with the harmonic number h, where h = 1 is the
fundamental, h = 2 is the second harmonic, etc., the equations describing the property
of the higher order echo are as follows. The higher order echo will have azimuthal

symmetry described by

116

= os
re aos
| = a1 ia fe a
sara 4 - se Fo f- a
oO — i ew. ak
url a Va i
[i141 = Applied ee a
Ow m=2 :
OQ Ps 0.01 = en ro
| Soa
fr € a 2
cm a f
an) pphed WA
=] —
5 & 0.001 m4 F
ag
2 7
m=?
echo
0.0001
6.0001 0.001 0.01

Apphed Pulsel and Pulse2
(units of central potential)

Figure 6.12: Response of the m; = 2 pulse, mz = 4 pulse, and m3 = 2 echo response
as a function of pulse 1 amplitude (and pulse 2 amplitude divided by 2). For low
amplitude pulses, the echo amplitude is proportional to the amplitude of pulse 1 mul-
tiplied by the amplitude of pulse 2. For larger pulses, the echo amplitude saturates.
The applied amplitudes are varied together, with both amplitudes increased by a
factor of \/2 from one point to the next. The applied amplitude of the m2 = 4 pulse
always has twice the applied amplitude and half of the pulse duration of the m; = 2
pulse. The dashed line has a slope of 2, indicating points that are proportional to the
pulse 1 amplitude multiplied by the pulse 2 amplitude.

117

m3 = home —_ hym4

and will occur at time

i- hgma
~ hgmg — hym,

where h is the harmonic number for pulse 1, and hz is the harmonic number for pulse

2. The echo will have an amplitude dependence
Ag ~ Am Abr

where A3 is the amplitude of the higher order echo, and A; and Ag represent the
amplitudes of pulse 1 and pulse 2. For large amplitudes, these higher order echoes
are visible in the results from the Nonlinear Fluid Code.

Figure 6.13 is a plot of the perturbed electric field at the wall for a set of higher
amplitude applied pulses. The gap from the plasma edge to the wall was changed
for this computation, with a = 0.95, so the higher m values could be more easily
observed. The first one cycle m, = 2 pulse occurs at time t = 0 with an amplitude of
4.95 x 10-%¢9. The second pulse has m2 = 4 symmetry, an amplitude of 9.90 x 107* dp,
and occurs at time 7 = 17.8t,,.. An echo occurs with m3 = 2 symmetry at time
t = 35.6t,,, and reaches a maximum amplitude of 0.384% of the unperturbed electric
field at the wall.

An additional echo is also visible in this trace located between the second pulse at
7 =17.8t,, and the m3 = 2 echo at time t = 35.6t,,. This echo has m = 6 azimuthal
symmetry and occurs at time t = 23.7t.,. This response is a higher order echo caused
by the interaction of the fundamental of the first m,; = 2 pulse, h; = 1, and the
second harmonic of the mz = 4 pulse, hy = 2. This interaction produces an echo
with symmetry m3 = hem2 — him, = 2:4—1-2 = 6, that should occur at time

_ heme _ __ 2-4 _ 4,
t= jgme ami? = B4-137 = 37 = 23. Ter.

118

1% ' q T T T T T T T T T T T T

Perturbed Electric Field ( % )

15 30 45
Time (Central Rotation Periods)

Figure 6.13: Perturbed electric field at the wall versus time for two high amplitude
applied pulses. A one cycle m; = 2 pulse is applied with amplitude 4.95 x 1073¢, at
time t = 0, and a one cycle mz = 4 pulse is applied with amplitude 9.90 x 10734)
at time 7 = 17.8t,,. An m3 = 2 echo occurs at time t = 35.6t,, while a higher order
m3 = 6 echo occurs at time t = 23.7t,,.

Further higher order echoes are also observable. The Fourier components of the
envelope of the perturbed electric field at the wall is plotted on a semilog scale in
figure 6.14.

In this plot, the initial m; = 2 pulse, m2 = 4 pulse, and corresponding m3 = 2
echo have the highest amplitude response in the perturbed electric field at the wall.
Additionally, an hy = 3 and h; = 2 harmonic are visible from the first m; = 2 pulse
at time t = 0, and an hy = 2 harmonic is visible from the second mz = 4 pulse at
T=1728.

These harmonics can interact to produce further echoes. The m3 = 6 echo at time

119

{a4 spl pulse

vessel e é ped ed moe

Log of Perturbed Electric Field
(each tick is a factor of 2)

Time (Central Rotation Periods)

Figure 6.14: Semilog plot of a computation demonstrating higher order echoes. An
m, = 2 pulse is applied at time t = 0, and an m, = 4 pulse is applied at time
7 = 17.8t,.,. Harmonics of these two pulses, and corresponding echoes are observed.

t = 23.7t,, is produced from the initial m, = 2 pulse fundamental h, = 1 interacting
with the hg = 2 harmonic of the second pulse. The m3 = 10 echo at time t = 21.4¢,,
is caused by the initial m, = 2 pulse fundamental h, = 1 interacting with the hg = 3
harmonic of the second pulse.

Echoes of echoes can also occur. The mg = 8 response at time t = 26.8¢,, is
caused by the initial m, = 2 applied pulse fundamental h; = 1 interacting with the
m = 10 echo at time t = 21.4t,,. The m3 = 4 response at time t = 35.6t,, is caused
by the initial m, = 2 applied pulse fundamental h, = 1 interacting with the m = 6
echo at time t = 23.7t¢,.

The idea of echoes of echoes explains the response occurring near time t = 75t¢,r

120
back in figure 6.2. This computation has been computed further in time to 90 central

rotation periods, and replotted in figure 6.15.

0.1 v T v T T T T T T T T T T T
oo L 4
BS
Neue” r 4
aS)
Y ;
-)
g :

A mrt DIR

oO 0
an
er)
oO
= L 4
a L ]
oO
poy

- i 1 L 1 { 1 n 1 1 | l 1 L n

Oo 30 60 90

Time (Central Rotation Periods)

Figure 6.15: Plot of the perturbed electric field at the wall versus time for a one cycle
m1 = 2 applied pulse with amplitude 1.0 x 10~%¢) at time t = 0, and a one cycle
my = 6 applied pulse with amplitude 3.0 x 1074, at time 7 = 25t,,. An echo occurs
at time t = 37.5t., with m3 = 4 azimuthal symmetry. A higher order m = 2 echo
occurs at time t = 75.0t¢,.

The applied one cycle m, = 2 pulse occurs at time t = 0 with an amplitude of
1.0 x 1073, and the one cycle mz = 6 applied pulse occurs at time 7 = 25t,, with an
amplitude of 3.0 x 107345. The mg = 4 echo occurs at time t = 37.5t,,. The mg = 4
echo then interacts with the m; = 2 initial pulse to produce an m = 4 echo of an
echo at time t = 75.0t,,.

The results from the Nonlinear Fluid Code exhibit higher order effects that were

not contained in the quasi-linear argument presented in section 5.2. At higher ampli-

121
tudes, the amplitude of the echo will saturate. Higher order echoes due to harmonics
of the applied pulses are also evident and echoes of echoes occur. The next section

will briefly discuss the effects of using a different density profile.

6.3.5 Slow Decay Density Profile

The echo has been clearly demonstrated using the Fast Decay Profile, no(r) =
no(0){1 — 3(£)? + 3(£)* — (£)°} with a = 0.8. However, this profile may not be
a good representative for density profiles found in a real experiment. The Slow Decay
Profile, no(r) = no(0){1 — 6(£)4 +8(£)® — 3(£)8} with a = 0.8 has a linear decay rate
much closer to the value found in the experiment, and is expected to more closely
resemble an experimental density profile.

For a one cycle pulse with m; = 2 azimuthal symmetry applied at time t = 0 with
an amplitude of 2.14 x 10~*@ , and a second one cycle pulse with mz = 4 azimuthal
symmetry applied at time 7 = 23.7¢,, with an amplitude of 1.71 x 10~°¢,, an m3 = 2
echo occurs at time t = 47.4t,,. The response of the perturbed electric field at the
wall is plotted in figure 6.16.

As seen in this figure, an echo with mg = 2 azimuthal symmetry occurs at time
t = 47.4t,,, with a perturbed response in the electric field at the wall of 0.00798%
of the unperturbed electric field at the wall. Since the first pulse does not decay as
quickly for this profile, the echo is not as clearly observed as in the computations
using the Fast Decay Profile, but is still evident at time t = 47.4¢,,.. The first pulse
amplitude was chosen so that the response to the m, = 2 pulse decays exponentially
without observable nonlinear effects. The second mz = 4 pulse was chosen with an
applied amplitude 8 times higher than the applied amplitude of the first m,; = 2
pulse. The previous results involving the time and amplitude dependence of the echo
using the Fast Decay Profile also hold for the Slow Decay Profile. Since an echo is
also observed using this profile, it should be possible to produce an echo in any profile

that has perturbations that decay.

122

0.05 q T T T T T T T T T T T T T

Perturbed Electric Field ( % )

-0 . 05 1 1 1 L l L i 1 1 l 1 1 1 1
0 20 AO 60

Time (Central Rotation Periods)

Figure 6.16: Plot of the perturbed electric field at the wall versus time for a com-
putation demonstrating the echo while using the Slow Decay Profile. A one cycle
m4 = 2 pulse with amplitude 2.14 x 10~*¢, is applied at time t = 0, and a one cycle
mg = 4 pulse with amplitude 1.71 x 1073¢, is applied at time 7 = 23.7t.,. An m3 = 2
echo occurs at time t = 47.4t,, and reaches a maximum amplitude of 0.00798% of the
unperturbed electric field at the wall.

6.4 Conclusion

A pure electron plasma under the fluid approximation can exhibit echoes. For an
echo to occur, an applied pulse with m, azimuthal symmetry must be followed by
a second applied pulse with a higher azimuthal symmetry m2. For the first pulse
occurring at time t = 0, and the second pulse at time t = 7, an echo will occur at
time t = —“—r with m3 = mz — m,; azimuthal symmetry.

mea-my

For small applied pulses, the amplitude of the echo is proportional to the amplitude

123

of the first pulse multiplied by the amplitude of the second pulse. The echo is also
proportional to the temporal interpulse spacing of the two pulses, 7. For larger applied
pulses, the amplitude of the echo will saturate. Also at higher amplitudes, echoes due
to harmonics of the two applied pulses and echoes of echoes become visible.

The echo can occur for a density profile which has perturbations that decay by
a collisionless phase mixing process. The second applied pulse partially unmixes the
density perturbations of the first pulse, producing an echo later in time. Echoes were
found to occur for a density profile that should be similar to a profile found in the
experiment. This suggests that an experimental search for echoes may be successful.

Nonideal fluid effects (transport phenomena) such as diffusion and viscosity are
expected to dissipate the echo at later times. If echoes are found experimentally,
measurement of the echo amplitude as a function of interpulse spacing may provide

a diagnostic to measure these additional effects.

124

Chapter 7 Beat-Wave Decay Instability

7.1 Introduction

This chapter presents preliminary results from the Nonlinear Fluid Code describing
a decay instability that occurs at very large amplitudes of the applied pulse. This is
quite a different phenomena from the results studied in the previous chapters, and
presents another example where the Nonlinear Fluid Code can be used to study a
phenomena that occurs in non-neutral plasmas.

In Chapter 5 we observed a modulation of the decay envelope caused by non-
linear fluid trapping. Although nonlinear trapping of fluid occurred, these plasma
perturbations were not visible in images of the total density. Only by subtracting
out the average density could the plasma dynamics be observed. This chapter will
discuss the response of the plasma to very high amplitude applied pulses. For these
high amplitude pulses, the perturbation to the plasma is visible in images of the total
density.

Such high amplitude perturbations in a pure electron plasma are subject to a decay
instability. Mitchell and Driscoll [22] performed an experiment which demonstrated
that a perturbation with azimuthal symmetry mode number m and angular frequency
Wm could decay into a mode m — 1 with angular frequency w,,_1, provided the beat
frequency Wm — Wm_—1 exists in a frequency range associated with a resonant radius
that is within the plasma. This decay is related to the three oscillator parametric
instability, and a theoretical discussion of this decay mechanism may be found in [36].

Density and electric potential may be decomposed into a Fourier series with
n(r,0) = > mm(r)e, and o(r,0) = >> ¢,,(r)e””’. The component of density n(r, #)

with azimuthal symmetry m is given by

125

2n
| n(r, Oe db

and the component of the electric potential ¢(r,@) with azimuthal symmetry m is
given by
1 20

Omir) = 52 for Bead.

A normalized density amplitude component A,, for a mode with azimuthal sym-

metry m was defined by Mitchell and Driscoll as

a2 = Jo" r rl dr (7.1)
So” r|no(r)|° dr
where no(r) is the unperturbed density profile and r,, is the radius of the wall. They
found the m— 1 mode to grow at the expense of the mode m, with the rate of change
given by

amt = (Im + Dynm—14m)Am=1 (7.2)

where 7/_1 is the linear decay rate for the mode m — 1, and [ym m_1 is the coupling
rate between mode m and mode m — 1 which is independent of amplitude A,, for
small values of Am. When Dnm—142, > Ym—1, expression 7.2 leads to exponential
growth.

Figure 7.1 is a reproduction of the main experimental results from Mitchell and
Driscoll. The figure plots the daughter mode growth rate 7,,_, versus the normalized
parent mode amplitude A,,, and shows the growth is approximately proportional to
Az,

Currently, the beat-wave decay instability is not well understood. Theoretical

calculations of the mode coupling rates [';,; do not agree well with the experiment

[37]. In order to address this problem, the Nonlinear Fluid Code was used to study

126

0 TUTitd t qT u Tri?
10 1 qv Lm my Li T T |
Pw -1 oOo,”
10 oO a 4
b a,
- 403 9” ’ °
—2 7 qd 7 “ °
—& 10 ® A’ 7
~ Me a
ro) _3 a& v4
2 10 ra 6
me 392 OW 4 21
oc
va 10—4 s® a
z 7h gy
& 5 7 o 4
10 é UZ
10-8 1 L Lil i Lt Lil n os co oo ee 8
10-2 1071 109

Amplitude A,,

Figure 7.1: Experimental results of Mitchell and Driscoll on the beat-wave decay
instability in a pure electron plasma. The measured growth rates 7,,_, for daughter
mode m — 1 versus normalized amplitude A,, of the parent mode, for m = 2, 3, and
4. The growth rates are scaled by the central rotation time t., = 5.54s. The dashed
lines indicate 7,1 = Pimjm—-1A2,-

the beat-wave decay instability, and determine if the numerical results agree with the

experiment. The Nonlinear Fluid Code results are presented in the next section.

7.2. Fluid Code Results

7.2.1 Choice of Density Profile

The density profile was chosen to approximate the experimental profile measured by
Mitchell and Driscoll [22]. The normalized Fermi function density profile with the
parameter (3 = 8, given by equation 4.1 in Chapter 4, is a close approximation to the

density profile measured in the experiment. This density profile and the corresponding

127

angular velocity profile are shown in figure 7.2. In Chapter 4, this density profile was
named the Beat Profile, and was found to have an m = 2 resonant frequency of w2 =
1.012w9(0) corresponding to a resonant radius of r, = 0.703r,, and a linear decay rate
of Yz = 0.00302w (0). A small amplitude computation with an m = 3 perturbation
was performed using this profile which resulted in an m = 3 resonant frequency
of w3 = 1.858w9(0) corresponding to a resonant radius of r, = 0.635r,, and a linear
decay rate of 73 = 0.0446w9(0). The m = 1 resonant frequency is given by the angular
velocity at the wall, and is equal to w, = 0.250w9(0). The beat frequency between
the m = 3 resonance and the m = 2 resonance is 0.846w0(0), and occurs at radius
0.534w9(0). This profile has a central electric potential of ¢) = —0.147522£ n9(0)ri,,
and an unperturbed electric field at the wall of —0.125047=no(0)rw.

Unperturbed Normalized Density and Angular Velocity vs. Radius

angular velocity |

Normalized Profiles
ts

Oo
th

Q.2 0.4 0.6 0.8 1

Normalized Radius

Figure 7.2: Plot of the normalized density and normalized angular velocity for the
Beat Profile. The m = 2 and m = 3 resonant radius, and the radius for the beat
frequency between these two modes are indicated by the dotted lines. The m = 1
resonant frequency is given by the angular velocity at the wall.

128
7.2.2. Demonstration of the Decay Instability

This section presents results from the Nonlinear Fluid Code for a computation that
demonstrates the beat-wave decay instability. A high amplitude m = 3 perturbation
is observed to decay into an m = 2 perturbation through an interaction with the beat
frequency between these two perturbations.

The Nonlinear Fluid code was used to study the response of a pure electron plasma
to a large amplitude externally applied m = 3 pulse. The plasma is first seeded with
a weak one cycle burst of voltage with m = 2 azimuthal symmetry and an amplitude
of 1.0 x 10~¢9. Directly after the one cycle weak m = 2 pulse, a one cycle burst of
voltage with m = 3 azimuthal symmetry and a large amplitude of 0.25¢9 is applied
to the plasma. The electric field at the wall for the m = 3 mode reaches a maximum
amplitude of 5.68% of the unperturbed electric field at the wall.

A plot of the perturbed electric field at the wall versus time is plotted on a semilog
scale in figure 7.3. In this figure, the decay of the m = 3 mode and growth of the
m = 2 mode is clearly visible. Additional nonlinear effects, which may be bounce
oscillations, are also observed. The m = 2 mode is found to have a growth rate of
Yo = 0.038w (0).

The decay of the m = 3 mode and growth of the m = 2 mode is also visible in
images of the total density. Figure 7.4 plots the total density for six successive times.
A clear transition is seen from a mode with m = 3 azimuthal symmetry to a mode
with m = 2 azimuthal symmetry. Images (a)-(f) correspond to the times 2.47t,,,

6.92t.,, 11.4t.,, 15.8t.,, 20.3t.,, and 24.7t,,, respectively.

7.2.3 Growth Rate Versus Amplitude

This section examines the growth rate of the m = 2 daughter mode as a function of
the m = 3 parent mode amplitude. Five computations were performed, varying the
amplitude of the applied m = 3 pulse. Each computation is seeded with an initial one
cycle pulse with m = 2 symmetry and amplitude 1.0 x 107*¢p, directly followed by
the larger m = 3 pulse. The perturbed electric field at the wall for the computation

129

Applied m=3 Amplitude: 0.259,

Log of Perturbed Electric Field
(each tick is a factor of 2)

o iB
Time (Central Rotation Periods)

Figure 7.3: Plot of the envelope of the perturbed electric field at the wall versus time
on a semilog scale for the m = 1,2 and 3 Fourier components due to a large applied
pulse with m = 3 azimuthal symmetry. The beat-wave decay instability causes the
decay of the m = 3 mode and growth of the m = 2 mode.

using an applied m = 3 amplitude of 0.250¢, was shown in figure 7.3. The perturbed
electric field at the wall for the remaining four computations using an applied m = 3
amplitude of 0.177¢ 9, 0.125¢9, 0.0884¢,, and 0.0625¢, is shown in figure 7.5. The
normalized amplitude A3 is proportional to the applied wall potential.

The initial m = 3 density amplitude A3 is computed using equation 7.1, and the
corresponding growth rate of the m = 2 amplitude A, is estimated. To determine
an m = 2 growth rate from the curves in figure 7.5 is difficult at best, so the growth
rates are not very accurate. The results are summarized in table 7.1.

The first column is the amplitude of the one cycle applied m = 3 pulse, and the

130

Figure 7.4: Images of the total density at times 2.47ter, 6.92ter, 11.4ter, 15.8t er, 20.3ber,
and 24.7t,,, after excitation from a large pulse with m = 3 azimuthal symmetry. Due
to beat-wave resonance damping, the m = 3 mode decays into an m = 2 mode.

131

Ss [ye 2 eg
2 — 2 1.
_ [Win - In, 4%
o 5 f\ & § V
i 8 "\ ~\ m8
3 J AW 38 VM An
oS 2 : VV 3 8 ia vey
Ema ae ay
Bos Vi eos A]
a 3 a 3 rey
6 8 ROP 5 8 [KR
Bo m=2 20 m=2
el . ; 4 : .

: Applied m=3 Amplitude: 0,177¢, \ Applied m=3 Amplitude: 0.1256,

0 15 30 45 15 30 45
Time (Central Rotation Periods) Time (Central Rotation Periods)
ss} 3
a lee ra)
co)
‘o Ss we ‘co ‘S m=3
3 8 VAVAV, Ie A Hf \ 3 w NN Dy Lo
a 8 Wr m8
gs “yw! ‘3 * precio
2 4 2 4
3 3 a3 [
bale fy iz iil
Sak sp 68 ff
eT \IN ES ee
| \ Applied m=3 Amplitude:0,0884¢, v | Applied m=3 Amplitude: 0.0625¢,
0 30 rs 90 0 30 60 90
Time (Central Rotation Periods) Time (Central Rotation Periods)

Figure 7.5: Plot of the envelope of the perturbed electric field at the wall for the
m = 3 and m = 2 mode versus time on a semilog scale. The m = 2 mode has
an initial applied amplitude of 0.001¢), and the m = 3 mode has an amplitude of
0.1779, 0.125¢,, 0.0884¢), and 0.0625¢, in the four traces. The time scale has been
doubled in the lower two figures.

second column is the resulting maximum perturbed electric field at the wall for the
m = 3 mode. The third column A3 is the maximum m = 3 mode density amplitude
given by equation 7.1.

The m = 2 growth rate was computed using a least squares fit of a straight line
to the logarithm of the perturbed electric field over the time just after the end of
the m = 3 pulse to the time of the maximum amplitude of the m = 2 mode. As
seen from figure 7.5, the perturbed electric field at the wall for the m = 2 mode has
many features and the logarithm of the perturbed field is not well fit by a straight

line, so the computed growth rates are not very accurate. The fourth column of table

132

m=3 Amplitude | Perturbed field | Ag Yo
0.0625¢, 1.32% 0.0505 | 0.0045w(0) = 0.028
0.0884¢, 1.88% 0.0703 | 0.0078w(0) = 0.049
0.125¢ 2.69% 0.105 | 0.024w9(0) = 0.154
0.177¢, 3.87% 0.141 | 0.026w (0) = 0.17=
0.250¢ 5.68% 0.184 | 0.038w(0) = 0.24

Table 7.1: Growth rate of the m = 2 mode as a function of m = 3 mode initial
amplitude.

7.1 is the calculated m = 2 mode growth rate normalized to the central rotation
frequency and to the central rotation period, so an easy comparison can be made to
the paper by Mitchell and Driscoll [22]. As the amplitude of the applied m = 3 pulse
is increased, the growth rate of the m = 2 mode is found to increase, with the m = 2
mode becoming the dominant mode earlier in time at higher amplitudes.

The m = 2 growth rate y. versus m = 3 mode amplitude Ag is plotted on a log-log
scale in figure 7.6. The straight line in figure 7.6 has slope 2, indicating yz = ['3,2.A3.
A least squares fit to the points has a slope of 1.7 + 0.1, but this error does not
account for the large uncertainty in the estimated growth rates. Using the slope 2 fit,

the coupling constant is calculated to be ['3,. = 10.55 + 2.4.

7.2.4 Comparison to Experimental Results

Comparing the amplitudes of the m = 3 mode with the amplitudes used in the
experiment by Mitchell and Driscoll, it is found that the amplitudes used in the
computation are near the highest amplitude experimental runs. The experimental
coupling constant is I'3 2 = 1.32, while the numerical results give a coupling constant
and growth rates that are 8 times higher. The experimental data has a single point
that does not fall along the trend line of the other experimental points. This highest
amplitude experimental result has an amplitude of A; = 0.11, and a corresponding
m = 2 growth rate of y. = 0.09- which differs by less than a factor of two from the
computation with amplitude A; = 0.105 and a corresponding m = 2 growth rate of

= 0.15-. At these amplitudes, the coupling constant may not be independent

133

ag
= °
wo 01
s) :

0.01

0.01 01 i
Amplitude A;

Figure 7.6: Estimated growth rates 7, for daughter mode m = 2 normalized to the
central rotation time versus normalized amplitude A3 of parent mode m = 3. The
straight line indicates 7. = 3,23.

of amplitude. The computational growth rates and experimental growth rates are
plotted versus parent amplitude on the same scale in figure 7.7.

In the experiment it was found for high applied m = 3 amplitudes, the m = 2
mode grows faster than A3. However, for the amplitudes used in the computation,
the m = 2 mode was found to grow slower than A}.

Some differences exist between the experiment and the computation. First, the
exact density profile of the experiment was not used in the computation, but rather
a smooth Fermi function approximation of it. The growth rate is expected to be
sensitive to aro at the beat resonant radius. Second, the m = 2 mode was seeded
with a moderate amplitude in the computation, whereas in the experiment, the m = 2
mode must grow from noise. These factors may account for the differences in the
results between the computation and the experiment. Additionally, the m = 2 mode
growth is not very smooth in this high amplitude range, so the estimation of the
growth rates is not very accurate. It is not possible to use lower applied amplitudes

in the computation since the slower growth rates would require computing further in

134

10° 7 T Fr Pee? 7 i T T YTrTrT? T T YT Trier
1 -—2 A 7
& 10 7
Ly
cf Ae
2 1073 y
: c
me 392
GS 1074 s®
= “bs
5 A
wa 7
_ il A rT i LL ul L lL L Lita
10 6 { L [os Oo
10-8 1071 109

Amplitude A,,

Figure 7.7: Plot of the experimental and computational growth rates for the m = 2
daughter mode versus parent mode m = 3 amplitude. The experimental data of
Mitchell and Driscoll [22] is represented by triangles, and the computational data is
represented by squares. The upper right section of the graph is the same scale shown
in figure 7.6. Higher growth rates were observed in the computation.

time, which is not practical due to the necessary computational time and accumulation

of errors.

7.3 Conclusion

The experimentally observed beat-wave resonance damping was successfully observed
using the Nonlinear Fluid Code. The growth rates of the daughter wave are larger in
the computation than in the experiment, with the coupling constant 8 times higher
than the experimental value.

The growth rate is not easily computed for the five computations discussed in this

section. Additional nonlinear phenomena cause fluctuations in the envelope of the

135
m = 2 mode as a function of time which degrade the least squares fit for the growth
rate. For the amplitudes used in the computation, the coupling constant may not be
independent of amplitude.

The results presented in this chapter are only preliminary. More work is neces-
sary to understand beat-wave resonance damping, the source of the variations in the
computed growth rate, and the discrepancies with the experiment. Further work may
resolve the difference between the growth rates found from the experiment and the

growth rates found from the computation.

136

Chapter 8 Conclusion and Future Work

8.1 Summary of Results

A cylindrical pure electron plasma can under certain approximations be described
by the equations for a two-dimensional, incompressible, inviscid fluid. By applying
numerical techniques to these equations, a two-dimensional, cylindrical, nonlinear
fluid code was created and used to study the response of a pure electron plasma to
various externally applied pulses.

For small amplitude pulses, the response of the plasma, quantified by the per-
turbed electric field at the wall, was found to decay exponentially in time. The decay
rate is strongly dependent on the density profile used in the calculation. The de-
cay is caused by a collisionless phase mixing process, due to differential rotation of
the plasma at different radii. This mixing process is clearly illustrated in successive
images of the perturbed density. As expected, the results from the Nonlinear Fluid
Code agree with linearized calculations for these small amplitude pulses.

For moderate amplitude pulses, fluid elements are trapped within the potential
well of the resulting wave. The circulation of the trapped fluid elements in the well
causes a modulation of the decay envelope, a phenomenon which is illustrated in
images of the perturbed density. For a slowly applied pulse, a perturbation, which is
nearly stationary in the wave frame, similar to a BGK mode, is excited.

For two small amplitude pulses that are separated in time, the code also shows
that an echo may occur after the second pulse. The response to the first pulse will
decay by the collisionless phase mixing process and the response to the second pulse
will also decay by this process. A nonlinear interaction between the two applied pulses
can then generate an echo response. Based on a simplified theory, several aspects of
echo behavior have been predicted and verified with the fluid code. An echo can

only occur if the azimuthal mode number of the second pulse is greater than that

137

of the first pulse. A nonlinear interaction between the two applied pulses allows the
echo to form, with the second pulse partially unmixing the first pulse. The echo has
azimuthal mode number mz — m, and occurs at time t = m27/(mz — m1), where the
first pulse of m; azimuthal symmetry occurs at time t = 0, and the second pulse of mz
azimuthal symmetry occurs at time t = r. The amplitude of the echo is proportional
to the amplitude of the first pulse, the amplitude of the second pulse, and the time
between the two applied pulses, 7. For larger pulses, the echo amplitude saturates.
Echoes due to the interaction with various harmonics, and echoes of echoes can also
occur.

Very large perturbations are known to be subject to a decay instability based on
previous experimental work. A parent wave with azimuthal symmetry m can decay
to a daughter wave with m — 1 symmetry by an interaction with the beat frequency
between the two waves. The growth rates seen in the computation are higher than the
growth rates reported in the experiment, and the coupling constant for the daughter

wave is 8 times higher than the experimentally derived value.

8.2 Further Work

As with any area of research, for every question answered, many new questions are
generated. This section will highlight some of the questions which still remain for

further research.

8.2.1 Improvements to the Nonlinear Fluid Code

The Nonlinear Fluid Code presented in this thesis is very useful for studying nonlinear
phenomena in non-neutral plasmas. Improvements to the Nonlinear Fluid Code could
be made by adding additional physics to the code, and by changing some of the

computational methods. Some possible improvements are presented in this section.

138
Viscosity

Viscosity was neglected in the Nonlinear Fluid Code. However, viscosity may be an
important effect in real plasmas. It is expected that over long time scales viscosity
will destroy the fine scale structure due to mixing. This would reduce the amplitude
of the echo. A careful study including the effects of viscosity on the echo could provide

a diagnostic for measuring the viscosity in experimental plasmas.

Higher Order Effects

The pressure gradient, diamagnetic drift, polarization drift, and other thermal effects
were neglected in the Nonlinear Fluid Code. These effects should be included to
see whether they change the results significantly. For an experimental plasma with
density profiles similar to those studied in this thesis, the diamagnetic drift velocity
may not be negligible compared to the E x B drift velocity. The diamagnetic drift
does not seem necessary to describe the nonlinear phenomena studied here, but it may
cause effects that could be quantitatively measured. A modification to the Nonlinear

Fluid Code should be able to explore the differences due to thermal effects.

Fourth Order Runge-Kutta Time Integration

The time integration method could be improved by using a fourth order Runge-Kutta
time integration method. For the same size time step, the fourth order method should
be much more accurate than the second order leapfrog method used in the Nonlinear
Fluid Code. However, the necessary computational time is doubled for the same time

step.

Spectral Code

The accuracy of the Nonlinear Fluid code could be greatly improved by using spectral
methods to solve the 2-D fluid equations. The azimuthal direction can be represented
as a Fourier series, and the radial direction can be represented by a summation of

Chebyshev polynomials. Such a spectral code would not contain the same aliasing

139
errors that exist in the finite difference methods of the current Nonlinear Fluid Code.
However, a much smaller time step would be required to maintain computational sta-
bility, resulting in much longer computational time. With the use of faster computers,

a spectral code may be practical.

8.2.2 Further Research on Non-neutral Plasmas
Fluid Trapping

It was shown that fluid elements could be trapped in an excited wave of moderate
amplitude. This trapping has not been explored as a function of density profile.
A different density profile may better match the experimental results of Pillai [14],
including the asymptotic decay at much later times.

For a small amplitude pulse applied for 90 cycles, a perturbation that is nearly
stationary in the wave frame and resembles a BGK mode is excited. Future com-
putations may be capable of applying a pulse over a much longer duration, perhaps
exciting a true BGK mode that is stationary in the wave frame. It may also be
possible to find an analytical expression describing a BGK mode in a non-neutral

plasma.

Echoes

For short times between pulses, the amplitude of the echo is proportional to the in-
terpulse spacing. For longer times, the echo amplitude saturates. The amplitude
dependence for large interpulse spacing remains unexplored. Such computations re-
quire increased resolution to accurately compute further in time.

Results from the Nonlinear Fluid Code suggest that echoes may be observable
experimentally. For sufficiently large pulses, the echo amplitude can approach the
amplitude of the response to the applied pulses. Based on these results, a future

experimental search for echoes may be successful.

140

Beat-Wave Resonance Damping

Many questions remain regarding the beat-wave resonance damping results. What
causes the modulation in the growth of the m = 2 mode? Why are growth rates found
from the computation higher than growth rates found in the experiment? Transitions
other than the m = 3 to m = 2 have not yet been explored with the Nonlinear
Fluid Code. The dependence on the density profile is unknown, and the growth
rate is probably quite sensitive to ao at the beat-wave resonant radius. It should
be possible to construct a density profile where the m = 3 mode and the m = 2
mode are undamped. For such a profile, the m = 2 mode will grow unhindered by
the linear decay mechanism, and beat-wave resonance damping could be studied for

smaller applied pulses.

Other Phenomena

Many other nonlinear phenomena that occur in non-neutral plasmas can be studied
using the Nonlinear Fluid Code. The hollow beam instability, vortex mergers and
vortex crystals are a few examples of other phenomena that could be examined. It is
expected that this code will be an important new tool to be used in future research

on these plasmas.

141

Bibliography

[1]
[2|

[3]

[9]

R. G. Greaves, M. D. Tinkle, and C. M. Surko, Phys. Plasmas 1, 1439 (1994).
G. Gabrielse, Sci. Am. 267, 78 (1992).

J. J. Bollinger, D. J. Prestage, W. M. Itano, and D. J. Wineland, Phys. Rev.
Lett. 54, 1000-1003 (1985).

T. B. Mitchell, M. M. Schauer, and D. C. Barnes, Bull. Am. Phys. Soc. 41, 1603
(1996).

L. Brillouin, Phys. Rev. 67, 260 (1945).
H. F. Webster, J. Appl. Phys. 26, 1386 (1955).
G. G. McFarlane and H. G. Hay, Proc. Phys. Soc. 63, 409 (1953).

A. H. W. Beck, Space Charge Waves and Slow Electromagnetic Waves (Pergamon
Press, 1958).

R. H. Levy, Phys. Fluids 8, 1288 (1965).

[10] R. J. Briggs, J. D Daugherty, and R. H. Levy, Phys. Fluids 13, 421 (1970).

[11] L. Landau, J. Phys. (Moscow) 10, 25 (1946).

[12] J. H. Malmberg and J. S. deGrassie, Phys. Rev. Lett. 35, 577 (1975).

[13] C. F. Driscoll and K. S. Fine, Phys. Fluids B 2, 1359 (1990).

[14] N.S. Pillai and R. W. Gould, Phys. Rev. Lett. 73, 2849 (1994).

[15] N. R. Corngold, Phys. Plasmas 2, 620 (1995).

[16] R. L. Spencer and S. N. Rasband, Phys. Plasmas 4, 53 (1997).

142
[17] C. W. Roberson and C. F. Driscoll, Non-Neutral Plasma Physics, AIP Conf.
Proc. No. 175 (American Institute of Physics, New York, 1988).

[18] J. Fajans and D. H. E. Dubin, Non-Neutral Plasma Physics II, AIP Conf. Proc.
No. 331 (American Institute of Physics, Woodbury, 1995).

[19] N.S. Pillai, Ph.D. thesis, California Institute of Technology (1995).

[20] I. B. Bernstein, J. M. Greene, and M. D. Kruskal, Phys. Rev. 108, 546 (1957).
[21] R. W. Gould, Phys. Plasmas 2, 2151 (1995).

[22] T. B. Mitchell and C. F. Driscoll, Phys. Rev. Lett. 73, 2196 (1994).

[23] Charles K. Birdsall and A. Bruce Langdon, Plasma Physics via Computer Sim-
ulation (McGraw-Hill, 1985).

[24] William H. Press, Brian P. Flannery, Saul A. Teukolsky, and William T. Vetter-
ling, Numerical Recipes, The Art of Scientific Computing (Cambridge University
Press, 1986).

[25] N. A. Phillips, Quart. J. Roy. Meteorol. Soc. 82, 123-164 (1956).

[26] N. A. Phillips, The Atmosphere and the Sea in Motion (Rockefeller Institute
Press and the Oxford Univ. Press, New York, 1959) pp. 501-504.

[27] A. Arakawa, J. Comput. Phys. 1, 119 (1966).

[28] K. M. Case, Phys. Fluids 3, 143 (1960).

[29] W. Thompson, Nature 23, 45 (1880).

[30] E. L. Hahn, Phys. Rev. 80, 580 (1950).

[31] R. M. Hill and D. E. Kaplan, Phys. Rev. Lett. 14, 1062 (1965).
[32] R. W. Gould, Phys. Lett. 19, 477 (1965).

[33] I. D. Abella, N. A. Kurnit, and S. R. Hartman, Phys. Rev. 141, 391 (1966).

143
[34] R. W. Gould, T. M. O’Neil, and J. H. Malmberg, Phys. Rev. Lett. 19, 219 (1967).

[35] J. H. Malmberg, C. B. Wharton, R. W. Gould, and T. M. O’Neil, Phys. Fluids
11, 1147 (1968).

[36] J. D. Crawford and T. M. O’Neil, Phys. Fluids 30, 2076 (1987).

[37| T. M. O’Neil to R. W. Gould, private communication (1997).