## Abstract

High Reynolds number mathematical models for convected vortex impinging against a local hump and sound scattering into Tollmien–Schlichting eigenmodes are introduced to simulate basic mechanisms of the disturbance excitation typical of turbomachinery environments. The streamwise dimension of the transitional flow on the suction side of a blade is evaluated on the assumption that the transition length is of equal order with the extent of viscous/inviscid interaction controlling the boundary-layer response. The triple-deck theory gives a simple power law correlation to express a value of the Reynolds number based on the transition length in terms of the Reynolds number calculated with the blade cord. Precisely the same correlation stems from processing experimental data for both smooth and rough surfaces. The computation shows the explosive development of highly modulated wave packets and their rapid breakdown brought about by erratic short-scaled wiggles riding on the primary long-scaled oscillation cycles. The filling-up of distant parts of the wavenumber spectrum is at the heart of the signal distortion. This process heralds the start of deep transition terminating in fully developed turbulent flow well before reaching the upper stability branch. With the time-harmonic excitation broadly used in experiments, transition requires a much longer distance to complete. An agreement between theoretical predictions based on the assumption of indefinitely large Reynolds numbers and experimental findings from wind-tunnel observations at finite Reynolds numbers is encouraging.

## 1. Introduction

As known, gas-turbine-engine blades typically operate at Reynolds numbers *Re* in the range of 450 000<*Re*<850 000. For these low Reynolds numbers, laminar flow occupies a significant portion of the suction side of the blade surface (Abu-Ghannam & Shaw 1980; Gostelow & Ramachandran 1983; Feiereisen & Acharya 1986; Gostelow & Blunden 1989; Kuan & Wang 1990). As a consequence, the aerodynamic performance of turbomachinery blades strongly depends on the transition modes developing in the suction-side flow. The transition process recorded in various wind-tunnel tests shows instabilities in the Tollmien–Schlichting spectral interval of frequencies and wavenumbers. There exist similarities in the transition process between wake-induced turbulent patches on the blades and artificially triggered turbulent spots in the flat plate boundary layer. These similarities are exploited in experiments by Gostelow *et al*. (1999) and Gostelow & Thomas (2005) to simulate the transition phenomena in turbomachinery blading on a larger scale and in more detail. With the time-harmonic excitation broadly in use in wind-tunnel observations, transition requires a long distance to terminate in fully developed turbulent flow.

The earlier experimental data discussed above have been processed and summarized by Mayle (1991) in an attempt to establish a correlation between the length of transition and a reference Reynolds number. The correlation proved to be extremely simple. It shows that the Reynolds number *Re*_{lT} calculated with the transition length is a linear function of , where designates the momentum thickness Reynolds number at the end of transition. This empirical dependence calls for a firm theoretical footing to be derived from hydrodynamic stability theory. The triple-deck scheme based on the assumption that a typical Reynolds number takes on sufficiently large values offers a workable means to attack the problem and underpin the findings from wind-tunnel measurements. However, the results set forth below go far beyond their immediate motivation. They demonstrate that the powerful asymptotic approach that holds for laminar boundary layers, as it was initially designed and discussed in Stewartson (1974), Smith (1982) and Kluwick (1998), can be extended to transitional regimes. The applicability of theoretical predictions is not limited to turbomachinery flows, similar inferences equally pertain to problems of external aerodynamics, aircraft wing sections being a direct analogue for turbine/compressor blades.

External perturbing sources operating in the pulse mode generate strongly modulated wave packets in the Tollmien–Schlichting spectral range which cause transition to complete long before reaching the upper stability branch. The continuous filling-up of the wavenumber spectrum by the overlapping wings of neighbouring peaks gives rise to erratic short-scaled wiggles, corrupting the primary long-scaled oscillation cycles. In this regard, transition induced by the time-harmonic agency shows marked distinctions from the explosive transitioning to turbulence experienced by the wave packets. With references to earlier contributions, the spiking stage in late transition was thoroughly analysed by Bowles *et al*. (2003). New experimental data corroborate the general concepts even for rough surfaces (Pinson & Wang 2000).

## 2. Theory

