PRES22_0231.docx DOI: 10.3303/CET2294166 Paper Received: 15 April 2022; Revised: 10 May 2022; Accepted: 20 May 2022 Please cite this article as: Strouhal J., Juřena T., 2022, CFD Simulations of Packed-Bed Combustion: Accounting for Irregularities in Particle Orientation in the Estimation of Bed Radiative Properties, Chemical Engineering Transactions, 94, 997-1002 DOI:10.3303/CET2294166 CHEMICAL ENGINEERING TRANSACTIONS VOL. 94, 2022 A publication of The Italian Association of Chemical Engineering Online at www.cetjournal.it Guest Editors: Petar S. Varbanov, Yee Van Fan, Jiří J. Klemeš, Sandro Nižetić Copyright © 2022, AIDIC Servizi S.r.l. ISBN 978-88-95608-93-8; ISSN 2283-9216 CFD Simulations of Packed-Bed Combustion: Accounting for Irregularities in Particle Orientation in the Estimation of Bed Radiative Properties Jiří Strouhal, Tomáš Juřena Brno University of Technology, Institute of Process Engineering, Technická 2, 61669 Brno, Czech Republic Jiri.Strouhal@vutbr.cz Modelling of packed-bed combustion includes the problem of modelling the radiative heat transfer. Thermal radiation significantly contributes to heat transfer both across the bed and from the hot gases to the top of the fuel layer. Several methods can be utilized aside from the so-called effective bed conductivity, such as the P1 or the Discrete Ordinates Model. When using the porous media for the fuel representation, its effective radiative properties must be determined, namely the extinction coefficient Gómez et al. (2014). It is often determined by assuming regularly oriented and mono-sized particles Shin and Choi (2000). Expressions for its estimation differ in the quantity representing the fuel and bed properties, for example, the bed porosity or particle emissivity. However, the absorption and scattering are influenced not only by how densely the particles are packed but also by their mutual position and orientation. The estimated radiative properties can significantly influence the combustion simulation, both quantitatively (e.g., the thickness of the reaction front, and ratios of the produced gases) and qualitatively (ignition of the fuel or the flame extinction). This work compares expressions of the radiative properties, that are modified by introducing irregularities in the positions or orientation of the fuel particles. The model for the radiative properties is implemented as a part of the radiation model based on the Discrete Ordinates method, modified for the case of two Eulerian phases representing the fuel and gas phase. In the presented work, these models for radiative properties are compared in a case of a fixed-bed reactor described by Zhou et al. (2005) and compared to the reported experimental results. The case is modelled in two dimensions. The radiation model is also described, including its implementation in ANSYS Fluent. 1. Introduction CFD based modelling of the packed-bed combustion provides valuable insight into related processes, including the thermal decomposition of the fuel or for example, the development of solid and gaseous pollutants. In the recent review of Hoang et al. (2021), focusing on modelling techniques for municipal solid waste combustion, all three basic heat transfer mechanisms are reported to be included in packed-bed combustion models. The emphasis is put on models using the so-called effective thermal conductivity. This was widely used by authors of earlier works in the packed-bed combustion modelling (Zhou et al., 2005). Effective conductivity provided by Yagi and Kunii (1957) combines conduction in both the gas and the fuel, radiation between surfaces and inter-particle voids and the effect of convection and lateral mixing gas. An important note on including the radiative heat transfer using effective conductivity was made highlighted by Strieder (1997). One of the basic assumptions is relatively large interparticle distances and extension to the densely packed beds are not reliable using only the basic theory of the radiative heat transfer. This extension requires empirical validation and possibly evaluation of the radiative properties of the bulk. Other approaches can be however used to model radiation inside densely packed beds. If a continuum phase is used to represent the fuel, it presents a second participating medium in the radiative heat transfer, which’s radiative properties must be defined. The continuum representation of the fuel was used in the presented work. During the implementation of the radiation model, the choice of the model for the radiative properties – namely 997 the extinction coefficient – was found to influence, if for given conditions fuel can be ignited. Calculation of the extinction coefficient according to Shin and Choi (2000), usually used in packed-bed combustion modelling (Gómez et al., 2014), was implemented in our model. However, as the authors themselves point out, the model is derived from an ideal arrangement of fuel particles. In the study presented here, the original expression given by Shin and Choi (2000) is compared to two modifications aiming to better capture the particles’ shape and arrangement. Comparison is made in the case of an experimental fixed-bed reactor (Zhou et al., 2005). 2. Mathematical model The fuel is represented as a porous medium defined over the whole computational domain. The bed location is then specified by choosing a non-zero volume fraction of the solid phase. In this work, particles are modelled as hollow cylinders (representing straw). The shape and size of particles are introduced via parameters of the porous media, e.g., specific surface and mass and heat transfer coefficient. Particles are assumed to be thermally thin. Transport equations for the mass conservation of both phases were adopted mainly from the work of Zhou et al. (2005), among others including the heterogeneous reactions such as moisture evaporation, devolatilization and fixed carbon oxidation. Mass conservation was solved for each of the four main fuel components – moisture, volatiles, fixed carbon, and ash. Permeability and inertial loss coefficient for modelling the Pressure loss in the bed were estimated using Ergun’s equation. The heat transport in the solid phase is modelled using a modified version of the model of Fjellerup et al. (2003), which calculates the effective heat conductivity of the bed. The original effective conductivity included effects of the flow and radiation, which in this work were excluded from the effective conductivity and included via the convective heat exchange between solid and gas and by solving the radiative heat transfer separately. 2.1. Model description At a large scale, the fuel bed acts like an optically thick medium. As the combustion is simulated in two or three dimensions, the Discrete Ordinates Model (DOM) is selected as appropriate for modelling the radiative heat transfer across the bed with an arbitrary directional resolution. The radiative heat transfer is solved for a defined set of discrete directions spanning the total solid angle, each represented by a vector 𝑠𝑠𝑖𝑖 and weight 𝑤𝑤𝑖𝑖 and solving a radiative transport equation (RTE), as expressed with Eq(1) (Modest, 2003). 𝑠𝑠𝑖𝑖 ⋅ ∇ 𝐼𝐼(𝑠𝑠𝑖𝑖) + �𝜂𝜂 𝛼𝛼𝑔𝑔 + (1 − 𝜂𝜂)𝛽𝛽𝑠𝑠� 𝐼𝐼(𝑠𝑠𝑖𝑖) = 𝜎𝜎𝑅𝑅 𝜋𝜋 �𝛼𝛼𝑔𝑔𝜂𝜂 𝑇𝑇𝑔𝑔4 + 𝛼𝛼𝑠𝑠(1 − 𝜂𝜂) 𝑇𝑇𝑠𝑠4� + 𝜎𝜎𝑠𝑠 4 𝜋𝜋 �𝑤𝑤𝑗𝑗 𝐼𝐼�𝑠𝑠𝑗𝑗 ′� 𝛷𝛷𝑠𝑠�𝑠𝑠𝑖𝑖 ⋅ 𝑠𝑠𝑗𝑗 ′� 𝑗𝑗 (1) 𝛽𝛽𝑠𝑠 is the extinction coefficient of solid, 𝜎𝜎𝑅𝑅 is the Stefan-Boltzmann constant, 𝛼𝛼 and 𝜎𝜎 are the absorption and scattering coefficients and 𝛷𝛷𝑠𝑠 is the scattering function. 𝜂𝜂 is the bed external porosity, i.e., the volume fraction of the gas outside the fuel particles. The absorption coefficient of the gas was modelled using the domain-based weighted-sum-of-grey-gases model. Scattering in the gas is assumed to be negligible. The extinction coefficient 𝛽𝛽𝑠𝑠 represents the ratio of the radiative intensity 𝐼𝐼 of the radiation transmitted through a medium by a distance Δ𝑥𝑥. and the initial intensity 𝐼𝐼0 (Modest, 2003). The ratio 𝐼𝐼/𝐼𝐼0 can be also interpreted as a probability 𝑃𝑃 of finding an unobstructed path of length Δ𝑥𝑥 , passing through the porous for the given direction. Then, 𝛽𝛽𝑠𝑠 can be expressed with Eq(2). 𝛽𝛽𝑠𝑠 = −𝑙𝑙𝑙𝑙 𝑃𝑃/Δ𝑥𝑥 (2) Using the same approach as in the work of Gómez et al. (2014), the absorption and scattering coefficients can be estimated using the surface emissivity, as expressed with Eq(3) and Eq(4), assuming opaque particles. 𝛼𝛼𝑆𝑆 = 𝜀𝜀𝑝𝑝 𝛽𝛽𝑠𝑠 (3) 𝜎𝜎𝑆𝑆 = �1 − 𝜀𝜀𝑝𝑝� 𝛽𝛽𝑠𝑠 (4) A simple model of the extinction coefficient was proposed by Shin and Choi (2000). They assumed, that the radiation in each direction is partially blocked by a layer of ideally arranged particles. Although it is not directly mentioned, the spherical shape of particles was assumed, as explicitly stated in other parts of their combustion model. Eq(5) is the proposed expression for the extinction coefficient. 𝛽𝛽𝑠𝑠 = − 2 ln 𝜂𝜂 /𝑟𝑟𝑝𝑝 (5) Assuming that the particles were arranged as assumed by the authors, the actual attenuation of radiation (in the direction normal to the particle layer) would be given by the surface area of the particles’ shadow instead of the volumetric fraction, leading to the expression for 𝛽𝛽𝑠𝑠 to be − ln(4𝜂𝜂/𝜋𝜋)/𝑑𝑑𝑝𝑝 for cylindrical and −ln (1.5 𝜂𝜂)/𝑑𝑑𝑝𝑝 998 for spherical particles, as illustrated in Figure 1. Assuming more particles than just a single layer leads to a further increase in the value of 𝛽𝛽𝑠𝑠. Figure 1: Particles’ arrangement assumed by Shin and Choi (2000) for estimation of 𝛽𝛽𝑠𝑠 . The grey area represents the actual transmitted radiation, dotted lines illustrate the amount of radiation as given by the proposed expression for 𝛽𝛽𝑠𝑠. Strieder (1997) in his work focused on solving radiative heat transfer inside a bed consisting of mono-sized cylinders. The probability of finding an unobstructed path is expressed with Eq(6). Each particle is assigned to a certain plane with a given orientation in space, to which the particle’s axis is perpendicular. The summation goes over all these particle groups, denoted by index 𝑗𝑗. 𝑙𝑙𝑗𝑗 is the number of particles in the group per unit area and 𝜃𝜃𝑗𝑗 is the angle between the plane and the direction of radiation. Mutual overlap of the particles is allowed. ln 𝑃𝑃 = �− 𝑑𝑑𝑝𝑝𝑙𝑙𝑗𝑗Δ𝑥𝑥 cos (𝜃𝜃𝑗𝑗) 𝑗𝑗 (6) Note, the original expression for 𝑃𝑃 used by Strieder (1997) contains a second term under the summation, related to a probability that the starting point of the path does not lie inside any particle. In our case, we investigate the absorption of radiation starting only from points inside voids between particles. The second term (probability of choosing the starting point inside some particle) is then zero by definition. The porosity 𝜂𝜂 can be viewed as the probability of randomly choosing a point lying inside bed voids. As such, 𝜂𝜂 can be expressed with Eq(7). ln 𝜂𝜂 = �−𝜋𝜋 𝑑𝑑𝑝𝑝2/4 𝑙𝑙𝑗𝑗 𝑗𝑗 (7) For a continuous variation in the particles’ orientation, the summation used in Eq(5) and Eq(6) is replaced by integration from 0 to 𝜋𝜋/2 and 𝑙𝑙𝑗𝑗 has to be substituted with 𝜉𝜉(𝜃𝜃) sin 𝜃𝜃. For two cases of particle orientation, the expressions for the extinction coefficient are listed in Table 1. Table 1: Three expressions for the extinction coefficient: fully random orientation (RP); original model of Shin and Choi (2000) (SC); all particles perpendicular to the radiation direction (PP). RP SC PP 𝛽𝛽𝑠𝑠 [1/m] − 2 𝜋𝜋 ln 𝜂𝜂 𝑑𝑑𝑝𝑝 − ln 𝜂𝜂 𝑑𝑑𝑝𝑝 − 4 𝜋𝜋 ln 𝜂𝜂 𝑑𝑑𝑝𝑝 2.2. Model implementation The combustion model was implemented via ANSYS Fluent. The porous zone model with physical velocity formulation was used to model the pressure loss and increase of the flow velocity inside the bed. However, combining this model with the DOM does not allow for using the Non-Equilibrium Model for modelling the porous medium temperature. The fuel enthalpy along with several other quantities related to the fuel was included using the so-called User Defined Scalars (UDS), enabling to define of a new transport equation for an arbitrary physical quantity (ANSYS, Inc. 2021). The porous zone model was then used mainly for capturing the reduction of the gas phase volume fraction and the pressure loss and the porous media density, thermal conductivity and specific heat capacity were set to zero to avoid its influence on the energy balance. The effect of the porous solid on radiation is included by adjusting the transport equation during each iteration. Similar (but not identical) adjustments were used by Gómez et al. (2013). First, the absorption coefficient of the gas is replaced by the volumetric average of 𝛼𝛼𝑔𝑔 and 𝛼𝛼𝑠𝑠 as expressed with Eq(8), using the macro DEFINE_WSGGM_ABS_COEFF. Also, the scattering coefficient is replaced with (1 − 𝜂𝜂) 𝜎𝜎𝑠𝑠. 999 𝛼𝛼 = 𝜂𝜂 𝛼𝛼𝑔𝑔 + (1 − 𝜂𝜂) 𝛼𝛼𝑠𝑠 (8) After this adjustment, the emission term in Eq(1) is still not in the correct form. To account for the difference between the gas and fuel temperature, a source term given by Eq(13) is introduced into the RTE using macro DEFINE_DOM_SOURCE. 𝑆𝑆𝑅𝑅𝑅𝑅𝑅𝑅 = (1 − 𝜂𝜂) 𝛼𝛼𝑠𝑠 𝜎𝜎𝑅𝑅 𝜋𝜋 �𝑇𝑇𝑠𝑠4 − 𝑇𝑇𝑔𝑔4� (13) Fluent automatically assigns all the absorbed/emitted radiative energy (calculated by the modified DOM) to the gas phase. The last adjustments include a user-defined source term Δ𝑆𝑆𝑟𝑟𝑟𝑟𝑟𝑟 , added to the solid and subtracted from the gas enthalpy transport equation. This source term is expressed in Eq(14). Δ𝑆𝑆𝑟𝑟𝑟𝑟𝑟𝑟 = (1 − 𝜂𝜂) 𝛼𝛼𝑠𝑠 ��𝑤𝑤𝑖𝑖 𝐼𝐼(𝑠𝑠𝑖𝑖) 𝑖𝑖 − 4 𝜎𝜎𝑅𝑅𝑇𝑇𝑠𝑠4� (14) 3. Simulation setup The experimental reactor described by Zhou et al. (2005) was simulated. A two-dimensional section of the reactor was used as the model geometry. Only the main dimensions reported by the authors, the air inlet and the flue gas outlet were derived from the scheme of the reactor in the original publication. The schematic representation of the used geometry is shown in Figure 2. The domain was discretized using approximately 9,200 quadrilateral cells with a maximum size set to 5 mm. 12 prismatic layers were used near the walls, with the boundary cell height set to 0.5 mm. Figure 2: The model geometry, adapted from the description of the experimental reactor. Bed dimensions and operating conditions were kept constant. Zhou et al. (2005) provided data on the fuel proximate and ultimate analysis. The LHV was estimated using a correlation provided by van Loo and Koppejan (2008). As in the work of Juřena and Hájek (2011), the composition of the devolatilization products was simplified using a representative compound, with properties derived from mass and energy balance, fuel composition and the LHV. According to the experimental observations of Zhou et al. (2005), the (external) porosity was kept at a constant value of 0.58 throughout the simulation, with straw particles having an outer diameter of 4.2 mm. The wall emissivity was set up according to Gómez et al. (2014), assuming a partially oxidized wall surface. Inlet boundary conditions and case initialization were set up according to the experiment described by Zhou et al. (2005). In their work, the straw emissivity is equal to the emissivity of straw char (𝜀𝜀𝑝𝑝 = 0.9), which was used in our work. The authors mention that by decreasing the amount of char in the solid residuum, the emissivity may significantly decrease. Other values are used in works on packed-bed combustion modelling. A value of 0.8 was reported as a typical value for solid fuels by Yang et al. (2007) and used for simulations of MSW combustion. The inlet air temperature was set to 20°C, its mass flow rate 3.35∙10-3 kg/s, with a composition of 75.94 % N2, 23 % O2, 1.02 % H2O and 0.04 % CO2 (by mass). 4. Results and discussion The data obtained from the simulated 200 s period are compared to those reported by Zhou et al. (2005). In the experiments, the straw was fed into the reactor and ignited with a gas burner. After the thermocouple below the bed top detected a temperature rise, the burner was removed. In the simulation done by the same authors, ignition was done by setting the temperature at the bed top to 800°C for 2 s. In the work presented here, the fixed temperature was set across the whole volume above the bed top and in the first layer of fuel containing cells. The temperature was kept at a minimum of 800°C for 6 s to keep the fuel burning after switching the ignition source off. The first observable results have shown, that in the case of randomly oriented particles (RO) the flame is extinguished very shortly after the ignition source is turned off, which can be seen on the graph of the maximum temperature during the calculation for all three cases (Figure 3). 1000 The relative imbalance of the total mass and energy was monitored throughout the calculation. The definition of the relative imbalance can be found in the work of Juřena and Hájek (2011). The values at 250 s and the maximum reported values are listed in Table 2. Table 2: Comparison of the maximum mass and energy relative imbalances recorded throughout the simulated period 250 s. RO – randomly oriented particles, SC – extinction coefficient is given by Shin and Choi (2000), PP – particles perpendicular to the radiative flux direction. RO SC PP mass maximum 0.0036 % 0.00024 % 0.00024 % final 0.0013 % 0.00023 % 0.00024 % energy maximum 0.040 % 0.39 % 0.57 % final 0.030 % 0.26 % 0.39 % Figure 3: Temperature profiles during the first minute of simulated combustion, obtained using three models for 𝛽𝛽𝑠𝑠. Temperatures in the range of 727°C – 1,027°C were measured during the experiments (Zhou et al., 2005). Nine points representing the thermocouples in the experimental reactor were placed at the reactor axis inside the bed. These were labelled from the bed bottom to the top as T2-T9. The period between the same given gas temperature is recorded at the furthest thermocouple (T6). The periods were almost identical for both the SC (218 s) and PP case (224 s) and the temperatures reported by experimenters (≈ 230 s). The average temperature recorder during the 250 s was 1,197°C (PP) and 1,173°C (SC), being significantly higher compared to the maximum temperatures observed during experiments (727°C–1,027°C). As both the PP and SC models differ from the experiment almost identically, the cause of the higher temperature might be expected to relate to using too high an estimation of the straw LHV. Predicted molar fractions of CO2 and CO at the bed top were 12.04 % and 4.7 % for the PP and 11.7 % and 5.11 % for the SC case. Predicted contents were lower compared to the experiment (approximately 16 % of CO2 and 13 % of CO). 5. Conclusions Simulation of three minutes of the combustion process demonstrated the difference in assuming a particular arrangement of fuel particles even when all the particle properties and the bed porosity are identical. For the case of the experimental reactor and straw bed, the main results are the following: • No significant difference was observed between the original expression for the extinction coefficient proposed by Shin and Choi (2000) and the expression derived from the randomly placed particles oriented perpendicularly to the radiation direction was observed. • Using the model with randomly orientated straw particles led to the flame extinction within 2 s after the artificial heat source representing the auxiliary gas burner was turned off. The value of the extinction coefficient for this particle arrangement The observed difference points to the fact, that the bed structure cannot be automatically considered as an insignificant parameter compared to other parameters, although in practice and namely in the case 0 200 400 600 800 1000 1200 1400 1600 1800 0 10 20 30 40 50 M ax im um g as te m pe ra tu re [° C ] Time [s] PP RO SC 1001 of heterogeneous fuels, e.g., municipal solid waste, the structure cannot be always predicted in such a detailed level to obtain the accurate distribution of particle orientations. Nomenclature dp – fuel particle diameter, m I – radiative intensity, W/(m2∙sr) n – number of particles per unit volume, 1/m3 P – probability, - s⃗ – directional vector, - SRTE – energy source term in RTE, W/(m3∙sr) T – temperature, °C w – weighting factor, - α – absorption coefficient, 1/m ΔSrad – enthalpy source term, W/m3 Δx – path distance, m ε – emissivity, - η – external porosity, - θ – angle, rad ξ – number of particles per unit area and unit solid angle, 1/(m2∙sr) σ – scattering coefficient, 1/m σR – Stefan-Boltzmann constant, W/(m2∙K4) Φs – scattering function, - Ω – solid angle, sr Subscripts g – gas s – solid p – particle Acknowledgements This research is realized within the project Quality Internal Grants of BUT (KInG BUT), Reg. No. CZ.02.2.69/0.0/0.0/19_073/0016948, which is financed from the OP RDE (Operational Program, Development and Education) References ANSYS, Inc., 2021, Ansys® Fluent, Release 21.2, Help System, Fluent User’s Guide, ANSYS, Inc. Fjellerup J.S., Henriksen U.B., Jensen A., Jensen P.A., Glarborg P., 2003, Heat Transfer in a Fixed Bed of Straw Char, Energy & Fuels 17 (5), 1251–58. Gómez M.A., Patiño D., Comesaña R., Porteiro J., Álvarez Feijoo M.A.., Míguez J.L., 2013, CFD Simulation of a Solar Radiation Absorber, International Journal of Heat and Mass Transfer 57 (1), 231–40. Gómez M.A., Porteiro J., Patiño D., Míguez J.L., 2014, CFD Modelling of Thermal Conversion and Packed Bed Compaction in Biomass Combustion, Fuel 117, 716–32. Hoang Q.N., Vanierschot M., Blondeau J., Croymans T., Pittoors R., Caneghem J., 2021, Review of numerical studies on thermal treatment of municipal solid waste in packed bed combustion, Fuel Communications 7, Juřena T., Hájek J., 2011, Energy Considerations in CFD Modelling of Biomass Combustion in an Experimental Fixed-bed Reactor, Chemical Engineering Transactions, 803–808 Modest M.F., 2003, Radiative Heat Transfer. 2nd ed. Amsterdam, The Netherlands, Boston, USA: Academic Press. Shin D., Choi S., 2000, The Combustion of Simulated Waste Particles in a Fixed Bed, Combustion and Flame 121 (1), 167–80. Van Loo S., Koppejan J., 2008, The Handbook of Biomass Combustion and Co-firing, Earthscan, London, UK. William S., 1997, Radiation Heat Transport in Disordered Media, Advances in Water Resources, Advances in Heat Transfer in Porous Media, 20 (2), 171–87. Yagi S., Kunii D., 1957, Studies on Effective Thermal Conductivities in Packed Beds, AIChE Journal 3, Vol 3, 373–81. Zhou H., Jensen A., Glarborg P., Jensen P.A., Kavaliauskas A., 2005, Numerical modeling of straw combustion in a fixed bed, Fuel 84, 389–403. 1002 PRES22_0343.pdf CFD Simulations of Packed-Bed Combustion: Accounting for Irregularities in Particle Orientation in the Estimation of Bed Radiative Properties