Studies on Structures of Turbulent Swirling Coal Flames by Large-Eddy Simulation

1. Introduction

Large-eddy simulation (LES) of combustion, including pulverized-coal combustion, is a powerful tool to study the instantaneous flame structures and can give more accurate statistical results than Reynolds-averaged Navier-Stokes (RANS) modeling,. For LES of coal combustion, many researchers did LES of coal jet flames. Kurose et al. [1] did LES of a pulverized-coal-air jet flame using a Smagorinsky SGS stress model and a conservative scalar-fast chemistry gas combustion model, together with coal-particle devolatilization and char combustion models and Lagrangian particle tracking approach. The LES obtained instantaneous particle dispersion and instantaneous temperature map were obtained. The LES statistically averaged gas velocity, temperature, species concentration and particle number density were also reported. However, no vortex structure was shown, and it is not clear about the particle-turbulence interaction in the coal flame. Also, no experimental validation of the velocity, temperature and species concentration fields was reported. Zhou, H [2] reported the so-called “DEM-LES” of coal combustion in a bubbling fluidized bed. For gas-particle flows the SGS kinetic energy stress model together with a DEM particle model is used. The coal devolatilization and char combustion models are taken into account. It seems that no turbulent combustion model was taken into account. As the modeling results, the particle temperature and heating rate and gas species concentration were reported. No instantaneous flow and flame structures were given and there were no experimental validation. Besides, just as in his case of “DEM-LES” of isothermal gas-particle flows, the grid size is too coarse, it looks like a U-RANS modeling, but not LES. Yamamoto et al. [3] studied the coal jet flame by LES using the code Open-FOAM-1.3. The Smagorinsky SGS stress model, DT radiative heat transfer model, Lagrangian particle tracking approach, and a very simple turbulent combustion model based on the concept of the equilibrium concentration and Kolmogorov time scale are used. The LES predicted time-averaged gas temperature and coal burn-out rate were compared with the measurement results. The agreement is good. However, no radial profiles are given. The flame lift-off height predicted by LES is in much better agreement with experiments than that predicted by RANS modeling using the k-ε model. The temperature maps are shown, where the time-averaged temperature maps obtained using LES and RANS are different from each other. However, there is no quantitative comparison between the LES and RANS modeling temperature distribution. Also, no flow vortex structures are reported and no experimental validation of the predicted species concentration field is given. Muto et al.[4] used LES to study the effect of oxygen concentration on NO formation in a coal jet flame. Smith Philip et al. [5] did LES of a coal jet flame using an Eulerian particle-phase modeling, called < Direct Quadrature Method of Moments (DQMOM). The results show that Eulerian-Eulerian LES can capture the detailed flame structures. Wan et al. [6] studied pilot-assisted pulverized-coal flame jet by LES, focusing on the effect of devolatilizatio model. Wen et al. [7] did LES of a coal jet flame adopting a velocity- scalar joint filtered density function model. The statistical results were validated by experiments.Alternative, an extended flamelet/progress variable (EFPV) model for LES of a coal jet flame was reported by Wen et al. [8] and Xing et al. [9].The results show that it is better than the EBU combustion model.