Reasoning from the experimental findings of Gostelow *et al*. (1999) and Gostelow & Thomas (2005) that instabilities develop in the Tollmien–Schlichting spectral range, we begin with the triple-deck scheme set forth in Stewartson (1974). This approach leans upon asymptotic expansions in a small parameter,(2.1)with *Re* being the characteristic Reynolds number. When considering experimental data, *Re* is as a rule evaluated with the distance *L*^{*} counted along a section from its nose. The empirical correlation by Mayle (1991) gives the transition length Reynolds number, , for the blade suction side in terms of the momentum thickness Reynolds number at the end of transition, , by(2.2)Thus, we have three Reynolds numbers: *Re* enters the asymptotic development, whereas processing of experimental findings is based on and . Some comments are due at this point. First, for a flat plate, the momentum thickness Reynolds number, , can be expressed through the displacement thickness Reynolds number, , as(2.3)An analogous relationship holds for any airfoil section, with the numerical factor on the right-hand side depending on its shape. Second, within the framework of the high Reynolds number approach, the transition length is assumed to be much shorter than the characteristic distance . Asymptotically, as *ϵ* tends to zero, (2.2) becomes(2.4)with the coefficient determined by using (2.3). Insofar asfor the Blasius boundary layer on a flat plate, we finally have an expression(2.5)defining the transition length through the characteristic distance . An analogous dependence is typical of an arbitrary blade or wing section. Asymptotically, (2.4) and (2.5) are tantamount to (2.2) and can be exploited as empirical correlations; however, the results from wind-tunnel measurements collapse without appreciable scatter onto a straight line drawn by Mayle (1991) from (2.2). For this reason, the momentum thickness Reynolds number seems to be preferable when processing experimental data.

Let and be the viscosity and temperature, respectively. With the subscript labelling the conditions at the upper edge of the boundary layer, the Chapman–Rubesin linear law frequently used in theory readsThe triple-deck time , where is the velocity of the oncoming stream, and spatial coordinates , are scaled in terms of Reynolds number through the small parameter *ϵ* introduced in (2.1) and include the ratio of the wall temperature to the external temperature , as well as the Mach number and the non-dimensional skin friction , in the incompressible Blasius flow. For an arbitrary boundary layer, the skin friction can be expressed through *λ* by means of . The above scaling defines a characteristic Strouhal number for free interaction. According to Stewartson (1974), Smith (1982) and Kluwick (1998), the independent variables in the viscous near-wall sublayer of the subsonic boundary layer are(2.6a)(2.6b)(2.6c)Related to these normalized coordinates, the velocity field,(2.7a)(2.7b)and the pressure variations,(2.8)obey the system of Prandtl equations,(2.9a)(2.9b)for an incompressible fluid, independent of the Mach number. The self-induced pressure comes from viscous/inviscid interaction rather than being given beforehand by a potential flow solution. In the absence of perturbing agencies, the pressure rises and drops are connected with the instantaneous displacement thickness −*A*(*t*, *x*) by the interaction law,(2.10)The limit condition(2.11)at the upper reaches of the viscous near-wall sublayer and the no-slip condition(2.12)at the body surface make the mathematical formulation complete.

In order to extend the above consideration, we may drop the Chapman–Rubesin viscosity law and introduce the independent variables through(2.13a)(2.13b)(2.13c)instead of (2.6*a*)–(2.6*c*). The non-dimensional velocity field and the pressure variations,(2.14a)(2.14b)(2.14c)come in place of (2.7*a*), (2.7*b*) and (2.8), respectively. Here and the ratio depend on a particular boundary layer in question. The scaling in terms of *ϵ* as well as -dependence remain intact; however, the normalization that accounts for the temperature factor varies with the type of the boundary layer under consideration. Nonlinear stability of boundary layers for disturbances of various sizes was discussed by Smith (1979*a*).

The length of the nonlinear interaction process governed by the system of Prandtl equations can be estimated as with defined in equations (2.6*b*) and (2.13*b*). Hence, it follows thatin view of (2.1). The Reynolds number *Re*_{li} calculated with this characteristic length is expressed by virtue of(2.15)in terms of the Reynolds number. Thus, the transition length in (2.4) exhibits the same scaling as the interaction length in (2.15). In fact, can be identified with because evidently . The result obtained goes far beyond the scope of the triple-deck theory which has been originally invented as applied to laminar boundary layers (Stewartson 1974; Smith 1982; Kluwick 1998). The identity of the two lengths implies that the asymptotic approach extends to cover the regime of laminar/turbulent transition (except when the intermittency becomes close to 1).

