Inverse Modeling of Thermal Decomposition of Flame-Retardant PET Fiber with Model-Free Coupled with Particle Swarm Optimization Algorithm.

The thermal decomposition model of flame-retardant polyethylene terephthalate (FRPET) fiber is essential for predicting its fire behavior and do relevant fire simulation. In this work, the thermal decomposition character of FRPET is investigated via thermogravimetric analysis at four heating rates. Two kinetic schemes are proposed based on the analysis of experimental data and model-free methods. The model-free methods (Friedman and advanced Vyazovkin methods) are employed to determine possible search range for particle swarm optimization algorithm with constriction factor (CFPSO). Thus, this coupled method could evaluate the kinetic parameters for two proposed schemes without initial guess. Both models could reasonably predict the experimental data with obtained parameters, and the second two-step consecutive model shows better performance. The performance of CFPSO on the second model is further compared with improved generalized simulated annealing algorithm, and CFPSO was found to be more effective. Furthermore, global sensitivity analysis was conducted via the Sobol method to investigate the influence of kinetic parameters for the second model on predicted results. The most influential parameters are ln A and E α of the second reaction and reaction order n of the third reaction.

Introduction

Polyethylene terephthalate (PET) has been widely used due to its merits such as high tensile strength, low cost, and light weight.[1] However, this material is relatively flammable due to its organic nature, and this undermines its applicability.[2] Thus, PET is usually added with flame retardants (FRs) to suppress the combustion process and meet the regulatory criteria when used in high-risk situations (e.g., upholstered furniture and mass transport). To assess the fire risk and mitigate the fire hazard in this situation, the knowledge of thermal decomposition character of flame-retardant PET (FRPET) is essential. This is because thermal decomposition is the first step of combustion—the combustible volatiles released in the thermal decomposition process could feed the flame zone and further facilitate the fire growth. The product yield in this process could be described using a kinetic model. Therefore, to further understand and predict the fire behavior of FRPET, the kinetic modeling of its thermal decomposition is important.[3,4]

The modeling of thermal decomposition usually starts with the kinetic parameter estimation. Those parameters could be obtained with thermogravimetric (TG) tests coupled with the model-free or model-fitting method. Model-free methods could evaluate kinetic parameters without knowing the reaction scheme.[5] However, those obtained kinetic parameters usually vary with the extent of conversion because of the complex reaction, which might be less useful in building the model. On the other hand, the model-fitting approach uses different algorithms to minimize the error between the experimental data and predictions to find the kinetic parameters. This approach, therefore, is a type of an inverse modeling method, and the kinetic parameters in this approach are invariant for a single reaction. There have been many algorithms used in the model-fitting method (i.e., Levenberg–Marquart algorithm,[6] Genetic algorithm,[7] Shuffled Complex Evolution algorithm,[8] etc.). Most of those studies need proper initial guess (from previous studies or from model-free methods) to determine the fittest kinetic parameters.[8,9] Few studies noticed that the Friedman method and advanced Vyazovkin method could provide proper search range for the model-fitting method, and thus, the model-free and model-fitting methods could be coupled to find the fittest kinetic parameters without initial guess.

Model-free and model-fitting methods have been employed in many reported literatures on decomposition of PET. Jenekhe et al. evaluated the kinetic parameters of PET in non-isothermal decomposition using the Flynn–Wall method.[10] Cooney et al. utilized Kissinger, Freeman-Carroll, and other model-free methods to evaluate the kinetic parameters of decomposition of PET in air.[11] They also established that at least three reaction stages occurred in decomposition of PET. On the other hand, Yang et al.[12] employed model-fitting methods on decomposition of PET in nitrogen and obtained an activation energy of 242 kJ/mol. Saha used the model-fitting method to evaluate the kinetic triplet of PET in non-isothermal and isothermal decomposition processes and found that nth order model could better predict the experimental data.[13] Martín-Gullón et al.[14] and Moltó et al.,[1] respectively, modeled the decomposition process of PET in a nitrogen atmosphere and air with the model-fitting method. Although model-free and model-fitting methods were frequently used for kinetic estimation and kinetic modeling of pyrolysis of pure PET, rare studies employed the model-free coupled with model-fitting method to model the pyrolysis process of FRPET.

