# The speed of breaking waves controls sea surface drag

Alex Ayet1,2<EMAIL_ADDRESS>Bertrand Chapron1  Peter Sutherland1
Gabriel G. Katul3  1Ifremer, CNRS, IRD, Univ. Brest/ Laboratoire
d’Océanographie Physique et Spatiale (LOPS), IUEM, Brest, France 2LMD/IPSL,
CNRS, École Normale Supérieure, PSL Research University, Paris, France
3Nicholas School of the Environment and the Department of Civil and
Environmental Engineering, Duke University, Durham, NC 27708-0328, USA

###### Abstract

The coupling between wind-waves and atmospheric surface layer turbulence sets
surface drag. This coupling is however usually represented through a roughness
length. Originally suggested on purely dimensional grounds, this roughness
length does not directly correspond to a measurable physical quantity of the
wind-and-wave system. Here, to go beyond this representation, we formalize
ideas underlying the Beaufort scale by quantifying the velocity of breaking
short waves that are the most coupled to near-surface wind. This velocity
increases with wind speed, reflecting the fact that stronger winds can be
visually identified by longer (and faster) breakers becoming predominant on
the sea surface. A phenomenological turbulence model further shows that this
velocity is associated with breaking waves that impede the most the formation
of turbulent eddies. Scales of such eddies are then constrained inside a so-
called roughness sub-layer. Unlike previous theoretical developments, the
proposed breaker velocity is a directly measurable quantity, which could be
used to characterize the coupling between wind and waves using remote sensing
techniques. This work provides a physical framework for new formulations of
air-sea momentum exchange in which the effects of surface currents and slicks
on surface drag can also be incorporated. Finally, it provides a long-sought
physical explanation for the Beaufort scale: a universal link between wave
breaking, wind speed and surface drag.

###### keywords:

##  1 Introduction

The interaction between near-surface turbulence (up to a few meters from the
surface) and wind-waves (waves shorter than ten meters) determines air-sea
fluxes of momentum, heat and gases. Yet, this link is not readily expressed in
current parameterizations of surface drag, which rely on a roughness length
(e.g. Kitaigorodskii, 1973; Edson et al., 2013). Originally suggested on
purely dimensional grounds (Charnock, 1955), this roughness length does not
directly correspond to a measurable physical quantity of the wind-and-wave
system (Kraus, 1967). More importantly, it does not reflect the fact that the
wind blows over a multi-scale and moving surface (each wave moving with its
own phase speed). As such, it is difficult to capture (experimentally and in
parameterizations) the sensitivity of air-sea fluxes to external parameters
such as ocean surface currents, which alter the speed of wind-waves and hence
surface drag.

However, seafarers have known for centuries how to estimate wind velocity by
observing the local sea surface. A windy sea indeed demonstrates an apparent
visual organization associated with the occurrence and intensity of breaking
waves. This led George Simpson, in 1906, to derive a scale for the surface
wind speed, labeled the Beaufort scale after the original classification of
the sea surface by Frances Beaufort in 1831. Such an estimation reflects the
existence of physical quantities, related to short-scale breaking waves, that
could be used to characterize the wind-and-waves system and hence air-sea
fluxes. This work presents one such quantity, aimed at characterizing the
momentum flux at the air-water interface.

In fact, it has long been known that the statistics of short breaking waves
are tightly linked to the local properties of near-surface turbulence, surface
drag, and wind speed, both from measurements (Melville & Rapp, 1985; Katsaros
& Ataktürk, 1992) and from theoretical studies (Melville, 1977; Phillips,
1985; Newell & Zakharov, 1992; Kudryavtsev et al., 2014). On the one hand,
breaking waves are steep disturbances of the sea-surface that influence its
overall roughness (Melville, 1977; Kudryavtsev et al., 2014). On the other
hand, wave breaking and the formation of whitecaps play a fundamental role in
the energy transfer between wave scales, and hence in the slope distribution
of short waves (Phillips, 1985), which is ultimately related to surface
roughness and drag (Munk, 1955). Including this link in parameterizations is
difficult due to the challenges associated with modeling a surface that is no
longer connected (i.e. disrupted by wave breaking), and for which the
governing equations of motion are not fully known (Newell & Zakharov, 1992).

Recently, technological advances resulted in high resolution images of the sea
surface, which can now be used to advance our understanding of short-scale
wave breaking statistics (Sutherland & Melville, 2013, 2015). In parallel,
phenomenological models of turbulence close to a boundary have emerged (Gioia
et al., 2010): those models establish links between the spectral properties
and mean properties of the turbulent flow (Katul & Manes, 2014; Ayet & Katul,
2020), and have been used to describe their modification due to changes in
atmospheric stability (Katul et al., 2011), the presence of roughness elements
(Gioia & Bombardelli, 2001; Bonetti et al., 2017) and, recently, the presence
of long wind-waves (Ayet et al., 2020a).

In this work, both advances are used to show that the complex interaction
between short, breaking wind-waves and wind can be expressed as function of a
representative velocity of short-scale breaking fronts. The wave properties
associated with the predicted phase speed correspond to the Beaufort scale.
This velocity is derived independently (i) from high resolution measurements
(in §2) and (ii) by extending a phenomenological model of turbulence (in §3).
Finally, in §4 we show that this representative speed is associated with
breaking waves that efficiently impede the generation of turbulent eddies.
Unlike previous theoretical developments, this velocity is a directly
measurable trace of the dynamical wind-wave coupling, possible from remote
sensing platforms.

##  2 Analysis of wave-breaking measurements

Figure 1: (a) Momentum flux $u_{*}^{2}$ versus $U_{10}$ from the COARE 3.5
parameterization (solid line) and the corresponding open-ocean measurements
(dots and bins, from Edson et al., 2013). The relation derived for a smooth-
wall boundary layer is included for reference (dashed line). The open-ocean
measurements used in this study (Sutherland & Melville, 2013, 2015) are shown
as diamonds. In §3, only the solid line is used in the turbulence
phenomenological model to obtain the representative breaker speed. (b)
Atmospheric stability versus 10-m wind speed for the open-ocean measurements
used in this study. Atmospheric stability was quantified by the dimensionless
stability parameter $z/L$, where $z$ is the height of the measurement, and $L$
the Obukhov length. Only measurements for which $-0.2\leq z/L\leq 0.2$ have
been used in §2.

###  2.1 Data

Open-ocean measurements of wave-breaking statistics from Sutherland & Melville
(2013, 2015) were used. The measurements were obtained during three field
campaigns on-board R/V _FLIP_ : the RaDyO experiment, which took place 120 km
south of the Island of Hawai’i in 2009 (Dickey et al., 2012), and the HiRes
and SoCal experiments, which took place 25km off the coast of Northern
California and in the Southern California Bight, respectively. As shown in
figure 1a (red diamonds), the 10-meter mean wind $U_{10}$ (averaged over
20-minute intervals) ranged from 2 to 16 m s-1 with few measurements for
$U_{10}>10$ m s-1. Figure 1a also shows the measured friction velocity $u_{*}$
defined by $u_{*}^{2}=\tau/\rho_{a}$, where $\tau$ is the measured turbulent
stress at $10$-m height and $\rho_{a}$ is the mean air density. These
measurements are consistent with other experiments (summarized in Edson et
al., 2013, dots and bins) used to calibrate the widely used COARE
parameterization (solid line). The presence of wind-waves at equilibrium with
the wind results in a bulk drag coefficient $C_{d}=u_{*}^{2}/U_{10}^{2}$ that
is higher than a $C_{d}$ associated with a flat-solid surface (dashed line,
and see also e.g. Ayet et al., 2020a, their figure 3). Except for very low
winds, atmospheric stability was near neutral during the measurements (see
figure 1b). For this reason, atmospheric stability effects are not discussed
further and low $U_{10}$ runs are filtered out so as to only maintain near-
neutral atmospheric stability conditions. In the measurements, no strong
correlation was found between variations in atmospheric stability conditions
for low winds and variations in $u_{*}$ for a given $U_{10}$.

To each point in figure 1a corresponds a $20$-minute averaged spectrum of
wave-breaking statistics, measured using a stereo pair of long-wave infrared
cameras. More precisely, the stereo images permit the estimation of the length
and speed of the breaker fronts (Sutherland & Melville, 2013). These
measurements can be used to compute an effective $20$-minute averaged spectrum
$\Lambda(c)$ of the breaker front length per unit area of sea surface per unit
increment of breaking front velocity $c$, which we analyse in the following.
In contrast to measurements in the visible range (Kleiss & Melville, 2010),
the use of infrared cameras allows the measurement of non-air-entraining
micro-breakers with crest speeds as slow as the gravity-capillary phase speed
minimum (corresponding to a wavelength of 1.7 $\times 10^{-2}$ m).

