1、- R, , *- , W NASA TECHNICAL NOTE NASA TN D-4563 e-/ taI z c 4 w 4 z CALCULATION OF AXISYMMETRIC SUPERSONIC FLOW PAST BLUNT BODIES WITH SONIC CORNERS, INCLUDING A PROGRAM DESCRIPTION AND LISTING by Jerry C. South, Jr. Lungley Resemch Center LungZey Stution, Hdmpton, Vu. . . 8c,det, the calculations
2、for small bluntness ratios are good approximations to the experimental results for both the shock shape and pressures. On the contrary, when 0, 5 Oc,det, the limit Rn/Rb - 0 does involve discontinu ous alteration of an important boundary condition; that is, a pointed cone (R,/Rb = 0) has an oblique
3、shock wave attached at the vertex, and hence the surface entropy is less than that corresponding to the detached, normal shock conditions which hold for any finite (4Rb 0). Even for a pointed cone, the passage from an attached obliqueblunting R shock to a detached normal shock at Oc = oC,det implies
4、 a discontinuous jump in the surface entropy, and hence a corresponding jump in both the vertex and sonic corner pressures. In spite of this, it seems that the drag of a pointed, finite cone is a continuous function of the cone angle while passing through Oc = oc,det. Busemann (ref. 14) argued lSome
5、 local nonuniform behavior is anticipated near the stagnation point as Rn/Rb -c 0, since the apex of the pointed cone (with detached shock) is a singular point. The behavior of the flow near the apex is analogous to that of a wedge in incompressible potential flow, where fluid properties vary like s
6、IJ, 0 8c,son, and increases monotonically and continuously through detachment up to Oc = 90. This offers a heu ristic picture, then, of the dependence of CD on 8c for a pointed finite cone. It remains to consider the effects of small, but finite, bluntness. A rigorous discussion of the interaction b
7、etween small blunting and the transonic singularities that arise in the pointed case for , 0.05). For Bc = 51, Rn/Rb = 0.25, forty-five 21n an exact analysis, the chief concern would be the point of origin of the limiting characteristic (ref. 1, pp. 202-203), downstream from which changes in the bod
8、y shape cannot affect the nose solution. In the one-strip approximation, however, it is the surface sonic point that is crucial. 13 Provided by IHSNot for ResaleNo reproduction or networking permitted without license from IHS-,-,-halving-mode iterations (producing about fourteen correct figures3 in
9、S(O)/Rn) were required for an accuracy criterion E = In a study of blunted cones of “small7*angle (i.e., with sonic point always on the spherical nose), Traugott (ref. 19) described a range of cone angles in which the integral method should fail, for a given y and M,. The upper limit is the critical
10、 angle just described; that is, the cone angle Bc for which the sonic point of the complete smooth sphere would occur at the sphere-cone junction. The lower limit is the cone angle for which the asymptotic, pointed-cone pressure is equal to the sonic pressure obtained by isentropic expansion from th
11、e normal shock conditions. The cone angles in between con stitute what Traugott called the “second sonic point“ region; wherein the flow along the surface will become supersonic on the spherical nose, but must become subsonic again if asymptotic conical pressure is to be reached downstream. While Tr
12、augott found the-method to fail upon approaching this range of angles from below, the present program fails upon approaching from above. Unfortunately, there are not available any extensive tabulations of the sonic-point location on a sphere for various combinations of y and M,. In fact, it is that
13、informa tion for the one-strip integral method that is needed here, for the purpose of avoiding the critical condition when using this program. A more practical limit, which can be charted profusely, is the sonic condition for the pointed cone, BC,SOn(y,M,). It was found that this condition Lies fai
14、rly close to the “smooth sphere sonic-point“ condition and will serve as a guide in warning the program user of expected sensitivity and possible nonconver gence. Figure 21 is a chart of the pointed sonic cone condition; the calculations were performed by using the algebraic equations which approxim
15、ate conical flow by setting the right-hand sides of equations (2) and (3) equal to zero (ref. 17). Convex or Concave Spherical Caps The same difficulties are to be expected for the spherical cap (convex to the stream) when the corner is near the smooth sphere sonic point. For the concave spherical c
16、ap, no extensive investigation was undertaken to determine the limits of applicability of the present program. An obvious geometrical condition to be avoided for this shape is the intersection of the surface normals within the shock layer - that is, since 3A reliable estimate for the number of halvi
17、ng-mode iterations Nhm required to produce Nd decimal figures in S(O)/Rn is 14 Provided by IHSNot for ResaleNo reproduction or networking permitted without license from IHS-,-,-K = -dQ/ds = -1 for the concave spherical cap, the condition K6 9 1 might occur. It is possible that this could happen for
18、low values of M, and/or values of Rn/Rb approaching 1.0 (from above). M, Too Small As M, decreases, the shock wave moves farther away from the body and the sonic point on the shock wave moves up, away from the axis. For a given body shape, a value of M, 1 will be reached where a considerable portion
19、 of the subsonic flow region lies beyond the normal drawn from the sonic corner of the body. It would seem that such solutions would be either nonexistent or meaningless in the one-strip integral approximation, since the calculation ends when the sonic corner is reached, and hence a part of the subs
20、onic flow would be ignored. Such a situation is, of course, both physically and mathematically objectionable. What appears to happen in the one-strip approximation, at least in all cases studied so far, is the following: the shock-wave sonic point always occurs ahead of the surface normal drawn from
21、 the sonic corner of the body, which implies proper closure of the subsonic region. But since the shock wave has less curvature near the axis as M, decreases, it must turn more rapidly as the sonic corner is approached to achieve the sonic shock angle; the factor a (plulvl )Iap in the denominator of
22、 dp/ds becomes small near the corner, resulting in the rapid increase in magnitude of dp/ds. Ultimately a value of M, is reached where a (plul ap passes through zero at the corner, and hence the shock has infinite curvature there. The disk provides a good example of this apdifficulty. In table I are
23、 given values of a (plulvl ,*/ as M, decreases from hyper sonic values; it is clear that two-fold smooth (continuous curvature) shock-wave solutions do not exist below M, 2.3. TABLE 1.- APPROACH TO SHOCK CURVATURE SINGULARITY FOR A DISK M, s(o)/Rb 10.0 0.470 6.0 .494 4.0 .541 3.0 .607 2.4 .703 2.3 .
24、732 a (p1u1v1)*/ 2.004 1.812 1.464 .997 .284 .024 15 Provided by IHSNot for ResaleNo reproduction or networking permitted without license from IHS-,-,-CONCLUDING =MARKS The present program has the capability for calculating good approximations to the shock-wave shape and surface pressure distributio
25、n for blunt bodies with sonic corners. Typical results for the disk, concave or convex spherical caps, and sphere-cones agree well with experimental data. The flow past a pointed cone was approximated by calculations for a sphere-cone with small bluntness ratio. It was found that results were genera
26、lly good for cone angles larger than the detachment angle; but results were not satisfactory as the cone angle was decreased to and below the detachment angle. The reason for this appears to be the inability of the method to account for a vanishingly thin entropy layer. For the blunted cone, inflect
27、ion points occur in the shock wave when the sonic con dition for a pointed cone is reached. This condition can be reached, for example, by decreasing the cone angle. Further decreases in cone angle cause the nose and afterbody solutions to become nearly independent of each other. It was found that r
28、easonable estimates could be made for angle-of-attack pressure distributions in the symmetry plane of a blunt cone. The windward and leeward distri butions were simulated by adding and subtracting the angle of attack from the cone angle, respectively . A critical condition for the method is evident
29、in the application to sphere-cones and spherical caps. It occurs when the combination of parameters (cone angle or surface inclination angle at the sonic corner, Mach number, and specific heat ratio) is such that the natural sonic-point location for the complete sphere lies at the sphere-cone juncti
30、on or the spherical cap corner. Depending on the bluntness ratio, extreme sensitivity and nonconvergence is encountered somewhat before this condition is reached. For a given body shape, there is a minimum Mach number below which no two-fold smooth (continuous curvature) shock solutions exist. The a
31、ccuracy of the solutions near this minimum value is questionable, particularly near the sonic corner. For the disk, the minimum Mach number was about 2.3 (in air). Langley Research Center, National Aeronautics and Space Administration, Langley Station, Hampton, Va. , February 13, 1968, 129-01-03-06-
32、23. 16 Provided by IHSNot for ResaleNo reproduction or networking permitted without license from IHS-,-,-APPENDIX A PROGRAM SYMBOLS The FORTRAN symbols appearing in the input, headings, and output of this program (appendix C) are given in the left-hand column, and their meanings defined in terms of
33、standard notation or the symbols defined in the section “Symbols“ are given in the right-hand column. These program symbols are as follows: Input: LC body shape trigger, 1, 2, 3, or 5 for disk, convex spherical cap, blunted cone, or concave spherical cap, respectively MM loop counter, number of case
34、s (combinations of y, M, OC, and R) for a given body shape (LC) GAMMA Y AMINF M, THETAC Bc, cone half-angle in degrees for blunt cone (LC = 3); dummy input value for other shapes (other LC) R sonic-point radius parameter, 1.0 for disk (LC = 1) and ro*IR, for other shapes (other LC) Case heading (oth
35、er than symbols defined for input): DEL0 6(0), initial standoff distance in appropriate length scale (as explained in appendix B) DS initial integration step size SI value of s at sphere-cone junction (LC = 3); dummy value for other shapes 17 Provided by IHSNot for ResaleNo reproduction or networkin
36、g permitted without license from IHS-,-,-APPENDIX A NIS number of integration steps in final run NIR number of integration runs to find correct 6(0) CD forebody pressure drag coefficient (referred to base area m0*2 for all shapes) .- for each integration step):Output (the output is printed across in
37、 a block of three rows -S BETAD THETAD xo RO uo PO RHO0 TO DELTA x1 R1 v1 u1 P1 18 S p, degrees 8, degrees XO 0 UO PO PO temperature ratio, T/T 6 X1 1 V1 U1 P1 Provided by IHSNot for ResaleNo reproduction or networking permitted without license from IHS-,-,-APPENDIX A DS integration step size DUODS
38、duo/dS POBAR SBAR s/s* Main program (some of the more important FORTRAN symbols used in the main program): DSC starting step size for each integration run LN overall print trigger, 1for no print and 0 for print (switched to zero when correct 6(0) is found) N step print trigger, 0 for no print and 1f
39、or print (only print when stable step is completed) L body shape trigger, equal to LC but switched to 4 for rescaling when LC = 3 (see description of body shape subroutine in appendix B) LL trigger for double-valued curvature at sphere-cone junction, 1 for sphere curvature and 2 for zero curvature L
40、LL step size trigger, 1 until step size is halved to achieve velocity bracket and 2 thereafter (see program statement 54) LLLL halving mode trigger, 1 until upper and lower bounds for 6(0) are found and 2 thereafter Bodv shaDe subroutine : THETA 8, radians AK K, surface curvature 19 Provided by IHSN
41、ot for ResaleNo reproduction or networking permitted without license from IHS-,-,-APPENDIX B PROGRAM OPERATION See main text section e titled “Limitations of Program“ before operati g this program. General Remarks The FORTRAN IV program listed in appendix C was originally developed for the IBM 7094
42、electronic data processing system, and subsequently it was modified for use in the Control Data 6600 computer system. Other than using control cards which are appropriate to an individual system, the changes necessary to operate the program on the IBM 7094 system are as follows: (1) Remove the two i
43、nstructions in statement 1, IF (EOF,5) 9999,1000 9999 STOP (2) In statement 2, change the thirteenth instruction to: IF (NIR.LT.31) GO TO 2000 The reason for change (2) is the shorter single-precision word length of the IBM system (eight decimal digits). More than 27 halving iterations for 6(0) will
44、 not change the eight decimal digits carried in the calculations. About four runs are allotted for finding upper and lower bounds for 6(0). A large number of comment cards were included to highlight the various sections of the program. A fourth-order Runge-Kutta integration scheme was built into the
45、 pro gram so that a “library“ routine would be unnecessary. Also, a rule-of-thumb step size variation scheme is used which is very simple yet adequate for these calculations. It tests on the shock-angle derivative dp/ds and seeks a step size such that the “pre dicted“ and “corrected“ values for d/3/
46、ds at the interval midpoint are within a few per cent of each other. Numerical stability and accuracy are not a problem in these calculations. Initial estimates for DEL0.- For each of the four shapes included, the initial esti mate for 6(0) is automatically computed, as follows: 20 Provided by IHSNo
47、t for ResaleNo reproduction or networking permitted without license from IHS-,-,- - APPENDIX B -LC 1 2 3 5 Shape I DELO I Source Disk 1.03/im Ref. 21 Convex sphere cap .667/E1(0) - 13 Ref. 21 Sphere-cone .667R/p1(0) - 3 Concave sphere cap 1.03R/m Accuracy criterion.- The accuracy criterion for locat
48、ing the sonic corner is listed as that is, IROST/R - 11 5 (see the fourth instruction from the end of state ment 60). This can be relaxed, if necessary, to achieve convergence closer to the “critical condition“ described in the main text in the section “Limitations of Program.“ Modifications for oth
49、er body shapes.- The program can be modified to include other shapes in the body shape subroutine (BSR) by adding more statements, triggers, and, if there are more than two segments to the shape, more junction locations SII, SUI, etc. Care must be taken to insure that an integration step coincides with each junction loca tion, and that