Numerical simulation of hypersonic intake starting characteristics is presented. Three dimensional RANS equations are solved alongwith SST turbulence model using commercial computational fluid dynamics (CFD) software. Wall pressure distribution and intake performance parameters are found to match well with experimental data for different free stream Mach number in the range of 3-8. The unstarting of the intake is traced from the sudden drop of mass capture ratio. Wall condition (adiabatic or isothermal) is seen to have pronounced effect in estimating the performance parameters in the intake. The computed unstarting Mach number is seen to be higher for adiabatic condition compared to isothermal condition. For unstarting case, large separation bubble is seen near the entrance of the intake, which is responsible for expulsion of the shock system out of the intake.

Keywords:    Numerical simulation ,  hypersonic intake,  computational fluid dynamics

The performance of a ramjet/scramjet powered hypersonic vehicle is determined by its inlet efficiency as the engine depends very much on the quantity and quality (uniformity and total pressure) of the flow required for its smooth performance. Hypersonic intakes are designed as mixed compression intake which is a combination of internal and external compression. A schematic of the flow pattern in the mixed compression intake is shown in Fig 1. The bow shock of the vehicle forebody compresses the air followed by the number of compression at the central body which coalesces at the cowl lip at design Mach number. The flow is turned inward to the axial direction in the internal compression zone by the lip geometry. The interaction of the reflected shock with the existing boundary layer on the ramp surface might lead to formation of a separation zone. The extent of separation will depend on the strength of the reflected shock and the condition of the boundary layer on the ramp surface.

One major problem of hypersonic intakes with internal compression is the ‘unstart problem’ which describes the phenomenon that supersonic internal flow is not reached in the internal compression region if the area ratio between the throat and capture area is too low when the flight Mach number is increased. The internal flow remains subcritical and the intake is chocked. Generally unstart of the intake is observed through expulsion of the shock system and massive spillage, leading to degraded pressure recovery and large flow distortion at the exit and hence there may be catastrophic effect in the vehicle performance. The unstart of the intake could occur due to several reasons, e.g. over-contraction, variation of flight conditions, perturbations in combustor operation, back pressure, angle of attack, etc., or due to a combined effect of these factors. Interaction of the boundary layer with shock reflections and subsequent thickening of the boundary layer inside the internal duct, are believed to be the prime cause of a separation leading to a complex oscillatory flow structure and expulsion of the shock and the unstart of the intake. Usually, to start supersonic inlets at any flight condition, variable intake geometry or bleed bypass is used. But in a hypersonic flow situation which contains high enthalpy flow with high total temperature (~1800 K), any complex mechanical control system may cause severe structural and cooling problems.

The prediction of intake unstart and the mitigation plan to reduce its occurrence or its effect is very much essential for hypersonic intake design. Experimental and numerical research is in progress to understand the causes of hypersonic intake unstart and means to avoid it. Schmitz and Bissinger¹. studied experimentally two fixed geometry hypersonic intakes and reported starting of the intake at M = 4.3 and stable operation up to M = 6. Various performance parameters of two geometries were also compared. Schneider and Koschel². studied both experimentally and numerically the start and throttling behavior of a supersonic intake system of 9 different configurations with geometric variation at different inlet Mach number and exit throttling conditions. It was shown that by proper geometry selection, the size of the separation bubble at the ramp surface could be minimized without applying boundary layer bleed and high intake performance could be achieved. Goonko, et al. ³ reported experimental studies of three dimensional inlet at Mach number range 4 to 8 with ramp and side-swept compression wedges and depicted complex system of multiple shock waves, expansion waves and vortex structures leading to significant nonuniformity of the flow field at the inlet entrance and exit.

