Fluid Dynamics 3 - 2015/2016

Jens Eggers
Preliminaries
Course information
• Lecturer: Prof. Jens Eggers, Room SM2.3

• Timetable: Weeks 1-12, Tuesday 1.00 (room 1.18 Queen’s building), Thursday 1.00
(room 1.18 Queen’s bulding) and Friday 10.00 (Physics 3.21). Course made up of
32 lectures.

• Drop-in sessions (formerly office hours): Monday 12:00 in my room, 2.3

• Prerequisites: Mechanics 1, APDE2 and Calc2. Need ideas from vector calculus,
complex function theory, separation solutions and PDE’s.

• Homework: Questions from 10 worksheets will be set and marked during the course.
Homework sheets will be handed out each Friday, starting in the first week. Solutions
to be returned the following Friday in the box marked “Fluids 3”.

• Web: Standard unit description includes detailed course information. Lecture notes
will be posted on Blackboard, along with homework and solution sheets.

• Lectures: There is no need to take notes during the lecture, as all material relevant
for the exam will be put on the web. The main purpose of the lecture is the live
development of the material, and a chance for you to ask questions!

Recommended texts
1. A.R. Paterson, A First Course in Fluid Dynamics, Cambridge University Press.
(The recommended text to complement this course - costs ≈ £50 from Amazon;
there are 6 copies in Queen’s building Library and 3 copies in the Physics Library)

2. D.J. Acheson, Elementary Fluid Dynamics. Oxford University Press

3. L.D. Landau and E.M. Lifshitz, Fluid Mechanics. Butterworth Heinemann

Films
There is a very good series of educational films on Fluid Mechanics available on YouTube,
produced by the National Committee for Fluid Mechanics Films in the US in the 1960’s.
Each film is also accompanied by a set of notes. I recommend them highly, and will point
out the appropriate ones throughout this course.
 The following 3 sections are useful for the course. For the purposes of an examination, I
would expect you to know the definition of div, grad, curl and the Laplacian in Cartesians
and grad and the Laplacian in plane polars. Other definitions would be provided.

Revision of vector operations

Let u = (u1 , u2 , u3 ), v = (v1 , v2 , v3 ) be Cartesian vectors. Let φ(r) be a scalar function
and f (r) = (f1 (r), f2 (r), f3 (r)) a vector field of position r = (x, y, z) ≡ (x1 , x2 , x3 ). Then

• The dot product is u · v = u1 v1 + u2 v2 + u3 v3

• The cross (or vector) product is u×v = (u2 v3 −u3 v2 )r̂+(u3 v1 −u1 v3 )ŷ+(u1 v2 −u2 v1 )ẑ
 
∂φ ∂φ ∂φ
• The gradient is ∇φ = , ,
∂x1 ∂x2 ∂x3
∂f1 ∂f2 ∂f3
• The divergence is ∇ · f = + +
∂x1 ∂x2 ∂x3
     
∂f3 ∂f2 ∂f1 ∂f3 ∂f2 ∂f1
• The curl is ∇ × f = − r̂ + − ŷ + − ẑ
∂x2 ∂x3 ∂x3 ∂x1 ∂x1 ∂x2
∂ 2φ ∂ 2φ ∂ 2φ
• The Laplacian is ∇2 φ = ∇ · ∇φ = + +
∂x21 ∂x22 ∂x23
     
∂f1 ∂f1 ∂f1 ∂f2 ∂f2 ∂f2 ∂f3 ∂f3 ∂f3
• (f ·∇)f = f1 ∂x1
+ f2 ∂x2
+ f3 ∂x 3
r̂+ f 1 ∂x1 + f 2 ∂x2 + f 3 ∂x3 ŷ+ f 1 ∂x1 + f 2 ∂x2 + f 3 ∂x3 ẑ

Formulae in cylindrical polar coordinates

Coordinate system is r = (r, θ, z) where the relationship to Cartesians is x = r cos θ,
y = r sin θ. The unit vectors are r̂ = r̂ cos θ + ŷ sin θ, θ̂ = −r̂ sin θ + ŷ cos θ and ẑ. In the
following, f = (fr , fθ , fz ) ≡ fr r̂ + fθ θ̂ + fz ẑ.
∂φ 1 ∂φ ∂φ
• The gradient is ∇φ = r̂ + θ̂ + ẑ
∂r r ∂θ ∂z
1 ∂(rfr ) 1 ∂fθ ∂fz
• The divergence is ∇ · f = + +
r ∂r r ∂θ ∂z

1 r̂ rθ̂ ẑ
• The curl is ∇ × f = ∂/∂r ∂/∂θ ∂/∂z .
r
fr rfθ fz

∂ 2 φ 1 ∂φ 1 ∂ 2φ ∂ 2φ
• The Laplacian is ∇2 φ = + + + 2
∂r2 r ∂r r2 ∂θ2 ∂z
 
fθ2
• (f · ∇)f = fr ∂f fθ ∂fr ∂fr
r̂ + fr ∂f + frθ ∂f + fz ∂f f r fθ


∂r
r
+ r ∂θ
+ f z ∂z
− r ∂r
θ
∂θ
θ
∂z
θ
+ r
θ̂ +
fr ∂f + frθ ∂f + fz ∂f
z z z


∂r ∂θ ∂z
ẑ
Formulae in spherical polar coordinates
Coordinate system is r = (r, θ, ϕ) where the relationship to Cartesians is x = r sin θ cos ϕ,
y = r sin θ sin ϕ, z = r cos θ.
The unit vectors are r̂ = r̂ sin θ cos ϕ + ŷ sin θ sin ϕ + ẑ cos θ, θ̂ = r̂ cos θ cos ϕ +
ŷ cos θ sin ϕ − ẑ sin θ and ϕ̂ = −r̂ sin ϕ + ŷ cos ϕ.
In the following, f = (fr , fθ , fϕ ) ≡ fr r̂ + fθ θ̂ + fϕ ϕ̂.
∂φ 1 ∂φ 1 ∂φ
• The gradient is ∇φ = r̂ + θ̂ + ϕ̂
∂r r ∂θ r sin θ ∂ϕ
1 ∂(r2 fr ) 1 ∂(sin θfθ ) 1 ∂fϕ
• The divergence is ∇ · f = 2
+ +
r ∂r r sin θ ∂θ r sin θ ∂ϕ

1 r̂ rθ̂ r sin θϕ̂

• The curl is ∇ × f = 2
∂/∂r ∂/∂θ ∂/∂ϕ .
r sin θ
fr rfθ r sin θfϕ

∂ 2φ
   
2 1 ∂ 2 ∂φ 1 ∂ ∂φ 1
• The Laplacian is ∇ φ = 2 r + 2 sin θ + 2 2
r ∂r ∂r r sin θ ∂θ ∂θ r sin θ ∂ϕ2
fθ2 +fϕ2 fϕ2 cot θ
   
fϕ ∂fr fϕ ∂fθ
• (f ·∇)f = fr ∂f ∂r
r
+ fθ ∂fr
r ∂θ
+ r sin θ ∂ϕ
− r
r̂+ fr ∂f∂r
θ
+ frθ ∂f
∂θ
θ
+ r sin θ ∂ϕ
+ frrfθ − r
θ̂+
 
fr ∂f
∂r
ϕ
+ fθ ∂fϕ
r ∂θ
+ fϕ ∂fϕ
r sin θ ∂ϕ
+ fr fϕ
r
+ fθ fϕ cot θ
r
ϕ̂
1 Introduction & Basic ideas
1.1 What is a fluid ?
In this course we will treat the laws governing the motion of liquids and gases. A liquid
or gas is characterized by the fact that there is no preferred rest state for the parts it
is composed of. If a hole is made in a water bottle, the water will flow out. If a drop
is placed on a solid surface, it will spread. By contrast, a solid retains a memory of its
original state. If one deforms a piece of metal, it will relax back to its original state once
the force is no longer applied. Collectively, we will call liquids or gases, which share this
fluid property, “fluids”.
Another important property of liquids and gases is that they are “featureless”. Viewed
from a particular point in space, all directions are equivalent. Solids, on the other hand,
often have an internal (lattice) structure. As a result, it makes a difference in which
direction they are deformed relative to their internal structure.
Fluid dynamics is an example of ‘continuum’ mechanics:

Definition 1.1.1 (Continuum) A continuum is any medium whose state at a given

instant can be described in terms of a set of continuous functions of position r = (x, y, z).
E.g. Density, ρ(r, t), velocity, u(r, t), temperature T (r, t). (These functions usually also
depend upon time, t).

We know that this description fails if we observe matter on small enough length scales:
e.g. typical molecule size (∼ 10−9 m) or typical mean free path in a gas (∼ 10−7 m). The
miracle is that on a scale only slightly larger than that, all microscopic features can be
ignored, and we end up with a “universal” description of all things fluid. We will also
ignore all effects of incompressibility, so for the purposes of this course, liquids and gases
are the same thing. We will almost always speak of a fluid, but which can mean either a
liquid or a gas for the purposes of this course.
The motion of liquids and gases is governed by the same underlying prin-
ciples.

1.2 What we would like to do

In this lecture course, we will first develop an equation of motion for the velocity field
u(r, t), which gives the fluid velocity at any instant in time t, everywhere in space. This
equation (or set of equations) will necessarily have the form of a partial differential equa-
tion. It will be based on Newton’s equations of motion, but for a continuum of
particles, distributed over space. One effect we will neglect is the friction be-
tween fluid particles. The mathematical idealization of this situation is called
an “Ideal Fluid”.
With the equations in hand, it is down to our ability to deal with the mathematical
complexities of solving a partial differential equation (PDE) to solve physical problems.
Some examples of problems dealt with rather successfully using the concept of ideal fluids
are the following:
 Figure 1: A jet of water from a bottle.

Figure 1 shows a jet of water from a bottle. Both the efflux of the water and the
trajectory of the resulting jet are well described by ideal fluid theory, which we will
describe.

Figure 2: An airplane is held up by the lift generated by the wings.

Another spectacular success is the theory of flight. The ideal flow of air around a
wing is able to describe the lift necessary for flight, and much more. What is much more
difficult is the theory of drag. Inviscid theory suggests that there should be no energy
cost to flight at all!
 Figure 3: Water waves on the surface of a lake.

The last topic we will be covering in this course is the huge area of water waves, see
Fig 3. This includes waves from the scale of millimetres up to huge tsunami waves. The
absence of any solid boundary results in very little friction, so the ideal theory works very
well.

1.3 Lagrangian and Eulerian descriptions of the flow

We now begin to develop a dynamical description of fluid flow, which will lead us to
formulate a PDE for fluid motion, known as the Euler equation. Before we can do
that, we must understand the motion of fluids a little better. The description of motion
is called kinematics. In this chapter, we will deal with kinematics. I encourage you to
look at the film “Fluid Mechanics (Eulerian and Lagrangian description) parts 1-3” on
YouTube.
There are two very different ways of describing fluid motion, known as the Eulerian
and Lagrangian description. Ultimately, they are equivalent, as they describe the same
thing. However, they serve different purposes, so we need them both.

Definition 1.3.1 (Eulerian description of the flow) . This is what the stationary
observer sees. Choose a fixed point, r to measure, for e.g. the velocity u(r, t). This pro-
vides a spatial distribution of the flow at each instant in time. This is the way continuum
equations are usually formulated, and our equation of motion will indeed be an equation
for the Eulerian field u(r, t). If the flow is steady, then u does not depend on time, t:
u = u(r).

Why do we need anything else? The reason is that to make contact with Newton’s
equations, we need to describe the flow as a moving particle would see it. This is the

Definition 1.3.2 (Lagrangian description of the flow) The observer moves with the
fluid. Choose a fluid particle (for example, we can place a small drop of ink in the fluid),
and follow it through the fluid. Measuring its velocity at a given time, t gives its ‘La-
grangian velocity’. Now we describe the whole velocity field this way, by labeling all mate-
rial points. A convenient way of doing so is to choose an initial time t0 , and to label all
fluid particle by their position r = a at that time. Then at time t > t0 , the particle is at

r = r(a, t), where r(a, t0 ) = a.

Of course, r(a, t) is very interesting in its own right. For example, it describes the
course of a balloon, launched at time t0 and at position a into the atmosphere. The
Lagrangian velocity is defined as

dr
v(a, t) = . (1)
dt af ixed

By definition, it is the velocity of a particle going with the flow. This is precisely how
velocity is defined in Newtonian Mechanics. As we said earlier, the Eulerian and the
Lagrangian velocity fields contain the same information; the relation between the two is:

v(a, t) = u(r(a, t), t). (2)

An example: Consider logs flowing along a narrowing section of river. A fixed observer
measures the velocity by observing the velocity of logs at a given point in space. By
observing many logs at different positions, he will be able to obtain the entire Eulerian
velocity field u(r, t). Unless the flow conditions are changing, this field will be time-
independent, u(r, t) = u(r). Now imagine each log being ridden by a moving observer,
each of whom reports his velocity as time goes by. If all observers are labeled by their
position a at some reference time t0 , this will produce the Lagrangian velocity field v(a, t).
The logs travel with the fluid and will see the flow accelerating, as the river bed becomes
narrower. Thus v(a, t) is manifestly time-dependent although the Eulerian field is not!

Observer moving with the flow sees the logs

accelerate

Stationary observer sees logs passing at constant speed.

Figure 4: Stationary and co-moving observers of a flow.

Definition 1.3.3 (Two-dimensional flow) A flow is two-dimensional if it is indepen-

dent of one of its components (in some fixed frame of reference). E.g. u = (u, v, 0).

Example 1.3.4 (Hyperbolic flow) Let us consider a very simple model flow for the
river shown above. Consider the two-dimensional, stationary Eulerian velocity field u =
(u, v, 0), defined by
u = kx, v = −ky. (3)
Figure 5: Local velocities of the flow (3). The red lines are streamlines, and denote the
bank of the river.

As can be seen in Fig. 5, the flow speeds up as the river contracts, which we can think of
as being confined by the river banks, shown as the red lines.

1.4 Particle paths and streamlines

Definition 1.4.1 (Pathlines) The particle paths (pathlines) are the paths followed by
individual particles. They are determined by the solution of the differential equation:
dr
= u(r, t), with initial condition r(t0 ) = a = (a1 , a2 , a3 ) (4)
dt
where u is assumed given.
The system of equations (4) specifies a unique curve. For some (simple) u, can be
integrated using elementary methods, but in general not. Let u = (u(r, t), v(r, t), w(r, t)),
then in components, (4) is
dx

= u(x, y, z, t) 

dt 


dy

= v(x, y, z, t)
dt 

dz


= w(x, y, z, t) 

dt
with x(t0 ) = a1 , y(t0 ) = a2 , z(t0 ) = a3 .

Example 1.4.2 (Pathlines of a steady flow) Find the pathlines for the flow field
(3).
Then the equations become
ẋ = kx, ẏ = −ky.
The initial positions a = (a1 , a2 ) at time t0 = 0, are

x(t) = a1 ekt , y(t) = a2 e−kt ,

which describes the path of particles in the flow.

To find the shape of the curve followed by a particle in Example 1.4.2, we eliminate
time between x(t) and y(t), to find

xy = a1 a2 , (5)

which is the equation of a hyperbola. An example is shown in Fig. 5.

Using the relation (2) between Eulerian and Lagrangian fields, we find

v(a, t) = ka1 ekt , −ka2 e−kt .

In particular, the Lagrangian velocity field is indeed time-dependent, while the Eulerian
field was steady.
Now let us consider a case where the Eulerian field is also time dependent.

Example 1.4.3 (Pathlines of an unsteady flow)

Calculate the Pathlines of the two-dimensional flow

u = (1, t). (6)

The equations of motion are

ẋ = 1, ẏ = t,
with solution
t2
x = x0 + t, .
y = y0 +
2
For example the path of a particle with initial position a = (0, 1) would be (x, y) =
(t, 1 + t2 /2).
What is the shape of a particle path? Eliminating time, we find

(x − x0 )2
y = y0 + , (7)
2
which is a parabola, but with different coefficients, depending on the initial particle posi-
tion.
Now we come to the definition of a streamline, which describes the geometry of a flow
at a given time t0 . As a result, the pattern of streamlines will in general change in time.

Definition 1.4.4 (Streamline)

A streamline of a flow u(r, t) at a given instant in time t0 , is a curve which is every-
where parallel to u(r, t0 ).
 Thus, along a streamline,
1 dr
= u(r, t0 ), (8)
λ(s) ds
where s parameterizes the streamline r(s), and λ(s) is an arbitrary function. This illus-
trates that the magnitude of u does not matter for the definition of the streamline, only
its direction. To find all possible streamlines, we have to solve (8) with different initial
conditions r(s = 0) = r0 . Eliminating λ(s)ds from (8), we find

dx dy dz
= = (= λ(s)ds), t = t0 . (9)
u v w
The function λ(s) does not change the shape of the streamline, but only the parameteri-
zation. To complete the picture, we also add arrows pointing in the direction of the flow
along a streamline (see Fig. 6 below).
The streamline pattern will in general be different for different times t0 , since the veloc-
ity field u(r, t0 ) is different. For a steady flow, the pattern is of course time-independent
as well. Streamlines can be visualized by taking a short-time exposure of illuminated
particles in a flow. Taking an arbitrary starting point, one always continues in the local
flow direction.

Remark 1.4.5 (Streamlines and pathlines of steady flows)

For steady flows, streamlines and pathlines coincide. Since u does not depend on
time, we can equate the left hand sides of the defining equations (8) and (4) to obtain

dr 1 dr
= .
dt λ(s) ds

Thus if λ(s) is chosen such that

ds 1
= ,
dt λ(s)
and thus Z s
t= λ(s)ds,
s0

the two definitions become identical. In other words, time is a particular parameterization
of the streamline.

Example 1.4.6 (Streamlines of a steady flow)

As in Example 1.4.2, we take the flow u = (kx, −ky).
 y

0 x

Figure 6: The streamlines of the flow (3). Arrows are drawn for the case k > 0.

Since this is a steady flow, we can eliminate t from the particle path to obtain (5),
which must also be a streamline. Alternatively, we have
dx dy
=−
kx kv
according to (9). Integrating, we find ln |x| + ln |y| = C, or

|x||y| = eC ,

which is (5), the equation for a hyperbola. Any such hyperbola is a possible equation for
the river bank, which must be a streamline of the flow, see Fig. 5.

Example 1.4.7 (Streamlines of an unsteady flow)

y t=0 y
t=1

0 x 0 x

Figure 7: The streamlines of the flow (10). The pattern changes in time.

As in Example 1.4.3, we take u = (1, t).

Thus the equation (9) for the streamline becomes
dy
dx = ,
t
with solution
y = tx + y0 . (10)
Thus streamlines of the flow are straight lines, but whose slope depends on time, see
Fig. 7. The streamlines are completely different from the pathlines (7) of the same flow,
which are parabolae.

1.5 The Lagrangian derivative

(a.k.a. the convective derivative, or the material derivative).

