0% found this document useful (0 votes)
5 views32 pages

Tensor Inversion and Vector Calculations

Uploaded by

劉德信
Copyright
© All Rights Reserved
We take content rights seriously. If you suspect this is your content, claim it here.
Available Formats
Download as PDF, TXT or read online on Scribd
0% found this document useful (0 votes)
5 views32 pages

Tensor Inversion and Vector Calculations

Uploaded by

劉德信
Copyright
© All Rights Reserved
We take content rights seriously. If you suspect this is your content, claim it here.
Available Formats
Download as PDF, TXT or read online on Scribd

Chapter 1

Solutions

Solution 1.3. To find the inverse of Pij = a ij + bzi zj , we look for a tensor of the same form: Qij =
c ij + dzi zj where the two constants c and d have to be determined. We then have

Qij Pjk = (c ij + dzi zj )(a jk + bzj zk ) (1.1)


= ac ij jk + adzi zj jk + bc ij zj zk + bdzi zj zj zk (1.2)
= ac ik + (ad + bc + bd)zi zk (using zj zj = 1) (1.3)

We want to have Qij Pjk = ik , so we require

ac = 1, (1.4)
ad + bc + bd = 0. (1.5)

Equation (1.4) gives c = 1/a, from which equation (1.5) becomes

b b
d(a + b) + = 0 i.e. d= . (1.6)
a a(a + b)

The inverse of Pij is therefore expressed as

1 b
Qij = ij zi zj (1.7)
a a(a + b)

Solution 1.4. In index notation, the product of the two tensors Aij = ui vj and Bjk = wj zk is written:

Aij Bjk = ui vj wj zk = (vj wj )ui zk = ↵ui zk (1.8)

with
↵ = vj wj = v · w (1.9)

Solution 1.7. In Gibbs notation: We wish to solve the equation

↵v + v ⇥ a = b. (1.10)

We first take the cross-product of equation (1.10) with a:

↵v ⇥ a + (v ⇥ a) ⇥ a = b ⇥ a. (1.11)

Applying the result of Problem 1.5, the double cross-product simplifies to

(v ⇥ a) ⇥ a = a ⇥ (v ⇥ a) = a2 v + (a · v)a. (1.12)

Using equations (1.10) and (1.12), equation (1.11) becomes

1
2 1 Solutions

↵(b ↵v) + (a · v)a a2 v = b ⇥ a. (1.13)

To eliminate (a · v), we take the dot product of equation (1.10) with a, yielding ↵(a · v) = a · b, where
we have used (v ⇥ a) · a = 0. Substituting this into equation (1.13) and solving for v, we obtain

a ⇥ b + ↵1 (a · b)a + ↵b
v= (1.14)
↵ 2 + a2

In index notation: We repeat same steps as above. The original equation (1.10) is written

↵vi + "ijk vj ak = bi . (1.15)

Take the cross-product with a:

↵"lmi am vi + "lmi "ijk am vj ak = "lmi am bi . (1.16)

Using the epsilon-delta identity, the second term simplifies to

"lmi "ijk am vj ak = ak ak vl ak vk al , (1.17)

and we also have "lmi am vi = ↵vl bl using equation (1.15). Equation (1.16) then becomes

↵(↵vl bl ) + a k a k v l ak vk al = "lmi am bi . (1.18)

Now, taking the dot product of (1.15) with a gives ↵ak vk = ak bk . We then solve for v as

"lmi am bi + ↵1 ai bi al + ↵bl
vl = (1.19)
↵ 2 + ak ak

Solution 1.15. Let us first show that for any symmetric tensor S and antisymmetric tensor T we have
S : T = 0. Indeed:

Sij Tji = Sji Tji (by symmetry of S) (1.20)


= Sji Tij (by antisymmetry of T) (1.21)
= Sij Tji (as i and j are dummy indices). (1.22)

Therefore 2Sij Tji = 0, i.e. S : T = 0.


Now for an arbitrary tensor A we have: A = sym(A) + skew(A). Then,

A : E = sym(A) : E + skew(A) : E = sym(A) : E (1.23)

since skew(A) is antisymmetric and E is symmetric. Similarly,

A : W = sym(A) : W + skew(A) : W = skew(A) : W (1.24)

since sym(A) is symmetric and W is antisymmetric.

Solution 1.18. Let us first simplify the relation between and ":
1 Solutions 3

ij = Cijrs "rs (1.25)


=[ ij rs + µ( ir js + is jr )]"rs (1.26)
= ij "rr + µ("ij + "ji ). (1.27)

Taking i = j = 1, this gives


11 = "rr + 2µ"11 . (1.28)
We need to find an equation for "rr : this can be done by taking the trace of equation (1.27):

rr = 3 "rr + 2µ"rr = (3 + 2µ)"rr . (1.29)

Recalling that rr = 11 + 22 + 33 = 11 for a uniaxial extension test, this gives


11
"rr = . (1.30)
3 + 2µ
Combining equations (1.28) and (1.30), we can solve for "11 as


"11 = 11 (1.31)
µ(3 + 2µ)
Chapter 2
Solutions

Solution 2.2.
(a) In index notation:

@ @vk @ 2 vk
r · (r ⇥ v) = "ijk = "ijk . (2.1)
@xi @xj @xi @xj

Now we note that "ijk is antisymmetric with respect to indices (i, j), whereas @ 2 vk /@xi @xj is symmetric.
Therefore, their dot product has to be zero (refer to Problem 1.15):

r · (r ⇥ v) = 0 (2.2)

(b) In index notation:


