Numerical Analysis of the Combustion Dynamics of Passively Controlled Jets Issuing from Polygonal Nozzles

: In the present work, the combustion of vitiated hydrogen jets issuing from differently shaped nozzles is modelled using the LES method. We investigate the impact of nozzle cross-sectional geometries (circular, square, triangular, hexagonal and hexagram) and the jet Reynolds numbers (Re = 18,000, 20,000 and 23,600) on the ﬂame lift-off height, its structure, the locations of the temperature maxima and species distributions. The triangular nozzle yields the highest mixing rate and therefore the fastest decay of axial velocity and the fastest growth of the average temperature along the ﬂame axis. It was found that for the largest Re, the zone of intense mixing and the reaction zone occur in distinct regions, while for the lower Re, these regions combine into an indistinguishable zone. Finally, it is shown that the lift-off height of the ﬂames and the mean temperature ﬁeld are non-linearly correlated with Re and strongly dependent on the nozzle shape. We use the LES method for modelling the of vitiated hydrogen the co-ﬂow. We investigate the impact of the nozzles’ cross-sections geometry (circular, square, triangular, hexagonal and hexagram) and jet Reynolds numbers (Re = 18,000, 20,000 and 23,600) on L h , the shape of the ﬂame and the of method for pressure-velocity coupling on half-staggered meshes. The convective terms in Equations (4) and (5) are discretized using the second order total variation diminishing scheme with van Leer’s limiter. Time integration was performed with the use of a predictor-corrector scheme that combines second order Adams–Bashforth/Adams–Moulton methods.


Introduction
In the present work, we focus on the application of passive control techniques to reactive jets. They are used in practical devices, e.g., burners or combustion chambers, to control jet dynamics and, indirectly, combustion characteristics. The passive methods rely in particular on the geometrical shape of jet nozzles and aim to meet the conditions required for a given configuration, e.g., an attached or lifted flame, the most uniform spatial temperature distribution, or pollution reduction. In this work, we consider hydrogen jet flames issuing from four polygonal nozzles, namely hexagonal, hexagram, square and triangular. Referring to the base flow configuration (a circular flame), we analyse the impact of nozzle geometry on the flame dynamics including the flame structure and lift-off heights altered by inlet flow conditions. This study is mainly motivated by the observation that non-circular jets, e.g., emanating from rectangular or elliptical nozzles, have a tendency to enhance the mixing process between the jet and the surrounding flow. In the case of non-reactive flows, this aspect has been extensively studied for many years, both experimentally [1,2] and numerically using LES (Large Eddy Simulation) [3,4] and DNS (Direct Numerical Simulation) [5,6] approaches. A comprehensive review concerning the theoretical issues of non-circular jets, in particular rectangular and elliptic, along with a description of their applications can be found in [7].
The research on the impact of fuel nozzle geometries on reactive flows was mainly conducted experimentally. Historically, Gutmark et al. in [8,9] studied the dynamics of jet burners with corners. In the latter work, it was found that their effect on jet dynamics was similar to cold flows. The small-scale structures emanating from sharp corners were shown to interfere with the large-scale structures originating from the flat sides. These combined effects were found to be beneficial for combustion. A comprehensive review of earlier works focusing on passive combustion enhancement can be found in [10]. Iyogun et al. [11] investigated the flames issuing from five types of nozzles, i.e., circular pipe, contracting circular, triangular, rectangular and square. Compared to the circular nozzle, it was found maxima, among others. It is shown that depending on the jet Reynolds number, the flames are characterised by distinct behaviour close to the nozzle.

