Analytical study on ignition time of PMMA exposed to time-decreasing thermal radiation using critical mass flux

This contribution addresses an analytical model to predict the ignition time of PMMA (Polymethyl methacrylate) subjected to a time-decreasing incident heat flux. Surface temperature, transient mass flux and ignition time of PMMA are thoroughly studied based on the exact solutions of in-depth temperature. Critical mass flux is utilized as the ignition criteria. An approximation methodology is suggested to simplify the unsolvable high order equations and deduce the explicit expressions of ignition time. A numerical model is employed to validate the capability of the developed model. The results show that no ignition occurs when the decreasing rate of heat flux increases larger than a critical value. The agreement of the transient mass flux between analytical and numerical models is good at high decreasing rate but turns worse as the decreasing rate declines. However, this enhanced discrepancy affects the ignition time prediction slightly. The inverse of the square root of the ignition time is linearly correlated with the decreasing rate of heat flux, and it becomes significantly sensitive to the decreasing rate when the decreasing rate approaching its critical value. Meanwhile, the value of critical mass flux has appreciable influence on the ignition time prediction.

While under time-increasing HF, Lamorlette 18 discussed an approach to determine the analytical solutions, including power law and polynomial HFs. Didomizio 10 experimentally and numerically studied the ignition of wood under fourth-order HF. Vermesi 12,13 investigated the ignition of PMMA and wood under parabolic HF by the FPA. Also, a numerical solver Gpyro 21 was employed to implement the simulation. Yang 14 and Ji 15 designed a linearly increasing HF in their tests to explore the ignition of wood species. An integral model was developed to analyze the experimental data. Later, Zhai 17 extended Yang and Ji's works to exponential HF. In all these studies, critical temperature and critical mass flux were commonly used for analytical and numerical models, respectively. Reszka 9 studied the ignition time of solid under time-dependent HF in forest fire, and found that the ignition time can be correlated with the total absorbed energy before ignition. Nonetheless, critical mass flux which is believed to be a more reasonable ignition criterion is barely studied under the time-varying HF in these models. Bilbao 24 examined the ignition behaviors of wood under time-decreasing HF in his experiment by turning off the power supply of the heater once the predetermined initial HF was reached. A numerical model was introduced in that study to estimate the ignition time. However, no analytical work was provided and the detailed information of the decreasing heat flux was avoided by using an average heat flux before ignition.
In this study, the ignition of PMMA under linearly and quadratically decreasing HFs is investigated analytically. Critical mass flux is utilized to obtain an explicit expression for ignition time with consideration of pyrolysis within solid. Also, transient mass flux is studied in this work. A previously developed numerical model 23 is employed to validate the proposed model. The capability of the numerical model has been validated through experimental measurements, including the surface temperatures and mass loss rates of several charring and non-charring polymers under constant incident heat flux, and good agreement is found. Furthermore, parametric study is implemented to analyze the dependence of ignition time prediction on ignition criterion.

