Chapter 2 :

Conservation Laws of Fluid Motion and

Boundary Conditions

Ibrahim Sezai
Department of Mechanical Engineering
Eastern Mediterranean University

Fall 2014-2015

Governing Equations of Fluid Flow and Heat Transfer

The governing equations of fluid flow represent

mathematical statements of the conservation laws
of physics.

The mass of fluid is conserved

The rate of change of momentum equals the sum of
the forces on a fluid particle (Newton’s second law)
The rate of change of energy is equal to the sum of
the rate of heat addition to and the rate of work done
on a fluid particle (first law of thermodynamics).

ME555 : Computational Fluid Dynamics 2 I. Sezai – Eastern Mediterranean University

1
The six faces are labelled
N, S, E, W, T, B

The center of the element is

located at position (x, y, z)

ρ = ρ (x, y, z, t) p=p
(x, y, z, t)
T = T (x, y, z, t) u=u
(x, y, z, t) Fig. 2-1

Fluid element for

c
o
n
s
e
r
v
a
t
i
o
n
Mass Conservation in Three Dimensions
l
a
 Rate of increase  Net rate of flow  w
 of mass in
 
 of mass into 
fluid element fluid element  s
    .
Fluid properties at faces are approximated by means of the two terms
of the Taylor series
The pressure at the W and E faces, can be expressed as

t ( x ytz)   
p(  x 
 y z)  u  (u)  xp

 y z2  u  (u) 1 x  y z (left)
 p  x1  x
1

p xand
 x
2
1

x 2  x 2

ME555 : Computational Fluid Dynamics 3 I. Sezai – Eastern Mediterranean University

  
 v  (v)
(v) 
1
 y  x z  v 
1
 y  x z
(front)




 
 w  (
(w) 
w) 1
z  x y    w 
1
z  x y
(bottom)





 2


Rates of change following a fluid particle and for a fluid element

The total or substantial derivative of with respect to time following a

fluid particle is

D     dx   dy
   
  dz Dt t  x dt
 y dt  z dt
A fluid particle follows the flow, so
dx / dt  u
dy / dt 
v dz / dt
w

Hence the
substantive
derivative of is
given by

D     
 u v w 
 u grad Dt t x y z
t
ME555 : Computational Fluid Dynamics 5 I. Sezai – Eastern Mediterranean University

D/Dt defines the rate of change of property per unit mass.

The rate of change of property per unit volume for a fluid particle is
ρD/Dt, hence
       u  grad 
D 

(2-8)
Dt  t 
 
ρ = mass per unit volume.
Lhs of the mass conservation equation (2-4) is

 div( u)
t
The generalization of these terms for an arbitrary conserved
property is
 ( )
 div(  u)

(2-9)
t
 Net rate of flow of 
 Rate of increase  
of  per unit volume   out of fluid element 
   per unit volume 

 
ME555 : Computational Fluid Dynamics 6 I. Sezai – Eastern Mediterranean University

3
Rewriting eq. (2-9) 0 (due to conservation
of mass)   
 ( )  
 div(  u)     u  grad      div( u)  
D
(2-10)
t   t     t  Dt

 Rate of increase  Net rate of flow of    Rate of increase

of  of   out of   of  for a 

 fluid element  
  fluid element 
  fluid particle 

     
of  for momentum
Relevant entries x-momentum u  Du
and energy
 (u)equations as defined in Eqn (2.10)
 div( uu)
are: Dt t
y-momentum v  Dv  (v)  div( vu)
Dt t
z-momentum w  Dw  (w)  div(
Dt wu)
t
energy E  DE  (E)  div(
Dt Eu)
ME555 : Computational Fluid Dynamics 7 t I. Sezai – Eastern Mediterranean University

The rates of increase of x-, y-, and z-momentum per unit volume are
Du Dv Dw
  
Dt Dt Dt
We distinguish two types of forces on fluid particles:
• surface forces - pressure forces
- viscous forces
• body forces - gravity
forces
- centrifugal
forces

