Single particle conversion of woody biomass using fully-resolved and Euler–Lagrange coarse-graining approaches

The conversion of woody biomass is studied by means of a layer-based model for thermally-thick biomass particles (Thunman et al. 2002, Ström et al. 2013). The model implementation is successfully validated against experiments that study particle conversion in a drop tube reactor. After this validation step, this work focuses on the well-known problem of grid dependence of two-phase numerical simulations using the standard Euler– Lagrange (EL) framework. This issue is addressed and quantified by comparing EL data that models the particle boundary layers to corresponding simulations which fully resolve these boundary layers (fully-resolved, FR, simulations). A comparison methodology for the conceptually different FR and EL approaches by extracting the heat transfer coefficient from the detailed FR simulations is proposed and confirms that the EL results are strongly grid-dependent. This issue is overcome by applying a set of coarse-graining methods for the EL framework. Two coarse-graining methods are evaluated, a previously suggested diffusion-based method (DBM) and a new approach based on moving averages referred to as MAM. It is shown that both DBM and MAM can successfully recover the detailed FR data for pure particle heating for a case where the grid size is half the particle diameter, i.e. when the standard EL method fails. Both coarse-graining methods also give improved results for an EL simulation that considers the more complex combined physics of particle heating, drying and devolatilisation, given that the CG model parameters that scale the corresponding CG interaction volumes are sufficiently large. Based on the available FR data, recommended model parameter ranges for DBM and MAM are provided as a function of normalised boundary layer thickness. The novel MAM approach is shown to be significantly more efficient than the DBM and therefore suitable for future EL simulations with multiple particles.


