A Novel Method of Kinetic Analysis and Its Application to Pulverized Coal Combustion under Different Oxygen Concentrations

Currently, many efforts have been made to improve the approach to build kinetic models. Based on mathematical algorithms, a novel method (named DIM method) of kinetic analysis was introduced in detail. A formula combining differential and integral was deduced and applied to the determination of the mechanism function f (α). Subsequently, multivariable linear regression was conducted to simultaneously obtain the apparent activation energy E, pre-exponential factor A, and oxygen concentration exponent n. In the application of pulverized coal combustion under different oxygen concentrations (3%, 5%, 10%, 15%, and 21%), E, A, and n were calculated as 258,164 J/mol, 6.660 × 1017 s−1, and 3.326, respectively, and the mechanism function f (α) was determined as the Avrami-Erofeev equation. A validation was performed under a 7% oxygen concentration, which shows that the DIM method has a higher accuracy. This work can provide a reference for the study of kinetic analysis.


Introduction
Coal is one of the main fossil fuels [1] and is widely used in industrial boilers, so the efficient use of coal is crucial to the economy.Moderate or intense low-oxygen dilution combustion can increase the efficiency and reduce the NOx emission of coal combustion [2,3].Furthermore, improving efficiency necessitates better understanding of the mechanisms of the coal combustion process.During novel method and model researches, it is an effective solution to make a theoretical analysis based on the engineering application or experimental data, which was used in the work of Brodny et al. [4][5][6], Pach et al. [7], and Adamczyk et al. [8].Usually, the kinetic triplet (mechanism function f (α), apparent activation energy E, and pre-exponential factor A) is vital to characterizing the process of coal combustion.Particularly, the mechanism function is the functional relationship between the reaction rate and conversion ratio, which reveals the chemical reaction mechanism and directly determines the tendency of the thermogravimetry (TG) curve.The mechanism function was firstly proposed in the 1920s because the reaction order in the homogeneous phase reaction was incapable of characterizing the complex heterogeneous reaction, and it was rapidly developed in the 1930-1940s based on numerous experiments [9].The exploration of the mechanism function laid a foundation for the kinetic model.
In the studies of the coal network model by Fletcher et al. [10], Shurtz et al. [11], and Niksaa et al. [12], multiple chemical and empirical parameters were prerequisites.In addition, due to the heterogeneous characteristics of coal constituents, the carbon-oxygen reaction kinetics in char combustion, for example, are believed to be influenced by ash, maceral group properties of coal, particle size, diffusion phenomena, and flame temperature [13], which indeed impact the complexity of network models.Thus, developing a global kinetic model with the proper mathematical form is a reasonable approach for modelling practical combustion and reproducing the major features of kinetic data.
For a long period of time, many scholars have been trying to explore the reaction mechanisms and the kinetic analysis methods to reproduce the basic characteristics of coal combustion.Wang et al. [14] performed a kinetic analysis of raw coal and pyrolysis carbocoal using TG data.The results showed that the first order kinetic model could fit the combustion process of raw coal and pyrolysis carbocoal well at different heating rates.They also studied the method for obtaining the kinetic model of coal combustion.With nine mechanism functions, the Achar differential and Coats-Redfern integral methods were introduced to obtain E and A from the slope and the intercept of corresponding plots, respectively.With the same mechanism function, if the two E values obtained through differential or integral methods were very close and the correlation coefficient R 2 value was greater than 0.98, the supposed mechanism function could be selected as the suitable one.
Niu et al. [15] conducted combustion experiments of bituminous and anthracitic coals in O 2 /N 2 and O 2 /CO 2 mixture atmospheres, respectively.Based on the Coats-Redfern method, they used an undetermined reaction order n in an assumed mechanism function, to determine E and A. R 2 acting as the selection criterion helped to pick the most fitting n from six options.
The volume model (VM) was used to depict the process of volatile combustion.The random pore model (RPM) presented by Bhatia et al. [16,17], in which the pore structure and its evolution were taken into account, could fit char combustion process well.Based on the combination of VM and RPM, Shao et al. [18] developed a new kinetic model (NEWM) for pulverized coal combustion.The results using NEWM showed greater agreement with experimental data and smaller deviation compared to separate application.
Liu et al. [19] regarded the process of coal pyrolysis as a first order reaction.After a linear least square fit, they found that three independent first order reactions fitted the experimental data well.The method employed the piecewise kinetic model to describe the coal pyrolysis in partial temperature regions.The similar dynamic analysis using a thermal analysis curve was also conducted in reference [20].
Irfan et al. [21] studied the kinetics of coal, palm shells, and their blends using a thermogravimetric analyzer (TGA) under different gas mixture compositions (N 2 /CO 2 /O 2 ).Doyle's and Coats-Redfern's approximations were employed to determine E due to their simplicity in manipulating the variables of equations into linear forms.
Generally, the mechanism function f (α) is obtained firstly and followed by apparent activation energy E, and pre-exponential factor A in the kinetic triplet solution using data of TG.When determining f (α), the empirical method [25,26] and discrimination method [27,28] are mainly used.In the empirical method, f (α) is determined referring to corresponding references in which they have similar mechanisms of reaction.The discrimination method uses various approaches to select f(α) from a batch of mechanism functions.With regard to E and A, they are almost calculated from the slope and intercept of corresponding plots.The kinetic triplet can help describe the reaction process and establish a related kinetic model.In consideration of expecting to obtain a reliable kinetic method with a high accuracy, this paper aims at a novel method combining differential and integral, and multivariable linear regression (named DIM method).In the application of the DIM method, an apparent kinetic model including f (α) E, A, and n can be obtained in pulverized coal combustion under different O 2 concentrations.Furthermore, validation shows the accuracy of the DIM method by comparison.