Simple physical reasoning provides arguments substantiating this inference. The triple-deck scaling of time and spatial coordinates is typical of viscous/inviscid interaction. On the other hand, the same type of interaction is known to control the boundary-layer eigenmodes in the Tollmien–Schlichting spectral range centred about the lower branch of the neutral stability curve (Smith 1979*b*; Zhuk & Ryzhov 1980). Thus, the scaling introduced above specifies the period and wavelength of oscillation cycles within the Tollmien–Schlichting wave packets. The zone of transitional flow was recorded in numerous wind-tunnel tests to be filled with vortex spots (wave packets) that amplify exponentially fast and turn into turbulent spots coalescing with each other, and thereby sharply increasing intermittency at later stages of transition. Between these spots the flow field remains laminar. So long as there is enough space between neighbouring spots they develop as if they were travelling in laminar surroundings. This is the reason why the transition and free interaction lengths have equal scaling in terms of the Reynolds number. Vigorous growth of the pulsation amplitudes within the spots does not change the scaling in (2.15), but can trigger erratic short wavelength wiggles to destroy the original well-organized vortex structures.

## 3. Numerical evidence

The last statement needs to be reinforced computationally. To this end, let us consider two mechanisms of the vortex-spot production.

### (a) Wave-packet excitation by a convected vortex

The first mechanism comes about by virtue of the potential vortex/surface roughness interaction. Strong similarities in the transition process between wake-induced turbulent patches on turbomachinery blades and artificially triggered turbulent spots in the flat plate boundary layer underlie this theoretical model. As mentioned, these similarities were exploited by Gostelow *et al*. (1999) and Gostelow & Thomas (2005) in their wind-tunnel tests to simulate the transition phenomena. Certainly, natural transition on turbomachinery blades is three-dimensional. However, in Gostelow & Thomas (2005), wakes were introduced into the flow by a tapered rod, of average radius 4.5 mm, cantilevered from a disc having axis 0.5 m upstream of the leading edge of the thin plate. The rod was mounted at a radius of 270 mm on the disc and was rotated at an angular velocity of 60 rpm, resulting in the introduction of two dissimilar wakes per second. Centreline phase-averaged velocity traces were obtained starting from the aforementioned similarities between the wake-induced turbulent patches and artificially excited turbulent spots. In fact, the covected-vortex model relates to conditions typical of the experimental set-ups rather than to the general case for which no wind-tunnel data are available.

Figure 1 gives an idea of the process and makes clear the triple-deck pattern of the velocity field with pertinent scalings involved. The full treatment of the problem set forth in Ryzhov & Timofeev (1995) will not be duplicated here. We point out briefly the most important features associated with the convected vortex. Asymptotically, as , the flow in the near-wall sublayer obeys the system of Prandtl equations (2.9*a*) and (2.9*b*). However, an additional term is included in the expression (2.10) for the excess pressure to take into account the vortex-induced contribution. The limit condition (2.11) with the instantaneous displacement thickness entering the right-hand side also contains a correction for the vortex-induced pressure. This inhomogeneous limit condition at the upper edge of the near-wall sublayer should be supplemented with the inhomogeneous limit condition at the entry, as , to the sublayer. The derivation of the inflow condition rests upon a detailed analysis of the vortex motion in a region upstream of the local roughness. Matching with the Stokes sublayer located close to the wall in this region shows that the instantaneous displacement thickness tends to zero as , while the vortex-induced velocity,(3.1)provides the time-dependent contribution. In consequence,as . The no-slip condition (2.12) applies at the body surface with a small hump or dent,(3.2)centred around the origin. A different mathematical model of the Tollmien–Schlichting wave generation by free-stream turbulence is considered in Duck *et al*. (1996).

Following Burggraf & Duck (1982) and Duck (1987), the boundary-value problem posed is recast by introducing the shear stress as a new desired function. Computing the Fourier spectral distributionand then inverting it to the physical spaceplay a dominant role in an iteration procedure devised for an equationthat arises upon differentiating (2.9*b*) with respect to *y*. Our main concern is with the onset of erratic wiggles in accompanied with filling up the spectrum with *t* increasing. A code with FFT (fast Fourier transform) algorithm has been tested in Ryzhov & Savenkov (1989) and Ryzhov & Timofeev (1995) when analysing the short-scaled large-amplitude distortions appearing in a few central soliton-like cycles of the wave packets which can be identified with the vortex spots in laminar/turbulent transition. The geometry(3.3)assumed in the computation specifies a hump mounted on an otherwise flat plate; accordingly, the height of the obstacle in (3.2) being *a*=1.