Introduction
The global energy demand is expected to increase in the next few decades as a result of the rising world population and economic growth.Fossil fuels are the most dominant global energy source, but the need for alternative sources increases due to the limited resources and negative environmental impact of fossil hydrocarbons, particularly coal.Solid biomass [1] is considered as a promising renewable alternative to coal due to the possibility of retrofitting existing coal-fired power plants.This offers the opportunity of making use of existing supply chains, infrastructure and proven power generation technology, while solid fuel particle conversion.A first attempt to develop thermophysical sub-models for biomass was to adapt reliable coal sub-models to biomass [14][15][16][17][18][19][20].Typical predictions show reasonable results for small particle diameters but a mismatch between simulations and experiments occurs for larger particle diameters [21].A major reason for these discrepancies is attributed to the fact that pulverised coal particle diameters are typically of (small) micrometre size, which satisfies the thermally-thin particle assumption, i.e. uniform intra-particle temperatures, that is invoked for most models.Conversely, typical pulverised biomass particle diameters are of (large) micrometre to millimetre size, which results in non-negligible intra-particle temperature gradients, i.e. thermally-thick particles, during the conversion process [22,23].To overcome this issue, several researchers have developed onedimensional thermophysical sub-models specifically for thermally-thick biomass particles.Those sub-models can generally be classified into two types [24], namely interface-based models [25][26][27][28][29][30][31][32] and mesh-based models [33][34][35][36].The mesh-based models solve the spatially discretised conservation equations on mesh points along the radius of the particle, which requires intra-particle discretisation that becomes costly for large particle ensembles.The underlying concept of interface-based models is to divide the particle into a set of layers, with different particle conversion processes like heating, drying, devolatilisation and char conversion occurring at the interface between two neighbouring layers.Since such interface-based sub-models are relatively easy to implement, have lower computational cost and ensure a higher model fidelity compared with simpler approaches, a layer-based particle conversion sub-model for woody biomass [25,29] is used in the present work.
In the CFD modelling of multiphase gas-solid systems, the interaction between the gas phase and the solid particles needs to be specified as illustrated in Fig. 1.Various coupling methods exist, but the most detailed approach to describe gas-solid interactions is the fully-resolved (FR) particle method [37][38][39][40].In FR simulations, the governing equations (for total mass, momentum, energy and species) are solved in the Eulerian framework and all relevant gradients are resolved on the computational grid, including the particle boundary layers.In the FR approach the exchange of momentum, mass and heat between gas and solid is considered through a set of boundary conditions at the interface between both phases, while -different from Fig. 1 -the governing equations do not contain any two-phase source terms.While this approach is highly accurate, it is associated with very high computational costs and therefore restricted to single particles and small particle groups [41][42][43].Conversely, the two-way coupled Euler-Lagrange (EL) approach provides a good trade-off between accuracy and computational cost for dilute two-phase flows.In the EL approach, the bulk gas phase outside the particle boundary layers is treated as a continuum and solved analogously to the FR approach.However, the particle boundary layers are modelled and the solid particles are treated as dispersed point-particles in the Lagrangian framework.To model the crucial exchange processes across the particle boundary layers the point-particles need to retrieve information on the mass, momentum, energy and composition of the surrounding gas phase.This information is typically obtained from the local Eulerian cell () in which the particles reside.The reverse coupling from the particles to the gas phase is realised through a set of source terms in the Eulerian governing equations, cf.Fig. 1.However, for a reliable EL method contradictory requirements arise if small Eulerian cells are required to resolve the relevant gas phase physics (e.g.small turbulent vortices), while large Lagrangian particles are present.A typical assumption of the EL method is that the particle diameters are much smaller than the Eulerian grid size, i.e.   ∕ ≪ 1.For flows laden with small evaporating fuel droplets detailed indicators for the EL resolution requirements exist [44], but due to the larger size of biomass particles and the inherent differences between droplet evaporation and solid particle conversion these recommendations are not directly applicable to biomass.For biomass, the ratio   ∕ is typically of order one or even larger, which invalidates the standard EL approach.In this situation, exaggerated coupling effects between gas and solid phase can occur, leading to unphysical results and simulation instabilities [45,46].To solve this issue, coarse-graining (CG) methods have been developed, see e.g.[47][48][49][50][51][52][53][54].Sun et al. [48] briefly summarise the existing CG methods, referred to as particle centroid method, divided particle volume method, statistical kernel method and two-grid method.All these methods follow a different approach but have the same aim: The distribution of the particle source terms to more cells than just the local Eulerian cell and the consideration of a fluid volume larger than the local cell to obtain farfield gas phase information (e.g.temperature) to trigger the particle processes (e.g.pyrolysis).Zhang et al. [51] have demonstrated that the interplay between the Eulerian and the Lagrangian framework strongly influences the prediction of solid fuel particle conversion due to the non-linear coupling processes.However, studies of EL approaches coupled to CG methods (EL-CG) in the past have shown that dependencies of grid size vs. particle diameter are present, but in most CG methods an unknown free parameter is required, which needs further investigation and quantification.
In the present study, we consider the conversion of single woody biomass particles using a layer-based particle conversion sub-model with two complementary modelling approaches for gas-solid coupling.At first, a FR laminar flow simulation is validated against experimental data to obtain well-resolved comparison data for the evaluation of the EL framework.Subsequently, the FR data is compared to the results from the standard EL method and a set of EL-CG approaches.The objectives of the present study are to • validate the thermally-thick layer-based biomass sub-model for particle conversion against recent experiments, • propose a methodology to meaningfully compare FR and EL approaches, • propose a simple and flexible CG method for enhanced EL-CG simulations, • compare the high-fidelity FR data with data from the standard (grid dependent) and coarse-graining (grid independent) EL approaches.
The remainder of the paper first introduces the modelling approach in Section 2, followed by experimental validation in Section 3 and the description of the computational setup in Section 4. Results are shown and discussed in Section 5, and the major conclusions are drawn in Section 6.

Gas phase
In the present section the governing equations of the gas phase are written in a general format including two-phase exchange terms.In the FR approach these source terms are set to zero and replaced by boundary conditions, while they are non-zero for EL simulations, as will be described later.The transport equations for total mass, momentum, sensible enthalpy and chemical species in the gas phase read