DIM Method
The study of chemical reaction kinetics dates from the late 19th century.Van't Hoff proposed the reaction order in 1884.The Arrhenius rate constant formula was presented in 1889.In 1891, Wilhelmy found the fact that the conversion rate is proportional to the amount of residual sugar in acidic conditions and established the original dynamic equation.Generally speaking, the homogeneous phase reaction kinetics equation under an isothermal condition was basically completed in the 19th century and shown as Equation (1) [9]: where c represents the concentration of the product; t is time; k(T) stands for the rate constant; T is temperature; and f represents the mechanism function.
In the early 20th century, kinetics research derived from a homogenous reaction and isothermal condition to a heterogeneous reaction and non-isothermal conditions.Basically, a slight corresponding adjustment was made to fit the conditions of a heterogeneous reaction and non-isothermal environment.
The conversion ratio α is expressed as Equation (2).
where w 0 , w t , and w ∞ are the initial weight, the weight at the examined time t, and the final weight at the end of the process, respectively.Since the heating rate β = dT/dt, dt = dT/β.Then, dT/β is used to replace dt in Equation (1).Furthermore, with a substitution of the conversion ratio α for product concentration c, the kinetic equation matching a heterogeneous reaction and non-isothermal conditions is described as the following expression: Then, Arrhenius expression k(T) = A exp(−E/RT), which is commonly used to express the reaction rate as a function of temperature, is taken into Equation (3).Equations ( 4) and ( 5) are known as the differential and integral kinetic equations of the first kind, respectively [29].
where A, E, and R denote the pre-exponential factor, activation energy, and universal gas constant, respectively.Usually, the purpose of the kinetic study is to obtain the kinetic triplet f (α), E, and A.
At the beginning of reaction, the reaction rate can be ignored due to relatively low initial temperature T 0 .Finally, the expression is written as Equation ( 6), with the temperature range from 0 to T.
u 2 du denotes the temperature integral, and u = E RT .Because of no analytical solution for P(u) in mathematics, a few scholars developed some kinetic analysis methods based on differential approximate solutions, for example, the common Coats-Redfern method [30].Inevitable approximation deviations exist in the temperature integral of the first-kind kinetic equations.Besides, many complex reactions exist at different temperature ranges of the coal combustion.When conducting kinetic analysis in separated temperature ranges, omitting T 0 can contribute to big deviations.Therefore, the second-kind kinetic equation is introduced [31].The differential expression is shown in Equation (7).
After the separation of variables, the integral is made at both sides of the Equation (7).Then, we have the deduction below:

