# Stochastic characterisation of unstable lean hydrogen-air annular premixed flames

Luís Fernando Figueira da Silva

► To cite this version:

Luís Fernando Figueira da Silva. Stochastic characterisation of unstable lean hydrogen-air annular premixed flames. Combustion Theory and Modelling, 2023, 28 (3), pp.1-27. 10.1080/13647830.2023.2284910. hal-03728087v4

HAL Id: hal-03728087

https://hal.science/hal-03728087v4

Submitted on 31 Oct 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
# Stochastic characterization of unstable lean hydrogen-air annular premixed flames

L. F. Figueira da Silvaª

ªInstitut Pprime, CNRS, ISAE-ENSMA, Université de Poitiers
1 ave. Clément Ader, 86961 Futuroscope Chasseneuil France

## ARTICLE HISTORY

Compiled October 31, 2023

## ABSTRACT

This study is devoted to the characterization of lean premixed hydrogen-air laminar flames stabilized on an annular burner. The flames studied are chosen to lie close to the lean flammability limit, spanning a range of Damköhler and Reynolds numbers, and correspond to situations where instabilities arise. The experimental characterization is achieved by examining the statistical moments of the flame hydroxyl radical (OH*) chemiluminescence and schlieren. The standard deviation and the skewness, that are needed to fully characterize the measured probability density functions, are shown to exhibit a non monotonic behavior with the distance from the burner surface, first increasing, then decreasing downstream the flame tip. An exponentially modified Gaussian distribution model is proposed to describe the skewed schlieren fluctuations probability density function. This model is closed by developing relations for the probability density function parameters based on the Damköhler and Reynolds numbers, and then used to describe the OH* stochastic behavior. The relation between the distribution parameters and the transport equation of a reactive scalar gradient is addressed also.

## KEYWORDS

Skewed probability density function; lean hydrogen flames; premixed combustion; annular burner; unstable flames

## 1. Introduction

Since it can be synthesized entirely from renewable sources, the hydrogen energy carrier has the advantage of not producing greenhouse carbon compounds or soot particles during its combustion with air. In comparison with near-stoichiometric mixtures, the use of fuel-lean mixtures in systems destined for power generation is interesting because the associated flame temperatures impose lesser thermal stresses on combustion chambers walls, of the lower fuel consumption, and of the smaller nitrogen oxides emissions. However, the hydrogen-air combustion process can be intrinsically unstable in fuel-lean situations, which should be anticipated to ensure a safe and reliable operation of new practical systems which could be designed to operate using ultra-lean mixtures. Therefore, this work investigates the behavior of unstable hydrogen-air flames, close to the lean flammability limit, stabilized on a newly developed annular burner.

The stable behavior of premixed laminar flames of hydrogen and air is well known

CONTACT L. F. Figueira da Silva. Email: luis-fernando.figueira-da-silva@ensma.fr
and has been reported in numerous studies, as can be seen in a recent review of the literature [1]. This stable behavior is observed when the Lewis number controlling the combustion is sufficiently close to one, which is the case for near-stoichiometric mixtures, or when the intrinsically unstable modes are unable to develop. When the concentration of hydrogen in the mixture is much lower than the stoichiometry, thermo-diffusive instabilities may arise, which characteristics have been shown to be controlled by the Lewis number of the limiting reactant [2–4]. These thermo-diffusive instabilities lead to flame surface wrinkling that eventually evolves to a cellular structure in unconstrained flames.

The presence of such instabilities in lean hydrogen mixtures has been observed in several flame configurations, notably in spherical flames propagation studies [5–8], and in flames stabilized in planar burners [9–12]. In particular, the first of these configurations has enabled to identify the effects of the mixture equivalence ratio, temperature and pressure on the growth of such instabilities. More specifically, the cellular structure that results from the unstable behavior, and manifests as flame surface wrinkling, was found to arise earlier during the spherical flame expansion as the equivalence ratio decreases [8]. However, the inherent transient nature of the spherical flames makes it difficult to obtain large samples of the observed phenomena, and thus to perform stochastic analysis. On the contrary, on studies performed on planar burner configurations the long time available for the observations permitted to characterize the dynamic response of the flame front instabilities. This characterization evidenced, in particular, the effect of the mixture hydrogen mass fraction on the flame front inherent wrinkling and on the corresponding dynamic modes energy content [12]. The heat transfer to the burner surface in this planar configuration seems to preclude flame stabilization in conditions close to the lean flammability limit, though.

Considering the two-dimensional, freely propagating flame configuration, recent numerical studies have addressed both the early, linear, unstable flame development, and also the long term intrinsic instability behavior [13–15]. More precisely, a stochastic analysis of the flame consumption speed and characteristic thickness evidenced the corresponding large variations along the unstable flame front [13]. This canonical configuration has enabled to identify and characterize the origin of the flame instabilities. However, even if these studies provide valuable information on the thermo-diffusive instability influence on the flame front, these computed quantities are not readily measurable. Furthermore, practical combustors may offer limits on the flame surface growth, due to a finite premixed flame height, for instance. As a consequence, the instabilities development could be constrained along the flame.

In order to describe the unstable combustion phenomena of lean premixed flames, models of high physico-chemical fidelity are often deployed, in particular those that can be solved by computational fluid dynamics techniques of the balance equations using detailed chemical kinetics and molecular transport descriptions [16–18]. However, these models are still associated with computation times and costs that are too high to be directly usable in real-time observation and control of industrial processes. To predict and control these unstable behaviors in real time, it is necessary to develop models capable of anticipating them from sensor data that may be installed in the practical systems [19, 20]. More specifically, it is interesting to combine the detailed results from this type of high-fidelity model with those from measurements that can be deployed in an industrial combustion process environment, thus eventually enabling the construction of reduced-order models [21–25]. In particular, stochastic models, based on transported [26, 27] or prescribed [28, 29] probability density functions of the combustion process properties, are often used to describe hydrocarbon-air com-
bustion in situations where inherent flame instabilities are not relevant. However, the
application of these techniques to inherently unstable premixed hydrogen-air combus-
tion could lead to modifications of these models. This is why the knowledge of the
different stochastic moments of the flame properties probability density functions is
required.

The objective of this study is to study unstable lean premixed hydrogen-air flames stabilized in an annular burner, characterized by a large mixture injection perimeter to slit ratio. This configuration has been chosen to enable the stabilization of flames in situations where circumferential unstable modes develop, but longitudinal ones are constrained by the flame height. More specifically, the aim is to propose a stochastic model of the measured flow-field properties in situations where the flow and the chemical time scales – and thus the Reynolds and Damköhler numbers – exhibit large variations, which is the case when the mixtures are close to the lean flammability limit. Furthermore, the goal is to examine the stochastic properties of schlieren and hydroxyl radical (OH*) chemiluminescence, thus enabling to model the observed skewed probability density functions, and to propose closure relations for the corresponding moments. Note that OH* chemiluminescence could be used in practical devices to infer hydrogen-air flame properties, in particular those associated to the unstable behavior.

On what follows, the experimental approach is first described. Then, the original
experimental results obtained for a range of Damköhler and Reynolds numbers are
discussed. The proposed stochastic model for the skewed probability density functions
of the schlieren and OH* chemiluminescence is then developed, together with closure
assumptions that describe the measured data.

## 2. Methodology

In this section is first described the developed burner, the associated experimental
methods and stochastic characterization, followed by the relevant non-dimensional
numbers definition.

### 2.1. *Experimental methodology*

Figure 1 shows a scheme of the the annular slit burner “ANÔ” developed for this
study. Note that this configuration bears some resemblance to one recently used to
characterize high frequency combustion instabilities [30]. These instabilities are not
found to arise in the present study, though. The slit outer and inner diameters are
$D_o = 57.99$ mm and $D_i = 53.45$ mm, which yields a slit width of $h = (D_o - D_i)/2 =
2.27$ mm, and the slit length is $L = 30$ mm. The measurement uncertainty on these
diameters is 0.2 %. The central part of the burner is hollow so as to enable the free flow
of ambient air. These dimensions have been chosen so as to enable the stabilization of
lean premixed flames close to the lean flammability limit, which is found to occur for a
hydrogen-air mixture equivalence ratio close to $\phi = 0.26$. Furthermore, the perimeter
to slit width ratio, $\approx 75$, should enable unstable circumferential modes to evolve, but
the development of unstable longitudinal modes is expected to be constrained by the
smaller flame height to slit width ratio, which will be seen to lie between 1.5 and 4. The
burner is made of stainless steel, and its uncooled surface temperature is monitored
with thermocouples at the external and internal surfaces. During the experiments the
maximum measured temperature is always smaller than 350 K.

Premixed hydrogen-air enters the slit radially. Fuel-air mixture is effected in a 1.5 m
<figure><img src="image_2.png" /><figcaption>**Figure 1.** Annular slit burner “ANÔ” scheme, with examples of mean OH\* chemiluminescence (letf) and schlieren (right) images obtained for Reynolds and Damköhler number values of 1000 and 16.10⁻³. The dashed rectangle (40×27.5 mm²) represents the area used for averaging the instantaneous schlieren results, whereas the dash-dot line is the burner symmetry axis. Slit outer and inner diameters: D₀ = 57.99 mm and Dᵢ = 53.45 mm.</figcaption></figure>

length tube upstream to the burner. Hydrogen (99.999 % purity) and air (compressed, oil and moisture free) flow rates are measured with mass flow controllers (Brooks models SLA5850 and SLA5851), which have a calibration uncertainty of 1%. As a consequence, the relative standard uncertainty associated to the hydrogen-air mixture equivalence ratio is 1.4%, whereas the Reynolds number uncertainty is 2%. The flow-meter ranges, i.e., 8, 20, 50 and 100 nlmpm are chosen to cover Reynolds and Damköhler numbers values that lie between 100 and 2000, and 2.10⁻³ and 100.10⁻³, respectively. These non-dimensional numbers are defined in section 2.2.

