# Numerical Simulation of Reactive Flows in Overexpanded Supersonic Nozzle with Film Cooling

Authors: Mohamed Sellam; Amer Chpoun
Journal: International Journal of Aerospace Engineering
 (2015)

---

## Abstract

Reignition phenomena occurring in a supersonic nozzle flow may present a crucial safety issue for rocket propulsion systems. These phenomena concern mainly rocket engines which use H2 gas (GH2) in the film cooling device, particularly when the nozzle operates under over expanded flow conditions at sea level or at low altitudes. Consequently, the induced wall thermal loads can lead to the nozzle geometry alteration, which in turn, leads to the appearance of strong side loads that may be detrimental to the rocket engine structural integrity. It is therefore necessary to understand both aerodynamic and chemical mechanisms that are at the origin of these processes. This paper is a numerical contribution which reports results from CFD analysis carried out for supersonic reactive flows in a planar nozzle cooled with GH2 film. Like the experimental observations, CFD simulations showed their ability to highlight these phenomena for the same nozzle flow conditions. Induced thermal load are also analyzed in terms of cooling efficiency and the results already give an idea on their magnitude. It was also shown that slightly increasing the film injection pressure can avoid the reignition phenomena by moving the separation shock towards the nozzle exit section.


---

## Body

## 1. Introduction

One of the major challenges that the aerospace industry continues to face is the continued increase in launchers payload. For example, the latest version of the Vulcain-II engine of the European Ariane 5 ECA launcher is able to put into geostationary orbit a payload of about 10 tons. The goal is however to reach 12 tons in the near future. There are two challenges here: the increase in payload and the performance consolidation in terms of reliability.These challenges promote development of nozzles with higher performances, which are substantially achieved by increasing the nozzle expansion area ratio or by developing new innovative nozzle concepts. The rocket engine nozzles, with high expansion area ratio, are generally optimized for operating at high altitudes. At sea level and at low altitudes, the nozzle operates in overexpanded flow conditions; that is, the ambient pressure is higher than the nozzle exit pressure. The resulting adaptation shock may lead to flow separation, unsteadiness, and shock interaction. The ensuing side loads may be detrimental for both nozzle and other engine components. In addition, these nozzles are designed to expand and accelerate combustion gases at high temperature. To avoid thermal loads, designers adopt several nozzle cooling methods. The most effective one uses the film cooling technique. For example, for Vulcain-II rocket engine, the cooling system is designed in two parts: a dump cooling for the first expansion part of the nozzle and GH2 film cooling for the second part [1, 2].The present study is a contribution to the works initiated by the CNES during the last ten years in the field of innovative nozzle researches. The main goal is to develop new supersonic nozzle concepts for cryogenic rocket engines. One issue that attracted the interest of scientists is that related to the risk of reignition that can occur in the main flow, resulting from combustion of GH2 used as film cooling. This phenomenon may happen when the GH2 film is mixed with the air engulfed into the separated region along the mixing shear layer at high temperature. In this occurrence, the film cooling efficiency is altered and the induced thermal loads may be critical for these nozzles.To understand these phenomena, a campaign of experimental measurements was conducted at the cryogenic LH2-LO2 ONERA/Mascotte test bench facility [3]. The results yet obtained for lowoxidizer/fuelratios (o/f), on the basis of PLIF and Kulite pressure measurements, already constitute a very useful database for understanding the involved physicochemical phenomena. The first series of measurements were performed under o/f ratios of up to 3 and combustion chamber pressure up to 40 bars. The first results obtained during this measurements campaign highlighted a reignition phenomenon occurring in the nozzle main flow, inside the separation area.Given the complexity of these phenomena, numerical simulations performed by different computing codes can give a better interpretation of the experimental results. The present numerical work is a further contribution in this regard and consists in simulating the turbulent reactive flows in the same operating conditions as those performed experimentally.This work was carried out within the framework and with the support of the French research group ATAC (Aérodynamique des Tuyères et Arrières Corps (Nozzles and Afterbodies Aerodynamics)).


## 2. Objectives

The aim of this paper is to perform 2D RANS turbulent reacting flows simulations, for three test cases carried out experimentally in the ONERA/Mascotte test bench. For this purpose, two chemical kinetic schemes and two turbulence models were tested to investigate their relevance in the reignition process inside the nozzle.Basically, the test bench consists of a subscale planar nozzle connected to a cylindrical combustion chamber fed with a cryogenic mixture LOx-LH2 (Figure 1(a)). This hardware version has been developed especially to investigate the flow separation in overexpanded regimes. The combustion chamber operates under controlled total pressure andoxidizer/fuelratios (o/f). The main geometrical characteristics of the nozzle are given in Figure 1(b). Pure hydrogen gas (GH2), serving as film cooling, is injected tangentially into the upper wall of the nozzle through a slightly supersonic injector. Five Kulite pressure transducers are taped along this wall. A schematic longitudinal cross section of the nozzle is shown in Figure 1(b).Figure 1
Mascotte test bench assembly (a) and sketch of longitudinal section of the planar nozzle (b).
(a)
From [22]
(b)The first results [3] obtained by PLIF visualizations performed in this test bench show an increase in concentrations of OH radical near the upper wall of the nozzle’s divergent, in the separation zone. This suggests a reactivation of the combustion inside this area (Figure 2). During these experiments, it was also revealed that there is unsteadiness in the radical OH emission inside the separation zone. The authors related this phenomenon to the pressure fluctuations observed in the combustion chamber. An example of PLIF images obtained by this technique is depicted in Figure 2 [3]. However, no detailed explanations have been given so far about the occurrence of this phenomenon.Figure 2
PLIF images of OH radicals from [3].Parallel to these experiments, steady RANS calculations have been carried out at ONERA/DEFA using the CEDRE Code. In these calculations the reduced chemical kinetic model of Eklund et al. [4] was used for chemical rates. However no process of reignition in the flow has been pointed out from the numerical results [3].From these experimental and numerical conclusions further calculations using new chemical kinetics based on the well-known Evans and Schexnayder model [5] were recommended. This last model was slightly modified by adding two additional reactions, without excessively penalizing the computation time. Furthermore, two RANS turbulence models have also been tested in these calculations. The main objective of this study was to investigate relevance of chemical kinetics and turbulence models, in the context of reignition phenomenon.


## 3. Numerical Code and Equation Formulation

### 3.1. Numerical Code

The calculations presented in this paper were performed using the FASTRAN code. This code was specifically designed for compressible flow studies at high Mach numbers and based on solving the multispecies Reynolds-Averaged Navier-Stokes (RANS) equations with a finite volume formulation. The code offers two upwind differencing schemes with a variety of high order slope limiters to calculate the convective terms in the transport equations. Both explicit and implicit iterative and noniterative time integration schemes are available for steady-state and time accurate flow simulations. The convective fluxes are calculated by means of either flux difference splitting scheme (Roe) or flux vector splitting scheme (van Leer). The code also offers a choice of several turbulence models (k-ε, k-ω, k-ω-SST-Menter, Spalart Allmaras, and Baldwin Lomax). The following sections recall the main fluid physical modelling and the flow field numerical formulations methodology used in this code for solving reactive Navier-Stokes equations.


