## Abstract

The use of numerical simulations to study the development of boundary-layer disturbances is illustrated for a number of different incompressible flow configurations. These include cases where the disturbances are generated by, or interact with, flow-control devices in the form of compliant panels, suction slots and microelectromechanical systems actuators. The velocity–vorticity system of governing equations used for the simulations is reviewed, along with the numerical discretization.

## 1. Introduction

Theoretical studies of disturbance evolution in incompressible boundary layers frequently require an assumption of laminar flow at an arbitrarily large Reynolds number in order for any progress to be made in obtaining approximate analytic solutions of an appropriately chosen set of governing equations. However, in corresponding physical experiments, it is usually only possible to maintain laminar flow over a very limited range of Reynolds numbers. It may or may not turn out to be the case that the largest Reynolds number that can be attained in a physical experiment is sufficiently large to ensure that the approximations made in the analytic theory can be applied with any confidence.

Similar difficulties in trying to match mathematical theory and experiment can also occur if computer experiments are considered instead of physical experiments. Resolution requirements, limited computational resources and the numerical stability properties of feasible algorithms all conspire to limit the range of Reynolds numbers that can be explored by means of computer simulation. Nevertheless, carefully conducted computer experiments do offer some new possibilities that allow them to play a complementary role to that of physical experiments. Particularly important are the facts that in a computer experiment, and unlike in any corresponding physical experiment, it is possible to: (i) study purely two-dimensional forms of disturbances, with no contamination at all arising from the effects of three-dimensionality, (ii) turn various effects of nonlinearity on or off as required, by performing a suitable linearization, (iii) turn on or off the effects of any spatial inhomogeneity in the basic boundary-layer flow, such as non-parallelism owing to boundary-layer growth. When the effects of either three-dimensionality, nonlinearity or spatial inhomogeneity are artificially excluded in a computer simulation, it can also be possible to study laminar disturbance development at much higher Reynolds numbers than could be attained in a physical experiment. Thus, computer simulations can sometimes be used to explore an idealized environment where the development of boundary-layer disturbances may be more readily compared with the behaviour predicted by mathematical analysis.

In this paper, I will describe a selection of numerical simulation results that have been obtained for the evolution of disturbances in a variety of boundary-layer flow configurations. Taken altogether, these simulation results are intended to illustrate how computer experiments can be used to provide an additional, idealized source of data that complements and extends the restricted data available from physical experiments.

## 2. Formulation of the disturbance evolution equations

The common feature of the numerical simulations that we will review is that they were all undertaken using a new velocity–vorticity formulation of the Navier–Stokes equations (Davies & Carpenter 2001). This formulation is particularly suitable for studying disturbance evolution in incompressible boundary layers. Thus, we begin by giving a brief account of the new formulation, followed by a summary of the associated numerical methods. The formulation will be stated for the simple case of a fixed Cartesian coordinate system, but this can readily be adapted to other cases such as cylindrical polar coordinates and a rotating frame of reference.

We will consider the governing equations for the flow perturbation variables only, by first decomposing the total velocity and vorticity fields into the form:(2.1)where the superscript *B* is used to distinguish the boundary-layer mean flow, which is assumed to be already known. The boundary layer is taken to be located entirely above a surface that is positioned at *z*=*η*(*x*,*y*,*t*). The components of the perturbation flow variables * u*=(

*u*,

*v*,

*w*),

*=(*

**ω***ω*

_{x},

*ω*

_{y},

*ω*

_{z}) are divided into two sets; namely, the primary variables {

*ω*

_{x},

*ω*

_{y},

*w*} and the secondary variables {

*ω*

_{z},

*u*,

*v*}. This terminology is adopted because the secondary variables can be explicitly defined in terms of the primary variables by performing an integration along the wall normal

*z*-direction, across the boundary layer, in the following manner:(2.2)(2.3)(2.4)