For better visualization of oscillation properties, a quantity(3.4)is meant by in what follows. This quantity equals for the steady flow set in prior to the vortex motion that induces in line with (3.1), an additional component of velocity. The hump fixed by (3.3) introduces moderate disturbances in the velocity field of the oncoming steady stream without separation accompanied by a bubble with reversed flow (for separation to occur, the height of the hump should be larger). The resulting boundary-value problem provides the simplest model of gas-turbine-engine environments, where vortical disturbances interact with small imperfections of the turbine-blade surface (especially covered by the thermal barrier coating). Recent wind-tunnel tests by Pinson & Wang (2000) lend credence to the applicability of the transition length correlation to rough surfaces.

Figure 2*a–c* shows the origin and rapid development of a wave packet induced by a convected vortex interacting with the hump shaped by (3.2) and (3.3). A value determines a fairly strong intensity of the vortex. At *t*=0, the excess shear stress is associated predominantly with the steady perturbations from the hump, the formation of unsteady oscillations downstream is barely discernible. The explosive character of the highly modulated wave packet propagation becomes clear at *t*=1.5 when a narrow soliton-like cycle appears in the central part of the signal. Figure 2*c* illustrates the onset at *t*=2.5 of short-scaled erratic wiggles that distort a portion of two central cycles of the wave packet with positive . The soliton-like portions of the same cycles with negative remain almost intact at this spiky stage typical of transition.

Figure 3*a–c* displays the spectral content of disturbances shown in figure 2*a–c*. The first of them presents the wavenumber spectrum of almost steady perturbations triggered by the hump, with the peak in located approximately at . The second plot exhibits the spectrum featuring the wave packet at *t*=1.5. The local maximum at becomes much smaller compared to the local peak at about that corresponds to the wavenumber of the fastest growing linear eigenmode. The second local maximum at forks into two small mounds divided by a shallow local minimum in . The third local smoothed-out maximum arises at . One more local maximum is scarcely visible close to . Thus, the filling-up of the wavenumber spectrum basically follows a weakly nonlinear scenario in spite of the fact that the signal in figure 2*b* is essentially nonlinear. To the contrary, the filling-up of the distant parts of the spectrum in figure 3*c* is at the heart of the short-scaled wiggles corrupting two central cycles of the signal in figure 2*c*. Here, the first two peaks in become almost equal, whereas the other local maxima attain the magnitudes comparable to these peaks. The signal as a whole nears to transitioning to a turbulent spot. This process embracing all the vortex spots in the boundary layer eventually brings up the intermittency close to 1.

### (b) Acoustic wave/surface roughness interaction

The second mechanism of the vortex-spot production comes about by virtue of sound scattering into Tollmien–Schlichting eigenmodes by small streamwise variations in surface geometry. This mechanism is even more important as applied to aircraft wings rather than turbomachinery blades. Time-harmonic oscillations were studied by Ruban (1984), Goldstein (1985) and Goldstein & Hultgren (1989). The closest resemblance to the wave-packet excitation by a convected vortex shows the propagation of a sound pulse impinging against a local obstacle located on an otherwise smooth surface (Ryzhov & Timofeev 1995). However, there exists a difference between the two problems, because a fluid in the potential flow domain cannot be regarded anymore to be incompressible even if velocities are small. The Mach number and temperature ratio dependencies should be incorporated into a set of similarity parameters in a manner specified in (2.6*a*)–(2.6*c*), (2.7*a*), (2.7*b*) and (2.8), or (2.13*a*)–(2.13*c*) and (2.14*a*)–(2.14*c*). The system of incompressible Prandtl equations (2.9*a*) and (2.9*b*) still controls the flow field in the viscous near-wall sublayer, where the excess pressure related to the instantaneous displacement thickness obeys (2.10) with an additional term accounting for sound scattering at subsonic velocities. What is more, the form *f*(*t*) of the sound pulse can be defined in such a way as to render it identical to a function *p*_{v}(*t*) that was used above to prescribe the vortex-induced pressure. This analogy enables the wave-packet excitation by a convected vortex to be experimentally simulated by means of acoustic generators.

