3
GENERALIZED GOVERNING EQUATIONS IN MULTIPHASE SYSTEMS: LOCAL INSTANCE FORMULATIONS
3.1 Introduction
A multiphase system, which is distinguished from a singlephase system by the presence of one or more interfaces separating the phases, can be considered as a field that is divided into singlephase regions by those interfaces, or moving boundaries, between phases. A straightforward method of describing a multiphase system is to use the standard local instance differential equations for each singlephase subregion, with appropriate jump conditions to match the solution of these differential equations at the interfaces. This method requires tracking the location of interfaces, and for this reason it is referred to as the interface tracking method. Since the locations of the interfaces are unknown a priori, the governing equations are solved separately for each phase; then iterations are performed at the anticipated location of the phase interface until convergence of certain interfacial balances is obtained. These interfacial balances must fulfill the requirement of mass, momentum, and energy conservation across the interface. With this approach, the position of the phase interface may be clearly identified and used to separate the solutions of the governing equations on both sides of the interface. This method is most appropriate for dealing with separatedphase cases (see Table 1.10), such as melting, solidification, film evaporation, film condensation, and single bubble (or droplet) growth or collapse. The local instance formulation of the generalized governing equations will be presented in this chapter. For the multiphase systems that have extremely complicated and deformable interfaces, it is not possible to solve the local instance fluid flow. Various averaging techniques can be applied to obtain the mean values of properties and eliminate the need to explicitly track interfaces and/or the local instance fluctuations of properties (see Chapter 4).

