## Abstract

The linear stability of confined, periodic, parallel fluid flows is examined. The flow fields considered consist of a steady pressure gradient-driven velocity field combined with a purely oscillatory component generated by either an oscillatory pressure gradient or by harmonically oscillating bounding surfaces. Plane channel and circular pipe geometries are studied and all possible combinations of the steady and oscillatory flow components investigated. Neutral stability curves and critical conditions for instability are computed for a selection of steady–unsteady velocity ratios, channel half-widths and pipe radii. The results obtained confirm previous investigations into the effects of small amounts of periodic modulation on the linear stability of the underlying steady flow, but provide much more comprehensive information on the linear stability regions of unsteady parallel flows in channels and pipes.

## 1. Introduction

The stability characteristics of time-periodic modulated steady shear layers are of considerable theoretical and practical interest. Time-modulation may help control the stability of steady boundary layers and other flows of industrial importance. Oscillatory plane Couette and Poiseuille flows are model problems for the investigation of more complicated time-periodic shear flows and Davis (1976) provides an early review of the stability theory relating to such flows. There have also been experimental investigations into the stability of unsteady flow in a circular pipe (Clamen & Minton 1977) for which physiological considerations were a strong motivating factor (Nerem & Seed 1972).

Plane Couette, plane Poiseuille and Hagen–Poiseuille flows are all well-known unidirectional steady solutions of the Navier–Stokes equations. Another familiar solution is the classical Stokes layer which is the paradigm purely oscillatory flow that is generated when an infinite flat plate oscillates harmonically in a semi-infinite region of viscous fluid. Inducing a small amount of oscillation in the above steady flows, either by superimposing an oscillatory pressure gradient or by allowing the bounding surfaces to oscillate parallel to the existing motion, creates a time-modulated unidirectional flow. In either case, an approximation to the classical Stokes layer is formed at the flow boundaries for sufficiently high frequencies of oscillation. The effects of a small component of oscillation on the linear stability of the underlying steady flow have been comprehensively investigated and a quick summary of the main results follows.