Instead of applying this analogy, we put in the computationwith . Thus, a function in the case under consideration is positively defined for any positive *t*. The sound-pulse-induced velocity comes from (3.1). Since with the values of parameters indicated, the wave packet emitted in the course of interaction is fairly weak and it needs to amplify before entering an essentially nonlinear stage. As obvious from figure 4*a*, the amplitude of pulsation cycles at *t*=2.5 still remains rather small, hardly amounting up to one-third of the size of steady perturbation provoked by the hump in the absence of the pulse. Note that in the immediate vicinity of the origin, the excess shear stress calculated from (3.4) only slightly deviates from the one obtained at *t*=0 in figure 2*a* for the vortex/hump interaction. The formation of travelling disturbances does not appreciably perturb the steady flow upstream. At *t*=4.0, the wave packet in figure 4*b* separates from the domain of steady flow around the hump and acquires tangibly nonlinear properties clearly embodied in a narrow soliton-like cycle at the centre. The further growth of pulsations leads to the occurrence of short-scaled wiggles that are randomly distributed over the upper portions of several large-amplitude cycles. As evident from figure 4*c* computed at *t*=5.5, the lower portions of the cycles with negative preserve their well-defined soliton-like shapes.

Figure 5*a–c* demonstrates the spectra of the sound-generated vortex spot. At *t*=2.5, there is a strong interaction of steady modes shown in figure 3*a* with the most amplified Tollmien–Schlichting eigenmode having . The mode interaction results in a toothed form of the wavenumber spectrum in figure 5*a*. The large-amplitude peak in at approximately dominates the spectrum in figure 5*b*. Here, the steady-state part of the spectrum gets much weaker compared to that generated by the downstream moving signal despite the fact that the steady and travelling disturbances in figure 4*b* are of nearly equal sizes. The filling-up of the distant parts of the spectrum gives rise to short-scaled wiggles. The main peak in figure 5*c* still relates to the fastest growing linear eigenmode though the other local maxima for modes of multiplicity *j* stand out clearly. The signal is close to transitioning to a turbulent spot. It is worth emphasizing that even though the amplification process depends on the precise form of the function , this dependence, as the computation suggests, is fairly weak for the hump in (3.3), which may be assigned to obstacles of moderate size.

## 4. Correlation

The formation of spikes in the computation presented above is the cornerstone underlying the assumption about the correlation of the viscous/inviscid interaction region predicted theoretically with the extent of transition flow observed experimentally in wind-tunnel tests with compressor/turbine blades. The computation points to a very rapid breakdown of central oscillation cycles in the wave packets excited by a convected vortex or sound pulse. The process is accompanied by unsteady flow separation within rapidly moving local bubbles which are seen in figures 2*a–c* and 4*a–c* as narrow zones, where total vorticity .

A strong support for the break-up of the self-focusing spikes comes from Smith (1988), who rigorously proved that a singular behaviour can occur at finite times in any unsteady interacting boundary layer. In fact, it is just this singularity that governs the spike development in central cycles of the wave packets. The shortening streamwise length scales associated with such a singularity cause the normal-to-wall velocity to increase to larger values than are admissible in (2.7*b*) and (2.14*b*). As a result, normal pressure gradients come into operation, rendering the interaction law (2.10) invalid (Li *et al*. 1998). This marks the start of late or deep transition terminating in fully developed turbulence. A thorough investigation of deep transition and unsteady separation, including a comparison between asymptotics and numerics, can be found in Bowles *et al*. (2003). The narrow eruptive spires are clearly observable in recent careful experiments and direct numerical simulations by Peacock *et al*. (2005).

All these analytical, numerical and experimental findings lend credence to our view that viscous/inviscid interaction controls transition up to the position where intermittency gets close to 1 owing to the spreading and ensuing amalgamation of spiking wave packets. Asymptotically, (2.4) yields an expression for the transition lengths on the blade suction side in terms of . However, the results from wind-tunnel measurements correlate better if the momentum thickness Reynolds number is used in place of . Accordingly, (2.2) supercedes (2.4), the only difference between these two equations is in numerical factors on their right-hand sides. Computationally, the values of the numerical factors vary depending on the maximum intensity of the convected vortex or sound pulse. With small , a value of the coefficient in the relationship giving as a function of does not appreciably deviate from that indicated in (2.4). The wave packet in figure 2 experiences an earlier breakdown with a fairly large used as a convected vortex intensity, whereas short-scaled wiggles destroy the positive oscillation cycles of the disturbance in figure 4 in compliance with (2.4), since specifies this case of sound scattering.