We know how to measure the time derivative of a physical quantity associated with
the fluid, for example that of the temperature T (r, t), at a fixed point in space (the
∂T
Eulerian derivative); it is . This quantity will describe the change of temperature at
∂t
a fixed location, for example air temperature in Bristol. However, this quantity will not
be a measure of how a mass of air heats up or becomes colder. The reason is that air is
swept away by the prevailing flow field u(r, t). In other words, to describe the change of
temperature of a piece of air, we need to consider the rate of change of T (r, t), following
a fluid particle with trajectory r(a, t).
Definition 1.5.1 (Lagrangian time derivative)
DT
The Lagrangian derivative, denoted as , is the time derivative of a field T (r, t)
Dt
along the trajectory of a fluid particle.
Thus using the chain rule, we have
DT dT ∂T
≡ = + (ṙ · ∇)T.
Dt dt ∂t
But according to (1) and (2), ṙ = v(a, t) = u(r, t), and so
DT ∂T
= + (u · ∇)T. (11)
Dt ∂t
A particularly important example is the change in velocity (the acceleration) of a fluid
particle, which we need to apply Newton’s equations to fluid motion. Of course, the
velocity is a vector quantity, which means we have to apply (11) to each component:
 
Du Du1 Du2 Du3
= , , .
Dt Dt Dt Dt
The the acceleration of a fluid particle becomes
Du ∂u
= + (u · ∇)u. (12)
Dt ∂t
Note: The Lagrangian derivative (i.e. following fluid particles) is given in terms of
Eulerian (i.e. fixed point) measurements. It is vital to understand exactly how to com-
pute expressions like (u · ∇)u for a given velocity.
Example 1.5.2 (advective derivative) Consider an accelerating fluid flow, such as
the log flowing through a narrowing channel. Suppose

u = (U + kx, U − ky, 0).

Then
∂u
• The flow is steady since = 0.
∂t
• The advective term is
 
∂ ∂
(u · ∇)u = (U + kx) + (U − ky) (U + kx, U − ky, 0).
∂x ∂y

Hence the acceleration of the log is

Du
= (k(U + kx), −k(U − ky), 0) .
Dt

1.6 Mass conservation

One of the fundamental laws of continuum mechanics is the law of mass conservation.
That is, fluid is neither created or destroyed.

Consider an arbitrary finite volume, V , which is fixed in a fixed frame of reference. V

is bounded by the surface S and n represents a unit normal on S outward from V .

A fluid occupies the space of which V is a subset. The fluid has velocity u(r, t) and
density ρ(r, t). Fluid can flow in and out of V , and its density can change (within V ).

ds

V S

Figure 8: The mass inside a volume V can change only because because mass flows in or
out through the surface S.

The mass contained in V is Z

ρ(r, t)dV.
V
 u
n

dS
Figure 9: The flux through a surface element is determined by the inner product u · n.

The flux of mass out of V (n points outward) is calculated using the following ob-
servations, see Fig. 9. In time dt the fluid on the surface element dS is transported a
distance |u|dt in the direction of u. Projected onto the direction n perpendicular to S,
this corresponds to a distance u · ndt, and in turn to a volume (u · n)dtdS. Then the
total mass flowing in or out of dS is ρ(u · n)dtdS, and the rate of mass transport is this
quantity divided by dt. This means the total flux of mass out of V is
Z
ρu · ndS.
S

Definition 1.6.1 (mass flux density) The quantity j = ρu is called the mass flux den-
sity.

Now the rate of change of mass in V must equal the rate of change of mass in/out of
V through S, which means that
Z Z
∂
ρ(r, t)dV = − j · ndS (13)
∂t V S

Since V does not change with t and using the divergence theorem for the right hand side
of (13) we get Z Z
∂ρ
dV = − ∇ · jdV
V ∂t V
or Z  
∂ρ
+ ∇ · j dV = 0.
V ∂t
This is true for any fixed V , so we must have

∂ρ
+∇·j=0 (14)
∂t
at every point in the fluid. This is called the mass conservation equation or the continuity
equation.

Remark 1.6.2 (conservation laws) The structure of (14) is very general, and applies
to any conserved quantity; an equation of the form (14) is also called a conservation law.
Its basic ingredients are a conserved density (e.g. mass density or energy density), and a
flux (i.e. a mass flux or energy flux).
1.7 Incompressibility
Definition 1.7.1 (Incompressibility) A fluid is said to be incompressible if the density
Dρ
of each fluid ‘particle’ is constant, (i.e. = 0).
Dt
From (14) we have

∂ρ ∂ρ Dρ
0= + ∇ · (ρu) = + u · ∇ρ + ρ∇ · u = + ρ∇ · u
∂t ∂t Dt
So (ρ > 0) an incompressible fluid satisfies

∇ · u = 0. (15)

No fluid is completely incompressible, but even gases are often sufficiently incompress-
ible for (15) to apply to their motion. Incompressibility is a valid approximation if

• The flow speed is much less than the velocity of sound;

• Timescales are much larger than (sound frequency)−1 .

Of course this leaves out everything to do with sound waves, which are due entirely to
compressible effects.

Remark 1.7.2 (Constant density) A particular case of an incompressible fluid is one

Dρ
in which ρ = const everywhere in space and time, and hence = 0. However, this need
Dt
not be the case. For example, air in the atmosphere has different temperatures, resulting
in different densities; as a result, hot air will rise. But in a frame moving with the air,
the density may still be the same, and incompressibility satisfied.

1.8 Streamfunction for two-dimensional, incompressible flows

In most circumstances, the incompressibility condition (15) is awkward to satisfy. Some
progress can be made by using the following

Theorem 1.8.1 (Vector potential) If ∇·u = 0, then there exists a vector field A(r, t)
s.t.
u = ∇ × A. (16)

The vector field A is called the vector potential. If u is represented like that, u is clearly
incompressible. However, A is far from unique: if the gradient of any scalar function is
added to it, the same u results. Thus to find a unique A, constraints must be applied
to it, and we are back to square one! However, the situation is different if the flow is
two-dimensional. In that case, A must point in the direction perpendicular to the plane,
and we can write A = ψ(x, y, t)ẑ, where ψ is a single scalar quantity. We will consider
the case of a Cartesian and of a polar coordinate system in the plane separately.
Cartesian coordinates : Here A = ψ(x, y, t)ẑ, so ψ is written as a function of the two
Cartesian coordinates x, y in the plane. A single scalar function ψ corresponds exactly to
the expected number of degrees of freedom of the flow: there are two components of the
velocity, and one constraint ∇ · u = 0. Then
 
∂ψ ∂ψ
u=∇×A= ,− ,0 ,
∂y ∂x

