Understanding Continuum in Fluid Mechanics
Understanding Continuum in Fluid Mechanics
[Link]
Peeter Joot — [Link]@[Link]
Revision [Link] commit 16e75b52d6a80e9861f94495cfd62fc4b44aa129
Apr/25/2012 [Link]
Keywords: PHY454H1S, PHY454H1, strain, displacement vector, stress, constitutive relation,
Navier-Stokes, incompressible fluid, capillary constant, Surface tension, Laplace pressure, channel flow,
pressure gradient, shear flow, Poiseuille flow, film flow, flux, traction vector, no-slip condition, interface,
hydrostatics, Buoyancy force, mass conservation, unit normal, unit tangent, surface, Reynold’s number,
scaling, Non-dimensionalisation, Lagrangian view, Eulerian view, similarity variable, Boundary layer,
Bernoulli equation, Blassius problem, Singular perturbation, oscillatory unstability, marginal unstability,
neutral stability, thermal stability, Rayleigh number, Prandtl number, Rayleigh-Benard problem
Contents
1 Motivation. 2
2 Vector displacements. 2
4 Conservation of mass 3
5 Constitutive relation 3
7 No slip condition 4
10 Hydrostatics. 7
12 Bernoulli equation. 8
1
13 Surface tension. 9
13.1 Surfaces, normals and tangents. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9
13.2 Laplace pressure. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10
13.3 Surface tension gradients. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10
13.4 Surface tension for a spherical bubble. . . . . . . . . . . . . . . . . . . . . . . . . . . . 10
13.5 A sample problem. The meniscus curve. . . . . . . . . . . . . . . . . . . . . . . . . . . 11
16 Boundary layers. 14
16.1 Impulsive flow. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14
16.2 Oscillatory flow. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14
16.3 Blassius problem (boundary layer thickness in flow over plate). . . . . . . . . . . . . 15
18 Stability. 16
18.1 Thermal stability: Rayleigh-Benard problem. . . . . . . . . . . . . . . . . . . . . . . . 17
1. Motivation.
Review of key ideas and equations from the fluid dynamics portion of the class.
2. Vector displacements.
Those portions of the theory of elasticity that we did cover have the appearance of providing
some logical context for the derivation of the Navier-Stokes equation. Our starting point is almost
identical, but we now look at displacements that vary with time, forming
1
Ωk = ∂i u j eijk (4a)
2
2
ωij = −Ωk eijk , (4b)
or
ω = ∇×u (5a)
1
Ω= (ω) a e a . (5b)
2
We were then able to put our displacement differential into a partial vector form
We are able to identity the divergence of the displacement as the relative change in volume per
unit time in terms of the strain tensor trace (in the basis for which the strain is diagonal at a given
point)
dV 0 − dV
= ∇ · u. (7)
dVδt
4. Conservation of mass
dρ ∂ρ
= + u · ∇ρ = −ρ∇ · u. (9)
dt ∂t
An important consequence of this is that for incompressible fluids (the only types of fluids
considered in this course) the divergence of the displacement ∇ · u = 0.
5. Constitutive relation
We consider only Newtonian fluids, for which the stress is linearly related to the strain. We
will model fluids as disjoint sets of hydrostatic materials for which the constitutive relation was
previously found to be
3
6. Conservation of momentum (Navier-Stokes)
dui ∂σij
ρ = + ρ fi , (11)
dt ∂x j
where f i are the components of the external (body) forces per unit volume acting on the fluid.
Utilizing the constitutive relation and explicitly evaluating the stress tensor divergence ∂σij /∂x j
we find
du ∂u
ρ =ρ + ρ(u · ∇)u = −∇ p + µ∇2 u + µ∇(∇ · u) + ρf. (12)
dt ∂t
Since we treat only incompressible fluids in this course we can decompose this into a pair of
equations
∂u
ρ + ρ(u · ∇)u = −∇ p + µ∇2 u + ρf. (13a)
∂t
∇·u = 0 (13b)
7. No slip condition
We’ll find in general that we have to solve for our boundary value conditions. One of the
important constraints that we have to do so will be a requirement (experimentally motivated) that
our velocities match at an interface. This was illustrated with a rocker tank video in class.
This is the no-slip condition, and includes a requirement that the fluid velocity at the boundary
of a non-moving surface is zero, and that the fluid velocity on the boundary of a moving surface
matches the rate of the surface itself.
For fluids A and B separated at an interface with unit normal n̂ and unit tangent τ̂ we wrote
the no-slip condition as
u A · τ̂ = u B · τ̂ (14a)
u A · n̂ = u B · n̂. (14b)
For the problems we attempt, it will often be enough to consider only the tangential component
of the velocity.
As well as matching velocities, we have a force balance requirement at any interface. This will
be expressed in terms of the traction vector
T = ei σij n j = σ · n̂ (15)
where n̂ = n j ei is the normal pointing from the interface into the fluid (so the traction vector
represents the force of the interface on the fluid). When that interface is another fluid, we are able
to calculate the force of one fluid on the other.
In addition the the constraints provided by the no-slip condition, we’ll often have to constrain
our solutions according to the equality of the tangential components of the traction vector
4
τi (σij n j ) A
= τi (σij n j ) B , (16)
We’ll sometimes also have to consider, especially when solving for the pressure, the force bal-
ance for the normal component of the traction vector at the interface too
ni (σij n j ) A
= ni (σij n j ) B . (17)
As well as having a messy non-linear PDE to start with, our boundary value constraints can
be very complicated, making the subject rich and tricky.
8.1. Flux
A number of problems we did asked for the flux rate. A slightly more sensible physical quan-
tity is the mass flux, which adds the density into the mix
dm dV
Z Z Z
=ρ =ρ (u · n̂)dA (18)
dt dt
• shear flow
d2 u G
2
=− (23)
dz µ
with solution
5
G 2
u(z) = (h − z2 ) + A(z + h) + B. (24)
2µ
Application of the no-slip velocity matching constraint gives us in short order
G 2 2 1
u(z) = (h − z ) + v 1 − (z + h) . (25)
2µ 2h
With v = 0 this is the channel flow solution, and with G = 0 this is the shearing flow solution.
Having solved for the velocity at any height, we can also solve for the mass or volume flux
through a slice of the channel. For the mass flux ρQ per unit time (given volume flux Q)
dm dV
Z Z Z
= ρ = ρ(∆A) u · τ̂, (26)
dt dt
we find
2Gh3
ρQ = ρ(∆y) + hv . (27)
3µ
We can also calculate the force of the boundaries on the fluid. For example, the force per unit
volume of the boundary at z = ±h on the fluid is found by calculating the tangential component
of the traction vector taken with normal n̂ = ∓ẑ. That tangent vector is found to be
∂u
σ · (±n̂) = − pẑ ± 2µei eij δj3 = − pẑ ± x̂µ
. (28)
∂z
The tangential component is the x̂ component evaluated at z = ±h, so for the lower and upper
interfaces we have
vµ
(σ · n̂) · x̂|z=−h = − G (−h) − (29)
2h
vµ
(σ · −n̂) · x̂|z=+h = − G (+h) + , (30)
2h
so the force per unit area that the boundary applies to the fluid is
vµ
force per unit length of lower interface on fluid = L Gh − (31)
2h
vµ
force per unit length of upper interface on fluid = L − Gh + . (32)
2h
Does the sign of the velocity term make sense? Let’s consider the case where we have a zero
pressure gradient and look at the lower interface. This is the force of the interface on the fluid, so
the force of the fluid on the interface would have the opposite sign
vµ
. (33)
2h
This does seem reasonable. Our fluid flowing along with a positive velocity is imparting a
force on what it is flowing over in the same direction.
6
10. Hydrostatics.
We covered hydrostatics as a separate topic, where it was argued that the pressure p in a fluid,
given atmospheric pressure p a and height from the surface was
p = p a + ρgh. (34)
As noted below in the surface tension problem, this is also a consequence of Navier-Stokes for
u = 0 (following from 0 = −∇ p + ρg).
We noted that replacing the a mass of water with something of equal density would not change
the non-dynamics of the situation. We then went on to define Buoyancy force, the difference in
weight of the equivalent volume of fluid and the weight of the object.
It was noted that mass conservation provides a relationship between the flow rates through
apertures in a closed pipe, since we must have
ρ1 A1 v1 = ρ2 A2 v2 , (35)
and therefore for incompressible fluids
A1 v1 = A2 v2 . (36)
So if A1 > A2 we must have v1 < v2 .
∆v = gt (37)
1 2
∆z = gt + v0 t (38)
2
Solving for v f = v0 + ∆v, we find
s
2g∆z
v f = v0 1+ . (39)
v20
Mass conservation gives us
7
Figure 1: Tap flow measurement.
or
2g∆z −1/4
s
v0
r (∆z) = R = R 1+ 2 . (41)
vf v0
For the image above I measured a flow rate of about 250 ml in 10 seconds. With that, plus
the measured radii at 0 and 6cm, I calculated that the average fluid velocity was 0.9m/s, vs a free
fall rate increase of 1.3m/s. Not the best match in the world, but that’s to be expected since the
velocity has been considered uniform throughout the stream profile, which would not actually be
the case. A proper treatment would also have to treat viscosity and surface tension.
In figure (2) is a plot of the measured radial distance compared to what was computed with
41. The blue line is the measured width of the stream as measured, the red is a polynomial curve
fitted to the raw data, and the green is the computed curve above.
g = −∇χ (42)
and utilizing the vector identity
1 2
(u · ∇)u = ∇ × (∇ × u) + ∇ u , (43)
2
we are able to show that the steady state, irrotational, non-viscous Navier-Stokes equation
takes the form
p 1 2
∇ + χ + u = 0. (44)
ρ 2
8
Figure 2: Comparison of measured stream radii and calculated.
or
p 1
+ χ + u2 = constant (45)
ρ 2
This is the Bernoulli equation, and the constants introduce the concept of streamline.
FIXME: I think this could probably be used to get a better idea what the tap stream radius is,
than the method used above. Consider the streamline along the outermost surface. That way you
don’t have to assume that the flow is at the average velocity rate uniformly throughout the stream.
Try this later.
φ = z − h( x, t) = 0. (46)
Computing the gradient we find
∂h
∇φ = −x̂ ẑ. (47)
∂x
Recalling that the gradient is normal to the surface we can compute the unit normal and unit
tangent vectors
1 ∂h
n̂ = r 2 − ∂x , 1 (48a)
∂h
1 + ∂x
1 ∂h
τ̂ = r 2 1, ∂x (48b)
∂h
1 + ∂x
9
13.2. Laplace pressure.
We covered some aspects of this topic in class. [1] covers this topic in typical fairly hard to com-
prehend) detail, but there’s lots of valuable info there. §2.4.9-2.4.10 of [2] has small section that’s
a bit easier to understand, with less detail. Recommended in that text is the “Surface Tension in
Fluid Mechanics” movie which can be found on youtube in three parts [Link]
[Link] HKw, [Link] which is very interesting and en-
tertaining to watch.
It was argued in class that the traction vector differences at the surfaces between a pair of fluids
have the form
σ
t2 − t1 = − n̂ − ∇ I σ (49)
2R
where ∇ I = ∇ − n̂(n̂ · ∇) is the tangential (interfacial) gradient, σ is the surface tension, a
force per unit length value, and R is the radius of curvature.
In static equilibrium where t = − pn̂ (since σ = 0 if u = 0), then dotting with n̂ we must then
have
σ
p2 − p1 = (50)
2R
2σ
∆p = , (52)
R
is only for spherical objects that have a single radius of curvature. This formula can in fact be
derived with a simple physical argument, stating that the force generated by the surface tension σ
along the equator of a bubble (as in 3), in a fluid would be balanced by the difference in pressure
times the area of that equatorial cross section. That is
10
Figure 3: Spherical bubble in liquid.
0 = −∇ p1 + ρg. (55)
With g = − gẑ we have
11
∂p1
0=− − ρg, (56)
∂z
or
ρgz 1 1
constant − = + . (59)
σ R1 R2
One of our axis of curvature directions is directly along the y axis so that curvature is zero
1/R1 = 0. We can fix the constant by noting that at x = ∞, z = 0, we have no curvature 1/R2 = 0.
This gives
constant − 0 = 0 + 0. (60)
That leaves just the second curvature to determine. For a curve z = z( x ) our absolute curva-
ture, according to [3] is
1 |z00 |
= . (61)
R2 (1 + (z0 )2 )3/2
Now we have to fix the sign. I didn’t recall any sort of notion of a signed radius of curvature,
but there’s a blurb about it on the curvature article above, including a nice illustration of signed
radius of curvatures can be found in this wikipedia radius of curvature figure for a Lemniscate.
Following that definition for a curve such as z( x ) = (1 − x )2 we’d have a positive curvature,
but the text explicitly points out that the curvatures are will be set positive if pointing into the
medium. For us to point the normal into the medium as in the figure, we have to invert the sign,
so our equation to solve for z is given by
ρgz z00
− =− . (62)
σ (1 + (z0 )2 )3/2
The text introduces the capillary constant
p
a= 2σ/gρ. (63)
Using that capillary constant a to tidy up a bit and multiplying by a z0 integrating factor we
have
2zz0 z00 z0
− = − , (64)
a2 (1 + (z0 )2 )3/2
we can integrate to find
z2 1
A− = . (65)
a 2 (1 + (z0 )2 )1/2
12
Again for x = ∞ we have z = 0, z0 = 0, so A = 1. Rearranging we have
−1/2
1
Z Z
dx = dz −1 . (66)
(1 − z2 /a2 )2
Integrating this with Mathematica I get
√
a 2a − 4a 2 − 2z2 sgn( a − z )
p a
x − x0 = 2a2 − z2 sgn( a − z) + √ ln . (67)
2 z
− cos θ
z 0 (0) = = − cot θ, (68)
sin θ
but we don’t have an explicit function for z.
u → Uu0 (69)
2 0
p → ρU p (70)
L
t → t0 (71)
U
1
∇ → ∇0 (72)
L
we can put Navier-Stokes in dimensionless form
∂u0 1
0
+ (u0 · ∇0 )u0 = ∇0 p0 + ∇0 u0 . (73)
∂t R
Here R is Reynold’s number
LU
R= (74)
ν
A relatively high or low Reynold’s number will effect whether viscous or inertial effects dom-
inate
The importance of examining where one of these effects can dominate was clear in the Blas-
sius problem, where doing so allowed for an analytic solution that would not have been possible
otherwise.
13
15. Eulerian and Lagrangian.
We defined
14
This was a bit more obvious as aq
boundary layer illustration since we see the exponential drop
2ν
off with every distance multiple of Ω.
∂u ∂u dU ∂2 u
u +v =U +ν 2 (84a)
∂x ∂y dx ∂y
∂p
=0 (84b)
∂y
∂u ∂v
+ =0 (84c)
∂x ∂y
With boundary conditions
U ( x, 0) = 0 (85)
U ( x, ∞) = U ( x ) = U0 (86)
V ( x, 0) = 0 (87)
f 000 + f f 00 = 0. (92)
This was done much more clearly in [4] and I worked this problem myself with a hybrid
approach (non-dimensionalising as done in class).
FIXME: the end result of this is a plot (a nice one can be found in [5]). That plot ends up being
one that’s done in terms of the similarity variable η. It’s not clear to me how this translates into an
actual velocity profile. Should plot these out myself to get a feel for things.
15
17. Singular perturbation theory.
1 2
(u · ∇)u = −∇ p + ∇ u (93)
Re
where the inverse of Reynold’s number
UL
Re = (94)
ν
can potentially get very small. That introduces an ill-conditioning into the problems that can
make life more interesting.
We looked at a couple of simple LDE systems that had this sort of ill conditioning. One of them
was
du
+ u = x,
e (95)
dx
for which the exact solution was found to be
1 du
e + u = ey. (98)
e dy
This gives us
du
+ u ≈ 0, (99)
dy
for which we find
18. Stability.
16
18.1. Thermal stability: Rayleigh-Benard problem.
We considered the Rayleigh-Benard problem, looking at thermal effects in a cavity. Assuming
perturbations of the form
u = ubase + δu = 0 + δu
p = ps + δp
ρ = ρs + δρ
∂δp ∂2 δρ
∇2 = −g 2 , (104)
∂z ∂z
from which we apply back to 103 and take just the z component to find
∂ 2 1 ∂δp δρ
− ν∇ δw = − − g (105)
∂t ρs ∂z ρs
With an assumption that density change and temperature are linearly related
and operating with the Laplacian we end up with a relation that follows from the momentum
balance equation
2
∂2
∂ 2 2 ∂
− ν∇ ∇ δw = gα + 2 δT. (107)
∂t ∂x2 ∂y
We also applied our perturbation to the energy balance equation
∂T
+ (u · ∇) T = κ ∇2 T (108)
∂t
We determined that the base state temperature obeyed
∂2
κ Ts = 0, (109)
∂z2
with solution
∆T
Ts = T0 − z. (110)
d
This and application of the perturbation gave us
17
∂δT
+ δu · ∇ Ts = κ ∇2 δT. (111)
∂t
We used this to non-dimensionalize with
x, y, z with d
t with d2 /ν
(112)
δw with κ/d
δT with ∆T
And found (primes dropped)
2
∂2
2 ∂ 2 ∂
∇ − ∇ δw = R + 2 δT (113a)
∂t ∂x2 ∂y
∂
Pr − ∇2 δT = δw, (113b)
∂t
where we’ve introduced the Rayleigh number and Prandtl number’s
gα∆Td3
R= , (114a)
νκ
ν
Pr = (114b)
κ
We were able to construct some approximate solutions for a problem similar to these equations
using an assumed solution form
Using these we are able to show that our PDEs are similar to that of
w = D2 w = D4 w = 0, (117)
z=0,1
where D = ∂/∂z. Using the trig solutions that fall out of this we were able to find the constraint
2
n2 π 2 + k 2 + σ n2 π 2 + k 2 −Rk2
0= , (118)
−1 n 2 π 2 + k 2 + Pr σ
which for σ = 0, this gives us the critical value for the Rayleigh number
( k 2 + n2 π 2 )3
R= , (119)
k2
which is the boundary for thermal stability or instability.
The end result was a lot of manipulation for which we didn’t do any sort of applied problems.
It looks like a theory that requires a lot of study to do anything useful with, so my expectation is
that it won’t be covered in detail on the exam. Having some problems to know why we spent two
days on it in class would have been nice.
18
References
[1] L.D. Landau and E.M. Lifshitz. A Course in Theoretical Physics-Fluid Mechanics. Pergamon Press
Ltd., 1987. 13.2, 13.5
[3] Wikipedia. Curvature — wikipedia, the free encyclopedia [online]. 2012. [Online; ac-
cessed 25-April-2012]. Available from: [Link]
Curvature&oldid=488021394. 13.5
[4] D.J. Acheson. Elementary fluid dynamics. Oxford University Press, USA, 1990. 16.3
[5] Wikipedia. Blasius boundary layer — wikipedia, the free encyclopedia [online]. 2012. [Online;
accessed 28-March-2012]. Available from: [Link]
Blasius_boundary_layer&oldid=480776115. 16.3
19