CFD analysis of the combustion in the BERL burner fueled with a hydrogen-natural gas mixture

The regulatory restrictions, currently acting, impose a significant reduction of the Greenhouse Gas (GHG) emissions. After the coal-to-gas transition of the last decades, the fossil fuel-to-renewables switching is the current perspective. However, the variability of energy production related to Renewable Energy Sources requires the fundamental contribution of thermal power plants in order to guaranty the grid stability. Moving toward a low-carbon society, the industry is looking at a reduction of high carbon content fuels, pointing to Natural Gas (NG) and more recently to hydrogen-NG mixtures. In this scenario, a preliminary study of the BERL swirled stabilized burner is carried out in order to understand the impact of blending natural gas with hydrogen on the flame morphology and CO emissions. Preliminary 3D CFD simulations have been run with the purpose to assess the best combination of combustion model (Non Premixed and Partially Premixed Falmelets), turbulence model (Realizable k ɛ and the Reynolds Stress equation model) and chemical kinetic mechanism (GriMech3.0, GriMech 1.2 and Frassoldati). The numerical results of the BERL burner fueled with natural gas have been compared with experimental data in terms of flow patterns, radial temperature profiles, O2, CO and CO2 concentrations. Finally, a 30% hydrogen in natural gas mixture has been considered, keeping fixed the thermal power output of the burner and the global equivalence ratio.


Introduction
In the last decades, natural gas covered a wide range of energy demand and it has proved essential in the coal-to-gas transition. To figure out, in the last years the substitution of coal with natural gas averted 95 Mt/year of CO 2 emissions. But it seems that this is no longer sufficient, indeed, the European Commission in 2014 proposed a cut of Greenhouse Gas (GHG) emissions by 40% with respect to the 1990 levels until 2030 [1].
For instance, according to the global decarbonization of the power generation, the gas turbine industry is currently engaged in the transition to 100% hydrogen combustion until 2030. At the same time, many researchers are occupied in developing reliable technologies for a hydrogen-based energy infrastructure, since hydrogen represents a carbon-free energy vector for energy-storing and combustion, as stated by the International Energy Agency (IEA) [2,3].
Even though Gas-fired power generation increased slowly in the last year, it accounts for 23% of the overall electric power production. Bearing this in mind, gas turbines fed with natural gas (NG) are characterized by a life-cycle greenhouse gas emissions of averagely 500 g eq /kWh of CO 2 , which is half of Coal power plants, but higher than Nuclear (averagely 65 g eq /kWh) and Renewable energy power plants (averagely 40 g eq /kWh) [4,5].
In the field of combustion and power generation, many studies have been carried out taking into account the hydrogen combustion, which has been applied to systems such as Rankine cycle [21], Double Reheat Rankine cycle [20] and cooktop burners [19].
In the transition from fossil fuels to pure hydrogen combustion, one of the main objectives for the industry is to guarantee the necessary fuel flexibility for operating with different NG-H 2 mixtures [2,6]. Even though, hydrogen is a clean burning fuel and has a high-energy content, however it is prone to flashback and combustion instabilities, which could further degrade, for instance, the performance of gas turbines working under lean conditions. Furthermore, many researchers showed that rising the hydrogen concentration increases the flame speed and the flammability limits of the mixture [7]. Up to now, for gas turbine burners a mixture of 70% methane and 30% hydrogen is suggested as a good compromise between the reduction of the global emission of CO 2 and the absence of significant changes to materials, designs, control and protection systems.
With the purpose to improve the understanding of the basics of NG-H 2 blends combustion mechanism for turbulent flow, a 300 kW burner tested at the Sandia Burner Engineering Research Laboratory (BERL) has been investigated by means of numerical simulations. The BERL burner, which is a NG-air swirled burner with a bluff body, has been widely investigated in the past with remarkable experimental results (benchmark) [8,9]. The burner is characterized by an axisymmetric geometry; indeed a large number of quasi-3D simulations have been performed in the last decades in order to reproduce the experimental results [10,15]. It consists of an axial air flow inlet with a swirl component and 24 holes arranged in the circumferential direction, which supply natural gas in the radial direction. In order to preserve the correct fuel injection penetration, herein a 3D periodic domain, corresponding to 3/24 of the entire geometry, has been investigated. Moreover, the choice to use such a periodicity is motivated by the need to make the vortical structures fully developed in the tangential direction. In this study, steady-state RANS simulations have been performed in the ANSYS Fluent environment, aiming at the assessment of the set up for simulating the NG-air combustion in the BERL combustor. Several combustion models (Non Premixed and Partially Premixed) have been used in combination with different chemical kinetic mechanisms (GriMech3.0 [17], GriMech 1.2 and Frassoldati [11] models) and two turbulence models (k − Realizable and the Reynolds Stress models). Finally, a NG-H 2 fuel blend (70% NG and 30% H 2 ) has been simulated. The air and fuel mass flow rates have been tuned with the aim to keep constant the burner thermal output, the equivalence ratio and the swirl number. In order to evaluate the effect of blending the Natural Gas with hydrogen, the CO emissions have been compared.

