References

' Pedras, M. H. J. and de Lemos, M. J. S.,
“Macroscopic turbulence modeling for incompressible
flow through undeformable porous media”,
International Journal of Heat and Mass Transfer, 44,
yp. 1081-1093, 2001.
* Pedras, M. H. J. and de Lemos, M. J. S., “Simulation
of Turbulent Flow in Porous Media Using a Spatially
Periodic Array and a Low Re Two-Equation Closure”,
Numerical Heat Transfer, 39, pp. 35-39, 2001.
> Pedras, M. H. J. and de Lemos, M. J. S., “On the
Definition of Turbulent Kinetic Energy for Flow in
Porous Media”, Comm. Heat Mass Transfer, 27,
pp.211-220, 2000.
4+ Nakayama, A. and Kuwahara, F., “A Macroscopic
Turbulence Model for Flow in a Porous Medium”,
Journal of Fluids Engineering, 121, pp. 427-433, 1999.
5 Kuwahara, F., Kameyama, Y., Yamashita, S. and
Nakayama, A., “Numerical Modeling of Turbulent
Flow in Porous Media Using a Spatially Periodic
Array”, Journal of Porous Media, 1, pp. 47-55, 1998.
® Antohe, B.V. and Lage, J.L., “A general two-equation
macroscopic turbulence-model for incompressible flow
in porous media”, Int. J. Heat Mass Transfer, 40, pp.
3013-3024, 1997.
7 Kaviany, M., Principles of Heat Transfer in Porous
Media, Springer-Verlag, New York, NY, 1991.
® Anderson, T. B. and Jackson, R., “A Fluid Mechanical
Description of Fluidized Beds”, /&EC Fi undamentals,
6, 1967.
° Gray, W. G. and Lee, P. C. Y., “On the Theorems for
Local Volume Averaging of Multiphase Systems”, /nt.
J. Multiphase Flow, 3, pp. 333-340, 1977.

© Gough, P. S. and Zwarts, F. J., “Modeling
Heterogeneous Two-Phase Reacting Flow’, A/AA
Journal, 17, pp. 17-25,1979.

"l Kuo, K. K., Principles of Combustion, A Wiley-
Interscience Publication, John Wiley & Sons, New
York, NY, 1986.

2 Slattery, J.C., “Flow of viscoelastic fluids through
porous media”, A./.Ch.E. J., 13, pp. 1066- , 1967.

'3 Slattery, J. C., momentum energy and mass transfer
in continua, second edition, Robert E. Krieger
Publishing Company, Inc., 1981.

' Kaplan, W., Advanced Calculus, Edition 3",
Addison-Wesley Publishing Company, Reading, MA,
pp. 257-259, 1984.

'S Ghosal, S. and Moin, P., “Basic Equations for the
Large Eddy Simulation of Turbulent Flow in Complex
Geometry", J. Comput. Phys., 118, pp. 24-37, 1995.

'© Ghosal, S., “Mathematical and Physical Constraints
on Large-Eddy Simulation of Turbulence”, AIAA
Journal, 37, 1999.

"7 Vailyev, O. V. , Lund, T.S. and Moin, P., “A General
Class of Commutative Filters for LES in Complex
Geometries,” Journal of Compuational Physics, 146,
pp. 82-104, 1998.

'S DesJardin, P.E and Gritzo, L.A., “On the
Development of a Clutter Model for Turbulent Flows”,
SAND Report (in preparation), 2002.

19 Tennekes, H. and Lumley, J.L., A First Course in
Turbulence, MIT Press, Cambridge, MA, 1972.

2 Pope, S.B., Turbulent Flows, Cambridge University
Press, New York, NY, 2000.

”| Brinkman, H. C., "A Calculation of the Viscous
Force Exerted by a Flowing Fluid on a Dense Swarm of
Particles", Appl. Sci. Res., Al, pp. 27-34, 1947.

” Brinkman, H. C., "On the Permeability of Media
Consisting of Closely Packed Porous Particles”, Appl.
Sci. Res., Al, pp. 81-86, 1948.

3 Blevins, R. D., Applied Fluid Dynamics Handbook,
Van Nostrand Reinhold Company Inc., New York, NY,
pp.333-334, 1984.

** Hoerner, S. F., Fluid-Dynamic Drag, published by
author, New York, NY, 1957.

25 Forchheimer, P. Z., Z. Ver. Deutsch Ing., vol. 45, pp-
1782, 1901.

6 Eroun,, S., “Fluid flow through packed columns”,
Chem. Eng., Prog., 48, 89-94, 1952.

27 Holen, J., Brostom, M. and Magnussen, B.F., “Finite
difference calculation of pool fires”, In Proceedings of
23"4 Symp. (Int.) on Combustion, 1677-1683, The
Combustion Institute, Pittsburgh, PA, 1990.

28 Jones, W.P. and Launder, B.E., “The prediction of
laminarization with a two-equation model of
turbulence”, Int. J. Heat Mass Transfer, 15, 301-314,
1972

>? Patankar, S.V., Numerical Heat Transfer and Fluid
Flow, Hemisphere Publishing Co, New York, NY,
1980.

