An investigation into the use of CFD to model the co-firing of Jatropha curcas seed cake with coal

ABSTRACT Jatropha curcas seed cake is a potential candidate for co-firing with coal. Combustion modelling using Ansys Fluent 14.0 was carried out to assess the combustion and co-firing characteristics of untorrefied and torrefied Jatropha curcas seed cake. The effect of torrefaction on the devolatilisation characteristics, flame properties and consequently NOx pollutant formation was established. Compared to the torrefied biomass, the untorrefied seed cake devolatilised earlier, had a more dispersed flame and higher NO formation. The higher reactivity of the biomass was shown to have a positive effect on the devolatilisation rate of the less reactive coal under co-firing simulations.


Introduction
Historically, progress in the field of combustion technologyparticularly with traditional fossil fuels such as coalhas been dependent on empirical data. This data can be sourced from either large-scale plants or laboratory/pilot-scale experiments. A limited number of measurements (such as those of effluent samples) can be obtained from large-scale plants, and they are also restrictive in terms of the range of process parameters. Smaller-scale experimental rigs solve these two issues, but have the disadvantage of questionable applicability when the results are extrapolated. Combustion modelling provides an intermediate stage in the process of designing and improving of combustion systems, whereby the two modes of empirical data can be linked (Eaton et al. 1999). For a model's predictions to be useful, they should be validated against experimental data. However, the usage of modelling allows confident predictions to be made using a small pool of experimental data, thus representing a cost saving (Barnes 2014).
CFD expertise in biomass and coal/biomass blends has been lagging relative to that of coal. Biomass combustion/ co-firing in itself has its share of issues, and so has its computational modelling. The thermochemical structure of biomass adds to the complexity of the model (Stroh et al. 2015). Also, the size and shape of the particles potentially has a more significant effect in the modelling of biomass, as opposed to coal; this is a consequence of internal heat transfer within the particles (Gubba et al. 2012(Gubba et al. , 2011. The biomass used in this study is derived from Jatropha curcas, a species of shrub which has been recognised as a potential source of renewable energy in the form of biodiesel; it is found natively and in a growing number of plantations in sub-tropical/tropical regions in Central America, Africa, India, China and Southeast Asia (Zhengguo Edrisi et al. 2015;Li et al. 2010). The plant is resistant to droughts and most pests (Openshaw 2000;Pandey et al. 2012). It grows and propagates rapidly (Sujatha, Reddy, and Mahasi 2008) and unlike palm oil and sugar cane, enjoys a lack of competition with the food industry due to its toxicity (Wever, Heeres, and Broekhuis 2012). The authors have previously explored the use of torrefaction as a pre-treatment process to improve thermochemical properties of Jatropha seed cake and enhance its viability as a feedstock for co-firing with coal (Madanayake et al. 2016) The overarching aim of this research work is to ascertain the suitability of Jatropha curcas seed cake for co-firing. Combustion of solid fuels is a complex process. This is even more so when two fuels, i.e. biomass and coal, are fired together. The fundamental differences between biomass and coal in properties such as reactivity, HHV, proximate and ultimate analysis, and particle characteristics can have a profound effect on their combustion characteristics. This is further compounded by the possibility of interaction effects between the two fuels. This study focused on the use of CFD as a tool to analyse the behaviour of Jatropha curcas seed cake when co-fired with a typical bituminous coal.

Biomass
bags before shipping to Nottingham, UK. The oil extraction was carried out by Statfold Seed Oil Ltd. uing twin-screw, cold-pressing expellers. 48.5 kg of solid residue (seed cake) was obtained from 55 kg of seeds, corresponding to an oil yield of 12%. This is lower than the oil yield of 28-40% reported in the literature (Ginwal, Rawat, and Srivastava 2004;Pandey et al. 2012;Singh et al. 2008). The seed cake was composed of two distinct componentshard, dark, rodlike structures in a soft, oily, loose soil-like matrix. The former was milled using a Retsch PM100 ball mill and then sieved to < 1 mm with a Retsch AS200 sieve shaker. The milled and sieved biomass were torrefied in a horizontal tube furnace under an inert (N 2 ) atmosphere at a temperature of 250°C for a period of 30 min.

Biomass characterisation
Thermogravimetric analysis (TGA) was carried out using a TA Instruments Q500 instrument. The TGA programme to determine the proximate composition was derived from British Standards guidelines (BS EN 14774-2:2009, BS EN 15148:2009, BS EN 14775:2009. The non-isothermal devolatilisation segment of the TGA run (constant heating rate of 20°C/min to 905°C under N 2 ) was also used to determine the devolatilisation kinetics. The devolatilisation reaction can be described by the following the expression (Gil et al. 2010): where x is the degree of conversion, t is time in s, f(x) is the reaction model, and k(T) is the rate constant which is a function of the absolute temperature T (K). The degree of conversion is defined as: where m is the mass of the sample at time t, and m 0 and m f are the mass at the beginning and end of the devolatilisation segment in question, respectively. The rate constant is governed by the Arrhenius law: where A is the pre-exponential factor, E is the activation energy and R is the universal gas constant of 8.314 J K −1 mol −1 . Hence, The first-order chemical reaction assumption is typically used in devolatilisation studies (Gil et al. 2010) (Gil et al. 2010): The method developed by Coats and Redfern (1964) is a widely-used model-fitting method used to calculate the kinetic parameters of devolatilisation. This method offers an approximate solution: It can be assumed that 2RT≪E (Damartzis et al. 2011). Hence, Equation 6 can be simplified to: ln A straight line is fitted to a plot of ln(-ln(1-x)/T 2 ) against 1/T. From the gradient and y-intercept of the fitted line, the activation energy and pre-exponential factor, respectively, can be calculated. The higher heating value (HHV) of the samples was measured using a IKA C5000 bomb calorimeter. The CHON content was measured using a Thermo Scientific FlashEA 1112 elemental analyser.

Drop tube furnace testing
The drop tube furnace (DTF) testing aimed to replicate the high heating rates found in industrial applications, for instance in pulverised fuel boilers (Ulloa et al. 2005). Figure 1 is a schematic representation of the DTF setup. The DTF consists of a vertically-aligned ceramic work tube which has an inner diameter of 50 mm and a length of 1.55 m. A water-cooled feeder probe with an internal diameter of 3 mm is present at the top of the work tube. At the base of the work tube is a water-cooled collector probe connected to a cyclone collector where the combusted material is deposited. The probe separation was 22 cm, which corresponds to a residence time of approximately 200 ms as calculated in the DTF's operating guidelines. The main gas inlet component flow rates (N 2 and O 2 in this case) and the feeder probe (carrier gas) gas flow rate (N 2 at 1 l/min) are controlled by rotameters. There are three heating units, all controlled via individual Eurotherm controllers linked to thermocouples present along the work tube: • Pre-heatera spirally-wound silicon carbide element around the feeder probe which brings the temperature of incoming gases close to the desired operating temperature. • Main heaterfour vertical silicon carbide elements around the work tube. • Trim heatertwo smaller U-shaped elements at the base of the work tube, which compensate for the water-cooling of the collector probe.
The main gas component flow rates and temperature settings are set as per the guidelines in the standard operating procedure, and depend on the desired operating temperature, residence time, and volumetric gas makeup. For the respective settings of 1300°C, 200 ms and N 2 (90%)/O 2 (10%) chosen for the runs, the N 2 and O 2 flow rates were set to 11.1 l/min and 1.23 l/min, respectively.
A weighed biomass sample (3.0 g) was manually fed through the open end of the feeder probe, which has a funnel attached to it. Care was taken to make sure that the feeding is slow, controlled and at an even rate. The time taken to feed the entire sample was measured, allowing the average feed rate to be calculated. After the run was complete, the combusted sample (char) was collected from the removable sample pot attached to the cyclone collector, and stored in a glass vial for subsequent analysis by TGA.
Samples of both untorrefied and torrefied material were combusted in the DTF. Triplicate runs were carried out in each case. After passing through the DTF, the proximate analysis of the samples was carried out using TGA. In order to determine the fraction of the original volatile matter content which had been lost (volatile matter conversion, VMC), it is necessary to know the total weight loss of the sample as it passes through the DTF. However, direct measurement of the weight loss in not possible since the collection efficiency of the DTF is low; particles can stick to the feeder and collector probes. Hence an indirect method known as the ash tracer method was used (Shuangning et al. 2006). This is widelyused in DTF testing, and relies on the assumption that the ash content does not undergo any change within the DTF. Using the ash tracer method, ΔW, the mass loss of the sample as a percentage of the original mass, is given by: where A and A' are the dry ash contents of the original sample and the char (after passing through the DTF), respectively. V c , the VM content of the char as a percentage of the original mass, was then calculated as follows: where V' is the dry VM content of the char measured by proximate analysis. The VMC is then calculated: where V is the dry VM content of the original sample measured by proximate analysis.

CFD meshing
A 2D axi-symmetric mesh was used for the CFD modelling. A 2D mesh in axi-symmetric mode cuts down on computational time compared to a full 3D mesh. Since the domain is truly axi-symmetric, 2D modelling of a "slice" of this cylinder is sufficient as this can then be rotated around the axis. Figure 1 shows the section of the DTF that is represented in the mesh, marked out by the dashed line. The initial geometry was created in Ansys DesignModeler. Ansys Meshing was then used to generate the mesh from the geometry. The "mapped meshing" and "fixed advanced sizing" features were used to obtain a grid of uniformly-sized quadrilateral faces (Figure 2). This method was chosen to minimise the orthogonal skew and maximise the orthogonal quality. These two mesh statistics are indicators of a highquality mesh (ANSYS Inc. 2011b). A summary of the boundary conditions is shown in Table 1, while the names of the boundaries are shown in Figure 2. A mesh-independence study was carried out to ascertain that the mesh size (how fine the mesh is) does not have a significant effect on the results, and to determine the optimum mesh size. By varying the maximum face size, six meshes were generated, as shown in Table 2. This table also shows the minimum orthogonal quality and the maximum orthogonal skewness of the generated meshes obtained from the "mesh statistics" feature, the former being closer to 1 and the latter being closer to 0 indicates a mesh with good quality (ANSYS Inc. 2011b).
An identical case was run on all six meshes, and three measured outputs (monitors) were chosen for comparisonthe temperature integrals at two x-positions along the mesh (0.15 m and 0.25 m), and the particle volatile matter conversion. The VMC is particularly important since it is used to validate the model against empirical data from the DTF. Also, there was some concern that the number of particle streams could have an effect on the VMC; since a surface injection is used, the number of particle streams is dependent on the number of faces on the injection surface, and this varies with the mesh density ( Table 2). The computational cost was measured using the average time per iteration reported by Fluent; a PC with an Intel Core-i5 3.1 GHz processor and 4GB of RAM was used for the CFD work.

CFD model setup
The combustion process inside the DTF was modelled as a two-phase flow. The continuous phase (or fluid phase) consists of the gases inside the work tube, i.e. the primary (N 2 ) and secondary (N 2 and O 2 mixture) inlet flows, volatiles evolved, reaction intermediates and products. The discrete phase (or solid phase) consists of the solid fuel (biomass and/or coal) particles. Fluent uses an Euler-Lagrangian approach, where the continuous phase is modelled using Navier-Stokes equations and the discrete phase is modelled by tracking the particles through the flow field. The system was modelled as steady state, since the DTF setup aims to achieve a steady, laminar-regime flow with minimum fluctuations.

Turbulence
The most commonly-used approach to modelling turbulence in co-firing applications is the k-ε model and its two variants, i.e. RNG and realisable (Tabet and Gökalp 2015;Wang and Yan 2008). The realisable k-ε model was implemented in this case, since this model has been favoured in recent studies due to its pragmatic compromise between accuracy and computational demand, and better solution convergence (Black et al. 2013;Zhengqi Li et al. 2009;Park et al. 2015;Rezeau et al. 2012). The underlying assumption of this model is that the turbulent viscosity is anisotropic or the ratio between the Reynolds stress and the mean rate of deformations differs depending on the direction or axis (Versteeg 2007).

Chemistry and turbulence-chemistry interaction
The chemistry model to be used is dependent on the turbulence-chemistry interaction model. The turbulence-chemistry interaction was modelled using the Eddy-Dissipation model (EDM) (Magnussen and Hjertager 1977). EDM is a commonly-used interaction models for biomass combustion studies, and is appropriate for one-step or two-step reaction chemistries. The EDM assumes that reaction rates are limited by turbulent mixing (as is the case with fast-burning fuels) and hence avoids computationally-expensive Arrhenius chemical kinetics calculations (Magnussen and Hjertager 1977).
The volumetric combustion chemistry, i.e. the reactions occurring in the gaseous phase, was modelled using the Species Transport model. This model involves solving conservation equations for each chemical species, i. Equation 8 is the general form of the conservation equation: where Y i is the local mass fraction of a species i, S i is its rate of addition from the dispersed phase and J i is its diffusion flux (modelled by Fick's Law). R i is the net rate of production of i by chemical reaction, and is calculated by the EDM (ANSYS Inc. 2011a). A two-step mechanism was used for the gaseous phase reaction, as follows: where the values of x, y and z are given in Table 4. Hence, the volatiles in the gaseous phase (which are represented by C x H y O z ) undergo oxidation via an intermediate CO species.

Discrete phase
The Discrete Phase Model (DPM) treats the fuel particles using a Lagrangian approach. The DPM calculates the particles' motion through the continuous phase, as well as their heat and mass exchange with the continuous phase.
The particle injections were defined as shown in Table 4. A surface injection from the fuel inlet was specified; a particle stream is emitted from each face on this surface. The velocity magnitude was calculated from the carrier gas flow rate (2.36 m/s), and the average flow rate was calculated from the time taken to feed a measured quantity of sample. The continuous and discrete phase calculations were coupled, i.e. the discrete phase affects the continuous phase (two-way coupling), with 10 continuous phase iterations being performed for every DPM iteration. This number was chosen by trial and error, based on solution convergence and stability.

Particle motion
The particle trajectory is governed by the following force balance: where u p is the particle component velocity, u is the fluid velocity, F D is the drag force, g x is the gravitational force component, ρ p is the particle density and ρ is the fluid density. The drag force is calculated by: where µ is the fluid viscosity, d p is the particle diameter, C D is the drag coefficient and Re is the Reynold's number. Using a spherical particle shape assumption might deviate from reality particularly for large biomass particles. A sphere gives a minimum in terms of the surface-area-to-volume ratio, which impacts significantly the particle motion (Ariyaratne et al. 2015). Thus, the drag coefficient is calculated using the nonspherical drag law developed by Haider and Levenspiel (1989). To account for the non-spherical nature of the particle, this model uses a shape factor, defined as the ratio of the surface area of a sphere having the same volume of the particle, to the actual surface area of the particle.

Particle heat and mass exchange
For the case of a combusting particle, the heat and mass exchange of the particle are governed by four laws, which are activated in succession: 2.5.5.1. Inert heating. This models the heating of the particle until it reaches its vaporisation temperature. A simple heat balance is used to calculate the particle temperature by taking into account the convective and radiative heat transfer.
2.5.5.2. Devolatilisation. The devolatilisation of the particle was modelled in this case by the single first-order reaction: where m p is the particle mass, m p,0 is the initial particle mass, and f v,0 is the initial volatile fraction. k is the rate constant, which is given by: where A is the pre-exponential factor, E is the activation energy, R is the universal gas constant, and T is the particle temperature. The kinetic parameters A and E were calculated from TGA runs and entered into the model. The values of A and E are given in Table 4.
2.5.5.3. Surface combustion. Surface combustion of the char is initiated after complete devolatilisation has taken place. A simple diffusion-limited model was used in this case, since accurate surface combustion kinetics were not available. This model makes the assumption that the surface combustion occurs at a rate determined by the diffusion of the gaseous oxidant to the surface of the particle (Gera et al. 2002). The diffusion-limited model has been used in previous biomass co-firing studies (Yin et al. 2010(Yin et al. , 2004. Complete devolatilisation was not expected in the present case due to the low particle residence time in the DTF. Hence, accurate modelling of the char combustion was not considered vital.
2.5.5.4. Inert cooling. This is similar to the inert heating law, in that the heat balance is used to calculate the particle temperature after complete devolatilisation.

Radiation
The discrete ordinates (DO) model was used to model the radiative heat transfer. This is a widely-used radiation model in biomass combustion and co-firing studies since it is more accurate and is suitable for all optic thickness (Álvarez et al. 2014;Bonefacic, Frankovic, and Kazagic 2015;Gubba et al. 2012;Stroh et al. 2015). The DO model solves the radiative transfer equation for a number of discrete solid angles (angular discretisation). Each solid angle is associated with a vector direction fixed in the global Cartesian coordinate system. The radiative transfer equation accounts for the absorption, emitting and scattering of radiation in the different direction vectors (ANSYS Inc. 2011a).

NOx
The formation of NO was modelled by using the Thermal model. The underlying mechanism of this model is temperature dependent and the model calculates the NO formed from the cracking of atmospheric N 2 at very high temperatures (ANSYS Inc. 2011a). It is based on the following equilibrium reactions: Other NO formation models that are available include the Prompt NO model and the Fuel NO model. However, these were not activated since they require inputs which were not known for this case and would introduce further uncertainty. For instance, the Fuel NO model requires the ratio to which the N content was split with the volatiles and the char.

Mesh independence
Excluding the coarsest mesh, a reasonable level of mesh-independence was seen. The maximum variation among the VMC and the temperature integrals were 4.2% and 1.4%, respectively. It is also noteworthy that although the variation of the monitors was small, the number of cells in the mesh (and hence the computational cost) increases exponentially as the maximum face size is decreased (Table 2). Taking into consideration these factors, mesh 3 was chosen to be used for the rest of the study.

Drop tube furnace testing
The VMCs for each DTF run, along with the calculation steps, are shown in Table 3. It was observed that the mean VMC of the torrefied material was higher than that of the untorrefied samples. However, it should be noted that there is a degree of overlap between the untorrefied and torrefied samples. A larger coefficient of variation (CV) was obtained in the untorrefied case. A significant degree of uncertainty would be expected considering the several variables that affect the operation of the DTF, in addition to the unique challenges posed by biomass testing such as the low collection efficiency and inhomogeneity within the biomass itself (especially in the untorrefied samples).
The empirical VMCs thus obtained would be used to validate the CFD model. However, taking the above-mentioned sources of error should be taken into account when interpreting the results from the CFD model.

Sensitivity analyses
The basic fuel properties such as the proximate analysis, ultimate analysis, HHV and devolatilisation kinetics were measured and input into the CFD model. However, there were several gaps in the required inputs where empirical measurements could not be made. Furthermore, the novelty of the use of Jatropha curcas seed cake as a fuel means that these unknown properties are unavailable in literature as well, as limited characterisation work has been done to date. The estimates used as inputs in these cases introduce uncertainties in the results; therefore, sensitivity analyses were conducted to quantify these effects. In each case, the uncertain property in question was varied while all other model parameters were fixed; the effect of this variation on the volatile matter conversion % (VMC) was studied.

Devolatilisation kinetics and particle size
The most fundamental uncertainties were with respect to the particle size and the devolatilisation kinetics. Although the samples used in the DTF had been sieved to less than 1 mm, an accurate distribution within this size fraction could not be determined. A tendency of the particles to agglomerate resulted in ineffective sieving at fine mesh sizes, and hence a sieve analysis could not be performed. This issue also meant that an optical particle analysis using a Retsch Camsizer (Retsch Technology 2016) also did not yield reliable data as clumping of particles occurred. Figure 3 is the particle size distribution for the untorrefied type A, obtained from the Camsizer. The normal-like distribution which would be expected is observed in the 0-300 µm range, with a centre of approximately 200 µm. A small peak is observed at approximately 400 µm, followed by a much larger peak centred at~600-700 µm. It can be postulated that these latter peaks are the result of the agglomeration of two or more 0-300 µm particles. The fact that the latter peaks occur approximately at multiples of 200 µm further supports this. Hence, although the Camsizer results are far from conclusive, they support the possibility that the majority of the individual particles are in the 200 ± 50 µm range.
The devolatilisation kinetics were measured using TGA data. This method is acceptable to compare the kinetics between different materials as highlighted in other studies utilising pinewood and coals of different rank (Jones et al. 2005), as well as willow coppice and eucalyptus and their torrefied counterparts (McNamee et al. 2015). However, data obtained from the TGA cannot be reliably applied to the DTF as there is a vast difference in the heating rate between the two instrumentsa heating rate of 15°C/min was used in the TGA runs, while that in the DTF is in the order of 10 4 -10 5°C /min. Hence, the approach taken was to carry out the CFD runs with the devolatilisation kinetics obtained from the TGA, along with constant particle diameter injections (this was specified in the DPM settings). The VMC calculated by the model was compared with that obtained from the DTF, and the impact of particle size on the VMC was gauged. Figure 4a,b show the effect of the particle diameter on the VMC for the untorrefied and torrefied cases, respectively. The VMC range obtained from the DTF is 59.7%-77.5% and 68.2%-77.9% for the untorrefied and torrefied samples, respectively (Table 3). For the untorrefied case, this corresponds to a particle diameter of < 40 µm. The corresponding diameter is higher for the torrefied case, at approximately 80-90 µm. A possible reason for this apparent increase is that particles have undergone agglomeration after torrefaction. Another factor is the possibility of particle swelling during the torrefaction process. This phenomenon has been observed before in a different type of biomass, miscanthus which underwent swelling at torrefaction temperatures up to 250°C (Dufour et al. 2012). However, the uncertainty of the TGA kinetics would be a definite and significant contributor to this difference in particle size.
Overall, the particle size range obtained is at the low end of the possible spectrum of 0 to 1000 µm. Although an accurate particle size distribution could not be measured empirically, as explained earlier in this section, the tentative empirical measurements indicate a mean particle size higher than 150 µm. The lower particle diameter predicted by the CFD runs is a reflection on the unreliable applicability of the TGA kinetics to the DTF due to the substantial difference in heating rates.
Single first order reaction (SFOR) devolatilisation kinetics parameters do not take into account heat transfer effects to and within the particle; they are based on the assumption that the particle is at a uniform temperature. TGA runs which are performed at relatively low heating rates allow larger particle sizes to be used while satisfying this assumption, but this assumption is no longer valid within the DTF. As the reaction kinetics are dependent on the heating rate, the TGA kinetics cannot be directly applied to the DTF to make accurate predictions of absolute values. However, considering the lack of high-heating rate kinetic data, the TGA kinetics were deemed suitable enough to make a comparative analysis between untorrefied biomass, torrefied biomass and coal as long as the kinetics for all the cases were derived under similar TGA conditions. Previous studies have also utilised TGA kinetics in combustion modelling cases (Stroh et al. 2015;Zhang et al. 2013).
For the rest of the validation CFD runs, it was assumed that the biomass injection can be modelled by a single particle diameter of 30 µm and 85 µm (based on the mean VMC obtained from the DTF) for the untorrefied and torrefied cases, respectively, while using the TGA kinetics

Density
The particle density affects its velocity through the fluid phase and hence its residence time in the high-temperature zone. This in turn would change the VMC.
Although the bulk density of Jatropha curcas has been reported in literature as 670 kg/m 3 (Kavalek et al. 2013), this is not indicative of the particle density which is what has to be input in to Fluent . The particle density would depend on the pressing process used to extract the oil out of the seeds. Since pelletisation would also require the material to be pressed, the typical density of biomass pellets was considered to be representative of the particle density of the seed cake. The density of biomass pellets can be in the range of 1000-1400 kg/m 3 (Bergström et al. 2008;Stelte et al. 2011); since the type A seed cake easily sinks in water, a density range of 1200-1500 kg/m 3 was used for the sensitivity analysis ( Figure 5). The VMC shows a mild sensitivity to the densitythe VMC decreased by 17% as the density was increased by 25% -in the torrefied case, with negligible change in the untorrefied case. Figure 6 illustrates the sensitivity of the VMC to the shape factor. The shape factor is a parameter used to account for the non-spherical nature of particles. A shape factor of 1 would be the simplest case, where the particle is modelled as a perfect sphere. A shape factor of 0.5 means that the particle has a surface area twice that of a sphere of the same volume. The shape factor has a negligible effect on the VMC in the case of the untorrefied biomass. For the torrefied case, the VMC increased by 8% when the shape factor is decreased by 60%, indicating a mild sensitivity. Coal particles can typically be modelled as simple perfect spheres, i.e. shape factor of 1 (Ma  ). However, the fibrous nature and irregular shape of biomass particles means that a lower shape factor would be more suitable. The shape factor of the untorrefied biomass was assumed to be 0.7, which was in line with typical values in the literature (Karampinis et al. 2012). Since it has been reported that the torrefaction process causes particles to become more spherical (Arias et al. 2008), a shape factor of 1 was assumed for the torrefied case.

Emissivity
Emissivity is a measure of the how effectively radiative heat is transferred to the particle. The baseline value of 0.9 is the default value used in Fluent. Although 0.9 is commonly used in biomass studies, values as low as 0.75 have also been utilised (Gupta, Yang, and Roy 2003;Ma et al. 2009;Pokothoane 2010). Figure 7 shows the variation of the VMC with the emissivity. There is a linear positive relationship between the emissivity and the VMC, as higher heat transfer to the particle would result in a higher particle temperature and hence more devolatilisation. The sensitivity is moderate in the case of the torrefied biomass -VMC drops by 11% when the emissivity is decreased by 22% from the baseline of 0.9 to 0.7. There is a relatively negligible effect in the untorrefied case.

Overview of sensitivity analyses
The particle size and kinetics were considered the most significant uncertainties. The TGA kinetics were fixed and the effect of varying the particle size was analysed. Based on the VMCs thus predicted from the model,and the empirical VMCs obtained from the DTF, representative particle diameters of 30 µm and 85 µm were chosen for the untorrefied and torrefied cases, respectively. These assumptions were considered sufficient for a comparative study without a detailed quantitative aspect. The density and shape factor had a mild effect on the VMC in the case of the torrefied biomass, while the emissivity had a more significant effect. However, the sensitivity of the untorrefied biomass to these three factors was relatively negligible. 3.3.6. Stochastic tracking Stochastic tracking (ST) uses the discrete random walk (DRW) model to model the effects of turbulence on the particle trajectories, i.e. the turbulent dispersion. ST was not performed on the sensitivity analyses and co-firing cases in order to reduce computational time and improve convergence. However, two simulations with ST enabled were carried out for the baseline untorrefied and torrefied cases, in order to determine if there is a change in solution. 10 stochastic "tries" were performed for each particle track that is injected at the fuel inlet, i.e. 10 particle trajectories per non-ST track. This results in a total of 20 different particle trajectories per symmetric half of the 2D mesh. Figure 8 shows the particle tracks coloured by particle temperature. In both cases, the majority of the particles (approximately 90%, as calculated from the particle mass flow rates at the inlet and outlet) exit the furnace through the collector probe. This indicates that the entrained flow conditions remain valid to a significant degree even with very small particles in the order of 20 µm.
The VMCs obtained with ST enabled were 67.2% and 71.9% for the untorrefied and torrefied cases, respectively. The corresponding values without ST were 74.8% and 77.4%. Despite the difference in VMCs between non-ST and ST cases, they all fall within the corresponding ranges of empirical VMCs obtained from the DTF runs. The particle size would have an effect on the turbulent dispersion modelled by ST. Since the particle sizes used in the CFD are not an accurate representation of the actual particle size used in the DTF runs (as they were adjusted to fit the TGA kinetics), it is possible that an additional degree of error is introduced by enabling ST. Furthermore, the present CFD study is comparative in nature. Considering these factors, using ST could not be justified for the subsequent CFD runs, and the simpler non-ST approach was chosen. Table 4 presents the material/particle properties that were used in the main study. Experimental measurements were made of the fundamental properties of the untorrefied and torrefied biomass, while the selection of uncertain or estimated properties was has been justified earlier. The properties of a typical coal (a high-volatile bituminous coal) was taken from literature (Biagini et al. 2002). This source was chosen since the devolatilisation kinetics were measured in similar conditions to those of the Jatropha seed cake (20°C/min in an N 2 flow). Figure 9c,d shows the temperature distribution for the untorrefied and torrefied cases, respectively. Both cases have a similar VMC (30 µm and 85 µm particle diameter for the untorrefied and torrefied cases, respectively). The "flame" can be considered the areas where the temperature is greater than the 1490-1580 K band (as the DTF temperature was set to 1573 K). For comparison, Figure 9a shows the temperature distribution with no combustion (inert particle) whereas Figure 9b shows the temperature distribution for coal. A distinct flame was not observed in the case of the coal (Figure 9b). This is reflected in the very low VMC of 7%. As expected from the kinetic parameters, the coal is significantly less reactive than both untorrefied and torrefied Jatropha seed cake. Meanwhile, the distinct shape of the flame is visible in both biomass cases (c & d). However, there is a difference seen in the size and position of the flame. In the untorrefied case, the flame starts closer to the top of the furnace, and is more widely distributed. The observations seen in Figure 9 concur with the findings of Li et al. (2013) who reported that the biomass flame volume increased as the torrefaction degree decreased or in other words, as the volatile content increased. This was attributed to the higher volatile content improving the diffusion of the combustible gases thereby resulting in a larger flame volume and diameter; the volatiles and their preliminary combustion products also enhanced gas phase diffusion from near-flame zone to a wider zone. This difference can also be attributed to the devolatilisation characteristics of the two materials. The untorrefied biomass has a lower activation energy. Hence, the reaction can commence at a lower particle temperature (closer to the top of the furnace). However, the pre-exponential factor is higher in the torrefied biomass by an order of magnitude. The pre-exponential factor determines the overall magnitude of the reaction rate. Hence, although the devolatilisation is delayed in the case of the torrefied biomass, it is more intense when it starts. This allows the torrefied case to achieve a VMC similar to that of the untorrefied biomass, even though the activation energy is higher. This effect is graphically illustrated in Figure 10a,b, which are contour plots of the devolatilisation rate for the untorrefied and torrefied cases, respectively. Although the onset of devolatilisation is delayed, a higher devolatilisation rate of 5.1x10 −9 kg/s is obtained in the torrefied case, compared to 3.6 × 10 9 kg/s for the untorrefied biomass.

Individual biomass/coal combustion
However, the delayed onset of devolatilisation in the torrefied case means that more volatiles combustion occurs within the collector probe. Since the collector probe is cooled, the temperature of the flame is controlled. This explains the lower-than-expected peak temperatures observed with the torrefied biomassthe peak temperatures recorded within the modelled volume was 1785°C and 1720°C for the untorrefied and torrefied cases, respectively. The lower temperatures obtained with the torrefied biomass are hence not a result of a lower energy release from the combustion reaction itselfthe heats of reaction calculated over the modelled volume was higher for the torrefied case at 127 W, compared to 104 W for the untorrefied biomass (as expected because of the torrefied biomass' higher HHV).
In a practical sense, this difference in devolatilisation characteristics and flame temperatures is important since it can affect the formation of NOx pollutants. The maximum NO formation rates were 1.3x10 −4 mol m −3 s −1 and 3.4 x10 −5 for the untorrefied and torrefied cases, respectively. The NO flow rates at the outlet were 4.9x10 −11 kg/s and 1.3x10 −11 kg/s for the untorrefied and torrefied cases, respectively. These results demonstrate the link between the temperature distribution and the NO formation whereby higher NO levels are formed with larger, higher-temperature flames, which is in agreement with the findings reported in the literature (Li et al. 2013).
Meanwhile, comparing with other fuel types, it has been published that torrefaction was able to reduce the NOx emissions of crop-derived wastes biomass such as corn-based distiller's dried grains with solubles, olive residue and rice husk (Ren et al. 2017). These reported results agree well with the outcomes of this work.

Coal co-firing
Coal co-firing cases were run with both untorrefied and torrefied seed cake. Two DPM injections were defined, one with the properties of the biomass and the other with those of the coal. The biomass % by mass was varied through 25%, 50% and 75% by changing the flow rates of each injection, while keeping the overall flow rate constant. Figure 11 shows the variation of the coal VMC with the co-firing ratio. As the biomass % is increased, the VMC of the coal increases as well. Although the increase is marginal in the case of the torrefied biomass, an increase by 51% is observed when the coal is co-fired with untorrefied biomass. The more reactive biomass undergoes combustion, and causes the temperature to rise (the flame region). As the coal particles enter this flame region, their temperature increases and this results in the coal undergoing faster devolatilisation. There was no regular trend in the change in VMC of the biomass under co-firing conditions, and the magnitude of the change was not significant. This effect is graphically illustrated in Figure 12, which shows the temperature distribution for the three co-firing ratios for the untorrefied case. As the biomass % increases, the flame region grows in size and intensity. Hence, the coal particles are subjected to both increasingly higher temperatures and longer residence times within the flame, thus experiencing significantly faster devolatilisation.
The temperature distributions for the torrefied biomass cofiring is shown in Figure 13. Here, it can be seen that due to the delayed devolatilisation of the torrefied biomass, the region of the flame is relatively small, even with the highest biomass fraction. Hence, the boost given to the coal devolatilisation is much less substantial than in the untorrefied cofiring case. This is reflected in the less significant increase in the coal VMC seen in Figure 11.
Since there was a significant difference in the volatile conversion between the biomass and coal particles, a direct comparison in the NO emission is not meaningful. Hence, the NO flow rate at the outlet was normalised by the total heat of reaction summed over the entire volume. This gives an indication of the efficiency of the co-firing in terms of NO reduction. Figure 14a,b show the normalised NO emissions when the coal is co-fired with the untorrefied and torrefied biomass, respectively. Although the uncertainties in the model do not justify detailed conclusions being drawn, the overall trends indicate that co-firing reduces the NO emissions when compared to combusting only coal. The reduction in NO emissions when co-firing agrees with the results of another CFD study which reported that NO reduced proportionally with the biomass mass fraction in the mixed fuel during co-combustion of wheat straw with bituminuous coal (Ghenai and Janajreh 2010). Likewise, co-firing cardoon with lignite was demonstrated to reduce NOx emissions by up to 10% in a numerical study by Karampinis et al. (2012). Li et al. (2012) also reported that NOx emissions decreased significantly with increasing biomass substitution in coal co-firing. Nevertheless, it should be noted that only a small fraction of the coal has reacted in this work. If the coal and biomass combustion continues to a greater extent (for instance with a longer residence time), the results from the model are liable to change.

Limitations and further work
The limitations of this work revolve around three main areas as follows: 1. sole focus on the devolatilisation stage of combustion; 2. application of low-heating rate TGA kinetics to high-heating rate DTF combustion; and 3. assumptions made in certain model parameters such as particle size distribution. As elaborated upon earlier, complete devolatilisation had to occur before surface combustion of char was initiated. Complete devolatilisation was not expected to occur during the short residence times experienced in the DTF. The authors were aware of this experimental limitation at the outset, and this study was designed to be the first phase of a larger scope of work. The devolatilisation characteristics were deemed an important (albeit not the only) primary indicator of how the overall combustion behaviour differs between coal and torrefied and untorrefied biomass. Another limitation stemming from the short residence times in the DTF is the lack of direct applicability of TGAderived devolatilisation kinetics. The single first order reaction kinetics derived from the TGA data do not consider the heat transfer effects to and within the particle. These heat transfer effects become critical considering the substantial difference in heating rates between the TGA and DTF -15°C /min in the TGA versus 104-105°C/min in the DTF. However, since this study was comparative in naturewith a more qualitative than quantitative biasthis limitation was recognised and justified. Furthermore, this approach has been utilised in other independent works as well (Stroh et al. 2015;Zhang et al. 2013).
The accuracy of the model depends on the model inputs; although the fundamental fuel characteristics could be measured, there were several uncertainties in the model inputs resulting from limited experimental data and the lack of such data in the literature (as Jatropha curcas seed cake is still a relatively novel source of energy). The model was adapted to fit the DTF data. Again, it was recognized that these shortcomings are justified in the light of the objective of this study being an analysis of the effect of torrefaction and co-firing from a qualitative/semi-quantitative standpoint, and to establish an understanding of the mechanisms. To this end, sensitivity analyses were also carried out to ascertain how significant the effect of the assumed parameters had on the results.
The combustion modelling can be strengthened by gathering further empirical characterisation data for the Jatropha curcas seed cake to reduce the uncertainties in the model inputs. These include data such as high-heating rate devolatilisation kinetics, accurate particle size distributions, and nitrogen partitioning between volatiles and char. Furthermore, longer residence times can be modelled and validated using DTF to ensure more complete devolatilisation as well as char oxidation occurring. A stronger, more complete model would enable more meaningful quantitative data to be extracted and provide a more complete understanding of the co-firing process.

Concluding remarks
A CFD model was developed to investigate the combustion behaviour of the Jatropha curcas seed cake in a DTF, where conditions which are close to real-world applications, i.e. high temperature (up to 1300°C) and heating rate (in the order of 10 4 -10 5°C /min). However, a model's accuracy is dependent on the accuracy of the inputs given to the model. Several properties of the biomass had to be estimated in the light of limited experimental resources and a dearth of detailed characterisation of Jatropha curcas seed cake in the literature. The most significant uncertainties were the lack of accurate particle size data, and the devolatilisation kinetics obtained from low-heating rate TGA. An approach was adopted where the TGA kinetics were used but the particle size used in the model was reduced until the volatile conversion predicted by the model matched that measured by from the DTF runs. Sensitivity analyses were also carried out to ascertain the potential inaccuracy introduced by other estimated inputs. Significant differences in the temperature distribution were observed between the untorrefied and torrefied combustion. The untorrefied biomass resulted in a more dispersed flame. This was reflected in the significant differences in NO emissions that were observed, as NO production is heavily temperaturedependent at the high temperatures present in the furnace.
The co-firing cases demonstrated an interaction between the coal and biomass in terms of reactivity. The higher reactivity of the biomass has the effect of boosting the reaction rate of the coal particles as well due to the extra heat generated by the biomass combustion. This effect was most profoundly seen with the untorrefied biomass.
In the light of the model input uncertainties, detailed quantitative analysis of the combustion modelling results cannot be justified. However, this study has highlighted some important distinctions between coal and biomass combustion, untorrefied and torrefied biomass combustion, and the interaction effects between the two fuels that occur during co-firing. The importance of using CFD as a tool to investigate and optimise the cofiring of torrefied Jatropha curcas seed cake was established.
From this baseline, the model can be further developed in the future to include more accurate model inputs (for example, using kinetic data obtained under high heating rates, use of the Fuel NO and Prompt NO mechanisms as well) and longer residence times so that the coal and biomass undergo a greater degree of conversion.