Case study
The geometry studied in this work reproduces the swirl-stabilized burner tested at the Burner Engineering Research Laboratory (BERL) located at the Sandia National Laboratories. The BERL combustor is an industrial natural gas burner which has been widely investigated in the past to study the scaling laws for industrial burners and how they affect their emission of pollutants [8]. Thus, various thermal loads, between 30 kW to 12 MW, have been investigated under staged and unstaged conditions, both with different thermal extractions, i.e., refractorylined or water cooled combustion chamber walls. This multiple set up allowed the researchers to study the performance, the flame stability and the NO x emission reduction capabilities of the burner. Furthermore, a lot of data about flow, turbulence and thermochemical patterns were collected.
The burner consists of a swirler and a bluff body with the aim to anchor the flame. A swirler upstream of the burner chamber can provide the air flow with a wide range of swirl numbers (S ) inside the air annulus. The fuel, namely natural gas (NG), comes from an axial pipe inside the annulus with 24 holes positioned in the circumferential direction, see Fig.1. For staging conditions, the fuel can be injected also from holes positioned in the bottom surface. The geometry under examination refers to the 300 kW experimental test case. Many laboratories have been constructed for testing different thermal loads and each of this has its dimensions based on the external annulus diameter (D o ), which for this laboratory is equal to 87 mm, see Fig.1. In this work an unstaged operating condition with 15% excess air has been investigated. The air mass flow rate is equal to 425.95 kg/h at 312 K, whereas the fuel mass flow rate is 21.91 kg/h at 308 K.
In the last decade, the BERL burner has been investigated via quasi-3D axisymmetric numerical simulations [10,15]. Herein, the results of 3D numerical simulations will be shown. In Fig.2, the numerical domain, which reproduces the actual scale swirled burner at Sandia Laboratories, is represented. A 1/8 of the entire geometry, namely a geometry with a rotational periodic boundary condition of 45 • , has been employed in order to reduce the computational costs. Notwithstanding the periodicity, the development of turbulence structure due to the interaction between a pair of holes in the tangential direction is guaranteed [14].

Numerical models and methods
The aim of this work is to perform a preliminary investigation of an industrial burner under different fuel compositions (100% NG and a 30% hydrogen in NG mixture). Herein the simulations have been performed with the commercial code ANSYS Fluent R . Before investigating the BERL with NG-H 2 mixture, a numerical campaign has been carried out with the purpose to select the more accurate and computationally cost effective combination of combustion model and turbulence closure. The application of Flamelets has been investigate with both Non-Premixed and Partially Premixed models. In addition, another important point is the interaction between the turbulence and the combustion models. Indeed, this is a challenge when using a RANS approach. Herein, two turbulence closures have been employed: the Realizable k − and the Reynolds Stress equation model (RS).

RANS method
As previously said, a numerical campaign has been performed in order to assess the CFD model. To do this, Reynolds Averaged Navier-Stokes equations have been employed. The continuity (eq.1) and the momentum (eq.2) equations can be written in Cartesian form as [12]: where u is the velocity, p the pressure, ρ the density and the subscripts i, j, k the three Cartesian coordinates. The RANS simulations have been run with a pseudo-transient approach which is a form of implicit under-relaxation for steady-state case. It helps in stabilizing the case and, at the same time, it gives faster convergence. The under-relaxation is controlled through a pseudo time step size [12].

Turbulence closure
Two turbulence models have been considered for the turbulence closure of the RANS simulations. The k − turbulence model and its modifications (RNG and Realizable) are widely used in simulating combustion processes where the wall treatment is not a primary concern. Looking at the Realizable one, it is characterized by a different formulation for the turbulent viscosity and the transport equation for the dissipation rate [12]. Thus, it enhances the performance for the flows involving rotation, separation and recirculation. Another model employed in this work is the Reynolds Stress equation model (RS), which is the most elaborate and computationally expensive model, since it solves the transport equations for the Reynolds stresses together with the dissipation rate. In 3D cases, it has 7 additional transport equations, whereas eddy-viscosity models need only 2 additional transport equations [12,14], being based on the Boussinesq assumption.