Flame imaging is effected combining schlieren and OH\* chemiluminescence techniques. The latter uses a Princeton Instruments PI-MAX4:1024 EMB low repetition rate (10 Hz) electron multiplying back-illuminated CCD (emICCD, Gen II Intensifier, model RB-FG with P43 phosphor) camera equipped with an objective (CERCO 2085 94 F/4.1, 75% transmissivity in the UV range) and a 10 nm band-pass filter with 65% FWMH transmissivity centered at 311 nm. Note that the OH (A²Σ⁺ − X²Π) transition from the electronic excited state is noted OH\*, and that quantification of the measured results, albeit possible, has not been attempted [31]. Since obtaining OH\* images that encompass the whole annular flame is of interest, the measurement spatial resolution is 13.5 px/mm. This leads to the slit width being resolved along the direction parallel to the ICCD with 31 px. Due to the low OH\* luminosity at the lean equivalence ratios of interest, the camera is operated at a gain of 6,500, and a gate of 300 µs is used. The raw, 1024 × 1024 px², images are filtered using a symmetric 5 × 5 two-dimensional median filtering. Image statistics are obtained from a set of 500 images, which represent an averaging time of 0.5 s. The resulting minimum signal to noise ratio is 10, which corresponds to the flame obtained with the leanest mixture studied. An example of a mean OH\* chemiluminescence field is shown in the left side of Fig. 1. In this figure both the Bunsen-like flame structure above the slit exit, and part of the annular OH\* distribution can be seen.
A dual-lens schlieren arrangement is used to obtain time-resolved flame images. After traversing the flame, collimated light from a diode source is focused by a 100 mm diameter, 750 mm focal length, plano-convex lens to a 45° knife edge. The deviated light is then focused by 95 mm diameter, 150 mm focal distance, bi-convex lens to the Photron FASTCAM Mini AX200 CMOS detector. The schlieren images spatial resolution is 10 px/mm, as a consequence the burner annular slit is resolved with 23 px along the direction parallel to the camera sensor. Throughout the experiments this camera is operated at 6,400 frames per second, with a shutter speed of 4 µs. For each flow and mixture condition a set of 1000 schlieren images is obtained. The same median filter used to process the OH* images is applied to the 1024×1024 px² schlieren images. A typical mean schlieren image can be seen in the right side of Fig. 1. The interfaces between the fresh unburned mixture and the combustion products, which lies above the slit, and between those and ambient air can be seen.

The set of instantaneous images is then used to obtain the statistical moments of order 1 to 3 of the experimental images, which should enable to infer the flame surface statistical behavior. These moments are defined as:

$$ \mu_1 = \mu \equiv E[Z] = \int_{-\infty}^{+\infty} z f(z) dz, \qquad (1) $$

$$ \mu_2 = \sigma \equiv \sqrt{E[(Z - \mu)^2]}, \qquad (2) $$

$$ \mu_3 \equiv E \left[ \left( \frac{Z - \mu}{\sigma} \right)^3 \right], \qquad (3) $$

where $\mu$ and $\sigma$ are the expected value and the standard deviation of the random variable $Z$, which has a density $f(z)$, and $\mu_3$ is the skewness. The skewness, that quantifies the departure of a density from a symmetric distribution, will be shown to be needed to fully describe the measured probability density functions. Note that this moment is scaled with respect to the standard deviation, and that the gaussian distribution skewness is zero. In the present case where the variable $Z$ represents the $N$ sampled images, $z_i$, Eq. (1) is rewritten as:

$$ \mu = \frac{1}{N} \sum_{i=1}^{N} z_i. \qquad (4) $$

The different image processing and statistics algorithms used here are those implemented in Matlab.

Aiming to account for differential illumination and for the average flame structure, each schlieren image is background-corrected by dividing the local intensity by the mean. Thus, if $S(x, y, t)$ is the instantaneous schlieren signal, and $\mu_S(x, y)$ its mean, the stochastic behavior of $s(x, y, t)$, defined as

$$ S(x, y, t) = (s(x, y, t) + 1)\mu_S(x, y), \qquad (5) $$

is of interest. Note that $s(x, y, t)$ is a zero-mean quantity at each $(x, y)$ position, and that the coordinate axis origin, given in Fig. 1, lie at the burner surface along the
symmetry axis.

## 2.2. Non-dimensional numbers

The studied flames are characterized using the Reynolds, Damköhler, Zel'dovich and Lewis numbers.

The Reynolds number definition used is based on the mixture bulk velocity, $\bar{U}$, and burner slit hydraulic diameter, $D_h$, and may be written as:

$$Re = \frac{\bar{U}D_h}{\nu} = \frac{4\dot{\nabla}_t}{\nu\pi(D_o + D_i)}, \quad (6)$$

where $\dot{\nabla}_t$ is the total flow rate and $\nu$ is the mixture kinematic viscosity.

The Damköhler number is defined as the ratio of a characteristic fluid time scale upstream to the annular flame, $\tau_f$, and the unstretched freely propagating planar flame time scale, $\tau_c$. The latter is obtained via Cantera [32] freely propagating one-dimensional adiabatic flame computations, being determined as the ratio between the flame thickness based on the maximum temperature gradient and the unstretched one-dimensional flame speed ($\delta_{L,max}/S_L^0$). These flame computations used the Alekseev chemical kinetics mechanism [33] supplemented by a model for the OH* chemistry [31]. Note that the computed flame speed is known to under-predict the existing experimental data for the lean mixtures of interest [33], which is also the case of other available chemical kinetics mechanisms [1]. This under-prediction could be attributed to the inherent inability of the one-dimensional model to reproduce the measured flames cellular structure, which would increase the burning rate [34]. Despite these shortcomings, these planar flame scales are used for the purpose of determining non-dimensional numbers. The following fits to the computed data are used:

$$S_L^0 = -1107 \phi^3 + 2016 \phi^2 - 788 \phi + 89.49 \ [cm/s], \quad (7)$$

$$\tau_c = 4.863 \cdot 10^5 \exp(-54.83 \phi) + .4165 \exp(-13.61 \phi) \ [s], \quad (8)$$

for which $R^2 > .999$ in the $0.28 < \phi < 0.6$ range. Within this range the flame speed increases from 2.6 to 103 cm/s, whereas the flame thickness decreases from 2.4 to 0.2 mm as the chemical time scale decreases from 104 to 0.19 ms. Note that this flame model is used with the aim of providing reference time and length scales only.

The fluid time scale is taken as being the inverse of the flow velocity gradient upstream to the flame that, under the hypothesis of fully developed channel flow at the burner annulus (with $h \ll D_i$), may be written as $\tau_f = \pi h^2 (D_o + D_i)/12\dot{\nabla}_t$. As a consequence, the Damköhler number reads as:

$$Da = \frac{\tau_f}{\tau_c} = \frac{\pi h^2 (D_o + D_i)}{12 \dot{\nabla}_t \tau_c} = \frac{\pi (D_o - D_i)(D_o^2 - D_i^2)}{48 \dot{\nabla}_t \tau_c}. \quad (9)$$

In order to compare the one-dimensional unstretched laminar premixed flame thickness, $\delta_L^0$, with the characteristic length scale of the flow velocity gradient, $\delta_f$, a non
dimensional thickness is used:

$$ \delta^* = \frac{\delta_L^0}{\delta_f} = \frac{S_L^0 \tau_c}{h/2}. \qquad (10) $$

Note that these non dimensional quantities are related via:

$$ DaRe \left(\frac{\delta^*}{2}\right)^2 = \frac{1}{3} \frac{S_L^0 \tau_c}{\nu/S_L^0} = \frac{1}{3} \frac{\delta_{L,max}}{\delta_{L,T} Pr}, \qquad (11) $$

where Pr is the Prandt number, the flame thickness based on the maximum temperature gradient and the thermal flame thickness are $\delta_{L,max}$ and $\delta_{L,T}$, respectively. Classically, the flame thermal thickness is $\delta_{L,T} = \alpha/S_L^0$, where $\alpha$ is the mixture thermal diffusion coefficient. If the approximation $\delta_{L,max}/\delta_{L,T} = 2(T_b/T_b)^{0.7}$ is used [35], Eq. (11) may be written as:

$$ DaRe \left(\frac{\delta^*}{2}\right)^2 Pr = \frac{2}{3} \left(\frac{T_b}{T_u}\right)^{0.7}, \qquad (12) $$

where the temperatures of the fresh and the burnt gases are $T_u$ and $T_b$, respectively. Given this relationship, and since Da and $\delta^*$ exhibit opposite trends with the hydrogen-air mixture equivalence ratio, the choice is made here to characterize the studied flames using Da and Re. Note that if the thermal and maximum temperature gradient flame thickness are supposed to be identical, the right hand side of Eq. (12) reduces to 1/3.

The hydrogen-air mixture effective Lewis number is computed via a formulation [36] that enables a continuous variation between fuel lean, and fuel rich limits and uses the Zel'dovich to that end:

$$ Le = Le_1 + H(Le_2 - Le_1), \qquad (13) $$

where $Le_2 = 2.33$ and $Le_1 = 0.33$ are the fuel-rich and oxidizer-rich values for hydrogen-air mixtures, and $H$ is given by:

$$ H = n \frac{G(m, n - 1, A)}{G(m, n, A)}, \\ G(m, n, A) = \int_{0}^{\infty} \xi^{m} (\xi + A)^{n} \exp(-\xi) d\xi, \qquad (14) $$