### 3.2. Thermodynamics Gas Properties

The thermal equation of state for a mixing or reacting gas is given by Dalton’s Law of partial pressures such that(1)p=∑sρsMwsRuT,where p is the static pressure, ρs is the species or mixture density, Ru is the universal gas constant, T is the static temperature, and Mws is the molecular weight of the species s.The caloric equation of state relates the total energy to other gas dynamic variables and gas properties. For calorically perfect gas(2)Et=ρcvT+12ρujuj,where Et is the total energy per unit volume and ρ, cv, T, and uj are the gas density, the heat capacity of the gas mixture at constant volume, the static temperature, and the mass averaged velocity, respectively.The form of the caloric equation of state for a mixing or reacting gas depends on the database used for describing the molecular properties of the species. Two databases are available: the first is based on molecular (or spectroscopic) data for chemical species and the second is based on fifth-order polynomial curve fits for each chemical species [6]. Using molecular properties, (2) can be written as(3)Et=∑sρscv,tr⁡,sT+Δhf,Tr,s0+Ev+12ρujuj,where ρs is the species density, Cv,tr⁡,s is the translational-rotational heat capacity for species s at constant volume, Δhf,Tr,s0 is the heat of formation at reference temperature Tr and pressure for species s, and Ev is the molecular vibrational energy per volume.Using of polynomial curve fits for properties gives two forms of caloric equation of state.For thermal equilibrium(4)Et=∑sρshs-RuTMws+12ρujuj.For thermal nonequilibrium(5)ET=∑sρscv,tr⁡,sT+Δhf,Tr,s0-Trcv,tr⁡,s-RuMws+Eint+12ρujuj,where hs is the sensible enthalpy per unit mass for species s defined as(6)hs=∫TrTcpsdT+Δhf,Tr,s0,where cps is calculated from fifth-order polynomial curve fits, for each chemical species, from Gordon database [6].The viscosity of each fluid speciess, in the case of reactive flows, is calculated by the relationship of Bird et al. [7](7)μs=2.66693·10-6Mw,sTσΩμ,where σ is the characteristic molecular diameter and Ωμ is the viscosity collision integral. The characteristic molecular diameter is based on Lennard-Jones potentials [8]. For the mixture, the viscosity is calculated by the semiempirical relationship of Wilke [9](8)μ=∑sXsμsΦs,where Xs is the mole fraction of specie s and Φs is given by(9)Φs=∑rXr1+μsμrMw,rMw,s1/4281+Mw,sMw,r-1.For the diffusivity vector Js, the mass diffusivity can be represented by either Fick’s law [10] or a binary diffusion model [11]. Consider(10)Js=-ρDs∇Ys,where Ds is the diffusion coefficient and Ys is the species mass fraction. Ds is given by(11)Ds=μρSc,where Sc is the Schmidt number.


### 3.3. Chemical Production

The reactive flow calculation is obtained by solving the flow conservation equations in which one integrates a source termωs expressing the mixture chemical composition variation resulting from chemical reactions.In the approach used for reacting flows, the general finite rate reaction is written as(12)∑s=1nsνsr′Ms⟺∑s=1nsνsr′′Ms,where νsr′ and νsr′′ are the stoichiometric coefficients of the reaction and Ms represents an arbitrary molecule in the reaction. According to Kuo [11], the source term for species s is given by(13)ωs=Mwsνsr′′-νsr′∑s=1nsβsrCs·Kfr∏s=1nsCsαsr′-Kbr∏s=1nsCsαsr′′,where βsr is the coefficient of efficiency of the third body for the reaction r, Cs is the species concentration, and Kfr and Kbr are forward and backward reaction rates of a reaction r, respectively. The concentration powers αsr′ and αsr′′ are identical to νsr′ and νsr′′, respectively, for most applications, particularly for chemical kinetic reaction governed by Arrhenius rates of reaction(14)Kfr=αfrTβfr·e-Ear/(R·T),where αfr, βfr, and Ear/R must be specified for each reaction under investigation.


### 3.4. Flow Field Numerical Method

Basically, the conservation equations with appropriate closure models are expressed in vector form as(15)∂Q∂t+∇·F⃑C-∇·F⃑D=S˙.In this expression, FC and FD represent the convective and diffusive fluxes, respectively, such as (16)Q=Eintρ1ρ2⋮⋮ρnsρuρvρwEt,FC=Eintuρ1uρ2u⋮⋮ρnsuρu2+pρuvρuwEt+pu,S˙=ω˙intω˙1ω˙2⋮⋮ω˙ns0000,FD=∑s=1nsev,sρsUs+kint∂Tint∂x-ρ1U1-ρ2U2⋮⋮-ρnsUsτxx-23ρkτxyτxzuτxx+vτxy+wτxz+k+kt∂T∂x+kint∂Tint∂x+∑hsρsUs.In this system of equations, u,v, and w are the velocity components (w=0 for 2D calculations), ρ is the mixture density, and τij is the shear stress tensor. ρs and Us are, respectively, the species density and species velocity. One can note that, for calorically perfect gas, only one “species” is tracked such that ns=1, ρs=ρ, and Us=0. Note that, for thermal equilibrium calculation, all terms relating the contribution of vibrational internal energy (Eint,Tint,evs,…) are no longer required.


## 3.1. Numerical Code

The calculations presented in this paper were performed using the FASTRAN code. This code was specifically designed for compressible flow studies at high Mach numbers and based on solving the multispecies Reynolds-Averaged Navier-Stokes (RANS) equations with a finite volume formulation. The code offers two upwind differencing schemes with a variety of high order slope limiters to calculate the convective terms in the transport equations. Both explicit and implicit iterative and noniterative time integration schemes are available for steady-state and time accurate flow simulations. The convective fluxes are calculated by means of either flux difference splitting scheme (Roe) or flux vector splitting scheme (van Leer). The code also offers a choice of several turbulence models (k-ε, k-ω, k-ω-SST-Menter, Spalart Allmaras, and Baldwin Lomax). The following sections recall the main fluid physical modelling and the flow field numerical formulations methodology used in this code for solving reactive Navier-Stokes equations.


## 3.2. Thermodynamics Gas Properties

