Differential Diffusion Effects in Numerical Simulations on Smoke Toxicity Evaluation

With the increasing frequency of fire caused by construction materials, smoke toxicity evaluation plays a key role in related fields. Numerical simulation has become a popular method to predict the toxicity of smoke. A computational study of differential diffusion effects on smoke toxicity evaluation is proposed in this study. Further, an effective Lewis number model derived from the Reynolds-averaged form of the Navier-stokes (RANS) transport equations is proposed in turbulent flames. The accuracy of the study is illustrated for a polyurethane foam fire in a 1/5 scale vertical shaft. The temperature and the concentrations of smoke composition are mainly discussed. From the comparison of the calculations with the direct numerical simulations (DNS) data it is observed that the temperature and mass fractions of species agree well with the DNS data when differential diffusion effects are taken into account. On the other hand, these numerical results are overestimated if differential diffusion effects are neglected. The FED values indicate that differential diffusion has a strong influence on smoke toxicity evaluation when using N-Gas


Introduction
The number of high-rise buildings is increasing rapidly with social and economic development and urbanization in recent years.Meanwhile, the demand for construction materials has been greatly increased because of the large number of architectures.Unfortunately, most of the construction materials will bring serious fire risks.According to previous researches, the frequency of fire caused by building materials is increasing year after year [1,2].Toxic smoke was identified to be the most hazardous factor to human beings in the fire caused by construction materials [3,4,5].Apparently, smoke toxicity evaluation plays a key role in related fields.
At present, there are several methods [6] for testing the toxicity of flue gas: (1) small scale test method.This method does not have a unified standard.In China, the smoke production principle follows the German standard [7]; (2) animal infected test method.Animal infected test was mainly used in the past.Usually there are several measurements, LC50 (Lethal Concentration), IC50 (Incapacitation Concentration), EC50 (Effete Concentration) and so on, and LC50 is the most commonly used [8,9]; (3) component analysis method.Component analysis is mainly based on different characteristics of the fire smoke components.Most of the methods can only analysis the intermittent sampling process, not the whole combustion process [10].The final toxicity values can be obtained by comprehensive evaluation after the detection of composition and concentration of gas by the above method.The evaluation models of flue gas toxicity mainly have N-Gas, FED, FEC and TGAS model [11].N-Gas model is considered to be one of the most precise and scientific one [12].There is an increasing tendency to use numerical simulation method instead of experimental method to test the toxicity of gas to prevent the use of animal testing for all except some very limited applications [13].Therefore, it is very necessary to predict the composition and concentration of gas accurately.
Nowadays, CFD (Computational Fluid Dynamics) simulation has become a very popular and convenient method to research building fires, especially to predict gas concentrations.As one of the CFD software, FDS (Fire Dynamics Simulator) is widely used to solve fire cases.In FDS, large eddy simulation is a more practical method to simulate building fires.But large eddy simulation adopts the assumptions of equal thermal and mass diffusivities, leading to the unity Lewis number for all combustion products of materials [14].Unfortunately, species found in industrial applications rarely have unity Lewis numbers [15].Actually, species and heat locally redistribute (this phenomenon is referred to as differential diffusion [16]) when Lewis number is not the unity [17].Moreover, the species concentrations and temperatures are quite different when differential diffusion effects are taken into account [18].
In this paper, FDS is chosen to predict the gas composition and concentration in a fire and a study on effects of differential diffusion on smoke toxicity is implemented by modifying FDS code.The following content are mainly divided into four parts.First, governing equations and computational methodology are introduced.The open source program is compiled after modifying the thermal and mass diffusivities formulas.In addition an effective Lewis number model derived from the Reynolds-averaged form of the Navier-stokes (RANS) transport equations is proposed.Second, a 1/5 scale vertical shaft model [19] is used as the physical model.And polyurethane is chosen to be the combustion material in the study.To verify the work, temperature and species concentrations will be compared with the data calculated by direct numerical simulations (DNS) which is considered as "Numerical experiment".Next, N-Gas model is used to evaluate the toxicity of smoke and values of FED are calculated to investigate the influence of differential diffusion on smoke toxicity evaluation.Finally, some conclusions obtained from the study are presented.