Consequently, the model-free coupled with model-fitting method are proposed in this work. It is applied to model the decomposition process based on two decomposition mechanisms, which are developed based on the experimental data and results of model-free methods. The model-free methods also determined possible search range for the model-fitting method. Then, particle swarm optimization algorithm using constriction factor (CFPSO) as the model-fitting method was used to find the fittest kinetic parameters. The obtained kinetic parameters for two proposed mechanisms are validated by reconstruction of the experimental data. To show the performance of CFPSO on the second model, the results of this method are compared with results from improved generalized simulated annealing algorithm. Furthermore, to analyze the influence of input kinetic parameters on output of the second model, global sensitivity analysis (SA) was performed via a Sobol method.

Experimental and Methods

Samples

The FRPET samples were cut from bus seat assembly cover materials (supplied by Zhongtong Bus Holding Co., Ltd.). Those cover materials followed the requirement of the FR performance of GB 38262-2019 standard—flammability of interior materials for buses. The samples were sliced into pieces (less than 1 mm) and were then dried for 24 h to eliminate the moisture.

TG Analysis

TG analysis experiments were conducted in SDT Q600 (TA instruments). 4–5 mg samples are used because smaller sample weight could reduce the thermal lag effect. The samples were heated from room temperature to 800 °C, and the heating rates were 5, 10, 15, and 20 °C/min, respectively. All tests were conducted under nitrogen flow (100 mL/min) using alumina crucibles without lid.

Kinetic Modeling

The reaction rate of the solid material in the non-isothermal experiment for a single reaction is usually modeled based on the following kinetic equation[15]where β = dT/dt is the heating rate, T is the absolute temperature, R is the gas constant, A and Eα are the pre-exponential factor and activation energy, and f(α) is the function that represents the reaction model. For polymer thermal decomposition, f(α) is usually represented by the nth order reaction model.[7,16,17] Thus, f(α) = (1 – α) was used in this work and α is the extent of conversion, which is given bywhere mi/mf are the initial mass/final mass, and m is the sample mass in the experimental process. Based on this equation, the TG curve can be transformed into the α–T curves.

Model-Free Isoconversional Methods

Model-free isoconversional methods could evaluate the dependence of activation energy on conversion. In this study, the Friedman method and advanced Vyazovkin method were used to obtain both activation energy Eα and pre-exponential factor A.

Friedman Method

Friedman proposed the following expression by taking the natural logarithm of 1(18)

In 3, the values of dα/dT can be evaluated numerically by differentiating the experimental data. Therefore, by plotting βα/dT) against 1/T for different heating rates at given α, the activation energy Eα and Af(α) can be obtained from the slope and intercept, respectively.

Advanced Vyazovkin Method

The Friedman method is sensitive to instantaneous experimental noise and might introduce inaccuracy.[5] Therefore, the activation energy Eα and Af(α) were also evaluated with the advanced Vyazovkin method.[19,20] This method could avoid inaccuracy in the Friedman method, and it is free of the approximations of temperature integral used in other integral methods like Flynn–Wall–Ozawa methods and Kissinger–Akahira–Sunose methods. According to Vyazovkin, for a set of n experiments with different heating programs, the activation energy at certain α can be obtained by finding the Eα value which minimizes the following functionwhere i and j indicate ith and jth experiment, and J(Eα) is defined as

This integration was evaluated numerically from α – Δα to α using modified Simpson’s rule, and the minimization of Φ(Eα) is solved with the Brent method.[21] This process is repeated for each certain α, then dependence of Eα on α can be obtained.

After that, the Af(α) for each α was evaluated using the following equation proposed by Lina[22]

This equation is based on the assumption that for small interval Δα, only one reaction occurs, and therefore, the Af(α) and Eα can be treated as constant, where is the average of J(Eα) for all experiments based on optimized activation energy.

For both the Friedman method and advanced Vyazovkin method, analysis is performed between α = 0.02 – 0.9, and Δα is set to 0.02. For particular α, the experimental data is interpolated using quadratic interpolation with scipy.interpolate.interp1d.[23]