Haberle and Gulhan⁴⁵ studied experimentally the effect of different bleed dimension and internal contraction ratios on starting characteristics of 2-D and 3-D scramjet inlets at hypersonic Mach number. Das and Prasad⁶. have conducted both experimental and numerical investigation of mixed compression intake flow field at Mach 2.2 with different cowl deflections and showed that small angle at the cowl lip leads to start of the intake and improve its performance. Experimental and numerical studies of Dirk, et al.⁷demonstrate that bleed can reduce boundary layer separation and improve the inlet total pressure recovery. Molder, et al ⁸ have studied the effect of mass extraction and overboard spillage as a means to start the intake and defined a starting index based on this study. They have also shown that for M>4, the starting operation regime can be increased by mass extraction. Reinhold, et al⁹ studied a multiple strut based 2-D hypersonic ramjet inlet flowfield using a Parabolised Navier Stokes code and presented good comparison of simulation results against experiment. Lind, et al ¹⁰ have studied the effect of forebody shock and cowl lip shock interaction demonstrating that, the shock-shock interaction can lead to very high pressure and temperature region around the cowl lip resulting in flow instability. It was observed that free stream variations have strong influence in making the flow unsteady. Brenneis, et al. ¹¹. have studied 2-D inlet at freestream Mach number 7.4 and shown that there is drastic change in the behavior of the flowfield with adiabatic wall compared to fixed temperature wall. Barber, et al. ¹². have carried out intercode comparisons to predict the starting characteristics of hypersonic intakes. It has been shown that starting characteristics is strongly dependent on viscous flow effects and choice of turbulence model has significant impact on the prediction of shock wave boundary layer interaction (SWBLI). Predicted transient calculations exhibit a time-lag effect related to SWBLI effects. Donde, et al ¹³. carried out numerical simulations of a starting problem in a variable geometry hypersonic intake with a movable cowl. Dynamic meshes have been used for depicting motion of the cowl. It was shown that the cowl needs to be rotated through 15.7° and then be brought back to the original position for restarting of the intake after an ‘unstart’. Yu, et al. ¹⁴. studied hypersonic inlet numerically for different freestream conditions and the backpressure using a RANS code with RNG k-ε turbulence model. A new method of pattern classification of inlet start/unstart is developed using the ‘numerical/experimental’ database by the support vector machine-recursive feature elimination algorithm. The control can sense inlet start/unstart through this pattern classification data. In order to develop these databases numerical tools played an important tool.

It is clear from the above discussion that the problem of starting of hypersonic intake has not been fully understood and requires further investigation. In this work, starting problem of hypersonic air-intake pertaining to experimental condition of Emami and Trexler ¹⁵ is simulated by solving 3-D RANS equations alongwith SST turbulence model using a commercial code ¹⁶. The computed results were compared with experimental values and insights were obtained about the complex process through analysis of various flow variables.

Figure 2 presents the intake geometry of Emami and Trexler ¹⁶. alongwith dimensions. In the geometry, the x-axis is taken along the longitudinal direction, while y and z axes are taken along the height and width of the model respectively. The intake is a 2-D geometry with ramp in the lower surface, which makes an angle of 11 degrees to the x-axis. The cowl is attached to the upper surface, and is tilted downwards by 3 degrees about x-axis. The throat height (h) is 0.01 m. All other dimensional details are provided in Fig 2. Side fencing of the intake starts ahead of the cowl lip. It is provided to prevent the flow compressed by ramp from spilling sidewise into ambient. 3-D computational domain is created to account for the simulation which is shown in Fig 3 along with all the boundaries and is extended upto a length of 0.34m in the longitudinal direction.

Since the intake is symmetric about mid-plane, the computational domain upto half of the intake width (25.4 mm) is simulated. Hexahedral grids of sizes 0.4 million (coarse grid) and 0.68 million (fine grid) are generated for grid variation study. In the forebody and in the intake interior, the number of cells in X,Y and Z directions are 338X49X30 and 259X32X30 for fine and coarse grid respectively. The grids are made very fine near the intake entry and near the throat for capturing the flow crisply. The axial distribution of Y⁺ along the length of the lower surface is presented in Fig. 4. In most of the length Y⁺ is less than 20; although in the downstream of the throat Y⁺ value is about 40. The minimum distance at the first cell is kept 0.05 mm and the sizes were exponentially stretched along the height. There are 5 to 10 cells in boundary layers at different regions.