###  2.2 Analysis and results

Figure 2: Analysis of the wave breaking measurements. (a) Average and standard
deviation of $c^{3}\Lambda(c)$ for $U_{10}$-bins, for different breaker
speeds, vs $U_{10}$ (the $U_{10}$-values correspond to the middle of the bin).
Dots and error bars correspond to data binned in $U_{10}$ intervals of 3 m
s-1, and shadings to intervals of 1 m s-1. For $U_{10}$-bins of $3$ m s-1, (b)
binned $c^{3}\overline{\Lambda}(c,U_{10})$ vs breaker speed for different
$U_{10}$. and (c) gradient of the binned $c^{3}\overline{\Lambda}(c,U_{10})$
with respect to changes in $U_{10}$ vs breaker speed, for different $U_{10}$.
The vertical dashed lines correspond, for each bin, to the breaker speed whose
gradient $\partial[c^{3}\overline{\Lambda}(c,U_{10})]/\partial U_{10}$ is the
largest, with colours consistent with (a).

The statistics of the breaker front length $\Lambda(c)$ were originally
introduced as a measurable quantity to characterise the wind-and-waves
equilibrium (Phillips, 1985). Here, we argue that the third moment of this
distribution ($c^{3}\Lambda(c)$) is the appropriate quantity to characterise
air-sea momentum fluxes. Originally related to mass fluxes at the air-sea
interface (Deike et al., 2017; Deike & Melville, 2018), $c^{3}\Lambda(c)$ can
also be linked to the roughness elements generated by breaking waves: near
breaking crests, local energy fluxes disrupt the sea surface and generate a
spectrum of roughness elements (i.e. roughness elements at different
horizontal scales). Newell & Zakharov (1992) suggested that these roughness
elements are at a Kolmogorov-type equilibrium. This equilibrium requires that
the production rate, the energy transfer rate across scales, and the
dissipation rate be balanced and scale-independent (i.e. a conservative
cascade). Based on this idea, Kudryavtsev et al. (2014) showed that this
energy flux across scales should be proportional to $c^{3}$ (where $c$ is the
phase speed of the breaker front generating the roughness elements). The
authors finally demonstrated that this spectrum of roughness elements modifies
near-surface turbulent stress and that the intensity of the modification
depends on $c^{3}\Lambda(c)$. Hence $c^{3}\Lambda(c)$ is physically related to
the effect of breaking waves on near-surface turbulence.

To analyse this effect further, data were first binned in $U_{10}$-intervals
of 3 m s-1. Dots and error bars in figure 2a show the resulting
$U_{10}$-dependent average $c^{3}\overline{\Lambda}(c,U_{10})$ and standard
deviation for each bin and different breaker speeds, normalised by the
$c$-dependent maximum of $c^{3}\overline{\Lambda}(c,U_{10})$ with respect to
$U_{10}$. Examples of the dependence of $c^{3}\overline{\Lambda}(c,U_{10})$
with $c$ are also shown in figure 2b for reference. Figure 2a shows that, for
a given breaker speed, the increase in $c^{3}\overline{\Lambda}(c,U_{10})$ is
not constant with increasing wind speed, but seems be strong around a
particular wind speed before reaching a saturation value: this is mostly
apparent for the blue line, corresponding to a breaker speed of $0.85$ m s-1.
To further quantify this effect, i.e. how _changes_ in $U_{10}$ leave imprints
on the local wind-waves breaking field, what is analysed in the following is
the _gradient_ of $c^{3}\overline{\Lambda}(c,U_{10})$ with respect to $U_{10}$
as shown in figure 2c for different wind speeds.

Figure 2c reveals that, for a given $U_{10}$ this gradient
($\partial[c^{3}\overline{\Lambda}(c,U_{10})]/\partial U_{10}$) is found
strongest around a speed $c_{r}$, termed representative breaker speed in the
following (vertical lines in figure 2b). This speed increases with $U_{10}$:
for $0.5<U_{10}<3.5$ (m s-1), it is $c_{r}=0.85$ m s-1 (blue vertical line in
figure 2c) while for $12.5<U_{10}<15.5$ (m s-1), it is $c_{r}=1.75$ m s-1
(orange vertical line in figure 2c). Note that, even though for
$15.5<U_{10}<18.5$ (m s-1) such a wave speed can still be defined
($c_{r}=2.35$ m s-1, yellow line in figure 2c), the low number of measurement
points contained in the bin induces noise in the gradient of
$c^{3}\overline{\Lambda}(c,U_{10})$.

Hence, for each $U_{10}$ interval, a representative breaker speed $c_{r}$ can
be defined, corresponding to the short-scale waves whose breaking statistics
($c^{3}\Lambda(c)$) are the most sensitive to changes in mean wind speed
(solid blue line in figure 3). Sensitivity to the bin size was tested using
10-m wind speed binning intervals of a width of 1 and 2 m s-1. The resulting
variability of the average value of $c^{3}\Lambda(c)$ in each bin is shown as
shadings in figure 2c, and as dashed and dotted lines in figure 3. The
outcomes are qualitatively consistent with the analysis above.

Finally, a sensitivity analysis on the choice of the moment $n$ of the wave
breaking distribution $c^{n}\overline{\Lambda}(c,U_{10})$ was also performed.
Indeed, even though there is a physical reason to the choice $n=3$, other
exponents could have been used to extract a representative breaker speed from
the data. As explained in appendix A, the exponent $n=3$ used above is the one
for which the representative speed $c_{r}$ is the most unambiguously defined.
This choice is also justified by the agreement with the turbulence model (see
§3), and by the tri-dimensional view of eddies filling the space between waves
presented in §4.

###  2.3 A measurable kinematic relation

A key property of the representative breaker speed $c_{r}$ is its kinematic
link with the friction velocity

$c_{r}=2.5u_{*}.$ (1)

While an approximate form of this relation can be obtained from a best fit to
the data (in figure 3), the turbulence model presented in the next section
provides an exact derivation of (1).

Such a relation was already discussed in Phillips (1977, p. 141) when
considering the interaction between small scale waves and surface drift
induced by viscous wind stress: Phillips (1977) argued that only waves with
phase speeds smaller than five times the surface drift would be strongly
impacted by it. The magnitude of the surface drift was further reported to be
$q=0.55u_{*}$ m s-1 (Wu, 1975). Hence, wave scales such that $c=5\times
0.55u_{*}\sim 2.5u_{*}$ are the smallest waves that do not significantly
interact with surface drift. Interactions with surface drift tend to reduce
the critical steepness needed for wave breaking. This leads to (i) a reduction
of the amplitude of drift-affected waves with respect to a no-drift situation,
and hence to (ii) a reduction of the strength of the breaking event, i.e. its
capacity to be visually detected.

Consequence (ii) can explain the existence of a maximum in
$\partial[c^{3}\overline{\Lambda}(c,U_{10})]/\partial U_{10}$ (evident in
figure 2c). With increasing winds, short wave-breaking occurrence is expected
to increase. However, this increase is truncated for drift-affected waves.
This leads to the Bell-shaped curves of figure 2c, i.e. to the fact that a
representative surface wave speed can be visually defined for a given wind:
the Beaufort scale. As for consequence (i), it implies that amplitudes of
drift-affected waves are reduced, and so is their contribution to the
roughness elements affecting atmospheric turbulence, such as discussed in the
next section.

To summarize, breaker data analysis suggests that it is possible to relate a
change in $U_{10}$ to a change in the properties of breaking crests with a
particular wave front velocity. This reflects both local wind stress
variations and changes in wind-induced surface drift.

##  3 A roughness sublayer from a phenomenological turbulence model

In §2, wave measurements were used to show the existence of a representative
breaker speed linked to the local properties of surface wind. In this section,
we follow the opposite reasoning by starting from atmospheric theory to find
the representative breaker speed. We ask the following question: assuming that
there exists a layer close to the surface where turbulence is entirely
constrained by the intermittent roughness elements associated with breaking
waves, what should the defining parameters of such roughness elements be? This
question is addressed using a a phenomenological model of turbulence. Its
original formulation, for flow above fixed roughness elements, is presented in
§3.1. In §3.2, the model is extended to a windy sea, providing an independent
justification for the results of the previous section, and in particular to
(1).