* Gritzo, L.A., Nicolette, V.F., Tieszen, S.R., Moya,
J.L. and Holen, J., “Heat transfer to the fuel surface in
large pool fires”, In Chan, S.H., editor, Transport
Phenomena in Combustion, pp. 701-712, Taylor &
Francis, 1995.

3! Gritzo, L.A. and Nicolette, V-F., “Coupling of large
fire phenomenon with object geometry and object
thermal response”, Journal of Fire Sciences, 15, 427-
442, 1997.

® CFD-ACETM, CFD-ACE Structured Flow Solver
Manual, Version 5, CFD Research Corporation,
Huntsville, Alabama, October, 1998.



For uniform filter size (ie. A, =Cconst), then

aG/dx, =-0G/ox’, and using the chain rule of
differentiation then Eq. (7) can be expressed as:

ra we J {4

The first term on the right hand size of Eq. (8) is the
two-phase average of the gradient of the vector quantity
and the second term can be further expressed in terms
of a surface integral using the divergence theorem.
Rearranging terms leads to an expression for the phase
averaged divergence of a vector quantity in terms of the
divergence of a phase averaged quantity and a surface
integral to account for microscopic clutter effects.

5B,@)_aB (x) xx 4°
= Fe *fo 7 pena 9)

OB (x/) a ,
ox, aie @ x

Ox; 0:

An important assumption in the derivation of Eq. (9) is
the filter width is assumed to be constant. That is, the
filter width is assumed to be invariant with space and
time. Alternatively, the development of a non-uniform
filter width could be pursued for which

aG/ox, # -dG/dx;, . However, the introduction of

such a filter introduces commutation error that has been
discussed in the Large Eddy Simulation (LES)
literature'®'*"” as well as an extra surface integral error
term'®, The former of these two sources of error has
been somewhat mitigated using filter functions that
satisfy commutivity up to some order of accuracy",
However, the second source of error is more difficult to
treat and leads to a violation of mass, momentum and
energy conservation when these quantities are
exchanged between phases'*. These issues affect
practical application. Often the filter width is assumed
to be proportional to the local grid size, eg.

A, =CA,, where A, is the local grid cell size and C

is a constant that is greater or equal to one. Therefore
in a practical engineering problem where the CFD grid
cell size changes in space, then the assumption of

A, =const is violated! In order to avoid this

contradiction, the filter width is assumed to be
independent from the grid. This choice also has the
added benefit to allowing for grid convergence (see
requirement 4 of the following discussion). The
disadvantages of a filter width that is independent of the
grid are two fold. First the clutter model must be able
to capture the relevant physics over a wide range of
filter to clutter size ratios (see requirement | of the

following discussion). Secondly, the coding of the
clutter model is further complicated with the
introduction of an integral of the clutter source terms to
account for the effects of intermittency (see
requirement 3 of below discussion).

Thus far, only phase averaging concepts have been
briefly reviewed. In addition, it is also desirable to time
average the governing equations of mass, momentum
and energy transport. Time averaging removes any
remaining high frequency turbulence signal and is
compatible with the Reynolds Averaged Navier Stokes
(RANS) transport equations that are commonly used for
solving practical engineering applications. | Time
averaging issues are not reviewed in this paper but can
be found in several excellent introductory texts on
turbulence!*””,

Phase and Time Averaged Transport
Equations

Applying the phase averaging relations given by Eqs.
(5) and (9), in addition to time averaging, to the
transport equations of momentum, turbulent kinetic
energy and turbulent dissipation results in a set of
phase-and-time averaged transport equations containing
unknown second order correlations and surface integral
terms. Both of these terms require explicit subgrid
scale modeling. In order to develop models for these
terms, a set of guiding principals are employed.

Guiding Principals for Clutter Subgrid Scale Model

The following four principals provide guidance for the
development of the SGS clutter model. These
requirements are summarized below with a diagram
illustrating each constraint for a two-dimensional
domain.

1. The clutter model(s) should be able to
accommodate a spectrum of clutter sizes ranging
from a single object (i.e. bluff body drag) to many
objects that are randomly placed (i.e. clutter that
resembles porous media) since the filter width of the
jing is assumed to be constant.

walla, <1

@®O

Range in clutier Sizes



Once Fc is calculated then the source term Szis

determined using Eq. (24) and summarized below for a
Cartesian coordinate system.
S,, =Nefetee =

Nop = = \P
n a, Pre Min tea Ki)
fos6cosylucos@cosy + vcosOsiny + wsin @]} (33)

+c A, \i)+(i)) Mi))
in @cosy|usin @cosy + vsin Asiny — wcos8}
+C,,A, (@,,) Sinylusiny - vcosv

prea “Al ii) + )+@a)f NG il

foxdsinyucondeony + veososny + wsino} (G4)
+C,A:)+(@,)) Miz)
$indsiny[usin@ cosy + vsin@siny — weosO}}
+C,,Ay (7) Ecosvlesiny veosy]}

NF ve

a) May)

Sp, = Mert = Ree prea CA a) + (
(35)

7 a,
$in6[ucosO cosy + casings wsino}}

+ CAL) +(Es)) lea)
Los6[usin @ cosy + vsin@siny — weosa

Conservation of Turbulent Kinetic Energy and
Dissipation Rate

Details on the derivation of the phase averaged k and &€
equations can be found in Ref. [1] and the modeled

form of them are summarized below for use with the
clutter model.

A} t0)+-2 (ey)
[eure ou)