The stability of modulated plane Couette flow has been previously examined by Kelly & Cheers (1970) and von Kerczek (1976). Both studies considered the flows generated when the two channel walls oscillate out of phase with each other, for varying frequencies of oscillation and modulation amplitudes. It was found that decay rates of disturbances are increased for a given parameter range so that time-modulation seems to have a stabilizing effect. Rather than re-visiting this perturbed plane Couette flow, the aim of the present work is to investigate the linear stability of both plane Poiseuille and Hagen–Poiseuille flows when combined with a purely oscillatory component. Moreover, the size of the oscillatory flow component is not restricted to be small when compared with the steady part, and the full range of parallel pressure gradient-driven flows in a plane channel and a circular pipe, extending from completely steady flows through to oscillatory flow with zero mean component, is examined. As this unsteady flow has two identifiable components, our intention is to map out the regions of the basic flow parameter space where the resulting flow is linearly stable. If the maximum steady velocity is denoted by *U*_{1} and the amplitude of the oscillatory velocity by *U*_{2}, then a useful parameter for describing these combined flows is *Γ*=*U*_{1}/*U*_{2} (and this is formally defined in equation (2.2) below). It is also convenient to associate a Reynolds number to each flow component, and so we define
1.1
as the Reynolds number for the steady and oscillatory components of the flow, respectively. In these expressions, *ω* is the frequency of the imposed oscillations, *d* is the channel half-width or pipe radius as appropriate and *ν* denotes the kinematic viscosity of the fluid. We note that a complete set of dimensionless parameters for the basic flow is provided by *Γ*, *R*_{1} and *R*_{2}.

Studies relevant to the current investigation include those of Grosch & Salwen (1968) based on numerical time integration of a Galerkin approximation to the governing equations; Herbert (1972), who employed energy methods; Hall (1975), who used perturbation methods based on Floquet theory and von Kerczek (1982), where numerical techniques also based on Floquet theory were implemented. These authors looked at the stability of oscillatory plane Poiseuille flow from the standpoint that the oscillatory component could be thought of as a perturbation of the underlying steady flow. Grosch & Salwen (1968) solved the system of disturbance equations for oscillatory plane Poiseuille flow and made stability calculations for several values of the Reynolds number, wavenumber, oscillation frequency and modulation amplitude. For relatively large values of *Γ* time-modulation appeared to have a stabilizing effect but when *Γ* is less than roughly 10 destabilization occurred. Herbert (1972) provided further qualitative information by examining the disturbance energy balance in modulated plane Poiseuille flow. He concluded that the oscillatory component of the basic flow inhibited energy transfer in some parameter regions and thus had a stabilizing effect on perturbations.

Using a Floquet theory approach to determine stability characteristics, von Kerczek (1982) re-examined oscillatory plane Poiseuille flow and found the modulation to have a stabilizing effect for values of *Γ* greater than about 4—a result that differs quite significantly from the findings of Grosch & Salwen (1968). Additional numerical computations by von Kerczek (1982) have shown that for very large and small frequencies of oscillation, the unsteady component has a slightly destabilizing effect. This conclusion is in qualitative agreement with that drawn by Hall (1975), who predicted that very high frequency modulations should promote disturbance growth.

At the time of the above studies, purely oscillatory flow and Stokes layers were thought to be strictly linearly stable for all Reynolds numbers (Yang & Yih 1977). Some semi-analytical calculations of Hall (1978), based upon Floquet theory, found only stable disturbance modes for Reynolds numbers *R*_{2} less than about 160 in the case of a Stokes layer in a semi-infinite fluid. In a similar study of an oscillating plane wall beneath a stationary parallel plate, von Kerczek & Davis (1974) also failed to detect unstable behaviour for Reynolds numbers up to approximately 400. Both of these investigations were unable to explore larger Reynolds numbers owing to the limitations of the computational resources available at the time. More recently, spectral methods have been used to solve the governing system of equations in a direct numerical simulation (DNS) study by Akhavan *et al.* (1991*a*). They considered both linear and nonlinear instabilities in oscillatory channels but again could only find linearly stable disturbance modes for *R*_{2} up to 500.

Analytical work by Cowley (1987) and Hall (2003) lent some support to the conclusions of these various theoretical studies. They suggested that there are no linearly unstable Floquet disturbances in the inviscid limit , and the quasi-steady calculations by Cowley (1987) also indicated that there might well be no linearly unstable modes for the flat Stokes layer irrespective of the size of the Reynolds number. However, this scenario is somewhat at odds with a welter of measurements of transition from a laminar to a turbulent state at Reynolds numbers as low as 150, but generally around 250 (Hino *et al.* 1976; Clamen & Minton 1977; Akhavan *et al.* 1991*b*; Lodahl *et al.* 1998). The theoretical and experimental results could then only be consistent should transition occur via some strictly nonlinear instability mechanism, or alternatively via a linear transient growth followed by nonlinear development, as suggested by Wu & Cowley (1995). However, both of these scenarios deny the existence of a purely linearly unstable perturbation to the underlying governing equations.

Contrary to the above evidence, some more recent theoretical investigations by Blennerhassett & Bassom (2002, 2006) (hereafter referred to as BB02 and BB06, respectively) have identified linearly unstable Floquet modes in Stokes layers and in oscillating channel and pipe flows. BB02 called upon the methods of Hall (1978) in their analysis of the semi-infinite Stokes layer and were able to examine the properties of disturbances at Reynolds numbers far greater than previously considered. Instability was found to first occur for *R*_{2} approximately 708. The later work in BB06 also determined regions of instability for oscillatory flow in channel and circular pipe geometries. In this second investigation, the linearized system of governing disturbance equations was numerically solved using methods based upon the pseudospectral techniques of Fornberg (1996) and Trefethen (2000). As would be expected, for sufficiently large channel widths and pipe diameters the disturbance characteristics were essentially equivalent to those of the semi-infinite Stokes layer. Further, these results have been independently confirmed by the DNS calculations of Luo & Wu (2010).

Following these discoveries on the stability of oscillatory flows with zero mean, it would appear to be worthwhile to re-visit the earlier studies on modulated shear layers to gain further understanding of the stability properties of these unsteady flows. Floquet theory and the numerical methods used by BB06 have been employed to determine the stability attributes for oscillatory plane Poiseuille flow and for oscillatory Hagen–Poiseuille flow, for varying frequencies of oscillation and velocity ratios *Γ*.

DNS methods were also used to obtain solutions of the linearized Navier–Stokes equations for the case of disturbances in oscillatory plane Poiseuille flow, thereby providing some independent verification of the results from the Floquet analysis. The DNS approach used was based on that originally developed by Davies & Carpenter (1997) but modified to suit the current basic flow. This scheme was chosen partly because the related velocity–vorticity formulation of Davies & Carpenter (2001) has already been successfully applied to the unbounded Stokes layer by Thomas *et al.* (2010), who obtained excellent agreement with BB02 for neutral stability conditions. In the Davies & Carpenter (1997, 2001) technique, a Chebyshev spectral method is used for the discretization of the normal-to-wall derivatives, while time-marching is achieved using a semi-implicit finite difference procedure.

There are four sections in the remainder of this paper. In the following §2, the basic flow profile for the plane channel flow is derived and the two numerical methods outlined. First, we describe the Floquet analysis in which the disturbance equations are solved numerically using pseudospectral methods and this is followed by a brief account of the DNS formulation of Davies & Carpenter (1997, 2001). Neutral stability curves and critical conditions are presented in §3 for a selection of channel half-widths and a wide range of velocity ratio *Γ*. The extension of the investigation to oscillatory pipe flow is tackled in §4 and, by varying *Γ*, we can describe the stability of a whole family of flows varying between purely oscillatory pipe flow at one end of the spectrum to steady Hagen–Poiseuille flow at the other. The paper closes with some discussion and concluding remarks.

## 2. Disturbances in oscillatory channel flow

Consider a two-dimensional viscous incompressible fluid that is bounded by two parallel plates separated by a distance 2*d*. A flow is generated by the superposition of two forcing terms: a mean pressure gradient that results in a steady flow with maximum velocity *U*_{1} and the oscillation of the boundary walls that generates an unsteady velocity component. The two bounding plates move sinusoidally along the *x*-direction with velocities , where *ω* is the frequency of oscillations. A dimensionless basic velocity profile is obtained by scaling all lengths on the Stokes layer thickness , velocities on the sum of the maximum of the steady and unsteady velocity profiles, *U*_{1}+*U*_{2}, and by defining the non-dimensional time *τ*=*ωt*. If the *y*-direction is oriented normal to the plate surfaces, then the undisturbed basic non-dimensional flow is given by
2.1a
where
2.1b
and
2.1c
Here ℜ denotes the real part and is the scaled channel half-width, while
2.2
are non-dimensional parameters that specify the relative magnitudes of the steady and oscillatory velocity parts. Varying the quantity *Γ* enables us to explore the stability of an entire family of flows from the purely oscillatory channel flow of BB06 (*Γ*=0) to steady plane Poiseuille flow ().

### (a) Floquet method

The linear stability of the basic profile (2.1) is determined by considering a disturbed flow of the form
2.3
where *ε*≪1 and *Ψ* denotes the stream function of a two-dimensional disturbance. Since Squire’s theorem has been extended to unsteady flows (Conrad & Criminale 1965; von Kerczek & Davis 1974), the above perturbation is sufficient for locating the critical conditions for linear instability. The disturbance stream function *Ψ* is expressed in the exponential form
2.4
where *a* is the wavenumber, *ψ*(*y*,*τ*) is taken to be 2*π*-periodic in *τ* and c.c. denotes the complex conjugate. In this representation, any exponential growth or decay of *Ψ* is incorporated in the Floquet exponent *μ*. Thus, the governing equation for *ψ*, when linearized in *ε*, may be reduced to
2.5a
subject to the boundary conditions
2.5b
where the prime denotes differentiation with respect to *y* and
2.6
Here, the Reynolds number *R* is defined by
2.7
In this formulation, *a* is real and *μ* generally complex, with the imaginary part *μ*_{i} of *μ* taken to lie within the interval [−0.5,0.5].

In the Floquet approach, the unknown stream function is decomposed into harmonics
2.8
so that equating coefficients of harmonics in equation (2.5) results in the infinite system of ordinary differential equations
2.9
where *u*_{1} and *u*_{−1} are as defined in equation (2.1c).

The system of equation (2.9) was solved numerically using the pseudospectral techniques described by Fornberg (1996) and Trefethen (2000). The differential operators appearing in equation (2.9) were replaced by their respective pseudospectral matrix approximations and each *ψ*_{n}(*y*) was represented as a vector *ψ*_{n} of its function values on a Chebyshev mesh over the interval −*h*≤*y*≤*h*. The introduction of the matrix operators,
2.10a
2.10b
induces a rearrangement of the general governing equation (2.9)
2.11
where **I** is the identity matrix and denotes the complex conjugate of **M**.

A finite system of equations was obtained by truncating the Fourier series (2.8) for *ψ* and setting *ψ*_{n}≡0 for all |*n*|>*N*. It was found empirically that *N* was required to be greater than roughly 0.8*aR*/(1+*Γ*) in order to ensure that the computed neutral stability curves were independent of the value of *N*. If *N* was too small, the above method would tend to overestimate the growth rates corresponding to the eigenvalues *μ*, thus incorrectly reducing the value of *R* for critical instability. In order to obtain values for *μ* that were independent of the spatial discretization, it was also necessary to consider a minimum of *K*=100 data points in the pseudospectral representation of the differential operators on the interval −*h*≤*y*≤*h*.

