16.522, Space Propulsion Prof. Manuel Martinez-Sanchez Lecture 11-12: SIMPLIFIED ANALYSIS OF ARCJET OPERATION 1. Introduction These notes aim at providing order-of-magnitude results and at illuminating the mechanisms involved. Numerical precision will be sacrificed in the interest of physical clarity. We look first at the arc in a cooled constrictor, with no flow, expand the analysis to the case with flow, and then use the results to extract performance parameters for arcjets. 2. Basic Physical Assumptions The gas conductivity model will be of the form σ = o aT? T e () ? ? ? T < T e ( ) T > T e () T e ≈ 6000? 7000K ( ) a ≈ 0.8Si / m / K (1) The termal conductivity k of the gas will be modelled as a constant (with possibly a different value outside the arc). This is a fairly drastic simplification, since in H 2 and N 2 k(T) exhibits very large peaks in the dissociation range (2000-5000K) and in the ionization range (12000-16000K). Because k always multiplies a temperature gradient, the combination d is relevant, and so the proper choice of Φ T()= kdT k to be used is the averaged value k = 1 T 2 ? T 1 kdT T 1 T 2 ∫ (2) over the range of temperatures intended. The arc gas is modelled as ideal, even though its molecular mass shifts strongly and its enthalpy increases rapidly in the dissociation and ionization ranges. In particular, c p = ?h ?T ? ? ? ? p has strong peaks, similar to those of k(T), and, once again, we should use temperature-averaged values for it. The arc is assumed quasi-cylindrical, with axial symmetry and with gradients which are much stronger in the radial than in the axial direction (similar to boundary layers). The flow region comprises three sub-domains: (a) The arc itself, for r , corresponding to < R a x() T >T e . This is the only part carrying current. (b) The outer gas, not ionized and with T <T e . 16.522, Space Propulsion Lecture 11-12 Prof. Manuel Martinez-Sanchez Page 1 of 18 (c) For the case with coaxial flow, a thin transition layer between (a) and (b) may be necessary for accuracy, but will be ignored in our analysis. 3. Constricted Arc With No Flow The typical arrangement is a strongly water-cooled cylindrical enclosure, made of mutually insulated copper segments, with the arc burning along its centerline (Fig. 1). Fig. 1. Constricted Arc Except for the near-electrode regions, the arc properties are constant along its length. In a cross-section, the axial electric field E = E x is independent of radius as well, and E r is small. The Ohmic dissipation rate is r j . r E per unit volume, or σE 2 , since . Here r j =σ r E σ varies strongly inside the arc, from zero at r=R a to a maximum σ c at the centerline; as a rough approximation, we take 1 2 σ c as a representative average, and so the amount of heat deposited ohmically per unit length is 1 2 πR a 2 σ c E 2 . This heat must be conducted to 16.522, Space Propulsion Lecture 11-12 Prof. Manuel Martinez-Sanchez Page 2 of 18 the arc’s periphery, and so it must equal 2πR a ()k ?T ?r ? ? ? ? r=R a . Representing the temperature gradient R a by (roughly) ? ?T ?r ? ? ? ? R a ? 2 T c ? T e R a , we obtain / π R a / 2 1 2 σ c E 2 = 2/ π R a k c 2 T c ? T e R a or E = 2 2k c T c ?T e ( ) σ c 1 R a (3) and since σ c = aT c ? T e (), E = 22 k c a 1 R a (4) This important result indicates that the arc field, and hence its voltage, is inversely proportionally to its radius: the dissipation must increase if the arc is constrained more tightly, which improves its cooling. But note that R a itself is not yet known, since it is only R, the constrictor diameter that is prescribed. The total arc current is . Once again, using I = 2πr σE() o R ∫ dr σ ? 1 2 σ c , we obtain I = πR a 2 σ c 2 E (5) and substituting (4) here, I = πR a 2 aT c ? T e () 2 22 k a a 1 R a I = π 2ak c T c ? T e ( )R a (6) Note also that, multiplying (4) and (6) together we obtain EI = 4πk c T c ? T e ( ) (7) which is another way to express the heat balance. Its main message is that the arc centerline temperature T c increases linearly with arc power per unit length, EI. 16.522, Space Propulsion Lecture 11-12 Prof. Manuel Martinez-Sanchez Page 3 of 18 Everything covered so far will also apply later to the arc in a flow. The difference is in how the heat is evacuated from the arc periphery. With no flow, this must be accomplished by heat conduction through the buffer gas. If k out denotes its thermal conductivity (probably much lower than k c ), and if we ignore cylindrical effects, we must have equality of heat flux (per unit area) on both sides of the arc’s edge: k c 2 T c ? T e R a = k out T e ? T w R ? R a (8) where T w is the temperature of the constrictor’s wall, controlled externally. Substituting T c -T e from (8) into (6), I = πR a 2ak c () 1 2 k out k c R a R ? R a T e ? T w () This is a quadratic equation for the arc radius Ra. To simplify algebra, introduce non- dimensional quantities: I * = I I ref ; I ref =πR 2ak c T e ?T w () (9) λ = k out 2k c (10) and r a = R a R (1) and so I * = λ r a 2 1? r a (1) Solving for r a , r a = 2 1+ 1+ 4λ I * (12) which approaches 1 from below as I * becomes large. We can now obtain other quantities of interest. From (8), T c ?T e T e ?T w = 1 2 k out k c R a R? R a = λ r a 1?r a = 2λ 1+ 4λ I * ?1 16.522, Space Propulsion Lecture 11-12 Prof. Manuel Martinez-Sanchez Page 4 of 18 or, rearranging, T c ?T e T e ?T w = 1+ 1+ 4λ I * 2 I * (13) which shows how Tc eventually increases linearly with I * , but its variation is faster ≈ I * ( ) at low current. The field itself follows from (4). Define a non-dimensional field E * = E E ref ; E ref = 2 2k c a 1 R (14) and then E * = 1 r a = 1+ 1+ 4λ I * 2 (15) which indicates a decreasing field (and voltage) as the current increases. These results are summarized in Fig. 2, calculated with λ =1/ 4. The negative slope of the line E * = fI * ( ) is typical of arc discharges, and creates some difficulties in their operation. We note first that the increase of E * as I * decreases does not continue indefinitely; below some current level, the thermal power input to the electrodes (particularly the cathode) is insufficient to sustain the electron emission required, and the discharge transitions to a different mode, probably an “anomalous glow discharge”. 16.522, Space Propulsion Lecture 11-12 Prof. Manuel Martinez-Sanchez Page 5 of 18 Figure 2 16.522, Space Propulsion Lecture 11-12 Prof. Manuel Martinez-Sanchez Page 6 of 18 This happens at current levels in the micro- to milliampere, and the complete V-I curve of an arc then appears as in Fig. 3: Figure 3. Complete V-I Arc Characteristic Assuming we wish to operate at point P on this line, if we simply connect the arc to a constant-voltage source, such as a battery, we obtain an unstable arrangement. This is because, if current, say, increases slightly above Ip, the arc would now demand less than the equilibirum voltage V p which the source delivers, and so I would increae even further, and run away (to the supply limit, or to destruction). Conversely, any negative current fluctuation would cause a rapid snap back to the stable operating point Q (at very low current). One solution is to insert a series resistance R B (“ballast”) and increase the source open-circuit voltage to, say V B (Fig. 3). Since the voltage available to the arc is now V B ? R B I , the supply line cuts the arc line form above (Fig. 3), and repeating the argument we notice stable operation at P. But, of course, we dissipate in the ballast the power I p V B ?V p ( ) , which leads to inefficient operation. Notice, however, that AC arcs can be efficiently ballasted by a series inductor, the familiar “choke” in fluorescent fixtures. Alternatively, one can use a current-regulated power supply, provided its regulating speed and authority are high enough. For space applications, this takes the form of high frequency solid state switching regulators, capable of stabilizing the arc with minimal losses (efficiency>90%). A series inductor can help with the high frequency part of the fluctuation spectrum. 16.522, Space Propulsion Lecture 11-12 Prof. Manuel Martinez-Sanchez Page 7 of 18 4. Constricted Arc in a Parallel Flow This is the desired configuration for an arcjet (Fig. 4). Here, the heat which arrives through internal conduction at the arc periphery is removed by convection in the outer flow. The gas which comes into contact with the arc is heated as it moves axially, and some of it reaches the ionization threshold T e and becomes part of the arc itself. This means the arc grows with distance, just as the thermal boundary layer adjacent to a heated surface does. Fig. 4. Constricted Arc in Arcjet Let ρu() e be the axial mass flux at the arc’s edge (r=R a ). The amount of “new” arc flow in a distance dx is then ρu() e 2πR a dR a dx dx . The energy absorbed in heating this gas from T out (the buffer gas temperature) to T e must be supplied by the arc power dissipation. Using (7), we then have ρu() e 2πR a dR a dx c p T e ? T out ()= 4πk c T c ? T e () (16) In order to estimate ρu() e , we now make the further assumption that the quantity ρu 2 is independent of r at a given x. This is motivated by the observation that in a parallel, 16.522, Space Propulsion Lecture 11-12 Prof. Manuel Martinez-Sanchez Page 8 of 18 inviscid flow, p +ρu 2 remains constant along the streamlines, while p is itself independent of r. Thus ρu 2 develops gradually at the expense of p, and this should lead to a radius-independent ρu 2 . Our flow is not exactly parallel or inviscid, but the approximation (verified from numerical 2-D solutions) is good enough for the present purposes. We then have, ρu = p ρu 2 () = p ρu 2 ( ) R g T (17) showing that most of the mass flow must occur in the cool, outside gas, since the numerator in (17) is independent of r. We then have, ρu () e ρu() out = T out T e (18) and the problem is now to calculate the gas flux ρu( ) out in the buffer gas. The simplest possible approximation is to state that all of the gas flow is carried by this uniform buffer flow: ρu() out ? Y m π R 2 ? R a 2 () (19) This yields ρu() e ? Y m π R 2 ? R a 2 () T ou t T e (20) and, substituting into (16), Y m π R 2 ? R a 2 () T ou t T e R a dR a dx c p T e ?T out ()? 2k c T c ?T e () The quantity T c -T e depends on current and arc radius through Eq. (6). Substituting, 16.522, Space Propulsion Lecture 11-12 Prof. Manuel Martinez-Sanchez Page 9 of 18 Y m π R 2 ? R a 2 () T ou t T e R a dR a dx c p T e ?T out ()? 2k c I π 2ak c R a or R a 2 R 2 ? R a 2 dR a dx = I 2k c a T e T out Y m c p T e ? T out () (21) We again non-dimensionalize using as in (9), (11), plus a non-dimensional distance I * ,r a x * = x x ref ; x ref = 1 2π Y m c p k c T ou t T e (2) with the result r a 2 1?r a 2 dr a dx * = T e ?T w T e ?T out T w T out ? ? ? ? ? ? I * (23) Since T e is much greater than either T w or T out (which are similar) the temperature function on the right of (23) is close to 1, and will be ignored. Integrating (23), from gives r a o()= r ao x * = 1 I * ln[ 1+r a 1?r a ?r a ] r ao r a (24) which is an inverse expression for radius vs. distance. For r a <<1 the bracketed quantity in (24) is approximated by Taylor expansion as r a 3 3 which indicates rapid arc growth near its upstream (cathodic) end. Therefore it is allowable to use r ao ?o in (24), with only minor effect on length. The remaining important question is the determination of the pressure (either P, the total pressure from upstream, or the pressure P ce at the constrictor exit) required for a flow Y m , given the current I and geometrical data. For this, recall the assumption that the outside layer carries all the flow, and is undisturbed, because the arc heat has not yet penetrated to it. It then is a subsonic ideal gas flow in a contracting area π R 2 ? R a 2 ( ) , and will reach sonic conditions at the constrictor exit, provided the initial nozzle divergence 16.522, Space Propulsion Lecture 11-12 Prof. Manuel Martinez-Sanchez Page 10 of 18 is sufficient. In terms of the total pressure P t , and temperature T t (the upstream conditions), we must then have Y m =Γ P t π R 2 ? R a, L () R g T t (25) with Γ= γ 2 γ +1 ? ? ? ? ? γ +1 2 λ?1() ? 2 3 (26) and with R a,L representing the arc radius at x-L, the constrictor exit. A non-dimensional mass flow rate is defined as Y m * = Y m Y m ref ; Y m ref =Γ P t πR 2 R g T t (27) and (25) reads then (28) Y m * =1?r a, L 2 These equations can be used in a variety of ways. Two of them are: (a) “Design” mode. Given P t , T t , I, Y m , R, etc., find the required constrictor length L (or given L, find R). (b) “Analysis” mode. Given T t , I, Y m , L, R, etc., find the pressure P t . This then will determine the thrust. Example For H 2 gas, say a=0.8Si/m/K,T e =7000K and (eyeballing from the chart handed out in class), k ? 7w / m / K. Use also c p ? 5×10 4 J / Kg. If the constrictor dimensions are R=1.25mm, L=8mm, and the buffer gas temperature is T t ? T w ? 500K , find the chamber pressure P t at operating conditions Y m = 0.1g / s, I =100A. From (9), I ref =π ×1.25×10 ?3 2×0.8× 7 7000?500( )= 85.4Amp I * = 100 85.4 = 1.171 From (22) 16.522, Space Propulsion Lecture 11-12 Prof. Manuel Martinez-Sanchez Page 11 of 18 x ref = 1 2π 10 ?4 ×5×10 4 7 500 7000 =0.0304m x * = 8 30.4 = 0.263 From (24) then, 0.263 = 1 1.171 ln 1+ r a, L 1? r a, L ? r a,L ? ? ? ? ? ? which can be solved to r a, L =0.805. Then, from (28) Y m * = 1? 0.805() 2 = 0.352 ; Y m ref = Y m Y m * = 0.284g / s and from (25), (27), Γ γ = 1.4()= 0.685, and Y m ref = 0.685 P t π 1.25×10 ?3 () 2 8.31 0.002 ×500 = 2.33×10 ?9 P t P t = 0.284×10 ?3 2.33×10 ?9 = 1.218×10 5 N / m 2 = 1.202atm] ? ? ? From the paper AIAA-90-2531, the experimental value under these conditions is 1.16 atm, a good agreement. Repeating for other conditions would yield: Y m g/ s() I(A) R a R (calc.) P t (atm) (calc.) P t (atm) (data) 0.1 100 0.805 1.20 (1.16) 0.2 100 0.686 1.60 (1.80) 0.3 100 0.616 2.05 (2.55) 0.1 60 0.718 0.88 (0.97) 0.1 150 0.870 1.74 (1.28) Table 1 As expected, the agreement is only approximate, but trends are well predicted. Notice how increasing the flow reduces the arc radius (higher buffer layer density, less radius growth for similar heated mass). The increase in pressure with current is also indicated, but the model exaggerates the effect. 16.522, Space Propulsion Lecture 11-12 Prof. Manuel Martinez-Sanchez Page 12 of 18 5. VOLTAGE AND POWER CALCULATION The axial field is still related to arc radius by E * = 1 r a . The voltage accumulated in the constrictor is then V const. = Edx = E ref o L ∫ x ref E * o L * ∫ dx * changing integration variable to r a , dx * =dr a / dr a dx * ? ? ? ? , V const. * ≡ V const . E ref x ref = 1 r a o r a, L ∫ r a 2 I * 1?r a 2 () dr a which integrates to V const. * = 1 2I * ln 1 1?r a, L 2 ? ? ? ? ? ? (29) One difficulty in modeling arcjets is the question of precisely where the arc attaches to the nozzle. This is still largely an unresolved issue, but it is experimentally observed often that the attachment point is well downstream of the constrictor exit, in the supersonic part of the flow. We will construct a crude model only for the voltage drop in this portion. The assumption will be that the arc growth beyond the constrictor end is simply due to the general flow expansion: R a x > L()? R a, L R N x() R c (30) In addition, the axial field is still governed by Eq. (4). Figure 5 Then 16.522, Space Propulsion Lecture 11-12 Prof. Manuel Martinez-Sanchez Page 13 of 18 ?V nozzle = 22 k c a dx R a x() L x>L ∫ and since dx = dR N tan? , using (3) we can integrate ?V nozzle ? 22 k c a 1 r a, L tan? ln R N,att. R c ? ? ? ? ? ? which can be evaluated if an empirical value for R N, att. (or for x att ) is available. In addition, the near-anode region of the arc, where electrons must traverse the normally cool wall boundary layer, is very resistive (although Ohmic heating does elevate the electron temperature there). As a consequence, an additional anodic drop occurs, which, once again, is beyond this model’s capabilities. Approximate empirical values are ?V anode ? 25 V water - cooled copper anode 0 -10 V Hot anode ? ? ? (3) Finally, a similar, but smaller cathodic drop also occurs near the cathode tip, this time due to the necessity for the electrons emitted from the cathode to gather sufficient energy to start the ionization process. This cathode drop is of the order of 1/2 to 1 times the gas ionization potential: (34) ?V cathode ? 6?12 V Example. We continue with the conditions of the Example in Sec. 4: H 2 , R c =1.25mm, L=8mm, I=100A Y m = 0.1g / s, T w =500K. We had found T a,L = 0.805, I * = 1.171 for these conditions. Suppose now we have a ? = 20 o nozzle, and that attachment is observed 4mm downstream of the constrictor: R att R c = 1+ 4tan 20 o 1.25 = 2.17 We also have E ref = 2 k c a 1 R c = 6690Volt, , and x ref = 0.0304m. Using now Eqs. (29)- (34), we estimate V const. = 6690× 0.0304 2×1.171 ln 1 1? 0.805() 2 = 90.7Volt 16.522, Space Propulsion Lecture 11-12 Prof. Manuel Martinez-Sanchez Page 14 of 18 V nozzle == 2 2 × 7 0.8 1 0.805tan20 o ln 2.17()= 31.5Volt ?V anode ? 25Volt ?V cathode ? 8Volt Total: V = 154.7 Volt The experimental value (AIAA 90-2531) is 148 Volt. Again, we can repeat this process for a few other conditions, as we did in Table 1, with the results: Table 2 Y m (g/ s) I(A) V CONSTR. V( ) V NOZZ. V( ) ?V a+c V( ) V(V) V exper. (V) 0.1 100 90.7 31.5 33 154.7 (148) 0.2 100 110.5 37.0 33 180.5 (180) 0.3 100 124.3 41.2 33 198.5 (203) 0.1 60 104.8 35.3 33 173 (160) 0.1 150 81.9 29.1 33 144 (142) The power absorbed by the thruster is simply now the product of the given current and the total calculated voltage P (35) = VI 6. THRUST CALCULATION Following the arc attachment, near the constrictor exit, the gas is very non- uniform in temperature, having “inherited” the arc’s termal distrubtion. From that point, each streamline will ideally expand in the nozzle to the nozzle exit pressure, P e . Since the starting point of the expansion had a uniform pressure and a uniform Mach number (despite the non-uniform T), and the exit condition also has uniform pressure, we must also have a uniform exit Mach number, M e . Each streamtube expands in accordance to the same Mach number all accross, and so the overall area ratio (exit/constrictor exit) is, as usual 16.522, Space Propulsion Lecture 11-12 Prof. Manuel Martinez-Sanchez Page 15 of 18 A e A ce = 1 M e 1+ γ ? 1 2 M e 2 γ +1 2 ? ? ? ? ? ? ? ? γ +1 2 γ ?1() (36) and the exit pressure is related to the stagnation pressure as P e P t = 1 1+ γ ?1 2 M e 2 ? ? ? ? γ /γ ?1 (37) We can then easily calculate the thrust by considering the exit plane: F = ρ e u e 2 + P e ( ) A e = P e A e 1+γM e 2 ( ) (38) The usually defined thrust coefficient is then c F = F P t A ce = P e P t A e A ce 1+γM e 2 ( ) (39) and, as can be verified from (36) and (37), this is identical to what would be obtained for an ordinary rocket for the same M e or area ratio. It is interesting, however, to examine the specific impulse c as well, or the “characteristic velocity”, c * = c / c F . The mass flow rate can be calculated at the (sonic) constrictor exit: Y m = ρudA ∫ () ce = P ce γ R g dA ce T ce ∫ (40) or, in terms of averaged quantities ψ = 1 A ψdA ∫ ? ? ? ? Y m = P ce A ce γ R g 1 T ce (41) For a given P t and A ce , F=c F P t A ce is as it would be if the temperature were uniform in the constrictor. But (41) indicates that Y m is affected by the T ce non-uniformity, and so will then c = F / Y m . For comparison, consider a uniformly heated gas stream that would use the same power per unit mass, and hence would have the same average temperature T ce . This hypothetical rocket would have a specific impulse c , and 16.522, Space Propulsion Lecture 11-12 Prof. Manuel Martinez-Sanchez Page 16 of 18 c c = Y m Y m = 1/ T ce() 1 T ce = 1 T ce () 1 T e ? ? ? ? ? ? (42) It can be proven in general that the quantity in the denominator of (42) is always greater than unity, unless T ce is uniform, and so our arcjet will deliver lower specific impulse than if the heating were uniform. This is a manifestation of a very general principle in propulsion: the highest thrust per unit power (or the least power per unit thrust) is achieved when the exit stream is uniform. As an example, suppose the inner 1/2 of the constrictor exit area is at 10,000K the rest at 1000K. We then have T = 10,000 +1000 2 = 5,500K, 1 T = 1 10,000 + 1 1000 2 = 0.02081. This gives c c = 1 74.16 × 0.02081 = 0.648 which is a substantial thrust loss. Performance example. We return to the hydrogen arcjet at the 0.1g/s, 100 Amp condition. The area ratio is now given as A e A ce =100 . We use γ =1.35 to represent the partially heated gas, and from (36), M e = 6.338. From this, and the previously calculated P t =1.20atm (Table 1), we find P e =0.000324atm. The thrust, from (38) is then F=0.890N The experimental value (AIAA 90-2531) at the power P = 100A×154.8V = 15.5KW is F=0.75l N. Our overprediction is probably related to ignoring friction losses and boundary layer blockage, not entirely negligible at the small Reynolds numbers of these devices. Once the power and the thrust are calculated, the efficiency follows from η= 1 2 F 2 Y m P = 1 2 0.890 2 10 ?4 () 1.55×10 4 () = 0.256 and the specific impulse is 16.522, Space Propulsion Lecture 11-12 Prof. Manuel Martinez-Sanchez Page 17 of 18 Isp= 1 9.81 0.890 10 ?4 = 907sec As before, we can generalize to other conditions, as follows Y m g/ s( ) I(A) F(N) (Calc.) F(N) (measured) η (Calc.) Isp(s) (Calc.) 0.1 100 0.890 (0.75) 0.256 907 0.2 100 1.196 (1.07) 0.195 607 0.3 100 1.52 (1.37) 0.194 516 0.1 60 0.652 (0.64) 0.205 664 0.1 150 1.29 (0.85) 0.385 1315 Table 3 The high current case appears anomalous (following from our overprediction of P t in Table 1), but all others are quite reasonably predicted. 16.522, Space Propulsion Lecture 11-12 Prof. Manuel Martinez-Sanchez Page 18 of 18