- Email: [email protected]
A hybrid numerical scheme for aeroheating computation of hypersonic reentry vehicles
International Journal of Heat and Mass Transfer 116 (2018) 432–444
Contents lists available at ScienceDirect
International Journal of Heat and Mass Transfer journal homepage: www.elsevier.com/locate/ijhmt
A hybrid numerical scheme for aeroheating computation of hypersonic reentry vehicles Zhen-Xun Gao a,⇑, Hai-Chao Xue b, Zhi-Chao Zhang a, Hong-Peng Liu a, Chun-Hian Lee a a b
National Laboratory for Computational Fluid Dynamics, School of Aeronautic Science & Engineering, Beihang University, Beijing 100191, China China Academy of Aerospace Aerodynamics, Beijing 100074, China
a r t i c l e
i n f o
Article history: Received 24 August 2016 Received in revised form 20 July 2017 Accepted 21 July 2017 Available online 17 October 2017 Keywords: Aeroheating Numerical scheme Hypersonic Blunt body
a b s t r a c t Hypersonic aeroheating simulations were conducted for three-dimensional (3D) reentry capsule for the Mars Science Laboratory (MSL) using the AUSM and the Roe schemes. The study confirms that AUSM family of schemes generate unphysical oscillations in the wall heat flux distributions due to the unstable simulation of the shock wave. And the wall heat flux is underpredicted by the Roe scheme due to excess dissipation, even though captures the shock wave accurately. In order to overcome these limitations, we develop a numerical scheme that can capture the shock wave stably and still maintain an accurate prediction of the heat flux. We propose a novel hybrid scheme that combines both the AUSM+ and Roe schemes. Two different hybrid techniques are proposed for switching between the two schemes, one based on the shock wave detection and other utilizing a criterion based on the distance from the wall. The implementation of the two hybrid schemes for the 3D MSL reentry capsule was studied for prediction accuracy compared to experimental data. It is found that the second hybrid scheme (based on the criterion for switching based on distance from the wall) can eliminate the unphysical oscillations in the wall heat flux distributions and improve the aeroheating prediction accuracy. Ó 2017 Elsevier Ltd. All rights reserved.
1. Introduction Space technology has been developing rapidly in recent years, and the problem of thermal protection for the hypersonic reentry vehicles is becoming increasingly prominent. High-precision prediction for the aerothermal environment is one of the crucial technologies for the design of a thermal protection system. With the advent of numerical methods and computational capability, computational fluid dynamics (CFD) gradually becomes an important tool for the hypersonic aeroheating prediction. However, it is still very difficult to achieve a robust and accurate wall heat flux prediction for hypersonic reentry vehicles using CFD. Previous research revealed that the grid-point distribution and numerical schemes are the two important factors influencing the accuracy of the heat flux prediction in CFD [1–11] for a hypersonic laminar flow. As early as 1990s, a series of numerical experiments were conducted by Hoffmann et al. [1] based on the studies of Klopfer and Yee [2] and Blottner and Larson [3], and the results showed that the computed wall heat flux depends highly on the ⇑ Corresponding author. E-mail addresses: [email protected] (Z.-X. Gao), xuehc31201@hotmail. com (H.-C. Xue), [email protected] (Z.-C. Zhang), [email protected] (H.-P. Liu), [email protected] (C.-H. Lee). http://dx.doi.org/10.1016/j.ijheatmasstransfer.2017.07.100 0017-9310/Ó 2017 Elsevier Ltd. All rights reserved.
density of the grid points near the wall. Subsequently, some researchers studied in detail the effects of the first grid point spacing off the wall on the accuracy of the computed heat flux, resulting in several grid generation criterions proposed to estimate an appropriate first grid point spacing off the wall for an accurate heat flux computation [5,6]. However, for a high Mach number flow over a blunt body where a strong detached shock wave exists, the accuracy of the computed wall heat flux depends not only on the near-wall grid spacing, but also the grid quality in the vicinity of the shock wave. By numerical experiments, Henderson and Menart [7] showed that, in order to achieve a reasonable aeroheating prediction, the shape of the mesh needs to be in accordance to the shock wave and grid points should also be clustered normal to the shock wave. A low-quality grid-point distribution near the shock wave is easy to lead to an oscillating resolution of the shock wave which would cause an anomalous distribution of the wall heat flux. While for the numerical schemes, Kitamura et al. [8,9] argued that a numerical scheme suitable for hypersonic wall heat flux computation should be equipped with the following three properties: shock stability/robustness, conservation of total enthalpy and accurate resolution of boundary layer. In Ref. [9], 15 kinds of numerical schemes were evaluated according to the above three properties, and it was found that no numerical scheme investigated in [9] possessed all the three properties, while the
433
Z.-X. Gao et al. / International Journal of Heat and Mass Transfer 116 (2018) 432–444
Nomenclature Ma
q T p uj H k
l sij e E R
a
free-stream Mach number density (kg/m3) temperature (K) pressure (pa) velocity (m/s) total enthalpy (J/kg) molecular thermal conductivity (W/(m K)) molecular viscosity (N s/m2) shear stress (N/m2) internal energy (J/kg) total energy (J/kg) gas constant (J/(kg K)) the angle of attack (deg)
AUSM family schemes [12–16] and Roe [17] scheme appeared to be promising. However, these two kinds of numerical schemes have different characteristics for the aeroheating prediction. The AUSM family schemes do not have tunable factors in their scheme design and generally have low numerical dissipation, and thus they could resolve the near-wall temperature gradient more precisely within the boundary layer. Hence, with the same grid spacing near the wall, AUSM family schemes could give a more accurate heat flux result than Roe scheme [10,11]. However, for a flow with strong shock wave, if the shape of the mesh is not in accordance to the shock wave or the grid points are not dense enough across the shock wave, the AUSM family schemes may resolve the shock wave with obvious oscillations [10], which would further generate an unphysical distribution of heat flux on the wall. While for the Roe scheme, an entropy correction function with a tunable factor is incorporated, which could help to control the added numerical dissipation and guarantee a stable shock wave resolution even with a low-quality mesh near the shock wave, but would also bring in additional numerical dissipation inside the boundary layer and thus decrease the prediction accuracy of the wall heat flux. Generally, a small value of the tunable factor in the entropy correction function could lead to relative accurate heat flux predictions for the Roe scheme. However, for a three-dimensional (3D) complicated case with high Mach number, this tunable factor needs to be taken as a large value to ensure a stable resolution of the shock wave and then the numerical resolution inside the boundary layer by the Roe scheme would be significantly reduced. In practical hypersonic aeroheating computations, for a twodimensional (2D) problem, it is still feasible to keep the mesh aligned with the shape of the shock wave and then maintain a stable resolution of the shock wave, so the AUSM family schemes could be adopted to enhance the accuracy of the wall heat flux computation. While for 3D complex flows, it would be difficult to adjust the computational mesh to achieve an accordance of the grid shape and shock wave. Even if it could be realized, the mesh will have to be re-generated for different free-stream conditions because the shape of the shock wave would change, thereby increasing the work for the mesh generation. To remove the above restriction on the grid generation for stable resolution of the shock wave, the Roe scheme is a better choice than the AUSM family schemes. However, the high-dissipation characteristic of Roe scheme in boundary layer is likely to introduce errors in the wall heat flux prediction. Thus Roe scheme generally needs more dense grid points near the wall than the AUSM schemes in order to obtain reasonable heat flux, which would severely limit the computational time step and thus dramatically increase the computational consumption especially for 3D complicated cases. The above dilemma indicates that more efforts should be put on the develop-
local heat flux (W/m2) Reynolds number radius of cylinder (m) damping function amplification factor of damping function the distance of the local grid point from the wall (m) empirical parameter the characteristic length (m)
Qw Re R fD Cam d d0 L
Subscripts 1 free-stream w at the wall
ment of a reliable numerical scheme which could guarantee a stable resolution of strong shock wave without strict restriction on the grid generation while in the meantime maintain an accurate computation of the wall heat flux for hypersonic vehicles. The objective of the present paper is to develop a hybrid numerical scheme combining the Roe and AUSM schemes, which would adopt the Roe scheme around the shock wave to guarantee a stable resolution of the shock wave without strict restriction on the grid generation while utilize the AUSM scheme within the boundary layer to maintain low numerical dissipation and high accuracy of the wall heat flux computation. First, we introduce the governing equations and numerical methodology used in the present study. Secondly, hypersonic aeroheating computations for a 3D reentry capsule from the Mars Science Laboratory (MSL) are performed using the Roe and AUSM family schemes respectively, and the aeroheating accuracy is examined independently for each of the two schemes. Then, new hybrid schemes combining the Roe and AUSM family schemes are proposed following two different techniques, and the performance of the new hybrid schemes is evaluated. Finally, by applying the hybrid schemes to the 3D reentry capsule case, improvements in aeroheating prediction accuracy of the hybrid schemes are examined in detail. 2. Governing equations and numerical methodology The compressible Navier-Stokes (N-S) equations can be written as follow: @ðqÞ @t
þ @x@ j ðquj Þ ¼ 0 @s
@ðqui Þ @t
þ @x@ j ðqui uj þ dij pÞ ¼ @xijj @ðqEÞ @T þ @x@ ðqHuj Þ ¼ @x@ sij ui þ k @x @t j
j
ð1Þ
j
where sij is the shear stress given as:
sij ¼ l
@ui @uj 2 @ui l þ dij 3 @xj @xj @xi
ð2Þ
E is the total energy defined as:
1 E ¼ e þ ui ui 2
ð3Þ
The equation of state is
p ¼ qRT
ð4Þ
All the numerical simulations are carried out using a finite difference code, ACANS (Aerodynamic, Combustion and Aerothermodynamic Numerical Simulation), developed by the authors. The reliability of ACANS code has been validated by a variety of simu-
Z.-X. Gao et al. / International Journal of Heat and Mass Transfer 116 (2018) 432–444
lations [18–20]. In this study, the inviscid flux vectors are discretized by the Roe scheme [17] and AUSM family schemes including AUSMDV [12], AUSM+ [13], AUSMPW [14], AUSMPW+ [15], AUSM+UP [16]. The entropy correction function in the Roe scheme is the Yee type [21], and the value of its tunable factor is taken as 0.125. The viscous flux vectors are discretized by the central-difference scheme. The massively parallel computation is realized using the message passing interface (MPI) based on the space-decomposed algorithm. In this study, a laminar perfect-gas model is used. 3. Hypersonic aeroheating computation of a 3D reentry capsule In this section, a 3D reentry capsule from the MSL [22,23] is employed for performing the aeroheating numerical simulations applying Roe and AUSM family schemes, respectively. This reentry capsule was tested in AEDC (the US Air Force’s Arnold Engineering Development Center) Tunnel 9 and the experimental condition of Run 3021 in [22] is selected in the present study: Ma1 = 9.47, p1 = 167.9 Pa, T1 = 54.8 K, a = 0°. The wall boundary condition is assumed as isothermal with Tw = 300 K. Only the forebody of the MSL capsule is kept in the simulations, as shown in Fig. 1(a). Meanwhile, Fig. 1(b) illustrates the 3D structured mesh generated for the forebody of MSL capsule with a total 243,700 grid points. There are 100 grid points in the direction normal to the wall and the stretching factor near the wall is 1.1. In the streamwise and spanwise directions, the grid points are clustered around the nose and shoulder of the capsule with the stretching factor no more than 1.2. Moreover, the grid points are clustered near the shock wave region in the normal direction, but the shape of mesh is not aligned with the shock wave. Three meshes with different distances of the first grid point off the wall are adopted to perform a grid convergence study: (a) 1 105 m, (b) 5 106 m, (c) 1 106 m. Fig. 2 shows the comparisons among the wall heat flux by AUSM+ scheme, using the above three meshes. It can be seen that the distribution of heat flux obtained using grid (b) nearly coincides with that by grid (c). Hence, for the following computations of the MSL capsule, grid (b) with the distance of the first grid point off the wall 5 106 m is adopted. All the computations for the MSL capsule case are performed with CFL number of 10 for Roe scheme while 5 for AUSM family schemes. Figs. 3 and 4 show the convergence histories of the Roe and AUSM+ schemes for the MSL capsule case. Both the residual and the stagnation-point heat flux are used to examine the convergence of the computations. It is seen from Fig. 3 that, for Roe scheme, the residual drops more than 6 orders of magnitude and meanwhile the stagnation-point heat flux value also achieves a constant, indicating a good convergence characteristic. While for the results of AUSM+ scheme shown in Fig. 4, the residual drops
(a) Geometry
160
AUSM+
140
1e-5m 5e-6m 1e-6m
120
Qw(kw/m2)
434
Exp
100 80 60 40 20 0
-1
-0.5
0
0.5
1
x/R Fig. 2. Grid convergence study of AUSM+ scheme for MSL capsule case.
only about 4 orders of magnitude and the stagnation-point heat flux does not converge to a invariant value but has a slight oscillation. Fig. 5 shows the contours of the computed heat flux on the wall of the MSL capsule, using the AUSM+ and AUSMPW schemes, from which it can be found that apparent unphysical oscillations of the heat flux occur for the two schemes, demonstrating the oscillatory convergence behavior for AUSM+ scheme. The computed heat flux distributions on the centerline by different schemes are given in Fig. 6, and it is seen that all the five AUSM family schemes: AUSMDV, AUSM+, AUSMPW, AUSMPW+, AUSM+UP used in this paper give wall heat flux results with unphysical oscillations, while the oscillations by the AUSM+ and AUSMPW schemes are relatively low among the five AUSM family schemes. Additionally, the computed heat flux by the Roe scheme is lower than the experimental data using the current mesh, which can be found in Fig. 6. Fig. 7 illustrates the pressure contours on the symmetric face given by the five AUSM family schemes, in order to show how the shock wave is captured. It is clearly seen that oscillating results arise for all the five AUSM schemes, which further influences the prediction of the wall heat flux and induce unphysical heat flux distributions. The wall heat flux computed by the Roe scheme is shown in Fig. 8. The distribution of the wall heat flux looks very smooth as shown in Fig. 8(a) and meanwhile Fig. 8(b) indicates a fairly stable resolution of the shock wave by the Roe scheme. Moreover, it is also seen from Fig. 6 that the computed heat flux values on the shoulder are even higher than the value on the stagnation point, showing a different pattern from that given by the measured data.
(b) Generated structured mesh
Fig. 1. 3D MSL capsule model.
435
Z.-X. Gao et al. / International Journal of Heat and Mass Transfer 116 (2018) 432–444
100
104
10-1 103
Qw(kw/m2)
Residual
10-2 10-3
102
10-4 101 10-5 10-6
0
2000
4000
6000
8000
100
10000
0
2000
4000
6000
8000
Iteration
Iteration
(a) Residual
(b) Stagnation-point heat flux
10000
100
104
10-1
103
Qw(kw/m2)
Residual
Fig. 3. Convergence histories by Roe scheme for MSL capsule case.
10-2
101
10-3
10-4
102
0
5000
10000
15000
20000
25000
100
0
5000
10000
15000
20000
Iteration
Iteration
(a) Residual
(b) Stagnation-point heat flux
25000
Fig. 4. Convergence histories by AUSM+ scheme for MSL capsule case.
2
2
Qw(kw/m )
Q w(kw/m ) 0
18
36
54
71
89
107
125
(a) AUSM+ scheme
0
20
39
59
78
98
117
137
(b) AUSMPW scheme
Fig. 5. Computational wall heat flux contours of AUSM+ and AUSMPW schemes.
In addition, observing the heat flux values around the stagnationpoint by the AUSM+ and AUSMPW schemes in Fig. 6 and comparing those with the measured stagnation-point heat flux value, it is found that the AUSM+ and AUSMPW schemes predict more accu-
rate stagnation-point heat flux, though unphysical oscillations exist in the distributions of the wall heat flux. The above results clearly reflect the features of the Roe and AUSM family schemes on the hypersonic aeroheating computation,
436
Z.-X. Gao et al. / International Journal of Heat and Mass Transfer 116 (2018) 432–444
300
tor, which refers to the method adopted in the construction of the hybrid high-order schemes [24], and the other one is using the distance from the wall to control the switch between the two kinds of schemes. Next, the accuracy of the aeroheating prediction of the five AUSM schemes is examined at first using a 2D cylinder case, in order to select an appropriate AUSM scheme for the construction of the hybrid scheme. Then, new hybrid schemes will be developed based on the above two techniques (mentioned above), and verified using the 2D cylinder case.
Roe
250
AUSMDV AUSM+ AUSMPW AUSMPW+
Qw(kw/m2)
200
AUSM+UP Exp
150 100
4.1. Comparative study on the aeroheating prediction accuracy of AUSM family schemes and Roe scheme
50 0
-1
-0.5
0
0.5
1
x/R Fig. 6. Comparison of the computed wall heat flux on the centerline with the experimental data.
which have been described in Section 1. The Roe scheme could ensure a stable resolution of the shock wave without strict restriction on the grid generation around the shock wave, and thus the computed wall heat flux shows a smooth distribution. However, the obtained wall heat flux results deviate from the measured data significantly. While for the AUSM family schemes, though they could predict more accurate stagnation-point wall heat flux values, unphysical oscillations exist among the distributions of the computed wall heat flux because of the unstable shock wave resolutions. 4. Development of new hybrid schemes based on the Roe and AUSM schemes Based on the characteristics of the Roe and AUSM family schemes, hybrid schemes are proposed here which would apply the Roe scheme only near the shock wave to ensure the stability of the strong shock wave resolution, while adopt the AUSM family schemes inside the boundary layer to improve the computational accuracy of the wall heat flux. To achieve the above aims, two different ideas are proposed: the first one is based on the shock detec-
P/pa
P/pa
22000 19810 17620 15430 13240 11050 8860 6670 4480 2290 100
(a)AUSMDV
A 2D cylinder experimental case by Holden et al. [25] is selected to perform the comparative study on the aeroheating prediction accuracy of the AUSM family schemes. The radius of the cylinder R is 38.1 mm. Two sets of experimental data are chosen in this study, and the test conditions [25,26] are shown in Table 1. A structured mesh with 101 91 grid points is generated for this cylinder case and shown in Fig. 9(a). The stretching factor near the wall is 1.1 in the normal direction to the wall and the distribution of the grid points in the streamwise direction is uniform. In the same way, the grid points are clustered near the shock wave in the normal direction without any adjustment for the mesh to be aligned with the shock wave. A grid convergence study is performed using three meshes with the first gird spacing off the wall 1 105 m, 1 106 m and 2 107 m, respectively, and the wall heat flux results for AUSM+ scheme and Roe scheme are shown in Fig. 9 (b) and (c). It is seen that the AUSM+ scheme could give a gridconverged and reliable heat flux result using the mesh with the first grid point spacing off the wall 1 106 m, while the result by the Roe scheme agrees with the experimental data until the first grid point spacing is refined to 2 107 m, implying excess dissipation for the Roe scheme. Based on the grid convergence study for the AUSM+ scheme, the mesh with the first grid spacing off the wall 1 106 m is selected to conduct the following simulations. Moreover, the CFL number is taken as 10 for the simulations of the 2D cylinder case, and the convergence histories of the residual and stagnation-point heat flux for Roe and AUSM+ schemes are illustrated in Figs. 10 and 11. Similar to above 3D capsule case, the residual drops more than 6 orders of magnitude for Roe scheme
P/pa
22000 19810 17620 15430 13240 11050 8860 6670 4480 2290 100
(b)AUSM+
P/pa
22000 19810 17620 15430 13240 11050 8860 6670 4480 2290 100
(c)AUSMPW
P/pa
22000 19810 17620 15430 13240 11050 8860 6670 4480 2290 100
(d)AUSMPW+
Fig. 7. Pressure contours on the symmetric face computed by the five AUSM schemes.
22000 19810 17620 15430 13240 11050 8860 6670 4480 2290 100
(e)AUSM+-UP
437
Z.-X. Gao et al. / International Journal of Heat and Mass Transfer 116 (2018) 432–444
2
Qw(kw/m )
0
13
25
38
51
64
76
89
P/pa
22000 19810 17620 15430 13240 11050 8860 6670 4480 2290 100
(a) Wall heat flux contours
(b) Pressure contours
Fig. 8. Computational results by Roe scheme.
Table 1 Free-stream conditions for the cylinder case. Ma1
p1/pa
T1/K
Tw/K
Re/m
8 16.34
855 82.95
125.07 52
294 294.4
4.72 106 3.74 106
vided in Ref. [25] is a little dispersed, it can be still seen that the heat flux distribution computed by the AUSM+ scheme has the best agreement with the measured data. Therefore, the AUSM+ scheme is chosen among the AUSM family schemes for the hybrid scheme together with the Roe scheme. In addition, to examine the influence of the tunable parameter d in the entropy correction of Roe scheme, Fig. 13 gives the computed heat flux distributions by Roe scheme with different values of d. It is seen that a small d could get higher wall heat flux, and with d = 0 the Roe scheme can also give accurate heat flux. However, as mentioned in the Introduction, d needs to be taken as a relative large value for 3D complex flows with strong shock waves, and thus the prediction accuracy of wall heat flux for Roe scheme would be reduced. It is worth noting that for this 2D cylinder case, unphysical oscillations do not occur in the wall heat flux results computed by the AUSM family schemes even when the shape of the mesh is not adjusted to align with the shock wave. It indicates that not all the hypersonic aeroheating simulations using a mesh whose shape is not consistent with the shock wave would generate unphysical oscillations on the wall heat flux for the AUSM family schemes. However, it is likely to result in unphysical oscillations
while only 4 orders for AUSM+ scheme. Both schemes results in a converged wall heat flux prediction. Fig. 12(a) shows the computed heat flux distributions on the wall by different AUSM schemes under the condition of Ma1 = 8. Comparing different numerical heat flux results with the experimental data, it can be found that the heat flux distribution by the AUSM+ scheme agrees with the experimental data best, while the stagnation-point heat flux values by other AUSM schemes are higher than the measured data. Moreover, Fig. 12(a) also presents the result by the Roe scheme using the same structured grid, and the computed heat flux is lower than the experimental data, which again verifies the lower aeroheating prediction accuracy of the Roe scheme. In Fig. 12(b), the wall heat flux under the condition of Ma1 = 16.34 are presented. Although the experimental data pro-
7
7
AUSM+ 1e-5m 1e-6m 2e-7m Exp
5 4 3 2 1 0
Roe
6
Qw/(105W m-2)
Qw/(105W m-2)
6
1e-5m 1e-6m 2e-7m Exp
5 4 3 2 1
0
15
30
45
60
75
90
0
0
15
30
/( )
45
60
75
/( )
(a) Structured
(b) Grid convergence study with
(c) Grid convergence study with
mesh
Ma∞=16.34 for AUSM+
Ma∞=16.34 for Roe
Fig. 9. 2D structural mesh and grid convergence study of AUSM + scheme for the cylinder.
90
438
Z.-X. Gao et al. / International Journal of Heat and Mass Transfer 116 (2018) 432–444
10
1
10
3
100 -1
Qw(105W/m2)
Residual
10
10-2 10-3 10
-4
102
101
10-5 10-6
0
2000
4000
6000
8000
10000
100
12000
0
2000
4000
6000
8000
10000
Iteration
Iteration
(a) Residual
(b) Stagnation-point heat flux
12000
Fig. 10. Convergence histories by Roe scheme for cylinder with Ma1 = 16.34
0
10
-1
10
-2
10
-3
10
Qw(105W/m2)
Residual
10
3
102
101
10-4
10
-5
0
5000
10000
15000
20000
10
25000
0
0
5000
10000
15000
20000
Iteration
Iteration
(a) Residual
(b) Stagnation-point heat flux
25000
Fig. 11. Convergence histories by AUSM+ scheme for cylinder with Ma1 = 16.34.
7
7
6
6 Roe AUSMDV AUSM+ AUSMPW AUSMPW+ AUSM+UP Exp
4 3
4 3
2
2
1
1
0
10
20
30
40
50
60
70
80
ROE AUSMDV AUSM+ AUSMPW AUSMPW+ AUSM+UP Exp
5
Qw/(105W m-2)
Qw/(105W m-2)
5
90
0
15
30
45
/( )
/( )
(a) Ma = 8
( b) Ma =16.34
Fig. 12. Wall heat flux distributions computed by different schemes.
60
75
90
439
Z.-X. Gao et al. / International Journal of Heat and Mass Transfer 116 (2018) 432–444
7
on the concept of smooth measurement factor (SMF) originally proposed in the construction of the WENO scheme [31]. For a grid point j, three stencils, each of which includes three grid points, are designed using the four adjacent grid points on both sides of j, as shown in Fig. 14. Each stencil has its own value of SMF. ¼ p=ðq1 V 21 Þ is employed to idenThe dimensionless pressure p tify the shock wave location. After mathematic deductions [30], the expressions for the SMFs on the three stencils for grid point j are
6 =0 =0.06 =0.125 Exp
Qw/(105W m-2)
5 4
8 j2 4p j1 þ 3p j Þ2 þ 13 j2 2p j1 þ p j Þ2 > SMF 1 j ¼ 14 ðp ðp > 12 < j1 p jþ1 Þ2 þ 13 j1 2p j þ p jþ1 Þ2 SMF 2j ¼ 14 ðp ðp 12 > > : j 4p jþ1 þ p jþ2 Þ2 þ 13 j 2p jþ1 þ p jþ2 Þ2 SMF 3j ¼ 14 ð3p ðp 12
3 2 1 0
dy=1e-6 m
0
15
30
60
45
75
90
/( ) Fig. 13. Wall heat flux distributions computed by Roe schemes with different d values.
ð5Þ
Generally, the SMFs defined by the above equations would approach to a small value in the smooth region, while they would be much larger in the discontinuous region. In present study, the SMF is introduced only to detect the discontinuity of the flow field, so it is reasonable to choose the maximum one among the three SMFs for point grid j. Furthermore, for 3D simulations, Eq. (5) can still be adopted to compute the three SMFs in each dimension, and then, for each SMF, the following treatment is adopted:
SMF m ¼
qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi SMF 2mi þ SMF 2mj þ SMF 2mk
ð6Þ
where the maximum value among the three SMFs is chosen for grid point (i, j, k). Then a damping function fD, whose value would be 0 near the shock wave but 1 in other smooth region, is designed based on the SMFs above. The hyperbolic tangent function is introduced by Zhang et al. [30] to define the damping function as follow:
f D ¼ 1 tanh ½cam maxðSMF 1 ; SMF 2 ; SMF 3 Þ
Fig. 14. Stencil selection for the construction of SMF at grid point j.
using the AUSM family schemes for 3D complex hypersonic flows with strong shock waves. 4.2. Construction of new hybrid schemes (1) Hybrid scheme based on the shock wave detector We evaluated the shock wave detectors proposed by Visbal and Gaitonde [24], Hill and Pullin [27], Larsson et al. [28], Priozzoli [29] and Zhang et al. [30], and eventually the one established by Zhang et al. [30] is selected. This shock wave detector proposed is based
Ma
where the cam is an amplification factor for the discontinuity in flow field, which could make the hyperbolic tangent value of the SMFs near the shock wave close to 1. According to the theoretical analysis by Zhang et al. [30], the value of cam is taken as 95. Here, the 2D cylinder is employed to examine the effectiveness of the above shock wave detector. Fig. 15 shows the results of the shock wave detector in the flowfield around the 2D cylinder with Ma1 = 8 and 16.34, where the top half is the Mach number contours, and the bottom half is the fD contours. It can be seen that the shock wave detector achieves a high-fidelity capture of the
Ma
8.0 6.9
16.0
13.7 11.4 9.1
5.7 4.9 3.8
6.9 4.6 2.3
2.8 1.6 0.0
0.0
fD
fD 1.0
1.0
0.9 0.7
0.8 0.6 0.5
0.6 0.4 0.3
0.3 0.2 0.1
0.1 0.0
(a) Ma =8
ð7Þ
0.0
(b) Ma =16.34
Fig. 15. Contours of the damping function and Mach number for 2D cylinder.
Z.-X. Gao et al. / International Journal of Heat and Mass Transfer 116 (2018) 432–444
6
7
5
6 Hybrid-1st Hybrid-2ed- 0 Hybrid-2nd-2 0 Hybrid-2nd-5 0 Hybrid-2nd-10 0 Hybrid-2nd-0.2 0 Hybrid-2nd-0.5 0 AUSM+ Exp
4 3 2
4 3 2
1 0
Hybrid-1st Hybrid-2nd- 0 Hybrid-2nd-2 0 Hybrid-2nd-5 0 Hybrid-2nd-10 0 Hybrid-2nd-0.2 0 Hybrid-2nd-0.5 0 AUSM+ Exp
5
Qw/(105W m-2)
Qw/(105W m-2)
440
1 0
10
20
30
40
50
60
70
80
0
90
0
20
40
/( )
60
80
/( )
(a) Ma =8
(b) Ma =16.34
Fig. 16. Wall heat flux distributions computed by the hybrid schemes.
10
1
10
3
100 -1
Qw(105W/m2)
Residual
10
10-2 10-3 10
-4
102
101
10-5 10-6
0
3000
6000
9000
12000
100
0
3000
6000
9000
Iteration
Iteration
(a) Residual
(b) Stagnation-point heat flux
12000
Fig. 17. Convergence histories by the second hybrid scheme for cylinder with Ma1 = 16.34.
shock wave and the fD value is close to 0 near the shock wave but 1 in other region. Finally, the numerical flux of the hybrid scheme is constructed based on the above shock wave detector:
e F iþ1 ¼ ð1 f D Þ e F AUSMþ F Roe þ f De iþ1 iþ1 2
2
2
ð8Þ
(2) Hybrid scheme based on the distance from the wall The above hybrid scheme based on the shock wave detector applies to the Roe scheme around the shock wave and the AUSM + scheme in other region, as can be seen in Fig. 15. Yet another way is to construct the hybrid scheme is to limit the utilization of the AUSM+ scheme only within the boundary layer and adopt the Roe scheme in other region including the shock wave region. This can be realized by designing a control function based on the distance from the wall. Here, a simple function form as follow is employed:
Hðd d0 Þ ¼
1 d > d0 0
d < d0
ð9Þ
where d0, an empirical parameter, would control the switch from the AUSM+ scheme to the Roe scheme. Using this function, the numerical flux of this hybrid scheme could be given as
e F iþ1 ¼ Hðd d0 Þ e F Roe þ ½1 Hðd d0 Þ e F AUSMþ iþ1 iþ1 2
2
2
ð10Þ
As for the parameter d0, its value should be a reflection of the boundary layer thickness around the vehicle. Because there is no available expression for the boundary layer thickness in hypersonic laminar flows, the following formula is used to calculate the boundary layer thickness for the incompressible flat-plate flow to give an estimation for d0:
pffiffiffiffiffiffiffiffi d0 ¼ L= ReL
ð11Þ
where ReL is the Reynolds number computed using the free-stream parameters and the characteristic length. Eq. (11) may not accurately estimate the boundary layer thickness for a hypersonic flow over blunt reentry vehicles. However, if the aeroheating computation using this hybrid scheme is not sensitive to the value of d0, which will be examined in the following numerical experiments, the estimation in Eq. (11) for d0 could be acceptable.
441
Z.-X. Gao et al. / International Journal of Heat and Mass Transfer 116 (2018) 432–444
Qw(Kw/m2)
0
23
46
69
91
114
137
160
P/pa
22000 19810 17620 15430 13240 11050 8860 6670 4480 2290 100
(a) wall heat flux contour
(b) pressure contour
Fig. 18. Computational results by the first hybrid scheme.
Qw(Kw/m2)
0
17
34
51
(a)
67
94
118
Qw(Kw/m2)
135
0
19
39
58
(b) 0.5
0
77
96
116
135
0
Fig. 19. Computational results by the second hybrid scheme with different distances from the wall.
160
=0
140
Hybrid-2nd-0.5
120
Hybrid-2ndRoe
0
0
Qw(kw/m2)
Exp
100 80 60 40 20 0
-1
-0.5
0
0.5
1
x/R Fig. 20. Comparison of the wall heat flux distributions on the centerline with experimental data.
4.3. Application of the hybrid schemes in the 2D cylinder case The 2D cylinder case is used to examine the effectiveness of the above two hybrid schemes. In the following numerical experi-
ments, the hybrid scheme using the shock wave detector is noted as the first hybrid scheme and the hybrid scheme based on the distance from the wall is called the second hybrid scheme. The computed heat flux distributions by the two hybrid schemes are shown in Fig. 16 (a: Ma1 = 8; b: Ma1 = 16.34), and the results obtained by the AUSM+ scheme are also given in Fig. 16. It can be seen that the heat flux results by the hybrid schemes agree well with the experimental data for the two cases. In addition, the results obtained by the hybrid scheme are very close to those using the AUSM+ scheme independently, which indicates that, for this 2D cylinder case, in spite of the introduction of numerical dissipation by the Roe scheme into the flowfield, the adoption of the AUSM+ scheme in the boundary layer could still maintain almost the same aeroheating prediction accuracy as applying the AUSM+ scheme in the whole flowfield. Furthermore, Fig. 16 also depicts the results for the hybrid scheme based on the distance from the wall with different values of d0, where the baseline d0 value is computed using Eq. (11). It is found that there is little influence on the wall heat flux prediction when d0 changes from 0.2 times of the baseline value to 10 times, which proves that the aeroheating simulation using the second hybrid scheme is not sensitive to the value of d0. Moreover, the convergence histories of the residual and stagnation-point heat flux for the second hybrid scheme are shown in Fig. 17. Compared to those for the AUSM+ scheme in Fig. 11, the hybrid scheme obviously improves the convergence characteristic
442
Z.-X. Gao et al. / International Journal of Heat and Mass Transfer 116 (2018) 432–444
and meanwhile could get a converged stagnation-point heat flux with almost same iteration numbers. For this 2D cylinder case, there is no unphysical oscillation existing in the wall heat flux distributions even solely applying the AUSM+ scheme, so this case indicates that the hybrid schemes could maintain a close computational accuracy of the wall heat flux compared to the AUSM+ scheme. However, whether the two hybrid schemes could eliminate the unphysical oscillations of the wall heat flux distributions induced by the unstable resolution of the shock wave should be further verified.
10
5. Aeroheating simulations on 3D MSL capsule using the hybrid schemes Here, the two hybrid schemes are adopted to simulate the 3D MSL capsule case in Section 3. Fig. 18 illustrates the results obtained by the hybrid scheme based on the shock wave detector. From the Fig. 18(a), it is seen that the obtained wall heat flux distribution still shows apparent unphysical oscillations. It indicates that the hybrid scheme based on the shock wave detector could not improve the stability of the shock wave resolution, which
0
4
10
10-1
103
10
Qw(kw/m2)
Residual
10
-2
-3
2
10
10-4
101 10-5 10
-6
0
0
5000
10000
15000
10
0
5000
10000
Iteration
Iteration
(a) Residual
(b) Stagnation-point heat flux
Fig. 21. Convergence histories by the second hybrid scheme for MSL capsule case.
Qw(kw/m2)
0
19
37
56
74
93
111
130
Qw(kw/m2)
0
19
39
(a) =4° Qw(kw/m2)
0
19
39
58
77
(c) =12°
58
77
96
116
135
114
137
160
(b) =8°
96
116
135
Qw(kw/m2)
0
23
46
69
91
(d) =24°
Fig. 22. Wall heat flux contours using the second hybrid scheme under different angles of attack.
15000
443
Z.-X. Gao et al. / International Journal of Heat and Mass Transfer 116 (2018) 432–444
can also be found from Fig. 18(b). One possible reason is that this hybrid scheme only activates the Roe scheme around the shock wave while use the AUSM+ scheme in other region, and meanwhile the switch between these two schemes occurs near the shock wave, which will not eliminate the oscillating solutions around the shock wave. Then, the hybrid scheme based on the distance from the wall is adopted to perform the numerical simulations on the 3D MSL capsule model, and the parameter d0 is estimated according to Eq. (11), whose value is approximately 1.35 104 m. Fig. 19(a) shows the wall heat flux contours computed by this hybrid scheme, and it can be seen that the heat flux oscillations are totally eliminated, proving the effectiveness of this hybrid scheme. In addition, Fig. 19(b) illustrates the wall heat flux result with the same hybrid scheme but the empirical parameter is changed to 0.5d0. The difference between Fig. 19(a) and (b) is limited, indicating that the wall heat flux prediction of this hybrid scheme is not sensitive to d0 in Eq. (10) for 3D case. Fig. 20 compares the wall heat flux distributions on the centerline computed by this hybrid scheme with the measured data, and good agreements are achieved for both d0 and 0.5d0. Compared to the results by solely using the Roe or AUSM+ schemes, significant improvements on the aeroheating prediction are observed by this hybrid scheme. The convergence histories for the MSL capsule case using the second hybrid scheme are given
160
in Fig. 21. By comparing Figs. 21 and 4, we demonstrate that the second hybrid scheme improves the convergence characteristic relative to the AUSM+ scheme. In addition to the case with a = 0°, a series of aeroheating simulations with a = 4°, 8°, 12° and 24° are conducted applying the hybrid scheme based on the distance from the wall. In these cases, the shape of the shock wave would change corresponding to different angles of attack, while keeping the same computational mesh with the one used for a = 0° case. The computed wall heat flux contours are shown in Fig. 22 for the four cases, and it is found that this hybrid scheme can maintain smooth distributions of the wall heat flux without any unphysical oscillation for all the four cases with different angles of attack. Fig. 23 compares the wall heat flux distributions on the centerline for the four cases with the corresponding experimental data, and reasonable agreements can be seen for all the four cases. The above applications of the hybrid schemes in the MSL capsule case demonstrate that the hybrid scheme based on the distance from the wall could effectively eliminate the unphysical oscillations in the distributions of the wall heat flux which are generated when solely using the AUSM+ scheme, and enhance the aeroheating accuracy compared to that of the Roe scheme. Moreover, for different freestream conditions, this hybrid scheme can also guarantee the accuracy of the wall heat flux prediction using
160
=4
140
100 80 60
80 60 40
20
20 -0.5
160
0
0.5
0 -1
1
0
x/R
(a) =4°
(b) =8° 160
=12
0.5
1
=24
140 Hybrid-2nd Exp
120
Qw(kw/m2)
80 60
100 80 60
40
40
20
20 -1
-0.5
0
0.5
Hybrid-2nd Exp
120
100
0
-0.5
x/R
140
Qw(kw/m2)
100
40
0 -1
Hybrid-2nd Exp
120
Qw(kw/m2)
Qw(kw/m2)
140
Hybrid-2nd Exp
120
=8
1
0 -1
-0.5
0
x/R
x/R
(c) =12°
(d) =24°
0.5
Fig. 23. Comparisons of the wall heat flux distributions on the centerline with experimental data under different angles of attack.
1
444
Z.-X. Gao et al. / International Journal of Heat and Mass Transfer 116 (2018) 432–444
the same mesh. The above advantages make the proposed hybrid scheme appropriate for the aeroheating prediction of the hypersonic reentry vehicles. 6. Conclusion Hypersonic aeroheating computations for a 3D capsule from the Mars Science Laboratory (MSL) are performed using the Roe and AUSM family schemes, and it is found that the AUSM family schemes may generate unphysical oscillations in the wall heat flux distributions due to their unstable resolutions of the shock wave, while the wall heat flux prediction accuracy by the Roe scheme is not adequate though the obtained heat flux distribution is smooth. In order to develop a numerical scheme that could capture the shock wave stably and meanwhile maintain a high-fidelity simulation of the wall heat flux, new hybrid schemes combining the Roe and AUSM+ schemes are developed in the present paper. Two different ideas are proposed to construct two hybrid schemes. The first hybrid scheme is to use the shock wave detector to detect the location of the shock wave, which would activate the Roe scheme near the shock wave and apply the AUSM+ scheme in other region. The second one is to utilize the distance from the wall to control the switch between the two schemes, and apply AUSM+ only in the boundary layer while the Roe scheme in other region. The applications of these two hybrid schemes for the 3D MSL capsule case show that the first hybrid scheme based on the shock wave detector does not show improvement in the unphysical oscillations of the wall heat flux, while the second hybrid scheme removes the unphysical oscillations and maintains accurate simulation of the wall heat flux. Additional computations of the 3D MSL capsule case with different angles of attack further demonstrate that the proposed hybrid scheme based on the distance from the wall could effectively eliminate the unphysical oscillations in the distributions of the wall heat flux which are generated when solely using the AUSM+ scheme, and remarkably enhance the aeroheating accuracy compared to that of the Roe scheme. The proposed hybrid scheme based on the distance from the wall is shown to be appropriate for the aeroheating prediction of the hypersonic reentry vehicles. Acknowledgments This work was supported by the National Natural Science Foundation of China through the Grant 91641123 and the Academic Excellence Foundation of BUAA for PhD students. Conflict of interest All the authors (Zhenxun Gao, Haichao Xue, Zhichao Zhang, Hongpeng Liu, Chun-Hian Lee) of the article ‘‘A Hybrid Numerical Scheme for Aeroheating Computation of Hypersonic Reentry Vehicles” state that there is no financial/personal interest or belief that could affect our objectivity, and all authors have seen and approved the final version of the manuscript being submitted. We warrant that the article is the authors’ original work. References [1] K.A. Hoffmann, M.S. Siddiqui, S.T. Chiang, Difficulties associated with the heat flux computations of high speed flows by the Navier-Stokes equations, in: AIAA paper, 1991, AIAA-1991-0457.
[2] G.H. Klopfer, H.C. Yee, Viscous hypersonic shock-on-shock interaction on blunt cowl lips, in: AIAA paper, 1988, AIAA-1988-0233. [3] F.G. Blottner, D.E. Larson, Navier-Stokes Code NS3D for Blunt Bodies, Sandia Report, 1988, SAND-1988-050411 UC-32. [4] M.S. Siddiqui, K.A. Hoffmann, A comparative study of the Navier-Stokes solvers with emphasis on the heat transfer computations of high speed flows, in: AIAA paper, 1992, AIAA-1992-0835. [5] P. Papadopoulos, E. Venkatapaphy, D. Prabhu, Current grid generation strategies and future requirements in hypersonic vehicle design, analysis and testing, Appl. Math. Model. 23 (9) (1999) 705–735, http://dx.doi.org/10.1016/ S0307-904X(99)00007-4. [6] I.S. Men’ shov, Y. Nakamura, Numerical simulations and experimental comparison for high-speed nonequilibrium air flows, Fluid Dyn. Res. 27 (5) (2000) 305–334, http://dx.doi.org/10.1016/S0169-5983(00)00010-1. [7] S.J. Henderson, J.A. Menart, Grid Study on Blunt Bodies with the Carbuncle Phenomenon, in: AIAA paper, 2007, AIAA-2007-3904. [8] K. Kitamura, P. Roe, F. Ismail, An evaluation of Euler fluxes for hypersonic flow computations, AIAA J. 47 (1) (2009) 44–53, http://dx.doi.org/10.2514/1.33735. [9] K. Kitamura, E. Shima, Y. Nakamura, P. Roe, Evaluation of Euler fluxes for hypersonic heating computations, AIAA J. 48 (4) (2010) 763–776, http://dx.doi. org/10.2514/1.41605. [10] J.H. Lee, O.H. Rho, Numerical analysis of hypersonic viscous flow around a blunt body using Roe’s FDs and AUSM + schemes, in: AIAA paper, 1997, AIAA1997-2054. [11] J.H. Lee, O.H. Rho, Accuracy of AUSM + scheme in hypersonic blunt body flow calculation, in: AIAA paper,1998, AIAA-1998-1538. [12] Y. Wada, M.-S. Liou, An accurate and robust flux splitting scheme for shock and contact discontinuities, SIAM J. Sci. Stat. Comput. 18 (1997) 633–657, http:// dx.doi.org/10.1137/S1064827595287626. [13] M.-S. Liou, A Sequel to AUSM: AUSM+, J. Comput. Phys. 129(2) (1996) 364– 382. 10.1006/jcph.1996.0256. [14] M.-S. Liou, Ten Years in the Making—AUSM Family, in: AIAA paper, 2001, AIAA2001-2521. [15] S.S. Kim, C. Kim, O.H. Rho, S.K. Hong, Methods for the accurate computations of hypersonic flows I. AUSMPW + Scheme, J. Comput. Phys. 174 (1) (2001) 38–80, http://dx.doi.org/10.1006/jcph.2001.6873. [16] M.-S. Liou, A Further Development of the AUSM + Scheme Towards Robust and Accurate Solutions for All Speeds, in: AIAA Paper, 2003, AIAA-2003-4116. [17] P.L. Roe, Approximate Riemann solvers, parameter vectors and difference schemes, J. Comput. Phys. 43 (1981) 357–372, http://dx.doi.org/10.1016/00219991(81)90128-5. [18] Z.X. Gao, C.W. Jiang, C.H. Lee, Improvement and application of wall function boundary condition for high-speed compressible flows, Sci. China Technol. Sci. 56 (10) (2013) 2501–2515, http://dx.doi.org/10.1007/s11431-013-5349-4. [19] Z.X. Gao, C.H. Lee, A numerical study of turbulent combustion characteristics in a combustion chamber of a scramjet engine, Sci. China Technol. Sci. 53 (8) (2013) 2111–2121, http://dx.doi.org/10.1007/s11431-010-3088-3. [20] Z.X. Gao, C.W. Jiang, S.W. Pan, C.H. Lee, Combustion Heat-release effects on supersonic compressible turbulent boundary layers, AIAA J. 53 (7) (2015) 1949–1968, http://dx.doi.org/10.2514/1.J053585. [21] H.C. Yee, G.H. Kolpfer, J.L. Montague, High-resolution shock capturing schemes for inviscid and viscous hypersonic flows, 1988, NASA TM 100097. [22] R.H. Brian, Turbulent aeroheating testing of mars science laboratory entry vehicle, J. Spacecraft Rockers 45 (3) (2008) 417–427, http://dx.doi.org/ 10.2514/1.31798. [23] R.H. Brian, Turbulent Aeroheating Testing of Mars Science Laboratory Entry Vehicle in Perfect-Gas Nitrogen, in: AIAA Paper, 2007, AIAA-2007-1208. [24] M.R. Visbal, D.V. Gaitonde, Shock capturing using compact differencing based methods, in: AIAA Paper, 2005, AIAA-2005-1265. [25] M.S. Holden, A.R. Wieting, J.R. Mpselle, Studies of aerothermal loads generated in regions of shock/shock interaction in hypersonic flow, in: AIAA Paper, 1988, AIAA-1988-0477. [26] K.M. Reery, S.T. Imlay, Blunt-Body Flow Simulations, in: AIAA Paper, 1988, AIAA-1988-2904. [27] D.J. Hill, D.I. Pullin, Hybrid tuned center-difference-WENO method for large eddy simulations in the presence of strong shocks, J. Comput. Phys. 194 (2004) 435–450, http://dx.doi.org/10.1016/j.jcp.2003.07.032. [28] J. Larsson, S.K. Lele, P. Moin, Effect of numerical dissipation on the predicted spectra for compressible turbulence, Center for turbulence research, annual research briefs (2007) 47–57. [29] S. Priozzoli, Numerical methods for high-speed flows, Ann. Rev. Fluid Mech. 43 (2011) 163–194, http://dx.doi.org/10.1146/annurev-fluid-122109-160718. [30] Z.C. Zhang, Z.X. Gao, C.W. Jiang, C.H. Lee, A RANS model correction on unphysical over-prediction of turbulent quantities across shockwave, Int. J. Heat Mass Transf. 106 (2017) 1107–1119, http://dx.doi.org/10.1016/j. ijheatmasstransfer.2016.10.087. [31] G.S. Jiang, C.W. Shu, Efficient implementation of weighted ENO schemes, J. Comput. Phys. 126 (130) (1996) 202–228, http://dx.doi.org/10.1006/ jcph.1996.0130.
Contents lists available at ScienceDirect
International Journal of Heat and Mass Transfer journal homepage: www.elsevier.com/locate/ijhmt
A hybrid numerical scheme for aeroheating computation of hypersonic reentry vehicles Zhen-Xun Gao a,⇑, Hai-Chao Xue b, Zhi-Chao Zhang a, Hong-Peng Liu a, Chun-Hian Lee a a b
National Laboratory for Computational Fluid Dynamics, School of Aeronautic Science & Engineering, Beihang University, Beijing 100191, China China Academy of Aerospace Aerodynamics, Beijing 100074, China
a r t i c l e
i n f o
Article history: Received 24 August 2016 Received in revised form 20 July 2017 Accepted 21 July 2017 Available online 17 October 2017 Keywords: Aeroheating Numerical scheme Hypersonic Blunt body
a b s t r a c t Hypersonic aeroheating simulations were conducted for three-dimensional (3D) reentry capsule for the Mars Science Laboratory (MSL) using the AUSM and the Roe schemes. The study confirms that AUSM family of schemes generate unphysical oscillations in the wall heat flux distributions due to the unstable simulation of the shock wave. And the wall heat flux is underpredicted by the Roe scheme due to excess dissipation, even though captures the shock wave accurately. In order to overcome these limitations, we develop a numerical scheme that can capture the shock wave stably and still maintain an accurate prediction of the heat flux. We propose a novel hybrid scheme that combines both the AUSM+ and Roe schemes. Two different hybrid techniques are proposed for switching between the two schemes, one based on the shock wave detection and other utilizing a criterion based on the distance from the wall. The implementation of the two hybrid schemes for the 3D MSL reentry capsule was studied for prediction accuracy compared to experimental data. It is found that the second hybrid scheme (based on the criterion for switching based on distance from the wall) can eliminate the unphysical oscillations in the wall heat flux distributions and improve the aeroheating prediction accuracy. Ó 2017 Elsevier Ltd. All rights reserved.
1. Introduction Space technology has been developing rapidly in recent years, and the problem of thermal protection for the hypersonic reentry vehicles is becoming increasingly prominent. High-precision prediction for the aerothermal environment is one of the crucial technologies for the design of a thermal protection system. With the advent of numerical methods and computational capability, computational fluid dynamics (CFD) gradually becomes an important tool for the hypersonic aeroheating prediction. However, it is still very difficult to achieve a robust and accurate wall heat flux prediction for hypersonic reentry vehicles using CFD. Previous research revealed that the grid-point distribution and numerical schemes are the two important factors influencing the accuracy of the heat flux prediction in CFD [1–11] for a hypersonic laminar flow. As early as 1990s, a series of numerical experiments were conducted by Hoffmann et al. [1] based on the studies of Klopfer and Yee [2] and Blottner and Larson [3], and the results showed that the computed wall heat flux depends highly on the ⇑ Corresponding author. E-mail addresses: [email protected] (Z.-X. Gao), xuehc31201@hotmail. com (H.-C. Xue), [email protected] (Z.-C. Zhang), [email protected] (H.-P. Liu), [email protected] (C.-H. Lee). http://dx.doi.org/10.1016/j.ijheatmasstransfer.2017.07.100 0017-9310/Ó 2017 Elsevier Ltd. All rights reserved.
density of the grid points near the wall. Subsequently, some researchers studied in detail the effects of the first grid point spacing off the wall on the accuracy of the computed heat flux, resulting in several grid generation criterions proposed to estimate an appropriate first grid point spacing off the wall for an accurate heat flux computation [5,6]. However, for a high Mach number flow over a blunt body where a strong detached shock wave exists, the accuracy of the computed wall heat flux depends not only on the near-wall grid spacing, but also the grid quality in the vicinity of the shock wave. By numerical experiments, Henderson and Menart [7] showed that, in order to achieve a reasonable aeroheating prediction, the shape of the mesh needs to be in accordance to the shock wave and grid points should also be clustered normal to the shock wave. A low-quality grid-point distribution near the shock wave is easy to lead to an oscillating resolution of the shock wave which would cause an anomalous distribution of the wall heat flux. While for the numerical schemes, Kitamura et al. [8,9] argued that a numerical scheme suitable for hypersonic wall heat flux computation should be equipped with the following three properties: shock stability/robustness, conservation of total enthalpy and accurate resolution of boundary layer. In Ref. [9], 15 kinds of numerical schemes were evaluated according to the above three properties, and it was found that no numerical scheme investigated in [9] possessed all the three properties, while the
433
Z.-X. Gao et al. / International Journal of Heat and Mass Transfer 116 (2018) 432–444
Nomenclature Ma
q T p uj H k
l sij e E R
a
free-stream Mach number density (kg/m3) temperature (K) pressure (pa) velocity (m/s) total enthalpy (J/kg) molecular thermal conductivity (W/(m K)) molecular viscosity (N s/m2) shear stress (N/m2) internal energy (J/kg) total energy (J/kg) gas constant (J/(kg K)) the angle of attack (deg)
AUSM family schemes [12–16] and Roe [17] scheme appeared to be promising. However, these two kinds of numerical schemes have different characteristics for the aeroheating prediction. The AUSM family schemes do not have tunable factors in their scheme design and generally have low numerical dissipation, and thus they could resolve the near-wall temperature gradient more precisely within the boundary layer. Hence, with the same grid spacing near the wall, AUSM family schemes could give a more accurate heat flux result than Roe scheme [10,11]. However, for a flow with strong shock wave, if the shape of the mesh is not in accordance to the shock wave or the grid points are not dense enough across the shock wave, the AUSM family schemes may resolve the shock wave with obvious oscillations [10], which would further generate an unphysical distribution of heat flux on the wall. While for the Roe scheme, an entropy correction function with a tunable factor is incorporated, which could help to control the added numerical dissipation and guarantee a stable shock wave resolution even with a low-quality mesh near the shock wave, but would also bring in additional numerical dissipation inside the boundary layer and thus decrease the prediction accuracy of the wall heat flux. Generally, a small value of the tunable factor in the entropy correction function could lead to relative accurate heat flux predictions for the Roe scheme. However, for a three-dimensional (3D) complicated case with high Mach number, this tunable factor needs to be taken as a large value to ensure a stable resolution of the shock wave and then the numerical resolution inside the boundary layer by the Roe scheme would be significantly reduced. In practical hypersonic aeroheating computations, for a twodimensional (2D) problem, it is still feasible to keep the mesh aligned with the shape of the shock wave and then maintain a stable resolution of the shock wave, so the AUSM family schemes could be adopted to enhance the accuracy of the wall heat flux computation. While for 3D complex flows, it would be difficult to adjust the computational mesh to achieve an accordance of the grid shape and shock wave. Even if it could be realized, the mesh will have to be re-generated for different free-stream conditions because the shape of the shock wave would change, thereby increasing the work for the mesh generation. To remove the above restriction on the grid generation for stable resolution of the shock wave, the Roe scheme is a better choice than the AUSM family schemes. However, the high-dissipation characteristic of Roe scheme in boundary layer is likely to introduce errors in the wall heat flux prediction. Thus Roe scheme generally needs more dense grid points near the wall than the AUSM schemes in order to obtain reasonable heat flux, which would severely limit the computational time step and thus dramatically increase the computational consumption especially for 3D complicated cases. The above dilemma indicates that more efforts should be put on the develop-
local heat flux (W/m2) Reynolds number radius of cylinder (m) damping function amplification factor of damping function the distance of the local grid point from the wall (m) empirical parameter the characteristic length (m)
Qw Re R fD Cam d d0 L
Subscripts 1 free-stream w at the wall
ment of a reliable numerical scheme which could guarantee a stable resolution of strong shock wave without strict restriction on the grid generation while in the meantime maintain an accurate computation of the wall heat flux for hypersonic vehicles. The objective of the present paper is to develop a hybrid numerical scheme combining the Roe and AUSM schemes, which would adopt the Roe scheme around the shock wave to guarantee a stable resolution of the shock wave without strict restriction on the grid generation while utilize the AUSM scheme within the boundary layer to maintain low numerical dissipation and high accuracy of the wall heat flux computation. First, we introduce the governing equations and numerical methodology used in the present study. Secondly, hypersonic aeroheating computations for a 3D reentry capsule from the Mars Science Laboratory (MSL) are performed using the Roe and AUSM family schemes respectively, and the aeroheating accuracy is examined independently for each of the two schemes. Then, new hybrid schemes combining the Roe and AUSM family schemes are proposed following two different techniques, and the performance of the new hybrid schemes is evaluated. Finally, by applying the hybrid schemes to the 3D reentry capsule case, improvements in aeroheating prediction accuracy of the hybrid schemes are examined in detail. 2. Governing equations and numerical methodology The compressible Navier-Stokes (N-S) equations can be written as follow: @ðqÞ @t
þ @x@ j ðquj Þ ¼ 0 @s
@ðqui Þ @t
þ @x@ j ðqui uj þ dij pÞ ¼ @xijj @ðqEÞ @T þ @x@ ðqHuj Þ ¼ @x@ sij ui þ k @x @t j
j
ð1Þ
j
where sij is the shear stress given as:
sij ¼ l
@ui @uj 2 @ui l þ dij 3 @xj @xj @xi
ð2Þ
E is the total energy defined as:
1 E ¼ e þ ui ui 2
ð3Þ
The equation of state is
p ¼ qRT
ð4Þ
All the numerical simulations are carried out using a finite difference code, ACANS (Aerodynamic, Combustion and Aerothermodynamic Numerical Simulation), developed by the authors. The reliability of ACANS code has been validated by a variety of simu-
Z.-X. Gao et al. / International Journal of Heat and Mass Transfer 116 (2018) 432–444
lations [18–20]. In this study, the inviscid flux vectors are discretized by the Roe scheme [17] and AUSM family schemes including AUSMDV [12], AUSM+ [13], AUSMPW [14], AUSMPW+ [15], AUSM+UP [16]. The entropy correction function in the Roe scheme is the Yee type [21], and the value of its tunable factor is taken as 0.125. The viscous flux vectors are discretized by the central-difference scheme. The massively parallel computation is realized using the message passing interface (MPI) based on the space-decomposed algorithm. In this study, a laminar perfect-gas model is used. 3. Hypersonic aeroheating computation of a 3D reentry capsule In this section, a 3D reentry capsule from the MSL [22,23] is employed for performing the aeroheating numerical simulations applying Roe and AUSM family schemes, respectively. This reentry capsule was tested in AEDC (the US Air Force’s Arnold Engineering Development Center) Tunnel 9 and the experimental condition of Run 3021 in [22] is selected in the present study: Ma1 = 9.47, p1 = 167.9 Pa, T1 = 54.8 K, a = 0°. The wall boundary condition is assumed as isothermal with Tw = 300 K. Only the forebody of the MSL capsule is kept in the simulations, as shown in Fig. 1(a). Meanwhile, Fig. 1(b) illustrates the 3D structured mesh generated for the forebody of MSL capsule with a total 243,700 grid points. There are 100 grid points in the direction normal to the wall and the stretching factor near the wall is 1.1. In the streamwise and spanwise directions, the grid points are clustered around the nose and shoulder of the capsule with the stretching factor no more than 1.2. Moreover, the grid points are clustered near the shock wave region in the normal direction, but the shape of mesh is not aligned with the shock wave. Three meshes with different distances of the first grid point off the wall are adopted to perform a grid convergence study: (a) 1 105 m, (b) 5 106 m, (c) 1 106 m. Fig. 2 shows the comparisons among the wall heat flux by AUSM+ scheme, using the above three meshes. It can be seen that the distribution of heat flux obtained using grid (b) nearly coincides with that by grid (c). Hence, for the following computations of the MSL capsule, grid (b) with the distance of the first grid point off the wall 5 106 m is adopted. All the computations for the MSL capsule case are performed with CFL number of 10 for Roe scheme while 5 for AUSM family schemes. Figs. 3 and 4 show the convergence histories of the Roe and AUSM+ schemes for the MSL capsule case. Both the residual and the stagnation-point heat flux are used to examine the convergence of the computations. It is seen from Fig. 3 that, for Roe scheme, the residual drops more than 6 orders of magnitude and meanwhile the stagnation-point heat flux value also achieves a constant, indicating a good convergence characteristic. While for the results of AUSM+ scheme shown in Fig. 4, the residual drops
(a) Geometry
160
AUSM+
140
1e-5m 5e-6m 1e-6m
120
Qw(kw/m2)
434
Exp
100 80 60 40 20 0
-1
-0.5
0
0.5
1
x/R Fig. 2. Grid convergence study of AUSM+ scheme for MSL capsule case.
only about 4 orders of magnitude and the stagnation-point heat flux does not converge to a invariant value but has a slight oscillation. Fig. 5 shows the contours of the computed heat flux on the wall of the MSL capsule, using the AUSM+ and AUSMPW schemes, from which it can be found that apparent unphysical oscillations of the heat flux occur for the two schemes, demonstrating the oscillatory convergence behavior for AUSM+ scheme. The computed heat flux distributions on the centerline by different schemes are given in Fig. 6, and it is seen that all the five AUSM family schemes: AUSMDV, AUSM+, AUSMPW, AUSMPW+, AUSM+UP used in this paper give wall heat flux results with unphysical oscillations, while the oscillations by the AUSM+ and AUSMPW schemes are relatively low among the five AUSM family schemes. Additionally, the computed heat flux by the Roe scheme is lower than the experimental data using the current mesh, which can be found in Fig. 6. Fig. 7 illustrates the pressure contours on the symmetric face given by the five AUSM family schemes, in order to show how the shock wave is captured. It is clearly seen that oscillating results arise for all the five AUSM schemes, which further influences the prediction of the wall heat flux and induce unphysical heat flux distributions. The wall heat flux computed by the Roe scheme is shown in Fig. 8. The distribution of the wall heat flux looks very smooth as shown in Fig. 8(a) and meanwhile Fig. 8(b) indicates a fairly stable resolution of the shock wave by the Roe scheme. Moreover, it is also seen from Fig. 6 that the computed heat flux values on the shoulder are even higher than the value on the stagnation point, showing a different pattern from that given by the measured data.
(b) Generated structured mesh
Fig. 1. 3D MSL capsule model.
435
Z.-X. Gao et al. / International Journal of Heat and Mass Transfer 116 (2018) 432–444
100
104
10-1 103
Qw(kw/m2)
Residual
10-2 10-3
102
10-4 101 10-5 10-6
0
2000
4000
6000
8000
100
10000
0
2000
4000
6000
8000
Iteration
Iteration
(a) Residual
(b) Stagnation-point heat flux
10000
100
104
10-1
103
Qw(kw/m2)
Residual
Fig. 3. Convergence histories by Roe scheme for MSL capsule case.
10-2
101
10-3
10-4
102
0
5000
10000
15000
20000
25000
100
0
5000
10000
15000
20000
Iteration
Iteration
(a) Residual
(b) Stagnation-point heat flux
25000
Fig. 4. Convergence histories by AUSM+ scheme for MSL capsule case.
2
2
Qw(kw/m )
Q w(kw/m ) 0
18
36
54
71
89
107
125
(a) AUSM+ scheme
0
20
39
59
78
98
117
137
(b) AUSMPW scheme
Fig. 5. Computational wall heat flux contours of AUSM+ and AUSMPW schemes.
In addition, observing the heat flux values around the stagnationpoint by the AUSM+ and AUSMPW schemes in Fig. 6 and comparing those with the measured stagnation-point heat flux value, it is found that the AUSM+ and AUSMPW schemes predict more accu-
rate stagnation-point heat flux, though unphysical oscillations exist in the distributions of the wall heat flux. The above results clearly reflect the features of the Roe and AUSM family schemes on the hypersonic aeroheating computation,
436
Z.-X. Gao et al. / International Journal of Heat and Mass Transfer 116 (2018) 432–444
300
tor, which refers to the method adopted in the construction of the hybrid high-order schemes [24], and the other one is using the distance from the wall to control the switch between the two kinds of schemes. Next, the accuracy of the aeroheating prediction of the five AUSM schemes is examined at first using a 2D cylinder case, in order to select an appropriate AUSM scheme for the construction of the hybrid scheme. Then, new hybrid schemes will be developed based on the above two techniques (mentioned above), and verified using the 2D cylinder case.
Roe
250
AUSMDV AUSM+ AUSMPW AUSMPW+
Qw(kw/m2)
200
AUSM+UP Exp
150 100
4.1. Comparative study on the aeroheating prediction accuracy of AUSM family schemes and Roe scheme
50 0
-1
-0.5
0
0.5
1
x/R Fig. 6. Comparison of the computed wall heat flux on the centerline with the experimental data.
which have been described in Section 1. The Roe scheme could ensure a stable resolution of the shock wave without strict restriction on the grid generation around the shock wave, and thus the computed wall heat flux shows a smooth distribution. However, the obtained wall heat flux results deviate from the measured data significantly. While for the AUSM family schemes, though they could predict more accurate stagnation-point wall heat flux values, unphysical oscillations exist among the distributions of the computed wall heat flux because of the unstable shock wave resolutions. 4. Development of new hybrid schemes based on the Roe and AUSM schemes Based on the characteristics of the Roe and AUSM family schemes, hybrid schemes are proposed here which would apply the Roe scheme only near the shock wave to ensure the stability of the strong shock wave resolution, while adopt the AUSM family schemes inside the boundary layer to improve the computational accuracy of the wall heat flux. To achieve the above aims, two different ideas are proposed: the first one is based on the shock detec-
P/pa
P/pa
22000 19810 17620 15430 13240 11050 8860 6670 4480 2290 100
(a)AUSMDV
A 2D cylinder experimental case by Holden et al. [25] is selected to perform the comparative study on the aeroheating prediction accuracy of the AUSM family schemes. The radius of the cylinder R is 38.1 mm. Two sets of experimental data are chosen in this study, and the test conditions [25,26] are shown in Table 1. A structured mesh with 101 91 grid points is generated for this cylinder case and shown in Fig. 9(a). The stretching factor near the wall is 1.1 in the normal direction to the wall and the distribution of the grid points in the streamwise direction is uniform. In the same way, the grid points are clustered near the shock wave in the normal direction without any adjustment for the mesh to be aligned with the shock wave. A grid convergence study is performed using three meshes with the first gird spacing off the wall 1 105 m, 1 106 m and 2 107 m, respectively, and the wall heat flux results for AUSM+ scheme and Roe scheme are shown in Fig. 9 (b) and (c). It is seen that the AUSM+ scheme could give a gridconverged and reliable heat flux result using the mesh with the first grid point spacing off the wall 1 106 m, while the result by the Roe scheme agrees with the experimental data until the first grid point spacing is refined to 2 107 m, implying excess dissipation for the Roe scheme. Based on the grid convergence study for the AUSM+ scheme, the mesh with the first grid spacing off the wall 1 106 m is selected to conduct the following simulations. Moreover, the CFL number is taken as 10 for the simulations of the 2D cylinder case, and the convergence histories of the residual and stagnation-point heat flux for Roe and AUSM+ schemes are illustrated in Figs. 10 and 11. Similar to above 3D capsule case, the residual drops more than 6 orders of magnitude for Roe scheme
P/pa
22000 19810 17620 15430 13240 11050 8860 6670 4480 2290 100
(b)AUSM+
P/pa
22000 19810 17620 15430 13240 11050 8860 6670 4480 2290 100
(c)AUSMPW
P/pa
22000 19810 17620 15430 13240 11050 8860 6670 4480 2290 100
(d)AUSMPW+
Fig. 7. Pressure contours on the symmetric face computed by the five AUSM schemes.
22000 19810 17620 15430 13240 11050 8860 6670 4480 2290 100
(e)AUSM+-UP
437
Z.-X. Gao et al. / International Journal of Heat and Mass Transfer 116 (2018) 432–444
2
Qw(kw/m )
0
13
25
38
51
64
76
89
P/pa
22000 19810 17620 15430 13240 11050 8860 6670 4480 2290 100
(a) Wall heat flux contours
(b) Pressure contours
Fig. 8. Computational results by Roe scheme.
Table 1 Free-stream conditions for the cylinder case. Ma1
p1/pa
T1/K
Tw/K
Re/m
8 16.34
855 82.95
125.07 52
294 294.4
4.72 106 3.74 106
vided in Ref. [25] is a little dispersed, it can be still seen that the heat flux distribution computed by the AUSM+ scheme has the best agreement with the measured data. Therefore, the AUSM+ scheme is chosen among the AUSM family schemes for the hybrid scheme together with the Roe scheme. In addition, to examine the influence of the tunable parameter d in the entropy correction of Roe scheme, Fig. 13 gives the computed heat flux distributions by Roe scheme with different values of d. It is seen that a small d could get higher wall heat flux, and with d = 0 the Roe scheme can also give accurate heat flux. However, as mentioned in the Introduction, d needs to be taken as a relative large value for 3D complex flows with strong shock waves, and thus the prediction accuracy of wall heat flux for Roe scheme would be reduced. It is worth noting that for this 2D cylinder case, unphysical oscillations do not occur in the wall heat flux results computed by the AUSM family schemes even when the shape of the mesh is not adjusted to align with the shock wave. It indicates that not all the hypersonic aeroheating simulations using a mesh whose shape is not consistent with the shock wave would generate unphysical oscillations on the wall heat flux for the AUSM family schemes. However, it is likely to result in unphysical oscillations
while only 4 orders for AUSM+ scheme. Both schemes results in a converged wall heat flux prediction. Fig. 12(a) shows the computed heat flux distributions on the wall by different AUSM schemes under the condition of Ma1 = 8. Comparing different numerical heat flux results with the experimental data, it can be found that the heat flux distribution by the AUSM+ scheme agrees with the experimental data best, while the stagnation-point heat flux values by other AUSM schemes are higher than the measured data. Moreover, Fig. 12(a) also presents the result by the Roe scheme using the same structured grid, and the computed heat flux is lower than the experimental data, which again verifies the lower aeroheating prediction accuracy of the Roe scheme. In Fig. 12(b), the wall heat flux under the condition of Ma1 = 16.34 are presented. Although the experimental data pro-
7
7
AUSM+ 1e-5m 1e-6m 2e-7m Exp
5 4 3 2 1 0
Roe
6
Qw/(105W m-2)
Qw/(105W m-2)
6
1e-5m 1e-6m 2e-7m Exp
5 4 3 2 1
0
15
30
45
60
75
90
0
0
15
30
/( )
45
60
75
/( )
(a) Structured
(b) Grid convergence study with
(c) Grid convergence study with
mesh
Ma∞=16.34 for AUSM+
Ma∞=16.34 for Roe
Fig. 9. 2D structural mesh and grid convergence study of AUSM + scheme for the cylinder.
90
438
Z.-X. Gao et al. / International Journal of Heat and Mass Transfer 116 (2018) 432–444
10
1
10
3
100 -1
Qw(105W/m2)
Residual
10
10-2 10-3 10
-4
102
101
10-5 10-6
0
2000
4000
6000
8000
10000
100
12000
0
2000
4000
6000
8000
10000
Iteration
Iteration
(a) Residual
(b) Stagnation-point heat flux
12000
Fig. 10. Convergence histories by Roe scheme for cylinder with Ma1 = 16.34
0
10
-1
10
-2
10
-3
10
Qw(105W/m2)
Residual
10
3
102
101
10-4
10
-5
0
5000
10000
15000
20000
10
25000
0
0
5000
10000
15000
20000
Iteration
Iteration
(a) Residual
(b) Stagnation-point heat flux
25000
Fig. 11. Convergence histories by AUSM+ scheme for cylinder with Ma1 = 16.34.
7
7
6
6 Roe AUSMDV AUSM+ AUSMPW AUSMPW+ AUSM+UP Exp
4 3
4 3
2
2
1
1
0
10
20
30
40
50
60
70
80
ROE AUSMDV AUSM+ AUSMPW AUSMPW+ AUSM+UP Exp
5
Qw/(105W m-2)
Qw/(105W m-2)
5
90
0
15
30
45
/( )
/( )
(a) Ma = 8
( b) Ma =16.34
Fig. 12. Wall heat flux distributions computed by different schemes.
60
75
90
439
Z.-X. Gao et al. / International Journal of Heat and Mass Transfer 116 (2018) 432–444
7
on the concept of smooth measurement factor (SMF) originally proposed in the construction of the WENO scheme [31]. For a grid point j, three stencils, each of which includes three grid points, are designed using the four adjacent grid points on both sides of j, as shown in Fig. 14. Each stencil has its own value of SMF. ¼ p=ðq1 V 21 Þ is employed to idenThe dimensionless pressure p tify the shock wave location. After mathematic deductions [30], the expressions for the SMFs on the three stencils for grid point j are
6 =0 =0.06 =0.125 Exp
Qw/(105W m-2)
5 4
8 j2 4p j1 þ 3p j Þ2 þ 13 j2 2p j1 þ p j Þ2 > SMF 1 j ¼ 14 ðp ðp > 12 < j1 p jþ1 Þ2 þ 13 j1 2p j þ p jþ1 Þ2 SMF 2j ¼ 14 ðp ðp 12 > > : j 4p jþ1 þ p jþ2 Þ2 þ 13 j 2p jþ1 þ p jþ2 Þ2 SMF 3j ¼ 14 ð3p ðp 12
3 2 1 0
dy=1e-6 m
0
15
30
60
45
75
90
/( ) Fig. 13. Wall heat flux distributions computed by Roe schemes with different d values.
ð5Þ
Generally, the SMFs defined by the above equations would approach to a small value in the smooth region, while they would be much larger in the discontinuous region. In present study, the SMF is introduced only to detect the discontinuity of the flow field, so it is reasonable to choose the maximum one among the three SMFs for point grid j. Furthermore, for 3D simulations, Eq. (5) can still be adopted to compute the three SMFs in each dimension, and then, for each SMF, the following treatment is adopted:
SMF m ¼
qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi SMF 2mi þ SMF 2mj þ SMF 2mk
ð6Þ
where the maximum value among the three SMFs is chosen for grid point (i, j, k). Then a damping function fD, whose value would be 0 near the shock wave but 1 in other smooth region, is designed based on the SMFs above. The hyperbolic tangent function is introduced by Zhang et al. [30] to define the damping function as follow:
f D ¼ 1 tanh ½cam maxðSMF 1 ; SMF 2 ; SMF 3 Þ
Fig. 14. Stencil selection for the construction of SMF at grid point j.
using the AUSM family schemes for 3D complex hypersonic flows with strong shock waves. 4.2. Construction of new hybrid schemes (1) Hybrid scheme based on the shock wave detector We evaluated the shock wave detectors proposed by Visbal and Gaitonde [24], Hill and Pullin [27], Larsson et al. [28], Priozzoli [29] and Zhang et al. [30], and eventually the one established by Zhang et al. [30] is selected. This shock wave detector proposed is based
Ma
where the cam is an amplification factor for the discontinuity in flow field, which could make the hyperbolic tangent value of the SMFs near the shock wave close to 1. According to the theoretical analysis by Zhang et al. [30], the value of cam is taken as 95. Here, the 2D cylinder is employed to examine the effectiveness of the above shock wave detector. Fig. 15 shows the results of the shock wave detector in the flowfield around the 2D cylinder with Ma1 = 8 and 16.34, where the top half is the Mach number contours, and the bottom half is the fD contours. It can be seen that the shock wave detector achieves a high-fidelity capture of the
Ma
8.0 6.9
16.0
13.7 11.4 9.1
5.7 4.9 3.8
6.9 4.6 2.3
2.8 1.6 0.0
0.0
fD
fD 1.0
1.0
0.9 0.7
0.8 0.6 0.5
0.6 0.4 0.3
0.3 0.2 0.1
0.1 0.0
(a) Ma =8
ð7Þ
0.0
(b) Ma =16.34
Fig. 15. Contours of the damping function and Mach number for 2D cylinder.
Z.-X. Gao et al. / International Journal of Heat and Mass Transfer 116 (2018) 432–444
6
7
5
6 Hybrid-1st Hybrid-2ed- 0 Hybrid-2nd-2 0 Hybrid-2nd-5 0 Hybrid-2nd-10 0 Hybrid-2nd-0.2 0 Hybrid-2nd-0.5 0 AUSM+ Exp
4 3 2
4 3 2
1 0
Hybrid-1st Hybrid-2nd- 0 Hybrid-2nd-2 0 Hybrid-2nd-5 0 Hybrid-2nd-10 0 Hybrid-2nd-0.2 0 Hybrid-2nd-0.5 0 AUSM+ Exp
5
Qw/(105W m-2)
Qw/(105W m-2)
440
1 0
10
20
30
40
50
60
70
80
0
90
0
20
40
/( )
60
80
/( )
(a) Ma =8
(b) Ma =16.34
Fig. 16. Wall heat flux distributions computed by the hybrid schemes.
10
1
10
3
100 -1
Qw(105W/m2)
Residual
10
10-2 10-3 10
-4
102
101
10-5 10-6
0
3000
6000
9000
12000
100
0
3000
6000
9000
Iteration
Iteration
(a) Residual
(b) Stagnation-point heat flux
12000
Fig. 17. Convergence histories by the second hybrid scheme for cylinder with Ma1 = 16.34.
shock wave and the fD value is close to 0 near the shock wave but 1 in other region. Finally, the numerical flux of the hybrid scheme is constructed based on the above shock wave detector:
e F iþ1 ¼ ð1 f D Þ e F AUSMþ F Roe þ f De iþ1 iþ1 2
2
2
ð8Þ
(2) Hybrid scheme based on the distance from the wall The above hybrid scheme based on the shock wave detector applies to the Roe scheme around the shock wave and the AUSM + scheme in other region, as can be seen in Fig. 15. Yet another way is to construct the hybrid scheme is to limit the utilization of the AUSM+ scheme only within the boundary layer and adopt the Roe scheme in other region including the shock wave region. This can be realized by designing a control function based on the distance from the wall. Here, a simple function form as follow is employed:
Hðd d0 Þ ¼
1 d > d0 0
d < d0
ð9Þ
where d0, an empirical parameter, would control the switch from the AUSM+ scheme to the Roe scheme. Using this function, the numerical flux of this hybrid scheme could be given as
e F iþ1 ¼ Hðd d0 Þ e F Roe þ ½1 Hðd d0 Þ e F AUSMþ iþ1 iþ1 2
2
2
ð10Þ
As for the parameter d0, its value should be a reflection of the boundary layer thickness around the vehicle. Because there is no available expression for the boundary layer thickness in hypersonic laminar flows, the following formula is used to calculate the boundary layer thickness for the incompressible flat-plate flow to give an estimation for d0:
pffiffiffiffiffiffiffiffi d0 ¼ L= ReL
ð11Þ
where ReL is the Reynolds number computed using the free-stream parameters and the characteristic length. Eq. (11) may not accurately estimate the boundary layer thickness for a hypersonic flow over blunt reentry vehicles. However, if the aeroheating computation using this hybrid scheme is not sensitive to the value of d0, which will be examined in the following numerical experiments, the estimation in Eq. (11) for d0 could be acceptable.
441
Z.-X. Gao et al. / International Journal of Heat and Mass Transfer 116 (2018) 432–444
Qw(Kw/m2)
0
23
46
69
91
114
137
160
P/pa
22000 19810 17620 15430 13240 11050 8860 6670 4480 2290 100
(a) wall heat flux contour
(b) pressure contour
Fig. 18. Computational results by the first hybrid scheme.
Qw(Kw/m2)
0
17
34
51
(a)
67
94
118
Qw(Kw/m2)
135
0
19
39
58
(b) 0.5
0
77
96
116
135
0
Fig. 19. Computational results by the second hybrid scheme with different distances from the wall.
160
=0
140
Hybrid-2nd-0.5
120
Hybrid-2ndRoe
0
0
Qw(kw/m2)
Exp
100 80 60 40 20 0
-1
-0.5
0
0.5
1
x/R Fig. 20. Comparison of the wall heat flux distributions on the centerline with experimental data.
4.3. Application of the hybrid schemes in the 2D cylinder case The 2D cylinder case is used to examine the effectiveness of the above two hybrid schemes. In the following numerical experi-
ments, the hybrid scheme using the shock wave detector is noted as the first hybrid scheme and the hybrid scheme based on the distance from the wall is called the second hybrid scheme. The computed heat flux distributions by the two hybrid schemes are shown in Fig. 16 (a: Ma1 = 8; b: Ma1 = 16.34), and the results obtained by the AUSM+ scheme are also given in Fig. 16. It can be seen that the heat flux results by the hybrid schemes agree well with the experimental data for the two cases. In addition, the results obtained by the hybrid scheme are very close to those using the AUSM+ scheme independently, which indicates that, for this 2D cylinder case, in spite of the introduction of numerical dissipation by the Roe scheme into the flowfield, the adoption of the AUSM+ scheme in the boundary layer could still maintain almost the same aeroheating prediction accuracy as applying the AUSM+ scheme in the whole flowfield. Furthermore, Fig. 16 also depicts the results for the hybrid scheme based on the distance from the wall with different values of d0, where the baseline d0 value is computed using Eq. (11). It is found that there is little influence on the wall heat flux prediction when d0 changes from 0.2 times of the baseline value to 10 times, which proves that the aeroheating simulation using the second hybrid scheme is not sensitive to the value of d0. Moreover, the convergence histories of the residual and stagnation-point heat flux for the second hybrid scheme are shown in Fig. 17. Compared to those for the AUSM+ scheme in Fig. 11, the hybrid scheme obviously improves the convergence characteristic
442
Z.-X. Gao et al. / International Journal of Heat and Mass Transfer 116 (2018) 432–444
and meanwhile could get a converged stagnation-point heat flux with almost same iteration numbers. For this 2D cylinder case, there is no unphysical oscillation existing in the wall heat flux distributions even solely applying the AUSM+ scheme, so this case indicates that the hybrid schemes could maintain a close computational accuracy of the wall heat flux compared to the AUSM+ scheme. However, whether the two hybrid schemes could eliminate the unphysical oscillations of the wall heat flux distributions induced by the unstable resolution of the shock wave should be further verified.
10
5. Aeroheating simulations on 3D MSL capsule using the hybrid schemes Here, the two hybrid schemes are adopted to simulate the 3D MSL capsule case in Section 3. Fig. 18 illustrates the results obtained by the hybrid scheme based on the shock wave detector. From the Fig. 18(a), it is seen that the obtained wall heat flux distribution still shows apparent unphysical oscillations. It indicates that the hybrid scheme based on the shock wave detector could not improve the stability of the shock wave resolution, which
0
4
10
10-1
103
10
Qw(kw/m2)
Residual
10
-2
-3
2
10
10-4
101 10-5 10
-6
0
0
5000
10000
15000
10
0
5000
10000
Iteration
Iteration
(a) Residual
(b) Stagnation-point heat flux
Fig. 21. Convergence histories by the second hybrid scheme for MSL capsule case.
Qw(kw/m2)
0
19
37
56
74
93
111
130
Qw(kw/m2)
0
19
39
(a) =4° Qw(kw/m2)
0
19
39
58
77
(c) =12°
58
77
96
116
135
114
137
160
(b) =8°
96
116
135
Qw(kw/m2)
0
23
46
69
91
(d) =24°
Fig. 22. Wall heat flux contours using the second hybrid scheme under different angles of attack.
15000
443
Z.-X. Gao et al. / International Journal of Heat and Mass Transfer 116 (2018) 432–444
can also be found from Fig. 18(b). One possible reason is that this hybrid scheme only activates the Roe scheme around the shock wave while use the AUSM+ scheme in other region, and meanwhile the switch between these two schemes occurs near the shock wave, which will not eliminate the oscillating solutions around the shock wave. Then, the hybrid scheme based on the distance from the wall is adopted to perform the numerical simulations on the 3D MSL capsule model, and the parameter d0 is estimated according to Eq. (11), whose value is approximately 1.35 104 m. Fig. 19(a) shows the wall heat flux contours computed by this hybrid scheme, and it can be seen that the heat flux oscillations are totally eliminated, proving the effectiveness of this hybrid scheme. In addition, Fig. 19(b) illustrates the wall heat flux result with the same hybrid scheme but the empirical parameter is changed to 0.5d0. The difference between Fig. 19(a) and (b) is limited, indicating that the wall heat flux prediction of this hybrid scheme is not sensitive to d0 in Eq. (10) for 3D case. Fig. 20 compares the wall heat flux distributions on the centerline computed by this hybrid scheme with the measured data, and good agreements are achieved for both d0 and 0.5d0. Compared to the results by solely using the Roe or AUSM+ schemes, significant improvements on the aeroheating prediction are observed by this hybrid scheme. The convergence histories for the MSL capsule case using the second hybrid scheme are given
160
in Fig. 21. By comparing Figs. 21 and 4, we demonstrate that the second hybrid scheme improves the convergence characteristic relative to the AUSM+ scheme. In addition to the case with a = 0°, a series of aeroheating simulations with a = 4°, 8°, 12° and 24° are conducted applying the hybrid scheme based on the distance from the wall. In these cases, the shape of the shock wave would change corresponding to different angles of attack, while keeping the same computational mesh with the one used for a = 0° case. The computed wall heat flux contours are shown in Fig. 22 for the four cases, and it is found that this hybrid scheme can maintain smooth distributions of the wall heat flux without any unphysical oscillation for all the four cases with different angles of attack. Fig. 23 compares the wall heat flux distributions on the centerline for the four cases with the corresponding experimental data, and reasonable agreements can be seen for all the four cases. The above applications of the hybrid schemes in the MSL capsule case demonstrate that the hybrid scheme based on the distance from the wall could effectively eliminate the unphysical oscillations in the distributions of the wall heat flux which are generated when solely using the AUSM+ scheme, and enhance the aeroheating accuracy compared to that of the Roe scheme. Moreover, for different freestream conditions, this hybrid scheme can also guarantee the accuracy of the wall heat flux prediction using
160
=4
140
100 80 60
80 60 40
20
20 -0.5
160
0
0.5
0 -1
1
0
x/R
(a) =4°
(b) =8° 160
=12
0.5
1
=24
140 Hybrid-2nd Exp
120
Qw(kw/m2)
80 60
100 80 60
40
40
20
20 -1
-0.5
0
0.5
Hybrid-2nd Exp
120
100
0
-0.5
x/R
140
Qw(kw/m2)
100
40
0 -1
Hybrid-2nd Exp
120
Qw(kw/m2)
Qw(kw/m2)
140
Hybrid-2nd Exp
120
=8
1
0 -1
-0.5
0
x/R
x/R
(c) =12°
(d) =24°
0.5
Fig. 23. Comparisons of the wall heat flux distributions on the centerline with experimental data under different angles of attack.
1
444
Z.-X. Gao et al. / International Journal of Heat and Mass Transfer 116 (2018) 432–444
the same mesh. The above advantages make the proposed hybrid scheme appropriate for the aeroheating prediction of the hypersonic reentry vehicles. 6. Conclusion Hypersonic aeroheating computations for a 3D capsule from the Mars Science Laboratory (MSL) are performed using the Roe and AUSM family schemes, and it is found that the AUSM family schemes may generate unphysical oscillations in the wall heat flux distributions due to their unstable resolutions of the shock wave, while the wall heat flux prediction accuracy by the Roe scheme is not adequate though the obtained heat flux distribution is smooth. In order to develop a numerical scheme that could capture the shock wave stably and meanwhile maintain a high-fidelity simulation of the wall heat flux, new hybrid schemes combining the Roe and AUSM+ schemes are developed in the present paper. Two different ideas are proposed to construct two hybrid schemes. The first hybrid scheme is to use the shock wave detector to detect the location of the shock wave, which would activate the Roe scheme near the shock wave and apply the AUSM+ scheme in other region. The second one is to utilize the distance from the wall to control the switch between the two schemes, and apply AUSM+ only in the boundary layer while the Roe scheme in other region. The applications of these two hybrid schemes for the 3D MSL capsule case show that the first hybrid scheme based on the shock wave detector does not show improvement in the unphysical oscillations of the wall heat flux, while the second hybrid scheme removes the unphysical oscillations and maintains accurate simulation of the wall heat flux. Additional computations of the 3D MSL capsule case with different angles of attack further demonstrate that the proposed hybrid scheme based on the distance from the wall could effectively eliminate the unphysical oscillations in the distributions of the wall heat flux which are generated when solely using the AUSM+ scheme, and remarkably enhance the aeroheating accuracy compared to that of the Roe scheme. The proposed hybrid scheme based on the distance from the wall is shown to be appropriate for the aeroheating prediction of the hypersonic reentry vehicles. Acknowledgments This work was supported by the National Natural Science Foundation of China through the Grant 91641123 and the Academic Excellence Foundation of BUAA for PhD students. Conflict of interest All the authors (Zhenxun Gao, Haichao Xue, Zhichao Zhang, Hongpeng Liu, Chun-Hian Lee) of the article ‘‘A Hybrid Numerical Scheme for Aeroheating Computation of Hypersonic Reentry Vehicles” state that there is no financial/personal interest or belief that could affect our objectivity, and all authors have seen and approved the final version of the manuscript being submitted. We warrant that the article is the authors’ original work. References [1] K.A. Hoffmann, M.S. Siddiqui, S.T. Chiang, Difficulties associated with the heat flux computations of high speed flows by the Navier-Stokes equations, in: AIAA paper, 1991, AIAA-1991-0457.
[2] G.H. Klopfer, H.C. Yee, Viscous hypersonic shock-on-shock interaction on blunt cowl lips, in: AIAA paper, 1988, AIAA-1988-0233. [3] F.G. Blottner, D.E. Larson, Navier-Stokes Code NS3D for Blunt Bodies, Sandia Report, 1988, SAND-1988-050411 UC-32. [4] M.S. Siddiqui, K.A. Hoffmann, A comparative study of the Navier-Stokes solvers with emphasis on the heat transfer computations of high speed flows, in: AIAA paper, 1992, AIAA-1992-0835. [5] P. Papadopoulos, E. Venkatapaphy, D. Prabhu, Current grid generation strategies and future requirements in hypersonic vehicle design, analysis and testing, Appl. Math. Model. 23 (9) (1999) 705–735, http://dx.doi.org/10.1016/ S0307-904X(99)00007-4. [6] I.S. Men’ shov, Y. Nakamura, Numerical simulations and experimental comparison for high-speed nonequilibrium air flows, Fluid Dyn. Res. 27 (5) (2000) 305–334, http://dx.doi.org/10.1016/S0169-5983(00)00010-1. [7] S.J. Henderson, J.A. Menart, Grid Study on Blunt Bodies with the Carbuncle Phenomenon, in: AIAA paper, 2007, AIAA-2007-3904. [8] K. Kitamura, P. Roe, F. Ismail, An evaluation of Euler fluxes for hypersonic flow computations, AIAA J. 47 (1) (2009) 44–53, http://dx.doi.org/10.2514/1.33735. [9] K. Kitamura, E. Shima, Y. Nakamura, P. Roe, Evaluation of Euler fluxes for hypersonic heating computations, AIAA J. 48 (4) (2010) 763–776, http://dx.doi. org/10.2514/1.41605. [10] J.H. Lee, O.H. Rho, Numerical analysis of hypersonic viscous flow around a blunt body using Roe’s FDs and AUSM + schemes, in: AIAA paper, 1997, AIAA1997-2054. [11] J.H. Lee, O.H. Rho, Accuracy of AUSM + scheme in hypersonic blunt body flow calculation, in: AIAA paper,1998, AIAA-1998-1538. [12] Y. Wada, M.-S. Liou, An accurate and robust flux splitting scheme for shock and contact discontinuities, SIAM J. Sci. Stat. Comput. 18 (1997) 633–657, http:// dx.doi.org/10.1137/S1064827595287626. [13] M.-S. Liou, A Sequel to AUSM: AUSM+, J. Comput. Phys. 129(2) (1996) 364– 382. 10.1006/jcph.1996.0256. [14] M.-S. Liou, Ten Years in the Making—AUSM Family, in: AIAA paper, 2001, AIAA2001-2521. [15] S.S. Kim, C. Kim, O.H. Rho, S.K. Hong, Methods for the accurate computations of hypersonic flows I. AUSMPW + Scheme, J. Comput. Phys. 174 (1) (2001) 38–80, http://dx.doi.org/10.1006/jcph.2001.6873. [16] M.-S. Liou, A Further Development of the AUSM + Scheme Towards Robust and Accurate Solutions for All Speeds, in: AIAA Paper, 2003, AIAA-2003-4116. [17] P.L. Roe, Approximate Riemann solvers, parameter vectors and difference schemes, J. Comput. Phys. 43 (1981) 357–372, http://dx.doi.org/10.1016/00219991(81)90128-5. [18] Z.X. Gao, C.W. Jiang, C.H. Lee, Improvement and application of wall function boundary condition for high-speed compressible flows, Sci. China Technol. Sci. 56 (10) (2013) 2501–2515, http://dx.doi.org/10.1007/s11431-013-5349-4. [19] Z.X. Gao, C.H. Lee, A numerical study of turbulent combustion characteristics in a combustion chamber of a scramjet engine, Sci. China Technol. Sci. 53 (8) (2013) 2111–2121, http://dx.doi.org/10.1007/s11431-010-3088-3. [20] Z.X. Gao, C.W. Jiang, S.W. Pan, C.H. Lee, Combustion Heat-release effects on supersonic compressible turbulent boundary layers, AIAA J. 53 (7) (2015) 1949–1968, http://dx.doi.org/10.2514/1.J053585. [21] H.C. Yee, G.H. Kolpfer, J.L. Montague, High-resolution shock capturing schemes for inviscid and viscous hypersonic flows, 1988, NASA TM 100097. [22] R.H. Brian, Turbulent aeroheating testing of mars science laboratory entry vehicle, J. Spacecraft Rockers 45 (3) (2008) 417–427, http://dx.doi.org/ 10.2514/1.31798. [23] R.H. Brian, Turbulent Aeroheating Testing of Mars Science Laboratory Entry Vehicle in Perfect-Gas Nitrogen, in: AIAA Paper, 2007, AIAA-2007-1208. [24] M.R. Visbal, D.V. Gaitonde, Shock capturing using compact differencing based methods, in: AIAA Paper, 2005, AIAA-2005-1265. [25] M.S. Holden, A.R. Wieting, J.R. Mpselle, Studies of aerothermal loads generated in regions of shock/shock interaction in hypersonic flow, in: AIAA Paper, 1988, AIAA-1988-0477. [26] K.M. Reery, S.T. Imlay, Blunt-Body Flow Simulations, in: AIAA Paper, 1988, AIAA-1988-2904. [27] D.J. Hill, D.I. Pullin, Hybrid tuned center-difference-WENO method for large eddy simulations in the presence of strong shocks, J. Comput. Phys. 194 (2004) 435–450, http://dx.doi.org/10.1016/j.jcp.2003.07.032. [28] J. Larsson, S.K. Lele, P. Moin, Effect of numerical dissipation on the predicted spectra for compressible turbulence, Center for turbulence research, annual research briefs (2007) 47–57. [29] S. Priozzoli, Numerical methods for high-speed flows, Ann. Rev. Fluid Mech. 43 (2011) 163–194, http://dx.doi.org/10.1146/annurev-fluid-122109-160718. [30] Z.C. Zhang, Z.X. Gao, C.W. Jiang, C.H. Lee, A RANS model correction on unphysical over-prediction of turbulent quantities across shockwave, Int. J. Heat Mass Transf. 106 (2017) 1107–1119, http://dx.doi.org/10.1016/j. ijheatmasstransfer.2016.10.087. [31] G.S. Jiang, C.W. Shu, Efficient implementation of weighted ENO schemes, J. Comput. Phys. 126 (130) (1996) 202–228, http://dx.doi.org/10.1006/ jcph.1996.0130.