The system of equation (2.11) was cast as an algebraic eigenvalue problem
2.12
where **A** is a sparse matrix and the vector ** ϕ** is given by
2.13
The eigenvalues

*μ*and eigenvectors

**were obtained using the eigenvalue routine**

*ϕ*`eigs`in Matlab.

### (b) A velocity–vorticity formulation

The results of our Floquet computations are described in §3 below but, in order to provide an independent check of our findings, a DNS was also implemented. The DNS approach used here was based upon the velocity–vorticity formulation of Davies & Carpenter (2001) and the numerical methods of Davies & Carpenter (1997) for a channel flow geometry. As for the Floquet approach only an overview of the salient points of the formulation and methods is presented here, but the interested reader is referred to the above papers for further detailed descriptions.

The total velocity and vorticity fields for the two-dimensional system can be written as
2.14
where *U*_{B} is the base flow (equation (2.1)), *u*_{x} and *u*_{y} are the streamwise and normal velocity perturbations and *ζ* is the vorticity perturbation. Assuming the perturbation quantities are periodic, with wavenumber *a* in the *x*-direction, the governing equations are then the vorticity transport and Poisson relationships
2.15a
and
2.15b
while the streamwise velocity is given by
2.16
which has been derived from the definition of the vorticity perturbation *ζ*. (Note that for simplicity, we have abused the notation slightly and replaced the real perturbation quantity *ζ*(*x*,*y*,*τ*) by *ζ*(*y*,*τ*)e^{iax} to obtain the above equations. Similar replacements have been used for *u*_{x} and *u*_{y}.) As has been kindly pointed out by a referee, it is possible to eliminate *ζ* from equations (2.15*a*,*b*) leaving a single fourth-order equation for *u*_{y}. The merits and drawbacks of this approach have been discussed elsewhere (for example, Davies & Carpenter (1997, 2001)), so are not repeated here.

