Laminar-to-turbulent transitions over an ablating reentry capsule

www.elsevier.com/locate/actaastro

PII:

Acta Astronautica Vol. 47, No. 10, pp. 745–751, 2000 ? 2000 Elsevier Science Ltd. All rights reserved Printed in Great Britain S0094-5765(00)00117-X 0094-5765/00/$ - see front matter

LAMINAR-TO-TURBULENT TRANSITIONS OVER AN ABLATING REENTRY CAPSULE KIMIYA KOMURASAKI Department of Aeronautics and Astronautics, University of Tokyo, Tokyo 113-8656, Japan and

GRAHAM V. CANDLER Department of Aerospace Engineering and Mechanics, University of Minnesota, Minneapolis, MN 55455, U.S.A (Received 27 April 1999; revised version received 23 November 1999)

Abstract—The aerodynamic mechanism of early transition phenomena over an ablating reentry capsule has been analytically examined. A two-equation turbulence model (k– model) was coupled with Reynolds averaged Navier–Stokes equations. Low-Reynolds-number eects on the solid wall were taken into account by modifying the Chien’s correction. As a result, transition occurred at the lower Reynolds number with higher ablation rate. The predicted transition-point Reynolds number was 3 × 104 at the surface-mass-injection rate of 100 g=sm2 . The principal mechanism of this early transition is thought as follows; the viscosity damping eects are reduced and re-laminarization is prevented in the downstream of the capsule surface, due to the turbulence on the surface and due to the pushing out of near-surface stream-lines from the surface by successive mass injection. ? 2000 Elsevier Science Ltd. All rights reserved

1. INTRODUCTION

A sample return mission MUSES-C has been proposed by the Institute of Space and Astronautical Sciences of Japan. A spacecraft will rendezvous with a near-earth asteroid Nereus and bring some samples of the asteroid material back to the Earth by a reentry capsule [1]. The spacecraft is scheduled to be launched in 2002 and to return in 2006. The atmospheric reentry speed of the capsule will exceed 12 km=s, and the corresponding reentry Mach number is 42. The laminar shock layer analyses conducted on this mission indicate that the maximum convective heating rate on the capsule surface would be approximately 10 MW=m2 and the radiative rate about 2 MW=m2 . In order to survive such an extremely severe heating, the reentry capsule will be shielded with a carbon-phenolic ablator. From a material test of the carbon-phenolic ablator using an arc-heater, the recession rate has been estimated at 100 g=sm2 at the heating rate of 10 MW=m2 [2]. This mass ux corresponds to several % of the free-stream mass

ux. Such a large amount of mass injection from the ablator might induce turbulence over the surface even near the stagnation point of the capsule, resulting in a higher heating rate than predicted by laminar ow analyses [3–5].

Several experimental studies on the laminarto-turbulent transitions over a body with surface mass injection have been done in the United States so far, and the results have shown transitions [6,7]. Generally speaking, boundary layer transitions hardly occur at a low Re of the order of 104 , which is a typical MUSES-C reentry condition. However, Kaattari’s experiment predicted the transition at Re=3 × 104 on a blunt body with surface-mass-injection. It implies that transitions might occur during the MUSES-C reentry. On the ablator surface, the gas injection would be inhomogeneous and unsteady. This surface condition could initiate turbulence in a ow, and eventually cause such a low-Re transition or an “early” transition [3]. The objective of this study is to analytically predict the transition Re and the resulting heating rate distribution over the MUSES-C reentry capsule. 2. NUMERICAL MODEL†

2.1. Reynolds averaged Navier–Stokes equation

Reynolds averaged Navier–Stokes equations coupled with the standard k– turbulent model [8]

745

†Nomencluture is given in Appendix at end of paper.

746

K. Komurasaki and G. V. Candler

in a cylindrical polar coordinate are expressed as @(FI − FV ) 1 @r(GI − GV ) @U + + = W; (1) @t @x r @r where



   u        U =  E ;    k  



 u  uu + p     u   FI =   (E + p)u  ;    ku  u



0 T   +  xx xx − 2=3k  T   xr + xr Fv = T (xx + −2=3k)u + (xr +T )v+ @x T xx xr   (M + T ) @x k (M + T =1:3) @x 



M = 1:462×10−6 T 1=2 =(1 + 112=T ) kg=ms: (7)

   ;   

