Biomass pyrolysis TGA assessment with an international round robin

The large variations found in literature for the activation energy values of main biomass compounds (cellulose, hemicellulose and lignin) in pyrolysis TGA raise concerns regarding the reliability of both the experimental and the modelling side of the performed works. In this work, an international round robin has been conducted by 7 partners who performed TGA pyrolysis experiments of pure cellulose and beech wood at several heating rates. Deviations of around 20 – 30 kJ/mol were obtained in the activation energies of cellulose, hemicellulose and conversions up to 0.9 with beech wood when considering all experiments. The following method was employed to derive reliable kinetics: to first ensure that pure cellulose pyrolysis experiments from literature can be ac- curately reproduced, and then to conduct experiments at different heating rates and evaluate them with isoconversional methods to detect experiments that are outliers and to validate the reliability of the derived kinetics and employed reaction models with a fitting routine. The deviations in the activation energy values for the cases that followed this method, after disregarding other cases, were of 10 kJ/mol or lower, except for lignin and very high conversions. This method is therefore proposed in order to improve the consistency of data acquisition and kinetic analysis of TGA for biomass pyrolysis in literature, reducing the reported variability.


Introduction
Biomass is currently the main renewable energy source and it is expected to play a key role to reach the target formulated in the special IPCC report to limit global warming to 1.5°C [1,2]. Biomass combustion for bioheat production is a consolidated technology and it is currently the main bioenergy use. Besides, biomass has the potential to play a more relevant role in the production of power, liquid fuels or chemicals. Thermo-chemical processes applied to lignocellulosic biomass are very promising for this purpose. Pyrolysis itself is a promising conversion process that can be used to produce liquid bio-oil [3,4] and biochar [5], and is a main intrinsic sub-process in other thermo-chemical conversion processes such as gasification or combustion.
Biomass mass loss behaviour is commonly determined by thermogravimetric analysis (TGA). TGA experiments with small samples and low heating rates can be conducted in a pure kinetic regime, i.e. without heat and mass transport limitations. For cellulose pyrolysis, it was concluded by Antal et al. [6] that mass loss at low heating rates can be described with a single first order reaction with a high activation energy (191-253 kJ/mol). Gronli et al. [7] conducted a round robin at 8 European labs with commercial cellulose Avicel PH-105, showing some although limited scattering in the results, with a temperature of peak conversion at 327 ± 5°C and activation energy of 244 ± 10 kJ/mol at 5°C/min.
For lignocellulosic biomass pyrolysis, the derivation of kinetics is more challenging. Mass loss description with a single reaction is inaccurate, and the most common approach is to employ three reactions with a parallel reaction scheme, where each component represents cellulose, hemicellulose and lignin, respectively [8]. The main peak in conversion rate corresponds to cellulose, the shoulder at lower temperatures to hemicellulose and the tail at high temperatures to lignin. The kinetic parameters for each component are determined with model-based (model-fitting) methods, where the reaction model (commonly first and n th order) are selected before the fitting. Very different activation energies are however reported for each component, and especially for lignin [9]. This leads to a concern about the reliability of TGA data [10]. It has been pointed out that experiments should be conducted at different heating rates in order to minimize the influence of compensation effects [11,12]. Besides, it has been reported that the selection of higher order reaction models for lignin leads to better predictions [9,13], as it was also the case for the use of a distributed activation energy model (DAEM) [14].
Isoconversional (model-free) methods are also applied for deriving kinetics removing the need for a reaction model assumption. Integral isoconversional methods as Kissinger-Akahira-Sunose (KAS), Flynn-Wall-Ozawa (FWO) or Vyazovkin can provide the activation energies along conversion when experiments at different heating rates are conducted. The obtained kinetic data with these methods cannot be straightforwardly employed in a particle and reactor model, as the kinetic parameters change along conversion, but their complementary application with fitting (modelbased) methods can increase the consistency of the obtained kinetic data though fitting, especially the selected reaction models, as suggested by Khawam and Flanagan [15] and applied by Anca-Couce et al. [9] for biomass pyrolysis. Anca-Couce et al. [16] recommended that for a reliable determination of biomass pyrolysis kinetics the reference experiments with pure commercial cellulose from Gronli et al [7] should be first reproduced.
Subsequently, experiments with biomass should be conducted at different heating rates and assessed with isoconversional methods, in order to verify the reliability of the experiments as well as the employed reaction models and obtained activation energies in a fitting routine.
In a thorough review conducted by Anca-Couce [16] in 2016, where concern about the variation of kinetics values in literature was raised, activation energy values for hemicellulose, cellulose and lignin pseudo-components in the parallel reaction scheme were reported from biomass pyrolysis experiments conducted at several heating rates. Cellulose values ranged from 190 kJ/mol to 250 kJ/mol, while the respective ranges for hemicellulose and lignin were of 70 -215 kJ/mol and 20 -190 kJ/mol, respectively. This review was updated in the present work to include investigations conducted after 2016. In Fig. 1, the activation energies are reported from the kinetic analyses conducted for pure biomass components (cellulose, hemicellulose and lignin) using various methods [17][18][19][20][21][22][23][24][25][26][27][28][29]. A very wide range of activation energies has been reported for each component, namely 114 -288 kJ/mol for cellulose, 34 -179 kJ/mol for hemicellulose and 7 -226 kJ/mol for lignin. Fig. 2 includes the results from kinetic analyses conducted for various biomass species employing three (or four in some cases) pseudo-components [30][31][32][33][34][35][36][37][38][39][40][41][42][43][44][45][46][47]. Very large variations are also present for each pseudo-component. The ranges for the activation energies of each pseudo-component were 72 -244 kJ/mol for cellulose, 58 -200 kJ/ mol for hemicellulose and 16 -205 kJ/mol for lignin. The initial sample mass employed in these works was on average of 9.2 ± 5.0 mg, which is a relatively high value that can cause thermal lag in certain cases [6,7]. The initial sample mass was however not directly correlated with the obtained activation energies. These large variations observed for the activation energy values of pure and pseudo components in very recent works raise concerns regarding the reliability of both the experimental and the modelling side of the analyses performed in the literature.
The objective of this work is to investigate the reproducibility of TGA biomass pyrolysis experiments and the deviations that can be expected when mass loss kinetics are derived from the same sample with different TGA devices. An international round robin has been conducted for this purpose with 7 European partners. Experiments have been first conducted with commercial cellulose to reproduce the results from a previous round robin conducted by Gronli et al. [7]. Subsequently, experiments were conducted at different heating rates with homogenized beech wood and the results were analysed with fitting and isoconversional methods for kinetics derivation. The reported activation energies for woody biomass pyrolysis in literature vary in a broad range and this variability is attributed to a certain extent to the inherent variability in biomass composition. This uncertainty in composition is removed in this study, which is focused on investigating the deviations that can be expected when TGA pyrolysis experiments are conducted with different devices by experienced users and the kinetics are derived using reliable methods. The results of this work will lead to relevant conclusions regarding the accuracy that can be expected in determining biomass pyrolysis mass loss kinetics.
A. Anca-Couce, et al. Fuel 276 (2020) 118002 with cellulose at a constant heating rate of 5 K/min and with beech wood at 4 constant heating rates in the range from 1 to 20 K/min (1, 5, 10 and 20 K/min for all partners, except #5 which conducted experiments at 2 K/min instead of 1 K/min). The samples were centrally distributed to all the partners involved in the round robin. Commercial Avicel® PH-101 cellulose (CAS Number: 9004-34-6) was purchased from Sigma-Aldrich. Beech wood chips were homogenized according to ISO 14,780 and milled to particle sizes < 200 μm. The same homogenized sample was used by every participant in the round robin. The proximate, ultimate and elemental analysis of beech is presented in Table 2.
The partners were requested to conduct experiments with their usual procedure, employing an initial mass sample as low as possible, ideally of 3 mg (following the ASTM E1641 -16 standard [48]). Some partners employed a higher initial mass, as shown in Table 1, following their commonly employed methods. In case #3 an open weighing system was employed, which inherently limits the precision of the system, and a compromise has to be found between background noise and measurement signal. In cases #6 and #7 it was employed the minimum quantity that is required to cover the whole surface of the crucible, in order to maximize the exposed surface with an even distribution of the sample. The ASTM E1641 -16 standard is similar to ISO 11,358 -2 but differs mainly in its mathematical treatment and it was employed as a general guideline in the present study. The TGA experiments were conducted from room temperature up to 500°C, however the results were analysed and presented only for the interval between 150°C and 500°C. In order to ensure complete moisture evaporation until the aforementioned lower limit, a holding time of 10 -15 min was employed from the partners at temperatures around 110°C. The char yield is reported at 500°C, considering the initial mass  [30][31][32][33][34][35][36][37][38][39][40][41][42][43][44][45][46][47]. The values denoted with a * are mean values obtained from experiments performed under various heating rates or with different biomass species. For cases with 4 components, the 4th is either extractives [45,46] or a second hemicellulose component [44]. Boxplots indicate the median and interquartile ranges.  as the one at a temperature of 150°C.

Kinetic analysis
The kinetic analysis is done is the range from 150 to 500°C. Conversion (α) is defined in Equation (1) and it is calculated as a function of the current mass (m) as well as the initial and final masses at 150 and 500°C, respectively. The raw data is treated to provide 200 points for conversion in this range for further kinetic analysis. Conversion is employed instead of mass loss to do not account for the differences in char yield, as the current study focuses on mass loss kinetics without considering differences in product composition. These conversion values for all data sets are provided in the supplementary information. Besides, standard deviations in this manuscript are calculated with the "n -1" method to estimate them in a conservative way.
Kinetics are calculated with the model fitting method considering n th order reactions to determinate the reaction rate (dα/dt), as shown in Equation (2), where A is the pre-exponential factor, E is the activation energy, R the gas constant, T the temperature and n the reaction order. Cellulose pyrolysis is described with one first order reaction (n = 1). Beech wood pyrolysis is described with 3 reactions, including one first order reaction representing cellulose and two n th order reactions representing hemicellulose and lignin. The model fitting routine derives the kinetic parameters from the derivative thermogravimetry (DTG) curves, where dα/dt is plotted as a function of temperature, and employs a least minimum squares method (nlinfit routine from Matlab [49]). The error in the fit for each experiment is normalized by the maximum of the DTG curve, as described elsewhere [9,50].
Besides, the isoconversional KAS method has been applied at defined conversion intervals. The temperatures T α,i at which a certain conversion α is achieved at each heating rate (ΔT/Δt) i are first calculated. The activation energy at that conversion E α is then calculated based on Equation (3). For that, an Arrhenius plot of the left side of Equation (3) versus 1/T α,i is derived, and the activation energy is calculated from its slope which is equal to -E α /R. It was previously shown that other integral isoconversional methods as FWO or Vyazovkin lead to very similar results for biomass pyrolysis [9].

Cellulose pyrolysis
The thermogravimetry (TG) and derivative thermogravimetry (DTG) results of cellulose pyrolysis experiments at 5 K/min are shown in Fig. 3. For TG the conversion α is shown in the range 150 -500°C and for DTG the reaction rate dα/dt (calculated from the percentage of conversion) is shown for the same temperature interval. The averaged heating rate, char yield and temperature of the peak in reaction rate (DTG) are shown in Table 3. The shape of the curves is similar for all cases, but some deviations are present. The maximum reaction rate was achieved on average at a temperature of 328.3 ± 9.2°C, while it was of 327 ± 5°C in a previous round robin with the same commercial cellulose, conducted by Gronli et al. [7]. The obtained averaged heating rates (5.12 ± 0.26 K/min) are close to the target of 5 K/min.
Regarding the char yields, significantly higher deviations were obtained among the partners. The obtained values from cases #1, #2, #4 and #5 (5.1 ± 2.5%), employing initial mass samples of 3 mg, are similar to the ones of the round robin by Gronli et al. [7] (7.2 ± 2.4% with an initial mass sample of 4.1 ± 1.3 mg). However, the obtained  char yields are higher for cases #6, #7 and especially for the case #3 (9.5%, 12.2%, 21.4%, respectively), which leads to higher values of the mean char yield and standard deviation. These 3 cases were also the ones employing a higher initial mass for the experiments (10, 12 and 100 mg for cases #6, #7 and #3, respectively). The correlation shown in Fig. 4 between initial mass sample and char yield may be due to for bigger samples, where the retention time of the volatiles in contact with the sample is increased, secondary charring reactions are promoted leading to higher char yields [16,51]. The char yield should also depend on the configuration of sample holder and inert flow.
The DTG experiments were modelled with a first order reaction. The obtained kinetic parameters and error in the fit for all cases are shown in Table 3. Moreover, the DTG curves for experiments and model are shown in Fig. 5 for the cases #1 (left) and #3 (right), considered as representative cases. The obtained activation energies and (log 10 ) preexponential factors were on average of 222.0 ± 22.1 kJ/mol and 17.03 ± 1.81 s −1 , respectively. In the previous round robin values of 244 ± 10 kJ/mol and 19 ± 1.1 s −1 were respectively achieved [7]. The averaged error in the fitting was of 2.9 ± 1.2%. For most of the cases the error was very low, around 2%, as shown in Fig. 5 -left for the case #1. This error took mainly place at temperatures around 350°C due to the tail of the DTG curve, which cannot be modelled with a single reaction. This tail was more pronounced for the cases #3 and #7, which leads to a higher error (see the DTG curve in Fig. 5 -right for case #3) and can be attributed to a more relevant char devolatilization in these cases with a higher char yield. The obtained char yields for cases #7 and specially #3 were higher than the range reported in the round robin by Gronli et al. [7], as seen in Fig. 4. In this previous round robin, the error in the fit was generally lower because it was calculated for the TG curves and in a narrower temperature range [7].
The obtained results show a good reproducibility of the previously conducted round robin by Gronli et al. [7]. The peak of reaction rate is at the same temperature, although with a higher variability. The case number #5 is the only one clearly outside the general trend, with the peak being observed at 15°C less than the average. The obtained activation energies are in a similar range, with the exception of case #5 and to a lower extent case #4, where lower values are obtained. It is also remarkable that for cases #3 and #7 a significant conversion takes place at temperatures higher than 370°C, which leads to a higher error in the fitting as previously commented. This may be related to the higher values of char yield observed for these two partners (see Table 3), probably caused by the higher initial cellulose mass (100 and 12 mg for the case #3 and #7, respectively) leading to char devolatilization in this region.
It can be concluded from the cellulose results that previous results from literature could be well reproduced for most of the cases, although with a slightly higher variability. The main discrepancies were present in the char yield, which is significantly affected by the initial mass sample. Mass loss in cellulose pyrolysis is already well understood and the current results show that deviations are to be expected among experiments conducted from different partners, but they are limited, especially if small initial sample masses are employed. The next section will analyse the results of the round robin with beech wood at different heating rates, which is the main novelty of this work.

Beech wood pyrolysis
The TG and DTG curves for beech wood pyrolysis at 1, 5, 10 and 20 K/min are shown in Fig. 6 and Fig. 7, respectively. The averaged heating rate, temperature of the peak and temperature of the shoulder in DTG curves and char yields are shown in Table 4. The temperature of the shoulder is calculated from the local minimum of the derivative of the DTG curve. For the experiments at 1 K/min, the data of #4 is analysed only in the range from 150 to 450°C, as the heating rate of the experiment significantly increases after 450°C, and data of #5 is not present as the experiment was conducted at 2 K/min.
As for cellulose, the shape of the curves is similar for all cases, but some deviations are present. The obtained peak and shoulder Fig. 4. Char yield as a function of initial sample mass for commercial Avicel® PH-101 cellulose pyrolysis at a heating rate of 5 K/min for partners #1 to #7.   [8]). The obtained deviations among experiments with different devices are of a similar order as in the previously presented cellulose pyrolysis round robin. Besides, the deviations from average presented the same tendency for most cases, e.g. peak of DTG at lower temperatures than average for #1, #3 and #5 and higher temperatures than average for #2 and #6.
The obtained averaged heating rates were very close to the targets for most of the cases, although some deviations were present. In case #4, the heating rate increased at the end of the experiment. The experiment at 1 K/min was analyzed only until 450°C, as previously stated, because at higher temperatures the heating rate became too high. However, the most critical experiment of #4 is the one at 20 K/ min, as the increase in heating rate takes place already when the conversion is not high (i.e. a significant fraction of biomass has not yet pyrolyzed, as opposite to other experiments when the increase in heating rate takes place when conversion is already very high). In the experiments of #3 the heating rate is higher at the begin of the experiment for the cases at 10 and 20 K/min. In all other cases the deviations were minor.
Regarding the char yields, the obtained deviations are lower than for cellulose. The highest char yields were obtained for cases #3, with the highest initial sample mass, and #5. The results are surprising for case #5 as the initial sample mass was of 3 mg. Cases #6 and #7, with an initial sample mass of 10 -12 mg had generally higher char yields than the other cases with an initial sample mass of 3 mg (#1, #2 and #4). Besides, the char yield was slightly lower at the lowest heating rate in some cases, which points out a higher uncertainty in the determination at these very low heating rates.

Isoconversional KAS method results
The isoconversional KAS method has been applied to determine the activation energies along conversion for all cases. The Arrhenius plots are shown in Fig. 8. For a certain conversion, one point is obtained from each experiment at a different heating rate. Therefore, 4 points are presented for each conversion, since experiments were conducted at 4 heating rates. The slope of the Arrhenius plots at a certain conversion determines the activation energy. The quality of the linear fit is good for most of the cases. There are however two exceptions: the experiment from #2 at 1 K/min for high conversions (see at α = 0.9 in Fig. 8) and the experiment from #5 at 20 K/min.
The obtained activation energies from the KAS method and Rsquared values are presented in Table 5 and Table 6, respectively. For #2 and #5 the values are reported for all heating rates and also without the experiments considered as outliers (at 1 K/min for #2 and at 20 K/ min for #5). Including all experiments, the obtained R-squared values were low at all conversions for #5 (values around 0.9) and very low at conversions higher than 0.8 for #2. These 2 outliers were therefore easily detected with the isoconversional KAS method, looking at the Rsquared values as well as the Arrhenius plots in Fig. 8. Without the outliers, the obtained R-squared values are very high for all cases (see Table 6), obtaining values above 0.99 for the seven data sets for conversions between 0.1 and 0.85 (except for #5 without outliers at conversions of 0.1 and 0.15, with R-squared values of 0.97 and 0.98). The R-squared values are a bit worse for some cases at low conversions (α = 0.05) and they are much worse at high conversions (especially at α = 0.95).
The obtained activation energies with the KAS method along conversion are presented in Fig. 9. The values are shown on the top for all experiments and on the bottom for the experiments without outliers in cases #2 and #5. The cases with outliers had (together with case #3) higher activation energies than the other 4 cases. Including them would lead to an erroneous higher determination of activation energies with this method. Excluding the experiments considered as outliers led to lower activation energies in these cases, following the general trend for case #2 and with lower values for case #5. The predicted activation energies on average (without outliers) are a bit above 150 kJ/mol at low conversions and they increase along conversion for all cases, achieving values above 170 kJ/mol at a conversion of 0.4 and up to around 180 kJ/mol at a conversion of 0.85 (see Table 5). The obtained standard deviations until this conversion are in the range 20 -25 kJ/ mol (25 -30 for all experiments including outliers). At yet higher conversions, the activation energies increase significantly, but the obtained deviations among all cases are much higher. The obtained results from the KAS method will be critically discussed in the next subsection.

Model fitting results
The results of the model fitting method are shown in Table 7. Three components are considered for the fit: cellulose, hemicellulose and lignin. In Fig. 10, the DTG curves for experiments and model fitting are shown for case #1 and for all heating rates. In Fig. 11, the DTG curves for other cases at 5 K/min are shown. The main component is cellulose (c ≈ 0.57), which has been modelled with a first order reaction and describes the main peak of the DTG curve. Hemicellulose (c ≈ 0.33) and lignin (c ≈ 0.1) describe respectively the shoulder at lower temperatures and the tail at high temperatures with an n th order reaction. The error in the fitting is generally low, with values around 2% for all cases. Similar results are obtained for all cases and the only exception is case #7, where the best fit was obtained when the first order reaction describes lignin and cellulose is described with a reaction with an order lower than one. The error in the fitting was also a bit higher in this case (3.1%). It should be mentioned that the same initial values were employed in the fitting routine for the fitted parameters in all cases, and it was checked that with different values a better fit was not obtained.
The obtained activation energies for cellulose with the fitting routine are around 180 kJ/mol. For hemicellulose lower activation energies are obtained, around 150 kJ/mol, with an order of reaction close to two. The obtained standard deviations from all cases are moderate, especially when not including case #7 (15 and 8 kJ/mol for cellulose and hemicellulose, respectively). For lignin the obtained activation energies (over 300 kJ/mol) and reaction orders (around 7) are high, and the variability is much higher than for the other compounds.

Comparison of isoconversional KAS and model fitting results
In Table 8 the activation energy values obtained from the fitting and the isoconversional KAS method are compared. For cellulose, the derived activation energy with the fitting method is compared to the value from the KAS method obtained at the conversion level with the maximum reaction rate in the DTG curve, which is around 0.7. A very good agreement is obtained between both activation energies and the highest difference is 5 kJ/mol for case #4. For hemicellulose, the shoulder of the DTG curve is selected as a representative point for the comparison. The shoulder is observed at conversions around 0.25 and at this point, where the curvature of the DTG curves changes, the cellulose and hemicellulose compounds have a similar relevance (see Fig. 10 and Fig. 11). Therefore, the activation energy from the KAS method at the conversion where the shoulder takes place is compared to the averaged activation energy of cellulose and hemicellulose from the fitting method. In this case a very good agreement is as well achieved, and the highest deviation is 7 kJ/mol for case #6. Therefore, the isoconversional KAS method validates the obtained activation energies from the fitting method obtained for cellulose and hemicellulose. For lignin, the comparison is more challenging. The lignin compound in the fitting method describes the tail at high temperatures in the DTG curve. This tail starts at conversions around 0.9. The activation energies from the KAS method at this stage increase for all cases, but also the R-squared values are lower. At a conversion of 0.95 the R-squared values are already quite low for most of the cases (see Table 6). For the comparison with the fitting method, the activation energies of the KAS method at a conversion of 0.92 are selected. At this stage the R-squared values are generally acceptable (see Table 8). From the comparison it can be seen that the results from the isoconversional method suggest that there is an increase in activation energy during the last stage of conversion, which was also obtained from the fitting method for lignin, although discrepancies in the obtained values are present.

Discussion
The conduction of TGA for biomass pyrolysis and its evaluation is inherently difficult. Cellulose mass loss at 5 K/min can be well described with a single first order reaction with a high activation energy (above 200 kJ/mol). The results of this round robin show however that standard deviations of 9°C in the position of the maximum in the DTG curve at 5 K/min and 22 kJ/mol in activation energy were obtained. Such deviations are therefore to be expected when conducting experiments with different devices.
For beech wood pyrolysis, the obtained standard deviations in the positions of the peak and shoulder of the DTG curves (around 10°C at 5 K/min) are of a similar order as for pure cellulose. These deviations arise therefore only from the use of different devices with different operation modes. The variability in TGA data for biomass pyrolysis in literature is commonly attributed to the inherent heterogeneity in biomass composition. This work shows that the use of different devices and methods with the same samples further adds variability in the results.
The determination of mass loss kinetics from wood pyrolysis is more challenging than for pure cellulose. As already discussed in the introduction, experiments at different heating rates are required. Isoconversional methods are a suitable tool to assess the reliability of wood pyrolysis TGA experiments. The integral isoconversional KAS method has been applied in this work and it was shown that outliers can be detected when low R-squared values are present. Two outliers (at 1 K/min for case #2 and at 20 K/min for case #5) were detected from the 7 data sets with different TGA devices evaluated and therefore removed from further analysis. High R-squared values, over 0.99 during most of the conversions (except at the very begin and specially at the end of the conversions) are to be expected from reliable experiments at different heating rates, as it is the case for the experiments without outliers in this work.
The kinetics were determined using a fitting method which employs three parallel reactions (cellulose, hemicellulose and lignin), and validated with the isoconversional KAS method. A variable reaction order was assumed for hemicellulose and lignin, and a first order for cellulose. The obtained activation energy for cellulose is high (around 180 kJ/mol), but lower than for pure cellulose. The activation energy for hemicellulose (around 150 kJ/mol) is lower than for cellulose. In both cases, the values obtained from the fitting are supported by the isoconversional KAS method with a good accuracy. For lignin, very high activation energies are obtained with both fitting and isoconversional methods (above 200 kJ/mol), but the deviations for different cases are large.
The obtained deviations in activation energies for cellulose and hemicellulose with the fitting method and for the KAS method for conversions up to 0.9 are around 20 kJ/mol (without outliers). These deviations are of similar order as the results with pure cellulose. A closer inspection of the results shows that the deviations are especially higher for case #7, obtaining much lower activation energies than the average and also significant differences in reaction orders. For the other Table 4 Average heating rate (K/min or°C/min -equivalent units), temperature of the peak and shoulder in DTG curve (°C) and char yield (%) for beech wood pyrolysis experiments at heating rates of 1, 5, 10 and 20 K/min (* result at 2 K/ min, not considered for the mean). Partner  cases, the differences are smaller, but it can yet be detected that lower activation energies are obtained for case #5 and higher for case #3. These 3 cases (#3, #5 and #7) were also the ones for which higher discrepancies were found for the cellulose experiments. The conversion took place at significantly lower temperatures for case #5, with the peak in the DTG curve at 15°C less than the average. This behaviour was also seen in the beech wood experiments for this case (see Fig. 6 and Fig. 7). For cases #3 and #7, a remarkably higher char yield was obtained for cellulose, as well as a higher error in the fitting due to char devolatilization at temperatures higher than 370°C. It is not clear if this influences also significantly the beech wood experiments, but anyhow in both cases the deviations were higher for cellulose and also the highest initial mass samples were employed. If these cases for which deviations were already present for cellulose are not considered for the beech wood analysis, the obtained standard deviations are much lower. The averaged kinetic parameters for cases #1, #2 (without 1 K/ min), #4 and #6 are shown in Table 9 and the obtained activation energies with the KAS method in Fig. 12. These selected cases are the ones which could reproduce with a good accuracy the experimental results with pure cellulose and for case #2 the experiment at 1 K/min was not considered after an analysis of the results applying the isoconversional KAS method. The previously described trends remain the same, but the obtained deviations are quite lower, with 10 or less kJ/ mol for the conversions up to 0.9 with beech wood, as well as for the cellulose and hemicellulose compounds in beech wood and for pure cellulose. Besides, the averaged activation energy of the cellulose component in the selected cases (183.9 ± 10.3 kJ/mol) is higher than when considering all cases (174.6 ± 25.3 kJ/mol) and very close to the values commonly reported in literature from detailed analysis (in the vicinity of 47 kcal/ mol − 197 kJ/mol -according to Burnham et al. [52]).
It also interesting to discuss possible causes leading to results outside the general trends for cases #3, #5 and #7, which are the finally not selected cases because they could not reproduce the results with pure cellulose. For cases #3 and #7 higher initial masses were employed than other cases (see Table 1). It is therefore advised to employ initial masses lower than 10 mg (the exact value probably depends on the device), as higher values led to deviations in this work. This is consistent with previous literature works with cellulose, which showed that thermal lag can only be completely avoided with samples of a few mg at low heating rates (up to 10 K/min) and even smaller samples at higher heating rates [52,53]. Despite it, initial sample masses of 10 mg and higher are still commonly employed in literature for determining biomass pyrolysis kinetics, as shown in the review of recent works presented in the introduction, where the initial sample mass was on average of 9.2 ± 5.0 mg. For case #5, with an initial sample mass of 3 mg, it was found after the experiments of the round robin that differences in the peak temperature (15°C lower than the average for pure cellulose) was due to: 1) the temperature calibration, as a recalibration considering the nonlinear behaviour of the furnace lead to a difference Table 5 Activation energies (kJ/mol) along conversion from KAS method for all beech wood pyrolysis experiments -data of cases #2 and #5 presented as well without outliers (#2 at 1 K/min and #5 at 20 K/min).   of 3.8°C in the peak temperature for cellulose (the furnace was originally calibrated with melting points of 5 pure metals including Ag and Au -melting points of 962 and 1064°C -, which led to deviations at lower temperatures), and 2) over geometry and purge gas flow direction, as the purge gas flows downwards in the employed device (as opposite to the most common devices) promoting the decomposition; and this effect can be diminished by using pan lids, leading to a difference of 7.2°C in the peak temperature for cellulose. Furthermore, deviations in the target heating rate were present for case #3 at 10 and 20 K/min. Additionally, for case #4 the heating rates were higher than the target. A detailed analysis shows that the experiment at 20 K/min of case #4 is slightly outside the general trend, which is probably due to the higher heating rates obtained in that experiment when a significant conversion was still taking place. Slightly higher activation energies would be obtained for case #4 without the experiment at 20 K/min, but the current deviations are considered acceptable.
Other characteristics of the employed TGA instruments and methods besides the initial sample mass have shown to not be critical for deriving reliable kinetics. The finally selected cases as the most reliable (#1, #2 without outlier, #4 and #6) span the whole range of tested possibilities regarding TGA configuration type (horizontal and vertical), N 2 purge flow, sensitivity and temperature accuracy (from the lowest to highest values), sample holder (ceramic, platinum and alumina) and calibration method (Curie point or melting point), as shown in Table 1.
The results of this work show that a certain variability is present in the experimental results employing several TGA devices with beech wood and pure cellulose pyrolysis, and it is of a similar order regarding temperature differences than in the round robin conducted by Gronli et al. [7] with pure cellulose. However, these deviations can lead to larger discrepancies in the calculated activation energies for beech wood pyrolysis, where several components decompose in parallel. This result should not discourage the use of TGA to determine kinetics, but this needs to be conducted with care, reproducing first TGA results of pure cellulose and employing isoconversional methods as suggested in this work to validate the reliability of the experiments. Besides, kinetics should be derived with reliable methods, as in this work combing fitting and isoconversional methods, which is not always the case in literature and adds further variability to the reported values.
Finally, it is not the target of this work to claim that the conducted kinetic analysis for beech wood pyrolysis is the optimal one. A parallel scheme was selected as it is the most common option in literature and it provides reasonable results for mass loss, although it cannot predict the yield of different products as char or bio-oil. An n th order reaction for hemicellulose and first order reaction for cellulose provided good results. The determination of the lignin parameters is more complicated A. Anca-Couce, et al. Fuel 276 (2020) 118002 and the employed n th order reaction is one possible alternative. With the current approach the lignin component is modelling only the tail at high temperatures. Lignin pyrolysis actually covers a wider temperature range and it starts reacting at low temperatures [54]. In fact, the proportion of lignin in hardwood is around 22% [55], which is higher than the ≈ 10% proportion of the lignin component in the fittings in this work. It is already known that the proportions of each component in a parallel reaction scheme derived from the fitting method with TGA data does not completely resemble the biomass composition, due to interaction of compounds or different char yields. The high deviation in this case for lignin is also due to lignin reactions taking place at lower temperatures which are not covered by this component in the current fitting. Additionally, a very high reaction order is obtained in this work for lignin. This high reaction order, and probably also the very high activation energy for lignin in some cases, is mainly an artefact so that it is possible to model a compound with a relatively flat DTG curve, which in principle can be obtained for one single reaction with a low activation energy (which is not meaningful according to the KAS method) or a high reaction order (which was the case of this work). This high reaction order does not represent a real chemical phenomenon and it rather masks the complexity of several consecutive chemical reactions being modelled with one single reaction. These issues (proportion and reaction order of lignin) are however not a main limitation of this work, as it was an advantage that the lignin compound was covering only the conversion at higher temperatures, to have an easier comparison with isoconversional methods at high conversions. Anyhow, the obtained high values for the activation energy of lignin should be taken with caution, specially when there are big differences with the values from isoconversional methods. Furthermore, this work supports that for lignin a first order reaction or low activation energies [9,54] (which are employed/derived in many works, and are commonly cause and consequence respectively), are not supported by isoconversional methods, at least for the high conversion region. However, other models may be more suited to describe lignin than the one employed here (e.g. DAEM [14]). Other reaction models can also be employed for cellulose and hemicellulose, but the obtained activation energies should be supported by isoconversional methods, as in this work.

Conclusions
The determination of mass loss kinetics for biomass pyrolysis is still an unresolved topic, due to the broad range of values reported in literature. An international round robin of TGA pyrolysis experiments with pure cellulose and beech wood has been conducted by 7 partners. Cellulose pyrolysis has been modelled with one first order reaction and beech wood pyrolysis with 3 parallel reactions employing a fitting routine. The isoconversional KAS method has been employed to support the kinetic analysis for beech wood. It was shown that certain deviations are obtained in DTG curves for all cases, of around 10°C in the position of the peak at 5 K/min and 20 -30 kJ/mol in the predicted activation energies for cellulose, hemicellulose and conversions up to 0.9 for beech wood. Higher deviations are obtained for yet higher conversions and for lignin.
The following method [9,16] for determining biomass pyrolysis kinetics has been employed and is hereby recommended in order to increase the reliability of kinetics derived from biomass pyrolysis: • Pure cellulose pyrolysis experiments conducted at 5 K/min from Gronli et al. [7] should first be reproduced, in order to assess the employed TGA device and methods. This study shows that a good reproducibility can be generally obtained and the cases with higher deviations were the ones who led to higher discrepancies in the determination of kinetics for beech wood. Relevant examples were two cases for which an initial mass higher than 10 mg was employed, leading to a significantly higher char yield. As often stated in literature, lower initial mass samples are recommended to avoid Table 8 Comparison of activation energies in kJ/mol with fitting method for cellulose (CEL), hemicellulose (HC) and lignin (LIG) and KAS isoconversional method at selected points derived from all experiments without outliers (#2 at 1 K/min and #5 at 20 K/min).  Table 9 Averaged kinetic parameters from model fitting method and comparison of activation energies to from fitting and isoconversional KAS method for selected cases #1, #2 (without 1 K/min), #4 and #6. • Pyrolysis experiments with biomass should be conducted at different heating rates and isoconversional methods should be employed in order to validate the reliability of the experiments, as well as of the derived kinetic parameters from a fitting method. This study shows that high R-squared values (over 0.99), except at very low and high conversions, are to be expected from isoconversional methods. Outliers who clearly did not fulfil this criterium were detected for two separate cases and not considered for further analysis. Finally, isoconversional methods should support the activation energies for cellulose and hemicellulose derived with a fitting method, while higher quantitative deviations are expected for lignin.
Following strictly this method, the data of three out of seven cases in this work was not selected for the final evaluation, and for a fourth one the data at one heating rate was disregarded. The deviations in the values of activation energy for these selected cases following the previous method were of 10 kJ/mol or lower (less than half than when considering all experiments), except for lignin. An activation energy of around 180 kJ/mol was obtained for the cellulose component in beech wood, which was a bit lower than for pure cellulose. A value of 150 kJ/ mol was derived for hemicellulose while the in the case of lignin the value was higher than 200 kJ/mol. This method does not completely guarantee that optimal kinetics are derived, but at least ensures that the obtained kinetics are chemically meaningful and can help in the effort to reduce the variability in biomass kinetics in literature, which can be attributed to a significant extent to the lack of consistency in data acquisition and kinetic analysis of TGA experiments.