The usual no-slip conditions *u*_{x}=*u*_{y}=0 apply at the walls *y*=±*h*. However, in order to create an initial perturbation to the basic flow, a time-dependent forcing *η* is introduced at time *τ*=0 to give the wall conditions
2.17 *a*,*b*
with
2.18
(Notice that large values of the parameter *σ* creates an almost impulsive disturbance in the wall position.) The substitution of equation (2.17*a*) into the definition (2.16) leads to the integral constraint on the vorticity
2.19
which replaces the no-slip condition (2.17*a*).

Numerical approximations to the governing equations (2.15) and relevant boundary constraints (2.17*b*) and (2.19) were obtained by expanding *ζ* and *u*_{y} in terms of even Chebyshev polynomials
2.20
where *T*_{k} denotes the *k*th-Chebyshev polynomial of the first kind. Here, *ξ* is the computational normal coordinate defined as
2.21
which maps the physical domain −*h*≤*y*≤*h* onto the computational domain −1≤*ξ*≤1. Note that *u*_{x} does not directly appear in the governing equations (2.15*a*,*b*); if need be it can be calculated using the definition (2.16).

The governing equations (2.15*a*,*b*) were then integrated twice with respect to the mapped coordinate *ξ*, which removed *y*-derivatives and replaced them with *ξ*-integral operators. This yields the system
2.22a
and
2.22b
where the integral operator **I** is defined by
2.23
for *k*=2,…,*M*. The substitution of equation (2.20) for the variables *ζ* and *u*_{y} into the integrated governing equations (2.23*a*,*b*) followed by a comparison of coefficients of *T*_{2(k−1)} generates a system of 2(*M*−1) ordinary differential equations for the 2*M* unknown variables *ζ*_{k} and *u*_{yk} with 1≤*k*≤*M*. The remaining two equations are furnished by the requirements (2.17*b*) and (2.19), which may be written as
2.24
where
2.25
A time-stepping procedure was implemented using a predictor–corrector method for the convective terms in the vorticity equation with the remaining terms treated implicitly. Further details are available in Davies & Carpenter (1997, 2001) and Thomas (2007).

The results of the DNS calculations are functions of the number of Chebyshev polynomials, *M*, retained in equation (2.20) and the time-step, Δ*τ*, in the time-marching scheme. After some experimentation, it was found that using *M*≥100 and gave results that were independent of these two artificial parameters. Note that the dependence of Δ*τ* on *R* mimics the dependence of the number of harmonics *N* retained in the eigenvalue technique as both quantities relate to the accuracy of the temporal representation of the disturbance.