As regard to LES of swirling coal flame, Edge et al. [10] studied two kinds of IFRF oxygen-coal burners using LES and RANS modeling with k-ε and k-ε-RNG turbulence models. The commercial code ANSYS FLUENT version 12 was used. In LES the Smagorinsky-Lilly SGS stress model, an EDC combustion model, and coal devolatilization and char oxidation models were selected. 1.6 million grid nodes with grid sizes of 1–2 mm were taken for the computation domain of near-burner region. The velocity vectors superposed with temperature distribution obtained by LES and RANS modeling were obtained. The authors like to compare and judge the LES and RANS modeling results, but from the predicted temperature (no experimental validation) and wall heat fluxes it is difficult to say that the LES results are better than the RANS modeling results. Gharebaghi et al. [11] did LES of oxygen-coal combustion in a furnace. The commercial code ANSYS-FLUENT was used by selecting the WALE SGS stress model, the EBU (EDM) volatile and CO combustion model, the DO radiative heat transfer model, together with Lagrangian particle tracking approach, coal-particle devolatilization and char oxidation models. The results give the temperature maps in the near-burner region by both LES and RANS modeling using the k-ε model. The comparison of predicted temperature change along the vertical direction from the burner centreline using both LES and RANS modeling with the measurement results was given. Both modeling results remarkably over-predict the temperature and LES results are not better than RANS results. The author’s opinion is that although LES can more accurately model the turbulence, but for this temperature over-prediction it may be caused by the inaccuracy in radiative heat transfer model and char combustion model. No experimental validation of the predicted gas velocity and species concentration was reported. Also no flow vortex structures are reported. Rabacal et al.[12] studied the coal particle history during combustion in a large-scale furnace using large eddy simulation. A specific feature of the devolatilization model was discussed. Watanabe et al. [13] did large-eddy simulation of coal combustion in a pulverized coal combustion furnace with a complex burner. The simulated gas temperature and NO concentrations are in general agreement with the experiments. Rieth et al [14] introduced the flamelet model into large eddy simulation of realistic coal furnaces. The results of the LES are focussied on instantaneous particle and gas phase data to gain additional insight into the coal conversion process inside the furnace. Shen et al. [15] studied coal combustion in a 660 MW utility boiler furnace using LES. However, the adopted grid sizes for such a large-size furnace are questionable for the requirement of LES. Olenik et al. [16] studied swirl-stabilised pulverised coal combustion in the semi-industrial (2.5 MW) IFRF furnace by means of large eddy simulation with a EBU combustion model and compared to results from the experimental results and RANS predictions. It was found that some temperature underprediction was observed. Franchetti et al. [17] studied a 100 kW(th) swirling oxy-coal furnace using LES. It was reported that in the downstream region the LES overestimated the combustion rates.

It is seen that in most of above-cited studies the focus was paid to the time-averaged LES results and the application of LES in practical coal furnaces. The instantaneous structures of swirling coal flames remain to be further studied. In this paper, swirling coal flames are studied by LES. The specific features of instantaneous strucutures of swirling coal flames are analyzed and discussed.

2. Filtered Governing Equations for LES of Coal Combustion

The filtered gas governing equations for Eulerian-Lagrangian LES of coal combustion can be given as:

$$\frac{\partial \rho }{\partial t}+\frac{\partial }{\partial x_i}\left(\rho {\overline{u}}_i\right)=\overline{S}\begin{array}{cccc}& & & \end{array}\overline{S}={\displaystyle \sum {\overline{n}}_k{\dot{m}}_k}$$

(1)

$$\frac{\partial }{\partial t}\left(\rho {\overline{u}}_i\right)+\frac{\partial }{\partial x_j}\left(\rho {\overline{u}}_i{\overline{u}}_j\right)=\frac{\partial }{\partial x_j}\left(\mu \frac{\partial {\overline{u}}_i}{\partial x_j}\right)-\frac{\partial \overline{p}}{\partial x_i}-\frac{\partial {\tau }_{ij}}{\partial x_j}+\frac{{\rho }_k}{{\tau }_{rk}}({\overline{u}}_{ki}-{\overline{u}}_i)+{\overline{u}}_{i}S$$

(2)

$$\frac{\partial }{\partial t}(\rho \overline{f})+\frac{\partial }{\partial x_j}(\rho {\overline{u}}_j\overline{f})=\frac{\partial }{\partial x_j}(\frac{\mu }{Sc}\frac{\partial \overline{f}}{\partial x_j})-\frac{\partial m_{jsgs}}{\partial x_j}+{\alpha }_{s}S$$

(3)

