Versteeg H., Malalasekra W. An Introduction to Computational Fluid Dynamics: The Finite Volume Method

Подождите немного. Документ загружается.

APPENDIX E 457

div(

V)dV =−

= 0 (E.3)

To develop the one-dimensional mass conservation equation we consider a

control volume with sufﬁciently small length ∆x so that we can use a one-

term Taylor-series expansion to write expressions for velocity and cross-

sectional area:

= u and u

= u +∆x (E.4a)

= A and A

= A +∆x (E.4b)

Substitution of equations (E.4a) and (E.4b) into (E.3) yields

−

uA +

u + ∆xA+∆x = 0 (E.5)

Multiplying out the brackets and neglecting terms of order ∆x

, we obtain

u ∆x +

A ∆x = 0 (E.6)

Dividing the equation by the control volume length ∆x and combining the

two terms on the left-hand side yields

= 0 (E.7)

The results (E.3) and (E.7) are familiar from introductory ﬂuid mechanics

texts; the former is the usual continuity equation for a macroscopic control

volume and the latter is the differential form.

x-momentum conservation:

+ div(

uV) =− +div(

grad(u)) + S

(2.37a)

For steady, incompressible, frictionless ﬂow this reduces to:

div(

uV) =− (E.8)

Once more we integrate both sides of the equation over the control volume

and apply Gauss’s theorem to obtain integrals over the bounding surfaces.

For the left-hand side we ﬁnd:

div(

uV)dV =

uV . ndA +

uV . ndA

=−

A +

u +∆x

A +∆x (E.9)



∆A



∆V

∂

(

∂

uA)



∆V

ANIN_Z05.qxd 29/12/2006 10:08 AM Page 457

458 APPENDIX E

Using (E.5) this can be rewritten as follows:

−

A +

u +∆x

A +∆x

uA −u + u +∆x =

uA ∆x (E.10)

The control volume integral of the right-hand side of equation (E.8) and

subsequent application of Gauss’s divergence theorem (2.41) yields:

− dV =− pi . ndA − pi . ndA − pi . ndA

= pA − p +∆xA+∆x − p + k ∆x −∆x (E.11)

The factor k(0 < k < 1) in the last term accounts for the fact that the

pressure on side wall area ∆A

takes a value between p on inlet area A

and

p + (dp/dx)∆x on outlet area A

. Moreover, the x-component of the force

due to this pressure applies in the negative x-direction and has a magnitude

that is equal to the product of this intermediate pressure p + k(dp/dx)∆x

and the projection of side wall area ∆A

into the x-direction, which is

− A

=−(dA/dx)∆x. Multiplying out the brackets of (E.11) and neglect-

ing terms of order ∆x

, we obtain

pA − p +∆xA+∆x + p + k ∆x ∆x

=−A ∆x (E.12)

We place equation (E.10) on the left-hand side of the integrated x-momentum

equation and equation (E.12) on the right-hand side:

uA ∆x =−A ∆x (E.13)

Dividing the result by ∆x yields the one-dimensional momentum equation

for frictionless ﬂow:

uA =−A (E.14)



∆A



∂



∆V

∂

ANIN_Z05.qxd 29/12/2006 10:08 AM Page 458

As shown in Figures 11.15 and 16 (Chapter 11), using co-ordinates

along

the line joining P and A, and

along the face of the control volume (i.e.

along line joining points a and b), the term grad

can be written with respect

to x, y co-ordinates or n,

co-ordinates as

grad

= i + j = n + e

(F.1)

where n and e

are unit vectors along normal and tangential directions:

grad

. n A

= i + j . nA

= . − . A

(F.2)

Using the chain rule we have

=+ (F.3)

=+ (F.4)

It is usual practice to write

∂φ

∂ξ

∂φ

∂

x =

∂φ

∂η

etc.

Equations (F.3) and (F.4) can be rewritten as

(F.5)

(F.6)

Equations (F.5) and (F.6) can be solved to obtain

and

. The solution is

= (

−

) =− (F.7)

= (−

) =− + (F.8)

where J is the Jacobian given by

J = (x

− x

) =− (F.9)

∂

∂ξ

∂

∂η

∂

∂η

∂

∂ξ

∂

∂ξ

∂φ

∂η

∂

∂η

∂φ

∂ξ

∂

∂ξ

∂φ

∂η

∂

∂η

∂φ

∂ξ

∂

∂η

∂φ

∂

∂η

∂φ

∂

∂φ

∂η

∂

∂ξ

∂φ

∂

∂ξ

∂φ

∂

∂φ

∂ξ

∆x

∂φ

∂

∆y

∂φ

∂

∂φ

∂

∂φ

∂

∂φ

∂η

∂φ

∂

∂φ

∂

∂φ

∂

Appendix F Alternative derivation for the

term (n . grad

φφ

) in Chapter 11

ANIN_Z06.qxd 29/12/2006 10:09 AM Page 459

In equations (F.7) and (F.8)

∂φ

∂ξ

is the gradient of

along the line joining

P and A, and

∂φ

∂η