= 2(6(i,))

ul 7) Me) - pe)

IW) (k) C,,S,, +C.,5o, "Ola, )aav

Ay i,

(36)

{rr Wap 2 oem? a)
~Cyplid, 7 ayaa # cod (37)

+ JeaAete,S, +6, falda) )av

Ari,

The unique feature of the closed modeled terms in Eqs.
(36) and (37) is again the use of the averaged source
terms over the filtering volume to account for the

effects of non-periodic clutter. The constants C,, and

C,, require calibration for the limits of porous media

and bluff body drag,
associated with the porous media limit, C,, , is set toa

respectively. The constant

value of 0.28 to be consistent with the results of Pedras
etal. ', The second constant, C,, , will eventually be

set based on the results of detailed measurements and
CED predictions, to be discussed in the results section.

Numerical Implementation

The clutter models discussed in the previous section
are implemented into a general-purpose fire simulation
code, VULCAN, which is based on the KAMELEON-
Fire code”’. VULCAN uses a RANS based model suite
including a k-e turbulence model’® Numerical
discretization is on a staggered, block-structured grid
with second-order upwind differencing for the
convective terms using a version of the SIMPLE
algorithm”. Previous studies have included using
VULCAN for pool fire simulations*TM*!, The source
terms appearing in Eqs. (12), (13), (17) and (33)-(35)
for momentum transport are treated in a semi-implicit
nature to maintain numerical stability. The source
terms appearing in Egs. (36) (36)and (37) are
implemented in an explicit manner.

Results

Assessment of Linear Blending of

Ky pand Ky ‘gAssumption

An assessment of using a linear blending of constitutive
models from the porous media and bluff body drag is
investigated for a collection of spheres where the clutter
is considered isotropic. Figure 4 (a) and (b) show
results of the drag force per unit volume using a
modified Darcy-Forchheimer Law and also a bluff body


dispersion and turbulent fluctuations ~—_ using
Brinkman”, 44g, and turbulent eddy”, £4, , viscosity

models, respectively. The source terms, Sp, and Sy ,
sper 'y. P; B;

represent the gas phase exchange at the surface of the
solid for the limits of a porous media and bluff body
drag, respectively. The exact functional form of these
terms will be discussed shortly for the isotropic and
anisotropic clutter classes. These two source terms are
blended depending on the size of the clutter relative to
the filtering volume. Recall, this blending is employed
to satisfy first constraint condition for of the clutter
model formulation previously defined. The linear
weighting used is not unique and other weightings
could have been chosen. An example of using a linear
blending of constitutive models for the isotropic clutter
class is presented in the results to assess the utility of
this approach. In addition to the blending, these source
terms are also averaged over the filtering volume. The
motivation for introducing this averaging is to satisfy
the third constraint condition for the clutter model to
account for non-periodic clutter.

Isotropic Clutter Class

The first type of clutter class considered is for isotropic
media that is represented by a collection of spheres.
Spheres are chosen to represent this class due to the
abundant amount of drag data available in the literature
over a wide range of Reynolds numbers. For this class,
the effects of the clutter are modeled as a combination
of a modified Darcy-Forchheimer Law for porous
media® and a bluff body drag relations for the linear
superposition of a collection of spheres.

5, =, u/K +Copq VK] 02)

8, =@/9)K@ Ku )l-ey/D ay

In Eq. (12), KE 119(u,)/@p/ax, ) is defined as the
permeability and is determined using the empirical
relation K =¢°D? Hea - | where C, is set
equal to 147 based on the work of Kuwahara et al. * .
The value of C, is set to either 1.75 as suggested by
Kuwahara et al.’ or to 2.3 based on the classical result

of Ergun’®, Calculations presented in the results shows
that the total drag from a collection of spheres shows
the overall drag is not sensitive to the value of C,. In
the bluff body limit, the clutter is treated as a collection
of spheres for which the total number of spheres in the
filtering volume can be related to porosity and ratio of
clutter to filtering size using the relation:

N, =6(1-@XD/A, J /7. In Eq. (13), Cp is the

coefficient of drag of each sphere given by the
following empirical relation.

ce fet +Re,"/6)/Re, for Re, $1000
_

(14)
0.424 for Re, > 1000

Anisotropic Clutter Class

In order to satisfy the second constraint for the clutter
model, an anisotropic clutter class is addressed. The
orientation of the flow relative to the clutter will greatly
affect the drag due to the change in projected area as
well as any lift force components that are normal to the
flow direction. The geometry chosen to represent the
anisotropic clutter class is a collection of cylinders for
which the force on each cylinder can be decomposed

into directions along the length of cylinder, &. and

normal to cylinder, @,, and @,,, as illustrated Figure

y Net drag and lift forces on cylinder, Fe

Figure 2: Decomposition of net force using a
coordinate system attached to the cylinder and a
fixed Cartesian coordinate system.

The relationship between the transformed coordinate
system, Canoe) , and the fixed Cartesian

coordinate system, (..é,,2, ) is given by following

rotation of basis coordinate transformation,