source terms
- Coriolis
forces
- electromagne
tic force
The pressure, a normal stress, is denoted by p.
Viscous stresses are denoted by τ.
ME555 : Computational Fluid Dynamics 8 I. Sezai – Eastern Mediterranean University 4
 Fig. 2-3 Stress components on three
faces of fluid element.
The suffices i and j in τij indicate that
the stress component acts in the j-
direction on a surface normal to the i-
direction.

First we consider the x-

components of the forces
due to pressure p and stress
components τxx, τyx and τzx
shown in Fig. 2-4.

Fig. 2-4 Stress components in

the x-direction.

ME555 : Computational Fluid Dynamics 9 I. Sezai – Eastern Mediterranean University

The net force in the x-direction is the sum of

the force components acting in that direction
on the fluid element.
On the pair of faces (E, W) we have
 p  
x   p  xx 1
x 
y z
   
 p  
1 x p   xx 1
x
 y z
    
 xx  p
    x y z
(2-12a)

x x 
 
The net force in the x-direction on
the pair of faces (N, S) is
   
 
 yx 1  1  y x
yx
 yxz  
  y  x z   x y
2 z  yx   yx (2-12b)
y y y
  

The net force in the x-direction on
the pair of faces (T, B) is
  zx 1
 

 zx 1   z  x

 zxy  
  z  x y   x
y2 
 zzx   (2-12c)
2
z z
zx
z 5
  

ME555 : Computational Fluid Dynamics 10 I. Sezai – Eastern Mediterranean University
 The total force per unit volume on the fluid due to these surface
stresses is equal to the sum of (2-12a), (2-12b), (2-12c) divided by the
volume δxδyδz:
( p  xx )   yx  zx
 

(2.13)
x y z
To find x-component of the momentum equation:
 Rate of change of   Total force in x-direction   Total force in x-direction

 x-momentum of    on the element due to
 fluid particle 
   on the element due to 
surface   body forces
stresses

      

Eqn.( 2.11) Eqn.( 2.13) SMx

Du  ( p  xx )   yx  zx  S
 
(2.14a)
Dt x y z Mx

SMx = Body force on the element per unit volume in x-direction

SMz = –ρg (body force due to gravity per unit volume in z-direction)
ME555 : Computational Fluid Dynamics 11 I. Sezai – Eastern Mediterranean University

Similarly, y-component of the momentum equation is

  (p yy ) 
 Dv  xy
  zy
S
(2.14b)
Dt x y z My

and, z-component of the momentum equation is

Dw    (p zz )
  xz
 yz
 S
(2.14c)
Dt x y z Mz

ME555 : Computational Fluid Dynamics 12 I. Sezai – Eastern Mediterranean University

6
 Energy Equation in Three Dimensions
The energy equation is derived from the first law of thermodynamics
which states that
 Net rate of   Net rate of
 Rate of increase
 to 

work 
of  heat added  done on
 energy of fluid particle
  fluid pacirtel   fluid particle 

 

DE

Dt

Work Done by Surface Forces = (Fsurface forces)(V)

V = velocity component in the direction of the force.

ME555 : Computational Fluid Dynamics 13 I. Sezai – Eastern Mediterranean University

The surface forces given by (2.12a-c) all act in the x-direction.

The net rate of work done by these forces acting in x-direction is

u( p   xx )  (uyx ) (uzx ) 

   x y z

(2.16a)

x y z 
 

Similarly, work done by surface stresses in y and z-direction are

 (v xy )    v( p  yy ) 
(v zy ) 
  (2.16b)  x  y z
 x y z 
 

 xz(w)   (wyz ) w(p  zz )

   x y z

(2.16c) 7
 
x y z
 
Summing (2.16a-c) yields the total rate of work done on the fluid
particle by surface stresses:

(uxx ) (u yx )  (uzx ) (v xy ) (v yy )

div(pu)     
(v zy )
x y z x
y z
(w xz )  (wyz ) (w zz )
  
x y z

where

(up) (vp) (wp)

div( pu)   
x y z