theoretical Analysis
Considering a thermally thick PMMA imposed to a time-decreasing HF, the one-dimensional heat transfer in solid is illustrated in Fig. 1. Relatively low HF, lower than 80 kW/m 2 4 , is focused in this study and only surface absorption of radiation, corresponding to opaque materials, is utilized instead of in-depth absorption. This is because Jiang 7 , Vermesi 12 , Beaulieu 5 , Bal 3,41 , Delichatsios 4 and Boulet 40 found surface absorption dominates the heat absorption process even for translucent materials under low HF, whereas in-depth absorption plays a more important role as HF increases. Surface heat loss including convection and reradiation is neglected in this study for simplification. The time-decreasing HF is expressed as: ″ is the initial HF, a and b are rational constants. Defining a relative temperature: where T and T 0 are the transient and initial temperatures, respectively. The energy conservation equation in solid, initial and boundary conditions can be written as: , k is thermal conductivity, ρ is the density, C p is the specific heat and x is the spatial variable. Equation (3) can be further decomposed as: 1 2 θ θ θ = − Figure 1. Schematic of heat and mass transfer in PMMA exposed to a time-decreasing HF.
θ 1 and θ 2 are the solutions of constant and time-increasing HF scenarios, respectively. θ 2 was derived in ref. 19 . θ 1 is equivalent to θ 2 when = ″  a q 0 , = b 0. Consequently, the solution of θ can be obtained as: , Eq. (6) can be rewritten as: (7) is the transient in-depth temperature in solid. Some researchers tried to integrate the pyrolysis reaction into analytical models by invoking some approximations. For instance, Lautenberger 8 replaced the pyrolysis rate by a power law function, and Snegirev 42 used the Frank-Kamenetskii decomposition to simplify the Arrhenius function. However, all these works were conducted under constant HF. In this section, we extend Lautenberger's 8 study from constant to time-decreasing boundary condition. Both linearly and quadratically decreasing HFs are examined.
Linearly decreasing heat flux. The total transient mass flux of a thermally thick solid subjected to an incident HF can be expressed as 19 : where  m″ is mass flux, Z is the pre-exponential factor and T a is the activation temperature. As the approximation results in ref. 8 , the pyrolysis rate can be replaced by a power law function when the pyrolysis temperature is between 250 °C and 450 °C and the T a ranges from 10000 K to 30000 K: where A and B are constants, T 1 and T 2 are characteristic temperatures. For linearly decreasing HF, b = 1, Eq. (7) can be rewritten as: Noticing that Rearranging Eqs (12)(13)(14) and (9), the total mass flux at linearly decreasing HF condition, Eq. (8) can be expressed as: www.nature.com/scientificreports www.nature.com/scientificreports/ In order to further simplify Eq. (15), the same exponential approximation in ref. 19 where the excellent accuracy of this approximation was verified, is utilized as follow:  (16) and (17), Eq. (15) can be rewritten as:

Combining Eqs
Unfortunately, Eq. (18) cannot be integrated analytically. Another approximation is introduced to simplify Eq. (18): The reliability of this assumption is validated in Fig. 2. Although the approximate curves deviate more greatly as ξ increases, they agree well with the exact solutions near the high temperature surface where the pyrolysis reaction contributes the most of the mass flux. Accordingly, the transient mass flux can be expressed as:  20), it can be derived that:  (18) and (19).
⁎ is the equivalent ignition temperature under linearly decreasing HF. Equation (22) is an implicit expression because θ ig ,1 ⁎ actually is a function of t max,1 . However, a characteristic value of t max,1 = 20 s is used in Eq.
(23) to obtain an explicit result. When t max,1 changes from 10 s to 50 s, ( ) varies from 0.95 to 0.92 and this small discrepancy can be neglected. Solving Eq. (22), the ignition range of a can be calculated as: The ignition time at the critical condition can be expressed as: If a is larger than a cri,1 , no ignition is observed. When m m cri Consequently, the ignition time can be calculated based on Cardano's formula for cubic equation: (27) provides an exact solution for ignition time, it is a relatively complicated correlation. An approximate method is introduced here to simplify the ignition time expression. Noticing that when < a a cri,1 , the ignition time is slightly less than the value obtained under constant HF. A reasonable approximation is to substitute the t ig,1 in brackets of Eq. (26) for the average value of t ig,const and t ig,cri,1 . Combining Eq. (25) and the correlation under constant hear flux: (26) can be approximately rewritten as: Under quadratically decreasing HF, the mass flux can be expressed according to Eqs (7-9), (13) and (31): Based on Eqs (16) and (33), Eq. (32) can be simplified as: Being similar to the simplification process from Eqs (18) to (19), another approximation is invoked to simplify Eq. (34): The reliability of this assumption is validated in Fig. 3 and the agreement is good especially in the vicinity of the high temperature surface. After integrating Eq. (35), the mass flux can be consequently expressed as: The ignition range of a can be calculated as:

numerical Model
A previously developed numerical model 23 was used to validate the established analytical model and experimental results. Both surface and in-depth absorptions and their combination were considered in the numerical model. The energy and mass conservation equations and decomposition rate in solid are: v where r and τ denote the reflectivity, absorptivity of top surface, S v is the rate of volatiles generation in solid, ΔH v is the heat of decomposition, C g is specific heat of gas, m  ″ is the mass flux in the controlled volume, Z is pre-exponential factor, E is activation energy and R is ideal gas constant. 0 τ = denotes in-depth absorption and 1 r τ = − means surface absorption. When τ ∈ −r (0, 1 ), both exist. The initial and boundary conditions are: where L denotes the thickness of sample. Although temperature-dependent thermal parameters were considered in this numerical model, constant values at ambient temperature were used for comparison with analytical model. More detailed information about the numerical model can be found in ref. 23 . The simulation time, spatial and time steps are 500 s, 0.05 mm and 0.2 s, respectively. In order to verify the thermally thick assumption, 30 and 50 mm thick sample cases were also simulated. Little discrepancy of the in-depth temperature was found between the two cases. Therefore, 30 mm thickness was used in this study to save the simulation time.

Results and Discussion
PMMA is selected for computation to validate the proposed analytical model in this section. The thermophysics, chemical kinetics and other relevant parameters of PMMA are listed in Table 1.
Surface temperature. The exact solutions of surface temperature can be expressed based on Eqs (7), (13), (14) and (31) for linearly and quadratically decreasing HFs: The comparison of surface temperature of PMMA between the developed analytical model and simulation results is illustrated in Fig. 4. a cri is crucial for ignition and the ranges of a used in Fig. 4 are in the vicinity of a cri . The calculated values of a cri through critical mass loss rate under the initial HFs of 50 and 70 kW/m 2 are listed in Table 2. For fixed  ″ m cri , a cri,1 and a cri,2 are proportional to ″ q 0 3  and q 0 5  ″ , respectively, based on Eqs (24) and (40), and thus a cri increases with increasing initial HF. When a a cri ≤ , ignition occurs in air and the established flame would affect the surface temperature. However, the purpose of this section is to verify the capability of the analytical www.nature.com/scientificreports www.nature.com/scientificreports/ model. Therefore, the heating is assumed to be conducted in anaerobic environment and no ignition is examined in Fig. 4.
As expected, the surface temperature gets higher as a decreases. When pyrolysis is considered in numerical model. The exact analytical solutions fit the simulation results relatively well, but still some discrepancy is observed in Fig. 4. Especially in the range near the peak value, the analytical model deviates from the numerical results most. This divergence is caused by the gasification heat absorbed during pyrolysis in the numerical model. A small portion of the absorbed energy is consumed in the thermal degradation to evaporate the condensed solid. However, in the analytical model all the energy is assumed to heat the thermally inert material. Consequently, the analytical curves are slightly higher than the numerical ones. In the early and late stages, the lower surface temperature leading to a weaker pyrolysis rate results in a better agreement. The maximum discrepancy of surface temperature caused by the thermal degradation in Fig. 4 between the analytical and numerical models is less than 20 K. Transient mass flux. Transient mass flux reflects the generation rate of volatiles and determines the subsequent ignition. Figure 5 demonstrates the comparison of mass flux of PMMA under linearly and quadratically decreasing HFs without considering the ignition effect. Also, a cri listed in Table 1 can be calculated through Eqs (25) and (41). The value of a used in analytical and numerical models in Fig. 5 spans the ignition and no-ignition Ref.

48
Activation temperature, T a (K)  domains. The resultant peak values higher and lower than the critical mass flux, 2.42 g/m 2 s, imply the ignition and no-ignition scenarios, respectively. The developed analytical model fits the simulation results reasonably well except the peak values in low a cases. This disagreement is also attributed to the temperature related pyrolysis rate. Higher temperature resulting from lower a vaule would enhance the thermal degradation and consume more absorbed energy during the gasification. While with higher a values, the agreement is good. Although the disagreement of peak values is significant in the ignition cases, the prediction of ignition time is barely affected since the agreement at critical mass flux, 2.42 g/m 2 s, is good during the increasing phase of the curves in Fig. 5. Predictably, this distinction of peak value would induce the maximal error of ignition time in the critical case, a = a cri , which will also be demonstrated later.
Ignition time. Figure 6 shows the predicted ignition time of PMMA under linearly and quadratically decreasing HFs through analytical and numerical models. Under linearly decreasing HF, Fig. 6 (a), Eq. (27) provides good agreement with the numerical simulation by solving the cubic equation even though several approximations are invoked, including Eqs (9), (16), (17), (19) and (23). Without solving the cubic and quintic equations, the approximate analytical solutions, Eqs (30) and (43), provide linear correlations between − . t ig 0 5 and a for both linearly and quadratically decreasing HFs. In Fig. 6, the straight lines agree well with the numerical simulation until a cri . When a approaches a cri , t ig 0 5 − . decreases sharply. However, the ignition time at a cri can be attained through Eqs (25) and (41).