Non-Premixed and Partially Premixed Flamelets
The Non-premixed model considers the combustion problem as a mixing problem. It involves the conservation of the scalar quantity known as mixture fraction, f , defined as follows: where Y stands for the generic atomic element mass fraction. The termochemistry si preprocessed and the interaction between turbulence and chemistry is collected in a Probability Density Function (PDF). The chemistry can be modeled as being in chemical equilibrium (Equilibrium model), near chemical equilibrium (Steady Laminar Flamelet model), or significantly far from chemical equilibrium (Unsteady Laminar Flamelet model) [12].
Other than that, the Partially Premixed model has been considered, since the two flow actually encounter each other before the exit of the burner throat. This model is a combination of the Non-Premixed and the Premixed models, since other than using a preprocessed PDF function, it solves the transport equation of the mean reaction progress variable, c, which is employed for premixed flames. Behind the flame front the progress variable is c = 1, whereas, c = 0 indicates that the mixture is composed by unburnt species, within the flame the progress variable is 0 < c < 1. Moreover, differently from the Non-Premixed model, the latter needs a laminar flame speed profile. ANSYS Fluent R uses fitted curves for laminar flame calculation, anyway they are not always trustworthy, thus in this work, the laminar flame speed profiles for the NG and the NG-H 2 mixture have been corrected with the ones calculated using CANTERA tool and presented by Gimeno-Escobedo [19]. In this work, the Zimont model has been applied for calculating the turbulent flame speed. Finally, both models allow the resolution of Unsteady Laminar Flamelet, which can be used for calculating chemical products derived from slow chemistry, such as NO x .

Grid refinement study
Before starting, a grid refinement study has been carried out on three grids by means of Reynolds Averaged Navier-Stokes simulations. An unstructured mesh has been chosen and three different grids have been generated by means of the ANSYS Workbench R , see Table 1. A density size meshing technique has been adopted with the purpose to increase the geometry refinement in the mixing zone, where the two turbulent flows encounter each other. The refinement zone is made of 3 levels, see Fig.3, so as to include the recirculation zones. In these zones there are 3 points along the streamwise direction where the flow properties have been recorded and compared with the experimental data, see Fig.3. For the grid refinement study, the temperature profiles of the three grids have been compared, see Fig.4, along the axial coordinate (at fixed radial distance y=0.4 m) and along the radial distance (at fixed axial distance x=0.027 m). Unfortunately the experimental results are available only along the radial distance. But in spite of this, the three grids show a similar trend. The grid refinement study has been performed with the Steady State Flamelet model, the k − Realizable turbulence model and the GriMech 3.0 kinetic mechanism [17]. The coarse one seems the less accurate as it deviates from the experimental data in the radial direction and it predicts smoother temperature peaks along the axial distance. Eventually, the fine grid has been chosen as, even if the time required is higher than the medium grid, the fine grid could resolve better the mixing process between fuel and air.

Boundary conditions
In this work, the operating conditions, i.e., flow conditions, turbulence characteristics, chemical fuel composition and the external boundary conditions (wall temperatures) have been inherited by the experimental reports [8]. It is worth to highlight that the BERL combustor has been tested at ambient pressure condition. In particular, the air velocity at the inlet has been imposed by importing velocity profiles for the axial and the tangential velocity components derived from the experimental measurements, for details [8]. The air inflow is characterized by a swirl number S =0.56 with turbulence intensity of 17% and an hydraulic diameter of 29 mm. The fuel, whose composition is summed up in Table 2, is characterized by a molecular weight of 16.52 g/mol, with an air-fuel ratio (A/F) equal to 16.58. The fuel is marked out by a turbulence intensity of 5% and an hydraulic diameter equal to 1.8 mm. The burner works under a global equivalence ratio condition of 0.872 (lean) and the Reynolds number (Re) is equal to 3.20 × 10 5 (D o = 87 mm, v = 60 m/s, ρ mix = 1.10 kg/m 3 and µ mix = 1.82 × 10 5 kg/(m s)) .

Reduced mechanism
In the last decades, several kinetic mechanisms have been developed for the methane combustion, which is characterized by peculiar properties, such as high ignition temperature and low flame speed. A renowned one is the GriMech 3.0, which consists of detailed chemistry, i.e., 325 reactions and 53 species [13,17]. Less expensive kinetic mechanisms have been developed, for instance some of them, already available in the literature, are: the GriMech 1.2 (22 species and 104 reactions) and the one proposed by Frassoldati (10 species and 6 reactions) [11]. Both of them have been used in this work. It is worth to highlight that the GriMech 3.0 includes all the reactions regarding the hydrogen combustion since the hydrogen-oxygen system is actually a subsystem for the oxidation of hydrocarbons [13]. Indeed, it is used as a benchmark for reduced mechanisms developed accounting for methane-hydrogen mixtures [19].