3-D Reynolds averaged unsteady Navier stokes equations with SST turbulence model are solved. A density based solver with 2^nd order spatially accurate Roe-Flux difference splitting scheme ¹⁷. is used for spatial discretisation and 2^nd order implicit Euler Scheme for temporal discretisation is used in the present solution.

3.1 Governing Equation

The appropriate system of equations governed the turbulent compressible gas may be written as
Continuity equation:
$\frac{\partial \rho }{\mathrm{\partial t}}+\frac{\partial }{\partial {\chi }_{\kappa }}(\rho {\mu }_{\kappa })=0$     _{k = 1,2,3}

Momentum equation :
$\frac{\partial }{\partial t}(\rho u_i)+\frac{\partial }{\partial x_{\mathit{\kappa }}}(\rho u_i{u}_k)+\frac{\partial \rho }{\partial x_{\mathrm{i}}}=\frac{\partial ({\mathit{\tau }}_{\mathrm{ik}})}{\partial x_{\mathrm{i}}}$     _{,i,k =1,2,3}

Energy equation :
$\frac{\mathit{\partial }}{\mathit{\partial }t}\mathit{(}\mathit{\rho }E\mathit{)}\mathit{+}\frac{\mathit{\partial }}{\mathit{\partial }{x}_k}\mathit{(}\mathit{\rho }{u}_{k}H\mathit{)}\mathit{=}\mathit{-}\frac{\mathit{\partial }}{\mathit{\partial }{x}_k}\mathit{(}{u}_j{\mathit{\tau }}_{\mathrm{jk}}\mathit{)}\mathit{+}\frac{\mathit{\partial }{q}_k}{\mathit{\partial }{x}_k}$ ,
j,k = 1,2,3

where, ρ,u i,p, E and Hare the density, velocity components, pressure and total energy and enthalpy respectively Turbulent shear stress is defined as
μ= μ_l + μ_t is the total viscosity; μ_l , μ_t being the laminar and turbulent viscosity
Laminar viscosity (μ_l ) is calculated from Sutherland law as
${\mu }_{\iota }={\mu }_{\mathrm{ref}}{\left(\frac{\mathrm{T}}{{\mathrm{T}}_{\mathrm{ref}}}\right)}^{\raisebox{1ex}{$3$}\!\left/ \!\raisebox{-1ex}{$2$}\right.}\left(\frac{{\mathrm{T}}_{\mathrm{ref}}}{+}\right)$
where, T is the temperature and μ_ref , T_ref and S are known coefficient.
In eddy viscosity models, the stress tensor is expressed as a function of turbulent viscosity (μ_t ). Based on dimensional analysis, few variables (k, ε, ω) are defined as given below,
Turbulent kinetic energy k,
$\mathrm{k}=\overline{{{\mathrm{u}}_{\mathrm{i}}}^{\prime }{{\mathrm{u}}_{\mathrm{i}}}^{\prime }}/2$

Turbulent dissipation rate ε,
$\mathit{ϵ}\mathit{\equiv }\mathit{\nu }\overline{\frac{\mathit{\partial }{u_i}^{\mathit{\prime }}}{\mathit{\partial }{x}_j}\left(\frac{\mathit{\partial }{u}_i}{\mathit{\partial }{u}_j},+,\frac{{{\mathit{\partial }u}^{\mathit{\prime }}}_j}{{\mathit{\partial }x}_i}\right)}$

Specific dissipation rate ω,
ω = ε /k
The turbulent viscosity μ_t is calculated as
${\mathit{\mu }}_t\mathit{=}{c}_{\mathit{\mu }}\frac{\mathit{\rho }{k}^{\mathit{2}}}{\mathit{ϵ}}$
The heat flux q_k is calculated as $q_k\mathit{}\mathit{=}\mathit{-}\mathit{\lambda }\frac{\mathit{\partial }T}{\mathit{\partial }{x}_k}$, λ is the thermal conductivity

3.2 k-ω Turbulence Model

The turbulent viscosity is calculated as function of k and ω ¹⁸