## 3. Results

Figure 1*a* displays the time histories (obtained using the DNS formulation) of impulsively excited disturbances at Reynolds numbers of 520, 529 and 540. Here, the wavenumber *a*=0.3, the channel half-width *h*=6 and velocity flow profile is defined with *Γ*=0.5. The vorticity perturbation at the wall, that is the real part of *ζ*(−*h*,*τ*), is plotted as a function of *τ*/2*π* along with the corresponding envelopes, ±|*ζ*(−*h*,*τ*)|. The temporal development is given over four cycles of oscillation, with the early transients of the time history omitted in order to demonstrate the development of the Floquet mode structure as time increases. During each cycle, the perturbation tends to grow and develops two distinct maximums, located near *τ*=*π*/2 and *π*, before decreasing towards the end of the period. This process is then repeated but the magnitudes of these maximums either increase or diminish depending on the net temporal growth or decay of the disturbance. The results pertaining to *R*=520 clearly demonstrate an overall temporal decay of the disturbance while net temporal growth is evident when *R*=540. To the accuracy of the plot, neither net growth nor decay is observed for the intermediate *R*=529, so this value is very close to neutral stability conditions for the wavenumber *a*=0.3.

The temporal growth *μ*_{r} of a disturbance was determined by calculating
3.1
for a range of values of *T*. At a given wavenumber *a* neutral conditions, *R*_{N}(*a*), were then obtained via linear interpolation on a pair of points either side of the neutral curve with . In figure 1*b*, the results of the DNS calculations (for the purely oscillatory channel flow with half-width *h*=6) are shown as open circles for selected wavenumbers, while the solid line in this figure corresponds to the smoothed neutral curve, as defined in BB06, that is obtained using the Floquet theory outlined in §2*a*. To the accuracy of the plot, the DNS prediction of conditions for neutral stability are indistinguishable from those obtained through Floquet analysis. Hence, the DNS method described in §2 provides confirmation, additional to that of Luo & Wu (2010), of the growing Floquet modes calculated by BB06 in oscillatory channels.

Our two techniques were applied to the case of plane Poiseuille flow (corresponding to ) in order to provide further confidence in the numerical codes. However, as the Reynolds number for the purely steady plane Poiseuille flow *R*_{1} is related to the combined Reynolds number according to equations (1.1) and (2.7), some care must be exercised when interpreting the values of *R*_{N} that are predicted for neutral conditions. Further, as the length scaling based on the Stokes layer thickness is not relevant to the purely steady flow situation, the appropriate wavenumber is now *ah* for the Poiseuille flow. With these transformations noted, part of the neutral stability curve for steady Poiseuille flow, in a channel of half-width *h*=6, is shown in figure 1*c*. Here, the solid line indicates the results obtained from the Floquet technique, while the open circle marks a prediction from the DNS approach. The results are consistent and suggest a critical Reynolds number of *R*=481.02, which corresponds to *R*_{2}=5772.24 and a critical wavenumber of *a*=0.17 (or *ah*=1.02); these conditions are close to the well-known classical results for plane Poiseuille flow (Drazin & Reid 1995). Thus, at both extremes of the possible values of *Γ*, the DNS and Floquet techniques produce results that are reassuringly in very good agreement with previously available data.

Following the nomenclature of BB06, we shall describe a disturbance as an even mode if the stream function *Ψ* is an even function of *y*. Owing to the symmetries in equations (2.5*a*,*b*), it is clear that both even and odd perturbations are possible, and the Floquet technique of §2*a* was usually able to determine the value of *μ* for both modes at any values of the parameters considered. However, as required by the techniques of Davies & Carpenter (2001), both vorticity and wall-normal velocity perturbations are expanded as even Chebyshev polynomials (2.20), so only even disturbance modes could be examined via the DNS technique. Thus, for a given flow, the Floquet technique was used to determine which of the even or odd disturbance mode was the more unstable and then the DNS technique was used to confirm the results of the Floquet approach for the even mode. For the range of values of channel half-width considered here, the even mode was always more unstable than the odd disturbance mode.