The thermal equation of state for a mixing or reacting gas is given by Dalton’s Law of partial pressures such that(1)p=∑sρsMwsRuT,where p is the static pressure, ρs is the species or mixture density, Ru is the universal gas constant, T is the static temperature, and Mws is the molecular weight of the species s.The caloric equation of state relates the total energy to other gas dynamic variables and gas properties. For calorically perfect gas(2)Et=ρcvT+12ρujuj,where Et is the total energy per unit volume and ρ, cv, T, and uj are the gas density, the heat capacity of the gas mixture at constant volume, the static temperature, and the mass averaged velocity, respectively.The form of the caloric equation of state for a mixing or reacting gas depends on the database used for describing the molecular properties of the species. Two databases are available: the first is based on molecular (or spectroscopic) data for chemical species and the second is based on fifth-order polynomial curve fits for each chemical species [6]. Using molecular properties, (2) can be written as(3)Et=∑sρscv,tr⁡,sT+Δhf,Tr,s0+Ev+12ρujuj,where ρs is the species density, Cv,tr⁡,s is the translational-rotational heat capacity for species s at constant volume, Δhf,Tr,s0 is the heat of formation at reference temperature Tr and pressure for species s, and Ev is the molecular vibrational energy per volume.Using of polynomial curve fits for properties gives two forms of caloric equation of state.For thermal equilibrium(4)Et=∑sρshs-RuTMws+12ρujuj.For thermal nonequilibrium(5)ET=∑sρscv,tr⁡,sT+Δhf,Tr,s0-Trcv,tr⁡,s-RuMws+Eint+12ρujuj,where hs is the sensible enthalpy per unit mass for species s defined as(6)hs=∫TrTcpsdT+Δhf,Tr,s0,where cps is calculated from fifth-order polynomial curve fits, for each chemical species, from Gordon database [6].The viscosity of each fluid speciess, in the case of reactive flows, is calculated by the relationship of Bird et al. [7](7)μs=2.66693·10-6Mw,sTσΩμ,where σ is the characteristic molecular diameter and Ωμ is the viscosity collision integral. The characteristic molecular diameter is based on Lennard-Jones potentials [8]. For the mixture, the viscosity is calculated by the semiempirical relationship of Wilke [9](8)μ=∑sXsμsΦs,where Xs is the mole fraction of specie s and Φs is given by(9)Φs=∑rXr1+μsμrMw,rMw,s1/4281+Mw,sMw,r-1.For the diffusivity vector Js, the mass diffusivity can be represented by either Fick’s law [10] or a binary diffusion model [11]. Consider(10)Js=-ρDs∇Ys,where Ds is the diffusion coefficient and Ys is the species mass fraction. Ds is given by(11)Ds=μρSc,where Sc is the Schmidt number.


## 3.3. Chemical Production

The reactive flow calculation is obtained by solving the flow conservation equations in which one integrates a source termωs expressing the mixture chemical composition variation resulting from chemical reactions.In the approach used for reacting flows, the general finite rate reaction is written as(12)∑s=1nsνsr′Ms⟺∑s=1nsνsr′′Ms,where νsr′ and νsr′′ are the stoichiometric coefficients of the reaction and Ms represents an arbitrary molecule in the reaction. According to Kuo [11], the source term for species s is given by(13)ωs=Mwsνsr′′-νsr′∑s=1nsβsrCs·Kfr∏s=1nsCsαsr′-Kbr∏s=1nsCsαsr′′,where βsr is the coefficient of efficiency of the third body for the reaction r, Cs is the species concentration, and Kfr and Kbr are forward and backward reaction rates of a reaction r, respectively. The concentration powers αsr′ and αsr′′ are identical to νsr′ and νsr′′, respectively, for most applications, particularly for chemical kinetic reaction governed by Arrhenius rates of reaction(14)Kfr=αfrTβfr·e-Ear/(R·T),where αfr, βfr, and Ear/R must be specified for each reaction under investigation.


## 3.4. Flow Field Numerical Method

Basically, the conservation equations with appropriate closure models are expressed in vector form as(15)∂Q∂t+∇·F⃑C-∇·F⃑D=S˙.In this expression, FC and FD represent the convective and diffusive fluxes, respectively, such as (16)Q=Eintρ1ρ2⋮⋮ρnsρuρvρwEt,FC=Eintuρ1uρ2u⋮⋮ρnsuρu2+pρuvρuwEt+pu,S˙=ω˙intω˙1ω˙2⋮⋮ω˙ns0000,FD=∑s=1nsev,sρsUs+kint∂Tint∂x-ρ1U1-ρ2U2⋮⋮-ρnsUsτxx-23ρkτxyτxzuτxx+vτxy+wτxz+k+kt∂T∂x+kint∂Tint∂x+∑hsρsUs.In this system of equations, u,v, and w are the velocity components (w=0 for 2D calculations), ρ is the mixture density, and τij is the shear stress tensor. ρs and Us are, respectively, the species density and species velocity. One can note that, for calorically perfect gas, only one “species” is tracked such that ns=1, ρs=ρ, and Us=0. Note that, for thermal equilibrium calculation, all terms relating the contribution of vibrational internal energy (Eint,Tint,evs,…) are no longer required.


## 4. Chemical Kinetic Models

To investigate the pertinence of the chemical reactions associated with the reacting mixture issued from the combustion chamber, it was expedient to test the most suitable kinetic schemes for reactive H2-O2 flow. Two kinetic schemes were selected for this study: the modified Evans-Schexnayder model and Eklund’s kinetic model, commonly used by ONERA and CNES.


### 4.1. Modified Evans-Schexnayder

This model is initially based on 7 species [O2, H2, OH, H2O, N2, and O] and 8 chemical reactions’ scheme [12]. In this system, N2 operates as the third body and does not dissociate. This mechanism has been widely used for simulation in supersonic and hypersonic flows, particularly in the case of combustion initiation around obstacles or in scramjets [13–15]. This model was proved to be less expensive in terms of computation time; however it presents weakness in modeling the self-ignition delay (induction time) and in estimating the reaction heat release. This is mainly due to the absence of hydroperoxyl radical (HO2) in this scheme. Indeed, studies have shown that fast three body recombination reactions, involving the radical HO2, have been identified as major contributor in the heat release process during the combustion of hydrogen with air [16].To overcome this deficit, two reactions taken from the model of Rogers and Chinitz [14], involving this radical, have been added to the original Evans model, adding no substantial computation time(17)H+O2+M⟺HO2+M(18)H+HO2⟺OH+OHAnother insufficiency attributed to this model is its autoignition delay, which is relatively long, especially for reactions at low temperature (≈1000 K). This problem is related to the absence of hydrogen peroxide (H2O2) in the model.Adding more reactions involving this species to correct this deficiency substantially complicates the model. An alternative solution would be to increase the production rate of Reaction 7 of the original model(19)O2+H⟺OH+OIndeed, this reaction has been identified as important in the case of inflammation at low temperatures [14]. Initially, the forward rate equation for this reaction is expressed as [5](20)kf=2.2·1014exp⁡-8455T.This value was obtained with an accuracy of 50% for a temperature range of 300 to 2000 K. By multiplying the coefficient αf,r by 3, the rate of hydroperoxyl radical production OH is increased which leads to reduction in the ignition delay [12]. The corresponding ignition delay becomes comparable to those obtained by more complex chemical kinetic models with more reactions. Finally, the modified Evans-Schexnayder kinetic model, with ten chemical reactions, is given in Table 1.Table 1
Modified Evans-Schexnayder reaction model,kf in cm3/mole·s. M is the third body with an efficiency = 1 for all species and θi = Ea/R.
Number
Reaction
α
i
β
i
θ
i
(K)
1
H2 + M          ⇔  H+ H + M
5.5 × 1018
−1.0
51987
2
O2 + M          ⇔   O + O + M
7.2 × 1018
−1.0
59340
3
H2O + M       ⇔   OH + H + M
5.2 × 1021
−1.5
59386
4
OH + M⇔   O + H + M
8.5 × 1018
−1.0
50830
5
H2O +O         ⇔   OH + OH
5.8 × 1013
0.0
9059
6
H2O + H        ⇔   OH + H2
8.4 × 1013
0.0
10116
7
O2 + H           ⇔   OH + O
6.6 × 1014
0.0
8455
8
H2 + O           ⇔   OH + H
5.5 × 1013
0.0
5586
9
H + O2 + M   ⇔   HO2 + M
2.3 × 1016
0.0
−403
10
H   + HO2⇔   OH + OH
2.4 × 1014
0.0
950Figure3 depicts the results obtained for validation of this kinetic model, in the case of H2-O2 combustion. The results are presented in terms of pressure and temperature rise, H2 consumption, and OH and H2O formation, from one-dimensional combustion simulation. In Figures 3(a) and 3(b), the initial and the modified Evans-Schexnayder model results are compared to those from more complex kinetic schemes of Rogers and Chinitz [17], Drummond [18], and Bitker and Scullin [19], respectively. The results clearly highlight the relevance of the added specific reactions on the ignition time delay [12].Figure 3
Results of one-dimensional combustion tube simulation obtained with some reaction models compared to initial and modified Evans-Schexnayder reaction model [12].
(a)
Results from the initial model
(b)
Results from the modified model