Chapter 3 Generalized Governing Equations: Local Instance Formulations 177
Two common techniques are used for describing fluid flow: Lagrangian and Eulerian. The Lagrangian approach requires that the properties of a particular element of the fluid particles be tracked as it traverses the flow. This approach is similar to what we used in particle and rigidbody dynamics. The location of this fluid element is described by its coordinates (x, y, z), which are functions of time. The fluid element can be identified by tracking it from its initial location (x0, y0, z0) at time t = 0, and the velocity of this element at an arbitrary time t is expressed as V = V ( x0 , y0 , z0 , t ) . In order to describe a fluid flow using the Lagrangian approach, the sensors that monitor fluid properties would have to move at the same velocity as the fluid element; this is an impractical and often impossible requirement to meet, especially for such complex cases as threedimensional transient flow. Therefore, the Lagrangian approach is rarely used in description of fluid flow. The Eulerian approach, on the other hand, observes the flow properties from a fixed location relative to a reference frame, which can be stationary or more generally moves at its own velocity. The Eulerian approach gives the values of the fluid variable at a given point (x, y, z) at a given time t. For example, the velocity can be expressed as V = V ( x, y , z , t ) , where x, y, and z are independent of t. The Eulerian approach requires that the fluid properties be measured at spatial locations that are fixed relative to the reference frame in the fluid field, so the sensors are not required to move with the individual particles. Since the Eulerian approach is consistent with conventional experimental observation techniques, it is widely used in fluid mechanics and will also be adopted in this textbook. The conservation laws for multiphase flow and heat transfer can be expressed for a fixedmass or for a control volume. While the former is a fixed collection of particles with constant mass (fixedmass), the latter is a defined region relative to the reference frame in space (fixedvolume). From a thermodynamic point of view, the fixedmass and the control volume can be considered as closed and open systems, respectively. Writing the governing equations for a fixedmass requires tracking the motion of the particles, i.e., the Lagrangian approach. As we have already mentioned, although the Lagrangian description is applicable to some fluid mechanics problems, it is not a very practical way to describe multiphase systems, with a few exceptions such as liquid droplet tracking in a thermal plasma. The governing equations written for a control volume, on the other hand, express the relationship between the change of properties inside the control volume and the property of the flow into or out of the control volume. The governing equations obtained by writing the conservation laws for a control volume are consistent with the Eulerian approach. However, all of the fundamental laws of mechanics, including conservation of mass, momentum, and energy, are formulated for a collection of particles with fixed identity; that is to say, they are Lagrangian in nature. Therefore, it is necessary to apply the fundamental laws to the fixedmass first, and then to convert them into expressions for the control volume.
178 Transport Phenomena in Multiphase Systems
Section 3.2 presents local instance macroscopic (integral) formulations of the governing equations of the multiphase system equations. The macroscopic (integral) formulation is obtained by performing mass, momentum, energy, entropy and species mass balances over a control volume that includes different phases as well as interfaces. The microscopic (differential) formulation of the governing equations of the multiphase systems is presented in Section 3.3. The local instance microscopic (differential) formulations are obtained by simplifying the integral formulations for control volumes with only one phase. The local instance differential equations must be supplemented by the jump conditions at the interface. The classification of PDEs and boundary conditions, as well as a rarefied vapor selfdiffusion model, are also discussed. Section 3.3 is closed by discussion of a rarefied vapor selfdiffusion model and application of the differential formulations to combustion.
3.2 Macroscopic (Integral) Formulation
A fixedmass system describes an amount of matter that can move, flow and interact with the surroundings, but the control volume approach depicts a region or volume of interest in a flow field, which is not unique and depends on the user. So conservation laws for a fixedmass system need to be transformed to apply to a control volume. A mathematical relation that allows one to mathematically and physically link the conservation laws for a control volume with that of a fixedmass system will be derived. Figure 3.1 shows the flow field under consideration. At time t, the control volume shown by the solid line coincides with a singlephase fixedmass system depicted by the dashed line. At time t+dt, a portion of the fixedmass system moves outside of the boundaries of the control volume. It can be seen from Fig. 3.1 that region I is occupied by the system at time t only, region II is common to the system at both t and t+dt, and region III is occupied by the system at t+dt
Figure 3.1 Relation between a fixed mass system and a control volume.
Chapter 3 Generalized Governing Equations: Local Instance Formulations 179
only. For a system with a fixedmass, the change of the general property Φk, which has a specific value (per unit mass) φk, can be written as (Welty, 1978): Φ k t + dt −Φ k t d Φk (3.1) = dt system dt where k denotes the kth phase in the multiphase system. Considering that the control mass occupies regions I and II at time t, one obtains Φ t = Φ I + Φ II . At time t + dt , the fixedmass system occupies both regions II and III, i.e.,
Φ t + dt = Φ II + Φ III . Therefore, eq. (3.1) can be
t + dt
rewritten as
d Φk dt =
system
Φ k II
−Φ k II
t
dt
+
Φ k III dt
t +Δt
−
Φk I dt
t
(3.2)
The first term on the righthand side of eq. (3.2) may be written as Φ k II t + dt −Φ k II t d Φ k II = dt dt and represents the rate of change of property Φk within the control volume, because region II becomes coincident with the control volume as dt → 0 . For the general case of variable Φk within a control volume, it is appropriate to write the time derivative of Φk as d Φ k II § d Φ k · ∂ =¨ ¸ = ∂t ³V ρ k φk dV dt © dt ¹CV where V is the volume of the control volume. The second and third terms on the righthand side of eq. (3.2), respectively, represent the property Φk leaving and entering the control volume due to mass flow across its boundary. If the absolute velocity for the kth phase is Vk, and the reference frame moves with a constant velocity Vref, the relative velocity of the kth phase is Vk , rel = Vk − Vref . For the control volume’s entire surface area, A, the rate of movement of property Φk due to mass flow may be written as (see Fig. 3.2) Φ k III t +Δt Φ k I t − = ³ ρ k (Vk ,rel ⋅ n k )φk dA A dt dt where nk is the normal direction of the control volume. Therefore, eq. (3.2) can be rewritten as d Φk ∂ = ³ ρ k φk dV + ³ ρ k (Vk ,rel ⋅ n k )φk dA (3.3) A dt system ∂t V
which is the final form of the transformation formula that relates the change of property for a fixedmass system to that of the control volume. It states that the rate of change of a property Φk for a fixedmass system is equal to the rate of change of Φk in the control volume (the first term on the righthand side) plus the
180 Transport Phenomena in Multiphase Systems
net rate of efflux of Φk by mass flow into or out of the control volume (the second term on the righthand side). It should be pointed out that the control volume moves with the reference frame, which moves with a constant velocity, Vref. The coordinate system is attached to and moves with the reference frame. In other words, the coordinate system is stationary relative to the reference frame. The reference frame is inertial, so Newton’s second law is valid in the coordinate system that moves with the reference frame. Equation (3.3) will be used to obtain the macroscopic (integral) formulation of the basic laws for a control volume. In the development of eq. (3.3), it was assumed that the control volume is located in a region occupied only by the kth phase. For multiphase flow and heat transfer problems, the control volume may include multiple phases separated by several interfaces. Figure 3.3 shows a fixed control volume V with surface area A that contains three phases (solid, liquid, and gas): (1) solid with volume Vs (t ) bounded by As (t ) , AIsg (t ) (interfacial area between solid and gas), and
AIsA (t ) (interfacial area between solid and liquid); (2) liquid with volume VA (t )
bounded by AA (t ) , AIAg (t ) (interfacial area between solid and liquid), and
AIsA (t ) ; and (3) gas with volume Vg (t ) bounded by Ag (t ) , AIsg (t ) and
AIAg (t ) . Since the interfaces that separate different phases move, the volume and surface area of each phase are functions of time. The volume and area of each phase satisfy Vs (t ) + VA (t ) + Vg (t ) = V and As (t ) + AA (t ) + Ag (t ) = A, because the control volume is fixed in space. For a control volume that contains three phases, eq. (3.3) becomes
z
Vk,rel
0 y x
Figure 3.2 Control volume in a flow field.
Chapter 3 Generalized Governing Equations: Local Instance Formulations 181
dΦ dt
=
system
∂ ∂ ∂ ³Vs (t ) ρ sφs dV + ∂t ³VA (t ) ρAφA dV + ∂t ³Vg (t ) ρ gφg dV ∂t +³
As ( t )
ρ s (Vs ,rel ⋅ n s )φs dA + ³ ρ g (Vg , rel ⋅ n g )φg dA
AA ( t )
ρ A (VA ,rel ⋅ n A )φA dA
(3.4)
+³
Ag ( t )
If the control volume contains only two phases, as shown in Fig. 3.3(b), eq. (3.4) may be reduced to dΦ ∂ ∂ = ³ ρ1φ1dV + ³ ρ 2φ2 dV V1 ( t ) dt system ∂t ∂t V2 ( t ) (3.5) + ³ ρ1 ( V1, rel ⋅ n1 )φ1dA + ³ ρ 2 ( V2, rel ⋅ n 2 )φ2 dA
A1 ( t ) A2 ( t )
Therefore, for control volumes containing only one phase, the local instance integral formulation can be obtained by using eq. (3.3). For the control volume that contains more then one phase separated by interfaces, the local instance integral formulation can be obtained by using eq. (3.4) or (3.5). In subsection 3.2.1 to 3.2.5, the integral form of the governing equations for a control volume that contains only one phase will be derived first, followed by an extension of the integral form of the governing equations to a control volume that contains multiple phases.
(a) Three phases
(b) Two phases
Figure 3.3 Control volumes that enclose multiple phases separated by interfaces. where the subscripts 1 and 2 can be s, A , or g, depending on the phases in the control volume.
182 Transport Phenomena in Multiphase Systems
3.2.1 Conservation of Mass
The law of the conservation of mass dictates that mass may be neither created nor destroyed. For a control volume that contains only one phase k, conservation of mass can be obtained by setting the general and specific property forms to Φ k = mk and φk = 1 in eq. (3.3), i.e.,
dmk dt =
system
∂ ρ k dV + ³ ρ k (Vk ,rel ⋅ n k )dA A ∂t ³V
(3.6)
Since the mass of a fixedmass system is constant by definition, and the fixedmass system contains only one phase, the resulting formulation of conservation of mass is ∂ ρ k dV + ³ ρ k (Vk ,rel ⋅ n k )dA = 0 (3.7) A ∂t ³V For a control volume containing multiple phases separated by interfaces, the conservation of mass can be similarly obtained by using eq. (3.4) or (3.5), i.e., Π ª∂ º (3.8) ¦ « ∂t ³Vk (t ) ρk dV + ³Ak (t ) ρk (Vk ,rel ⋅ n k )dA» = 0 ¼ k =1 ¬ where Π is the number of phases.
3.2.2 Momentum Equation
Newton’s second law states that, in an inertial reference frame, the time rate of momentum change of a fixed mass system is equal to the net force acting on the system, and it takes place in the direction of the net force. Mathematically, Newton’s second law of motion for fixedmass in a reference frame that moves at a constant velocity Vref (see Fig. 3.2) is written as d (mV ) (3.9) ¦ Fk = dt k ,rel where the lefthand side is the net force vector acting on the fixedmass system, and the righthand side is the rate of momentum change. For control volumes that contain only one phase k, the integral form of Newton’s second law can be obtained by using eq. (3.3). With the applicable value of Φk and φk in eq. (3.3) defined as Φ k = mk Vk , rel and φk = Vk , rel , one obtains ∂ (3.10) ¦ F k = ∂t ³V ρk Vk ,rel dV + ³A ρk (Vk ,rel ⋅ n k )Vk ,rel dA which is in the vector form and is valid in all three directions. Forces acting on the control volume include body forces and contact forces that act on its surface. For example, for a multicomponent system that contains N components, if the body force per unit volume acting on the ith species in the kth phase is Xk,i, the total body force acting on the control volume is
Chapter 3 Generalized Governing Equations: Local Instance Formulations 183
where ρ k is the density of the k phase that contains N components. The stress tensor acts on the surface of a fluid control volume, and includes both normal and shear stresses. The net force may be written as ªN º Fk = ³ « ¦ ρ k ,i X k ,i »dV + ³ ′ , rel ⋅ n k dA (3.11) ¦ V Ak ¬ i =1 ¼ where ′ is the total stress tensor and n k is the local normal unit vector on k surface A. The dot product of a tensor of rank two, ′ , and a vector, n k , is a k vector that represents the force acting on the surface of the control volume per unit area. Combining eqs. (3.10) and (3.11), we obtain the momentum equation for the control volume: ªN º ³V «¦ ρk ,i Xk ,i »dV + ³A ′k ,rel ⋅ n k dA ¬ i =1 ¼ (3.12) ∂ = ³ ρ k Vk , rel dV + ³ ρ k (Vk , rel ⋅ n k )Vk , rel dA A ∂t V where the two terms on the lefthand side represent, respectively, the body force and stress on the control volume, and the two terms on the righthand side represent, respectively, the rate of momentum change in the control volume and the rate of the momentum flow into or out of the control volume. Vk , rel is the bulk velocity of the kth phase that contains N components. When the control volume includes multiple phases, integrations must be performed for each subvolume. In that case, the momentum equation becomes Π ª º ªN º ¦ « ³Vk (t ) «¦ ρk ,i Xk ,i »dV + ³Ak (t ) ′k ,rel ⋅ n k dA» k =1 ¬ ¬ i =1 ¼ ¼ (3.13) Π ª∂ º = ¦ « ³ ρ k Vk , rel dV + ³ ρ k (Vk , rel ⋅ n k )Vk ,rel dA» Vk ( t ) Ak ( t ) ¼ k =1 ¬ ∂t Equations (3.12) and (3.13) are momentum equations in a coordinate system that is attached to and moves with an inertial reference frame. For a fixed coordinate system that does not move with the reference frame (while the control volume still moves with the reference frame at velocity Vref), one can substitute the general variables Φ k = mk Vk and φk = Vk into eq. (3.3) to obtain the momentum equation.
th
ªN º th ³V «¦ ρk ,i Xk ,i »dV , where ρk ,i is the mass concentration of the i species in the ¬ i =1 ¼ kth phase. If the body force per unit mass is the same for different species and phase (as is the case with gravity), the body force term is reduced to ³ ρ k XdV ,
V
184 Transport Phenomena in Multiphase Systems
3.2.3 Energy Equation
The first law of thermodynamics for a fixedmass system (closed system) states that ˆ dEk δ Qk δ Wk = − (3.14) dt system dt dt where Qk is positive if heat is transferred into the system, and Wk is positive if the work is done by the system to the surrounding. The mass of system can store ˆ energy internally in a number of different forms. Therefore, the total energy E is due to internal energy, E, knietic, potential, electomagnitic, surface tension and other forms. For the development of equations in Chapter 3 and 4, the contributions of internal and kinetic energies are considered. Other contributions can also be added as a source term or boundary condition based on the physical model. For a control volume including only one phase k, the energy equation for the control volume can be obtained by setting the general property as Φ k = Ek + mVk2,rel / 2 and φk = ek + Vk2, rel / 2 in eq. (3.3), i.e.,
· § Vk2,rel · ¸ dV + ³A ρ k (Vk ,rel ⋅ n k ) ¨ ek + ¸ dA (3.15) ¸ ¨ 2¸ system ¹ © ¹ where Ek and ek are, respectively, internal energy and specific internal energy,
dEk dt =
§ V2 ∂ ρ k ¨ ek + k ,rel 2 ∂t ³V ¨ ©
and Vk2,rel = Vk , rel ⋅ Vk ,rel . In arriving at eq. (3.15), it is assumed that potential energy is negligible. Substituting eq. (3.15) into eq. (3.14), an integral form of the energy equation may be written as § § V2 · V2 · δ Qk δ Wk ∂ − = ³ ρ k ¨ ek + k ,rel ¸ dV + ³ ρ k (Vk , rel ⋅ n k ) ¨ ek + k ,rel ¸ dA (3.16) ¨ A 2¸ 2¸ ∂t V ¨ dt dt © ¹ © ¹ The heat flow, δ qk / dt , is attributed to conduction (and radiation, for the case of participating medium) across the boundary and/or to internal generation, i.e., δ Qk ′′′ = ³ −q′′ ⋅ n k dA + ³ qk dV k A V dt where q′′ is the heat flux vector at the control volume surface, which can be k caused by temperature or concentration gradients, as indicated in eq. (1.64). The dot product of the heat flux vector q′′ with the unit normal vector n k gives k ′′′ heat conducted out of the control volume. qk is the internal heat generation per unit volume, and it can be caused by chemical reaction, electrical heating, etc. The work rate can include work done by normal and shear stress on the surface of the control volume, δ Wk , s / dt , as well as work done by the body force, δ Wk ,b / dt . The surface work is evaluated by taking the scalar product of
Chapter 3 Generalized Governing Equations: Local Instance Formulations 185
the force acting on the surface of the control volume per unit area, n k ⋅ the velocity, Vk ,rel , over the entire surface area of the control volume:
k
, and
= − ³ (n k ⋅ ′ ) ⋅ Vk , rel dA k A dt The work done by the body force is δ Wk ,b ªN º = − ³ « ¦ ρ k ,i X k ,i » ⋅ Vk , rel dV V dt ¬ i =1 ¼ Finally, the integral form of the energy equation (3.16) becomes § § Vk2,rel · Vk2, rel · ∂ ρ k ¨ ek + ¸ dV + ³A ρ k ( Vk , rel ⋅ n k ) ¨ ek + ¸ dA ¨ ∂t ³V ¨ 2¸ 2¸ © ¹ © ¹ (3.17) ªN º ′′′ = − ³ q′′ ⋅ n k dA + ³ qk dV + ³ (n k ⋅ ′ , rel ) ⋅ Vk , rel dA + ³ « ¦ ρ k ,i X k ,i » ⋅ Vk , rel dV k Ak V A V ¬ i =1 ¼ where the lefthand side represents the rate of total energy change in the control volume and the rate of total energy flow into or out of the control volume. The four terms in the righthand side represent heat transfer across the boundary of the control volume, internal heat generation, work done by the stress at the boundary of control volume, and the work done by the body force. For control volumes including multiple phases, the energy equation becomes Πª § § Vk2,rel · Vk2, rel · º ∂ ¦ « ∂t ³Vk (t ) ρk ¨ ek + 2 ¸ dV + ³Ak (t ) ρk (Vk ,rel ⋅ n k ) ¨ ek + 2 ¸ dA» ¨ ¸ ¨ ¸ k =1 « © ¹ © ¹» ¬ ¼
′′′ = ¦ ª − ³ q′′ ⋅ n k dA + ³ qk dV + ³ (n k ⋅ ′ , rel ) ⋅ Vk , rel dA k « Ak ( t ) k Vk ( t ) Ak ( t ) k =1 ¬ +³
Π
δ Wk , s
(3.18)
º ªN º « ¦ ρ k ,i X k ,i » ⋅ Vk , rel dV » Vk ( t ) ¬ i =1 ¼ ¼ where the second integral in the bracket on the lefthand side is advection of energy due to mass flow.
3.2.4 The Second Law of Thermodynamics
The second law of thermodynamics states that the entropy generation in a closed system (fixedmass) must be greater than or equal to zero. The entropy change for a system with fixedmass and contains only one phase k can be obtained by setting Φ k = Sk , φk = sk in eq. (3.3), i.e., dSk dt =
system
d ρ k sk dV + ³ ρ k (Vk ,rel ⋅ n k ) sk dA A dt ³V
(3.19)
The change of entropy in a closed system results from heat transfer and/or entropy generation:
186 Transport Phenomena in Multiphase Systems
dSk dt
system
=³
A
−q′′ ⋅ n k q′′′ k ′′′ dA + ³ k dV + ³ sgen ,k dV VT V Tk k
(3.20)
where, on the righthand side of eq. (3.20), the first term represents the change of entropy due to heat transfer across the boundary of the control volume, and the second term represents the change of entropy due to internal heat generation in the control volume. The last term represents entropy generation, which should always be greater than or equal to zero, i.e., (3.21) gen ³ s′′′ ,k dV ≥ 0
V
Combining eqs. (3.19) and (3.20) and applying eq. (3.21), one obtains the integral form of the second law of thermodynamics: d ρ k sk dV + ³ ρ k (Vk ,rel ⋅ n k ) sk dA A dt ³V (3.22) q′′ ⋅ n q′′′ + ³ k k dA − ³ k dV = ³ s′′′ ,k dV ≥ 0 A VT V gen Tk k If the control volume includes Π phases, the second law of thermodynamics must be obtained by integrating over the two phases separately, Π ªd q′′ ⋅ n q′′′ º ¦ « dt ³Vk (t ) ρk sk dV + ³Ak (t ) ρk (Vk ,rel ⋅ n k )sk dA + ³Ak (t ) kT k dA − ³Vk (t ) Tk dV » k =1 ¬ k k ¼
= ¦³
k =1
Π
Vk ( t ) gen , k
s′′′ dV + ³
AI ( t ) gen , I
s′′
dA ≥ 0
(3.23)
The entropy generation for a control volume including Π phases consists of entropy generation in each phase, plus that in the interfaces. The second law of thermodynamics requires that each of these entropy generations be greater than or equal to zero.
3.2.5 Species
The continuity equation states that the total mass for a closed system is constant. For a system containing one phase but more than one component, the total mass of the system is composed of different species. If the concentrations of each of these species are not uniform, mass transfer occurs in a way that makes the concentrations more uniform. Therefore, it is necessary to track the individual components by applying the principle of conservation of species mass. For a system of multiple components, each component can have its own mass density and velocity. The conservation of mass for the ith species in the kth phase for a fixedmass system is obtained by applying eq. (3.3) for the ith species with Φ k = mk ,i and φk = ρ k ,i / ρ k , i.e., dmk ,i dt
system
=
∂ ρ k ,i dV + ³ ρ k ,i ( Vk ,i ,rel ⋅ n k )dA A ∂t ³V
(3.24)
Chapter 3 Generalized Governing Equations: Local Instance Formulations 187
If there is no chemical reaction, the total mass of the ith species for a closed system remains constant. Chemical reactions, on the other hand, will result in the production or consumption of the ith species, which can be modeled as a mass source or sink for the ith species in the kth phase, i.e., dmk ,i ′′′ (3.25) = mk ,i dV dt system ³V
' where mk",i is the mass production rate for the ith species in the kth phase, which can be either positive (mass production) or negative (mass consumption). For the ' cases without chemical reaction, mk",i = 0 . The conservation of species mass for a control volume with phase k is therefore ∂ ′′′ ρ k ,i dV + ³ ρ k ,i (Vk ,i ,rel ⋅ n k ) dA = ³ mk ,i dV (3.26) A V ∂t ³V which is valid for each component in the control volume. If the total number of species is N, summation of the conservation of mass of all species results in ∂ ªN º ªN º ªN º (3.27) ³V «¦ ρk ,i » dV + ³A «¦ ( ρk ,i Vk ,i ,rel ) ⋅ n k » dA = ³V «¦ mk′′′,i » dV ∂t ¬ i =1 ¼ ¬ i =1 ¼ ¬ i =1 ¼ The summation of the densities of the individual species is equal to the bulk density of the multicomponent substance:
ρ k = ¦ ρ k ,i
i =1
N
(3.28)
The bulk velocity of the multicomponent substance is the massaveraged velocity of the velocity of each individual species:
ρ k Vk ,rel = ¦ ( ρ k ,i Vk ,i ,rel )
i =1
N
(3.29)
The righthand side of eq. (3.27) must be zero, because the total mass of species produced is equal to the total mass of species consumed, i.e.,
¦ m′′′
i =1
N
k ,i
=0
(3.30)
Substituting eqs. (3.28) – (3.30) into eq. (3.27), the continuity equation (3.7) is obtained. This means that we can write conservation of species mass for N species, but only N1 of these equations are independent. In practice, one can use these N1 equations for conservation of species mass in combination with the continuity equation, eq. (3.7), to describe multicomponent systems. In each equation of the conservation of species mass, there are two new unknowns: mass density of the ith species ρk,i, and the species mass velocity Vk ,i , rel . This yields more unknowns than the number of equations. So in order to have a properlyposed equation set, it is necessary to reduce the number of unknowns in each equation from two to one. The second term on the lefthand side of eq. (3.26) represents the species i mass flow across the surface of the control volume, which results from convection by bulk flow and diffusion relative to the bulk convection.
188 Transport Phenomena in Multiphase Systems
where J k ,i is the diffusive mass flux vector of species i in the kth phase, which includes mass fluxes due to ordinary diffusion driven by the concentration gradient, and thermal (Soret) diffusion (see Table 1.7). Substituting eq. (3.31) into eq. (3.26), one obtains an expression for the conservation of species mass that contains only one new additional variable, ρk,i: ∂ ′′′ ρ k ,i dV + ³ ρ k ,i (Vk ,rel ⋅ n k )dA = − ³ J k ,i ⋅ n k dA + ³ mk ,i dV (3.32) A A V ∂t ³V This analysis is based on the assumption that the control volume contains only one phase. For a control volume containing Π phases, the conservation of species mass is Π ª∂ º ¦ « ∂t ³Vk (t ) ρk ,i dV + ³Ak (t ) ρk ,i (Vk ,rel ⋅ n k )dA» ¼ k =1 ¬ (3.33) Π ª− º ′′′ = ¦ ³ (J k ,i ⋅ n k ) dA + ³ mk ,i dV « Ak ( t ) » Vk ( t ) ¼ k =1 ¬
³
A
ρ k ,i (Vk ,i ,rel ⋅ n k )dA = ³ ρ k ,i (Vk , rel ⋅ n k )dA + ³ J k ,i ⋅ n k dA
A A
(3.31)
3.3 Microscopic (Differential) Formulation
The microscopic (differential) formulations to be presented here include phase equations and jump conditions. The former apply within a particular phase, and the latter are valid at the interface that separates two phases. The phase equations for a particular phase should be the same as those for a singlephase system. Most textbooks (e.g., White, 1991; Incropera and DeWitt, 2001; Bejan, 2004; Kays et al., 2004) obtain the governing equations for a singlephase system by performing mass, momentum, and energy balances for a microscopic control volume. We will obtain the phase equations by analyzing the integral equations for a finite control volume that includes only one phase. Jump conditions at the interface will be obtained by analyzing the integral equations for a control volume that includes two phases separated by an interface. In order to obtain the differential formulations, it is necessary to apply the divergence theorem from vector calculus. For a general vector quantity , which is continuously differentiable, and for a control volume V, enclosed by a piecewise smooth control surface A, the divergence theorem states that (3.34) ³ ⋅ ndA = ³ ∇ ⋅ dV
A V
The desired differential equations may be obtained by applying this relationship to the integral form of the basic laws. Furthermore, since the control volume shape and size are fixed in time, the Leibniz’s rule for the specific general quantity φk is also valid: ∂ ( ρ k φk ) d (3.35) ³V ρkφk dV = ³V ∂t dV dt
Chapter 3 Generalized Governing Equations: Local Instance Formulations 189
3.3.1 Conservation of Mass
The surface integral in eq. (3.7) may be converted to a volume integral by applying eq. (3.34) as follows: (3.36) ³ ρk (Vk ,rel ⋅ nk )dA = ³ ∇ ⋅ ρk Vk ,rel dV
A V
Substituting the eq. (3.36) into eq. (3.7) and considering eq. (3.35), the entire lefthand side of eq. (3.7) is included in a single volume integral, i.e., § ∂ρ k · (3.37) ³V ¨ ∂t + ∇ ⋅ ρk Vk ,rel ¸ dV = 0 © ¹ The only condition that ensures that eq. (3.37) is true for any arbitrary shape and size of the control volume is that the integrand must be equal to zero, i.e., ∂ρ k + ∇ ⋅ ρ k Vk ,rel = 0 (3.38) ∂t which is the differential form of the law of the conservation of mass. Equation (3.38) can also be rewritten as Dρk + ρ k ∇ ⋅ Vk , rel = 0 (3.39) Dt where D/Dt is the substantial derivative (or material derivative) defined by D∂ ≡ + Vk , rel ⋅ ∇ (3.40) Dt ∂t For a stationary reference frame, Vk , rel = Vk , and eq. (3.39) becomes
Dρk + ρ k ∇ ⋅ Vk = 0 (3.41) Dt For incompressible flow, in which the density of the fluid is constant ( ρ k ≡ const ), eq. (3.41) simplifies as ∇ ⋅ Vk = 0 (3.42) For steadystate compressible flow, the continuity equation can be obtained by simplifying eq. (3.38), i.e., (3.43) ∇ ⋅ ρ k Vk , rel = 0 Example 3.1 For a stationary reference frame, give the continuity equation in the Cartesian coordinate system for both compressible and incompressible fluid flow. Solution: The ∇ operator in the Cartesian coordinate system is [see eq. (C.6) in Appendix C] ∂ ∂ ∂ ∇=i + j +k (3.44) ∂x ∂y ∂z The velocity vector in a Cartesian system is Vk = iuk + jvk + kwk (3.45)
190 Transport Phenomena in Multiphase Systems
where uk , vk , and wk are the velocity components in the x, y, and zdirections, respectively. The continuity equation for compressible fluid flow in a stationary reference frame can be obtained by substituting eqs. (3.44) and (3.45) into eq. (3.41), i.e., § ∂u Dρk ∂v ∂w · + ρk ¨ k + k + k ¸ = 0 (3.46) Dt ∂y ∂z ¹ © ∂x where D ρ k ∂ρ k ∂ρ ∂ρ ∂ρ = + uk k + vk k + wk k (3.47) ∂t ∂x ∂y ∂z Dt For incompressible flow, eq. (3.46) simplifies as ∂uk ∂vk ∂wk + + =0 (3.48) ∂x ∂y ∂z
3.3.2 Momentum Equation
The integral form of Newton’s second law for a control volume that includes only one phase is expressed by eq. (3.12). The surface integral terms in eq. (3.12) can be rewritten using the divergence theorem: (3.49) ³ ′k ,rel ⋅ n k dA = ³ ∇ ⋅ ′k ,rel dV
A V
³
A
ρ k ( Vk , rel ⋅ n k )Vk ,rel dA = ³ ∇ ⋅ ρ k Vk ,rel Vk ,rel dV
V
(3.50)
Substituting eq. (3.49) and (3.50) into eq. (3.12), and considering eq. (3.35), the entire equation can be rewritten as a volume integral: N ∂ ª º ∇ ⋅ ′ ,rel + ¦ ρ k ,i X k ,i − ( ρ k Vk , rel ) − ∇ ⋅ ρ k Vk , rel Vk , rel » dV = 0 (3.51) k ³V « ∂t i =1 ¬ ¼ As was the case for the continuity equation, the integrand must equal zero to assure the general validity of eq. (3.51); so, one obtains the desired differential form of the momentum equation: N ∂ ( ρ k Vk , rel ) + ∇ ⋅ ρ k Vk , rel Vk , rel = ∇ ⋅ ′ , rel + ¦ ρ k ,i X k ,i (3.52) k ∂t i =1 The derivatives on the lefthand side of eq. (3.52) may be expanded to yield ª ∂V º ª ∂ρ º Vk ,rel « k + ∇ ⋅ ρ k Vk , rel » + ρ k « k ,rel + Vk , rel ⋅ ∇Vk , rel » ¬ ∂t ¼ ¬ ∂t ¼ (3.53) = ∇ ⋅ ′ ,rel + ¦ ρ k ,i X k ,i k
i =1 N
The first bracketed term on the left vanishes, as required by the continuity eq. (3.38). The second term may be written more simply in substantial derivative form, and the entire equation becomes
Chapter 3 Generalized Governing Equations: Local Instance Formulations 191
(3.54) Dt i =1 The stress tensor, ′ , rel , is the sum of an isotropic thermodynamic stress, k (3.55) Substituting eq. (3.55) into eq. (3.54), the momentum equation becomes N DVk , rel (3.56) ρk = −∇pk + ∇ ⋅ k , rel + ¦ ρ k ,i X k ,i Dt i =1 The viscous stress tensor measured in the reference frame, k , rel , can be determined by using Newton’s law of viscosity [see eq. (1.53)]: 2 μk (∇ ⋅ Vk ,rel )I (3.57) k , rel = 2 μ k D k , rel − 3 where μk is viscosity of the kth phase. Dk,rel is the rate of strain tensor in phase k, i.e., T 1 (3.58) Dk , rel = ª∇Vk , rel + ( ∇Vk , rel ) º « » ¬ ¼ 2 and I in eq. (3.57) is the unit tensor that satisfies a ⋅ I = I ⋅ a = a for any tensor a. The diagonal components of I are equal to one and all other components are zero: 1 i = j (i, j = 1, 2,3) I ij = ® (3.59) ¯0 i ≠ j If the fluid is incompressible ( ρ k =const ), the second term on the righthand side of eq. (3.55) will be zero according to eq. (3.42). The momentum equation (3.54) then becomes N DVk , rel (3.60) ρk = ¦ ρ k ,i X k ,i − ∇pk + ∇ ⋅ ( μ k ∇Vk , rel ) Dt i =1 where the lefthand side is the inertial term (mass per unit volume ρk times acceleration, DVk,rel/Dt). The three terms on the righthand side represent body force per unit volume, pressure force per unit volume, and viscous force per unit volume, respectively. For DVk , rel / Dt = 0 , we have Stokes’ flow or creep flow, and eq. (3.60) becomes elliptic and is similar to the steadystate conduction equation. In a Cartesian coordinate system, the vector form of the momentum equation, eq. (3.60), for incompressible and Newtonian fluid with constant viscosity can be written as three equations in the x, y, and zdirections: N § ∂ 2 uk ,rel ∂ 2uk ,rel ∂ 2uk ,rel · Duk ,rel ∂p (3.61) = ¦ ρ k , i X k ,i − k + μ k ¨ + + ρk ¸ ¨ ∂x 2 Dt ∂x ∂y 2 ∂z 2 ¸ i =1 © ¹
k , rel
ρk
DVk ,rel
= ∇ ⋅ ′ ,rel + ¦ ρ k ,i X k ,i k
N
− pk I , and the viscous stress tensor, ′ ,rel = − pk I + k
k , rel
[defined in eq. (1.50)], i.e.,
192 Transport Phenomena in Multiphase Systems
ρk ρk
Dvk , rel Dt Dwk ,rel
= ¦ ρ k ,iYk ,i −
i =1 N
N
§ ∂ 2 vk , rel ∂ 2 vk ,rel ∂ 2 vk ,rel ∂pk + μk ¨ + + ¨ ∂x 2 ∂y ∂y 2 ∂z 2 ©
· ¸ ¸ ¹
(3.62)
· (3.63) ¸ ¸ Dt i =1 ¹ where X k ,i , Yk ,i , and Z k ,i are the components of body force per unit volume acting on the ith species of the kth phase in the x, y, and z directions, respectively. For the case that the only body force is gravity, X k ,i = g , eq. (3.60) becomes DVk ,rel ρk = ρ k g − ∇pk + ∇ ⋅ ( μk ∇Vk ,rel ) (3.64) Dt For natural convection problem, it is often assumed that the fluid is incompressible except in the first term on the righthand side of eq. (3.64); this is referred to as the Boussinesq assumption. The density of a mixture is a function of temperature and mass fractions of species. It can be expanded using a Taylor’s series near the vicinity of a reference point ( Tk , ωk ,1 , ωk ,2 , "ωk , N ): = ¦ ρ k ,i Z k ,i
§ ∂ 2 wk ,rel ∂ 2 wk ,rel ∂ 2 wk , rel ∂pk − + μk ¨ + + ¨ ∂x 2 ∂z ∂y 2 ∂z 2 ©
ρk = ρk +
N ∂ρ k ∂ρ (Tk − Tk ) + ¦ k (ωk ,i −ωk ,i ) + " ∂Tk i =1 ∂ωk ,i N
(3.65)
≈ ρ k − ρ k β k (Tk − Tk ) − ρ k ¦ β mk (ωk ,i −ωk ,i )
i =1
where ρ k is density at the reference point, β k = −(∂ρ k / ∂T ) / ρ k is the coefficient of thermal expansion, and β mk = −(∂ρ k / ∂ωk ,i ) / ρ k is composition coefficient of volume expansion. Substituting eq. (3.65) into eq. (3.64), the momentum equation for natural convection is obtained DVk ,rel ρk = ( −∇pk + ρ k g ) − ρ k gβ k (Tk − Tk ) Dt (3.66) N − ρ k g ¦ β mk (ωk ,i −ωk ,i ) + ∇ ⋅ ( μk ∇Vk , rel )
i =1
where the second and third terms on the righthand side of eq. (3.66) describe the effect of buoyance force due to temperature and composition variation with in the system, respectively.
3.3.3 Energy Equation
The surface integrals on the righthand side of eq. (3.17) can be rewritten as volume integrals using eq. (3.34), i.e.,
Chapter 3 Generalized Governing Equations: Local Instance Formulations 193
§ Vk2, rel ³A ρk (Vk ,rel ⋅ n) ρk ¨ ek + 2 ¨ ©
ª · § Vk2, rel ¸ dA = ³V ∇ ⋅ « ρ k Vk , rel ρ k ¨ ek + ¸ ¨ 2 « ¹ © ¬
V
·º ¸ » dV ¸ ¹» ¼
(3.67) (3.68) (3.69)
³
A
−q′′ ⋅ n k dA = ³ −∇ ⋅ q′′dV k k
V
³n
A
k
⋅ ′ ,rel ⋅ Vk ,rel dA = ³ ∇ ⋅ ( ′ k ,rel ⋅Vk , rel )dV k
Substituting eqs. (3.35) yields ∂ ª § Vk2, rel ° ® « ρ k ¨ ek + ³V ∂t « ¨ 2 °¬© ¯
(3.67) – (3.69) into eq. (3.17) and considering eq.
ª ·º § Vk2,rel ¸ » + ∇ ⋅ « ρ k Vk ,rel ¨ ek + ¸ ¨ 2 « ¹» © ¼ ¬
N
·º ¸» ¸ ¹» ¼
(3.70)
½ § · ′′′ +∇ ⋅ q′′ − qk − ∇ ⋅ ( ′ k ,rel ⋅Vk ,rel ) − ¨ ¦ ρ k ,i X k ,i ¸ ⋅ Vk , rel ¾ dV = 0 k © i =1 ¹ ¿ For eq. (3.70) to be true for any arbitrary control volume, the integrand must be zero. The result is the general differential form of the energy equation: ª § Vk2,rel · º Vk2,rel · º ∂ª § « ρ k ¨ ek + ¸ » + ∇ ⋅ « ρ k Vk ,rel ¨ ek + ¸» ¨ 2 ¸» 2 ¸» ∂t « ¨ « ¹¼ © ¹¼ ¬© ¬ (3.71) N § · ′′′ = −∇ ⋅ q′′ + qk + ∇ ⋅ ( ′ k , rel ⋅Vk , rel ) + ¨ ¦ ρ k ,i X k ,i ¸ ⋅ Vk , rel k © i =1 ¹ The lefthand side of eq. (3.71) can be rewritten as ª∂ § § § V2 · V 2 ·º V 2 · ª ∂ρ º ek + k ,rel ¸ « k + ∇ ⋅ ( ρ k Vk ,rel ) » + ρ k « ¨ ek + k , rel ¸ + Vk , rel ⋅ ∇ ¨ ek + k , rel ¸ » ¨ ¨ ¸ ¨ ¸ ¨ ¸ 2 ¹ ¬ ∂t 2¹ 2 ¹» ¼ « ∂t © © © ¬ ¼ N § · ′′′ (3.72) = −∇ ⋅ q′′ + qk + ∇ ⋅ ( ′ k , rel ⋅Vk , rel ) + ¨ ¦ ρ k ,i X k ,i ¸ ⋅ Vk , rel k © i =1 ¹ According to the continuity eq. (3.38), the first bracketed term on the lefthand side of eq. (3.72) is zero. The second term on the lefthand side may be written more simply in substantial derivative form, i.e., Vk2,rel · D§ ′′′ ρk ¨ ek + ¸ = −∇ ⋅ q′′ + qk k Dt ¨ 2¸ © ¹ (3.73) §N · +∇ ⋅ ( ′ k , rel ⋅Vk ,rel ) + ¨ ¦ ρ k ,i X k ,i ¸ ⋅ Vk , rel © i =1 ¹ which is the total energy (including thermal and mechanical energies) balance equation. It would be more convenient to remove the mechanical energy terms from eq. (3.73). The mechanical energy balance equation can be obtained by forming a dot (scalar) product of the momentum equation (3.54) with the velocity vector Vk,rel, i.e.,
194 Transport Phenomena in Multiphase Systems
§N · (3.74) ⋅ Vk ,rel = ( ∇ ⋅ ′ ,rel ) ⋅ Vk , rel + ¨ ¦ ρ k ,i X k ,i ¸ ⋅ Vk , rel k Dt © i =1 ¹ which can be rearranged to obtain 2 D § Vk ,rel · §N · ρk ¨ ¸ = ∇ ⋅ ( ′ ,rel ⋅ Vk , rel ) − ∇Vk : k , rel + ¨ ¦ ρ k ,i X k ,i ¸ ⋅ Vk , rel (3.75) k ¨2¸ Dt © © i =1 ¹ ¹ Subtracting eq. (3.75) from eq. (3.73) yields the following thermal energy equation: De ′′′ ρ k k = −∇ ⋅ q′′ + qk − p∇ ⋅ Vk , rel + ∇Vk , rel : k ,rel (3.76) k Dt where the lefthand side represents the rate of gain of internal energy per unit volume. The terms in the righthand side are the rate of internal energy input by heat transfer per unit volume, the internal heat generation per unit volume, the reversible rate of internal energy increase per unit volume by compression, and the irreversible rate of internal energy increase per unit volume by viscous dissipation, respectively; the viscous dissipation is the contraction [see eq. (C.35)] of two tensors of rank two: ∇Vk , rel and k , rel , i.e.,
ρk
DVk , rel
∇Vk ,rel :
k , rel
=
∂uk ,rel
τ yz + τ zx + τ zy + τ zz ∂z ∂x ∂y ∂z Equation (3.76) is the energy equation expressed in terms of internal energy. In order to obtain an equation that contains enthalpy, the definition of enthalpy is employed: p hk = ek + (3.78)
+
∂x ∂vk , rel
τ xx +
∂uk ,rel
∂y ∂wk ,rel
τ xy +
∂uk , rel
∂z ∂wk , rel
τ xz +
∂vk , rel
∂x ∂wk , rel
τ yx +
∂vk , rel ∂y
τ yy
(3.77)
ρk
Substituting eq. (3.78) into eq. (3.76), and considering continuity equation (3.38), the energy equation in term of enthalpy and temperature is obtained: Dh Dp ′′′ ρ k k = −∇ ⋅ q′′ + + qk + ∇Vk , rel : k , rel (3.79) k Dt Dt For a multicomponent system, the enthalpy can be expressed as
hk = ¦ ωk ,i hk ,i
i =1 N
(3.80)
where the ωk ,i and hk ,i are the mass fraction and specific enthalpy of the ith component in the kth phase. Substituting eq. (3.80) into eq. (3.79), the energy equation for a multicomponent system becomes N Dp ªD º (3.81) ¦ « ρk Dt (ωk ,i hk ,i )» = −∇ ⋅ q′′ + Dtk + ∇Vk ,rel : k ,rel k ¼ i =1 ¬ The heat flux vector q′′ can be obtained from eq. (1.64) k
Chapter 3 Generalized Governing Equations: Local Instance Formulations 195
T xk ,i xk , j Dk ,i § J k ,i J k , j · (3.82) − ¨ ¸ ρ k ,i Dk ,ij ¨ ρ k ,i ρ k , j ¸ i =1 i =1 © ¹ where the three terms on the righthand side represent the heat flux due to Fourier conduction, interdiffusional convection, and concentration gradient (Dufour effect). The kth phase is assumed to be isotropic, so its thermal conductivity is the same for any direction. Substituting eq. (3.82) into eq. (3.81), the energy equation becomes N N ªD ρk (ωi hk ,i )º = ∇ ⋅ (kk ∇Tk ) − ¦ ∇ ⋅ ( J k ,i hk ,i ) ¦ « Dt » ¼ i =1 ¬ i =1 (3.83) T Nxx ª § J k ,i J k , j · º Dpk k ,i k , j Dk , i −∇ ⋅ « cRuT ¦ − + ∇Vk , rel : k , rel ¨ ¸» + ρ k ,i Dk ,ij ¨ ρ k ,i ρ k , j ¸ » Dt i =1 « © ¹¼ ¬ where the second and third terms on the right hand side, which represent the contributions of interdiffusional convection and Dufour effect, are both negligible for most applications. For a pure substance, eq. (3.83) is reduced to: Dh Dp ρ k k = ∇ ⋅ (kk ∇Tk ) + k + ∇Vk ,rel : k , rel (3.84) Dt Dt The enthalpy of a pure substance (single component) can be expressed as a function of temperature and pressure, hk = f (Tk , pk ), i.e.,
q′′ = − kk ∇Tk + ¦ hk ,i J k ,i + cRuT ¦ k
N
N
Dhk § ∂hk · DTk § ∂hk · Dpk = +¨ ¸ Dt ¨ ∂T ¸ P Dt © ∂p ¹T Dt © ¹ Thermodynamic relations give us § ∂hk · 1 − β k Tk § ∂hk · ¨ ∂T ¸ = c p , ¨ ∂p ¸ = ρ © ¹P © ¹T k
(3.85)
(3.86)
where 1 § ∂ρ k · (3.87) ¨ ¸ ρ k © ∂Tk ¹ p is the coefficient of thermal expansion. Therefore, eq. (3.85) becomes Dhk DTk 1 − β k Tk Dpk = c pk + (3.88) Dt Dt ρk Dt Substituting eq. (3.88) into eq. (3.84), the energy equation becomes DTk Dpk ′′′ = ∇ ⋅ (kk ∇Tk ) + Tk β k + qk + ∇Vk , rel : k , rel ρ k c pk (3.89) Dt Dt which can be simplified for different substances as demonstrated below. For ideal gases, substituting the equation of state for ideal gas [ ρ k = pk /( Rg Tk ) ] into eq. (3.87) yields β k = 1 / Tk , and eq. (3.89) becomes
βk =
196 Transport Phenomena in Multiphase Systems
DTk Dpk ′′′ = ∇ ⋅ (kk ∇Tk ) + + qk + ∇Vk , rel : k , rel (3.90) Dt Dt For incompressible fluids, the constant density yields β k = 0 and eq. (3.89) becomes DTk ′′′ = ∇ ⋅ (kk ∇Tk ) + qk + ∇Vk , rel : k , rel ρ k c pk (3.91) Dt For solids, the density is constant and the velocity is zero, and the energy equation becomes ∂T ρ s c ps s = ∇ ⋅ (ks ∇Ts ) + qs′′′ (3.92) ∂t which is a heat conduction equation.
ρ k c pk
Example 3.2: Determine the viscous dissipation for incompressible flow in the Cartesian coordinate system. The reference velocity can be assumed to be zero. Solution: Since the reference velocity is zero, the viscous dissipation is ∇Vk : k . According to eq. (C.33),
ª ∂u ∂u ∂u º « ∂x ∂y ∂z » « » « ∂v ∂v ∂v » ∇V = « » « ∂x ∂y ∂z » « ∂w ∂w ∂w » « » ¬ ∂x ∂y ∂z ¼ The shear stress tensor for incompressible flow is ª § ∂u ∂v · ∂u § ∂u ∂w · º μ ¨ + ¸ μ ¨ + ¸» « 2μ ∂x © ∂z ∂x ¹ » © ∂y ∂x ¹ « « § ∂v ∂u · § ∂v ∂w · » ∂v 2μ μ ¨ + ¸» = 2 μ Dk = « μ ¨ + ¸ k ∂y « © ∂x ∂y ¹ © ∂z ∂y ¹ » « » ∂w » « μ § ∂w + ∂u · μ § ∂w + ∂v · 2μ ¨ ¸ « ¨ ∂x ∂z ¸ » ∂z ¹ © ∂y ∂z ¹ ¬© ¼ Therefore, the viscous dissipation is
ª ∂u « ∂x « « ∂v =« « ∂x « ∂w « ¬ ∂x ∂u ∂y ∂v ∂y ∂w ∂y
(3.93)
∇Vk :
k
§ ∂u ∂v · § ∂u ∂w · º ∂u º ª 2μ ∂u μ ¨ + ¸ μ ¨ + ¸» »« ∂x © ∂z ∂x ¹ » © ∂y ∂x ¹ ∂z » « « § ∂v ∂u · § ∂v ∂w · » ∂v » ∂v 2μ μ ¨ + ¸» » : «μ¨ + ¸ ∂z » « © ∂x ∂y ¹ ∂y © ∂z ∂y ¹ » « » ∂w » « § ∂w ∂u · § ∂w ∂v · ∂w » » μ¨ + μ¨ + ¸ 2μ » ∂z ¼ « © ∂x ∂z ¸ ∂z ¹ © ∂y ∂z ¹ ¬ ¼
Chapter 3 Generalized Governing Equations: Local Instance Formulations 197
The contraction of two tensors of rank two yields 2 ª§ ∂u · 2 ∂u ∂v º ª§ ∂u · 2 ∂u ∂w º § ∂u · ∇Vk : k = 2 μ ¨ ¸ + μ «¨ ¸ + » + μ «¨ ¸ + » © ∂x ¹ « » «© ∂z ¹ ∂z ∂x » ¬ ¼ ¬© ∂y ¹ ∂y ∂x ¼ ª§ ∂v · 2 ∂v ∂u º ª§ ∂v · 2 ∂v ∂w º § ∂v · + μ «¨ ¸ + » + 2 μ ¨ ¸ + μ «¨ ¸ + » © ∂y ¹ «© ∂x ¹ ∂x ∂y » «© ∂z ¹ ∂z ∂y » ¬ ¼ ¬ ¼ 2 2 2 ª§ ∂w · ∂w ∂v º ª§ ∂w · ∂w ∂u º § ∂w · + μ «¨ » + 2μ ¨ » + μ «¨ ¸+ ¸+ ¸ ∂x ∂z » ∂y ∂z » © ∂z ¹ «© ∂y ¹ «© ∂x ¹ ¬ ¼ ¬ ¼
2
(3.94) The final form of the viscous dissipation for incompressible flow is obtained by rearranging eq. (3.94), i.e., 2 2 ª § ∂u · 2 § ∂v · § ∂w · ∇Vk : k = μ « 2 ¨ ¸ + 2 ¨ ¸ + 2 ¨ ¸ © ∂z ¹ « © ∂x ¹ © ∂y ¹ ¬ (3.95) 2 2 2 º § ∂v ∂u · § ∂w ∂v · § ∂u ∂w · +¨ + ¸ + ¨ + ¸ +¨ + ¸» © ∂x ∂y ¹ © ∂y ∂z ¹ © ∂z ∂x ¹ » ¼ The viscous dissipation in cylindrical and spherical coordinate systems can be found by similar means using information in Appendix C.
3.3.4 The Second Law of Thermodynamics
To obtain the differential form of the second law of thermodynamics, the surface integrals in eq. (3.22) can be rewritten as volume integrals: (3.96) ³ ρk (Vk ,rel ⋅ n)sk dA = ³ ∇ ⋅ ( ρk Vk ,rel sk )dV
A V
³
A
§ q′′ · q′′ k ⋅ n k dA = ³ ∇ ⋅ ¨ k ¸ dV V Tk © Tk ¹
(3.97)
Substituting eqs. (3.96) and (3.97) into eq. (3.22) and considering eq. (3.35) yields ª∂ § q′′ · qk º ′′′ k gen ³V « ∂t ( ρk sk ) + ∇ ⋅ ( ρk Vk ,rel sk ) + ∇ ⋅ ¨ Tk ¸ − Tk » dV = ³V s′′′ ,k dV ≥ 0 (3.98) « » ©¹ ¬ ¼ In order for eq. (3.98) to be true for any arbitrary control volume, the integrand in eq. (3.98) should always be positive, i.e., § q′′ · q′′′ ∂ (3.99) ( ρk sk ) + ∇ ⋅ ( ρ k Vk ,rel sk ) + ∇ ⋅ ¨ k ¸ − k = s′′′ ,k ≥ 0 gen ∂t © Tk ¹ Tk Equation (3.99) can be rewritten as
198 Transport Phenomena in Multiphase Systems
§ q′′ · q′′′ ª ∂ρ º ª ∂s º sk « k + ∇ ⋅ ( ρ k Vk , rel ) » + ρ k « k + Vk , rel ⋅ ∇sk » + ∇ ⋅ ¨ k ¸ − k = s′′′ , k ≥ 0 gen ¬ ∂t ¼ ¬ ∂t ¼ © Tk ¹ Tk (3.100) Considering the continuity equation, eq. (3.38), and definition of the substantial derivative, eq. (3.40), eq. (3.100) can be reduced to § q′′ · q′′′ Ds ρ k k + ∇ ⋅ ¨ k ¸ − k = s′′′ ,k ≥ 0 (3.101) gen Dt © Tk ¹ Tk where the three terms on the lefthand side represent rate of change of entropy per unit volume, rate of change of entropy per unit volume by heat transfer and internal heat generation, respectively. Equation (3.101) means that the entropy generation per unit volume must not be negative at any time or location. ′′′ For a multicomponent system without internal heat generation ( qk = 0 ), Curtiss and Bird (1999; 2001) obtained the entropy flux vector and the entropy generation as N 1§ · s′′ = ¨ q′′ − ¦ hk ,i J k ,i ¸ (3.102) k k T© i =1 ¹
N N§ Ng · cR T § · ′′′ Tk s′′′ ,k = − ¨ q′′ − ¦ hk ,i J k ,i ¸ ⋅ ∇ ln Tk − ¦ ¨ J k ,i ⋅ u d k ,i ¸ − : ∇V − ¦ k ,i mk ,i gen k ¨ ¸ ρ k ,i i =1 i =1 © i =1 M k ,i © ¹ ¹ (3.103) where dk,i is the diffusional driving force [see eq. (1.110) – (1.114)], g k ,i is partial molar Gibbs free energy, q′′ is total heat flux, obtained by eq. (3.82), k including conduction, interdiffusional convection, and Dufour effect.
3.3.5 Species
The surface integrals in eq. (3.32) may be converted to volume integrals by applying eq. (3.34) as follows: (3.104) ³ ρk ,i (Vk ,rel ⋅ nk )dA = ³ ∇ ⋅ ρk ,i Vk ,rel dV
A V
³
A
J k ,i ⋅ n k dA = ³ ∇ ⋅ J k ,i dV
V
(3.105)
Substituting eqs. (3.104) and (3.105) into eq. (3.32) and considering eq. (3.35), the entire lefthand side of eq. (3.32) is included in a single volume integral, i.e., § ∂ρ k ,i · (3.106) ³V ¨ ∂t + ∇ ⋅ ρk ,i Vk ,rel + ∇ ⋅ J k ,i − mk′′′,i ¸ dV = 0 © ¹ The only condition that makes eq. (3.106) true regardless of the shape and size of the control volume is that the integrand must equal zero, i.e., ∂ρ k ,i ′′′ (3.107) + ∇ ⋅ ρ k ,i Vk , rel = −∇ ⋅ J k ,i + mk ,i , i = 1, 2,", N ∂t
Chapter 3 Generalized Governing Equations: Local Instance Formulations 199
The first term on the lefthand side is the rate of increase of mass of the species i per unit volume, and the second term is net rate of additions of mass of the ith species per unit volume by convection. The terms on the righthand side are the net rate of mass of ith species per unit volume by diffusion, and rate of production of species i by chemical reaction. Equation (3.107) is the equation of conservation of mass for species. If the total number of species in the kth phase is N, a total of N1 independent equations for conservation of species mass can be obtained. After defining the mass fraction of species i in phase k as
ωk ,i =
ρ k ,i ρk
(3.108)
eq. (3.107) can be rewritten as ª ∂ωk ,i º ª ∂ρ k º ′′′ « ∂t + ∇ ⋅ ρ k Vk , rel » ωk ,i + ρ k « ∂t + ∇ ⋅ ωk ,i Vk , rel » = −∇ ⋅ J k ,i + mk ,i (3.109) ¬ ¼ ¬ ¼ According to eq. (3.38), the first bracket on the lefthand side of eq. (3.109) is zero. The second bracket on the lefthand side is the substantial derivative of the mass fraction. Therefore, the conservation of species mass in terms of mass fraction becomes Dωk ,i ′′′ (3.110) = −∇ ⋅ J k ,i + mk ,i ρk Dt Assuming binary syetem of A and B, one can use Fick’s law in eq. (3.110) to yield Dωk , A ′′′ (3.111) = − ρ∇ ⋅ ( DAB ∇ωk , A ) + mk , A ρk Dt which is useful in determing the diffusion in dilute liquid solution at constcant ′′′ temperature and pressure. Equation (3.111), with mk ,i = 0 , is similar to energy equation (3.91) with no internal heat source and viscous dissipation and therefore it is used for analogy between heat and mass transfer analysis. In the proceding discussion to develop eq. (3.107), the mass fraction and mass flux were used. The species equation can also be developed in term of molar concentration (or molar fraction) and molar flux. By following a similar procedure, the species equation is ∂ck ,i ′′′ (3.112) = −∇ ⋅ n′′,i + nk ,i k ∂t where the molar flux relative to the stationary coordinate axes can be obtained from eq. (1.97), i.e., n′′ = ci V * + J * (3.113) i i Substituting eq. (3.113) into eq. (3.112), we have ∂ck ,i ′′′ (3.114) + ∇ ⋅ (ck ,i V * ) = −∇ ⋅ J * ,i + nk ,i k ∂t
200 Transport Phenomena in Multiphase Systems
where the first term on the lefthand side is rate of increase of mole of the species i per unit volume, and the second term is net rate of additions of mole of the ith species per unit volume by convection. The terms on the righthand side are net rate of mole of ith species per unit volume by diffusion, and the molar rate of prodiction of species i by chemical reaction. For a binary system of components A and B with constant pressure, eq. (3.114) reduces to ∂ck , A ′′′ (3.115) + ∇ ⋅ (ck , A V * ) = −c∇ ⋅ ( DAB ∇x A ) + nk , A ∂t where c is the mixture molar concentration and V * is molaraveraged velocity defined in Chapter 1. Equation (3.115) can be applied to low density gases with constant temperature and pressure. ′′′ ′′′ The production rate of the ith species in the kth phase, mk ,i (or nk ,i ), is still unknown at this point, but it can be obtained by analyzing the chemical reaction. If the number of chemical reactions taking place in the system is Nc, the mass production rate is (Kleijn, 1991; Mahajan, 1996) ′′′ mk ,i = ¦ aij M i ℜ j
j =1 Nc
(3.116)
where aij is the stoichiometric coefficient that describes the proportions of the mole numbers of reactants disappearing and mole numbers of products appearing as a result of the reaction process (see Section 2.3.3). The net reaction rate of the jth chemical reaction ℜ j is the difference between the forward and backward reactions, i.e.,
ℜ j = ℜ j+ − ℜ j−
(3.117)
If chemical reactions take place in gas mixture, the forward and backward reaction rates are
Nr §p · xi ¸ ℜ j + = rj + ( p, T )∏ ¨ i =1 © Ru T ¹ Np aij
(3.118)
aij
§p · xi ¸ ℜ j − = rj − ( p, T )∏ ¨ (3.119) i =1 © Ru T ¹ where Nr and Np are number of reactants and products, respectively. The two reaction rate constants, rj + and rj − , depend on the specific chemical reaction
under consideration. These constants are related by
¦ aij 1 § RuT · i=1 = (3.120) ¨ ¸ rj + ( p, T ) K j (T ) © p 0 ¹ where p0 is the standard pressure and Kj(T) is the thermodynamic equilibrium constant of the jth chemical reaction:
rj − ( p, T )
N
Chapter 3 Generalized Governing Equations: Local Instance Formulations 201
ª −ΔG 0 (T ) º j (3.121) K j (T ) = exp « » « RuT » ¬ ¼ where ΔG 0 (T ) is the standard Gibbs energy for the jth chemical reaction; its j
value depends on the specific chemical reaction considered. For the case in which only one chemical reaction ( N c = 1 ) leads to production of the ith component, the molar production rate can be obtained by a simplified from ′′′ ′′′ n nk ,i = kn ck ,i (3.122) ′′′ where the index n represents the order of the reaction and kn is a rate constant n (1/s ) depends on the temperature. As before, the mass flux, J k ,i , in eqs. (3.107) and (3.110) includes mass fluxes due to ordinary diffusion driven by the concentration gradient, pressure diffusion, body force diffusion, and thermal (Soret) diffusion [see eq. (1.116)]. For a binary mixture, the mass flux can be calculated by eq. (1.143). Equation (3.110) is valid for the case that chemical reaction occurs in the entire volume – referred to as homogeneous reacting system. For the case where the chemical reaction takes place on a surface – referred to as heterogeneous reaction – the source term in eq. (3.110) will not appear, and the rate of production will be accounted as a boundary condition. More discussion about modeling of the heterogeneous reacting system can be found in discussion of Chemical Vapor Deposition (CVD) in Section 7.3.2.
Example 3.3 Couette flow is found in many engineering applications, such as lubricant flow in a journal bearing. Consider the laminar viscous Couette flow where two plates are of infinite extent, steady state conditions apply, and there is an applied pressure gradient in the xdirection, dp/dx (see Fig. 3.4). The fluid is Newtonian and incompressible. Starting from the basic equations and using mathematical, physical arguments to make suitable simplifying assumptions, derive the differential equations governing the fluid velocity and temperature. The fluid properties can be assumed to be constant, but viscosity dissipation should be accounted for. Solution: This is an incompressible fluid flow and heat transfer problem. If the reference frame is stationary, the continuity equation is ∂u ∂v ∂w ++ =0 (3.123) ∂x ∂y ∂z The momentum equations in the x, y and zdirections are § ∂ 2u ∂ 2u ∂ 2u · § ∂u ∂u ∂u ∂u · ∂p ρ ¨ + u + v + w ¸ = ρ X − + μ ¨ 2 + 2 + 2 ¸ (3.124) ∂x ∂y ∂z ¹ ∂x ∂y ∂z ¹ © ∂t © ∂x
202 Transport Phenomena in Multiphase Systems
U2 T2
Plate 2
L y, v z, w
T1
Plate 1
x, u
Figure 3.4 Couette Flow.
U1=0
ρ¨
ρ¨
§ ∂ 2v ∂ 2v ∂ 2v · § ∂v ∂v ∂v ∂v · ∂p + u + v + w ¸ = ρY − + μ ¨ 2 + 2 + 2 ¸ (3.125) ∂x ∂y ∂z ¹ ∂y ∂y ∂z ¹ © ∂t © ∂x
§ ∂2w ∂2w ∂2w · § ∂w ∂w ∂w ∂w · ∂p +u +v + w ¸ = ρZ − + μ¨ 2 + 2 + 2 ¸ ∂x ∂y ∂z ¹ ∂z ∂y ∂z ¹ © ∂t © ∂x
(3.126) Assuming that there is no internal heat generation and the fluid is incompressible, the energy equation is ª § ∂u · 2 § ∂ 2T ∂ 2T ∂ 2T · § ∂T ∂T ∂T ∂T · ρcp ¨ +u +v +w ¸ = k ¨ 2 + 2 + 2 ¸ + μ «2 ¨ ¸ ∂x ∂y ∂z ¹ ∂y ∂z ¹ © ∂t « © ∂x ¹ © ∂x ¬ 2 2 2 2 2 § ∂v · § ∂w · § ∂v ∂u · § ∂w ∂v · § ∂u ∂w · º + 2¨ ¸ + 2¨ + ¸ +¨ + ¸ +¨ + ¸ +¨ ¸ » (3.127) © ∂z ¹ © ∂x ∂y ¹ © ∂y ∂z ¹ © ∂z ∂x ¹ » © ∂y ¹ ¼ Since the upper plate moves along the xdirection at a constant velocity, the velocity component in the xdirection is not a function of x, ∂u =0 (3.128) ∂x and the velocity component in the zdirection is zero, i.e., w=0 (3.129) Substituting eqs. (3.128) and (3.129) into eq. (3.123) yields ∂v =0 (3.130) ∂y Since there is no blowing or suction on the bottom plate (v=0 at y=0), eq. (3.130) gives us v=0 (3.131)
Chapter 3 Generalized Governing Equations: Local Instance Formulations 203
Since the velocity and temperature of the plates do not vary along the zdirection, the velocity and the temperature do not vary in the zdirection, i.e., ∂u ∂v ∂w ∂T = = = =0 (3.132) ∂z ∂z ∂z ∂z Substituting eqs. (3.128) – (3.132) into eqs. (3.124) – (3.127), one obtains ∂u ∂p ∂ 2u ρ = ρX − + μ 2 (3.133) ∂t ∂x ∂y ∂p 0 = ρY − (3.134) ∂y ∂p 0 = ρZ − (3.135) ∂z
§ ∂ 2T ∂ 2T · § ∂u · ∂T · § ∂T ρcp ¨ +u ¸ = k¨ 2 + 2 ¸+ μ¨ ¸ ∂x ¹ ∂y ¹ © ∂t © ∂y ¹ © ∂x It is further assumed that there is no body force, X = Y apply the steadystate condition ∂u ∂T = =0 ∂t ∂t eqs. (3.134) – (3.135) are reduced to ∂p =0 ∂y ∂p =0 ∂z The governing equations of the problem become dp ∂ 2u =μ 2 dx ∂y § ∂ 2T ∂ 2T · § ∂u · ∂T ρ c pu = k¨ 2 + 2 ¸+ μ¨ ¸ ∂x ∂y ¹ © ∂y ¹ © ∂x
2 2
(3.136)
= Z = 0 , and
(3.137)
(3.138) (3.139)
(3.140) (3.141)
Example 3.4 Simplify the energy equation for a multicomponent system, eq. (3.83), for the case in which the Dufour effect is negligible but the interdiffusion term needs to be considered. It is assumed that the reference frame is stationary ( Vref = 0 ) and the viscous dissipation is
negligible.
Solution: With the Dufour effect and viscous dissipation neglected, eq. (3.83) becomes
204 Transport Phenomena in Multiphase Systems
¦ «ρ ¬
i =1
N
ª
k
N Dp D (ωi hk ,i )º = ∇ ⋅ (kk ∇Tk ) − ¦ ∇ ⋅ ( J k ,i hk ,i ) + Dtk » Dt ¼ i =1
(3.142)
or
ª º ¦ « ∂t (ω ρ h ) + ∇ ⋅ ( V ω ρ h )» ¬ ¼
N i k k ,i k i k k ,i i =1 N
∂
(3.143) Dpk = ∇ ⋅ (kk ∇Tk ) − ¦ ∇ ⋅ ( J k ,i hk ,i ) + Dt i =1 The diffusive mass flux can be rewritten in terms of the velocity of one component, Vk ,i compared to the massaveraged velocity, Vk [see eq. (1.93)] J k ,i = ωk ,i ρ k ( Vk ,i − Vk ) (3.144) The advection and diffusive term can be combined as ωk ,i ρ k Vk + J k ,i = ωk ,i ρ k Vk ,i (3.145)
The mass averaged velocity, Vk , is a function of the component velocities
Vk = ¦ ωi Vk ,i
i =1 N
(3.146)
Using the component velocity, Vk ,i , defined in eq. (3.145), in the energy equation, eq. (3.143), yields: N Dp ª∂ º (3.147) ¦ « ∂t (ωi ρk hk ,i ) + ∇ ⋅ (ωi ρk Vk ,i hk ,i )» = ∇ ⋅ (kk ∇Tk ) + Dtk ¼ i =1 ¬ For a stationary reference frame, the species conservation equation, (3.107), becomes ∂ ′′′ (ωk ,i ρ k ) + ∇ ⋅ (ωk ,i ρ k Vk ) = −∇ ⋅ J k ,i + mk ,i (3.148) ∂t Substituting eq. (3.145) into eq. (3.148), the conservation of species equation can be written in terms of the velocity of the ith component ∂ (3.149) (ωi ρk ) + ∇ ⋅ (ωi ρk Vk ,i ) = mk′′′,i ∂t For convenience, let’s define a substantial derivative in terms of the velocity of the ith component Dk ,i ∂ (3.150) ( ) = ( ) + Vk ,i ⋅ ∇ ( ) Dt ∂t Using the species equation, eq. (3.149), the energy equation (3.147) can be rewritten as N N Dhº ª Dp ª ′′′ º ωi ρ k k ,i k ,i » = ∇ ⋅ (kk ∇Tk ) + k − ¦ ¬ mk ,i hk ,i ¼ (3.151) ¦« Dt ¼ Dt i =1 i =1 ¬
Chapter 3 Generalized Governing Equations: Local Instance Formulations 205
The specific enthalpy of the ith component is a function of temperature and pressure, hk ,i = hk ,i (Tk , pk ) . The substantial derivative of enthalpy is [see eq. (3.88)] Dk ,i hk ,i D T 1 − β k ,iTk Dk ,i pk (3.152) = c p , ki k ,i k + Dt Dt ρ k ,i Dt where 1 § ∂ρ k ,i · (3.153) ¨ ¸ ρ k ,i © ∂Tk ¹ p is the coefficient of thermal expansion for the ith component in the kth phase. Substituting eq. (3.152) into eq. (3.151) yields: N Dk , T D ,i p º ª ¦ «ωi ρk c pk ,i Dit k + (1 − β k ,iTk ) kDt k » i =1 ¬ ¼ (3.154) N Dpk ′′′ ¼ = ∇ ⋅ (kk ∇Tk ) + − ¦ ª mk ,i hk ,i º Dt i =1 ¬ Since the massaveraged velocity is a function of the component velocity, a simplification can be made on one of the pressure terms. N N Dpk N Dk ,i pk −¦ = Vk ⋅ ∇pk − ¦ Vk ,i ⋅ ∇pk = −¦ (1 − ωi ) Vk ,i ⋅ ∇pk (3.155) Dt i =1 Dt i =1 i =1 Substituting this relation into the energy equation, eq. (3.154), yields N Dk ,iT º ª ¦ «ωi ρk c pk ,i Dt k » = ∇ ⋅ (kk ∇Tk ) i =1 ¬ ¼ (3.156) N N Dk ,i pk º N ª ′′′ ¼ −¦ (1 − ωi ) Vk ,i ⋅ ∇pk + ¦ « β k ,iTk » − ¦ ª mk ,i hk ,i º Dt ¼ i =1 ¬ i =1 i =1 ¬ where the second and third terms on the right hand side of the energy equation above are due to the difference in pressure work from the mass average velocity, and the pressure work from each component. The last term is due to species generation or consumption. For a single component system, the energy equation reduces to DTk Dpk ρ k c pk = ∇ ⋅ (kk ∇Tk ) + β k Tk (3.157) Dt Dt
β k ,i =
3.3.6 Jump Conditions at the Interfaces
The phase equations introduced above can be applied within each phase and up to an interface. However, they are not valid across the interface, where sharp changes in various properties occur. Appropriate boundary conditions at the interface must be specified in order to solve the governing equations for heat,
206 Transport Phenomena in Multiphase Systems
mass, and momentum transfer in the two adjoining phases. The interface conditions will serve as boundary conditions for the transport equations in the adjacent phases. These will be formulated based on the integral formulation for a control volume that includes two phases separated by an interface, as shown in Fig. 3.3(b). Jump conditions at the interface can also be obtained by applying the basic laws (conservation of mass, momentum, energy, and the second law of thermodynamics) at the interface. For the sake of simplicity, we will omit – for now – the effects of surface tension, disjoining pressure and interfacial thermal resistance; these subjects will be taken up in Chapter 5, when we address interfacial phenomena.
Conservation of Mass
The integral formulation for conservation of mass for a control volume that includes two phases is given by eq. (3.8), with Π = 2 . It can be rewritten for Fig. 3.3(b) by adding and subtracting the interface term as 2 ª∂ º ¦ « ∂t ³Vk (t ) ρk dV + ³Ak (t )+ AI (t ) ρk (Vk ,rel ⋅ n k )dA» ¼ k =1 ¬ (3.158) 2 ª º − ³ « ¦ ρ k (Vk , rel ⋅ n k ) » dA = 0 AI ( t ) ¬ k =1 ¼ In Fig. 3.3(b), the surface areas of V1 (t ) and V2 (t ) are A1 (t ) + AI (t ) and A2 (t ) + AI (t ) , respectively. Using eqs. (3.34) and (3.35), eq. (3.158) becomes
ª2 º § ∂ρ k · (3.159) ¨ ∂t + ∇ ⋅ ρ k Vk ,rel ¸ dV − ³AI ( t ) « ¦ ρ k (Vk ,rel ⋅ n k ) » dA = 0 Vk ( t ) © ¹ k =1 ¬ k =1 ¼ Substituting eq. (3.38) into eq. (3.159), one obtains ª2 º (3.160) ³AI (t ) «¦ ρk (Vk ,rel ⋅ n k ) » dA = 0 ¬ k =1 ¼ The only condition that ensures the general validity of eq. (3.160) is that the integrand equals zero, i.e., ρ1 (V1,rel ⋅ n1 ) + ρ 2 (V2,rel ⋅ n 2 ) = 0 (3.161)
¦³
2
If the reference frame velocity, Vref , equals the velocity of the interface,
VI , eq. (3.161) becomes ρ1 (V1 − VI ) ⋅ n1 + ρ 2 (V2 − VI ) ⋅ n 2 = 0 (3.162) th Introducing the mass flux of the k phase at the interface, ′′ mk = ρ k ( Vk − VI ) ⋅ n k (k = 1, 2) (3.163) the mass balance at the interface can be expressed as ′′ ′′ m1 + m2 = 0 (3.164) Since the thickness of the interface is considered to be zero, unit normals at an interfacial point in two adjacent phases are related by (see Fig. 3.5)
Chapter 3 Generalized Governing Equations: Local Instance Formulations 207
Figure 3.5 Unit normal vectors at interface.
n1 = −n 2 (3.165) Substituting eq. (3.165) into eq. (3.162) yields ρ1 (V1 − VI ) ⋅ n1 = ρ 2 (V2 − VI ) ⋅ n1 (3.166) If there is no phase change at the interface (a waterair twophase system, for example), the mass flow rate at the interface should be zero, i.e., ′′ ′′ m1 = − m2 = 0 and the mass balance equation at the interface simply reduces to V1 = V2 = VI (3.167)
Momentum Balance
The momentum integral equation for a control volume containing two phases separated by an interface is given by eq. (3.13) with Π = 2 . It can be rewritten as 2 ª º §N · ª2 º ′ ,rel ⋅ n k dA» − ³ « ¦ ′ , rel ⋅ n k » dA ρ k ,i X k ,i ¸ dV + ³ ¦ « ³Vk (t ) ¨ ¦ k k Ak ( t ) + AI ( t ) k =1 ¬ © i =1 ¹ ¼ ¼ AI ( t ) ¬ k =1
2 ª∂ º = ¦ « ³ ρ k Vk , rel dV + ³ ρ (Vk , rel ⋅ n k )Vk ,rel dA» Vk ( t ) Ak ( t ) + AI ( t ) k ¼ k =1 ¬ ∂t
(3.168) ª2 º − ³ « ¦ ρ k (Vk ,rel ⋅ n k ) Vk , rel » dA AI ( t ) ¬ k =1 ¼ The surface integrals over Ak (t ) + AI (t ) in eq. (3.168) can be converted to volume integrals over Vk (t ) using eq. (3.34):
N ∂ ª º ∇ ⋅ ′ + ¦ ρ k ,i X k ,i − ( ρ k Vk ) − ∇ ⋅ ρ k Vk ,rel Vk » dV ¦ ³Vk (t ) « k i =1 ∂t k =1 ¬ ¼
2
ª2 º = ³ « ¦ ª ′ ,rel ⋅ n k − ρ k (Vk ,rel ⋅ n k )Vk ,rel º » dA ¬k ¼ AI ( t ) ¬ k =1 ¼ Substituting eq. (3.52) into eq. (3.169), one obtains ª2 º ¬ ¼ ³AI (t ) «¦ ª ′k ,rel ⋅ nk − ρk (Vk ,rel ⋅ n k )Vk ,rel º » dA = 0 ¬ k =1 ¼
(3.169)
(3.170)
208 Transport Phenomena in Multiphase Systems
As stated above, the only condition that assures the general validity of eq. (3.170) is that the integrand equals zero, so ′ ′ (3.171) 1 ⋅ n1 + 2 ⋅ n 2 = ρ1 ( V rel ⋅ n1 ) V rel + ρ 2 ( V2, rel ⋅ n 2 ) V2, rel 1, 1, If the reference frame velocity equals the velocity of the interface, eq. (3.171) becomes ′ ′ 1 ⋅ n1 + 2 ⋅ n 2 = ρ1[( V − VI ) ⋅ n1 ]( V − VI ) + ρ 2 [( V2 − VI ) ⋅ n 2 ]( V2 − VI ) 1 1 (3.172) Substituting eqs. (3.163) and (3.165) into eq. (3.172), the momentum balance becomes (3.173) ( 1′ − ′2 ) ⋅ n1 = m1′′(V1 − V2 ) Substituting eq. (3.55) into eq. (3.173), one obtains ′′ ( p2 − p1 )n1 + ( 1 − 2 ) ⋅ n1 = m1 (V1 − V2 ) (3.174)
Energy Balance
The total energy (including mechanical and thermal energy) balance at the interface for a control volume that includes two phases is expressed by eq. (3.18) with Π = 2 . It can be rewritten as 2 2 § § V2 · V2 · d ¦ dt ³Vk (t ) ρk ¨ ek + k2,rel ¸ dV + ¦ ³Ak (t )+ AI (t ) ρk (Vk ,rel ⋅ n k ) ¨ ek + k2,rel ¸ dA ¨ ¸ ¨ ¸ k k =1 © ¹ © ¹ 2 ª2 § V ·º − ³ « ¦ ρ k (Vk , rel ⋅ n k ) ¨ ek + k , rel ¸ » dA ¨ AI ( t ) 2 ¸» « k =1 © ¹¼ ¬ 2 2 2 § · ′′′ = −¦ ³ q′′ ⋅ n k dA + ³ ¨ ¦ q′′ ⋅ n k ¸ dA + ¦ ³ qk dV k k Ak ( t ) + AI ( t ) AI ( t ) V (t ) k =1 k =1 k © k =1 ¹ 2 ª2 º +¦ ³ (n k ⋅ ′ , rel ) ⋅ Vk , rel dA − ³ « ¦ (n k ⋅ ′ , rel ) ⋅ Vk , rel » dA (3.175) k k Ak ( t ) + AI ( t ) AI ( t ) k =1 ¬ k =1 ¼ N 2 § · + ¦ ³ ¨ ¦ ρ k ,i X k ,i ¸ ⋅ Vk , rel dV Vk ( t ) k =1 © i =1 ¹ The surface integrals for Ak (t ) + AI (t ) in eq. (3.175) can be converted to the volume integral for Vk (t ) by using eq. (3.34), i.e.,
¦³
k =1
2
∂ ª § Vk2,rel ° ® « ρ k ¨ ek + ¨ Vk ( t ) ∂t 2 °«© ¯¬
·º ¸» ¸ ¹» ¼
ª § V2 +∇ ⋅ « ρ k Vk ,rel ¨ ek + k ,rel ¨ 2 « © ¬
½ ·º §N · ° ′′′ ¸ » + ∇ ⋅ q′′ − qk − ∇ ⋅ ( ′ k ⋅Vk , rel ) − ¨ ¦ ρ k ,i X k ,i ¸ ⋅ Vk , rel ¾ dV k ¸ © i =1 ¹ » ° ¹¼ ¿
Chapter 3 Generalized Governing Equations: Local Instance Formulations 209
º½ · ° ¸ + q′′ ⋅ n k − (n k ⋅ ′ , rel ) ⋅ Vk , rel » ¾ dA k k ¸ »° ¹ ¼¿ (3.176) Substituting eq. (3.71) into eq. (3.176), one obtains 2ª º½ § V2 · ° ° ρ k (Vk ,rel ⋅ n k ) ¨ ek + k ,rel ¸ + q′′ ⋅ n k − (n k ⋅ ′ ,rel ) ⋅ Vk , rel » ¾ dA = 0 k k ³AI (t ) ®¦ « ¨ ¸ 2¹ »° ° k =1 « © ¼¿ ¯¬ (3.177) In order to satisfy eq. (3.177) for any shape and size of the interface, the integrand in eq. (3.177) must be zero. Thus, the energy balance at the interface becomes 2ª º § Vk2, rel · (3.178) « ρ k (Vk ,rel ⋅ n k ) ¨ ek + ¸ + q′′ ⋅ n k − (n k ⋅ ′ , rel ) ⋅ Vk , rel » = 0 ¦ k k ¨ ¸ 2¹ k =1 « » © ¬ ¼ If the reference velocity is set as the velocity of the interface, eq. (3.178) becomes 2ª º § Vk2,rel · «[ ρ k ( Vk − VI ) ⋅ n k ] ¨ ek + ¸ + q′′ ⋅ n k − (n k ⋅ ′ , rel ) ⋅ Vk , rel » = 0 (3.179) ¦ k k ¨ 2¸ k =1 « » © ¹ ¬ ¼ Substituting eq. (3.163) into eq. (3.179), the energy balance at the interface becomes 2ª º § Vk2,rel · " (3.180) ¦ « mk ¨ ek + 2 ¸ + q′′ ⋅ nk − (nk ⋅ ′k ,rel ) ⋅ Vk ,rel » = 0 k ¸ ¨ k =1 « » ¹ ¬© ¼ where the three terms in the bracket represent the contributions of advection, heat transfer, and work done by normal and shear stresses.
The Second Law of Thermodynamics at an Interface
2ª § Vk2,rel ° = ³ ®¦ « ρ k (Vk ,rel ⋅ n k ) ¨ ek + ¨ AI ( t ) 2 °« © ¯ k =1 ¬
The entropy inequality for a control volume that includes two phases separated by an interface is eq. (3.23) with Π = 2 , which can be rewritten as 2 2 ∂ ª2 º ρ k sk dV + ¦ ³ ρ k (Vk , rel ⋅ n k ) sk dA − ³ « ¦ ρ k (Vk ,rel ⋅ n k ) sk » dA ¦ ∂t ³Vk (t ) Ak ( t ) + AI ( t ) AI ( t ) k =1 k =1 ¬ k =1 ¼ +¦ ³
k =1 2 2 2 § 2 q′′ ⋅ n · ′′′ qk q′′ ⋅ n k k dA − ³ ¨ ¦ k k ¸ dA − ¦ ³ dV Ak ( t ) AI ( t ) Vk ( t ) T Tk k =1 k © k =1 Tk ¹
= ¦³
k =1
Vk ( t ) gen , k
s′′′ dV + ³
AI ( t ) gen , I
s′′
dA ≥ 0
(3.181)
The surface integrals over Ak (t ) + AI (t ) in eq. (3.181) can be converted to volume integrals over Vk (t ) , i.e.,
210 Transport Phenomena in Multiphase Systems
¦³
k =1
2
ª∂ ′′′ § q′′ · qk º k « ( ρ k sk ) + ∇ ⋅ ( ρ k Vk , rel sk ) + ∇ ⋅ ¨ ¸ − » dV Vk ( t ) ∂t « © Tk ¹ Tk » ¬ ¼
+³
2ª q′′ ⋅ n k º ½ ° ° k ®¦ « − ρ k (Vk ,rel ⋅ n k ) sk − » ¾ dA AI ( t ) Tk ¼ ° ° k =1 ¬ ¯ ¿
2 k =1 Vk ( t ) gen , k
= ¦³
s′′′ dV + ³
AI ( t ) gen , I
s′′
dA ≥ 0
(3.182)
The entropy generations in V1 (t ) and V2 (t ) are greater than or equal to zero, as required by eq. (3.99). The entropy generation at the interface should also be greater than or equal to zero in order to satisfy eq. (3.182), i.e., 2ª q′′ ⋅ n º ½ ° ° ′′ (3.183) sgen , I dA = ³ ®¦ « − ρ k (Vk , rel ⋅ n k ) sk − k k » ¾ dA ≥ 0 ³AI (t ) AI ( t ) Tk ¼ ° ° k =1 ¬ ¯ ¿ The integrand in eq. (3.183) must be greater than or equal to zero in order to satisfy eq. (3.183), i.e., 2 ª q′′ ⋅ n º ′′ sgen , I = ¦ « − ρ k (Vk , rel ⋅ n k ) sk − k k » ≥ 0 (3.184) Tk ¼ k =1 ¬ If the reference frame moves with the interface, eq. (3.184) becomes 2 ª q′′ ⋅ n º ′′ sgen , I = ¦ « −[ ρ k (Vk − VI ) ⋅ n k ]sk − k k » ≥ 0 (3.185) Tk ¼ k =1 ¬ Substituting eq. (3.163) into eq. (3.185), the interfacial entropy generation becomes 2 ª q′′ ⋅ n º ′′ ′′ sgen , I = ¦ « −mk sk − k k » ≥ 0 (3.186) Tk ¼ k =1 ¬
Species
The conservation of species mass for a control volume that contains two phases separated by an interface is eq. (3.33) with Π = 2 . It can be rewritten as 2 2 ∂ ª2 º ¦ ∂t ³Vk (t ) ρk ,i dV + ¦ ³Ak (t )+ AI (t ) ρk ,i (Vk ,rel ⋅ n k )dA − ³AI (t ) «¦ ρk ,i (Vk ,rel ⋅ n k ) » dA k =1 k =1 ¬ k =1 ¼ = −¦ ³
k =1 2 Ak ( t ) + AI ( t )
( J k ,i ⋅ n k )dA + ³
AI ( t )
¦J
k =1
2
k ,i
⋅ n k dA + ¦ ³
k =1
2
Vk ( t )
′′′ mk ,i dV
(3.187)
The surface integrals over Ak (t ) + AI (t ) in eq. (3.187) can be converted to a volume integral over Vk (t ) by using eq. (3.34), i.e.,
¦³
k =1
2
V
§ ∂ρ k ,i · ′′′ + ∇ ⋅ ρ k ,i Vk , rel + ∇ ⋅ J k ,i − mk ,i ¸ dV ¨ © ∂t ¹
Chapter 3 Generalized Governing Equations: Local Instance Formulations 211
2 ½ (3.188) ®¦ ª ρ k ,i ( Vk ,rel ⋅ n k ) + J k ,i ⋅ n k º ¾ dA ¬ ¼ AI ( t ) ¯ k =1 ¿ Substituting eq. (3.107) into eq. (3.188), the conservation of species mass at the interface becomes 2 ½ ª º (3.189) ³AI (t ) ®¦ ¬ ρk ,i (Vk ,rel ⋅ n k ) + J k ,i ⋅ n k ¼ ¾ dA = 0 ¯ k =1 ¿ In order to satisfy eq. (3.189), the integrand in eq. (3.189) must be zero, i.e., =³
¦ ªρ ¬
k =1
2
k ,i
(Vk , rel ⋅ n k ) + J k ,i ⋅ n k º = 0 ¼
(3.190)
If the reference velocity is equal to the velocity of the interface, the conservation of mass of species at the interface becomes
¦ ªρ ¬
k =1
2
k ,i
(Vk − VI ) ⋅ n k + J k ,i ⋅ n k º = 0 ¼
(3.191)
Equation (3.191) can also be rewritten in terms of mass flux at the interface and mass fraction of the ith species, i.e.,
¦ ªω ¬
k =1
2
k ,i
′′ mk + J k ,i ⋅ n k º = 0 ¼
(3.192)
where the two terms in the bracket represent species transport due to phase change and diffusion across the phases due to gradients in concentration, temperature, etc.
Supplementary Conditions at an Interface
The jump conditions specified above must be supplemented by additional boundary conditions because they cannot provide sufficient matching conditions at the interface to uniquely define the problem. The supplementary conditions are similar to the constitutive equations for the governing equations. These interfacial constitutive laws satisfy the restrictions imposed by the entropy inequality. In order to obtain more straightforward jump conditions at the interface, eq. (3.184) can be combined with the equation of energy balance at the interface – which is derived from eq. (3.180) – to obtain (Delhaye, 1974; Hetsroni, 1982) 2ª º §1 1· 1 ′′ s′′ , I = ¦ «(q′′ ⋅ n k + mk sk Tk ) ¨ − ¸ − ( k ⋅ n k ) ⋅ (Vkt − V t ) » gen k k =1 « » © TI Tk ¹ TI ¬ ¼ (3.193) 2 º m′′ ª 1 1 + ¦ k « g k + (Vk − VP ) 2 − ( k ⋅ n k ) ⋅ n k » 2 ρk k =1 TI ¬ ¼
212 Transport Phenomena in Multiphase Systems
where Vp is the velocity of an interfacial fluid particle. The normal component of Vp is the interfacial displacement velocity, VI , and the tangential component is V t . Thus
Vp = ( VI ⋅ n k )n k − V t
(3.194)
For most cases, the discontinuity of temperatures, tangential velocities, and chemical potentials at the interface has a relatively insignificant effect on the system. As a limiting case, it is assumed that the interfacial transfers are reversible; therefore, the entropy generation of the interface s′′ , I is equal to gen zero regardless of the values of mass flux, viscous stress tensor, and heat fluxes. To satisfy zero entropy generation at the interface, the following conditions must be satisfied: T1 = T2 = TI (3.195)
V1t = V2t = V t
2 2
(3.196)
V − VI V − VI ½ 1 ½ 1 ° ° g1 − g 2 = ® 2 −1 ¾ − ® ( 2 ⋅ n 2 )n 2 − ( 1 ⋅ n1 )n1 ¾ (3.197) ρ1 2 2 ¿ ° ° ¯ ρ2 ¯ ¿ These conditions are in addition to the boundary conditions required by conservations of mass, momentum, and energy [see eqs. (3.166), (3.173), and (3.180)]. Equation (3.197) means that the difference in chemical potential compensates for the mechanical effects of the relative kinetic energy difference and of the normal stress. Substituting eqs. (3.163) and (3.164) into eq. (3.197), the phase change boundary condition becomes ½ 1 1 1 ½ 1 1 ′′ g1 − g 2 = mk ® 2 − 2 ¾ − ® ( 2 ⋅ n 2 )n 2 − ( 1 ⋅ n1 )n1 ¾ (3.198) ρ1 ¿ ¯ ρ 2 ρ1 2 ¯ ρ2 ¿ The boundary condition at the interface specified by eq. (3.195) is identical to the thermodynamic equilibrium condition at the interface specified by eq. ′′ (2.146). This means that the phase change process ( mk ≠ 0 ), a nonequilibrium process, does not alter the temperatures of the two phases or the interface. Thermodynamic equilibrium at the interface also requires that the chemical potentials for two phases in equilibrium be equal [see eq. (2.148)]. For a singlecomponent system, the chemical potential is the specific Gibbs energy. Therefore, the Gibbs free energies of two phases at thermodynamic equilibrium ′′ should be equal. Since the phase change process ( mk ≠ 0 ) is not an equilibrium process, the specific Gibbs free energies of two phases are no longer equal, as indicated by eq. (3.198). Example 3.5 A container with an open top is exposed to superheated water vapor at temperature Tv, with the bottom and the side walls insulated. At time t=0, the bottom surface of the vessel comes into contact with the air at temperature T∞, which is below the vapor temperature Tv. The sides and top of the vessel remain insulated for t > 0 [see Fig. 3.6(a) and (b)].
Chapter 3 Generalized Governing Equations: Local Instance Formulations 213
Formulate the governing equations for the condensation problem. If the vapor is saturated (Tv = Tsat), find the instantaneous thickness of the condensate. Solution: This is a onedimensional condensation problem because the sides of the vessel are insulated. The thickness of the bottom wall is very small and therefore the temperature drop across the bottom wall can be neglected. In addition, it is assumed that the thermal properties of both liquid and vapor are constants. The interfacial velocity is related to the instantaneous liquid thickness by dδ wI = (3.199) dt The continuity equation at the interface can be obtained by using eq. (3.162) in the zdirection, i.e., ρA ( wA − wI ) = ρ v ( wv − wI )
(a) Insulated (t < 0) .
(b) Bottom cooled ( t > 0 ) Figure 3.6 Thinwalled vessel containing twophase mixture.
214 Transport Phenomena in Multiphase Systems
The condensate is stationary, so wA = 0 . Therefore, the vapor phase velocity is §ρ · dδ (3.200) wv = ¨ A − 1¸ © ρv ¹ dt For the liquid phase, the governing energy equation and the corresponding initial and boundary conditions are ∂TA ∂ 2T (3.201) = α A 2A ∂t ∂z TA ( z,0) = Tsat (3.202) ∂T kA A = h[TA (0, t ) − T∞ ] (3.203) ∂z For the vapor phase, the governing energy equation and the corresponding initial and boundary conditions are · d δ ∂Tv ∂Tv § ρA ∂ 2T + ¨ − 1¸ = α v 2v (3.204) ∂t © ρv ∂z ¹ dt ∂z Tv ( z ,0) = Tv (3.205) Tv (∞, t ) = Tv (3.206) The liquid and vapor temperatures at the interface must be the same, a condition required by eq. (3.195), i.e., Tv (δ , t ) = TA (δ , t ) = Tsat (3.207) The energy balance at the interface can be obtained by using eq. (3.180). Since the velocity in the liquid is zero and the velocity in the vapor phase is not very high, the kinetic energy and the viscous dissipation in eq. (3.180) can be neglected. The resulting energy balance equation at the interface is ∂T (δ , t ) ∂T (δ , t ) dδ − kv v + kA A = ρA hAv (3.208) ∂z ∂z dt The analytical solution of this problem will be very difficult. If the vapor is saturated, the governing equations for the vapor phase will not be necessary because the vapor temperature will be uniformly equal to the saturation temperature. This is a onedomain problem because only the liquid domain needs to be investigated. The problem can be further simplified by introducing the quasisteady state assumption. This step eliminates the transient term in the energy equation for the liquid phase because condensation is a very slow process and the change of internal energy in the liquid can be neglected. The energy equation in the liquid phase is thus reduced to ∂ 2T 0 = α A 2A (3.209) ∂z which is subject to the boundary conditions specified by eqs. (3.203) and (3.207). The energy balance at the interface is reduced to
Chapter 3 Generalized Governing Equations: Local Instance Formulations 215
∂TA (δ , t ) dδ (3.210) = ρA hAv ∂z dt The temperature in the liquid film can be obtained by integrating eq. (3.209) twice: T ( z ) = c1 z + c2 (3.211) where c1 and c2 are integral constants. They can be determined using eqs. (3.202) and (3.203), i.e., kA c1 = h ( c2 − T∞ ) (3.212) kA Tsat = c1δ + c2 (3.213) which can be solved to yield h(Tsat − T∞ ) h(Tsat − T∞ ) c1 = c2 = Tsat − δ (t ) (3.214) kA + hδ (t ) kA + hδ (t ) Substituting eq. (3.214) into eq. (3.211) yields the temperature distribution in the liquid film h(Tsat − T∞ ) TA ( z , t ) = Tsat + [ z − δ (t )] (3.215) kA + hδ (t ) The temperature distribution in the vessel can be illustrated by Fig. 3.7. At this point, eq. (3.215) can be substituted into eq. (3.210). An ordinary differential equation for the liquid thickness is obtained as a result: d δ kA h(Tsat − T∞ ) = (3.216) ρA hAv dt kA + hδ (t ) which is subject to the following initial condition: δ (0) = 0 (3.217) Rearranging eq. (3.216) so that δ appears only on the lefthand side, and t appears only on the righthand side, gives
Figure 3.7 Temperature distribution for a twophase mixture in a thinwalled vessel.
216 Transport Phenomena in Multiphase Systems
ªh º h(Tsat − T∞ ) dt (3.218) « δ + 1» d δ = ρA hAv ¬ kA ¼ Integrating eq. (3.218) with the boundary condition, eq. (3.217) yields h(Tsat − T∞ ) h2 t =0 δ +δ − (3.219) 2 kA ρ A hAv which is a quadratic equation in the form of ax 2 + bx + c = 0 . This allows for one positive and one negative solution of δ , and the latter does not physically make sense. Therefore, the liquid film thickness is ª 2h 2 (Tsat − T∞ )t º 2 hδ (t ) = −1 + «1 + » kA kA ρA hAv ¬ ¼
1
(3.220)
Example 3.6 In a countercurrent condenser shown in Fig. 3.8, the liquid flows downward with a mass flow rate of mA , while the vapor flows upward ′′ with a mass flow rate of mv . The heat flux at the external wall is qw . Derive the jump conditions using mass and energy balances at the liquidvapor interface. Solution: The mass balance at the interface can be obtained by analyzing the mass balances of the liquid and vapor phases in the control volume shown by the dashed line. The mass balances for the liquid and vapor phases are, respectively dm ′′ mA ( z ) + mδ dAδ = mA ( z ) + A dz (3.221) dz Vapor z
mA ( z ) ′′ qw mA ( z + dz ) mv ( z )
z
′′ mδ m v ( z + dz ) ′′ qw
z+dz
Liquid
Figure 3.8 Countercurrent condensation.
Liquid
Chapter 3 Generalized Governing Equations: Local Instance Formulations 217
dmv dz (3.222) dz where dA is the area of the liquid vapor interface. Equations (3.221) and (3.222) can be simplified to yield dm dAδ dm dAδ ′′ mδ = A (3.223) =− v dz dz dz dz ′′ which means that the rate of condensation mδ is related to the increase of the liquid mass flow rate and decrease of the vapor mass flow rate. Similarly, the energy balances for the liquid and vapor phases are respectively d ′′ ′′ ′′ mA hA + mδ hAδ dAδ + qAδ dAδ = mA hA + (mA hA )dz + qw dAw (3.224) dz d ′′ ′′ (3.225) mv hv + mδ hvδ dAδ + qvδ dAδ = mv hv + (mv hv )dz dz ′′ ′′ where qAδ and qvδ are conduction heat flux at the interface in the liquid and vapor phases, respectively. hvδ and hAδ are the enthalpy for saturated vapor and liquid, respectively. Equations (3.224) and (3.225) can be simplified as d ′′ ′′ ′′ (3.226) mδ hAδ dAδ + qAδ dAδ = (mA hA )dz + qw dAw dz d ′′ ′′ (3.227) mδ hvδ dAδ + qvδ dAδ = ( mv hv )dz dz Subtracting eq. (3.226) from eq. (3.227), one obtains d d ′′ ′′ ′′ ′′ mδ ( hvδ − hAδ )dAδ + (qvδ − qAδ )dAδ = (mv hv )dz − (mA hA )dz − qw dAw dz dz (3.228) The overall energy balance for the control volume including both phases is d d ′′ mA hA + mv hv + (mv hv ) = mA hA + (mA hA )dz + mv hv + qw dAw (3.229) dz dz i.e., d d ′′ (mv hv )dz = (mA hA )dz + qw dAw (3.230) dz dz Substituting eq. (3.230) into eq. (3.228), one obtains ′′ ′′ ′′ (3.231) −qvδ + qAδ = mδ hAv where hAv = hvδ − hAδ is latent heat of vaporization. Equation (3.231) is consistent with eq. (3.180) when kinetic energy and viscous dissipation at the interface in eq. (3.180) are neglected. ′′ mv ( z ) − mδ dAδ = mv ( z ) +
218 Transport Phenomena in Multiphase Systems
3.3.7 Classification of PDEs and Boundary Conditions
A general transport equation, whether it is mass, momentum, energy, or species, can be written as: ∂ (3.232) ( ρ k Φ k ) + ∇ ⋅ ( ρ k Vk Φ k ) = ∇ ⋅ ( Γ k ∇Φ k ) + F ( x, t , Φ k ,…) ∂t where the dependent variable, Φ , is one for mass, any component of velocity ( uk , vk , wk ) , enthalpy hk , or the species mass fraction ωk ,i . The gradient operator is the partial derivative of the dependent variable with respect to all the spatial directions x for a given coordinate system. A partial differential equation (PDE) is an equation of a function and its partial derivatives. In general, a PDE is classified by its linearity or nonlinearity, and by its order. Its order is considered by its highest derivative. A general second order PDE for two independent variables, η and ζ , is ∂ 2Φ ∂ 2Φ ∂ 2Φ ∂Φ ∂Φ (3.233) + 2B +C 2 + D +E + FΦ + G = 0 ∂η 2 ∂η∂ζ ∂ζ ∂η ∂ζ This equation is linear if all the coefficients (A, B, C, D, E, F and G) are a function of η and ζ , a constant, or zero. If they are a function of Φ or any of its derivatives, then the PDE is nonlinear. A partial differential equation is called quasilinear if it is linear in the highest derivatives. In eq. (3.233), this means that A, B and C are a function of η and ζ , a constant, or zero. If a differential equation is quasilinear, it can be classified as an elliptic ( AC − B 2 > 0 ), parabolic ( AC − B 2 = 0 ) or hyperbolic ( AC − B 2 < 0 ) equation. The classification of a PDE describes how disturbances or changes propagate through a domain. If you are interested in a single point on the domain, you need to know two things: what region in the domain affects that point, and, if you do something at that point, what region in the domain it will affect. What region affects a point and what region that point affects are called the zone of dependence and the zone of influence, respectively. A representation of the zone of dependence and the zone of influence of elliptic, parabolic, and hyperbolic equations are displayed in Fig. 3.9. A multidimensional PDE that varies with time may be classified by looking at only the relation of two independent variables at a time. For instance look at how the PDE would be classified only in the x and y directions. Then examine the PDE in x and t. Finally examine the PDE in y and t. The zone of dependence is where the zone of dependence for the three cases intersects. The zone of influence is where the zone of influence for the three cases intersects. An example of an elliptic equation is the Laplace equation. This equation governs physical problems, such as twodimensional steadystate heat conduction, electrical potential, the free stream characteristics in boundary layer equations, and the pressure field in a porous medium with Darcy’s assumption. The independent variables, η and ζ , are spatial coordinates in these problems, A
Chapter 3 Generalized Governing Equations: Local Instance Formulations 219
Elliptic
Parabolic
Hyperbolic
η
η
η
ζ
(a)
ζ η
ζ
η
η
ζ η
T2 T1 T2 T2 (b)
ζ η η
ζ
T2
ζ
(c)
ζ
Flow across a cylinder
Isotherms in steadystate heat conduction
ζ
Wake formed by a boat
Figure 3.9 Propagation of disturbances for different types of PDEs displayed as the (a) Zone of Dependence, the (b) Zone of influence, and (c) a physical example
and Φ is the potential field, representing temperature, voltage, the stream function, or pressure, respectively, for the given cases. Some examples of parabolic problems are unsteady heat conduction, boundary layer equations for momentum, energy and species as well as small vibrations in an elastic beam. For the unsteady heat equation, η is a spatial coordinate, and ζ is time. For the boundary layer equations, η is a spatial coordinate normal to the free stream velocity, and ζ is a spatial coordinate parallel to the free stream. A hyperbolic function is often called the wave equation. The most common examples of the wave equation are the propagation of shock waves and sound
220 Transport Phenomena in Multiphase Systems
Neumann [ ∂Φ / ∂η = f (ζ ) , or ∂Φ / ∂ζ = f (η ) ], and the mixed [ a (∂Φ / ∂η ) + bΦ = f (ζ ) , or a(∂Φ / ∂ζ ) + bΦ = f (η ) ] type. The number of boundary conditions that must be applied for each independent variable is equal to the highest order of that variable. For example, if the PDE is second order in both η and ζ , then two boundary conditions are needed for each independent variable for a total of four boundary conditions. If ζ is time and first order, than an initial Dirichlettype boundary condition is needed, and if it is second order, then both a Dirichlet and Neumann boundary condition are needed. There are implications to experimental measurements and numerical analysis based on the classification of the governing PDE. For example, experimental measurements in an incompressible flow field with moderate to low Reynolds numbers are very difficult, because the governing equations are highly elliptic. The elliptic nature means that disturbances downstream greatly affect the upstream flow field. Therefore, any measurements that require a device in the flow field may change the nature of that flow field, and are inherently inaccurate. Numerical simulations are very reliable for these cases when the flow is in the laminar regime. In the laminar flow regime, the full NavierStokes equations can be directly solved as an elliptic problem with no approximations. In the compressible flow regime with high Reynolds numbers, the characteristic of the flow is parabolic or hyperbolic, depending on the Mach number. In these cases, disturbances produced downstream do not affect the upstream flow field; therefore measurement devices in the flow will give an accurate depiction of the flow field without the device. In numerical simulations, it is computationally efficient to solve a parabolic flow field, because disturbances propagate in one direction. However, a flow field that is truly parabolic usually has a high Reynolds number and therefore is turbulent, which means it is threedimensional and transient. Turbulence modeling either involves a very fine mesh, which is very computationally expensive, or an averaging technique, which requires a closure problem. The classification of a particular problem will help reduce the time it takes to get reliable results.
waves. The interface between phases is also hyperbolic. Other less intuitive examples of the use of a hyperbolic function is heat conduction in micro and nanoscales. It becomes relevant in these situations because heat actually conducts at a finite rate, and the limits of these rates are relevant at such small scales. For any problem to be well defined, there are boundary/initial conditions that must be applied. There are three basic boundary conditions for second order PDEs. These boundary conditions are the Dirichlet [ Φ = f (η , ζ ) ], the
Chapter 3 Generalized Governing Equations: Local Instance Formulations 221
3.3.8 Rarefied Vapor SelfDiffusion Model
The discussions so far are limited to the case where the density of the fluid is sufficiently high to permit the continuum assumption. There are cases in which the assumption of continuum is not valid ( Kn ≥ 0.01 , see Section 1.4.1). For example, in the early stage of heat pipe startup from the frozen state, the vapor pressure and density are very small in the heat pipe core. Because of the low density, the vapor in the rarefied state is somewhat different from the conventional continuum state. Also, the vapor density gradient is very large along the axial direction of the heat pipe. The vapor flow along the axial direction is mainly caused by the density gradient via vapor molecular diffusion. The lowdensity vapor state that has partly lost its continuum characteristics is referred to as rarefied vapor. Neglecting the presence of noncondensable gases, the rarefied vapor flow can be simulated by a selfdiffusion model. The term selfdiffusion here means the interdiffusion of particles of the same mass due to a gradient in density. The governing equations for startup of a heat pipe from the frozen state were derived by applying the principles of the conservation of mass and energy in a differential cylindrical control volume in conjunction with the definition of mass flux (Cao and Faghri, 1993). The mass selfdiffusion equation is ∂ρ ∂ § ∂ρ · 1 ∂ § ∂ρ · (3.234) − ¨ Dv ¸− ¨ rDv ¸=0 ∂t ∂z © ∂z ¹ r ∂z © ∂r ¹ and the energy equation is ∂ ( ρ cvT ) 1 ∂ ∂ 1∂§ ∂T · ∂ § ∂T · ′′ (mr′′c p rT ) + (mz c pT ) = + ¨ rkv ¸ + ¨ Dv ¸ ∂t ∂r ∂r ¹ ∂z © ∂z ¹ r ∂r r ∂z © (3.235) ′′ where the mass fluxes mr and m′′ are z ∂ρ mr′′ = ρ v = − Dv (3.236) ∂r ∂ρ m′′ = ρ w = − Dv (3.237) z ∂z where Dv is the selfdiffusion coefficient, and kv is the vapor molecular conductivity. The evaluation of lowdensity properties such as Dv is carried out using the kinetic theory of gases. The coefficient of selfdiffusion is obtained from the relation based on the ChapmanEnskog kinetic theory (Hirschfelder et al., 1966): T 3 / MW (3.238) pσ 2 Ω D (kbT / ε ) where p is the pressure in atmospheres, σ is the collision diameter in Å, ε is the maximum energy of attraction between a pair of molecules, kb is the Boltzmann constant, and Ω D is the collision integral for mass diffusion. For Dv = 2.628 × 10−7
222 Transport Phenomena in Multiphase Systems
sodium, both σ and ε / kb can be found from the table for constants of the LennardJones potential model (Edwards et al., 1979), which gives σ = 3.567 Å, and ε / kb = 1375 K. The value of the collision integral Ω D can also be found from the same reference, which is listed as a function of ε / kb . The results of numerical simulation of heat pipe startup from the frozen state will be discussed in Section 6.9.2.
3.3.9 An Extension: Combustion
Combustion is an exothermic chemical reaction process between fuel and oxidant. If combustion involves a liquid fuel, the liquid fuel does not actually burn as a liquid; it is vaporized first and diffuses away from the liquidvapor surface. Meanwhile, the gaseous oxidant diffuses toward the liquidvapor interface. Under the right conditions, the mass fluxes of vapor fuel and gaseous oxidant meet and the chemical reaction occurs at a certain location known as the flame (Lock, 1994; Avedisian, 1997, 2000). The flame is usually a very thin region with a color dictated by the temperature of the combustion. The temperature and mass concentration distributions during a combustion process can be represented as shown in Fig. 3.10. The initial temperature of a liquid fuel is T0, and the temperature of the liquid–vapor interface, TI, is at the dew point temperature of the fuel. The temperature reaches a maximum at the location of the flame, and decreases with increasing x until it reaches T. The
Liquidvapor interface
T,ω
T∞
ωo
TI
ωf
ωp
T0
0
sI
sf
x
Figure 3.10 Combustion near a planar surface.
Chapter 3 Generalized Governing Equations: Local Instance Formulations 223
mass fraction of fuel, f, is maximal at the liquidvapor interface and decreases as location of the flame is approached. The mass fraction of the oxidant, o, on the other hand, is maximal at infinity and decreases as location of the flame is approached. The mass fractions of both vapor fuel and oxidant reach their minimum values at the location of the flame. The mass fractions of the products of combustion, p, are at their maximum at the location of the flame and decrease as x either increases or decreases. The transport phenomena involved in combustion include heat and mass transfer in both the liquid fuel and the gaseous mixture. The governing equations for the gaseous phase will be discussed below. Firstly, the gaseous mixture must satisfy the continuity equation, i.e., Dρ + ρ∇ ⋅ V = 0 (3.239) Dt and the gas flow is governed by the momentum equation, DV ρ = ∇ ⋅ + ρg (3.240) Dt where the shear stress can be obtained from eq. (3.57) and gravity is the the only body force considered. The energy equation for combustion is Dh ρ = ∇ ⋅ k ∇T (3.241) Dt which is a simplified version of eq. (3.79), with internal heat generation and viscous dissipation neglected. The effect of the substantial derivative of pressure on the energy is also neglected. The specific enthalpy for the mixture is related to the specific enthalpy of each component in the mixture by h = ¦ ωi hi
i =1 N
(3.242)
It is a common practice in thermodynamic and heat transfer analyses to consider the change in enthalpy during a chemical reaction process rather than the absolute values of enthalpy. For a process that does not involve chemical reaction, we can choose any reference state for an individual substance and define the enthalpy at that reference state as zero. For example, the reference state for water is often chosen at the triple point, while the reference state for an ideal gas is often chosen as zero K. However, when chemical reaction is involved in a process, as is the case with combustion, the composition of the system changes during the chemical reaction; therefore, the reference states for all reactants and products must be the same. One convenient option in such a situation is to segregate the enthalpy of any component into two parts: (1) the enthalpy due to its chemical composition at the standard reference state (at 25 ÛC and 1 atm), and (2) the sensible enthalpy due to any temperature deviation from the standard reference state. Therefore, the enthalpy for the ith component in the mixture can be expressed as hi = hiD + hiT (3.243)
224 Transport Phenomena in Multiphase Systems
Table 3.1 Enthalpy of formation ho for selected substances at 25 ÛC and 1 atm Substance Acetylene Benzene Carbon Carbon monoxide Carbon dioxide Ethyl alcohol Ethyl alcohol Ethylene Ethane Hydrogen Methane Nitrogen Nitrogen nDodecane nOctane nOctane Oxygen Propane Water Water vapor Formula (phase) C2H2 (g) C6H6 (g) C (s) CO (g) CO2 (g) C2H5OH (g) C2H5OH ( A ) C2H4 (g) C2H6 (g) H2 (g) CH4 (g) N2 (g) N (g) C12H26 ( A ) C8H18 (g) C8H18 ( A ) O2 (g) C3H8 (g) H2O ( A ) H2O (g) Enthalpy of formation (kJ/kmol) 226,730 82,930 0 110,530 393,520 235,310 277,690 52,280 84,680 0 74,850 0 472,650 291,010 208,450 249,950 0 103,850 241,820 285,830
where hiD is the enthalpy of formation for the ith component, i.e., the enthalpy due to its chemical composition at the standard reference state. The enthalpy of formation for selected substances is shown in Table. 3.1. The sensible enthalpy, hiT , is related to temperature by hiT = ³ D c pi dT
T N T
(3.244)
Substituting eqs. (3.243) and (3.244) into eq. (3.242), the enthalpy of the mixture becomes h = ¦ ωi hiD + ³ D c p dT
T i =1 T
(3.245)
where cp is the average specific heat of the mixture, defined as c p = ¦ ωi c pi
i =1 N
(3.246)
Substituting eq. (3.245) into eq. (3.241), the energy equation for combustion becomes N D (c pT ) Dωi ρ = ∇ ⋅ k ∇T − ρ ¦ hiD (3.247) Dt Dt i =1 If the fuel is consumed at a rate of m′′′ per unit volume, hc – the heat of f combustion – is defined by hc = −
¦h m′′′
f i =1
ρ
N
i
Dωi Dt
(3.248)
Chapter 3 Generalized Governing Equations: Local Instance Formulations 225
Substituting eq. (3.243) into eq. (3.248) gives ρ N D Dωi ρ N T Dωi hc = − ¦ hi Dt − m′′′ ¦ hi Dt m′′′ i =1 f f i =1
(3.249)
The contribution of the second term on the righthand side of eq. (3.249) is negligible, since hc hiT . Dropping the second term from the righthand side of eq. (3.249) and substituting the result into eq. (3.242), the energy equation becomes DT ρcp (3.250) = ∇ ⋅ k ∇T + m′′′hc f Dt The mass fraction of each component (fuel, oxidant, and product) is dominated by Dωi ρ = −∇ ⋅ J i + mi′′′ (3.251) Dt where the subscript i can be f (fuel), o (oxidant), or p (product). The ratio of the rates of oxygen and fuel consumption is defined as the oxygen/fuel ratio: m′′′ γ= o (3.252) m′′′ f The above analysis applies to combustion occurring on a planar surface. For many applications, combustion of the liquid fuel is usually preceded by breaking up a fuel jet into liquid droplets so that combustion occurs around a spherical liquid droplet. Combustion of a falling liquid droplet, as shown in Fig. 3.11, is analyzed here. The liquid fuel droplet vaporizes at the dew point, Td, which is the saturation temperature corresponding to the partial pressure of the fuel vapor in the mixture near the liquidvapor interface. To simplify the analysis, it is assumed that the temperature in the liquid fuel droplet is uniformly equal to the saturation temperature of the fuel at the total system pressure, i.e., the mass fraction of the fuel at the liquidvapor interface equals one. It is further assumed that the shapes of both the liquid fuel droplet and the flame are spherical, which allows for application of a onedimensional symmetric model (Lock, 1994), which is presented here. If the combustion process is assumed to be in a quasisteady state (neglecting the transient term in the governing equation), the energy equation for combustion becomes 1d 2 1 d§ dT · (3.253) (r ρ c p uT ) = 2 ¨ r 2 k f ¸ + m′′′hc r 2 dr r dr © dr ¹ where u is the velocity of the gaseous mixture in the radial direction. The mass fraction of the oxidant can be obtained by simplifying eq. (3.251), i.e., d ωo º 1d 2 1 dª2 ′′′ (3.254) (r ρ uωo ) = 2 2 « r ρ Do dr » − mo r dr r dr ¬ ¼
226 Transport Phenomena in Multiphase Systems
T, ω
Tf
T∞ Tsat
ωf
ωp
ω0
rI
rf
r
Figure 3.11 Combustion near a liquid fuel droplet.
where Do is the mass diffusivity of the oxidant in the gaseous mixture. Multiplying eq. (3.254) by hc / γ and combining the resulting equation with eq. (3.253) yields § h · º 1 d ª 2 dT d § h ·º 1 dª2 rk + r 2 ρ Do ¨ ωo c ¸ » (3.255) « r ρ c p u ¨ T + ωo c ¸ » = 2 ¨ 2 ¸ » r dr « r dr ¬ dr dr © γ ¹ ¼ γ c p ¹¼ « ¬ © where eq. (3.252) was used to eliminate the rates of fuel and oxygen consumption. Assuming that the Lewis number [ Le = k /( ρ c p Do ) ] equals one and the specific heat is constant, eq. (3.255) can be simplified as follows: 1d 2 1 d§ dT * · r ρ c p uT * ) = 2 ¨ r 2 k ( ¸ r 2 dr r dr © dr ¹ where T * = T + ωo hc γ cp (3.257) (3.256)
is a modified temperature in the combustion process. Equation (3.256) can be rearranged to the following form: d ª 2§ dT * · º r ¨ ρ c p uT * − k (3.258) « ¸» = 0 dr ¬ © dr ¹ ¼ The continuity equation requires that 4π r 2 ρ u = 4π rI2 ( ρ u ) I = m (3.259) Since it has been assumed that the mass fraction of the fuel in the mixture at the liquidvapor interface equals one, the mass flow rate m reflects the mass
Chapter 3 Generalized Governing Equations: Local Instance Formulations 227
flow rate of the fuel. The mass flux is often used in combustion analysis, and it is defined as m m′′ = (3.260) 4π r 2 Equations (3.258) and (3.259) can be rewritten in terms of mass flux, i.e., d ª 2§ dT * · º (3.261) r ¨ m′′c pT * − k « ¸» = 0 dr ¬ © dr ¹ ¼ 4π r 2 m′′ = 4π rI2 mI′′ = m (3.262) Integrating eq. (3.261) from the liquid fuel droplet surface ( r = rI ) to an arbitrary radius ( r > rI ) yields § dT * · 2 § dT * · r 2 ¨ m′′c pT * − k = rI ¨ m′′c pT * − k ¸ ¸ dr ¹ dr ¹ I © © Substituting eq. (3.262) into eq. (3.263), one obtains § dT * · dT * r 2k = rI2 mI′′c p (T * − TI* ) + rI2 ¨ k ¸ dr © dr ¹ r = r
I
(3.263)
(3.264)
The mass flow rate at the surface of the liquid fuel droplet, mI′′ , is the same as the fuel burning rate, m′′ , because the mass fraction at the surface of the f droplet is one. Introducing the excess modified temperature, θ = T * − TI* eq. (3.264) becomes ª dθ 1 § dθ · º r 2k = rI2 m′′ c p «θ + f ¨k ¸» dr m′′ c p © dr ¹ r = rI » « f ¬ ¼ Introducing a new dependent variable, 1 § dθ · ϕ =θ + ¨k ¸ m′′ c p © dr ¹ r = rI f eq. (3.266) becomes rI2 m′′ c p 1 1 f dϕ = dr ϕ k r2 Integrating eq. (3.268) over the interval of (r , ∞) , one obtains ln
2 ϕ∞ rI m′′ c p 1 f = ϕ k r
(3.265)
(3.266)
(3.267)
(3.268)
(3.269)
The fuel burning rate, Gf, in eq. (3.269) is related to the heat transfer at the liquid droplet surface as follows § dθ · (3.270) m′′ hAv = qI ′′ = k ¨ f ¸ © dr ¹ r = rI
228 Transport Phenomena in Multiphase Systems
Substituting eq. (3.270) into eq. (3.267), the new dependent variable ϕ becomes h ϕ = θ + Av (3.271) cp Substituting eq. (3.271) into eq. (3.269) and letting r = rI , an equation for the fuel burning rate is obtained: § c pθ ∞ · k ln ¨1 + (3.272) m′′ = ¸ f rI c p © hAv ¹ Equation (3.272) can also be used to determine the transient liquid fuel droplet size, because the fuel burning rate is related to the size of the droplet by d ( ρ f 4π rI3 / 3) / dt dr (3.273) = −ρ f I m′′ = − f 2 4π rI dt where ρ f is the density of the liquid fuel. Substituting eq. (3.273) into eq. (3.272), one obtains § c pθ ∞ · dr k ln ¨1 + rI I = − ¸ ρ f cp © dt hAv ¹
(3.274)
Integrating eq. (3.273) and considering the initial condition rI = ri at t = 0, the liquid droplet radius becomes 2kt § c pθ ∞ · (3.275) ln ¨ 1 + rI2 = ri 2 − ¸ ρ f cp © hAv ¹ Equation (3.275) can be used to estimate the time needed to completely burn the liquid fuel droplet: ρ f c p ri 2 (3.276) tf = § c pθ ∞ · 2k ln ¨ 1 + ¸ hAv ¹ © Example 3.7 Liquid hydrocarbon fuel is suspended in air in the form of small droplet, and the oxygen/fuel ratio is γ = 3.5 . Before the liquid fuel is sprayed into the air, the mass fraction of oxygen is ωo ,∞ = 0.21 , and the temperature of the air is T∞ = 25 o C . The latent heat of evaporation for the fuel is hAv = 360 kJ/kg , and the enthalpy of combustion is hc = 4.4 × 104 kJ/kg . The specific heat and thermal conductivity of the gaseous mixture are c p = 1.15 kJ/kgK and k = 0.07 W/mK,
respectively. The initial temperature of the liquid droplet can be assumed to be equal to the saturation temperature of the fuel, 187 °C. The density
Chapter 3 Generalized Governing Equations: Local Instance Formulations 229
of the liquid fuel is ρ f = 700 kg/m3 . Estimate the time needed to completely burn a liquid fuel droplet with a diameter of 40 m.
Solution: The excess modified temperature at the locations far away from the fuel droplet can be obtained by substituting eq. (3.257) into eq. (3.265): h * θ ∞ = T∞ − TI* = (T∞ − TI ) + (ωo,∞ − ωo, I ) c γ cp
4.4 × 104 = 2133.65 o C 3.5 × 1.15 The time required to completely burn the liquid fuel droplet is then obtained by using eq. (3.276), i.e., ρ f c p ri 2 700 × 1.15 × 103 × (20 × 10−6 ) 2 tf = = = 1.12ms § c pθ ∞ · § 1.15 × 2133.65 · 2k ln ¨1 + ¸ ¸ 2 × 0.07 × ln ¨1 + 360 © ¹ hAv ¹ © = (25 − 187) + (0.21 − 0) ×
References
Avedisian, C.T., 1997, “Soot Formation in Spherically Symmetric Droplet Combustion,” Physical and Chemical Aspects of Combustion, edited by Irvin Glassman, I., Dryer, F.L., and Sawyer, R. F., pp. 135160, Gordon and Breach Publishers. Avedisian, C.T., 2000, “Recent Advances in Soot Formation from Spherical Droplet Flames at Atmospheric Pressure,” Journal of Propulsation and Power, Vol. 16, pp. 628656. Bejan, A., 2004, Convection Heat Transfer, 3rd ed., John Wiley & Sons, New York. Cao, Y., and Faghri, A., 1993, “Simulation of the Early Startup Period of High Temperature Heat Pipes From the Frozen State by a Rarefied Vapor SelfDiffusion Model,” ASME Journal of Heat Transfer, Vol. 115, pp. 239246. Curtiss, C.F., and Bird, R.B., 1999, “Multicomponent Diffusion,” Industrial and Engineering Chemistry Research, Vol. 38, pp. 21152522. Curtiss, C.F., and Bird, R.B., 2001, “Errata,” Industrial and Engineering Chemistry Research, Vol. 40, p. 1791. Delhaye, J.M., 1974, “Jump Conditions and Entropy Sources in Two Phase Systems. Local Instant Formulation,” International Journal of Multiphase Flow, Vol. 1, pp. 395409.
230 Transport Phenomena in Multiphase Systems
Edwards, D.K., Denny, V.E., and Mills, A.F., 1979, Transfer Process, Hemisphere, New York. Hetsroni, G., 1982, Handbook of Multiphase Systems, Hemisphere Publishing Co., Washington, DC. Hirschfelder, J.O., Curtiss, C.F., and Bird, R.B., 1966, Molecular Theory of Gases and Liquids, Wiley, New York. Incropera, F.P., and DeWitt, D.P., 2001, Fundamentals of Heat and Mass Transfer, 5th ed., John Wiley & Sons, New York. Kays, W.M., Crawford, M.E., and Weigand, B., 2004, Convective Heat Transfer, 4th ed., McGrawHill, New York, NY. Kleijn, C.R., 1991, “A Mathematical Model of the Hydrodynamics and Gas Phase Reaction in Silicon LPCVD in a Single Wafer Reactor,” Journal of Electrochemical Society, Vol. 138, pp. 21902200. Lock, G.S.H., 1994, Latent Heat Transfer, Oxford Science Publications, Oxford University, Oxford, UK. Mahajan, R.L., 1996, “Transport Phenomena in Chemical VaporDeposition Systems,” Advances in Heat Transfer, Academic Press, San Diego, CA. Welty, J.R., 1978, Engineering Heat Transfer, John Wiley & Sons, New York. White, F.M., 1991, Viscous Fluid Flow, 2nd ed., McGrawHill, New York.
Problems
3.1. The local instance continuity equation can also be obtained by performing a mass balance for the differential control volume shown in Fig. P3.1. Derive the continuity equation and compare your result with eq. (3.41). 3.2. Derive the local instance continuity equation in a cylindrical coordinate system using the control volume shown in Fig. P3.2. 3.3. The continuity equation for incompressible flow is eq. (3.42). Is it valid for transient flow? Why or why not? 3.4. In order to take water without stopping the train, a narrow trough of several thousand feet long can be placed in the midway between the rails of a railroad track. As a locomotive passes over the narrow trough, the water can be forced up into the tender through a scoop by the speed of the locomotive (see Fig. P33). Assuming the crosssectional area of the scoop is A, find the force acting on the train by the water using momentum balance in (a) a coordinate system that is attached to and moves with the locomotive, and (b) a fixed coordinate system.
Chapter 3 Generalized Governing Equations: Local Instance Formulations 231
Figure P3.1
Figure P3.2
Locomotive tender
Train Speed, U
trough
Figure P3.3
3.5. Determine viscous dissipation for compressible flow in the Cartesian coordinate system. 3.6. Write down the governing equations for twodimensional, compressible flow with the effect of viscosity being accounted for. Assume that the constant properties assumptions apply and that the only body force is gravity. 3.7. Fourier’s law of conduction [see eq. (1.60)] assumes that heat flux responds immediately to a temperature gradient. However, for cases with a very high temperature gradient or extremely short duration, heat is propagated at a finite speed, which can be accounted for by a modified heat flux model:
q′′ + τ ∂q′′ = − k ∇T ∂t
where τ is the thermal relaxation time. Derive the energy equation for a pure substance using the above model and eq. (3.79). Internal heat generation and viscous dissipation can be neglected.
232 Transport Phenomena in Multiphase Systems
3.8. The generalized mass balance on an interface is given in eq. (3.166). Rewrite it in a threedimensional Cartesian coordinate system. 3.9. The enthalpy for the ith component in the kth phase is hk ,i is function of temperature, Tk, and pressure, pk. Show that the substantial derivative in terms of the velocity of the ith component, Vk,i, can be obtained from eq. (3.152). 3.10. Formulate the hydrodynamics for developing flow in a duct formed by two infinite flat plates. Assume the flow to be laminar, incompressible, 2D and steady. It can also be assumed that the inlet velocity at the entrance of the duct is uniform, as shown in Fig. P3.4. Using the integral method, find out the hydrodynamic entrance and the skin friction coefficient in terms of Reynold’s number. Is the analysis valid for a circular pipe?
u
Inlet Boundary layer region D y ,v
0
x, u Hydrodynamic entrance length
Figure P3.4
Fully developed flow
3.11. Formulate the hydrodynamics when flow is fully developed in a duct formed by two infinite flat plates and a round pipe. Assume the flow to be laminar, incompressible, twodimensional and steady. Because flow is fully developed, it can be assumed that the velocity in the duct does not change in the xdirection. Determine the skin friction coefficient. Clearly describe the methods used to solve the momentum and continuity equations. 3.12. Obtain the Nusselt number for constant wall temperature and constant wall heat flux cases for fully developed flow in terms of temperature and velocity profiles in a circular pipe. Assume laminar, incompressible, twodimensional and steady state. 3.13. Develop analytical expressions for the Nusselt number for constant wall temperature and constant heat flux cases in a circular pipe with fully developed velocity profile but developing temperature profile. Assume
Chapter 3 Generalized Governing Equations: Local Instance Formulations 233
laminar, incompressible, twodimensional and steady state with uniform inlet velocity. 3.14. Consider steady state, twodimensional, fully developed thermal and hydrodynamic laminar flow with constant properties in a circular pipe. Let there be heat transfer to or from the fluid at a constant rate per unit of tube length. Determine the Nusselt number if the effect of frictional heating due to viscosity (viscous dissipation) is included in the analysis. How does frictional heating affect the Nusselt number? Identify the pertinent nondimensional parameters appearing in this problem. U0
L
y x Figure P3.5 3.15. Develop expressions for the temperature and velocity profiles for the Couette flow configuration shown in the Fig. P3.5. The lower plate is stationary while the upper plate is drawn with a constant velocity of u0. Viscous dissipation is important in the problem. The upper and lower plates are held at constant temperatures. Assume the flow to be laminar, twodimensional and steady. 3.16. Develop formulation for the case of heating with no evaporation for the flow configuration shown below over a rotating disk (Fig. P3.6) Assume that the flow is laminar, incompressible, steady, axisymmetric and the disk has a constant wall temperature. The flow is introduced from a central collar that directs the liquid radially outward with velocity u0 over a gap height of h0. is the liquid film thickness and r0 is the initial radius of the disk at the outlet of the collar, h0 is the collar height and T is the thickness of the thermal boundary layer. Obtain Nusselt number Nu = hr0 / k , and film thickness δ / h0 as a function of the Reynolds number and Roosby
2 numbers [Re = u0 r0 /ν and Ro = u0 /(ω 2 r02 )] using the integral method.
234 Transport Phenomena in Multiphase Systems
Liquid free surface
u0
h0
liquid
Figure P3.6 Schematics of the thin film over a rotating disk. 3.17. Equation (3.166) is valid for any interface between two phases. Rewrite eq. (3.166) so that it is applicable to the melting problem. Suppose the liquid phase is the first phase and solid phase is the second phase. 3.18. The momentum balance equation (3.173) is given in vector/tensor notation. For a twodimensional interface, the velocity of kth phase, Vk, can be decomposed into normal and tangential components, i.e., Vk = n kVkn + t kVkt , where nk and tk are unit vectors in the normal and tangential directions, respectively. Write momentum balances in normal and tangential directions for a twodimensional interface. Assume both phases are Newtonian fluid. 3.19. The jump conditions at an interface can also be obtained by applying the various conservation laws to a control volume that includes the interface. Apply the conservation of mass principle to the control volume shown in Fig. P3.7, and show that when the thickness of the control volume Δx goes to zero, the conservation of mass is ρ1 (u1 − uI ) = ρ 2 (u2 − u I ) . 3.20. For the control volume shown in Fig. P3.7, show that the momentum balance at the interface is p1 − p2 = ρ 2 (u2 − u I )u2 − ρ1 (u1 − uI )u1 when the thickness of the control volume Δx goes to zero. 3.21. A tube with radius Ro and zero thickness passes through a liquid phase change material at its melting point, Tm (see Fig. P3.8). A cold fluid with a temperature of Tc, which is below Tm, flows through the inside of the tube. Solidification will occur on the outer surface of the tube. Formulate the problem by giving the governing equations and corresponding boundary conditions in the cylindrical coordinate system. If the problem can be assumed to be quasisteady state (see Example 3.5), solve for the instantaneous location of the solidliquid interface.
Chapter 3 Generalized Governing Equations: Local Instance Formulations 235
Liquid R R0 r Tc Solid
Figure P3.7 Figure P3.8
3.22. The mass and energy balance in a countercurrent condenser as shown in Fig. 3.8 was analyzed in Example 3.6. Redo the problem for the case of cocurrent condensation, i.e., both liquid and vapor flow downward. 3.23. Redo example 3.6 for the case of countercurrent evaporation and discuss the effect of the heat transfer direction on the energy balance at the interface. 3.24. The time it takes to completely burn a fuel droplet of specified size as estimated by eq. (3.276) is valid for conductioncontrolled combustion processes. For a liquid droplet that has sufficient inertia, the heat transfer between vapor and liquid is dominated by convection, i.e., eq. (3.270) is replaced by m′′ hAv = qI ′′ = hθ ∞ , where h is the heat transfer coefficient f between the gaseous mixture and liquid droplet, and can be obtained by NuD = hD / k = C Re1/ 2 . What is the time it takes to completely burn the D liquid fuel droplet? Compare your result with eq. (3.276). 3.25. Redo Example 3.7 using the result from Problem 3.24 and compare your result with the result in Example 3.7. The relative velocity of the fuel droplet is V p = 5m/s and the heat transfer correlation is Nu D = 0.6 Re1 / 2 . D The viscosity of the gaseous mixture is 4 × 10−5 m 2 /s . 3.26. A fuel droplet with an initial diameter of Dp moves at a velocity vp in an oxidant gas. The gravitational acceleration vector is g and the local oxidant gas velocity vector is vg. The fuel droplet evaporates as it moves in the combustor. If the drag coefficient of the fuel droplet is CD, write the momentum equation for the fuel droplet.
236 Transport Phenomena in Multiphase Systems
3.27. If the temperature of the fuel droplet in Problem 3.26 can be assumed to be uniform at any time, and the latent heat of evaporation of the fuel is hsA , derive the energy equation for the fuel droplet. 3.28. A miniature tube with diameter d and length 2L is bent into a Ushaped tube with the two ends sealed (see Fig. P3.9). Two evaporator sections are near the two closed ends, and each of them has a length of Lh. The condenser section, with a length of 2Lc, is located at the bottom of the Ushaped tube. The wall temperatures at the heating and cooling sections are Te and Tc, respectively. A liquid slug with length 2Lp is located at the bottom of the Ushaped tube. The location of the liquid slug is represented by displacement, xp, which is zero when the liquid slug is exactly in the middle of the Ushaped miniature tube. Suppose the initial value of xp0 is greater than zero, part of the vapor plug in the left side of the heat pipe is in contact with the condenser section, and condensation occurring in the left vapor plug will result in a decrease in the pressure of the left vapor plug, pv1. On the other hand, part of the right evaporator is in contact with the liquid slug and boiling may occur at the contact area of the right evaporator and the liquid slug, which causes increasing vapor pressure of the right vapor plug, pv2. The liquid plug will be pushed back to the left side due to the pressure difference between the two vapor plugs, Δp = pv1 − pv 2 . When xp becomes zero, there is no evaporation or condensation in either vapor plug, but the liquid slug keeps moving due to its inertia. When part of the liquid slug enters the left evaporator, part of the right vapor plug will be in contact with the condenser. At this point, the boiling in the left vapor plug and the condensation in the right plug will change the sign of p, and this will result in the liquid slug moving to the right side. The oscillation of the liquid slug can be maintained by alternative boiling and condensation in the two vapor plugs. Assuming the pressure drop in the bend is negligible, specify the governing equations that describe the oscillatory flow in this Ushaped tube.
Figure P3.9
Chapter 3 Generalized Governing Equations: Local Instance Formulations 237