The governing equations for the disturbances can then be written in the form:(2.5)(2.6)(2.7)where *Re* is the Reynolds number and the convective quantity * N*=(

*N*

_{x},

*N*

_{y},

*N*

_{z}) is given by(2.8)

The two vorticity transport equations and the single Poisson equation stated immediately above specify a velocity–vorticity system of three governing equations for the three unknown primary variables {*ω*_{x},*ω*_{y},*w*}. The secondary variables, which are only required to evaluate the components of the quantity * N*, can all be eliminated because each can be expressed solely in terms of the primary variables. It may be shown that the stated velocity–vorticity system remains fully equivalent to the primitive variables formulation of the Navier–Stokes equations, provided that fairly weak conditions can be imposed on the primary variables for

*z*→∞. (Davies & Carpenter (2001) contains a more detailed discussion of these conditions.)

An obvious advantage of the velocity–vorticity system in equations (2.5)–(2.7) is that there is a reduction in the number of variables and governing equations, compared with other formulations of the Navier–Stokes equations for three-dimensional incompressible flow. Another major advantage is that the no-slip conditions can be imposed in an elegant and mathematically consistent fashion. If the perturbation fluid velocity is specified as on the bounding surface located at *z*=*η*(*x*, *y*, *t*), then the no-slip conditions can be written, using the definitions of the secondary variables, in the form(2.9)(2.10)

Each of these two integral conditions can be associated with the vorticity transport equation for the corresponding component of the perturbation vorticity. They thus provide a very natural means of imposing constraints on the vorticity evolution that arise directly from the no-slip conditions. Moreover, the no-slip conditions are employed once, and only once, to provide such constraints. The remaining condition to be applied at the surface *z*=*η*(*x*, *y*, *t*); namely, the normal velocity-matching condition(2.11)may be imposed on the solution of the Poisson equation (2.7) in a relatively straightforward manner.

The coupling between the no-slip conditions and the vorticity evolution that is facilitated by the integral conditions in (2.9) and (2.10) remains in force even when a linearized version of the velocity–vorticity system is employed. This is not the case for some other possible velocity–vorticity formulations where artificial vorticity boundary conditions are introduced (e.g. Meitz & Fasel 2000). Such artificial conditions do not guarantee that there is any well-posed connection between the vorticity and the no-slip conditions if nonlinear terms are excluded.

For each of the numerical simulations that are reported in the present paper, at least one simplification was made to the full velocity–vorticity system in equations (2.5)–(2.7). This accords with the remarks made in the introduction about the possible advantages of conducting computer experiments in which the effects of three-dimensionality, nonlinearity or spatial inhomogeneity are artificially eliminated. To study purely two-dimensional disturbances, the spanwise perturbation velocity *v* and spatial variations along the spanwise direction must both be neglected. This gives a system comprising a single vorticity transport equation and a Poisson equation that, taken together, determine the evolution of the primary variables {*ω*_{y}, *w*}. The streamwise velocity component *u*, which is the only non-vanishing secondary variable, remains explicitly defined by the integral relation (2.2), while the streamwise no-slip condition is ensured by the corresponding integral constraint (2.9). Some further details for the two-dimensional case can be found in Bowles *et al*. (2003). Studies that neglect nonlinear perturbation effects may be performed by simply dropping the product * ω*×

*in the calculation of the convective quantity*

**u***. The effects of spatial inhomogeneity in the mean flow can also be quite simply removed. For example, in the case of an essentially unidirectional flow configuration, such as Blasius flow over a flat plate, an approximation of the form*

**N**

**U**^{B}≏(

*U*(

*z*),0,0) may be used to eliminate terms that would otherwise be associated with the non-parallel nature of the mean flow.

## 3. Numerical methods

