Fractal Geometry of a Numerically Simulated Flame Surface

Fractal Geometry of a Numerically Simulated Flame Surface. Lance R. Collins. Ind. Eng. Chem. Res. , 1995, 34 (8), pp 2588–2597. DOI: 10.1021/ie00047...
0 downloads 0 Views 2MB Size
Ind. Eng. Chem. Res. 1995,34, 2588-2597

2588

Fractal Geometry of a Numerically Simulated Flame Surface Lance R. Collins Department of Chemical Engineering, The Pennsylvania State University, 118-B Fenske Laboratory, University Park, Pennsylvania 16802

Two-dimensional, direct numerical simulations of premixed turbulent combustion in the socalled flamelet regime were done in an earlier study. The flame “surfaces” passed through stationary isotropic turbulence that was maintained by forcing at either large scales (LS)or small scales (SS) producing inertial ranges with power-law exponents of -3 and -5/3, respectively. In the present study, the impact of the inertial range scaling on the fractal dimension of the flame surface is considered. It is shown that the SS simulations, having the more classical -513 inertial range, yield fractal dimensions that are in closer agreement with experimental measurements in the literature than the LS simulations with the -3 inertial range. Indeed, the latter produce consistently lower fractal dimensions than the SS simulations and experiment. A theoretical expression for the fractal dimension of a n isocontour of a passive scalar developed by Vascillicos and Hunt does not compare favorably with the simulation results for the flame surface. We propose a n alternative empirical expression (DK= d - 0.42 - 0.15q, where d is the dimension of the system and q is the inertial range exponent) that agrees more closely with the simulation results, particularly at high turbulence intensities.

I. Introduction

IOJ