Solid phase
The thermally-thick spherical particles are described with the layerbased model originally proposed by Thunman et al. [25] and further developed by Ström et al. [29].The particles are divided into four layers, namely moist/wet wood (1), dry wood (2), char (3) and ash (4), see Fig. 2. The intra-particle temperature gradient is estimated assuming that each inner layer (subscript <>) has its own set of mass and temperature and each inner boundary (subscript <>) has its own individual temperature.At these boundaries, the physical conversion processes of drying, devolatilisation and char conversion take place.The outer boundary represents the particle surface (  =  ,4 ).All temperatures are calculated through a heat balance inside the particle, the outer surface temperature of which is additionally affected by the gas phase.The volume and surface areas for each layer are predicted from the mass balance, e.g. the progress of devolatilisation.The evaluation of the heat and mass balance can be simplified to a one-dimensional problem along the radial direction of the particle.The set of governing equations for all layers is given as where  is the th layer, −∇⋅ is the conductive heat transfer term,  is the heat of conversion on the boundaries and   is the heat source/sink term due to released gas within different layers (e.g.water vapour or volatiles).Eq. ( 6) may take different forms depending on the location of its evaluation [51].For layers, −∇ ⋅  is calculated as and for boundaries it is expressed as Here,   is the layer heat conductivity,   = 4 2  is the boundary surface area with radius .The temperature gradient along the radius is approximated for spherical particle as [29] Eqs. ( 7) and ( 8) describe the intra-particle heat balance between layers and boundaries.The interphase mass and heat balance between the solid phase and the gas are discussed later in Section 2.3.A shrinkage model following [25] is applied to consider the volume change in each layer due to mass conversion as (10) where  , is the change of mass in each layer and   is the shrinkage factor.
A drying model [55] is applied at boundary 1, between layers 1 and 2. The rate of drying is calculated from the heat flux to this boundary.The temperature of the drying front is limited by the boiling temperature of water at atmospheric pressure.A two-stage devolatilisation model is employed to model wood pyrolysis [56], which takes place at boundary 2 and is considered to be heat-neutral ( = 0) [23].One-step or two-step competing rates models commonly describe the devolatilisation process for solid fuels.The one-step model assumes a fixed char-to-volatile yield ratio, which only applies to systems with a single heating rate.The two-competing rates model overcomes this limitation and allows for varying char and volatile yields under different heating conditions.Both models neglect the details of the physical process of biomass degradation.However, as shown for example in [57], the leading order effects in the conversion of woody biomass can reliably be described by a two-stage model.In the primary stage, dry wood is decomposed into the three fractions light gas, tar and char via three competing release steps.In the secondary stage, tar further decomposes into light gas and char.Tar is treated as a gaseous species.The Arrhenius reaction rate coefficients of the devolatilisation model are listed in Table 1.The selected rates are validated against experiments on woody biomass conversion in Section 3. Char burn-out is neglected in the present work.

Gas/solid interphase
The interphase mass and heat balance differ depending on whether the FR or EL framework is employed.

Gas/solid interphase for FR simulations
When employing the FR approach the particle conversion processes are treated as boundary conditions for the gas phase and, therefore the interphase source terms in Eqs. ( 1)-( 4) are set to zero [38]: At the particle surface (subscript <>) two conditions need to be fulfilled for every surface cell layer in order to obtain consistent particle surface temperatures.These are for the temperature value itself and for the temperature gradient Here, subscript <> denotes gas phase properties at particle position.From Eqs. ( 12) and ( 13), the mean surface temperature ( ,4 ) averaged over the fully-resolved particle surface is calculated and serves as an input value for the intra-particle temperature change in Eqs.(7) and (8).The mass released by drying and devolatilisation is imposed uniformly and normal to the particle surface direction with a Robin boundary condition [61] as where   =  ,4 is the surface area of the particle and  , is the gas diffusivity at particle surface.

Gas/solid interphase for EL simulations
When conducting EL simulations, the particle boundary layers are not resolved, hence, correlation models for mass and heat transfer are required.The source terms Ṡ, , Ṡ, , Ṡℎ  , and Ṡ, on the RHS of Eqs. ( 1)-( 4) denote the gas/solid coupling terms and are calculated for a single particle as where   is the volume of the local computational cell,   is the particle mass,  p, is the particle specific heat capacity,  is the temperature and  con is the convective heat transfer time scale. moist is the mass fraction of moisture (i.e.= H 2 O) content in the particle,  ,vol is the mass fraction of species  in the volatile composition,  ,drying is the mass release of moisture and  ,devol is the mass release of volatile content in the particle.The convective heat transfer time scale follows Ranz-Marshall [62]   = 1 6 Pr Nu with Here,  is the dynamic viscosity and Nu is the non-dimensional Nusselt number.Subscript < > denotes gas properties at film temperature evaluated with the 1/3-law as The boundary layer heat balance affecting the particle surface is given as with the Euler-Lagrange heat transfer coefficient defined as