${\mu }_{\mathrm{t}}=\mathrm{f}\left(\frac{{\rho }^{\mathrm{k}}}{\omega }\right)$

Turbulent kinetic energy (k) equation:

$\frac{\partial }{\partial \mathrm{t}}(\rho \mathrm{k})+\frac{\partial }{\partial {\mathrm{x}}_{\mathrm{i}}}(\rho \mathrm{k}{\mathrm{u}}_{\mathrm{i}})=\frac{\partial }{\partial {\mathrm{x}}_{\mathrm{j}}}\left({\Gamma }_{\mathrm{k}}\frac{\partial \mathrm{k}}{\partial {\mathrm{x}}_{\mathrm{j}}}\right)+{\mathrm{G}}_{\mathrm{k}}-{\mathrm{Y}}_{\mathrm{k}}$ (2.4.13)
Specific dissipation rate (ω) equation:
$\frac{\mathit{\partial }}{\mathit{\partial }t}\mathit{(}\mathit{\rho }\mathit{\omega }\mathit{)}\mathit{+}\frac{\mathit{\partial }}{\mathit{\partial }{x}_i}\mathit{(}\mathit{\rho }\mathit{\omega }{u}_i\mathit{)}\mathit{=}\frac{\mathit{\partial }}{\mathit{\partial }{x}_j}\left({\mathit{\Gamma }}_{\mathit{\omega }}\frac{\mathit{\partial }\mathit{\omega }}{\mathit{\partial }{x}_j}\right)\mathit{+}{G}_{\mathit{\omega }}\mathit{-}{Y}_{\mathit{\omega }}$
where are the production, dissipation and diffusion terms of the k and ω equations respectively

3.3 SST Turbulence Model

To retain the robust and accurate formulation of Wilcox’s k-ω model in the near wall region, and to take advantage of the freestream independence of the k-ε model in the outerpart of the boundary layer, Menter ¹⁹. blended both the models through a switching function. k-ε model was transformed into Wilcox’s k-ω formulation and was multiplied by (1–F 1) and added to original k-ω model multiplied by F1. The blending function F1will be one in the near wall region and zero away from the surface.
In the second step, the definition of eddy viscosity was modified in the following way to account for the transport of the principal turbulent shear stress ($\mathit{\tau }\mathit{}\mathit{=}\overline{\mathit{-}\mathit{\rho }{u}^{\mathit{\prime }}{v}^{\mathit{\prime }}}$ )

${\nu }_{\mathrm{t}}=\frac{a_{\mathit{1}}k}{\max ({\mathrm{a}}_1\omega ;\Omega {\mathrm{F}}_2)}$

3.4 Boundary Condition

At the inlet of the intake, freestream static pressure and total temperature (T₀) condition is imposed. By varying the inlet total pressure (P₀), different free stream Mach numbers have been simulated. At the intake exit boundary, supersonic outflow boundary condition has been employed. For the ambient outlet also, the supersonic outflow condition is imposed. Adiabatic, no slip wall condition has been imposed for top, bottom, side and cowl walls. To find the effect of wall temperature on intake performance, simulations were also performed with 300K temperatures on the wall. All simulations is performed at zero angle of attack.

4.1 Steady State Simulation for Inlet with Free Stream Mach number 4

Figure 5 shows Mach number distribution in the intake mid plane for free stream Mach number 4 and zero angle of attack. Zoomed view of Mach number around intake cowl lip is also shown in the figure. It can be seen that oblique shock generated from the lower ramp is above the cowl lip indicating flow spillage. The shock from cowl lip impinges inside the intake on the ramp surface. Further, the cowl shock gets reflected from both the ramp and upper wall multiple times to create a train of oblique shocks inside the intake which further compresses the flow. The axial distribution of lower ramp wall pressure are compared with experimental data ¹⁵. in Fig. 6. The pressure and distances are nondimensionalised with free stream pressure (pinf) and intake throat height (h) respectively. Though, in general a good agreement is obtained between prediction and experiment, there is a difference in the position of the cowl shock impingement point. The effect of grid on the results is studied by computing the flow on fine grid. The comparison of ramp pressure distribution for both the grids is also plotted in Fig. 6.

 