is the gradient of

along the line joining a and b (the

face). Other terms like

∂

∂ξ

∂

∂ξ

∂

∂η

and

∂

∂η

are simple geo-

metrical quantities. Expressions for those are given below.

By substituting (F.7) and (F.8) into (F.2) we have

grad

. n A

= . − . A

=−∆y −− + ∆x

=∆y +∆x −∆y +∆x (F.10)

The two ﬂux gradients involved in this expression may be calculated from

= and =

where ∆

= d

is the distance between points A and P and ∆

= d

is the

distance between vertices a and b (or in this case equal to A

). The ﬁrst term

involving

∂φ

∂ξ

is called the direct gradient and the second term involving

∂φ

∂η

is called the cross-diffusion term.

The other geometrical quantities in equation (F.10) may be obtained from

Using Mathur and Murthy’s (1997; see the references in Chapter 11)

approach, expression (F.9) can be written using normal unit vector n and

two other unit vectors e

and e

in the directions of

and

respectively.

Recall that

n = i − j = i − j (F.11)

The unit vector along the

-direction may be deﬁned as

= i + j (F.12)

The unit vector along the

-direction may be deﬁned as

= i + j (F.13)

− y

∆

− x

∆

− y

∆

− x

∆

− x

∆

− y

∆

∆x

∆y

− y

∆

∂

∂η

− x

∆

∂

∂η

− y

∆

∂

∂ξ

− x

∆

∂

∂ξ

−

∆

∂φ

∂η

−

∆

∂φ

∂ξ

∂φ

∂η

∂

∂ξ

∂

∂ξ

∂φ

∂ξ

∂

∂η

∂

∂η

∂

∂ξ

∂φ

∂η

∂

∂η

∂φ

∂ξ

∂

∂ξ

∂φ

∂η

∂

∂η

∂φ

∂ξ

∆x

∂φ

∂

∆y

∂φ

∂

460 APPENDIX F

ANIN_Z06.qxd 29/12/2006 10:09 AM Page 460

APPENDIX F 461

The Jacobian in the expanded form is

J =−

=−

J = n . e

(F.14)

Using the above simpliﬁcation for J the ﬁrst term of equation (F.10) may be

written as

∆y +∆x = ( y

− y

)

+ (x

− x

)

= [n . nA

]

= A

(F.15)

Similarly the cross-diffusion term of equation (F.10) may be written as

−∆y +∆x =− ( y

− y

)

+ (x

− x

)

=− [e

. e

∆

]

=− (F.16)

The complete diffusion contribution in vector form is

n . grad

=− (F.17)

∂φ

∂η

. e

n . e

∂φ

∂ξ

n . nA

n . e

∂φ

∂η

. e

n . e

∂φ

∂η

n . e

∂φ

∂η

− x

)

∆

( y

− y

)

∆

n . e

∂φ

∂ξ

∂

∂ξ

∂

∂ξ

∂φ

∂ξ

n . n

n . e

∂φ

∂ξ

n . e

∂φ

∂ξ

− x

)

∆

( y

− y

)

∆

n . e

∂φ

∂ξ

∂

∂η

∂

∂η

− y

)

∆

− x

)

∆

− y

)

∆

− x

)

∆

∂

∂ξ

∂

∂η

∂

∂η

∂

∂ξ

ANIN_Z06.qxd 29/12/2006 10:09 AM Page 461

G.1 Some applications

In this appendix we give three examples of simple CFD applications. An

interested reader with access to commercial or any other CFD software

could repeat these exercises to gain some hands-on experience. The bound-

ary conditions and problem speciﬁcation are brieﬂy described and specimen

results are presented to illustrate how CFD and modelling of combustion can

be applied to practical situations. The examples are presented without ﬁner

details of the calculations, which can be found in the cited references.

G.2 Flow in a sudden pipe contraction

The problem considered

This problem was selected to illustrate the application of CFD to a bench-

mark problem with a set of well-documented data for the comparison of pre-

dictions with experiments. The problem considered here is laminar pipe ﬂow

in a sudden contraction shown in Figure G.1. Durst and Loy (1985) have

provided the experimental data for a range of Reynolds numbers. The ﬂow

with a Reynolds number (Re =

UD/

) of 372 was considered for CFD

modelling, where U is the average velocity in the pipe with diameter D.

Appendix G Some examples

Figure G.1 Experimental

conﬁguration of Durst and

Loy (1985)

CFD simulation

The geometry was modelled with a two-dimensional axisymmetric grid of

100 × 60. The velocity proﬁle for fully developed laminar ﬂow was imposed

at the inlet, and the no-slip condition was applied at wall boundaries. At the

exit plane, all derivatives in the axial direction were set to zero. The CFD

calculation was carried out using the SIMPLER algorithm and the hybrid

differencing scheme.

ANIN_Z07.qxd 29/12/2006 04:50PM Page 462

APPENDIX G 463

Specimen results

Since the ﬂow is laminar, the governing equations are exact (i.e. no turbu-

lence modelling involved here). The predicted streamlines of the ﬂow are