$m$ and $n$ are the reaction orders with respect to the deficient and abundant reactants, respectively, here it is assumed $m = n = 1$, for simplicity. The normalized final concentration of the abundant reactant is either $A = \beta(\phi - 1)/Le_2$ or $A = \beta(\phi - 1)/Le_1$, and $\beta = E_a(T_b - T_u)/(RT_b^2)$ is the Zel'dovich number.

In order to determine the Zel'dovich number, $\beta$, a recently proposed methodology is employed [37], that permits to continuously compute $\beta$ throughout the whole equivalence ratio range. This methodology optimally fits an ad-hoc premixed flame heat release rate expression

$$ \dot{q}(\theta) \propto \{1 - \exp[\beta'(1 - 1/\theta)]\}^{\beta''} \exp\left\{\frac{\beta(\theta - 1)}{(1 + \gamma\theta)/(1 + \gamma)}\right\}, \qquad (15) $$
to the values computed of the above-mentioned detailed chemical kinetics calculations. In this expression $\theta = (T - T_u)/(T_b - T_u)$ is the temperature based reaction progress variable, and $\gamma = (T_b - T_u)/T_u$ is the flame thermal expansion, whereas $\beta'$ and $\beta''$ control, respectively, the initial and the final $\dot{q}$ decrease towards equilibrium as $\theta \to 1$. The resulting Zel'dovich number correlation that is used in Eq. (14) reads:

$$ \beta = 5.432 \cdot 10^4 \exp(-30.22 \phi) + 6.166 \exp(-2.657 \phi), \quad (16) $$

for which $R^2 > .999$ in the (0.28, 0.6) equivalence ratio range. The Zel'dovich number is found to decrease from $\approx 10$ to $\approx 1.3$ in this range, whereas the Lewis number increases from 0.38 to 1.0.

Finally, it is worth mentioning that, given the above definitions, the ratio between the mixture bulk velocity at the slit and the laminar flame speed may be written as:

$$ \frac{\bar{U}}{S_L^0} \propto \frac{Re \delta^*}{Pr Le}. \quad (17) $$

### 3. Results and discussion

In this section the studied flow conditions are first specified in terms of the hydrogen and air flow rates and the relevant non dimensional parameters. The overall flame structure, as unveiled by the schlieren and OH* distributions is then presented, followed by an analysis of the corresponding stochastic properties. A model of the probability density function of these measured properties is then developed, together with the required closure hypothesis. Finally, the relation between this new model and the stochastic properties of the gradient of a reactive scalar is assessed.

#### 3.1. *Studied flow conditions*

In this study are examined the results of 12 different flow conditions, given in Tab. 1 together with selected laminar one-dimensional unstretched flame properties and the corresponding non-dimensional numbers defined in section 2.2. These flow and mixture studied have been chosen to span a Reynolds number range of (400, 1000) and a Damköhler number range from $2.10^{-3}$ to $16.10^{-3}$, which correspond to lean fuel/air mixtures close to the lean flammability limit. In this table a given letter corresponds to a Reynolds number, whereas the integers indicate the Damköhler number. The burner is operated in 1 atm and the combustible mixture fed at 294 K.

The values of $\delta^*$ reported in this table indicate that the laminar flame thickness ($\delta_{L,max}$) may be either smaller or larger than the slit half-width. Indeed, for the leanest and richest fuel/air mixtures studied, the computed flame thickness is 2.8 and 0.76 mm, respectively, which correspond to $\delta^*$ of 2.5 and 0.67. Furthermore, the leanest mixture Zel'dovich number is of the order of 14, whereas for the richest mixture $\beta < 4$. Note that this leads to a significant change of the computed chemical heat release rate distribution within the premixed flame, not shown here for the sake of brevity. Indeed, for the smaller values of $\beta$, the heat is released for intermediate temperature values, whereas for the larger values this heat release occurs closer to the maximum flame temperature.

As shown in Tab. 1, the effective Lewis number of the studied cases are all significantly smaller than one. Therefore, thermo-diffusive instability is expected to influence
<table><caption>**Table 1.** Experimental flow conditions, laminar flame properties and non dimensional numbers.</caption><thead><tr><th>Case</th><th>∇<sub>air</sub><br>(nlpm)</th><th>∇<sub>H<sub>2</sub></sub><br>(nlpm)</th><th>φ</th><th>S<sub>L</sub><sup>0</sup><br>(cm/s)</th><th>Re</th><th>Da<br>×10<sup>3</sup></th><th>δ*</th><th>β</th><th>Le</th></tr></thead><tbody><tr><td>A1</td><td>35</td><td>4.15</td><td>.282</td><td>2.79</td><td>405</td><td>2.27</td><td>2.49</td><td>13.7</td><td>.38</td></tr><tr><td>A2</td><td>35</td><td>4.49</td><td>.305</td><td>5.32</td><td>407</td><td>7.00</td><td>1.53</td><td>8.08</td><td>.43</td></tr><tr><td>A3</td><td>35</td><td>4.74</td><td>.322</td><td>7.88</td><td>409</td><td>14.7</td><td>1.07</td><td>5.82</td><td>.47</td></tr><tr><td>B1</td><td>52</td><td>6.35</td><td>.291</td><td>3.58</td><td>601</td><td>2.33</td><td>2.09</td><td>11.2</td><td>.40</td></tr><tr><td>B2</td><td>52</td><td>6.87</td><td>.314</td><td>6.62</td><td>606</td><td>7.09</td><td>1.26</td><td>6.73</td><td>.45</td></tr><tr><td>B3</td><td>52</td><td>7.28</td><td>.333</td><td>9.80</td><td>610</td><td>15.0</td><td>0.87</td><td>4.85</td><td>.50</td></tr><tr><td>C1</td><td>70</td><td>8.65</td><td>.294</td><td>3.95</td><td>810</td><td>2.05</td><td>1.95</td><td>10.3</td><td>.40</td></tr><tr><td>C2</td><td>70</td><td>9.46</td><td>.322</td><td>7.76</td><td>818</td><td>7.14</td><td>1.09</td><td>5.89</td><td>.47</td></tr><tr><td>C3</td><td>70</td><td>10.1</td><td>.343</td><td>11.8</td><td>825</td><td>15.8</td><td>0.74</td><td>4.17</td><td>.53</td></tr><tr><td>D1</td><td>85</td><td>10.65</td><td>.298</td><td>4.42</td><td>985</td><td>2.06</td><td>1.78</td><td>9.42</td><td>.41</td></tr><tr><td>D2</td><td>85</td><td>11.7</td><td>.328</td><td>8.78</td><td>995</td><td>7.45</td><td>0.97</td><td>5.31</td><td>.48</td></tr><tr><td>D3</td><td>85</td><td>12.5</td><td>.350</td><td>13.2</td><td>1004</td><td>15.9</td><td>0.67</td><td>3.82</td><td>.55</td></tr></tbody></table>

the flame development above the annular slit.

These studied flow conditions lie closer to the burner blow-off limit than to the flashback one. Blow-off is considered to occur when the flame anchoring ceases at either of the annular surfaces – internal or external, whereas the flashback is annotated once the flame fully enters the annular space, since partial entrances are not visible. Note that the schlieren images are used to characterize these limit situations. The blow-off limit is nearly independent of the Reynolds number, being observed for $\phi \approx 0.27$, and thus $Da \approx 1.5.10^{-3}$. The studied conditions are such that blow-off is not observed, and the flame wrinkling will be shown to be purely combustion driven. For the considered Reynolds number range the critical flashback Damköhler number is of the order of $80.10^{-3}$.

## 3.2. Annular flame structure

Typical instantaneous, average and background-corrected schlieren images are shown in Figs. 2a and 2b. The corresponding flow conditions, given in Tab. 1, were chosen to represent the smallest (A1) and largest (D3) Damköhler and Reynolds number values studied. A first assessment of the annular flame structure may be obtained by examining the instantaneous schlieren images, $S(x, y, t)$, shown at the top images. At the left and right parts of each image the Bunsen-like conical flame surface is clearly seen above the burner surface, as is the interface between the burned gases and ambient air that surrounds this flame. The center part of these images is occupied by collected light intensity fluctuations which could stem from a wrinkled flame surface. In the average images, $\mu_S(x, y)$, given at the center of these figures, these wrinkled structures are absent for case D3, thus confirming their instantaneous nature. Concerning case A1, a reminiscence of the wrinkles is seen, one which is associated to the use of a fixed sampling time for all experiments, and to the slower fluctuating motion of this case, when compared to D3. The flame cone and interface between the burned gases and ambient air can be clearly seen also. In these average images the flame cone tip seems to be open, which is a classical feature of low Lewis number flames, such as those. The bottom images of these figures display the instantaneous background corrected signal, $s(x, y, t)$, defined in Eq. (5). It is possible to verify that the background correction eliminates the conical flame structure above the slit, as well
<figure><img src="image_3.png" /><figcaption>**Figure 2.** Schlieren images: instantaneous $S(x, y, t)$ (top), average $\mu_S(x, y)$ (center) and background-corrected instantaneous $s(x, y, t)$ (bottom) on a 72.5 × 27.5 mm$^2$ field of view above the burner surface.</figcaption></figure>

as the interface between the burned gases and ambient air. Furthermore, these wrinkles seem to be vertically aligned close to the burner surface, displaying a seemingly two-dimensional structure along the flame, but exhibit a more complex structure further downstream, in particular for the low-$Da$, high-$Re$ situations. Downstream the flame cone these wrinkles progressively disappear, which suggests that they are combustion-driven. The quantitative analysis of the flame wrinkling length scale is deferred to section 3.3.

To further characterize the annular flow structure, Fig. 3 shows the average OH* signal, $\mu_{OH}$, corresponding to cases A1 (top) and D3 (bottom), given in Tab. 1, and which schlieren results are depicted in Figs. 2a and 2b, respectively. Note that the color scales used in these images are different, since the maximum value of $\mu_{OH}$ in case D3 is approximately five times larger than that of case A1, which is explained by the higher equivalence ratio of the former mixture. In this figure the flame cone can clearly be seen at the left and right sides, and the cone tip opening that has been hypothesized to exist in the schlieren averaged images analysis is confirmed. Since these are line-of-sight integrated images the OH* intensity is larger at the part of the cone closer to the burner symmetry axis than at the outer parts. Along the center part of these figures the OH* signal is fainter, which is also a feature of this line-of-sight measurement. Due to the low signal to noise ratio of the smaller equivalence ratio cases, no attempts were made to perform image deconvolution, that would have enabled to describe the OH* signals along the burner symmetry plane.

Further comparing the OH* (Fig. 3) and schlieren (Figs. 2a and 2b) distributions, it is clear that the wrinkling is associated to the combustion process, i.e., to the flame surface. Therefore, the analysis of the schlieren and OH* stochastic behavior should provide a description of these flame front properties.
<figure><img src="image_4.png" /><figcaption>**Figure 3.** Average OH* images ($\mu_{OH}(x, y)$): case A1 (top) and case D3 (bottom) on a 64 × 23.5 mm² field of view above the burner surface. The maximum OH* value of case D3 is 5.2× that of case A1.</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>**Figure 4.** Vertical evolution of the schlieren standard deviation ($\sigma_S$, top) and skewness ($\mu_{3,S}$, bottom) for the different studied cases.</figcaption></figure>

### 3.3. Flame stochastic behavior

It may be expected for the statistical properties of the schlieren and OH* images to change with *Da* and *Re*. Accordingly, the standard deviation and the skewness of the background-corrected schlieren images, $\sigma_S$ and $\mu_{3,S}$, given in Fig. 4 as a function of the vertical coordinate, enables to quantify the flame surface fluctuations. In this figure are compared, for a given Damköhler number, the results obtained for each Reynolds number. These results correspond to the averaging, along the direction parallel to the burner surface, *x*, of a given statistical moment computed for the set of 1000 images, i.e.,

$$ \bar{\mu}_{i,S}(y) = \frac{1}{\Delta x} \int_{-\Delta x/2}^{\Delta x/2} \mu_{i,S}(x, y) \, dx, \qquad (18) $$

where *i* = 2 corresponds to the standard deviation, and *i* = 3 to the skewness. Since only $\bar{\mu}_{i,S}(y)$ are analyzed henceforth, the bar is omitted in the remaining of the text. In order to restrict the focus to the central portion of the annular flame, this averaging
<table><caption>**Table 2.** Notable experimental values.</caption><thead><tr><th>Case</th><th>y<sub>σ<sub>S</sub></sub><br>(mm)</th><th>σ<sub>S,max</sub><br>×10<sup>3</sup></th><th>τ<sub>P</sub>/σ<sub>P</sub><sup>2</sup><sub>max</sub><br>×10<sup>3</sup></th><th>σ<sub>P,max</sub><br>×10<sup>3</sup></th><th>τ<sub>P,max</sub><br>×10<sup>3</sup></th></tr></thead><tbody><tr><td>A1</td><td>4.77</td><td>126</td><td>7.96</td><td>99.4</td><td>78.3</td></tr><tr><td>A2</td><td>3.67</td><td>126</td><td>9.03</td><td>95.4</td><td>82.3</td></tr><tr><td>A3</td><td>3.18</td><td>127</td><td>8.51</td><td>97.5</td><td>80.8</td></tr><tr><td>B1</td><td>5.96</td><td>193</td><td>9.00</td><td>127</td><td>145</td></tr><tr><td>B2</td><td>4.77</td><td>183</td><td>10.1</td><td>118</td><td>140</td></tr><tr><td>B3</td><td>3.97</td><td>179</td><td>10.0</td><td>117</td><td>136</td></tr><tr><td>C1</td><td>6.75</td><td>231</td><td>8.87</td><td>143</td><td>181</td></tr><tr><td>C2</td><td>5.16</td><td>219</td><td>10.2</td><td>131</td><td>175</td></tr><tr><td>C3</td><td>4.77</td><td>212</td><td>10.9</td><td>125</td><td>171</td></tr><tr><td>D1</td><td>6.75</td><td>246</td><td>9.32</td><td>146</td><td>198</td></tr><tr><td>D2</td><td>5.66</td><td>241</td><td>10.7</td><td>137</td><td>199</td></tr><tr><td>D3</td><td>4.77</td><td>237</td><td>12.0</td><td>129</td><td>199</td></tr></tbody></table>

is performed along $\Delta x = 40$ mm only (400 px). The sample size is thus 40,000 at a given vertical position above the burner surface.

Figure 4 shows that a similar behavior is observed in all cases, i.e., the standard deviation first increases, then decreases above the burner surface. The non-zero value of the standard deviation at the burner surface ($\sigma_S(y = 0) = \sigma_{S,0} \approx 0.09$) and initial increase along the vertical position seem to be independent both of *Da* and of *Re*. However, the maximum standard deviation value and the corresponding vertical position are found to increase with *Re*, being less sensitive to *Da* changes. In particular, for the smaller *Da* values, when the Reynolds number increases from 400 to 1000, the maximum standard deviation goes from $\sigma_{S,max} \approx 0.13$ to $\approx 0.25$. Table 2 gives these maximum values for the different studied cases. Furthermore, this non monotonic standard deviation behavior with respect to *y* suggests that this maximum position could be used as a measure of the flame height, which will be examined below.

In order to expand the statistical analysis of the flame front fluctuations, Fig. 4 also depicts the corresponding skewness. This statistical moment also exhibits a non-monotonic behavior above the burner surface, first increasing, then decreasing with height. Furthermore, similarly to the standard deviation, the maximum skewness value is found to increase with the Reynolds number, but is nearly independent of *Da*. In particular, the skewness increases from $\approx 0.5$, for *Re* = 400 to $\approx 1.0$, when *Re* = 1000. These non-zero values of the skewness indicate that the probability density function of *s* should be an asymmetrical distribution, which will be examined below.

To further the analysis, the transversely averaged first and second moments of the OH* images, i.e., $\mu_{OH}$ and $\sigma_{OH}$, are given in Fig. 5. The choice of depicting the average OH* evolution is motivated by the classical use of chemiluminescence to characterize combustion heat release, however no attempts to quantify the concentration has been made [31]. In this figure these moments have been arbitrarily divided by the global mean maximum, that is found to occur for case D3, which thus enables to quantitatively compare the cases among themselves. Similarly to what has been observed above for the schlieren moments, both $\mu_{OH}$ and $\sigma_{OH}$ increase with the Reynolds and the Damköhler numbers, the latter seems however more pronounced here than for $\sigma_S$. The overall evolution of $\mu_{OH}$ and $\sigma_{OH}$ with the vertical position, $y$, is quite similar. Indeed, an increase near the burner surface which is steeper than that of $\sigma_S$ is first seen. Further downstream the moments decrease to zero, which is consistent with the ten-
<figure><img src="image_6.png" /><figcaption>Figure 5. Vertical evolution of the OH* mean deviation ($\mu_{OH}$, top) and standard deviation ($\sigma_{OH}$, bottom) for the different studied cases.</figcaption></figure>

dency of reaching chemical equilibrium downstream the flame, but a different behavior is found to occur for intermediate $y$ values. Indeed, concerning $\mu_{OH}$, a slow decrease is seen to occur, the extent of which increases with $Re$ but decreases with $Da$. On the contrary, in this intermediate $y$ range, $\sigma_{OH}$ slowly increases. This slow increase could be a manifestation of the increased schlieren fluctuations, since, by comparing Figs. 4 and 5, the decreasing behavior is seen to initiate immediately downstream the maximum $\sigma_S$ position. Therefore, the $y$ position at which $\sigma_{S,max}$ occurs, $y_{\sigma_S}$, will be used as a measure of flame height. The $y_{\sigma_S}$ values obtained for the different studied cases is also given in Tab. 2. These values may be fitted by

$$y_{\sigma_S} = a_y \frac{Re^{3/8}}{Da^{2/10}}, \qquad (19)$$

with $a_y = 0.155$ mm for a $R^2 = 0.96$ with 95 % confidence bounds.

Another relevant feature of the wrinkled flame front is the characteristic length scale of these wrinkles. The average length scale of the schlieren fluctuating signal seen in Figs. 2a and 2b, $l_S$, is given in Fig. 6 as a function of the non-dimensional vertical position, $y/y_{\sigma_S}$. This length scale is defined as

$$l_S(y) = \frac{\Delta x}{N} \sum_{i=1}^{N} \left( \int_{-\Delta x/2}^{-\Delta x/2} \delta(s_b(x, y, t_i)) dx \right)^{-1}, \qquad (20)$$

where $\delta(x)$ is the Dirac delta function, and $s_b(x, y, t_i)$ is the binarized schlieren signal. This enables to determine the number of positive and negative $s(x, y, t_i)$ regions for each $y$. The sum is performed throughout the $N$ images, and $\Delta x$ is the width of the considered region. More precisely, this length scale has been determined by the following procedure. First each instantaneous image is binarized, thus enabling to identify regions with opposite $s$ signals to finite regions of zeros or ones. Then, at each vertical position, $y$, a given binarized image is treated as a square wave, which
<figure><img src="image_7.png" /><figcaption>**Figure 6.** Mean transversal length scale of the schlieren fluctuations $l_s$ as a function of the height above the burner surface, $y/y_{\sigma_S}$. Damköhler numbers: $2.10^{-3}$ (top), $7.10^{-3}$ (center), $15.10^{-3}$ (bottom).</figcaption></figure>

pulse widths are determined, then averaged at a given height for all images, thus leading to $l_s$. Given this procedure, Fig. 6 shows that the value of $l_s$ remains nearly constant for $0 < y/y_{\sigma_S} < 1.5$, and increases further downstream. The values of $l_s$ within this plateau region are independent of the Reynolds number, and decrease with the Damköhler number. Indeed, as $Da$ increases from $\approx 2.10^{-3}$ to $\approx 15.10^{-3}$, $l_s$ goes from 1.3 to 1.0 mm. This suggests that the processes that control $l_s$ remain invariant in that $y/y_{\sigma_S}$ range, which should be the case in regions of the flow where $s$ is controlled by the premixed flame, i.e., by the balance between chemical and fluid time scales. Furthermore, an increase of $l_s$ associated to viscous damping of the fluctuations was to be expected downstream the flame, i.e. for $y/y_{\sigma_S} > 1.5$ since, as shown above, $\sigma_S$ decreases in this region. These results show that the initial flame wrinkling is purely combustion driven, and that downstream the flame tip viscous dissipation acts to dampen these wrinkles. This is indeed the case of the flame structure results depicted in Fig. 2, where both the instantaneous, $S(x, y, t)$, and the fluctuating, $s(x, y, t)$ schlieren exhibit quite regular structures along the transversal direction close to the burner surface ($y_{\sigma_S} < 1$). Downstream to the average flame tip ($y_{\sigma_S} > 1$), these structures are progressively dampened.

To further characterize the fluctuating behavior of the schlieren images, Figs. 7a and 7b present the histograms (in blue) of the probability density function of $s$, $P(s)$, at different heights, for cases A1 and D3, respectively. The modeled PDF, given by the red line, will be discussed in section 3.4. These heights have been chosen to correspond to $y/y_{\sigma_S}$ values of (0.5, 1.0, 1.5 and 2.0) with the purpose of illustrating the $P(s)$ evolution around the maximum $\sigma_S$ region. These similar normalized heights correspond to different $y$ values, which are indicated above each PDF. Following what has been done above for the stochastic moments, this PDF has been constructed by accumulating the $s$ values in a transverse region spanning $-20 \le x \le 20$ mm. Aiming to ease the comparison between each case, that correspond to rather different standard deviation values, the horizontal axis has been normalized by the corresponding maximum value, $\sigma_{S,max}$, which are given in Tab. 2. Both these figures show that $P(s)$ spans to smaller $s/\sigma_{S,max}$ away from $y/y_{\sigma_S} = 1.0$, i.e., as the value of $\sigma_S$ decreases
<figure><img src="image_8.png" /><figcaption>**Figure 7.** Probability density function of the schlieren fluctuations $s(x, y, t)/\sigma_{S,max}$ at different heights above the burner surface, $y/y_{\sigma_S} = (0.5, 1, 1.5, 2)$, the corresponding y values are given above each image. Blue histogram: experiment, red line: model.</figcaption></figure>

upstream or downstream this maximum position. However, these cases correspond to rather different $P(s)$ shapes. Indeed, the case A1 exhibits a nearly symmetrical PDF with respect to $s = 0$, whereas the PDF of case D3 is significantly skewed. Such behavior is connected to the different skewness values seen in Fig. 4, i.e., the maximum $\mu_{3,S}$ values are $\approx 0.45$ and $\approx 1.2$, respectively.

A model for this probability density function, $P(s)$, will now be developed. However, the discussion of the OH* probability density function, which also exhibits a significantly skewed shape, is deferred to after the presentation of different $P(s)$ moments models.

### 3.4. Proposed stochastic model

A stochastic model of the schlieren fluctuating signal $s(x, y, t)$ is now developed. This development proceeds by considering first the corresponding probability density function, $P(s)$, then a relation between the parameters of this PDF and finally, an evolution equation for these parameters. Relations between the obtained expressions and those that ensue from the transport of a reactive scalar function are discussed also.

#### 3.4.1. Prescribed probability density function

Given the nature of the observed stochastic process of $s(x, y, t)$, which involves a decay downstream to a region where forcing, $F(s)$, occurs, it is proposed that at each height, $y' = y/y_{\sigma_S}$, the PDF of $s$, $P(s)$, is described by:

$$ \tau_P \frac{dP}{ds} + P = F(s). \qquad (21) $$

In the absence of forcing, i.e., $F(s) \equiv 0$, it is evident that $P(s)$ is an exponential PDF with decay rate $\tau_P$. When forcing is present and $\tau_P \to 0$, one obtains $P(s) \to F(s)$. Here it is assumed that the forcing is a Gaussian PDF with zero mean and standard
deviation $\sigma_P$, i.e.,

$$F(s) = N(0, \sigma_P) = \exp \left( -\frac{1}{2} \left( \frac{s}{\sigma_P} \right)^2 \right). \qquad (22)$$

This choice leads to a closed-form solution of Eq. (21) as an exponentially modified Gaussian function:

$$P(s) = \frac{1}{2\tau_P} \exp \left( \frac{1}{2} \left( \frac{\sigma_P}{\tau_P} \right)^2 - \frac{s - \mu_P}{\tau_P} \right) \times \operatorname{erfc} \left( \frac{1}{\sqrt{2}} \left( \frac{\sigma_P}{\tau_P} - \frac{s - \mu_P}{\sigma_P} \right) \right), \qquad (23)$$

where erfc is the complementary error function. This PDF has been used to describe the shape of chromatographic measurements, for instance [38, 39]. The PDF parameters $\mu_P$, $\tau_P$ and $\sigma_P$ are related to the moments of $P(s)$ following

$$\begin{align}
\mu_P &= \mu_1 - \mu_2 \left(\frac{\mu_3}{2}\right)^{1/3}, \nonumber \\
\sigma_P^2 &= \mu_2^2 \left(1 - \left(\frac{\mu_3}{2}\right)^{1/3}\right), \nonumber \\
\tau_P &= \mu_2 \left(\frac{\mu_3}{2}\right)^{1/3}, \tag{24}
\end{align}$$

where the sample mean, standard deviation and skewness ($\mu_1$, $\mu_2$ and $\mu_3$) are given by Eqs. (1)-(3). Accordingly, provided these moments are known, the PDF is completely determined. Note that the definition of $s$ trivially leads to $\mu_1 \equiv 0$, thus two moments are left to be specified.

A first obvious choice is to use the measured values of $\mu_2$ and $\mu_3$ depicted in Fig. 4 to determine the PDF parameters in Eqs. (24). The resulting $P(s)$ is given in Figs. 7a and 7b as the solid red lines. The agreement between the measured PDF and the prescribed exponentially modified function is excellent, and applies also to all the studied cases and for $y$ values different than those depicted. This excellent agreement permits to say that this PDF model fully explains the measured fluctuating schlieren data. However, the closure of the PDF transport does require that the second and third order PDF moments be specified, which is addressed now.

### 3.4.2. PDF parameters closure

In order to provide closed expressions for the exponentially modified PDF, parameters $\sigma_P$ and $\tau_P$, given in Eqs. (24), are assumed to be functions of the Reynolds and Damköhler numbers, defined in Eqs. (6) and (9), respectively.

A first modeling hypothesis is that $\tau_P$ describes $P(s)$ relaxation and, thus, could be supposed to be dependent on the Reynolds number only. Since this relaxation should reflect the balance between inertia and viscous forces, it is then assumed that $\tau_P \propto Re$. A second modeling hypothesis is that the PDF forcing, supposed to follow a normal distribution (Eq. (22)), includes both chemistry and transport effects, such that $\sigma_P = \sigma_P(Re, Da)$ is a function to be determined.
<figure><img src="image_9.png" /><figcaption>**Figure 8.** Ratio between the exponentially modified PDF parameters, $\tau_P/\sigma_P^2$, a function of the height above the burner surface, $y/y_{\sigma_S}$, for the different studied cases.</figcaption></figure>

However, before developing these closure functions, it is interesting to note that $\sigma_P$ and $\tau_P$ exhibit a peculiar relationship. This can indeed be observed in Fig. 8, where the evolution of $\tau_P/\sigma_P^2$, as determined by using $\mu_{2,S}$ and $\mu_{3,S}$ in Eq. (24) is given as a function of the non dimensional vertical coordinate, $y/y_{\sigma_S}$. This figure permits to verify that a nearly constant value of $\tau_P/\sigma_P^2$ arises, for each (Re, Da), in the range $0.25 < y/y_{\sigma_S} < 2.0$, i.e., along the flame region. A rather clear separation arises among the different cases studies, which should enhance the corresponding discrimination in the (Re, Da) plane, also. Furthermore, this behavior suggests that the relation between the parameters values at the particular position where the PDF variance is maximum, i.e., $y/y_{\sigma_S} = 1$, could be used to represent the same relation elsewhere, thus enabling to express $\tau_P$ as a function of $\sigma_P$.

Accordingly, these parameters at the maximum standard deviation position are expressed as functions of (Re, Da),

$$ \tau_{P,max} = f_{\tau}(Re) = a_{\tau}Re, \qquad (25) $$

$$ \frac{\tau_P}{\sigma_P^2} = \frac{\tau_P}{\sigma_P^2}\bigg|_{max} = f_{\tau,\sigma}(Re, Da) = a_{\tau,\sigma}Re^{b_{\tau,\sigma}}Da^{c_{\tau,\sigma}}, \qquad (26) $$

where $\tau_P/\sigma_P^2|_{max}$ is the value at $y/y_{\sigma_S} = 1$, and the parameters $a_{\tau}$, $a_{\tau,\sigma}$, $b_{\tau,\sigma}$ and $c_{\tau,\sigma}$, obtained by fitting the measured data given in Tab. 2 are

$$ a_{\tau} = 0.21.10^{-3}; \quad R^2 = 0.95, \qquad (27) $$

$$ a_{\tau,\sigma} = 3.2, \quad b_{\tau,\sigma} = \frac{1}{4}, \quad c_{\tau,\sigma} = \frac{1}{10}; \quad R^2 = 0.88, \qquad (28) $$

with 95 % confidence bounds.
<figure><img src="image_10.png" /><figcaption>Figure 9. Exponentially modified PDF parameter $\sigma'_P$ a function of the height above the burner surface, $y/y_{\sigma_S}$, for the different studied cases. Also show are the Eq. (31) model results for two values of $n$.</figcaption></figure>

This first closure for the PDF parameters assumes that a relation exists between $\tau_P$ and $\sigma_P$ which is a known function of the Reynolds and Damköhler numbers. As a consequence, three parameter exponentially modified Gaussian PDF is reduced to a single parameter function, $\sigma_P$. Thus, reworking Eqs. (24), one may express this parameter as a function of $\mu_{2,S}$, the measured PDF standard deviation:

$$ \sigma_S^2 = \mu_{2,S}^2 = \sigma_P^2 + (f_{\tau,\sigma} \sigma_P^2)^2. \qquad (29) $$

Note that the results of Fig. 8 show that the order of magnitude of $f_{\tau,\sigma}$ is $\approx O(10)$, which implies that

$$ \sigma_S = f_{\tau,\sigma} \sigma_P^2 \qquad (30) $$

is a very good approximation for most of the $y/y_{\sigma_S}$ range, except when $y/y_{\sigma_S} \to 0$.

A final closure step involves determining an expression for $\sigma_P$. Accordingly, this parameter normalized by the corresponding maximum and minimum values, i.e., $\sigma'_P = (\sigma_P - \sigma_{P,0})/(\sigma_{P,max} - \sigma_{P,0})$ is shown in Fig. 9 as a function of the non dimensional vertical coordinate, $y/y_{\sigma_S}$ for the different studied cases. It is possible to observe that all measured results collapse in the range $0 \le y/y_{\sigma_S} \le 1.5$, exhibiting discrepancies for larger $y/y_{\sigma_S}$ values. The observed behavior is that of an initial increase, followed by a saturation and, then, a progressive relaxation. This behavior could be described by prescribing $\sigma'_P$ as a function of the non dimensional coordinate, $y' = y/y_{\sigma_S}$, such as:

$$ \sigma'_{P}(y') = y' \exp \left( \frac{(1 - (y')^n)}{n} \right), \qquad (31) $$