### 4.2. Eklund Model

This reaction simplified scheme, proposed by Eklund et al. [4] and implemented on both CEDRE and CPS codes, has been widely used by ONERA and CNES [3, 20] for nozzle reactive flow studies. This scheme consists of 7 reversible reactions and 6 chemical species [O2, H2, OH, H2O, O, and H]. As can be seen in Table 2, this scheme does not involve any third body reaction which can present an advantage in terms of computing time.Table 2
Reduced Eklund reaction model [4]. kf   in cm3/mole·s.
Number
Reaction
α
f
i
β
f
i
θ
f
i (K)
α
b
i
β
b
i
θ
b
i (K)
1
H2 + O2⇔  2OH
1.7 × 107
0
24240
27810
0.363
14620
2
H + O2⇔  O + OH
1.42 × 108
0
8258
4.17 × 105
0.3905
−205
3
OH + H2⇔  H2O + H
3.16 × 101
1.8
1526
895.3
1.572
9346
4
O + H2⇔  OH + H
2.07 × 105
0.0
6924
1.153 × 108
−0.02726
5769
5
OH + OH⇔  H2O + O
5.5 × 1008
0.0
1511
2.797 × 1010
−2.2004
10490
6
H + OH⇔   H2O
1.105 × 1012
−2.0
0.0
2.365 × 106
1.313
50910
7
H + H⇔  H2
3.265 × 1007
−1.0
0.0
3.08 × 106
0.664
47201


## 4.1. Modified Evans-Schexnayder

This model is initially based on 7 species [O2, H2, OH, H2O, N2, and O] and 8 chemical reactions’ scheme [12]. In this system, N2 operates as the third body and does not dissociate. This mechanism has been widely used for simulation in supersonic and hypersonic flows, particularly in the case of combustion initiation around obstacles or in scramjets [13–15]. This model was proved to be less expensive in terms of computation time; however it presents weakness in modeling the self-ignition delay (induction time) and in estimating the reaction heat release. This is mainly due to the absence of hydroperoxyl radical (HO2) in this scheme. Indeed, studies have shown that fast three body recombination reactions, involving the radical HO2, have been identified as major contributor in the heat release process during the combustion of hydrogen with air [16].To overcome this deficit, two reactions taken from the model of Rogers and Chinitz [14], involving this radical, have been added to the original Evans model, adding no substantial computation time(17)H+O2+M⟺HO2+M(18)H+HO2⟺OH+OHAnother insufficiency attributed to this model is its autoignition delay, which is relatively long, especially for reactions at low temperature (≈1000 K). This problem is related to the absence of hydrogen peroxide (H2O2) in the model.Adding more reactions involving this species to correct this deficiency substantially complicates the model. An alternative solution would be to increase the production rate of Reaction 7 of the original model(19)O2+H⟺OH+OIndeed, this reaction has been identified as important in the case of inflammation at low temperatures [14]. Initially, the forward rate equation for this reaction is expressed as [5](20)kf=2.2·1014exp⁡-8455T.This value was obtained with an accuracy of 50% for a temperature range of 300 to 2000 K. By multiplying the coefficient αf,r by 3, the rate of hydroperoxyl radical production OH is increased which leads to reduction in the ignition delay [12]. The corresponding ignition delay becomes comparable to those obtained by more complex chemical kinetic models with more reactions. Finally, the modified Evans-Schexnayder kinetic model, with ten chemical reactions, is given in Table 1.Table 1
Modified Evans-Schexnayder reaction model,kf in cm3/mole·s. M is the third body with an efficiency = 1 for all species and θi = Ea/R.
Number
Reaction
α
i
β
i
θ
i
(K)
1
H2 + M          ⇔  H+ H + M
5.5 × 1018
−1.0
51987
2
O2 + M          ⇔   O + O + M
7.2 × 1018
−1.0
59340
3
H2O + M       ⇔   OH + H + M
5.2 × 1021
−1.5
59386
4
OH + M⇔   O + H + M
8.5 × 1018
−1.0
50830
5
H2O +O         ⇔   OH + OH
5.8 × 1013
0.0
9059
6
H2O + H        ⇔   OH + H2
8.4 × 1013
0.0
10116
7
O2 + H           ⇔   OH + O
6.6 × 1014
0.0
8455
8
H2 + O           ⇔   OH + H
5.5 × 1013
0.0
5586
9
H + O2 + M   ⇔   HO2 + M
2.3 × 1016
0.0
−403
10
H   + HO2⇔   OH + OH
2.4 × 1014
0.0
950Figure3 depicts the results obtained for validation of this kinetic model, in the case of H2-O2 combustion. The results are presented in terms of pressure and temperature rise, H2 consumption, and OH and H2O formation, from one-dimensional combustion simulation. In Figures 3(a) and 3(b), the initial and the modified Evans-Schexnayder model results are compared to those from more complex kinetic schemes of Rogers and Chinitz [17], Drummond [18], and Bitker and Scullin [19], respectively. The results clearly highlight the relevance of the added specific reactions on the ignition time delay [12].Figure 3
Results of one-dimensional combustion tube simulation obtained with some reaction models compared to initial and modified Evans-Schexnayder reaction model [12].
(a)
Results from the initial model
(b)
Results from the modified model


## 4.2. Eklund Model