Biomass properties
The biomass employed in this work is torrefied wood, used in the experimental investigations by Tolvanen et al. [63,64].The proximate and ultimate analyses are given in Table 2. Following [51], the released light gases have a presumed composition as listed in Table 3. Tar is represented by C 6 H 8 O with thermophysical properties taken from benzene [25].The thermal properties for each solid layer are shown in Table 4.

Coarse-graining methods
As previously outlined in Section 1 EL methods, where the gas phase governing equations ( 1)-( 4) contain interphase exchange terms, are in need of coarse-graining schemes for situations where the particle diameter is of the order of the Eulerian grid size or larger.For such cases, a CG method must be employed to distribute the particle source terms to more cells than the local Eulerian cell and to interrogate a fluid volume larger than the local cell for far-field gas data.The vast majority of CG approaches is based on four basic methods, namely the particle centroid method (PCM), divided particle method (DPM), twogrid method (TGM) and statistical kernel method (SKM).The original SKM is a method that applies statistical kernel/weight functions to map input data to an alternative space.Due to their flexibility and simplicity of implementation we only focus on the PCM and variants of the SKM in the present work, since these methods can be applied to the most general case of unstructured grids, moving particles and parallel computing with relative ease.

Particle centroid method (PCM, i.e. standard EL)
The most common and widely used method is the PCM, or standard EL approach.In the PCM, two-phase coupling is restricted to the grid cell, where the particle centroid is located, see Fig. 3 (left).All gas phase properties ( gas ) are retrieved from this local Eulerian cell and the source terms ( source ) are assigned to the same cell.For ratios of   ∕ ≪ 1, when using reliable closures for the momentum, heat and mass transfer across the particle boundary layers, the PCM works well.However, for cases where the ratio   ∕ is close to unity or larger the error becomes unacceptable and causes unphysical behaviour, which may also lead to numerical instability.

Diffusion based method (DBM)
To overcome the limitations of the PCM, the diffusion based method (DBM) can be employed, see Fig. 3 (middle).The DBM belongs to the previously mentioned class of SKM, as it is based on the concept of a statistical kernel function to retrieve the reference gas field and to distribute the source terms.However, the DBM does not use an explicit kernel function, but solves a diffusion equation instead, which is equivalent to applying a kernel function [48].Before the particle  retrieves the far-field information from the gas phase, the unsteady diffusion equation is solved to obtain the reference gas phase properties  gas .Conversely, to distribute the two-phase source terms to a larger volume, these source terms are first calculated and then diffused to the neighbouring cells according to the unsteady diffusion equation The diffusion equation is solved for a time period 0 ≤  ≤  which is equivalent to a Gaussian kernel function with bandwidth

Moving average method (MAM)
In the present work, we propose a simple and flexible alternative CG approach, referred to as moving average method (MAM).Similar to the DBM from Section 2.5.2, the MAM approach can be considered as a variant of the SKM, employing a uniform kernel function, which is equivalent to applying the number average method.The MAM uses an averaging operation for retrieving  gas for the particle and for applying  source to the gas phase, see Fig. 3 (right).For the averaging operation, various geometrical shapes can be considered.In the present work, we choose a cubic distribution volume with edge length   .In the first step of MAM, the particle position is tracked and taken to be the centroid of the cubic distribution volume.Subsequently the gas phase properties are averaged according to with the gas phase properties   of each Eulerian cell inside the cubic volume and the total number of cells within the cube  cells .Conversely, the source terms from the underlying Eulerian grid are first summed up to a total source term for the entire cube  source and then equally distributed to all cells within the cube where   is the source term for each Eulerian cell within the cube.In a general EL-CG framework where particles are transported by a carrier flow, the location of the each moving particle (i.e. the centroid of the averaging volume) can be updated at every time step such that the averaging volume moves along with the particle, hence the term moving average method.Finally, we remark that the functionality and performance of MAM is independent of any given mesh structure.
The novel MAM approach constructs a look-up table of each cell's neighbours at the start of the simulation that remains unchanged and is efficiently read over the course of the simulation, to define the averaging volumes for MAM.

