2.1 Orifice flow revisited

For incompressible fluid flow (at Mach number <0.3), the orifice flow Q_orf is commonly described in terms of the isentropic flow rate derived from Bernoulli's equation,

(1)\begin{equation}{Q_{orf}} = {C_d} \times {A_{orf}} \times \sqrt {\textrm{ }2\Delta {P_{orf}}/\rho } ,\end{equation}

where ρ is the fluid density, A_orf is the orifice cross-section and C_d is the discharge coefficient, which is around 0.7 typically (Reference Borutzky, Barnard and ThomaBorutzky, Barnard, & Thoma, 2002; Reference LaffertyLafferty, 1998; Reference Oertel, Prandtl and BöhleOertel, Prandtl, & Böhle, 2008), but may exceed 1 due to the presence of a cross-flow (Reference Berger, Gourdain, Bauerheim and DevillezBerger, Gourdain, Bauerheim, & Devillez, 2019). Isentropic flow occurs under adiabatic conditions and, for pressure drops ΔP_orf < 200 Pa, the adiabatic cooling is expected to be less than 5 K (for diatomic gas) in horizontal flame gas sampling tubes under atmospheric conditions. Accordingly, ρ changes by less than 0.2 % and can be regarded as constant. Equation (1) predicts the orifice flow to be proportional to the square root of the pressure drop and to be independent of the gas viscosity. However, despite the vast number of publications about orifice flow, to the authors’ knowledge a high-temperature validation (T > 700 K) of (1) is not available in the literature.

2.2 The jet in cross-flow approach

From the fluid mechanics point of view, the flame gas sampling by horizontal dilution tubes represents a jet in cross-flow arrangement, wherein the gas sample flow is the laminar jet and the cold dilution gas flow is the transverse flow. Due to their wide range of practical applications, jets in cross-flow have been investigated intensively for many years (Reference Chang, Shao, Hu and ZhangChang, Shao, Hu, & Zhang, 2016; Reference MaheshMahesh, 2013; Reference Muppidi and MaheshMuppidi & Mahesh, 2005). The flow characteristics are described by several characteristic numbers:

• • the jet to cross-flow velocity ratio $R = {({\rho _j}/{\rho _\infty })^{0.5}} \times {u_j}/{u_\infty }$;

• • the jet to cross-flow density ratio $S = {\rho _j}/{\rho _\infty }$;

• • the jet to cross-flow momentum flux ratio J = R ²;

• • the jet Reynolds number $R{e_j} = {u_j}{D_j}/{\mu _j},{\mu _j} = {\eta _j}/{\rho _j}$;

where u is velocity, ρ is density, η is dynamic viscosity, μ is kinematic viscosity, D is diameter and the subscripts are j for the jet and ∞ for the cross-flow. Different from flame gas sampling, experimental investigations usually employ long jet pipes to ensure fully developed flow. High R, J and Re_j were found to favour the mixing of the jet into the cross-flow, but even in case of laminar jets (low Re_j and R), turbulent vortex structures can develop which induce a high mixing efficiency. Yet, extremely small jet-to-cross-flow momentum flux ratios, J < 5, may prevent the jet from penetrating deeply into the cross-flow and may cause the jet to adhere to the wall or wall boundary layer (Reference KaragozianKaragozian, 2014; Reference MaheshMahesh, 2013). As shown earlier (Reference Mätzing and StapfMätzing & Stapf, 2019), such a situation may arise under flame gas sampling conditions, when a high-temperature and high-viscosity jet enters into a much colder cross-flow with R = 0.14–0.77 and J = 0.02–0.6 (these values hold in the present study, too).

3.1 CFD simulations