$$\frac{\partial \left(\rho \overline{h}\right)}{\partial t}+\frac{\partial }{\partial x_j}\left(\rho \overline{h}{\overline{u}}_j\right)=\frac{\partial }{\partial x_j}\left(\frac{\mu }{\mathrm{Pr}}\frac{\partial \overline{h}}{\partial x_j}\right)-\frac{\partial q_{sgs,j}}{\partial x_j}+n_{k}Q{}_k+c_p\overline{T}S$$

(4)

$$\frac{\partial \rho {\overline{Y}}_s}{\partial t}+\frac{\partial }{\partial x_j}\left(\rho {\overline{u}}_j{\overline{Y}}_s\right)=\frac{\partial }{\partial x_j}\left(\frac{\mu }{Sc}\frac{\partial {\overline{Y}}_s}{\partial x_j}\right)-{\overline{w}}_s-w_{sgs}-\frac{\partial g_{jsgs}}{\partial x_j}+{\overline{\alpha }}_s\overline{S}$$

(5)

3. Sub-Grid Scale (SGS) Stress Models

For single-phase flows, frequently used SGS stress models are the Smagorinsky eddy-viscosity model [18], Germano´s dynamic eddy-viscosity model [19] and Kim´s SGS energy-equation model [20]. The SGS stress in the momentum equation is defined by:

$${\tau }_{ij}\equiv \rho \overline{u_i{u}_j}-\rho {\overline{u}}_i{\overline{u}}_j$$

The Smagorinsky eddy viscosity model is given by

where C_s=0.16 is the Smagorinsky constant, Δ is the filtered scale. The shortcoming of the Smagorinsky model is its too large dissipation rate of the kinetic energy. Germano proposed a dynamic Smagorinsky eddy viscosity model, in which the coefficient C_s is not a constant, but is determined by the following filtration process as

$$M_{ij}=2\rho C_s{\overline{\Delta }}^2\left|\widehat{\overline{S}}\right|{\widehat{\overline{S}}}_{ij}-2\rho C_s{\widehat{\overline{\Delta }}}^2\left|\widehat{\overline{S}}\right|{\widehat{\overline{S}}}_{ij}{L}_{ij}=2\rho C_s{\widehat{\overline{\Delta }}}^2\left|\widehat{\overline{S}}\right|{\widehat{\overline{S}}}_{ij}-2\rho C_s{\overline{\Delta }}^2\left|\widehat{\overline{S}}\right|{\widehat{\overline{S}}}_{ij}$$

It was indicated by many investigators that in many cases the dynamic eddy-viscosity model gives better results than that given by the standard Smagorinsky model. Kim proposed a SGS kinetic energy model as

$${\tau }_{ij}=-2\rho v_t(\stackrel{-}{S_{ij}}-\frac{1}{3}\stackrel{-}{S_{kk}}{\delta }_{ij})+\frac{2}{3}\rho k^{sgs}{\delta }_{ij}$$

$$\frac{\partial \overline{\rho }{k}^{sgs}}{\partial t}+\frac{\partial }{\partial x_i}\left(\overline{\rho }{\stackrel{-}{{\displaystyle u}}}_i{k}^{sgs}\right)=P^{sgs}-D^{sgs}+\frac{\partial }{\partial x_i}\left(\frac{\overline{\rho }{\nu }_t}{{\mathrm{Pr}}_t}\frac{\partial k^{sgs}}{\partial x_i}\right)+\stackrel{.}{W_s}$$

$$k^{sgs}=\frac{1}{2}(\overline{u_i{u}_j}-{\overline{u}}_i{\overline{u}}_j)v_t=C_v{(k^{sgs})}^{1/2}\stackrel{-}{\Delta }$$

The merit of the SGS kinetic energy model is that it avoids the negative value of eddy viscosity; hence it is better than the eddy-viscosity models. In application it was reported that the results of SGS kinetic energy model are near to those obtained using the dynamic eddy-viscosity model. In this paper, a combined

Smagorinsky-Lilly’s SGS stress model [21] is used as