RT
A kinetic integral formula containing an analytical solution is finally evolved as Equation ( 8), which is the analytical solution rather than the approximate solution to the temperature interval.
Using βdt to denote dT in Equation ( 7) and through transposition, Equation ( 9) emerges. dα/dt Equation ( 10) is obtained by the logarithms at both sides of Equation ( 9).ln Transposing the Equation ( 8) and taking the logarithms at both sides leads to Equation (11).
The right parts of Equations ( 10) and ( 11) are the same.By substitution, we have Equation (12).
Expression of E is consequently deduced as Equation ( 13).
The aforementioned method to deduce Equation ( 13) combines differential and integral, which has the ability of compensating for the weakness when a differential or integral method is separately applied.Crucially, effective kinetic analysis is available in a range of temperatures for complex reactions, thus approximate deviation contributed by inappropriate division of temperature can be reduced.
In Equation (11), the plots of ln βG(α) T−T 0 against 1/T are generated and the least square fit method is adopted to test data with correlation coefficient R 2 as the criterion.The comparatively bigger R 2 , which is usually required to be bigger than 0.98 [15,21,32], indicates a better fitting effect.E is calculated from the slope of the linear plot of ln βG(α) T−T 0 against 1/T.The another way to calculate E (E expection) is through Equation ( 13) within the temperature range.If the E values by the two different ways are approximate, it is deemed that the assumed mechanism function is of the most probability.
After determining the mechanism function f (α), multivariable linear regression is used to obtain kinetic parameters.The reaction process with changing O 2 concentration is taken as an example.
The apparent kinetic model for different O 2 concentrations is assumed as: in which, n and f o 2 denote the O 2 concentration exponent and O 2 concentration, respectively.After the replacement with f (α), logarithm is made at both sides of Equation ( 14): In order to calculate the values of A, E, and n, a multivariable linear regression model is established as: Two input variables x 1 , x 2 and an output y exist in this model.In Equation (15), the left side is regarded as y; lnA, E, n, −1 RT and ln( f o 2 ) at the right side are respectively corresponding to β 0 , β 1 , β 2 , x 1 , and x 2 .
y, x 1 , and x 2 in Equation ( 16) are replaced with m sets of TG experiment data at different O 2 concentrations.Equation ( 17) can be obtained: in which, m and k are the numbers of observation data and explanatory variables, respectively.Equation ( 17) can be expressed in matrix as Equation ( 18): βj denotes the estimated value of β j , which can be calculated by Equation (19): X is the inverse matrix of the matrix X.
ŷi as the estimated value of y i is written by Equation ( 20): After acquiring the regression coefficients β 0 , β 1 , and β 2 , residual error is calculated to analyze whether it accords with a normal distribution, which is used as the validation index of the regression model.The residual error is defined as: e = y − ŷ (21) The confidence interval of β j is described as: where β Dj and β Uj represent the lower and upper boundaries of the confidence interval, respectively; t γ/2,m−k−1 is the critical value of t distributions with a significant level γ and m − k − 1 degree of freedom.The estimated value of variance of βj is: The estimated values of covariance of βi and βj are: Multivariable coefficient of determination R 2 is defined as: In Equation (31), in which, y is the mean value of y.
A hypothesis test is subsequently conducted to test whether a linear relationship exists between at least one regression variable and output.
The above-mentioned method involving determining f (α) and using multivariable linear regression to get kinetic parameters is the DIM method.