With obtained Af(α), the pre-exponential factor A can be determined by substituting f(α) = (1 – α) with different reaction order. In this study, the reaction order n was set to [0, 5] as in Ding’s work.[24] Thus, by substituting the f(α) with n = 0 and n = 5 for each Af(α), the range of A can be obtained. Finally, those obtained kinetic parameters can be used to determine the search range of model-fitting methods.

Non-linear Model-Fitting Methods

All model-fitting methods involve minimizing the error between predictions and the experimental data.[5] However, the classical linear model-fitting method shows worse performance when dealing with complex reactions. Thus, the non-linear model-fitting methods are employed in this study.

Before doing model fitting, the decomposition mechanism of FRPET must be determined. This mechanism is based on the analysis of the results of TG test and model-free methods. With the determined decomposition mechanism, a system of ordinary differential equations (ODEs) can be developed to describe the reaction rate of FRPET in the decomposition process. This part will be detailed discussed later.

When the ODE system is developed, the unknown parameters in the ODE system can be estimated by fitting the experimental data. In this process, the error between the experimental data and solutions of the ODE system are minimized with optimization tools. Therefore, the optimal parameters which could best reproduce the experimental data were determined. However, this inverse problem is highly non-linear with high dimensional search space. To solve this problem, particle swarm optimization (PSO) algorithm as an efficient stochastic global optimization algorithm is employed in this study.[25]

Particle Swarm Optimization Method Using Constriction Factors

CFPSO is proposed by Clerc,[26,27] who has established that using the constriction factor could ensure the convergence of PSO. According to Clerc, the velocity and position of each particle are updated, according to the following equationwhere k and i indicate the kth iteration and ith particle. r1 and r2 are two random values between [0, 1] following the uniform distribution, while p and p are the best position and global best position of particles for each iteration. The constriction factor K is given bywhere ϕ = c1 + c2, ϕ > 4. Typically, ϕ is set to 4.1 and c1 = c2 = 2.05. Thus, K = 0.729.

Improved Generalized Simulated Annealing Algorithm

The improved generalized simulated annealing (IGSA) algorithm is derived from Tsallis’s work.[28,29] IGSA algorithm combined the local search strategy and generalized simulated annealing (GSA) algorithm. The GSA algorithm is developed by Tsallis by generalizing classical simulated annealing algorithm and fast simulated annealing algorithm, according to Tsallis statistics.[29]

The original GSA algorithm used distorted Cauchy–Lorenz visiting distribution, which is governed by parameter qv.where t is the artificial time and D is the dimension of search space. This visiting distribution is used to generate a trial jump distance Δx(t) of variable x(t) under artificial temperature T(t). The artificial temperature is decreased according to

Then, a generalized Metropolis algorithm is used for the acceptance probabilitywhere β = 1/KT is the Lagrange parameter and ΔE is the energy spectrum. The details of GSA can be found in refs (29) and (30).

Finally, the original GSA is improved with Broyden–Fletcher–Goldfarb–Shanno (BFGS) algorithm,[31] which is a large-scale bound-constrained local search strategy. The IGSA has been proved to show good performance,[28,32] therefore it is employed in this work.

Objective Function

The objective function is used for measuring the difference between the experimental data and predictions, as mentioned before. Many objective functions have been employed in previous studies for inverse modeling problems.[33] Bustamante Valencia tested and compared different objective functions, and then, he developed a new objective function considering the phase difference and distance error between curves. This objective function can be expressed as followswhere and are vectors of experimental and estimated mass loss rate (MLR) as a function of temperature.

Global Sensitivity Analysis

The model proposed in this study has multiple input variables. Therefore, we want to identify which input parameter has more influence on the result of objective function via SA. SA methods could be classified into the local SA method and global SA method. Local SA studies consider from small input variation on the model output, while global SA considers the whole variation of the inputs and tries to apportion the output uncertainty to input uncertainty.[34,35]