As already noted, the correlation (2.2) was first suggested by Mayle (1991) as a consequence of processing experimental data from different sources (Abu-Ghannam & Shaw 1980; Gostelow & Ramachandran 1983; Feiereisen & Acharya 1986; Gostelow & Blunden 1989; Kuan & Wang 1990). His correlation was extended by Pinson & Wang (2000) to rough surfaces characteristic of gas turbine blades in the course of their long operation life. A remarkable fact is that the asymptotic theory with the concept of the interacting boundary layer at the basis predicts precisely the same dependence. A solid line in figure 6 may be equally treated as coming from all experimental findings or deriving from asymptotic analysis. Experimental points collapse without appreciable scatter onto the solid line, independent of the type of blade surface, be it smooth or rough. Both correlations, found theoretically as well as experimentally, are identical in pointing to the same 5/4-power law in spite of the fact that the Reynolds number is assumed to tend to infinity in the triple-deck approach, whereas it takes finite values in wind-tunnel tests.

## 5. Discussion

It is shown that predictions from the asymptotic triple-deck scheme are in excellent agreement with the experimental correlation for the length of the transition zone in turbine/compressor blade flows, thereby creating a firm theoretical footing for the preliminary design strategy. The theoretical dependence of the transition length on the reference Reynolds number agrees with an empirical behaviour derived by Mayle (1991) from smooth-surface measurements reported in different studies. New experimental data by Pinson & Wang (2000) for rough surfaces are in line with these earlier findings. An assumption underlying theoretical inferences is that the length of the nonlinear viscous/inviscid interaction and the transition length are of equal order in magnitude. Then the scaling inherent in the high Reynolds number asymptotic approach automatically provides an estimate for the size of the transition region. Simple physical reasoning is used to substantiate the basic assumption. Vigorous growth of the pulsation amplitudes triggering erratic short-scaled wiggles lends credence to physical arguments.

Theoretical conclusions are reinforced computationally. Two mechanisms of the vortex-spot excitation are examined using the Prandtl equations with the self-induced pressure. In both cases, we start with the external disturbance/surface roughness interaction and then pursue the nonlinear development of the emitted wave packet into spiking and deep transition. A convected vortex furnishes the simplest possible mathematical model of a vortical patch impinging against a tiny surface roughness. Sound-pulse scattering into Tollmien–Schlichting eigenmodes represents another mechanism typical of the vortex-spot production. The computation reveals the boundary-layer receptivity to both types of external perturbing agencies. Wave packets emitted amplify in an explosive way characterized by the fast filling-up of distant parts of the wavenumber spectrum. Short wavelength disturbances randomly riding on the primary oscillation cycles in figures 2*c* and 4*c* are the outcome of the filling-up process. Thus, the final phase of transition does not appreciably change its length, primarily determined by viscous/inviscid interaction and a singularity arising in the velocity field (Smith 1988). The formation of eruptive spires is seen in dye visualization of material spikes by Peacock *et al*. (2005).

A comment is due in conclusion. The triple-deck space scaling defined in (2.6*b*) and (2.13*b*) depends on the non-dimensional skin friction , the Mach number and the ratio , which can vary with thermal conditions on the blade surface. Nevertheless, the experimental correlations (2.2), (2.4) and (2.5) involve solely the Reynolds numbers defined in a different way. It would be extremely advisable to put to a test experimentally all the theoretical predictions ensuing from the similarity law cast in (2.6*a*)–(2.6*c*) or (2.13*a*)–(2.13*c*), the more so that the temperature dependence is rather strong and can exert a great impact on design. As mentioned, no experimental data exist so far for the general case of natural transition on turbomachinery blades. In this regard, results from DNS (direct numerical simulation) would be of considerable importance.

## Acknowledgments

Thanks are due to both referees and the adjudicator for their helpful comments. The work was sponsored by the Air Force Office of Scientific Research, Air Force Material Commands, USAF under grant no. F49620-03-1-0043 and monitored by Dr A. Nachman.

## Footnotes

- Received July 20, 2005.
- Accepted December 21, 2005.

- © 2006 The Royal Society