Available online at www.sciencedirect.com

ScienceDirect Procedía

Engineering

Procedía Engineering 62 (2013) 13 - 27 =

www.elsevier.com/locate/procedia

The 9th Asia-Oceania Symposium on Fire Science and Technology

Mechanism of buoyant turbulent diffusion flames

John L. de Ris*

FM Global Research Division, Norwood, MA, USA

Abstract

The paper addresses the basic combustion mechanism of buoyant turbulent diffusion flames with emphasis on flame generated turbulence. The vorticity equation and Kelvin's theorem show that the buoyant generation of vorticity and its dissipation within the flamelets drives the combustion of fuel and oxidant in buoyant turbulent flames, characteristic of fires. It is an example of a Rayleigh-Taylor instability mechanism characterized by a constant Rayleigh number. The mechanism explains why measured volumetric heat release rates and radiant fractions are constants and independent of the overall fire size. It conforms to the fundamental requirement that the volumetric heat release be proportional to p4/3 g2/3. The paper presents the cascades of vorticity as well as rotational and translational turbulent kinetic energies. A set of modified LES equations using the EDC concept are suggested for simulating buoyant turbulent diffusion flames that should yield mesh-size independent solutions.

© 2013 International Association forFire SafetyScience.Published by Elsevier Ltd.AllRightsReserved Selectionandpeer-review under responsibility of the Asian-Oceania Association of Fire Science and Technology

Keywords: Buoyant turbulent diffusion flames; Rayleigh-Taylor Instability; LES; Vorticity

Nomenclature Greek.

a f Flame area per unit volume [1/m] a Thermal diffusivity [m2/s]

CP Specific heat [kJ/kgK] ka Wave-number 1/ A [1/m]

D Dissipation [kg/m3s2], Diffusivity [m2/s] A Flamelet size [m]

g Gravitational acceleration [m/s2] P Density [kg/m3]

AHC Heat of combustion [kJ/g] X Fraction [--]

h Specific enthalpy [kJ/g] V Kinematic viscosity [m2/s]

k Turbulent kinetic energy [m2/s2] T Time constant [s] or Stress [kg/ms2]

L Length [m] O Vorticity [1/s ]

m m Mass generation rate [kg/m3s] O Stress tensor [kg/m2s2]

P Volumetric production of TKE [kg/ms3] A Eddy or Mesh size [m]

Prf Turbulent Prandtl number [--] Subscripts

P Pressure [kg/ms2] F Fuel

Q Heat release rate [kW] f Flame

qch Volumetric heat release rate [kW/m3] Mix Mixing

R Radius [m] O Oxidant

S, S Stoich.oxidant to fuel mass ratio [--], Strain rate [s-1] P Product

* Corresponding author. Tel.: +1 508 698 0648; fax: +1 781 255 4025. E-mail address: john.deris@verizon.net

1877-7058 © 2013 International Association for Fire Safety Science. Published by Elsevier Ltd. All Rights Reserved Selection and peer-review under responsibility of the Asian-Oceania Association of Fire Science and Technology doi: 10.1016/j.proeng.2013.08.040

ELSEVIER

T ft Torque per unit volume [kg/s2m] R Radiant, Rotational

T CA Dispersion of vorticity [kg/m3s2] 0 Ambient

T Temperature [K], Transfer TKE [kg/ms3] T Translational

t Time [s] B Buoyant

u, v Velocity [m/s] ch Chemical

x, z Space coordinates [m] A Eddy or Mesh size

Y Mass fraction [-- ] t Turbulent

CO Vorticity

V Viscous

1. Introduction

In 1956 Morton, Taylor and Turner [1] successfully modeled buoyant turbulent plumes by applying G. I. Taylor's entrapment concept that air being drawn into a thermal plume at a particular height is proportional to the characteristic plume velocity at that height. Turbulent shear forces draw in the surrounding air by transferring some of the plume momentum. The proportionality is an empirical constant that is independent of the size of the fire, but does depend on the type of turbulent flow —whether it is axisymmetric or planar, or forced or buoyant; and hence is not truly a fundamental quantity. In 1963 Thomas [2] used the entrainment concept to model visible flame heights of fires by showing that the amount of air entrained into a fire plume up to the flame height is proportional to the stoichiometric air requirement of the burning fuel. Soon thereafter, Steward [3] provided a more detailed combustion model for line fires by assuming that only a certain fraction of the entrained air participates in the combustion. In agreement with experiment the model predicted flame heights increase with the 2/3 power of fire heat release per unit length. He similarly modeled axisymmetric flames [4] showing flame heights increasing with 2/5th power of the fire heat release rate. Then in 1979 McCaffrey [5] observed that the flames burn continuously up to about 40% of the flame height and then intermittently up to the time-averaged flame height. The fire heat release occurs principally in the lower continuous flaming region. In 1983, Heskestad [6] correlated of visible flame heights measured by many investigators. In 1985 both Delichatsios et al. [7] and Zukoski et al. [8] measured the rate of entrainment of air into fire plume verses height. They showed that the total entrainment up to the visible flame height is as much as ten times the stoichiometric air requirement of the fire. The above understanding of buoyant turbulent diffusion flames was based, almost entirely, on the concept of entrainment. There was little, if any, understanding of the subsequent processes as to how, the air, following entrainment, mixes with the fuel and ultimately combines with the fuel in a diffusion flame.

In contrast, at that time, there was much better understanding of momentum-dominated turbulent diffusion flames. The turbulent kinetic energy supplied externally at large-scale cascades downward to smaller and smaller scale eddies until it eventually dissipates by the action of viscosity at the Kolmogorov micro-scale where most of the heat is released. At the same time progress was made in understanding radiation from fires. In 1985 Markstein [9] showed that the radiant emission from buoyant turbulent diffusion flames is a constant fraction of the fire heat release rate (see Fig. 1a). In the same study,