This reaction simplified scheme, proposed by Eklund et al. [4] and implemented on both CEDRE and CPS codes, has been widely used by ONERA and CNES [3, 20] for nozzle reactive flow studies. This scheme consists of 7 reversible reactions and 6 chemical species [O2, H2, OH, H2O, O, and H]. As can be seen in Table 2, this scheme does not involve any third body reaction which can present an advantage in terms of computing time.Table 2
Reduced Eklund reaction model [4]. kf   in cm3/mole·s.
Number
Reaction
α
f
i
β
f
i
θ
f
i (K)
α
b
i
β
b
i
θ
b
i (K)
1
H2 + O2⇔  2OH
1.7 × 107
0
24240
27810
0.363
14620
2
H + O2⇔  O + OH
1.42 × 108
0
8258
4.17 × 105
0.3905
−205
3
OH + H2⇔  H2O + H
3.16 × 101
1.8
1526
895.3
1.572
9346
4
O + H2⇔  OH + H
2.07 × 105
0.0
6924
1.153 × 108
−0.02726
5769
5
OH + OH⇔  H2O + O
5.5 × 1008
0.0
1511
2.797 × 1010
−2.2004
10490
6
H + OH⇔   H2O
1.105 × 1012
−2.0
0.0
2.365 × 106
1.313
50910
7
H + H⇔  H2
3.265 × 1007
−1.0
0.0
3.08 × 106
0.664
47201


## 5. Results and Discussions

### 5.1. Initial Conditions and Implementation of Calculations

The calculations were performed over a computational domain which includes the nozzle, the injector, and the outside experimental environment. The calculation was initiated at the nozzle inlet using the same initial data as the experimental test conditions described below. No-slip conditions along the nozzle walls were assumed. For the outlet condition, 1 bar fixed pressure is applied at the downstream exit section. Adiabatic no-slip conditions are imposed for the rest of the surrounding block boundaries.Combustion temperature and species mass fractions for the cryogenic LH2-LO2 combustion products, at desired operating pressure chamber and o/f ratio, are obtained from separate calculations using the CEA thermochemical code [21].In order to perform 2D CFD simulations of the nozzle’s flow field, finite volume grids have been constructed using algebraic grid generator software. Multiblock structured grids have been used in this calculation. The computational domain includes the convergent-divergent parts of the nozzle, the injector, the zone downstream, and the area located on both sides of the nozzle. The dimensions of the downstream block are more than70h in both x and y directions, h=0.028 m being the nozzle height at its exit section. A part of the meshing used in the computational domain, representing the nozzle and part of the outer domain, is illustrated in Figure 4(a).Figure 4
Grids system used for CFD calculation: (a) nozzle with part of the surrounding area, (b) nozzle with a zoom of the region near the injection slot.The right part of Figure4(b) depicts the mesh system used for the injection slot and the nearby region. It also depicts the mesh refinement near the wall. To achieve the CFD calculation, up to 40000 cells were required. The grid refinement near the wall, in the y direction, such as the first cell-height leads to y+=1, where y+ is the dimensionless distance defined as(21)y+=yντwρwwith τw, ρw, and ν being the shear stress at the wall, the density, and the kinematic viscosity, respectively. An average of 25 cells in the normal direction of the nozzle wall was used to resolve the boundary layer thickness. Furthermore, a preliminary study of the mesh refinement sensitivity was beforehand performed before achieving the final calculations.The fluid is considered to be an ideal gas mixture in thermal equilibrium. The mixture heat capacity at constant pressurecp(T) is calculated as a function of temperature from the calorific capacity of each species, expressed in polynomial form. The Schmidt number Sc and the turbulent Prandtl number are set equal to 0.9 for all calculations.Boundary conditions for turbulence model are also needed. For thek-ε turbulence model, the turbulent kinetic energy k and the turbulent dissipation rate ε required at the inlets for both nozzle and injector are calculated on the basis of an estimated turbulence intensity of 5% for this nozzle [3].Moreover, the Spalart Allmaras turbulence model requires the kinematic turbulence viscosityνt of the mixture. For Reynolds numbers Re>105, which is the case here, a value of the turbulent viscosity ratio, namely, the ratio of turbulent viscosity to laminar viscosity (β=νt/ν), up to 100 is a reasonable estimate.In this study, the fluxes are evaluated at each time step using the Roe scheme associated with the second-order spatial accuracy MINMOD flux limiter. The time integration is performed by a fully implicit scheme with a local time step based on CFL number ramped from 0.1 to 0.5 within the first 1000 iteration cycles. Convergence to10-6 residuals was usually attained after 20,000 to 25,000 iterations.Finally the code is based on the finite rate mechanism in calculating the chemical composition. Each chemical reaction obeys the law of mass action; hence the backward rate coefficient is calculated on the basis of the equilibrium constant for each reversible reaction.The results obtained by the numerical simulations are presented in terms of flow Mach number field, temperature, and OH mass fraction fields. In addition to the PLIF pictures, five experimental wall pressure values are available. These pressures are measured along the cooled upper nozzle wall and will be compared to the calculated values in the same plots. Each test case is simulated using two chemical kinetic schemes, namely, “the modified Evans-Schexnaiyder model” and “the Eklund model.” Two turbulence models are also considered in this calculation, the standard two-equationk-ε and the one-equation Spalart-Allmaras model.The operating conditions of the tested cases and the initial data for the numerical simulations inputs are summarized in Table3.Table 3
Operating test cases conditions and initial inputs for CFD calculation.
Test cases
1
2
3
4*
Measured quantities
LO2 mass flow
(g/s)
44.9
46.1
67.2
x
LH2 mass flow
(g/s)
24.2
14.2
35.9
x
o/f
1.86
3.25
1.87
6
Conditions at the H2 cooling injector
P
0
j cooling (H2)
(bar)
3.1
3.2
4.3
4.3
T
0
J cooling (H2)
(K)
291.35
294.25
295.65
294
Conditions at the combustion chamber
P
0 (combustion chamber)
(bar)
25.9
22.2
36.5
25.9
T
0 (calculated)
(K)
1693
2577
1945
3372
Calculated species mass fractions
Y
H
2
O
%
0.73229
0.85903
0.73366
0.87603
Y
OH
%
0.
0.00190
0.
0.06548
Y
H
2
%
0.26771
0.13840
0.26634
0.03777
Y
O
%
0.
0.00001
0.
0.00580
Y
H
%
0.
0.00065
0.
0.00325
Y
O
2
%
0.
0.00001
0.
0.01168
∗This case is performed only numerically.


### 5.2. Aerodynamic Field