a, 8] [eosacosy cosdsiny sind V2.) (45)
a, b=748, b=| -siny cosy 0 Re
a, é,| |+sin@cosy —sin@siny cose ||é,

where the direction angles6 and y are defined using
the following relations.


2. The clutter model(s) should provide a mechanism
to represent isotropic (i.e. rotationally invariant) as
well as anisotropic (rotationally dependent) clutter
with minimal amount of user input.

Isotropic clutter Anisotropic clutter

3. The clutter model(s) should allow for
intermittency within the filtering volume so that
non-periodic SGS clutter may be addressed.

Periodic clutter

CLO

So] lo

4. The clutter model(s) should allow for grid
convergence studies to be performed. This implies
that the numerical grid should be chosen
independent of the phase averaging filtering volume.

Non-periodic clutter

Imposing these constraints guide the final form for the
closure terms needed for the phase and time averaged

transport equations. The descriptions of these models
are broken into two classes of clutter models. The first
is the isotropic class where the clutter is treated as a
collection of spheres. The second is the anisotropic
class for which the clutter is treated as a collection of
cylinders.

Conservation of Linear Momentum

Phase and time averaging the Navier-Stokes equations
for momentum transport and assuming constant density
and molecular viscosity results in the following,

02 6())- 2 ba,)|-
_a(z)), 2%) _ able, Me

Ox, ax,

-J P6, -7, kar

co

where the (--) notation refers to a time-and-phase

The term (i, \(z,) - (ua)
represents the effects of the subgrid fluctuations in time
and space and can be decomposed into three terms
using the double decomposition approach of Pedras et
al’. This term physically represents the effects of i)
turbulent Reynolds stresses due to macroscopic phase
averaged velocity, ii) dispersion due to spatial
fluctuations of macroscopic time averaged velocity and
iii) dispersion due to microscopic time averaged
velocity. The first two of these effects are modeled
using a turbulent eddy viscosity? and Brinkman”)?
models, respectively. The dispersion of momentum due
to microscopic time-averaged velocity is simply
neglected due for lack of an established model. The
closure for the surface integral term is based on a linear
weighting of constitutive models from the porous media
literature’’ and simple bluff body lift and drag
models”*”*, A weighting function,  , is defined as the

ratio of the length scale of the clutter to the size of the

averaged quantity.

filtering volume, ie. W=l-/A,, and is used to

develop a modeled transport equation for Eq. (10) as
follows.

A] 2 6.))-;2-rane))
ao) abr, (a)4@,))_

Tox, Ox, (11)
«alot o0-vnsf ED. 5]

The first term on the right hand side of Eq. (11) are
models to represent the effects of microscopic

~ gi flt-w, +s,


flow. Approximations have to be introduced for
cylinders exposed to an arbitrary orientation relative to
the mean flow. This need is addressed by assuming

that the total force on the cylinder, Fo is the result of

the superposition of the lift and drag along and normal
to the cylinder, shown in Fig. 3, in addition to a drag
force from gas flowing across the cylinder for which

C,, is set to a value of 1.2. Making these

approximations, the total net force on the cylinder may
be expressed as:

Fo= Fe+ Fo, + Fe,
~APIAG «(af GA

+ Pius, Ya)

+ PICA 7). IF Gs) Ga)

where Ay, A,, and A,,are the reference (projected)

(23)

areas of the cylinder and are equal to ADL, DL and
DL, respectively. The net force on the cylinder or
group of cylinders is related to the momentum source

term, Ss . of Eq. (11) using the relation,

S,=NcFe [Ay (24)
where Ne is the total number of cylinders that are
prescribed within the filtering volume.

In order to determine these forces, the angle @ is

required for the coefficients of lift and drag and is
obtained using the relation:

G)-(G)+ Ga) (25)
Kak) + (Ha) }

where (i) (ii) and (a2) are the components of the

a@=cos"

time and phase-averaged velocity in the 2.22)

coordinate system have the following definitions.

(a) = Me
(i,,)= Wey (26)
(i,2) Sue.

The magnitudes of the velocities in Eq. (26) can be

expressed in terms of velocity components in a
Cartesian coordinate system using the transformation.

uy u
uy e=Tyv
Ws, w (en)
cos@cosy cos@siny sin@ |{u
=| —sinw cosy 0 v
-sinOcosy —sin@siny cos@||w

Substituting Eqs. (26) and (27) into Eq. (25) and
after significant algebraic manipulation then @ can be
expressed in terms of known scalar velocities, u, Vv,

w and the angles 6 and y .
= COS" Cyn ! amon) (28)
where,

aun = |ulcosly — 8 )+ cosy + 4)]}+ vfsin(y —8)+ sin(y + 4)}+ 2wsin 9|

and
Cjnm = Ole? + v2 + 2w?+ (2 =v? cos(2y)+ 2uv sinay)]}? +

In obtaining the numerator and denominator of Eq.
(28), the magnitude of the various velocity components
and sums of components expressed in terms of a
Cartesian coordinate system is required. These terms
are summarized below for completeness.

a,)|=Ileos & - 8)+ cos Y +8)
+ vfsiny - @)+ sin w +0)] (29)

+2wsin @)|/2
Ki | = |veosy- usiny| (30)