###  3.1 A turbulence model above fixed roughness elements

The phenomenological model, originally presented in Katul et al. (2011)
following ideas of Gioia et al. (2010), describes turbulence in the Surface
Boundary Layer (SBL), where the flow is horizontally homogeneous and
stationary with no subsidence, and hence $u_{*}^{2}$ is height-independent.
Turbulence is modeled through a turbulence kinetic energy balance and so-
called energy-containing eddies (red circles in figures 4b,c), which are
leading order contributors to the vertical momentum flux at a given height
(dashed lines in figures 4b,c).

Those eddies are defined by their streamwise extension or height $s_{e}$ and
their turnover velocity, proportional to $u_{*}$ (e.g. Hunt et al., 1988). For
neutral conditions, the mean wind shear can then be expressed as

$\frac{\partial U(z)}{\partial z}=\frac{u_{*}}{\kappa s_{e}(z)},$ (2)

where $U$ is the mean wind speed in the streamwise ($x$) direction, and $z$ is
height, whose origin is the mean sea surface height, and $\kappa\sim 0.4$ is
the von Kármán constant. The derivation of (2) from the model hypotheses can
be found in appendix B.

Figure 3: Representative breaker speed of the wind-wave local coupling as a
function of $10$-m wind speed. Blue lines are the results from analysis of
wave-breaking statistics $c^{3}\overline{\Lambda}(c,U_{10})$, binned in
different wind intervals. Solid black line is the result from the
phenomenological model of turbulence. The wavelengths ($\lambda=2\pi c^{2}/g$)
correspond to those of deep-water linear waves travelling with phase speed
$c$, and have been included for reference only.

The phenomenological model defined in (2) thus links the mean wind shear to
the friction velocity $u_{*}$ and the size of the energy-containing eddies
$s_{e}$. For flow past a rough wall, Gioia & Bombardelli (2001) and Bonetti et
al. (2017) proposed that the surface boundary can be decomposed into two
sublayers, depicted in figure 4a. The highest sublayer is a logarithmic layer,
in which the size of the energy-containing eddies scales with distance from
the wall ($s_{e}(z)=z$), corresponding to attached eddies (as proposed in
Townsend, 1976, and depicted in figure 4b). From (2), the mean wind profile is
hence logarithmic

$U_{\mathrm{log}}(z)=\frac{u_{*}}{\kappa}\log(z/z_{0}),\;\;\mathrm{for}\;z\geq
h_{r},$ (3)

where $z_{0}$ is the roughness height. Below the logarithmic sublayer lies the
roughness sublayer of height $h_{r}$. In this layer, the size of the energy-
containing eddies scales with the height of the roughness elements that is
proportional to $h_{r}$. To ensure continuity of $\partial U/\partial z$ at
the interface between both layers, the scaling constant is considered to be
equal to one, i.e. $s_{e}(z)=h_{r}$ for $z\leq h_{r}$ (see figure 4c and
Bonetti et al., 2017). Equation (2) then yields a linear wind profile

$U_{\mathrm{lin}}(z)=\frac{u_{*}}{\kappa}\frac{z}{h_{r}},\;\;\mathrm{for}\;0\leq
z\leq h_{r}.$ (4)

The height of the roughness elements $h_{r}$ is different from the roughness
height $z_{0}$ (which is derived from an extrapolation of the logarithmic wind
profile to $U=0$) and, by definition, $h_{r}\geq z_{0}$ (Raupach et al.,
1991).

Finally, a relation between $z_{0}$ and $h_{r}$ can be obtained by requiring
continuity of $U$ at $z=h_{r}$, and reads

$z_{0}=h_{r}\exp(-1).$ (5)

In the aerodynamically smooth regime, $h_{r}$ is the height of the viscous
sublayer, $\exp(1)z_{0}^{v}$, where $z_{0}^{v}=\gamma\nu/u_{*}$ is the viscous
roughness height (and $\nu$ and $\gamma\sim 0.11$ are the kinematic viscosity
of the air and the roughness Reynolds number of smooth flows, respectively).
For this regime, the resulting relation between $u_{*}^{2}$ and $U_{10}$ is
plotted in figure 1a (dashed line).

Figure 4: Phenomenological model describing the surface boundary layer (SBL)
in the presence of breaking waves. (a) The SBL is divided into two sublayers:
(b) the logarithmic sublayer, where energy-containing eddies scale with $z$
(attached eddies); (c) the roughness sublayer, where their size is independent
of height and set to $h_{r}$, the _effective height_ modeling the intermittent
effect of breaking waves. In the roughness sublayer above waves, no mean wind
profile is prescribed, contrary to the case above fixed roughness elements
where it is linear (Bonetti et al., 2017).

###  3.2 Model extension to an ocean surface with breaking waves

An ocean surface with breaking waves of multiple scales and velocities is now
considered. We assume that the SBL can still be decomposed into two sublayers.
In the logarithmic sublayer, the conventional picture of attached eddies
holds. In the roughness sublayer, we seek the measurable properties of the
roughness elements that set the geometry of the energy-containing eddies. To
do so, we adopt a different interpretation of the roughness sublayer than
above fixed roughness elements (Bonetti et al., 2017).

First, we do not assume a linear wind profile in the roughness sublayer.
Instead, we only require the bulk mean wind shear across the roughness
sublayer to follow

$\frac{U(h_{r})}{h_{r}}=\frac{u_{*}}{\kappa}\frac{1}{h_{r}}.$ (6)

This is a weaker assumption than requiring the validity of (2) at all heights
(and hence the linear wind profile (4)). This assumption hence covers the
occurrence of more complex wind profiles that can occur close to breaking
waves, within the roughness sublayer. In fact, the mean wind speed in the
roughness sublayer does not have a well-defined physical interpretation. In
the logarithmic sublayer, $U(z)$ can be interpreted as the convection velocity
of energy-containing eddies of size $s_{e}=z$ (see e.g. Phillips, 1957) while,
in the roughness sublayer, the convection velocity could be ill-defined due to
transient roughness elements (Kraus, 1967). Certainly, it is more reasonable
to assume that the convection velocity of eddies at all heights below $h_{r}$
is $U(h_{r})$, i.e. the convection velocity of the roughness-sublayer eddies
of size $s_{e}=h_{r}$.

The second difference with respect to flow over fixed roughness elements is
the interpretation of $h_{r}$. We interpret $h_{r}$ as being an _effective_
height of the roughness elements, modeling the bulk effect of intermittent
wind-wave breaking events on the SBL. It does not correspond to a measurable
height (e.g. a monochromatic wave), but is instead only measurable through its
dynamical link with $c_{r}$.

This link results from the dynamical relation between the bulk wind shear and
the properties of energy-containing eddies in the roughness sublayer. Energy-
containing eddies can indeed be seen as resulting from the interaction between
mean wind shear and turbulence: rearranging (2), the streamwise extension of
energy-containing eddies (which defines their properties) can be expressed as
$s_{e}\propto w_{e}T_{e}$, i.e. as a product of an eddy turnover velocity
$w_{e}\propto u_{*}$ and time $T_{e}\propto(dU/dz)^{-1}$. For roughness-
sublayer eddies, the turnover velocity further reads
$T_{e}=(U(h_{r})/h_{r})^{-1}$, indicating that those are controlled by the
_bulk_ mean wind shear over the sublayer (defined in (6)). Roughness-sublayer
eddies can thus be viewed as originating from an instability of a similar form
than above forest canopies (Raupach et al., 1996), for which eddies have a
size $s_{e}$ proportional to the scale of the wind shear $h_{r}$ (this scale
is also called the “vorticity thickness” in e.g. Harman & Finnigan, 2007).
Now, if the wind shear in the roughness sublayer is set by intermittent yet
intense wave breaking events, fluid is entrained at the speed $c_{r}$ of the
breaker fronts on top of the roughness sublayer (as suggested by Melville,
1977). Hence $U(h_{r})$, the magnitude of the bulk wind shear across the
roughness sublayer should be equal to the wave speed, i.e.

$U(h_{r})=c_{r}=u_{*}/\kappa,$ (7)

providing a theoretical justification to (1).