To understand the flow topology, we plotted the aerodynamic flow fields in terms of Mach number of the different tested cases. As depicted in Figure5, for each test case, both turbulence models used in the computation predict a similar flow pattern. Indeed, in such nozzle regimes (overexpanded conditions), a separation shock is created inside the nozzle extension to adapt the flow to the outside pressure. This configuration is known as the Free Separation Shock (FSS). As sketched in Figure 5(a), in all configurations, the separation shock impacts the opposite wall and leads the boundary layer to separate similarly as in the so-called Shock Wave Boundary Layer Interaction (SWBLI). The other point the flow topology reveals is that relating to the shock position. As expected, this position depends on both nozzle pressure ratio NPR, for a given turbulence model, and turbulence model, for a given NPR, as illustrated in the comparisons given in Figure 5. This point will be emphasized with more details in wall pressure discussion below.Figure 5
Calculated Mach number fields for Cases 1, 2, and 3.
(a)
Case 1-k-ε-Evans (NPR = 25.9)
(b)
Case 1-Spal-Evans (NPR = 25.9)
(c)
Case 2-k-ε-Evans (NPR = 22.2)
(d)
Case 2-Spal-Evans (NPR = 22.2)
(e)
Case 3-k-ε-Evans (NPR = 36.5)
(f)
Case 3-Spal-Evans (NPR = 36.5)Moreover, it should be noted that an increase in the GH2 injection pressure, used in the film cooling, as in the calculated Case 4 pushes back the separation shock downstream the injection slot towards the nozzle exit section. This result is depicted in Figure 6 in terms of calculated pressure profiles for Case 1 (NPR = 25.9, P0j = 3 bars) and for Case 4 (NPR = 25.9, P0j = 4.3 bars). As shown, the relative displacement of shock, determined with respect to the length l of the nozzle part concerned with the fim cooling (Δx/l), is about 36.3%. A similar trend has been reported in previous works [23], where similar film injection is used to control the flow separation in overexpanded supersonic nozzles.Figure 6
Comparison between Case 4 and Case 1k-ε and Evans.
(a)
(b)


### 5.3. Wall Pressure

The calculated pressure profiles are compared to the experimental pressure measured by means of Kulite transducers. Figure7 depicts the pressure distribution along the upper wall of the nozzle’s divergent, downstream of the injection slot.Figure 7
Calculated and measured wall pressure in the zone downstream the injection slot.
(a)
Calculated and measured wall pressure  (Case 1)
(b)
Calculated and measured wall pressure (Case 2)
(c)
Comparison between calculated Cases 1 and 4The results in Figures7(a) and 7(b) show that the turbulence models used in this work were able to predict with good accuracy the pressure rise in the separation zone, for Test Case 1 (NPR = 25.9) and Test Case 2 (NPR = 22.2). The numerical results also show that the reignition phenomenon significantly affects the flow field in the separation zone. Due to the combustion process, the separation area becomes larger than the corresponding frozen cases and the separation shock moves upstream. This result highlights the chemical reactions influence on the separation shock displacement and therefore on the pressure rise within the recirculation zone. This result can be connected to the boundary layer, which becomes thicker in reactive flow and easier to detach in response to an adverse pressure gradient. This finding is quite consistent with previous works such as in [24]. Note that Test Case 3 is not represented here in terms of comparison between the calculated and measured pressure profiles. Indeed, for this particular test case, large pressure fluctuations in the combustion chamber were recorded during the course of experiments. This led to great separation shock instabilities and large discrepancy in measurements.Additional Test Case 4 was carried out only by numerical means, with anoxidizer/fuelratio (o/f=6), close to that of real rocket engines, functioning with LO2-LH2. The nozzle pressure ratio NPR is fixed the same way as for Test Case 1; that is, NPR = 25.9, while the total pressure of the H2 film cooling is set to P0j=4.3 bars. As can be seen in Figure 7(c), augmenting the film cooling total pressure leads the separation shock to move downstream as well as the reignition attachment point consequently. This possibility can be used, as will be mentioned below, as means to avoid the reignition of the mixture inside the nozzle.


### 5.4. Chemical Production

To highlight the reignition phenomenon in the separation zone, we plotted the flow OH mass fraction. The results obtained show that OH radical concentrations are at levels high enough to readily suggest a reignition of the H2-Air mixture, in the separated region. The examples illustrated in Figures 8(a) and 8(b) suggest that the entrained fresh oxygen, from the outside into the separation zone, leads to the inflammation GH2 used in the film cooling across the mixing shear layer. A like-diffusion flame appears and extends downstream along the shear layer. A schematic of the flow configuration, with an attached flame to the nozzle upper wall in the separation area, is depicted in Figure 9. It is worth noting that these results were observed for both turbulence models, with modified Evans-Schexnayder chemical kinetics model but not with Eklund model. This can be explained by the temperature levels in the nozzle main flow, across the shear layer, which are not high enough to reactivate the chemical reactions with the above mentioned reaction model. On the other hand, when the nozzle operates at pressure ratios close to the adaptation, as in Case 3 (NPR = 36.5), no reignition process has been observed in the flow inside the nozzle. Indeed, as mentioned above, the separation zone size is too small that it does not allow the full completion of the combustion process within this area. In this case, the combustion takes place in the outer caisson, far away from the nozzle exit section. This situation is depicted in Figure 10 for both OH concentrations and temperature plots.Figure 8
Calculated OH radicals and H2 concentrations: Evans and k-ε calculations.
(a)
Case 1
(b)
Case 4Figure 9
Flow configuration as deduced from CFD calculations.Figure 10
Calculated OH concentration and flow temperature Case 3, Evans andk-ε calculation. The combustion takes place outside the nozzle.
(a)
(b)


### 5.5. Thermal Loads

The film cooling efficiency is evaluated by a widely used parameterη, describing the difference between the film cooled wall temperature, at each point, and the combustion chamber temperature, relative to the difference between the coolant temperature and the combustion chamber temperature. The cooling efficiency is expressed by the nondimensional efficiency η defined as(22)η=Tw-TcTf-Tc,where Tw is the calculated wall temperature, Tc the calculated combustion chamber temperature, and Tf the injected film temperature.The cooling efficiency along the wall depends, as before, on whether the nozzle flow is fully expanded or overexpanded. Figures11(a) and 11(b) represent an example of wall temperature calculated in the Case 1 (NPR = 25.9) and Case 2 (NPR = 22.2), respectively, where the nozzle operates in overexpanded conditions. The combustion of injected GH2 leads to an increase in wall temperature in the region where the flame is attached. The wall temperature plots are characterized by marked peaks up to 1600 K, compared to frozen flow calculations. These results represent a further confirmation of the hypotheses raised previously about the reignition phenomenon. As can be seen, these peaks are shifted according to the observations made above about the separation zone. In this case the cooling efficiency falls to values close to 0.05 and 0.4, respectively, for Case 1 and Case 2 (Figures 12(a) and 12(b)), where the temperature of the attached flame is close to the adiabatic temperature of the combustion chamber, particularly for Case 1. Meanwhile Case 3 shows a typical behavior of film cooling applications as it was observed in several works [2]. In this case the temperature of the film cooling grows gradually by mixing with the hot main flow gas (Figure 11(c)) and the corresponding efficiency (η>0.8) is within the range of the expected values (Figure 12(c)).Figure 11
Wall temperature profiles with reignition, Cases 1 (a) and 2 (b), and without reignition for Case 3 (c).
(a)
(b)
(c)Figure 12
Cooling film efficiency as function of distance from the injection slot for (a) Case 1, (b) Case 2, and (c) Case 3.
(a)
(b)
(c)


## 5.1. Initial Conditions and Implementation of Calculations