where $n$ is a parameter that controls the $\sigma'_{P}(y')$ relaxation as $y'$ increases. This proposed expression for the vertical evolution of the PDF forcing standard deviation
corresponds to a variable coefficient ordinary differential equation, i.e.,

$$y' \frac{d\sigma'_{P}}{dy'} - \sigma'_{P} = -(y')^{n}\sigma'_{P}. \qquad (32)$$

The relation between this equation and a reactive scalar transport equation will be explored in section 3.5.

The results of such model for two arbitrary values of $n$ (2 and 4) are also given in Fig. 9. These values are found to yield to $\sigma'_{P}(y')$ that encompass the measured values, thus leading to a full closure if $\sigma_{P,max}$ and $\sigma_{P,0}$ are known. More precisely, for $y' > 1$, the smallest and largest Reynolds number $\sigma_{P}(y')$ curves are bounded by Eq. (29) curves for $n = 4$ and 2, respectively. This inverse behavior suggests a simple ad-hoc model for the relationship between the Reynolds number and this decay parameter,

$$nRe = C_n, \qquad (33)$$

where the value $C_n = 1600$ has been used to close this relation. The maximum value, $\sigma_{P,max}$, is taken to be the function of $(Re, Da)$ determined by relating Eqs. (25) and (26),

$$\sigma_{P,max}^2 = \frac{a_{\tau}}{a_{\tau,\sigma}} \frac{Re^{3/4}}{Da^{1/10}}. \qquad (34)$$