This paper is concerned with the coupling between turbulent fluctuations and the mean rate of reaction in premixed flames at relatively large Reynolds numbers. It is well-known that the rate of reactant consumption can be increased severalfold due to fluid mechanical enhancements resulting from turbulence in the reactant mixture. Naturally an effect of this magnitude has been of practical interest t o designers of combustion systems and of fundamental interest t o combustion scientists who have studied the mechanisms of the interaction over the past several decades. Perhaps no single paper has had more influence on the modern conceptual picture of premixed flames than the pioneering study by Damkohler (19401,who identified a combustion regime referred to as flamelet combustion. This regime is most easily described by considering the diagram shown in Figure 1. The axes of the diagram represent the two parameters that control turbulent combustion, namely the dimensionless turbulence intensity U’/ULand the dimensionless turbulence length scale Z/&. The solid lines in Figure 1 represent further subdivisions of the parameter space based on unit values of the Reynolds number (U’Zh), Damkohler number ((Z/U’)/(ddu~)),and Karlovitz number (6;/v2). The intrinsic rate of reaction in flames is sufficiently high that most practical devices operate in the flamelet regime (shown with a cross hatch in the diagram). Flamelet combustion is characterized by very thin reaction fronts (relative to the scales of the turbulence) that are bent, distorted, and wrinkled by turbulence. Remarkably, the intrinsic reactions within the thin fronts remain largely unaffected by the outer flow. Damkohler recognized that, in this limit, the dominant effect of the turbulence is t o increase the surface area of the flame without significantly altering the chemistry within the thin front. On the basis of simple scaling arguments he postulated that the turbulent flame speed would be given by the following expression:

IO’

10’

IO0

lo-’ 10.‘

IO0

IO’

IO‘

IOJ

r/ 4 Figure 1. Parameter space for premixed combustion as suggested by Williams (1988).Cross-hatched region is the “flamelet” regime considered in this study.

where ut is the turbulent flame speed, U L is the laminar flame speed for an unwrinkled surface, and As and A, refer to the surface area and cross-sectional area of the flame, respectively. The conceptual picture that Damkohler proposed for flamelet combustion remains a central feature of most modern theories of turbulent combustion 50 years later. One reason for its longevity is that it implies a decoupling of the very complex chemistry within the flamelet from the surrounding turbulent flow field. Indeed eq 1 is based on the postulate that calculating the turbulent flame speed is equivalent t o estimating the surface area of the flame per unit cross section. In spite of the incredible simplification implied by eq 1, the issue of predicting turbulent flame speeds over a wide range of

Q888-5885I95/2634-2588$09.Q~lQ0 1995 American Chemical Society

Ind. Eng. Chem. Res., Vol. 34, No. 8, 1995 2689 conditions remains fundamentally unsolved. Indeed, there is a complication to this simple formula for ut. A flame surface exposed to turbulent fluctuations must wrinkle a t a variety of different length scales. That is, large-scale turbulent fluctuations will distort the surface at large scales while simultaneously small-scale eddies will create finer wrinkles. The surface area generated by wrinkles of all scales will contribute to the overall rate of reaction; thus a predictive theory of flamelet combustion must account for wrinkles over a multitude of length scales. There appear to be two fundamentally different approaches to accounting for the multiscaled nature of the surface wrinkling process. The first uses a concept called fractals, introduced by Mandelbrot (1975, 1983) to explain a variety of non-Euclidean surfaces that arise in natural phenomena. A Euclidean surface (embedded in a three-dimensional (3D) space) by definition has a topological dimension of 2; however, if the surface is a fractal, it may fill a greater fraction of the space and thereby reside in a dimension that is greater than 2 (e.g., isosurfaces of scalar concentration in a turbulent fluid have a fractal dimension of 2%). Fractal surfaces have been observed in a variety of natural phenomena including clouds, coastlines, diffusion-limited aggregates, and, most importantly to this paper, flame surfaces (Mantzaras et al., 1989; North and Santavicca, 1990). Models for predicting turbulent flame speeds based on a fractal characterization of the surface have been developed and have had some successes (Gouldin, 1987; Kerstein, 1988; Peters, 1988). However, a major difficulty with the approach is determining the physical cutoffs (inner and outer) that define the scales over which the surface is fractal. Experimentally the inner and outer cutoffs are not found to be sharp demarcations as most theories imply, but rather smooth transitions (North and Santavicca, 1990). Nevertheless, fractal flamelet approaches have provided important physical insights into the turbulence-flame interaction, as will be shown shortly. The alternative approach is based on defining an appropriate Fourier spectrum for the flame surface. This approach has been described in detail elsewhere (Collins, 1995; Dandekar and Collins, 1995; Peters, 1992; Bray and Libby, 1986; Bray et al., 1984). Although this paper is primarily concerned with fractal descriptions of flames, it will be demonstrated that the fractal dimension is closely related t o the power-law exponent of the inertial range of the scalar spectrum, suggesting that the two approaches may be connected. The present analysis considers results from direct numerical simulations of flame surfaces passing through preexisting turbulence (Collins, 1995). The simulations are based on the flamelet concept described above. The dynamics of the flame surfaces are determined by solving a scalar equation for a quantity G(x,t),isocontours of which constitute flame surfaces. The transport equation for G(x,t), often referred to as the Huygen equation, accounts for convection due to the local fluid velocity and propagation resulting from reaction. The simulations resolve all features of the surface for length scales ranging over 2 orders of magnitude. The flame speeds that result from the simulations are in agreement with experimental measurements in the literature (Collins, 1995). The trends in the simulations indicate a sensitivity of the flame speed to the shape of the energy spectrum. The present study will explore the

relationship between the fractal dimension of the surface and the turbulent energy spectrum. Section I1 describes the simulations in detail and summarizes the results from the earlier study. Section I11 shows the technique used t o measure the fractal dimension and how it was implemented to analyze the present results. Results of eight different simulations are presented in section IV along with comparisons with experimental measurements in the literature. Section V analyzes the relationship between the fractal dimension and the scalar spectrum, followed by a summary of the results in section VI. 11. Direct Numerical Simulation of Flame

Surfaces The term direct numerical simulation (hereafter DNS) refers to a broad category of numerical solutions of the exact equations of motion for fluid systems at low to moderate Reynolds numbers. In the past 15 years, DNS has revolutionized the fundamental understanding of turbulent systems by providing a simple, inexpensive procedure for obtaining precise data on complex turbulent flows. DNS employs highly accurate numerical algorithms to solve the equations of motion on a discretized grid with as much spatial resolution as is computationally feasible. Increasing the resolution enables a simulation to achieve a higher Reynolds number, thereby closing the gap between simulations and most experimental conditions. The simulations presented here extend the Reynolds number envelope in three ways: (i)the simulations use a highly efficient Galerkin procedure based on a finite Fourier series representation of the velocity and scalar fields (Canuto et al., 1988), (ii) the simulations are performed in two dimensions using 5122points (256 wavenumbers), and (iii) a hyperviscosity formulation is used to limit the dissipation region (Maltrud and Vallis, 1991, 1993; Maltrud, 1992). The range of length and time scales that the present simulations resolve is sufficient t o identify fractal characteristics of the surface. Also, despite the well-known differences of the dynamics of 2D (two-dimensional) versus 3D turbulence (Maltrud and Vallis, 1991, 1993; Maltrud, 19921, it was shown in a previous study (Collins, 1995) that the response of the flame surface in the 2D system was similar to that in the 3D system (Dandekar and Collins, 1995) for a particular set of parameters. Thus, the 2D simulation increases the spatial resolution of the simulation while providing results that are meaningful for 3D systems. Because the flame is considered passive, the velocity field is assumed to be incompressible. It is therefore convenient to express the velocity field in terms of a stream function w as shown below:

u,=*

a

and u,,= - 9 &

A transport equation for the z-component of fluid vorticity defined as

is then given by

2590 Ind. Eng. Chem. Res., Vol. 34,No. 8, 1995 Table 1. Parameter Values, Wbulent Flame Speeds, and Fractal Dimensions for Eight Runs

where

J(5,W) represents the convective term for 5, F is a stochastic forcing function used to achieve stationary turbulence, D is a low-wavenumber filter required in forced, 2D simulations to eliminate the accumulation of energy in the first couple of modes (the accumulation results from the inverse cascade of energy in two dimensions, as suggested by Kraichnan (1967,1971)and Maltrud and Vallis (1993)), and n is an integer hyperviscosity coefficient. Hyperviscosities reduce the size of the dissipation range relative to the inertial range, thereby increasing the effective Reynolds number of the flow. In the calculations shown in the present study, the value of n is frxed a t 4 (note n = 2 corresponds to Navier-Stokes turbulence). The boundary conditions in the x- and y-directions are periodic for all variables. Stochastic forcing introduces energy over a narrow band of wavenumbers so as to mimic the source terms that are present in engineering flows due to mean flow gradients, while not significantly affecting the intrinsic turbulence processes (Eswaran and Pope, 1988). The advantage of forced turbulence is that it is statistically stationary, which is the simplest possible fluid mechanical environment for observing the dynamics of flame propagation. Two different forcing ranges are used in the present study. The large-scale forcing (hereafter the “LS” case) adds energy to the wave band defined by 10 Ik I 14. As shown by Kraichnan (1967,1971), largescale forcing leads to a forward entrophy cascade and an energy spectrum with a -3 exponent. A parallel study using small-scale forcing (hereafter the “SS”case) with a range of 80 Ik I84 has also been conducted and contrasted with the first. In the SS case, the inertial range is dominated by an inverse energy cascade that results in a spectrum with a more classical -5/3 exponent (Kraichnan, 1967, 1971). Comparisons between the flame response associated with the LS and SS simulations will demonstrate the importance of the shape of the energy spectrum (i.e., its inertial range exponent) on the topological shape of the flame surface (i.e., its fractal dimension). The flame surface is described by a scalar field variable, G(x,t),defined such that isocontours of G(x,t) represent individual flame fronts. G(x,t)is initialized with a uniform mean gradient that is maintained throughout the simulation. It is therefore convenient to express G(x,t)in terms of a periodic scalar function g(x,t),as shown below:

g ( x , t )= G(x,t)- x

(6)

Note that the flame is assumed to propagate on average in the x-direction (without loss of generality). The governing equation for the periodic function g(x,t)is as follows (Collins, 1995):

@ + J(g,v)at

where

and

* ay

= -uLO

(7)

run

velocity

y =VIUL

UtIUL

DK

1 2 3 4

ss ss ss ss

2.0 1.0 0.5 0.25 2.0 1.0 0.5 0.25

3.39 1.97 1.49 1.22 2.21 1.75 1.37 1.16

1.23 1.15 1.09 1.06 1.13 1.12

5

LS

6

LS

7 8

LS

LS

1.06 1.03

4 eXp(+$K) UL

= U L0

4

3 + exp(+5zK)

(9)

is the functional form of the laminar flame velocity, which depends on the local curvature of the flame surface defined by

K = -

&+& !%L?2+ 0

0”

Note that effects due to strain have been neglected in the above formulation. Equations 4, 5, and 7-10 constitute the governing equations for the velocity field and flame surface. The equations were solved using a pseudospectral Galerkin procedure (Canuto et al., 1988) for the spatial derivatives and a fourth order Runge-Kutta for the time integration. Simulations were run for several eddy turnover times until the flame surface reached its fully developed state. Table 1 shows the parameter values for the eight runs that were performed. Typical flame surfaces are shown in Figure 2 for the SS (a) and LS (b) cases, respectively. It is apparent that the SS flames are wrinkled over a wider range of scales than the equivalent LS flame surfaces, although the values of y are the same in the two simulations. Collins (1995) argued that the differences in the flame surfaces demonstrated the significance of the spectral distribution of turbulent energy on flame propagation. Thus, although the conceptual picture proposed by Damkohler (1940) may be accurate, eq 1 appears to be an oversimplification of the process of surface wrinkling because it does not realistically account for the multitude of length and time scales present in turbulent flows. Theories that account for this important feature of flame surfaces are required to provide a more fundamental basis for turbulent flame prediction. 111. Fractal Dimension

Mandelbrot, in his pioneering work on fractals (Mandelbrot, 1975, 19831, applied a geometrical construct used mostly in the mathematics community to a range of so-called “complex”physical problems. Mandelbrot observed that many surface geometries in nature were more accurately described by a fractional dimension (fractal) than by the classical Euclidean dimension. Fractal objects can be thought of as fundamentally rough in a self-similar fashion over a wide range of length scales. It is apparent, based on the physical

Ind. Eng. Chem. Res., Vol. 34,No. 8, 1995 2591 60

...............................................................................................

+

........... 1.......................

1 : ...........

L...

0 4 .a3

!

:...........

........... ........... ........... ........... L

I

................................................................................

3.00

1.50

6.00

4.50

X

6.0

Figure 3. Schematic of the computational grid used to calculate

, $.....

t

................................................................

i

............................ . . . . . . G . . . . . .

__ .......... _,

L

1

.............L...........

...................................

i

l

0

a3

I 50

300

4 50

600

X

Figure 2. Fully developed turbulent flame surfaces corresponding to G(x,t)= 2x13 and G(x,t)= 4x/3 for (a) SS and (b) LS simulations a t y = 2 (runs 1and 5, respectively). Notice that the SS surface is wrinkled over a wider range of length scales than the relatively smooth LS surface.

arguments presented earlier, that flame surfaces would be likely candidates for fractal analysis, since they are wrinkled over a range of scales. In this section we shall briefly describe the relevant features of fractal surfaces and the mathematical tools that can be used to describe them. Consider a surface that is embedded in a three dimensional space. Suppose we define a unit of measure E and the number of +units required to cover the surface is N(E). If the surface is Euclidean we would expect the following relationship to hold:

where the dimension DK = 2. (The subscript K refers to the Kolmogorov capacity. See section V for a precise definition of this terminology.) A fractal object is defined as one that has an exponent that is greater than its topological dimension. For example flame surfaces have been experimentally observed to have dimensions that lie in the range 2 I DK -= Z1/3 (North and Santavicca, 1990; Mantzaras et al., 1989), as has been observed with passive scalar surfaces as well (Prasad and Sreenivasan, 1990a,b;Sreenivasan and Meneveau, 1986; Sreenivasan et al., 1989). We can therefore conclude that flame surfaces are so wrinkled that they “fill” more than the two dimensions implied by their topological dimension. For flame “lines” embedded in a two-dimensional Euclidean space the arguments are

the fractal dimension. The gray boxes represent boxes that touch the flame surface (shown here as a line). To the right is an expanded view of an individual box. The black points are the places the function G(x,t)is evaluated. If the bounding points bracket the value in question, then the flame is assumed to pass through the box. A count of the number of boxes required to cover the flame, N(E),is recorded for each box size, c.

identical, except that the fractal dimensions are reduced by 1, i.e., 1 IDK < 11/3 for a 2D system. The method we use to measure the fractal dimension is based on the box counting technique suggested by Prasad and Sreenivasan (1990a,b). The 512 x 512 discretized values of G(x,t) represent 5122 boxes as shown schematically in Figure 3. For a particular flame surface described by G(x,t)= C, where C is a constant, one can determine the number of boxes through which the flame passes as a function of the box size. An individual box is assumed to cover the flame surface if the boundary points bracket the value C. The analysis is repeated for boxes of increasing size up to the system box size. The fractal dimension, DK,is then defined as the negative of the slope of a plot of In N(E)versus In(€). A surface is referred to as “fractal”if this slope is greater than the topological dimension (in this case, unity).

IV. Results From Simulations The fractal dimension of each of the eight runs are computed from the logarithmic plots shown in Figures 4 and 5. The results are summarized in Table 1. As is clear from the figures, the simulated surfaces yield relatively straight lines and well-defined fractal dimensions based on the algorithm summarized in section 111. One observation that is immediately apparent is the fractal dimensions of the SS runs are consistently greater than the equivalent LS runs. Figure 6 shows a plot of the fractal dimensions observed in the SS and LS simulations compared with experimental measurements by North and Santavicca (1990). The solid line is the empirical relationship derived by North and Santavicca to correlate their data. It appears that the SS runs more closely resemble the values from experiment while the LS results are consistently below the measured values. The discrepancy between the two sets of runs is consistent with the results from our earlier study which demonstrated that the flame velocity and scalar spectrum were influenced by the shape of the energy spectrum and not just the magnitude of the turbulence intensity (Collins, 1995). Furthermore, the closer agreement of the SS runs with the experiments can be attributed to the similarity in the inertial ranges.

2692 Ind. Eng. Chem. Res., Vol. 34, No. 8, 1995

-

1.25

1.20-

d

-

1.15

1107

t

I

“.Q

I 2

I 6

I 4

I

R

.

0

IC@{/ 00

In D

5

y

:

I

I

I

I

I

1

os

10

15

20

25

TO

Y

Figure 4. Plot of In N ( E )versus In E for the four SS runs (runs 1-41. The lines passing through each set of points are linear fits to the data. Note the data (and the relevant linear curve) corresponding to y = 0.5, 1.0, and 2.0 are shifted vertically by 1, 2, and 3 respectively for visual clarity.

I+ -

1

v=2

1

I

Figure 6. Fractal dimension of SS (solid circle) and LS (solid diamond) flame surfaces as a function of the parameter y . The open diamonds are experimental measurements from North and Santavicca (1990), and the solid line is their proposed correlating curve.

power-law exponent of the inertial range of the spectrum that defines that surface. The numerical simulation results described in section Tv provide us with an opportunity to explore this relationship, since both the spectrum and fractal dimension can be readily computed from the simulation database. In this section, we consider the relationship between the scalar inertialrange exponent and the fractal dimension of the surface, taking into consideration the previously reported results for flame surfaces. First, it will be necessary to clearly distinguish among various “fractal dimensions” defined in the literature. In particular VH identified two mathematically distinguishable dimensions. The first is called the Kolmogorov capacity and is equivalent to the fractal dimension we have been discussing in the previous sections. The capacity of a surface is the value DK shown in eq 11in the limit as E 0. Hereafter we shall adopt the terminology of VH and refer to a surface as a “K-fractal” if the dimension DK is greater than the topological dimension. The second definition of a fractal is based on the Hausdorf dimension, which is defined by the following expression:

-

I

I

2

4

I

I

6

R

In e

Figure 6. Plot of In N ( E )versus In E for the four LS runs (runs 5-81, The lines passing through each set of points are linear fits to the data. Note the data (and the relevant linear curve) corresponding to y = 0.5, 1.0, and 2.0 are shifted vertically by 1, 2,and 3 respectively for visual clarity.

North and Santavicca’s experiments of freely propagating flames in three dimensions had an energy spectrum with a k-5/3 inertial range (see Videto and Santavicca (1990) for details on the energy spectrum), which is closer to the SS runs in this study.

V. Relationship Between the Fractal Dimension and the Scalar Spectrum Vascillicos and Hunt (1991; hereafter VH) and a related paper by Hunt and Vassilicos (1991) considered the topological characteristics of surfaces and how those features were related to the spectrum of the defining function. Their analysis lead to a linear relationship between the fractal dimension of a surface and the

(12)

-

where Ei is constrained by Cis E . The limit of the above expression as 6 0 is either 0 or 00 depending on the value of D. The critical value of D is referred t o as the Hausdorf dimension and is designated as DH.A noninteger value of DH implies a surface that is an “Hfractal”. The Hausdorf dimension is the dimension identified by Mandelbrot in his work that originally popularized the term “fractal”. K-fractals are a broader class of fractals than H-fractals. K-fractals imply local self-similarity while H-fractals require global selfsimilarity; thus it is possible for a surface t o have a noninteger value of DK even though its Hausdorf dimension DH is equal t o its topological dimension. It should be noted that the box-counting technique described above calculates the Kolmogorov capacity. Experimental investigations have relied on the box-

Ind. Eng. Chem. Res., Vol. 34, No. 8,1995 2693 counting techniques, and thus it is known that many surfaces observed in turbulent flows are K-fractals; however, there is no experimental evidence to confirm the hypothesis that H-fractals exist in turbulent flows. VH went through an extensive analysis to develop a relationship between the Kolmogorov capacity DK and the inertial-range exponent of the scalar spectrum. Their analysis was restricted to a scalar function that was discontinuous at the interface (i.e.,g(x-b,t) = 0 and g(x+S,t) = 1,where S is an infinitesimal vector pointing in the direction normal t o the surface). Furthermore, they stated that the scalar function must be homogeneous and isotropic for their analysis to be valid. By considering the correlation between two points in the fluid g(x)g(x+r),they predicted that a scalar spectrum with an inertial range power law of the form

r(K)

K-Q

(13)

has a Kolmogorov capacity given by

D,= 1+ d - 9

(14)

where d is the topological dimension (d = 2 in the present simulations) and r(k)is the classical scalar spectrum defined by

r(K)= -Ld-' d-1

lg(x)g(x+r)e-ik'r d r

(2nId-'

According to eq 14, the fractal dimension for the SS runs (g = 5/3) is 4/3 and for the LS runs (q = 3) is unity (note, the fractal dimension cannot be less than the topological dimension). There are several discrepancies between this result and the fractal dimensions found in the numerical simulations. Firstly, although the trend shown in eq 14 is correct, the fractal dimension for the SS case is overpredicted while that for the LS case is underpredicted. Secondly, the simulated flame surfaces and experimental measurements of flame surfaces (North and Santavicca, 1990) demonstrated a sensitivity of the fractal dimension to the parameter y that is not reflected in eq 14. (It should be noted that the agreement between eq 14 and the SS simulations is best for the largest value of y . We shall return to this point momentarily.) We now consider each question. Effect of the Inertial Range Exponent. In order to gain additional insight into the relationship between the spectrum r(k)and the fractal dimension of its isocontours, a numerical experiment was performed using a spectrum of the following generic form:

(15)

where

r(k)is

the spectrum normalized such that G'r(k)dk = 1 leading t o

The parameter km defines the location of the maximum in the spectrum and q is the inertial-range power-law exponent at large wavenumbers (Le., fork z=k,, T(k) k-4). Based on the above spectrum, it is possible to define the complete spectral function &(k)by random-

-

izing the phasing as shown below (Eswaran and Pope, 1988):

where 8 is a random number uniformly distributed between 0 and 1and m(k)is the number of points within the shell k - (AM21 Ik k (AM2). The function defined by eq 17 will have a spectrum that is identical to eq 16 and will be homogeneous and isotropic. In computing the fractal dimension of the random scalar functions defined above, several realizations were used to assure statistical convergence. The study considered the fractal dimensions of spectra with values of q = 1,2,5/3, and 3 and values of Km = 1,8,16, and 32. The results are shown collectively in Figure 7. It is interesting to note that a linear region is observed for all simulations at low values of In E ; however at larger values of E there is curvature. Furthermore, the degree of curvature appears to increase with increasing km. Unfortunately, the presence of curvature implies some ambiguity in determining the fractal dimension, but since all curves have nearly the same slope a t small values of In E , we have taken that slope as the fractal dimension. Figure 8 shows the K-fractal dimension as a function of the power-law exponent q along with the curve predicted by the VH analysis. The error bars reflect the uncertainty resulting from varying the value of k,. There appears to be a linear relationship between the fractal dimension and q as predicted by VH;however the slope of the line is quite different than their prediction. From this numerical experiment, we obtain the following empirical relationship for the fractal dimension:

+

D, = d - 0.42 - 0.159

(18)

(shown with a straight line). It should be noted that the VH theory (eq 14) is in agreement with eq 18 for a classical k-5/3 spectrum; otherwise the disparity is significant. Effect of y. The question of how the parameter y affects the fractal dimension of the flame is not resolved by the above numerical simulations or by eq 18. For example, the dimension predicted by eq 18 for the LS runs (setting q = -3) is reasonably close to the values observed in the direct numerical simulations a t y = 2 and y = 1(runs 5 and 6);however a t smaller values of y the fractal dimension is observed to decrease. Likewise, the fractal dimensions for the SS runs show a similar trend with y (runs 1-4). The experimental measurements by North and Santavicca (1990) are consistent with this finding, as noted above. Furthermore, their model of the fractal dimension predicts a smooth transition from the Euclidean dimension at small y to the passive scalar value of 2'/3 in the limit y 00, In this subsection, we attempt t o explain this transition in terms of changes in the scalar spectrum. Our previous DNS study of flame propagation (Collins, 1995) demonstrated that the scalar spectrum was sensitive to changes in the value of y . It is therefore useful to consider the spectrum of the flame surfaces a t different values of y . Figure 9 shows a semilog plot of k5l3r(k)vs k for the energy and SS scalar at three values of y (2, 1, and 0.5). A classical inertial range would appear as a horizontal line in these coordinates. The energy spectrum (solid line) shows a very well defined inertial range up to a wavenumber of 80, which

-

2594 Ind. Eng. Chem. Res., Vol. 34, No. 8, 1995

18 16

-2 -

14

12

e

IO

8

6

0

2

6

4

8

0

2

4

6

8

In E

In E 16 14

12

-

10

-

R

t

6 4

I 0

I

I

1

I

I

I

2

4

6

8

1

2

I 4

In E

1

I

6

R

In E

Figure 7. Plots of In N ( E )versus In E for numerically “generated”surfaces with inertial range exponents of (a) 1.0, (b) 1.67, (c) 2.0, and (d) 3.0. The plus, circle, triangle, and diamond correspond to different values of k,. The linear fit is based on the lowest three values of In e . 2.(

18

1.6

O :i

14

I .2

1 in I

I

1.5

2.0 9

I 25

I

Figure 8. Fractal dimension of “artificially generated” surfaces as a function of the inertial range slope. The solid line is the empirical formula D K = 1.57 - 0.15q, and the dashed line is the formula derived by Vascillicos and Hunt (1991).

corresponds t o the beginning of the range of forcing. In contrast, the scalar spectrum at y = 2 has only a small inertial range (through wavenumber 201, and those at y = 1 and 0.5 have no distinct inertial range at all. Apparently there is a systematic decrease in the size of the inertial range as y is decreased. Physically, the

I

I

1

I

20

40

bo

I

k

so

I I on

1

120

Figure 9. Scalar spectrum for SS runs ( y = 2, 1,0.5) multiplied by k5I3.A classical inertial range would appear as a horizontal line in these coordinates. The energy spectrum is also shown for comparison (solid line). Notice the well-defined inertial range in the energy spectrum.

effect of lowering y is to increase the rate of flame front propagation (for k e d U‘),thereby increasing the surface “smoothing” that results from propagation, causing a net reduction in the size of the inertial range. This is analogous to the effect that the diffusivity has on the

Ind. Eng. Chem. Res., Vol. 34, No. 8, 1995 2595 20

IS

10

-

Y,

e

05

00

-05

25

30

35

40

45

50

55

Ink

Figure 10. Linear plot of the function @(k) (see eq 20 for the definition) for SS runs with y = 2, 1, 0.5. Least squares fits of each function (solid lines) yield slopes of 0.99, 0.99, and 0.95, respectively.

turbulent transport of a passive scalar. In that case, increasing the diffusivity of a passive scalar increases the dissipation range of the scalar spectrum a t the expense of the inertial range. However, it should be noted that the effects due to diffusion are fundamentally different from those due to propagation because the former is a linear process while the latter is highly nonlinear. The analogy is only intended to be qualitative. A second feature of the scalar spectrum apparent from Figure 9 is the exponential dependence of the dissipation range. On the basis of this result and previous studies (Collins, 1995; Peters, 1992), we propose the following modified flame spectrum:

where K D is the characteristic wavenumber of the dissipation range and p is the power of the exponential tail. For turbulent flames at high Reynolds number and large y , we expect that k d k m >> 1 and that K D increases with increasing y . Several theoretical values for the exponent p have been proposed in the literature (Peters, 1992; Collins, 1995). We elected to determine its value from the numerical simulations by plotting

@(k) = ln{

-hrq]} vs In K

(20)

A straight line at large In k will emerge with a slope of p if eq 19 is accurate. Figure 10 shows the data from the same three SS runs ( y = 2,1,and 0.5). The values of p as determined from a linear regression are 0.99, 0.99, and 0.95, respectively. We conclude that the exponent is approximately unity. We extend the analysis performed by VH to the more general spectrum shown in eq 19, withp = 1. In this case, the results are found t o be sensitive t o the range of length scales of interest. That is, the fractal dimension valid for length scales in the inertial range differs

from the dimension found for the dissipation scales. For the spectrum shown in eq 19, length scales within the inertial range (i.e., l / k ~