The Sobol method is a variance-based global SA method which decomposes the output variance into parts attributed to input variables and combinations of variables.[35] Most of the Sobol method used the first-order effect index Si and total effect index STi to measure the main effects of inputs on the outputs and the contributions from inputs to the outputs including interactions among inputs, respectively.[36] Since the Sobol method has been widely used,[37,38] it was employed in this work.

According to Saltelli et al.,[36] given a model of the form Y = f(X1, X2,...,X), it can be decomposed into 14. In this study, Y is the result of objective function and Xi is the input kinetic parameter.

By assuming f(x) is square-integrable, the function could be squared and integrated, then we could get variance (VAR) for Ywhere

Dividing both sides of 15 by Var(Y), we could obtainwhere Si is the first-order effect index and therefore it can be evaluated with 18.

While the total effect index STi can be evaluated with 19

The inputs X are generated with Sobol sequence, which is a low-discrepancy sequence. This sequence shows better performance when used in integration with higher dimensions and therefore was used.[39]

Implementation

Both model-free and model-fitting numerical methods and global SA method used in this study were implemented in Python.

The details of implementation of the model-free method can be found in our previous work.[40]

For CFPSO and IGSA model-fitting methods, the NumPy, SciPy, and Matplotlib module were mainly used.[23,41,42] To speed up the calculation, the Numba is used.[43] The ODE system in the model-fitting method is solved with odeint from SciPy. This module uses LSODA algorithm and could automatically select algorithm to deal with non-stiff and stiff problems.[44] The population size and iteration number of CFPSO are set to 2500 and 5000 in this study. The IGSA comes from the scipy.optimize module. The max iteration number is set to 5000.

For global SA, the SALib module was used in this study.[45] The sampling number N of the Sobol method is set to 10,000.

Results and Discussion

The normalized TG curves of FRPET at different heating rates with pure PET at 10 °C/min are shown in Figure . For FRPET, the TG curves show that at least two decomposition stages are involved for all heating rates. In the first stage, the sample loses about 3% of mass. While in the second decomposition stage, about 79% of mass was lost, with the residue mass of 15.5%. Since the first stage is not shown in pure PET, this stage is probably related to the decomposition of additives.

Figure 1

TG curves of FRPET at 5, 10, 15, and 20 °C/min and pure PET at 10 °C/min.[14]

The differentiated TG (DTG) curves of FRPET are shown in Figure . It shows two peaks in the decomposition process. The peaks of DTG curves would shift to higher temperature with the increasing heating rate. The first decomposition stage is between 150 and 250 °C, while the second decomposition stage occurs between 330 and 510 °C. The temperature range of the second stage is similar to previous findings on pure PET.[14,46] Since the first stage is not shown in pure PET, this stage might be related to the decomposition of FR additives like ammonium polyphosphate. This is because phosphorus-based FR additive has been widely used for many years in PET textiles.[47] Meanwhile, the second stage is corresponding to decomposition of pure PET. Bednas et al.[48] established that FR do not greatly influence the mechanism of pyrolysis of PET. Therefore, the reaction mechanism might be split into two parts—decomposition of FR and decomposition of pure PET, respectively. For the simplicity of kinetic analysis, the following parallel decomposition mechanism is proposed in this work

Figure 2

DTG curves of FRPET at 5, 10, 15, and 20 °C/min.

Kinetic Analysis

Figure shows the plotted lines obtained from the Friedman method for α from 0.05 to 0.85. It can be seen that plotted points fitted line very well and therefore show high R2 values. The curves of activation energy Eα versus extent of conversion α obtained from Friedman and advanced Vyazovkin methods are presented in Figure . The trends of Eα from two methods are similar. The Eα increases at first for α < 0.04, which is between 62 and 140 kJ/mol. The Eα values obtained from the Friedman method are different from values evaluated with the advanced Vyazovkin method. This is because the reaction rate is quite small at this stage and can be affected easily by noise. For 0.04 < α < 0.7, the Eα fluctuates around 200 kJ/mol. The fluctuation at the start of this stage might be related to the stop of the first decomposition reaction and the start of the second decomposition reaction of FRPET. At high conversion (α > 0.7), the kinetic values increased slowly to about 250 kJ/mol. Therefore, the reaction interval can be split into three parts: [0, 0.04], [0.04, 0.7], and [0.7, 1]. Therefore, the Eα – α curve revealed three decomposition stages. This is not shown in DTG curves of FRPET (Figure ), indicating two overlapping reactions occurred in the second mass loss stage of FRPET, which is related to the decomposition of pure PET. Buxbaum[49] and Martín-Gullón et al.[14] have shown that pure PET followed a two consecutive reaction mechanism. Therefore, another two-step decomposition model for FRPET can be developed

