III Hydrodynamic Stability - Absolute and convective instabilities

2Absolute and convective instabilities

III Hydrodynamic Stability

2 Absolute and convective instabilities

So far, we have only been considering perturbations that are Fourier modes,

, p

, ρ

) = [

u(z), ˆp(z), ˆρ(z)]e

i(k·x−ωt

This gives rise to a dispersion relation

(

k, ω

) = 0. This is an eigenvalue problem

for

(

), and the

th mode is unstable if the solution

(

) has positive imaginary

part.

When we focus on Fourier modes, they are necessarily non-local. In reality,

perturbations tend to be local. We perturb the fluid at some point, and the

perturbation spreads out. In this case, we might be interested in how the

perturbations spread.

To understand this, we need the notion of group velocity. For simplicity,

suppose we have a sum of two Fourier modes of slightly different wavenumbers

± k

∆

. The corresponding frequencies are then

± ω

∆

, and for small

∆

, we

may approximate

∆

= k

∆

∂ω

∂k



k=k

We can then look at how our wave propagates:

η = cos[(k

+ k

∆

)x − (ω

+ ω

∆

)t] + cos[k

− k

∆

)x − (ω

− ω

∆

)t]

= 2 cos(k

∆

x − ω

∆

t) cos(k

x − ω

= 2 cos

∆

w −

∂ω

∂k



k=k

cos(k

x − ω

Since

∆

is small, we know the first term has a long wavelength, and determines

the “overall shape” of the wave.

As time evolves, these “packets” propagate with group velocity

∂ω

∂k

. This

is also the speed at which the energy in the wave packets propagate.

In general, there are 4 key characteristics of interest for these waves:

– The energy in a wave packet propagates at a group velocity.

– They can disperse (different wavelengths have different speeds)

– They can be advected by a streamwise flow.

– They can be unstable, i.e. their amplitude can grow in time and space.

In general, if we have an unstable system, then we expect the perturbation to

grow with time, but also “move away” from the original source of perturbation.

We can consider two possibilities — we can have convective instability, where

the perturbation “goes away”; and absolute instability, where the perturbation

doesn’t go away.

We can imagine the evolution of a convective instability as looking like this:

whereas an absolute instability would look like this:

Note that even in the absolute case, the perturbation may still have non-zero

group velocity. It’s just that the perturbations grow more quickly than the group

velocity.

To make this more precise, we consider the response of the system to an

impulse. We can understand the dispersion relation as saying that in phase space,

a quantity χ (e.g. velocity, pressure, density) must satisfy

D(k, ω) ˜χ(k, ω) = 0,

where

˜χ(k, ω) =

∞

−∞

∞

−∞

χ(x, t)e

−i(kx−ωt)

dx dt.

Indeed, this equation says we can have a non-zero (

k, ω

) mode iff they satisfy

(

k, ω

) = 0. The point of writing it this way is that this is now linear in

˜χ

, and

so it applies to any χ, not necessarily a Fourier mode.

Going back to position space, we can think of this as saying



−i

∂

∂x

, i

∂

∂t



χ(x, t) = 0.

This equation allows us to understand how the system responds to some external

forcing F (x, t). We simply replace the above equation by



−i

∂

∂x

, i

∂

∂t



χ(x, t) = F (x, t).

Usually, we want to solve this in Fourier space, so that

D(k, ω) ˜χ(k, ω) =

F (k, ω).

In particular, the Green’s function (or impulse response) is given by the response

to the impulse

ξ,τ

(

x, t

) =

(

x − ξ

)

(

t − τ

). We may wlog assume

= 0,

and just call this

. The solution is, by definition, the Green’s function

(

x, t

satisfying



−i

∂

∂x

, i

∂

∂t



G(x, t) = δ(x)δ(t).

Given the Green’s function, the response to an arbitrary forcing

is “just” given

χ(x, t) =

G(x − ξ, t − τ)F (ξ, τ) dξ dτ.

Thus, the Green’s function essentially controls the all behaviour of the system.

With the Green’s function, we may now make some definitions.

Definition (Linear stability). The base flow of a system is linearly stable if

lim

t→∞

G(x, t) = 0

along all rays

= C.

A flow is unstable if it is not stable.

Definition

(Linearly convectively unstable)

An unstable flow is linearly con-

vectively unstable if lim

t→∞

G(x, t) = 0 along the ray

= 0.

Definition

(Linearly absolutely unstable)

An unstable flow is linearly absolutely

unstable if lim

t→∞

G(x, t) 6= 0 along the ray

= 0.

The first case is what is known as an amplifier, where the instability grows

but travels away at the same time. In the second case, we have an oscillator .

Even for a general

, it is easy to solve for

˜χ

in Fourier space. Indeed, we

simply have

˜χ(k, ω) =

F (k, ω)

D(k, ω)

To recover χ from this, we use the Fourier inversion formula

χ(x, t) =

(2π)

˜χ(k, ω)e

i(kx−ωt)

dk dω.

Note that here we put some general contours for

and

instead of integrating

along the real axis, so this is not the genuine Fourier inversion formula. However,

we notice that as we deform our contour, when we pass through a singularity,

we pick up a multiple of

i(kx−ωt)

. Moreover, since the singularities occur when

(

k, ω

) = 0, it follows that these extra terms we pick up are in fact solutions to

the homogeneous equation

(

−i∂

, i∂

)

= 0. Thus, no matter which contour

we pick, we do get a genuine solution to our problem.

So how should we pick a contour? The key concept is causality — the

response must come after the impulse. Thus, if

(

x, t

) = 0 for

t <

0, then we

also need

(

x, t

) = 0 in that region. To understand this, we perform only the

temporal part of the Fourier inversion, so that

˜χ(k, t) =

2π

F (k, ω)

D(k, ω; R)

−iωt

dω.

To perform this contour integral, we close our contour either upwards or down-

wards, compute residues, and then take the limit as our semi-circle tends to

infinity. For this to work, it must be the case that the contribution by the

circular part of the contour vanishes in the limit, and this determines whether

we should close upwards or downwards.

If we close upwards, then

will have positive imaginary part. So for

−iωt

not to blow up, we need

t <

0. Similarly, we close downwards when

t >

0. Thus,

if we want

to vanish whenever

t <

0, we should pick our contour so that it

it lies above all the singularities, so that it picks up no residue when we close

upwards. This determines the choice of contour.

But the key question is which of these are causal. When

t <

0, we require

(

k, t

)

0. By Jordan’s lemma, for

t <

0, when performing the

integral, we

should close the contour upwards when when perform the integral; when

t >

we close it downwards. Thus for

(

k, t

) = 0 when

t <

0, we must pick

to lie

above all singularities of ˜χ(k, ω).

Assume that

has a single simple zero for each

. Let

(

) be the corre-

sponding value of ω. Then by complex analysis, we obtain

˜χ(k, t) = −i

F [k, ω(k)]

∂

∂ω

[k, ω(k)]

−iω(k)t

We finally want to take the inverse Fourier transform with respect to x:

χ(x, t) =

2π

F [k, ω(k)]

∂

∂ω

[k, ω(k)]

−iω(k)t

dk.

We are interested in the case

(

)

(

), i.e.

(

k, ω

) = 1. So the central

question is the evaluation of the integral

G(x, t) = −

2π

exp(i(kx − ω(k)t)

∂

∂ω

[k, ω(k)]

dk.

Recall that our objective is to determine the behaviour of

t → ∞

with

fixed. Since we are interested in the large

behaviour instead of obtaining

exact values at finite

, we may use what is known as the method of steepest

descent.

The method of steepest descent is a very general technique to approximate

integrals of the form

H(t) =

−i

2π

f(k) exp (tρ(k)) dk,

in the limit t → ∞. In our case, we take

f(k) =

∂

∂ω

[k, ω(k)]

ρ (k) = i (kV − ω(k)) .

In an integral of this form, there are different factors that may affect the

limiting behaviour when

t → ∞

. First is that as

gets very large, the fast

oscillations in

causes the integral to cancel except at stationary phase

∂ρ

∂k

= 0.

On the other hand, since we’re taking the exponential of

, we’d expect the

largest contribution to the integral to come from the

where

(

) is the greatest.

The idea of the method of steepest descent is to deform the contour so that

we integrate along paths of stationary phase, i.e. paths with constant

, so that

we don’t have to worry about the oscillating phase.

To do so, observe that the Cauchy–Riemann equations tell us

∇ρ

· ∇ρ

= 0,

where in

∇

we are taking the ordinary derivative with respect to

and

viewing the real and imaginary parts as separate variables.

Since the gradient of a function is the normal to the contours, this tells us

the curves with constant

are those that point in the direction of

. In other

words, the stationary phase curves are exactly the curves of steepest descent of

Often, the function

has some stationary points, i.e. points

∗

satisfying

(

∗

) = 0. Generically, we expect this to be a second-order zero, and thus

looks like (k − k

∗

)

near k

∗

. We can plot the contours of ρ

near k

∗

as follows:

∗

where the arrows denote the direction of steepest descent of

. Note that

since the real part satisfies Laplace’s equation, such a stationary point must

be a saddle, instead of a local maximum or minimum, even if the zero is not

second-order.

We now see that if our start and end points lie on opposite sides of the ridge,

i.e. one is below the horizontal line and the other is above, then the only way to

do so while staying on a path of stationary phase is to go through the stationary

point.

Along such a contour, we would expect the greatest contribution to the

integral to occur when ρ

is the greatest, i.e. at k

∗

. We can expand ρ about k

∗

ρ(k) ∼ ρ(k

∗

) +

∂

∂k

∗

)(k −k

∗

)

We can then approximate

H(t) ∼

−i

2π

f(k)e

tρ(k)

dk,

where we are just integrating over a tiny portion of our contour near

∗

. Putting

in our series expansion of ρ, we can write this as

H(t) ∼

−i

2π

f(k

∗

tρ(k

∗

)

exp



∂

∂k

∗

)(k −k

∗

)



dk.

Recall that we picked our path to be the path of steepest descent on both sides

of the ridge. So we can paramretrize our path by K such that

(iK)

∂

∂k

∗

)(k −k

∗

)

where K is purely real. So our approximation becomes

H(t) ∼

f(k

∗

tρ(k

∗

)

2π

tρ

∗

)

−ε

−K

dK.

Since

−K

falls of so quickly as

gets away from 0, we may approximate that

integral by an integral over the whole real line, which we know gives

√

. So our

final approximation is

H(t) ∼

f(k

∗

tρ(k

∗

)

2πtρ

∗

)

We then look at the maxima of

(

) along these paths and see which has the

greatest contribution.

Now let’s apply this to our situation. Our ρ was given by

ρ(k) = i(kV − ω(k)).

So k

∗

is given by solving

∂ω

∂k

∗

) = V. (∗)

Thus, what we have shown was that the greatest contribution to the Green’s

function along the

direction comes from the the modes whose group velocity

is V ! Note that in general, this k

∗

is complex.

Thus, the conclusion is that given any

, we should find the

∗

such that (

∗

)

is satisfied. The temporal growth rate along

= V is then

σ(V ) = ω

∗

) − k