The comparison between the proposed closure models and the measured PDF is given in Figs. 10a and 10b. The solid red lines correspond to $P(s)$ computed with the experimental values of $\mu_{2,S}$ and $\mu_{3,S}$, which were shown in section 3.3 to reproduce the experimental histograms exactly, and are thus identified here to the experimental PDF. Two models are considered. The first, Mod1, assumes the relationship between $\tau_P$ and $\sigma_P$ given by Eq. (26), thus using the measured $\sigma_P$ together with Eq. (29) to describe $P(s)$. The second model, Mod2, further assumes that $\sigma'_{P}(y')$ is given by Eq. (31), thus requiring the values of $y_{\sigma_S}$, $\sigma_{P,max}$, and $\sigma_{P,0}$ also. In these figures one may verify that the agreement between the first model and the experimental PDF is excellent. Such an agreement is also observed for the different cases studied, which leads to the conclusion that the closure given by Eq. (26) and the knowledge of the standard deviation of the forcing function ($\sigma_P$) suffice to determine the PDF. Concerning now the second model, where $\sigma'_{P}(y')$ is prescribed, the agreement with the experimental PDF is also very good for case A1, and for case D3 close to $y' = 1$. However, for this latter case, discrepancies are evidence both close to the burned surface and further downstream. In particular, the modeled standard deviation is larger than the experimental one, which leads to wider $P(s)$. Nevertheless, the model adequately describes the overall variation of the PDF shape for all the studied experimental conditions, which is not shown here for the sake of brevity.