$$\begin{array}{c}{\tau }_{ij}-\frac{1}{3}{\tau }_{kk}{\delta }_{ij}=-2{\mu }_t{\overline{S}}_{ij},{\overline{S}}_{ij}\equiv \frac{1}{2}\left(\frac{\partial {\overline{u}}_i}{\partial x_j}+\frac{\partial {\overline{u}}_j}{\partial x_i}\right),{\mu }_t=\rho L_s^2\left|\overline{S}\right|,\left|\overline{S}\right|\equiv \sqrt{2{\overline{S}}_{ij}{\overline{S}}_{ij}}，\\ L_s=\min (\kappa d,C_s{V}^{1/3})。\end{array}$$

(6)

where：${\delta }_{ij}$-Kronecker Delta，V -cell volume，$L_s$-SGS scale，$S_{ij}$-stress deformation rate，${\mu }_t$-SGS viscosity coefficient， d -closest distance to the walls，$\kappa$-Von Kármán constant，$C_s$-empirical constant, taken as $C_s$＝0.1。

The SGS mass flux and heat flux are closed by gradient modeling as

$$g_{sgs,j}=\rho \left(\overline{u_j{Y}_s}-{\overline{u}}_j{\overline{Y}}_s\right)=\frac{{\mu }_t}{{\sigma }_Y}\frac{\partial {\overline{Y}}_s}{\partial x_j}$$

(7)

$$q_{sgs,j}=\rho \left(\overline{u_{j}T}-{\overline{u}}_j\overline{T}\right)=\frac{{\mu }_t}{{\sigma }_T}\frac{\partial \overline{T}}{\partial x_j}$$

(8)

where ${\sigma }_T$ and ${\sigma }_Y$ are model constants

4. Combustion Models for LES of Coal Combustion

For gas combustion (CO and volatile combustion) during coal combustion, two models were adopted. One is the Magnusen’s modified eddy-break-up (EBU) model [22]

$${\overline{w}}_s{w}_{sgs}=c_1\rho \left|\overline{S}\right|\min \left\{\overline{Y_f},\frac{\overline{Y_{ox}}}{\beta },c_2\frac{\overline{Y_P}}{1+\beta }\right\}$$

(9)

where $c_1$=4.0，$c_2$=0.5

The gas-phase combustion includes two reactions:

$$C{H}_4+2O_2\to C{O}_2+2H_{2}O\left(\mathrm{reaction}1\right)$$

(10)

$$2CO+O_2\to 2C{O}_2\left(\mathrm{reaction}2\right)$$

(11)

Alternatively, a presumed-PDF fast chemistry gas combustion model [23] is used. In the latter case a mixture fraction equation should be solved

$$\frac{\partial }{\partial t}(\rho \overline{f})+\frac{\partial }{\partial x_j}(\rho {\overline{u}}_j\overline{f})=\frac{\partial }{\partial x_j}(\frac{\mu }{Sc}\frac{\partial \overline{f}}{\partial x_j})-\frac{\partial m_{jsgs}}{\partial x_j}$$

(12)

where $m_{jsgs}=\rho \left(\overline{u_{j}f}-{\overline{u}}_j\overline{f}\right)=\frac{{\mu }_t}{{\sigma }_f}\frac{\partial \overline{f}}{\partial x_j}$ and $\overline{ff}-\overline{f}\overline{f}$$=\frac{1}{2}{L}_s{}^2{\left|\nabla \overline{f}\right|}^2$

For radiative heat transfer the P1 model is used as：

$$-{\overline{q}}_r=aG-4a\sigma {\overline{T}}^4$$

(13)

For particle combustion, neglecting moisture evaporation and taking the sub-models of two-equation devolatilization and diffusion-kinetic char oxidation model with total formation of CO are used [22].

5. Particle-Phase Equations

The Lagrangian particle continuity, momentum and energy equations are [23]

$${\displaystyle {\int }_A{n}_k}{v}_{kn}dA=N_k=const$$

(14)

