JID:AESCTE AID:4278 /FLA
[m5G; v1.224; Prn:9/11/2017; 8:35] P.1 (1-9)
Aerospace Science and Technology ••• (••••) •••–•••
1
Contents lists available at ScienceDirect
2
67 68
3
Aerospace Science and Technology
4 5
69 70 71
6
72
www.elsevier.com/locate/aescte
7
73
8
74
9
75
10
76
11 12 13 14 15 16 17
Assessment of a hybrid RANS/LES simulation method and URANS method in depicting the unsteady motions of flow structures in a scramjet combustor Yiyu Han ∗ , Yuanyuan He, Ye Tian, Fuyu Zhong, Jialing Le
18 19
24 25 26 27 28 29 30 31 32 33
79 80 81 82 83 85 86
21 23
78
84
Science and Technology on Scramjet Laboratory of Hypervelocity Aerodynamics Institute, CARDC, China
20 22
77
87
a r t i c l e
i n f o
Article history: Received 26 September 2017 Received in revised form 25 October 2017 Accepted 1 November 2017 Available online xxxx Keywords: Scramjet combustor Unsteady motion Hybrid RANS/LES simulation Zonal detached eddy simulation (ZDES) Unsteady Reynolds averaged Navier–Stokes (URANS)
34 35 36 37
a b s t r a c t
88
A scramjet combustor is simulated using both a hybrid Reynolds averaged Navier–Stokes (RANS)/large eddy simulation (LES) method, namely zonal detached eddy simulation (ZDES), and unsteady RANS (URANS) method to conduct an assessment of their abilities in depicting the unsteady motions of flow structures in the scramjet combustor. Although ZDES and URANS present a similar ability with regard to calculating the mean variables, they present distinct differences concerning the ability to capture the unsteady motions of the flow structures. The ZDES case successfully captures the unsteady motions of the shock structures. The resolved flow is comprised of a series of oscillating cycles that mixes strong cycles and weak cycles, and one or more strong cycles are followed by one or more weak cycles. In strong cycles, the oscillating amplitude of the pressure is large and the duration of cycle is long. The range of the motion of shock structures is large, and a Mach disc is formed when the shock structures move upstream. In weak cycles, the oscillating amplitude of the pressure is small and the duration is short, and no Mach disc appears throughout a cycle. All of these coincide fairly well with the experimental shadowgraph images. In contrast, the results of URANS distinctly disagree with those of the ZDES case and the experiments. © 2017 Elsevier Masson SAS. All rights reserved.
89 90 91 92 93 94 95 96 97 98 99 100 101 102 103
38
104
39
105
40 41
106
1. Introduction
42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58
As a key component of hypersonic air-breathing vehicle, scramjet has gradually drawn a great deal of interest. Researches on scramjets are developed mainly through three ways: flight experiments, shock tunnel experiments and numerical simulations. Flight experiments, such as the HyShot 2 flight experiment [1], provide the most realistic flow data. However, they have two drawbacks. One is the limitation of the high costs and long cycle of a whole fight experiment. The other is that flow data acquired from these experiments is always restricted to pressure distribution, temperature distribution, etc. and no flow structures can be acquired directly from the experiments. Shock tunnel experiments (see for example, [2–4]) can provide some visible images of flow structures thanks to optical measurements such as schlieren, shadowgraph, PIV, PLIF, etc. Nevertheless, these methods also have their own limitations. For example, schlieren images represent the integration of the gradient of density. They are typically used for the description
59 60 61 62 63 64 65 66
*
Corresponding author at: 6 south section, second ring road, Mianyang, Sichuan, 621000, China. E-mail address:
[email protected] (Y. Han). https://doi.org/10.1016/j.ast.2017.11.003 1270-9638/© 2017 Elsevier Masson SAS. All rights reserved.
of shock wave structures and hardly can any other flow details be extracted. Besides, the range of the flow field photographed is restricted by the sight window, so usually only part of the whole flow field can be acquired. By contrast, numerical simulations can provide intuitionistic description of the whole flow field, and thus are beneficial for an in-depth understanding of the flow. Still, with regard to numerical simulations there are also some challenges, of which the most critical one is to select the most appropriate method for a certain type of flow from various numerical simulation methods, in order to reach a compromise between the accuracy and the computational cost. The most widely used numerical simulation method for turbulent flows in practical engineering applications is RANS method, owing to its low cost and its ability to yield satisfying accuracy in attached flow and mildly separated flow. However, it has significant deficiency in predicting massively separated flow, which may exist in scramjet combustors. A more questionable challenge for RANS method is its ability to capture the dynamic characteristics of unsteady flow in scramjet combustors since even URANS method may yield steady results for many unsteady flows [5]. LES method resolves the large-scale turbulent structures and only models the smallest scales, and thus is supposed to yield higher fidelity results for such flows. However, enormous computational resources are
107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132
JID:AESCTE
AID:4278 /FLA
[m5G; v1.224; Prn:9/11/2017; 8:35] P.2 (1-9)
Y. Han et al. / Aerospace Science and Technology ••• (••••) •••–•••
2
1 2
67
Nomenclature
68
3 4 5 6 7 8
69
RANS URANS LES ZDES MUSCL
9 10 11
LAM
Reynolds averaged Navier–Stokes unsteady RANS large eddy simulation zonal detached eddy simulation monotonic upstream-centered scheme for conservation laws hybrid optimized low-dissipation and adaptive MUSCL scheme
PIV
particle image velocimetry
70
PLIF
planar laser-induced fluorescence
71
LU-SGS
lower-upper symmetric Gauss–Seidel
72
SST
shear stress transport
73
MPLD
optimized low-dissipation scheme with monotonicity
74
preserving procedure PSD
power spectral density
12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66
75 76 77 78
needed to use LES in most wall-bounded flows at high Reynolds number in practical engineering applications due to the excessive computational cost for resolving the turbulent structures in the inner layer of the boundary layer [6–9]. For example, Bisek [10] pointed out that when LES was used to simulate the HIFiRE-6 flow path, the resolution of a grid with as high as 640 million cells was inadequate in some regions of the domain. Hence for these flows, it is preferable to use RANS in wall boundary layers or at least in the inner regions of boundary layers, and use LES elsewhere. This is the principle of hybrid RANS/LES methods. Increasing applications of hybrid RANS/LES methods on scramjet combustors have been reported during the last decade. References [11–22] carried out simulations on scramjet combustors using hybrid RANS/LES methods, and revealed some flow structures that RANS can’t capture, and thus to a certain extent helped to better understand flows in scramjet combustors. Nevertheless, none of these papers presented any analysis on the dynamic evolution of the flow structures such as the unsteady motions of the shock structures and separation bubbles. Choi et al. [23] carried out a numerical analysis of the flow dynamics of a scramjet combustor using detached eddy simulation, a version of hybrid RANS/LES methods. However, their simulations were based on two-dimensional grids which were obviously improper and without experimental data for comparison. MaloMolina et al. [24] calculated an axisymmetric scramjet combustor with uniform inflow condition using LES and showed the dynamic development of small scale structures. However, they computed only a one-eighth circumferential sector to reduce the computational costs which was also questionable for LES, and no experimental data was presented for comparison either. Hence, the credibility of these results is doubtful. Boles et al. [25] conducted a hybrid RANS/LES simulation on a simplified inlet-isolator configuration at unstart condition. Koo et al. [26] calculated the same configuration using LES. Both of the simulations presented distinct deviations from the experiments with regard to quantitative variables: Boles et al. [25] presented larger separation region compared with the experiment and the shock propagation speed from simulations in Koo et al. [26] was 3–4 times larger than that from the experiment. Nevertheless, both of the simulations successfully captured the dynamic characteristics of the upstream propagation of shock-train and separation bubbles, which qualitatively agreed with the experiments. This suggests that hybrid RANS/LES simulation or LES method can help to understand the dynamic characteristics of the flow structures in the internal flow path, at least qualitatively. In this paper, we aim at a further assessment of hybrid RANS/LES simulation and URANS in depicting the unsteady motions of the flow structures in scramjet combustors. This is achieved by simulating a scramjet combustor in which flow structures evolve dynamically and comparing the numerical results with wind tunnel experimental data. The rest of the paper is structured as follows. In Section 2, we briefly introduce the wind tunnel experiments for this study, including the facility and scramjet config-
79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95
Fig. 1. Photo of the CARDC’s supersonic combustion facility.
96 97
uration. In Section 3, numerical methods, including time marching and spatial discretization schemes and inflow turbulent boundary conditions are shown, as well as the URANS and hybrid RANS/LES methods. Then we give the settings of the numerical simulations. In Section 4, results of numerical simulations are presented and discussed in detail, combined with the experimental data. Finally some concluding remarks are presented in Section 5.
98 99 100 101 102 103 104 105
2. Experimental configuration
106 107
The experiments were carried out at China Aerodynamics Research and Development Center (CARDC)’s direct-connect supersonic combustion facility [3,4] as shown in Fig. 1. A H2 /O2 vitiated air heater was used to generate high enthalpy airflow supplied into the combustor. The vitiated air had the mole fraction of 21% O2 , 12% H2 O, and 67% N2 and a mass flow rate of 2.68 kg/s. It entered the nozzle and got accelerated by the nozzle. The entrance of the isolator of the scramjet was directly connected to the exit of the nozzle. The nominal flow conditions at the exit of the nozzle were: Ma = 2.0, P = 104136 Pa and T = 555.3 K. The schematic illustration of the scramjet adopted in this study is presented in Fig. 2. The scramjet was comprised of an isolator and a combustor. The cross-sectional area was 30 × 150 mm2 at the isolator exit. The isolator contained two sections, namely a straight section of 220 mm length and an expansion section of 80 mm length whose upper-wall had a 1.4 degree divergent angle. For clarity, we denote the location of the entrance of the isolator as x = 0. The combustor also contained two sections, namely a cavity which extended from x = 300 mm to x = 418.8 mm and a four part expansion section which extended from x = 418.8 mm to x = 1070 mm. The first part ranged from x = 418.8 mm to x = 523.8 mm, with the wall divergent angle being 1.4 degree. The second part ranged from x = 523.8 mm to x = 668.4 mm, with the wall divergent angle being 2.0 degree. The third part ranged from x = 668.4 mm to x = 810 mm, with the wall divergent an-
108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132
JID:AESCTE AID:4278 /FLA
[m5G; v1.224; Prn:9/11/2017; 8:35] P.3 (1-9)
Y. Han et al. / Aerospace Science and Technology ••• (••••) •••–•••
3
1
67
2
68
3
69
4
70
5
71
6
72
7
73
8
74
9
75
10
76
11
77
12
78
13 15 16 17 18 19 20 21 22 23 24 25
gle being 8.0 degree. The fourth part ranged from x = 810 mm to x = 1080 mm, with the wall divergent angle being 15.0 degree. A fuel-injector was placed inside the cavity, but in this study no fuel was injected into the combustor. A high-speed digital camera was utilized for taking shadowgraph images. The shadowgraph images were collected at 3000 frame/s. The sight window was located at the cavity, as shown in Fig. 2. 3. Numerical methods and settings
26 27
3.1. Numerical methods
28 29 30 31 32 33 34 35 36 37 38 39 40 41
79
Fig. 2. Schematic illustration of the scramjet combustor in this study.
14
The simulations presented in this study have been completed using the in-house CFD code AHL3D. AHL3D solves the compressible Navier Stokes equations using a three-dimensional, parallel, cell-averaged finite volume flow solver on multi-block, structured meshes. It has been validated by numerous applications (see for example, [27,28]). In this study, the URANS method uses Menter’s k–ω SST model [29] for turbulence modeling:
∂ ∂t
=
42 43
ρk ρω
45
+
∂ ∂ xi +
44
∂ ∂ xi
ρ ku i ρωu i
(μ + σk μt ) ∂∂xki (μ + σω μt ) ∂∂ xωi
τi j ∂∂ ux ij
− ρk
3/2
/˜l
γ ρ ∂ ui 2 νt τi j ∂ x j − βω ρω + 2σω2 (1 − F 1 ) ω max(∇ k · ∇ ω, 0)
46
(1)
80 81 82 83 84 85
Fig. 3. The computational grid with every fourth point plotted in each direction.
lRANS if d w < dinterface w lZDES = min(lRANS , lLES ) otherwise
49
Here d w is the distance to the wall of the local cell, dinterface is a w prescribed location at which the RANS mode transits to LES mode. To the authors’ knowledge, there are not any strict criteria for set. Nevertheless, when setting dinterface two principles ting dinterface w w should be considered: on the one hand, as suggested by Deck [7], ‘for a better skin friction prediction, it is shown that the RANS/LES interface should be high enough in the boundary layer’; but on the other hand, in order to resolve the large scale structures in the outer layer of the boundary layer, the RANS/LES interface can’t be too high. Thus, dinterface can be located at the upper edge of w the inner layer, which is about 0.2δ away from the wall. In this paper, dinterface is set to 1 mm, which is close to 20% of the boundw ary layer thickness at x = 220 mm (this thickness is denoted as δ0 hereinafter). The second-order-accurate dual-time stepping method [32] is adopted for time marching, with LU-SGS [33] for implicit approach. LAM scheme [8], which is developed in house, is used for spatial discretization. It is a blend of MPLD scheme [34] and adaptive MUSCL scheme [35] based on a smoothness indicator. It ensures stability in large gradient or strong discontinuous regions while keeps a high-resolution, low-dissipation performance elsewhere.
52 53 54 55 56 57 58 59 60 61 62 63 64 65 66
89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113
Here k is the turbulent kinetic energy and ω is the specific dissipation rate. The eddy viscosity νt can be calculated using k and ω :
50 51
87 88
(4)
47 48
86
νt =
a1 k max(a1 ω, F 2 )
(2)
Detailed meanings of other variables in Equations (1) and (2) can be found in [30] and here we only focus on the turbulent length scale l of the destruction term D k = ρ k3/2 / l in the k-transport equation. In Menter’s original k–ω SST model, ˜l = ˜lRANS = k1/2 /(βk ω). If lRANS is replaced by the length scale ˜lLES = C S G S , the modified SST model can be used as a subgrid-scale model [30]. Hence, a hybrid RANS/LES method can be fabricated by hybridizing the two different length scales using a blending function f blend :
l = f blend lRANS + (1 − f blend ) lLES
(3)
In this paper, we adopt ZDES, a version of hybrid RANS/LES methods proposed by Deck et al. [31]:
3.2. Numerical settings
114 115
To avoid confusion, hereinafter, ZDES simulation and URANS simulation are denoted as the ZDES case and the URANS case, respectively. The computational domain for both the ZDES case and the URANS case extends from x = 220 mm to x = 1072 mm in streamwise direction. In order to reduce the cost, URANS is used instead of ZDES at x > 668 mm for the ZDES case so that the grid in this domain can be coarsened. The grid for the ZDES case consists of 33.2 million cells. The topology of the grid is illustrated in Fig. 3. In spanwise direction, the length of computational domain is 20 mm, with periodic condition at the spanwise boundary. The cell size in this direction is 0.25 mm which is close to 0.05δ0 . The cell size in streamwise direction near the cavity is shown in Fig. 4. From this we can see that x/δ0 is basically around or less than 0.1 near the cavity. Thus, it can be concluded that the cell sizes in streamwise direction and in spanwise direction meet the typical demand of hybrid RANS/LES simulations (see for example, [36]). For the URANS
116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132
JID:AESCTE
AID:4278 /FLA
[m5G; v1.224; Prn:9/11/2017; 8:35] P.4 (1-9)
Y. Han et al. / Aerospace Science and Technology ••• (••••) •••–•••
4
1
67
2
68
3
69
4
70
5
71
6
72
7
73
8
74
9
75
10
76
11
77
12
78
13
79
14 15
80
Fig. 4. The cell sizes in streamwise direction near the cavity, normalized by the boundary layer thickness at x = 220 mm.
81
16
82
17
83
18
84
19
85
20
86
Fig. 6. Mean pressure distribution at the down wall.
21
87
22
88
23
89
24
90
25
91
26
92
27
93
28
94
29
95
30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58
96
Fig. 5. y + at the wall near the cavity.
97 98
case, two-dimensional grid is adopted and it is a slice extracted from the grid of the ZDES case. For the URANS case, the inflow condition is fixed by the profile extracted from the result of a previous steady RANS simulation. For the ZDES case, unsteady fluctuations need to be added on the mean profile at the inlet. The fluctuations are created through a simultaneous separate simulation of a flat plate boundary layer using the recycling/rescaling method proposed by Edwards et al. [37]. At the walls of the combustor, a constant temperature of 300 K is imposed.
99 100 101 102 103 104 105 106 107 108 109
4. Results and discussions
110 111
Both the ZDES case and the URANS case were restarted from the flow field from a previous steady RANS simulation, and then run with a time step of 4E-7s. The ZDES case was run for approximately 20 ms and the URANS case was run for approximately 36 ms. From the steady RANS simulation we can deduce that the cell sizes in wall-normal direction at the wall are within one wall unit ( y + < 1) as shown in Fig. 5, which further enhances the suitability of the computational grid. We first show the mean pressure distribution at the lower-wall of the URANS case and the ZDES case in Fig. 6. The mean pressure distribution agrees quite well with each other. This suggests that ZDES and URANS present a similar ability with regard to calculating the mean variables. Next, we focus on the ability of ZDES and URANS to capture the unsteady motions of flow structures.
112 113 114 115 116 117 118 119 120 121 122 123 124
59 60 61 62 63 64 65 66
125
4.1. Pressure signals and spectra
126
We note here in advance that the experimental shadowgraph images show that there are shock wave structures oscillating around x = 385 mm (details are to be shown in Section 4.2). Thus, a probe is disposed at this location for both the ZDES case and the URANS case (probe A in Fig. 6). Fig. 7 illustrates the calculated
128
127 129 130
Fig. 7. Pressure signal at x = 385 mm at the down wall. a) the URANS case; b) the ZDES case.
131 132
JID:AESCTE AID:4278 /FLA
[m5G; v1.224; Prn:9/11/2017; 8:35] P.5 (1-9)
Y. Han et al. / Aerospace Science and Technology ••• (••••) •••–•••
5
1
67
2
68
3
69
4
70
5
71
6
72
7
73
8
74
9
75
10
76
11
77
12
78
13
79
14
80
15
81
16
82
17
83
18
84
19
85
20
86
21 22 23
87
Fig. 8. Power spectral density (PSD) function of the pressure signal at x = 385 mm at the down wall. For the URANS case, PSD function is calculated using the pressure signal from t = 3 ms to t = 26 ms.
88 89
24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44
90 91
unsteady pressure signal of this probe and shows that in both of the two cases the pressure signal presents a character of quasiperiodicity. To further affirm this, the PSD function of the pressure signal of this probe is illustrated in Fig. 8. The spectra display a peak at about 320 Hz for both of the two cases, indicating the quasi-periodicity. Despite the similarities, there are significant differences between the two cases with regard to the evolution of the pressure signal. In the URANS case, the amplitude of the pressure signal is regularly and gradually decreasing as time goes on, and the pressure signal finally approximates a steady status. This implies that for the URANS case the unsteady motions of flow structures around the probe are likely to decay gradually and finally vanish. However, the pressure signal in the ZDES case maintains large scale oscillation all the time. Besides, it seems that the oscillating duration and the oscillating amplitude vary irregularly in different cycles for the ZDES case (see for example, cycle II and cycle IV in Fig. 7b). This phenomenon is to be further discussed in the subsequent section.
92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110
45 46
111
4.2. Unsteady motions of the flow structures
112
47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66
113
In this subsection, we get interested in the evolution of flow structures captured by ZDES and URANS, combined with the experimental shadowgraph images. First we focus on the ZDES case. Two representative oscillating cycles, namely cycle II and cycle IV in Fig. 7b are selected for detailed discussions. In Fig. 9, the evolution of flow structures in cycle II is illustrated using numerical shadowgraph images which are acquired by calculating ∇ 2 ρ , the Laplacian of density. Six typical instantaneous images are shown corresponding to states T1–T6 in cycle II. In Fig. 10, the evolution of flow structures in cycle IV is illustrated with four typical instantaneous images corresponding to states t1–t4 in cycle IV. Before we start to discuss the evolution of the flow structures, two things need to be clarified here. The first is that experimental shadowgraph images show that the flow structures oscillate around x = 385 mm and in one test there are dozens of cycles whose durations can vary from 3 ms to 5 ms (not shown in this paper). It seems meaningless to compare the exact time of each state in Fig. 7b with that of experimental data. Thus, we focus on
114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129
Fig. 9. The evolution of flow structures in cycle II of Fig. 7b, illustrated by numerical shadowgraph images from ZDES Calculation (downside) and experimental shadowgraph images (upside).
130 131 132
JID:AESCTE
AID:4278 /FLA
[m5G; v1.224; Prn:9/11/2017; 8:35] P.6 (1-9)
Y. Han et al. / Aerospace Science and Technology ••• (••••) •••–•••
6
1
67
2
68
3
69
4
70
5
71
6
72
7
73
8
74
9
75
10
76
11
77
12
78
13 14
79
Fig. 9. (continued)
80
15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66
81
whether ZDES can capture the unsteady motions of the flow structures which exist in the experiments, and in Fig. 9 and Fig. 10 shadowgraph images taken from two cycles of the experimental data are shown for comparison with the numerical shadowgraph images. The second is that a wave (wave A shown in Fig. 9a) coming from the isolator exists in each experimental shadowgraph image. This wave originates from a small crevice on the upper-wall of the experimental isolator model, so it doesn’t exist in numerical shadowgraph images. Besides, in each experimental shadowgraph image, the edge of the cavity is redrawn using dash lines. The small region which is purely black under the dash lines is caused by the fuel-injector and the gasket between the cavity model and the sight window. Nevertheless, these matters don’t affect the main characteristics of the flow field. From Fig. 9 and Fig. 10, we can see that the evolution of the main flow structures, including the wave structures around x = 385 mm and the wave structures at the step of the cavity in experimental shadowgraph images are fairly captured by ZDES calculation. Fig. 9a shows that in state T1, the start of cycle II, no shock structures reach upstream of probe A (x = 385 mm) near the lower-wall, so the pressure signal of probe A is in the trough. The wave structures at the step of the cavity (wave structures B shown in Fig. 9a) are approximately parallel to the upper-wall of the cavity, indicating that the cavity is open. In Fig. 11a, the streamlines starting from the boundary layer of the upper-wall of the isolator develop directly to the downstream of the cavity without reattaching to the upper-wall of the cavity, thus also demonstrating that the cavity is open in this state. From Fig. 11a we can also see that no distinct reverse flow exists near the lower-wall. Besides, the streamwise velocity is high almost throughout the flow path in the region above the boundary layer, which we term as the mainflow hereinafter. As the shear layer shedding from the step of the cavity moves towards the upper-wall of the cavity, the relevant wave structures (wave structures B shown in Fig. 9b) deflect to the upper-wall too, and get reflected at point R (shown in Fig. 9b). Then the reflected shock wave (shock wave C shown in Fig. 9b) interacts with the lower-wall boundary layer as an incident shock, leading to an adverse pressure gradient and inducing a separation zone (shown more clearly in Fig. 11b). Then at the leading edge of the separation zone, a reflected shock wave (shock wave E shown in Fig. 9b) is formed and raises the pressure signal of probe A. These flow structures reproduce the phenomena of typical incident oblique shock/ turbulent boundary layer interaction. From Fig. 11b we can see that in this state the streamlines nearly reattach to the upperwall of the cavity, indicating that the cavity is transiting to a closed cavity. Besides, the streamwise velocity near the lower-wall downstream of the separation zone decreases, compared with that in state T1.
82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118
Fig. 10. The evolution of flow structures in cycle IV of Fig. 7b, illustrated by numerical shadowgraph images from ZDES calculation (downside) and experimental shadowgraph images (upside).
119 120 121 122
As wave structures B deflect more to the upper-wall, the reflected shock wave C gets stronger. Meanwhile, the separation zone S moves upstream, along with the reflected shock wave E. The pressure signal of probe A keeps rising, as shown in Fig. 7b. In state T3, shock wave C and shock wave E interact with each other and form a Mach disc, as shown in Fig. 9c. As time goes on, these flow structures keep moving upstream and the state comes to T4. In this state, in addition to the shock structures occurring in state T3, a bow shock wave (shock wave G shown in Fig. 9d) occurs near the exit of the cavity. The low-speed area downstream of sepa-
123 124 125 126 127 128 129 130 131 132
JID:AESCTE AID:4278 /FLA
[m5G; v1.224; Prn:9/11/2017; 8:35] P.7 (1-9)
Y. Han et al. / Aerospace Science and Technology ••• (••••) •••–•••
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66
ration zone S becomes larger than that of state T3 owing to these two shock structures, as shown in Fig. 11d. From this figure we can also see that the streamlines completely reattach to the upper-wall of the cavity, indicating that the cavity is a closed cavity in this state. After state T4, wave structures B start to deflect against the upper-wall and the reflected shock wave C starts to move downstream, along with shock wave E and separation zone S. From Fig. 9f we can see that in state T5 the structure of Mach disc disappears and the flow structures are similar to those of state T2. The pressure signal of probe A starts to decrease, as shown in Fig. 7b. The cavity also gradually turns from a closed cavity to an open cavity, as shown in Fig. 11e. In state T6, separation zone S and the ambient shock structures disappear, and the streamwise velocity of the mainflow at the downstream becomes large again. The cavity has almost turned into an open cavity. After state T6, oscillating cycle II is completed. From Fig. 10 we can see that cycle IV also consists of a complete oscillating cycle of the motion of shock structures around x = 385 mm and the change of the status of the cavity. However, cycle IV does not contain a state in which a Mach disc occurs, compared with cycle II. Besides, the most upstream position that the shock structures can move to in cycle IV is less than that in cycle II. These are consistent with the character of the pressure signal: both the oscillating amplitude and the duration of cycle IV are less than those of cycle II, as illustrated in Fig. 7b. We denote the oscillating cycles in which the evolution of flow structures is similar to that of cycle II and in which the evolution of flow structures is similar to that of cycle IV as strong cycles and weak cycles, respectively. Experimental shadowgraph images show that the flow is comprised of numbers of oscillating cycles that mix the strong cycles and weak cycles: one or more strong cycles are followed by one or more weak cycles. (They can’t be entirety shown in this paper, though.) The precise number of continuous appearance for strong cycles or weak cycles appears to be irregular. The ZDES case successfully reproduces the phenomena, as indicated by Fig. 7b: cycles I–V are weak cycle, strong cycle, strong cycle, weak cycle and strong cycle, respectively. Interestingly, the behavior of mixing is somewhat analogous to a mixed oscillatory pattern of hypersonic inlet buzz which mixes the “big buzz” with “little buzz” as shown in [38], although these two phenomena are apparently arisen from different flows. Next we focus on the URANS case. Fig. 12 and Fig. 13 show the unsteady motions of flow structures in cycle 2 of Fig. 7a. Note that some unphysical structures are shown in Fig. 12 (illustrated by unphysical structures U in Fig. 12a). They are caused by the extremely large aspect ratio of the grid cells at their locations, as illustrated in Fig. 3. The analogous unphysical structures exist in the ZDES case too, though they are not highlighted in Fig. 9 and Fig. 10. Nevertheless, all of these unphysical structures don’t affect the analysis of the flow field. From Fig. 12 we can see that during cycle 2 of Fig. 7a, the shape of the shock structures around x = 385 mm (illustrated by shock structures S in Fig. 12a) is basically similar for each state, and the shock structures exist throughout the cycle. The unsteady motion of the shock structures is on a small scale, and the most upstream location that the shock structures can affect is about x = 380 mm. From the streamlines shown in Fig. 13 we can see that the cavity keeps open throughout the whole cycle, contradicting with the experiments and the ZDES case, in which the cavity transforms between the open state and the closed state in an oscillating cycle. Moreover, as shown in Section 4.1, the amplitude of the pressure signal is gradually decreasing as time goes on and finally approximates a steady status in the URANS case, indicating that the unsteady motions of flow structures around x = 385 mm are likely to decay gradually as time goes on and finally vanish for the
7
67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117
Fig. 11. The evolution of flow structures in cycle II of Fig. 7b, illustrated by contours of velocity u from ZDES calculation. The grey region of each subfigure shows the region where u < 0.
118 119 120 121
URANS case. This is further confirmed by the numerical shadowgraph images in cycle 8 of Fig. 7a, as shown in Fig. 14. In this figure, state t1 and state t2 refer to the states when the pressure signal reaches the trough and crest respectively, as shown in Fig. 7a. The shock structures S in these two states are almost of the same shape, and the difference of the locations of the shock structures between these two states is rather tiny. These indicate that the unsteady motion of the shock structures has almost vanished in this cycle. It thus can be concluded from the above that the results of the URANS case distinctly disagree with those of the ZDES case
122 123 124 125 126 127 128 129 130 131 132
JID:AESCTE
8
AID:4278 /FLA
[m5G; v1.224; Prn:9/11/2017; 8:35] P.8 (1-9)
Y. Han et al. / Aerospace Science and Technology ••• (••••) •••–•••
1
67
2
68
3
69
4
70
5
71
6
72
7
73
8
74
9
75
10
76
11
77
12
78
13
79
14
80
15
81
16
82
17
83
18
84
19
85
20
86
21
87
22
88
23
89
24
90
25
91
26
92
27
93
28
94
29
95
30
96
31
97
32
98
33
99
34
100
35 36 37
101
Fig. 12. The evolution of flow structures in cycle 2 of Fig. 7a, illustrated by numerical shadowgraph images from URANS calculation.
102 103
38
104
39
105
40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66
106
and the experiments, indicating that the reliability of URANS in capturing unsteady motions of flow structures in scramjet is questionable. Next we’d like to conduct analysis on why the results of the URANS case are problematic. As shown in Menter et al. [5], turbulent flows can be divided into three categories: ‘global unstable flow’, ‘locally unstable flow’ and ‘stable flow’. It’s cumbersome to describe the definitions of these three types of flow and readers can refer to [5] for details. URANS simulation can provide an unsteady output for ‘global unstable flow’ among these three types of flows. For the other two types of flow, URANS simulation is in all probability unable to provide a reliable unsteady output, unless there are other proper influencing factors, such as highly unsteady boundary conditions. In this study, there are two most significant flow structures which contribute to the flow unsteadiness in the scramjet combustor: one is the shear layer shedding from the step of the cavity and the other is the shock wave/turbulent boundary layer interaction. With respect to these two flows, the former belongs to the category of ‘locally unstable flow’, and the latter belongs to the category of ‘stable flow’ [5], so URANS simulation is unable to provide reliable unsteady results for these two flows. Besides, the deficiency of URANS in simulating shock wave/turbulent boundary layer interaction is also stated in Georgiadis et al. [9]. Thus, it’s not surprising that the results of the URANS case in this study are problematic with regard to capturing the unsteady motions of flow structures in the scramjet combustor.
107 108 109 110
Fig. 13. The evolution of flow structures in cycle 2 of Fig. 7a, illustrated by contours of velocity u from ZDES calculation. The grey region of each subfigure shows the region where u < 0.
111 112 113 114
5. Conclusions
115 116
In this paper, a scramjet combustor is simulated using both URANS method and ZDES method, to conduct an assessment of their abilities in depicting the unsteady motions of flow structures in the scramjet combustor. Although ZDES and URANS present a similar ability with regard to calculating the mean variables, they present distinct differences concerning the ability to capture the unsteady motions of flow structures. The numerical shadowgraph images of the ZDES case successfully capture the complicated unsteady motions of the shock structures around x = 385 mm: the resolved flow is comprised of numbers of oscillating cycles that mixes the strong cycles and weak cycles, and one or more strong cycles are followed by one or more weak cycles. Besides, the transiting of the open state and closed state of the cavity is also captured by ZDES. These coincide with the experiments fairly well. In contrast, the results of URANS distinctly disagree with those of the ZDES case and the experiments.
117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132
JID:AESCTE AID:4278 /FLA
[m5G; v1.224; Prn:9/11/2017; 8:35] P.9 (1-9)
Y. Han et al. / Aerospace Science and Technology ••• (••••) •••–•••
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23
Fig. 14. The evolution of flow structures in cycle 8 of Fig. 7a, illustrated by numerical shadowgraph images from URANS calculation.
24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56
Conflict of interest statement The authors declare there are no conflicts of interest regarding the publication of this paper. Acknowledgements The project was supported by National Natural Science Foundation of China [grant numbers91641107, 51706237]. References [1] M.K. Smart, N.E. Hass, Flight data analysis of the HyShot 2 scramjet flight experiment, AIAA J. 44 (10) (2006) 2366–2375. [2] Y. Tian, S. Yang, J. Le, et al., Investigation of combustion and flame stabilization modes in a hydrogen fueled scramjet combustor, Int. J. Hydrog. Energy 41 (2016) 19218–19230. [3] Y. Tian, S. Yang, J. Le, Study on flame stabilization of a hydrogen and kerosene fueled combustor, Aerosp. Sci. Technol. 59 (2016) 183–188. [4] Y. Tian, S. Yang, J. Le, et al., Investigation of the effect of fuel injector locations on ignition and flame stabilization in a kerosene fueled scramjet combustor, Aerosp. Sci. Technol. 70 (2017) 310–316. [5] F.R. Menter, J. Schütze, M. Gritskevich, Global vs. zonal approaches in hybrid RANS-LES turbulence modelling, in: Progress in Hybrid RANS-LES Modelling, Springer, Berlin, Heidelberg, 2012, pp. 15–28. [6] J. Larsson, S. Laurence, I. Bermejo-Moreno, et al., Incipient thermal choking and stable shock-train formation in the heat-release region of a scramjet combustor. Part II: large eddy simulations, Combust. Flame 162 (2015) 907–920. [7] S. Deck, N. Renard, R. Laraufie, et al., Zonal detached eddy simulation (ZDES) of a spatially developing flat plate turbulent boundary layer over the Reynolds number range 3150 ≤ Reθ ≤ 14000, Phys. Fluids 26 (2014) 025116. [8] Y. Han, G. Ding, Y. He, et al., Assessment of IDDES method acting as wall modeled LES in the simulation of spatially developing supersonic flat plate boundary layers, Eng. Appl. Comput. Fluid 12 (1) (2018) 89–103. [9] N.J. Georgiadis, D.P. Rizzetta, C. Fureby, Large-eddy simulation: current capabilities, recommended practices, and future research, AIAA J. 48 (8) (2010) 1772–1784.
9
[10] N.J. Bisek, High-fidelity Simulations of the HIFiRE-6 Flow Path, 2016, AIAA 2016-1115. [11] E.B. Tylczak, D.M. Peterson, G.V. Candler, Hybrid RANS/LES Simulation of Injection and Mixing in the CUBRC Combustion Duct, 2011, AIAA 2011-3216. [12] J.R. Edwards, A. Potturi, J.A. Fulton, Large-Eddy/Reynolds-Averaged Navier– Stokes Simulations of Scramjet Combustor Flow Fields, 2012, AIAA 2012-4262. [13] J.A. Fulton, J.R. Edwards, H.A. Hassan, et al., Large-Eddy/Reynolds-Averaged Navier–Stokes Simulations of a Dual-Mode Scramjet Combustor, 2012, AIAA 2012-0115. [14] P.A.T. Cocks, C. Bruno, J.M. Donohue, et al., IDDES of a Dual-Mode Ethylene Fueled Cavity Flameholder with an Isolator Shock Train, 2013, AIAA 2013-0116. [15] J. Liang, C. Gong, M. Sun, Hybrid RANS/LES Simulation of the Fuel Injection and Mixing Process in a Supersonic Combustor with Cavity Flame-Holder, 2013, AIAA 2013-3698. [16] K. Makowka, T. Sattelmayer, M. Gurtner, et al., Semi-zonal hybrid RANS/LES turbulence modeling with RANS sensor based interfacing applied to supersonic flows, 2014, AIAA 2014-3087. [17] A.S. Potturi, J.R. Edwards, Hybrid large-eddy/Reynolds-averaged Navier–Stokes simulations of flow through a model scramjet, AIAA J. 52 (7) (2014) 1417–1429. [18] D.M. Peterson, E. Hassan, Simulating Turbulence and Mixing in Supersonic Combustors Using Hybrid RANS/LES, 2015, AIAA 2015-0636. [19] J. Shin, H. Sung, Numerical Simulation of Hydrogen Combustion in Model SCRAMJET Combustor Using IDDES Framework, 2015, AIAA 2015-3625. [20] E. Hassan, D.M. Peterson, J. Liu, et al., Reacting Hybrid Reynolds-Averaged Navier–Stokes/Large Eddy Simulation of a Supersonic Cavity Flameholder, 2016, AIAA 2016-4566. [21] R.A. Baurle, Hybrid Reynolds-averaged/large-eddy simulation of a cavity flameholder: modeling sensitivities, AIAA J. 55 (2) (2017) 524–543. [22] K. Makowka, N.C. Dröske, J. von Wolfersdorf, et al., Hybrid RANS/LES of a supersonic combustor, Aerosp. Sci. Technol. 69 (2017) 563–573. [23] J. Choi, V. Yang, F. Ma, et al., DES Combustion Modeling of a Scramjet Combustor, 2006, AIAA 2006-5097. [24] F.J. Malo-Molina, D.V. Gaitonde, H.B. Ebrahimi, et al., Three-dimensional analysis of a supersonic combustor coupled to innovative inward-turning inlets, AIAA J. 48 (3) (2010) 572–582. [25] J. Boles, J. Choi, J. Edwards, et al., Simulations of High-Speed Internal Flows Using LES/RANS Models, 2009, AIAA 2009-1324. [26] H. Koo, V. Raman, Large-eddy simulation of a supersonic inlet-isolator, AIAA J. 50 (7) (2012) 1596–1613. [27] J. Le, S. Yang, W. Liu, et al., Massively Parallel Simulations of Kerosene-Fueled Scramjet, 2005, AIAA 2005-3318. [28] Y. Tian, S. Yang, J. Le, Numerical study on effect of air throttling on combustion mode formation and transition in a dual-mode scramjet combustor, Aerosp. Sci. Technol. 52 (2016) 173–180. [29] F.R. Menter, Two-equation eddy-viscosity turbulence models for engineering applications, AIAA J. 32 (8) (1994) 1598–1605. [30] M. Strelets, Detached Eddy Simulation of Massively Separated Flows, 2001, AIAA 2001-0879. [31] S. Deck, Recent improvements in the zonal detached eddy simulation (ZDES) formulation, Theor. Comput. Fluid Dyn. 26 (2012) 523–550. [32] T.H. Pulliam, Time Accuracy and the Use of Implicit Methods, 1993, AIAA 1993-3360. [33] A. Jameson, S. Yoon, L.U. Implicit, Scheme with Multiple Grid for the Euler Equations, 1986, AIAA 1986-0105. [34] J. Fang, Z.R. Li, L.P. Lu, An optimized low-dissipation monotonicity-preserving scheme for numerical simulations of high-speed turbulent flows, J. Sci. Comput. 56 (2013) 67–95. [35] G. Billet, O. Louedin, Adaptive limiters for improving the accuracy of the MUSCL approach for unsteady flows, J. Comput. Phys. 170 (2001) 161–183. [36] M. Shur, P.R. Spalart, M. Strelets, et al., A rapid and accurate switch from RANS to LES in boundary layers using an overlap region, Flow Turbul. Combust. 86 (2011) 179–206. [37] J.R. Edwards, J. Choi, J.A. Boles, Hybrid LES/RANS simulation of a Mach5 compression–corner interaction, AIAA J. 46 (4) (2008) 977–991. [38] J. Chang, L. Wang, W. Bao, et al., Novel oscillatory patterns of hypersonic inlet buzz, J. Propuls. Power 28 (6) (2012) 1214–1221.
67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122
57
123
58
124
59
125
60
126
61
127
62
128
63
129
64
130
65
131
66
132