*2.2. CFD Modeling*

The mesh was constructed for the boiler using 3,529,358 cells. This mesh was selected after assessing the sensitivity by comparing the degree of numerical di ffusion between coarser (1.23 million cells) and finer (5.47 million cells) meshes. In these meshes, the cell fineness was varied only in the burner zone to have a respective average volume of 0.0114, 0.0034, and 0.0022 m<sup>3</sup>/cell. Compared to the finer mesh, the selected one exhibited very small deviations in the key performance parameters (0.12% in carbon conversion, 1.9% in the exit NO concentration, and 0.2% in the heat absorption on the furnace wall). In the detailed comparison for the flow pattern along the vertical centerline, the deviations from the results of the finer mesh were 4.7% in the velocity profile and 0.71% in the temperature profile. By contrast, the coarser mesh had deviations of 17.58% and 1.11%, respectively. The computation time for the selected mesh to reach a converged solution was 69% of that for the finer mesh. The details of the mesh sensitivity test together with the iteration strategy have been reported elsewhere [19].

The CFD simulations were performed using ANSYS Fluent (version 17.2) [20] with submodels for reactions, turbulence, and radiation, most of which have been reviewed by Sankar et al. [5] as common models applied to a coal-fired boiler. The reaction submodels are summarized in Table 3. Coal particles were tracked using the discrete phase method for 10 size fractions ranging between 5.9 μm and 204 μm, with a total of 49,560 particles. FLASHCHAIN [21] was used to determine the input parameters for coal devolatilization in a drop tube furnace at 1200 ◦C which was close to the heat transfer condition of pulverized particles in the boiler. This code predicted the reaction dynamics based on the semi-empirical coal network mode to provide the product yields (tar, CO, CO2, H2O, CH4, C2H4, C2H6, C3H6, HCN, H2S, and solid char) and reaction kinetics as listed in the table. Because the heating rate was much faster than that in the proximate analysis, the total volatiles yield (58.56% daf) was larger than the volatile matter content (41.67% daf) in Table 1. Because of sharing the same mass source from devolatilization and having similar reaction rates, C2+ hydrocarbons were

simplified to an imaginary species of CxHy. The yield and composition of tar were modified to include the H and O fractions of char so that the char could be modelled as a pure carbon solid and ash. For simplicity, it also incorporated the two minor volatiles (HCN and H2S). Char conversion by O2, CO2, and H2O was solved using the unreacted core shrinking model (UCSM) [22] which is suitable for high-temperature reactions in the industrial-scale furnace and the release of the UBC in fly ash. This model considers the three competing rates of the heterogeneous reaction on the char core surface and the gas di ffusions onto the particle and through the ash layer. With the UCSM, the decrease in the reactivity toward the end of char conversion can be simulated, which was required to predict the UBC in fly ash. The multiple volatile products and the UCSM were implemented into the CFD code using user defined functions (UDF).


**Table 3.** Summary of submodels adopted for boiler simulation.

Gaseous reactions (R4)−(R10) were based on the global reaction scheme of Jones and Lindstedt [23] for hydrocarbon and the tar oxidation rate for (R4) was taken from Smoot and Smith [24]. In the reaction rate, the turbulence-chemistry interaction was considered using the kinetic rate/eddy dissipation rate model [25] with the realizable k-ε model employed for turbulence.

Regarding the heat transfer, radiation was solved using the discrete ordinate method with the weighted sum of the gray gases model for gas absorption [28]. In the boundary condition, the entire furnace wall (evaporator) was set to have the average water/steam temperature of 652.15 K with an overall heat transfer coe fficient of 250 <sup>W</sup>/m2K and inner wall emissivity of 0.7. The heat transfer coe fficient was taken as the average of measured values 3.51–4.37 m2K/kW of thermal resistance equivalent to 229–285 <sup>W</sup>/m2K depending on the coal types [29]. The tube bundles from primary SH to the economizer were simplified as porous zones, with local source terms calculated for flow resistance, convection, and radiation based on the detailed geometry and steam conditions [30]. In brief, the inertial resistance in the lateral and transverse directions was calculated using the Jakob correlation [31]:

$$
\Delta \mathbf{p} = \frac{2f' \mathbf{G}\_{\text{max}} \, ^2 \text{N}\_L}{\rho} \left(\frac{\mu\_s}{\mu}\right)^{0.14} \tag{1}
$$

$$f' = [0.044 + \frac{0.08\left(\frac{S\_L}{D}\right)}{\left(\frac{S\_T}{D} - 1\right)^{0.43 + 1.13D/S\_L}}]Re^{-0.15} \tag{2}$$

Similarly, the heat transfer rate in each cell was calculated as the sum of convection (.*qconv*) and radiation (.*qrad*) for the specific tube surface area per unit volume (*As*):

$$\dot{q}\_{\rm conv} = f\_{\rm conv} A\_s \frac{k \,\mathrm{Nu}}{D} (T\_{\rm gus} - T\_s) \left[ \mathrm{W/m^3} \right] \tag{3}$$

$$\dot{q}\_{rad} = f\_{rad} \varepsilon \sigma A\_s \left( T\_{\text{gas}}{}^4 - T\_s{}^4 \right) \left[ \text{W/} \text{m}^3 \right] \tag{4}$$

The Nusselt number (*Nu*) for convection was determined from the Zukauskas correlation [32]:

$$Nu = 0.40 Re\_D^{0.6} Pr^{0.36} \tag{5}$$

The above equations include two correction factors (*fconv* and *frad*) that were introduced to consider the slagging/fouling factors and tuned to match the design values. The formulations for the flow resistance and heat transfer rates were implemented using UDFs.

NO reactions were calculated by post-processing of the CFD results. Thermal NO reactions were based on the extended Zeldovich mechanism with the rate constants taken from Hanson and Salimian [26] and the radical concentrations acquired from the respective partial equilibrium assumptions. In the fuel NO mechanism [27], the partitioning of fuel-N between volatile-N and char-N was determined using FLASHCHAIN. The release of N intermediates during the devolatilization of the low rank coal was assumed to be HCN 5:NH3 1, whereas the char-N was converted directly to NO. The reduction of NO to N2 on the active char surface was considered with a microscopic surface area that was assumed to be 100 m<sup>2</sup>/g. Prompt NO was ignored because its concentration is known to be minor [3].