Though, there is no significant change in pressure distribution up to x = 0.28 m, the pressure peaks have become sharper with fine grid. Beyond x = 0.28m, shock reflections were captured more crisply in the fine grid. The difference in surface pressure between the experimental and numerical values near x/h=20 is due to the inability of the turbulence model to predict the pressure ahead of the separation bubble. Incorporation of unsteady correction term is suggested ²⁰. in recent literature to improve the prediction. Recent version of k-ω model²¹. has also shown improved performance in predicting the high speed flow. However, in the present study, SST turbulence model in its original form is used. The velocity vector plot inside the intake is shown in Fig. 7 where two shock induced separation bubbles are seen. To locate the starting of separation bubble, the axial wall shear stress distribution is presented in Fig. 8. The wall shear stress is seen to be nearly zero at x/h= 22.6 and 26.8 indicating the location of separation bubble inside the intake.

 

4.2 Simulation for different Free Stream Inlet Mach Number–unstarting of Intake

In order to bring out the condition at which the intake unstarting occurs, steady state simulations were also carried out for Mach 8, 7, 6, 5, 4, 3.5, 3.3, 3.2 and 3 conditions respectively by varying the total pressure at inlet. Initially, all the simulations were carried out under adiabatic wall condition. As the test duration is very short ¹⁵., it is expected that wall temperature will not change much during test duration. Hence,

simulations with wall temperature 300 K were also carried out. For determining the unstart condition of the intake, simulations are done first at a higher Mach number at which the intake starts and the Mach number is then gradually reduced to lower values. The converged solution corresponding to higher Mach number is fed as the initial flow field values for lower inlet free stream Mach numbers. Computed values of Mass capture ratio and mass averaged total pressure ratios for different free stream Mach numbers for adiabatic and isothermal wall conditions are compared with experimental results in Figs 9 and 10 respectively.

 

It is observed that for computed values of mass flow ratios and pressure recoveries match well with the isothermal calculation, while adiabatic wall condition shows a lower value. The difference in intake performance due to change in wall temperature has also been reported by Brenneis, et al. ¹¹It is observed that there is sudden drop in mass capture ratio and total pressure recovery around Mach number 3 whereas this drop is observed at Mach number 4 for adiabatic wall condition. This drop in values can be considered as the unstarting of the intake. The drop in mass capture ratiois happening because of strong flow separation due to shock-boundary layer interaction which is responsible for most of the flow spillage. For adiabatic wall condition the flow spillage is more compared to isothermal wall condition. This fact will be evident from the comparison of Mach number distribution at symmetry plane for adiabatic and isothermal case for M∞=3.5 in Fig 11. Large separation bubble near the intake entrance for the adiabatic case is responsible for expulsion of shock system outside the intake and causing the intake to unstart. It can be noticed that in adiabatic case, an oblique shock is generated ahead of flow separation in the cowl entry. This causes a large flow spillage and drop in mass capture. No separation bubble is seen for the isothermal case. The higher gas temperature for adiabatic wall, (900 K, the recovery temperature) reduces density and momentum of flow; whereas for isothermal wall condition, the boundary layer flow is relatively at a lower temperature causing higher momentum. As a result, for adiabatic wall condition the boundary layer flow is expected to separate at a higher Mach number because of lower momentum compared to isothermal wall condition.

 

Starting characteristics of a 2-D hypersonic intake with side fencing is presented. 3-D RANS equations are solved with SST turbulence model using commercial CFD software. The steady state simulations have been done for different Mach numbers ranging from 3.0 to 8.0 for both adiabatic and isothermal wall condition (Tw = 300K). Computed wall pressure distribution and intake performance parameters (mass flow ratio and total pressure recovery) match reasonably well with experimental results. The free stream Mach number for which intake unstarting occurs is found out from the sudden drop in the mass capture ratio. It is observed that wall boundary condition for temperature (adiabatic or isothermal) has a pronounced effect in determining the starting Mach number. Computed free stream Mach number for which unstarting occur is higher for adiabatic condition compared to isothermal condition. Heated boundary layer for adiabatic condition is seen to cause large separation bubble at the intake entrance causing flow unstarting; while flow separation bubble is not observed for isothermal condition for same free stream Mach number.