\@.2)| =[alsin &w -6)-snw+e)]
+ vicos y + @)- cos y -0)] (31)

+2weos @)|/2

(a) +(@a) =e 4° 20")
+(u -v* Jcos2y) + 2uvsiny))2}"

32)


AIAA 2002-0984

On the Development of a Subgrid Scale Clutter Model

P.E. DesJardin, J. M. Nelsen, L. A. Gritzo
Sandia National Laboratories
Albuquerque, NM 87185

T. A. Ghee
NAVAIR
Patuxent River, MD 20670

40" AIAA Aerospace Sciences Meeting and Exhibit
January 14-17, 2002
Reno, NV

For permission to copy or to republish, contact the copyright owner named on the first page.
For AIAA-held copyright, write tolAIAA Permissions Department,
1801 Alexander Bell Drive, Suite 500, Reston, VA, 20191-4344.



Periodic porous
z= _
oO

Flow direction

Figure 6:
problem.

Schematic of periodic porous media

Figure 7 (a) and (b) show predictions of normalized
turbulent kinetic energy and its dissipation rate as a
function of downstream locations using the
macroscopic clutter model. The initial kinetic energy
and dissipation rate are chosen to be 10 and 30 times
the steady state values, respectively, in order to match
the previous studies of the same problem from
Nakayama et al.* and Pedras et al'. Predictions using
the clutter model agree well with these previous studies
in this limiting where the clutter model degenerates to
just the porous media contribution.
= -
9 @ Nakayama ot al. 1999

B—_-Pedras ot al. 2001
Present Result

<be/<k>

(a)

2 @ = —Nayayama et al., 1999,
@ —_Pedras et al., 2004
Present Results

‘0!

<g>/<g_>

ld,
(b)
Figure 7: Normalized turbulent kinetic energy and

its dissipation rate verses x/A; using clutter model in
the limit of porous media.

Single Clutter Element in Turbulent Flow

The second problem under consideration is a cylinder in
well-characterized homogeneous turbulent _ flow.
Experimental measurements are taken at the U.S. Naval
Academy (USNA) using their low-speed (0 to 20 m/s)
tunnel as shown in Fig. 8, for which high turbulence
intensities could be generated at the inlet. Hot wire

anemometry is used to measure instantaneous stream-
wise velocity to obtain mean and RMS profiles before
and after the cylinder.

PVC grid to generate well characterized
homogenous turbulent flow conditions

(b)

Figure 8:: Pictures of (a) USNA wind tunnel facility
and (b) blow up region of cylinder representing
clutter object under highly turbulent flow
conditions.


drag relations using Eqs. (12) and (13), respectively.
The results show that the porous media and bluff body
drag have a similar functional dependence on Reynolds
number with the two predicting almost the same force
per unit volume for either low porosity at high
Reynolds number or for high porosity at low Reynolds
number. This suggests that the linear blending of these
two limits is a reasonable approach to span the possible
range of clutter sizes for the isotropic clutter class.
Further work is required to see if the same approach is
also appropriate for the anisotropic clutter.

Porosity = 0.1

wt
wo
wt —— Bluff body drag
. Porous flow drag, Coy £1.75
o®” Porous flow drag, C,,, =2.3
Ew ‘
ow
es
2 ww
S
LL
10
107
10° 7 A
10° 0 10) ie 10° 10° 10°
Re,
(a)
Porosity = 0.9
10!
10
ie —— Bluff body drag
. ~ Porous flow drag, C,,,
@ 10 - Porous flow drag, C,
Ee
ow
a
2
9)
wo
10°
10°
10"
10" 10" 10° 10° 10° 10° 10°
Re,
(b)
Figure 4: Force/volume using constitutive models

for a collection of spheres from porous media and
bluff body drag relations using a porosity of (a) 0.1
and (b) 0.9.

Flow direction —>

Verification of Anisotropic Clutter Class for Single
Cylinder in a Cross Flow

The anisotropic clutter class is developed based on
empirical relations from Blevins” for lift and drag on
an inclined cylinder from the relations normal and
aligned to a cylinder shown in Figure 3. In the limiting

of v=w=0, u=-U,, w=Oand
6 = Q then the drag and lift force in Eqs. (33) and (35)
may be directly compared with the result from

HoernerTM for drag and lift on a cylinder at an angle of
attack equal to @, illustrated in the sketch below.

case

C,, = 1.1sin’? acosa

Cy = 1.1sin® @ +0.02

——S

Figure 5: Sketch of cylinder at angle of attack, a ,
and coefficients of lift and drag from Hoerner [24].

Substituting in the drag and lift relations from Figure 3
into Eqs. (33) and (35) for this case results in the

following for Fe.

-pu.