It should be stressed that the argument leading to (7) relies on two
conditions. First, it requires that at the top of the roughness sublayer
$h_{r}$, $U(h_{r})$ is non-zero and is physical. Hence $z_{0}$ could not have
been chosen as a roughness sublayer height, since it is the extrapolated
height where the wind speed $U(z_{0})$ is zero for a logarithmic wind profile.
Second, a logarithmic profile for $z>h_{r}$ is required to obtain (7). This
second point can be questioned since, for $z>h_{r}$, the presence of waves and
wave growth induces wave-induced motions, which cause a deviation of the mean
wind from a logarithmic profile (e.g. Miles, 1965). In appendix C, the
phenomenological model is extended to account for these motions above the
roughness sublayer (following Ayet et al., 2020a). In this case, equation (5)
is more complex, but the resulting $c_{r}$ is of the same order of magnitude
than (7) (see figure 9c).

The extension of the phenomenological model also revealed that as wind
increases and wave-induced motions become more energetic, the ratio
$s_{e}/h_{r}$ increases (up to 2.2 for 16 m s-1, while it is one in the
absence of wave-induced motions, figure 9b). This finding is consistent with
detailed laboratory experiments of flow over granular roughness elements
(Manes et al., 2007) that reported a log region commencing at approximately
1.6 times the mean grain height, coincident with the turbulent flow statistics
attaining spatial uniformity (i.e. the effect of the roughness heterogeneity
is no longer felt). The ratio $s_{e}/h_{r}$ being larger than one also
supports the interpretation of the roughness elements associated to breaking
waves as being spatially and temporally intermittent, with an average size
($h_{r}$) smaller than the scale ($s_{e}$) of the eddies they generate.
Overall, these findings reveal the existence of a direct interaction between
wave growth and roughness-sublayer eddies, through the presence of wave-
induced motions (which are perturbations of the airflow, coherent with wind-
waves, at heights $z>h_{r}$).

To summarise, the phenomenological model used here does not contain _any_ a
priori spectral information on the wave field, and solely requires
specification of the bulk information on turbulence contained in open-ocean
measurements of $u_{*}$ and $U_{10}$ (e.g. the measurements presented in
figure 1a, solid line). Yet, the resulting wave velocity is similar to the one
obtained from analysis of independent wave-breaking statistics. This is shown
in figure 3, where the wave velocity $c_{r}$ corresponding to the velocity
obtained from (7) (black line) is of the same order of magnitude than that of
obtained from data (blue lines), with a similar trend as a function of
$U_{10}$. For very low winds, model and data trends disagree (see inset of
figure 3), but the resulting breaker phase speed is close to the capillary-
gravity transition, and hence also to the measurement limit of the wave-
breaking data.

This analysis shows that there is a correspondence between the properties of
the wave breaking field and those of the energy-containing eddies. This
positive result is used below to propose a three-dimensional interpretation of
the phenomenological model in the presence of roughness elements associated to
wave breaking.

##  4 A three-dimensional view of wave-breaking-constrained turbulence

Figure 5: Three-dimensional interpretation of the phenomenological model
(presented on figure 4) over a rough sea surface. (a) Breaker fronts of speed
$c$ and spanwise extension $\Lambda dc$ constrain the spanwise extension
$l_{\Lambda}$ of vortices advected at a speed $U=c$ (purple circles). The
streamwise scale of energy-containing eddies $s_{e}$ (red circle) results from
the spiraling motion of a fluid parcel (red line) due to advection at the mean
wind speed and rotation at a turnover velocity $u_{*}$. Over their liftetime
(longer than the eddie lifetime), breaker fronts cover an approximate
streamwise distance $g^{-1}c^{2}$, the dotted area. Note the difference in the
geometry of the individual breaker front and of the effective roughness
element (black line in the $(x,z)$ transect). (b) Required $l_{\Lambda}$ for
logarithmic-layer energy-containing eddies, i.e. for $s_{e}=z$ and $U\propto
u_{*}\ln(z/z_{0})$.

In the previous section, a roughness height was determined and was related to
the streamwise scale of energy-containing eddies $s_{e}$. Theoretical
(Csanady, 1985; Eifler, 1993), experimental (Reul et al., 1999) and numerical
(Suzuki et al., 2011) work further revealed that the interaction between
breaking waves and turbulence is a three-dimensional processes in which the
spanwise length of the breaking crests plays an essential role. In the
following, the phenomenological model of turbulence is extended so as to
relate this spanwise length to the streamwise scale of energy-containing
eddies. The model is sketched in figure 5a and is described below. The
spanwise direction is labelled as $y$.

Defining $l_{\Lambda}$ to be a spanwise extension of surface roughness
elements (figure 5a), $l_{\Lambda}$ is then proportional to the length of the
individual crests associated with breaking fronts of speed $c$, and the sum of
the latter, per unit area of ocean surface, is $\Lambda(c)dc$. Roughness
elements are assumed to trigger spanwise atmospheric eddies of similar
extension with a turnover velocity proportional to $u_{*}$. These eddies are
advected at the mean wind speed $U$ (compare the two purple circles in figure
5a). The turnover time of the eddies is $T_{\Lambda}\propto
l_{\Lambda}/u_{*}$, and the streamwise distance, $s_{e}$, traveled by a
representative fluid parcel during a turnover time $s_{e}=UT_{\lambda}$, is

$s_{e}(z)\propto U(z)l_{\Lambda}(z)/u_{*}.$ (8)

As sketched in figure 5a, the distance $s_{e}$ can be represented as the
streamwise Lagrangian distance separating the ascending and descending
branches of the eddy, in between which a representative fluid parcel would
undergo a spiraling motion. It is also the streamwise extension of the energy-
containing eddies defined in the phenomenological turbulence model (red circle
in figure 5a). As for the phenomenological model, this representation is to be
understood in an ensemble mean sense, unlike similar but instantaneous
representations of turbulence over flat walls (Kline et al., 1967; Blackwelder
& Eckelmann, 1979) and waves (Csanady, 1985; Eifler, 1993) that would result
from conditional sampling of the flow.

This three-dimensional interpretation relates, through (8), the streamwise
extension of the energy-containing eddy $s_{e}$ to the spanwise extension of
the roughness elements constraining its size $l_{\Lambda}$. For a logarithmic
layer, for which $s_{e}=z$ and $U\propto u_{*}\ln(z/z_{0})$, equation (8) then
yields the height-dependent $l_{\Lambda}(z)$ required to obtain attached
energy-containing eddies at each height $z$ satisfying the law of the wall. As
shown in figure 5b, after some elementary algebra, $l_{\Lambda}(z)$ has a
minimum, reached at a height corresponding to the roughness sublayer height
$h_{r}$ defined in (5). For a wavy surface, an increase of $l_{\Lambda}$ with
decreasing height is not physical (as happens for $z<h_{r}$). This would
indeed imply an increase of the spanwise extension of breaking fronts
(proportional to $l_{\Lambda}$) as height decreases. Noting that, for
decreasing heights, the speed (and wavelength) of the representative breakers
decreases (from (7)), this would imply an increase in the breaking-front
average length with decreasing wavelength, which is not realistic (see Belcher
& Vassilicos, 1997; Reul & Chapron, 2003, which, among others, highlight the
self-similarity of breaking fronts). Hence, breaking fronts can only imprint
their spanwise scale on energy-containing eddies for $z\geq h_{r}$.

Subsequently, $h_{r}$ represents the smallest height at which roughness
elements associated with breaking fronts can set the spanwise and streamwise
extension of energy-containing eddies. Below this height, the scale of energy-
containing eddies is $h_{r}$. Furthermore, $h_{r}$ is proportional to the
spanwise size of the roughness elements associated to breaking fronts
($l_{\Lambda}$). This interpretation offers an additional argument for the
importance of $h_{r}$ in the characterization of the wind-over-waves coupling.

From the wave-breaking data, the scale $h_{r}$ emerged from the variations of
$c^{3}\overline{\Lambda}(c,U_{10})dc$ with $U_{10}$. Assuming foam that
patches have a lifetime proportional to the underlying wave period (Reul &
Chapron, 2003), $g^{-1}c^{2}\Lambda(c)dc$ is related to the the fraction of
sea-surface turned over by breaking fronts, weighted by their lifetime. Hence
$c^{3}\Lambda(c)dc$ also contains information about both the lifetime of
roughness elements and their momentum (proportional to $c$). The fact that a
change in the properties of roughness-sublayer energy containing eddies is
coupled to a change in $c_{r}^{3}\Lambda(c_{r})dc$ thus indicates that not
only the size of the breaking fronts ($\propto l_{\Lambda}$) but also their
momentum ($\propto c_{r}$) and lifetime ($\propto c_{r}/g$) are needed for the
description of near-surface turbulence properties.