Methodology
In this section, the practical CFD tool FDS was used.FDS is a free and open source code developed by the National Institute of Standards and Technology in the USA (NIST).It solves the Navier-Stokes equations for low Mach number flows with an emphasis on smoke and heat transport from fires [14].The effective species Lewis numbers are obtained derived from the species and energy transport equations.And an effective Lewis number model which depends on turbulent Reynolds numbers is proposed.

Transport equation
The transport equation for species i can be written as where ρ is the mixture density, u is the velocity vector, Y i , D i and m i are the mass fraction, the mixture-averaged diffusivity and the production rate of species i, respectively.
Neglecting the temporal variations in pressure, the Low Mach number transport equation approximation for energy can be written as [14] where h s is the sensible enthalpy, T is the temperature, α is the mixture thermal diffusivity, h f is the heat fluxes for diffusive and radiative, q c c c is the heat release rate.And q c c c can be expressed as where is the heat of formation of species i.
Fig. 1 shows the contribution from each term in Eq. ( 2).Each contribution is normalized by the largest absolute value of all contributions.Obviously, the last term on the right side of Eq. ( 2) has negligible contribution in comparison to the other terms.The energy equation can be simplified as The use of the Reynolds-averaged form of the Navierstokes equation (RANS) in CFD models allowed us to obtain the information about the ensemble average of the transported quantity.Hence, based on Eq. ( 1) and Eq. ( 4), the Reynolds-averaged species and energy equations can be written as (after applying the Favre filter) where D M and D T denote the Reynolds-averaged and the Favre average, respectively.D M and D T are the eddy diffusivity for species and thermal.From Eq. ( 5) and Eq. ( 6), a formula for the effective species Lewis numbers can be obtained , where D mass and D heat represents the mass diffusivity and the thermal diffusivity, respectively.In LES, the mass diffusivity and the thermal diffusivity are taken from the molecular and turbulent viscosities divided by the turbulent Schmidt number Sc t and turbulent Prandtl number Pr t , respectively.And Pr t and Sc t are assumed to be constant for a given scenario and the default value are 0.5 for both [14] In the present study, in order to calculate the unclosed terms in Eq. ( 7), the definition of the mixture-average binary Fickian diffusion [20] is used.The species mixture-averaged diffusivity can be written as where N is the total number of species, W is the mean molecular weight of the mixture, X j and W j are the mole fraction and the molecular weight of the jth species, respectively.D ji is the binary diffusion coefficient of species j in species i and is defined as following: where k B is the Boltzmann constant, m ji is the reduced molecular mass for the (j,i) species pair, P is the pressure, σ ji is the reduced collision diameter, and Ω (1,1)* is the collision integral.The detailed information of the binary diffusion coefficient can be found in the Ref. [21].
In Eq. ( 9), the round-off is accumulated in roughly the same way in both the numerator and denominator, and thus the formula is undefined if the mixture is exactly a pure species.So a small number ε is added to the actual mole fraction to overcome this difficulty.We assume in Eq. ( 9) that where ε equals to 10 -12 .
As for the mixture thermal diffusivity, it can be computed from the thermal conductivity which is calculated by a combination averaging formula: where λ i is the thermal conductivity of the ith species.All the thermodynamic properties and the transport properties used to calculate the diffusivity can be obtained from the JANAF tables [22].The method mentioned above which is used to solve molecular diffusion problems has been successfully applied in different researches such as Ref. [20].
Using the Deardorff's model [14] to calculate the remaining two unclosed terms in Eq. ( 7), the eddy diffusivity can be written as where ∆ is the integral length scale, k sgs is the kinetic energy, and C υ is the model constant which is set to the literature value 0.1.
From the previously literature, it is well known that diffusion of temperature and species is dominated by turbulence level [23].Hence, a choice of the effective Lewis numbers depend on the turbulent Reynolds number is made.From Eq. ( 7), the effective Lewis number can be obtained for each species where Le i is the species molecular Lewis number.
Using the definition of the mixture Prandtl number Pr and the turbulent Reynolds number Re T , the ratio of the eddy to the thermal diffusivities and to the mass diffusivities can be written as where , Pc is the root-mean-square velocity fluctuation which can be written as is the kinetic energy, l is the characteristic length, and υ is the kinematic viscosity [24].
So, the species effective Lewis number becomes where It is difficult to evaluate the parameters C 1 and C 2 exactly, but an estimate of them can be provided.In FDS, the default value of turbulent Prandtl number and turbulent Schmidt number are 0.5 for both and Pr is set to 0.7 [14].From the discussion above, the two parameters have the same value 0.14 approximately.Hence, the species effective Lewis number can be written as In Eq. ( 19), Le i and Re T are the inputs.Le i can be obtained from literature survey, Re T is a known quantity from the simulations and C T here is used as a fitting coefficient and should be modified in specific cases.
A 1/5 scale vertical shaft model which was first chosen to research about the smoke in vertical shaft by Marshall [25] was considered.In the scaling model, there are five floors, each floor is 0.66 meters.The doors on the ground floor and the top floor are open.The height of the door is 0.3 m and 0.44 m, respectively.
A geometrically simplified model (2D) is illustrated in Fig. 2. The polyurethane foam (polyurethane is defined by the chemical formula C 25 H 42 O 6 N 2 in this study [26]) is chosen as the material and is set in the corner of the fire compartment.And there are 25 thermocouples in the shaft, they are set every 0.13 m and the first one is set at 0.08 m.
To demonstrate that the simulations are sufficiently well-resolved, mesh sensitivity study is presented for all simulations.In the mesh resolution study, a relatively coarse mesh is used first, and then gradually refines the mesh until there are no appreciable differences in our results.In the case of large eddy simulations, three different grids are used, grid 1 (1 cm), grid 2 (2.5 cm) and grid 3 (5 cm).
Fig. 3 shows the average temperature in the vertical shaft simulated with LES.Grid 1 and grid 2 depart immediately from grid 3, indicating that 5cm is insufficient in this case, but are in close agreement with each other.The turbulent nature of the flow makes it inherently unstable, and there are no appreciable differences between grid 1 and grid 2 demonstrating that the simulation is being well-resolved with grid 2. And 5 mm was selected for DNS simulation after doing the same mesh sensitivity study.In this section, the fire described above is simulated with the effective and the unity species Lewis numbers.The temperature and the main species mass fractions are compared with the DNS data in validation.
In order to observe the comparison results, mean results for temperature and species mass fractions on the center line of the shaft are presented in Fig. 4. The DNS data are indicated with symbols while the results with effective Lewis numbers or unity Lewis numbers are presented with solid and dashed lines, respectively.
Obviously, the temperature profile is well captured by the simulations with the effective Lewis numbers and the same applies for the species mass fractions profiles.In Fig. 4, it is observed that the temperatures and mass fractions of CO 2 , CO and H 2 O are overestimated while the mass fractions of C 25 H 42 O 6 N 2 and O 2 are underestimated when differential diffusion effects are not taken into account.Moreover, the peak temperature simulated with the effective Lewis numbers is at the position of about z = 268 cm, coinciding with the maximum values for CO 2 , CO and H 2 O mass fractions.On the other hand, if differential diffusion effects are neglected, the temperature distribution is not well captured which directly leads to the discrepancies in the species mass fractions.From the discussion above, species concentration is influenced by temperature, especially for CO 2 , CO and H 2 O.In reference [27], the conditional means of species concentrations with regard to temperature <Yi│T> for average of relatively steady flame are calculated.Fig. 5 shows the conditional mean mass fractions of CO 2 , CO and H 2 O as a function of temperature at the location of z = 268 cm where temperature reached the peak value.The DNS data are indicated with solid lines while the results with effective Lewis numbers or unity Lewis numbers are presented with dashed and dotted lines, respectively.
As mentioned earlier, the simulation results can be more accurately captured by LES method with an appropriate change in Lewis numbers of different species.Compare the data of these three groups, it is observed that the concentration of combustion products will reach a higher value at the same temperature when simulated with the unity Lewis number.In the case of differential diffusion effects are taken into account, there are also