$$\frac{d{u}_{ki}}{d{t}_k}=\frac{1}{{\tau }_{rk}}(u_i-u_{ki})$$

(15)

$$\frac{d{T}_k}{d{t}_k}=[Q_h-Q_k-Q_{rk}+{\dot{m}}_k(c_{p}T-c_k{T}_k)]/(m_k{c}_k)$$

(16)

6. Simulated Pulverized-Coal Swirl Combustor and Numerical Method

The LES of swirling pulverized-coal combustion was made for a combustor, shown in Figure 1.

Figure 1. A coal combustors with swirl burner.

Figure 1. A coal combustors with swirl burner.

The swirl combustor was measured at Shanghai Jiaotong University (Figure 1a) [24]. In this combustor, non-swirling primary air laden with coal particles was issued from the inner annular inlet. The secondary swirl air was issued from the outer annular inlet. The flow parameters of Cases 1 and 2 are given in Table 1. The flow parameters of Case 3 are the same as those of Case 2 with only a difference in adding a inlet perturbation using a vortex method.

Table 1. Flow Parameters for Cases 1 and 2.

	Primary air flow laden with coal particles			Secondary swirling air flow
	Axial velocity（m/s）	Temperature（K）	Coal feeding rate（kg/h）	Axial velocity（m/s）	Temperature（K）	Swirl number
Case 1	23.1	353	14.0	38.3	573	1.0
Case 2	28	353	16.4	45.1	573	1.03

Table 1. Flow Parameters for Cases 1 and 2.

	Primary air flow laden with coal particles			Secondary swirling air flow
	Axial velocity（m/s）	Temperature（K）	Coal feeding rate（kg/h）	Axial velocity（m/s）	Temperature（K）	Swirl number
Case 1	23.1	353	14.0	38.3	573	1.0
Case 2	28	353	16.4	45.1	573	1.03

The proximate and element analyses of coal are given in Table 2 and Table 3

Table 2. Coal proximate analysis (Weight％).

moisture	volatile	fixed carbon	ash
6.3	35.8	53.7	4.2

Table 2. Coal proximate analysis (Weight％).

moisture	volatile	fixed carbon	ash
6.3	35.8	53.7	4.2

Table 3. Coal element analysis（Weight％).

carbon	hydrogen	nitrogen	sulfur	oxygen
72.6	5.05	1.29	1.55	15.31

Table 3. Coal element analysis（Weight％).

carbon	hydrogen	nitrogen	sulfur	oxygen
72.6	5.05	1.29	1.55	15.31

The particle size distribution are shown in Table 4

Table 4. Particle size distribution.

<200$\mu m$％	<75$\mu m$％	<40$\mu m$％	<25$\mu m$％	<10$\mu m$％
100	80	60	40	15

Table 4. Particle size distribution.

<200$\mu m$％	<75$\mu m$％	<40$\mu m$％	<25$\mu m$％	<10$\mu m$％
100	80	60	40	15

Kinetic constants of coal devolatilization:Bv1=2e+5, Ev1=1.046e+8j/kgmol，Bv2=1.3e+7，Ev2=1.674e+8 j/kg-mol。Kinetic constants of char reaction: Bh=0.86kg/m².s (N/m²)，Eh=1.02$\times$10⁸J/kg-mol。Gross Calorific Value:29.29MJ/kg。

The grid sizes are taken as 0.1mm to 1mm. The total grid number is 1943214. The contribution of the resolved part of the total kinetic energy is $\ge$75% in most of the computational region. The effect of grid size has been tested and it is found that the statistical results are independent of the grid sizes smaller than the adopted ones. The time step is taken as 1e-3s, which is chosen based on the computation test. Within each time step the convergence can be reached after 30 to 70 iterations. The filtered governing equations are discretized by a finite-volume scheme. For the numerical procedure, the pressure-implicit with splitting operators (PISO) algorithm is used for p-v corrections; the second-order implicit difference scheme is used for the time-dependent term, the convection and diffusion terms. The stochastic component of the flow field at the inlet boundary is taken into account by superposing random perturbations on each velocity component. The boundary condition at the exit is based on pressure outlet boundary conditions. At the near-wall grid nodes the wall function approximation is employed..

