EFFECT OF STABILITY: MONIN-OBUKHOV SIMILARITY THEORY
SCALING AND UNDERLYING HYPOTHESIS So far, we have described wall-bounded turbulent flow in neutral conditions. However, in the ASL the presence of buoyancy forces alters the scaling of the previous section by introducing an additional scale, the Obukhov scale [Obukhov, 1946, Monin and Obukhov, 1954] where g is gravity acceleration, and θ0 and θ′ are a reference potential temperature and its turbulent deviation respectively. This length accounts for a non-zero turbulent heat flux w′θ′, which is positive (upwards) for an unstable atmosphere, and negative (downwards) for a stable atmosphere. 8 Scaling the height with the Obukhov length defines the stability parameter ζ = z/L, positive for a stable atmosphere and negative for an unstable atmosphere. The absolute value of this dimensionless length is characteristic of (i.e. proportional to) the height of the « dynamical sublayer », below which the effect of stratification is negligible [Obukhov, 1971].
This scale, forms the basis of the Monin-Obukhov Similarity theory [MOST, see the review by Foken, 2006]. Using the Buckingham Π theorem, the following dimensionless relations were posited
where φm and φh are the so-called MOST universal functions for momentum and heat fluxes respectively, equal to one for neutral conditions. It is interesting to highlight that these relations can be interpreted as a modification of the diffusion properties of turbulence in the presence of heat fluxes. In fact, the momentum Km and heat Kh turbulent diffusion coefficients diffusion coefficients. Equation (1.8) shows that (i) φm accounts for the stability-dependent deviation of the turbulent diffusion coefficient from its neutral value in the logarithmic sublayer (κz)2(dU /dz) 9 (ii) φh accounts for the influence of stability on the dissimilarity between turbulent transport of momentum and heat (see the recent analysis of Li et al. [2012, 2015] and the review of Li ).
This interpretation indirectly highlights two of the fundamental assumptions of MOST. The first assumption is the existence of linear gradient-flux relationship (i.e. the fact that turbulent fluxes are related to the local properties of their averaged counterparts, Eq. (1.8)). It is well known that this assumption may not be valid in the presence of large scale convection or in the roughness sublayer, where non-local (counter-gradient) fluxes are present [e.g. Kaimal and Finnigan, 1994]. The second assumption is the existence of a logarithmic-law scaling, which relies on several other conditions. As shown in Eq. (1.5) (and the related discussion) those are the ABL hypothesis and Townsend’s view of energy-containing eddies scaling with z (Fig. 1.2). The ABL hypothesis can be violated for inhomogeneous and non-stationary ABLs. According to Wyngaard , for a typical height z = 10 m corresponding to standard measurements, the horizontal length scale on which the ABL has to be homogeneous should be larger than a few kilometers, while the time scale of variation should be greater than 5 minutes. The ASL constant-stress assumption can also be invalidated by the presence of pressure gradients or of dispersive fluxes, which change the momentum balance (1.2).
With respect to the validity of Townsend’s view of energy-containing eddies, let us first mention that, while the presence of stratification should obviously alter the length scale of the eddies, it is still proportional to z in MOST, consistent with Townsend’s scaling. Nonetheless, additional non-local scales can enter the problem for several reasons. First, strong inner-outer layer interactions, in particular in the presence of large convective structures, introduce outer-layer or environmental parameters in the similarity theory [e.g. Li et al., 2018, Fodor et al., 2019]. Additional scales can also be introduced due to the presence of roughness elements [e.g. Garratt and Hicks, 1990, Zilitinkevich et al., 2006]. Note that in the limiting free convective conditions [of strong instablity, Tennekes, 1970] and very stable conditions [z-less conditions, Wyngaard, 1973], which overlap some of the cases mentioned above, MOST scaling is also invalidated. Here and in the rest of the thesis, the discussion will be restricted to mildly stable and unstable conditions, defined for |z/L| ≤ 1.
THE O’KEYPS EQUATION Empirical forms of φm and φh have been derived in numerous experiments, in particular by Businger et al.  using data from the Kansas experiment. Those are not presented here, and can be found in Appendix A of Sec. 2.2. Note however that the general behavior is an increase of φm and φh with ζ, reminiscent of a decrease of turbulent intensity for stable conditions and an increase for unstable conditions.
To match the Businger empirical functions, several authors proposed the so-called O’KEYPS equation (after Obukhov, Kaimal, Elliot, Yamamoto, Panofsky, Sellers, see Panofsky  and Businger ). The equation reads φm4 z φm3 = 1, (1.9) −γ L where γ is a « universal » constant which varies, depending on the experiments, between 5 and 18.
The O’KEYPS equation was derived in Ellison  to fit two limiting regimes mentioned above: the free convective regime for ζ ≪ −1 (where wind shear can be neglected) and the stable regime for ζ ≫ 1. Alternative heuristic derivations where presented in Panofsky  and Sellers , considering that shear- and buoyancy-induced turbulence contribute independently to the total turbulent diffusion, in which case γ appears as a proportionality coefficient between these two processes. Using these hypothesis, these studies stress that the O’KEYPS equation is valid only for an unstable atmosphere. Further, these derivations show that γ is proportional to the turbulent Prandtl number Pr. Hence the stability-independence of γ implies that of Pr, which can be questioned in ABLs [Li, 2019]. The « universality » of γ was discussed in Katul et al. , who related γ to turbulent transport and to properties of energy-containing eddies, showing its sensitivity to stability and environmental conditions.
TKE BALANCE We have mentioned that MOST accounts for the change in turbulent intensity due to buoyancy. The origin of this modification can be traced back to the TKE equation (1.5) (page 30) which, in the presence of stratification, reads dU g (1.10) −u′w′ + w′θ′ = ǫ. dz θ0
The additional term with respect to Eq. (1.5) is production or destruction of TKE by buoyancy forces. This additional term explains the changes in turbulence properties when stability is included. It defines the flux Richardson number Rif Rif = u′w′ ∂U∂z .
The flux Richardson number can be related to the gradient Richardson number Ri= (g/θ0) [∂zθ/(∂zU)2] and to ζ using the flux gradient relations (1.7), and reads Ri f = Kh Ri = ζφ−1. (1.11) This gives another interpretation of ζ in terms of the energetic properties of turbulence, as a the ratio between energy production/destruction by buoyancy and mechanical energy production.
The TKE balance can be rewritten in dimensionless form by multiplying Eq. (1.10) by κz/u∗3, and reads
In the above TKE balance, transport terms have again been neglected. Measurements for mildly unstable conditions, shown in Fig. 1.4, indeed show that the two components of turbulent transport (transport by turbulent fluctuations and by pressure fluctuations) almost compensate each other. For mildly unstable conditions, a first order correction can nonetheless be included in this budget, by replacing ζ by (1+β2)ζ, where β2 ∼ 1 is a proportionality coefficient relating the dimensionless transport to ζ [e.g. Katul et al., 2011].
A PHENOMENOLOGICAL SPECTRAL LINK
In the previous section, the description of wall-bounded turbulence has been related to statistical eddies, transporting momentum and energy both in physical space (the downward turbulent momentum flux) and in spectral space (the 3D isotropic cascade). In particular, the central role of « attached eddies », which interact with the mean wind shear at a given height, has been emphasized. Relying on this phenomenology, several authors have proposed a model, termed « spectral link », relating the averaged (or bulk) properties of turbulence to its spectral properties. The original model of Gioia et al.  and its extension for a stratified flow by Katul et al.  are first reviewed to introduce their use in the modeling of the mean velocity and MOST functions. They rely on the geometry of the eddies, without necessity of discussing spectral budgets. The latter are then discussed, along the lines of Katul et al.  and Katul and BULK MODEL At the core of the concept of an attached eddy is the concept of a momentum-transporting eddy, i.e. which supports u′w′ at a given height. This was made explicit by Gioia et al.  by expressing the momentum flux across a horizontal surface at a height z as the product of an energy-containing eddy turnover velocity ve, depending on the streamwise scale of the eddy sh, and a streamwise momentum difference transported by the eddy across a vertical distance sv . To first order, the momentum difference u(z + sv ) −u(z −sv ) reads 2∂U∂z sv , which yields the following expression for the momentum flux
If the mean wind speed follows a logarithmic law (Eq. (1.4)), the scaling ǫ = u∗3/(κz) is recovered by requiring sv (z) = sh(z) = z in Eq. (1.14) [Gioia et al., 2010]. This can be seen as indicative that eddies considered in this model are related to Townsend’s attached eddies, since both have sizes which scale height from the ground (see Fig. 1.2a). Note that, for consistency, the proportionality factor in Eq. (1.14) has to be in this case κ−4.
For a stratified atmosphere, the mean wind profile deviates from the logarithmic scaling (as mentioned in the previous section). In addition, Katul et al.  argued that, while the height of the energy-containing eddies sv should still be equal to z (since they are attached to the surface, as in Fig. 1.2a), their horizontal extension, and hence energy, can vary due to buoyancy forces. The authors hence introduced the so-called anisotropy factor fa(ζ) = sh(ζ)/z, larger than one for an unstable atmosphere (in which eddies are more energetic) and smaller than one for a stable atmosphere (in which eddies are “squeezed” due to the restoring buoyancy forces).
For unstable conditions, the anisotropy factor can be experimentally determined as a change in the ratio between the spanwise and streamwise integral length scales of vertical velocity Salesky et al. . For stable conditions, Li et al.  showed that a more appropriate estimate is the Ozmidov length scale (the size of the largest eddy unaffected by buoyancy). In both cases, using empirical estimates for fa(ζ), the authors were able to predict experimentally-consistent values of φm, giving support for the capacity of the phenomenological model to predict mean flow variables. THE SPECTRAL LINK The underlying idea of the phenomenological model is to relate the bulk properties of the flow (here TKE dissipation) to its spectral properties, i.e. the energy and size of energy-containing eddies. In fact, Gioia et al.  had suggested that the eddy turnover velocity ve(s) should be computed as an integral of the TKE spectrum, truncated at the energy-containing eddy wavelength. Katul and Manes  went one step further by using a co-spectral budgets, i.e. budgets of co-spectrum amplitude at a given eddy scale, to derive Eq. (1.13).
which is consistent with the posited stability-dependent anisotropy factor fa(ζ)−1. Equations (1.19) and (1.21) establish a « spectral link » relating spectral properties of the flow to its bulk characteristics. They show that the latter are in fact not only determined by properties of the eddies whose size is commensurate to z, but by a weighted contribution of eddies of all scales.
With respect to the bulk analysis of Gioia et al.  and Katul et al. , the co-spectral budget unravels the physical hypothesis and shortcomings of the phenomenological model. The first assumption concerns the terms retained in the co-spectral budget. While, as for its bulk counterpart, the spectral budget is expected to be satisfied in flat-wall ASLs, the flux-transport term (the spectral analogue of the vertical transport term in the bulk budget) has been ignored. It can be expected that this term becomes important for some other type of flows, e.g. in the roughness sublayer, even if measurements for its quantification are lacking. Setting aside that assumption, the spectral link relies on two models: for the vertical velocity spectrum (Eq. (1.20), which sets the production of co-spectrum, and for the return-to-isotropy term (Eq. (1.18)). Both models are expected to be valid for mildly stable and unstable conditions, in a variety of ASL flows [see the references in Katul et al., 2013, Katul and Manes, 2014]. Nonetheless, it is important to stress that the LRR-IP model can fail for strongly inhomogeneous flows, in which the timescale of return to isotropy can exceed the modeled T or τ˜.
As a final note, let us mention that, herein, we have discussed the spectra only in the wavenumber domain. Conversion to the frequency domain requires a dispersion relation for turbulent motions. The simplest conversion relies on Taylors’ hypothesis of frozen turbulence [Taylor, 1938] which assumes that the advection time of turbulent motions of wavenumber k is smaller than their development time (defined as the time for turbulent structures to change significantly in their frame of reference, see Phillips ). For turbulent structures advected at the mean wind speed U, the conversion is then ω = kU . (1.22)
However, more complex models have been used [see the review Wallace, 2014], including e.g. the turbulent advection of eddies of a given scale by larger eddies [Wilczek et al., 2015]. Finally, the presence of surface roughnesses can invalidate Taylor’s hypothesis. Squire et al.  and references therein used measurements to demonstrate that the advection of small scale structures can be modified, but it is uncertain to which extent. The signature is especially important on the vertical velocity spectrum. This should be carefully considered when using the spectral link in complex flows, since the latter relies on spectra expressed in the wavenumber domain, which are often measured in the frequency domain.
FINAL REMARKS The phenomenological model relates the bulk turbulent momentum flux to the spectrum of the vertical velocity, whose properties can be represented by means of an idealized attached eddy. This statistical eddy has dimensions which should be the expression of the complex instantaneous features of turbulence, e.g. the coherent structures discussed in Sec. 1.1. This model allows recovering some of the « universal » MOST functions, and gives hints on their non-universal behavior. It goes beyond MOST scaling arguments, which gives support to its applicability to establish other spectral links.
As an important example of such a link, the phenomenological model has been used to predict the properties of rough pipe flows, using its bulk formulation [Gioia and Bombardelli, 2001] and the co-spectral budget [Bonetti et al., 2017]. The effect of roughness elements has been included by assuming that, in the roughness sublayer, the size of the energy-containing eddy sh (or the spectral peak wavenumber kp) is constant with height and equal to the size of the roughness elements. This is consistent with the measurements presented by Raupach et al. , and their interpretation as indicative of energy-containing eddies being mixing-layer eddies, whose size is set by the height of the roughness elements (see Fig. 1.3e).
The model is hence an interesting candidate for a simplified description of the complex properties of turbulence over rough surfaces, of which the windy sea is an example of. With this in mind, we first review below the properties of turbulence on top of waves, from the theoretical modeling perspective.
THE DYNAMICAL INTERACTION BETWEEN NEAR-SURFACE TURBULENCE AND WAVES
Understanding the properties of turbulent motions close to the sea surface is of uttermost importance for air-sea interactions, and is a complex problem due to the presence of ocean surface waves [see, e.g., the recent reviews on wind-wave interactions in Jones et al., 2001, Janssen, 2004, Sullivan and McWilliams, 2010, LeMone et al., 2019]. Field observations [e.g. Edson et al., 2013] indicate that, in the presence of surface waves and for sufficiently strong winds, the turbulent momentum flux on top of the so-called wave boundary layer (WBL), of height of about 10 m, is increased with respect to a flat surface, and has a one-to-one dependence on the mean wind speed (for averaging periods of about 30 minutes). Hence the disturbances generated by surface waves, whose amplitude is coupled to atmospheric motions, result in an overall change in the properties of turbulence in the WBL, for a prescribed mean 10m-wind [see experiments of Edson and Fairall, 1998, Sjöblom and Smedman, 2002]. In the following we review the dynamical properties of the WBL, with emphasis on theoretical (and, when possible, analytical) models for the interaction between atmospheric turbulence and waves. What follows is by no means an exhaustive or historical review, but rather one with a specific focus theoretical studies of the WBL.
Flow over surface waves, while sharing some similarities with flow with hills and wavy boundaries, exhibits some unique features due to the intrinsic properties of surface waves [see the review by Belcher and Hunt, 1998, where both flows are compared].
First, surface waves are moving undulations of the sea surface, which follow a dispersion relation. More precisely, the phase speed c = ω/k of a monochromatic wave depends on its wavenumber k and frequency ω, which are linked as ω2 = gk + Tswk3, where g is the gravity acceleration and Tsw is the dynamical surface water tension. The phase speed has a minimum for waves of wavelength (λ = 2π/k) of 1.6 cm, which marks the transition between capillary (smaller) and gravity (larger) waves. The phase speed of gravity waves increases with their size, as opposed to capillary waves. This first feature implies that the the impact of waves on atmospheric turbulence should depend not only on their geometry, but also on relative velocity of the wave with respect to the airflow, and hence on the scale of the wave [see for instance Kitaigorodskii, 1973, p. 27 to 36, where the surface is modelled as a linear superposition of moving roughness elements]. The relative velocity of the wave is termed wave age, c/u∗, with u∗ the friction velocity (the square root of the momentum flux on top of the WBL).
Second, the airflow bottom boundary condition is non-uniform, and depends on the wave steepness. As explained by Kraus , the difference in dynamic viscosity between the air and the water results in the air-sea interface moving at the wave orbital velocity.12 The wave orbital velocity varies periodically in the reference frame of the wave, and its magnitude, akc, depends on the wave steepness ak (a is the wave amplitude). As waves grow and decay locally under the action of the wind, their steepness changes [Longuet-Higgins, 1987]. Hence the airflow bottom boundary condition is non-uniform, and is dynamically coupled to the wind and wave fields.
Third, the growth and decay of waves under the action of the wind is a complex mechanism. It cannot be simply explained by Kelvin-Helmoltz type of instability [as originally proposed by Lord Kelvin, Thomson, 1871]. Indeed, the suction on top of wave crests required to overcome the gravity restoring force can only be strong enough for winds of about 6-7 m s−1, for which waves are already present [Ursell, 1956]. In addition, wave breaking is an essential component of the sea surface, which is associated with intense energy dissipation events from the wave field towards ocean currents, but also to intense ejection events and momentum flux spikes in the atmosphere [Banner and Melville, 1976, Kawamura and Toba, 1988, Melville, 1996].
Finally, a realistic sea surface is described by a broadband wave spectrum: it is a multiscale interface which can be, to some extent, viewed as a sum of individual sinusoidal components, each moving at a different phase speed. This implies three additional difficulties: (i) the roughness properties of waves of different sizes (quantified by e.g. the roughness Reynolds number) span a wide range of regimes, and hence the overall effect of the multiscale sea surface is non trivial [Kitaigorodskii, 1973]; (ii) the response of a given wave component to the airflow depends on its non-linear interaction with other wave components, which redistributes energy among different wavenumbers [Phillips, 1985]; (iii) the response of the turbulent airflow to a multiscale surface can be qualitatively different from the superposition of responses to individual wave components [Deardorff, 1967].
From a practical point of view, it is interesting to represent the impact of those various features on the atmosphere through a roughness height z0 (the level at which the mean wind speed, following a logarithmic law, cancels). A very influential scaling was proposed by Charnock  for sufficiently high winds: z0 ∝ u∗2/g, where g is the gravity acceleration. This scaling can be interpreted in light of the observations of Francis  and Munk , which indicate that the wind stress (equivalent to z0) should essentially depend on short wind-waves (smaller than 1m). As argued by Phillips  (p. 194), the amplitude of short wind-waves depends on wind-induced surface drift (related to u∗) and wave breaking, and hence scales with u∗2/g, explaining the Charnock relation.
Table of contents :
1 Turbulence and waves: a literature review
1.1 Wall-bounded turbulence with stratification
1.1.1 General considerations for neutral conditions
1.1.2 Effect of stability: Monin-Obukhov Similarity Theory
1.2 A phenomenological spectral link
1.3 The dynamical interaction between near-surface turbulence and waves
1.3.1 Momentum balance in the wave boundary layer
1.3.2 Theoretical models of wave-induced stress
1.3.3 Turbulence in the wave boundary layer
1.3.4 Concluding remarks
1.4 Objectives of the present work
2 Towards a « spectral link » for the vertical velocity spectrum?
2.2 Article: « Scalewise return-to-isotropy in stratified boundary layer flows »
2.3 Article: « Scaling laws for the Length Scale of Energy-containing Eddies in a Sheared and Thermally Stratified Atmospheric Surface Layer »
3 Geometrical impact of wind-waves on energy-containing eddies
3.2 Article: « On the impact of long wind-waves on near-surface turbulence and momentum fluxes »
4 Dynamical impact of wind-waves on energy-containing eddies
4.2 Article: « Revisiting Beaufort scale: the dynamical coupling between turbulence and breaking waves »
4.2.1 Main text
4.2.2 Supplementary material
Annex A: Some steps for the derivation of the spectral budget
Annex B: Comparison of several return-to-isotropy models
Annex C: Details on the preliminary numerical simulation