##  5 Concluding discussion

Figure 6: Relation between the two wind-dependent quantities proposed in this
work and the Beaufort scale (red text and curve, with speeds corresponding to
$U_{10}$ intervals): a representative speed of breaking fronts (horizontal
axis, converted to a wavelength for reference), and the height of the
roughness sublayer, equivalent to the smallest size of energy-containing
eddies (the vertical axis). Note that the classification of breaker types on
the horizontal axis (from Katsaros & Ataktürk, 1992) corresponds to the
qualitative sea-surface description related to the Beaufort scale (red text).

In this work, we identified a surface speed $c_{r}$ corresponding to small
scale breakers that are most coupled to the atmosphere. Infrared measurements
of wave breaking statistics have indeed revealed that variations of wave
breaking statistics ($c^{3}\overline{\Lambda}(c,U_{10})$) with ten-meter wind
are significant only for a particular breaker speed $c_{r}$, which increases
with friction velocity (equation (1)). Following the ideas of the Beaufort
scale, breakers travelling at this speed can be interpreted as being the most
visually related to a particular near-surface wind speed. This is shown in
figure 6, where the visual characteristics of the surface associated to the
Beaufort scale (red text) are compared to the breaker speed identified in this
analysis and its classification following Katsaros & Ataktürk (1992, italic
text in the figure).

We then extended a phenomenological turbulence model to interpret this
measurable relation as a macroscopic property of the near-surface turbulent
air flow. Without specifying any information on the wave spectrum, we found
that the breaker speed $c_{r}$ is associated with the sea-surface roughness
elements that affect the most energy production of near-surface turbulence
eddies. This creates a roughness sublayer in which the scale of near-surface
eddies is constrained by the breaker speed $c_{r}$ (§3), but also (as shown in
§4) by its spanwise length and lifetime (i.e. the dotted area in figure 5).
The height $h_{r}$ of the roughness soublayer is small, ranging from $10^{-4}$
m to $10^{-2}$m (see figure 6).

At first glance, the roughness sublayer height could be interpreted as (i) the
height at which airflow separation events induced by waves of phase speed
$c_{r}$ affect near-surface turbulence (Kraus, 1967; Melville, 1977), or (ii)
the height at which inviscid instabilities above a free surface are generated
(Miles, 1957). While both mechanisms certainly play a role in the generation
of roughness-sublayer eddies, they are associated with wind-wave generation,
which is not what we aim to describe through roughness-sublayer energy-
containing eddies. The latter should instead be viewed as wave-coherent
motions (as defined in Stewart, 1961), which are an indirect imprint of wind-
wave growth on atmospheric turbulence. The imprint is indirect since: (i)
wind-wave growth occurs due to the coupling between waves slower than $c_{r}$
and turbulence at heights greater than $h_{r}$ (Kudryavtsev et al., 2014; Ayet
et al., 2020a, see also appendix C) (ii) this coupling and wave growth then
sets the statistics of breaking waves (Phillips, 1985), (iii) this finally
affects roughness-sublayer eddies.

While the existence of a roughness sublayer over a windy sea has already been
suggested elsewhere (e.g. Kraus, 1967; Kitaigorodskii, 1973; Csanady, 1985),
the properties of roughness-sublayer eddies proposed herein are not supported
by empirical evidence, but based on an analogy with other types of surfaces
(Raupach, 1988; Bonetti et al., 2017). Importantly, energy-containing eddies
have been related to the spectral properties of turbulence in prior work
(Katul & Manes, 2014). This enables the new hypotheses on wind-wave
interactions formulated in this work to be tested using (i) direct numerical
simulations (Yang & Shen, 2009), which could leverage the spectral link
between energy-containing eddies and the turbulent pressure-strain
correlations (Ayet et al., 2020b) or (ii) in situ measurements (Ortiz-Suslow &
Wang, 2019), to directly analyze the shape of the turbulence spectrum (related
to energy-containing eddies, see Ayet & Katul, 2020).

Finally, the proposed breaker speed $c_{r}$ is a more general descriptor of
the wind-and-waves dynamical coupling than the roughness height $z_{0}$. The
speed $c_{r}$ can be determined from direct measurements without a priori
specification of turbulence or wave properties. We anticipate that $c_{r}$ is
sensitive to variations in parameters affecting the wind-and-waves coupling
(besides variations of $U_{10}$ already captured by $z_{0}$). Indeed, $c_{r}$
results from the selective attenuation of the amplitude of short waves by
surface drift (as discussed in §2.3). The amplitude of those waves is
generally described by a wave action budget balancing wind input, breaking
wave dissipation and non-linear wave-wave interactions (Phillips, 1985). On
the other hand, the size of roughness-sublayer energy-containing eddies is
also dependent on the presence of wave-induced motions (appendix C). Hence,
the dynamical wind-wave coupling expressed by the representative velocity
should be sensitive to modifications of these different processes by
environmental conditions, e.g. the presence of slicks, surface currents and
modulating longer waves that can alter surface drift (Phillips & Banner,
1974), the wave action budget (Kudryavtsev et al., 2005) and wave-induced
motions (Ayet et al., 2020a). Additionally, atmospheric stability effects, not
considered in the present work, can be integrated into the proposed
phenomenological framework, as modifying the scale of energy-containing eddies
(Katul et al., 2011; Li et al., 2012). Assessing the sensitivity of $c_{r}$ to
these parameters will be the focus of future work.

More generally, $c_{r}$ captures the effect of a multiscale moving ocean
surface on atmospheric turbulence and its modification due to external
parameters. It is a promising measurable candidate for the characterization of
air-sea interactions both by in situ or remote sensing methods. Measuring
$c_{r}$ is challenging and would imply capturing a characteristic velocity
(associated to a breaker scale) that changes depending on wind speed (and
probably on the external parameters mentioned above). As such, this can guide
methods for the interpretation of remote sensing observations, sensitive to
breakers and foam-coverage properties, as well as help design future
satellite-borne instruments, especially to directly retrieve sea surface
Doppler estimates.

## Declaration of Interests

The authors report no conflict of interest.

## Acknowledgments

The code used to generate the figures is available on Github [LINK WILL BE
PROVIDED UPON ACCEPTANCE]. AA was supported by DGA grant No D0456JE075, the
French Brittany Regional Council and ANR grant ANR-10-IEED-0006-26 . PS was
supported by funding from the European Research Council (ERC) under the
European Union’s Horizon 2020 research and innovation programme (grant
agreement No 805186). GGK was supported by the US National Science Foundation
(Grants NSF-AGS-164438, NSF-IOS-175489, and NSF-AGS-2028644).

##  Appendix A Sensitivity analysis to the choice of $n=3$ in $c^{n}\Lambda$

Figure 7: Width of the peak of $\partial(c^{n}\overline{\Lambda})/\partial
U_{10}$ for different values of the exponent $n$ and different bin sizes. The
width is defined as $w_{p}(U_{10},n)=|\log(c_{u})-\log(c_{d})|$, where $c_{u}$
and $c_{d}$ are respectively the smallest and largest breaker front speeds for
which $\partial(c^{n}\overline{\Lambda})/\partial U_{10}$ drops below 95% of
its maximal value. (a) Mean and (b) standard deviation of $w_{p}$ over all
$U_{10}$-bins.

Figure 2c shows that the gradient of $c^{3}\overline{\Lambda}(c,U_{10})$ with
respect to $U_{10}$ has a $U_{10}$-dependent peak associated with the
characteristic speed $c_{r}(U_{10})$. In §2, we also provided theoretical
arguments justifying the choice of the exponent $n=3$ when analysing
variations of $c^{n}\overline{\Lambda}(c,U_{10})$. To provide further support
for these arguments, sensitivity tests are performed on the exponent $n$. The
idea underlying the analysis of this appendix is to find the exponent for
which the peak in $\partial(c^{n}\overline{\Lambda})/\partial U_{10}$ is best
defined.