Fig. 1. (a) Radiative power vs. heat release rate of buoyant turbulent diffusion flame; (b) heat release rate of a pool fire vs. photographically measured flame volumes.

he showed the radiant fraction is correlated by the smoke-point flame height of the laminar flame of the fuel. Additionally, in a study of laminar diffusion flames, Markstein et al. [10] showed that the chemical time for soot formation in diffusion

flames is proportional to the smoke-point flame height of the fuel. This means that the laminar smoke-point flame height of a fUel, which is governed by the flame chemistry, is the principal fUel property controlling the radiant emission from hydrocarbon diffusion flames. During the same period Delichatsios et al. [11] went further and established that the radiant emission from momentum-dominated diffusion flames is correlated by the ratio of the Kolmogorov micro-scale flow time to the chemical time for soot formation. These were all very important advances. As discussed below, the understanding of radiation is helpful in understanding buoyant turbulent combustion.

Theories for buoyant combustion, of that era (i.e. 1980's), presumed a close coupling of combustion with entrainment.

Those theories suggested that the volumetric heat release rate q"h ~ Q/if ~ Q"15 decreased with the minus one-fifth power

of the fire heat-release-rate (or equivalently, the minus one-sixth power of the flame volume). To test this supposition Orloff and de Ris [12] took photographs of pool fire flames to infer their volumes. The average volumes, shown in Fig. 1b, showed to our surprise, that the volumetric heat release rate is constant within the range of the data.

32. 6 kV/m 1

= 1.215 1.525

0.5 1.0

Height Above Burner, m

Axisymetric_Buoyant Flame

2a corr

0.2 0.4 0.6

z =(z - zo )/ lr

Free-Burning Planar Flame

Propylene 32. 6 kV/m

Free-Burning Planar Flame

Johnson Fit n = 2.349, m -e = -0.1434 s =

0.5 1.0

Height Above Burner, tn

= 1.902 1.703

32.6 kV/* 49.5 Mf/.

67.7 kW/m 87.1 kV/«

107.7 kV/m

129.8 M/* 153.1 bM/m 178.0 kV/m

Planar Flame Against Wall Johnson Fit n = 1.061, m = 0.925 e = 0.0428 s = 0.966

Planar Flame Against Wall

153.2 k¥/m

0,2 0.4 0.6 0.6

z =(z - z o ) / lr

0.2 0.4 0.6

z =(z -Zo ) / Lr

Fig. 2. Vertical distributions of radiant emission from buoyant turbulent diffusion flames showing the radiant emission per unit volume independent of fire size: (2a) axisymmetric flames - data & correlation; (2b) free-burning planar flames from slot burner - data & correlation; (2c) planar flame from slot burner against cold wall- correlation & non-correlation.

Figure 2 shows Markstein et al.'s [13] subsequent scanning radiometer measurements of buoyant turbulent diffusion flames, which show that the radiant emission per unit volume is independent of the overall fire size. Since the radiant fraction is also constant for these flames, their volumetric heat release rate must similarly be independent of fire size. More specifically, correlation of the effective radiative flame heights of axisymmetric flames, LR ~ Q1/3, shown in Fig. 2a, scale with the 1/3 power of heat release rate, this implies that q"'h ~ Q>/~ const. Similarly, correlation of the effective radiative flame heights of planar (i.e. line fire) flames, IR ~ QQ,in , shown in Fig.'s 2b and 2c scale with the % power of heat release per unit length, q'"h ~ QQ'/(.R ~ const. Finally, the bottom right hand Fig. 2c shows that the earlier predictions, If ~ Q,2/3

based on entrainment- visible flame height, do not correlate the wall fire data. The three figures within Fig. 2, showing successful correlations (i.e. 2a_corr, 2b_corr and 2c_corr) provide best curve-fit parameters for Johnson SB

function,SB (n,s,e,m) = {nsj42k (Z -e)(s-Z + e)}exp^-0.5(m + nln[(Z -e)/(s -Z + e)^[14]. 2. Fluid dynamics

There is increasing interest in vortical structures in fires. In 1993 Cetagen et al's [15] experiments established that the puffing phenomenon of pool fires comes from the buoyant flow rising from the pool fire surface causing the formation and periodic shedding of toroidal vortices. Tieszen et al. [16] describes many observations of fluid dynamic phenomena occurring in large-scale pool fires -especially fires in the presence of wind. The baroclinic generation of large toroidal vortices is particularly prominent for large-scale pool fires.

Instead, the present paper focuses on the mechanism of combustion occurring at the turbulent micro-scale, which controls the radiant emission and volumetric heat release rate. Vorticity plays the key role. The vorticity equation describes its generation, transport and dissipation. Kelvin's theorem describes the conservation of vorticity as it moves with the fluid. Kelvin's theorem explains the difference between momentum-controlled and buoyancy-controlled turbulent diffusion flames. The cascade of vorticity provides the key concept of this paper. More specifically, the generation of vorticity increases with decreasing scale and peaks at the micro-scale near where the generation and viscous dissipation come into balance. For scales greater than the micro-scale, the generation of vorticity exceeds its dissipation making it a dynamic process. The micro-scale occurs at a constant Rayleigh number, gApf Xf Ipfvfaf , having both characteristic length, Xf and time scale, t f independent of fire size. The mechanism is an example of Rayleigh-Taylor Instability whose disturbances generally originate at some small scale, grow, and eventually encompass the entire flow. The importance of flame generated turbulence increases with fire size, or to be more specific, increases with the ratio of overall flame height, IR , to Xf .

2.1. Vorticity generation:

The equation for momentum provides

where p is the total pressure, g is the downward acceleration and O is the viscous stress tensor that in indicial notation is

PD = p Ik +(v*v)v = -vp + Pg+O

dTji d \ 2 „ dvk fdVi dVj)

—— =---/dSii—- +1! —L + —J-

dx. dx. 3 1 dx, dx. dx.

J 1 k V 1 i /

Recalling that Vx^A = yVx A + Vy x A , the curl of Eq. (1) is

Fig. 3. (a) Baroclinic generation of vorticity in a fire; (b) Kelvin's theorem for the evolution of vorticity as it is carried by the flow. because both the curl of a gradient and curl of a constant vector are zero. Upon rearranging, the equation becomes

p V x-= p V x

+ (V «V )v

= Vpx | g - I + V x O

The first term on the right is the baroclinic generation of vorticity. Fig. 3 a shows how gravity [or other acceleration] exerts a torque by acting on local density gradients causing a local rotation, or generation of vorticity, ( = V x V. The net effect is to push the fluid up where the density is least. For buoyant flames, the local acceleration Dv I Dt in Eq. (4) will most often be in an upward direction, often exceeding the absolute magnitude of |g| . The acceleration is most intense when the local density is small. In contrast, the gravitational acceleration g is in the downward direction, so we expect the difference of these two acceleration vectors, [g - Dv/Dt] to be in a net downward direction with magnitude typically greater than |g| . Assuming constant viscosity, ¡i , the final term on the right describes diffusion of vorticity [after noting that the curl of a gradient is zero]

V x o = ¡v (+¡Vx-Vfvvy

One can express the left hand side of Eq. (4) in terms of vorticity using the following two vector identities

(v*V)v = -2 V (v *v )- v x V x V = 2 V (v *v )- v x ( V x v x ( = v(() - ((V*v) + ((*V)v - (v*V)((

Remembering that the divergence of a curl is zero, one obtains the general vorticity equation^

+ (v «V ) + (o (V«V )- (c5 «V ) V

= Vpx| g - D I + ^V (

The sum of the first two terms D( Dt = d( dt + (V «V )cc is the well-known "Lagrangian" convective derivative. The third term ( (V«V ) is the product of the vorticity and the volumetric expansion rate. The

'By inserting the convective derivative from Eq. (4) into Eq. (7), one can also arrive at the alternative vorticity equation

p dC + (v «V )(( + cC (V«V )- (c «V )V = Vpx fV^ 1 + pV x f O

I p J I p J

divergence V«v = (q'S'ensiUe/p0CpT0) is proportional to the sensible heat release rate (assuming constant specific heat) [17]. The release of heat causes the vorticity to spread out over a larger volume; but in itself does not change the overall amount of vorticity. The fourth term, (co«V )v describes how velocity gradients, when aligned with the vorticity, stretches it and

concentrates it into a smaller volume; while velocity gradients orthogonal to the vorticity cause it to tilt and be reoriented. This latter term is responsible for the unique behavior of vorticity including how it facilitates dispersion of energy to both larger and smaller scales. The next section shows why the terms on the left hand side, above, do not alter the total vorticity carried by the flow.

2.2. Kelvin's theorem

Kelvin's Theorem describes the evolution of vorticity as it travels with the fluid. Fig. 3b, above, shows a thin vortex tube surrounding the vortex line 3 . Consider the circulation, r (t) around a closed path, C (t) surrounding this vortex tube.

C (t) is attached to fluid elements as they move with the fluid and A (t) is the area enclosed by C (t) . Then by Stokes

Theorem, the circulation r (t) is also equal to the vorticity going through the area A (t)

r (t) = (j)c v ■ dl = $ 3-ndA (8)

The rate of change of r (t) = (jj^ v ■ d 1 is

d r (t)/dt = (jj^ {))(Dv/Dt )d 1 + <jjc () v •(( 1/Dt) = £ ()V x(Dv/Dt )ndA + j^v^Vj) (9)

Stokes Theorem has been applied to the first integral on the right, while the second integral is zero, because the integrand v-(DljDt) = v.(d 1-V)v = vJd£i (dvJIdxt) = ± V|v|2 •d £ ,when evaluated around the closed path, C (t) , returns to its initial value. Upon recalling the derivation from Eqs. (4) through (7), the right hand side of Eq. (9) then becomes

drdT = Jaw {d3+3 (Vv> («^[n^ (t) = JaW {vpx [g - D]+»v23 ^•ndA (t) (10)

This says that there is no change in circulation for an inviscid fluid having constant density, because the right hand side of the above Eq. (10) would be zero. In such cases, the total vorticity going through area A (t) remains invariant as it moves in

time. Viscous diffusion in a flow field does not create vorticity. It only dissipates it and/or moves it around. For purely buoyant turbulent diffusion flames, without an external supply of vorticity, all the vorticity must come from this baroclinic term. In contrast, in the absence of gravity, the acceleration Vpx Dv/Dt would be small, so that there is little, if any, baroclinic production of vorticity and most of the vorticity present in the flow must come from external sources, such as, a jet orifice supplying fuel to a momentum dominated turbulent diffusion flame. The supplied vorticity eventually dissipates by viscosity, after being concentrated by the turbulent flow. In contrast, a buoyant turbulent diffusion flame receives virtually all its vorticity from the flames.

2.3. Cascades of vorticity and turbulent kinetic energy

Here we study the cascade mechanism connecting the fluid mechanics of fires taking place at large scale with the molecular diffusion processes occurring at small scale where the dissipation of vorticity and turbulent kinetic energy occurs and determines: (1) soot formation and oxidation, (2) flame radiation, as well as (3) incomplete combustion in turbulent diffusion flames. A cascade mechanism is an integral part of Large Eddy Simulation (LES) models that explicitly describe large-scale processes and boundary conditions while modeling the smaller, subgrid-scale (SGS) processes. The premise is that the problem as a whole is governed at large scale by its boundary conditions, whereas the processes taking place in the smaller subgrid eddies are less sensitive to the specific problem details and thus can be modeled. Such cascades are a natural outcome of the fact that fluid eddies interact most strongly with eddies of similar size or nearly similar sizes while interacting only weakly with eddies having distinctly different sizes. Tennekes and Lumley [18] as well as Ertesvag and

Magnussen [19] conveniently modeled the cascade process by considering a stepwise transfer of kinetic energy between eddies close in size. Viscous effects become important only at small scale where the kinetic energy is consumed.

LES modeling techniques are well developed for momentum-dominated flows for which the turbulent kinetic energy originates at large scale and dissipates relatively unchanged as it flows down a cascade of eddy sizes until it is consumed at Kolmogorov scale by molecular processes. The fact that it is unchanged in the inertial region makes the cascade particularly simple for momentum-dominated flows. The cascade process is quite different for buoyancy dominated turbulent diffusion flames. The next three sections describe how the presence of buoyancy adds both vorticity and turbulent kinetic energy at each cascade stage that is larger than the micro-scale where dissipation dominates.

2.4. Vorticity cascade

The cascade process is distinctively different for buoyancy dominated turbulent diffusion flames for which the baroclinic torque generates additional vorticity at each successive cascade step. Consider a particular scale A with wave number ka = 1/A . For some readers it is helpful [but not necessary] to think of the sizes of successive steps decreasing by

factors of two with associated cascade interval covering the wave-number range kJV2 < ka < -\/2ka . Equation (7) says the baroclinic production of vorticity by the flames within the ka interval is of the order

PA = Vp x [g - Dv/Dt] - ApfgKA (11)

Here the gradient Vp for scale A is approximated by ApfJA = (p„- pf ))a . We employ the full difference pm - pf because buoyant diffusion flames generally have large intermittency with clumps of ambient air distributed throughout the flaming volume. For simplicity, the vector [g - Dv/dt] is replaced by g here for clarity in understanding the physics. The approximation should not be necessary for an LES formulation, for which the acceleration Dv/Dt is available from within the resolved-scale equations. Notice that the baroclinic production increases as ka increases [i.e. as the scale A decreases]. Equation (7) also says that the viscous dissipation for the interval ka is of the order

D(A = pvV(a ~ pv(ak_2 (12)

Figure 4a, below, and Eqs. (11) and (12) show that the dissipation of vorticity, D increases more strongly than the

production P with increasing ka . The last two terms on the left of Eq. (7), namely ( (V«v) and ((«V)v , describe,

respectively, the dispersion (i.e. expansion) of vorticity toward larger and the dispersion (i.e. concentration) to smaller cascade intervals. Both terms, are of the order

T(A ~ p(a I^aI ~ p(A (13)

since the characteristic strain rate, = \dv/d*|A is of the order of mA . In the absence of viscous effects, the dispersion, T , both to and from interval ka should be proportional to the production, P . Most of what is produced spills out to the neighbouring intervals. That is, from Eqs. (11) and (13)

T(A ~ p(A ~ ApfgKA (14)

This suggests that (A ~ JApfgKA /p which when substituted into Eq. (12) yields the viscous dissipation

D„a ~ pvjApfgIp K52 (15)

One sees that the viscous dissipation increases very strongly with ka , indeed to the 5/2 power. The viscous dissipation is thus unimportant at all scales except near or beyond the final stage where diffusion processes can consume all the available vorticity. See Fig. 4a. This allows one to identify the final scale Xf = 1/Kf . Production exceeds dissipation for A > Xf

while dissipation exceeds production for À < Xf . Solving for Xf by equating the vorticity production, Eq.(11), and dissipation, Eq. (15), one has

ÀPfghf/ pv2 = const.

We note here that the viscous, thermal and species diffusivities are all of similar magnitude and all involved in the combustion process. Substituting the thermal diffusivity, a , for one of the kinematic viscosities, v , the result says that the combustion of buoyant turbulent diffusion flames is controlled by a Rayleigh-Taylor instability mechanism at the micro-

Ra = ÀpfgXj I pv j a j = const.

The Rayleigh number depends only on local properties and thus is independent of the overall size of the fire. In particular, the characteristic flamelet size for the Rayleigh-Taylor instability mechanism

h ~ \?fafl(Af/P)]"3

is independent of the overall scale of the fire. The corresponding characteristic flow time

Tf ~ hfl iff ~ ff afl (APfg/p)2

is similarly independent of the size of the fire. This explains why the mean volumetric radiant emission, qR is independent on the size of the fire [10]. It also explains why the volumetric heat release rate, qfh is independent of fire size,

<ich ~ PhTf ~ Ph {ÀPfg/p)) yfy

h = Xch ÀHC/ (1 + S )

Fig. 4. (a) Vorticity cascade; (b) rotational and translational turbulent kinetic energy cascades.

is the chemical energy released with the production of unit mass of products. Finally, recalling that both vf and af vary inversely with pressure Eq. (20) implies that the volumetric rate of heat release scales with gravity and pressure according to qr"h ~ g2 3 p4 3 . Appendix A shows that is a requirement of the Navier Stokes equations for both laminar and turbulent buoyant diffusion flames.

Figure 4a shows that both the vorticity, œÀ , and production of vorticity, P , peak at the Rayleigh wave number,

ka = Kf (i.e. A = Xf ), at which molecular diffusion first becomes dominant. This peak yields a single Rayleigh number

characterizing the most unstable wavelength. For very small fires, would be the approximate size of the initial flame

cells. For very large fires, one expects to observe a predominance of cells of this size embedded in larger scale flow. This peaking of vorticity explains why a Rayleigh-Taylor instability mechanism operating at the turbulent micro-scale governs the combustion of purely buoyant turbulent diffusion flames independent of the overall fire size. This is in contrast to purely momentum dominated diffusion flames, which (by virtue of Kelvin's Theorem) receive their vorticity from external and/or large-scale sources external to the cascade. Practical applications may include combinations of these buoyant and momentum dominated mechanisms.

2.5. Turbulent kinetic energy cascade

To gain a deeper understanding and assist the development of LES models we examine cascades of both rotational and translational kinetic energies. The two cascades turn out to be quite similar, but differ in the manner buoyancy generates the turbulent kinetic energy. Consider a clump of fluid having size A . It can possess both rotational and translational kinetic energy. Horizontal density gradients acting on buoyancy produces a torque that can generate rotational kinetic energy; whereas the vertical density gradients can create the forces that generate translational turbulent kinetic energy.

For the rotational case, consider a cylindrical volume having radius R = A in the flaming zone. The baroclinic torque per unit volume, T exerted by the flames is T - Vpx[g - Dv/Dt]R2 ~ ApfgR . Upon multiplying this by the angular velocity, a it becomes the power per unit volume PBR = Ta ~ ApfgaR . However aR is the circumferential velocity which in turn is proportional to the square root of the kinetic energy, aR ~ kA2. Thus the baroclinic production of rotational kinetic energy per unit volume is

Pbr ~ Apfgk^2 (22)

For the translational case, when lower density eddies rise up and exchange places with higher density eddies above, they lose potential energy and release turbulent kinetic energy. Conversely, when higher density eddies, lying below, rise up and exchange places with lower density above, they require additional potential energy and consume turbulent kinetic energy. Consider first the usual situation of a hot buoyant eddy rising in the flaming zone where the surrounding density increases with height. The upward buoyant force exerted on the clump of size A per unit volume is of the order A (d p/dz)g ~ Apfg .

By multiplying the force by vertical velocity, it becomes the rate of energy release. For an isotropic turbulent flow having turbulent kinetic energy, kA the production of turbulent kinetic energy in a flaming zone is

Pbt ~ Apfgk'l2 (23)

Thus PBR and PBT are of the same order. The viscous dissipation of turbulent kinetic energy is

Dv ~ pVjkKl (24)

The viscous dissipation, being proportional to the square of ka becomes important only for large ka . The turbulent dissipation (i.e. transfer) of rotational kinetic energy to smaller and larger cascade intervals is of the order

Dt ~ pk'^KA (25)

Equating the turbulent dissipation, Eq. (25), with the baroclinic production, Eqs. (22) and (23) one has

kA ~ (pfg/p) (26)

Substituting this into Eq. (22) the baroclinic production becomes

PBR PBT ~ APfg (APfg/p) KA1/2 (27)

Similarly substituting this into Eq. (24) the viscous dissipation becomes

Dv ~ Vf ApfgK (28)

The production and dissipation intersect at the Rayleigh scale Xf = 1/ka = A where the kinetic energy is consumed. See Fig. 4b. Solving for hf by equating the last two equations one obtains the result

A Pfghf I pv2 = const (29)

This is the same as Eq. (16) and confirms the Rayleigh-Taylor instability result.

2.6. Non-flaming buoyant generation of turbulent kinetic energy

In volumes where there is no combustion, that is, where there is either no fuel or no oxidant, one cannot assume a maximum density difference equal to Apf = p0 - pf . Instead the sum of the baroclinic and translational generation of turbulent kinetic energy in for an eddy of size A is

Pb = Pbr + Pbt ~ k1 [|Vp X g| - Vp.g] Ka-1 (30)

It is noted that when dpdz is positive, that is in the opposite direction to g kinetic energy is released as eddies of lesser density rise up and exchange places with eddies of greater density above. Conversely, when dp dz is negative and eddies of greater density rise up and must exchange places with eddies of lesser density above kinetic energy is consumed. Equation (30) shows that both the vertical and lateral density gradients contribute to the production of turbulent kinetic energy. The corresponding viscous dissipation is

Dv ~ pvkK (31)

while the turbulent dissipation is

Dt ~ pkfK (32)

Equating Eqs. (30) and (32) the turbulent kinetic energy is kA ~ (1/p)[|Vp x g| - Vp«g] kA2 . Substituting this into Eq. (31) one obtains DV ~v [Vp«g +1 Vp x g| ] which is independent of A . Finally, substituting kA into Eq. (30) and identifying the eddy size, A = hd, at which the buoyant production equals the viscous dissipation, one obtains the modified Grashof number

[|Vpxg-Vp.g]h4/pv2 = 1 Dt ~ pkfK, (33)

This result should be applicable quite broadly to LES solutions for buoyant turbulent free convection problems.

2.7. Eddy dissipation concept

The eddy dissipation concept (EDC) of Magnussen et al. [20] is well suited for the lazy turbulent diffusion flames characteristic of fires. It provides a simple expression for mixing controlled combustion of turbulent diffusion flames assuming infinitely fast chemistry. The reaction proceeds as soon as the fuel and oxidant are mixed. The reaction rate is proportional to whichever reactant is more limited and the reaction stops whenever the more limited reactant is depleted. The volumetric production rate of stoichiometric products is equal to (1 + S) times the volumetric consumption rate of fuel and is equal to (1 + S)S times the volumetric consumption rate of oxidant. That is

m P = - (1 + S )m F = -[(1 + S )/S ] m 0 (34)

The EDC concept expresses the volumetric production of fuel as

mF(t) = pdYF/dt = -(p/TMlx)Min[7^ (t),Yo (t)/S] (35)

The corresponding volumetric heat release rate is

q'l = X* \_AHC/(1 + S)]mP(t) = (ph(1 + S)/tMx)Min[Y„ (t),Yo (t)/S] (36)

where from Eq. (21) h = xch AHCj(1 + S) is the chemical energy released with the production of unit mass of products.

Comparing Eqs. (20) and (36) one sees that the mixing time TMix proportional to the characteristic flow time,

- __2-|1/3

TJ = ~vfOf (Apjg/p) for buoyant turbulent diffusion flames. This should be true for CFD computations having

mesh sizes, A , both much greater and much less than the characteristic flamelet size, Xj . Clearly, for A » Xj the average heat release rates should be similar. The argument becomes more subtle in the limit, A « Xj , corresponding to Direct Numerical Simulation (DNS) calculations. In this limit the area of the diffusion flame per unit volume is

aJ ~ 1/Xf (37)

The production of products per unit flame area, m"P , is proportional to the density times the square root of the product of diffusivity times the characteristic strain rate, 1jtj , of the flames. That is,

mP ~ p(1 + S)jDj7~f (38)

Upon multiplying by the flame area per unit volume the corresponding volumetric heat release rate is

qh = hm P ~ hm Pa J ~ ph (1 + S)D/ (X2fzf ) = p h/ Tf . (39)

Thus for both A » Xj and A « Xj the volumetric heat release rate is proportional to 1jt f . Experience, suggests the results are not sensitive to the particular value of 1/TMix .

2.8. Large eddy simulation (LES)

Large Eddy Simulation (LES) shows promise for modeling buoyant turbulent diffusion flames. LES solutions explicitly solve the fully resolved equations and boundary conditions for the problem at large scale, while modeling the smaller unresolved subgrid (SGS) processes. The premise is that the problem as a whole is governed at large scale by its boundary conditions, whereas the processes taking place in the smaller subgrid eddies are less sensitive to the specific problem details and thus can be modeled. The subgrid equations describe how the turbulent kinetic energy cascades down to smaller and smaller scales and ultimately dissipates by viscosity. Cascades are a natural outcome of the fact that fluid eddies interact most strongly with eddies of similar or nearly similar sizes and interact only weakly with eddies having distinctly different sizes.

The physics in current LES fire models [21-24], is suitable for momentum dominated turbulent diffusion flames, which receive all their kinetic energy (and vorticity) from sources external to the flames. However, the equations for the subgrid need to include the buoyancy physics important for buoyant turbulent diffusion flames. Omitting subgrid buoyancy physics will make the solution change with mesh size, since the amount of included physics changes. By including the buoyancy physics in the subgrid the solutions should become more mesh-size independent and allow for larger mesh-sizes without significant loss of accuracy.

2.9. Suggested LES equations

Conservation equations for typical LES fire models solve the Navier-Stokes equations using fully implicit solution schemes. The partial derivatives in the conservation equations are discretized using the finite volume (FV) method. The LES algorithm uses the eddy viscosity model for the subgrid scale (SGS) stress tensor, and the gradient diffusion model with unity turbulent Lewis number to model the SGS scalar and enthalpy flux. The combustion model employs the eddy dissipation concept (EDC) with conservation of energy providing the sensible enthalpy.

The Favre filtered conservation equations for mass, momentum, state, species and enthalpy are:

dp + d_pû± = o dt dx,.

dpû. dpùûj = d p d

dx dx,.

p(v +vt )

N + dU - 2 dû* 5..

dx. v j

dx 3 dx.

Po = pRT

P = P - Po

dpY + dpû Y = _d_

dt dx, dx,

v, ïdY p| a +— I—'-1 Pr I dx.

for ' = F, O

dphs dpùjhs d

p| a +

Pr, I dx

-in . ■ m

-9R + Ich

here hs is the sensible enthalpy and vt = ck À k^2 is the turbulent viscosity. The above conservation equations are all conventional.

The suggested modification of the one-equation SGS kinetic energy model, such as employed in OpenFOAM [22], is

dpkÀ + dpûjkÀ__d_

dt dx, dx,

-( \ dkÀ p(v +vt ^^

= PSÀ + PBÀ- psà

dùi dùj 2 dùk

—'- + —J----- 5'

dx. dx 3 dx, ■ v J ' k

- p^Sj

dû' dx.

CBlÀpsgK2

flaming

Pbà = Pbr + Pbt ^ cb2ÀkÀ/2 [| Vp x g - Vp-g] non-flaming

s. = C„À-1k3'

The left hand side of Eq. (46) serves to smooth the turbulent kinetic energy kA among neighboring cells. Equation (47) describes the turbulent dissipation, PS A being transferred by turbulent shear stresses from the resolved scale to the subgrid; while Eq. (48) describes buoyant rotational [Eq. (22)] and translational [Eq. (23)] production of turbulent kinetic energy in the subgrid for both flaming and non-flaming volumes. Finally, Eq. (49) expresses the dissipation of turbulent kinetic

energy, sA at the particular mesh size, A . Both the buoyant production, PBA and the production by the resolved equations, PSA should decrease with A1/2 as the mesh size A approaches the characteristic Rayleigh scale, Xf .

The important difference between momentum-dominated and buoyancy dominated flames occurs when one evaluates the inverse of the micro-scale flow time, af = 1/zf ,which is proportional to the volumetric heat release rate and determines the

radiant fraction. The inverse flow time from Eq. (19) expressed in terms of the kinetic energy from Eq. (26) for eddies of size A is

(Apfg/p)2 1/3 kA 1/3 i 1/2 A kA A 4/3 -- JVf af

Jvf af af a2 _ A2

which is independent of the mesh size A since kA is proportional to A . Conventional treatments based on the Kolmogorov

i \1/2

cascade for momentum-dominated situations give the inverse micro-scale flow time as mf ~ (s/vf ) or

1 s 1/2 " k3/2 " 1/2 r kA/2A~|

Tf vf _ _ vi _

Notice that Eqs. (50) and (51) have different exponents on their Reynolds numbers, ReA = / . Otherwise, the equations are the same. Now what happens as one changes the mesh size when calculating buoyant flamers? From Eq. (26) we see that the turbulent kinetic energy, kA, is proportional to the mesh size, A ; so that the inverse flow time given by Eq. (50) for buoyant situations is independent of mesh size, while the inverse flow time using the Kolmogorov formula, (51) is proportional to A1/4 and thus depends on the mesh size. The latter difference is not large, but may prevent establishing mesh size independence.

The final suggestion is to use Eq. (50) when for evaluating the mixing rate1/zMx for the eddy dissipation concept described in Section 2.7. More specifically,

1 1 {APfg/P)2 1/3 kA

TMix Tf ylVf af _V^^ ^^A 2 _

This choice is exact for purely buoyant flames and reasonably accurate for flames with added supply of external vorticity. 3. Conclusions and suggestions for further research

The buoyant generation of vorticity and its dissipation within the flamelets drives the combustion of fuel and oxidant in buoyant turbulent diffusion flames, characteristic of fires. The generation of vorticity peaks the Rayleigh micro-scale near where the vorticity generation comes into balance with its dissipation. It is an example of a Rayleigh-Taylor Instability mechanism whose disturbances generally originate at small-scale, grow and eventually encompass the entire flow. The mechanism is characterized by a nearly constant Rayleigh number independent of fuel type, which explains why experiments show that both radiant fractions and volumetric heat release rates of buoyant turbulent diffusion flames are independent of the overall size of the fire. The peak characterizes the most unstable wavelength. This peak yields a single Rayleigh number characterizing the most unstable wavelength, Xf . For very small fires, the mechanism suggests that this

most unstable wavelength should characterize the initial size of the flame cells at the laminar to wiggly laminar transition. For very large fires, the mechanism suggests one should observe a predominance of flame cells of this most unstable size embedded in the larger scale flow field.

Kelvin's theorem explains the contrasting roles vorticity plays in momentum-dominated and buoyancy-dominated turbulent diffusion flames. Vorticity for momentum-dominated flames comes principally from external sources, whereas for buoyancy-dominated flames it comes from buoyancy driven baroclinic torque. Flame-generated turbulence becomes of increasing importance as the size of the fire increases, or to be more specific, as the ratio of overall flame height, IR to flamelet size, Xf increases.

Cascades of vorticity as well as turbulent kinetic energy provide further understanding of buoyant flames. Cascades are a natural outcome of the fact that fluid eddies interact most strongly with eddies of similar size or nearly similar sizes while interacting only weakly with eddies having distinctly different sizes. Cascade mechanisms form the basis for the subgrid-scale (SGS) models used in Large Eddy Simulations (LES).The paper suggests a modified set of LES equations that promise to provide mesh-size independent solutions. The LES equations are applicable to both buoyant and momentum controlled diffusion flames as well as when both momentum and buoyant effects are present. They also apply to situations with and without combustion. The suggested LES equations require testing to establish their mesh-size independence.

The proposed concept that buoyant turbulent diffusion flames as an example of a Rayleigh-Taylor instability mechanism requires much further study and confirmation. In general, initial disturbances of Rayleigh-Taylor Instability mechanisms originate at some small scale, grow, and eventually encompass the entire flow. How might this apply to hazardous fire situations? Other areas of possible future research include performing spectral analyses of measurements of buoyant turbulent diffusion flames to identify their characteristic time, t f and length, hf scales. Do t f and hf have their predicted

pressure and gravity dependencies? One might similarly perform spectral analyses of DNS or finely resolved LES solutions of buoyant turbulent diffusion flames to examine the cascades of vorticity and turbulent kinetic verses eddy size. For mathematically inclined researchers it would be interesting to perform a classical Rayleigh stability analysis by analytically solving the perturbation equations for the initial size and growth rate of buoyant diffusion flame. A possible configuration could be that of cellular and turbulent ceiling fires [25, 26].

In addition understanding of fire protection physics require systematic efforts in many other areas, flame radiation and heat transfer, incomplete combustion and emission of toxic products, pyrolysis of fuels, spray-plume interaction, water transport on burning objects, flame quenching, etc.

Acknowledgements

The author expresses gratitude for the ideas and encouragement from his colleagues, Dr.'s Sergey Dorofeev, Yi Wang and Prateep Chatterjee of FM Global and also to Prof. Liu Nai'an of the State Key Laboratory of Fire Science, whose personal care has made this research possible. The author is also grateful for the financial support of FM Global and the Chinese Academy of Sciences.

Appendix A. Some fundamental scaling requirements

The earlier study on the "Pressure Modeling of Fires" [27] showed that the Navier-Stokes and Shvab-Zeldovtch equations for buoyant diffusion flames equations require that the volumetric heat release increase with the 4/3 power of ambient pressure. The requirement is quite general. It applies to both laminar and turbulent buoyant diffusion flames. The study also extended modeling principles to flame spread. For the concept to be valid one has to assume either that the flame radiation is negligible or that the radiant fraction of the heat release rate is constant. The latter assumption is at least approximately correct. Experiments confirmed the prediction by varying the ambient pressure over thirty atmospheres. Virtually the same analysis technique shows that the volumetric heat release rates increase with the 2/3 power of the gravitational acceleration. The scaling requirement provides independent support for the Rayleigh-Taylor Instability mechanism that predicts the volumetric heat release should the same p4/3g 2/3 proportionality. The older entrainment based flame height predictions do not satisfy this fundamental scaling requirement for the volumetric heat release rates of buoyant flames.

The analysis technique is easy to apply. Take for example, Eqs. (40) through (50) that provide the complete set of governing LES equations for buoyant turbulent diffusion flames. One increases the ambient pressure by the factor "a" and increases the acceleration of gravity by the factor "A" while scaling all other quantities according to Eq. (53) below. After performing the scaling operations, all the governing equations should increase by the same factor, but otherwise continue to be satisfied, but with an altered set of time and space scales. This analysis presumes the boundary conditions scale in a similar manner with the same set of altered time and space scales. If this is the case then the two problems are an exact replica of each other except for their time and space scales. The underlying physics is the same. Mathematically we have shown that the governing equations are invariant under the transformation,

Po. ^ aPo

g. ^ Ag p. ^ ap,

t. ^ a"1'3A-2l3t x„ ^ a-2l3A-mxl u, ^ a-1'3All3ul m,^ a4,3A2,3< C, ^ a^A^ql

Ra, ^ Ra Re, ^ Re Gr, ^ Gr Pr, ^ Pr T,,Y ,,h, ^T, Y,h M, ^ M Vf, ^ a~Vf

Ps, ^ a2,3A4,3PS Pg» ^ a A e, ^ a-1,3A4,3s p, ^ a1'3A2'3p Xf, ^ a-2,^A-I,3A/ Tf, ^ a"1'3A-2,3rf

1/3 A_2/3-

a ^ a-2l3A~113A

References

[1] Morton, B. R., Taylor, G. I., Turner J. S., 1956. Turbulent Gravitational Convection from Maintained and Instantaneous Sources, Proceedings of the

Royal Society A 234, p. 1.

[2] Thomas, P.H., 1963. The Size of Flames of Natural Fires, Proceedings of the Combustion Institute 9, p. 844.

[3] Steward, F. R., 1964. Linear Flame Heights for Various Fuels, Combustion and Flame 8, p 171.

[4] Steward, F. R., 1970. Prediction of the Height of Turbulent Diffusion Buoyant Flames, Combustion Science and Technology 2, p.203.

[5] McCaffrey, B. J.: Purely Buoyant Diffusion Flames: Some Experimental Results, Technical. Report NBSIR79-1910, 1979.

[6] Heskestad, G. 1983. Luminous Flame Heights of Turbulent Diffusion Flames, Fire Safety Journal 5, p. 103.

[7] Delichatsios, M. A., Orloff, L., 1985. Entrainment Measurements in Turbulent Buoyant Jet Flames, Proceedings of the Combustion Institute 20, p. 367.

[8] Zukoski, E. E., Cetagen, B. M., Kubota, T., 1985. Visible Structure of Buoyant Diffusion Flames, Proceedings of the Combustion Institute 20, p. 361.

[9] Markstein, G. H., 1985. Relationship between Smoke Point and Radiant Emission from Buoyant Turbulent and Laminar Diffusion Flames, Proceedings of the Combustion Institute 20, p. 1055.

[10] Markstein, G. H., de Ris, J. L., 1985. Radiant Emission and Absorption by Laminar Ethylene and Propylene Diffusion Flames, Proceedings of the Combustion Institute 20, p. 1637.

[11] Delichatsios, M. A., Orloff, L., 1985. Effects of Turbulence on Flame Radiation from Diffusion Flames, Proceedings of the Combustion Institute 20, p. 1271.

[12] Orloff, L., de Ris, J. L., 1982. Froude Modeling of Pool Fires, Proceedings of the Combustion Institute 19, p. 885.

[13] Markstein, G. H., de Ris, J.L., 1991. Wall-fire Radiant Emission. Part 1: Slot-burner Flames, Comparison with Jet Flames, Proceedings of the Combustion Institute 23, p. 1685.

[14] Hahn, C. J., Shapiro, S. S., 1968. Statistical Models in Engineering, John Wiley & Sons. Ltd, New York.

[15] Cetegen, B. M., Ahmed, T. A., 1993. Experiments on the Periodic Instability of Buoyant Plumes and Pool Fires, Combustion and Flame 93, p. 157

[16] Tieszen, S. R., Nicolette, V. F., Gritzo, L. A., Moya, J. L., Holen, J. K., Murray, D. Vortical Structures in Pool Fires: Observation, Speculation, and Simulation, Tech. Rept. SAND—96-2607, Sandia National Labs, 1996.

[17] Baum, H. R., McCaffrey, B. J., 1989. "Fire Induced Flow Field—Theory and Experiment," Fire Safety Science - Proceedings of the Second International Symposium, International Association of Fire Safety Science, pp. 519.

[18] Tennekes, H., Lumley, J. L., 1972. A First Course in Turbulence, M. I. T. Press, Cambridge, MA.

[19] Ertesvag, I. S., Magnussen, B. F., 2000. The Eddy Dissipation Turbulence Energy Cascade Model, Combustion Science and Technology 159, p. 213.

[20] Magnussen, B. F., Hjertager, B. H., 1977. On Mathematical Modeling of Turbulent Combustion with Special Emphasis on Soot Formation and Combustion, Proceedings of the Combustion Institute 16, p. 719.

[21] Wang, Y., Chatterjee, P., de Ris, J. L., 2011. Large Eddy Simulation of Fire Plumes, Proceedings of the Combustion Institute 33, p. 2473.

[22] OpenFOAM User Documentation, http://www.opencfd.co.uk/openfoam/doc/user.html.

[23] McGrattan, K., Hostikka, S., Floyd, J., Baum, H. R., Mell, W., McDermott, R., Fire Dynamics Simulator (Version 5), Technical Reference Guide, Volume 1: Mathematical Model. NIST Special Publication 1018-5, NIST, Gaithersburg, MD, Sept. 2010.

[24] Chatterjee, P., de Ris, J. L., Wang, Y., Dorofeev, S. B., 2011, A Model for Soot Radiation in Buoyant Diffusion Flames, Proceedings of the Combustion Institute 33, p. 2665.

[25] Orloff, L., de Ris, J., 1972. Cellular and Turbulent Ceiling Fires, Combustion and Flame 13, p. 389.

[26] Orloff, L., de Ris, J., 1971. Modeling of Ceiling Fires, Proceedings of the Combustion Institute 13, p. 979.

[27] de Ris, J., Kanury, A. M., Yuen, M. C., 1973, Pressure Modelling of Fires, Proceedings of the Combustion Institute 14, p. 1033.