Figure 3

Friedman plots of FRPET.

Figure 4

Evaluated activation energy and search range for CFPSO.

The pre-exponential factor A is evaluated using equations and 6 for Friedman and advanced Vyazovkin methods using obtained activation energy, and then, they were converted to the logarithm form. The evaluated results are shown in Figure . As can be seen, the evaluated ln A based on the different reaction order are similar at first. Then, with the increasing α, the differences of ln A increased. In the final stage, the higher reaction order n leads to higher ln A. Obviously, the ln A with the different reaction order 0 < n < 5 is between evaluated ln A curves.

Figure 5

Evaluated ln A and search range for CFPSO.

Comparing Figures and 5, the trends of Eα and ln A curves are similar, this can be explained by the kinetic compensate effect. However, this technique was not suitable for analysis of the consecutive reaction mechanism and therefore was not used.

Kinetic Model

Based on the decomposition mechanism discussed above, the ODE system can be derived. The MLRs for each component of parallel model 1 can be expressed as

For two-step consecutive model 2, the MLR for each reaction is as followswhere k is the stoichiometric yield and mFR,0 and mPET,0 represent the normalized initial mass fraction of FR and PET.

For both models, mFR,0 + mPET,0 = 1, thus only one initial mass fraction needs to be determined. Therefore, for model 1, nine parameters are needed for calculating the reaction rate: two kinetic triplets (A, E, and n), two stoichiometric yields k, and one initial mass fraction mFR,0, respectively. While for model 2, another four parameters (Aint, Eint, nint, and kint) are needed. As a consequence, there are 9 and 13 unknown parameters for model 1 and model 2, respectively.

As mentioned before, the results obtained from the model-free method could provide guidance for search range of kinetic parameters in the model-fitting method. In this way, those two methods were coupled to estimate the kinetic parameters which minimize the objective function.

The search range of Eα and ln A is set to 80% of smaller values and 120% of the larger values evaluated from model-free methods. It is shown as filled area in Figures and 4. Based on the analysis of the experimental data, model 1 contains two reaction stages (0 < α < 0.04 and 0.04 < α), while model 2 shows three decomposition stages (0 < α < 0.04, 0.04 < α < 0.7, and 0.7 < α, respectively). Consequently, the search range of Eα and ln A for each reaction is based on the upper and lower bounds of filled area for each reaction stage of model 1 and model 2. However, for reaction order n, the search range is assumed to be 0–5, as mentioned before.

The initial mass of FR is set to 0–0.25. This is because most of FR additives used for PET are less than 25%.[47] However for stoichiometric yield k for each reaction, no related reference range can be found. Thus, k was assumed to be between [0, 1] in this study. Table summarizes the search range of kinetic parameters for two models. Although some studies show that the results of the model-free method can be used for initial guess, it was not used in this study and the initial guess is generated randomly in search space.

Table 2

Samples Obtained With Sobol Sequences for the Sobol Methoda