E= (A, [2sin’ a+ ,f0.083c0s' a—0.035cos' ay

+A, [1.2sin? acosar]- A,|0.083cosersinar—0.035cos" asinalp, |

Comparing the leading order terms of this result to the
lift and drag relations in Figure 5 verifies that the clutter
model recovers the independent lift and drag results of
HoernerTM for a cylinder at angle of attack.

Unidirectional Turbulent Flow in Porous Media
Two limiting cases are under consideration for

assessing the clutter model. The first is unidirectional
flow in periodic porous media as illustrated in Figure 6.



sin6 =é¢é,
3 (16)

y

cos Osiny = é
cos @ cosy = é, ¢ é,

In the porous media limit, the modified Darcy-
Forchheimer Law is again used except now the
permeability is expressed in terms of a second order

tensor resulting in the following expression for Sp
i

5, <8 [Ky + CoAT NE, (z,)

If the clutter is assumed to be orthotropic then the
permeability tensor can be expressed in terms of a
permeability components aligned with the cylinders,

K,, and components perpendicular to the cylinder,

yan

K,, and K,, inthe 81-412) coordinate system.

K, 0 0
K'@peu-8a)= Oo K, 0 (18)
0

0 Kk,

These coordinate directions represent the principal
directions of the permeability tensor. The superscript *
denotes that the tensor should be applied to a velocity

field expressed in terms of the &.é., »€, 5) coordinate

system. The empirically based models for K,, and
K_,, are same as used in Eq. (12),

K,=K,=9D' \ca-6"] a9

while the relation for K, comes from the use of

Carman-Kozey theory for flow along an array of
cylinders ’.

= spp ein tl-9))-3+ 40-9) 0-99] 2

In order to make use of the permeability tensor in a
calculation it must be expressed in terms of a Cartesian
coordinate system using the following transformation,

(21)

K=2'K

Is

where Tis defined in Eq. (15) and T
is a unitary matrix. Carrying out the matrix

To
=T since T

multiplication in Eq. (21) results in the following for
K,, | from Eq. (17),

K, Ky" Ky"

- “1 ~ -
K,, '= Ky Ky ' Ky, '
Ky! Ky! Ky

where,

Ky! -oyl +o ow) | So + Sey

Ky K

W028) oy yee. Ginoy

42

1 CosoF. Gino} (22)
Ky =
Ky Ky
- . ey  Ginay 1
K,,'=K,," =cosysiny| Gosoy
a Ky K, Ka Ka
1
Ky! =K,) = 6) —-—
au 13 =cosysinOcos. [x z|

1 |
= 8cos6| —-—
siny sin 6 cos E ]

Ky Ky

For anisotropic clutter class in the limit of bluff body
drag, empirical correlations for the coefficient of lift
and drag are available in the handbook by Blevins® for
long cylinders at an incline as illustrated in the sketch
below.

Flow direction

aaa “4
uy + U2

Coefficient of lift and drag for tangential force:
C, =0.083cosa@—0.035cos’ a

Coefficient of lift and drag for normal force:
C,, =1.2sin? @

Figure 3: Sketch of inclined cylinder and lift and
drag correlations from Blevins [23].

Note, that the relations for lift and drag in Fig. 3 are
only valid under the conditions that there is no cross


more value is the intrinsic average, (B) , defined as the

V,,, for

‘ee
which constitutive and thermodynamic properties are
well defined. The intrinsic average is defined as:

(B) = Blo and is the variable of interest to solve for

after phase averaging the transport equations for mass,
momentum and energy. The variable, @, is the void
fraction and is defined as the volume of gas divided by

local average of f over the gas phase volume,

the averaging volume, ie. 9=V,/V, =V,/A , :

Temporal__and__ Spatial__Derivatives _for__Phase
Averaging

Applying the spatial filtering function to the transport
equations requires expressing the phase averaged time
and spatial partial derivatives in terms of temporal and
spatial derivatives of phase averaged quantities. In the
following development, the volume of the gas phase is
assumed to change as a function of both time and space
to accommodate for the complexities that may be
included in future efforts (e.g. decomposing or burning
clutter). Simplifications are then imposed to limit the
scope to the focus of this work on turbulence modeling
on rigid, stationary clutter.

e de!

Relations for the temporal derivative are obtained by
taking the partial derivative of B .

(2)

Since the volume of the gas is assumed to be a function
of time, Leibnitz rule'* has to be applied to commute
the temporal derivative term inside the volume
averaging operator as follows:

oP _ G| ==: BE) yy
a yrw

; ot
oats oom (a,
f

(3)

+faG

Ag. 1)

A

In Eq. (3), the first term on the right hand size is simply
the definition of the phase averaged time rate in change

of B, ie. Of/dt. The second term is a surface

integral that accounts for material crossing the
boundary of the phase averaging volume as a function

of time where w, is the velocity of the surface (pointing
into the solid clutter, see Figure 1) of the gas phase

volume, V,, on the surface, A,. The gas phase

8
surface can be further broken into the contribution

associated with the clutter, A,, and rest of the surface
area associated with the phase averaging volume,

Aye: Since G goes to zero at the boundary on the

surface A,_.then the contribution of the surface

integral in those regions is equal to zero resulting in the
following relation for the phase averaged time
derivative.

_~
oP _9 x \p (et) aa’
Fa fe £ |B(x!)w, (x7) dA, (4)
a at gy | Ay ‘ ‘

Assuming the clutter surface is rigid and non-porous
then the normal velocity to the clutter surface is

identically zero (i.e. W, = 0) allowing for the temporal
derivative to commute with the filter operation resulting
in the simple relation.
—~
op _
ot

B 6)

Ye