The simulations were performed for free-standing and embedded dilution tubes, as outlined in figure 1(a,b). The free-standing dilution tube served as a base case. The experimental set-up was used in our laboratory previously (Reference Frenzel, Krause and TrimisFrenzel, Krause, & Trimis, 2017) and it is very similar to other common set-ups (Reference Zhao, Yang, Johnston, Wang, Wexler, Balthasar and KraftZhao, Yang, Johnston, et al., 2003). The uncooled, free-standing ceramic (Al₂O₃) dilution tube was 9 mm inner diameter, with a wall thickness of 0.5 mm; the pinhole diameter was D_orf = 0.3 mm. The input nitrogen dilution flow was 30 s.l.p.m. $(1\ {\rm s.l.p.m.} = 1.66 \times 10^{-5}\ {\rm Nm}^{3}\ {\rm sec}^{-1})$, which corresponds to a Reynolds number Re = 5000. At the measurement locations ΔP1 and ΔP2 in figure 1(c), low differential pressure sensors were installed and the pressure was assumed to drop linearly along the tube axis. This assumption was confirmed during the simulations. Note that in this set-up, the outlets to the measurement device (DMA) and to the pump are in parallel, at slightly different vertical positions. Other set-ups may have different arrangements, e.g. outlets at equal vertical positions and also outlets which are oriented perpendicular to each other. Such a configuration was simulated earlier (Reference Mätzing and StapfMätzing & Stapf, 2019) with no difference in major findings about the flow properties. Hence, the current results are taken to be valid for most of the common geometric configurations and tube materials.

Figure 1.

Sketches of free-standing (a) and embedded (b) dilution tubes used in flame studies. Geometry of the free-standing dilution tube and wall boundary conditions (c). Tube wall temperature is 298 K under cold (calibration) conditions and 400 K for embedded tubes. For free-standing tubes, under hot (flame) conditions, the tube wall temperature profile is indicated above. Initial turbulence intensity of 5 % at N₂ inlet (see § 4.1).

The solid tube wall is not modelled explicitly; therefore the orifice is represented by a small tube which is 0.5 mm long, i.e. the wall thickness. Occasional clogging of the orifice can pose experimental problems, as will be discussed below, but it is not the object of this study.

To simulate embedded dilution tubes, the same geometric configuration is used, but the wall temperature is set to a constant appropriate value, e.g. 400 K, which results from the intense cooling by the stagnation plate (Reference Saggese, Cuoci, Frassoldati, Ferrario, Camacho, Wang and FaravelliSaggese et al., 2016).

The simulations were performed for a single phase gas flow. Under sooting conditions, in particular, when the soot forming (inception) zone is studied, the soot volume fraction is typically below 10⁻⁷ and, the particle size is well below 100 nm, often below 10 nm. Thermophoretical effects can eventually accelerate the particle motion relative to the gas in the cooling region in front of the orifice (Reference Saggese, Cuoci, Frassoldati, Ferrario, Camacho, Wang and FaravelliSaggese et al., 2016). In and behind the orifice, the nanoparticle trajectories can be expected to follow the gas flow and, particle–flow interactions are negligible due to the very small volume fraction (Reference Corson, Mulholland and ZachariahCorson, Mulholland, & Zachariah, 2018; Reference FriedlanderFriedlander, 2000). Hence, as first approximation, two-phase flow modelling is not really mandatory to capture the important flow properties and processes in the dilution tube.

The CFD simulation was performed in three dimensions using Ansys Fluent. Two grids of different cell number and cell geometry were generated:

1. (i) an unstructured grid with 18.5 × 10⁶ tetrahedral cells, using Fluent 17.2 for RANS and LES;

2. (ii) an unstructured grid with 3.2 × 10⁶ polyhedral cells, using Fluent 14.5 for RANS only.

The calculations were started with a precursor RANS calculation on grid (ii), in order to estimate the integral length scale ${l_0} = {k^{0.5}}/({C_\mu }\ \omega )$, wherein k is the turbulent kinetic energy, ω is the specific dissipation rate and C_μ = 0.09 (Reference GerasimovGerasimov, 2016). Since approximately 5–10 cells across l ₀ (∼270 μm) are required to resolve ∼80 % of the turbulent energy spectrum, grid (i) was generated for the LES calculations with roughly six times smaller cell volume. The two grids were also used for a study of grid independence.