We first quantify the width around the peak of the curves shown in figure 2c
by computing $|\log(c_{u})-\log(c_{d})|$, where $c_{u}$ and $c_{d}$ are
respectively the smallest and largest breaker front speeds for which
$\partial(c^{n}\overline{\Lambda})/\partial U_{10}$ drops below 95% of its
peak value. The natural logarithm is used since, for different $U_{10}$
values, the width of the peak is comparable only in order of magnitudes (i.e.
in logarithmic coordinates, as shown in figure 2c). The $U_{10}$-averaged
value of this width, shown in figure 7a as a function of the exponent $n$,
reveal that, besides $n=3$, $n=0$ is also an exponent for which the peak is
narrow. The variance of this width with across all measured $U_{10}$ values
(figure 7b), confirms that this statement is statistically correct.

Figure 8: (a) $\partial\overline{\Lambda}/\partial U_{10}$ as a function of
the breaker speed for 3 m s-1 $U_{10}$ bins. This figure shows that, for a low
value of the exponent (here $n=0$), there exists a second peak in the gradient
of $c^{n}\overline{\Lambda}$ which is negative. (b) Difference between the
value of the negative peak and of the positive peak, normalized by the value
of the positive peak. This difference is averaged aver all $U_{10}$-bins (each
corresponding to one green/blue curve in (a)), and plotted as a function of
the exponent $n$. This figure shows that $n\sim 3$ is the exponent for which
the negative peak is the smallest with respect to the positive peak.

However, for $n<3$, $\partial(c^{n}\overline{\Lambda})/\partial U_{10}$
exhibits a second, negative peak, at lower values of $c$ (see figure 8a for an
example for $n=0$). This implies that there is no longer a unique
characteristic breaker speed that can be defined from the analysis. This
observation can be quantified by computing the difference between the absolute
values of the positive and the negative peaks, normalized to the value of the
positive peak (figure 8b). For low $n$, the difference is negative, implying
that the negative peak has a higher value than the positive peak. For $n$
close to $3$, the negative peak has a very small value with respect to the
positive peak (the solid blue curve in figure 8b has its maximum for $n~{}3$):
this supports the choice $n=3$ and discards the choice $n=0$.

Note that in figure 8b we did not show results from the data binned in
$U_{10}$-intervals of 1 m s-1. This is because, for this choice of bins, the
value of the positive peak is very weak for $n$ close to $0$, and hence the
relative difference becomes too large to be compared to the two curves shown
in Fig 8b. Besides, binning data in $U_{10}$-intervals of 1 m s-1 results in a
low number of points per bin size.

##  Appendix B Derivation of the phenomenological turbulence model

The derivation of (2) is offered starting from a phenomenological model of
wall-bounded turbulence (Gioia et al., 2010; Katul et al., 2011). The
calculations follow closely those of Ayet et al. (2020a).

We model an idealized SBL, defined as the lowest part of the atmospheric
boundary layer (adjacent to the surface) where the flow has high Reynolds
number and is horizontally homogeneous and stationary with no subsidence. All
averaged atmospheric quantities are invariant with respect to the streamwise
and spanwise directions, and depend only on height ($z$) from the surface. For
near-neutral atmospheric stability conditions ($|z/L|<0.1$, where $L$ is the
Obukhov length), the Turbulence Kinetic Energy (TKE) budget equation is then a
balance between mechanical production and dissipation

$u_{*}^{2}\frac{\partial U}{\partial z}=\epsilon,$ (9)

where $U$ is the mean wind speed (in the streamwise direction), $u_{*}$ is the
friction velocity, and $\epsilon$ is the TKE dissipation rate.

To obtain a link between the mean wind speed and friction velocity, a closure
for $\epsilon$ is needed. To this end, it was hypothesized (Gioia et al.,
2010; Katul et al., 2011) that turbulent structures of half streamwise and
vertical extension $s_{e}$ are attached to the surface (red circles in figures
4b,c). The turnover velocity $w_{e}$ of these structures was first related to
$\epsilon$ using the Kolmogorov 4/5 law for the third-order velocity structure
function (Monin & Yaglom, 1975),

$w_{e}(s_{e})\propto(\epsilon s_{e})^{1/3}.$ (10)

Then, the momentum flux $u_{*}^{2}$ through a given surface at a height $z$
(dashed lines in figures 4b,c) was expressed as a function of the vertical
transport of momentum by the turbulent structures, reading

$u_{*}^{2}(z)\propto w_{e}(s_{e})[U(z+s_{e})-U(z-s_{e})]\propto
w_{e}(s_{e})~{}\frac{\displaystyle\partial U}{\displaystyle\partial
z}~{}2s_{e}.$ (11)

The TKE dissipation was finally obtained by inserting (11) in (10). Using this
outcome in (9), the scaling $w_{e}\propto u_{*}$ and (2) are easily obtained.

##  Appendix C Wave-induced motions and their link with the phenomenological
model

Using open-ocean measurements of $u_{*}$ and $U_{10}$ and their COARE fit
(solid line in figure 1a), a wind-dependent $z_{0}$ can be inferred and a
corresponding $h_{r}$ can be computed from (5). As shown in figure 9a (black
solid line), $h_{r}$ ranges from $10^{-4}$ m to $10^{-2}$m. The log-log inset
further reveals that both $h_{r}$ (black solid line) and $z_{0}$ (dashed line)
follow a power law whose exponent is negative for winds lower than $3$ m s-1,
and tend towards $u_{*}^{2.5}$ for wind speeds larger than $8$ m s-1. These
two regimes are reminiscent of the dependence of $z_{0}$ with $u_{*}$, which
is usually parameterized as (Charnock, 1955)

$z_{0}=z_{0}^{v}+\zeta_{c}\frac{u_{*}^{2}}{g}.$ (12)

where $\zeta_{c}$ is the Charnock coefficient, which depends on $U_{10}$ or
$u_{*}$ (see e.g. Edson et al., 2013).

In wind-over-waves models (e.g. Kukulka et al., 2007; Kukulka & Hara, 2008a,
b; Kudryavtsev et al., 2014), waves generate _wave-induced motions_ in the SBL
due to: (i) airflow separation events associated to breaking waves (Reul et
al., 1999), and (ii) the contribution of waves with low steepness that deform
the overlying airflow (a process called wave-induced stress in the following,
see e.g. Plant, 1982; Belcher & Hunt, 1993). For a given wind-wave field, both
mechanisms generate stresses which decay exponentially with height, and
dashed-dotted and dotted lines in figure 9a correspond to their vertical
integral length scales divided by 100. When scaled, the integral length scales
of both mechanisms seems to have the same order of magnitude of $h_{r}$, but
do not exhibit the correct power-law exponent (as shown in the log-log inset).
Characterising as a whole the wind-over-wave coupling, $h_{r}$ does not model
any of those mechanisms, but corresponds to the height of an _effective
roughness element_ representing the bulk effect of individual wind-waves on
the energy-containing eddies at the bottom of the SBL (in the roughness
sublayer).

It is however well known that above the roughness sublayer, the stresses
caused by air-flow separation events and wave-induced stresses result in a
deviation of the wind profile from its logarithmic form (e.g. Miles, 1965;
Makin & Kudryavtsev, 1999; Hara & Belcher, 2004; Kudryavtsev et al., 2014).
This deviation results from TKE dissipation being modified by the presence of
wave-induced motions (see Ayet et al., 2020a), leading to a mean wind gradient
of the form

$\frac{dU}{dz}(z)=\frac{u_{*}[1-\alpha_{c}(z)]^{3/4}}{\kappa
s_{e}(z)},\;\;\mathrm{for}\;z\geq h_{r}$ (13)

where
$\alpha_{c}(z)=(u_{*}^{2}-\overline{u^{\prime}w^{\prime}}(z))/u_{*}^{2}<1$ is
the coupling coefficient that accounts for attenuation of the turbulent
momentum flux $\overline{u^{\prime}w^{\prime}}(z)$ as height decreases, due to
motions being increasingly coherent with the waves (i.e. the presence of
additional stresses, which induce a deviation of
$\overline{u^{\prime}w^{\prime}}$ from its value $u_{*}^{2}$ at the top of the
SBL).

In the roughness sublayer, as mentioned in §3, it is assumed that turbulent
motions are indistinguishable from wave-induced motions, and hence

$\frac{dU}{dz}(z)=\frac{u_{*}}{\kappa s_{e}},\;\;\mathrm{for}\;z\leq h_{r}.$
(14)

The deviation of the mean wind profile affects the wind speed at the bottom of
logarithmic layer and hence the properties of the roughness sublayer. Indeed,
the condition of continuity for the mean wind shear at the top of the
roughness sublayer leads to the following relation

$s_{e}=h_{r}[1-\alpha_{c}(h_{r})]^{-3/4}.$ (15)