∗i

This is the growth rate we would experience if we moved at velocity

. However,

it is often more useful to consider the absolute growth rate. In this case, we

should try to maximize

. Suppose this is achieved at

max

, possibly

complex. Then we have

∂ω

∂k

max

) = 0.

But this means that

∂ω

∂k

is purely real. Thus, this maximally unstable mode

can be realized along some physically existent V = c

We can now say

– If ω

i,max

< 0, then the flow is linearly stable.

– If ω

i,max

> 0, then the flow is linearly unstable. In this case,

◦ If ω

0,i

< 0, then the flow is convectively unstable.

◦ If ω

0,i

> 0, then the flow is absolutely unstable.

Here

is defined by first solving

∂ω

∂k

(

) = 0, thus determining the

that

leads to a zero group velocity, and then setting

(

). These quantities

are known as the absolute frequency and absolute growth rate respectively.

Example.

We can consider a “model dispersion relation” given by the linear

complex Ginzburg–Landau equation



∂

∂t

+ U

∂

∂x



χ − µχ − (1 + ic

)

∂

∂x

χ = 0.

This has advection (

), dispersion (

) and instability (

). Indeed, if we replace

∂

∂t

↔ −iω and

∂

∂x

↔ ik, then we have

i(−ω + Uk)χ − µχ + (1 + ic

χ = 0.

This gives

ω = Uk + c

+ i(µ − k

We see that we have temporal instability for

|k| <

√

, where we assume

µ >

This has

= U + c

On the other hand, if we force ω ∈ R, then we have spatial instability. Solving

− i) + Uk + (iµ − ω) = 0,

the quadratic equation gives us two branches

−U ±

− 4(iµ − ω)(c

− i)

2(c

− i)

To understand whether this instability is convective or absolute, we can complete

the square to obtain

ω = ω

(k −k

)

where

= 2(c

− i)

2(i − c

)

= −

4(1 + c

)

+ i



µ −

4(1 + c

)



These k

are then the absolute wavenumber and absolute frequency.

Note that after completing the square, it is clera that

is a double root at

. Of course, there is not nothing special about this example, since

was

defined to solve

∂ω

∂k

) = 0!