A full description of a numerical scheme for discretizing the velocity–vorticity system of governing equations is included in Davies & Carpenter (2001) for the particular case of disturbances developing in the three-dimensional von Kármán boundary layer over a rotating disc. The details of the discretization do not need to be altered in a very significant manner when a Cartesian coordinate system is used instead of the cylindrical polar coordinates that are appropriate for the geometry of the rotating disc. The main features of the discretization can be summarized as follows. (i) The streamwise (or radial) variation is discretized using finite differences, which are typically of at least fourth-order accuracy. (ii) Fourier expansions are deployed for the variation in the spanwise (or azimuthal) direction. (iii) Chebyshev expansions are used for the discretization in the wall normal direction across the boundary layer. An algebraic coordinate transformation of the form *Z*=*l*/(*l*+*z*), where *l* is a suitable mapping constant, is deployed to map the semi-infinite physical domain onto a finite computational domain. (iv) To ensure numerical stability without any overly restrictive limit on the size of the time-step, the temporal discretization of the vorticity transport equations is taken to be implicit for viscous terms that involve second derivatives in the wall normal direction, but all other terms (including Coriolis terms, if appropriate) may be treated explicitly. (v) The integral constraints on the vorticity, which are used to impose the no-slip conditions, can also be treated in an explicit fashion. This makes it possible for the solution of the discretized vorticity transport equation and the discretized Poisson equation to be fully decoupled within the time-stepping procedure. (vi) The wall normal discretization is formulated so as to involve only pentadiagonal matrix operations. This facilitates the direct solution of the discretized transport equations using a highly efficient Thomas algorithm. The Poisson equation can be solved by combining the Thomas algorithm with either an iterative streamwise line-marching procedure or a direct method that involves a fast sine-transform along the streamwise direction. (vii) A pseudo-spectral transform technique is used to compute the nonlinear and other product terms in the transport equations that appear via the convective quantity * N*. It should be noted that in all of the numerical simulations that we describe in the present paper, nonlinear perturbation effects are only retained for two-dimensional forms of boundary-layer disturbance. Thus, it was only necessary to implement the pseudo-spectral calculation of product terms for the Chebyshev expansions that are used to discretize the wall normal variation.

## 4. Numerical simulation results

### (a) Disturbances in boundary layers over compliant surfaces

The growth of two-dimensional linear Tollmien–Schlichting instability waves provides a possible first step in the process of laminar–turbulent transition for many boundary-layers flows over solid bodies. Coating the surface of the solid with a compliant material can lead to a stabilization of Tollmien–Schlichting waves (Carpenter *et al*. 2000, 2001; Davies 2003). Thus, compliant surfaces may be employed to delay the transition to turbulence in boundary layers, with a consequent reduction in skin friction drag.

We have conducted a range of numerical experiments to study the behaviour of Tollmien–Schlichting waves and other forms of disturbance as they propagate in boundary layers over compliant surfaces. Figure 1 shows some typical results from a two-dimensional linearized simulation that was conducted to study the spatial development of a Tollmien–Schlichting wave in a Blasius boundary layer. A parallel flow approximation has also been applied, but mainly to simplify the interpretation of the simulation results. The wave is very strongly stabilized as it propagates over a compliant panel that is located between the streamwise positions *x*=80 and 320. The wall is taken to be rigid at all other streamwise locations. It can be seen that the vorticity at the wall drops in a dramatic fashion as the Tollmien–Schlichting wave passes over the rigid compliant wall join at the upstream end of the compliant panel. In fact, the Tollmien–Schlichting wave is so strongly stabilized that, in effect, it is removed almost as soon as it passes over the join. What remains is a relatively long wavelength, surface-based wave that is generated by the localized wall motion forced by the incoming Tollmien–Schlichting wave. This surface wave then propagates downstream along the compliant panel until it reaches the compliant rigid wall join. At the downstream end of the compliant panel, a second surface-based wave is generated, but this propagates along the compliant panel in an upstream direction. The wall motion in the vicinity of the downstream end of the compliant panel also acts as a generator to create a new Tollmien–Schlichting wave that then propagates away in the downstream rigid walled section. It should be noted that despite the somewhat complex dynamics that are found over the compliant panel, there is still a strong reduction in the overall disturbance growth. In the absence of the compliant panel, the incoming Tollmien–Schlichting wave would have grown by a factor of approximately 30 as it propagated over the full range of streamwise positions that are shown in figure 1.