Validation
As a crucial validation of our numerical methodology, we compare our most accurate-and computationally expensive-modelling approach for conversion of the woody biomass considered in this work to experimental data.The modelling approach with the highest fidelity (and cost) is achieved by combining the layer-based biomass model for intra-particle conversion from Section 2.2 with the FR closures to accurately describe the boundary layer and far-field gas phase processes outlined in Section 2.3.1.The experimental data used for validation are the measurements conducted by Tolvanen et al. using a drop tube reactor (DTR) at different temperatures with solid fuel particles in an inert atmosphere (N 2 ) [63,64].The experimental campaign aimed to investigate the pyrolysis behaviour of coal, peat and various types of woody biomass, while we focus on the woody biomass data only.The experimental configuration comprises of a particle-feeding silo, a water-cooled movable injector, a steel tube with heating elements and a particle-collecting system.From the experiments it is known that the biomass particles are non-spherical and highly elongated.However, the experimentalists also showed that a simplified model, based on a spherical equivalent diameter, is able to capture the rate of mass loss [64].Therefore, the experimentally determined spherical equivalent diameter is used in this work to set up a fully spherical FR case for comparison with the experiments (see Fig. 5 in Section 4 for more details on the computational configuration).Fig. 4 shows the dry/ash-free mass loss of a torrefied biomass particle with   = 160 μm over time at different gas temperatures inside the DTF.The initial and boundary conditions have been extracted from previous experiments [63,64,71].The biomass particle needs  ≈ 0.1 s to heat up for the two high gas temperature cases (  = 1123∕1173 K) before starting its pyrolysis process.The conversion process takes ≈ 0.15 s until reaching the final state with only char and ash remaining in the particle.Due to the inert gas atmosphere, char conversion processes are suppressed, which explains that the mass loss does not reach a final value of one.However, a higher final volatile yield than the proximate analysis is due to the low heating rates in the chemical kinetic studies [63].The comparison of the FR predictions with the experimental data shows an excellent agreement both of the temporal evolution and the final mass release for the two high temperature cases.At lower gas temperatures (  = 973 K), the FR simulation predicts the conversion process to start slightly later and the whole process takes around 0.5 s, which is two times longer than with higher gas temperatures.A comparison with the experimental evidence shows differences in the time range 0.1 <  < 0.3 s, where the numerical model predicts a slightly later start of the conversion process.At late times the experimental data and FR predictions agree again very well.The differences between the numerical and experimental data at early times for   = 937 K may be attributed to the simplification of the particle shape (spherical) in the simulation, whereas elongated particles are present in the experiments.Moreover, there is some uncertainty in the activation energies reported in Table 1 (particularly No. 1 and 2) that govern the mass loss profiles show in Fig. 4, where the fits seem to work well for higher temperatures/heating rates, but are not ideal for the lower temperature/heating rate case   = 937 K.However, further experimental and numerical analyses are required to come to a definitive conclusion on these remaining deviations.
Overall, the predictions from the FR simulation approach using the layer-based biomass sub-model for intra-particle conversion show a good agreement with the experimental evidence.The FR approach correctly predicts the dependence of the mass loss rate and final particle conversion ratio on the ambient gas temperature in the DTF.Therefore, the detailed FR data is taken as a reference for the subsequent evaluation of the EL approach and EL-CG methods.