7. Simulation Results and Discussion

The simulated time-averaged results and their assessment by experiments can be found in Ref. [25]. The following are simulated instantaneous results for structures of swirling coal flames, which has not reported before. Figure 2, Figure 3 and Figure 4 show the instantaneous vorticity isolines of swirling coal flames, The large-eddy coherent structures can be seen for all cases. These structures are produced near the inlet and rapidly damped in the downstream region. It is seen that the vortex structures of Cases 2 and 3 are stronger than those of Case 1, indicating that stronger turbulence exist in Cases 2 and 3.

Figure 2. Instantaneous vorticity isolines (1/s)（Case1）.

Figure 2. Instantaneous vorticity isolines (1/s)（Case1）.

Figure 3. Instantaneous vorticity isolines (1/s)（Case 2）

Figure 3. Instantaneous vorticity isolines (1/s)（Case 2）

Figure 4. Instantaneous vorticity isolines (1/s)（Case 3）

Figure 4. Instantaneous vorticity isolines (1/s)（Case 3）

Figure 5, Figure 6 and Figure 7 give the instantaneous velocity vectors of Cases 1 to 3. It can also be seen that the vortex structures of Cases 2 and 3 are stronger than that of Case 1.

Figure 5. Instantaneous velocity vectors（Case 1）

Figure 5. Instantaneous velocity vectors（Case 1）

Figure 6. Instantaneous velocity vectors（Case 2）

Figure 6. Instantaneous velocity vectors（Case 2）

Figure 7. Instantaneous velocity vectors（Case 3）

Figure 7. Instantaneous velocity vectors（Case 3）

Figure 8, Figure 9 and Figure 10 show instantaneous particle concentration isolines. In Cases 1 to 3, the particles are concentrated in the near-axis zones. In other regions the particle concentration is low. The particle distribution of Case 2 is more uniform than that of Case 1, and that of Case 3 is more uniform than that of Case 2. This is because of different inlet velocities and coal feeding rates and adding inlet purturbation in Case 3. The near-wall particle concentrtion is higher for Case 3 than that of other cases.

Figure 8. Instantaneous particle concentration isolines（Case1）

Figure 8. Instantaneous particle concentration isolines（Case1）

Figure 9. Instantaneous particle concentration isolines（Case 2）

Figure 9. Instantaneous particle concentration isolines（Case 2）

Figure 10. Instantaneous particle concentration isolines（Case 3）

Figure 10. Instantaneous particle concentration isolines（Case 3）

Figure 11, Figure 12 and Figure 13 give the instantaneous temperature isolines. The flames concentrated in the near-inlet highly shear region and then spreaded in the downstream region. The flame length of Case 1 is longer than that of Case 2 and Case 3 due to different inlet velocities and better mixing of coal particles in Cases 2 and 3.

Figure 11. Instantaneous temperature isolines (K）（Case 1）

Figure 11. Instantaneous temperature isolines (K）（Case 1）

Figure 12. Instantaneous temperature isolines (K）（Case 2）

Figure 12. Instantaneous temperature isolines (K）（Case 2）

Figure 13. Instantaneous temperature isolines (K）(Case 3）

Figure 13. Instantaneous temperature isolines (K）(Case 3）

Conclusions

(1) The structures of swirling coal flames were studied in detail by LES

(2) For all cases coherent vortex structures exist in flames, indicated by vorticity isolines and velocity vectors.

(3) For all cases flames are located in the regions of high particle concentration and high Vorticities.

(4) The annular inlet of primary air flow laden with coal particles leads to higher vorticity, stronger mixing and higher temperature near the inlet, hence intensifies combustion than the central inlet.

(5) Larger secondary swirling air flow enhances mixing and intensifies combustion.