Recent direct numerical simulations of turbulent reacting ows show that the chemical reactions affect the transition process by either damping or amplifying the disturbances [10]. However, this eect has not been modeled yet. 2.3. Low-Reynolds-number eects

 0   0     (p −  )=r ;  W =  0   R   ( ∇)u −  R (=k)(1:44( ∇)u − 1:92 

E = p=( − 1) + 0:5(u2 + v2 ) + k:

(2)

2.2. Transport coecients

In the rst-order turbulent closure model, the Reynolds shear stress R is expressed through the eddy viscosity, following Boussinesq’s assumption: Rij = Tij − (2=3)kij ;

(3)

Tij

is the stress component proportional to where the eddy viscosity as Tij = T (@i uj + @j ui − (2=3)∇ · uij ) = (T =M )ij : (4) In the standard k– model, the eddy viscosity is expressed as T = 0:09k 2 =:

(5)

Thermal conductivity coecient can be de ned through the Prandtl numbers. Assuming the laminar and turbulent Prandtl numbers for air ows constant at 0.72 and 0.9, respectively, total thermal conductivity is expressed as  = Cp (M =0:72 + T =0:9):

In this study, neither chemical reactions nor real gas eects are taken into account, because the boundary transition phenomena is strongly governed by the boundary layer con guration and would not be aected by the location of bow shock nor by the chemical composition. Because of this simpli cation, calculated post-shock temperature goes up to 80; 000 K, which is several times higher than actual post-shock temperature. The molecular viscosity coecient is obtained using an approximation for non-reacting air with a chemical composition frozen at standard conditions [9].

(6)

In order to predict the boundary ow transition, the limiting behavior of the uctuating velocities approaching a solid boundary has to be considered. This eect is often expressed in following corrections: (1) “non-isotropic dissipation” of the turbulent energy, and (2) a “damping” eect on eddy viscosity, as proposed by Jones and Launder [11], Launder and Sharma [12], Chien [13] and others. If the body surface is smooth and non-ablative, the uctuating velocity satis es the no-slip boundary condition and also satis es conservation of mass. Therefore, from expanding the uctuating velocity in Taylor series near a solid boundary, the tangential and normal components u0 and v0 must behave as follows: u0 = A(s; t)y + O(y2 );

v0 = B(s; t)y2 + O(y3 ); (8)

respectively, as y → 0. The behavior of k and  is deduced from these asymptotic variations of uctuation velocities as k = u02 + v02 = A2 y2 + O(y3 ); q √  = (u0 2y + v0 2y ) = A2 + O(y):

(9)

k decreases rapidly as y → 0; while  remains nite at y = 0. This is called as “non-isotropic dissipation” of turbulent energy. To achieve this asymptotic consistency, the “wall” dissipation, 0 is added to the turbulent energy equation in most of the low-Reynolds-number corrections. In Chien’s correction, 0 is expressed as 0 = 2k=y2 :

(10)

Laminar-to-turbulent Transitions

747

Table 1. Free stream conditions

Fig. 1. Schematic of velocity uctuations on a charring ablator.

As for the “damping” eect on eddy viscosity, Eq. (5) is modi ed using a damping function, f as T = 0:09f k 2 =

(11)

so that T → 0 asy+ → 0 . This eect preserves boundary layer ows from transitions until Re ∼ 106 . In the Chien’s correction, f depends on a wall variable, y+ . f = 1 − exp(−0:0115y+ ):

(12)

However, in the case of a charring ablator such as a carbon-phenolic one, there can exist nite velocity uctuations even at y = 0, because the pyrolysis gas ejects at subsonic speed through the porous char remnant and allows disturbances to propagate from the exterior ow into the ablator as schematically shown in Fig. 1. Assuming the limiting case when surface uctuation velocities perfectly couple with the ones of exterior ow, the wall dissipation, 0 is neglected in this study. As for the viscosity damping eect, eddy does not vanish on the ablator surface when the injected gas is turbulent. In order to specify the eddy scale on the ablator surface, we de ne a new wall variable using a non-zero wall variable yw+ on the wall as y+ = y(u =) + yw+ :

(13)

In this study, two cases are discussed: (1) the eddy scale is larger than boundary layer thickness so that the damping eect is negligible: yw+ = ∞

(14)

and (2) the eddy scale is given by the Park’s wall mixing length model [3] as yw+ = (d=K)(u =);

(15)

d =  e vw :

(16)

76 km altitude

Density Temperature Mach number Reynolds number

3:5 × 10−5 kg=m3 195 K 42 6500

60 km altitude

Density Temperature Mach number Reynolds number

30:5 × 10−5 kg=m3 256 K 38 45,000

Here vw is the mean injection velocity on the surface. A time constant, e is an average time interval of uctuation of injected ow. A reasonable value of the time constant can be taken to be e =2 × 10−4 s for carbon–carbon composite [3]. Although some more damping functions are added to the dissipation rate equation in the conventional low-Reynolds-number corrections, they are neglected here because their eects are very small for turbulent-surface conditions. 2.4. Free-stream conditions

The free-stream Reynolds number is a key parameter for the ampli cation or suppression of eddy viscosity in the boundary layer. Reentry conditions at the ight altitude of 60 and 76 km were tested, corresponding to the free-stream Reynolds numbers of 45,000 and 6500, respectively. The conditions are listed in Table 1. (In the MUSES-C mission, peak heating is predicted to be at approximately 60 km of altitude.) The free stream is assumed turbulence free. 2.5. Ablator surface conditions

The surface is assumed isothermal, and the surface temperature is assumed at 4000 K, which would be the maximum temperature for a carbon-phenolic-type ablator. The local mass-injection rates are assumed in proportion to the local heating rates as ˙ w vw = q:

(17)

The injection rate can be changed by tuning the coecient . The gas is injected normal to the surface with a given turbulent kinetic energy. The turbulence intensity on the surface is speci ed to 0:03 − 0:3. The boundary condition of dissipation rate, W was given so that (0:09k=2 )W = (M )W , while W doesnot give any distinguishable eect on the results. The eect of surface roughness has not been considered in this study, though it would be important when the surface mass injection is small.

748

K. Komurasaki and G. V. Candler

Fig. 2. Plot of a 120 × 30 grid.

3. NUMERICAL SCHEME

Flux vector splitting (FVS) method is often used for calculating hypersonic blunt-body ow- elds, because of its robustness and suitability for use in implicit schemes [14], while it comes at a price of reduced accuracy due to numerical dissipation. The numerical dissipation is unfavorable especially to the turbulence calculation, because the turbulence calculation is essentially a calculation of turbulence-dissipation processes. On the other hand, ux dierence splitting (FDS) method is very accurate, though operation count and complexity are increased for complete linearization of ux formulas for implicit schemes. Thereby, the advection upstream splitting method (AUSM) [15] is employed in this study. This scheme combines the eciency of FVS and the accuracy of FDS, and enables us to eciently capture the hypersonic shock and to accurately evaluate boundary layer ows. The Gauss–Seidel line relaxation method [16] is adopted to implicitly solve the nonlinear equations. Figure 2 is a plot of the grid used in this study. The master equations were discretized with a third-order upwind scheme using the MUSCL-type extrapolation. A grid convergence was obtained for heating rate distributions and for transition positions. A grid convergence has been con rmed: the grid points of 120 in the wall normal direction is found enough to have a convergence in heating rate distribution and the grid points of 30 in the wall trangential direction is enough to have a convergence in transition points.

Fig. 3. Heating rate distributions on an ablator surface. + = ∞; T 0 = 0:1. (a) Re = 45; 000, (b) Re = 6; 500. yW

4. RESULTS

4.1.

+ yW

= ∞ case

The calculated heating rate distributions over the capsule surface are shown in Fig. 3 . In the case of Re = 45; 000, the stagnation-point heating is decreased with the increase in surface injection rate due to the heat-blockage eect. Every plot shows a peak heating in the downstream region due to the laminar-to-turbulent transition, except for no-injection case. At the injection rate of 80 g=sm2 , the heating rate pro le peaks at s=12 cm, and the maximum rate was approximately 25% higher than the stagnation-point one. In the case of Re = 6500, the ow stayed laminar for any injection rate. Figure 4 shows the pro les of mean velocity and turbulence kinetic energy at Re =45; 000. The plots are enlarged in the surface normal direction. Turbulence kinetic energy is ampli ed in the stream-wise direction in the middle of boundary layer. Figure 5 shows the distributions of the eddy viscosity and molecular viscosity. These plots are also enlarged in the surface normal direction. The eddy viscosity surpassed the molecular viscosity in the downstream region. The peak eddy viscosity along wall normal lines is plotted in the stream-wise direction in Fig. 6. It was normalized by the local molecular

Laminar-to-turbulent Transitions

749

viscosity. The nonlinear ampli cation of the viscosity indicates the transition. The transition point shifts upstream with an increase in the surface-mass-injection rate. + 4.2. yW = (d=K)(u =) case

Fig. 4. Pro les of mean velocity and turbulent + kinetic energy. Re = 45; 000, yW =∞, T 0 = 0:1, 2 (W vW )SP =80 g=s m . Enlarged in the wall normal direction by a factor of 5.

Figure 7 shows the distributions of heating rate and peak eddy viscosity along the surface contour + = (d=K)(u =). Comparing with in the case of yW + the case of yW =∞, the transition is delayed due to the viscosity damping eect. The ow stayed laminar at the injection rate less than 80 g=sm2 . Figure 8 shows the heating rate distributions for various injection rates. The transition occurs earlier with larger injection rate. The distance from the stagnation point to the transition point is plotted for the cases with various surface turbulence intensities in Fig. 9. The gure shows a weak dependence on the surface turbulence intensity. This is due to the fact that the kinetic energy of the injection gas is much smaller than that of main ows and the turbulent energy in the vicinity of surface is determined through the diusion of turbulent energy from the main ow. Owing to the same reason, the calculation is insen-

Fig. 5. Pro les of eddy viscosity and molecular viscosity. + Re = 45; 000, yW =∞, T 0 = 0:1, (W vW )SP =80 g=sm2 . Enlarged in the wall normal direction by a factor of 5. Fig. 7. Distributions of heating rate and eddy viscosity. Re = 45; 000, (W vW )SP =80 g=sm2 , T 0 = 0:1.

+ Fig. 6. Eddy viscosity distribution. Re=45; 000, yW =∞, T 0 = 0:1.

Fig. 8. Heating rate distributions on an ablator surface. + = (d=K)(u =), T 0 = 0:1. Re = 45; 000, yW

750

K. Komurasaki and G. V. Candler

Fig. 9. Surface turbulence intensity and transition distance from stagnation point. Re = 45; 000, + yW = (d=K)(u =).

Fig. 11. Transition point Re and injection rate ratio.

into the middle of boundary layer by successive mass injection in the downstream region of the surface. At a small injection rate, the streamline stays in the region where the viscosity damping eect is strong, while at a large injection rate, the streamline is pushed into the region where local Re is relatively high and damping eect is weak, resulting in a transition. 4.3. Injection rate and transition point Re

Figure 11 shows a comparison between the Kaattari’s transition experiment and the present prediction. The injection rate normalized by the free-stream mass ux is on the abscissa and the transition point Re is on the ordinate. + =∞, calculated results show In the case of yW earlier transition than the experimental observations. Transitions occur even at the quite small injection rate. That seems quite unrealistic. On the + = (d=K)(u =), tranother hand, in the case of yW sitions do not occur at the injection rate less than 80 g=sm2 . The calculated transition envelope is well agreed with the experimental results. Fig. 10. A streamline of the gas injected near stagnation point and contours of (T =M ). Re = 45; 000, + yW = (d=K)(u =), T 0 = 0:1, (W W )SP =80 g=sm2 . Enlarged in the wall normal direction by a factor of 20.

sitive to the dissipation rate condition on the surface, W . If T 0 = 0 on the surface, the ow stays laminar at any injection rate because there is no seed of turbulence. Figure 10 shows a streamline of the gas injected near the stagnation point. The seed of turbulence given on the surface is carried along the streamline and is gradually ampli ed by the work of the main ow against the Reynolds stress. The streamline is separated from the surface and is swelled

5. DISCUSSIONS

The mechanism for the early transition induced by surface-mass-injection would consist in following two points: (1) reduction of the viscosity damping eect due to the surface turbulent condition, and (2) pushing out of the injected-gas stream-lines toward the middle of boundary layer by successive mass injection in the downstream region of the body. Both of them assist the ampli cation of eddy viscosity in the boundary ows. Since the predicted transition point Re shows a good agreement with the Kaattari’s experimental results, it can be concluded that this numerical model well represents the early transition mechanism due to surface-mass-injection.

Laminar-to-turbulent Transitions

In the MUSES-C reentry, the increase in heating rate due to the transition would be small so far as the ablation rate is less than 100 g/sm2 . REFERENCES

1. Park, C. and Abe, T. Research on the heatshield for MUSES-C Earth Reentry. 7th Joint Thermophysics and Heat Transfer Conference, AIAA Paper 98-2852, Albuquerque NW, 1998. 2. Ishii, N., Inatani, Y., Hiraki, K. and Yamada, T., Thermal protection system of Muses-C Reentry Capsule. 21st International Symposium on Space Technology and Science, Paper 98-d-04, Omiya, Japan, 1998. 3. Park, C., Injection-induced turbulence in stagnation-point boundary layer. AIAA Journal, 1984, 22, 219 –225. 4. Park, C. and Balakrishnan, A., Ablation of Galileo probe heat-shield models in a ballistic range. AIAA Journal, 1985, 23, 301–308 5. Laurien, E., Numerical investigation of laminarturbulent transition on reentry capsules. Journal of Spacecraft Rockets, 1996, 33, 313–318. 6. Demetriades, A., Laderman, A. J. and Von Seggern, L., Eect of mass addition on the boundary layer of a hemisphere at Mach 6. Journal of Spacecraft and Rockets, 1976, 13,508–509 7. Kaattari, G. E., Eects of mass addition on blunt-body boundary-layer transition and heat transfer. NASA Technical Paper 1139, 1979. 8. Launder, B. E., Spalding, B., Mathematical Models of Turbulence. Academic Press, New York, 1972. 9. Hansen, C. F., Approximation for the thermodynamic and transport properties of high temperature air. NASA TR-R-50, 1959. 10. Martin, M. P. and Candler, G. V., Eect of chemical reactions on decaying isotropic turbulence. Physics of Fluids, 1998, 10, 1715 –1724. 11. Jones, W. P., Launder, B. E., The prediction of laminalization with a two-equation model of turbulence. International Journal of Heat Mass Transfer, 1972, 15, 301–314. 12. Launder, B. E. and Sharma, B. I., Application of the energy dissipation model of turbulence to the calculation of ow near a spinning disc. Letters in Heat and Mass Transfer, 1974, 1, 131–138. 13. Chien, K-Y., Predictions of channel and boundary layer ows with a low-reynolds number turbulence model. AIAA Journal, 1982, 20, 33–38. 14. Candler, G. V. and MacCormack, R. W., The computation of hypersonic ionized ows in chemical and thermal nonequilibrium. Journal of Thermophysics and Heat Transfer, 1991, 5, 266 –273.

751

15. Liou, M. S. and Steen, C. J., A new ux splitting scheme. Journal of Computational Physics, 1993, 107, 23–39. 16. MacCormack, R.W., Current status of the numerical solutions of the Navier–Stokes Equations. AIAA Paper 95-0032, Reno NV, 1998.

APPENDIX Nomenclature Cp d E F; G k K p Pr q˙ R Re Retr s t T T0 u;  u 0 ; 0 u U W x; r;  y y+ 

 M T    R

Speci c heat at constant pressure Park’s mixing length on the turbulent wall Total energy per unit mass Flux vectors in x and r directions, respectively Turbulence kinetic energy per unit mass Von Karman constant, 0.41 Static pressure Prandtl number Heating rate on the capsule surface Nose radius, 0:2 m Reynolds number=∞ V∞ R=∞ Transition point=∞ V∞ Str =∞ Distance from stagnation point along the body contour Time Temperature p Turbulence intensity= 2k=3vw3 Mean velocity components Fluctuating velocity pcomponents Friction velocity= w = Vector of conservative variables Source vector Cylindrical polar coordinates Wall normal distance Wall variable, y+ = y(u =) Turbulence dissipation rate Speci c heat ratio Thermal conductivity Molecular viscosity Eddy viscosity Kinematic molecular viscosity Density Shear stress Reynolds shear stress

Subscripts I SP T tr V W ∞

Inviscid Stagnation point Turbulent Transition Viscous Wall condition Free-stream condition