The general form of the conservation equations, which are solved in CFD calculations, is

(2)\begin{gather}\frac{{\partial \rho }}{{\partial t}} + \boldsymbol{\nabla }\boldsymbol{\cdot }(\rho \boldsymbol{u}) = {S_m}\quad \textrm{mass}\;\textrm{conservation}\;\textrm{or}\;\textrm{continuity},\end{gather}

(3)\begin{gather}\frac{\partial }{{\partial t}}(\rho \boldsymbol{u}) + \boldsymbol{\nabla }\boldsymbol{\cdot }(\rho \boldsymbol{uu}) = {-} \boldsymbol{\nabla }P + \boldsymbol{\nabla }\boldsymbol{\cdot }({\overline{\overline \tau } } ) + {\boldsymbol{F}_{\boldsymbol{b}}}\quad \textrm{momentum},\end{gather}

(4)\begin{gather}\rho \frac{{\textrm{D}E}}{{\textrm{D}t}} = {-} \boldsymbol{\nabla }(P\boldsymbol{u}) + \boldsymbol{\nabla }\boldsymbol{\cdot }({\boldsymbol{u}}\overline{\overline \tau }) + \boldsymbol{\nabla }\boldsymbol{\cdot }({\lambda _G}\boldsymbol{\nabla }T) + {S_E}\quad \textrm{energy},\end{gather}

(5)\begin{gather}\frac{{\partial \rho {y_i}}}{{\partial t}} + \boldsymbol{\nabla }\boldsymbol{\cdot }(\rho \boldsymbol{u}{y_i} + {\boldsymbol{J}_{\boldsymbol{i}}}) = {R_i} \quad \textrm{species}\ \textrm{transport,}\end{gather}

wherein t (s) is time, P (Pa) is pressure, τ (Pa) is stress tensor, E (J) is energy, λ_G (W (m K)⁻¹) is heat conductivity of the gas, T (K) is temperature, y_i (−) is the mass fraction of species i, J_i (kg (m² s)⁻¹) is the diffusion flux of species i and R_i (kg (m³ s)⁻¹) is a source term of species i. The source term for the mass, S_m, is zero if no condensation or evaporation processes are considered. Similarly, R_i may be zero in the absence of chemical reactions. Here F_b is the source term for momentum which can be chosen to include gravity. The source term for energy, S_E, includes heat transfer processes as well as reaction heats and heat exchange due to phase transitions, if applicable.

In the RANS calculations, the turbulence model SST k–ω was used with corrections for curvature and low Reynolds numbers Re. Near wall treatment was done by high resolution inflation layers. The spatial resolution at the walls was better than 170 μm in mesh (ii) and better than 30 μm in mesh (i). The dimensionless wall distance y ⁺ was always less than unity in the laminar sublayer for both meshes which means that transport properties at the wall boundary layer were fully captured (Reference Zirwes, Häber, Zhang, Kosaka, Dreizler, Steinhausen and BockhornZirwes et al., 2021).

The pressure based solver was used to solve the equation system with the SIMPLE algorithm for pressure–velocity coupling which is acceptable for transient calculations, too, according to the Ansys Fluent User's Guide. Second-order upwind schemes were used for discretization. For good convergence, the calculations were often continued until residuals dropped below 10⁻⁷.

The LES calculations were performed using mesh (i) according to the ANSYS guidelines (Reference GerasimovGerasimov, 2016). The spatial resolution was better than one fifth of the turbulence length scale which allows us to resolve roughly 80 % of the turbulent kinetic energy. For subgrid-scale modelling, the dynamic Smagorinsky model was used with wall Prandtl number set to 0.85. The time step was set to 15 μs, i.e. the average travel time across a cell.

