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

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

Conduction, convection and radiation are the three modes of heat transfer.

In previous chapters we have seen how problems involving multi-dimensional

conductive and convective heat transfer can be solved using CFD.

Conduction-only problems are simply dealt with by solving the heat conduc-

tion equation, i.e. diffusion equation (4.1). When ﬂuid ﬂow is involved the

resulting convective heat transfer problem is solved by tackling the enthalpy

equation (2.27) alongside the Navier–Stokes equations (2.32a–c) and the

continuity equation (2.4). The boundary conditions for the enthalpy equation

take care of heat transfer into and out of the computational domain across its

boundaries. The internal distribution of heat source and sink processes and

the transport of heat by means of diffusion and convection determine the

enthalpy distribution due to ﬂuid ﬂow.

The third mechanism of heat transfer – thermal radiation – is caused

by energy emission in the form of electromagnetic waves or streams of

photons. Radiative energy sources emit a broad-band spectral distribu-

tion with maximum energy content at a wavelength that is determined by

the source temperature. Most engineering systems emit thermal radiation at

infrared wavelengths (0.7–100 µm). The peak wavelength of sources at room

temperature will be around 10 µm, whereas an increasing fraction of the

radiation will be emitted at visible wavelengths when the source temperature

exceeds 1000–1500 K.

Engineering problems with significant radiation effects

Radiative heat transfer is often neglected in CFD calculations, because the

majority of engineering problems are dominated by high rates of convective

heat transfer. There are, however, several practically important categories

of problems where radiative heat transfer should be considered. Three pro-

minent examples are:

1 Manufacturing processes involving lasers or other high-energy beams:

these involve exposure of materials to radiation with very high energy

density. The radiative heat input will dominate these problems.

2 Combustion equipment: there is usually vigorous convection in

combustors, but the chemical reactions generate operating temperatures

that are sufﬁciently high for radiative heat ﬂuxes to be of similar order

of magnitude.

3 Naturally ventilated spaces in buildings: the temperatures are low so

radiative heat ﬂuxes are modest but comparable in size with convective

heat ﬂuxes since buoyancy-driven ﬂow velocities are often small.

Chapter thirteen Numerical calculation of

radiative heat transfer

Introduction13.1

ANIN_C13.qxd 29/12/2006 04:49PM Page 417

All these problems involve combined-mode or conjugate heat transfer: coupled

radiation and conductive heat transfer or, if ﬂuid ﬂow occurs, coupled radi-

ation, conductive and convective heat transfer. Laser-based material proces-

sing falls within the complex category of free boundary problems with phase

change. We do not discuss these further. In this section we focus on coupled

radiation/CFD problems in enclosures or ﬂow domains with ﬁxed bound-

aries. Radiation acts to redistribute energy through interactions within the

ﬂuid, through radiative exchanges between boundary surfaces, and through

interactions between the surfaces and the ﬂuid.

Definitions

We give a selection of basic facts and deﬁnitions of radiative heat transfer

needed to support our review of the main solution algorithms.

The emitted radiative heat ﬂux is a strong function of the temperature

of the substance. Materials at higher temperature will emit more radiation.

The rate of heat ﬂow per unit surface area emitted by a radiating surface is

called its emissive power E (units W/m

). For a so-called black body the

emissive power is related to its temperature by E

, where

= 5.67 × 10

−8

(W/m

) is Stefan–Boltzmann’s constant and T is absolute temperature (K).

The incident radiative heat ﬂux at a certain location varies as a

function of orientation of the receiver relative to the radiation source. The

intensity I (units W/m

.sr) is the rate of heat ﬂow received per unit area

perpendicular to the rays and per unit solid angle (steradians = sr), and is a

quantity with a magnitude that will vary with direction.

The black-body emissive power E

can also be expressed as a (non-

directional) black-body intensity using I

= E

. The black-body

intensity is used to compute the emitted intensity from surfaces and ﬂuids.

Surface properties

The emissive power of an ideal black surface is, of course, given by E

where T

is the surface temperature. Real surfaces usually emit less radiative

heat. The ratio of the heat ﬂux emitted by a real surface and a black surface

at the same temperature is called the surface emissivity

. Hence, the emis-

sive power of a real surface with surface emissivity

is given by E

εσ

The (non-directional) emitted intensity I

of a real surface is the product

of its surface emissivity and the black-body intensity: I

εσ

Figure 13.1 shows that surfaces interact with incident radiation in three

ways. The incident radiation can be (i) absorbed, (ii) reﬂected or (iii) trans-