@
r · (v ⇥ w) = ("ijk vj wk ) (2.3)
@xi

@vj @wk
= "ijk w k + vj (2.4)
@xi @xi
@vj @wk
= wk "kij vj "jik using "ijk = "kij = "jik . (2.5)
@xi @xi
In Gibbs notation, this is written

r · (v ⇥ w) = w · (r ⇥ v) v · (r ⇥ w) (2.6)

Solution 2.6.
(a) In index notation:
@2
r2 r n
= r n
. (2.7)
@xi xi
Let us first calculate the gradient of r n
:

@ n n 1 @r
r = nr by the chain rule. (2.8)
@xi @xi

Now recall that r = (xk xk )1/2 . Applying the chain rule one more time,

@ 1 @ 1 xi
(xk xk )1/2 = (xk xk ) = ( ki xk + xk ki ) = . (2.9)
@xi 2(xj xj )1/2 @xi 2r r

So equation (2.8) becomes


@ n n 2
r = nr xi . (2.10)
@xi

5
6 2 Solutions

Now take the divergence of equation (2.10):

@2 n @ n 2
⇥ n 4 n 2
⇤ n 2
r = n (r xi ) = n ( n 2)r xi xi + r ii = n( n + 1)r .
@xi xi @xi
(2.11)

Therefore:
r2 r n
= n( n + 1)r n 2
(2.12)
Note that this result could also have been obtained by using the expression for the Laplacian in spherical
coordinates:
✓ ◆ ✓ ◆
1 @ @f 1 @ @f 1 @2f
r2 f (r, ✓, ) = 2 r2 + 2 sin ✓ + 2 2 . (2.13)
r @r @r r sin ✓ @✓ @✓ r sin ✓ @ 2
(b) In index notation:
@ @aj @xj
[r(a · x)]i = (aj xj ) = xj + aj = 0 + aj ji = ai , (2.14)
@xi @xi @xi
or in Gibbs notation:
r(a · x) = a (2.15)
(c) In index notation:
@
[rf (x)]i = (Ajk xj xk + bj xj + c) (2.16)
@xi
@Ajk @xj @xk @bj @xj @c
= xj xk + Ajk xk + Ajk xj + x j + bj + (2.17)
@xi @xi @xi @xi @xi @xi
= 0 + Ajk ji xk + Ajk xj ki + 0 + bj ji + 0 (2.18)
= Aik xk + Aji xj + bi , (2.19)

or in Gibbs notation:
rf (x) = (A + AT ) · x + b (2.20)

Solution 2.11. Apply Stokes’ theorem:


I Z
f (r)x · t d` = n · [r ⇥ (f (r)x)] dS. (2.21)
C S

In index notation,
@
n · [r ⇥ (f (r)x)] = ni "ijk (f (r)xk ) (2.22)
@xj

@r @xk
= ni "ijk f 0 (r) xk + f (r) (2.23)
@xj @xj
h x x i
j k
= ni "ijk f 0 (r) + f (r) jk (refer to Problem 2.6). (2.24)
r
Now, "ijk is antisymmetric with respect to indices (j, k), whereas the tensor between brackets is sym-
metric. Therefore, their dot product has to be zero (refer to Problem 1.15), from which:
2 Solutions 7
I
f (r)x · t d` = 0 (2.25)
C

Solution 2.12. To calculate the second-order tensor T, we choose an orthonormal basis (e1 , e2 , e3 ) such
that e3 = p, and we decompose the unit vector u into a component uk parallel to p = e3 , and a
component u? normal to it (in the plane spanned by e1 and e2 ).
The component parallel to p is obtained as the projection of u along e3 :

uk = (u · e3 )e3 = cos ✓ e3 . (2.26)

In the plane normal to p, the vector u? describes a circle of radius sin ✓, which can be parametrized by
an angle 2 [0, 2⇡) such that:
u? = sin ✓(cos e1 + sin e2 ). (2.27)
In the basis (e1 , e2 , e3 ), the unit vector u therefore has components:
2 3
sin ✓ cos
u = 4 sin ✓ sin 5 , (2.28)
cos ✓

and the dyadic product uu is


2 3
sin2 ✓ cos2 sin2 ✓ sin cos sin ✓ cos ✓ cos
uu = 4 sin2 ✓ sin cos sin2 ✓ sin2 sin ✓ cos ✓ sin 5 . (2.29)
sin ✓ cos ✓ cos sin ✓ cos ✓ sin cos2 ✓

Also note that the element of length d` along the circle is d` = sin ✓ d . The tensor T can therefore be
calculated as Z 2⇡
T= uu sin ✓ d . (2.30)
0
Integrating each component with respect to , we obtain:
2 3
⇡ sin3 ✓ 0 0
T=4 0 ⇡ sin3 ✓ 0 5. (2.31)
2
0 0 2⇡ cos ✓ sin ✓

In order for this tensor to be proportional to the identity tensor, i.e. T = ↵I, we require that all three
diagonal components be equal, which leads to the following condition on the angle ✓:

sin2 ✓ = 2 cos2 ✓ sin ✓, (2.32)

which simplifies to
tan2 ✓ = 2 (2.33)

leading to a magic angle of ✓ ⇡ 54.7o . The corresponding value of ↵ is ⇡ sin3 ✓ ⇡ 1.71.

Solution 2.13. For an arbitrary body of volume V , we can write the moment inertia with respect to the
axis of direction n as Z Z
I= ⇠ 2 dV = ⇠i ⇠i dV, (2.34)
V V
where the vector ⇠ denotes the component of the position vector x that is normal to n. We note that
⇠ = x (x · n)n or in index notation: ⇠i = ( ij ni nj )xj . Therefore:
8 2 Solutions
Z
I= ( ij ni nj )xj ( ik ni nk )xk dV (2.35)
V
Z
=( ij ni nj )( ik ni nk ) xi xk dV (2.36)
V
Z
=( jk nj nk ) xi xk dV, (2.37)
V

where we have used the fact that unit vector n isRfixed and therefore can be taken outside of the integral.
To find I, we need to evaluate the tensor Ajk = V xj xk dV for a particular shape (cube and sphere).

(a) Cube: Because Ajk does not depend on n, we can choose any coordinate system to evaluate it. Let
us simply choose the coordinate system with axes parallel to the cube faces. In this basis, the tensor Ajk
is written 2 3
Z Z Z a/2 x2 xy xz
Ajk = 4 xy y 2 yz 5 dx dy dz. (2.38)
a/2 xz yz z 2

Performing the integrals, we obtain:


2 3
1 0 0
a5 4 a5
Ajk = 0 1 05 = jk . (2.39)
12 12
0 0 1

Noting that I = ( jk nj nk )Ajk and ( jk nj nk ) jk = kk nk nk = 2, we obtain

a5
Icube = (2.40)
6

(b) Sphere: By symmetry, and having gained insight from the case of the cube, we expect the tensor Ajk
to be of the form Ajk = jk where is a constant to be determined. To find , we can calculate the
trace of Ajk : Z Z
Ajj = xj xj dV = r2 dV = jj =3 . (2.41)
V V
In spherical coordinates, we have
Z Z 2⇡ Z ⇡ Z a/2
⇡a5
r2 dV = r4 dr d✓ sin ✓d = . (2.42)
V =0 ✓=0 r=0 40

Therefore, = ⇡a5 /120, and I = 2 gives

⇡a5
Isphere = (2.43)
60
Chapter 3
Solutions

Solution 3.3. (a) At t = 0, we have u1 = x1 and u2 = 1. The equation for the streamline is

dx1 dx2
= , (3.1)
u1 u2
which leads to
dx1
= dx2 , (3.2)
x1
This integrates to
ln |x1 | = x2 + C. (3.3)
To determine the integration constant C, we specify that the streamline must pass through the point
(x1 , x2 ) = (1, 1), which leads to C = 1, and the equation for the desired streamline becomes

ln |x1 | = x2 1 (3.4)

(b) Pathline: The equations for particle trajectories or pathlines are

dx1 dx2
= u1 = x1 (1 + t), and = u2 = 1, (3.5)
dt dt
with initial conditions x1 (0) = 1, x2 (0) = 1. The equation for x1 (t) is solved by separation of variables:

dx1
= (1 + t) dt, (3.6)
x1
which integrates to
t2
ln |x1 | = t + + C1 , (3.7)
2
and the initial condition on x1 gives C1 = 0. Similarly, the equation for x2 (t) integrates to

x2 (t) = t + C2 , (3.8)

with C2 = 1 from the initial condition. Equations (3.7) and (3.8) can be combined to eliminate t and
obtain the equation for the pathline:
x2 1
ln |x1 | = 2 (3.9)
2

9
10 3 Solutions

(c) Streakline: We still solve equations (3.5), with general solutions given by (3.7) and (3.8). However,
we now determine the constants of integration by using the condition that the fluid particle must have
passed through the point (1, 1) at a prior time t̂:

x1 (t̂) = 1, x2 (t̂) = 1. (3.10)

This can be used to find C1 = t̂ t̂2 /2 and C2 = 1 t̂, and therefore

t2 t̂2
ln |x1 (t; t̂)| = t + t̂ , and x2 (t; t̂) = 1 + t t̂. (3.11)
2 2
Equations (3.11) give the position at time t of a fluid particle that passed through (1, 1) at time t̂. To
determine the equation for the streakline at time t = 0, we set t = 0 in equations (3.11), yielding:

t̂2
ln |x1 (t̂)| = t̂ , and x2 (t̂) = 1 t̂. (3.12)
2
These can be combined to eliminate t̂, leading to the equation of the streamline:

1
ln |x1 | = (x2 1)(x2 3) (3.13)
2

Solution 3.6. (e) he velocity field is given by uj = xi ij where


2 3
2 0 0
ij =
4 0 1 15. (3.14)
0 1 1

We first calculate the rate of strain tensor as


2 3
2 0 0
1
Eij = ( ij + ji ) = 4 0 1 0 5 . (3.15)
2
0 0 1

Since it is diagonal, we see that the principal axes are the coordinate axes, with compression along x1
with rate 1 , and extension along x2 and x3 with equal rates 2 = 3 = 1.
We also calculate the rate-of-rotation tensor:
2 3
0 0 0
1 40 0 15.
Wij = ( ij ji ) = (3.16)
2
0 1 0

The corresponding angular velocity vector (or dual vector of the tensor Wij ) is obtained as

1
⌦k = ✏kij Wij = [1, 0, 0]. (3.17)
2
Therefore we see that there is also rotation of fluid elements in the counter-clockwise direction around
the x1 axis.
Based on these observations, we expect fluid particles to move towards the (x2 , x3 ) plane, which
moving away from the x1 axis, while rotating around the x1 axis. In other works, the streamlines should
look like expanding spirals around the x1 axis:
3 Solutions 11
Chapter 3
Solutions
Chapter 3
20 Chapter 3 Solutions
0
Solutions
Chapter 3
x3 20
Solutions 0

40
Chapter 3 Solution 3.7. We consider the sink flow:
0.5
60
20 Solutions
0
1
x
Solution 3.7. We consider the sinkuflow:
i (xj ) =
xi
,
20 40 60
r3
x2 80 100 1 xi
where 0 and r the )1/2 . ui (xj ) =
Solution 3.7. We>consider = (x k xkflow:
sink r3
Solution 3.7. We consider the sink flow:
Solution 3.7. We consider the sink flow: (a) We calculate
where the
> velocity
0 and r =gradient
(xk xkin index
)1/2 .
i j
notation:
xi
u (x ) = xi ,
xi 3
ui (xj ) = , @u @ ⌘ rui (xj ) = r3 ,
⇣ x (3.18)
r3 (a) We calculate the
j velocity gradient
=
j in index notation:
Solution 3.7. whereWe consider 0 and the
> where r => sink xflow:
(x0kand k)r =
1/2 @x
. (xi x )1/2@x ✓.j i
r3
where > 0 and r = (xk xk )1/2 . k k @u
1=x@x @ ⇣ x@j ⌘ ◆
j 3
(a) We calculate
(a) We thecalculate
velocity the gradient ui (x
velocity
=
inj )index
=@xi rnotation:
gradient
i
in
+ xj notation:
@x
, i index r3 (r )
(a) We calculate the velocity gradient in index notation: 3 3@x
r ✓ i @x i
✓ 1 ⌘@x ◆j @
@uj @ ⇣ xj ⌘ @uj ⇣
@= j xj =
@u ⌘ ij @ ⇣ xxi xj j + xj (r 3
where > 03 and r = (xk xk ) .
= 1/2
= = 3
(3.19) r 3 @x i (refer @x toi Pro
@xi @xi r @x @x
@x i r 3 r 3@xi ✓r r35 ◆
i

1 @xj the @velocity
◆ ✓i ✓ ij ◆ xi xj ◆
(a)
= We calculate + xj (r ) 3 gradient in index notation: = 1 @x= j ( = 1@x̂ @x
ij (3.20) x̂ j) 3
i jr 3 3+ xj r 5 (r
@ 3 (ref
)
r3 @xi @xi = 3 @x
r3 + xj r3 @x (ri ) @xi
✓ ◆
@uj @ ⇣ xj ✓ ⌘ r i ✓ i @x ◆
xi xj =x̂ ◆ ( ij xiNoting
xx̂ji x̂j ) that the v
=
ij
3 where we
(refer have
=
to introduced
Problem 2.6) a unit vector
= xi xj 3 = ij3x/|x|.
(3.21)
r 3 (refer to
r3 r5
already @x symmetric,
i @x
wei=find r3 thatijthe rate-of-strain
3 r◆ (refer and r 5vorticity
to Problem tensors2.6

= ( 3x̂ x̂ ) where we have 1 @x introduced
j
r3 @a unit r5 3vector(3.22) x̂ = x/|x|. Noting
ij i j
r3 =
already symmetric, 3 @x we find
+ xj = that (r () ijrate-of-strain x̂i x̂j )
r= i (Eij = @xx̂ii x̂rj3)the
( ij x̂i x̂j ) and Wij
and vortici
where we have introduced a unit vector x̂ = x/|x|. Noting that the velocity ✓ gradient r3 x in ij ◆
equation r 3 (3.22) is
where we have introduced ij ixajunit vector x̂ = x/|x|. Noting that
already symmetric, we find that the rate-of-strain and vorticity tensors = are given 3by: 5 (refer to Problem
Erate-of-strain 2.6)
where we already
have symmetric,
introduced a unit wevector
find that
x̂ =ofthe ij = Noting ( ijand
3 that
x̂ i x̂j )
vorticity an
teg
thethe velocity
r 3 r
To find the principal axes and rates x/|x|.
extension r(i.e. eigenvectors
already thesymmetric, and we
i x̂components
find
ofij E=that
(ij0ijinthe rate-of-strain
a x̂basis. and(3.23)vorticity tensors are given
Eij = ( ij 3x̂ j) =W i x̂j ) Rather than choosing an arbitrary b
r3 To find the
r 3 principal
significantly), let us choose a basis in whichrthe axes E and
ij = rates of ( extension x̂ x̂
3 first unit vector e1 is
ij i j(i.e.
) the
andeig
shown the components of Eeij in ax̂).
basis. Rather than choosing an a
To find the principal axes andwhere
rates ofwe have introduced
extension (i.e. the below a (in vector
unit
eigenvectors other
and E words,
eigenvalues
= 1 = of
( Noting
Eij ), x̂we i x̂that
jwrite the
) which and velocity ij gradient
To significantly),
find the principal
x̂ =
letaxes
ij x/|x|.
us choose
and
r 3 ij
ratesa basis
of in
extension theWfirst
(i.e. the
=unit
0 ve
eigenve
the components of Eij in a basis. Rather than choosing
already symmetric, we findshown O an arbitrary basis (which
that the rate-of-strain would complicate
and vorticity the
tensors are given by:
algebra significantly), let us choose a basis in which xthe first unit below
thecomponents vector of (in other
eE1ijis in a words,
alignedbasis. ethe
withRather1 =radial x̂).
than choosing an arbitr
vector x as shown below (in other words, To e1 find
= x̂).theesignificantly),
principal
1 = x̂
Oaxes and let usrates
choose of extension
a basis (i.e.
in which the eigenvectors
the first unitand vectoreige
the components shown of xE
below
E in
ij = a
(in basis.
other ( Rather
words,
ij x̂ than
i x̂
e j )=choosing
and
x̂). anW arbitrary
ij = 0 basis (whi
e2 ij r3 1
significantly), e3 O let usechoose
1 = x̂ a basis in which the first unit vector e is aligned w
1
shown below x
(in e2 words, e = x̂).
other
To find the principal axes and rates of extension 1 (i.e. the eigenvectors and eigenvalues
e1 = ex̂3
the components O of Eij in a basis. Rather than choosing an arbitrary basis (which comp
e2
significantly),xlet us choose e3
a basis in which the first unit vector e1 is aligned with the
shown belowe(in 1 = other
x̂ words, e1 = x̂).
O e2
x e3
e1 = x̂
e2
e3
ui (xj ) = ,
r3 @u
Solution 3.7. We consider where
the sink > flow: 0 and r = (xk xk ) . 1/2

where > 0 and r = (x x ) 1/2


. @x
Solution 3.7. We 12 consider the sink flow: k k
xi
Solution 3.7. We consider the sink flow:
(a) We calculate the
ui (x velocity gradien
2j ) = 3 3
,
(a) We calculate the velocity gradient in index notation: xi 2 0 r0
ui (xj ) = , x@u
= 4 0 we 1, 0= 5intr
i j
⇣ x. ⌘ E ru3iji (x=j )where .
where > 0 and @ujr = (xk x@k )1/2 j r3 r 3@xhave
0 0 1
@xi r33 Solutions already symmetric,
12 = i
where
Solution 3.7. 0 and the
We consider
> where r(a)
=>sink
(x0kand
We xflow:1/2 @x
k)r =
calculate . (x i
the k xvelocity . gradient in index notation:
1/2
k) ✓ ◆
Since Eij is diagonal in this basis, 1 x@x wej know @ that the eigenvectors = a
In this basis, we have 2(a)3We calculate the velocity
2 3 = i +the xj notation: (r 3
⌘ )of thewhere w
1 (a) read
We the
calculate
1 0 0 gradient
corresponding
the u
velocity
i (xin j index
eigenvalues
) =
gradient
@u notation:
3
rj r= in
3@xi
, on
index @ ⇣diagonal
@xijx matrix.
ofj E ✓/r > 0,@x ◆ 3 e2 where and e3 we already boths
have
=ij4 with
0 0 0extension rate whereas are
3
x̂i = 4 0 5 and x̂i x̂ 5, ⇣ x2@x ⌘ i ij @ ⇣xx (3.24)
x i⌘ r =
@u @
@u ✓
i j
already symme ◆ th
)1/[Link]
0 j0 that= the@x direction=j 3of e23and
j e are(refer indeterminate in
where >00 and r = (xk xk/r
3 = j
r
where
r 5 1 3@xwe
3 To
j
have
find to introduced
the
@ Problem
principa
3
@xi
form an orthonormal @x 3
i i r Locally, =@xi there r is stretching + xj along (r the )find
therefore the components of the rate-of-strain tensor obtained in equation (3.23) ✓[Link] ✓ already
r3 @x ◆ the symmetric,
i @ components @xi ◆ we = of
ra
(a) We calculate the velocity gradient in
the orthoradial directions. index notation:
=1 @x = ( 1@
x̂ @x ✓
x̂ )
j ◆ 3 r3
2 3 = r 3 + x 3 (r 3+
j ij i j
) xjsignificantly),
x x (r ) let u
= j r @x i ij @x i
3◆ i5 below
j
(refer tot
2 0 0
(b)
@ujWe 5now @assume ⇣ x ⌘that r3 @x theiwhere ✓ @xi israligned
dumbbell shown To findve
(in ot
Eij = where 4 0 we j ✓ ◆ we have
(3.25)
3
x x rwith the
introduced a axis
unit of
1 have= introduced 3a unit vector
0 .
xi xx̂j =3symmetric,
ij x/|x|. 3 Noting
i j that the velocit
ij = already OTo find
(refer the
the
to prin
comp
Prob
0 direction
0 i 1 of@x e (radialr thatdirection):
r3 @x we find that the
already symmetric, ✓we1i=find the rate-of-strain
3 = r◆ (refer ( and r 5vorticity
to x̂ Problem
x̂ ) tensors 2.6) are gi
1 @xj r @ r 3 r
3 5 3 ij i
xthe component
j significa
() ij, e To x̂ifind
x̂j ) theesignificantly),
principal axes al
1 = x̂ shown b
Since Eij is diagonal in this basis, we know that the eigenvectors = arer 3the@x basis xj = (r
+ vectors 2 3 ), and
where we have = introduced
i (Eij = @xx̂iiax̂(e 31 , e
runit
) ( vectorthe x̂ components

x̂ =
) x/|x|.
shown
and Noting
of EO
Wbelow
E =ijthat
in0= (a
we can directly read the corresponding eigenvalues on the diagonal✓of the rmatrix. 3 ij ◆ Therefore,r
j
3 e1 = x̂
ij i is j e2 ij ij
already
3 where we have introduced symmetric, x wex
a find
unit that
vector the rate-of-strain and
Noting vorticity
that the tev
tosignificantly), let us choose
a principal axis of Eij with extension rate 2 /r > 0, whereas x̂ = x/|x|. O
= e2 and3 e3 3are both principal
(refer axes with
ij i j x
Problem e2.6)
3
/r3 . Note that thewhere we already
have symmetric,
introduced a unit
r wevector
find r 5that the rate-of-strain and vorticity
Notingshown that the
below x and
velocity
(in tensors
gradien
other
compression rate directions of e and e are indeterminate in
x̂ =the
To plane
x/|x|.
find
To find the principal axes and rates of extension (i.e. the eigenvectors eand
2 3 normal
the to
principal axes 1wor
rates =e
e1 , as long as they form an orthonormalalready
[Link]
Locally,
symmetric, there weis stretching
find that thealong the
rate-of-strain radial Edirection,
and =vorticity and ( tensors x̂ ej1)given
are
x̂ = x̂ and by:
components = of 3E(ij ijin a x̂basis. i x̂j ) Rather the components
than choosing
ij
rO 3 ofan ij
ij in a basis.
Earbitrary
i e
basis2 R (w
compression along the orthoradial directions.
significantly), let us r choose a basisEinij =
significantly),
which the
3
( first
ij xlet x̂ius
unit ) eand
x̂vector
j choose 2
e ais e W
basis
aligne
3 ij i
r e3 1
(b) We now assume that thewhere we have
dumbbell shownwithbelow
introduced
is aligned theaTounit
findvector
(in
axis of the
other principal
words,
maximum
E x̂ij = e1axes =( x̂).
extension,
= x/|x|. andwhich
shown
Noting
ij
rates
x̂ i below
that
x̂ j ) ofthe
is inextension
and e(in
3 rates of extension (i.e. the eigenvectors
= other
velocity
the
1 x̂
W (i.e.
ij words,0the eigenve
gradient
= ein1 equ=
direction of e1 (radial direction): OTo find the
the principal
components axes
already symmetric, we find that the rate-of-strainijand vorticityof Eand
r in a basis.
O Rather
tensors than choosing
e2 are given by: an arbitr
xthe components
significantly), of Elet ij inusachoosebasis. Rather axbasis than in which choosing
e3 first the first an unit
arbitraryvector ba
To find theesignificantly),
principal
1 = x̂ R shownaxes and
let
below usrates
choose
(in of
other extension
a basis
words, in(i.e.
e which
= the x̂).eigenvectors
the unit andvectoreigenvalue
e 1 is a
e1and 1= x̂
the components
shown
e2 of E
below
E O
ij in
ij = a
(in basis.
other
r3
( Rather
words,
ij x̂than
i x̂
e j
1 )=choosing
x̂). anW arbitrary
ij = 0 basis (which com
e2
significantly),
e3 O let us choose
x a basis in which the first unit vector e1 is aligned with th
x and e3
To find theshown below
principal (in
axes otherrates
e1words,
= of x̂ extension
e1 = x̂). (i.e. the eigenvectors and eigenvalues of Ei
e1 = x̂e
the components O of Eij in a basis. 2 Rather than choosing an arbitrary basis (which complicate
e2 e3 in which the first unit vector e1 is aligned with the radial
significantly),xlet us choosee3
a basis
shown belowe(in 1 =other
x̂ words, e1 = x̂).
O e2
x e3
The equation for the end-to-endevector
1 = x̂Ri is
e2
dRi @ui 1
e3 Rj + Ri = 0. (3.26)
dt @xj 2

Noting that Rj = Rx̂j (since Rj and xj are both in the radial direction), we have

@ui 2
Rj = Rx̂j ( ij 3x̂i x̂j ) = R(x̂i 3x̂i ) = Ri , (3.27)
@xj r3 r3 r3

which is not surprising since Rj is aligned with the principal axis of extension and is therefore an eigen-
vector of the symmetric velocity gradient with eigenvalue 2 /r3 . Therefore, equation (3.26) becomes

dRi 2 1
Ri + Ri = 0. (3.28)
dt r3 2
This can be projected along the radial direction, leading to the following equation for the end-to-end
distance R:
1 dR 2 1
= 3 (3.29)
R dt r 2
We see that the dumbbell will either stretch or contract depending on the sign of the right-hand side,
which is a function of the distance r = |x| from the origin:
3 Solutions 13
p
• If r < p
3
4 , the dumbbell stretches (dR/dt > 0);
• If r > p
3
4 , the dumbbell contracts (dR/dt < 0);
• If r = 3 4 , the length of the dumbbell does not change (dR/dt = 0).

(c) In the previous question, we assumed that the dumbbell was fixed in space, so the distance r from the
origin was a constant. Now if the dumbbell is free to move, r will change with time. The equation for the
position of the dumbbell is
dxi xi
= , (3.30)
dt r3
where xi = rx̂i . Projecting this equation along the radial direction, we get

dr
= , (3.31)
dt r2
which integrates to r3 = r03 3 t (where r0 = r(0)). Equation (3.29) for the length of the dumbbell then
becomes ✓ ◆
dR 2 1
= dt. (3.32)
R r03 3 t 2
This can be integrated to give (for t > 0):

2/3
t/2 3
R(t) = R0 e 1 t (3.33)
r03

where R0 = R(0). The condition for the dumbbell to stretch over the time interval [0, tf ] is R(tf ) > R0 ,
which gives
2/3
tf /2 3
e 1 tf >0 (3.34)
r03

Solution 3.9. (a) We consider the velocity field:


xi
ui (x) = (3.35)
r3
where x is the position vector (pointing from origin O) and r = |x|. Since at every point in space the
velocity is in the radial direction, the streamlines are straight lines emanating from the origin (source
flow):

(b) We now consider the deformation of a ring of dye centered at the origin, see figure below. Based on the
velocity field above, we expect the radius R(t) to increase with time. The cross-section of the ring, which
14 3 Solutions

is initially circular with radius b0 , will also deform. The velocity gradient for this flow was calculated in
Problem 3.7, where we find that there is compression along the e1 direction (with rate = 2/R3 ), and
stretching in the orthogonal direction e2 (with rate = 1/R3 ). Let us denote by b(t) and c(t) the drop
dimensions along these two directions, where b(0) = c(0) = b0 .

We first solve for R(t) by writing that point A shown above, with position R(t)e1 , moves with the fluid:

d(Re1 ) 1
= u(Re1 ) = 2 e1 . (3.36)
dt R
Projecting this along the direction e1 , we get

dR 1
= 2, (3.37)
dt R
which integrates to
q
R(t) = 3
R03 + 3t (3.38)

Next, we use the fact that b(t) ⌧ R0 (t) to write

d(be1 ) 2b
⇡ be1 · ru = e1 (3.39)
dt R3
where we have used the fact that e1 is an eigenvector of ru with eigenvalue 2/R3 . Therefore,

db 2
= dt, (3.40)
b R03 + 3t

which integrates to
✓ ◆ 2/3
3t
b(t) = b0 1+ 3 (3.41)
R0
A similar analysis in the direction of e3 yields
✓ ◆1/3
3t
c(t) = b0 1+ 3 (3.42)
R0
Chapter 4
Solutions

Solution 4.2. In index notation:


Z Z
d d @f
f ni dS = dV (Gauss theorem) (4.1)
dt S(t) dt V (t) @xi
Z ✓ ◆ Z
@ @f @f
= dV + w j nj dS (Leibnitz theorem) (4.2)
V (t) @t @xi S(t) @xi
Z ✓ ◆ Z
@ @f @f
= dV + w j nj dS (equality of mixed partials) (4.3)
V (t) @x i @t S(t) @x i
Z Z
@f @f
= ni dS + w j nj dS (Gauss theorem) (4.4)
S(t) @t S(t) @x i

or, in Gibbs notation:


Z Z Z
d @f
f n dS = n dS + (w · n)rf dS (4.5)
dt S(t) S(t) @t S(t)

which can be interpreted as a version of Leibnitz’s theorem for surface integrals.

Solution 4.4. (a) The continuity equation r · v = 0, in cylindrical coordinates (r, ✓, x), gives

1 @ 1 @v✓ @vx
(rvr ) + + =0 (4.6)
r @r r @✓ @x
with vr = ↵r and v✓ = 0. This provides an equation for vx :

@vx 1 @
= (↵r2 ) = 2↵. (4.7)
@x r @r
Using the boundary condition vx (0) = V0 , this integrates to

vx (x) = 2↵x + V0 (4.8)

(b) For the velocity vector to be tangent to the wall at the nozzle boundary, we must have

dR vr
= at r = R(x), (4.9)
dx vx
which expresses that the channel wall and velocity have the same slope and or equivalently that the wall
is a streamline. This yields

dR ↵R
= , (4.10)
dx 2↵x + V0

15
16 4 Solutions

dR 1 dx
i.e., = . (4.11)
R 2 x + V0 /2↵

This integrates to ✓ ◆ ✓ ◆
R 1 x + V0 /2↵
ln = ln , (4.12)
R0 2 V0 /2↵
from which we find:
✓ ◆1/2
V0 /↵
R(x) = R0 (4.13)
2x + V0 /↵

The shape R(x) of the channel should not depend on the velocity of the flow, and therefore V0 /↵ should
be a constant independent of V0 . This implies that V0 and ↵ must be proportional.
(c) The mass flow rate across cross-section S at streamwise position x is expressed as:
Z
Qm (x) = ⇢vx (x) dS = ⇢vx (x)⇡R(x)2 (4.14)
S

At the nozzle entrance (x = 0), we find:

Qm (0) = ⇡⇢V0 R02 (4.15)

At the nozzle exit (x = L),

V0 /↵
Qm (L) = ⇡⇢(V0 + 2↵L)R02 = ⇡⇢V0 R02 (4.16)
2L + V0 /↵

We find that Qm (0) = Qm (L), which is a consequence of conservation of mass in steady flow.

Solution 4.10. In index notation, and using the Newtonian viscous stress ⌧ = µ(rv + rvT ):
✓ ◆
@⌧ji @ @vi @vj @ 2 vi @vj
=µ + =µ +µ . (4.17)
@xj @xj @xj @xi @xj @xj @xj @xi

The last term is also written ✓ ◆


@vj @ @vj
= = 0, (4.18)
@xj @xi @xi @xj
where we have used the incompressibility condition @vj /@xj = 0. Equation (5.1) then becomes:

r · ⌧ = µr2 v (4.19)

Solution 4.11. Using the expression for the Newtonian viscous stress, we have

⌧ : rv = [ (r · v)I + 2µE] : rv (4.20)


= (r · v)(I : rv) + 2µE : rv. (4.21)

We note that I : rv = ij (@vi /@xj ) = @vi /@xi = r · v, and = 2µ/3 from the Stokes assumption.
Also, decomposing rv into its symmetric and antisymmetric parts, we have

E : rv = E : (E + W) = E : E + 0 (4.22)
4 Solutions 17

where we have used the result that the scalar product of symmetric tensor E with antisymmetric tensor
W is zero (refer to Problem 1.15). Therefore, equation (5.5) for the viscous dissipation becomes

2
= µ(r · v)2 + 2µE : E (4.23)
3

Solution 4.12. (a) To determine whether the flow is incompressible, we calculate the divergence of the
velocity field in spherical coordinates:

1 @ 2 1 @ 1 @v
r·v = 2
(r vr ) + (sin ✓ v✓ ) + (4.24)
r @r r sin ✓ @✓ r sin ✓ @
 ✓ 3
◆  ✓ ◆
1 @ 2 a 1 @ 2 a3
= 2 U cos ✓ r + U sin ✓ 1 + 3 (4.25)
r @r r r sin ✓ @✓ 2r
✓ ◆ ✓ ◆
U cos ✓ a3 U a3
= 2r + 2 1 + 3 2 sin ✓ cos ✓ (4.26)
r2 r r sin ✓ 2r
=0 (4.27)

Since r · v = 0 , the flow is incompressible.

(b) For an incompressible flow, the viscous stress is ⌧ = µ(rv + rvT ). The non-zero components of
the velocity gradient in spherical coordinates are
2 3
@vr @v✓ 2 3
6 0 7
6 1 @v@r v @r 7 3 2 cos ✓ sin ✓ 0
1 @v✓ vr 3U a
rv = 6 6
r ✓
+ 0 7=
7
4 sin ✓ cos ✓ 0 5 . (4.28)
4 r @✓ r r @✓ r v 5 2r4
r v ✓ 0 0 cos ✓
0 0 + cot ✓
r r
We note that rv is symmetric, and therefore the viscous stress tensor is given by
2 3
2 cos ✓ sin ✓ 0
3µU a3 4
⌧ = sin ✓ cos ✓ 0 5 (4.29)
r4
0 0 cos ✓

(c) For an incompressible flow, the viscous dissipation is given by

= 2µE : E = 2µtr(E · ET ). (4.30)

Using the expression for the rate-of-strain tensor obtained in (b), we get
2 3 2 3
✓ ◆
3 2 2 cos ✓ sin ✓ 0 2 cos ✓ sin ✓ 0
3U a 4 sin ✓ 5 ⇥ 4 sin ✓
E · ET = cos ✓ 0 cos ✓ 0 5 (4.31)
2r4
0 0 cos ✓ 0 0 cos ✓
2 2
3
3 cos ✓ + 1 cos ✓ sin ✓ 0
9U 2 a6 4 5.
= cos ✓ sin ✓ 1 0 (4.32)
4r8
0 0 cos2 ✓

Taking the trace and multiplying by 2µ, we get


18 4 Solutions

9µU 2 a6
= (2 cos2 ✓ + 1) (4.33)
r8

(d) Integrating over the volume outside of the sphere of radius a gives
Z Z 1 Z ⇡ Z 2⇡
dV = (r, ✓) r2 sin ✓ dr d✓ d . (4.34)
V r=a ✓=0 =0

All calculations done, we find:


Z
dV = 12⇡µU 2 a (4.35)
V
Chapter 5
Solutions

Solution 5.3. Circular Couette flow


(a) The momentum equation in cylindrical coordinates is written
 ✓ ◆
@ur @ur 1 @ur u✓ @ur @p
⇢ + ur + u✓ + uz =
@t @r r @✓ r @z @r
 ✓ ◆
@ 1 @ 1 @ 2 ur @ 2 ur 2 @u✓
+µ (rur ) + 2 2
+ ,
@r r @r r @✓ @z 2 r2 @✓
(5.1)
 ✓ ◆
@u✓ @u✓ 1 @u✓ ur @u✓ 1 @p
⇢ + ur + u✓ + + uz =
@t @r r @✓ r @z r @✓
 ✓ ◆
@ 1 @ 1 @ 2 u✓ @ 2 ur 2 @ur
+µ (ru✓ ) + 2 + + 2 ,
@r r @r r @✓2 @z 2 r @✓
(5.2)

@uz @uz u✓ @uz @uz @p
⇢ + ur + + uz =
@t @r r @✓ @z @z
 ✓ ◆
1 @ @uz 1 @ 2 uz @ 2 uz
+µ r + 2 + . (5.3)
r @r @r r @✓2 @z 2

Noting that the flow is azimuthal (u(r, ✓, z) = u✓ (r)e✓ ) and axisymmetric (@/@✓ ⌘ 0), these equations
simplify to

u2✓ @p
⇢ = , (5.4)
r @r ✓ ◆
d 1 d
0=µ (ru✓ ) , (5.5)
dr r dr
@p
0= . (5.6)
@z
In particular, we find that @p/@✓ = @p/@z = 0, and therefore p = p(r).
(b) In order to determine the velocity, integrate equation (5.5):

1 d d
(ru✓ ) = A i.e. (ru✓ ) = Ar. (5.7)
r dr dr
Integrate once more to yield

r2 B
ru✓ = A +B i.e. u✓ (r) = Cr + . (5.8)
2 r
To find the integration constants C and B, we apply the no-slip condition at the inner and outer cylinders:
u✓ (r = Ri ) = !i Ri , and u✓ (r = Ro ) = !o Ro . This yields a system of equations for the integration
constants:

19
20 5 Solutions

B B
CRi + = !i Ri , CRo + = ! o Ro , (5.9)
Ri Ro
the solution of which is easily obtained as

!i Ri2 !o Ro2 (!o !i )Ri2 Ro2


C= , B= . (5.10)
Ri2 Ro2 Ri2 Ro2

Substituting these expressions into equation (5.8) provides the velocity field:
✓ ◆ ✓ ◆
!i Ri2 !o Ro2 !o !i Ri2 Ro2
u✓ (r) = r+ (5.11)
Ri2 Ro2 Ri2 Ro2 r

(c) To determine the pressure distribution, we use equation (5.18):


✓ ◆2 ✓ ◆
dp u2 ⇢ B B2 2CB
2
=⇢ ✓ = Cr + =⇢ C r+ 3 + . (5.12)
dr r r r r r

This integrates to ✓ ◆
C 2 r2 B2
p(r) = ⇢ + 2CB ln r + D, (5.13)
2 2r2
where D is an arbitrary constant. Subtituting the values of C and D found above, we obtain

⇢ r2 R2 R2
p(r) = (!o Ro2 !i Ri2 )2 (!o !i )2 i 2 o
(Ro2 Ri2 )2 2 2r (5.14)

2(!o !i )(!o Ro2 2 2 2
!i Ri )Ri Ro ln r + D.

Solution 5.4. Flow down an inclined plane


(a) The momentum equation for an incompressible flow in two dimensions is written in components as
✓ ◆ ✓ 2 ◆
@vx @vx @vx @p @ vx @ 2 vx
⇢ + vx + vy = ⇢fx +µ + , (5.15)
@t @x @y @x @x2 @y 2
✓ ◆ ✓ 2 ◆
@vy @vy @vy @p @ vy @ 2 vy
⇢ + vx + vy = ⇢fy +µ + . (5.16)
@t @x @y @y @x2 @y 2

In the present problem, the flow is steady (@/@t = 0), and unidirectional in the x direction (vy = 0
and @/@x = 0). We also note that (fx , fy ) = (g sin ✓, g cos ✓). Using these conditions, the equations
simplify to

d2 vx
⇢g sin ✓ + µ =0, (5.17)
dy 2
dp
⇢g cos ✓ + =0. (5.18)
dy

(b) Integrating equation (5.18) with respect to y gives

p(y) = ⇢g cos ✓ y + C. (5.19)

The integration constant C is obtained using the fact that p(h) = pa (atmospheric pressure), yielding
dy 2 µ dx
This can be integrated twice:
⇢g dh 2
u(y) = y + Ay + B.
2µ dx
We use the boundary conditions of part (b) to solve for the integration constants A
5 Solutions 21
no-slip condition u(0) = 0 gives B = 0. The condition on the shear stress at the free s
p(y) = pa ⇢g cos ✓ (y h) (5.20)
d2 u ⇢g dh 6 Solutions
S S ⇢g dh
2
(h 0 ) = h0 + A = i.e. A = h0 .
dystress ⌧ are µ dx
(c) The components of the normal n and of the viscous µ µ µ dx
2 33 2
The velocity profile is2 @v
therefore
x @vxgiven
+ by:
@vy dvx
0
6 7 6 candy be7
n = [0, 1] , ⌧ = µ ru + ruT = µ 4 @v @x@v @y @v @x 5 =This
4 dv 5integrated
. (5.21) twice:
x
+
y
2
y x
0
✓ 2 ◆
@y @x @y S dy ⇢g dh y u(y)
u(y) = y+ yh0
µ
The condition of zero-traction t = n · ⌧ = 0 at y = h therefore is written
µ dx 2
We use the boundary conditions of part
5 Solutions
Note that the slope no-slip
dvx of the free surface remains condition u(0) = 0 gives B = 0.21T
unknown.
(y = h) = 0 (5.22)
dy (5.20)
(e) The volumetric flow rate perp(y)
unit=width
pa ⇢g incos
the✓ (y h)
z direction is d2 u ⇢g dh
(h0 ) = h
At the wall (y = 0), the boundary condition is no slip: Z nh0and of the viscous  ⌧ are
Z hstress ✓ dy 2 ◆ µ dx
(c) The components of the normal 2
0
S ⇢g dh y
vx (y = 0)Q ==0 u(y) dy2 = The y +
@v velocity
yh
3 is2therefore
@v profile
@vy(5.23)
3 dy.
dvx0 given by:
0 2 0 x µ x µ dx
+ 0 2
6 7 6 7
n = [0, 1] , ⌧ = µ ru + ruT = µ 4 @v @x@v @y @v @x 5 = 4 dv dy 5 . (5.21)
All calculations
(d) First integrate equation (5.17) once: done, x
+
y
2
y x
0
@y @x @y dy u(y) =
2
dvx ⇢g sin ✓ Sh0 ⇢g dh 3
The condition of zero-traction ty =
= +n A.· ⌧ = 0 at y Q
== h therefore is writtenh
(5.24)
dy µ 2µ 3µ dx 0
dv Note that the slope of the free surface rem
The integration constant A is obtained using boundary condition (5.22): x (y = h) = 0 (5.22)
(f) At steady state, we must have Q = dy0 (zero flow rate across any vertical cross-secti
(e) The volumetric flow rate per unit widt
dvx ⇢g sin ✓ ⇢g sin ✓
us
(h)to
= solve for the
h + slope
A = 0,of thei.e. free
A =surface: h. (5.25)
dy At the wall (yµ = 0), the boundary condition is no µ slip: Z h0
Now, we integrate equation (5.24) once more: Q= u(y) (5.23)
dy =
Sh20 ⇢g
vx (ydh
= 0) = 0 dh 3S 0
h30 = 0, i.e. =
⇢g sin ✓ 2 ⇢g sin 2µ
✓ 3µ dx dx 2⇢gh0
(d)vFirst
x (y)integrate
= equation
y +(5.17) once:
hy + B. All calculations done,
(5.26)
2µ µ
Substituting this into equation (6.22) dvx yields⇢g the
sin ✓final expression for the velocity profile
(5.24)Q
Applying the no-slip condition at y = 0 shows that B = 0. Therefore the = velocity profile
y + [Link] given by
dy µ
✓ ◆ 3 ⇣ y
Sh0 y (5.22): ⌘
The integration constant
⇢g sinA✓ is obtained
y 2 using boundary condition 0
vx (y) = yh u(y) (f)
= At steady3state, 2 must have Q = 0 (
(5.27)
µ dv 2⇢g sin ✓
2
4µ h0 h⇢g sinwe

h + A us
= 0,to solve
i.e. Afor
= the slope h. of the free surface
x
(h) = (5.25)
dy µ µ
The velocity profile is parabolicxas shown below:
Now,y we integrate equation (5.24) once more: 1 Sh2 ⇢g d
0
S ⇢g sin ✓ 2 ⇢g sin ✓ 2µ 3µ d
vx (y) = y + hy + B. (5.26)
u(y) 2µ 0 µ
Substituting
Applying the no-slip condition at y = 0 shows that B this
= 0. Therefore the into equation
velocity profile is(6.22) yield
given by
1-
✓ ◆
⇢g sin ✓ y2
vx (y) = yh (5.27)
u(y
µ2- 2

x3-
Solution 5.7. Pipe flow of two immiscible liquids
y
01 (a) The momentum
8 6 equation
4 in cylindrical
2 coordinates
0 isSwritten
u(y)
22 5 Solutions

(e) We obtain the flow rate Q (per unit length in the z direction) by integrating the velocity profile:
Z h Z h ✓ ◆
⇢g sin ✓ y2
Q= vx (y) dy = yh dy. (5.28)
0 0 µ 2

All calculations done,


⇢gh3 sin ✓
Q= (5.29)

Solution 5.7. Pipe flow of two immiscible liquids


(a) The momentum equation in cylindrical coordinates is written
 ✓ ◆
@vr @vr 1 @vr v✓ @vr @p
⇢ + vr + v✓ + vz =
@t @r r @✓ r @z @r
 ✓ ◆
@ 1 @ 1 @ 2 vr @ 2 vr 2 @v✓
+µ (rvr ) + 2 2
+ ,
@r r @r r @✓ @z 2 r2 @✓
(5.30)
 ✓ ◆
@v✓ @v✓ 1 @v✓ vr @v✓ 1 @p
⇢ + vr + v✓ + + vz =
@t @r r @✓ r @z r @✓
 ✓ ◆
@ 1 @ 1 @ 2 v✓ @ 2 vr 2 @vr
+µ (rv✓ ) + 2 2
+ + 2 ,
@r r @r r @✓ @z 2 r @✓
(5.31)

@vz @vz v✓ @vz @vz @p
⇢ + vr + + vz =
@t @r r @✓ @z @z
 ✓ ◆
1 @ @vz 1 @ 2 vz @ 2 vz
+µ r + 2 + . (5.32)
r @r @r r @✓2 @z 2

For a velocity field of the form v(r, ✓, z, t) = v(r)ez , these simplify to


✓ ◆
@p @p @p µ d dv
= = 0, and = r . (5.33)
@r @✓ @z r dr dr

Therefore p = p(z), and we let dp/dz = G, which must be a constant. The two velocity profiles v1 (r)
and v2 (r) therefore both satisfy:
✓ ◆
1 d dv G
r = (5.34)
r dr dr µ
Multiply equation (5.34) by r and integrate:

dv G r2
r = + A. (5.35)
dr µ 2
Now divide equation (5.35) and integrate once more:
G 2
v(r) = r + A ln r + B (5.36)

Therefore the two velocity profiles are of the form
5 Solutions 23

G 2 G 2
v1 (r) = r + A1 ln r + B1 v2 (r) = r + A2 ln r + B2 (5.37)
4µ 4µ

where the four constants A1 , B1 , A2 and B2 have yet to be determined.


(b) At the pipe wall, the no-slip condition applies: v2 (r = R) = 0
(c) At the interface between the two liquids, the boundary conditions are continuity of the velocity and
shear stress ⌧rz = µ(dv/dr):

dv1 dv2
v1 [r = (1 ")R] = v2 [r = (1 ")R], and µ1 [r = (1 ")R] = µ2 [r = (1 ")R]. (5.38)
dr dr
We use these various boundary conditions to solve for the unknown integration constants:
• First, we require that v1 (r) remain finite at the centerline r = 0, which implies A1 = 0 .
• Next, the continuity of the shear stress at r = (1 ")R gives
✓ ◆ ✓ ◆
2G 2G A2
µ1 (1 ")R + 0 = µ2 (1 ")R + , (5.39)
4µ1 4µ2 (1 ")R

which shows that A2 = 0 .


• The no-slip condition at r = R gives

GR2 GR2
v2 (R) = + B2 = 0, i.e. B2 = (5.40)
4µ2 4µ2

• Finally, continuity of the velocity at r = (1 ")R is written

G G GR2
(1 ")2 R2 + B1 = (1 ")2 R2 , (5.41)
4µ1 4µ2 4µ2
from which
G 2 G 2
B1 = R [(1 ")2 1] R (1 ")2 (5.42)
4µ2 4µ1

(d) We now that assume that " ⌧ 1, in which case (1 ")2 ⇡ 1 2" (where we neglect terms of order
"2 ). Then,
G 2 G 2
B1 ⇡ R " R (1 2"), (5.43)
2µ2 4µ1
and the approximate expressions for the velocity profiles are
 ✓ ◆
G µ1 G 2
v1 (r) ⇡ r2 R 2
2"R 2
1 and v2 = (r R2 ) (5.44)
4µ1 µ2 4µ2

(e) From equation (5.44), it is obvious that the presence of liquid 2 on the periphery of the pipe increases
the speed of liquid 1 and therefore its flow rate if (µ1 /µ2 1) > 0 (recall that G < 0 for a positive
flow rate). This implies that µ2 < µ1 , i.e. the fluid on the periphery should be less viscous. This makes
physical sense: the less viscous liquid on the periphery reduces the viscous drag near the wall, where it is
the highest owing the large shear rate there. This could be made more quantitative by actually calculating
the flow rates.
• Finally, continuity of the velocity at r = (1 ")R is written
• Finally, continuity of the velocity at r = (1 ")R is written

G G GR2 G G GR2
(1 ")2 R2 + B1 = (1 ")2 R2 , (1 ")2(6.38)
R 2 + B1 = (1 ")2 R2 , (6.38)
4µ1 4µ2 4µ2 4µ1 4µ2 4µ2
from which from which
G 2 2 G 2 G 2 G 2
B1 = R [(1 ") 1] R (1 ")2 B1 = (6.39)
R [(1 ")2 1] R (1 ")2 (6.39)
24 4µ2 4µ1 4µ2 4µ1 5 Solutions

(d) We now that assume that


The"⌧ 1, in which
velocity casein(1the two
profiles 2
")(d)⇡cases
We1 now (where
µ2"2 that we µ
1 and
< µassume neglect
that ⌧
2 >" µterms
1, inofwhich
1 have order
the following
case (1 shapes:
")2 ⇡ 1 2" (where we neglect terms of order
"2 ). Then, "2 ). Then,
G 2 G 2 G G 2
B1 ⇡ R " R (1 2"), B1 ⇡ (6.40)R2 " R (1 2"), (6.40)
2µ2 1.0
4µ1 2µ2 4µ1 1.0

and the approximate expressions for the velocity profiles are


and the approximate expressions for the velocity profiles are
0.5
 ✓ ◆ 
0.5
✓ ◆
G µ1 G G µ1 G 2
v1 (r) ⇡ r2 R 2
2"R 2
1 and v2 v=1 (r) ⇡(r2 Rr2 ) R2 (6.41)
2"R 2
1 and v2 = (r R2 ) (6.41)
4µ1 µ2 4µ2 4µ1 µ2 4µ2
0.0 0.0

(e) From equation (6.41), it is obvious that the presence of (e)


liquid
From2 on the periphery
equation (6.41), itofisthe pipe increases
obvious that the presence of liquid 2 on the periphery of the pipe increases
the speed of liquid 1 and therefore its flow rate if (µ1 /µ2the speed
-0.5
1) > of 0 (recall
liquid 1that < 0 for its
andGtherefore a positive
flow rate if (µ1 /µ2 1) > 0 (recall that G < 0 for a positive
-0.5

flow rate). This implies that µ2 < µ1 , i.e. the fluid on theflow
periphery should
rate). This be less
implies viscous.
that µ2 < µ This makes
1 , i.e. the fluid on the periphery should be less viscous. This makes
physical sense this is allows to reduce the viscous drag near
-1.0
the wall,
physical sense where it allows
this is is the highest owing
to reduce the
the viscous
-1.0
drag near the wall, where it is the highest owing the
large shear rate there. This can be made more quantitative by
0.0
actually
large
0.1
calculating
shear rate
0.2 0.3
the can
there. This
0.4
flowberates.
0.5
made more quantitative by actually calculating the flow rates.
0.6 0.0 0.5 1.0 1.5
The velocity profiles in the two cases µ2 < µ1 and µ2 > µThe 1 arevelocity
as follows:
profiles in the two cases µ2 < µ1 and µ2 > µ1 are as follows:

Solution 5.8. Flow of paint down a wall


(a) The x component of the momentum equation is written
✓ ◆
@ux @ux @ux @ux @p @⌧xx @⌧yx @⌧zx
⇢ + ux + uy + uz = ⇢g + + + . (5.45)
@t @x @y @z @x @x @y @z

Assuming steady parallel flow of the form u = ux (y)ex , this simplifies to

d⌧
+ ⇢g = 0 where ⌧ = ⌧yx (y). (5.46)
dy

Using the boundary condition ⌧ (H) = 0, this integrates to

⌧ (y) = ⇢g(H y) (5.47)

The shear stress profile is linear: it is maximum near at the wall and zero at the free surface:

(b) To solve for the velocity, we need to distinguish several cases:


• Case I: ⌧0 ⇢gH
In this case, ⌧ (y)  ⌧0 for all 0  y  H, i.e. the shear stress is everywhere less than the yield stress.
This implies that = dux /dy = 0, which integrates to ux (y) = A. Since the no-slip condition applies
5 Solutions 25

at the wall, we must have A = 0, therefore there is not flow:

ux (y) = 0 for all y. (5.48)

• Case II: 0  ⌧ < ⇢gH.


In this case, we have ⌧ (y) ⌧0 when 0  y  a = H ⌧0 /⇢g (near the wall), but ⌧ (y)  ⌧0 when
a  y  H (near the free surface), and both intervals should be treated separately:

When 0  y  a, we use the relation ⌧ = µ + ⌧0 , i.e.


dux
⇢gH ⇢gy = µ + ⌧0 , (5.49)
dy
which is also written
dux ⇢g
= (a y). (5.50)
dy µ
This integrates to ✓ ◆
⇢g y2
ux (y) = ay + B, (5.51)
µ 2
and applying the no-slip condition at y = 0 gives: ux (0) = B = 0.
Now, when a  y  H, we have = 0 and therefore ux = C. The continuity of the velocity at y = a
requires ✓ ◆
⇢g a2 ⇢ga2
C= a2 = . (5.52)
µ 2 2µ

To summarize, the velocity profile when ⌧0 < ⇢gH is given by:


8 ✓ 2

> ⇢g ay y
>
< if 0  y  a,
ux (y) = µ 2 (5.53)
>
> ⇢ga2
: if a  y  H .

where a = H ⌧0 /⇢g. Fluid deformation only occurs near the wall where the shear stress exceeds the
yield stress (see figure below for a sketch of the velocity profile).

(c) To obtain the flow rate, we integrate the velocity across the film thickness. We again distinguish two
cases:
26 5 Solutions

• Case I: ⌧0 ⇢gH (or ↵ 1). In this case, ux (y) = 0, therefore Q = 0 .


• Case II: 0  ⌧ < ⇢gH (or ↵  1). We calculate the flow rate as
Z H Z a ✓ ◆ Z H
⇢g y2 ⇢ga2
Q= ux (y) dy = ay dy + dy. (5.54)
0 0 µ 2 a 2µ

Performing the integrals and using that a = H(1 ↵), this gives after simplifications

⇢gH 3 ⇣ ↵⌘
Q= (1 ↵)2 1 + (5.55)
3µ 2
Problem Set 6 - Solutions

Problem 6.1
Problem 6.2
Problem 6.3
Problem 6.4
Problem 6.5
Problem 7.6 (continued)

You might also like