I really ought to say something about Bers’ method at this point.

In practice, it is not feasible to apply a

function of perturbation and see

how it grows. Instead, we can try to understand convective versus absolute

instabilities using a periodic and switched on forcing

F (x, t) = δ(x)H(t)e

−iω

where H is the Heaviside step function. The response in spectral space is then

˜χ(k, ω) =

F (k, ω)

D(k, ω)

D(k, ω)(ω − ω

)

There is a new (simple) pole precisely at

on the real axis. We can invert

for t to obtain

˜χ(k, t) =

2π

−iωt

D(k, ω)(ω − ω

)

dω =

−iω

D(k, ω

)

−iω(k)t

(ω(k) − ω

)

∂

∂ω

[k, ω(k)]

We can then invert for x to obtain

χ(x, t) =

−iω

2π

ikx

D(k, ω)

| {z }

(x,t)

2π

i[kx−ω(k)t]

[ω(k) − ω

]

∂

∂ω

[k, ω(k)]

| {z }

(x,t)

The second term is associated with the switch-on transients, and is very similar

to what we get in the previous steepest descent.

If the flow is absolutely unstable, then the transients dominate and invade

the entire domain. But if the flow is convectively unstable, then this term is

swept away, and all that is left is

(

x, t

). This gives us a way to distinguish