For the purely oscillatory channel flow, the eigenvalues *μ* occur as complex conjugates *μ*_{r}±*iμ*_{i}, and correspond to left- and right-travelling waves with the same growth rate. This occurs as a result of the time-symmetry property
which is applicable when *Γ*=0. When *Γ*≠0, this symmetry is destroyed, natural upstream and downstream directions are defined, and left- and right-travelling disturbances no longer have the same growth rate. Figure 2*a* displays the growth rates *μ*_{r} (obtained using Floquet theory) of even mode disturbances as a function of the non-dimensional velocity ratio parameter *Γ*. Here, the Reynolds number has been fixed at 1100, the wavenumber *a*=0.25, the channel half-width *h*=6 and a logarithmic scale has been used for the *Γ*-axis. The solid and dashed lines denote two even Floquet modes emanating from a Stokes layer mode at *Γ*=0. For the range of *Γ* shown here, the effect of the mean flow is to decrease the growth rate of the upstream travelling wave and to increase the growth rate of the downstream wave.

Neutral stability curves in the (*a*,*R*)-plane for even disturbance modes to the basic flow corresponding to the values of the parameter *Γh*=0, 0.5, 2, 5, 10, 30, 50 and are shown in figure 2*b*. The channel half-width is again taken to be equal to 6, and lines again refer to results obtained using Floquet analysis while open circles (°) denote the corresponding DNS solutions. The two methods are again in excellent agreement. The upper most neutral stability curve (solid line) corresponds to the purely oscillatory flow and as the parameter *Γh* increases monotonically the critical conditions move to smaller wavenumbers across the plot, asymptoting towards the neutral curve for the plane Poiseuille flow, the dotted line located on the left-hand side of the plot with critical conditions of (*a*_{c},*R*_{c})=(0.17,481.02).