Experimental Results
The experiment of pulverized coal combustion under different O 2 concentrations was conducted by TGA.During the experiment, the TGA/SDTA851 e (Mettler Toledo, Zurich, Switzerland) analyzer was adopted to measure the mass loss percentage and mass loss rate of samples.Table 1 shows the proximate and ultimate analyses of Yanzhou bituminous coal.In the experiment, the sample was initially weighted to (10 ± 0.1) mg and five samples of an equal weight were prepared.The temperature increased linearly from room temperature to 1273 K at the heating rate of 10 K•min −1 .The reaction gas is the mixture of O 2 and N 2 with various O 2 concentrations (O 2 = 3%, 5%, 10%, 15%, and 21%), and the flux is kept at 30 mL•min −1 .
TG and DTG (differential thermogravimetry) curves of Yanzhou pulverized coal with different O 2 concentrations are described as Figures 1 and 2, respectively.
Energies 2018, 11, x FOR PEER REVIEW 7 of 15 In the experiment, the sample was initially weighted to (10 ± 0.1) mg and five samples of an equal weight were prepared.The temperature increased linearly from room temperature to 1273 K at the heating rate of 10 K•min −1 .The reaction gas is the mixture of O2 and N2 with various O2 concentrations (O2 = 3%, 5%, 10%, 15%, and 21%), and the flux is kept at 30 ml•min −1 .
TG and DTG (differential thermogravimetry) curves of Yanzhou pulverized coal with different O2 concentrations are described as Figures 1 and 2, respectively.

The DIM method Application and Discussion
The DIM method presented in Section 2.1 was applied to build the kinetic model of pulverized coal combustion under different O2 concentrations.The conversion ratio calculated by TG data was inserted into forty-one common kinetic mechanism functions summarized as Table 2 [28,33,34].

The DIM Method Application and Discussion
The DIM method presented in Section 2.1 was applied to build the kinetic model of pulverized coal combustion under different O 2 concentrations.The conversion ratio α calculated by TG data was inserted into forty-one common kinetic mechanism functions summarized as Table 2 [28,33,34].

Number
Model Mampel power law α (1 − α) As presented in Section 2.1, five sets of experimental data at different O 2 oncentrations were used to plot ln βG(α) T−T 0 against 1/T (Equation ( 11)), and then 5 × 41 sets of E values were obtained.Furthermore, 5 × 41 sets of E values were also calculated by using Equation (13).At each O 2 oncentration, 41 pairs of E values from the two above-mentioned algorithms were compared for finding the closet E values.Finally, mechanism function No. 19 has the closest E values for the five sets of experimental data, and therefore it is the most probable mechanism funciton.For simplicity, 41 pairs of E values at a 10% O 2 concentration were presented as an example to show the comparison of E values in Table 3.
Table 3.The comparison of E values with the forty-one mechanism functions at a 10% oxygen concentration.