mitted. The fraction absorbed is denoted by

, the fraction reﬂected by

and the fraction transmitted by

. Of course, the sum of the absorbed,

reﬂected and transmitted fractions equals unity:

= 1. In many

coupled CFD/radiation problems the surfaces will be (i) opaque, so

= 0,

and (ii) diffusely reﬂecting, i.e. the reﬂected radiation leaves the surface in

all directions irrespective of the incident angle of the incoming radiation.

Real surfaces may have different properties, e.g. non-zero transmissivity

or specular reﬂection (equal angle between surface normal and incident and

reﬂected radiation). All the surface properties are dependent on the type of

material, the surface roughness and the presence of surface contaminants,

as well as the temperature, and surface properties may also depend on the

direction and wavelength of the incident radiation.

418 CHAPTER 13 CALCULATION OF RADIATIVE HEAT TRANSFER

ANIN_C13.qxd 29/12/2006 04:49PM Page 418

13.1 INTRODUCTION 419

Fluid properties

We have noted before that air is transparent to the infrared radiation typically

associated with room temperature. The species present in clean atmospheric

air do not absorb readily at these wavelengths, and this ﬂuid medium will

not participate in the radiative heat exchanges. Products of combustion reac-

tions, on the other hand, contain large amounts of carbon dioxide and water

vapour, which are both strong absorbers/emitters in the infrared part of the

spectrum. Moreover, scattering may occur under certain conditions if a com-

bustion reaction produces particulate matter in the form of soot, or if the fuel

is in the form of solid particles. Thus, in combusting ﬂows incident radiation

is absorbed, transmitted and/or scattered and radiation may be emitted (see

Figure 13.2). In such cases the ﬂuid is termed a participating medium.

Figure 13.1 Radiation processes

at a surface

Figure 13.2 Schematic of

a pencil of ray in a medium

subjected to emission,

absorption, scattering and

transmission processes.

The strength of the interactions between a participating ﬂuid medium

and radiation can be measured in terms of its absorption coefﬁcient

and its scattering coefﬁcient

, both of which have units m

−1

. The sum of

these two properties is called the extinction coefﬁcient

. The

emitted intensity I

of a participating ﬂuid medium is the product of the

absorption coefﬁcient and the black-body intensity: I

κσ

ANIN_C13.qxd 29/12/2006 04:49PM Page 419

where T

is the ﬂuid temperature. The distribution of emitted intensity by

a point radiation source in a participating ﬂuid medium is uniform in all

directions, but the scattered intensity is generally not. The distribution of

the latter can described by the so-called scattering phase function Φ(s

, s),

which is deﬁned as the fraction of radiation incident upon the medium

along direction vector s

, which is subsequently scattered in direction s. The

absorption coefﬁcient, scattering coefﬁcient and the scattering phase function

are properties of the participating medium. The absorption coefﬁcient

will depend on the ﬂuid temperature and species concentrations and may

also depend on pressure. The scattering coefﬁcient

and the scattering

phase function Φ(s

, s) depend on the size, concentration, shape and

material characteristics of suspended particulate matter in the ﬂuid. All these

radiative properties generally have a complex wavelength dependence.

The nature of coupled CFD/radiation problems

Radiation can be considered as an electromagnetic wave phenomenon.

Therefore, the propagation speed of thermal radiation is the speed of light,

which, at 3 × 10

m/s, is at least 10

times as fast as the speed of sound of

common ﬂuids encountered in engineering problems. This large separation

of the velocity scales means that radiative heat exchanges are always in

quasi-steady state. This means that radiation propagates sufﬁciently fast to

effectively adjust itself immediately to variations in ﬂow conditions and/or

boundary conditions.

There is no direct coupling between radiation and the ﬂow ﬁeld, since

radiation or radiation properties of ﬂuids and boundaries do not depend

directly on the ﬂuid velocity. However, the ﬂow ﬁeld inﬂuences the spatial

distributions of temperature and species concentration. These determine the

intensity of radiation emitted by boundary surfaces and participating ﬂuids

as well as the radiation properties of a participating medium. This ensures

that there is strong indirect coupling between the ﬂow ﬁeld and radiation

environment, which is particularly important in combustion problems where

product species are at high temperature. We now review the consequences

for CFD of these coupling effects:

1 If the ﬂuid is absorbing/emitting and/or scattering there will be an

additional radiation source term in the energy equation

2 Radiation effects will cause changes in the boundary conditions of the

energy equation

Radiation source term in energy equation