Computational Approach
In the present study, the results were obtained using two different numerical LES solvers. First, we used the second order accurate ANSYS Fluent solver [20] for modelling the flow inside nozzles and generating the outlet velocity signals. In the next step, they were applied as inlet boundary conditions for the main computations performed using a SAILORin-house high-order compact difference LES solver. We used two different eddy viscosity models. For simulations in the Fluent solver, the Wall-Adapting Local Eddyviscosity (WALE) subgrid model of Nicoud and Ducros [21] was used since it is well suited for wall-bounded flows. In the SAILOR code, the turbulent viscosity was calculated using the Vreman model [22] since it is characterized by excellent accuracy, especially for jet type flows [23].
The SAILOR code solves the spatially and Favre filtered set of the governing equations assuming the low Mach number approximation. The continuity and Navier-Stokes equations complemented with the equation of state in the framework of the LES method are given as: ∂ρ ∂t The bars and tildes denote the spatially and Favre filtered quantities. The symbols p 0 and R stand for the thermodynamic pressure and the mixture gas constant, respectively. The variables: u, ρ, p and T denote the velocity, density, hydrodynamic pressure and temperature, respectively. The termsτ ij and τ sgs ij represent the viscous and the subgrid stress tensors, respectively. The former is defined asτ ij = µ ∂ũ i and the latter as τ sgs ij = 2µ sgs S ij in which S ij is the Favre filtered rate of strain and µ sgs is the subgrid scale viscosity, modelled using Vreman's model [22]. The species mass fractions (Y α ) and the enthalpy (h) transport equations are given by: ∂ρh ∂t where the subscript α denotes the species index from one to N-species. The symbols σ and σ sgs refer to the Prandtl/Schmidt numbers and their turbulent analogues. They are assumed to be equal to 0.7 and 0.9, respectively. The termsω α (Y, h) are the production/destruction terms of the species, and they are computed using the laminar closure assumingω α (Y, h) ≈ ω α Y, h . It is worth noting that in [24], this approach was termed the Implicit LES (ILES), but it should not be confused with the popular ILES method with implicit modelling of the subgrid terms. The terms in Equations (1) and (2) and the diffusion terms in Equations (4) and (5) are discretized using the sixth order compact difference method combined with the projection method for pressure-velocity coupling on half-staggered meshes. The convective terms in Equations (4) and (5) are discretized using the second order total variation diminishing scheme with van Leer's limiter. Time integration was performed with the use of a predictorcorrector scheme that combines second order Adams-Bashforth/Adams-Moulton methods.
The applied ILES method for the chemical source terms was proven to yield correct results for sufficiently dense meshes that ensure appropriate resolution of the scalar fields as shown in [25,26] using the SAILOR code. The chemical reaction rates are modelled using the detailed mechanism of Mueller et al. [27] consisting of nine species and 21 reactions for the hydrogen oxidation. A detailed description of the numerical solver can be found in [28,29]. The applied code has been used and verified in numerous studies of non-reacting [30][31][32] and reacting flows [25,[33][34][35].

Computational Configuration
As was mentioned before, the simulations were performed using a two-stage approach. In the first stage, the ANSYS Fluent LES solver was used to model the flow within each investigated nozzle. A view of the computational mesh applied on the nozzles is presented in the upper part of Figure 1. Each nozzle is comprised of a straight rod with a length of 25D e and starts with a circular cross-section (up to 5D e ), which smoothly converts into the final shape over a distance of 5D e (see Cross-section 2-2 in Figure 1) and then continues up to a distance of 15D e . The computational meshes generated inside the pipes were block-structured and precisely matched the inner wall shape of each nozzle, i.e., circular, hexagonal, square, hexagram and triangular, as can be seen in the lower part of Figure 1. The computational meshes were discretized to achieve the near-wall cell height y + < 1 for a proper resolution of the turbulent boundary layers developing in the pipes. The unsteady velocity signals were acquired at the nozzle outlet planes for a period of 300D e /U j . To avoid the impact of the outlet boundary conditions on the velocity signal imposed at the inlet boundary plane of the main computational domain, the instantaneous velocity components were extracted from the plane placed 2D e upstream, i.e., at y/D e = 23.
In the second part of the simulations, the velocity signals were superimposed onto the inlet of the computational domain shown in Figure 2. It can be seen that the nozzles are surrounded by a solid wall (black circle with a diameter of 1.64D e ) that is assumed here to mimic a real situation where the nozzles would be manufactured inside the long rods of a given diameter. The computational domain was discretized using two mesh resolutions, i.e., 192 × 264 × 192 and 240 × 312 × 240 in the "x", "y" and "z" directions. Additionally, in both meshes, the regions close to the inlet plane and towards the jet axis were refined to account for the high gradients occurring in the shear layer. The details of the mesh size used for each case are given in Table 1. Outside the nozzle wall, a uniform co-flowing stream was assumed, with the velocity equal to U c = 0.03U b . At the side boundaries, the vertical velocity (in the y-direction) was assumed to be equal to U c , and the remaining velocity components were assumed equal to zero, i.e., there was no flow through the side boundaries. The pressure at these boundaries was computed from the Neumann condition n · ∇p = 0 with n being the outward normal vector. With these settings, the side boundaries can be treated as moving walls. At the outlet plane, the velocity components were computed from a convective boundary condition preventing the back-flow, and the pressure was set constant to 1.01325 × 10 5 Pa. Except for the solid element surrounding the nozzle and the fact that in the present work, we consider polygonal nozzles, the following flow conditions closely correspond to a setup used by Cabra et al. [36] who analysed a lifted hydrogen diluted flame originating from a thin-wall circular nozzle. The fuel jet is a mixture of hydrogen and nitrogen with the molar concentrations equal to X H 2 = 0.25 and X N 2 = 0.75, respectively. The jet initial temperature is equal to T j = 305 K, and the jets are issued into the hot co-flow (T c f = 1045 K) composed of O 2 (X O 2 = 0.1474), H 2 O (X H 2 O = 0.0989) and N 2 (X N 2 = 0.7534). In the study of Cabra et al. [36], the Reynolds number was equal to Re = 23,600. In the present study, besides Re = 23,600, we additionally consider the cases with Re = 18,000, 20,000 that allow us to analyse the impact of the Reynolds number on the flame dynamics. Furthermore, for validation purposes, we analyse the non-reacting flow at Re = 23,600. All the details concerning the simulated cases along with their abbreviations that will be used further in the text are given in Table 1. Table 1. Parameters of the computational cases for non-reacting (C, cold) and reactive (H, hot) jets with different nozzle geometries, specified by the Reynolds number Re, and the size of the computational grids (F, fine). The cells sizes are expressed as: ∆x * = ∆x/D e , ∆y * = ∆y/D e and ∆z * = ∆z/D e with ∆x, ∆y and ∆z taken from the non-uniform grid at x = y = z = 0.

Case Name
Nozzle  Figure 3 shows a comparison of the present time-averaged solutions obtained for the non-reactive case with available experimental data along the axis. In these simulations, the jet and the co-flow temperature were equal to 305 K. The mean velocity profiles were normalised by the jet velocity U j measured at the nozzle centre point. The velocity fluctuations were normalised by U j for the circular jet and by the centreline velocity U cl for the square nozzle.

Model Validation
In Figure 3a, the LES results obtained for the circular nozzle on the coarse (Ci-23C) and fine mesh (Ci-23CF) are compared with the LES results obtained by da Silva and Métais [37] (Re = 1.0 × 10 4 ) and the experimental data of Crow and Champagne [38] (Re = 1.03 × 10 5 ), Zaman and Hussain [39] (Re = 3.2 × 10 4 ) and Aleyasin et al. [40] (Re = 1.0 × 10 4 (a) and 2.0 × 10 4 (b)), who analysed the jet dynamics in the range Re = 6.0 × 10 3 − 2.0 × 10 4 . As was previously reported by Ball et al. [41], they confirmed that in circular jets, the impact of the Reynolds number on the mean velocity decay and the potential core length was weak when Re ≥ 1.0 × 10 4 . In fact, their results for Re = 1.0 × 10 4 and 2.0 × 10 4 were very close to each other. The discrepancies between the results reported by different authors originated from the differences in the details of nozzles (e.g., nozzle length or contraction) and the inlet turbulence characteristics (e.g., turbulence intensity, turbulent length scales) rather than from the Reynolds number. In general, the agreement of the mean velocity profile with the literature results was satisfactory, especially with the data reported by Aleyasin et al. [40] for Re = 20,000. This is most likely because the turbulence intensity distribution at the nozzle exit was in this case very similar to that in the present work. Regarding the velocity fluctuations, up to a distance of y/D e = 7, the agreement between the present simulation results and experimental data was very good. Further downstream, the simulation results seemed to slightly underestimate the velocity fluctuation level, which was caused by the fact that in a fully turbulent flow region, a portion of the energy was contained in the unresolved sub-grid flow scales interacting with discretization errors [42].  Similar conclusions can be drawn considering the results from the simulations of the jets issuing from the square nozzle using either a coarse (Sq-23C) or a fine (Sq-23CF) mesh, shown in Figure 3b. In this case, the simulations results were validated against the experimental data of Quinn and Militzer [43] (Re = 1.85 × 10 5 ) and Xu et al. [44] (Re = 2 × 10 4 ). In terms of the time-averaged axial velocity profiles along the jets axis, the agreement between LES and experimental results of Xu et al. [44] was very good. Both the length of the potential core and velocity decaying rate along the jet axis agreed with the experiment very well. Regarding the velocity fluctuations, a good agreement was achieved up to y/D e = 12. Further downstream, they appeared to be slightly lower compared with the results of Xu et al. [44]. As in case of the circular jets, the small discrepancies with the literature data can be attributed to minor differences in the jet parameters at the inlet, i.e., the shape of the mean velocity profile, momentum thickness and turbulence intensity. An important fact is that for the jet issuing from both the circular and square nozzle, the computations performed on the fine and coarse mesh led to very similar results. Although in the reacting flows, the mesh density or the numerical method can have a larger impact on the solutions, as observed by Cocks et al. [45], only the coarse mesh was used for further investigations, without performing the grid dependence analysis. Such studies were carried out by Rosiak and Tyliszczak [34] and Tyliszczak [33] for similar flames (a vitiated hydrogen in a hot co-flow at Re = 4000-23,600) and using the same numerical code. There, it was shown that the applied high-order discretization method provided accurate and mesh-independent results even on meshes coarser than those used in the present research.

Impact of the Nozzle Geometry
To investigate how the nozzle shape affects the flow field, the time-averaged axial velocity and its fluctuations are presented in Figure 4a. Clearly, the different nozzle shapes affect both the near-field and far-field of the jet. The length of the potential core (L c ) taken as the axial distance from the nozzle exit to the point where U / U j = 0.95, given with respect to the nozzle shape, is as follows: triangular (y/D e = 3.82), hexagram (y/D e = 4.03) square (y/D e = 4.09), circular (y/D e = 4.52) and hexagonal (y/D e = 4.82). Interestingly, the longer the potential core, the faster the decaying rate of axial velocity in the intermediate-field. Consequently, from the location y/D e = 12 downstream, the velocity profiles coming from different nozzle shapes converged to the same value. The velocity fluctuations were very similar for all the cases up to y/D e = 3. After this point, the values obtained from triangular, hexagram and square nozzle showed a higher level of turbulence intensity in comparison to circular and hexagonal jets. The results from all cases converged to a similar value around y/D e = 10, and from this point, the fluctuations seemed to be at an approximately constant level.   Figure 4b shows the entrainment rate, which is a measure that allows quantifying the mixing of the fluid from the surroundings into the jet. The mixing efficiency plays a key role in combustion devices, i.e., burners, combustion chambers, etc. Enhanced mixing allows faster fuel conversion and hence makes the devices more compact. In this study, following Wygnanski and Fiedler [46], the entrainment rate is computed as: where Q is the volume flux, c is a path outside the jet at a given downstream location, u and n are the velocity and normal vectors, respectively, u r is the radial velocity component and θ is the azimuthal coordinate. The side boundaries in the present configuration are assumed as moving walls; therefore, the entrainment calculated along the external boundaries would be zero, as there is no flow through them. Therefore, the formula given in Equation (6) is applied only for a part of the domain surrounding the vicinity of the jet by a cylindrical boundary with the radius equal to r = 4D e . As can be seen in Figure 4b, the jets issuing from the circular and hexagonal nozzles result in the lowest level of entrainment. The biggest difference between the cases is around the location y/D e = 9, where the entrainment rate was equal to 0.0267, 0.0281, 0.0297, 0.0301 and 0.0304 for the hexagonal, circular, square, hexagram and triangular nozzles, respectively. Those differences can lead to substantial changes in the mixing process and affect the combustion dynamics, which will be investigated in the following parts of the paper.

Reacting cases
In this section, the effect of nozzle geometry on the flame dynamics is investigated at Re = 23,600 for all nozzle shapes and with Re = 18,000 and 20,000 for the square, triangular and circular nozzles as in these cases, the results differed the most.
Independently of the Reynolds number and the nozzle shape, the flames showed complex, turbulent behaviour. As an example, Figure 5 presents the isosurfaces of the Qparameter (Q = 10(U j /D e ) 2 ) coloured by the instantaneous temperature values obtained for the triangular nozzle cases Tr-18H, Tr-20H and Tr-23H. The Q-parameter is one of the quantities that is commonly used as a good indicator of vortical structures. It is defined as Q = 1/2(Ω ij Ω ij − S ij S ij ) where S ij and Ω ij are symmetrical and anti-symmetrical parts of the velocity gradient tensor, respectively. Here, it shows the formation of large-scale vortices near the nozzle exit and how these vortices evolve flowing downstream. It can be seen that depending on the Reynolds number, their strength and location where they break down into the smaller vortices change. The vortical structures near the nozzle are the least coherent for Re = 23,600. This is the manifestation of the interactions between the flow in the jet centre and the small-scale rib-like vortices formed in the jet corners. They lead to the destruction of the toroidal vortices, which for Re = 23,600 seem to be weaker compared to the cases with Re = 18,000, 20,000. One can also observe that the Reynolds number affects both the region occupied by the flame and the position at which the flame stabilises above the nozzle. For Re = 23,600, its anchoring point oscillates around y/D e = 12.5, whereas for Re = 18,000, 20,000, the flames are seen already in between y/D e ≈ 2-4. The instantaneous and time-averaged temperature contours for Re = 23,600 are presented in Figure 6. Following Stankovic et al. [47], the lift-off height is defined as the location of the first occurrence of the OH species with the mass fraction above the value 2 × 10 −4 . L h values are equal to 12.56, 12.78, 12.95, 13.54 and 13.59 for the triangular, hexagram, circular, hexagonal and square nozzles, respectively.
As by visual inspection of the temperature contours, it is hard to grasp much of a difference between the jets, we present the distributions of the time-averaged velocity and temperature with their fluctuations along the jet axis in Figure 7. One of the conclusions drawn from the non-reactive cases was that the use of the triangular nozzle results in the highest entrainment rate. This is also true for the reactive case. As can be seen in Figure 7a, the triangular jet shows the fastest decay of the axial velocity along the axis. The square and hexagram nozzles also show a stronger velocity decay rate in comparison to the circular and hexagonal nozzles. Regarding the time-averaged temperature along the axis (see Figure 7b), one can notice that the results for the triangular nozzle case again differ from other results. Enhanced mixing leads to a quicker transport of the co-flow gases to the jet region and an increase of the temperature in the jet axis. This effect can also be seen for temperature fluctuations. In all cases, there are two peaks: the first peak is related to the location where the inert mixing with co-flow occurs, whereas the second one is observed in the region where the chemical reactions take place. Both peaks for the triangular case are shifted towards the nozzle exit, which indicates faster mixing and shorter flame lift-off height. The symbols in Figure 7b denote the experimental data of Cabra et al. [36], which can be compared with the case Ci-23H. Although in the present configuration, the nozzle is surrounded by the thick wall, the agreement is very good. This concerns both the region of the fuel/oxidizer mixing in between y/D e = 5-15, where the temperature increases almost linearly, and the fully burning regime in the far-field y/D e > 20, where the flame behaviour is less dependent on the inlet conditions and the temperature is very high.   The profile of the H 2 O species shows that up to y/D e = 13, its level increases, though the combustion does not occur in this region yet. In such a case, the H 2 O species can originate only from the co-flowing stream, and therefore, it can be used as a measure of mixing in the jet axis. It turns out to be the most intense for the jet issuing from the triangular nozzle. Further downstream, from around y/D e = 15, a rapid increase of the H 2 O species is observed. The fact that it rises in an analogous way as the temperature values (see Figure 7b) confirms that at this particular location, the chemical reactions are initiated. Analysis of the OH mass fraction shows that a direct consequence of the intensified mixing is that the auto-ignition and the reaction zone occur closer to the inlet. Compared to the remaining cases, the peak of OH mass fraction for the case Tr-23H is shifted upstream from y/D e ≈ 19.9 to 18.2. Noticing that the nozzle shape changes the flame base location, one can presume that it also affects its downstream behaviour in the radial direction. Figure 9a shows the time-averaged velocity profiles along the x-direction at the streamwise location equal to y/D e = 20. As can be seen, the results are similar for all nozzle shapes, though the triangular and square nozzles seem to cause the flames to be slightly wider. Regarding the radial temperature distributions shown in Figure 9b, one can notice that when the triangular nozzle is used, the temperature in the axis is the highest and the flame is also the widest. For the circular, hexagonal and square jets, the temperature drops to the coflow temperature around a radial distance of x/D e = 4, whereas for the hexagram and triangular cases, the region of higher temperature extends up to x/D e = 5.3. Larger differences are found in the species distributions. The radial profiles of the mass fractions of OH and H 2 at the location y/D e = 20 are presented in Figure 10a,b, respectively. The concentration of OH seems to differ only near the jet axes, while for x/D e > 2, the profiles converge to virtually the same value. This is also the case for the H 2 species, though, in this case, the differences between particular solutions are more pronounced. This is especially noticeable near the axis where for the triangular nozzle (Tr-23H) hydrogen mass fraction is much smaller than in the remaining cases, almost five times compared to the hexagonal jet (He-23H). This is related to the following facts: (i) L h in the case Tr-23H is the shortest, and this causes hydrogen at the jet periphery to start to be consumed by the reaction closer to the nozzle; (ii) the velocity in the case Tr-23H decays the fastest (see Figure 7a) as a manifestation of the most intense mixing. This translates into the fastest rise of temperature in the mixing zone in the region of y/D e = 5 − 15 (see Figure 7b). The small level of H 2 at y/D e = 20 results from the combined effect of the chemical reactions that started earlier and the mixing with the co-flow that dilutes the fuel in the centre of the flame. The differences in the results obtained for Re = 23,600 are mainly quantitative and not as substantial as one could expect taking into account the evident differences in the nozzle shapes. This is consistent with the observations of Iyogun et al. [11], who reported that nozzle shape is important mainly at low Reynolds numbers. The results for Re = 18,000 and 20,000 are presented only for the circular, square and triangular nozzles that suffice to show how the flames respond to changes in Re and which of them are more prone to such changes. Figure 11 shows the instantaneous and time-averaged temperature distributions for all three Reynolds numbers. The first thing that can be noticed is that decreasing the Reynolds number from 23,600 to 20,000 leads to a drastic change in L h . It is reduced from y/D e = 12.95 to 3.50 (cases Ci-23H → Ci-20H), from 13.59 to 3.23 (Sq-23H → Sq-20H) and from 12.56 to 5.29 (Tr-23H → Tr-20H). This is somewhat surprising since in the case Tr-23H, the L h was the shortest, while for Tr-20H, it was the highest.
The analysis of the results obtained using the circular nozzle (Figure 11a) reveals that the reduction of the Reynolds number from 23,600 to 20,000 (Ci-23H → Ci-20H) resulted in a drastic change of L h . However, reducing the Reynolds number to 18,000 did not lead to any substantial differences in L h . On the other hand, reducing the Reynolds number for the jets issuing from the square and triangular nozzles (cases Sq-20H→Sq-18H and Tr-20H→Tr-18H) causes remarkable differences, as can be seen in Figure 11b,c. First of all, one can observe that in these cases, the flames are discontinuous. They are almost attached to the nozzle; however, further downstream, they are extinguished, and the new flame fronts are formed in close vicinity. Most likely, this is the effect of the small-scale mixing between the cold fuel stream and hot oxidiser. It could be seen in Figure 5 that deformed toroidal structures appear in this region. Close inspection of the results shows that they rotate in a way that the cold fuel from the jet core is transported outside the jet region causing local blow-off of the flame formed in the closest proximity of the nozzle. Note that for Re = 20,000, this phenomenon occurs also for the square nozzle, but not for the triangular one. This shows how much the lift-off height and flame stability in its initial phase are sensitive to the inlet conditions, both related to the flow parameters, as well as to the nozzle shapes, particularly their corners. The regions of an intense reaction and unsteady flame behaviour caused by corners were observed also in [8,19].   To further investigate the flame discontinuity, a temporal evolution of the temperature was measured at a single point P (coordinates: x/D e = 0.8, y/D e = 4.7, z/D e = 0) shown in the middle of Figure 11c as a green circle. Its location corresponds to the anchoring point of the flame above the gap that is seen in the time-averaged results for the case Tr-18H. The recorded time signals for the cases Tr-18H and Tr-20H are presented in Figure 12 where the red line denotes the co-flow temperature. As one can see, temperature variability is large, and depending on the Reynolds number, it rises more or less frequently above T c f . The distinct peaks denote the time instances when, due to the chemical reactions, the heat release is large. Evidently, more such points are observed for the case Tr-18H. For Tr-20H, they appear rarely, and in this case, the temperature is mostly below T c f . In effect, only in the case Tr-18H, the time-averaged temperature at the point P indicates the occurrence of ignition.
The impact of the Reynolds number on the jet behaviour is studied in the following in detail based on the time-averaged temperature and its fluctuating axial signals normalised by T c f and T cl , presented in Figure 13. The time-averaged temperature profiles for Re = 23,600 exhibit three different regions that can be distinguished. The first one, indicated in Figure 13 as (I), corresponds to the unmixed flow inside the potential core (y/D e < 5) where the temperature is constant. Then, in the region (II) localised in between 5 < y/D e < 15, the mixing of the fuel with the co-flow takes place, and the temperature rises almost linearly. Further downstream, the region (III) is found, inside which the slope of the temperature profiles suddenly changes. This is caused by the occurrence of the ignition spots in this part of the flow, as can be seen in Figure 11. The fact that this process is strongly unsteady does not change the qualitative behaviour of the mean temperature. Up to half of the region (III), the axial temperature is lower than T c f such that the convective transport of the cold fuel affects the mean temperature more than the chemical reactions. Depending on the nozzle, the combustion process ( T /T c f > 1) prevails only from y/D e ≈ 16. At the end of the region (III), the flame is fully developed. When the Reynolds number decreases, the region (I) still exists, but the regions (II) and (III) cannot be distinguished. The analysis of the results presented in Figure 11 shows that for Re = 18,000 and 20,000, the high temperature regions near the axis appear already at y/D e =6-7. This indicates that the temperature rise in the region y/D e > 7 has the same character as that found for Re = 23,600 in the second part of the region (III) (y/D e = 12-15). It can be said that depending on Re and the nozzle type, the convective fuel/oxidiser mixing influences the combustion process in a different way. Tr-18H Tr-20H The decrease in the Reynolds number changes not only the flame development mechanism. It also affects the rate of the temperature rise. As was already discussed, for Re = 23,600, the fastest growth of the temperature along the axis (see Figure 13a) is observed for the cases Tr-23H, Sq-23H, and Ci-23H, respectively. When the Reynolds number is reduced to 20,000, the situation changes, and in this case, the order is as follows: Tr-20H, Ci-20H and Sq-20H. Further reduction of the Reynolds number alters the sequence again, i.e., for Re = 18,000, the sequence is Ci-18H, Sq-18H and Tr-18H. It is somewhat surprising that for the square and triangular jets, the mean temperature growth is faster for Re = 20,000 than for Re = 18,000. The differences between the temperatures in the transient region (y/D e ≈ 7-20) are the largest for Re = 18,000. For instance, at y/D e = 15, we have T /T c f = 1.15, 0.95, 1.06 for Ci-18H, Sq-18H and Tr-18H, respectively.
The temperature fluctuations along the jet axis shown in Figure 13b reflect the differences observed in the mean temperature profiles. Two distinct peaks can be observed for Re = 23,600. The first one is located at y/D e ≈ 5, i.e., in between the regions (I) and (II). Note that before the first maximum is reached, the temperature fluctuations grow the fastest for the triangular cases, which seems to confirm their ability to improve the mixing process the most. The second maximum occurs at y/D e = 16, which can be related to the beginning of the zone inside the region (III) where T /T c f > 1. These peaks correspond to the borders between different mixing/reaction regions discussed above. For Re = 20,000 and 18,000, the fluctuations exhibit significantly different behaviour. Only single peaks are seen at locations y/D e = 6-8, and their levels depend on the Reynolds number and nozzle shape. The maxima change in the following orders: Tr-20H, Ci-20H, Sq-20H and Ci-18H, Tr-18H, Sq-18H. The former sequence coincides with the order at which the mean temperature grows, while the latter does not. This supports a common view that the instantaneous flame dynamics, the strength of which is often expressed in terms of temperature fluctuations, cannot be one-to-one related to the mean temperature field.

Conclusions
The present study focused on the LES-based analysis of the impact of nozzle geometry on the dynamics of a diffusion hydrogen flame. Five nozzle shapes (circular, square, triangular, hexagonal and hexagram) and three Reynolds numbers (18,000, 20,000 and 23,600) were investigated. The analysis was preceded by a detailed verification of the results for the non-reactive cases and also for the flame, though only partially in this case. In the non-reacting cases, it was shown that nozzle geometry affects the mean velocity field, its fluctuations and the entrainment rate, which is in line with the research of Zhou and Gutmark [8,15].
The flame lift-off height for all nozzle shapes for Re = 23,600 was located approximately at y/D e = 12.5 − 13.6, which was outside the region where the largest differences in the entrainment rate were found. Analysis of the H 2 O species concentration showed that in the jet issuing from the triangular nozzle, the mixing process was the most intense. This caused the movement of the reaction zone closer to the inlet and almost a 2D e shift in the peak of the OH mass fraction. For Re = 20,000, the lift-off height was significantly reduced, i.e., up to y/D e ≈ 3.5 for the circular and square nozzle and to y/D e = 5.29 for the triangular jet. Such flame behaviour agrees well with the observation of Iyogun et al. [11], who reported that nozzle shape is important mainly at low Reynolds number values. Its further decrease to 18,000 did not lead to significant differences in the case of the circular nozzle. However, for the triangular and square nozzles, the flames became almost attached to the nozzle exit, and the changes of their instantaneous topology were very dynamic and dominated by spontaneously occurring local extinctions. Analysis of the time-averaged temperature and its fluctuations along the axes led to the observation that for Re = 23,600, there are three distinct flow regions along the jet axis: the region, where the jet temperature is constant, the region where the mixing with the hot co-flow occurs and the third region, where chemical reactions take place. Reducing the Reynolds number to 20,000 and 18,000 resulted in a sudden change of the flow characteristics. The first region with the constant temperature (no-mixing zone) still existed; however, the zone where the mixing with the co-flow occurred and where chemical reactions took place were combined. It seems that this happens due to the large values of temperature fluctuations at the end of the first region. This causes a rapid temperature rise accompanied by a simultaneously occurring combustion process, making the second and third region hardly distinguishable.
In conclusion, the present study shows that nozzle geometry plays an important role and can affect the combustion process significantly. It was observed that the impact of nozzle shape on flame dynamics is strongly dependent on the Reynolds number, and a slight change of its value can lead to significant differences in the flame lift-off height, the occurrence of local extinctions, the temperature and the species distributions. In further research, we will analyse to what extent the present results are dependent on the turbulent combustion model that are already available in our code (conditional moment closure or Eulerian PDF model) and whether its impact on the flame dynamics is not larger than the nozzle shape.

Data Availability Statement:
The data presented in this study are available on request from the corresponding author.