ME555 : Computational Fluid Dynamics 15 I. Sezai – Eastern Mediterranean University

Energy Flux due to Heat Conduction

The heat flux vector has three
components qx, qy, qz

The net rate of heat transfer to

the CV due to heat flow in x-
direction is
 q
qq  x 1 x    q

qx
x 1
x y  z  
 x  y z
(2.18b-c)

  

ME555 : Computational Fluid Dynamics 16 I. Sezai – Eastern Mediterranean University

8
Similarly, the net rates of heat transfer to the fluid due to heat flows in
the y- and z-direction are
q q (2.18b-c)
– y x y z and  z x y z
y z
The net rate of heat added to CV per unit volume is the sum of (2.18a-
c) divided by δxδyδz
 qx  qy  qz
–  div q  (2.19)
x y z
T T T
q  k q  k q  k
x y z
x y z
This can be written in vector form as
q  k grad T

(2.20)
Combining (2.19) and (2.20) yields the rate of heat addition to the
CV due to heat conduction
div q  div(k grad T )
ME555 : Computational Fluid Dynamics 17 I. Sezai – Eastern Mediterranean University

Energy Equation
sum of the net rate of work done on the CV
b y su rfa ce stre sse s(2 .17 )
                 
  (uxx )  (uyx )  (uzx )  (vxy ) 
 div( pu)    
DE
   x y z x 
Dt  (v ) (v )  (w) (w )
 (w) (2.22)

 y z x y z


 div(k gradT )  SE

net rate of heat addition rate of increase of energy

to the fluid (2.21) due to sources
E  i  1 (u2  v2  w2 )
2      
kinetic energy
i  internal (thermal) energy
SE = source of energy per unit volume per unit time (i.e. effects of
gravitational potential energy changes)
ME555 : Computational Fluid Dynamics 18 I. Sezai – Eastern Mediterranean University

9
Multiplying
the x-momentum equation (2.14a) by u
the y-momentum equation (2.14a) by v
the z-momentum equation (2.14a) by w

and adding the results together