In a steady ﬂow without radiation there is a balance within each ﬂuid control

volume between the sources and sinks of energy on the one hand and diffu-

sive and convective ﬂuxes across the boundaries on the other hand (see

Figure 13.3). Depending on the radiative properties of the ﬂuid inside the

control volume, it may absorb, emit and/or scatter radiation. When these

effects are signiﬁcant there will be additional heat ﬂuxes across the control

volume boundaries due to radiation. The net radiative heat ﬂux into or out

of the control volume will appear in the energy equation as a source or sink,

respectively.

Thus, the steady ﬂow energy balance for the control volume in the two

cases is as follows:

420 CHAPTER 13 CALCULATION OF RADIATIVE HEAT TRANSFER

ANIN_C13.qxd 29/12/2006 04:49PM Page 420

13.1 INTRODUCTION 421

Figure 13.3 Heat ﬂuxes and

sources for a ﬂuid control volume

without and with radiation

Without radiation: S

dV = q

conv

. n dA + q

diff

. n dA (13.1a)

Net energy Net convective Net diffusive

source in CV = ﬂux across + ﬂux across

boundaries boundaries

With radiation: S

dV + S

h,rad

= q

conv

. n dA + q

diff

. n dA + q

rad

. n dA (13.1b)

Non-radiative Net radiative Net Net Net radiative

energy source + energy source = convective + diffusive + ﬂux across

in CV in CV ﬂux ﬂux boundaries

h,rad

= q

rad

. ndA = (q

−

− q

)dA (13.1c)

where D

h,rad

is the net radiative source per unit volume for the energy equation.



∆V



∆V



ANIN_C13.qxd 29/12/2006 04:49PM Page 421

422 CHAPTER 13 CALCULATION OF RADIATIVE HEAT TRANSFER

Figure 13.4 Angular notation

for equation (13.2), incoming

radiative ﬂux

In (13.1c) the net radiative heat ﬂux q

rad

. n at a point on the control

volume boundary has been written as the difference between the incident

radiative heat ﬂux q

−

and the outgoing heat ﬂux q

. The integral of the

difference q

−

− q

over the bounding surface of the control volume gives the

net radiative heat ﬂux into the control volume to obtain the net radiative heat

source for the energy equation. To obtain the incident heat ﬂux q

−

at a point

it is necessary to integrate the intensity, which varies with direction, and hence

we must consider all possible incoming ray directions. It can be shown (Modest,

2003) that this can be written as an integration over a unit hemisphere

(i.e. over a solid angle of 2

steradians) surrounding the point just outside the

control volume boundary (indicated as a semi-circle in Figure 13.3). So,

−

= I

−

(s) s.ndΩ= I

−

(

) cos(

) sin(

) d

(13.2)

The dot product s.nin the ﬁrst integral ensures that we consider the

component of the incident radiative heat ﬂux vector in the direction of

the outward surface normal n. We note that the swap of argument of the

incident intensity I

−

from ray direction vector s to angular co-ordinates

(

) is simply a matter of notation. The co-ordinate system is illustrated

in Figure 13.4.



2π

To obtain the outgoing radiative heat ﬂux q

we integrate over all pos-

sible ray directions pointing outwards at the point. This corresponds to an

integration over a unit hemisphere (indicated as a semi-circle in Figure 13.3)

surrounding the point just inside the control volume boundary. So,

= I

(s) s.ndΩ= I

(

) cos(

) sin(

) d

(13.3)

(s) and I

(

) represent the outgoing intensity of radiation in the direc-

tion of vector s corresponding to angular direction (

). The limits in the



2π



ANIN_C13.qxd 29/12/2006 04:49PM Page 422

13.1 INTRODUCTION 423

-integration have been adjusted to indicate the unit hemisphere just inside

the control volume with outward-pointing radiation rays.

The effect of radiation on boundary conditions

Without radiation there is a balance between the conductive and/or convec-

tive heat ﬂuxes at the interface that constitutes the imaginary boundary

between two materials. Radiation will cause an additional heat ﬂux towards

the interface due to incident radiation and an extra outgoing heat ﬂux asso-

ciated with emission of radiation.

The overall heat balance for the interface between the ﬂuid and the

boundary states that the heat ﬂux out of the interface must be equal to the

heat ﬂux into the interface. Thus, as illustrated in Figure 13.5,

Without radiation: q

= q

EXT

(13.4a)

With radiation: q

+ q

= q

EXT

+ q

−

⇒ q

= q

EXT

+ (q

−

− q

) (13.4b)