As discussed in §3, this interesting relation reveals that the properties of
roughness-sublayer energy-containing eddies are coupled with wave-induced
motions aloft.

In general, the dependence of $\alpha_{c}$ with height and wind speed is non
trivial and analytical expressions of $h_{r}$ dcannot be derived. In figure
9a, red line, we show $h_{r}$ using a wind-dependent
$\alpha_{c}(h_{r},U_{10})$ from the wind-over-waves model of Kudryavtsev et
al. (2014) (dashed line in figure 9b). In this model, $\alpha_{c}$ includes
both airflow separation events and wave-induced stress to yield a momentum
flux and short wind-wave spectrum consistent with measurements and the COARE
parameterization. Figure 9a reveals that the power-law dependence of $h_{r}$,
discussed above, is affected by the presence of wave-induced motions in the
logarithmic sublayer (see the inset of figure 9a). However, the order of
magnitude of $h_{r}$ is similar to that when wave-induced motions are
neglected (compare black and red lines). When translated in terms of critical
wind speed $c_{r}=U(h_{r})$, the disagreement is similar (figure 9c). Finally,
in figure 9b, we further show the resulting ratio between $s_{e}$ and $h_{r}$
(solid line), discussed in §3.

Figure 9: Properties of the roughness sublayer in the presence of wave-induced
motions. (a) Different sublayer heights are compared: the roughness sublayer
height from this study without (solid black line) and with (solid red line)
wave-induced motions, and the roughness height $z_{0}$ using the Charnock
parameterization (dashed line). Also shown are the vertical integral length
scales of airflow separation (dotted-dashed line) and wave-induced stress
(dotted line) from the wind-over-wave model of Kudryavtsev et al. (2014). (b)
The ratio of the size of roughness-sublayer eddies and the height of the
roughness sublayer when wave-induced motions are included in the logarithmic
sublayer (solid line, left axis). This ratio and the red lines in (a) and (c)
have been computed using a coupling coefficient at the top of the roughness
sublayer (dashed line, right axis) from the wind-over-wave model of
Kudryavtsev et al. (2014). (c) Representative breaker speed $c_{r}$ with
(black line) and without (red line) wave-induced motions (corresponding the
roughness heights of (a)).

For completeness, a simplified case for which the roughness sublayer height
can be computed explicitly is illustrated. To that end, a simplified height
dependence of wave-induced stress is assumed following Makin et al. (1995),
and is given as

$\alpha_{c}(z)=\alpha_{c}^{0}\mathcal{H}(h_{w}-z),$ (16)

where $\mathcal{H}$ is the Heaviside step function, and $h_{w}$ is the height
at which the effect of waves on momentum flux becomes negligible. This
expression assumes that the coupling coefficient is constant and equal to
$\alpha_{c}^{0}$ below a height $h_{w}$, and zero above.

Equation (13) reveals that wave-coherent motions induce a deviation from the
logarithmic profile. Nevertheless, by extrapolating a logarithmic profile, a
roughness height $z_{0}$ can still be defined from the wind at an arbitrary
height $z$, $U(z)=(u_{*}/\kappa)\log(z/z_{0})$. Upon using this expression,
the simplified height dependence of $\alpha_{c}$, and integrating (13) from an
arbitrary height $z$ down to $h_{r}$, the mean wind speed at the roughness
sublayer height reads

$U(h_{r})=\frac{u_{*}}{\kappa}\log\left(\frac{h_{r}^{(1-\alpha^{0}_{c})^{3/4}}}{z_{0}}h_{w}^{1-(1-\alpha^{0}_{c})^{3/4}}\right)$
(17)

By further using (15), the mean wind speed and the roughness sublayer height
can be expressed as

$U(h_{r})=(u_{*}/\kappa)(1-\alpha_{c}(h_{r}))^{-3/4},\quad
h_{r}=z_{0}^{(1-\alpha^{0}_{c})^{-3/4}}h_{w}^{1-(1-\alpha^{0}_{c})^{-3/4}}e^{(1-\alpha_{c}^{0})^{-3/4}}.$
$None$

The above expressions show that as $\alpha_{c}$ increases, the roughness
sublayer height increases from $z_{0}\exp(1)$, and so does the mean wind speed
on top of the roughness sublayer. The limit $\alpha_{c}\to 1$ is out of the
range of the model, since the mean wind shear can no longer be continuous at
$h_{r}$ in this regime.

