3768-26-T Report of BAMIRAC APPROXIMATE METHODS FOR CALCULATING THE STRUCTURE OF JETS FROM HIGHLY UNDEREXPANDED NOZZLES T. C. ADAMSON, JR. Associate Professor Department of Aeronautical and Astronautical Engineering June 1961 Infrared Laboratory 5tittUe so Scuece a<W 7ec4i,4o THE UNIVERSITY OF MICHIGAN Ann Arbor, Michigan

NOTICES Sponsorship. The work reported herein was conducted by BAMIRAC of the Institute of Science and Technology for the U. S. Air Force Cambridge Research Laboratories, Contract AF 19(604)-7350, implementing Advanced Research Projects Agency Order Number 30. Contracts and grants to The University of Michigan for the support of sponsored research by the Institute of Science and Technology are administered through the Office of the Vice-President for Research. ASTIA Availability. Qualified requesters may obtain copies of this document from: Armed Services Technical Information Agency Arlington Hall Station Arlington 12, Virginia Final Disposition. After this document has served its purpose, it may be destroyed. Please do not return it to the Institute of Science and Technology.

Institute of Science and Technology The University of Michigan PREFACE BAMIRAC, the-Ballistic Missile Radiation Analysis Center, is a facility established at the Institute of Science and Technology under Contract AF 19(604)-7350 with the Geophysics Research Directorate of the Air Force Cambridge Research Laboratories (AFRD). The contract is an implementation of Advanced Research Projects Agency Order Number 30. The objective of BAMIRAC is to function as a technical information center on phenomenology, theory, and technology pertaining to the fundamental phenomena associated with ballistic missiles and space vehicles which may be significant in any way to the formulation of defense measures against missile and space-vehicle systems. BAMIRAC collects and processes information concerned with electromagnetic and acoustic radiation emanating from or caused by ICBM's or IRBM's during their entire trajectory from launch to impact. The information includes field measurements, laboratory studies, and theoretical studies, and is available at the BAMIRAC reference library for use by representatives from all organizations presenting a properly authorized request. BAMIRAC conducts analyses in which experimental and theoretical results are evaluated and examined for correlations. Some theoretical and experimental investigations are carried out, and the results are combined with the technical information obtained from outside sources. In its capacity as a technical information center, BAMIRAC disseminates information by means of such technical media as: abstracts, indexes, bibliographies; technical reports and journal articles; technical meetings; Proceedings of the Anti-Missile Research Advisory Council; and technical guidance services. BAMIRAC is under the technical direction of the Infrared Laboratory. It draws also, however, upon the capabilities of the Computation Department of the Institute of Science and Technology, and upon those of the Aircraft-Propulsion Laboratory of the Department of Aeronautical and Astronautical Engineering within The University of Michigan's College of Engineering. iii

Institute of Science and Technology The University of Michigan ACKNOWLEDGMENT The author is indebted to E. K. Latvala of ETF, ARO, Inc., for the use of photographs of jet boundaries and shocks from which the data used in Figure 7 were obtained. IV

Institute of Science and Technology The University of Michigan CONTENTS Notices..................................... ii Preface...................................... iii Acknowledgment................................. iv List of Figures.................................. vi List of Symbols..................................vii Abstract...................................... 1 1. Introduction................................ 1 2. Jet Boundary.............................. 2 3. Intercepting Shock.............................. 14 4. First Mach Disc............................... 20 5. Discussion of Results.............................. 21 References.................................... 22 V

Institute of Science and Technology The University of Michigan FIGURES 1. Aerodynamic Structure of Jet from Highly Underexpanded Nozzle...... 2 2. Reflection of Characteristics at Jet Boundary................ 3 3. Notation Used in Jet-Boundary Calculations................. 4 4. Spherical Cross-Sectional Areas Used in Boundary Approximation..... 5 5. Comparison of Approximate Solutions Based on Fictitious-Jet Method with Solution by Method of Characteristics from Love and Grigsby................................... 7 6. Comparison of Approximate Solutions Based on Fictitious-Jet Method with Solution by Method of Characteristics from Love and Lee................................... 8 7. Comparison of Approximate Solutions Based on Fictitious-Jet Method with Experimental Values Obtained from Schlieren Photographs by Latvala........................ 9-11 8. Application of Approximation Based on Spherical Source-Flow to Jet-Boundary Computation.......................... 12 9. Comparison of Approximate Solutions Based on Circular-Arc Method with Solution by Method of Characteristics from Love and Grigsby...13 10. Comparison of Approximate Solution Based on Radial Source-Flow ("Exponential" Solution) with Solution by Method of Characteristics from Love and Grigsby.......................... 14 11. Formation of Envelope of Right-Running Characteristics Reflected from Boundary...............................15 12. Notation Used in Intercepting-Shock Calculation...............16 13. Geometric Interpretation of Circular-Arc Approximation to Boundary and Intercepting Shock...........................19 vi