No.
E 1 (Equation ( 11)) (J/mol) E 2 (Equation ( 13)) (J/mol) Based on the DIM method, the most suitable kinetic mechanism function for this experiment is the Avrami-Erofeev equation with n = 4 (No.19 in Table 2) and exhibited as: The calculation process of multivariable linear regression above is attained by computer programming.The estimated value of the regression coefficient β j , the confidence interval with its confidence level higher than 95% [35], and the R 2 are shown as the following: β 0 = 41.040, and the confidence interval is [39.961,42.118]; β 1 = 258,164, and the confidence interval is [250,853, 265,475]; β 2 = 3.326, and the confidence interval is [3.272, 3.380 Then, the regression model is: Therefore, A = e 41.040 = 6.660 × 10 17 s −1 , E = 258,164 J/ mol, and n = 3.326.
A, E, n, R = 8.314 J/(mol•K) and f (α) are used to substitute the corresponding values in Equation (14).After rearrangement, the consequently apparent kinetic model of coal based on this experiment data is described as Equation (36).

Validation
The DIM method was discussed above for selecting the mechanism function and calculating the kinetic parameters.To validate the DIM method, the calculated curve using Equation ( 36) was compared with the experimental curve at a 7% O 2 concentration ( f o 2 = 7%).The experiment of coal combustion with a 7% O 2 concentration was conducted in the same conditions as the experiment aforementioned in Section 2.2.In addition, the Achar differential method [36], Coats-Redfern integral method, and the method developed by Wang et al. [14] (called the Wang method below in this paper) were respectively applied and the comparisons were provided for the prediction effect with the DIM method.
The universal Achar differential is: The universal Coats-Redfern integral is: The ln dx dT f (x) and ln are separately regarded as y; the ln(A/β) and ln AR βE are separately regarded as b; − E R is k; and 1 T is x.Then, Equations ( 37) and ( 38) can be converted into linear equations.The plots of ln dx dT f (x) and ln G(x) T 2 vs. 1  T were drawn for each mechanism, and E values were determined by the slope after selection of the model with the biggest R 2 [15,33,37].
The Wang method includes both Achar and Coats-Redfern methods.After the linear regression, two E and A are respectively obtained from the slope and intercept of the fitting curves with regard to one mechanism function.The suitable mechanism function can be judged by the closest two E values with minimum relative difference, and the relatively greater R 2 .
With experimental data, the forty-one sets of E, A, and R 2 using Achar and Coats-Redfern are listed in Table 4  Based on the data in Table 4, the Achar method and Coats-Redfern method select the same mechanism function No. 19 shown in Table 2, but apply their own E and A. The Wang method chooses the mechanism function No. 6 shown in Table 2, but selects the E and A from the Coats-Redfern method.With the calculated E, A and the selected f (α), an experimental curve and four calculated curves for the kinetic model are generated in Figure 3. Four calculated curves include the DIM method curve, Achar method curve, Coats-Redfern method curve, and Wang method curve.Based on the data in Table 4, the Achar method and Coats-Redfern method select the same mechanism function No. 19 shown in Table 2, but apply their own E and A. The Wang method chooses the mechanism function No. 6 shown in Table 2, but selects the E and A from the Coats-Redfern method.With the calculated E, A and the selected (α), an experimental curve and four calculated curves for the kinetic model are generated in Figure 3. Four calculated curves include the DIM method curve, Achar method curve, Coats-Redfern method curve, and Wang method curve.Moreover, deviation is used to evaluate the accuracy of the model.The average deviation between the calculated and experiment curves is introduced as Equation (39) [38].Moreover, deviation is used to evaluate the accuracy of the model.The average deviation δ between the calculated and experiment curves is introduced as Equation (39) [38].
where δ is the average deviation; x i,exp , x i,calc , and max x i,exp respectively represent experimental values, calculated values, and the maximum experimental value; and N is the amount of experimental points.Based on Equation (39), the average deviations from different methods are calculated and shown in Table 5.Table 5 shows that, compared with other methods, the DIM method has the lowest average deviation (1.26%), showing a good agreement with experimental observations.

Conclusions
The DIM method of kinetic analysis was introduced, and its application to coal combustion under different oxygen concentrations was demonstrated and validated.

1.
A formula combining differential and integral was deduced through analytical approach, which can offset the defects of a single differential or integral method.2.

3.
The experimental TG curve with a 7% O 2 concentration was compared with four calculated curves generated by DIM, Wang, Achar, and Coats-Redfern methods, respectively.The DIM method showed a good accuracy with 1.26% average deviation.

Figure 3 .
Figure 3.The four kinds of calculated curves (DIM, Wang, Achar, and Coats-Redfern methods) and experimental TG curve at a 7% O2 concentration.

Figure 3 .
Figure 3.The four kinds of calculated curves (DIM, Wang, Achar, and Coats-Redfern methods) and experimental TG curve at a 7% O 2 concentration.

Table 1 .
Proximate and ultimate analyses of Yanzhou bituminous coal.
. The No.19, No.19 and No.6 mechanism functions shown in Table 2 were singled out with the Achar differential method, Coats-Redfern integral method, and Wang method, respectively.

Table 4 .
E, A, and R 2 of pulverized coal combustion with the forty-one mechanism functions.

Table 5 .
Average deviation comparison of different methods.