no.	YFR	ln A1	Eα1	n1	k1	ln A2	Eα2	n2	k2	ln A3	Eα3	n3	k3
1	0.07	15.40	139.64	3.77	0.30	29.10	227.23	3.16	0.32	28.86	281.50	0.14	0.04
2	0.13	20.60	49.45	0.02	0.55	34.79	201.09	1.91	0.07	36.92	249.37	1.39	0.29
3	0.00	10.20	109.58	2.52	0.06	46.18	148.80	4.41	0.56	53.03	185.11	3.89	0.79
4	0.16	12.80	64.48	4.39	0.18	43.33	188.01	1.29	0.19	40.94	265.44	2.02	0.91
5	0.04	23.20	124.61	1.89	0.68	31.95	240.30	3.79	0.69	57.06	201.17	4.52	0.41
6	0.22	18.00	94.55	3.14	0.43	26.25	214.16	0.04	0.94	49.00	169.04	3.27	0.16
7	0.10	7.60	154.67	0.64	0.92	37.64	161.87	2.54	0.44	32.89	233.30	0.77	0.66
8	0.05	8.90	56.97	3.46	0.74	27.68	155.33	1.60	0.25	38.93	193.14	1.08	0.60
9	0.18	19.30	117.09	0.96	0.24	39.06	207.62	4.10	0.75	55.05	257.40	3.58	0.10
10	0.11	24.49	87.03	4.71	0.98	44.75	233.77	0.35	1.00	46.99	289.54	4.83	0.35
11	0.24	14.10	147.16	2.21	0.49	33.37	181.48	2.85	0.50	30.87	225.27	2.33	0.85
12	0.08	11.50	102.06	0.33	0.37	30.52	168.40	0.98	0.87	42.96	241.34	4.21	0.97
13	0.21	21.90	162.19	2.83	0.86	41.91	220.69	3.48	0.38	26.84	177.07	1.71	0.47
14	0.02	16.70	72.00	1.58	0.12	36.22	246.84	2.23	0.13	34.90	209.20	0.46	0.22
15	0.14	6.30	132.13	4.08	0.61	24.83	194.55	4.73	0.63	51.02	273.47	2.96	0.72
16	0.11	6.95	75.76	1.11	0.40	46.89	223.96	3.32	0.10	56.05	205.19	3.11	0.88
17	0.23	17.35	135.88	3.61	0.89	35.50	171.67	0.82	0.59	39.94	269.45	0.61	0.38
18	0.04	22.55	105.82	2.36	0.15	29.81	197.82	4.57	0.84	31.88	237.32	1.86	0.13
19	0.17	12.15	165.95	4.86	0.64	41.20	250.11	2.07	0.35	47.99	173.06	4.36	0.63
20	0.01	14.75	90.79	2.99	0.52	38.35	210.89	3.95	0.97	27.85	285.52	2.49	0.50

Taking the first 20 sets of kinetic parameters as example.

The experimental data at 5, 10, and 20 °C/min were used to estimate the parameters for two proposed mechanisms with CFPSO, while the data at 15 °C/min were used for validation. The optimized results are listed in Table . The obtained activation energy for the main reaction (second reaction) of model 1 and model 2 are consistent with previous studies.[10,12] However, the differences between obtained kinetic values and other studies might be related to the choice of the model and material.

Table 1

Search Range of Parameters for Both Models and Optimized Parameters

parameters	search range of model 1	optimized values	search range of model 2	optimized values (CFPSO)	optimized values (IGSA)
YFR	[0, 0.25]	0.09	[0, 0.25]	0.13	0.18
ln A1 [ln(1/s)]	[5.06, 25.79]	20.65	[5.06, 25.79]	20.9	11.96
Eα1 (kJ/mol)	[47.45, 167.7]	85.18	[47.45, 167.7]	86.14	53.70
n1	[0, 5]	1.45	[0, 5]	1.83	1.10
k1	[0, 1]	0.74	[0, 1]	0.78	0.84
YPET		0.91		0.87	0.82
ln A2 [ln(1/s)]	[24.36, 58.95]	38.44	[24.36, 47.13]	34.93	34.30
Eα2 (kJ/mol)	[148.01, 293.89]	228.26	[148.01, 252.59]	206.2	203.36
n2	[0, 5]	0.98	[0, 5]	0.84	0.94
k2	[0, 1]	0.15	[0, 1]	0.54	0.46
ln A3 [ln(1/s)]			[26.72, 58.95]	29.25	48.46
Eα3 (kJ/mol1)			[165.36, 293.89]	165.36	237.27
n3			[0, 5]	1.82	5.0
k3			[0, 1]	0.16	0

Model Performance

Figures and 7 show the comparison between experimental and predicted TG and MLR data at the 5, 10, and 20 °C/min. The prediction of model 1 fits well with the experimental data at the first stage. However, this model cannot explain the mass loss in the final stage. Meanwhile, the prediction of model 2 agrees well with the experimental data for the whole decomposition process.