between

Note that in the forcing term, we get contributions at the singularities

(

Which one contributes depends on causality. One then checks that the correct

result is

(x, t) = iH(x)

i(k

(ω

)x−ω

∂D

∂k

(ω

), ω

]

− iH(−x)

i(k

−

(ω

)x−ω

∂D

∂k

−

(ω

), ω

]

Note that

(

) may have imaginary parts! If there is some

such that

−k

(

)

0, then we will see spatially growing waves in

x >

0. Similarly, if

there exists some

such that

−k

−

(

)

0, then we see spatially growing

waves in x < 0.

Note that we will see this effect only when we are convectively unstable.

Perhaps it is wise to apply these ideas to an actual fluid dynamics problem.

We revisit the broken line shear layer profile, scaled with the velocity jump, but

allow a non-zero mean U :

z = 1

z = −1

U = −1 + U

U = z + U

U = 1 + U

III

−α(z−1)

αz

+ Ce

−αz

α(z+1)

We do the same interface matching conditions, and after doing the same compu-

tations (or waving your hands with Galilean transforms), we get the dispersion

relation

4(ω − U

α)

= (2α − 1)

− e

−4α

It is now more natural to scale with

rather than ∆

2, and this involves

expressing everything in terms of the velocity ratio

∆U

. Then we can write

the dispersion relation as

D(k, ω; R) = 4(ω −k)

− R

[(2k −1)

− e

−4k

] = 0.

Note that under this scaling, the velocity for

z < −

1 is

= 1

−R

. In particular,

R <

1, then all of the fluid is flowing in the same direction, and we might

expect the flow to “carry away” the perturbations, and this is indeed true.

The absolute/convective boundary is given by the frequency at wavenumber

for zero group velocity:

∂k

) = 0.