shown in Figure G.2. The velocity proﬁles are shown in Figure G.3 for six

different cross-sections of the domain: three before, and three after, the con-

traction. The experimental data of Durst and Loy (1985) are also included

for comparison. It can be seen that the predictions agree well with the experi-

mental measurements. Further grid reﬁnement did not cause signiﬁcant

changes in the predictions, and therefore these results can be considered to

be grid independent. It should be noted that comparisons for locations other

than those shown in the ﬁgure and for other Reynolds numbers also agree

well with the experimental data. This simple example shows the capability of

CFD to predict practical ﬂow situations, with a good degree of accuracy.

Figure G.2 Predicted streamline

pattern

G.3 Modelling of a fire in a test room

The problem considered

In contrast to the previous benchmark problem we now study a case at the

other end of the spectrum of complexity. We compare CFD calculations

with experimental ﬁre tests carried out by the Lawrence Livermore National

Laboratory (LLNL) in the test room shown in Figure G.4. The details of the

experiments have been reported in Alvarez et al. (1984). The ﬁre was at the

centre of the ﬂoor and clean air was introduced along the ﬂoor of the test cell,

which is approximated in the model by a 0.12 m high and 2 m long slot for

air entry, located 0.1 m above the ﬂoor. The ﬁre sources in the experiments

were a burner, a spray and a pool of fuel in a tray. The products of combus-

tion were extracted near the top of the cell using an axial ﬂow fan through a

rectangular 0.65 m square duct placed 3.6 m above the ﬂoor, as shown in

Figure G.4. A total of 27 tests were reported by Alvarez et al. (1984), and the

one designated MOD08 has been selected for CFD modelling here. In this

test, a spray of isopropyl alcohol from an opposed-jet nozzle located at the

centre of the pan was used, and the fuel evaporated quickly to burn in a way

similar to a natural pool ﬁre. The fuel injection rate was 13.1 g/s with a

total heat release rate of 400 kW. These data were used to specify burner

conditions at the ﬁre source. The measured extraction rate, 400 l/s in the

steady state, was used to specify the outﬂow. The mass ﬂow rate of air into

the domain and the inlet and outlet velocities are calculated as part of the

solution. The walls, the ﬂoor and the ceiling of the compartment were of

0.1 m thick refractory. The estimated thermal conductivity, density and

speciﬁc heat were, respectively, 0.39 W/m.K, 1400 kg/m

and 1 kJ/kg.K

for the walls and 0.63 W/m.K, 1920 kg/m

and 1 kJ/kg.K for the ceiling

and the ﬂoor. The walls were assumed to be perfectly black for radiation

calculations.

ANIN_Z07.qxd 29/12/2006 04:50PM Page 463

464 APPENDIX G

Figure G.3 Comparison of predictions and experimental results at six different locations

ANIN_Z07.qxd 29/12/2006 04:50PM Page 464

APPENDIX G 465

CFD simulation

The simulation of the aerodynamics and combustion was carried out using

a three-dimensional CFD procedure based on the SIMPLE algorithm and

the hybrid differencing scheme for discretisation. Turbulence was modelled

with the k–

turbulence model with buoyancy terms, and combustion

modelling assumed fast chemistry (SCRS). The discrete transfer model

of thermal radiation (Lockwood and Shah, 1981) was used to calculate

radiative heat transfer. The wall temperatures were obtained from a one-

dimensional wall heat transfer model. A numerical grid of 14 × 13 × 12,

although not very ﬁne, was considered adequate to predict the overall prop-

erties of the ﬁre. Further details of the model can be found in Malalasekera

(1988) and Lockwood and Malalasekera (1988). Some specimen results are

presented below.

Specimen results

Figure G.5 shows the predicted steady state ﬂow pattern in the Y–Z plane at

X = 3.25 m. The buoyancy-generated ﬂow is clearly reproduced by the simu-

lation, which also shows the entrainment induced by the strong buoyancy

effects. The predicted temperature distribution in the Y–Z plane at X = 3.00 m

(Figure G.6) shows the hot gases around the central ﬂame and the forma-

tion of a hot layer at ceiling level. The ﬂame structure and tilt due to induced

air ﬂow are also clearly visible. Figure G.7 compares the room temperature

predictions with the experimental data of Alvarez et al. (1984). The experi-

mental temperatures were recorded using two thermocouple rakes (TR1,

east rake; and TR2, west rake) with 15 thermocouples each placed 1.5 m on

either side of the ﬁre and located in the central plane as shown in Figure G.4.

Figure G.4 Schematic diagram

of the Lawrence Livermore

National Laboratory (LLNL)

ﬁre test cell

ANIN_Z07.qxd 29/12/2006 04:50PM Page 465

466 APPENDIX G

Figure G.5 Predicted ﬁow inside

the compartment: velocity vector

plots in the Z–Y plane at

X = 3.25 m

Figure G.6 Predicted

temperature (K) ﬁeld in the

Y–Z plane at X = 3.00 m

ANIN_Z07.qxd 29/12/2006 04:50PM Page 466