Institute of Science and Technology The University of Michigan SYMBOLS A area M Mach number P pressure r radial coordinate r dimensionless radial coordinate, r/r R radius of curvature ds element of length along jet boundary x axial coordinate x dimensionless axial coordinate, x/r n a angle made by jet boundary with axis at nozzle edge, immediately after expansion y ratio of specific heats, CP/CV Q half angle made by nozzle wall with axis at exit A Mach angle, A = sin 1/M v Prandtl-Meyer angle (D angle made by jet boundary with axis Superscript * conditions at nozzle throat Subscripts b boundary of jet plume e conditions at jet boundary at nozzle exit immediately after expansion initial radius of equivalent radial flow corresponding to first point of jet boundary where f = a (used only with radius of curvature) n conditions at nozzle exit before expansion oo ambient conditions external to jet flow VII

I

APPROXIMATE METHODS FOR CALCULATING THE STRUCTURE OF JETS FROM HIGHLY UNDEREXPANDED NOZZLES ABSTRACT Approximate expressions are presented for the calculation of the jet-boundary and intercepting-shock locations of a jet flow from a highly underexpanded rocket nozzle. The gas is assumed to be an inviscid nonreacting perfect gas with constant specific heats, exhausting to an ambient atmosphere. One new method and two previously known methods of approximating the jet boundary are discussed, and numerical examples of each are compared with calculations based on the method of characteristics and with experimental results. A new method for approximating the intercepting-shock position is presented and compared with theoretical and experimental results. Although the accuracy varies, depending on the approximation employed, it is shown that with any of the methods presented, the initial part of the jet boundary and intercepting shock is well represented by the approximate solutions. The location of the first Mach disc is discussed briefly. 1 INTRODUCTION The hot exhaust gases issuing from a highly underexpanded rocket nozzle expand rapidly, forming a plume with a very complex structure. Although very little reaction takes place within the jet itself, the aerodynamic structure of the hot gases is complicated by internal shocks. Further, in the viscous mixing region along the periphery of the plume, considerable reaction does take place. Because the thickness of the mixing region is small compared to the diameter of the jet, at least in the region near the nozzle, it is possible to neglect viscous and reaction effects in considering the aerodynamic structure. Thus, the problem is generally divided into two parts. First, the jet is taken to be an ideal gas, and the internal structure and boundary are calculated. Second, a mixing problem with chemical reaction is considered, with one fluid being the external flow and the other being the jet fluid at those conditions occurring along the idealized boundary. Also, in order to further simplify the problem, first calculations usually involve the assumption that the external fluid is at rest, so that a constant-pressure boundary condition results. The effects of altitude may still be included by varying the external pressure, but effects due to the velocity of the vehicle are neglected. Presumably, methods devised for the ambient atmosphere case may be extended to cover the case of an external streanl with a relative velocity. 1

Institute of Science and Technology The University of Michigan In this report, only the aerodynamic structure is considered. An inviscid fluid with constant specific heats is assumed to expand from a rocket nozzle into an ambient atmosphere. Only approximate solutions to the problem of calculating the jet boundary and internal shock structure are presented. By approximate solutions are meant those solutions which do not make use of the exact method of characteristics. Such solutions are important because they will allow rapid, relatively accurate estimates to be made of the growth and development of the jet structure as altitude, fluid, and nozzle-exit conditions are varied. Further, because the solutions are functional rather than numerical, they will allow more insight as to which parameters are important in assessing a given variation in structure. Finally, it is clear that even with modern computing machines, the accuracy of the characteristic solution decreases the further downstream of the nozzle the solution is continued; hence simpler, less expensive approximations become attractive alternatives. In order to illustrate the terms used hereafter, a sketch of a typical underexpanded jet is given in Figure 1. _ Constant Pressure Flow Direction /g Immediately After/ in Sho c k Expansion at / Nozzle Lip 4/ Slip Line N/oi — - F Dotted Lines Show ^"-^_^.- ^SS^^11^^,~ ~-^ P~nanoinn T~an t M ach Disc \,, jNozzle~ o \ Expansion Fan at isio n at \\\ Nozzle Lip Flow Direction at \ End of Nozzle \ \ \ --------- -------- ---- --- - Center Line FIGURE 1. AERODYNAMIC STRUCTURE OF JET FROM HIGHLY UNDEREXPANDED NOZZLE 2 JET BOUNDARY Previous papers (References 1-4) describe the original work done on approximating the jet boundary. Hence only a short review of the results and improvements will be presented in this report. Since viscous effects are generally neglected in aerodynamic jet calculations, either exact or approximate, the boundary of the jet becomes the bounding streamline along which the pres2

Institute of Science and Technology The University of Michigan sure is constant. Since the turning of the boundary is done isentropically, the total pressure and total temperature are constant. Hence, the Mach number, velocity, and density are also constant along the jet boundary. Figure 2 illustrates the changes in Mach number, for example, which occur due to the interaction of the boundary and left- and right-running characteristics which, in this case, may be pictured as weak expansion and compression waves, respectively. Although variations in Mach number, pressure, etc., exist across the weak expansions, they are cancelled by identical negative variations across the weak compressions. The weak expansion is reflected as a weak compression because of the constant-pressure condition. Figure 2 illustrates physically the mechanisms which occur at the boundary corresponding to the results found in a characteristic solution. -,D/v& =X M=M + AM- = M1 M1 / M1 + AM Right-Running Characteristic / (Weak Compression) / Left-Running Characteristic (Weak Expansion) FIGURE 2. REFLECTION OF CHARACTERISTICS AT JET BOUNDARY The approximate methods so far developed make use of the condition of constant pressure but replace the actual mechanisms described above by two approximations. Thus the expansion process is assumed to be that due to an infinitesimal increase in area, in an isentropic channel flow, and the compression back to atmospheric pressure is that due to turning the supersonic flow through an infinitesimal change in angle. The two effects balance each other, resulting in a constant pressure. However, it is important to note that the isentropic expansion process in a channel implies a change in Mach number. This is true only insofar as the calculations are concerned. That is, the resulting boundary configuration is that of an isentropic channel flow wherein the differential change in pressure with respect to change in angle, at a given boundary point, is that which would be given by the linearized expression for the change in pressure due to a weak wave turning a supersonic flow through the given differential angle. Hence the approximate boundary shape is found from a fictitious isentropic nozzle fulfilling the above-mentioned boundary conditions on the pressure. The results of these calculations (Reference 1) is that along the boundary v + = constant = v + a (1) where v is the Prandtl-Meyer angle associated with the Mach number at a given point on the boundary, q is the angle made by the boundary with respect to the flow axis, v is the Prandtl3

Institute of Science and Technology The University of Michigan Meyer angle of the flow at the nozzle exit immediately after expanding to atmospheric pressure, and a is the initial inclination of the free-jet boundary with respect to the flow axis. The boundary is calculated using Equation 1 and the relation dr - = tan f (2) db as well as an equation relating v to r. Figure 3 indicates the notation used in Equation 2o As Latvala (Reference 2) has pointed out, the best results are obtained if, when the fictitious channel areas are calculated, spherical rather than planar areas are considered. Thus, with the notation given in Figure 4, the spherical area corresponding to a given point on the boundary is 27r2 + cos (3) * Further, if A, the area corresponding to M = 1, is assumed to be constant all through the fictitious flow (all along the calculated boundary), then A e 2 ( + cosa\ ( * *b 1 + cos A A where A is the spherical area at the exit of the nozzle (at the beginning of the free jet) where e f = a, and rb - (5) r nt _ Bnd Boundary ozz|le n I — - -- _ -_,1 -_ _ _ _ _ _ - ---- Center Line Xb x FIGURE 3. NOTATION USED IN JET-BOUNDARY CALCULATIONS. Ga = 6 + ve - vn. vn = PrandtlMeyer angle at nozzle exit. ve = Prandtl-Meyer angle after expansion. 4

Institute of Science and Technology The University of Michigan Boundary I/ \ Projection of Spherical / Surface in Meridian Plane FIGURE 4. SPHERICAL CROSS-SECTIONAL AREAS USED IN BOUNDARY APPROXIMATION is the dimensionless boundary radius. r is the actual radius at the nozzle exit. Thus, * n, A /A is the area ratio corresponding to Me, and A/A is the area ratio corresponding to some M along the boundary. A dimensionless axial distance xb is defined as b b - (6) r n In order to use Equation 2, it is necessary to calculate d(A/A )/drb. This involves the derivative dq/drb. However, from Equation 1, d dv _ dv dM d(A/A ) dr dr dM * dr db db dd(A/A) drb Now, A 2 )(y+l)/2(y-) (1 1M2)( )/2(1) 1*-^ ) *^(1+ L^1M) (7) and dv M -1 dM M(1 + M2) 5

Institute of Science and Technology The University of Michigan For details see any standard reference on compressible flow, e.g., Reference 5. Thus, d0/drb and the derivative of Equation 4 may be combined to give d(A/A) _ 2 A sin d -- 1+ 2 j (9) drb ^bA - 1 (1 + cos ) Then, integrating Equation 2, one finds for xb, A/A =lj |b sin A x, =? J --- ^ -- 1 + d (10) b 1 2(A/A ) tan \IM - 1 (1 + cos (A0 with A 1 + cos (4') b * \1 + cos a) a= a-(v- v) (1') and a= 6 + - v (11) e n being the necessary additional relations for the boundary calculations. Simpson's rule for numerical integration may be used to find xb. The integrand at any point is found as follows: (1) Find ve, v for given M and M and calculate a (Equation 11). (2) Choose an A/A (3) Knowing A/A, find M and v from compressible flow table for the proper y. (4) Calculate 0 (Equation 1'). (5) Calculate rb (Equation 4') (6) Calculate integrand. Latvala (Reference 2) has used difference equations rather than the above formulation, but it is believed that this numerical integration is nearly as easy to use and gives more accurate results. Calculations for a number of initial conditions are compared with experimental results and solutions calculated by the method of characteristics in Figures 5 to 7. The interceptingshock results presented on these graphs will be considered later. The curves labeled experimental on Figure 7, were obtained from schlieren photographs; some error certainly exists as a result of transferring the results from photograph to graph, although great pains were taken to make this error as small as possible. It should be noted that this approximate method does not give accurate results insofar as calculation of the maximum radius of the jet is con6

Institute of Science and Technology The University of Michigan cerned. However, for the initial part of the boundary, the comparison both with experimental results and the method of characteristics is quite satisfactory, the maximum error in radius being approximately 10%. There are other approximate methods which promise great simplification and better accuracy. However, at the present time these other methods are incomplete in that considerable experimental data are required for their use. For example, use of the so called circular-arc approximation (References 2 and 3) depends on the knowledge of the maximum radius of the jet, or some comparable piece of information. Hence, at this time, such methods are useful only for filling in the boundary between the nozzle exit and the maximum radius. For example, they would not be useful for predicting the growth of the free jet as the pressure ratio increases. Of course, the latter prediction is most desirable, and so considerable effort has been expended in attempting to analytically calculate the maximum radius without performing a compiete characteristic solution. Up to the present time, such efforts have been unsuccessful. However, it is informative to study the derivation of these other methods and to assess their accuracy in the event that the necessary additional information becomes available. In the following derivation, the turning of the flow along the jet boundary is assumed to be due to two separate mechanisms, in a manner similar to that considered in the first method. | Boundary.r = r/r / >-, = C Intercepting Shock r = r/rn Characteristics Solution ---- Approximate Solution 0 1 2 3 4 5 x = x/r n FIGURE 5. COMPARISON OF APPROXIMATE SOLUTIONS BASED ON FICTITIOUSJET METHOD WITH SOLUTION BY METHOD OF CHARACTERISTICS FROM LOVE AND GRIGSBY. (Reference 3.) Mn = 1; y = 1.40; Pn/Po = 10; 6 = 0. 7

Institute of Science and Technology The University of Michigan 10 Jet Boundary Intercepting Shock 7 _/ -r=r/rn /3t~ 0~~ — -tCharacteristics Solution / --- Approximate Solution 2 10 1 2 3 4 5 6 x = x/r n FIGURE 6. COMPARISON OF APPROXIMATE SOLUTIONS BASED ON FICTITIOUS-JET METHOD WITH SOLUTION BY METHOD OF CHARACTERISTICS FROM LOVE AND LEE. (Reference 6.) Mn = 2.5; y = 1.4; Pn/Poo = 1346; 6 = 15~. 8

Institute of Science and Technology The University of Michigan 10 r = r/ rn / 4 / / /.Experimental Je Boun / --- I Approximate Solution 0 1.0 2.0 3.0 4.0 5.0 6.0 x = x/r n FIGURE 7a. COMPARISON OF APPROXIMATE SOLUTIONS BASED ON FICTITIOUS-JET METHOD WITH EXPERIMENTAL VALUES OBTAINED FROM SCHLIEREN PHOTOGRAPHS BY LATVALA. (See Reference 2.) Mn = 1.5; y = 1.4; Pn/Po, = 305; 6 = 15~.

Institute of Science and Technology The University of Michigan 10 i/ / 8 Jet Boundary I. Intercepting Shock 3 —- Experimental 1 —-- Approximate Solutions -1.0 0 1.0 2.0 3.0 x = x/r n FIGURE 7b. COMPARISON OF APPROXIMATE SOLUTIONS BASED ON FICTITIOUS-JET METHOD WITH EXPERIMENTAL VALUES OBTAINED FROM SCHLIEREN PHOTOGRAPHS BY LATVALA. (See Reference 2.) Mn = 1.5; y = 1.4; Pn/Po = 2510; 6 = 15~. 10

Institute of Science and Technology The University of Michigan Jet Boundary r=r/r Intercepting Shock 2 7 —-- /r --- Experimental -- -- Approximate Solutions 0 1.0 2.0 3.0 4.0 5.0 x=x/r FIGURE 7c. COMPARISON OF APPROXIMATE SOLUTIONS BASED ON FICTITIOUSJET METHOD WITH EXPERIMENTAL VALUES OBTAINED FROM SCHLIEREN PHOGRAPHS BY LATVALA. (See Reference 2.) Mn = 2.97; y = 1.3; Pn/Poo = 30; 6 = 15~. Thus, at a given point on the boundary, the flow expands along a surface (represented by a line in the meridian plane) tangent to the boundary at this point, and turns through an angle sufficient to compress it back to the original pressure. Furthermore, the pressure decrease in the expansion process is taken to be that due to an increase in the radial length of a supersonic spherical source-flow. In a spherical source-flow, the Mach number is a function of radius alone. Hence the pressure gradient for the given expansion is calculated at the radius corresponding to the Mach number at the jet boundary. In essence, this method may be pictured as follows: a radial element of a spherical source flow, with length R, is placed tangent to the boundary so that its end point coincides with the point in question. Then the pressure decrease due to a differential increase in R is balanced by the pressure increase due to a differential decrease in 0 (see Figure 8 for notation). It is important to note that the two mechanisms need not occur in the order given. One could consider, also, a decrease in angle followed by an increase in the length of the tangential-line element. In view of the above remarks, the equation which holds along the jet boundary is dP =0 = dR + d (12) M=M M=M e e 11

Institute of Science and Technology The University of Michigan ddrb Boundary FIGURE 8. APPLICATION OF APPROXIMATION BASED ON SPHERICAL SOURCEFLOW TO JET-BOUNDARY COMPUTATION where 3P 2yPM2 aR ~ ~ ~ ~ ~ ~ ~!2 (13) 9R R(M2- 1) and 3P PM2 (14) -- Since M = M = constant, if a is defined as ~M 2 - 1 e a = e2 2tan (15) 2 2 tan Aji then Eauation 12 becomes dR -= aR (16) Within the order of the approximation, dr. dR = ds =- (17) sin Two possibilities exist in the interpretation of Equation 16. First, if one considers that R is a function of M alone, and due to the constant M is thus constant, then Equations 16 and 17 may be combined and integrated to yield rb - 1 = aR(cos - cos a) (18) This is the equation which gives a circular-arc approximation with radius of curvature Rb = aR (19) On the other hand, it may be argued that in reality R/R = f(M), where R is the radius corresponding to M = 1, and if one allows the sonic radius to vary in the same way as R so that 12

Institute of Science and Technology The University of Michigan R/R = constant, then R is not necessarily constant, and Equation 16 may be integrated with the following result: R Rea(a-) (20) R = R.eaThus (20) where R. is the radius corresponding to 0 = a. Thus, dR = -aR.ea, and with Equation 17 one may solve for rb. rb - 1 = R a(a- ) (a sin 0 + cos 0) - a(sin a + cos a (21) (1 + a2) so that, in this case, there is an exponential dependence of rb on 0. It should be noted that since dxb/ds = cos 0, one may easily find equations for xb for each of the above cases. Although the above two cases are encouraging in that they result in equations which agree well with experimental and characteristic-solution results, they both require additional information in that R or R. must be known. Efforts to derive relations defining R or R. in terms of the nozzle parameters have been unsuccessful. In order to compare each of the above results with known solutions, R and R. were calculated by substituting a known value for the maximum radius for a given nozzle and flow condition in Equations 18 and 21. Using the corresponding equations for xb in each case, the boundaries were calculated and compared with a characteristic solution for a sonic nozzle with pressure ratio P /P =10 and y = 1.4. The results are n 0o known in Figures 9 and 10. It is clear that the circular-arc approximation is quite good, and the "exponential" solution is even better. However, until a means of calculating R or R. for a given set of nozzle parameters is discovered, these equations have limited value. _- - Jet Boundary 2 = _______ __ 2_. S.. _ Intercepting Shock r = r/rn Characteristics Solution --- Circular-Arc Approximation 0 1 2 3 4 5 x=x/rn FIGURE 9. COMPARISON OF APPROXIMATE SOLUTIONS BASED ON CIRCULAR-ARC METHOD WITH SOLUTION BY METHOD OF CHARACTERISTICS FROM LOVE AND GRIGSBY. (Reference 3). Mn = 1; y = 1.40; Pn/Poo = 10; 6 = 0. 13

Institute of Science and Technology The University of Michigan 3 _ —--- Characteristics Solution r = r/r 0 Points from Approximate Solution 0 1 2 3 4 5 x = x/rn FIGURE 10. COMPARISON OF APPROXIMATE SOLUTION BASED ON RADIAL SOURCE-FLOW ("EXPONENTIAL" SOLUTION) WITH SOLUTION BY METHOD OF CHARACTERISTICS FROM LOVE AND GRIGSBY. (Reference 3). Mn = 1; = 1.4; Pn/Po = 10; 6 = 0. 3 INTERCEPTING SHOCK A sketch of the typical shock structure in a highly underexpanded jet is shown in Figure 1. The intercepting shock is one which is very weak immediately downstream of the nozzle, but which grows in strength in the downstream direction. Initially, this shock is tangent to the last ray of the expansion fan which exists at the nozzle lip. In a gross fashion one might explain the existence of the intercepting shock by the fact that as the outer part of the jet expands to atmospheric pressure, the inner core of the jet overexpands, and a shock is then necessary to relieve the pressure difference. More precisely, the shock is formed by the coalescence of the many weak compression waves emanating from the jet boundary. The boundary is curved inward as a result of the constant-pressure condition. It is well known that in a characteristic solution, the right-running characteristics from the jet boundary intercept each other and form an envelope (References 3 and 7) called the limiting line (or surface in this case), unless additional calculations involving the shock equations are included. The existence of a discontinuity in this region is implied mathematically by the fact that three simultaneous solutions exist downstream of the limiting line. Generally, this difficulty is overcome by assuming that the envelope, or limiting line, gives the intercept14

Institute of Science and Technology The University of Michigan ing-shock position and by choosing the proper solution so that the calculations may be continued downstream of the envelope. Evidently this assumption is quite good since the comparison of theory and experiment gives good results (Reference 3), i.e., when calculated results are superimposed on schlieren pictures of actual jets, the position of the intercepting shock agrees satisfactorily with the calculated envelope. Hence, in this approximate development, it is assumed that the intercepting shock exists very nearly at the envelope of the right-running characteristics initiated at the shock boundary. A consideration of a typical characteristic solution reveals that although there are small local variations, the right-running characteristics which form the limiting line are nearly straight lines. As the pressure ratio increases, this appears to be a better and better approximation, and so this assumption is made. Further, since the Mach number along the boundary of the jet is a constant, the Mach angle, pL, the angle which the right-running characteristics make with the boundary tangent, is a constant. Therefore, finding the position of the intercepting shock reduces to the calculation of an envelope of straight lines, each of which leaves a known curve with a known angle. In summary, the intercepting shock is the curve which is initially the last right-running characteristic of the expansion fan at the nozzle, and thereafter the envelope of the (assumed) straight right-running characteristics which leave the boundary at the relative angle Ji (Figure 11). Initially, the intercepting shock is of infinitesimal strength, but the shock strength grows as more and more right-running waves coalesce with it. The straight lines which form the envelope are a one-parameter family; for example, prescription of the axial distance to a point on the jet boundary allows one to calculate the Tangents to Boundary - Bounda - < ^^^^ ^ b o u n da —~~ —'-Boundary \ /- Last Right-Running Characteristic | \ D^ of Expansion Fan at Nozzle Envelope -Intercepting Shock -- - Right-Running Characteristics from Boundary. Assumed to be Straight Lines Nozzle Center Line FIGURE 11. FORMATION OF ENVELOPE OF RIGHT-RUNNING CHARACTERISTICS REFLECTED FROM BOUNDARY 15

Institute of Science and Technology The University of Michigan boundary radius and boundary angle. Then, since,i is known, a line can be constructed. The notation employed is given in Figure 12. The equation describing the line which includes the boundary point xb, rb is r - rb (x -. x) tan (- g) (22) Equation 22 may be written in the form f(x, r, c) = 0 where c is the parameter indicating which line is being considered. In order to find the envelope, then, one solves the above equation simultaneously with the relation =0 (23) ac (See, e.g., Reference 8.) In this case, one may choose c to be xb, the axial distance, as rb, Xb, 0, and a, the parameters appearing in Equation 22 are all functions of xb. Then Equation 23 becomes -drb 2 dx + tan ( - ) - (x- xb) sec2 ( ) d(-) =0 (24) dx b dxb However, as noted previously, dr = tan ~ (2) db Tangent to r = r/r \ Boundary ~n ~\ ^ Boundary / d/ /^ Right-Running Characteristic //<L/ ~" (x, r) hf\ (X b rb) No zz le x/ -- x = x/r FIGURE 12. NOTATION USED IN INTERCEPTING-SHOCK CALCULATION 16

Institute of Science and Technology The University of Michigan and since Ai = constant = y e Equation 24 becomes, with the use of some trigonometric identities, x - x -sin /iCcos =1 (1 + tan 1 tan m) (25) -b' (dP/dx b) t a Next it is necessary to calculate do/dxb, the rate of change of the angle which the boundary makes with the axis, with respect to the axial distance. Hence do/dxb depends on the boundary approximation being considered. The method first derived will be used in the following calculations. It is important to keep in mind the fact that, in this approximate method, the boundary Mach number of the fictitious jet is not constant. In order to emphasize this point, the subscript b will be employed to indicate this fictitious boundary Mach number in the following intercepting-shock calculations. From Equation 1, do - -b (26) dxb dxb where vb is the Prandtl-Meyer angle corresponding to Mb. Also, dvb dvb dMb d(A/A*) drb db dMbd(A/A* db db and employing Equations 2, 7, 8, and 9, d _ -2 tan sin l) (27) b~s b - xb rb Mb-1 I -21(1 + cos ) Substitution of Equation 27 into Equation 25 gives the equation for x. The parametric equations which define the intercepting shock are then x=xb +rb sn 1 (cot +tan L M ), si (28) bM.~~~~ bb2 ~- 1 (1 + cos 0) and r = rb + (x - xb) tan ( - /A) (22') Several calculations of the intercepting shock were made for the same cases considered in the jet-boundary examples. They can be seen in Figures 5 to 7, where they are compared to theoretical and experimental results. It can be seen that the downstream distance to the region where the approximate results diverge considerably from either the characteristic 17

Institute of Science and Technology The University of Michigan solution or the experimental curves is roughly equal to the corresponding distance for the jetboundary curves. Before this divergence occurs, the accuracy is somewhat better for the intercepting shock than for the boundary approximations. It is interesting to derive the equations for the intercepting shock employing the circulararc approximation for the jet boundary. In Section 2 it was shown that in this case, substituting Equation 19 into Equation 18, rb - 1 = Rb(cos - cos a) (29) The corresponding equation for x is b = Rb(sin a - sin 0) (30) Hence, do _ 1 (31) dxb Rb cos and substituting Equation 30 into Equation 25 one obtains, x =xb + Rbsin / cos ( - ) (32) Equations 32 and 22 then define the intercepting-shock position. In Figure (9) a comparison of the circular-arc approximation with corresponding intercepting-shock calculation is shown, compared to a characteristic calculation. This comparison indicates that with the circular-arc boundary, the approximate shock radii are larger than shown by the characteristic calculations. Of course, this is a very low pressure ratio; for larger pressure ratios, the agreement should be much better. It is informative to consider the geometry of the boundary and shock structure when the circular-arc approximation is employed. Equations 32 and 22 may be written in the following form: r - r = Rb cos [cos (q- ) -) cos (a - Mi)] (33) x - x0 = Rb cos [sin (a - g) - sin (0 - L)] (34) where r and x0 are the values of r and x when a = a. Thus, r0 = 1 + Rb sin A sin (a - Aj) (35) x0 = Rb sin 11 cos (a - Mi) (36) 18

Institute of Science and Technology The University of Michigan It should be noted that in Equations 33 and 34, as in all equations dealing with the intercepting shock, the parameter 0 refers to the angle made by a tangent to the boundary at the point xb, rb corresponding to point x, r of the shock. The above relations indicate that the intercepting shock is another circular arc with radius R, such that Rs = Rb cos A (37) R has the same center as Rb, and the arc is tangent to the straight line leaving the nozzle lip at an angle a -A with respect to the axis (i.e., tangent to the straight-line approximation to the last ray of the Prandtl-Meyer expansion fan at the nozzle). Figure 13 illustrates the construction of the boundary and the shock in this case. It can be seen both from the equations and from Figure 13 that according to this approximation, the maximum radius of the shock occurs at the same point at which the maximum radius of the boundary occurs. This is not found in characteristic solutions, but is a result of the circular-arc approximation used for the boundary. This fact indicates that the approximate circular-arc equations, while quite good for the boundary, give less accurate results when used for the approximate shock location. However, the error should decrease as the pressure ratio increases, since the proper limit is given as the pressure ratio tends toward infinity. Tangent to Boundary at Nozzle / B1~3 / xf,,,,,\ Sh~ck i \ Rs - -________\ \- \ - l Center Line R FIGURE 13. GEOMETRIC INTERPRETATION OF CIRCULAR-ARC APPROXIMATION TO BOUNDARY AND INTERCEPTING SHOCK 19

Institute of Science and Technology The University of Michigan At the present time, no calculations have been made using the so-called exponential boundary equations (Equation 21) as a basis for the intercepting shock. This will be done in the future. In summary, it is clear that the accuracy of the calculation of the intercepting shock, with the approximate methods described above, depends to a very great extent on the curvature of the jet boundary. With the first boundary approximation made, employing a fictitious nozzle, the curvature is initially very large, so that the calculated shock lies below the proper position. In the circular-arc approximation, the boundary has a constant average curvature which is too small initially, so that the calculated intercepting shock lies above the proper position. Hence, it is desirable to find a closer approximation to the actual jet boundary. In view of the more accurate representation given by the "exponential" boundary approximation, it is believed that a more accurate shock position will also result when the shock relations are based on this boundary approximation. 4 FIRST MACH DISC A method for calculating the position of the first Mach disc in the jet plume from a highly underexpanded nozzle has been given by Adamson and Nicholls (Reference 1). In this method, it was assumed that the pressure behind the shock was atmospheric, and a superposition method was given for calculating the disc position for various nozzle Mach numbers and pressure ratios. The comparison between theory and experiment was very good. Since the original work was published, experiments in a different facility have indicated that the pressure behind the wave is indeed very nearly atmospheric, but that it could vary up to 1.3 times atmospheric. Hence, some spot checks were made on the theoretical calculations, using a back pressure of 1.3 atmospheres. It was found that very little change in the Mach-disc location resulted in the Mach-number and pressure-ratio range covered in the original report. (1 Mnozle < 3; 10 <P /P < 100.) - nozzle ooAnother check on the original calculation was made using the exact relation for the leading Mach line in the flow emanating from the nozzle as opposed to the straight-line approximation used originally. Again no change could be noted in the final results. In order to apply the above-mentioned method, it is necessary that the Mach-number distribution along the center line of the jet be known. While the superposition method reported seems to give good results in the given pressure-ratio range, it is clear that more accurate 20

Institute of Science and Technology The University of Michigan results would be given by employing the Mach-number distribution obtained from a characteristic solution for a given nozzle. This Mach-number distribution would apply for all pressure ratios for this nozzle since the inner-core flow, up to the enclosing shocks, is identical for all pressure ratios. 5 DISCUSSION of RESULTS Three approximate models, each with varying levels of approximation, have been discussed with regard to calculating the boundary and intercepting shock of a jet plume. Only the method first discussed (fictitious-jet model) is completely defined analytically at the present time, in that no additional information beyond the nozzle parameters and pressure ratio is necessary. However, since this model does not predict the maximum diameter reached by the jet, it may be applied only a few jet diameters downstream of the nozzle. Within this region of applicability both the jet-boundary and intercepting-shock locations are predicted satisfactorily. The second model considered (circular-arc approximation) gives excellent results for the jet boundary, presuming the maximum radius is known. However, due to the dependence of the equations derived for the intercepting shock on the curvature of the boundary, relatively poor results are obtained for the intercepting shock. The last model presented (equivalent radial flow with variable radius- "exponential" model) seems to give excellent results for the boundary for the conditions in the one case considered. However, more comparisons, especially at higher pressure ratios, must be made before further comment can be made. If this model is as successful on high pressure ratios, then, since the curvature obtained is a much better approximation to the actual boundary curvature than those given by the other two methods, the intercepting-shock calculations based on this model should give correspondingly better results. However, use of this model is still seriously hindered by the fact that as yet it has not been possible to derive an analytic expression for Ri, the initial radius of curvature. Thus, as in the circular-arc model, additional information is necessary. 21

Institute of Science and Technology The University of Michigan REFERENCES 1. T. C. Adamson, Jr., and J. A. Nicholls, "On the Structure of Jets from Highly Underexpanded Nozzles into Still Air," J. Aeronaut. Sci., January 1959, Vol. 26, No. 1, pp. 16-24. 2. E. K. Latvala, Spreading of Rocket Exhaust Jets at High Altitudes, Report Number AEDC-TR-59-11, Arnold Engineering Development Center, ARO, Inc., Tullahoma, Tenn., June 1959 (UNCLASSIFIED). 3. E. S. Love and C. E. Grigsby, Some Studies of Axisymmetric Free Jets Exhausting from Sonic and Supersonic Nozzles into Still Air and into Supersonic Streams, Report Number NACA-RM-L54L31, National Advisory Committee on Aeronautics, Washington, D. C., May 1955 (UNCLASSIFIED). 4. W. T. Lord, On Axisymmetrical Gas Jets, with Application to Rocket Jet Flow Fields at High Altitudes, Report Number 2626, Royal Aircraft Establishment, Farnborough, England, July 1959 (UNCLASSIFIED). 5. H. W. Liepmann, and A. Roshko, Elements of Gas Dynamics, Wiley, New York, N. Y., 1957. 6. E. S. Love and L. P. Lee, Shape of Initial Portion of Boundary of Super-Sonic Axisymmetric Free Jets at Large Jet Pressure Ratios, Report Number NACA-TN4195, National Advisory Committee on Aeronautics, Washington, D. C., January 1958 (UNCLASSIFIED). 7. R. Courant and K. O. Friedrichs, Supersonic Flow and Shock Waves, Interscience Publishers, New York, N. Y., 1948, p. 390. 8. I. S. Sokolnikoff and E. S. Sokolnikoff, Higher Mathematics for Engineers and Physicists (2nd ed.), McGraw-Hill, New York, N. Y., 1941, p. 280. 22

AD Div. 12/1 UNCLASSIFIED AD Div. 12/1 UNCLASSIFIED Inst. of Science and Technology, U. of Mich.,Ann Arbor I. Title: BAMIRAC Inst. of Science and Technology, U. of Mich.,Ann Arbor I. Title: BAMIRAC APPROXIMATE METHODS FOR CALCULATING THE II. Adamson, T. C., Jr. APPROXIMATE METHODS FOR CALCULATING THE II. Adamson, T. C., Jr. STRUCTURE OF JETS FROM HIGHLY UNDER EX- III. Advanced Research Proj- STRUCTURE OF JETS FROM HIGHLY UNDER EX- III. Advanced Research ProjPANDED NOZZLES by T. C. Adamson, Jr. Rept. of ects Agency PANDED NOZZLES by T. C. Adamson, Jr. Rept. of ects Agency BAMIRAC. Jun 61. 22 p. incl. illus., 8 refs. IV. Air Force Cambridge Re- BAMIRAC. Jun 61. 22 p. incl. illus., 8 refs. IV. Air Force Cambridge Re(Rept. no. 3768-26-T) search Laboratories, (Rept. no. 3768-26-T) search Laboratories, (Contract AF 19(604)- 7350) Unclassified report Bedford, Mass. (Contract AF 19(604)- 7350) Unclassified report Bedford, Mass. Approximate expressions are presented for the calcu- V. Contract AF 19(604)-7350 Approximate expressions are presented for the calcu- V- Contract AF 19(604)-7350 lation of the jet-boundary are intercepting-shock lo- lation of the jet-boundary are intercepting-shock locations of a jet flow from a highly underexpanded rocket cations of a jet flow from a highly underexpanded rocket nozzle. The gas is assumed to be an inviscid nonre- nozzle. The gas is assumed to be an inviscid nonreacting perfect gas with constant specific heats, ex- acting perfect gas with constant specific heats, exhausting to an ambient atmosphere. One new method hausting to an ambient atmosphere. One new method and two previously known methods of approximating and two previously known methods of approximating the jet boundary are discussed, and numerical ex- the jet boundary are discussed, and numerical examples of each are compared with calculations based on amples of each are compared with calculations based on (over) Armed Services (over) Armed Services Technical Information Agency Technical Information Agency UNCLASSIFIED UNCLASSIFIED AD Div. 12/1 UNCLASSIFIED AD Div. 12/1 UNCLASSIFIED Inst. of Science and Technology, U. of Mich.,Ann Arbor I. Title: BAMIRAC Inst. of Science and Technology, U. of Mich.,Ann Arbor I. Title: BAMIRAC APPROXIMATE METHODS FOR CALCULATING THE II. Adamson, T. C., Jr. APPROXIMATE METHODS FOR CALCULATING THE II. Adamson, T. C., Jr. STRUCTURE OF JETS FROM HIGHLY UNDER EX- III. Advanced Research Proj- STRUCTURE OF JETS FROM HIGHLY UNDER EX- III. Advanced Research ProjPANDED NOZZLES by T. C. Adamson, Jr. Rept. of ects Agency PANDED NOZZLES by T. C. Adamson, Jr. Rept. of ects Agency BAMIRAC. Jun 61. 22 p. incl. illus., 8 refs. IV. Air Force Cambridge Re- BAMIRAC. Jun 61. 22 p. incl. illus., 8 refs. IV. Air Force Cambridge Re(Rept. no. 3768-26-T) search Laboratories, (Rept. no. 3768-26-T) search Laboratories, (Contract AF 19(604)- 7350) Unclassified report Bedford, Mass. (Contract AF 19(604)- 7350) Unclassified report Bedford, Mass. Approximate expressions are presented for the calcu- V. Contract AF 19(604)-7350 Approximate expressions are presented for the calculation of the jet-boundary are intercepting-shock lo- lation of the jet-boundary are intercepting-shock locations of a jet flow from a highly underexpanded rocket cations of a jet flow from a highly underexpanded rocket nozzle. The gas is assumed to be an inviscid nonre- nozzle. The gas is assumed to be an inviscid nonreacting perfect gas with constant specific heats, ex- acting perfect gas with constant specific heats, exhausting to an ambient atmosphere. One new method hausting to an ambient atmosphere. One new method and two previously known methods of approximating and two previously known methods of approximating the jet boundary are discussed, and numerical ex- the jet boundary are discussed, and numerical examples of each are compared with calculations based on amples of each are compared with calculations based on (over) Armed Services (over) Armed Services Technical Information Agency Technical Information Agency UNCLASSIFIED UNCLASSIFIED + + +

AD UNCLASSIFIED AD UNCLASSIFIED the method of characteristics and with experimental re- DESCRIPTORS the method of characteristics and with experimental suits. A new method, for approximating the intercepting- Guided missiles suits. A new method for approximating the interceptingshock position is presented and compared with theoret- Guidnfrared mishock sition is presth theor ical and experimental results. Although the accuracy Infrared ical and experimental results. Although the accuracy varies, depending on the approximation employed, it is Nozzlesd, iNozzles shown that with any of the methods presented, the initialshown that with any of the methods presented, the initial part of the jet boundary and intercepting shock is well part of the jet boundary and intercepting shock is well represented by the approximate solutions. The location represented by the approximate solutions. The location of the first Mach disc is discussed briefly. of the first Mach disc is discussed briefly. UNCLASSIFIED UNCLASSIFIED AD UNCLASSIFIED AD UNCLASSIFIED the method of characteristics and with experimental re- DESCRIPTORS the method of characteristics and with experimental re- DESCRIPTORS suits. A new method for approximating the intercepting- Guided missiles suits. A new method for approximating the interceptingshock position is presented and compared with theoret- Infrared shock position is presented and compared with theoretical and experimental results. Although the accuracy Nozzles ical and experimental results. Although the accuracy Infrared varies, depending on the approximation employed, it is varies, depending on the approximation employed, it is Nozzles shown that with any of the methods presented, the initial shown that with any of the methods presented, the initial part of the jet boundary and intercepting shock is well part of the jet boundary and intercepting shock is well represented by the approximate solutions. The location represented by the approximate solutions. The location of the first Mach disc is discussed briefly. of the first Mach disc is discussed briefly. UNCLASSIFIED UNCLASSIFIED

+ + + AD Div. 12/1 UNCLASSIFIED AD Div. 12/1 UNCLASSIFIED Inst. of Science and Technology, U. of Mich., Ann Arbor I. Title: BAMIRAC Inst. of Science and Technology, U. of Mich., Ann Arbor I. Title: BAMIRAC APPROXIMATE METHODS FOR CALCULATING THE II. Adamson, T. C., Jr. APPROXIMATE METHODS FOR CALCULATING THE II. Adamson, T. C., Jr. STRUCTURE OF JETS FROM HIGHLY UNDER EX- III. Advanced Research Proj- STRUCTURE OF JETS FROM HIGHLY UNDER EX- III. Advanced Research ProjPANDED NOZZLES by T. C. Adamson, Jr. Rept. of ects Agency PANDED NOZZLES by T. C. Adamson, Jr. Rept. of ects Agency BAMIRAC. Jun 61. 22 p. incl. illus., 8 refs. IV. Air Force Cambridge Re- BAMIRAC. Jun 61. 22 p. incl. illus., 8 refs. IV. Air Force Cambridge Re(Rept. no. 3768-26-T) search Laboratories, (Rept. no. 3768-26-T) search Laboratories, (Contract AF 19(604)- 7350) Unclassified report Bedford, Mass. (Contract AF 19(604)- 7350) Unclassified report Bedford, Mass. V. Contract A F 19(604)-73 50,,,.,.Cnrc 964-35 Approximate expressions are presented for the calcu- V Approximate expressions are presented for the calculation of the jet-boundary are intercepting-shock lo- lation of the jet-boundary are intercepting-shock locations of a jet flow from a highly underexpanded rocket cations of a jet flow from a highly underexpanded rocket nozzle. The gas is assumed to be an inviscid nonre- nozzle. The gas is assumed to be an inviscid nonreacting perfect gas with constant specific heats, ex- acting perfect gas with constant specific heats, exhausting to an ambient atmosphere. One new method hausting to an ambient atmosphere. One new method and two previously known methods of approximating and two previously known methods of approximating the jet boundary are discussed, and numerical ex- the jet boundary are discussed, and numerical examples of each are compared with calculations based on amples of each are compared with calculations based on (over) Armed Services (over) Armed Services Technical Information Agency Technical Information Agency UNCLASSIFIED I UNCLASSIFIED AD Div. 12/1 UNCLASSIFIED AD Div. 12/1 UNCLASSIFIED Inst. of Science and Technology, U. of Mich., Ann Arbor I. Title: BAMIRAC Inst. of Science and Technology, U. of Mich., Ann Arbor I Title: BAMIRAC APPROXIMATE METHODS FOR CALCULATING THE II. Adamson, T. C., Jr. APPROXIMATE METHODS FOR CALCULATING THE II. Adamson, T. C., Jr. STRUCTURE OF JETS FROM HIGHLY UNDER EX- III. Advanced Research Proj- STRUCTURE OF JETS FROM HIGHLY UNDER EX- III. Advanced Research ProjPANDED NOZZLES by T. C. Adamson, Jr. Rept. of ects Agency PANDED NOZZLES by T. C. Adamson, Jr. Rept. of ects Agency BAMIRAC. Jun 61. 22 p. incl. illus., 8 refs. IV. Air Force Cambridge Re- BAMIRAC. Jun 61. 22 p. incl. illus., 8 refs. IV. Air Force Cambridge Re(Rept. no. 3768-26-T) search Laboratories, (Rept. no. 3768-26-T) search Laboratories, (Contract AF 19(604)- 7350) Unclassified report Bedford, Mass. (Contract AF 19(604)- 7350) Unclassified report Bedford, Mass. Approximate expressions are presented for the calcu- V. Contract AF 19(604)-7350 Approximate expressions are presented for the calculation of the jet-boundary are intercepting-shock lo- lation of the jet-boundary are intercepting-shock locations of a jet flow from a highly underexpanded rocket cations of a jet flow from a highly underexpanded rocket nozzle. The gas is assumed to be an inviscid nonre- nozzle. The gas is assumed to be an inviscid nonreacting perfect gas with constant specific heats, ex- acting perfect gas with constant specific heats, exhausting to an ambient atmosphere. One new method hausting to an ambient atmosphere. One new method and two previously known methods of approximating and two previously known methods of approximating the jet boundary are discussed, and numerical ex- the jet boundary are discussed, and numerical examples of each are compared with calculations based on amples of each are compared with calculations based on (over) Armed Services (over) Armed Services Technical Information Agency Technical Information Agency UNCLASSIFIED UNCLASSIFIED ++ +

AD UNCLASSIFIED AD UNCLASSIFIED the method of characteristics and with experimental re- DESCRIPTORS the method of characteristics and with experimental re- DESCRIPTORS suits. A new method for approximating the intercepting- Guided missiles suits. A new method for approximating the intercepting- Guided missiles shock position is presented and compared with theoret- shock position is presented and compared with theoretical and experimental results. Although the accuracyNozzles ical and experimental results. Although the accuracy varies, depending on the approximation employed, it is varies, depending on the approximation employed, it is shown that with any of the methods presented, the initial shown that with any of the methods presented, the initial part of the jet boundary and intercepting shock is well part of the jet boundary and intercepting shock is well represented by the approximate solutions. The location represented by the approximate solutions. The location of the first Mach disc is discussed briefly. of the first Mach disc is discussed briefly. UNCLASSIFIED UNCLASSIFIED c AD UNCLASSIFIED AD UNCLASSIFIED ^ ) the method of characteristics and with experimental re- DESCRIPTORS the method of characteristics and with experimental re- DESCRIPTORS suits. A new method for approximating the intercepting- Guided missiles suits. A new method for approximating the intercepting- Guided missiles shock position is presented and compared with theoret- shock position is presented and compared with theoret- I ical and experimental results. Although the accuracy Nozzlrare ical and experimental results. Although the accuracy varies, depending on the approximation employed, it is varies, depending on the approximation employed, it is shown that with any of the methods presented, the initial shown that with any of the methods presented, the initial part of the jet boundary and intercepting shock is well part of the jet boundary and intercepting shock is well represented by the approximate solutions. The location represented by the approximate solutions. The location of the first Mach disc is discussed briefly. of the first Mach disc is discussed briefly. UNCLASSIFIED UNCLASSIFIED