J. Fluid Mech. (2008), vol. 603, pp. 271–304.
c
2008 Cambridge University Press
doi:10.1017/S0022112008001109 Printed in the United Kingdom
271
Convective instability and transient growth
in flow over a backward-facing step
H. M. BLACKBURN
1
,
D. BARKLEY
2
AND S. J. SHERWIN
3
1
Department of Mechanical and Aerospace Engineering, Monash University,
Victoria 3800, Australia
2
Mathematics Institute, University of Warwick, Coventry CV4 7AL, UK, and
Physique et M
´
ecanique des Milieux H
´
et
´
erog
`
enes, Ecole Sup
´
erieure de Physique et Chimie Industrielles
de Paris, (PMMH UMR 7636-CNRS-ESPCI-P6-P7), 10 rue Vauquelin, 75231 Paris, France
3
Department of Aeronautics, Imperial College London, SW7 2AZ, UK
(Received 30 July 2007 and in revised form 13 February 2008)
Transient energy growths of two- and three-dimensional optimal linear perturbations
to two-dimensional flow in a rectangular backward-facing-step geometry with
expansion ratio two are presented. Reynolds numbers based on the step height
and peak inflow speed are considered in the range 0–500, which is below the value for
the onset of three-dimensional asymptotic instability. As is well known, the flow has
a strong local convective instability, and the maximum linear transient energy growth
values computed here are of order 80×10
3
at Re = 500. The critical Reynolds number
below which there is no growth over any time interval is determined to be Re =57.7
in the two-dimensional case. The centroidal location of the energy distribution for
maximum transient growth is typically downstream of all the stagnation/reattachment
points of the steady base flow. Sub-optimal transient modes are also computed and
discussed. A direct study of weakly nonlinear effects demonstrates that nonlinearity
is stablizing at Re = 500. The optimal three-dimensional disturbances have spanwise
wavelength of order ten step heights. Though they have slightly larger growths than
two-dimensional cases, they are broadly similar in character. When the inflow of the
full nonlinear system is perturbed with white noise, narrowband random velocity
perturbations are observed in the downstream channel at locations corresponding
to maximum linear transient growth. The centre frequency of this response matches
that computed from the streamwise wavelength and mean advection speed of the
predicted optimal disturbance. Linkage between the response of the driven flow and
the optimal disturbance is further demonstrated by a partition of response energy
into velocity components.
1. Introduction
Flow over a backward-facing step is an important prototype for understanding the
effects of separation resulting from abrupt changes of geometry in an open flow setting.
The geometry is common in engineering applications and is used as an archetypical
separated flow in fundamental studies of flow control (e.g. Chun & Sung 1996), and of
turbulence in separated flows (e.g. Le, Moin & Kim 1997), which may further be linked
to the assessment of turbulence models (e.g. Lien & Leschziner 1994). The backward-
facing step geometry is also an important setting in which to understand instability
of a separated flow. However, the linear instability of the basic laminar flow in such
272 H. M. Blackburn, D. Barkley and S. J. Sherwin
Space
Time
(a)(b)(c)
Figure 1. Schematic of absolute and convective instabilities. An infinitesimal perturbation,
localized in space, can grow at a fixed location leading to an absolute instability (a) or decay
at a fixed points leading to a convective instability (b). In inhomogeneous, complex geometry
flow we can also observe local regions of convective instability surrounded by regions of stable
flow (c).
a geometry is not properly understood. While well-resolved numerical computations
by Barkley, Gomes & Henderson (2002) have determined to high accuracy both the
critical Reynolds number and the associated three-dimensional bifurcating mode for
the primary global instability for the case with expansion ratio of two, these results
have little direct relevance to experiment. Only through careful observation has it been
possible to see evidence of the intrinsic unstable three-dimensional mode (Beaudoin
et al. 2004). This is because the numerical stability computations determined one
type of stability threshold (of an asymptotic, or large time, global instabiliy) whereas
the flow is actually unstable at much lower Reynolds numbers to a different type of
instability (transient local convective instability). Moreover, the dynamics associated
with the two types of instability are very different for this flow. In the present work we
investigate directly the linear convective instability in this fundamental non-parallel
flow by means of transient-growth computations.
To understand the issues in a broader context as well as with respect to the work
presented in this paper, it is appropriate to review and contrast different concepts and
approaches in (linear) hydrodynamic stability analysis. In all types of linear stability
analysis one starts with a flow field U, the base flow, and considers the evolution
of infinitesimal perturbations u
to the base flow. The evolution of perturbations is
governed by linear equations (linearized about U). Generally speaking, if infinitesimal
perturbations grow in time, the base flow is said to be linearly unstable. However, one
must distinguish between absolute and convective instability (Huerre & Monkewitz
1985). If an infinitesimal perturbation to parallel shear flow, initially localized in
space, grows at that fixed spatial location (figure 1a), then the flow is absolutely
unstable. If on the other hand, the perturbation grows in magnitude but propagates
as it grows such that the perturbation ultimately decays at any fixed point in space
(figure 1b), then the instability is convective.
In practice, one is often interested in inhomogeneous flow geometries in which
there is a local region of convective instability surrounded upstream and downstream
by regions of stability (figure 1c). For illustration, we indicate a backward-facing-step
geometry, but many other inhomogeneous open flows, such as bluff-body wakes,
behave similarly. A localized perturbation initially grows, owing to local flow features
Transient growth in flow over a backward-facing step 273
near the step edge, and simultaneously advects downstream into a region of stability
where the perturbation decays.
At this point, we must distinguish between different current research directions in
hydrodynamic stability analysis of open flows. In one approach, arguably initiated by
the work of Orr and Sommerfeld, we think primarily in terms of velocity profiles U (y)
(streamwise velocity as a function of the cross-stream coordinate) and analyse the
stability of such profiles. The profiles may be analytic or may result from numerical
computation and likewise the stability analysis may be analytic or it may contain a
numerical component. Such an approach is called local analysis. In most practical
cases of interest, however, the base flow is not a simple profile depending on a
single coordinate. In problems in which the base flow does not vary too rapidly as a
function of streamwise coordinate (i.e. U(x,y)), we can legitimately consider a local
analysis of each profile (Huerre & Monkewitz 1990). Therefore through local sectional
analysis, we can identify regions in which the flow is locally stable or unstable, and
if unstable, whether the instability is locally convective or absolute. It is sometimes
possible to extend these local analyses to a global analysis and even in some cases
to a nonlinear analysis. However, rapid variations in flow geometry typically result
in base flows which are either far from parallel or which do not vary slowly with
streamwise coordinate, or both.
There is a second, largely distinct, approach to hydrodynamic stability analysis. In
this approach we use fully resolved computational stability analysis of the flow field
(see e.g. Barkley & Henderson 1996). We call this direct linear stability analysis in
analogy to the usage direct numerical simulation – DNS. We have the ability to fully
resolve in two or even three dimensions the base flow, e.g. U(x,y,z), and to perform a
global stability analysis with respect to perturbations in two or three dimensions, e.g.
u
= u
(x,y,z, t). We typically do not need to resort to any approximations beyond
the initial linearization other than perhaps certain inflow and outflow conditions.
In particular, we can consider cases with rapid streamwise variation of the flow.
By postulating modal instabilities of the form: u
(x,y,z, t)=
˜
u(x, y, z) exp(λt)orin
the case where the geometry has one direction of homogeneity in the z-direction
u
(x,y,z, t)=
˜
u(x, y) exp(iβz + λt), absolute instability analysis becomes a large-
scale eigenvalue problem for the global modal shape
˜
u and eigenvalue λ. There are
algorithms and numerical techniques which allow us to obtain leading (critical or near
critical) eigenvalues and eigenmodes for the resulting large problems (Tuckerman &
Barkley 2000). This approach has been found to be extremely effective for determining
global instabilities in many complex geometry flows, both open and closed (Barkley
& Henderson 1996; Blackburn 2002; Sherwin & Blackburn 2005) including weakly
nonlinear stability (Henderson & Barkley 1996; Tuckerman & Barkley 2000).
Direct linear stability analysis has not been routinely applied to convective
instabilities that commonly arise in problems with inflow and outflow conditions. One
reason is that such instabilities are not typically dominated by eigenmodal behaviour,
but rather by linear transient growth that can arise owing to the non-normality of the
eigenmodes. A large-scale eigenvalue analysis simply cannot detect such behaviour.
(However, for streamwise-periodic flow, it is possible to analyse convective instability
through direct linear stability analysis, see e.g. Schatz, Barkley & Swinney 1995.)
This brings us to a third area of hydrodynamic stability analysis known generally
as non-modal stability analysis or transient growth analysis (Butler & Farrell 1992;
Trefethen et al. 1993; Schmid & Henningson 2001; Schmid 2007). Here, the linear
growth of infinitesimal perturbations is examined over a prescribed finite time
interval and eigenmodal growth is not assumed. Much of the initial focus in this
274 H. M. Blackburn, D. Barkley and S. J. Sherwin
area has been on large linear transient amplification and the relationship of this
to subcritical transition to turbulence in plane shear flows (Farrell 1988; Butler &
Farrell 1992). However, as illustrated in figure 1 (c), the type of transient response
due to local convective instability in open flows is nothing other than transient
growth. While this relationship has been previously considered in the context of the
Ginzberg–Landau equation (Chomaz, Huerre & Redekopp 1990; Cossu & Chomaz
1997; Chomaz 2005), it has not been widely exploited in the type of large-scale direct
linear stability analyses that have been successful in promoting the understanding of
global instabilities in complex flows. This appoach has been employed in flow over a
backward-rounded step (Marquet, Sipp & Jacquin 2006). Also, with different emphasis
from the present approach, Ehrenstein & Gallaire (2005) have directly computed
modes in boundary-layer flow to analyse transient growth associated with convective
instability.
This paper has two related aims. The first, more specific, aim is to accurately
quantify the transient growth/convective instability in the flow over a backward-
facing-step flow with an expansion ratio of two. As stated at the outset, despite the
long-standing interest in this flow, there has never been a close connection between
linear stability analysis and experiments in the transition region for this flow. We
shall present results that should be observable experimentally.
The second aim is more general. We are of the opinion that large-scale direct
linear analysis provides the best, if not the only, route to understanding instability
in geometrically complex flows. This potency has already been demonstrated for
global instabilities. We believe the backward-facing-step flow studied in this paper
demonstrates the ability of direct linear analysis to also capture local convective
instability effects in flows with non-trivial geometry. Within the timestepper-
based approach, this merely requires a change of focus from computing the
eigensystem of the linearized Navier–Stokes operator to computing its singular value
decomposition (Barkley, Blackburn & Sherwin 2008). We also retain the ability to
perform a full complement of time-integration-based tasks within the same code-base,
in particular fully nonlinear simulations.
2. Problem formulation
In this section, we present the equations of interest, but largely avoid issues of
numerical implementation until §3. Since the numerical approach we use is based on
a primitive variables formulation of the Navier–Stokes equations, with inflow and
outflow boundary conditions, our exposition is directed towards this formulation.
Apart from such details, the mathematical description of the optimal growth problem
found here follows almost directly from the treatments given by Corbett & Bottaro
(2000), Luchini (2000) and Hœpffner, Brandt & Henningson (2005). The objective is
to compute the energy growth of an optimal linear disturbance to a flow over a given
time interval τ.
2.1. Geometry and governing equations
Figure 2 illustrates our flow geometry and coordinate system. From an inlet channel of
height h, a fully developed parabolic Poiseuille flow encounters a backward-facing step
of height h. The geometric expansion ratio between the upstream and downstream
channels is therefore two. We choose to fix the origin of our coordinate system at the
step edge. The geometry is assumed homogeneous (infinite) in the spanwise direction.
The other geometrical parameters, namely the inflow and outflow lengths, L
i
and L
o
,
Transient growth in flow over a backward-facing step 275
h
h
L
i
L
o
z
x
y
Figure 2. Flow geometry for the backward-facing step. The origin of the coordinate system
is at the step edge. The expansion ratio is two. The inflow and outflow lengths, L
i
and L
o
,are
not to scale.
are set to ensure that the numerical solutions are independent of these parameters. As
part of our convergence study in §3.4 we have found acceptable values to be L
i
=10h
and L
o
=50h for the range of flow conditions considered.
We work in units of the step height h and centreline velocity U
∞
of the incoming
parabolic flow profile. This defines the Reynolds number as
Re ≡ U
∞
h/ν,
where ν is the kinematic viscosity, and means that the time scale is h/U
∞
.
The fluid motion is governed by the incompressible Navier–Stokes equations,
written in non-dimensional form as
∂
t
u = −(u ·∇)u − ∇p + Re
−1
∇
2
u in Ω, (2.1a)
∇·u =0 inΩ, (2.1b)
where u(x,t)=(u, v, w)(x,y, z,t) is the velocity field, p(x,t) is the kinematic (or
modified) pressure field and Ω is the flow domain, such as illustrated in figure 2. In
the present work, all numerical computations will exploit the homogeneity in z and
require only a two-dimensional computational domain.
The boundary conditions are imposed on (2.1) as follows. First, we decompose the
domain boundary as ∂Ω = ∂Ω
i
∪ ∂Ω
w
∪ ∂Ω
o
, where ∂Ω
i
is the inflow boundary,
∂Ω
w
is the solid walls, i.e. the step edge and channel walls, and ∂Ω
o
is the outflow
boundary. At the inflow boundary we impose a parabolic profile, at solid walls we
impose no-slip conditions, and at the downstream outflow boundary we impose a
zero traction outflow boundary condition for velocity and pressure. Collectively, the
boundary conditions are thus
u(x,t)=(4y − 4y
2
, 0, 0) for x ∈ ∂Ω
i
, (2.2a)
u(x,t)=(0, 0, 0) for x ∈ ∂Ω
w
, (2.2b)
∂
x
u(x,t)=(0, 0, 0),p(x,t)=0 forx ∈ ∂Ω
o
. (2.2c)
2.2. Base flows and linear perturbations
The base flows for this problem are two-dimensional time-independent flows U(x,y)=
(U(x, y),V(x, y)), therefore U obeys the steady Navier–Stokes equations
0 = −(U ·∇)U − ∇P + Re
−1
∇
2
U in Ω, (2.3a)
∇·U =0 inΩ, (2.3b)
where P is defined in the associated base-flow pressure. These equations are subject
to boundary conditions (2.2)
Our interest is in the evolution of infinitesimal perturbations u
tothebaseflows.
The linearized Navier–Stokes equations governing these perturbations are found by