The calculations were performed over a computational domain which includes the nozzle, the injector, and the outside experimental environment. The calculation was initiated at the nozzle inlet using the same initial data as the experimental test conditions described below. No-slip conditions along the nozzle walls were assumed. For the outlet condition, 1 bar fixed pressure is applied at the downstream exit section. Adiabatic no-slip conditions are imposed for the rest of the surrounding block boundaries.Combustion temperature and species mass fractions for the cryogenic LH2-LO2 combustion products, at desired operating pressure chamber and o/f ratio, are obtained from separate calculations using the CEA thermochemical code [21].In order to perform 2D CFD simulations of the nozzle’s flow field, finite volume grids have been constructed using algebraic grid generator software. Multiblock structured grids have been used in this calculation. The computational domain includes the convergent-divergent parts of the nozzle, the injector, the zone downstream, and the area located on both sides of the nozzle. The dimensions of the downstream block are more than70h in both x and y directions, h=0.028 m being the nozzle height at its exit section. A part of the meshing used in the computational domain, representing the nozzle and part of the outer domain, is illustrated in Figure 4(a).Figure 4
Grids system used for CFD calculation: (a) nozzle with part of the surrounding area, (b) nozzle with a zoom of the region near the injection slot.The right part of Figure4(b) depicts the mesh system used for the injection slot and the nearby region. It also depicts the mesh refinement near the wall. To achieve the CFD calculation, up to 40000 cells were required. The grid refinement near the wall, in the y direction, such as the first cell-height leads to y+=1, where y+ is the dimensionless distance defined as(21)y+=yντwρwwith τw, ρw, and ν being the shear stress at the wall, the density, and the kinematic viscosity, respectively. An average of 25 cells in the normal direction of the nozzle wall was used to resolve the boundary layer thickness. Furthermore, a preliminary study of the mesh refinement sensitivity was beforehand performed before achieving the final calculations.The fluid is considered to be an ideal gas mixture in thermal equilibrium. The mixture heat capacity at constant pressurecp(T) is calculated as a function of temperature from the calorific capacity of each species, expressed in polynomial form. The Schmidt number Sc and the turbulent Prandtl number are set equal to 0.9 for all calculations.Boundary conditions for turbulence model are also needed. For thek-ε turbulence model, the turbulent kinetic energy k and the turbulent dissipation rate ε required at the inlets for both nozzle and injector are calculated on the basis of an estimated turbulence intensity of 5% for this nozzle [3].Moreover, the Spalart Allmaras turbulence model requires the kinematic turbulence viscosityνt of the mixture. For Reynolds numbers Re>105, which is the case here, a value of the turbulent viscosity ratio, namely, the ratio of turbulent viscosity to laminar viscosity (β=νt/ν), up to 100 is a reasonable estimate.In this study, the fluxes are evaluated at each time step using the Roe scheme associated with the second-order spatial accuracy MINMOD flux limiter. The time integration is performed by a fully implicit scheme with a local time step based on CFL number ramped from 0.1 to 0.5 within the first 1000 iteration cycles. Convergence to10-6 residuals was usually attained after 20,000 to 25,000 iterations.Finally the code is based on the finite rate mechanism in calculating the chemical composition. Each chemical reaction obeys the law of mass action; hence the backward rate coefficient is calculated on the basis of the equilibrium constant for each reversible reaction.The results obtained by the numerical simulations are presented in terms of flow Mach number field, temperature, and OH mass fraction fields. In addition to the PLIF pictures, five experimental wall pressure values are available. These pressures are measured along the cooled upper nozzle wall and will be compared to the calculated values in the same plots. Each test case is simulated using two chemical kinetic schemes, namely, “the modified Evans-Schexnaiyder model” and “the Eklund model.” Two turbulence models are also considered in this calculation, the standard two-equationk-ε and the one-equation Spalart-Allmaras model.The operating conditions of the tested cases and the initial data for the numerical simulations inputs are summarized in Table3.Table 3
Operating test cases conditions and initial inputs for CFD calculation.
Test cases
1
2
3
4*
Measured quantities
LO2 mass flow
(g/s)
44.9
46.1
67.2
x
LH2 mass flow
(g/s)
24.2
14.2
35.9
x
o/f
1.86
3.25
1.87
6
Conditions at the H2 cooling injector
P
0
j cooling (H2)
(bar)
3.1
3.2
4.3
4.3
T
0
J cooling (H2)
(K)
291.35
294.25
295.65
294
Conditions at the combustion chamber
P
0 (combustion chamber)
(bar)
25.9
22.2
36.5
25.9
T
0 (calculated)
(K)
1693
2577
1945
3372
Calculated species mass fractions
Y
H
2
O
%
0.73229
0.85903
0.73366
0.87603
Y
OH
%
0.
0.00190
0.
0.06548
Y
H
2
%
0.26771
0.13840
0.26634
0.03777
Y
O
%
0.
0.00001
0.
0.00580
Y
H
%
0.
0.00065
0.
0.00325
Y
O
2
%
0.
0.00001
0.
0.01168
∗This case is performed only numerically.


## 5.2. Aerodynamic Field

To understand the flow topology, we plotted the aerodynamic flow fields in terms of Mach number of the different tested cases. As depicted in Figure5, for each test case, both turbulence models used in the computation predict a similar flow pattern. Indeed, in such nozzle regimes (overexpanded conditions), a separation shock is created inside the nozzle extension to adapt the flow to the outside pressure. This configuration is known as the Free Separation Shock (FSS). As sketched in Figure 5(a), in all configurations, the separation shock impacts the opposite wall and leads the boundary layer to separate similarly as in the so-called Shock Wave Boundary Layer Interaction (SWBLI). The other point the flow topology reveals is that relating to the shock position. As expected, this position depends on both nozzle pressure ratio NPR, for a given turbulence model, and turbulence model, for a given NPR, as illustrated in the comparisons given in Figure 5. This point will be emphasized with more details in wall pressure discussion below.Figure 5
Calculated Mach number fields for Cases 1, 2, and 3.
(a)
Case 1-k-ε-Evans (NPR = 25.9)
(b)
Case 1-Spal-Evans (NPR = 25.9)
(c)
Case 2-k-ε-Evans (NPR = 22.2)
(d)
Case 2-Spal-Evans (NPR = 22.2)
(e)
Case 3-k-ε-Evans (NPR = 36.5)
(f)
Case 3-Spal-Evans (NPR = 36.5)Moreover, it should be noted that an increase in the GH2 injection pressure, used in the film cooling, as in the calculated Case 4 pushes back the separation shock downstream the injection slot towards the nozzle exit section. This result is depicted in Figure 6 in terms of calculated pressure profiles for Case 1 (NPR = 25.9, P0j = 3 bars) and for Case 4 (NPR = 25.9, P0j = 4.3 bars). As shown, the relative displacement of shock, determined with respect to the length l of the nozzle part concerned with the fim cooling (Δx/l), is about 36.3%. A similar trend has been reported in previous works [23], where similar film injection is used to control the flow separation in overexpanded supersonic nozzles.Figure 6
Comparison between Case 4 and Case 1k-ε and Evans.
(a)
(b)


## 5.3. Wall Pressure