Computational setup
With our main focus on comparing detailed FR simulation data to results from the standard EL method and EL-CG methods for large particles, two basic computational setups for FR and EL simulations are considered.The three-dimensional computational domains and grid designs for both approaches are shown in Fig. 5.To obtain a welldefined setup for the FR vs. EL comparison the gas phase is initialised as a quiescent medium with  ,init = 1500 K and  N 2,init = 1.Due to the inert gas environment, no homogeneous reactions occur such that ωℎ  = ω in Eqs. ( 3) and ( 4) are zero.A single biomass particle is located in the centre of the domain with  ,init = 118 μm and  ,init = 350 K.
The computational domain for the FR setup has a size of   =   =   = 40 ⋅   and consists of 113,280 cells, the size of which varies from 6.5 ≤    ≤ 113 μm, with the smallest cells located at the particle surface and increasing cell sizes with increasing distance from the particle to fully resolve the particle boundary layer.The mesh is considered fixed and does not change during the simulations.
For the EL simulations, the domain has either 1331 (standard EL approach) or 274,625 (EL-CG methods) cells which are uniform and cubic.The domain size varies for the various investigated cases as   =   =   =  ⋅   ⋅  with  = {0.5, 1, 2, 4, 6, 8, 10, 20} and the domain is discretised in -, -, -direction with  grid cells, which corresponds to  EL =  ⋅   μm.For the standard EL approach  = 11, whereas  = 65 for the EL-CG methods.It is noted that, for increased computational efficiency, the domain size of the EL domain with the finest resolution is smaller than the size of the FR reference domain.However, the relevant physics are captured within the bounds of the smallest EL domain already, such that the boundary conditions do not unduly affect the predictions.
All the simulations are performed with a low-Mach second-order finite volume multiphase solver based on the OpenFOAM library.The transport equations of the gas phase are solved with a constant timestep of   = 5 ⋅ 10 −7 s, while the biomass particle sub-model considers sub-steps of size   = 1 ⋅ 10 −7 s.

Comparison methodology to study grid independence -FR/EL
Since the FR and EL modelling frameworks are fundamentally different, a detailed comparison methodology to compare results from the two approaches is required.This comparison methodology also serves to demonstrate the grid dependence of the standard EL and EL-CG methods.The steps of the comparison methodology are given as follows: 1.A FR reference case is solved for a specific time interval.2. As reflected by Eq. ( 13) in FR simulations the gas-solid coupling directly occurs at the (resolved) particle surface, from which a nominal heat transfer coefficient can be calculated.3. Detailed information on the gas phase thermophysical properties (i.e.  ), particle surface temperature (  ), temperature of the first cell layer around the particle (  ), initial gas temperature ( ,init ), surface area of the particle (  ) and the distance between the particle surface and the centroid of the first gas cell layer are retrieved from the FR simulation.4. The nominal FR heat transfer coefficient ℎ FR is evaluated at every time-step as and stored in a table.5.The boundary layer heat balance in the EL approach is given in Eq. (22), where (the modelled) ℎ EL is now replaced by the nominal value from the FR simulation ℎ FR : 6.At every time-step, ℎ FR is retrieved from the previously generated table and used for Eq.(29).7.With the heat flux provided by the FR simulation, the only remaining unknown parameter is   , which is the local gas cell temperature where the particle centroid resides.
The local gas cell temperature highly depends on the grid resolution, which allows for the evaluation of grid (in-)dependence of the EL approach.A schematic representation of the comparison methodology is shown in Fig. 6.

Results and discussion
After validating the layer-based sub-model for woody biomass conversion within the detailed FR simulation framework against the experimental evidence in Section 3 and providing a detailed comparison strategy between FR and EL data in Section 4.2 we can now evaluate the performance and grid (in-)dependence of the standard EL approach (EL-PCM) and the two alternative EL-CG approaches (namely EL-DBM and EL-MAM, see Sections 2.5.2 and 2.5.3).The evaluation is carried out separately for pure heating of the biomass particle in Section 5.1 and for combined heating, drying and devolatilisation in Section 5.2.

Pure heating
Fig. 7 shows the time evolution of the particle surface temperature over time for the FR setup, as indicated by the red solid line.The dashed lines show the predicted particle surface temperature over time for the standard EL approach (EL-PCM) with varying grid sizes.
Considering the FR data, the particle surface temperature strongly increases from the beginning of the simulation until  = 0.2 ms due to the injection of the initially cold particle (  = 350 K) into the high