## References

  * Ayet et al. (2020a) Ayet, A, Chapron, Bertrand, Redelsperger, Jean-Luc, Lapeyre, Guillaume & Marié, L 2020a On the impact of long wind-waves on near-surface turbulence and momentum fluxes. Boundary-Layer Meteorology 174 (3), 465–491.
  * Ayet & Katul (2020) Ayet, A & Katul, GG 2020 Scaling laws for the length scale of energy-containing eddies in a sheared and thermally stratified atmospheric surface layer. Geophysical Research Letters 47 (23), e2020GL089997.
  * Ayet et al. (2020b) Ayet, A, Katul, GG, Bragg, AD & Redelsperger, JL 2020b Scalewise return to isotropy in stratified boundary layer flows. Journal of Geophysical Research: Atmospheres 125 (16), e2020JD032732.
  * Belcher & Hunt (1993) Belcher, SE & Hunt, JCR 1993 Turbulent shear flow over slowly moving waves. J Fluid Mech 251, 109–148.
  * Belcher & Vassilicos (1997) Belcher, SE & Vassilicos, JC 1997 Breaking waves and the equilibrium range of wind-wave spectra. J Fluid Mech 342, 377–401.
  * Blackwelder & Eckelmann (1979) Blackwelder, Ron F & Eckelmann, Helmut 1979 Streamwise vortices associated with the bursting phenomenon. J Fluid Mech 94 (3), 577–594.
  * Bonetti et al. (2017) Bonetti, Sara, Manoli, Gabriele, Manes, Costantino, Porporato, Amilcare & Katul, Gabriel G 2017 Manning’s formula and strickler’s scaling explained by a co-spectral budget model. J Fluid Mech 812, 1189–1212.
  * Charnock (1955) Charnock, H 1955 Wind stress on a water surface. Q J R Meteorol Soc 81 (350), 639–640.
  * Csanady (1985) Csanady, GT 1985 Air-sea momentum transfer by means of short-crested wavelets. J Phys Oceanogr 15 (11), 1486–1501.
  * Deike et al. (2017) Deike, Luc, Lenain, Luc & Melville, W Kendall 2017 Air entrainment by breaking waves. Geophysical Research Letters 44 (8), 3779–3787.
  * Deike & Melville (2018) Deike, Luc & Melville, W Kendall 2018 Gas transfer by breaking waves. Geophysical Research Letters 45 (19), 10–482.
  * Dickey et al. (2012) Dickey, T, Banner, ML, Bhandari, P, Boyd, T, Carvalho, L, Chang, G, Chao, Y, Czerski, H, Darecki, M, Dong, C & others 2012 Introduction to special section on recent advances in the study of optical variability in the near-surface and upper ocean. Journal of Geophysical Research: Oceans 117 (C7).
  * Edson et al. (2013) Edson, James B, Jampana, Venkata, Weller, Robert A, Bigorre, Sebastien P, Plueddemann, Albert J, Fairall, Christopher W, Miller, Scott D, Mahrt, Larry, Vickers, Dean & Hersbach, Hans 2013 On the exchange of momentum over the open ocean. J Phys Oceanogr 43 (8), 1589–1610.
  * Eifler (1993) Eifler, Walter 1993 A hypothesis on momentum and heat transfer near the sea-atmosphere interface and a related simple model. Journal of marine systems 4 (2-3), 133–153.
  * Gioia & Bombardelli (2001) Gioia, G & Bombardelli, FA 2001 Scaling and similarity in rough channel flows. Physical review letters 88 (1), 014501.
  * Gioia et al. (2010) Gioia, Gustavo, Guttenberg, Nicholas, Goldenfeld, Nigel & Chakraborty, Pinaki 2010 Spectral theory of the turbulent mean-velocity profile. Phys Rev Lett 105 (18), 184501.
  * Hara & Belcher (2004) Hara, Tetsu & Belcher, Stephen E 2004 Wind profile and drag coefficient over mature ocean surface wave spectra. J Phys Oceanogr 34 (11), 2345–2358.
  * Harman & Finnigan (2007) Harman, Ian N & Finnigan, John J 2007 A simple unified theory for flow in the canopy and roughness sublayer. Boundary-Layer Meteorol 123 (2), 339–363.
  * Hunt et al. (1988) Hunt, JCR, Leibovich, S & Richards, KJ 1988 Turbulent shear flows over low hills. Q J R Meteorol Soc 114 (484), 1435–1470.
  * Katsaros & Ataktürk (1992) Katsaros, Kristina B & Ataktürk, Serhad S 1992 Dependence of wave-breaking statistics on wind stress and wave development. In Breaking Waves, pp. 119–132. Springer.
  * Katul et al. (2011) Katul, Gabriel G, Konings, Alexandra G & Porporato, Amilcare 2011 Mean velocity profile in a sheared and thermally stratified atmospheric boundary layer. Phys Rev Lett 107 (26), 268502\.
  * Katul & Manes (2014) Katul, Gabriel G & Manes, Costantino 2014 Cospectral budget of turbulence explains the bulk properties of smooth pipe flow. Physical Review E 90 (6), 063008.
  * Kitaigorodskii (1973) Kitaigorodskii, Sergei Aleksandrovich 1973 The physics of air-sea interaction. Israel Program for Scientific Translations Jerusalem.
  * Kleiss & Melville (2010) Kleiss, Jessica M. & Melville, W. Kendall 2010 Observations of wave breaking kinematics in fetch-limited seas. J. Phys. Oceanogr. 40, 2575–2604.
  * Kline et al. (1967) Kline, Stephen J, Reynolds, William C, Schraub, FA & Runstadler, PW 1967 The structure of turbulent boundary layers. J Fluid Mech 30 (4), 741–773.
  * Kraus (1967) Kraus, EB 1967 Wind stress along the sea surface. Advances in Geophysics 12, 213–255.
  * Kudryavtsev et al. (2005) Kudryavtsev, Vladimir, Akimov, D, Johannessen, Johnny & Chapron, Bertrand 2005 On radar imaging of current features: 1\. model and comparison with observations. Journal of Geophysical Research: Oceans 110 (C7).
  * Kudryavtsev et al. (2014) Kudryavtsev, V., Chapron, Bertrand & Makin, V. 2014 Impact of wind waves on the air-sea fluxes: A coupled model. J Geophys Res Oceans 46 (2), 1022–1037.
  * Kukulka & Hara (2008a) Kukulka, Tobias & Hara, Tetsu 2008a The effect of breaking waves on a coupled model of wind and ocean surface waves. part i: Mature seas. Journal of physical oceanography 38 (10), 2145–2163.
  * Kukulka & Hara (2008b) Kukulka, Tobias & Hara, Tetsu 2008b The effect of breaking waves on a coupled model of wind and ocean surface waves. part ii: Growing seas. J Phys Oceanogr 38 (10), 2164–2184.
  * Kukulka et al. (2007) Kukulka, Tobias, Hara, Tetsu & Belcher, Stephen E 2007 A model of the air–sea momentum flux and breaking-wave distribution for strongly forced wind waves. J Phys Oceanogr 37 (7), 1811–1828.
  * Li et al. (2012) Li, Dan, Katul, Gabriel G & Bou-Zeid, Elie 2012 Mean velocity and temperature profiles in a sheared diabatic turbulent boundary layer. Phys Fluids 24 (10), 105105.
  * Makin & Kudryavtsev (1999) Makin, VK & Kudryavtsev, VN 1999 Coupled sea surface-atmosphere model: 1. wind over waves coupling. J Geophys Res Oceans 104 (C4), 7613–7623.
  * Makin et al. (1995) Makin, VK, Kudryavtsev, VN & Mastenbroek, C 1995 Drag of the sea surface. Boundary-Layer Meteorol 73 (1), 159–182.
  * Manes et al. (2007) Manes, Costantino, Pokrajac, Dubravka & McEwan, Ian 2007 Double-averaged open-channel flows with small relative submergence. Journal of Hydraulic Engineering 133 (8), 896–904.
  * Melville & Rapp (1985) Melville, WK & Rapp, Ronald J 1985 Momentum flux in breaking waves. Nature 317 (6037), 514–516.
  * Melville (1977) Melville, W K 1977 Wind stress and roughness length over breaking waves. Journal of Physical Oceanography 7 (5), 702–710.
  * Miles (1957) Miles, John W 1957 On the generation of surface waves by shear flows. Journal of Fluid Mechanics 3 (2), 185–204.
  * Miles (1965) Miles, John W 1965 A note on the interaction between surface waves and wind profiles. Journal of Fluid Mechanics 22 (4), 823–827.
  * Monin & Yaglom (1975) Monin, A. & Yaglom, A.M. 1975 Statistical Fluid Mechanics. MIT press, Cambridge, MA.
  * Munk (1955) Munk, WH 1955 Wind stress on water: An hypothesis. Quarterly Journal of the Royal Meteorological Society 81 (349), 320–332.
  * Newell & Zakharov (1992) Newell, Alan C & Zakharov, Vladimir E 1992 Rough sea foam. Physical Review Letters 69 (8), 1149.
  * Ortiz-Suslow & Wang (2019) Ortiz-Suslow, David G & Wang, Qing 2019 An evaluation of kolmogorov’s-5/3 power law observed within the turbulent airflow above the ocean. Geophys Res Lett .
  * Phillips (1985) Phillips, OM 1985 Spectral and statistical properties of the equilibrium range in wind-generated gravity waves. Journal of Fluid Mechanics 156, 505–531.
  * Phillips & Banner (1974) Phillips, OM & Banner, ML 1974 Wave breaking in the presence of wind drift and swell. Journal of Fluid Mechanics 66 (4), 625–640.
  * Phillips (1957) Phillips, Owen M 1957 On the generation of waves by turbulent wind. Journal of fluid mechanics 2 (5), 417–445.
  * Phillips (1977) Phillips, Owen M 1977 The dynamics of the upper ocean. Cambridge University Press.
  * Plant (1982) Plant, William J. 1982 A relationship between wind stress and wave slope. J Geophys Res Oceans 87 (C3), 1961–1967.
  * Raupach (1988) Raupach, MR 1988 Canopy transport processes. In Flow and transport in the natural environment: advances and applications, pp. 95–127. Springer.
  * Raupach et al. (1996) Raupach, M_R, Finnigan, JJ & Brunet, Y 1996 Coherent eddies and turbulence in vegetation canopies: the mixing-layer analogy. In Boundary-layer meteorology 25th anniversary volume, 1970–1995, pp. 351–382. Springer.
  * Raupach et al. (1991) Raupach, M. R., Antonia, R. A. & Rajagopalan, S. 1991 Rough-Wall Turbulent Boundary Layers. Applied Mechanics Reviews 44 (1), 1–25.
  * Reul et al. (1999) Reul, Nicolas, Branger, Hubert & Giovanangeli, J-P 1999 Air flow separation over unsteady breaking waves. Phys Fluids 11 (7), 1959–1961.
  * Reul & Chapron (2003) Reul, Nicolas & Chapron, Bertrand 2003 A model of sea-foam thickness distribution for passive microwave remote sensing applications. Journal of Geophysical Research: Oceans 108 (C10).
  * Stewart (1961) Stewart, RW 1961 The wave drag of wind over water. Journal of fluid mechanics 10 (2), 189–194.
  * Sutherland & Melville (2013) Sutherland, Peter & Melville, W Kendall 2013 Field measurements and scaling of ocean surface wave-breaking statistics. Geophys Res Lett 40 (12), 3074–3079.
  * Sutherland & Melville (2015) Sutherland, Peter & Melville, W Kendall 2015 Field measurements of surface and near-surface turbulence in the presence of breaking waves. Journal of Physical Oceanography 45 (4), 943–965.
  * Suzuki et al. (2011) Suzuki, Nobuhiro, Hara, Tetsu & Sullivan, Peter P 2011 Turbulent airflow at young sea states with frequent wave breaking events: Large-eddy simulation. Journal of the atmospheric sciences 68 (6), 1290–1305.
  * Townsend (1976) Townsend, Albert A 1976 The structure of turbulent shear flow, Vol. 2. Cambridge university press.
  * Wu (1975) Wu, Jin 1975 Wind-induced drift currents. Journal of Fluid Mechanics 68 (1), 49–70.
  * Yang & Shen (2009) Yang, Di & Shen, Lian 2009 Characteristics of coherent vortical structures in turbulent flows over progressive surface waves. Phys Fluids 21 (12), 125106.