The calculated pressure profiles are compared to the experimental pressure measured by means of Kulite transducers. Figure7 depicts the pressure distribution along the upper wall of the nozzle’s divergent, downstream of the injection slot.Figure 7
Calculated and measured wall pressure in the zone downstream the injection slot.
(a)
Calculated and measured wall pressure  (Case 1)
(b)
Calculated and measured wall pressure (Case 2)
(c)
Comparison between calculated Cases 1 and 4The results in Figures7(a) and 7(b) show that the turbulence models used in this work were able to predict with good accuracy the pressure rise in the separation zone, for Test Case 1 (NPR = 25.9) and Test Case 2 (NPR = 22.2). The numerical results also show that the reignition phenomenon significantly affects the flow field in the separation zone. Due to the combustion process, the separation area becomes larger than the corresponding frozen cases and the separation shock moves upstream. This result highlights the chemical reactions influence on the separation shock displacement and therefore on the pressure rise within the recirculation zone. This result can be connected to the boundary layer, which becomes thicker in reactive flow and easier to detach in response to an adverse pressure gradient. This finding is quite consistent with previous works such as in [24]. Note that Test Case 3 is not represented here in terms of comparison between the calculated and measured pressure profiles. Indeed, for this particular test case, large pressure fluctuations in the combustion chamber were recorded during the course of experiments. This led to great separation shock instabilities and large discrepancy in measurements.Additional Test Case 4 was carried out only by numerical means, with anoxidizer/fuelratio (o/f=6), close to that of real rocket engines, functioning with LO2-LH2. The nozzle pressure ratio NPR is fixed the same way as for Test Case 1; that is, NPR = 25.9, while the total pressure of the H2 film cooling is set to P0j=4.3 bars. As can be seen in Figure 7(c), augmenting the film cooling total pressure leads the separation shock to move downstream as well as the reignition attachment point consequently. This possibility can be used, as will be mentioned below, as means to avoid the reignition of the mixture inside the nozzle.


## 5.4. Chemical Production

To highlight the reignition phenomenon in the separation zone, we plotted the flow OH mass fraction. The results obtained show that OH radical concentrations are at levels high enough to readily suggest a reignition of the H2-Air mixture, in the separated region. The examples illustrated in Figures 8(a) and 8(b) suggest that the entrained fresh oxygen, from the outside into the separation zone, leads to the inflammation GH2 used in the film cooling across the mixing shear layer. A like-diffusion flame appears and extends downstream along the shear layer. A schematic of the flow configuration, with an attached flame to the nozzle upper wall in the separation area, is depicted in Figure 9. It is worth noting that these results were observed for both turbulence models, with modified Evans-Schexnayder chemical kinetics model but not with Eklund model. This can be explained by the temperature levels in the nozzle main flow, across the shear layer, which are not high enough to reactivate the chemical reactions with the above mentioned reaction model. On the other hand, when the nozzle operates at pressure ratios close to the adaptation, as in Case 3 (NPR = 36.5), no reignition process has been observed in the flow inside the nozzle. Indeed, as mentioned above, the separation zone size is too small that it does not allow the full completion of the combustion process within this area. In this case, the combustion takes place in the outer caisson, far away from the nozzle exit section. This situation is depicted in Figure 10 for both OH concentrations and temperature plots.Figure 8
Calculated OH radicals and H2 concentrations: Evans and k-ε calculations.
(a)
Case 1
(b)
Case 4Figure 9
Flow configuration as deduced from CFD calculations.Figure 10
Calculated OH concentration and flow temperature Case 3, Evans andk-ε calculation. The combustion takes place outside the nozzle.
(a)
(b)


## 5.5. Thermal Loads

The film cooling efficiency is evaluated by a widely used parameterη, describing the difference between the film cooled wall temperature, at each point, and the combustion chamber temperature, relative to the difference between the coolant temperature and the combustion chamber temperature. The cooling efficiency is expressed by the nondimensional efficiency η defined as(22)η=Tw-TcTf-Tc,where Tw is the calculated wall temperature, Tc the calculated combustion chamber temperature, and Tf the injected film temperature.The cooling efficiency along the wall depends, as before, on whether the nozzle flow is fully expanded or overexpanded. Figures11(a) and 11(b) represent an example of wall temperature calculated in the Case 1 (NPR = 25.9) and Case 2 (NPR = 22.2), respectively, where the nozzle operates in overexpanded conditions. The combustion of injected GH2 leads to an increase in wall temperature in the region where the flame is attached. The wall temperature plots are characterized by marked peaks up to 1600 K, compared to frozen flow calculations. These results represent a further confirmation of the hypotheses raised previously about the reignition phenomenon. As can be seen, these peaks are shifted according to the observations made above about the separation zone. In this case the cooling efficiency falls to values close to 0.05 and 0.4, respectively, for Case 1 and Case 2 (Figures 12(a) and 12(b)), where the temperature of the attached flame is close to the adiabatic temperature of the combustion chamber, particularly for Case 1. Meanwhile Case 3 shows a typical behavior of film cooling applications as it was observed in several works [2]. In this case the temperature of the film cooling grows gradually by mixing with the hot main flow gas (Figure 11(c)) and the corresponding efficiency (η>0.8) is within the range of the expected values (Figure 12(c)).Figure 11
Wall temperature profiles with reignition, Cases 1 (a) and 2 (b), and without reignition for Case 3 (c).
(a)
(b)
(c)Figure 12
Cooling film efficiency as function of distance from the injection slot for (a) Case 1, (b) Case 2, and (c) Case 3.
(a)
(b)
(c)


## 6. Conclusions

This numerical study was conducted with the goal of reproducing some experimental tests carried out previously. The main objective was to check to which extent the chemical kinetics and turbulence models may be relevant to highlight the reignition phenomenon, in the case of a subscale supersonic nozzle hot flow cooled by H2 film cooling device. The main findings of this study can be summarized as follows.Even though the real separated flow is naturally unsteady, the steady RANS approach, usingk-ε and Spalart-Allmaras models, combined with a modified kinetic model, was found to be able to reproduce the accurate averaged flow field and the reignition phenomenon observed experimentally.Globally, the calculated wall pressures are in good agreement with the experimental data for the whole tested cases and the effect of chemical reactions on these quantities has been pointed out.The computed steady reacting regions are slightly different from the instantaneous PLIF measurements which illustrate the complexity of this phenomenon. Furthermore, it is demonstrated that the steady attached flame occurs along the mixing shear layer and the combustion regime is mostly diffusion.When reignition occurs, the temperature at the nozzle wall is such that the film cooling becomes locally inoperative and the induced thermal loads is so high that the nozzle’s aerodynamic properties can be altered up to being detrimental for the nozzle structural integrity.As expected, the choice of a turbulence model does not globally impact the phenomena that depend on the chemical kinetics, as the reignition process; however, it significantly affects the phenomena related to the boundary layer. This is how Spalart-Allmaras model gives a better prediction of the shock position and the pressure rise in the separation zone, compared to the experimental results.



---

*Source: 252404-2015-04-05.xml*