Heating, drying and devolatilisation
After demonstrating that the CG-DBM and CG-MAM EL approaches are capable of recovering the detailed FR data for pure particle heating, we turn our attention to the more general case of heating, drying and devolatilisation of a biomass particle.Fig. 10 presents particle surface temperatures vs. time in a similar format as in the previous Section 5.1.The surface temperature profile of the FR case differs from the pure heating case because of the additional heat exchange due to the release of water vapour and volatile gases from the particle.The gases are released at particle surface temperature in the direction normal to the particle surface.Furthermore, the absorbed energy is used for the evaporation of moisture in the inner-most particle layer, that is included in the gases released from the particle surface and it also heats up the growing char layer.
The profiles predicted by the FR and (standard) EL approaches are identical for the benchmark case with  = ∞, again validating the comparison approach.Similar to the previous case of pure heating, refining the EL grid by decreasing the grid size also leads to increasingly larger deviations from the FR result for the present case of heating, drying and pyrolysis.Comparing the surface temperature at  = 80 ms between the cases  = ∞ and  = 0.5 ⋅   results in ≈ 50.5% relative error.
In analogy to Section 5.1 the standard EL case with  = 0.5 ⋅   is now taken as a reference and compared with results from the EL-DBM (Fig. 11) and EL-MAM (Fig. 12) approaches.As can be observed, both CG methods lead to significant improvements of the EL surface temperature predictions.This is achieved by either increasing the bandwidth of the DBM or increasing the edge length of the cube in MAM, both of which have a similar effect of considering a larger gas-solid interaction volume for CG.
With very large  MAM or  DBM , the FR data can faithfully be reproduced due to   → const., which is equivalent to the  = ∞ case in Fig. 10.Since it can be hypothesised that the required size of the CG interaction volume may be correlated with the thickness of the thermal boundary layer   around the particle, next Section 5.3 studies the relationship between   and suitable interaction volume sizes.

Boundary layer thickness and its relation to CG interaction volume
To establish correlations of boundary layer thickness   and CG interaction volume, we extract   from the FR data as the distance from the particle surface where   = 0.99 ⋅  ∞ for pure particle heating and   = 0.99 ⋅  ,max for particle heating, drying and devolatilisation.Correspondingly, the times when the (absolute) differences between the particle surface temperatures from FR and CG deviate by a threshold of 50 K are recorded for given values of the model parameters  MAM and  DBM from the cases shown in Figs. 8, 9, 11 and 12.Using this data, time ranges of validity of the tested parameter values can be evaluated.For example in Fig. 12 before the above-derived expressions for the CG model parameters can be employed.

Computational time using EL-CG methods
Apart from accuracy, an important aspect of employing EL-CG methods is the increased computational time associated with employing the coarse-graining approach.Using the CPU time for the general case of particle heating, drying and devolatilisation in the standard EL framework (EL-PCM) as a reference, the normalised CPU time for the CG methods is PCM = 1.0,DBM = 58.62 and MAM = 1.1514.These numbers reveal that the MAM is slightly more expensive than the PCM, whereas DBM is almost 60 times more costly.This is due to the fact that the solid-gas interaction volume in MAM is limited to the MAM cube (or alternative shape), whereas for the DBM, by default, diffusion is applied to the entire computational domain.At the beginning of the simulation, the MAM identifies the cell neighbourhood relations, finds the neighbouring cells of each control volume inside   and stores them in a table.During the simulation, these neighbours are simply retrieved from the table without further search steps for increased efficiency.Moreover, MAM is computed in the two principal steps of averaging and redistributing the average to the underlying grid cells only, whereas the DBM may involve the computation of a large number of diffusion time steps for large bandwidths.While the efficiency of the DBM could also be increased by choosing an appropriate diffusion range, such improvements have not been attempted here and solely the CPU time for the standard DBM is reported.