Space derivatives with uniform filter width (Ze. A, =

Analogous to the development for time derivatives, a
spatial derivative relation is obtained by taking the
divergence of a phase-averaged vector quantity.

0 4 re) nt %—x ) oo,
_ 8
ar, B,(x,) 2, J AiG i, 6)

Since the volume averaging is independent of the
location, the partial derivative can be brought inside the
integration and applied directly to the filtering function

G.

0 4 A)
3 B= J BEG

const )


studies, time and phase-averaged transport equations of
momentum, turbulent kinetic energy and turbulent
dissipation rate are derived and constitutive models
suggested for accounting for unknown SGS correlations
and surface integral terms. These constitutive models
are applicable to porous flow for which it is assumed
that the characteristic size of the rigid porous matrix is
much smaller than the size of the system 7. This effort
extends these approaches to clutter environments for
which the length scale of the clutter is not necessarily
small relative to the size of the system.

Discussion begins with first the mathematical
formulation of the two-phase system based on the use
of spatial filtering concepts. The filtering of the
governing equations results in unknown correlations
that require modeling. Constitutive models are
formulated based on a linear combination of porous
media and bluff body lift and drag relations. Results
are then presented in support of the modeling
assumptions and predictions of turbulent flow in
porous media with comparisons to data from the
literature. Experimental measurements of turbulent
flow over a single cylinder are presented with
comparisons to detailed CFD predictions. The
information from the combined experimental and CFD
predictions provide further calibration of the clutter
model constants as part of future efforts in the clutter
model development. Lastly, conclusions are drawn and
future efforts are summarized.

Mathematical Formulation of Two-Phase
Mathematical Formulation of Two-Phase
System

The following mathematical description is limited to a
summary of the development of stationary clutter. The
section starts with an introduction to the two-phase (i.e.
spatial) averaging and time-averaging. The reasons for
reviewing these mathematical formalities are to
familiarize the reader with the concepts of spatial
filtering and, more importantly, highlight the
restrictions imposed by this approach with respect to
the implementation of the clutter model into a CFD
code. In the next section, these averaging concepts are
applied to conservation equations for momentum,
turbulent kinetic energy and dissipation transport. This
averaging results in unknown second order correlations
and surface integral terms that represent subgrid
physics that must be closed with a clutter model.
Requirements for the clutter model are provided based
on both physical and computational requirements. The
end result of this section is to provide a closed set of
phase and time averaged transport equations for
momentum, turbulent kinetic energy, and turbulent
kinetic energy dissipation rate.

Two-Phase Averaging

The formal averaging for two-phase media was first
introduced by Anderson and Jackson* through the use
of a local filtering function and later refined by Gray er
al.” and Gough et al." as Presented in text book form
by Kuo '!. Alternatively, Slattery!" offers a different
derivation based on a spatially dependent volume of
integration. The presentation here follows the procedure
of Kuo. Figure 1 illustrates a typical phase averaging

volume showing the total volume of interest, V;., the

volume of the solid clutter, V.., and the volume of the

gas, V,.
oe "
ro.
Er c=
ai MK Ve

Figure 1: Illustration of phase averaging volume.

Phase averaged properties are obtained by first defining
a spatial filtering function, G(x, - xVA,|, with the

normalization _ property: ac —x VA, jav’= i
v

For volume averaging using a cubic volume of length
L, =(V,)" on a side then A, =L, andG is
Any, ,
Th-x +4, /2)-H{y—« A, /2)
ici
where H is the Heaviside function. Convoluting G
with the a gas property of interest, Bh » yields a gas
phase average quantity,

Bls,)= J Aled

defined as: Gu

dv’ (1)

which physically represents a spatially averaged
Property over the volume, V,,. It should be emphasized
that the volume of integration, V,,_ . represents all of the
regions occupied by the gas and so does not depend on
the location, x,, where the averaging takes place. Of


ON THE DEVELOPMENT OF A SUBGRID SCALE MODEL FOR
CLUTTER

Paul E. DesJardin, James M. Nelsen and Louis A. Gritzo
Fire Science and Technology

Sandia National Laboratories”
Albuquerque, NM 87185

Terry A. Ghee
NAVAIR
Patuxent River, MD 20670

Abstract