Figure 6

Experimental TG and predicted curves for model 1 and model 2 at 5, 10, and 20 °C/min.

Figure 7

Experimental DTG and predicted curves for model 1 and model 2 at 5, 10, and 20 °C/min.

The validation of optimized parameters at the heating rate of 15 °C/min for model 1 and model 2 is presented in Figures and 9. Obviously, both reconstructed MLR curves closely match the experimental data where R2 values are greater than 0.98. However, the two-step model shows better performance. This result is consistent with discussion before. The predicted MLR curve for pseudo components of FR, PET, and residue could help us understand the decomposition process of FRPET.

Figure 8

Predicted MLR based on optimized values for model 1 (lines) compared with the experimental data at 15 °C/min (symbols).

Figure 9

Predicted MLR based on optimized values for model 2 (lines) compared with the experimental data at 15 °C/min (symbols).

As a consequence, the kinetic parameters obtained from the model-free coupled with CFPSO method could reasonably predict the experimental data. It means that this newly developed method is an effective tool for kinetic inverse modeling.

Performance Comparison of CFPSO and IGSA

To further compare the performance of CFPSO, the IGSA was employed to evaluate kinetic parameters for model 2. The estimated parameters from IGSA are also listed in Table . We could see those parameters are quite different from parameters obtained from CFPSO for the first reaction and third reaction. Therefore, the predicted MLR curve based on those parameters from IGSA at 15 °C/min is shown in Figure . Compared with predicted results of CFPSO, the biggest difference is that no residue is generated in the third reaction and the most of residue are produced in the first decomposition stage of FRPET. This means that residue is generated from additives instead of PET. Obviously, this could not be true. Figure a shows that residue weight is about 10% after the decomposition of pure PET. Therefore, the obtained parameters from IGSA are not acceptable from this viewpoint.

Figure 10

Predicted MLR of model 2 (lines) based on parameters from IGSA compared with the experimental data at 15 °C/min (symbols).

The performance of IGSA is further compared with CFPSO using the convergence curve (Figure ). The fitness of IGSA would converge to its best fitness after about 750,000 function evaluations, while CFPSO would reach to its best value with 1,000,000 function evaluations. Since the optimal fitness of CFPSO is lower than IGSA, the results of IGSA are trapped in local minima. Adenson et al.[50] suggested to use algorithms which show better performance on obtaining the global minimum. Therefore, CFPSO is more effective than IGSA on this inverse modeling problem.

Figure 11

Convergence curve of CFPSO and IGSA for model 2.

Since model 2 shows better performance and contains more input parameters, the Sobol method was employed to conduct the global SA on those kinetic parameters. As mentioned before, the input parameters were generated using Sobol sequence, the samples of those input parameters are listed in Table .

The evaluated Sobol first-order indices and total order indices are shown in Figures and 13. The ln A2, Eα2, and n3 are top three most influential input parameters, indicating that objective function is mainly sensitive to those three parameters and they should be paid more attention. However, the first-order indices k1, k2, and YFR appear to have little influence. On the other hand, the total order indices for all parameters are higher than first-order indices, indicating higher-order interactions between each input parameter.

Figure 12

First-order indices of the Sobol method.

Figure 13

Total order indices of the Sobol method.

Conclusions

The decomposition kinetic parameters for two models of FRPET are estimated with a new method by coupling model-free and model-fitting methods. The model-free methods—Friedman and advanced Vyazovkin methods supplied guidance for the search range of model-fitting method, and therefore, initial guess is unnecessary. The CFPSO is employed as a model-fitting method to find the optimal kinetic parameters. With the possible decomposition mechanism based on analysis of experimental results, both developed model in this study could accurately predict experimental data, and the second two-step consecutive model shows better performance. The performance of CFPSO on the second model is compared with improved generalized simulated annealing algorithm. The CFPSO shows better performance on determining the global optimum on this problem. The global SA with the Sobol method shows that the top three influential kinetic parameters for the second model are ln A and Eα of the second reaction and reaction order n of the third reaction.