During post-processing, the Q criterion, the λ ₂ criterion as well as pressure and velocity fluctuations were used to detect vortices. In brief, the Q criterion defines a vortex in a region where the square of vorticity is larger than the square of strain rate meaning that rotational deformation is larger than stretching deformation. The λ ₂ criterion determines vortex regions from the eigenvalues of the symmetric and antisymmetric parts of the velocity tensor: if two of the eigenvalues are negative, a local pressure minimum exists which induces vortex formation. None of these methods is unique for vortex identification, as has been discussed and reviewed by several authors (Reference FröhlichFröhlich, 2006; Reference HolménHolmén, 2012; Reference Jiang, Machiraju and ThomsonJiang, Machiraju, & Thomson, 2005; Reference KolárKolár, 2011). For example the Q criterion may not distinguish between shear and rotation induced vortices, the λ ₂ criterion may not resolve closely neighboured vortices and some authors recommend scalar fluctuation monitoring to avoid problems related to numerical noise. But, as will be shown later, the methods used here gave very similar results showing little vortex formation in the dilution tube except close to and in the outlet pipes.

The mass or mole fractions of CO₂ were used to monitor and assess the mixing of the flame gas into the dilution gas.

3.2 Gas properties and flow rate evaluation

The calculations were performed for the following cold and hot conditions:

1. (i) ambient air conditions to simulate the experimental calibration procedure; air was taken to be dry and free of trace gases (mole fractions ${x_{{\textrm{O}_2}}} = 0.21$, ${x_{{\textrm{N}_2}}} = 0.79$);

2. (ii) flame conditions were simulated without the flame gas flowing around the dilution tube, just by defining the temperature to be 1000 K at the orifice (Reference Zhao, Yang, Wang, Johnston and WangZhao, Yang, Wang, et al., 2003) and using the tube wall temperature profile shown in figure 1, assuming a CH₄/O₂ oxy-fuel flame, equivalence ratio Φ = 2.4, C/O = 0.6, the flame gas concentrations were estimated from Reference Köhler, Oßwald, Xu, Kathrotia, Hasse and RiedelKöhler et al. (2016) and from Reference Stelzner, Weis, Habisreuther, Zarzalis and TrimisStelzner, Weis, Habisreuther, Zarzalis, & Trimis (2017). They are listed in table 1. The selected gas mixture has a relatively high viscosity η, $\eta _{flame,1000\;\textrm{K}} = 3.64 \times {10^{ - 5}}\;\textrm{kg}\;{\textrm{(m}\ \textrm{s})^{ - 1}} = 364\;\mathrm{\mu }\textrm{P}$, approximately 10 % higher than in usual fuel–air flames (see table 2). Its molar weight and density are half that of air.

Table 1.

Flame gas composition, x = mole fraction.

Table 2.

Gas viscosities η of air and flame gas between 298 K and 1000 K.

The orifice volume flows Q_orf (s.l.p.m.) were determined from the calculated mass flows (kg s⁻¹) and densities ρ (kg m⁻³) in the following way:

(6)\begin{equation}{Q_{orf}}(\textrm{s.l.p.m.)} = \frac{{\dot{m}(\textrm{kg}\;\textrm{se}{\textrm{c}^{ - 1}})}}{{\rho (\textrm{kg}\;{\textrm{m}^{ - 3}})}} \times 1000(\textrm{l}\;{\textrm{m}^{-3}}) \times 60(\textrm{sec}\; {\textrm{min}^{-1}}) \times \frac{{273.15\;\textrm{K}}}{{{T_{orf}}(\textrm{K})}} \times \frac{{{P_{orf}}(\textrm{Pa})}}{{101\;325\;\textrm{Pa}}},\end{equation}

where Q_orf, T_orf and P_orf all refer to the inlet (‘flame side’) conditions at the orifice entrance. The DR was determined from the dilution gas flow Q _{N2, 273.15 K} and Q_orf according to (7), and also from the calculated mole fractions of CO₂ at the orifice entrance and at the DMA outlet. Equation (7) reads

(7)\begin{equation}DR = 1 + \frac{{{Q_{{\textrm{N}_2},273.15\;\textrm{K }}}(\textrm{s.l.p.m.})}}{{{Q_{orf}}(\textrm{s.l.p.m.})}}.\end{equation}