Conclusions
In the present work, a layer-based sub-model for woody biomass conversion [25,29] is coupled to both FR and EL simulation frameworks for solid particle conversion studies.A single cold biomass particle is located in a hot, quiescent, inert gas environment.The inert gas consists of nitrogen to suppress volatile combustion and char oxidation, and the heating, drying and devolatilisation of the particle is studied.The detailed FR setup is successfully validated against experimental data on torrefied woody biomass conversion in a drop tube reactor.Subsequently the layer-based sub-model is solved in the standard EL framework, showing that the EL approach is capable of reproducing the FR data for infinitely large grid sizes, while the standard EL method fails for highly refined computational grids where the grid size is close to or smaller than the particle size.To resolve this issue, a previously suggested diffusion based method is employed for coarsegraining.Moreover, a new simple and efficient CG approach based on moving averages referred to as CG-MAM is proposed.The two CG methods recover the detailed FR data for both pure particle heating and heating, drying and devolatilisation if the CG model parameters that scale the CG interaction volumes are sufficiently large.Using the available FR data, recommended model parameter ranges are provided as a function of normalised boundary layer thickness.CG-MAM is shown to give similar results as CG-DBM, but at a much reduced computational cost, which makes it suitable for future simulations with multiple Lagrangian particles.We note that the cases investigated in the present work do not include chemical reactions or particleparticle interaction effects that occur in reacting dense particle clouds.Including chemical reactions and particle-particle effects requires detailed analyses of computationally-expensive FR simulations of particle groups in chemically reacting flow for a relevant range of conditions.Such FR simulations, see e.g.[42,43], and their intricate relationship with CG methods for EL approaches (that would need to account for a range of thermal boundary layer thicknesses, and relative flameparticle and particle-particle distances to only name some principal parameters) are beyond the scope of the present study, but should be considered by the research community in future work.

Fig. 4 .
Fig. 4. Comparison of the dry/ash-free mass loss vs. time  from biomass particles with an equivalent diameter   = 160 μm at different temperatures between FR simulations (solid line) and experiments [63,64] (square, triangle and plus) in a drop tube reactor.(For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Fig. 5 .
Fig. 5. Computational setup for single biomass particle conversion in the context of FR (left) and EL (right) simulations.

Fig. 6 .
Fig. 6.Schematic representation of the comparison method between the FR and EL simulation approaches.

Fig. 7 .Fig. 8 .Fig. 9 .
Fig. 7. Pure particle heating: Comparison of particle surface temperature   vs. time from the FR (solid line) and standard EL (dashed lines) simulation approaches, where the EL grid size has been varied.(For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Fig. 10 .Fig. 11 .Fig. 12 .
Fig. 10.Particle heating, drying and devolatilisation: Comparison of particle surface temperature   vs. time from the FR (solid line) and standard EL (dashed lines) simulation approaches, where the EL grid size has been varied.(For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.) ) > indicating the coordinate directions, subscript <> particle properties,  is the density,  is the velocity, ℎ  is the sensible enthalpy,  ℎ  is the energy release from homogeneous chemical reactions,   and   are the mass fraction and homogeneous chemical reaction rates of chemical species  and  denotes time.The spatial coordinate is expressed by ,  is the dynamic viscosity,  is the Kronecker delta and Pr = Sc = 0.7 are the the non-dimensional Prandtl and Schmidt numbers.

Table 1
Arrhenius coefficients for the two-stage devolatilisation process.

Table 4
Thermophysical properties of the solid layers: Initial density , thermal conductivity , specific heat capacity   and shrinkage factor . Wet wood  1, = 931.069(assume  moist = 0.001 with  p,moist = 1000 J/(kg K)) [72]> 5⋅  is only valid for 0 ≤  < 48.6 ms, before a deviation greater than 50 K from the FR data occurs.Making use of the FR data on boundary layer thickness as a function of time (  =  ()), these validated CG model parameters can be mapped to given intervals of boundary layer thickness.Finally, the following general expressions for  MAM and  DBM as functions of normalised boundary layer thickness   ∕  can be found using regression  DBM ≈ 18.679 ⋅ 10 −5 (  ∕  ) 0.113  MAM ∕  ≈ 2.802(  ∕  ) 0.146 .Equivalently, for the combined physics of particle heating, drying and devolatilisation, the recommended CG parameter values can be reduced to  DBM ≈ 8.472 ⋅ 10 −5 ⋅ 1.083 (  ∕  ) MAM ∕  ≈ 0.992 ⋅ 1.117 (  ∕  ) . the above expressions, the thickness of the thermal boundary layer was extracted from the FR data, which is generally not available in (inherently under-resolved) EL simulations.Therefore, in practical EL-CG simulations the boundary layer thickness must first be estimated, either based on the heat flux for quiescent environments or with the aid of the Nusselt number for convective cases[72]