Given the excellent agreement between the prescribed PDF and the schlieren experimental data, and that a this PDF requires a model for the PDF forcing standard deviation, the relation between this parameter and flowfield transported properties is now examined.
<figure><img src="image_11.png" /><figcaption>Figure 10. Comparison between experiment and models of the probability density function of the schlieren fluctuations $s(x, y, t)/\sigma_{S,max}$ at different heights above the burner surface, $y/y_{\sigma_S} = (0.5, 1, 1.5, 2)$, the corresponding $y$ values are given above each image.</figcaption></figure>

### 3.5. Relation to the reactive scalar transport

A relevant analysis point concerns the relation between the schlieren PDF forcing parameter, $\sigma_P$, and the standard deviation of a transported reactive scalar gradient.

In premixed combustion, the temperature is often described by a reactive scalar $(0 \le \theta = (T - T_u)/(T_b - T_u) \le 1)$ transport equation

$$ \rho \mathbf{v} \cdot \nabla \theta = \rho D \nabla^2 \theta + \rho \dot{\omega}, \qquad (35) $$

where $\rho$, $T$, $\mathbf{v}$, $D$ and $\dot{\omega}(\theta)$ are the mixture density, temperature, velocity, diffusion coefficient and chemical source term, respectively, whereas the initial and adiabatic flame temperature are $T_u$ and $T_b$. This equation is obtained from the energy transport equation under the classical hypothesis of constant thermodynamic pressure, low Mach number, heat transfer due temperature gradients only, unity Lewis number for all chemical species, and constant $\rho D$.

To mimic the schlieren separation in averaged and instantaneous components (Eq. (5)), it is assumed that the progress variable gradient may be written as $\nabla\theta = (1 + \vartheta)\nabla\mu_{\theta}$, where $\vartheta$ is a scalar perturbation to the average progress variable gradient.

Schlieren images result from the integration along the line-of-sight (direction $z$ here) through the refraction index field, $m(x, y, z, t)$, that is inversely proportional to the temperature, followed by the corresponding gradient component selection by the knife-edge on the focal plane [40]. The schlieren increment may thus be written as:

$$ dS(x, y, t) \propto \frac{1}{\theta} \nabla \theta \, dz \approx (1 + \vartheta) \frac{1}{\mu_{\theta}} \nabla \mu_{\theta} \, dz, \qquad (36) $$

where use has been made of the temperature gradient decomposition, and of the Taylor series expansion of the progress variable. The comparison of Eqs. (5) and (36) provides an analogy between the schlieren and temperature gradient fluctuations, $s$ and $\vartheta$ respectively. Even if this analogy is contingent on the used hypothesis, it enables to identify the processes controlling the fluctuations variations in the $y$ direction, and
thus of the schlieren standard deviation.

Accordingly, Eq. (35) may be rewritten as an instantaneous equation for $\vartheta$ that, after multiplication by $\vartheta$ and averaging, leads to an equation for $\mu_{\vartheta^2}$, the variance of $\vartheta$:

$$ \frac{1}{2} \frac{\rho D}{\mu_{\rho\dot{\omega}}} \nabla\mu_{\theta} \cdot \nabla\mu_{\vartheta^2} - \mu_{\vartheta^2} = - \frac{\mu_{\rho\dot{\omega}\vartheta}}{\mu_{\rho\dot{\omega}}}, \qquad (37) $$

where the average source term is noted $\mu_{\rho\dot{\omega}}$, and $\mu_{\rho\dot{\omega}\vartheta}$ is the correlation between the chemical source term and $\vartheta$. This equation shows that the $\mu_{\vartheta^2}$ production is controlled by this correlation, whereas the term $\frac{1}{2}\frac{\rho D}{\mu_{\rho\dot{\omega}}} \nabla\mu_{\theta}$ provides characteristic length scales. To compare with the model of the preceding section, Eq. (37) is integrated in the $(x, z)$ directions,