1. Schmitz, D.M. Bissinger, N.C. Design and testing of 2-D fixed –Geometry hypersonic intakes. AIAA Paper No. 98-1529, 1998.

2. Schneider, A. Koschel, W.W. Detailed analysis of a mixed compression hypersonic intake. ISABE Paper No. 99-7036, 1999.

3. Goonko, Y.P.; Latypov, A.F.; Mazhul, I.I.; Khartinov, A.M.; Yaroslavtsev, & Rostland, P. Structure of flow over a hypersonic inlet with side compression wedges. AIAA Journal, 2003, 41(3), 436-47

4. Haberle, J. & Gulhan, A., Experimental investigation of a two-dimensional and three dimensional scramjet inlet at Mach 7. J. Propulsion Power, 2008, 24(5), 1023-034.

5. Hebrale, J. & Gulhan, A. Internal flowfield investigation of hypersonic inlet at Mach 6 with bleed. J. Propulsion Power, 2007, 23(5), 1007-017.

6. Das, S. & Prasad J.K. Starting characteristics of a rectangular supersonic airintake with cowl deflection. Aeronautical Journal, 2010, 114(3), 177-89.

7. Dirk, S.; Andreas, H. & Ulrich, W. Reduction of shock induced boundary layer separation in hypersonic inlets using bleed. Aerospace Sci. Technol., 1998, 2(4), 231-39.

8. Molder, S.; Timofeev, E.V. Tahir, R.B. Flow starting in high compression hypersonic air inlets by mass spillage.AIAA Paper No. 2004-4130, 2004.

9. Reinhold, A.; Gerbsch & Agarwal, Ramesh K. Computation of hypersonic ramjet-inlet flowfields using an upwind parabolized Navier-Stokes code. AIAA Paper No. 88-2828, 1988.

10. Lind, C.A. & Lewis, M.J. The effect of shock/shock interactions on the design of hypersonic inlets. AlAA Paper No. 90-2217, 1990.

11. Brenneis, A. & Wanie, K.M. Navier-Stokes results for hypersonic inlet flows. AlAA Paper No. 91-2472, 1991.

12. Barber, T.J.; Hiett, D. & Fastenberg, S. CFD modeling of the hypersonic inlet starting problem. AIAA Paper No. 2006-0123, 2006.

13. Donde, P.; Marathe, A.G. & Sudhakar, K. Starting in hypersonic intakes. AIAA Paper No. 2006-4510.

14. Yu, D.; Chang, J.; Bao, W. & Zong, X. Optimal classification criterions of hypersonic inlet start/unstart. J. Propulsion Power, 2007, 23(2), 310-16.

15. Emami, S. & Trexler, C. A. Experimental investigation of inlet combustor – isolators for a dual mode scramjet at a Mach number of 4. NASA Technical Report No. NASA TP3502, May 1995.

16. Fluent 6.3 User’s guide, 2006, Fluent Inc, USA

17. Roe, P.L. Characteristic based schemes for the Euler equations. Ann. Rev. Fluid Mechanics, 1986, 18, 337-65.

18. Wilcox, D.C. Multiscale model for turbulent flows. AIAA Journal, 1988, 26(11), 1311-320.

19. Menter, F.R. Performance of popular turbulence models for attached and separated adverse pressure gradient flows. AIAA Journal, 1992, 30(8), 2066-072.

20. Pasha, A.A. & Sinha, K. Simulation of hypersonic shock/turbulent boundary-layer interactions using shock-unsteadiness model. J. Propulsion Power, 2012, 28(1), 46-60.

21. Wilcox, D.C. Formulation of the k −ω turbulence model revisited. AIAA Paper No. 2007-1408.