16.522, Space Propulsion Lecture 22
Prof. Manuel Martinez-Sanchez Page 1 of 8
16.522, Space Propulsion
Prof. Manuel Martinez-Sanchez
Lecture 22: A Simple Model For MPD Performance-onset
It is well known that rapidly pulsed current tends to concentrate near the surface of
copper conductors forming a “skin”. A similar effect occurs when current flows
through a highly conductive and rapidly moving plasma: current tends to concentrate
near the entrance and exit of the channel. The reason is the appearance of a strong
back EMF which tends to block current over most of the channel’s length. This is
most easily seen if we “unwrap” the annular chamber of an MPD thruster into a
rectangular 1-D channel.
Ampère’s law:
0
1
j= B?
μ
×
nullnullnull
(1)
In our case x=l,
x
?
?
?
null
so
y
z
0
dB
1
jj=+
dx
≡
μ
and calling
y
B-B≡ ,
0
1dB
j=-
dxμ
(2)
Ohm’s law (ignoring Hall effect) is
()
j= E+u × Bσ
nullnullnullnullnull
(3)
16.522, Space Propulsion Lecture 22
Prof. Manuel Martinez-Sanchez Page 2 of 8
or, using
zyx
V
E E = , B -B , u = u
H
≡≡
()j= E-uBσ
(4)
Combining (2) and (4),
()
0
dB
=- E-uB
dx
σμ (5)
The flow velocity u evolves along x according to the momentum equation (ignoring
pressure forces)
()
x
du dP
m+A=jBA=jBwH
dx dx
×
i nullnullnull
(6)
neglect for now
Substitute (2) into (6):
2
00
du 1 dB wH d B
m=- B wH=-
dx dx dx 2
??
??
μμ
??
i
(7)
Integrate:
22
0
0
00
BB
mu+wH =mu +wH
22μμ
ii
neglect
22
0
0
B-BwH
u=
2
m
μ
i
(8)
Putting this in Equation (5),
()
22
00
0
dB wH
=- E- B B -B
dx
2m
??
??
σμ
μ
??
i
(9)
If we approximate the conductivity σ as a constant, this can be integrated as
0
B
0
22
B
0
0
dB
x=
wH
E- B(B -B )
2m
σμ
μ
∫
i
(10)
This integral can actually be calculated analytically, but the resulting expression is
not very transparent. It is more useful to examine its behavior qualitatively. The
16.522, Space Propulsion Lecture 22
Prof. Manuel Martinez-Sanchez Page 3 of 8
denominator in the integrand is the driving field (applied field E, minus back emf,
uB). The field B
0
at x=0 is a measure of the current I, because integrating (2)
between x=0 and x=1 gives
L
00
0
00
BII
jdx = = B =
ww
μ
?
μ
∫
(11)
On the other hand, carrying (10) all the way to x=L, gives
0
B
0
22
0
0
0
dB
L=
wH
E- B(B -B )
2m
σμ
μ
∫
i
(12)
where, once I and m
i
are specified, only E remains as an unknown. This is then the
equation for voltage V=EH. Consider conditions where the maximum value of the
back emf uB reaches almost the level E. This means that the integrand will be very
large as long as this condition prevails, and it must indicate a large value of
0
Lσμ . By
the same token, Equation (5) says that B will remain flat when E-uB<<E, and from
(2), there will be little current in this region. Schematically:
We see here that two strong current concentrations develop, near x=0 and x=L. Let
us investigate when this situation will arise. From (12), the denominator is minimum
at a B value that maximizes
22
0
B(B -B ), namely,
22
0
B-3B=0
16.522, Space Propulsion Lecture 22
Prof. Manuel Martinez-Sanchez Page 4 of 8
MAX
0
1(uB)
B
B Β =
3
≡ (13)
and then
()
3
0
MIN
0
2BwH
E-uB =E-
33
2mμ
i
,
and, since by assumption, this difference is much less than E, we find
3
0
0
2BwH
E
33
2mμ
i
null (14)
or, in terms of voltage and current,
3
2
0
1H I
V
w
33
m
??
μ
??
??
i
null (15)
Returning now to Equation (5), we notice that near both x=0 and x=L, uB<<E, so
0
dB
-E
dx
σμnull , and so the thickness l of the thin current layers (where B varies
substantially) can be estimated as follows:
01 1
0e
00
B-B B
l ; l
EEσμ σμ
nullnull (16)
where
0
1
B
B=
3
.
Using (16),
()
0e22
00
33-1m
3m
l ; l
wHB wHBσσ
i
i
nullnull (17)
Remembering (from (8)) that the exit velocity is
22
00
e
0
BIwH 1 H
u= =
22w
mm
μ
μ
ii
(18)
these results can be written as
()0e0
3
ul 3-1
2
σμ null ;
0e0
3
ul
2
σμ null (19)
The non-dimensional group
0
ulσμ is called the Magnetic Reynolds Number (R
m
)
(based on length l). What we have seen is that this R
m
, when based on the current
16.522, Space Propulsion Lecture 22
Prof. Manuel Martinez-Sanchez Page 5 of 8
layer thickness, is of order unity. Since we started out by assuming conditions when
these layers are thin, i.e., l
e
, l
0
<<L, we can now state that this will occur when
()
m0e
RL uL>>1≡σμ (20)
This is indeed the condition for operation in the pure MPD regime.
Effects of Dissipation
The high-current inlet and exit layers are very dissipative. Their resistances can be
estimated as
0
0
4
H
3
R=
wlσ
;
e
e
4
H
3
R=
wlσ
(21)
(the 4/3 factor accounts for the “triangular” current distribution in the layers) and so
the Ohmically dissipated power is
22
00 ee
D=IR +IR (22)
where
()
0
001
00
wBw11
I = B -B = 1- =I 1-
33
????
????
μμ
????
(23)
and
e0
I=I-I=I 3 (24)
Substituting, we find
()
2 2
22 24
2 00
0
1
1-
HI I14 4 H
3
D=I 1- H =
33w
w3 3 -1 mw m
σμ??
??
?? ??
????
σ
μ
ii
2
22 24
2 00
e
4
H
HI I141H
3
D=I =
3w 39w
3wm m
σμ
??
??
σ
??
μ
ii
and, in total,
2
24
0
I4H
D=
w
93
m
μ??
??
??
i
(25)
Part of this dissipation goes to heating the gas, but the major portion is used in
ionizing and exciting (followed by radiation) the gaseous atoms. Let
'
i
eV 2 to 3null
16.522, Space Propulsion Lecture 22
Prof. Manuel Martinez-Sanchez Page 6 of 8
times (eV
i
) be the effective ionization energy per atom, and
e
α the degree of
ionization at the exit. We then have
ei
i
m
'
DeV
m
α
i
null (26)
or, from (25),
2
2
0i
e
i
'
Im 4H
w
eV 93
m
??
μ
??α
??
i
null (27)
This indicates very rapid increase of the ionization fraction as the current increases,
or as the flow is reduced.
Instability “onset”
For a given thruster, as
2
Im
i
is increased,
e
α increases rapidly. When it reaches
unity, the behavior of the plasma near the exit changes drastically. This is because
any extra dissipation cannot be absorbed into ionization anymore, and goes instead
directly into heating the plasma (or perhaps the electron component only). This
causes conductivity to increase whenever the current concentrates, which leads to
further current concentration. We have here the classical prescription for constriction
into an arc, and one can expect heavy arcing (with the corresponding damage to
electrodes) when
e
α approaches 1. This behavior has indeed been observed
repeatedly, and has been the focus of a lot of attention, because it limits the
practically achievable value of
2
Im
i
. Since (as we will see) efficiency increases with
2
Im
i
, this is a major hurdle in the path towards efficient MPD operation. It has been
dubbed “the onset condition and we are now in a position to see what it implies
quantitatively. Setting (27) to unity, we get
2
0 i
i
'
I eVH93
=
w4m
m
μ
i
(28)
and from the exit velocity expression (18),
ii
e
''
eV eV193
u = = 0.987
24m m
(29)
The velocity
i
i
'
2eV
m
at which the particle’s kinetic energy would be capable of ionizing
it is called the “Alfvèn critical speed”. Many years ago Alfvèn used this conversion of
kinetic to ionization energy to construct a model of the “condensation” of matter
expanding from the proto-Sun to form the existing planets. We see here that the
exhaust speed of an MPD thruster (or a PPT, which works on the same principles) is
16.522, Space Propulsion Lecture 22
Prof. Manuel Martinez-Sanchez Page 7 of 8
limited roughly to the Alfvèn critical speed of the gas used (if exit arcs are to be
avoided).
For various gases, assuming
ii
'
V=2V, we find
Gas Hydrogen Nitrogen Argon Lithium
M
i
(g/mol) 1 14 40 7
V
i
(volts) 13.6 14.6 15.8 5.4
(I
sp
)
MAX
(s) 5,120 1,420 870 1,220
Expressing
i
a
M
m=
N
(N
a
=Avogadro’s number), and using kA for I and g/s for m
i
,
Equation (28) can also be rewritten as
()
1
2
2 2
i
'
kA g/mol
IM w
15.4 V
g/s H
m
??
??
??
i
null (30)
volts
For Argon at 6 g/s and w/H=4, this predicts an “onset current”
*
I18kAnull .
Experimental values tend to cluster around
*
I20-23 kA~ , in reasonable agreement.
Much of the difference is simply due to the geometry (coaxial vs. rectangular). The
scaling of
*
I with
1
2
m
i
and with
1
-
2
Μ is also well documented experimentally.
Efficiency
Accounting only for the power lost to ohmic dissipation and to near-electrode voltage
drops (
cathode anode
V= V + V 10-20 Volts?? ? ~ ), we have
2
e
2
e
1
mu
2
=
1
mu +D+I V
2
η
?
i
i
(31)
From the equations derived before,
2
2
0
1
=
32 2I V
1+ +
93
I1H
m
2w
m
η
?
??
μ
??
??
i
i
3.05
2
2
0
2
2
224
00
I11H
m
22w
m
=
II11H 4H
m+ +IV
22w w
93
mm
??
μ
??
??
η
??
μμ
??
?? ?
??
??
??
i
i
i
ii
16.522, Space Propulsion Lecture 22
Prof. Manuel Martinez-Sanchez Page 8 of 8
()
2
23
0
1
=
8Vm
w
3.05+
HI
η
?
??
??
μ
??
i
(32)
Even if the voltage drops could be eliminated, this efficiency is no higher than
MAX
1
= = 0.328
3.05
η . The best that can be done in the presence of the voltage drop
is to approach “onset”, at which point one gets
*
1
32
4
i
i
0
'
1
=
mwH
3.05+2.88 V
eV
m
η
??
??
??
?
??
??
μ
??
i
(33)
For Argon, m=6 g/s
i
,
w
=4, V=10 Volts
H
? this gives
*
=0.259η .
Values of this order have been reported very often, and it has proven very difficult to
exceed =0.3η with Argon at least. Most of the inefficiency is seen to arise from the
strong dissipation in the inlet and exit layers. This is intrinsic to the constant area
geometry. Although it may not be obvious at this point, a convergent-divergent
geometry has the effect of weakening these dissipative layers, and can conceivably
be exploited to improve efficiency (and retard “onset”). There is some evidence for
this in experiments. The analysis showing the effect of a convergent-divergent
geometry (or, in a more limited form, of a purely divergent geometry), can be found
in Ref. 1.
References:
Ref. 1: Martinez-Sanchez, M. “Structure of Self-Field Accelerated Plasma Flows.” J. of Propulsion and
Plasma 7, no. 1 (Jan-Feb 1991): 56-64.