Results and discussion
exist some differences when compared results with DNS data.To claim the errors, the L2-norm errors of the species conditional mean mass fractions between DNS and LES are calculated.The error can be written as (take carbon dioxide as an example): where ) ( from LES at T = T i with the effective Lewis numbers.Fig. 6 shows the L2-norm errors between the DNS and LES (with unity or effective Lewis numbers).Obviously, the errors between DNS and LES with effective Lewis numbers are much smaller than that with unity Lewis numbers for species in each case.And all the errors are less than 10% which is acceptable to engineering application [28].When unity Lewis number is used, the errors of the reactants are bigger than 10%, especially for oxygen.In this sub-section, N-Gas model was used to evaluate the toxicity of the smoke.The 3-gas value can be written as [29]: where [X i ] represents the concentration of species  With the emergence of new building materials and the use of existing materials, smoke toxicity evaluation becomes one of the most important topics in related fields.
There is an increasing tendency to use numerical simulation methods to predict the composition and concentration of smoke which are needed for toxicity evaluation model.Large eddy simulation is considered to be one of the most potential methods to simulate building fires.Extending the ability of LES method to predict the smoke concentration accurately becomes very meaningful.
In the present study, the effective Lewis numbers of species derived from the species and energy transport equations have been presented to incorporate differential diffusion in CFD simulations.The study was verified by simulating a polyurethane fire in a 1/5 scale vertical shaft and the numerical results were compared with the data derived from DNS. Temperature and mass fractions of CO 2 , CO and H 2 O are overestimated while the mass fractions of C 25 H 42 O 6 N 2 and O 2 are underestimated when simulated with unity Lewis number.On the other hand, temperature and species mass fractions are in well agreement with DNS data when differential diffusion effects are taken into account.Values of FED are calculated based on the data simulated with the effective or unity Lewis numbers using N-Gas model.Results indicate that the prediction model will overestimate smoke toxicity when differential diffusion effects are neglected.Moreover, a model derived from the transport equations for the effective Lewis number which depends on turbulent Reynolds numbers is proposed.The model is found to be applicable in turbulent flames.More data are needed to fully evaluate the model, for example, flames with Lewis number of fuel less than unity shall to be compared to the model.Future work will be concentrated on the fully validation and the applications of the model.
This article provides a method for CFD modeling to capture the behavior of smoke in building fires more accurately with an appropriate change in the Lewis number of species.The modified CFD model can be used to predict the combustion products of materials in building fire and provide more accurate data for the gas toxicity evaluation.

)Figure 1 .
Figure 1.Contribution from each term in the energy equation as a function of temperature.