Numerical simulations have also been undertaken to examine the effects of surface compliance on the linearized development of highly unstable forms of disturbance in the three-dimensional von Kármán boundary layer over a rotating disc (Davies & Carpenter 2003). In addition to confirming that wall compliance can be stabilizing, the simulations results served to highlight the strongly non-parabolic nature of the rotating disc boundary layer. It was found that a compliant annulus inserted into an otherwise rigid disc was capable of stabilizing disturbances at locations that were radially inboard from the compliant part of the disc surface. Simulations conducted for an entirely rigid rotating disc also demonstrated that, contrary to expectation, the spatial inhomogeneity of the von Kármán flow was very significant in determining the long-term evolution of disturbances. By comparing the results of simulations conducted with and without a parallel flow type of approximation, it was found that the mean flow spatial inhomogeneity prevented the absolute instability from giving rise to a temporally growing global mode of disturbance.

Three-dimensional linearized numerical simulations have recently been conducted to study the effects of surface compliance on the spatial development of transient streaks that are associated with so-called Klebanoff modes of disturbance in transitional boundary layers (Ali 2003). These modes involve a spanwise modulation of the streamwise velocity. For the Blasius boundary layer, the numerical simulation results demonstrated that the growth of streaks can be diminished by a compliant surface, but also that these more weakly growing streaks may persist over a wider range of spanwise wavenumbers than is the case when the wall is rigid.

### (b) Disturbances generated by flow-control devices

Wall suction is perhaps one of the oldest known forms of boundary-layer control. Nevertheless, the disturbance field generated by a localized suction device still offers some intriguing behaviour for theoretical study (Ovenden 2001). We have conducted a series of numerical experiments for both Blasius boundary-layer flow and for an idealized flow with a uniform shear to examine how the disturbance field generated by a two-dimensional suction slot is altered as an increasing level of nonlinearity is introduced. More specifically, we studied what happens when the width of the suction slot is decreased while the streamwise distribution of the suction is simply rescaled to ensure that the total flux of the fluid drawn into the slot remains constant. Figure 2 illustrates the typical form of the jump in the pressure field that is found when the width of the suction slot is larger than or comparable to the boundary-layer thickness. (Within the context of the new velocity–vorticity formulation, the pressure can be defined as another secondary variable that is easily determined from the primary variables. It can be computed using an integral across the boundary layer that is obtained, in effect, by simply integrating the wall normal momentum equation. Further details may be found in Davies & Carpenter 2001.) The numerical simulation results for such relatively large suction slots suggest that various features of the disturbance field scale with the width of the slot in a relatively straightforward manner. For example, it appears that for a fixed (but not too large) value of total mass flux, the height of the pressure jump is inversely proportional to the quarter power of its width, where this width, in turn, is directly proportional to the width of the suction slot.

Suction slots that are much smaller than the boundary-layer thickness were found to display a rather different kind of behaviour. If the mass flux is set to be sufficiently large, then it becomes possible to find a region of reversed flow over part of the suction slot itself. The character of the reversed flow observed in the numerical simulations appears to be in broad agreement with the predictions that were made by Ovenden (2001) using an asymptotic analysis for indefinitely large Reynolds numbers.

The linearized version of the full velocity–vorticity system in equations (2.5)–(2.7) has also been used to conduct computer experiments concerned with the spatial and temporal development of small amplitude three-dimensional boundary-layer disturbances introduced by microelectromechanical systems flow control devices (Lockerby 2001; Carpenter *et al*. 2002). The numerical simulation results suggest that, for practicable devices, it is often very important to take full account of any interactive coupling that can occur between the device and the boundary-layer disturbances that it generates. For example, in Blasius flow, a certain sort of jet-type control device was found to be susceptible to a resonance that could lead to the generation of strong boundary-layer disturbances from low levels of background noise, even when the device was left in a supposedly inactive state.