Steady Laminar Flamlet
Firstly, a campaign of RANS simulations has been run in order to investigate the Non-Premixed Steady Laminar Flamlet model, which is the less time consuming model available in ANSYS Fluent R resolving Flamelets (par.3). Flamelets have been simulated with two detailed mechanisms, i.e., GriMech3.0, GriMech1.2 and a skeletal one proposed by Frassoldati [11,17]. For the sake of brevity, the numerical results shown in this paper have been compared with the experiments available in the literature [8] only at 27 mm from the burner outlet along the radial direction. It can be seen that the Non-Premixed Steady Laminar Flamelet model provides generally good results for the temperature profile and the O 2 concentration. On the other hand, it overestimates the axial velocity (u) and it is not able to predict the slow chemistry inherent to the CO. Moreover, the combination of Flamelet and the k − Realizable turbulence model underestimated the tangential velocity (w), see Fig.5. Finally, the Frassoldati kinetic model shows the largest differences in the chemical production calculation. After that, the Partially Premixed model with Steady laminar Flamelets has been set up in order to improve the numerical results. The laminar flame speed calculated by Fluent has been modified according to Gimeno-Escobedo [19]. The results are depicted in Fig.6, where it can be seen that the flow pattern is more accurate as well as the chemical products, which are smoother than the previous ones reducing their discrepancies from the experimental data. Since the tangential velocity values were almost far from the experimental data, the latter set up model has been used in combination with the Reynolds Stress equation (RS), see Fig. 7. For this case, the results with the GriMech 3.0 are exclusively reported. The results point out that the Reynolds Stress model provides good values for the tangential velocity resolving the two characteristic peaks. Furthermore, one remarkable result is that this set up allows the calculation of the mole fraction of CO with the lowest errors.

Hydrogen -Natural gas mixture
Once the suitable combination of combustion and turbulence models was found (i.e., Partially Premixed Combustion Model, with the GriMech 3.0 reaction mechanism, and the Reynolds

Stress Equation turbulence
Model), a simulation of the BERL burner fed with NG enriched with hydrogen has been performed. It is well-known that the hydrogen is characterized by physical properties completely different from the NG. This means that enriching NG with H 2 affects the density, the viscosity, the thermal conductivity, the thermal capacity and the lower heating value of the mixture. Another point is the choice of the operating condition for this new simulation. In order to fix the thermal load of the burner (300 kW), in this study, the air and fuel mass flow rates have been conveniently modified under the assumption that the number and the size of the fuel suppliers do not change. In order to guarantee the same global equivalence ratio (φ=0.872), the fuel inlet velocity has been increased due to the lower density of the mixture (ρ 30%H2 =0.456 kg/m 3 vs ρ NG =0.657 kg/m 3 ), but at the same time, the air mass flow rate has been reduced keeping the same swirl number (S =0.56) at the inlet.
The following results show the comparison of the flow patters, temperature and O 2 , CO 2 and CO mole fractions close to the burner throat (Fig.8). Moreover, the contours of the static temperature and the CO mole fraction for a fuel composition of 70% of NG + 30% of H 2 and 100% NG, (see Fig.9) have been compared with the aim to investigate the flame morphology and the CO emissions. Looking at Fig.9, the hot flame zone for the 70% of NG + 30% of H 2 case seems to be narrower than the baseline, with the core of the flame which moves far from the burner outlet. The flame zone has been highlighted by dashed lines in Fig.9.
Moreover, a lower quantity of CO has been calculated as expected (Fig.9), since the new mixture lacks of chemical compounds containing carbon (C). Eventually, the CO mass fraction has been calculated by means of the volume mass-weighted integral and it shows a reduction of 56% when the NG is enriched with 30% of H 2 . The remaining fraction of CO seems to be controlled by the low temperature, which surrounds the hot zone in the 30% H 2 case.

Conclusions
In this work, a preliminary study of the Natural Gas (NG) swirled stabilized BERL burner has been carried out with the purpose to select the suitable set up for investigating the NG-H 2 mixture in the ANSYS Fluent R environmental. To do this, a RANS simulation campaign has  been performed by varying combustion models (Non-Premixed and Partially Premixed), turbulence models (k − Realizable and the Reynolds Stress equation model (RS)) and chemical kinetic mechanisms (GriMech 3.0, GriMech 1.2 and Frassoldati model). The results show that the Partially Premixed model (with the introduction of validated laminar flame speed) solves the axial velocity component (u) and the chemical compounds better than the Non-Premixed one. Moreover, when the latter is coupled with the RS turbulence model the results are further improved in terms of tangential velocity (w) and CO mole fraction. Given the good results, a 70% NG-30% H 2 mixture has been simulated by modifying the air (at constant S ) and the fuel mass flow rates to fix the global equivalence ratio (φ) and the burner thermal output. The temperature contours show that the hot zone shrinks and moves far from the recirculation zone when NG is enriched with H 2 . Furthermore, the CO production is remarkably reduced (-57%), which is more than the reduction of carbon in the fuel composition. In the future, numerical simulations will be run to account for the NO and NO x production.