This gives us

= k

−

[2k

− 1 + e

−4k

Plugging this into the dispersion relations, we obtain

[2k

− 1 + e

−4k

] − [(2k

− 1)

− e

−4k

] = 0.

We solve for

, which is in general complex, and then solve for

to see if it

leads to

0,i

0. This has to be done numerically, and when we do this, we see

that the convective/absolute boundary occurs precisely at R = 1.

Gaster relation

Temporal and spatial instabilities are related close to a marginally stable state.

This is flow at a critical parameter

with critical (real) wavenumber and

frequency: D(k

, ω

; R

) = 0 with ω

c,i

= k

c,i

= 0.

We can Taylor expand the dispersion relation

ω = ω

∂ω

∂k

; R

)[k −k

We take the imaginary part

∂ω

∂k

, R

)(k

− k

) +

∂ω

∂k

, R

For the temporal mode, we have k

= 0, and so

(T )

∂ω

∂k

, R

)(k

− k

For the spatial mode, we have ω

= 0, and so

0 =

∂ω

∂k

, R

)[k

− k

] +

∂ω

∂k

, R

(S)

Remembering that c

∂ω

∂k

. So we find that

(T )

= −c

(S)

This gives us a relation between the growth rates of the temporal mode and

spatial mode when we are near the marginal stable state.

Bizarrely, this is often a good approximation when we are far from the

marginal state.

Global instabilities

So far, we have always been looking at flows that were parallel, i.e. the base

flow depends on

alone. However, in real life, flows tend to evolve downstream.

Thus, we want to consider base flows U = U(x, z).

Let

be the characteristic wavelength of the perturbation, and

be the

characteristic length of the change in

. If we assume

ε ∼



1, then we may

want to perform some local analysis.

To leading order in

, evolution is governed by frozen dispersion relation at

each

εx

. We can then extend notions of stability/convective/absolute to

local notions, e.g. a flow is locally convectively unstable if there is some

such

that ω

i,max

(X) > 0, but ω

0,i

(X) < 0 for all X.

However, we can certainly imagine some complicated interactions between

the different region. For example, a perturbation upstream may be swept

downstream by the flow, and then get “stuck” somewhere down there. In general,

we can define

Definition

(Global stability)

A flow is globally stable if

lim

t→∞

(

x, t

) = 0 for

all x.

A flow is globally unstable if there is some

such that

lim

t→∞

(

x, t

)

→ ∞

For a steady spatially developing flow, global modes are

χ(x, t) = φ(x)e

−iω

The complex global frequency

is determined analogously to before using

complex integration.

It can be established that

G,i

≤ ω

0,i,max

. This then gives a necessary

condition for global instability: there must be a region of local absolute instability

within the flow.

Sometimes this is a good predictor, i.e.

' R

. For example, with mixing

layers, we have

= 1

315 while

= 1

34. On the other hand, it is sometimes

poor. For example, for bluff-body wakes, we have Re

= 25 while Re

= 48.5.