Figure 3 .
Figure 3. Average temperature in the vertical shaft.

Figure 4 .
Figure 4. Temperature and species mass fractions at the center line of shaft.DNS data: symbols; effective Lewis numbers: solid lines; unity Lewis numbers: dashed lines.

Figure 5 .
Figure 5. Conditional mean species mass fraction profiles as a function of temperature.DNS data: symbols; effective Lewis numbers: solid lines; unity Lewis numbers: dashed lines.

Figure 6 .
Figure 6.L2-norm errors for species between DNS and LES with the effective or unity Lewis numbers.
X i (ppm); LC 50 (O 2 ) represents the lethal concentration of oxygen (%), and it is taken as 5.4% here; based on the empirical values, m = -18, b = 122000.From the comparison of the results of animal experiments with the N-Gas model it is shown that 50% animals will die when the value of FED is equal to 1. Fig. 7 shows the FED values at the center line of shaft.Apparently, the values of FED are different with the effective or unity Lewis numbers.When differential diffusion effects are not considered, values of FED are overestimated and the largest differences are observed at z ≥ 200 cm where fire smoke toxicity has reached the critical value of death.On the other hand, the values of FED are always less than 1 by the simulation with the effective Lewis numbers.

Figure 7 .
Figure 7. FED values at the center line of shaft.Effective Lewis numbers: solid lines; unity Lewis numbers: dashed lines.