which determines how the velocity field is determined in terms of the stream function. In
other words, we have represented the components u, v of a two-dimensional velocity field
as
∂ψ ∂ψ
u= , v=− . (17)
∂y ∂x
By construction, this representation guarantees that ∇ · u = 0, as can be checked directly.
In addition, ψ has a very convenient and useful physical interpretation. Let (x(s), y(s)
be a streamline, parameterized by s. Then the change of ψ along the streamline is
∂ψ ∂ψ ∂x ∂ψ ∂y ∂ψ ∂ψ
= + = λ(s)u + λ(s)v = λ(s) (−vu + uv) = 0,
∂s ∂x ∂s ∂y ∂s ∂x ∂y

having used (8) in the first step, and (17) in the second. It follows that ψ(x, y, t) = const
along a streamline of the flow. This motivates the following

Definition 1.8.2 (Streamfunction) We call the function ψ(x, y, t) which is constant

along streamlines the streamfunction of the flow.

ψ= constant

Figure 10: The streamfunction is constant along streamlines of the flow.

Note: For steady flows, the streamlines do not cross each other and fluid does not
cross the streamlines.

Example 1.8.3 (Vortex) Consider the following flow:

a
u= (y, −x) ≡ (u, v). (18)
x2 + y2
 y

Figure 11: Streamlines of the vortex (18)

First we must confirm that this flow is incompressible, i.e. that

∂u ∂v
+ = 0.
∂x ∂y
Now
∂u 2axy ∂v 2ayx
=− 2 , = 2 ,
∂x (x + y 2 )2 ∂y (x + y 2 )2
which indeed add up to zero.
The streamfunction must satisfy
∂ψ ay ∂ψ ax
= 2 , = 2 .
∂y x + y2 ∂x x + y2

From the first equation, ψ = a ln(x2 + y 2 )/2 + f (x), and then from the second equation
f 0 (x) = 0. Thus the streamfunction is
a
ln x2 + y 2 ,


ψ= (19)
2
only defined up to a constant, of course. Clearly, the streamlines ψ = const are circles,
as expected.

Polar coordinates : Clearly, a streamfunction such as (19) is represented more easily

in polar coordinates. In that case, we have A = ψ(r, θ, t)ẑ. Using the formula for the curl
in cylindrical polars, we find
1 ∂ψ ∂ψ
u=∇×A= r̂ − θ̂,
r ∂θ ∂r
so that the two velocity components are represented as
1 ∂ψ ∂ψ
ur = , uθ = − . (20)
r ∂θ ∂r
From the geometrical interpretation of ψ it is clear that streamlines are given by ψ(r, θ) =
const (of course this can also be checked by direct calculation). Thus

ψ(r, θ, t) = ψC (r cos θ, r sin θ, t),

where we have written ψC for the streamfunction in Cartesians.

Example 1.8.4 (Vortex in polar coordinates) From the above it follows that the stream-
function (19) of Example 1.8.3 can be written in polar coordinates as
a
ln x2 + y 2 = a ln r.


ψ=
2

Figure 12: Streamlines of a source.

Example 1.8.5 (Source in polar coordinates) In the case of a source, the flow is
purely radial, see Fig. 12, so that uθ = 0 and ur is independent of θ.

From (20) we have

1 ∂ψ 

ur = f (r) = 
r ∂θ
∂ψ
uθ = 0 = − . 

∂r
Thus we must have ψ = ψ(θ) but ∂ψ/∂θ independent of θ. So ψ = Aθ, where A is a
constant. It follows that f (r) = A/r.
Definition 1.8.6 (Source strength) The source strength is the flux of fluid from the
source point.

The flow is incompressible, so the flux at the origin equals the flux through any closed
boundary surrounding the origin. Found, most conveniently, by measuring the flux
through a circle radius r centred at the origin. The source strength is
I Z 2π Z 2π
A
m= u · nds = u · nrdθ = rdθ = 2πA,
circle 0 0 r
where n = r̂ in polars.
This means the mass flow through a circle surrounding the origin is independent of
the radius of the circle, and equals the source of mass flux at the origin. In fact, one can
show that the mass flux through any closed surface is m, using Green’s theorem (how
exactly?). In conclusion, the two-dimensional velocity field of a source is
m
u= r̂, (21)
2πr
and the stream function is
mθ
ψ= . (22)
2π
The velocity fields and streamfunctions derived here will be summarized in Appendix D.
Note: In Cartesians,
m
ψ= arctan(y/x).
2π

Another take on the above statement is provided by the following interesting property
of the streamfunction:

Theorem 1.8.7 (Flux across a line) Let r0 and r1 be two points in the plane, and ψ
the streamfunction. Then the volume flux crossing any curve from r0 to r1 is given by
ψ(r1 ) − ψ(r0 ).

Proof: The volume flux Q across the curve C joining the points r0 and r1 is
Z r1
Q= u · n ds,
r0

where ds is the arc length along the curve. We choose n in the mathematically posi-
tive, counterclockwise direction. Since the tangent vector in the direction of the curve is
t = ( dx, dy)/ ds, this is the orthogonal unit vector n = ( dy, − dx)/ ds, pointing in the
counterclockwise direction, going from r0 to r1 . Now
Z r1   Z r1 Z r1
∂ψ ∂ψ ( dy, − dx) ∂ψ ∂ψ
Q= ,− . ds = dy + dx = dψ = ψ(r1 ) − ψ(r0 ),
r0 ∂y ∂x ds r0 ∂y ∂x r0

as claimed. In particular, if r0 and r1 lie on the same streamline, the flux across the
streamline vanishes, which is obvious. Applied to the source (21), choose r0 = (r, 0) and
r1 = (r, 2π) in polar coordinates, so that ψ(r1 ) − ψ(r0 ) = m, according to (22), which is
indeed the source strength.
2 Flow dynamics for an incompressible inviscid flow
2.1 Forces on a fluid
Fluids move in response to the forces that act on each fluid particle. In order to apply
Newton’s second law in straightforward manner, we need to define a fluid “particle” (or
fluid element) in such a way that the mass is constant. We also want it to be big enough
(encompassing a large number of molecues), such that continuum theory applies. This
can be acccomplished for example by marking each molecule in an initial lump of fluid
(for example a cube). At later times, we think of the fluid particle as the volume δV
containing the original molecules. As the molecules move about, the fluid volume will
deform, but its mass remains constant.
There are two types of forces on a fluid element:

(i) Body forces are forces on a each parcel throughout the bulk of the fluid. In a gravi-
tational potential, the force is proportional to the volume δV of the fluid element.
In the typical case of a uniform gravitational field (as near the surface of the earth)
the force is δFb = ρδV g. We will denote the force density (per unit volume) by f (r).

(ii) Surface forces are transmitted across a surface element dS of a fluid parcel, so the
force is exerted by the fluid on the exterior of δV or vice versa.

When the fluid is a rest, the surface force must be in the direction of the normal, n,
since by definition a fluid cannot sustain any shear (it would flow). When a fluid is in
motion, tangential components of the force on dS can occur. These are associated with
viscosity which describes the effect that one layer of fluid in motion has on an adjacent
layer.

Definition 2.1.1 (Inviscid fluid) A fluid is said to be inviscid (also known as an ideal
fluid), when surface forces act in the normal direction only.

For an inviscid fluid, the ‘surface stress’ (i.e. Force per unit area) is in the direction
n normal to an infinitesimal surface element dS even when the fluid is in motion. In
other words, the force exerted by the exterior fluid on any surface element dS of the fluid
element δV is dFs = −pndS, where p(r, t) is the pressure. The pressure has units of force
per unit area. It is directed inwards because fluids are usually in a state of compression.

2.2 Equation of motion

We apply Newton’s law to a fluid element with volume δV and position r(a, t) as defined
above. The mass δm of this fluid element is constant by construction. Then Newton’s
equation of motion is
Du
δm = δF, (23)
Dt
where δF is the total force acting on the fluid particle. According to the above, δF is
composed of surface and body forces:
Z
δF = δFb + δFs = f (r, t)δV − pndS,
δS
where the integral is over the surface of the fluid element. Using Gauss’ theorem, the
surface integral can be written as
Z Z
pndS = ∇pdV ≈ ∇pδV,
δS δV

if the fluid particle is sufficently small so ∇p can be treated as constant over its volume.
Thus δF = (f (r, t) − ∇p)δV , and since we also have δm = ρ(r, t)δV , (23) becomes
Du
ρ = −∇p + f . (24)
Dt
Using the definition of the convective derivative (12), this can be written explicitly as
 
∂u
ρ + u · ∇u = −∇p + f . (25)
∂t
Either (25) or (24) are known as Euler’s equation (due to Euler 1756), on which the rest
of this course will be based. Notice that in deriving (25), we have made no assumption
about the fluid being incompressible.
What at first seems strange about Euler’s equation is that the pressure p on the right
is not specified. Indeed, a casual count shows that (25) is a set of three equations for 4
unknowns (three components of u and the pressure. We need another equation! In this
course, this will be the incompressibility condition (15): ∇ · u = 0. Together with (25),
we then have 4 equations for 4 unknowns.
A good analogy for this situation comes from point mechanics: Euler’s equation (25)
plays the role of Newton’s equation ma = F. However, only the body forces f (such as
gravity) are specified explicitely, p is unknown. The reason is that p plays the role of
a constraint force, with the constraint being ∇ · u = 0. We have to find p, such that
the constraint is satisfied at all times. Just as in mechanics, solving the constraint is the
most difficult bit. While Lagrangian mechanics teaches us how to solve the problem by
choosing an appropriate coordinate system, we don’t have this luxury in fluid mechanics:
it is not possible in general to find a representation of the velocity field u with two
independent components, which satisfies the constraint identically. The only exception,
as we have seen, is the case of two-dimensional flows, for which we have the stream
function representation.

2.3 The momentum flux

In Section 1.6 we derived the law for the conservation of mass. In fact, Euler’s equation
(24) is the statement of another conservation law, that of momentum conservation. To
obtain the momentum, we multiply the velocity u by the density ρ, to find the momentum
density per unit volume ρu. To be able to write conservation of momentum in the form
of the continuity equation (14), we write the equation for each Cartesian component
separately:
∂ ∂Mij
(ρui ) + = fi . (26)
∂t ∂xj
On the right-hand-side is the i’th component force density f , since according to Newton’s
law the time derivative of the momentum is the force. Indeed, if we define the flux of
x-momentum px by
(p )
fj x = M1,j ,
and there are no external forces, then it follows from (26) that
∂px
+ ∇ · f (px ) = 0. (27)
∂t
This describes the conservation of x-momentum, and has exactly the same form as (14).
We now want to show that (26) is equivalent to (25), with the momentum flux tensor
Mij = ρui uj + pδij . (28)
Indeed, inserting (28) into (26), we have
∂ρ ∂ui ∂ρ ∂ui ∂uj ∂p
ui + ρ + ui uj + ρ uj + ρui + δij =
∂t ∂t ∂xj ∂xj ∂xj ∂xj
 
∂ρ ∂uj ∂ui ∂ρ ∂ui ∂uj ∂p
− uj + ρ ui + ρ + ui uj + ρ uj + ρui + =
∂xj ∂xj ∂t ∂xj ∂xj ∂xj ∂xi
 
∂ui ∂ui ∂p
ρ + uj + = fi ,
∂t ∂xj ∂xi
which is precisely the Euler equation (25); in the second step we have used the continuity
equation (14).
For steady flows it is often useful to consider the momentum balance over a control
volume V . Since the flow is steady, the left hand side of (26) is zero. We also assume that
f = −ρ∇Φ, and that ρ is constant. Then
Z Z Z Z
∂Mij ∂Φ
0= dV + ρ dV = (Mij nj + ρΦni ) dS = (ui uj nj + (p + ρΦ)ni ) dS
V ∂xj V ∂xi S S

by the divergence theorem. So, back in vector form,

Z
ρu(u · n) + (p + ρΦ)ndS = 0 (29)
S

for any closed surface S inside or bounding the fluid. This is the momentum integral theorem,
of which we will make lots of use!
We now give a few simple examples of flows which are solutions of (25) in the absence of
an external force f = 0. In each case, it is crucial to also confirm that the incompressibility
condition ∇ · u = 0 is satisfied. Expressions for (u · ∇)u in various coordinate systems
are provided in chapter 0. We assume there is no external force.
Example 2.3.1 (Uniform flow) We have u = U. This must be a solution, as each
particle travels along straight lines, and thus no forces are involved.
Obviously, ∇ · u = 0, so the flow is incompressible. Turning to (25), the flow is steady,
∂u
and (u · ∇)u = 0, since u is constant. Thus ∇p = 0, and it follows that p = const,
∂t
as expected.

Example 2.3.2 (Source in three dimensions) The flow field must point radially out-
ward from the origin, and thus u = f (r)r̂. Then in spherical coordinates,
1 ∂(r2 f (r))
∇·u= .
r2 ∂r
For this to vanish, we must clearly have f (r) = A/r2 . The total flux from the origin must
equal the flux through a sphere of radius R. Since the outward normal is n = r̂, we have
Z
m= u · ndS = 4πR2 f (R) = 4πA.
SR

Thus the incompressible velocity field corresponding to a source of strength m is

m
u= r̂. (30)
4πr2
Now we test for the solution of Euler’s equation. Again the flow is steady, and in a
spherically symmetric situation Euler’s equation becomes

∂ur m2 ρ
ρ(u · ∇)u = ρur r̂ = − 2 5 r̂ = −∇p.
∂r 8π r
Thus the pressure only depends on the radial coordinate, with

∂p m2 ρ
= 2 5,
∂r 8π r
and so
m2 ρ
p = p0 − .
32π 2 r4
Thus with increasing r, the pressure increases. This is expected, as the flow slows down
in the radial direction, which means a force must be holding it back.

2.4 Hydrostatics
For a fluid at rest, u = 0 and so it follows from (24) that

∇p = f .

In the case of a graviational potential (which is conservative), the force is proportional to

the mass, mulitplied by the gradient of the potential ∇Φ. But since f is the force density
per unit volume, we obtain f = −ρ∇Φ. Then, if ρ = const, we obtain ∇(p + ρΦ) = 0.
Integrating, this leads to
p + ρΦ = 0. (31)
In the case of a constant gravitational force, f = −ρgẑ, where ẑ is pointing vertically
upwards. The corresponding potential is Φ = gz, so that

p + ρgz = const, in the fluid. (32)

In the case that ρ = ρ(z) (as is relevant for oceans or atmospheres),

Z z
∇p = −ρ(z)gẑ, ⇒ p(z) = p(z0 ) − g ρ(z 0 )dz 0 .
z0
 z

S n Stationary
fluid, density ρ
V

Figure 13: The force on a submerged body equals the weight of the displaced fluid.

Example 2.4.1 (Archimedes’ principle (∼ 250BC)) Let us apply these ideas to a

submerged body.
Z
Force on submerged body = − pndS
ZS
=− ∇pdV, Divergence Theorem
Z V

= ρgẑdV, ẑ unit normal in z-direction

V
= ρgV ẑ

Therefore Archimedes’ principle becomes: force on body = weight of fluid displaced is ρgV .

A body will float if its total weight is smaller than that of the fluid it displaces. According
to an anecdote, Archimedes used this principle to determine if a crown was made of pure
gold.

Figure 14: Left: the crown and a gold nugget used for reference have the same weight.
Right: The crown of lower density receives a greater upward force, so it becomes lighter
when submerged.
 patm
patm

H1 p(z)
1 z p2(z) H2

Figure 15: Forces act on both sides of a lock gate; everything is per unit length (into the
page).

Example 2.4.2 (Lock gates) Consider water on both sides of a lock gate.

To the left of the gate, we can use (32), and the atmosphere surrounding the fluid
surface is at constant pressure patm . Since the pressure must be continuous across the
interface (we will derive this condition more formally in the next section), we also have
p = patm on z = H1 , just inside the fluid. It follows that

p1 (z) = patm + ρg(H1 − z), 0 < z < H1 .

Similarly, to the right of the gate,


patm , z > H2 ,
p2 (z) =
patm + ρg(H2 − z), 0 < z < H2 ,

Thus the force exerted by the fluid on the gate (per unit width of the gate) is
Z H1
F1 = p1 (z)ndz, (33)
0

where n = x̂ is the vector pointing out of the fluid toward the gate. Thus the x-component
of the force from the left is
Z H1
F1 = F1 · x̂ = p1 dz = patm H1 + 21 ρgH12 ,
0

and the x-component of the force (per unit width) exerted by the fluid on the gate from
the right is Z H1
F2 = − p2 dz = −patm H1 − 12 ρgH22 ,
0
where now n = −x̂, and hence the - sign. Thus the net force in the x-direction is
1
F = F1 + F2 = ρg(H12 − H22 ),
2
which is toward the right (H1 > H2 ). In one of the problems, we will calculate the torque
on the gate.
2.5 Kinematic boundary condition
As with any PDE, it is important to understand the boundary conditions Euler’s equation
has to satisfy. The kinematic boundary condition expresses the conditions on the motion
of the fluid at the boundary. We will deal with solid boundaries first, and then generalize
to moving boundaries, for example those between two fluids or between fluid and air.

2.5.1 Solid boundary

Figure 16: A solid boundary with normal pointing into the fluid.

Let us consider the flow near a solid surface, whith normal n pointing into the fluid.
We begin by assuming that the solid be at rest. At any point on the surface, it is clear
that the fluid cannot penetrate the solid. This means that the normal component of the
velocity must be zero on the surface: u(r, t) · n = 0 for any r ∈ S. Instead, the flow can
only have a component parallel to the solid.
Now let us take the slightly more general situation that the solid is moving. At any
point on the surface, we can go into a frame of reference moving with the solid, such that
the point is at rest in this coordinate system. This leads us to the requirement

u(r, t) − u(solid) (r, t) · n = 0 for r ∈ S,

(34)

which is the boundary condition to be satisfied at solid surfaces.

Example 2.5.1 (A solid wall at rest) Let the wall occupy the plane z = 0, and the
fluid be at z > 0. Then n = ẑ, and the boundary condition is u · ẑ = 0, or uz = 0 for
z = 0.

Note that the other components of the velocity field, ux and uy , are arbitrary! The fluid
can “slip” over the solid without resistance, which is consistent with out notion of a
frictionless or ideal fluid.
 Figure 17: A solid boundary with normal pointing into the fluid.

Example 2.5.2 (A sphere of radius R moving through a fluid at velocity u) According

to (34), (u − u(solid) ) · n = 0 on the surface of the sphere.
At a given instant, let the centre of the sphere be at the centre of a spherical polar
coordinate system, and choose the z-axis such that u(solid) = U ẑ. Then n = r̂, so that
u · n = ur , and u(solid) · n = U (ẑ · r̂) = U cos θ. Thus the boundary condition reads
ur = U cos θ, for r = R.

2.5.2 Moving boundary

Now we formulate a generalization of the boundary condition of the preceeding section,
in which the boundary is allowed to deform in any way, for example the surface of the
ocean. To this end we describe the (time-dependent) boundary by the level line of some
function: let S(r, t) = 0 describe the equation of a surface between two fluids. As the
flow evolves, particles remain on the surface S if
DS ∂S
= + u · ∇S = 0 (35)
Dt ∂t
Namely, we want that a Lagrangian particle r(a, t), once on the surface, remains on
surface. This means that
d ∂S ∂S
0 = S(r(a, t), t) = + ṙ · ∇S = + u · ∇S,
dt ∂t ∂t
as claimed.

Oil
S=0
Water
Figure 18: The interface between a phase of oil and a phase of water is given by the level
curve S(r, t) = 0.
 Let us apply this to the interface between two liquids, as shown in Fig. 18. Then
DS ∂S

= + u(oil) · ∇S = 0, on S = 0 from above  
Dt ∂t
DS ∂S
= + u(water) · ∇S = 0, on S = 0 from below 

Dt ∂t
Now ∇S is normal to the surface S = const, so n = ∇S/|∇S| is the unit normal to the
surface S = 0. It follows that

u(oil) · n = u(water) · n. on S = 0 (36)

This expresses the intuitive notion that the normal component of the velocity on either
side of the interface must be equal.
Note that (36) expresses exactly the same thing as (34) for a solid boundary, except
that the velocity of the solid is now prescribed (and can only be a translation or rotation).

Example 2.5.3 (A collapsing/expanding bubble) An example we will consider in

section 3.4 below, is the collapse of a spherical bubble inside a fluid. Let us assume that
the radius R(t) of the bubble changes according to some time program, which we will
calculate later on. If we choose

S(r, t) = r − R(t), (37)

S(r, t) = 0 indeed describes the surface of the bubble.

Now computing (35), we find

∂S ∂S
0= + u · ∇S = −Ṙ + u · r̂ = −Ṙ + ur ,
∂t ∂r
where ur = u ·r̂ is the radial compoment of the velocity in spherical coordinates. In other
words, the boundary on the surface of the bubble is

ur = Ṙ, (38)

which is very intuitive, and could have been written down without any calculation!

2.6 Dynamic boundary condition

To derive a boundary condition the pressure has to satisfy on the boundary of a free surface
between two fluids (e.g. interface between oil/water) we use momentum conservation (26).

n (s) oil
Vε ε C

water

Figure 19: A “pillbox” argument at the interface between two phases.

 Let us integrate (26) throughout a small disc V of thickness  and let  → 0 (cf.
Fig. 19). We only consider the direction n(s) normal to the interface, in which u · n(s) is
continuous. Then Z
∂
(ρu · n(s) )dV → 0, as  → 0,
V ∂t
and therefore
Z Z
(s) ∂
h i h i
(s) (s) (s) (s) (s)
0 = lim ni Mij dV = lim ni Mij nj dS = ni Mij nj − ni Mij nj .
→0 V

∂xj →0 S

+ −

But this means that

ρ(u · n(s) )2 + pu · n(s) + − ρ(u · n(s) )2 + pu · n(s) − = 0.

   
(39)

Since u · n(s) is continuous across the interface, we are left with

p+ = p− , (40)

which says that the pressure is continuous across an interface. This is very intuitive, and
in fact we have already used (40) in example 2.4.2. In the atmosphere, the pressure has
the constant value patm , so (40) says that the fluid pressure at a liquid-air interface must
be p = patm .

Example 2.6.1 (Axisymmetric jet impinging on a wall) As a first example of the

application of the momentum integral theorem (29), let us consider the situation shown in
Fig. 20. What is the force exerted by the flow on the wall ? Ignore the effects of gravity.

Figure 20: A jet impinging perpendicularly on a plate. The dashed line marks the control
surface S.

Take a control surface S = Ss + Sj + Sf + Sw of cylindrical shape that is composed of the

free surface of the jet on one side, and the wall on the other. The surface is closed by a
circular surface through the jet far from impact (marked “jet”), and a cylindrical surface
where the fluid flows parallel to the wall (marked “film”). Assume uniform flow parallel
to wall in the film region far from impact and across the cross section of the jet u = −U ẑ.
In the preceeding section we have shown that the pressure on all boundaries exposed to
the atmosphere is patm . This includes the jet itself, which has assumed the pressure patm .
Then in the absence of body forces,
Z
(ρu(u · n) + pn) dS = 0. (41)
S=Sj +Ss +Sf +Sw

Now Z Z Z Z
pndS = (p − patm )ndS + patm ndS = − (p − patm )ẑdS
S S S Sw
Z Z
since patm ndS = ∇patm dV = 0 and p = patm everywhere else. Also,
S V
Z Z Z
2
ρu(u · n)dS = ρU ẑdS + ρUr2 r̂dS
S Sj Sf

since u · n = 0 on the wall and on the surface of the liquid. Thus (41) becomes
Z Z Z
2
(p − patm )ẑdS = ρU ẑdS + ρUr2 r̂dS.
Sw Sj Sf

Multiplying this by ẑ to only consider the z-component of the momentum balance one
obtains Z
(p − patm )dS = ρU 2 A,
Sw

where A is the cross section of the jet. But then the force on the wall in the ẑ-direction is
Z
Fz = − (p − patm )dS = −ρU 2 A.
Sw

2.7 Bernoulli’s equation for steady flows

The Euler equation for steady motion is

ρ(u · ∇)u = −∇p − ρ∇Φ,

assuming that the body forces possess a potential Φ. We can transform the term on the
right using the vector identity (see Appendix A)

u · ∇u = ∇( 12 u · u) − u × ω (42)

where ω = ∇ × u is the vorticity. Thus we obtain (taking ρ = const)

−u × ω = −∇(p/ρ + Φ + 21 u2 ). (43)

To make the left hand side disappear we take the dot-product with u on both sides to
find
u · (u × ω) = 0 = u · ∇(p/ρ + Φ + 21 u2 ).
 dr
The definition of a streamline is = λ(s)u, so
ds
d dxi ∂
(p/ρ + Φ + 12 u2 ) = (p/ρ + Φ + 12 u2 ) = λ(s)u · ∇(p/ρ + Φ + 12 u2 ) = 0.
ds ds ∂xi
Consequently, the quantity

H = p + ρΦ + ρu2 /2 = const, along any streamline r(s) in the flow. (44)

In the problems we show that the energy flux is jE = uH. Thus (44) can be seen
as a statement of energy conservation, where the Bernoulli quantity H is conserved for
a fluid element along a streamline. Note: A very nice discussion of the Bernoulli effect
and of pressure forces is presented in the film: “Pressure fields and acceleration”, found
on Youtube.

Example 2.7.1 (A three-dimensional source) According to Example 2.3.2, the ve-

locity file is
m
u= r̂.
4πr2
Evidently, stream lines are rays pointing radially outward from the center, and

1 2 m2
u = .
2 32π 2 r4
According to (44),
m2 ρ
p=− + const,
32π 2 r4
where the constant could in principle be different for each streamline. However, the
pressure must go to the same constant far from the source, so we reach the same conclusion
as in Example 2.3.2. However, Bernoulli’s theorem is particularly powerful for more
complicated flows, for which simple solutions are not available, as the next example shows.

Example 2.7.2 (Flow out of a tank) A tank of uniform cross section A0 has a small
hole, with area Ae , and at a height h above the base, see Fig. 21. The height of the fluid
above the hole is H. What is the flow speed out of the drain ?
 patm area
Ao

H Uo area
Ae
Ue

z streamlines
h

Figure 21: Water flowing out of a tank.

First, conservation of mass gives

ρA0 U0 = ρAe Ue , (45)
since the mass fluxes across two boundaries equal. Since the top fluid level decreases very
slowly, the flow is approximately steady. Second, we want to apply Bernoulli’s theorem
to a streamline that exits the hole. It is important to notice that any such streamline
must originate from the surface of the container (dashed lines), and ends at the hole.
Namely, all streamlines must be parallel to a stationary wall, and thus cannot originate
from the wall. By contrast, the top surface is moving downward (albeit slowly), and
thus the velocity field has a component pointing downward, in the local direction of the
streamline. The potential in (44) is Φ = gz. The pressure at the top surface and just
outside the hole is p = patm . Then
• At surface: patm + ρg(h + H) + 12 ρU02 = C
• At exit: patm + ρgh + 12 ρUe2 = C,n
and along the same streamline the constant C is the same. Eliminating C and patm and
using (45) we find
Ue2 [1 − (Ae /A0 )2 ] = 2gH.
It is reasonable to assume that Ae  A0 , so
p
Ue ≈ 2gH. (46)
This result corresponds exactly to the free fall velocity of a mass having been dropped
from a height H. The reason for this is clear: a fluid particle, having left the top surface,
moves along a (perhaps complicated) path dictated by the geometry of the vessel and
the incompressibility constraint ∇ · u = 0. However, since this motion is frictionless, the
shape of this path is irrelevant, and the final speed is dictated by energy conservation
alone.
With the main result (46) in hand, we can address two practical questions:
 1. What is the draining time td for an initial filling height H0 to reach 0 ?
The equation of motion for the surface is
dH Ae p
= −U0 ≈ − 2gH.
dt A0
and so Z 0
A0 0 dH
Z
dt A0 p
td = dH = − √ = 2H0 /g. (47)
H0 dH Ae H0 2gH Ae
2. What is the horizontal distance x travelled by the jet?
The fluid particles making up the jet are subject to gravitational acceleration only,
since the pressure along the jet has the constant value p = patm . Hence√ they obey
the same equations as a projectile fired horizontally
p with speed Ue = 2gH. The
time taken by a particle to hit the floor is t = 2h/g, and so
p p √
x = Ue t = 2gH 2h/g = 2 Hh. (48)

Finally we can fix a total filling height (measured from the floor) to be Hf , so that
Hf = H + h. Where do we have to place h to maximise x? To answer this question,
we maximise q
x = 2 (Hf − h)h
with respect to h. This gives h = 21 Hf , and so
xmax = Hf .

Figure 22: Experimental test of (48).

Experiment: We took data from the bottom hole, which is at h = 2.4cm. The length
of the trajectory x was measured as for the filling height Hf = H +h as given in the table:

p
Hf [cm] H [cm] 2 (Hf − h)h [cm] xexp [cm]
3.5 ± 0.1 1.1± 0.1 3.2± 0.2 2.2± 0.1
8.5 ± 0.1 6.1 ± 0.1 7.7± 0.2 6.6± 0.1
12. ± 0.1 9.6± 0.1 9.6± 0.2 8.1± 0.1
There are signifcant deviations, but they decrease for increasing Hf .

Example 2.7.3 (Flow through a slowly diverging channel) We can use Bernoulli’s
equation to calculate the pressure distribution in a pipe whose cross section is changing.

area
A2
area
A1
U1 streamline U2

Figure 23: A slowly widening channel.

From mass conservation: A1 U1 = A2 U2 , and from Bernoulli’s equation (ignoring grav-

ity):
1 2 1
ρU1 + p1 = const = ρU22 + p2 .
2 2
Hence, the pressure drop is
1 1
∆p = p2 − p1 = ρ(U12 − U22 ) = ρU12 (1 − A21 /A22 ) > 0 (49)
2 2
for a widening channel. As a channel widens, the flow slows down, and the pressure
increases. In the Venturi tube (see Fig. 24), this effect is used to measure the flow rate
Q = A1 U1 . Namely, since A1 ,A2 , and ∆p can be measured, we can express
s
2∆pA21 A22
Q=
ρ(A22 − A21 )
in terms of measurable quantities.

Figure 24: The Venturi tube measures the pressure difference between parts of a pipe
with different cross sections.
2.8 The vorticity equation
The vorticity is the key to many phenomena in fluid mechanics. It is often concentrated
into small regions of space, for example near a point (a point vortex), a line or a sheet (a
vortex sheet). Understnding the motion of the vorticity thus often holds the key to the
entire evolution of the flow. We start with Euler’s equation having rewritten the transport
term using (42), namely

∂u
− u × ω = −∇(p/ρ + Φ + 21 u2 ).
∂t
Here we have assumed that the density ρ is constant. Taking the curl of this equation
eliminates the right hand side:
∂ω
− ∇ × (u × ω) = 0.
∂t
We can use a vector identity:

∇ × (u × ω) = u(∇ · ω) − ω(∇ · u) + (ω · ∇)u − (u · ∇)ω = (ω · ∇)u − (u · ∇)ω

since the flow is incompressible and ∇ · ω = ∇ · (∇ × u) = 0. Hence

∂ω
+ (u · ∇)ω = (ω · ∇)u, (50)
∂t
or
Dω
= (ω · ∇)u (51)
Dt
in either form, (50) or (51) are called the vorticity equation. Its advantage is that we have
succeeded to eliminate the pressure from the Euler equation. However, we have payed a
heavy price: (50) contains both the velocity and the vorticity fields u, ω. Thus whenever
ω has been updated using (50), we need to reconstruct u from it.
The situation is more promising in two-dimensional flows, u = (u(x, y), v(x, y), 0) and
so, by definition ω = (0, 0, ω(x, y)) = ω(x, y)ẑ. It is then clear that the term

∂u
(ω · ∇)u = ω(x, y) =0
∂z
and so
Dω
=0 (52)
Dt
That is, vorticity is conserved as it moves with the flow.

Note 2.8.1 (Conservation of vorticity in two dimensions) If ω = 0 at time t = 0,

then ω = 0 for all time. Vorticity cannot be generated in a 2D flow.

2.9 Kelvin’s circulation theorem

In three dimensions conservation of vorticity (which corresponds to conservation of angular
momentum in mechanics) takes a somewhat more subtle form.
Definition 2.9.1 (Circulation) The circulation of a velocity field is defined to be
I
Γ= u · dl, (53)
C(t)

where C(t) is a closed loop which moves with the fluid.

Note 2.9.2 (relation to a patch of vorticity) By Stokes’ theorem,
I Z Z
Γ= u · dl = (∇ × u) · ndS = ω · ndS,
C(t) S(t) S(t)

where S(t) is a surface whose edges connect with C(t).

Theorem 2.9.3 (Kelvin’s theorem) The derivative of the circulation along a closed
loop convected by the flow is constant:
DΓ
= 0. (54)
Dt
Proof: Since C(t) moves with the fluid, we can think of it as being made up of Lagrangian
markers s, which move with the fluid. Thus C(t) is defined by the curve r(s, t), with (say)
a ≤ s ≤ b, and we can parameterize Γ as
Z b
∂r
Γ= u · ds.
a ∂s
Dr(s, t)
Since s is a Lagrangian marker, = u. Thus the total derivative of Γ is
Dt
Z b Z b Z b Z b
DΓ Du ∂r ∂ Dr(s, t) ∂r ∂u
= · ds+ u· ds = − ∇ (p/ρ + Φ)· ds+ u· ds, (55)
Dt a Dt ∂s a ∂s Dt a ∂s a ∂s
where we have used Euler’s equation (24) for constant ρ.
But both integrals on the right hand side of (55) vanish: The first,
Z b I
∂r
∇ (p/ρ + Φ) · ds = ∇ (p/ρ + Φ) · dl = 0
a ∂s C(t)

is the integral of a gradient field over a closed loop. The second,

Z b
1 b ∂ 2
Z
∂u
u· ds = u ds = 0,
a ∂s 2 a ∂s
vanishes again since s parameterizes a closed loop, and so u(a) = u(b). Hence we obtain
the result.

Figure 25: A circle being convected by the shear flow (56).

Example 2.9.4 (A shear flow) Let a closed loop of particles C(t) be defined by (see
Fig. 25):
r = a(cos s + αt sin s, sin s, 0), 0 ≤ s < 2π,
where each value of s corresponds to a different fluid particle, and a, α > 0.

On Sheet 1, Q. 1, we worked out how C(t) changed in time. We also showed that the
Eulerian velocity field is
u = (αy, 0, 0), (56)
which clearly is a solution of Euler’s equation.
Since s parameterizes the curve, we can calculate the circulation at any time as
Z Z
∂r
Γ= u · dl = u · ds.
C(t) 2π ∂s

Now
∂r
= a(− sin s + αt cos s, cos s, 0)
∂s
and
∂r
u(s) = = aα(sin s, 0, 0),
∂t
so that
Z Z
2 2
Γ = αa sin s(− sin s + αt cos s)ds = −αa sin2 sds = −παa2 .
2π 2π

This is constant, consistent with Kelvin’s theorem.

3 Irrotational flows: potential theory
Kelvin’s circulation theorem states that the circulation around any loop convected by
the flow cannot change, if the fluid is inviscid. In particular, if the circulation is zero,
it will remain zero. Now if there is no vorticity present in the flow at some initial time,
the circulation around any loop in the flow is zero. Hence Kelvin’s circulation theorem
guarantees that there will never be any circulation in the flow, and thus ω = ∇ × u = 0
throughout the domain. No vorticity can be created in a flow domain if the viscosity is
small.

Figure 26: Patches of vorticity remain confined by a loop around them.

However, we will see in a little while that problems arise with this argument when
applied near solid boundaries. The reason is that the circulation theorem applies to loops.
However, it is impossible to surround a point on a solid boundary by a loop that stays
in the flow. Thus no conclusions can be drawn on parts of the flow that originate from a
solid wall. However, even in this case Kelvin’s theorem has very important consequences.
Imagine a localised region of vorticity (shaded region) having been created by shedding
from a boundary, as seen in Fig. 26. Placing a loop Γ around it, which is convected by
the flow, both are convected together, and the vorticity will stay inside the loop. Then
Kelvin’s theorem guarantees that the total strength remaines constant. In other words,
the region of the flow which contains vorticity remains localised within a given region.

Example 3.0.5 (Circulation of a vortex) The most important example of an isolated

region of vorticity is that of a vortex, whose vorticity is confined to a point. We have
already looked at this in Example 1.8.3 and on Sheet 3, Q2.

In polar coordinates, the velocity field is u = A/rθ̂. Then ω = ∇ × u = 0 apart from at

r = 0. But I Z 2π
A
Γ= u · dl = rdθ = 2πA.
C 0 r
In terms of the total circulation, the velocity field becomes
Γ
u= θ̂
2πr
 We will concentrate on flows which either do not contain vorticity, or whose vorticity
is concentrated into singular regions. In that case, very powerful tools from potential
theory can be brought to bear on fluid problems, as we will see now.

3.1 The velocity potential

Theorem 3.1.1 Existence of the velocity potential Let ω = ∇ × u = 0 throughout a flow
domain D (apart from at isolated singularities). SUch a flow is called irrotational. Then
there exists a scalar field φ(r, t) such that

u = ∇φ. (57)

The field φ is called the velocity potential or simply the potential of u. If the velocity
field u is given, the potential is calculated from the path integral over u,
Z r
φ= u · dl,
r0

provided of course the flow is indeed potential. In particular, the result is independent of
the path taken, provided the domain is simply connected (see Section 4.2 below).
If the fluid is also incompressible then, ∇ · u = 0 and

(∇ · ∇)φ = ∇2 φ ≡ 4φ = 0. (58)

This is Laplace’s equation. When the fluid domain communicates with a (moving) solid
boundary, we impose the kinematic boundary condition (34), which can be written as
u · n = f (r), or n · ∇φ = f (r) for some given f .
The field equation (Laplace’s) and the boundary conditions are an example of a
Neumann boundary-value problem. An important feature of the equations are that they
are linear and so we can use superposition of solutions. This may seem odd, since the
Euler equation is a nonlinear equation, and does not obey the superposition principle.
However, as we will see very soon, the pressure depends on φ in a non-linear way (later),
which reflects the nonlinear character of the Euler equation.

Theorem 3.1.2 (Uniqueness of the velocity potential) Suppose an incompressible

irrotational fluid occupies a simply connected domain D, so u = ∇φ and that on the
boundaries S of D, u · n = f (r). Then the potential φ is unique up to a constant.

. Proof: Suppose φ were non-unique, so that both φ1 and φ2 satisfy

∇2 φi = 0, and n · ∇φi = f (r), for i = 1, 2

and φ1 6= φ2 . Let ψ = φ1 − φ2
Z Z Z
2
|∇ψ| dV ≡ ∇ψ · ∇ψdV = ∇ · (ψ∇ψ)dV, since ∇2 ψ = 0
D D ZD

= ψ∇ψ · ndS = 0, using Gauss’ theorem

S
since ∇ψ · n = ∇φ1 · n − ∇φ1 · n = 0 on S.
Since the integral over |∇ψ|2 vanishes, |∇ψ|2 must also vanish throughout the domain:
if |∇ψ|2 were non-zero somewhere, it would make a strictly positive contribution to the
integral. But sine there are only positive contributions to the integral, the total would
also be strictly positive, contradicting what we showed. Hence we conclude that |∇ψ| = 0
in D, and so ∇ψ = 0 in D. In a simply connected domain this implies that ψ = const in
D. But this means that φ1 and φ2 only differ by a constant, which is irrelevant, as only
terms in ∇φ appear. We will come back to the non-uniqueness that exists in multiply
connected domains below.

3.2 Some simple flows and their potentials

3.2.1 Uniform flow
If the velocity field is constant, u = U = const, a simple integration yields
φ = U · r.
If u is chosen to point in the x-direction, u = U x̂, then φ = U x. In two dimensions, the
streamfunction would then be ψ = U y.

3.2.2 A three-dimensional source

The flow field has been derived in Chapter 2, see (30). Integrating radially from infinity,
we find
m
φ=− .
4πr
It is shown readily that Laplace’s equation is indeed satisfied:
1 xi 3 x2
∂i2 = −∂i 3 = − 3 + 3 5i = 0.
r r r r
Taking the gradient, one finds the velocity field:
m mxi
ui = −∂i = ,
4πr 4πr3
and thus
mr m
u= 3
= r̂,
4πr 4πr2
where r̂ is the basis vector in spherical coordinates.

3.2.3 Stagnation point flow

y

stagnation point

x
 We begin with the two-dimensional case of a flow onto a flat plate, see Figure. We
consider the simplest case in which the flow is two-dimensional, u = (u, v, 0). There must
be a stagnation point u = 0 somewhere along the plate, which we place at the origin,
x = y = 0. Performing a Taylor expansion about the stagnation point, to lowest order
the x-component of the velocity will increase linearly as one moves away from the line of
symmetry: u = Ax. On account of incompressibility, we have
∂Ax ∂v
A= =− ,
∂x ∂y
and thus v = −Ay if one moves away from the plate against the flow. Integrating, one
finds the potential:
A
φ = (x2 − y 2 ) (59)
2
It is easy to confirm that this is a solution of Laplace’s equation, and thus satisfies the
fluid equations. More importantly, since (59) was derived from a local expansion, it is
valid locally near any two-dimensional stagnation point, regardless of the shape of the
surface, or of the exterior flow. The stream function corresponding to (59) is
ψ = Axy, (60)
as is confirmed directly from (17).

Figure 27: An axisymmetric flow around a so-called Rankine body, showing a stagnation
point near the tip of thre body to the left. The flow in the neighborhood of this point
will be described by (61).

The axisymmetric version of stagnation point flow is relevant if any axisymmetric body
is placed in an oncoming stream. The potential is in that case (exercise)
A 2
φ= (r − 2z 2 ). (61)
4

3.3 Bernoulli’s theorem for unsteady, irrotational flows

Within potential flow theory, we have the slightly confusing situation that we are seem-
ingly solving flow problems only by solving Laplace’s equation, which in turn comes from
the incompressibility condition ∇ · u = 0. We are not using Euler’s equation of motion
(25) to solve flow problems (at least this is almost true; we have been using (25) indi-
rectly inasmuch we used it to derive Kelvin’s theorem, which motivated the potential flow
assumption).
However, the Euler equation does come in, albeit in a sneaky way, as we will now
show. According to (43), the Euler equation can be written

∂u
− u × ω = −∇(p/ρ + Φ + 21 u2 ).
∂t
For irrotational flows, u = ∇φ and ω = 0, and so
 
∂φ p 1 2
∇ + + Φ + 2 u = 0.
∂t ρ

It follows that
∂φ p
+ + Φ + 12 u2 = C(t), (62)
∂t ρ
where C(t) is an arbitrary function of time.

Remark 3.3.1 (Unsteady Bernoulli equation)

1. Note that (62) applies throughout the fluid, not just on a streamline.

2. By defining a new velocity potential

Z t
φ̃ = φ − C(t0 )dt0 ,

the function C(t) can always be eliminated from the problem. The potential φ̃ leads
to the same flow field! (why?)

The Bernoulli equation (62) is best viewed as an equation for the pressure. Thus solving
a potential flow problem then becomes a two-step process:

1. Solve Laplance’s equation for φ, using the boundary conditions n · ∇φ = f (r) on

the boundary of D.

2. Given φ, compute u = ∇φ and then (62) to compute the pressure. The pressure
depends on u in a non-linear way, a reflection of the nonlinearity of the Euler
equation.

3.4 The collapse of a spherical cavity

As a first application of this program, consider a spherical cavity that has opened up in a
large expanse of liquid. Such holes or “cavitation bubbles” are produced quite commonly
near ships’ propellors, where the force of the propellor “rips apart” water molecules and
leaves a small hole behind. Thus there is little pressure inside the cavity, but a large
pressure in the surrounding fluid, which forces the hole to close again. Since the hole
is small, it is justified to neglect gravity, and to consider a constant positive pressure
p0 > 0 far from the cavity, while the pressure inside the cavity vanishes. For simplicity,
we suppose that at t = 0 we have a spherical bubble, initially at rest, with radius R0 in an
infinite fluid. The pressure difference causes the bubble to collapse, so R(t) < R0 , where
R(t) is radius at time t. We want to find R(t) and the time of collapse tc .
As the bubble collapses, to the outside fluid it acts like a sink at its centre. Thus we
try the corresponding potential
A(t)
φ(r, t) = ,
r
where A(t) is to be determined. We have seen that 4φ = 0. Now
∂φ A
ur = = − 2,
∂r r
and the boundary condition is
dR
ur = Ṙ ≡ , on r = R, ⇒ A = −R2 Ṙ,
dt
and hence
R2 Ṙ
φ(r, t) = −. (63)
r
To find R(t), we use the unsteady Bernoulli equation, so that
∂φ
ρ + p + 21 ρu2 = C(t) = p0 ,
∂t
since φ → 0 and u → 0 as r → ∞. Now, on surface of bubble p = 0 (by assumption) and

∂φ R2 R̈ + 2RṘ2
|u|2 = Ṙ2 , and =− .
∂t R
Substituting into Bernoulli gives

−RR̈ − 2Ṙ2 + Ṙ2 /2 = −RR̈ − 3Ṙ2 /2 = p0 /ρ,

which is an ODE for R(t). To integrate it, we multiply by 2R2 Ṙ, so that

d  3 2 p0 2p0 dR3
R Ṙ = 3R2 Ṙ3 + 2R̈ṘR3 = − 2R2 Ṙ = − .
dt ρ 3ρ dt
Integrating, we find
2 p0 3
R3 Ṙ2 = − R + C,
3ρ
where C = 2R03 p0 /(3ρ) is a constant, determined by the initial condition Ṙ = 0 when
R = R0 . Hence
2 p0 R03
 
2
Ṙ = −1 (64)
3 ρ R3
or r 1/2
2 p0 R03 − R3

dR
=− ,
dt 3ρ R3
where the negative sign must be chosen because Ṙ < 0 by assumption. Integrating up
gives
Z R(t) Z t 1/2 1/2
R3/2 dR

2 p0 2 p0
=− dt = −t .
R0 (R03 − R3 )1/2 0 3ρ 3ρ
If the bubble collapses at time t = tc , then R(tc ) = 0 implies (R = R0 u)
1 1/2
u3/2 du
Z 
tc 2 p0
= (65)
0 (1 − u3 )1/2 R0 3ρ

The integral on the left-had side of (65) can be evaluated to give

Z 1 √
u3/2 du
   
3 5 2
3 1/2
= √ Γ Γ ≈ 0.747.
0 (1 − u ) 2 π 6 3

Hence the collapse time is

 1/2
ρ
tc = 0.915 R0 . (66)
p0
We can also analyse the asymptotic behaviour of the collapse during its final moments,
when the radius goes to zero. Such singular events are usually described by power laws,
so we try
R = B(tc − t)α .
Now Ṙ2 ∝ (tc − t)2α−2 , and the right hand side of (64) scales as R−3 ∝ (tc − t)−3α . Thus
2α − 2 = −3α, and the power law exponent is found to be
2
α= .
5
Plugging this back into (64), we find
1/5
2p0 R03

R= (tc − t)2/5 ,
3ρ

in the limit of R approaching zero.

What is remarkable about the collapse is that the inertia of the surrounding liquid is
focused into a point by the converging motion. If one calculates the speed of the collapse,
one finds
Ṙ ∝ (tc − t)−3/5 ,
which is diverging as the radius goes to zero.
If the cavity is filled with gas, the gas can be compressed to very high pressures during
the last moments, reaching around 10,000 K in the interior of the bubble. This may
cause ionisation of the gas, which may start to glow. This phenomenon is known as
sonoluminescence, since the collapse of the bubble is started off by an external sound
field.
3.5 Flow past a sphere

Figure 28: Schematic of the flow around a sphere in a steady stream.

One of the most fundamental problems of fluid mechanics is to understand the flow around
obstacles, which stand in the way of a flow. Equivalently (think of Galilean invariance)
one can think of a body moving inside a fluid which is at rest (for example an aeroplane).
Let us begin studying this problem by considering the perhaps simplest possible body, a
sphere of radius R. The sphere is placed inside a uniform stream, which we choose in the
direction of the z or symmetry axis: U = U ẑ.
Far from the body, the flow will be uniform. The sphere introduces a perturbation to
the flow which decays as r → ∞, i.e. as one moves away from the body. According to
the superposition principle (the Laplace equation is a linear equation), this perturbation
is to be added to the uniform flow:
u = U + ∇φ, (67)
where φ is the potential of the perturbation. The simplest such flow we know is a source,
with potential φ ∝ 1/r. However, such a flow cannot describe the physical situation at
hand. Far from the sphere, the total mass flow through a closed control surface should be
zero, not finite as for a source. Another way of putting it is that the flow field of a source
decays too slowly (like 1/r2 ) away from the body.
This can be helped by taking the derivative of the source solution. Taking the deriva-
tive in the µ-direction, one arrives at the potential
1 µ·r
φ=µ·∇ =− 3 , (68)
r r
which is called the dipole flow (why?). The dipole is oriented in the µ-direction and has
strength |µ|. This flow is also a solution of Laplace’s equation, since
   
1 1
4 µ·∇ = µ · ∇4 = 0.
r r
In summary, from (67) and (68) our ansatz for the velocity is
µ j xj µi xi
ui = U δi3 + ∂i φ = U δi3 − ∂i 3
= U δi3 − 3 + 3µj xj 5 ,
r r r
or
1  r
u = U ẑ − µ − 3µ · r . (69)
r3 r2
 Now we have to satisfy the boundary condition on the surface of the sphere, which is
u · n = 0. The normal vector to the sphere is n = r/R. Now on the surface of the sphere
(r = R):
1  µ · n µ·n
u · n = U (ẑ · n) − 3 µ · n − 3 3 = U (ẑ · n) + 2 3 ,
R R R
3
UR
so choice µ = − ẑ will guarantee that this is zero for r = R. In summary, the velocity
2
field which satisfies all the boundary conditions is (z = ẑ · r):
 3
U R  z 
u = U ẑ + ẑ − 3 n . (70)
2 r r
Now we want to calculate the pressure field, and in particular the pressure on the
surface of the sphere. This will permit us to calculate the total force on the sphere. On
the surface we have
3U  z 
u= ẑ − n ,
2 R
so
9U 2 z2 9U 2 z2 9U 2
   
2 z
1 − cos2 θ ,


u = 1 − 2 ẑ · r̂ + 2 = 1− 2 =
4 R R 4 R 4
z
using that = cos θ.
r
According to Bernoulli’s equation
u2
p+ρ = C(t),
2
everywhere in the fluid. But the pressure is p = patm far from the sphere, and the flow is
uniform: u2 = U 2 , which implies that C = patm + ρU 2 /2. As a result, the pressure on the
surface of the sphere is
ρ(U 2 − u2 ) ρU 2
9 cos2 θ − 5 ,


p = patm + = patm + (71)
2 8

Figure 29: The pressure distribution around a sphere as function of the angle.

This result has the remarkable property that the pressure distribution is symmetric
around the equator of the sphere θ = π/2. Thus while the sphere is pushed by the flow in
the downstream direction, there is an equally high overpressure in the back. As a result,
without further calculation, we can conclude that the total drag force on the sphere is
zero. By definition, the drag on a body in a uniform stream is the force in the direction
of the stream. The force perpendicular to the stream is called the lift, which is of course
zero if the flow is axisymmetric, as is the case here.
Let us check these statements by doing the calcuation explicitly. Integrating over the
surface of the sphere, the total force is
Z Z 2π Z π
2
F=− pndS = −R pn sin θdθdφ,
r=R 0 0

where n = (sin θ cos φ, sin θ sin φ, cos θ) is the normal in spherical coordinates. Performing
the φ-integral over the x and y components, we obtain
Z 2π Z 2π
cos φdφ = sin φdφ = 0,
0 0

so the first two components of F vanish. This is clear because on account of the symmetry
of the problem, there can be no component normal to the direction of the flow, which is
in the z-direction. It remains to calculate the force in the flow direction (n · ẑ = cos θ):
Z 2π Z π 2 Z π
2 2 ρU
9 cos2 θ − 5 cos θ sin θdθ =


Fz = F · z = −R p cos θ sin θdθdφ = −2πR
0 0 8 0
Z 1  1
π 2 2 2

π 2 2 9 4 5 2
− ρR U 9x − 5 xdx = − ρR U x − x = 0,
4 −1 4 4 2 −1

as claimed. In conclusion, F = 0, all components of the force vanish!

3.6 d’Alambert’s paradox

At first sight, one might think that the fact that there is no force acting on a sphere in
a stream U at steady state is related to the high symmetry of the sphere. In fact, this is
not the case: we now show that the force F acting on a body of arbitrary shape vanishes
at steady state. At least as far as the component of F in the direction of U is concerned,
this is a statement about energy conservation, since F · U is the power. Since there is
no friction, this energy input cannot be dissipated; the only other possible option would
be for the energy to be transported off to infinity. The fact that this cannot occur is a
statement about how fast the perturbation introduced by the body into the stream decays
at infinity.
 Figure 30: The force on an arbitrary body in a uniform stream.

Let us write the total velocity field in the form

u = U + v,
where v = ∇φ is the perturbation introduced by the body. As argued before for the case
of the sphere, we will assume that
φ = O r−2 for r → ∞,

since the solution that decays like 1/r is a source or a sink. Now the force on the body is
calculated from integrating the pressure force over the body:
Z
F = − pndS;
S

the minus sign comes from the normal pointing outward, opposite the direction of the
∂φ
pressure force. To find p, we use the Bernoulli equation (62) with = 0, since the flow
∂t
is steady:
p + ρU 2 /2 + ρv2 /2 + ρU · v = patm + ρU 2 /2.
In other words,
Z Z Z
ρ 2
F = − patm ndS + v ndS + ρ (U · v)ndS.
S 2 S S

The first integral on the right is zero, as we have seen before. As a mathematical corollary,
we will show below that Z Z
2
v ndS = 2 v(v · n)dS. (72)
S S
But on the surface of the body the normal component of the total velocity vanishes:
n · (U + v) = 0, so the force becomes
Z Z
F = ρ [v(v · n) + (U · v)n] dS = F = −ρ [v(U · n) − (U · v)n] dS.
S S
In components, this means that
Z
Fi = −ρUj (vi nj − vj ni ) dS. (73)
S

From the symmetry of this expression it is clear immediately that the drag of the body,
which is defined as
D =F·U
is zero. However, the components of F normal to U are also zero. Let V be the volume
exterior to the body but bounded by a closed surface S∞ very far from the body (see
Figure). Then it follows from the divergence theorem that
Z   Z Z
∂vi ∂vj
− dV = − (vi nj − vj ni ) dS − (vi nj − vj ni ) dS;
V ∂xj ∂xi S S∞

note that n points into the volume V . If v decays faster than 1/r2 at infinity, the integral
over S∞ goes to zero as the radius of S∞ goes to infinity. As for the integrand of the
volume integral on the left, it clearly vanishes if i = j. If i 6= j, then the integrand is (up
to a sign) ωk where k 6= i, j. But since the flow is potential, this is again zero. It follows
that the surface integral appearing on the left hand side of (73) is zero, and thus the force
is altogether zero.
The corollary
We now demonstrate (72). From the divergence theorem we conclude that
Z Z Z
2 2
v ni dS = − ∂i v dV − v2 ni dS
S V S∞

The integrand of the second integral on the right behaves like 1/r6 , so its contribution
vanishes. To deal with the volume integral, we note that according to Appendix A,
∇v2 = 2(v · ∇)v + 2v × (∇ × v), and so ∇v2 = 2(v · ∇)v for an irrotational flow. Since
the fluid is incompressible, we also have ((v · ∇)v)i = vj ∂j vi = ∂j (vi vj ). Thus

∂i v2 = 2∂j (vi vj ),

and we can write

Z Z Z Z
2 2
v ni dS = − ∂i v dV = −2 ∂j (vi vj ) dV = 2 vi vj nj dS,
S V V S

using the divergence theorem in the last step. But this proves the corollary.
This is a truly remarkable result of potential flow theory, but it also presents a series
problem, as absence of drag is clearly not in accord with observation. If our argument were
completely airtight, airline companies would surely be doing some serious overcharging!
(maybe they do anyway).
3.7 Drag

Figure 31: A wake of almost stagnant fluid forms behind an obstacle in a strong flow.

The experimental picture of a sphere in a stream reveals what the problem is: the flow
around the sphere is far from rear-aft symmetric. While the flow in front of the sphere
is nicely attached to the body and well described by (70), the back is very different. The
so-called Reynolds number of the flow is
U Rρ
Re = ≈ 1.5 × 104 (74)
η
which is a dimensionless measure of the size of the viscosity η. Since the Reynolds number
is very large, one might believe that viscous effects are small and the potential flow as-
sumption represents a faithful representation. However, this is not the case. In particular,
the flow does not remain attached to the body, and rather there are fluid particles which
originate from the surface of the sphere and the enter the liquid. This is also the place
where vorticity enters the flow. Clearly, vortices have formed in so-called wake of the flow.
Apart from the vortices, the wake appears to be a relatively stagnant region of the
flow, which is thus at constant pressure p = pcav , provided it remains closed (which it
seems to be in the experimental Figure). To estimate the value of pcav , we assume that
the pressure is continuous across the point of separation.
It is also interesting to consider the streamline that separates from the body and
bounds the cavity. Since the pressure is constant in the cavity, the pressure on this
streamline is also constant, by continuity of the pressure. Such a streamline is called
a free streamline. According to the steady Bernoulli equation (44) the velocity along
this streamline is also constant. Therefore, there is a jump of the velocity at the free
streamline to the almost vanishing velocity in the wake. (Nothing prevents the velocity to
be discontinuous if there is no viscosity). Such a surface over which the velocity changes
abruptly is called a shear layer or a vortex sheet. Now pressure in the cavity will be the
pressure on the separating streamline, as it leaves the body.
 To do any calculation of the effect of the wake, one needs an estimate of the point
where the flow leaves the body, to produce the wake. A very simple criterion (which
can be justified in more detail using full viscous theory), is the concept of the “adverse
pressure gradient”. At the front of the body the pressure maximum, and then decreases
to a minimum value at the equator, where the speed is maximum; from there the pressure
decreases again. It is intuitive and can be motivated in greater detail from the viscous
flow theory, that the flow along the body is stable as long as the pressure decreases along
the path of a fluid particle. If however a fluid article has to push against an increasing
pressure (“adverse pressure gradient”), it prefers to leave the surface of the body.
Thus on the basis of the pressure distribution which we have calculated, one can give
an estimate of the place where separation is going to occur, which is at θ = π/2 in the case
of a sphere. Now we can attempt to recalculate the drag force on the sphere, assuming
separation at θs = π/2. The pressure in the cavity will be
5
pcav = patm − ρU 2 :
8
in our calculation, there is an underpressure in the cavity which pulls the sphere along.
We are now in a position to calculate the drag on the sphere, doing the front and the
back separately. The result cannot depend on patm which exerts equal forces on front and
back; therefore, we can put it to zero. Since n · ẑ = cosθ, the force from the front in the
flow direction is
1 2 2 2π π
Z Z Z
(f ront)
Fz =− p(θ) cos θdS = − ρU R sin θ cos θ(9 cos2 θ − 5)dθdφ =
f ront 8 0 π/2
Z 1
π 2 2 π
ρU R x(9x2 − 5)dx = − ρU 2 R2 ,
4 0 16
and from the back
Z 2π Z π/2 Z 1
5π 2 2 5π 2 2
Fz(back) = −R 2
sin θ cos θpcav dθ = ρU R xdx = ρU R .
0 0 4 0 8
Thus the total force or drag is
9π 2 2 9
D = Fz = ρU R = ρAU 2 ,
16 16
where A = πR2 is the projected area of the sphere.
For a general body, the drag coefficient Cd is defined by
ρ
D = U 2 ACd ,
2
and so Cd = 9/8 for our calculation of the sphere. We hasten to add that our calculation
is not particularly good, see the Figure. In the relevant range of Reynolds numbers,
CD ≈ 0.5. Our intention is to give a general idea of where the drag is coming from, and
how d’Alambert’s paradox arises.
Figure 32: The drag coefficient CD of a sphere as a function of Reynolds number.
4 Two-dimensional flows
If the flow is in the plane, and there are only two independent variables x, y, the flow prob-
lem is much simplified. Moreover, we will see that some of the physics is fundamentally
different from three dimensions. One fact we know of already is that no vorticity is created
in two dimensions, ω is simply convected with the flow. Of course, real flows are never
truly two-dimensional; however if the geometry extends very far in one direction (think
of the wing of an aeroplane), a two-dimensional description is a good approximation.

4.1 Flow past a cylinder

Figure 33: The streamlines around a cylinder in a uniform flow.

We want to find the flow around a stationary cylinder in a steady stream U = U x̂. This
is very closely analogous to the flow around a sphere; the effect of the sphere is modelled
by a (two-dimensional) dipole, which is the derivative of a source, whose velocity field is
given by (21). Thus the corresponding potential is
m
φ= ln r, (75)
2π
and a dipole in the µ-direction has potential
µ·r
φ = (µ · ∇) ln r = .
r2
The ansatz for the velocity field is thus
µ 2(µ · r)r
u = U x̂ + ∇φ = U x̂ + − .
r2 r4
Now the boundary condition for r = R is u · n = u · r/R = 0, and from the velocity field
we find (r = R):
x µ·x
u·n=U − .
R R3
Thus the boundary condition is satisfied if we choose

µ = R2 U x̂ = R2 U.

Thus the final answer for the velocity field is

R2
u=U+ [U − 2n(U · n)]. (76)
r2
Now we plug this into Bernoulli’s equation to compute the pressure:

u2 p U 2 patm
+ = + .
2 ρ 2 ρ
On the surface

u2 r=R
= (2U − 2n(U · n))2 = 4U 2 − 4(U · n)2 = 4U 2 (1 − cos2 θ) = 4U 2 sin2 θ,

so in other words
ρU 2
p = patm + (1 − 4 sin2 θ). (77)
2
The pressure is once more symmetric about θ = π/2, i.e. about the midsection of the
cylinder. It follows that the total force on the cylinder is again zero.

4.2 Non-uniqueness of the potential

One particularly important aspect of two-dimensional flow is that any solid body placed
in the flow domain will create a domain that is no longer simply connected.

Definition 4.2.1 (Simply connected domain) A closed curve C is reducible in a do-

main D if it can be shrunk to a point without ever leaving D. If every closed curve is
reducible then D is simply connected.

Example 4.2.2 (Simply connected domains)

1. If D is the interior of a circle, then it is simply connected.

2. If D is the exterior of a circle then it is not simply connected.

3. Let D be the domain a < r < b, 0 < θ < 2π in cylindrical polars. D is not simply
connected.
 i ii

Figure 34: On the left, and example of a simply connected domain: all curves in D can be
collapsed to a point. On the right, a domain which is not simply connected: some curves
can be collapsed to a point, but those going around the hole cannot.

We now place a point vortex (introduced at the beginning of chapter 3) at the centre,
which creates a flow
Γ
u= θ̂. (78)
2πr
The strength of the vortex is measured in terms of its circulation Γ. Since the flow is
irrotational, there exists a velocity potential s.t.
∂φ 1 ∂φ
ur = = 0, uθ = = Γ/(2πr).
∂r r ∂θ
The solution is
Γ
φ(r, θ) = (θ + A), (79)
2π
for any constant A.
The boundary conditions at the outer and inner walls of D are u · n = u · r̂ = ur = 0.
Clearly, (78) satisfies the boundary conditions, while the circulation Γ can assume any
value: the solution is not unique! This is related to another problem with the potential: for
example, φ(r, 0) 6= φ(r, 2π) so the potential is discontinuous along θ = 0 or, alternatively,
φ is multivalued. Evidently, whenever there is circulation
I
Γ= ∇φ · dl 6= 0,
C

a finite amount of potential is picked up as one goes around the curve. The potential of
such a flow is thus necessarily multivalued. This situation can arise only in a multiply
connected domain: if it were simply connected, the closed curve C could be reduced to a
point, and thus the integral would be zero. This can be achieved by placing a branch-cut
in D to make the domain simply connected. In simply-connected domains, φ is unique,
as seen in Fig. 35.
 b

a
Needs a branch cut to
make domain simply connected
and the potential single−valued

Figure 35: A branch cut avoids the non-uniqueness of the potential

Now assume any flow field bewteen the two cylinders, which must satisfy 4φ and
n · ∇φ = 0 on the boundaries. To any such solution, we can add the vortex (79), for any
value of the circulation Γ. The result will be another solution, which also satisfies the
boundary conditions. In other words, the flow between the cylinders can have any value
of the circulation. In general, we see that the flow in a non-simply connected domain is
determined only up to a flow with circulation Γ, where Γ can take any value.

4.3 Steady flow past a circular cylinder with circulation

Now we apply this idea to the flow past a circular cylinder derived in Section 4.1. The
potential for a uniform stream in the x-direction is φ = U x = U r cos θ, and so the total
potential for the steady flow past a cylinder of radius R is

R2
 
φ=U r+ cos θ (80)
r

(cylinder fixed at origin, flow speed U from left to right).

Suppose now that there is circulation round the cylinder, by adding a vortex rotating
in the clockwise direction:
R2
 
Γ
φ=U r+ cos θ − θ. (81)
r 2π

In the problem without circulation (cf. Fig.33), the stagnation points u = (0, 0) are at
the front and the back of the cylinder, i.e. at θ = 0 and θ = π along the surface r = R. To
find the stagnation points on r = R in the general case (81), we have to look for solutions
of uθ = 0 (ur = 0 by construction):

R2
 
1 ∂φ Γ
0 = uθ = = −U 1 + 2 sin θ − .
r ∂θ R 2πR

Thus stagnation points are at angles θ with

Γ
sin θ = − . (82)
4πU R
Equation (82) has two roots in 0 < θ < 2π (and hence two stagnation points) if Γ/U R <
4π (see Fig.36, left). If Γ/U R = 4π, the two stagnation points coincide at θ = 23 π
(middle). On the other hand, if Γ/U R > 4π there are no stagnation points on the surface
(right); rather, the stagnation point has moved into the fluid.

Figure 36: The streamlines around a circular cylinder with circulation Γ < 4πU R, Γ =
4πU R, and Γ > 4πU R.

To calculate the force on the cylinder, we use the unsteady Bernoulli equation,

p + 12 ρu2 = C,

noting that the flow is steady (C(t) = C, ∂φ/∂t = 0). With p → patm , u2 → U 2 as
r → ∞,
p = patm + 12 ρU 2 − 21 ρu2 .
On the surface of the cylinder,
Γ
uθ = −2U sin θ − (83)
2πR
and with n = (cos θ, sin θ), the total force is
Z 2π
F=− pnR dθ.
0

It follows that F · r̂ = Fx = 0, so there is no drag force, as before. However, the

transverse force (or lift) is now
Z 2π
p0 + 21 ρU 2 − 21 ρ(2U sin θ + Γ/2πR)2 sin θR dθ
 
Fy = F · ŷ = −
Z 2π0
1
= 2
ρ4U sin θ(Γ/2πR) sin θR dθ = ρU Γ. (84)
0

Hence there is an upward lift force on the cylinder (Magnus Effect - see Annalen der Physik
164, 129 (1853)). This can explain the forces acting on spinning objects (for example
footballs or ping-pong balls). A quantitative analysis of the resulting trajectories has been
performed by Dupeux et al., New Journal of Physics 12 093004 (2010).
Figure 37: On the left, an illustration of the direction of the Magnus force. On the right,
an application to table tennis.

Note that the above calculation refers to cylinder at rest, with a flow past it. It
is equivalent to a cylinder moving to the left in a stationary medium. Thus a cylinder
moving to the right will experience a downward force and deviate in the direction of
rotation, which is fairly intuitive, as shown in the Figure. Dupeux et al. investigate the
effect more quantitatively, and calculate the resulting spiral flight path. Among others,
they show that the Magnus effect is responsible for a famous goal by the Brazilian player
Roberto Carlos, see Fig. 38.

Figure 38: A spinning football describes an increasingly tight spiral.

The effect is also closely related to the lift forces on airplane wings, to which we will
come back in more detail. The big question in that case of course is what sets the value
of Γ!

4.4 The complex potential

We now introduce the formulation of potential flow using a complex formulation, where
each point in the x−y plane is represented by a complex number z. Only this formulation
reveals the full power of two-dimensional methods. The first remarkable result is that the
potential and the streamfunction are simply the real and imaginary part of the same
complex potential.
Definition 4.4.1 (Complex potential)
Let u = (u, v, 0) be an irrotational, incompressible flow, and z = x + iy. The the
complex potential w is defined as
w(z) = φ(x, y) + iψ(x, y). (85)

Theorem 4.4.2 (Complex potential)

The complex potential is an analytic function.

. Proof: We have that u = ∇φ, and the stream function ψ satisfies (17). Thus

∂φ ∂ψ 
u= = 
∂x ∂y

(86)
∂φ ∂ψ 
v= =− .  
∂y ∂x
But the equations (86) are the Cauchy-Riemann equations for w = φ + iψ, and so w is
analytic. In particular, if w(z) = w(x, y) is analytic then
dw ∂φ ∂ψ ∂φ ∂ψ
= +i = −i + ,
dz ∂x ∂x ∂y ∂y
i.e. w0 (z) is independent of the direction of the derivative.
Properties of the complex potential:
(i) ∇2 φ = ∇2 ψ = 0 follow directly from (86).
∂φ ∂ψ ∂φ ∂ψ
(ii) ∇φ · ∇ψ = + = 0 meaning equipotential curves (where φ = const)
∂x ∂x ∂y ∂y
are perpendicular to streamlines, (where ψ = const).
dw
(iii) = u − iv = qe−iχ where q = (u2 + v 2 )1/2 is the speed of the flow, χ is the angle
dz
dw
the flow makes to the x-axis; is known as the complex velocity (note the minus
dz
sign)!

Example 4.4.3 (Stagnation point (straining) flow)

According to (60), ψ = 2kxy, where k is a constant.

Then 
∂ψ ∂φ
u= = 2kx = 

∂y ∂x

, ⇒ φ = k(x2 − y 2 )
∂ψ ∂φ
v=− = −2ky =



∂x ∂y
Hence w(z) = φ + iψ = k(x2 − y 2 ) + i2kxy = kz 2 is the complex potential.

Example 4.4.4 (Line source)

In three dimensions, a two-dimensional source is a line of sources pushing fluid radially
outward. According to (21), the flow field is
m
u= r̂.
2πr
Now the potential and stream functions are
m m
φ= log r, ψ= θ. (87)
2π 2π
It is recognised immediately that these are the real and imaginary parts of the complex
potential
m
w= log z,
2π
since in polars z = reiθ where r = |z| and θ = arg(z) and so log(z) = log(r) + log(eiθ ) =
log(r) + iθ.

Example 4.4.5 (Line vortex)

We have seen in Section 4.2, that a vortex has the flow field
Γ
u= θ̂,
2πr
and the potential is (cf. (79))
Γ
φ= θ.
2π
But this is the same as the streamfunction (87) of a source! Thus we know directly that
the complex potential of a vortex must be the same as that of a source, up to a complex
rotation:
iΓ
w = − log z.
2π
In particular, we can conclude without further calculation that the stream function is
Γ
ψ = Im{w} = − ln r.
2π
This is a simple example for the application of a powerful theorem from complex analysis:
once one knows either the real or the imaginary part of an analytic function, the whole
complex function is determined completely. In three dimensions, this will be a line vortex
which is completely straight and perpendicular to the xy plane. In general, line vortices
are curved, but we disregard this effect here.

Example 4.4.6 (Flow around a cylinder)

We reconsider the problem treated in Section 4.1, but using complex notation.

We have seen that the flow can be described as a superposition of a uniform stream and
a dipole. First, a uniform stream in the x-direction has the complex potential

w(z) = U z,

since the complex velocity is u − iv = U . Second, a dipole is the derivative of a source,

which has the potential w ∝ ln z. Thus, we arrive at
µ
w(z) = U z −
2πz
for the potential flow around a cylinder. As before, we want to adjust µ such that the
boundary conditions on the surface |z| = R of the cylinder are satisfied. Namely, the
stream function
µy
ψ = ={w} = U y +
2π|z|2
must be a constant for |z| = R. Indeed, if we choose µ = −2πR2 U , ψ = 0 on the surface.
In conclusion, the required complex potential is
R2
 
w(z) = U z + . (88)
z
This agrees with the flow field (76) found earlier; namely
dw U R2 U R2 2 2


u − iv = =U− 2 =U− x − y − 2ixy ,
dz z |z|4
so that
U R2 2 2

2U R2
u=U+ y − x , v = − xy.
r4 r4
On the other hand, since U · x̂ = U and n · x̂ = x/r, (76) is equivalent to
U R2 x2 U R2
 
1 − 2 2 = U + 4 x2 + y 2 − x2


u=U+ 2
r r r
and
U R2  xy 
v= −2 2 .
r2 r

4.5 Interaction of vortices

Figure 39: A cyclone (point vortex) in the atmosphere.

A line vortex intersects a plane in one point, in two dimensions one can thus describe the
vortex by a point (the position z0 = x0 + iy0 of the vortex, and its strength Γ. According
to Kelvin’s theorem, the circulation is conserved, and thus the strength of the vortex is
constant in time. The only thing we need to keep track of is the position of the vortex.
In complex notation, the velocity field generated by a vortex located at z0 is
iΓ d log(z − z0 ) iΓ 1 iΓ z − z0
u − iv = − =− =−
2π dz 2π z − z0 2π |z − z0 |2
Now each vortex moves in the velocity field generated by all the other vortices, and so
dz1
= u − iv|allothervortices
dt
In the simplest case of two vortices of strength Γ1 and Γ2 , located at z1 and z2 , respectively,
the equations of motion are
dz1 iΓ2 z1 − z2 dz2 iΓ1 z2 − z1
= , = . (89)
dt 2π |z1 − z2 |2 dt 2π |z1 − z2 |2

Decomposing into real and imaginary parts, the first equation (89) is

dx1 Γ2 y1 − y2 dy1 Γ2 x1 − x2
=− , = ,
dt 2π (x1 − x2 )2 + (y1 − y2 )2 dt 2π (x1 − x2 )2 + (y1 − y2 )2

and correspondingly for the other vortex. The generalisation for N vortices is obvious
and involves the sum over all the vortices except for the one that is being convected. How
does the motion of two vortices look like ? Let us look at the two possible cases:

4.5.1 Two corotating vortices

Figure 40: Two vortices with the same sense of circulation rotate around a common
center, turning in opposite directions.

The vortex strengths are Γ1 and Γ2 , both being positive (the two vortices rotate coun-
terclockwise). Let the distance between the two vortices be h. The speeds of the two
vortices are then
Γ2 Γ1
U1 = , U2 = ,
2πh 2πh
respectively. From the construction it is clear that at each instant, the vortices move
in a direction perpendicular to the line connecting them. This means they have to go
around in circles, whose centre G is along the line connecting the vortices. The point G
is stationary in space. The some of the radii R1 , R2 must satisfy R1 + R2 = h. Since both
vortices make a full turn in the same time, we must have
R1 U1 Γ2
= = .
R2 U2 Γ1
Thus for example,
Γ2 h
R1 = .
Γ1 + Γ2
For two equal vortices of strength Γ1 = Γ2 = Γ, it follows that R1 = R2 = R. This
means the vortices are always at the opposing sides of a circle of radius R. The time to
complete one cycle (the period) is

2πR 8π 2 R2
T = = ,
U Γ
and so the whole configuration rotates at an angular velocity of
2π Γ
Ω= = .
T 4πR2

4.5.2 Two counterrotating vortices

Figure 41: Two counterrotating vortices are created by moving a paddle in water; they
travel in the same direction.

Now consider the case that the two vortices have opposite sign. This situation is produced
easily by the stroke of a paddle, see Fig. 41. In fact since the total vorticity of the system
should remain zero according to Kelvin’s theorem, typically Γ1 = −Γ2 . After the paddle
is withdrawn, the vortices continue to move according to the other vortex’ velocity field.
Let the distance between the two vortices again be h.
Figure 42: Two vortices with the opposite sense of circulation rotate around a common
center, turning in the same direction.

The vortices now move in the same direction, and rotate about the same centre. The
radii of this circular motion obeys R1 /R2 = −Γ2 /Γ1 , noting that the circulations have
Γ2 h
opposite signs. We now have R2 − R1 = h, and thus R1 = − . Thus when the
Γ1 + Γ 2
vortices are equal and opposite (Γ = Γ1 = −Γ2 ), the radii of both circles tend to infinity.
Γ
In this case the vortices move in a straight line, and translate uniformly at speed U = .
2πh
This is what’s seen in the experimental pictures of two vortices being created behind a
paddle.

4.6 The method of images

4.6.1 A vortex next to a wall

Figure 43: A vortex next to a wall moves parallel to it.

Consider the motion of a vortex located at some position z0 = a + ib in the upper half of
the complex plane. The x axis is the location of a solid wall. Thus we have to satisfy the
boundary condition of no flow through the wall: v = 0 for y = 0. Simply taking the flow
field of the vortex,
iΓ x − a − i(y − b)
u − iv = − , (90)
2π (x − a)2 + (y − b)2
does not satisfy this condition: the flow field must be modified by the presence of the wall.

Figure 44: Streamlines of a vortex next to a wall.

The method of images is a technique that permits to find solutions to Laplace’s equa-
tion which satisfy the right boundary condition on the solid wall. Near the vortex at
z0 , the solution should look like (90). The method of solution is indicated in the figure:
imagine replacing the wall by another image vortex that is placed at an equal distance
on the other side of the wall. Its position is where an image in a mirror would appear.
And equally like a mirror image, the sense of rotation of the image vortex is reversed.
By symmetry it is clear that the v component of the image vortex has opposite sign, so
on the line of symmetry (the locus of the wall), the two contributions cancel and v = 0.
The boundary condition on the wall is satisfied and we have solved our problem! In other
words, the flow problem we want to solve (one vortex in the presence of a wall) is exactly
the same as the flow problem of two vortices (the original vortex and its image) without
the wall.
Let’s check if our reasoning was correct. It is easiest to do the calculation using
complex notation. If the vortex is at z0 , its image is at z0 . Since the image also has the
opposite sense of rotation, there is a minus sign in front of it. The total potential (vortex
+ image) becomes
iΓ iΓ
w(z) = − log(z − z0 ) + log(z − z0 ). (91)
2π 2π
We have to verify that Im{w} = ψ = const for z = x real! Now since ln(z) = ln(z)
(why?) we have
iΓ iΓ
w(z) = log(z − z0 ) − log(z − z0 ).
2π 2π
But if z is real, z = z, and thus w(x) − w(x) = 0, and so ψ = 0 along the wall.
In fact, we have shown that the problem of a vortex near a wall is mathematically
equivalent to the problem of two counterrotating vortices of equal strength, which in the
previous section we showed to move at constant speed in the direction perpendicular to
the line connecting them, in other words, parallel to the wall. Let us check this fact as
well, using complex notation. The vortex in question moves in the flow field of its image,
i.e.
iΓ x − a − i(y + b)
u − iv = ,
2π (x − a)2 + (y + b)2
evaluated at z = a + ib. Thus
 
iΓ −2ib Γ
u − iv = = ,
2π 4b2 4πb

Γ
and the vortex moves with a horizontal speed of , which corresponds to the earlier
4πb
result with h = 2b.
This leads us to the following result, which permits to find the flow generated by any
combination of singularities next to a wall. Let the wall be the real axis of the complex
plain, and let f (z) describe all singularities (vortices, sources, etc.) present in the flow
iΓ µ
domain y > 0. For example, f (z) = ln(z − z0 ) for a vortex and f (z) = −
2π 2π(z − z0 )
for a dipole.

Theorem 4.6.1 (Method of images)

Let w(z) = f (z) would be the complex potential of the flow without the wall. Then

w(z) = f (z) + f (z) (92)

is the potential in the presence of the wall.

Note that in the the second member of (92), we take the complex conjugate of every-
thing except of z itself.
Proof: We have to show that ψ = Im{w(z)} = 0 if z = x is real. But if this is the case
then z = z, and so

w(x) − w(x) = f (x) + f (x) − (f (x) + f (x)) = 0,

which proves the theorem.

Example 4.6.2 (Source next to a wall)

 Figure 45: The streamlines of a source next to a wall.
m
Consider a source of strength m placed at z = ib above a wall. Then f (z) = ln(z − ib)
2π
m
and f (z) = ln(z + ib). Thus the solution we are seeking is
2π
m m
w(z) = (ln(z − ib) + ln(z + ib)) = ln(z 2 + b2 ). (93)
2π 2π
The streamfunction is ψ = ={w}, whose level lines are the streamlines, plotted in Fig. 45.

4.6.2 A vortex pair

Figure 46: Two counterrotating vortices next to a wall. On the left, the configuration; on
the right the experimentally recorded trajectory of one of the vortices.

A pair of counter-rotating vortices of strength Γ is produced by the wings of a plane

moving on a runway. To determine the flow field we have to use the method of images;
the the position of the vortices be (s, h) and (−s, h), as shown in Fig. 46. Let z0 = s + ih
be the position of the right vortex; then −z 0 = −s + ih is the position of the vortex on
the left. Now the complex potential of the vortices in the upper half plane is
iΓ
w= (− ln(z − z0 ) + ln(z + z 0 )) ≡ f (z).
2π
According to our recipe, the total potential in the presence of the wall at y = 0 is
iΓ
w = f (z) + f (z) = (− ln(z − z0 ) + ln(z + z 0 ) + ln(z − z 0 ) − ln(z + z0 )) .
2π
w3
The velocity at z0 due to the other three vortices is , where
dz z=z0
iΓ
w3 = (ln(z + z 0 ) + ln(z − z 0 ) − ln(z + z0 )) .
2π
Thus the velocity at z0 is
   
iΓ 1 1 1 iΓ 1 1 1
u − iv = + − = + −
2π z0 + z 0 z0 − z 0 2z0 2π 2s 2ih 2(s + ih)
s2 h2
 
Γ
= +i 2 .
4π h(s2 + h2 ) s(s + h2 )
Thus the equations of motion for the vortex at (s, h) are
Γ s2 Γ h2
ṡ = , ḣ = − . (94)
4π h(s2 + h2 ) 4π s(s2 + h2 )
Solving the equations is left as an exercise.

4.7 Mappings and transformations

The most powerful tool of complex analysis consists in mapping the flow domain into
another. As soon as I have found a mapping that is conformal, I have transformed the
problem from one geometry to another.

Figure 47: Via a complex mapping f (z) a corner can be mapped into the upper half plane.

Imagine we want to solve a flow problem in a given domain D. This means we are
looking for a function w(z) (the complex potential) with the following properties:
 • w(z) is analytic in the domain D, except perhaps for one or more singularities at z0
(and other places).

• Im{w(z)} = const on the boundaries of the domain.

To find such a complex potential, it is enough to find a mapping into a new domain D1
such that the flow problem can be solved in D1 . We have the following:

Theorem 4.7.1 (Complex mapping)

Let ζ = f (z) be a mapping that is analytic everywhere in D, and which maps onto a
new domain D1 . We seek a complex potential w(z) with prescribed singularities z0 ∈ D,
and boundary condition Im{w(z)} = const on ∂D.
Let the singularity at z0 be mapped to a new place ζ0 = f (z0 ) and the boundary of D is
mapped to the boundary of D1 . Let w1 (ζ) be the solution to the flow problem in D1 (which
has the singularity corresponding to z0 at ζ0 ). Then

w(z) = w1 (f (z)) (95)

is a solution to the original flow problem.

Proof: Since w1 (ζ) is analytic in D1 (apart from singularities), w(z) is analytic in D; at

z0 it has the prescribed singularity. This is the first requirement. In addition, Im{w(ζ)}
is constant on the boundary of D1 . But this means that Im{w(z)} = Im{w1 (f (z))} is
also constant, which meets the second requirement, and proves the theorem.

Example 4.7.2 (A source in a corner)

Figure 48: The flowlines of a source in a corner, given by the complex potential (96).
The function ζ = f (z) ≡ z 2 , z = ζ 1/2 maps the first quadrant onto the half-plane.
Namely, if z = reiθ , then ζ = r2 e2iθ = r2 eiφ . Thus the domain 0 ≤ θ < π/2 is mapped
onto 0 ≤ φ < π. The source at z0 maps to a source at ζ0 = z02 . Now the solution in the
ζ-domain is
m m
w1 (ζ) = log(ζ − ζ0 ) + log(ζ − ζ 0 )
2π 2π
as we have seen before, using the method of images. So
m
w(z) = w1 (z 2 ) = log(z 2 − z02 ) + log(z 2 − z 20 )


(96)
2π
m
= (log(z − z0 ) + log(z + z0 ) + log(z − z 0 ) + log(z + z 0 ))
2π
Thus the result can be written as a superposition of four sources of equal strength, located
at z0 , −z0 , z 0 , and −z 0 . The same result could have been obtained (albeit a bit more
laboriously), using the method of images. First, we mirror at the x-axis, then at the
y-axis. Fig. 48 is for the case z0 = 2 + i.

Example 4.7.3 (Rotation)

Figure 49: The flowlines of a uniform flow around a cylinder, which arrives under an angle
α.

Consider any flow problem, for example the flow around a cylinder, given by (88):

R2
 
w1 (ζ) = U ζ + .
ζ

Let ζ = ze−iα ≡ f (z), which is a rotation of the axes by α in the anticlockwise direction.
Namely, for z = reiθ , ζ = rei(θ−α) .
Then in the z-plane,

R2 iα
 
−iα −iα
w(z) = w1 (ze ) = U ze + e . (97)
z

The flowlines (equipotential lines of ψ = ={w}) are plotted in Fig. 49.

Note 4.7.4 (Calculation of velocities)

Velocities are most easily found from
dw dw1 dζ
= u − iv =
dz dζ dz

E.g. in above u − iv = (u1 − iv1 )e−iα .

4.8 The Joukowski mapping: circles to ellipses

A particularly useful application of the mapping idea concerns the flow around bodies. We
have solved the problem of the flow around a cylinder. Thus if we can find a conformal
mapping between the unit circle and any given shape, we have solve the flow problem
around this shape.

Figure 50: The Joukowski mapping (98) maps a circle to an ellipse.

Consider the (inverse) mapping

b2
z=ζ+ , (98)
ζ
called the Joukowski mapping, and illustrated in Fig. 50. Consider a circle of radius R,
whose surface in the ζ-plane is described by the polar representation ζ = Reiθ . Under the
Joukowski mapping,

b2 −iθ b2 b2
   
iθ
z = Re + e = R + cos θ + i R − sin θ = c cos θ + id sin θ
R R R
is an ellipse with axes length 2c and 2d.
 The semiaxes come out to be
b2 b2
c=R+ , d=R− . (99)
R R

Figure 51: The flowlines areound a circle are mapped to an ellipse by the Joukowski
mapping (98).

Now consider the uniform flow past an ellipse. To model an arbitrary angle α between
the direction of the flow and the semi-major axis, we consider the flow around a cylinder
that approaches the x-axis under an angle α:
R2 iα
 
−iα
w1 (ζ) = U ζe + e .
ζ
In principle, one can find w(z) using the inverse of the Joukowski mapping
1 √
ζ = f (z) = (z + z 2 − 4b2 ),
2
so that w(z) = w1 (f (z)).
However, the resulting expressions are often not so useful. For example, to find the
streamlines, it is much easier to find the streamlines of the w1 (ζ) in the ζ-plane, and then
to transform them using (98). This is how the Fig.51 was produced. If one wants to
calculate the velocity, one uses
dw dw1 1 U (e−iα − R2 eiα /ζ 2 )
u − iv = = = .
dz dζ dz/ dζ (1 − b2 /ζ 2 )
On the cylinder, ζ = Reiθ , so
U (e−iα − eiα e−2iθ ) 2iU sin(θ − α)
u − iv = −2iθ
= iθ . (100)
2 2
(1 − (b /R )e ) (e − (b2 /R2 )e−iθ )
This means there are stagnation points at θ = α and θ = −π + α. This point is where
a streamline leaves the surface. In other words, this streamline (plotted in red) has the
same value of the streamfunction ψ then the surface of the ellipse.
4.9 Lift
Now we want to fly! In principle, we know how to construct the flow around a wing of
arbitrary shape, we only have to find the transformation, starting from a circle. We have
seen already that the key ingredient is to have circulation around the wing. We have
calculated the lift in the case of a cylindrical cross-section, but what is it for an arbitrary
shape? To answer this question, we have the following theorem, which gives the force on
a body of arbitrary shape, using complex notation.

Theorem 4.9.1 (Blasius’ theorem)

Figure 52: If χ is the angle made by the tangent line with the x-axis, the tangent and
normal is given by (103).

Let C be a simple, closed curve, describing a body in a stready stream U . Let the
complex potential outside of the body be w(z). Then the force F = (Fx , Fy ) on the body is
given by 2
I 
iρ dw
Fx − iFy = dz. (101)
2 C dz

Proof: The complex velocity is given by

dw
= u − iv = qe−iχ , (102)
dz
where |u| = q is the speed of the flow, and χ is the angle the flow direction makes to the
horizontal. So, according to Bernoulli (no gravity):

p = p0 − 12 ρq 2 ,

where p0 = patm + ρU 2 /2 is a constant; U is the flow speed at infinity.

But on the surface of body, the flow is in the tangential direction. So χ is also the
angle the tangent line of the body makes with the real axis, and as shown in the above
Figure, we have t = (cos χ, sin χ), n = (sin χ, − cos χ), for the tangent and normal vectors,
respectively. If we interpret both components as real and imaginary parts of a complex
number, this is the same as
t = eiχ , n = −ieiχ , (103)
where t and n are complex numbers representing the tangent and the normal.
 To compute the total force F on the body, we have to do the integral
I
F = − pn ds,
C

where s is the arclength along C. On the other hand, a complex line element along C is
 
dx dy
dz ≡ dx + idy = +i ds = t ds = eiχ ds.
ds ds
Here we have used that by definition, the tangent vector is
 
dr dx dy
t= = , ,
ds ds ds
and (103) in the last step.
Now we can convert everything to complex notation, where it proves convenient to
define a complex force F = Fx − iFy . Then
I I
F = Fx − iFy = − pn ds = −i pe−iχ ds
C IC I
iρ
= −i p0 (dx − idy) + q 2 e−iχ ds
C 2 C
I I  2
iρ −iχ 2 iρ dw
= (qe ) dz = dz
2 C 2 C dz
(the term proportional to p0 vanishes because the integral of a total differential over a
closed loop is zero), and we have used (102) in the last step. This proves Blasius’ theorem.

Example 4.9.2 (Cylinder with circulation) Consider the problem of Section 4.3,
R2 iΓ
w(z) = U z + U − log z,
z 2π
so
dw R2 iΓ
=U −U 2 − .
dz z 2πz
Applying (101) to the circle C gives
2
R2
I 
iρ iΓ
Fx − iFy = U −U 2 − dz
2 C z 2πz
I  
iρ 2 iU Γ −2 −3 −4
= U − + terms {z , z , z } dz
2 C πz
I I
dz dz
Now use Cauchy’s residue theorem = 2πi, n
= 0, n 6= 1 and
C z C z

Fx − iFy = iρU Γ,

which agrees with (84).

Example (4.9.2) shows that Blasius’ theorem yields a much more general result, since
only the residue of the integral comes into play. This is the
Theorem 4.9.3 (The Kutta-Joukowski lift theorem) Consider a body of arbitrary
cross-section C, in a uniform stream U in the x-direction, which generates circulation of
strength Γ. Then the force on the body is given by:

Fx − iFy = iU ρΓ. (104)

Hence the drag (force in the direction of the stream) on an arbitrary body is zero and
the lift force is Flif t = −U ρΓ. Note that the lift force is by definition the force in the
direction normal to the flow.
Proof: Far away from the body,
iΓ
w(z) → U z − log z
2π
(origin inside C) so
dw iΓ
≈U− , z → ∞.
dz 2πz
Since there are no singularities in the flow between C and |z| = R, the integral must be
independent of any chosen closed path outside C:
I  2 I  2
dw dw
dz = lim dz.
C dz R→∞ |z|=R dz
So I  2 I  2 I
dw iΓ iΓU 1
dz = U− ,= − dz = 2U Γ,
C dz |z|=R 2πz π |z|=R z
using the same reasoning as for the cylinder. This proves the theorem.

4.10 Oblique flow past plates

.

Figure 53: The mapping of a circle to a plate, using (105).

We have developed a wonderful theory for lift. The only problem is that we don’t
know how to determine the value of Γ which determines the lift. Indeed, the wings of an
airplane do not have circular cross section: if cylinder is placed in a uniform stream (left
side of the Figure), there is no reason why the flow should turn either way, and thus there
is no circulation. The secret is to fashion the shape so as to induce lift. Namely, the tail
of wing is shaped to have a sharp corner, so as to force the flow to separate at that point.
As the simplest possible model for such a situation consider the flow around a flat
plate, which is placed in a uniform stream at an angle α. A flat plate is obtained from
an ellipse by letting d → 0. According to (99), this amounts to putting b = R in the
Joukowski transformation (98):
R2
z=ζ+ . (105)
ζ
Then the width of the plate is 2c = 4R. As is seen from the Figure, the flow separates
from the plate beyond the trailing edge, so that fluid particles very close to the plate have
to go around the sharp corner before leaving the plate.
This situation is reflected by the velocity on the surface of the plate, which is
U sin(θ − α)
u − iv = ,
sin θ
putting b = R in (100). Clearly, u diverges as as θ → 0, π, i.e. at the sharp edges of the
plate. This is a physically untenable situation; fluid particles are not likely to make such
a sharp turn without leaving the plate altogether, especially not at infinite speed! The
same conclusion can be drawn from the concept of “adverse pressure gradient”, introduced
earlier. The speed is very great at the trailing edge and decreases as one moves up the
plate. According to Bernoulli, this means that the pressure increases, in other words fluid
particles experience an adverse pressure gradient. Instead of following a path of adverse
pressure, fluid particles prefer to leave the body, producing a point of separation at the
trailing edge. This requirement is known as the Kutta condition.

Figure 54: Flow around a plate satisfying the Kutta condition (107).

As we have seen in our calculation of the flow around a cylinder with circulation, the
point of separation can be made to move by adding circulation:
R2 iα
 
−iα iΓ
w1 (ζ) = U ζe + e − log ζ.
ζ 2π
Then the complex velocity of the flow around a plate (b = R) is
−1
R2 iα R2
   
dw1 dζ −iα iΓ
u − iv = = U e − 2e − 1− 2 .
dζ dz ζ 2πζ ζ
Putting ζ = Reiθ , we find the velocity on the plate:
e−iθ
   
−iα+iθ iα−iθ

iΓ Γ 1
u − iv = U e −e − = U sin(θ − α) − . (106)
2πR 1 − e−2iθ 4πR sin θ
Let us focus on the trailing edge, corresponding to θ → 0. (assume for example that the
front of the plate is slightly ‘rounded’, so that separation is not as important there). Thus
we want u − iv to be finite for θ → 0, a requirement which is another incarnation of the
Kutta condition. This can be achieved by choosing Γ such that the square bracket on the
right of (106) vanishes for θ = 0. In other words,
Γ = −4πRU sin α, (107)
which is the Kutta condition for the plate. This is the value chosen for the above Figure,
where the flow is seen to leave the plate smoothly in the direction of the orientation of
the plate. Now the singularity cancels and the velocity becomes
U
u − iv = (sin(θ − α) + sin α) → U cos α, as θ → 0,
sin θ
a finite value!

Figure 55: The experimentally measured lift coefficient for a plate, compared to the
theoretical result (110).

The wonderful thing is that we have now determined the circulation uniquely, so we
can calculate the lift using Blasius’ theorem:
Flif t = −U ρΓ = 4πRρU 2 sin α, (108)
where 4R is the width of the plate. Once more, note that this is the force acting normal
to the flow, so it is at an angle α relative to the orientation of the plate. It is customary
to define a lift coefficient cL by

ρU 2
Flif t = cL Aw , (109)
2
where Aw = 4R is the area of the wing (per unit length), also called the chord. Thus we
have found that for the plate
(plate)
CL = 2π sin α ≈ 2πα (110)

In Fig.55, (110) is compared to experiment, and it works really well, as long as the angle
of attack α is small. If however α becomes too large the theory fails abruptly. The reason
is clear from the two photographs below. As long as α is small, the flow remains laminar
and attached to the wing. As α is too great, the flow separates and a completely different
type of description must be sought, see Fig. 56.

Figure 56: The flow around a wing. At small angles α, the flow separates smoothly. As
α becomes large, the flow separates from the wing and vortices are formed.

Finally, we comment on the presence of circulation around the wing, which is the
crucial ingredient needed for flying. In the first instance, there is nothing wrong with that
from the point of view of a potential flow description. However, where is this circulation
coming from when one imagines starting up the flow, with zero circulation? Since the
net circulation in a large circle around the wing must vanish initially, and the Kutta
condition requires a circulation Γ < 0 around the wing, this means that in the process
of the point of separation moving to the trailing edge, vortices of positive circulation
(rotating counterclockwise) are shed from the wing. Eventually the vortices are convected
downstream, and no longer matter for the problem.
5 Waves and free surface flows
Free surface flows are in some ways very different from what we have done before. In
all problems considered so far, the domain D in which to solve the problem is given (for
example some box or the exterior of an aeroplane wing). A free surface, on the other
hand, moves, so the domain varies in time. The key feature, however, is that the domain
does not vary according to some predetermined program. Instead, it moves in response
to the flow itself, that is in response to the flow solution (which of course itself depends
on the domain).
This leads to solutions which are quite different in character to the fixed-boundary solu-
tions we were considering so far (and generally speaking much more interesting solutions)!
Another nice feature of free surface flows is that they are easy to observe experimentally,
one just needs to track the motion of the free surface.

5.1 Nonlinear free-surface motion

We begin with the general equations for free-surface flow, assuming that the flow inside
the fluid is potential, and the fluid is incompressible. We assume for now that the flow is
2D.

Pressure = patm z ζ(x,t) = z

x
z=0

Density = ρ

z=−h

Figure 57: Schematics for wave motion on a fluid layer of thickness h above a solid wall.

Hence the fluid flow is described by

u = ∇φ, ∇2 φ = 0,

where u = (u, 0, w) and φ = φ(x, z, t). We choose z = 0 to coincide with the undisturbed
free surface, and the bottom of the fluid is at z = −h, see Fig. 57. Let the surface of the
water in motion be given by z = ζ(x, t).
The interesting part of the problem are its boundary conditions:

(i) At the lower boundary z = −h, the vertical velocity must vanish (kineamtic b.c.):

∂φ
0 = n · ∇φ = = w.
∂z
In other words,
∂φ
=0 on z = −h. (111)
∂z
 DS
(ii) On z = ζ(x, t), the kinematic boundary condition on a moving surface is = 0,
Dt
where the free surface is the zero set of S(r, t). If the surface can be represented by
a height function ζ(x, t) (i.e. we are not allowed overhangs) we can define a function
S(x, z, t) = z − ζ(x, t)
which is indeed zero at the free surface. Thus
 
∂ ∂ ∂ ∂ζ ∂φ ∂ζ ∂φ
0= +u +w (z − ζ(x, t)) = − − + on z = ζ,
∂t ∂x ∂z ∂t ∂x ∂x ∂z
so the kinematic boundary condition becomes
∂ζ ∂φ ∂ζ ∂φ
+ = , on z = ζ. (112)
∂t ∂x ∂x ∂z
(iii) The dynamic boundary condition at the free surface is that the pressure equals the
exterior athmospheric pressure: p = patm (const). on z = ζ. The unsteady Bernoulli
equation (62) for a constant force of gravity Φ = gz is
∂φ
p/ρ + 21 |∇φ|2 + + gz = C(t),
∂t
and so the third boundary condition becomes
∂φ
patm /ρ + 21 |∇φ|2 + + gζ = C(t), on z = ζ. (113)
∂t
This concludes our description of the fully nonlinear problem.

5.2 Linear gravity waves

∂ζ
By the small amplitude assumption, |ζ|  h and  1. If we assume that the flow is
∂x
driven by the motion of the free surface (no additional driving from below the surface), it
must be that the flow is also week: |φ|  1. Thus we will throw away all terms quadratic
in the two variables ζ and φ, i.e. containing ζ 2 , φ2 , or products ζφ. Note that (113)
contains the velocity at quadratic (nonlinear) order, which is a remnant of the nonlinear
character of the Euler equation, as it remains in Bernoulli’s equation. The kinematic
∂φ ∂ζ
boundary condition (112) contains , which is quadratic as well. Both terms will be
∂x ∂x
neglected in our linearised treatment.
But this represents only one type of nonlinearity. The most significant problem is that
in the formulation of the boundary conditions, the position of the free surface is part of
the solution, so don’t know where to apply the boundary conditions ! We can deal with
this difficulty by linearising about z = 0 using Taylor’s expansion about z = 0 for all
z-dependent terms, for example
 
∂φ ∂φ ∂ ∂φ
= +ζ + ....
∂t z=ζ ∂t z=0 ∂z ∂t z=0
The second term on the right is quadratic, and will be neglected. Thus we observe that
linearizing the boundary conditions implies that instead of at z = ζ(x, t), the boundary
condition is imposed at the equilibrium position z = 0. Thus we obtain
 ∂φ
(i) On z = −h, =0
∂z
∂ζ ∂φ
(ii) On z = 0, =
∂t ∂z
∂φ
(iii) On z = 0, + gζ = C(t) − patm /ρ.
∂t
The last condition can be prettyfied by redefining the potential according to
Z t
0
φ = φ − patm t/ρ + C(t0 )dt0 .

The extra terms are time dependent only, so they do not affect the flow velocities,
Laplace’s equation or the boundary conditions. So now (dropping primes) (iii) is
replaced with
∂φ
(iii) On z = 0, + gζ = 0.
∂t

The linearised pressure in the fluid is

p(x, z, t) ∂φ
+ + gz = C(t)
ρ ∂t
and because of φ → φ0 this transforms to
p(x, z, t) − patm ∂φ
=− − gz. (114)
ρ ∂t

Now we look for particular solutions to ∇2 φ = 0 with (i), (ii), (iii) motivated by observa-
tions, which describe a propagating wave.
Such a propagating wave of amplitude H is described by a solution of the form

ζ(x, t) = H cos(kx − ωt), (115)

where an arbitrary phase factor δ could have been added as well. The solution (115) is
periodic both in time (with period τ = 2π/ω) and space (with wavelength λ = 2π/k).
Most importantly, it moves to the right with speed c = ω/k, where c is called the phase
speed. To see that, let ξ = x − ct, so that

ζ = H cos(kξ)

is unchanged for ξ constant. But this means that

dξ dx
0= = − c,
dt dt
so c is the speed at which the point x is moving.
∂φ
Now we need the potential φ, which describes the velocity field. Since = −gζ
∂t
reasonable to write
φ(x, z, t) = sin(kx − ωt)Z(z)
for some Z(z) – this is a separation solution φ = X(x)Z(z). Then from 4φ = 0,
−k 2 sin(kx − ωt)Z(z) + sin(kx − ωt)Z 00 (z) = 0
and so Z 00 (z) − k 2 Z(z) = 0. Then
Z(z) = A cosh k(z + h) + B sinh k(z + h),
for constants A, B. Because of (i), need B = 0, so
φ(x, z, t) = A sin(kx − ωt) cosh k(z + h). (116)
We still have (ii) and (iii) to apply, but only A to find. The second condition will lead
to an equation for ω. First from (iii), inserting (115) and (116) into gζ = −∂φ/∂t|z=0 , we
find
gH cos(kx − ωt) = Aω cos(kx − ωt) cosh kh,
so that
A = gH/(ω cosh kh). (117)
In other words, the potential is
gH
φ(x, z, t) = sin(kx − ωt) cosh k(z + h). (118)
ω cosh kh
Then from (ii), inserting (115) and (118) into ∂φ/∂z|z=0 = ∂ζ/∂t results in
kA sin(kx − ωt) sinh kh = ωH sin(kx − ωt).
This implies  
gH
ωH = k sinh kh,
ω cosh kh
or
ω2h
= kh tanh kh, (119)
g
which is called the dispersion relation of the wave, plotted in Fig. 58.
3
2
ω h
g 2.5

1.5

0.5

0
-3 -2 -1 0 1 kh 2 3

Figure 58: The dispersion relation (119).

 The relation (119) is the principal result of our calculation. It establishes how the
frequency of a monochromatic wave is related to its wavelength. Very roughly, large λ ⇒
small k ⇒ small ω ⇒ large τ = 2π/ω. So long wavelengths imply long periods, and vice
versa.
If the dispersion relation between ω and k is non-linear, one speaks of a dispersive
problem: parts of a wave containing different frequencies travel at different speeds, and
thus disperse. If the relation is linear (as for exlectromagnetic waves), the phase speed is
a constant.
We remark that changing k → −k, (119) still holds. Then (115) becomes
ζ(x, t) = H cos(kx + ωt),
which is a wave travelling to the left, with speed c = −ω/k. The constant k = 2π/λ is
called the wavenumber (the number of wavelengths that can be fit into 2π).

5.2.1 Two special cases

Two important special cases of (119) exist:
1. Long wavelengths This is the case kh  1 (or λ/h  1), which means that the
wave length is long compared to the depth. Then tanh kh ≈ kh or
p
ω = ghk,
so the dispersion relation is a linear√function. The constant of proportionality is the
(constant) wave speed c = ω/k = gh. Evidently, the shallow water problem has
no dispersion.
2. Short wavelengths If on the other hand kh  1 (or λ/h  1), the water is deep
relative to the length of the wave. Then tanh kh ≈ 1 and the dispersion relation is
p
ω = gk.
p p
Thus the wave speed c = ω/k = g/k = gλ/2π does depend on the wavelength,
and there is dispersion. Since cosh k(z + h) ≈ ekz ekh /2 the potential (118) becomes
H kz
φ(x, z, t) = e sin(kx − ωt). (120)
c
The amplitude decays exponentially as one goes away from the surface.

5.3 Particle trajectories

Figure 59: The trajectories of particles underneath a wave.

The fact that a traveling wave moves at a speed c = ω/k does not mean that particles
inside the fluid also travel at this speed! Since the amplitude is small, we consider particles
at a position (x + ξ, z + η), where ξ and η are small departures from the mean position
(x, z) of the particle. Using (116), the flow velocities are

u = kA cos k(x+ξ −ct) cosh k(z +η+h), w = kA sin k(x+ξ −ct) sinh k(z +η+h). (121)

Now for small ξ and η, we can put ξ = η = 0 in (121) to obtain the equations of
motion
dξ dη
= kA cos k(x − ct) cosh k(z + h), = kA sin k(x − ct) sinh k(z + h),
dt dt
which we can integrate to yield the trajectory
A A
ξ=− sin k(x − ct) cosh k(z + h), η= cos k(x − ct) sinh k(z + h). (122)
c c
Thus we find that as long as the amplitude A of the wave is small, the motion (122) of
the particles is small as well. Hence we were justified in putting ξ = η = 0 as a first
approximation. In Section 5.4 below we will pursue the expansion to the next (quadratic)
order.
Taking the square of each of the components of (122), we find

ξ2 η2 A2
+ = , (123)
cosh2 k(z + h) sinh2 k(z + h) c2

which is the equation for an ellipse with major and minor semi-axes (in the x and z
directions)
A A
a = cosh k(z + h), a = sinh k(z + h).
c c
Using the expression (117) for A and the dispersion relation (119), A/c = H/ sinh kh, and
so
cosh k(z + h) sinh k(z + h)
a=H , b=H , (124)
sinh kh sinh kh
while the ratio between the two axes is
b
= tanh k(z + h). (125)
a
This means that if the point of observation is far from the bottom (k(z + h)  1),
a ≈ b and the trajectory is a circle. This is clearly visible in Fig. 59, where particle
trajectories are visible as bright lines. The trajectories near the free surface are close to
being circles, while closer to the lower boundary (where z + h becomes small), trajectories
are flattened in the vertical direction. The size of the ellipse decreases as well.

5.4 Stokes drift

A major conclusion of the previous section is that although the wave is propagating, there
is no particle transport, as particles move on closed orbits. G.G. Stokes discovered that
there is in fact a little bit of transport in the wave direction, since particle paths are not
exactly closed. To discover this, one has to go to next order in a perturbation expansion
in ξ, η.
To simplify the discussion, we consider the case (120) of infinite depth, for which the
equation of motion reads
dξ dη
= Acek(z+η) cos k(x + ξ − ct), = Acek(z+η) sin k(x + ξ − ct), (126)
dt dt
gHk gHk 2
where A = = = Hk is a dimensionless number. Expanding for small values
cω ω2
of ξ,η, one finds
dξ
= Ac ekz cos k(x − ct) + ηkekz cos k(x − ct) − ξkekz sin k(x − ct) ,
 
dt
and
dη
= Ac ekz sin k(x − ct) + ηkekz sin k(x − ct) + ξkekz cos k(x − ct) .
 
dt
Inserting the leading order expressions

A A kz
ξ = − ekz sin k(x − ct), η= e cos k(x − ct),
k k
this results in
dξ 2 dη 3
= Acekz cos k(x − ct) + A ce2kz , = Acekz sin k(x − ct) + O(A ).
dt dt
This is of course solved to give

A 2
ξ = − ekz sin k(x − ct) + A ce2kz t, η = Aekz cos k(x − ct). (127)
k
To leading order, the trajectories are circles, and no global transport takes place. The
only, but crucial difference is that at quadratic order, particles are pushed forward an
extra horizontal distance, so the circles do not close. On average, particles are pushed
along with the wave with a drift velocity of
2
US = A ce2kz , (128)

which is known as Stokes drift.

5.5 Oscillations in a container

Liquids readily slosh back and forth in a closed container (e.g. tea in a tea-cup). There
are many associated important practical problems: sloshing in road tankers, water on
decks of ships, resonance in harbours, etc.
 z
z= ζ(x,t)
z=0
x

z=−h

x=0 x=L
Figure 60: Standing waves in a container

As a model for such situations, consider a two-dimensional rectangular box with rigid
walls at x = 0, L and a bottom on z = −h, filled with fluid to z = 0, see Fig. 60. We
use the small-amplitude theory from before, so that linearised equations are ∇2 φ = 0 for
(x, y, z) in box,
∂φ
(i) = 0 on z = −h.
∂z
∂ζ ∂φ
(ii) Kinematic: = on z = 0.
∂t ∂z
∂φ
(iii) Dynamic: = −gζ on z = 0.
∂t
Also need no-flow conditions (kinematic) on walls:
∂φ
(iv) = 0 on x = 0, L;
∂x
However now we no longer have travelling waves, but are interested in describing a
phenomenon that is periodic in time, so we write in two dimensions:

ζ(x, t) = ζ(x) sin ωt

and assume we can separate variables in φ:

φ(x, z, t) = X(x)Z(z) cos ωt.

Then Laplace’s equation gives

X 00 (x)Z(z) cos ωt + Z 00 (z)X(x) cos ωt = 0

which implies
X 00 (x)Z(z) + Z 00 (z)X(x) = 0
This separates:
X 00 (x) Z 00 (z)
=− = −k 2
X(x) Z(z)
where −k 2 is the separation constant.
Solving for Z(z), with (i) gives Z(z) = A cosh k(z + h) as before. Combining (ii) and
(iii) to eliminate ζ (that’s g(ii) +∂/∂t(iii)) gives
∂ 2φ ∂φ
= −g , on z = 0
∂t2 ∂z
which implies
−ω 2 X(x)Z(0) cos ωt = −gX(x)Z 0 (0) cos ωt
and it follows that ω 2 /g = k tanh kh as before. (We expect this, as the vertical structure
of the fluid is independent of the vertical lateral walls)

Now solving for X(x) gives

X(x) = B cos(kx) + C sin(kx)
subject to (from (iv)), X 0 (0) = X 0 (L) = 0. Easy to show that must have C = 0 and
k = nπ/L and so
X(x) = B cos(nπx/L)
So pulling everything together we have
φ(x, z, t) = A0 cos(nπx/L) cosh k(z + h) cos ωt
for some constant A0 = AB while from (iii) the free surface is given by
ωA0 cosh kh
ζ(x, t) = cos(nπx/L) sin ωt
g
and here, k = nπ/L, so that (119) reads
ω 2 /g = (nπ/L) tanh(nπh/L) ≡ ωn2 /g, n = 0, 1, 2, . . .
and therefore defines a set of discrete wave frequencies at which these oscillations may
occur. Here, n is a mode number and tells you how many oscillations are occurring across
the box. Crucially, by considering waves in a finite box, we have selceted a discrete
spectrum of allowed frequencies.
We cannot have n = 0 as then ω = 0 and φ = 1, and ζ = 0. So there is no motion
in the fluid. (In fact, you can discount a sloshing mode with no x-dependence – a flat
surface oscillating up and down – as it would violate mass conservation in the tank).
The fundamental and higher order modes are illustrated in Fig. 60. p The higher√ n, the
higher the frequency. As n → ∞, tanh(nπh/L) → 1 and so ωn → gnπ/L = gk. For
large n, the discreteness of the spectrum hardly plays a role any more: any wavelength is
approximated closely by one of the discrete modes.

n=1 n=2 large n

Figure 61: Different modes in a container

 The Fundamental frequency is the lowest frequency, given here by n = 1.

ω12 L
T1 [s] h [cm] g
π tanh πh
L
4.15 0.7 0.17 0.093
2.7 1.4 0.41 0.19
2.1 2.3 0.68 0.31

The Table gives experimental data for measurements in a rectangular tank (see pic-
ture), whose length was L = 74 cm. One end of the tank was lifted, to excite the fun-
damental sloshing mode. The third and fourth columns compare experiment and theory,
as
ω12 L πh
= π tanh ,
g L
and T1 = 2π/ω1 . Remarkably, the period decreases with increasing depth h. The trend is
reproduced correctly, but there is a consistent discrepancy by a factor of two, perhaps due
to the fact that the tank is too long, so we didn’t really excite the findamental mode...
Appendix A: Vector calculus
We shall revise some vectors operations that you should have already met before this
course. We will be using suffix notation throughout this course

A.1 Suffix notation and summation convention

Suppose that we have two vectors u = (u1 , u2 , u3 ) and v = (v1 , v2 , v3 ). Then the
dot product is defined to be
3
X
u·v = ui vi or, more simply, write u · v = ui vi
i=1

(drop the summation symbol on the understanding that repeated suffices imply summa-
tion.
Defn A.1.2: The Kronecker delta is defined by
 
1, i=j
δij =
0, i 6= j

So in summation convention, δij aj = ai since

3
X
δij aj ≡ δij aj = ai
j=1

Examples:
1. δii = 3

2. δij ui vj = uj vj ≡ u · v.
Defn A.1.3: The antisymmetric symbol ijk is defined by
• 123 = 1

• ijk is zero if there are any repeated suffices. E.g. 113 = 0.

• Interchanging any two suffices reverses the sign: e.g. ijk = −jik = −kji

• Above implies invariant under cyclic rotation of suffices: ijk = jki = kij
With this definition all 27 permutations are defined. There are only 6 non-zero compo-
nents,
123 = 231 = 312 = 1, 213 = 132 = 321 = −1,
Defn A.1.4: The cross product is defined by

w = u × v = (u2 v3 − u3 v2 )r̂ + (u3 v1 − u1 v3 )ŷ + (u1 v2 − u2 v1 )ẑ

But can now be written in component form as

wi = ijk uj vk
where summation over j and k occurs. [Check].

Example: Consider the triple product,

w · (u × v) = wi ijk uj vk

where summation is over i, j, k so result is scalar. It follows that

w · (u × v) = jki wi uj vk = u · (v × w)
= kij wi uj vk = v · (w × u)
= −jik wi uj vk = −u · (w × v)
= −ikj wi uj vk = −w · (v × u)

Defn A.1.5: The double product of ijk is

ijk ilm = δjl δkm − δjm δkl

[Check]

Defn A.1.6: The vector triple product is defined by the result

w × (u × v) = (w · v)u − (w · u)v
Proof:

[w × (u × v)]i = ijk wj [u × v]k

= ijk wj klm ul vm
= kij klm wj ul vm
= (δil δjm − δjl δim )wj ul vm
= wj ui vj − wj uj vi = (w · v)ui − (w · u)vi

True for i = 1, 2, 3, hence result.

A.2 Differential operations

Here we consider operations on a function φ(r) and a vector field
 f (r) where r = (x1 , x2 , x3 ).
∂ ∂ ∂
One can regard ∇ as the vector operator , , . Without using any infor-
∂x1 ∂x2 ∂x3
mation (but always remembering the true meaning of the symbol!), one can also use the
much sleaker notation ∇ ≡ (∂1 , ∂2 , ∂3 ). I will usually do that whenever using summation
convention. Then
∂φ
• The gradient is ∇φ. So [∇φ]i = ≡ ∂i φ
∂xi
∂fi
• The divergence is ∇ · f = ≡ ∂i fi (in summation convention)
∂xi
 ∂fk
• The curl is ∇ × f where [∇ × f ]i = ijk ≡ ijk ∂j fk .
∂xj
Examples:
1. [∇(xj )]i = ∂j xi = δij

2. ∇ · r = ∂i xi = δii = 3

3. [∇ × r]i = ijk ∂k xj = ijk δkj = ijj = 0

4. ∇r = rr , where r2 ≡ r · r.
√ ∂i x2j xj ∂ i xj xj δij xi
Proof: [∇r]i = ∂i xj xj = q = = = .
2 x2j r r r

A.2.1 Useful vector identities

1. ∇ · (φf ) = ∂i (φfi ) = φ∂i fi + fi ∂i φ = φ∇ · f + f · ∇φ

2. ∇ × (φf ) = φ(∇ × f ) + (∇φ × f )

Proof: [∇ × (φf )]i = ijk ∂j (φfk ) = φijk ∂j fk + fk ijk ∂j φ = φ[∇ × f ]i + [∇φ × f ]i .

3. f × (∇ × f ) = ∇( 12 f · f ) − (f · ∇)f

Proof:

[f × (∇ × f )]i = ijk fj [∇ × f ]k = ijk fj klm ∂l fm = kij klm fj ∂l fm

1
= (δil δjm − δim δjl )fj ∂l fm = fj ∂i fj − fj ∂j fi = ∂i fj fj − (f · ∇)fi
2
Hence result.

A.3 Integral results

A.3.1 The divergence theorem
n

ds

V S

Very important. Consider a volume V bounded by a closed surface S with outward unit
normal n. Then for a vector field f (x)
Z Z
∇ · f dV = f · ndS
V S
Corollary: Let f (r) = aφ(r) where a is an arbitrary constant vector, and φ(r) a scalar
function. Since ∇ · (aφ(r)) = a · ∇φ(r), so the divergence theorem reduces to
Z Z
a· ∇φdV = a · φndS
V S

True for any a, so Z Z

∇φdV = φndS
V S

A.3.2 Stokes’ theorem

ds
S

C
dl

Let C be a closed curve bounding a surface S with unit normal n.

Then for a vector field f (r),
I Z
f · dl = (∇ × f ) · ndS
C S

where dl is a line element on C.

Corollary: Let f (r) = aφ(r) where a is an arbitrary constant vector. Then

I Z
aφ · dl = (∇ × aφ) · ndS
C S

and from §A.2.1(2), (∇ × aφ) = φ(∇ × a) + (∇φ × a) = (∇φ × a). Using the vector
triple product result, (∇φ × a) · n = −a · (∇φ × n) and then
I Z
a· φdl = −a · ∇φ × ndS
C S

Therefore I Z
φdl = − ∇φ × ndS
C S
Appendix B: Curvilinear coordinate systems
Many problems can be approached more simply by choosing a coordinate system that
fits a given geometry. Instead of writing the position vector r as a function of Cartesian
coordinates (x, y, z), r is now written as function of three new coordinates: r(q1 , q2 , q3 ).
The coordinate lines are swept out by varying one of the coordinates, keeping the other two
constant. We will deal only with the by far most important case of orthogonal coordinate
systems, in which the coordinate lines always intersect one another at right angles.
∂r
Evidently, is a vector which points in the direction of the i-th coordinate line,
∂qi
see the figure. If each of these vectors are normalized to unity, we obtain the local basis
system:
1 ∂r ∂r
q̂i = , hi ≡ . (129)
hi ∂qi ∂qi
The quantities hi (q1 , q2 , q3 ) are called scale factors or metric coefficients. The fact that
the coordinate system is orthogonal means that all q̂i , computed at a point (q1 , q2 , q3 ),
are orthogonal. However, the direction of q̂i of course changes as one goes along the
coordinate lines.

Example: The cylindrical polar coordinate system.

The position vector is
r = (r cos θ, r sin θ, z) ,
and the coordinates are (r, θ, z). Thus the scale factors become h1 = 1, h2 = r, and
h3 = 1, and the local basis is
∂r
r̂ = = (cos θ, sin θ, 0) ,
∂r
1 ∂r
θ̂ = = (− sin θ, cos θ, 0) ,
r ∂θ
∂r
ẑ = = (0, 0, 1).
∂z
It is clear that r̂, θ̂, ẑ are indeed mutially orthogonal.

Remark: As the coordinates are varied by δq1 , δq2 , and δq3 , respectively, r describes a
volume whose sides are orthogonal. The length of each side is hi δqi , and thus the volume
of the cuboid is d3 x = h1 h2 h3 dq1 dq2 dq3 . Thus integration in a curvilinear coordinate
system is achieved by the formula
Z Z
3
f (r)d x = f (q1 , q2 , q3 )h1 h2 h3 dq1 dq2 dq3 . (130)
V V
In cylindrical polars, the volume element is rdrdθdz.

To do vector calculus in the curvilinear coordinate system, we have to work out what
the ∇-operator is. For any scalar function f (r) (and suspending the summation conven-
tion for the moment), we have
3  
∂f X ∂xj ∂f ∂r
= = · ∇ f = hi (q̂i · ∇) f
∂qi j=1
∂q i ∂x j ∂q i
by the chain rule, and using (129). But up to the factor hi , the object on the right-
hand side is the qi component of ∇f . This means that in the local basis, ∇ has the
representation
q̂1 ∂ q̂2 ∂ q̂3 ∂
∇= + + . (131)
h1 ∂q1 h2 ∂q2 h3 ∂q3
For example in cylindrical polars, the gradient of a scalar function φ is
∂φ 1 ∂φ ∂φ
∇φ = r̂ + θ̂ + ẑ,
∂r r ∂θ ∂z
in agreement with chapter 0.
Now we are almost there. The only but crucial thing that remains to be considered
is the fact that both scale factors and unit vectors depend on the coordinates qi , so they
∂ q̂i
need to be differentiated. In particular, we want to represent in terms of the basis
∂qj
vectors. It is not difficult to write down general expressions for the derivatives of the basis
vectors if only the hi are known, see exercises. However, for a given coordinate system,
it is generally much easier to simply calculate the derivatives directly. For example, in
cylindrical polars
∂r̂ ∂ θ̂
= θ̂, = −r̂,
∂θ ∂θ
and all the other derivatives are zero.
Now it is easy to derive the other formulae in chapter 0. In cylindrical polars, the
Laplacian becomes
! !
∂ θ̂ ∂ ∂ ∂φ θ̂ ∂φ ∂φ
4φ = ∇ · ∇φ = r̂ + + ẑ r̂ + + ẑ =
∂r r ∂θ ∂z ∂r r ∂θ ∂z
!
∂ 2φ 1  1 ∂ 2φ ∂ 2φ
 
∂r̂ ∂φ 1 ∂ θ̂ ∂φ 
(r̂ · r̂) 2 + θ̂ · + 2 θ̂ · + θ̂ · θ̂ 2 2 + (ẑ · ẑ) 2 =
∂r r ∂θ ∂r r ∂θ ∂θ r ∂θ ∂z
∂ 2 φ 1 ∂φ 1 ∂ 2φ ∂ 2φ
+ + + 2.
∂r2 r ∂r r2 ∂θ2 ∂z

In local coordinates, a vector field is written as u = ur r̂ + uθ θ̂ + uz ẑ. Then the

divergence is
!
∂ θ̂ ∂ ∂   ∂u
r ur 1 ∂uθ ∂uz
∇ · u = r̂ + + ẑ · ur r̂ + uθ θ̂ + uz ẑ = + + + ,
∂r r ∂θ ∂z ∂r r r ∂θ ∂z

again in agreement with chpater 0.

Appendix C: The streamfunction
If ∇ · u = 0, then it follows that there exists a vector field A(r, t) s.t.

u = ∇ × A.

Conversely, it is clear from this representation that the flow is incompressible. The rep-
resentation is particularly useful in two dimensions, in which case A can be written in
terms of a single scalar function, the streamfunction. By definition, the streamfunction
is constant along streamlines.
(i) Cartesians: Consider Cartesian coordinates (x, y, z), and take the flow in the
(x, y)-plane. Then A = ψ(x, y, t)ẑ, and

∂ψ ∂ψ
u= , v=− .
∂y ∂x
Now confirm that ψ is indeed constant along streamlines, which are defined by
dr
= λu,
ds
and thus dx = λuds and dy = λvds. It follows that
∂ψ ∂ψ
dψ = dx + dy = vdx − udy = 0,
∂x ∂y
if indeed ψ is varied along streamlines. Thus ψ is constant along a streamlines, as adver-
tised.
Defn: We call the function ψ(x, y, t) the streamfunction of the flow.
ψ= constant

(ii) Cylindrical polar coordinates: Now consider flow in the same two dimensional
plane z = 0, which is represented in cylindrical polars (r, θ, z). The relation to Cartesian
coordinates is x = r cos θ, y = r sin θ. Since the stream function is defined completely by
the geometry of the streamlines, it must be the same in any coordinate system. In other
words if ψc (x, y, t) is the Cartesian version, and ψp (x, y, t) the streamfunction in polar
coordinates, then ψp (r, θ, t) = ψc (r cos θ, r sin θ, t).
As before, the vector potential is A = ψ(r, θ, t)ẑ, and thus from the curl in polar
coordinates
1 ∂ψ ∂ψ
ur = , uθ = − ,
r ∂θ ∂r
where u = ur r̂ + uθ θ̂.
Now verify that ψ(r, θ) is indeed constant along streamlines. In polar coordinates,
r = rr̂, and thus
dr dr dθ dr̂ dr dθ
= r̂ + r = r̂ + r θ̂,
ds ds ds dθ ds ds
 dr̂
where we have used that = θ̂. In other words, we find that dr = λur ds and rdθ =
dθ
λuθ ds. Now as before,
∂ψ ∂ψ
dψ = dr + dθ = −uθ dr + rur dθ = 0,
∂r ∂θ
if ψ is indeed varied along streamlines.
Appendix D Some simple flows and their potentials
D.1 Two dimensional flows
We will give the answer either in Cartesians or polars, depending on which is more con-
venient.
Reminder:
∂φ ∂ψ ∂φ ∂ψ
(i) in Cartesians, u = (u, v, 0) and u = = , v= =− where φ is the
∂x ∂y ∂y ∂x
velocity potential and ψ is the streamfunction.
∂φ 1 ∂ψ 1 ∂φ ∂ψ
(ii) In polars, u = (ur , uθ , 0) and ur = = , uθ = =− .
∂r r ∂θ r ∂θ ∂r
Type of flow flow field u potential φ streamfunction ψ complex pot. w
Uniform stream par- Ur̂ Ux Uy Uz
allel to x axis
U cos αr̂+
Uniform stream at U r cos(θ − α) U r sin(θ − α) U ze−iα
U sin αŷ
angle α to x axis
u = kx k 2 k 2
Stagnation point flow (x − y 2 ) kxy z
v = −ky 2 2
at origin
mr̂ m mθ m
Line source, strength log r ln z
2πr 2π 2π 2π
m, at r = 0

Γθ̂ Γθ Γ iΓ
Line vortex, circula- − log r − ln z
2πr 2π 2π 2π
tion Γ, at r = 0
µ  z  cos θ sin θ µ
Horizontal dipole, − ẑ − 2 r̂ −µ µ −
2πr2 r 2πr 2πr 2πz
strength µ, at r = 0

D.2 Axisymmetric flows

We use cylindrical or spherical polars, whichever is more convenient. In cylindrical coor-
dinates, (r, θ, z), u = (ur , uθ , uz ); in terms of the potential φ(r, z):

∂φ 1 ∂Ψ ∂φ 1 ∂φ
ur = =− , uz = = .
∂r r ∂z ∂z r ∂r
In spherical polar coordinates, (r, θ, ϕ), u = (ur , uθ , uϕ ); in terms of the potential φ(r, θ):

∂φ 1 ∂Ψ 1 ∂φ 1 ∂Ψ
ur = = 2 , uθ = =− .
∂r r sin θ ∂θ r ∂θ r sin θ ∂r
Type of flow flow field u potential φ coordinate system
Uniform stream U ẑ Uz cylindrical
aligned with axis
of symmetry
k
ur = r k 2
Stagnation point 2 (r − 2z 2 ) cylindrical
uz = −kz 4
flow at origin
m m
Point source, r̂ − spherical
4πr2 4πr
strength m, at
r=0
µ  z  cos θ
Dipole, strength µ, − ẑ − 3 r̂ −µ spherical
4πr3 r 4πr2
in ẑ-direction