parametric Study of critical Mass flux
Critical mass flux is a determinant input parameter in analytical and numerical ignition models for specified materials. The objective of this section is to investigate the influence of ignition criteria on ignition time.
Although the critical mass flux of PMMA used in Section 3.3 is 2.42 g/m 2 s, it is hard to be measured accurately in the tests and it is sensitive to the experimental condition. Therefore, utilization of a valid range of this parameter may be a more practicable choice. The reported critical mass flux ranges of PMMA by Bal 3 and Vermesi 12 are 1.82-3.75 and 1.9-3.2 g/m 2 s, respectively. However, a larger range, 1-8 g/m 2 s is used in this study to examine its   (23) and (39), and the values are listed in Table 3 for the focused two types of HFs. T ig ⁎ is independent of  ″ q 0 and increases with ″ m cri  because it is proportional to m cri (23) and (39). All these values are lower than T ig , 655.1 K. Furthermore, the values of a cri of PMMA computed through Eqs (24) and (40) using critical mass flux are given in Table 4. Not surprisingly, for fixed initial HF a cri decreases with increasing ″ m cri  , indicating that lower a cri is prerequisite to heat the solid to higher temperature and generate more volatiles. a cri also increases with initial HF because a q cri,1 0 3 ∝ ″  and ∝ ″  a q cri,2 0 5 for linearly and quadratically decreasing HFs, respectively. Figure 7 shows the influence of m cri ″  on − . t ig 0 5 of PMMA with initial HFs of 50 and 70 kW/m 2 . The agreement between analytical and numerical results in this figure is not very good but acceptable because several additional approximations are employed when considering thermal degradation reaction within solid. As a getting closed to a cri , the discrepancy gets larger. This is caused by the fact that the induced error by invoking the approximations approaches a climax when a = a cri , which has been interpreted in the end of Section 3.2. As shown in Fig. 7, the critical mass flux has significant effect on the ignition time.
conclusions An approximate analytical model is established in this study to investigate the ignition of solids exposed to a time-decreasing incident HF. Critical mass flux is used in the model to calculate the surface temperature, transient mass flux and ignition time based on the exactly obtained in-depth temperature. Linearly and quadratically decreasing HFs are focused. An approximation methodology is proposed to simplify the complicated expressions during the derivation. The reliability of the developed model is verified by comparison with a numerical model employing PMMA as the reference material.
For time-decreasing HF,   q q at b 0 ″ = ″ − , a must be lower than a critical value, a cri , to ensure the ignition occurrence. This a cri is proportional to ″  q 0 3 and ″  q 0 5 for linearly and quadratically decreasing HFs, respectively. A linear dependence is found between − .     www.nature.com/scientificreports www.nature.com/scientificreports/ tion time under constant HF. An equivalent ignition temperature, ⁎ T ig , which is identical to T ig , is found to keep the ignition time expression similar to the one in critical temperature case. T ig ⁎ is independent of q 0  ″ but is proportional to ″ m cri B 1/


. The accuracy of the proposed model is good for high values of a and decreases as a gets lower when predicting the transient mass flux. However, this increased inaccuracy does not compromise the capability of ignition time prediction. The error induced by invoking the approximations gets larger as a approaching a cri . According to the parametric study, critical mass flux has great effect on the ignition time prediction.
Although only linearly and quadratically decreasing HFs are studied in this work, the derivation procedure and the approximation method can be extended to other time-decreasing HF scenarios. The main criticism of this work is the negligence of the surface heat loss. The temperature-dependent heat loss in boundary condition would greatly complicate the derivation, and thus more future studies are needed to solve this issue. Meanwhile, in-depth absorption of thermal radiation within infrared translucent solids and the corresponding heat transfer and ignition problems also deserve more attention.