This numerical study is conducted using the opposed flow velocity and temperature and the solid fuel thickness as parameters to investigate flame ignition and the subsequent downward flame spread behavior.
Notably, in the present simulation, a finite-length fuel slab is used, the ignition/combustion is in a two-dimensional wind tunnel, and both gas and solid phase radiations are considered. The ambient oxygen concentration in this model is fixed at 0.233.
The numerical calculation is initiated from a prescribed external radiant heat flux, qex, input on the solid fuel surface at time t =0. The profile of the incident radiation flux is a Gaussian distribution and its non-dimensional form can be expressed as
(
2)
maxexp x
Q
qex = −α (2.63) where Qmaxis the non-dimensional peak value of the external radiant heat flux, α is a Gaussian distribution shape factor and the x is the non-dimensional distance along the solid fuel surface. The solid fuel absorbs the external heat to raise its temperature gradually and pyrolyze the fuel vapor mixing the air to form the flammable mixture. The ignition occurs as soon as the gas phase temperature raises high enough to enhance the chemical reaction and the flame stars to propagate downward subsequently.
The finite difference equation for each variable is obtained by integrating the differential equation over the relative computational cell, associating with specified interpolation applied on the interface between the variable of two adjacent grid points. The detailed derivation of finite difference equations is carried out in Appendix A. The resulting finite difference equations are:
continuity equation: The energy, fuel and oxidizer species equations can be expressed as following general form of
φ
φφP nbφnb C
P a S
a =
∑
+ (2.67) where φ represents temperature, fuel and oxidizer mass fractions. Thesummations are applied over E, W, N and S points, the values of dependent variables u, v, φ and p are evaluated at (n+1) iteration and the expressions of coefficients such as A’s, a’s, S’s and C’s are also given in Appendix A.
The appearance of Equs. (2.64)-(2.67) seems to be linear, but it is not because the coefficients are also the function of dependent variable itself.
The numerical scheme utilizes the SIMPLE (Semi-Implicit Method for Pressure-Linked Equations) algorithm developed by Patankar (1980).
Since the momentum equations can be solved only when the pressure field is given or somehow estimated. An equation for pressure is needed to resolve the problem. The thinking of SIMPLE procedure is that unless the correct pressure field is employed, the resulting velocity will not satisfy the continuity equation. Therefore, a pressure correction equation is derived by linking the continuity equation. The derivation is outlined next.
Because the pressure field is unknown in the beginning, a P* is guessed in momentum equation. An imperfect velocity field based on guessed P* will be denoted by u* and v*. This crossed velocity field is resulted from the following equations.
u The continuity equation generally will not be satisfied by employing u*
and v* into it, instead a net non-zero mass source mp will be generated way of improving the guessed P* such that the resulting crossed velocity will progressively get closer to satisfy the continuity equation. First, a correct pressure is proposed that
p'
p
p= + + (2.71) where p' is called the pressure correction. Suppose that the true velocity components u and v respond to this pressure change in the following ways:
) Then, substitute the velocity components into the continuity equation using the above velocity-correction formulas (Eqs. (2.72) and (2.73)), we can obtain, after rearrangement, the following difference equation for p':
P
Once p' is obtained from Eq. (2.74), the velocities and pressure will be updated by p' through Eqs. (2.71)-(2.73). In the computation, as soon as the numerical results meet the criteria of convergence, the value of mP
will come out to be practically zero for all control volumes. Therefore,
' =0
p at all grid points will be acceptable solutions of Eq. (2.74) and the crossed velocities and pressure will be the correct velocities and pressure.
The present radiation model incorporates the subroutine RADCAL developed by Grosshandler (1993) to determine the gas absorption coefficient. Since the radiation subroutine is complicated and consumes much computing time, it is executed once after ten iterations in each time step. The model is solved with a marching time step. At each time step,
the gas and solid phase equations are solved separately. The solid phase equations are coupled with the gas phase equations through a mass and energy balance linkage. An iteration procedure is performed until all variables in the gas and solid phases converge to their respective acceptable criteria. For each grid point a residual R can be calculated as
p
Obviously, when the discretization equation is completely satisfied, R should be zero. In practice, this is impossible to satisfy. Because of this, the suitable convergence criterion should be selected that the largest value of R in the computational domain is less than a convergence criterion, ε. The value of ε is set to be 0.01 in this work. When the conditions of
ε
max <
R satisfied simultaneously for each dependent variable, the iteration procedure is stopped. Thereafter, the procedure moves to the next time step. Computations are carried out using a non-uniform mesh distribution as shown in the figure 2.2. The calculation performed with a non-uniform mesh according to the formula listed as below:
042 The smallest grid is 0.01 cm wide. Most of the grid points are clustered in an external radiative heating region to capture drastic variations in the flame; the grids then expand upstream and downstream.
The tests of the independence of the grid-size were conducted in advance and the results are shown in Table 2. According to the grid-independence test, a non-dimensional time step of Δt = 10 (equivalent to a real time of 0.02 s) and non-uniform grid dimensions of 290×95 were found to optimize the balance among resolution, computational time and memory space requirements. The time step that was selected in this work is much
smaller than those used in previous studies, such as 0.0548 s in Lin and Chen (2000) and 0.05 s in Wu and Chen (2004). Hence, the present computation is expected to be more accurate and suitable for examining gas phase ignition. In sum, the computational procedure is shown in figure 2.3 and the whole produce is expressed briefly as follows.
1. Read initial conditions for both gas and solid phase.
2. Solve solid phase conservation equations.
3. Combine interface conditions.
4. Guess the pressure field P*
5. Solve the momentum equations to obtain u*, v*. 6. Solve pressure correction equation to obtain P'. 7. Calculate P by adding P' to P*
8. Calculate u, v from the velocity correction formulas
9. Solve the discretization equations for other variables (temperature and concentration), return to step 2, and repeat the whole procedures until a converged solution is obtained.
10. Update initial conditions for all variables and march next time step.
The computational time associated with each case was approximately two days on a 2.8GHz Intel Pentium 4 PC at National Chiao Tung University.
Chapter 3
Results and discussion
As mentioned previously, this dissertation consists of three parts. In this chapter, the computation results as well as the detailed discussion for each part are given in the following sections. These results are also compared with corresponding experimental measurements and numerical predictions.
The method and apparatus used to conduct the experiment are briefly described as follows. The wind tunnel used in the experiment is 70 cm long and has 10*10 cm2 rectangular cross section. The test specimens are PMMA slabs, which are 30 cm long, 3 cm wide and have thicknesses of 0.82 cm and 1.74 cm, respectively. Each specimen is mounted on the groove in the test section, whose sides are covered with asbestos plates to minimize side effects. Remind that the solid fuel is assumed to be homogeneous that its compositions are uniform and its thickness is remained constant by assuming the flame spreads relatively fast enough that the fuel surface near the flame base remains approximately flat. The air, heated by the heater, is drawn into the entire test section and then flows over the specimen. The specimen is ignited using an electrically heated Ni-Cr wire, which is placed above the PMMA surface. A 15V A.C.
current is passed through the Ni-Cr wire. The current is cut off when the flame is ignited. The laser holographic interferometry and K- type thermocouples are utilized to measure the local gas temperature in the test section and the PMMA surface temperature, separately. The thermocouples are separated by 5 cm. A thermocouple signal is recorded using a multi-channel Yokogawa DA-2500 analyzing recorder and the
flame spread rate is determined by dividing 5 cm by the time recorded to have elapsed between the peaks on each of the two temperature traces.
3.1 Numerical study for downward flame spread over a finite-length PMMA slab with radiation effect in a two-dimensional wind tunnel
In this topic, several aspects of the original combustion model of Wu et al. (2003) are modified to investigate the downward flame spread over a finite-length PMMA slab with radiation effect in a two-dimensional wind tunnel. Parametric studies were conducted by changing the opposed flow velocity and temperature and the solid fuel thickness, in the same ranges as in the experiment of Pan (1999), to enable the results to be compared fairly with those of Pan (1999) and Wu et al. (2003). Notably, in the simulation of Wu et al. (2003), the fuel slab extends infinitely in both directions, the ignition/combustion is in an open atmosphere and both gas and solid radiations are neglected. The discrepancies between the previous combustion models of Wu et al. (2003) and present work are listed in the Table 3. The ambient oxygen concentration in the present model is fixed at 0.233. Table 4 present the physical data used in this study and the non-dimensional parameters are in Table 5.
In the initial state, the solid fuel of 298K is heated by the hotter opposed flow of 313K, 333K and 353K, respectively, and its temperature raises gradually. As soon as the solid fuel temperature reaches a steady one, an external heat flux stars to heat the solid fuel surface. Figure 3.1.1 shows the steady state temperature distributions along the PMMA surface before the external radiative heat flux is incident on the fuel under three opposed flow velocities u∞ =40cm/s,70cm /sand 100cm /sand opposed flow temperatures from 313K to 353K. The solid and dashed lines represent the Pan’s experimental (1999) and present numerical results, respectively. The experimental error in the temperature measurement is
K 2 .
±5 . The predictions are quite consistent with the measurements.
From this figure, it can be found that the differences decline as the flow velocity increases at fixed flow temperature, but increases with the flow temperature at a fixed flow velocity. The former result follows for the fact that the upstream boundary layer may not have been fully developed in the experimental test section, especially in the low velocity regime.
However, the flow in the numerical simulation is assumed to be fully developed and the boundary layer thickness is fixed. For the latter result, the temperature distribution in the experimental tunnel is expected to be less uniform as opposed flow temperature increases at fixed flow velocity.
In addition, the errors of the instruments used to measure temperature such as the thermocouple and data recorder, also increase with the opposed flow temperature.
Figure 3.1.2 displays the time history of the flame profiles from ignition to subsequent flame spread for u∞ =40cm/s , Ti =313K and
cm 82 .
=0
τ . The right half of the figure is the fuel and oxidizer mass fraction distributions whereas the left half of the figure shows the temperature contours and flow velocity vector distributions. At t=0, the incident radiation flux of the Gaussian distribution starts to heat the solid fuel. The solid fuel absorbs the heat and raises its temperature gradually.
While the solid fuel reaches the pyrolysis temperature, the solid fuel pyrolyzes the fuel vapor mixing with the air to form the flammable mixture, as demonstrated in figure 3.1.2(a). However, the flame is not ignited due to a small amount of fuel vapor and the low gas phase temperature.
Thereafter, the concentration of the flammable mixture increases continuity and the gas phase temperature raises high enough to enhance the chemical reaction as well. In this interval, the fuel vapor and air are now well premixed by the convection and diffusion. According to Nakabe et
al. (1994) and Ferkul and T’ien (1994), the gas phase ignition as occurring as soon as the dimensional reaction rate (ω&F ) reaches 10−4 g cm3s. The starting point of gas phase ignition occurs at t =14.12s, as shown in figure 3.1.2(b). The thermal plume downstream is longer because of the high opposed flow. An extra source for the gas phase temperature rise is from the chemical reaction in addition to the solid conduction.
After that, the chemical reaction rate increases sharply and releases much heat, resulting in the flame grows quickly within a very short period and the flame size reaches a maximum value at t=14.14s, as illustrated in figure 3.1.2(c). As following, the flame size and temperature reduce simultaneously, as displayed in figure 3.1.2(d). The flame gradually transfers from a premixed flame to a diffusion one. From the fuel and oxidizer mass fraction distributions in the figure 3.1.2(d), it can be seen that the fuel vapor and oxidizer is premixed in the flame front, whereas the fuel vapor and oxidizer is mixed by the diffusion in the flame downstream region. In order to sustain the flame itself, the flame front begins to extend upstream (downward direction) to pyrolyze the solid fuel to generate more fuel vapor to form the flammable mixture for support itself, as depicted from figure 3.1.2(e) to figure 3.1.2(g). After that, the steady flame spreads downward with the flame front.
Figure 3.1.3 plots the ignition delay time versus the opposed flow temperature at three opposed flow velocities u∞ =40cm/s, 70cm/s and 100cm/s and solid fuel thicknesses of 0.82 cm and 1.74 cm. Notably, the solid fuel length herein is finite. The values expressed by solid symbols are ignition delays for τ =0.82cm, whereas those represented by hollow symbols are for τ =1.74cm. This figure reveals that the ignition delay time increases with the opposed flow velocity, because the thermal boundary layer is thinner in regime of higher opposed flow velocity, in
which more produced fuel vapors are carried downstream, hindering the accumulation of fuel vapor near the solid fuel surface, increasing the formation time of the flammable mixture, delaying the ignition. However, at a fixed opposed flow velocity, the ignition delay falls as the opposed flow temperature increases, because a flow at higher temperature can heat the solid fuel more effectively, generating more fuel vapor, and forming the flammable mixture sooner, shortening the ignition delay time.
Moreover, the ignition delay for τ =1.74 cm is longer than that for 82
.
=0
τ cm at a fixed flow velocity and temperature. A thicker solid fuel has a greater thermal inertia (the ability of a material to conduct and store heat) and requires more energy to reach the ignition temperature, increasing the ignition delay. The experimental observations of Pan (1999) and Chen (1999) and the predictions of Wu et al. (2003) and Wu and Chen (2003) have confirmed these findings.
Figure 3.1.4 plots the ignition delay time as a function of the opposed flow temperature with and without the radiation effects and the solid fuel lengths are finite and infinite, respectively. The ignition delay times with radiation (cases a,b,e) are longer than that without radiation (cases c,d,f), which fact can be explained as follows. Figure 3.1.5 displays the distributions of heat fluxes along the solid fuel surface at the instant just before ignition (t =13.72s). qex and qc are the external input radiant heat flux with a Gaussian distribution and the conductive heat flux from the gas, respectively. qgr and qsr represent the gas phase radiation feedback to the solid fuel and the radiation heat loss from the solid fuel to the ambient.
The sum of total heat fluxes is the net heat flux on the surface of the solid fuel, qnet. A positive value means that the solid fuel gains energy from the gas phase and a negative value represents the loss of heat from the surface of the solid fuel. It shows that the magnitude of gas phase radiation
changes insignificantly and can be negligible because the flame ignites in a high speed flow regime. Furthermore, this figure demonstrates that the net heat flux, qnet, near the origin is lower when the radiation effect is considered. The solid fuel must receive more energy for pyrolysis to produce fuel vapor, spending more time to form the flammable mixture.
Hence, the ignition occurs more slowly with radiation. The inset in figure 3.1.5 plots the distributions of solid fuel temperature and the density contours just before ignition. It reveals that the solid fuel temperature is lower with radiation, reducing the magnitude of fuel vapor. Furthermore, the differences between the ignition delay time in the figure 3.1.4 increase with the opposed flow temperature (cases a,c or cases b,d) and decrease with the opposed flow velocity (cases e,f). The former result follows for the fact that qsr is proportional to Ts4. The solid fuel surface temperature is higher when it is immersed in a hotter flow, so the heat lost by radiation from the solid fuel surface to the surroundings is increased. The later result is because that the convection effect is more effective at higher flow velocity. In other words, the radiation is less important in the high velocity regime. Additionally, whether the solid fuel length is finite or infinite does not significantly influence the ignition delay (cases a,b or cases c,d). This is because the ignition delay time is dominated by the opposed flow velocity and temperature and the solid fuel thickness, but not the solid fuel length.
Figure 3.1.6 plots the pyrolysis front positions as a function of time at cm
82 .
=0
τ , Ti =333K and u∞ =70cm/s. The pyrolysis front position is defined as the first upstream position of ρs =0.99. The steady flame spread rate can be obtained from the slope of a best fit line that passes through the pyrolysis front positions. Computations are carried out by using a non-uniform mesh distribution. Most of the grid points are
clustered around the external radiative heating region (x=0) and then they are expanded upstream and downstream. Hence, the flame displacement shows the stair-step pattern when the flame spreads far away from the origin. However, it does not alter the expected constant flame spread rate, indicated by that the straight line passes through these stair-step points.
The steady flame spread rate at various opposed flow velocities and temperatures can also be determined by this way, and they are presented in the next figures.
Figures 3.1.7 and 3.1.8 present the flame spread rate versus the opposed flow temperature at three opposed flow velocities for solid fuel thicknesses of 0.82 cm and 1.74 cm, respectively. The symbols represent the data measured by Pan (1999), whereas the solid and dashed lines indicate the results simulated herein and the prediction of Wu et al. (2003), respectively. There are three computed results in each line. The conditions, such as the opposed flow velocity and temperature, used in the computations are all the same as those used in experiments. These two figures indicate that the flame spread rate increases with the opposed flow temperature at a fixed opposed flow velocity. The flame spread rate falls as the opposed flow velocity increases at a fixed opposed flow temperature.
Figures 3.1.7 and 3.1.8 present the flame spread rate versus the opposed flow temperature at three opposed flow velocities for solid fuel thicknesses of 0.82 cm and 1.74 cm, respectively. The symbols represent the data measured by Pan (1999), whereas the solid and dashed lines indicate the results simulated herein and the prediction of Wu et al. (2003), respectively. There are three computed results in each line. The conditions, such as the opposed flow velocity and temperature, used in the computations are all the same as those used in experiments. These two figures indicate that the flame spread rate increases with the opposed flow temperature at a fixed opposed flow velocity. The flame spread rate falls as the opposed flow velocity increases at a fixed opposed flow temperature.