Figure 13.5 Boundary

conditions with and without

radiation

Both the incident radiative heat ﬂux q

−

and the outgoing radiative

heat ﬂux q

, which is the sum of emitted and reﬂected heat ﬂuxes, depend

on the intensity of the radiation. To obtain the incident heat ﬂux at a point

on the surface we use equation (13.2) with the surface normal n pointing into

the ﬂuid region.

The outgoing radiative heat ﬂux is given by

+ (1−

−

εσ

+ (1 −

−

(13.5)

The ﬁrst term on the right hand side is the emitted heat ﬂux and the second

term is the reﬂected heat ﬂux. The latter is just equal to the fraction (1 −

)

times the incident heat ﬂux q

−

in the case of diffuse reﬂection. The corres-

ponding non-directional intensity of the outgoing radiation is

= q

= [

εσ

+ (1 −

−

(13.6)

If we use a ﬁxed temperature boundary condition, we need to bear in mind

that the net radiative heat ﬂux q

−

− q

will vary across the boundary surfaces.

Radiation problems always involve combined-mode heat transfer, and surface

locations experiencing positive values of q

−

− q

are likely to heat up com-

pared with areas with zero or negative q

−

− q

. Similarly, if we use a heat ﬂux

boundary condition, a CFD computation requires us to specify the heat ﬂux

into the ﬂuid q

, but it is generally more convenient to specify the external

heat ﬂux q

EXT

. As we can see in the interface heat balance, q

= q

EXT

for a

ANIN_C13.qxd 29/12/2006 04:49PM Page 423

424 CHAPTER 13 CALCULATION OF RADIATIVE HEAT TRANSFER

ﬂow without radiation. In an environment with radiation, however, the heat

ﬂux going into the ﬂuid q

at a boundary point will be different from the

external heat ﬂux q

EXT

by an amount equal to the local net radiative heat

ﬂux q

−

− q

. These remarks serve to emphasise that completely accurate

speciﬁcation of surface boundary conditions is difﬁcult without a priori

knowledge of the radiation environment.

The incident intensity I

−

and the outgoing intensity I

, which are needed in

equations (13.2) and (13.3) for the radiative heat ﬂuxes, need to be computed

from the governing equation of radiative transfer. The general relationship

that governs the changes in intensity at a point along a radiation ray due to

emission, absorption and scattering in a ﬂuid medium is as follows (Modest,

2003):

(r) −

I(r, s) −

I(r, s) + I

−

) Φ(s

, s) dΩ

(13.7)

Rate of change of Emitted Absorbed Out-scattered In-scattered

intensity per unit = intensity − intensity − intensity + intensity

path length

where I(r, s) is the radiation intensity at a given location indicated by

position vector r, in the direction s within a small pencil of rays, travelling

through a participating medium (see Figure 13.2 and Modest, 2003). The

in-scattering integral on the right hand side of equation (13.7) accounts for

the effect of intensity I

−

) incident upon the point at r from all possible

directions s

. This involves integration over a unit sphere (i.e. a solid angle of

steradians) surrounding the point in the medium.

The absorbed and out-scattered intensity are often drawn together by

deﬁning the extinction coefﬁcient

as the sum of the absorption

coefﬁcient and the scattering coefﬁcient. Thus,

(r) −

I(r, s) + I

−

) Φ(s

, s)dΩ

(13.8)

This is the radiative transfer equation (RTE), which governs radiative heat

transfer along a ray path indicated by the direction vector s. If the incident

intensity ﬁeld for the in-scattering integral (last term on the right hand side

of the RTE) is somehow known, the RTE is a ﬁrst-order ordinary differen-

tial equation.

Many computational algorithms for radiative heat transfer proceed to

simplify the radiative transfer equation by introducing non-dimensional

optical co-ordinates:

= (

)ds′=

ds′ (13.9)



dI(r, s)



dI(r, s)

Governing

equations of

radiative heat transfer

13.2

ANIN_C13.qxd 29/12/2006 04:49PM Page 424

13.2 GOVERNING EQUATIONS OF RADIATIVE HEAT TRANSFER 425

and the single scattering albedo

deﬁned as

== (13.10)

The radiative transfer equation can be re-written as

=−I(

, s) + (1 −

(

) + I

−

) Φ(s

, s)dΩ

(13.11)

Finally, we introduce the source function, which is deﬁned as follows:

, s) = (1−

(

) + I

−

) Φ(s

, s)dΩ

(13.12)

This enables us to write the radiative transfer equation in a form with a very