### (c) Nonlinear two-dimensional disturbance development

It is well known that in physical experiments, the later stages of boundary-layer transition to turbulence are always highly three-dimensional. Nevertheless, there is still much of interest in theoretical studies conducted for the idealized case of purely two-dimensional forms of disturbance. Such studies may identify relatively simple two-dimensional versions of physical mechanisms that play a significant role in the transition process. However, the inherent difficulty to comparing the results of a two-dimensional theoretical analysis with the fully three-dimensional results obtained from physical experiments is often exacerbated by the fact that the theoretical analysis usually proceeds using an assumption of an indefinitely large Reynolds number. Thus, a useful preliminary step before attempting any direct comparisons between theory and physical experiments may be to conduct two-dimensional numerical simulations for various values of the Reynolds number. The aim is to ascertain whether or not the theoretically predicted behaviour retains any measure of agreement with that of simulated two-dimensional disturbances as the Reynolds number is lowered to values that are pertinent to physical experiments.

For example, we have compared theoretical predictions and numerical simulation results for the weakly nonlinear development of two-dimensional Tollmien–Schlichting waves for Blasius flow (Houten *et al*. 2000; Houten 2004). This work has been followed up by far more demanding computer experiments that simulated the strongly nonlinear development of such waves (Bowles *et al*. 2003). The simulations were conducted for a range of Reynolds numbers, based on the boundary-layer thickness, between *Re*=1000 and 100 000. Although there was no difficulty at all in using a non-parallel form for the base flow, most of the simulations were carried out using a parallel flow approximation. This was entirely for the purpose of simplicity. Complications in making comparisons between the relatively short-scaled structures seen in the simulations and theoretical predictions would otherwise have arisen because of the streamwise dependence of the stability properties for an inhomogeneous base flow.

Even with the imposed artificial restriction to two dimensions, the behaviour uncovered by the simulations for the later stages of the disturbance development was found to be extremely complex. However, it was still possible to use the predictions of high Reynolds number theory (Li *et al*. 1998) to correctly identify and describe many of the features that were found in the evolution of the simulated disturbance. In particular, the mechanisms highlighted by the two-dimensional theory for the so-called spiking stages could be clearly traced in the data supplied by the numerical simulations.

Recently, we have extended our two-dimensional simulations of strongly nonlinear Tollmien–Schlichting waves to include cases where the mean flow is prescribed as a Falkner–Skan boundary-layer flow, subjected to a parallel flow approximation. This was to assess the impact of streamwise pressure gradients on the disturbance development. Figure 3 is a typical early time snapshot of the perturbation vorticity for a case where there is an adverse pressure gradient. The elongated structure that extends almost vertically upwards from the wall corresponds to a region of negative vorticity that is being pulled away from the wall as part of a vortex-formation process. This behaviour is quite similar to developments that were identified in the simulations for Blasius flow.

## 5. Concluding comments

We have outlined various merits of using computer experiments to complement physical experiments in studying the development of disturbances in boundary layers. In particular, we have attempted to illustrate the benefits of conducting simulations for idealized cases where the effects of three-dimensionality, nonlinearity or mean flow inhomogeneity are systematically neglected. All of the simulations that we have reported were carried out using the new velocity–vorticity system of governing equations that is reviewed in §2.

Restrictions of space have meant that our descriptions of numerical simulation results have been kept as short as possible. We have also refrained from mentioning many essential details of the theoretical background and motivation for the simulations. Because these limitations might create a misleading impression, we conclude by reiterating our view that the most important role of computer experiments is to provoke, support or undermine specific explanatory theories. As with physical experiments, a computer experiment conducted for boundary-layer disturbance evolution does not, by itself, provide any substitute for a well-founded theoretical explanation.

## Footnotes

One contribution of 19 to a Theme ‘New developments and applications in rapid fluid flows’.

- © 2005 The Royal Society