D  1 (u 2 2v 2  w2       
  u  grad pu xx

 yx
 zx

Dt  x y z 

  xy  yy  zy 
v   
x y z 


  xz  yz
w  
 zz 
u  S

(2.23)

x
y 
z
M



ME555 : Computational Fluid Dynamics 19 I. Sezai – Eastern Mediterranean University
Subtracting (2.23) from (2.22)
Di u u
   p div u  div(k grad T )    (2.24)
xx
Dt  x yx  y
u v
  
v v w
 
zx xy yy zy xz
z x y z x
w w
yz  zz i S
y z

where Si =SE – u.SM

For an incompressible fluid  i = cT and div u = 0 (c = specific
heat)
DT  u u  u
c T
 div(k grad 
)  (2.25)
Dt xx
x yx
y zx
z
v v v w
   
xy yy zy xz
x y z x
w w
yz  S zz i
y z
ME555 : Computational Fluid Dynamics 20 I. Sezai – Eastern Mediterranean University

10
 hip/ and h  h  1 (u2  v2 
w2 )
o 2

Specific enthalphy Specific total enthalphy

Combining these two definitions with the one for specific energy E
h  i  p /   1 (u2  v2  w2 )  E  p / 
(2.26)
o 2

Substituting of (2.26) into (2.22) yields the (total) enthalphy equation

 (ho)
 div( h u) div(k grad T )
t o


 p  (u xx )  (uyx )  (u zx )
t x y z
 (v xy )  (v yy )  (v zy )
  
x y
z
 (wxz )  (wyz )
    (wzz )  S

(2.27)
x  y
Equations
z of State h

ME555 : Computational Fluid Dynamics 21 I. Sezai – Eastern Mediterranean University

Thermodynamic variables: ρ, p, i and T.
Relationships between the thermodynamic variables can be obtained
through the assumption of thermodynamic equilibrium.
Equations of state for pressure p and specific internal energy i:
p = p(ρ, T) and i = i(ρ, T)
For a perfect gas, equations of state are
p = ρRT and i = CvT
In the flow of compressible fluids the equations of state provide the
linkage between the energy equation and mass conservation and
momentum equations.
Liquids and gases flowing at low speeds behave as incompressible
fluids.
Without density variations there is no linkage between the energy
equation and the mass conservation and momentum equations.

ME555 : Computational Fluid Dynamics 22 I. Sezai – Eastern Mediterranean University

11
 Navier-Stokes Equations for a Newtonian Fluid
We need a suitable model for the viscous stresses τij.
Viscous stresses can be expressed as functions of the
local deformation rate (or strain rate).
In 3D flows the local rate of deformation is
composed of the linear deformation rate and the
volumetric deformation rate.
All gases and many liquids are isotropic.
The rate of linear deformation of a fluid element has
nine components in 3D, six of which are
independent in isotropic fluids.
They are denoted by the symbol eij.

ME555 : Computational Fluid Dynamics 23 I. Sezai – Eastern Mediterranean University

There are three linear elongating deformation components:

u v w
e  e  e 
xx yy zz
x y z

There are also shearing linear deformation components:

 v   u  w
e  11 
ee e  
   
 v w
 e yz e
zy
1
2

z y
 

The volumetric deformation is given by

1 DV 1 dV
   

u v w V Dt
   
V dt xx yy zz

x y z
ME555 : Computational Fluid Dynamics 24 I. Sezai – Eastern Mediterranean University

12
Strain Rate Tensor
We can combine linear strain rate and shear strain rate
into one symmetric second-order tensor called the strain-
rate tensor.

 u 1 u v 1   u  w 
 
   2  z x  
x 2  yx   
    
   1 v
xx xy xz
u v 1 v
 w 
  yx  yy  yz    2  x  y  y 
 2 z
ji

    y 
    
 zx zy zz   
1   w   u  1  w  v   w 
 2  x  z   2 y 
 z  z 


ME555 : Computational Fluid Dynamics 25 I. Sezai – Eastern Mediterranean University

In a Newtonian fluid the viscous stresses are

proportional to the rates of deformation.
The 3D form of Newton’s law of viscosity for
compressible flows involves two constants of
proportionality:
- The (first) dynamic viscosity, μ, to relate stresses
to linear deformations, ( ij )1  2ij
- The second viscosity, λ, to relate stresses to the
volumetric deformation;
( xx )2  ( yy )2  (zz )2  (exx  eyy  eyy )

ME555 : Computational Fluid Dynamics 26 I. Sezai – Eastern Mediterranean University

13
The nine viscous stress components, of which six are independent, are
u v w
xx   div
 u2  yy 2   div u  zz 2
 div u
x y z
 v  u  u w 
xy  yx 
  xz

zx z x
y x
   

 v w (2.31)
    
 

Not much is known about the second viscosity λ, because its effect is
small.
For gases a good working approximation is λ = –⅔μ
Liquids are incompressible so the mass conservation equation is
div u = 0
ME555 : Computational Fluid Dynamics 27 I. Sezai – Eastern Mediterranean University

Substitution of the above shear stresses (2.31) into (2.14a-c) yields the
Navier-Stokes equations
  
 Du   p   2 u   div u    u  v
(2.32a)
Dt x x x  
 y x 
 y
  

   v w  
   S
z z  x  Mx
  

Dv  p     u  v     v
      2 

div u

(2.32b)
Dt
  y  x    y 
x
y y
    

   v w
   S
z z  y  My
  

Dw p    u w      v
       
w 
14
 Often it is useful to rearrange the viscous stress terms as follows:

 
2
u     u  v      u  w  
  div u 

 x   x  y  y  x   z 
z  x 
   
 




 u

u div( grad

u)
x x  y y   z z 
     

 u
 

w 
v
 ( div u)

sMx
 x x  y  x  z   x x
      
 div( grad u)  sMx
The viscous stresses in the y- and z-momentum equations can be re-cast in a similar
manner.
To simplify the momentum equations:
‘hide’ the smaller contributions to the viscous stress terms in the momentum
source.
Defining a new source by
SM = SM + sM

(2.33)
theME555
Navier-Stokes equations
: Computational can be written
Fluid Dynamics 29 in the most useful
I. Sezai –form
Easternfor the
Mediterranean University
development of the finite volume method:
 Du   p  div( grad u)  S
(2.34a)
Dt x Mx

Dv p
  div( grad v)  S

(2.34b)
Dt y My

 div( grad w)  S  Dw   p
(2.34c)
Dt z Mz

Note that for incompressible fluids, sM is zero if μ is constant. For example

consider
sMx:
   u    v    w  
s Mx
       ( div u)
x  x  y  x  z  x  x
      

   u    v    w     u    v    w 
      
x x  y  x  z  x  x  x  x  y  x  
z
        
    

   u v w 
   0
  
 x  x y z
    15
0 from continuity equation
ME555 : Computational Fluid Dynamics 30 I. Sezai – Eastern Mediterranean University
If we use the Newtonian model for viscous stresses in the internal
energy equation (2.24) we obtain
Di
   p div u  div(k grad T )    S
(2.35)
Dt i

The dissipation function Φ is

     u   v 
2 2
v
2 2
 w  u
 2           
 x
    y    z  y  x  
  2 
 (div u) (2.36)
 
2 2

 u  w  
v w  
  
z x  
z y 
     

The dissipation function represents a source of internal energy due

to deformation work on the fluid particle.
ME555 : Computational Fluid Dynamics 31 I. Sezai – Eastern Mediterranean University

Conservative Form of the Governing Equations of Fluid Flow

The conservative or divergence form of the time-dependent 3-D flow

and energy equations of a compressible Newtonian fluid:

Mass  div( u)  0
(2.4)
t
(u) p
x-momentum  div(uu) 
div( grad u)  S
(2.37a)
t x
Mx

 (v) p
 div( grad v)  S  div( vu)  
y-momentum
(2.37b)
t y
My

 (w) p
z-momentum  div( wu) grad w)  S
div(
(2.37c)
t z
Mz

 (i)
internal energy  div( iu)   p div u  div(k grad T )   
S (2.38c)
t
i

equations of state p  p( ,T ) and i  i( ,T )

(2.28)
e.g. for a perfect gas:

(2.29)
p   RT and i  CvT
T able 2.1 16
A system of seven equations with seven unknowns  this system is mathematically
 Differential and Integral Forms of the General Transport Equations
Equations in Table 2.1 can usefully be written in the following form:

 ( )
 div(  u)  div( grad  )  S

t

(2.39)
 Rate of increase  Net rate of flow   Rate of increase  Rate of increase
 of  of fluid 
 element  of  out of 
     of  due to 
diffusion   of  due to  
fluid element sources
       

Rate of change term convective term diffusive term source term

Equation (2.39) is used as the starting point in finite volume method.

   1, u, v, w, i 
 
By setting   0, , k  we obtain equations in Table 2.1.
 S  0, (S M    p /  x),..., 
ME555 : Computational Fluid Dynamics 33 I. Sezai – Eastern Mediterranean University

Flow Equations in Vectorial Form

The x-, y- and z-momentum equations can be written in matrix form
as  u u u
  

 x y z

 u   uu uv uw   x  
 v v v
 x
  v     vu 

 vv  vw       
   
t       y  x y
 z    
y

 
w      wu  wv  w
w  
  z

   z
 
  w w
w
v  vv   
x y

z
 
v

 u v w
 x 
x

x  p 
 f   x    
u v w
x  x

        p 
    f 
 y y 
y
y   
y  y
  z
u v w
 p     z     f z  

   
 17

pI
 z f z z 
Then, the conservation of mass and momentum equations can be
written in vectorial form as
  ( v) 0

( v)  ( vv) ((v  v T )) 
pI  f
t
I : unit tensor
v: gradient of velocity vector v.
v is a tensor given by
 u u u 
 x y z 


 vv v 
v 
 x y z 


 w

w
w 

x

y
Similarly, the product vv appearing on the left hand side of
 z is a tensor given by
momentum equations
f: body forces
 uu
per unit volume


ME555 : Computational Fluid Dynamics uv35 I. Sezai – Eastern Mediterranean University

uw 
vv  v  v   vu vv
 vw


wu

wv

ww

and is obtained from the definition of the dyadic product of two

vectors, a and b, which can be expressed as matrix product of two
matrices, the first being a column matrix and the second a row matrix, 18
Momentum equations can also be written as

( v) ( vv)  (τ  pI)
t
where, the stress tensor τ for a Newtonian fluid is given
in terms of the deformation rates as

τ   (v  v T )
Neglecting dissipation effects, the energy equation can
be written in terms of temperature T as

(c T )  (c vTp ) (kTp )
t

ME555 : Computational Fluid Dynamics 37 I. Sezai – Eastern Mediterranean University

In the finite volume method Eqn. (2.39) is integrated over 3-D control
volume yielding
 ( )
 dV   div(u)dV   div( grad  )dV  S
dV
(2.40)
t 
CV CV CV CV

For a vector a Gauss’ divergence theorem states

 div adV n  adA

(2.41)
CV A

Applying Gauss’ divergence theorem, equation (2.40) can be written

as
  
 dV  n  (u)dA  n  ( grad )dA  
S dV

(2.42) 19
 CV  A A CV
Equation (2.42) can be expressed as follows:
 Net rate of   Rate of increase
 Rate of
   decrease of  due to    of  due to 
 Net
increase of  
convection across
rate of

diffusion across  creation of 
 
the boundaries  the boundaries
   
In steady state problems the rate of change term of (2.42) is equal to
zero.

 n  (u)dA  n  ( grad  )dA   S

dV

(2.43)
A A CV

Integrating (2.42) with respect to time

  
  dV t n  (u)dAdt
t  CV  t A
  n  ( grad )dAdt    S dVdt

Auxiliary Conditions
(2.44) for Viscous Fluid Flow Equations
t A t CV
Table 2-5
ME555 Boundary
: Computational Fluid conditions
Dynamics for39 compressible viscous
I. Sezai flow. University
– Eastern Mediterranean

Initial conditions for unsteady flows:

• Everywhere in the solution region ρ, u and T must be given at time t = 0
Boundary conditions for unsteady and steady flows:
• On solid walls u = uw (no-slip condition)
T = Tw (fixed temperature) or k∂T/∂n = –qw (fixed
heat flux)
• On fluid boundaries inlet: ρ, u and T must be known as a function of position outlet:
–p +μ∂un/∂n =Fn and –p +μ∂ut/∂n =Ft (stress continuity)

Suffices: n → normal direction

to boundary
t → tangential direction
F → given surface stress
For incompressible viscous flows:
Table 2.5 is applicable, except that there are no conditions on the density ρ.
ME555 : Computational Fluid Dynamics 40 I. Sezai – Eastern Mediterranean University

20
Outflow boundaries:
• High Re flows far from solid objects in an external flow
• Fully developed flow out of a duct.

For these boundaries:

Pressure = specified
∂un/∂n = 0
∂T/∂n = 0

Sources and sinks of mass are placed on the inlet and outlet
boundaries to ensure the correct mass flow into and out of domain.

ME555 : Computational Fluid Dynamics 41 I. Sezai – Eastern Mediterranean University

Boundary conditions for an internal flow problem.

ME555 : Computational Fluid Dynamics 42 I. Sezai – Eastern Mediterranean University

21
Boundary conditions for an external flow problem.

ME555 : Computational Fluid Dynamics 43 I. Sezai – Eastern Mediterranean University

Example to symmetry boundary conditions:


 0
r

ME555 : Computational Fluid Dynamics 44 I. Sezai – Eastern Mediterranean University

22
Example to cyclic boundary conditions:

Cyclic b.c.: 1
 2

1
2
1 2

ME555 : Computational Fluid Dynamics 45 I. Sezai – Eastern Mediterranean University

23