$$ \left\langle \frac{1}{2} \frac{\rho D}{\mu_{\rho\dot{\omega}}} \nabla \mu_{\theta} \right\rangle \frac{d}{dy} \langle \mu_{\vartheta^2} \rangle - \langle \mu_{\vartheta^2} \rangle \approx \left\langle \frac{1}{2} \frac{\rho D}{\mu_{\rho\dot{\omega}}} \nabla \mu_{\theta} \cdot \nabla \mu_{\vartheta^2} \right\rangle - \langle \mu_{\vartheta^2} \rangle = - \left\langle \frac{\mu_{\rho\dot{\omega}\vartheta}}{\mu_{\rho\dot{\omega}}} \right\rangle, \quad (38) $$

where $\langle \rangle$ represents this spatial integration.

Returning now to the ad-hoc model developed in section 3.4.2, the combination of Eqs. (30) and (32) leads to a model equation for the variance of the schlieren signal

$$ \frac{y'}{4} \frac{d\sigma_S^2}{dy'} - \sigma_S^2 = -(y')^n \sigma_S^2, \qquad (39) $$

which – given the analogy between $s$ and $\vartheta$ (Eqs. (5) and (36)) – suggests a closure expression to Eq. (38), i.e.,

$$ \left\langle \frac{2\rho D}{\mu_{\rho\dot{\omega}} y \sigma_S} \nabla \mu_{\theta} \right\rangle = y'. \qquad (40) $$

As shown in Fig. 11, this hypothesis is qualitatively supported by the vertical evolution of the average schlieren signal, at least in the range $0.5 \le y' \le 1.25$. The discrepancies that are observed below this range are due to the light beams that are deviated by the flame toward the burner surface, which effectively decreases the incident light level. The large variations among the different cases that are seen in this figure for $y' > 1.5$ arise in the burned gases region, as discussed in section 3.3 in connection to Figs. 4 and 5. Accounting for these variations is thus beyond the possibility of the variance model proposed.

The right hand side of Eq. (38) represents the correlation factor between $\rho\dot{\omega}$ and $\vartheta$, which extremum corresponds to the maximum $\sigma_S$, and that should decay to zero in the burned gases region faster than $\mu_{\rho\dot{\omega}}$ does. At this stage, a complete closure for this term should require assumptions regarding the velocity field, which has not been characterized. However, the simpler approach adopted, that is limited to the present burner configuration, is to write that

$$ \left\langle \frac{\mu_{\rho\dot{\omega}\vartheta}}{\mu_{\rho\dot{\omega}}} \right\rangle = f_y(y') \langle \mu_{\vartheta^2} \rangle. \qquad (41) $$

The closure of Eq. (38) is effected by prescribing $f_y(y') \propto (y')^n$, and by recognizing
<figure><img src="image_12.png" /><figcaption>Figure 11. Mean schlieren signal as a function of the height above the burner surface, y/y<sub>σS</sub>, for the different studied cases.</figcaption></figure>

that $\sigma_S^2 = \langle \mu_{\vartheta^2} \rangle$.

Given Eq. (30) and the PDF model (Eq. (23)), it is clear that this correlation between the source term and the fluctuating schlieren signal is directly related to the PDF forcing, $\sigma_P$. Improving the proposed ad-hoc model for this parameter could be done by learning from experimental data, for instance, which will be subject of future work. Finally, since the maximum variance value corresponds to $y' = 1$, Eqs. (30), (34) and (41) imply that

$$ \left\langle \frac{\mu_{\rho\omega\vartheta}}{\mu_{\rho\dot{\omega}}} \right\rangle_{\text{max}} \propto \frac{Re^{3/2}}{Da^{1/5}}. \qquad (42) $$

### 3.6. OH* stochastic behavior

Under the premise that the schlieren signal is closely related to temperature gradients, and that temperature is a controlling property in flames, it is interesting to examine the limitations of the prescribed PDF model for the description of the OH* fluctuations. Indeed, in the limit of very lean flames, it has been shown that all radicals follow a quasi-steady-state approximation [34, 41]. However, to the best of the author's knowledge, the relation between the schlieren and OH* has not been explored. Accordingly, Fig. 12 depicts the computed maximum values of $\nabla T$ and $\nabla T/T$ as a function of the OH* molar fraction maximum, $(X_{OH^*})_{max}$. These values have been obtained considering a steady adiabatic freely propagating planar one-dimensional flame for different equivalence ratio of the fresh unburned gases, $\phi \in (0.3, 0.6)$ – the smaller equivalence ratios corresponding to smaller $X_{OH^*}$ values. These temperature gradient expressions have been chosen because the schlieren can be shown to be proportional to $\nabla m/m$, where m is the refraction index of the medium, that is inversely proportional to the temperatures for gases [40]. This figure clearly shows that an increasing relationship exists between $(\nabla T)_{max}$, $(\nabla T/T)_{max}$ and $(X_{OH^*})_{max}$. More specifically, $(\nabla T/T)_{max} \propto (X_{OH^*})_{max}^{1/3}$, this 1/3 exponent corresponding to the depicted straight line in this log-log representation. Furthermore, for the equivalence ratios considered
<figure><img src="image_13.png" /><figcaption>**Figure 12.** Maximum values of $\nabla T$ (circles) and $\nabla T/T$ (squares) as a function of $(X_{OH^*})_{max}$ computed for a steady adiabatic freely propagating planar one-dimensional laminar flame. The straight line corresponds to a 1/3 derivative in this log-log representation.</figcaption></figure>

in the experimental part of this work, which correspond to $X_{OH^*} \in [5, 200] \cdot 10^{-14}$, $(\nabla T)_{max} \propto (\nabla T/T)_{max}$ is a good approximation. Note, also, that the two order of magnitude variation of OH* molar fraction in the studied range underscores the measurement difficulty associated to this species in these lean hydrogen-air flames. As a consequence, in this limit, it could be expected for the OH* stochastic behavior to closely follow that of the temperature gradient.

To verify this, Figs. 13a and 13b show the OH* fluctuations probability density function for cases A1 and D3, respectively, at different equally-spaced heights. Following what has been done for the schlieren PDF in section 3.3 to enable the comparison between these cases, the OH* fluctuation is divided by the corresponding standard deviation. In these figures the histograms (blue) correspond to the experimental PDF, whereas the curve (red) is the modeled result obtained by combining the measured standard deviation, $\sigma_{OH^*}(y)$, given in Fig. 5, with Eqs. (24)-(26). Both these figures show that the measured OH* PDF is significantly skewed, however the agreement between model and experiments is quite different. Indeed, the first of these figures, which correspond to an equivalence ratio of 0.28, show that this agreement is excellent, whereas in the second – obtained for $\phi = 0.35$ – discrepancies exist. More precisely, given $\sigma_{OH^*}(y)$, the OH* PDF skewness is fairly well described by the model relating the Reynolds and Damköhler numbers to the PDF forcing for the smaller equivalence ratios, whereas modeled skewness values that are larger than those measured arise for higher equivalence ratios. Such a behavior indicates that, even if the OH* fluctuations PDF seems to be described by Eq. (21), the functional dependency of the forcing and relaxation parameters on the relevant non dimensional numbers should be reexamined. The reasons for such a behavior change are unclear at present, but are subject of an on-going work, where computational fluid dynamics model results of the studied flames are expected to complement this experimental study. Note that, as the equivalence ratio increases from 0.28 (case A1) to 0.35 (case D3), the chemical pathways that control the combustion process have been shown to change [41]. More precisely, the chemical equilibrium hypothesis could hold for the OH* radical in the leanest mixtures, but not for the richer ones. Therefore, in the leaner cases the OH* fluctuations would be controlled by the temperature fluctuations, whereas a more complex dependency on other chemical species should arise for richer mixtures. Furthermore, the characteristic length scales of OH* and $\nabla T/T$ could be different as well, and could be at the origin
<figure><img src="image_14.png" /><figcaption>**Figure 13.** Comparison between experiment (blue histogram) and model (red curve) of the probability density function of the OH* at different equally-spaced heights above the burner surface, the corresponding y values are given above each image.</figcaption></figure>

of the observed discrepancies.

## 4. Conclusions

Hydrogen-air lean premixed laminar flames stabilized on an annular slit burner have been studied. The experimental conditions spanned over a range of Reynolds and Damköhler numbers close to the lean flammability limit. The studied flames, which were characterized via OH* chemiluminescence and time-resolved schlieren, were shown to exhibit unstable flame surface fluctuations along the direction transversal to the flow-field. The analysis of the transversely-averaged statistical moments of these schlieren images showed that both the schlieren variance and the skewness first increased with the height above the burner surface, then decreased downstream the flame tip. A similar behavior has also been evidenced for the OH* mean and variance.

The probability density function of the schlieren fluctuations was shown to be described by an exponentially modified Gaussian function, which is a skewed distribution that models the considered property decay when subject to a Gaussian forcing. Furthermore, this distribution parameters, that are linked to the measured statistical moments, were found to exhibit a peculiar relation, which is a function of the Reynolds and the Damköhler numbers. This relation, together with a hypothesis on the decay dependence on the Reynolds number, enabled to describe the schlieren fluctuations using the Gaussian forcing parameter only. To completely close the stochastic model with an ad-hoc, flow-field dependent, expression for this forcing was proposed. The connection between the corresponding equation and a reactive scalar transport equation enabled to relate the proposed closure to the scalar transport terms.

Finally, the extent to which the OH* chemiluminescence probability density function may be described by the proposed model was examined, showing that the model is adequate for the smaller Reynolds and Damköhler numbers values, but should be modified for larger ones. These results are expected to pave the way towards the development of hybrid reduced order models, in which the unmeasured flow-field dynamics may be learned from the experimental data.
Acknowledgements

This work has been performed while the author was on leave from Pontifícia Univers-
sidade Católica do Rio de Janeiro (PUC-Rio, Brazil). The author is grateful for the
experimental support provided by V. Montassier and A. Claverie, and to L. Pirateque
Henao and N. Lopes M. B. Junqueira for performing the premixed flame computations.

This work pertains to the French government program “Investissements d’Avenir”
(EUR INTREE, reference ANR-18-EURE-0010).

For the purpose of Open Access, a CC-BY public copyright license has been applied
by the authors to the present document and will be applied to all subsequent versions
up to the Author Accepted Manuscript arising from this submission.

Disclosure statement

The author has no conflicts to disclose.

Data Availability Statement

The data that support the findings of this study are available from the corresponding
author upon reasonable request.

References

[1] A.A. Konnov, A. Mohammad, V.R. Kishore, N.I. Kim, C. Prathap, and S. Kumar, *A comprehensive review of measurements and data analysis of laminar burning velocities for various fuel+air mixtures*, Prog. Energy Combust. Sci. 68 (2018), pp. 197–267.

[2] G. Sivashinsky, *Nonlinear analysis of hydrodynamic instability in laminar flames—I. derivation of basic equations*, Acta Astronaut. 4 (1977), pp. 1177–1206.

[3] G. Joulin and P. Clavin, *Linear stability analysis of nonadiabatic flames: Diffusional-thermal model*, Combust. Flame 35 (1979), pp. 139–153.

[4] P.D. Ronney, *Near-limit flame structures at low Lewis number*, Combust. Flame 82 (1990), pp. 1–14.

[5] F. Liu, X. Bao, J. Gu, and R. Chen, *Onset of cellular instabilities in spherically prop-
agating hydrogen-air premixed laminar flames*, Int. J. Hydrog. Energy 37 (2012), pp.
11458–11465.

[6] J. Goulier, A. Comandini, F. Halter, and N. Chaumeix, *Experimental study on turbulent expanding flames of lean hydrogen/air mixtures*, Proc. Combust. Inst. 36 (2017), pp. 2823–2832.

[7] R. Grosseuvres, A. Comandini, A. Bentaib, and N. Chaumeix, *Combustion properties of H₂/N₂/O₂/steam mixtures*, Proc. Combust. Inst. 37 (2019), pp. 1537–1546.

[8] F. Tinaut, M. Reyes, A. Melgar, and B. Giménez, *Optical characterization of hydrogen-air laminar combustion under cellularity conditions*, Int. J. Hydrog. Energy 44 (2019), pp. 12857–12871.

[9] W. Jin, J. Wang, Y. Nie, S. Yu, and Z. Huang, *Experimental study on flame instabilities of laminar premixed CH₄/H₂/air non-adiabatic flat flames*, Fuel 159 (2015), pp. 599–606.

[10] W. Jin, J. Wang, S. Yu, Y. Nie, Y. Xie, and Z. Huang, *Cellular instabilities of non-adiabatic laminar flat methane/hydrogen oxy-fuel flames highly diluted with CO₂*, Fuel 143 (2015), pp. 38–46.
[11] A. Kaewpradap and S. Kadowaki, Instability influenced by CO₂ and equivalence ratio in oxyhydrogen flames on flat burner, Combust. Sci. Technol. 189 (2017), pp. 438-452.

[12] W. Jin, J. Wang, W. Zhang, X. Cai, Y. Nie, and Z. Huang, Investigation of the heat loss effect on cellular flames via proper orthogonal decomposition, Combust. Sci. Technol. 190 (2018), pp. 803-822.

[13] T. Howarth and A. Aspden, An empirical characteristic scaling model for freely-propagating lean premixed hydrogen flames, Combust. Flame 237 (2022), p. 111805.

[14] L. Berger, A. Attili, and H. Pitsch, Intrinsic instabilities in premixed hydrogen flames: Parametric variation of pressure, equivalence ratio, and temperature. part 1 - dispersion relations in the linear regime, Combust. Flame 240 (2022), p. 111935.

[15] L. Berger, A. Attili, and H. Pitsch, Intrinsic instabilities in premixed hydrogen flames: parametric variation of pressure, equivalence ratio, and temperature. part 2 – non-linear regime and flame speed enhancement, Combust. Flame 240 (2022), p. 111936.

[16] D. Michaels and A.F. Ghoniem, Leading edge dynamics of lean premixed flames stabilized on a bluff body, Combust. Flame 191 (2018), pp. 39-52.

[17] L. da Costa Ramos, L.F. Figueira da Silva, F. Di Meglio, and V. Morgenthaler, Modelling of pulsating inverted conical flames: a numerical instability analysis, Combust. Theory Model. 26 (2022), pp. 260-288.

[18] F. Zhang, T. Zirwes, Y. Wang, Z. Chen, H. Bockhorn, D. Trimis, and D. Stapf, Dynamics of premixed hydrogen/air flames in unsteady flow, Phys. Fluids 34 (2022), p. 085121.

[19] J. Pathak, A. Wikner, R. Fussell, S. Chandra, B.R. Hunt, M. Girvan, and E. Ott, Hybrid forecasting of chaotic processes: Using machine learning in conjunction with a knowledge-based model, Chaos 28 (2018), p. 041101.

[20] L. da Costa Ramos, F. Di Meglio, V. Morgenthaler, L.F. Figueira da Silva, and P. Bernard, Numerical design of Luenberger observers for nonlinear systems, in 2020 59th IEEE Conference on Decision and Control (CDC), 12. IEEE, 2020, pp. 5435-5442.

[21] M. Raissi, P. Perdikaris, and G. Karniadakis, Physics-informed neural networks: A deep learning framework for solving forward and inverse problems involving nonlinear partial differential equations, J. Comput. Phys. 378 (2019), pp. 686-707.

[22] F. Chinesta, E. Cueto, E. Abisset-Chavanne, J.L. Duval, and F.E. Khaldi, Virtual, digital and hybrid twins: A new paradigm in data-based engineering and engineered data, Arch. Comput. Methods Eng. 27 (2020), pp. 105-134.

[23] S.L. Brunton, B.R. Noack, and P. Koumoutsakos, Machine learning for fluid mechanics, Annu. Rev. Fluid Mech. 52 (2020), pp. 477-508.

[24] M. Adams, X. Li, L. Boucinha, S.S. Kher, P. Banerjee, and J.L. Gonzalez, Hybrid digital twins: A primer on combining physics-based and data analytics approaches, IEEE Softw. 39 (2022), pp. 47-52.

[25] R. Geelen, S. Wright, and K. Willcox, Operator inference for non-intrusive model reduction with quadratic manifolds, Comput. Methods Appl. Mech. Eng. 403 (2023), p. 115717.

[26] V.A. Sabel'nikov and L.F. Figueira da Silva, Partially stirred reactor: study of the sensitivity of the monte-carlo simulation to the number of stochastic particles with the use of a semi-analytic, steady-state, solution to the pdf equation, Combust. Flame 129 (2002), pp. 164-178.

[27] E. Mendoza Orbegoso and L.F. Figueira da Silva, Study of stochastic mixing models for combustion in turbulent flows, Proc. Combust. Inst. 32 (2009), pp. 1595-1603.

[28] V.A. Sabel'nikov, B. Deshaies, and L.F. Figueira da Silva, Revisited flamelet model for nonpremixed combustion in supersonic turbulent flows, Combust. Flame 114 (1998), pp. 577-584.

[29] M. Pfitzner and P. Breda, An analytic probability density function for partially premixed flames with detailed chemistry, Phys. Fluids 33 (2021), p. 035117.

[30] D. Park, T. Lee, and K.T. Kim, Rotational symmetry-driven modal dynamics of high-frequency transverse instabilities in a lean-premixed multislit hydrogen combustor, Combust. Flame 245 (2022), p. 112356.

[31] M. Zhao, D. Buttsworth, and R. Choudhury, Experimental and numerical study of OH*
*chemiluminescence in hydrogen diffusion flames*, Combust. Flame 197 (2018), pp. 369–377.

[32] D.G. Goodwin, H.K. Moffat, and R.L. Speth, *Cantera: An object-oriented software toolkit for chemical kinetics, thermodynamics, and transport processes*, Caltech, Pasadena, CA 124 (2009).

[33] V.A. Alekseev, M. Christensen, and A.A. Konnov, *The effect of temperature on the adiabatic burning velocities of diluted hydrogen flames: A kinetic study using an updated mechanism*, Combust. Flame 162 (2015), pp. 1884–1898.

[34] A.L. Sánchez and F.A. Williams, *Recent advances in understanding of flammability characteristics of hydrogen*, Prog. Energy Combust. Sci. 41 (2014), pp. 1–55.

[35] R.J. Blint, *The relationship of the laminar flame width to flame speed*, Combust. Sci. Technol. 49 (1986), pp. 79–92.

[36] G. Joulin and T. Mitani, *Linear stability analysis of two-reactant flames*, Combust. Flame 40 (1981), pp. 235–246.

[37] M.C. de Jesus Vieira and L.F. Figueira da Silva, *Premixed flame heat release-based optimum global single-step chemistry for H₂, CH₄, and C₃H₈ mixtures with air*, J. Braz. Soc. Mech. Sci. Eng. 44 (2022), p. 131.

[38] I.G. McWilliam and H. Bolton, *Instrumental peak distortion. i. relaxation time effects*, Anal. Chem. 41 (1969), pp. 1755–1762.

[39] E. Grushka, *Characterization of exponentially modified gaussian peaks in chromatography*, Anal. Chem. 44 (1972), pp. 1733–1738.

[40] G.S. Settles, *Schlieren and Shadowgraph Techniques*, Springer Berlin Heidelberg, 2001.

[41] D. Fernández-Galisteo, A.L. Sánchez, A. Liñán, and F.A. Williams, *The hydrogen-air burning rate near the lean flammability limit*, Combust. Theory Model. 13 (2009), pp. 741–761.