The objective of this research is to develop a subgrid
scale model to account for the effects of small objects,
or “clutter”, that are difficult to resolve using a CFD
mesh. These objects may range in size from being
small relative to the grid (i.e. porous media li
scales, which are on the same order as the grid size (i.e.
a bluff body limit). Transport equations for
momentum, turbulent kinetic energy and turbulent
kinetic energy dissipation rate are derived based on
time and spatial filtering concepts, introducing
unknown correlations and surface integral terms that
require explicit subgrid scale (SGS) model closures.
The SGS models are formulated using a combination of
well-established constitutive relations borrowed from
the porous media literature and classical relations for
drag on bluff bodies. The modeling methodology is
exercised for two classes of problems. The first is flow
in a porous media for which the obstructions are
relatively dense. For this problem, predictions using
the SGS clutter model are compared to established
correlations taken from the porous media literature. The
second class of problems is a cylinder in cross-flow for
which both detailed experimental measurements and
detailed CFD calculations of velocity deficit and kinetic
energy profiles are available. Preliminary results
indicate that this modeling approach offers a promising
new approach to capture the macroscopic effects of
small-scale objects without using excessive mesh
resolution.

* Sandia is a multiprogram laboratory operated by Sandia
Corporation, a Lockheed Martin Company, for the United States
Department of Energy under Contract DE-AC04-94-AL85000.

‘This paper is declared a work of the U.S. Government and is
not subject to copyright protection in the United States.

Introduction

Capturing object features using a CFD grid that are
smaller than one one-hundredth the size of the
computational domain requires an excessively large
number of grid points and extremely small time steps
for numerical integration. Relevant applications
include large-scale problems such as forest fires down
to smaller scale phenomena such as heat transfer in a
porous media. An alternative approach to resolving
these features is to use a subgrid scale (SGS) model to
represent the macroscopic effects of these small
features using reasonably sized CFD grid cells. The
focus of this research is to develop engineering SGS
models of flow transport in cluttered environments.
This work is being carried out in a joint, iterative
computational/experimental approach for — which
detailed experimental and CFD predictions are used to
develop macroscopic models of clutter.

Subgrid clutter has two main physical effects on the
macroscopic flow field. The first is to provide a
momentum sink due to viscous and pressure drag forces
of gas flowing through clutter. The second is to either
increase or decrease the turbulent kinetic energy levels

depending on the local clutter size, 1, relative to the
characteristic length scale of turbulence, /,. In a
qualitative sense, if 1, < /, the turbulent kinetic energy
will decrease and if /. > /; the turbulent kinetic energy
will increase.

Previous investigations of turbulent flow in porous
media can be found in previous studies of Pedras et
al.!?3, Nakayama et al.*° and Antohe ef al.°. In these


One of the challenges in comparing the CFD
predictions to hot-wire data is that hot-wire only
measures the absolute value of velocity. This limitation
becomes especially problematic in regions of the flow
undergoing strong recirculation where the measured
time-averaged velocity using hot-wire can be much
larger than the actual velocity. Alternatively, one can
compare the average of the velocity squared for which
no ambiguity exists. Figure 9 show measurements
(shown with symbols) of time-averaged velocity
squared. This information is used in a momentum
balance calculation to find the coefficient of drag (Cp)
and is determined to be 0.97, which is 24% lower than
coefficient of drag of 1.28 determined with a turbulence
intensity level of 2% without the turbulence generator
in place. This result indicates that the effects of
upstream turbulence on the drag characteristics of
clutter are significant and should be included in future
refinement of the clutter model.

Also in Figure 9 are several detailed 2D and 3D CFD
predications using the CFD-ACE software package”
(shown with lines). The purpose of conducting the
detailed CFD calculations is to provide full flow field
details to construct phase-averaged quantities that
cannot be obtained experimentally. These details are
needed for calibration of the unknown clutter model

coefficient, C,, , in Eqs. (36) and (37) for the bluff
body drag limit.

As shown in Figure 9, good qualitative agreement is
obtained using 3D calculation while the 2D predictions
tend to under-predict the extent of the velocity deficit
profile for all of the turbulence models considered.
These differences may be attributed to the limitations
performing a 2D calculation when the flow has a strong
3D dimensional component of stream-wise vorticity
generation. The 3D prediction using the RNG k-e
model shows good quantitative agreement to the
experimental data just before and after the cylinder.
However, the agreement tends to degrade further
downstream. Similar results were obtained using a k-@
turbulence model. The current effort is to resolve these
differences so that the detailed CFD calculations can be
used to obtain phase-averaged stream-wise velocity and
turbulent kinetic energy profiles as a function of
downstream distance. This information will be
compared to predictions to macroscopic model and will
help determine the final modeling constants needed in
the bluff body limit of the clutter model.

nDeos wed

Vertical Traverse (mm)

vi (m/s?)
Figure 9: Comparisons of mean streamwise velocity
profiles at several down stream locations. Solid lines
are CFD predictions and symbols are experimental
measurements.

Conclusions

A preliminary clutter model has been formulated based
on the use of two-phase averaging concepts. The
implementation of phase averaging introduces
constraints in the requirements of the clutter model.
The most stringent of these constraints is that the model
must be able to represent a wide range of clutter length
scales ranging from a porous media limit to bluff body
drag. This requirement is satisfied using a linear
combination of constitutive models from the porous
media literature and relations for bluff body drag.
Based on this formulation, two simplified problems are
considered for calibration of the model constants. The
first is turbulent flow in a porous media for which a
preliminary version of the clutter model is able to
reproduce established trends taken from the literature.
The second problem is cylinder in a highly turbulent
cross flow for which there is little existing published
data. Results of experimental measurements show a
significant decrease in drag coefficient due to the highly
turbulent inflow that will be important in future

calibration of the model. Detailed numerical
predictions of the same problem show good qualitative
agreement to experimental measurements but

quantitative differences exist. Future work is to perform
this calibration so that the macroscopic model will be
able to reproduce velocity deficit and turbulent kinetic
energy production in the bluff body limit.

Acknowledgments

The authors would like to thank Dr. Dave Keyser of
Naval Air Systems Command for coordinating this
activity. This research is part of the Department of
Defense’s Next Generation Fire Suppression
Technology Program, funded by the DoD Strategic
Environmental Research and Development Program.