simple appearance:

+ I(

, s) = S(

, s) (13.13)

Boundary conditions

Equation (13.8) or (13.13) is solved by integration along ray paths starting

at the boundaries of the computational domain, which requires initial con-

ditions on these surfaces. For a point on a diffusely emitting and reﬂecting

opaque surface located at position vector r

we use expression (13.6) to

obtain the initial condition of a ray propagating along direction vector s:

I(r

, s) =

) + I

−

, s

)n.s

dΩ

(13.14)

This is the most widely used boundary condition in combustion-related

radiative transfer problems. Expressions for other possible boundary condi-

tions can be found in Modest (2003).

Incident intensity integrals

The RTE (13.8) and boundary condition (13.14) contain integrals of intensity

over a sphere and hemisphere, respectively. To evaluate these it is necessary

to know the integrands beforehand as a function of the direction of the incid-

ent radiation. Since the incident intensity integrals are unknown and must

be solved alongside the intensity along an RTE ray path, equation (13.8) is

not an ordinary differential equation but actually an integro-differential

equation. This type of equation is clearly very different in nature from the

transport equations we have seen before. This requires solution of RTEs

along all the relevant ray paths, so the transport of heat by radiation is always

three-dimensional. Moreover, radiative transfer is heavily dependent on the



(1 −

)

dI(

, s)



dI(

, s)

ANIN_C13.qxd 29/12/2006 04:49PM Page 425

geometry. Therefore, the calculation of radiative heat transfer requires an

approach that is quite different from the ﬁnite volume method we have dis-

cussed elsewhere in this book. We will discuss the most prominent general-

purpose solution algorithms in the next section.

Exact analytical solutions of the RTE are not available except in a small

number of idealised cases. In most practical problems further simpliﬁcation

of the transfer equation is not possible. In particular, the boundary conditions

are dictated by the problem geometry, and full three-dimensional effects and

all angular directions must be considered in solving the transfer equation.

Before we discuss the details of solution methods we examine the RTE and

boundary conditions to gain insight into some of the general features of

likely solution procedures. First, we note that the incident intensity integrals

in (13.2) and (13.8) are unknown at the start of a calculation, so an iterative

approach is required. Assumed values of the surface intensities are initially

used and, after solving the RTE for a sufﬁciently large number of ray paths,

the incident intensity integrals can be evaluated. This enables us to make

improved estimates of the boundary conditions and the in-scattering inte-

grals to carry out another round of RTE solutions. This process is iterated

until there is no further change in the solution.

Further complications arise in combusting systems, because the equation

of transfer contains terms involving the unknown temperature ﬁeld and

radiative properties, which are dependent on temperature and composition.

Combustion effects dominate the temperature ﬁeld, and the radiative prop-

erties of the medium should be determined from the concentration of species

of combustion. The interaction between radiative heat transfer and ﬂow

enters into the ﬂuid ﬂow calculation through the source term (see equation

(12.80) in Chapter 12) and wall heat transfer effects. Due to this coupling of

temperature and radiative properties of the medium the solution of the RTE

in combustion problems requires an outer iteration loop, which is executed

until the solution satisﬁes all the ﬂuid ﬂow equations and the RTE.

Finally, the radiation properties of combustion products are dependent

on the wavelength of the radiation and, for accurate simulations, it is neces-

sary to perform spectrally resolved radiation calculations or to model the

effect of this wavelength dependence.

The calculation of radiation is numerically challenging and resource

intensive, since algorithms have to compute radiation intensity as a function

of position (x, y, z) in the computational domain, angular direction (

)

and, in the most accurate calculations, radiation wavelength

. In our brief

review of techniques for radiation calculations we focus our attention on

computational methods for the dependence of radiation intensity on position

and direction and illustrate their application by means of three examples of

increasing complexity. We ﬁnish the chapter by making some brief remarks

on issues relating to wavelength-dependent radiation calculations.

Over the years many methods have been developed for the solution

of radiative heat transfer. These include various analytical approximation

techniques and a suite of numerical methods. Some early methods have

now largely been abandoned because of their limited applicability to general

situations. Other methods, such as the zone method, P–N methods, ﬂux

methods and ﬁnite element methods are not discussed here. Details can be

found in the reviews and texts by Saroﬁm (1986), Viskanta and Mengüç

426 CHAPTER 13 CALCULATION OF RADIATIVE HEAT TRANSFER

Solution methods13.3

ANIN_C13.qxd 29/12/2006 04:49PM Page 426