The critical Reynolds numbers *R*_{c} and wavenumbers *a*_{c}, for instability of even disturbance modes are illustrated in figure 3*a*,*b* as functions of *Γ*. Results are plotted using a logarithmic scale along the *Γ*-axis, for the choices *h*=4, 6, 8 and 10. Of course purely oscillatory flow is associated with the limit *Γ*→0 and steady Poiseuille flow with and for the velocity parameter range 0.1≤*Γ*≤10 the critical conditions differ significantly from these extremes. The most rapid variations in the critical Reynolds number occur around *Γ*≈1, while the critical wavenumbers decrease monotonically as *Γ* grows and asymptotes to the plane Poiseuille flow value *a*_{c}=1.02/*h*.

For the smaller channel half-widths (for example, *h*=4), the introduction of an increasing oscillatory component has a destabilizing effect on the purely steady shear flow; the critical Reynolds number is smallest near *Γ*=1 and then increases again and asymptotes to the critical value for the purely oscillatory channel. For larger channel half-widths, the effect of the modulation reverses and is stabilizing. For example, when *h*=10 the critical Reynolds number increases as the velocity ratio parameter *Γ* is lowered, taking a maximum at the sharp peak located near *Γ*=1. Thus, for sufficiently large channel half-widths, the oscillatory motion stabilizes the steady plane Poiseuille flow.

At this stage, it should be mentioned that most previous work discusses the stability characteristics of the current unsteady flows in terms of the frequency of the oscillatory component of the basic flow. Owing to the time scaling adopted here the dimensionless frequency of the oscillatory component is fixed, thereby forcing the channel half-width, *h* to play the role of proxy for the frequency imposed by the basic flow. With other physical quantities held fixed, small values of *h* correspond to a low frequency oscillatory component, while increasing the value of *h* is equivalent to increasing the frequency of the oscillatory component of the basic flow. In terms of frequency, the preceding paragraph can be summarized by reporting that for low frequency oscillations in the basic flow, increasing the magnitude of the oscillatory component destabilizes the basic flow, while with high frequencies present in the basic flow a stabilization is predicted for the range of amplitudes .

Further information on the linear stability of modulated unidirectional flows may be deduced by calculating the respective Reynolds numbers *R*_{1} and *R*_{2} relating to the plane Poiseuille flow and Stokes layer, respectively (recall equation (2.7)). Using the results in figure 3*a* for the combined critical Reynolds number *R*_{c}, the values for *R*_{1,c} and *R*_{2,c} may be determined separately and their dependence is depicted in figure 3*c*. Flows at (*R*_{1},*R*_{2}) combinations lying in the region to the left of the curves are stable and this plot clearly illustrates the influence that the channel half-width (frequency of oscillation) exerts on the stability characteristics of oscillatory plane Poiseuille flow. For relatively, small *h* (or low frequency of oscillation), *R*_{1,c} decreases and the modulation is seen to be destabilizing. In contrast, for larger *h* (higher frequency of oscillation), modulation appears to have the opposite effect. Indeed, when the channel half-width *h*=10, *R*_{1} takes a maximum of roughly 10^{4}, which is almost twice the steady flow value.

To focus more closely on the stability characteristics when a small oscillatory component is added to a steady flow (i.e. *Γ*≫1), figure 4 displays the relationship between *R*_{1} and *R*_{2} when *Γ*≥100 and *h*=20,40,60,80,100 and 120. The first four data lines (*h*≤80) indicate that modulation of the plane Poiseuille flow is stabilizing, while at the remaining larger channel half-widths (even higher frequencies in the basic flow), modulation promotes instability. This behaviour is consistent with the earlier investigations of Hall (1975) and von Kerczek (1982), who have, respectively, shown that destabilization occurs at very large *h* (or very high frequencies).

The neutral conditions *R*_{c} and *a*_{c} for odd disturbance modes are given as functions of *Γ* in figure 5. As *Γ* increases so the critical Reynolds number grows and the corresponding critical wavenumber diminishes. In all our simulations, odd mode disturbances were always found to occur at larger Reynolds numbers than their even counterparts implying that odd modal structures will inevitably be difficult to observe as they will be swamped by the faster growing even perturbations. We were unable to obtain satisfactory results for the odd mode disturbances beyond those given in figure 5. This was due to the necessity to retain a huge number of harmonics in the calculations, a difficulty previously noted (BB06) in connection with the Floquet approach. However, since it is the even modes that appear to be the more important in practice, this lack of reliable results for a wide parameter range is unlikely to be of great significance.

## 4. Axisymmetric disturbances in oscillatory pipe flow

The methods of the previous sections are now applied to fully developed oscillatory flow in circular pipes. Much of the formulation is preserved; of particular note is the fact that lengths are again scaled on the Stokes layer thickness with the consequence that *h* should now be interpreted as the dimensionless pipe radius rather than as the non-dimensional channel half-width. The governing equations are written with reference to the cylindrical polar coordinate system (*x*,*r*,*θ*), in which the *x*-direction points along the longitudinal axis of the pipe. This work is restricted to axisymmetric disturbances so that all basic flow and disturbance quantities are assumed to be independent of the azimuthal angle *θ*.

The non-dimensional basic flow generated when steady Hagen–Poiseuille flow, with maximum velocity *U*_{1}, is superimposed on a pipe that oscillates in the *x*-direction with a velocity , is given by
4.1a
where
4.1b
and
4.1c
with *J*_{0} denoting the Bessel function of order zero. As with the channel problem, all velocities are scaled on *U*_{1}+*U*_{2} and *γ*_{1}, *γ*_{2} and *Γ* are as defined in equation (2.2); when *Γ*=0 the basic flow is the zero time-mean, purely oscillatory flow in a circular pipe, with stability properties as determined by BB06, while as steady Hagen–Poiseuille flow is obtained.

The derivation of the governing equations for the Floquet eigenvalue problem closely follows that outlined in §2*a*. For axisymmetric disturbances, velocity perturbations are given by (*u*,*v*)=(*U*_{B},0)+*εr*^{−1}(*Ψ*_{r},− *Ψ*_{x}), where *Ψ*(*x*,*r*,*τ*) is the stream function, which may be decomposed as Here, *a* is again the wavenumber, *μ* the Floquet exponent and *ψ*(*r*,*τ*) is 2*π*-periodic in *τ*. The linearized governing equation for *ψ* is then
4.2
subject to the boundary conditions
4.3
where the operator and the subscript *r* denotes differentiation with respect to *r*. Note that the *U*_{B}′′(*y*,*τ*) coefficient appearing in equation (2.5a) has been replaced by in equation (4.2) and that the definitions of the differentiation matrices **V**,**L**,**M** and **P** are now
4.4a
4.4b
where
4.5
The change in geometry has meant that the matrix operator **P** no longer contains any curvature term and the matrix operator **M** is now more complicated owing to the presence of the Bessel functions in the unsteady component of the basic flow. The Fourier series (2.8) for *ψ*, with *y* replaced by *r*, was again truncated by setting *ψ*_{n}=0 for all |*n*|>*N* and the eigenvalue problem (2.12) again obtained. (Further details are given in BB06.) The eigenvalue routines within Matlab were again called upon to determine the eigenvalues and eigenvectors of the problem.

Figure 6*a* displays the critical Reynolds numbers for instability in pipe flow as a function of the velocity ratio parameter *Γ*, for several pipe radii. In every case, the critical Reynolds number grows as the size of the mean flow component increases away from zero. Such behaviour is to be expected, since it is thought that Hagen–Poiseuille flow, which is the flow in the limit, is strictly linearly stable for all *R*. It was found that for *Γ*∼*O*(1), the corresponding *R*_{c} was *O*(10^{3}) so that it became difficult to extend the calculations to larger *Γ* owing to the increasing number of harmonics needed in the Floquet representation of the perturbation. Thus, it proved impossible to confirm the existence, or otherwise, of a value *Γ*_{a} such that as *Γ*→*Γ*_{a}.

The analogue of figure 3*c* is shown here as figure 6*b*. The lack of a linear stability critical Reynolds number for Hagen–Poiseuille flow makes this figure very different from figure 3*c* and also from the equivalent graph from the experimental results of Lodahl *et al.* (1998).

## 5. Concluding remarks

In this work, the stability of a family of unidirectional flows in channels and pipes has been examined. The formulation of the stability problems presented has enabled the determination of the growth or decay of perturbations in a spectrum of motions from purely oscillatory flows with zero-mean through to flows with increasingly large steady flow components. The results of our calculations have confirmed previous studies of the effects of small components of oscillation superimposed on a steady flow and also extended the knowledge of these stability boundaries to flows comprising both mean and oscillatory components of significant sizes.

The Floquet analyses conducted in BB02 and BB06 concerning the stability of the classical flat Stokes layer and of the oscillatory flow, with zero mean, in a channel or circular pipe revealed the possibility of neutrally stable modes that are stationary relative to the basic flow. These stationary modes were induced by the complex conjugate symmetry of the eigenvalues *μ* and were found on several extremely narrow bands in the wavenumber space. Their presence gives rise to a neutral stability curve of fairly complicated structure, which is fully described in BB02 and BB06. Here, it is sufficient to note that the majority of the neutral curve in wavenumber–Reynolds number space was of a standard smooth shape but, around the critical conditions, this underlying structure was periodically punctuated by slender finger-like features that correspond to the stationary modes. However, even a small component of mean flow destroys the complex conjugate symmetry of the eigenvalues needed for the creation of stationary modes. Thus, in the presence of a mean flow, the only possibility is of isolated wavenumbers where there may be stationary neutral modes. We note that none of the results reported here suggested the existence of any stationary neutral modes when *Γ* was non-zero.

These observations may have implications for resolving the discrepancy between theoretical and experimental determinations of the critical conditions for the instability of the flat Stokes layer. As was outlined in §1, there has been a long history in the quest for an analytical prediction of the critical Reynolds number for instability. Even the first estimate provided by BB02 is considerably larger than the value suggested by most experiments. The BB02 calculation of *R*_{c}≈708 is at least two or three times the measured critical values, but there is very little consistency in the estimates of *R*_{c} provided by the various experiments. For purely oscillatory flows (*Γ*=0), Luo & Wu (2010) have clearly demonstrated that small amplitude waviness in the bounding surfaces can produce very large disturbance velocities in the flow. They attributed this to near-resonance of the wall-forced flow field with the weakly damped Floquet modes of the unforced problem. We speculate that the largest receptivity responses seen by Luo & Wu (2010) are associated with the presence of the stationary Floquet modes. Thus, these stationary modes could then be responsible for the sensitivity of the classical Stokes layer to imperfections in the experimental apparatus and hence cause instability to be triggered at Reynolds numbers much lower than predicted by linear theory. Such a suggestion could only be tested by conducting receptivity calculations in the spirit of that of Luo & Wu (2010) for oscillatory flows with both zero and non-zero mean components. For a flow consisting of non-zero mean and a weak oscillatory component, the case , Luo & Wu (2004) have already carried out the necessary analysis, including nonlinear effects. The case of *Γ*∼*O*(1) could be considered by combining the ideas in these two papers.

As a small non-zero mean flow component should be sufficient to disable the stationary modes, the results displayed in either figure 6*a* or *b* suggest that measurements of the critical Reynolds number for a pure Stokes layer might be less susceptible to experimental imperfections if a small mean flow were present. The critical Reynolds numbers reported (Lodahl *et al.* 1998) do not directly support this hypothesis, but the smallest steady flow Reynolds numbers used in these experiments are rather larger than those indicated in figure 6*b*. Thus, the possible role of the resonances found by Luo & Wu (2010) remains to be determined.

## Acknowledgements

This work was supported by the Australian Research Council by grant DP0880463. All referees are thanked for their suggestions which have led to a much improved paper.

- Received September 6, 2010.
- Accepted March 18, 2011.

- This journal is © 2011 The Royal Society