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 near the entrance and exit of the channel. The reason is the appearance of a strong through a highly conductive and rapidly moving plasma: current tends to concentrat back EMF which tends to block current over most of the channels length this is most easily seen if we" unwrap"the annular chamber of an MPd thruster into a ectangular 1-D channel B=0 了 H Amperes law j=-V×B (1) In our case and calling dB Ohms law(ignoring Hall effect)is j=G(E+uxB Pro f a spa m artie ssn Lecture 22 1 of 8
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 ∇ µ × G JG (1) In our case = l,x x ∂ ∇ ∂ G so y z 0 1 dB j j=+ dx ≡ µ and calling B -By ≡ , 0 1 dB j=- µ dx (2) Ohm’s law (ignoring Hall effect) is j = E + u × B σ ( ) G GGJG (3)
or, using j=σ(E-uB) 4 Combining(2)and (4) d (E The flow velocity u evolves along x according to the momentum equation(ignoring pressure forces neglect for now Substitute(2)into(6) u wh d B2 μodx Integrate mu+wH B WH B2-B2 Putting this in Equation(5), B(B:-B2 If we approximate the conductivity o as a constant this can be integrated as B(B2-B2) 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 2 of 8
16.522, Space Propulsion Lecture 22 Prof. Manuel Martinez-Sanchez Page 2 of 8 or, using z yx 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 = j B A = jBwH dx dx × i G JG (6) neglect for now Substitute (2) into (6): 2 0 0 du 1 dB wH d B m = - B wH = - dx dx dx 2 ⎛ ⎞ ⎜ ⎟ µ µ ⎝ ⎠ i (7) Integrate: 2 2 0 0 0 0 B B mu + wH = mu + wH 2 2 µ µ i i neglect 2 2 0 0 wH B -B u = 2 m µ i (8) Putting this in Equation (5), ( ) 2 2 0 0 0 dB wH =- E- B B -B dx 2 m ⎡ ⎤ σµ ⎢ ⎥ ⎣ ⎦ µ i (9) If we approximate the conductivity σ as a constant, this can be integrated as B0 0 2 2 B 0 0 dB x = wH E - B(B - B ) 2 m σµ µ ∫ 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
denominator in the integrand is the driving field(applied field E, minus back emf, uB). The field Bo at x=0 is a measure of the current I, because integrating(2) between x=0 and x=l gives I Bo B dI (11) On the other hand, carrying(10) all the way to x=L, gives d B (12) B where, once i and m 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 ooL. 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 TTTTTTTTTTTT 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 B(B:-B2), namely, B2-3B2=0 16.522, Space P pessan Lecture 22 Prof. Manuel martinez Page 3 of 8
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 B0 at x=0 is a measure of the current I, because integrating (2) between x=0 and x=1 gives L 0 0 0 0 0 I B I jdx = = B = w w µ ⇒ µ ∫ (11) On the other hand, carrying (10) all the way to x=L, gives B0 0 2 2 0 0 0 dB L = wH E - B(B - B ) 2 m σµ µ ∫ 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 σµ0L . 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 2 2 B(B - B ) 0 , namely, 2 2 B - 3B = 0 0
nd th (E-UBMIN=E WH 2B 33 and since by assumption, this difference is much less than E, we find H B 3√ or, in terms of voltage and current, 1H),I3 33 Returning now to Equation (5), we notice that near both x=0 and x=L, uB<<E,so d B dx - OHoE, and so the thickness i of the thin current layers(where B varies substantially) can be estimated as follows HDE where B Using(16) 3(3-1) aWHB2 (17) NHB Remembering(from( 8)) that the exit velocity is Bo 1 HH. these results can be written as cHou.lo The non-dimensional group foul is called the Magnetic Reynolds Number(rm) (based on length I). What we have seen is that this rm, when based on the current 16.522, Space Propulsion Lecture 22 Prof. Manuel martinez-s of 8
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 wH 2B E - uB = E - 3 3 2 m µ i , and, since by assumption, this difference is much less than E, we find 3 0 0 wH 2B E 3 3 2 m µ i (14) or, in terms of voltage and current, 3 2 0 1H I V 3 3 w m ⎛ ⎞ µ ⎜ ⎟ ⎝ ⎠ i (15) Returning now to Equation (5), we notice that near both x=0 and x=L, uB<<E, so 0 dB - E dx σµ , and so the thickness l of the thin current layers (where B varies substantially) can be estimated as follows: 0 1 1 0 e 0 0 B -B B l ; l σµ σµ E E (16) where 0 1 B B = 3 . Using (16), ( ) 0 e 2 2 0 0 3 3 -1 m 3m l ; l σ σ wHB wHB i i (17) Remembering (from (8)) that the exit velocity is 2 2 0 0 e 0 wH 1 H B I u= = 2 2w m m µ µ i i (18) these results can be written as 0 e0 ( ) 3 u l 3 -1 2 σµ ; 0 e0 3 u l 2 σµ (19) The non-dimensional group 0 σµ ul is called the Magnetic Reynolds Number (Rm) (based on length l). What we have seen is that this Rm, when based on the current
layer thickness, is of order unity. Since we started out by assuming conditions when these layers are thin i.e. le, lo>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 (21) owo (the 4/3 factor accounts for the"triangular"current distribution in the layers)and so the ohmically dissipated power is D=IRo+IR。 (22 ere Io B1) d I=I-5=/3 Substituting, we find 1 Dn=2|1- 1 H aHui √3(H)aI W3|√3-1mw =33wm D.=n13帅=41邮 and, in total D √3 (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 ev=2 to 3 16.522, Space P pessan Lecture 22 Prof. Manuel martinez Page 5 of 8
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., le, l0> 1 m 0e ( ) ≡ σµ (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 2 2 D=IR +IR 00 ee (22) where ( ) 0 0 01 0 0 w 11 wB I = B -B = 1- =I 1- 3 3 ⎛ ⎞⎛ ⎞ ⎜ ⎟⎜ ⎟ µ µ ⎝ ⎠⎝ ⎠ (23) and e 0 I =I-I =I 3 (24) Substituting, we find ( ) 2 22 24 2 2 0 0 0 1 1 - 14 4 H HI I 3 D =I 1- H = 3 33 3 3w w3 3 - 1 m w m ⎛ ⎞ σ µ ⎛ ⎞ ⎜ ⎟ ⎜ ⎟ ⎝ ⎠ ⎝ ⎠ σ µ i i 22 24 2 2 0 0 e 4 H 1 41 H 3 HI I D =I = 3 w 39 w 3wm m σ µ ⎛ ⎞ ⎜ ⎟ σ ⎝ ⎠ µ i i and, in total, 2 2 4 0 4 H I D = 9 3 w 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 3
times(evi be the effective ionization energy per atom, and a, the degree of ionization at the exit we then have (26) from(25) m. 4H HoI 9(W (27) This indicates very rapid increase of the ionization fraction as the current increases or as the flow is reduced Instabilityonset For a given thruster, as I/m is increased a 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 a, 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 I2/m. Since(as we will see)efficiency increases with I /m, 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 lantitatively Setting(27)to unity, we get H HoP2-19v3 ev and from the exit velocity expression(18), 1e=098 4 The velocity at which the particle's kinetic energy would be capable of ionizing it is called the" Alfven critical speed". Many years ago Alfven 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 6 of 8
16.522, Space Propulsion Lecture 22 Prof. Manuel Martinez-Sanchez Page 6 of 8 times (eVi) be the effective ionization energy per atom, and αe the degree of ionization at the exit. We then have e i i ' m D eV m α i (26) or, from (25), 2 2 i 0 e i ' m 4 H I eV 9 3 w m ⎛ ⎞ µ α ⎜ ⎟ ⎜ ⎟ ⎝ ⎠ i (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 I mi 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 I mi . Since (as we will see) efficiency increases with 2 I mi , 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 ' H 93 I eV = w 4m m µ i (28) and from the exit velocity expression (18), i i e i i ' ' 1 93 eV eV u = = 0.987 2 4m 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
limited roughly to the Alfven critical speed of the gas used (if exit arcs are to be For various gases, assuming V= 2V, we find Hydrogen Nitrogen argon Lithium (3)(s)512014208701220 M Expressing m=n(Na=Avogadro's number), and using kA for I and g/s for m Equation(28) can also be rewritten as PM(kA?(g/mol =154 (30) volts For Argon at 6 g/s and W/H=4, this predicts an "onset current"I= 18kA Experimental values tend to cluster around r-20-23 ka, in reasonable agreement. Much of the difference is simply due to the geometry(coaxial vs rectangular ).The scaling of I with m and with M /2 is also well documented experimentally Efficiency Accounting only for the power lost to ohmic dissipation and to near-electrode voltage drops(△V=△V cathod+△V 10-20 Volts ) we have 1 u 2 muz +d+IAv From the equations derived before 212w 1Hs2+4(H 22 9 2I△V 1+ 1HHoI' 3.05 16.522, Space P pessan Lecture 22 Prof. Manuel martinez Page 7 of 8
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 i i ' V = 2V , we find Gas Hydrogen Nitrogen Argon Lithium Mi (g/mol) 1 14 40 7 Vi (volts) 13.6 14.6 15.8 5.4 (Isp)MAX (s) 5,120 1,420 870 1,220 Expressing i a M m = N (Na=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 ' IM w kA g/mol 15.4 V g/s H m ⎛ ⎞ ⎜ ⎟ ⎜ ⎟ ⎝ ⎠ i (30) volts For Argon at 6 g/s and w/H=4, this predicts an “onset current” * I 18kA . Experimental values tend to cluster around * I 20 - 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 ( ∆∆ ∆ V = V + V 10 - 20 Volts cathode anode ∼ ), 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+ + 9 3 1 H I m 2 w m η ∆ ⎛ ⎞ µ ⎜ ⎟ ⎜ ⎟ ⎝ ⎠ i i 3.05 2 2 0 2 2 24 2 0 0 1 1H I m 2 2w = m 1 1H 4 H I I m + +I V 2 2w w m m 9 3 ⎛ ⎞ µ ⎜ ⎟ ⎜ ⎟ ⎝ ⎠ η ⎛ ⎞ µ µ ⎛ ⎞ ⎜ ⎟ ⎜ ⎟ ∆ ⎜ ⎟ ⎝ ⎠ ⎝ ⎠ i i i i i
n 8(△V)m 1(H Even if the voltage drops could be eliminated, this efficiency is no higher than 3.05 0.328. The best that can be done in the presence of the voltage drop is to approach onset, at which point one gets n 1 (33 3.05+2.:H For Argon, m=6 g/s 4, Av=10 Volts this gives n =0.259 H Values of this order have been reported very often, and it has proven very difficult to exceed n=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 16.522, Space P pessan Lecture 22 Prof. Manuel martinez Page 8 of 8
16.522, Space Propulsion Lecture 22 Prof. Manuel Martinez-Sanchez Page 8 of 8 ( ) 2 2 3 0 1 = 8 Vm w 3.05 + I H η ∆ ⎛ ⎞ ⎜ ⎟ µ ⎝ ⎠ 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 2 3 4 i i 0 ' 1 = w H m 3.05 + 2.88 V m eV η ⎛ ⎞ ⎛ ⎞ ⎜ ⎟ ∆ ⎜ ⎟ ⎜ ⎟ ⎜ ⎟ µ ⎝ ⎠ ⎝ ⎠ 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