1 Introduction
Flows in geophysical fluids span a remarkably broad range of scales. From a practical point of view, the interest concentrates on the large scales,
$O(10^{3}{-}10^{5}~\text{m})$
, most relevant for our experience of, e.g., the climate, the concentration of pollutants and the ocean productivity. Geophysical flows are, however, highly nonlinear, implying that the behaviour of the large scales is not independent of the small ones,
$O(10^{-3}{-}10^{2}~\text{m})$
. This is one of the main motivations underlying the study of small-scale turbulence in geophysical flows.
In this work, the vertical turbulent transport of a scalar, namely temperature, is investigated, using observations collected above a sloping sea floor. The specific interest for benthic boundary layers, sloping in particular, is connected to the suggestion (Munk Reference Munk1966) that the regions directly above the sea floor dominate the globally averaged vertical buoyancy transport in the ocean. Here, we address this issue from the point of view of the relation between the flux of temperature (
${\it\phi}_{{\it\theta}}$
), the gradient of temperature and the distance above the solid boundary. The temperature flux is related to the heat flux by a factor
${\it\rho}c$
, with
${\it\rho}$
the density and
$c$
the specific heat. By considering the temperature instead of the heat flux, we assume that the density and specific heat are approximately constant.
The non-dimensional equilibrium flux
$f_{0}$
(from (2.8) and (2.7)) as a function of the non-dimensional gradient
$g_{0}$
is shown using black lines for two values of the non-dimensional height above the bottom
$z$
(see legend). The loci
$(f_{max},g_{max})$
and
$(f_{min},g_{min})$
are also shown in red and blue respectively. The results refer to the parameter values of table 1.
The flux–gradient relation has received attention in particular because it may explain the intermittency of density gradients in stratified turbulence, i.e. the presence of layers. Layering is ubiquitous in the observations we present here, as in most observations from geophysical flows. Phillips (Reference Phillips1972) and Posmentier (Reference Posmentier1977) hypothesised that the time evolution of buoyancy in a stratified flow could be modelled, assuming horizontal homogeneity, by a diffusion equation of the form
$$\begin{eqnarray}b_{t}=({\it\phi}(b_{z}))_{z},\end{eqnarray}$$
where the subscripts indicate the derivation,
$b$
is the buoyancy,
${\it\phi}(b_{z})$
is the buoyancy flux, seen as a function of the buoyancy gradient,
$t$
is time and
$z$
is the vertical coordinate. If the flux function
${\it\phi}(b_{z})$
is non-monotonic, an instability can occur for
$\text{d}{\it\phi}(b_{z})/\text{d}b_{z}<0$
, leading to the unbounded growth of sharply stratified sheets divided by weakly stratified layers. This mechanism for layer formation was later recognised in the laboratory by Ruddick, McDougall & Turner (Reference Ruddick, Mcdougall and Turner1989) and Park, Whitehead & Gnanadeskian (Reference Park, Whitehead and Gnanadeskian1994).
Building upon this idea, Balmforth, Llewellyn Smith & Young (Reference Balmforth, Llewellyn Smith and Young1998) (hereafter B98) proposed a model for the evolution of a stratified turbulent flow, again assuming horizontal homogeneity, in terms of two partial differential equations for the turbulent kinetic energy density (
$e$
) and the buoyancy gradient (
$g=b_{z}$
). Despite its simplicity, the model is remarkably successful in describing the formation and evolution of layers in laboratory experiments. The fundamental advance in B98 was the recognition that
$\text{d}{\it\phi}(g)/\text{d}g$
must be greater than zero for large
$g$
, in order to stop the growth of the sharply stratified layers before they grow infinitely thin and strong, i.e. the
${\it\phi}(g)$
function must be ‘N-shaped’ (see, e.g., figure 1, black lines). B98 argued that this implies a forcing term in the kinetic energy equation proportional to a velocity scale (squared) divided by the eddy turnover time scale, leading to the correct N-shape of the equilibrium
${\it\phi}(g)$
.
Holford & Linden (Reference Holford and Linden1999) further tested the model predictions in the laboratory, suggesting that the N-shape of
${\it\phi}(g)$
was the consequence rather than the cause of layer formation. Martin & Rehmann (Reference Martin and Rehmann2006) also addressed the problem of layer formation in the laboratory, observing the development of an N-shaped
${\it\phi}(g)$
after layers formed. Wunsch & Kerstein (Reference Wunsch and Kerstein2001) reformulated the model in stochastic terms, also suggesting that the dynamics predicted by the model could be related to the zigzag instability (Billant & Chomaz Reference Billant and Chomaz2000).
Another important hypothesis in B98 is that
$e$
and
$b$
are transported by the turbulent flow according to a mixing-length argument, with the flux proportional to the gradient and to
$l\,U$
, where
$l$
is a length scale and
$U$
is a velocity scale. In B98,
$l$
was assumed to be the geometric average between a length scale proportional to the Ozmidov scale and a constant length scale
$d$
, set by the stirring device:
$$\begin{eqnarray}l=\frac{de^{1/2}}{(e+{\it\gamma}d^{2}g)^{1/2}},\end{eqnarray}$$
with
${\it\gamma}$
a non-dimensional constant. The B98 model implies that the mixing length, relevant for computing the turbulent flux or turbulent diffusivity, is set by the Ozmidov scale under strong stratification, while it is set by
$d$
under weak stratification. The constant
${\it\gamma}$
controls the range where
$l$
is influenced by the stirring device.
In this work, we will reinterpret the B98 model to describe the dependence of the buoyancy flux on the height above the bottom (HAB) and on the buoyancy gradient. The assumption of a mixing length proportional to the distance from the boundary is the approach in the classical problem of the ‘law-of-the-wall’ for a constant-density fluid, see, e.g., Landau & Lifshitz (Reference Landau and Lifshitz1987, Ch. 5) or the original discussion in Prandtl (Reference Prandtl and Durand1935, § 21). As in B98, we assume horizontal homogeneity but include the effect of the solid sea floor, in a highly idealised way which also ignores the slope of the boundary (the boundary normal and vertical directions thus coincide). In doing so, we reinterpret
$d$
in (1.2) as a different length scale, namely the HAB
$z$
, which is the external length scale limiting the mixing length close enough to the solid boundary. Differently from
$d$
,
$z$
is not a constant.
The predictions of the model will be compared with estimates of the flux obtained from high-resolution temperature observations above a sloping sea floor. The estimates are derived using the method suggested in Winters & D’Asaro (Reference Winters and D’Asaro1996). This method can be seen as an extension of the more common one by Osborn & Cox (Reference Osborn and Cox1972), and is based on the estimation of the small-scale buoyancy gradient, and on the adiabatic reordering of buoyancy profiles to remove inversions. The method proposed in Winters & D’Asaro (Reference Winters and D’Asaro1996) has three main advantages over the one by Osborn & Cox (Reference Osborn and Cox1972): (i) it provides a clear definition of the background stratification to be used for estimating the flux, (ii) it provides an estimate of the flux due to molecular mixing (irreversible) rather than due to advection (reversible) and (iii) it provides a natural framework for studying local variations of the flux as a function of the background buoyancy gradient and HAB.
The flux estimates will also be compared with those obtained using the classical method of reordering buoyancy overturns, first outlined by Thorpe (Reference Thorpe1977). This method has been used in several observational studies (see, e.g., Dillon Reference Dillon1982; Itsweire Reference Itsweire1984; Gregg Reference Gregg1987; Gargett Reference Gargett1989). Despite recent criticism (Mater & Venayagamoorthy Reference Mater and Venayagamoorthy2014; Chalamalla & Sarkar Reference Chalamalla and Sarkar2015), it can still provide a reliable order-of-magnitude estimate of turbulence parameters, after sufficient averaging (Chalamalla & Sarkar Reference Chalamalla and Sarkar2015; Mater et al. Reference Mater, Venayagamoorthy, Laurent and Moum2015). We also stress that the simulations of Chalamalla & Sarkar (Reference Chalamalla and Sarkar2015) are not directly relevant for the observations presented here, where a single overturn evolving from a laminar flow is never observed. In the present work, we will show that the flux estimates based on overturn reordering are correlated with those from the method of Winters & D’Asaro (Reference Winters and D’Asaro1996), despite a systematic low bias of the latter due to instrumental resolution limitations. Furthermore, we will present a mixing-length argument supporting the use of the estimates from the reordering of overturns in the vicinity of a solid boundary, and provide observational evidence in favour of this argument.
The article is organised as follows. In § 2, the extension of the model by B98 is briefly presented. In § 3.1, the temperature observations are briefly described. In §§ 3.2 and 3.3, the two methods used for estimating the flux are presented in some detail. In § 4, the flux estimates are presented and compared with the model. Discussion of the results and conclusions follow in § 5.
2 Theory
We develop an extension of the B98 model with the aim of applying it to stratified turbulence above a sea floor. While the most notable feature of the B98 model is layer development, our interest in it is also motivated by its focus on the flux–gradient relation, which we can estimate in the observations. Following the example of B98, the aim is to develop a model that captures the essential phenomenology in the simplest mathematical form. The starting point is equation (2.1a,b) of B98, which we rewrite here using a slightly modified notation:
$$\begin{eqnarray}\displaystyle & \displaystyle g_{t}=(le^{1/2}g)_{zz}, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle e_{t}={\it\beta}(le^{1/2}e_{z})_{z}-le^{1/2}g-{\it\alpha}l^{-1}e^{3/2}+{\mathcal{P}}. & \displaystyle\end{eqnarray}$$
$z$
is zero at the bottom and increases upwards,
${\it\alpha}$
and
${\it\beta}$
are dimensionless constants, and
${\mathcal{P}}$
is the energy production term. Equation (2.1) describes the coupled evolution of the buoyancy gradient and the kinetic energy in a stratified turbulent flow, realistically enough to reproduce all of the essential features of layer formation and evolution in the laboratory experiments discussed in B98. Horizontal homogeneity is assumed; the slope of the sea floor is thus ignored. The possible impact of this far-fetched assumption is discussed in § 5. Equation (2.1) is considered here only for the highly idealised case of the steady state (
$\partial /\partial t\equiv 0$
).
We depart from B98 in the definitions of the mixing length
$l$
and energy production term
${\mathcal{P}}$
. Along the same lines as B98, the mixing length is defined as
$$\begin{eqnarray}l={\it\delta}\frac{ze^{1/2}}{(e+{\it\gamma}z^{2}g)^{1/2}},\end{eqnarray}$$
with
${\it\delta}$
a non-dimensional constant. While
$l$
has the same form as (1.2), the variable
$z$
is used instead of the constant
$d$
. As suggested already by Thorpe (Reference Thorpe1977) for the case of the surface boundary layer, at small
$z$
the density overturns are limited by the presence of the boundary. Here, we argue that the mixing length is similarly limited by the presence of the sea floor, in analogy to the case of a wall-bounded boundary layer in a constant-density fluid (Landau & Lifshitz Reference Landau and Lifshitz1987). The geometric average in (2.2) is chosen for reasons of mathematical simplicity. We do not imply that this is the correct mixing length, but that (2.2) captures the essential characteristics of the system, namely that from a far-field value
$l\propto (e/g)^{1/2}$
, the mixing length interpolates to
$l\propto z$
close to the boundary. Identification of the exact dependence of
$l$
on
$z$
and
$g$
is beyond the scope of this work. Since the observations suggest that the influence of the sea floor extends throughout the mooring (100 m), the constant
${\it\gamma}$
is expected to be much smaller than 1, implying that the asymptotic value of
$l$
is approached only above the top of the mooring. For any realistic values of
$e$
and
$g$
, the constant
${\it\delta}$
is also expected to be smaller than
$1$
(the values used here are listed in table 1).
Values of the dimensional and non-dimensional parameters used in the numerical solution of the model.
In B98, the energy production term was defined as
${\mathcal{P}}={\it\alpha}U^{2}e^{1/2}l^{-1}$
, with
$U$
a constant velocity scale. It was shown that
${\mathcal{P}}\propto e^{1/2}l^{-1}$
is essential for having an N-shaped flux–gradient relation. In other words,
${\mathcal{P}}$
being inversely proportional to the eddy turnover time is essential for the model to produce layers of finite thickness. Rather than by a stirring device moving at a constant velocity proportional to
$U$
, here the forcing is due to the interaction of the internal wave field (mainly at the tidal frequency) with the sea floor. We hypothesise that the velocity profile in the region observed by the mooring (
$z=5{-}105~\text{m}$
) is approximated by the classical log layer for constant-density fluids (Landau & Lifshitz Reference Landau and Lifshitz1987). The hypothesis is not unreasonable, since the flow is highly turbulent (the Reynolds number is estimated to be at least
$10^{4}$
), but it neglects the effects of density stratification, internal wave reflection and rotation. According to this hypothesis, the production term can be written as
$$\begin{eqnarray}{\mathcal{P}}={\it\alpha}e^{1/2}l^{-1}u_{\ast }^{2}\log ^{2}\left(\frac{u_{\ast }}{{\it\nu}}z\right),\end{eqnarray}$$
with
$u_{\ast }$
the shear velocity, assumed constant, and
${\it\nu}$
the molecular viscosity. It should be noted that, as in the classical theory, the factor
$u_{\ast }/{\it\nu}$
in the argument of the logarithm is large if
$u_{\ast }$
is of the order of the velocity fluctuations
$u_{\ast }~(10^{-3}~\text{m}~\text{s}^{-1})$
. Our observations of velocity have insufficient resolution to evaluate how well a log layer approximates the real velocity profile. However, the small decrease of the mean velocity approaching the sea floor (not shown) is not incompatible with (2.3).
Non-dimensionalisation of the problem is carried out by taking
$$\begin{eqnarray}z={\it\gamma}^{-1/2}u_{\ast }N_{0}^{-1}\hat{z},\quad g=N_{0}^{2}{\hat{g}},\quad e=u_{\ast }^{2}\hat{e}\quad \text{and}\quad l={\it\delta}{\it\gamma}^{-1/2}u_{\ast }N_{0}^{-1}\hat{l},\end{eqnarray}$$
with hats indicating the non-dimensional quantities.
Combining the equations of the system (2.1) in steady state (
$\partial /\partial t\equiv 0$
) with (2.2) and (2.3), and applying the non-dimensionalisation, we obtain
$$\begin{eqnarray}\displaystyle & \displaystyle \left[\frac{zeg}{(e+z^{2}g)^{1/2}}\right]_{zz}=0, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle {\it\beta}r{\it\gamma}\left[\frac{zee_{z}}{(e+z^{2}g)^{1/2}}\right]_{z}-r\frac{zeg}{(e+z^{2}g)^{1/2}}+\frac{(e+z^{2}g)^{1/2}}{z}(\log ^{2}({\it\xi}z)-e)=0, & \displaystyle\end{eqnarray}$$
$r={\it\delta}^{2}/({\it\alpha}{\it\gamma})$
and
${\it\xi}=u_{\ast }^{2}/(N_{0}{\it\gamma}^{1/2}{\it\nu})$
. We could not derive an analytical solution of the system (2.5). Following B98, we rather seek equilibrium solutions of (2.5b
) neglecting the first term on the left-hand side. We thus solve an algebraic equation rather than a system of ordinary differential equations, by assuming that the first term on the left-hand side of (2.5b
) is much smaller than 1 (i.e.
${\it\beta}r{\it\gamma}\ll 1$
). From a physical point of view, this means that we neglect the vertical transport of kinetic energy in comparison to that of buoyancy. This is a somewhat artificial assumption, used in order to make analytical progress, and which should be considered only as a first step in the analysis of (2.5). The investigation of the complete system (2.5) is the topic of further investigation, going beyond the scope of this work.
We thus solve the following equation for
$e$
:
$$\begin{eqnarray}rz^{2}eg-(e+z^{2}g)(\log ^{2}({\it\xi}z)-e)=0,\end{eqnarray}$$
obtaining the (positive physically meaningful) solution for the relationship between the equilibrium values
$(g_{0},e_{0})$
:
$$\begin{eqnarray}e_{0}={\textstyle \frac{1}{2}}\{\log ^{2}({\it\xi}z)-z^{2}g_{0}(1+r)+[(\log ^{2}({\it\xi}z)-z^{2}g_{0}(1+r))^{2}+4z^{2}g_{0}\log ^{2}({\it\xi}z)]^{1/2}\}.\end{eqnarray}$$
This expression is very similar to (3.1) of B98, with additional terms
$z^{2}$
and
$\log ^{2}({\it\xi}z)$
due to our different definitions of
$l$
and
${\mathcal{P}}$
. In the unstratified case,
$g_{0}=0$
and
$e_{0}=\log ^{2}({\it\xi}z)$
, with the kinetic energy profile completely determined by the forcing.
Dependence of the non-dimensional equilibrium flux
$f_{0}$
on
$z$
and stratification
$g_{0}$
in the model. The flux is shown with grey contours in logarithmic scale. The red (blue) dashed line marks the position of the local maximum (local minimum) of the flux predicted by the model (
$\propto \!\log ^{2}({\it\xi}z)/z^{2}$
). The results refer to the parameter values of table 1.
The equilibrium value of the flux can be computed as
$$\begin{eqnarray}f_{0}=le_{0}^{1/2}g_{0}=\frac{ze_{0}g_{0}}{(e_{0}+z^{2}g_{0})^{1/2}}.\end{eqnarray}$$
Stationary points of the flux
$f_{0}$
are found by solving
$\partial f_{0}/\partial g_{0}=0$
:
$$\begin{eqnarray}\displaystyle & \displaystyle g_{max}=\frac{2}{3(1+r)^{2}}[2(r-1)-(1-14r+r^{2})^{1/2}]\frac{\log ^{2}({\it\xi}z)}{z^{2}}, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle g_{min}=\frac{2}{3(1+r)^{2}}[2(r-1)+(1-14r+r^{2})^{1/2}]\frac{\log ^{2}({\it\xi}z)}{z^{2}}. & \displaystyle\end{eqnarray}$$
$f_{max}$
and
$f_{min}$
respectively, are complicated and unenlightening, and for this reason are only evaluated numerically. As in B98, these stationary points of the flux are present only above a critical value,
$r_{c}=7+4\sqrt{3}$
. Since we are interested in this model because of the presence of the stationary points in the flux–gradient relation, we focus our analysis on the case of
$r$
well above
$r_{c}$
. Different values of
$r$
above the critical value give qualitatively similar results. The non-dimensional flux as a function of the gradient is shown in figure 1 for two different values of
$z$
and the parameters of table 1.
Figure 2 shows, in non-dimensional units, the flux
$f_{0}$
as a function of both
$g_{0}$
and
$z$
. From figures 1 and 2, the main predictions of the model can be summarised:
-
(i) the stationary points of the flux (
$f_{max}$
and
$f_{min}$
) are shifted to smaller values of the gradient (
$g_{max}$
and
$g_{min}$
) farther from the bottom; -
(ii) the values of the gradient at the stationary points,
$g_{max}$
and
$g_{min}$
, decrease with
$z$
as
$\log ^{2}({\it\xi}z)/z^{2}$
; -
(iii) the value of
$f_{max}$
increases for decreasing
$z$
; -
(iv) the flux tends to
$f_{0}=e_{0}g_{0}^{1/2}$
for
$z^{2}g_{0}\gg e_{0}$
(this asymptotic behaviour is the same as in the B98 model).
With these predictions in mind, we now move to the estimation of the flux from observations.
Description of the mooring from which the data are obtained.
3 Methods
3.1 Mooring description and context
The data used in this work were collected using a mooring deployed between spring and late summer 2013 on the slopes of Seamount Josephine (North-Eastern Atlantic Ocean, see table 2). The mooring was made of a 205 m long nylon-coated steel cable, attached approximately 5 m above a 500 kg ballast weight at the sea floor, and to an elliptical buoy at the top. A Teledyne/RDI 75 kHz acoustic Doppler current profiler (ADCP) was mounted in the top buoy of the mooring, 210 m above the sea floor. The pressure and tilt sensors of the ADCP recorded movements smaller than 0.2 m in the vertical direction and 3 m in the horizontal direction.
On the cable, 144 ‘NIOZ4’ thermistors were taped at intervals of 0.7 m. The instruments are an evolution of the model described in van Haren et al. (Reference van Haren, Laan, Buijsman, Gostiaux, Smit and Keijzer2009). The sampling rate of the thermistors was 1 Hz and the precision of the temperature estimate was better than 1 mK, with a noise level of approximately
$5\times 10^{-5}~\text{K}$
. Each thermistor recorded a time series approximately
$10^{7}~\text{s}$
long.
The top-mounted ADCP was pointed downwards, recording the horizontal velocity vector and echo intensity every 5 m vertically, between 15 and 185 m above the sea floor. The sampling interval of the ADCP was 600 s.
The calibration information for the raw temperature data was obtained from a calibration bath developed and operated at the Royal Netherlands Institute for Sea Research. Compensation was made in the temperature data for the drift in the response of the thermistor electronics, which can be visible in some instruments over periods longer than a few weeks.
At the mooring location, semidiurnal tidal motions are dominant in both the velocity and temperature signals (Cimatoribus, van Haren & Gostiaux Reference Cimatoribus, van Haren and Gostiaux2014; Cimatoribus & van Haren Reference Cimatoribus and van Haren2015). We divided the data into two subsets following Cimatoribus & van Haren (Reference Cimatoribus and van Haren2015), based on the phase of the tide. The cooling and warming phases are defined by the sign of the time derivative of the vertically averaged temperature signal (negative and positive respectively), band-pass filtered at the semidiurnal frequency. Velocities are mostly aligned to the isobath direction. During the cooling phase, the cross-isobath component of velocity is preferentially upslope; vice versa for the warming phase. Mean velocities are larger during the cooling phase. The vertical shear detected by the ADCP is small (due to limited vertical resolution).
The temperature profiles are often organised in layers, with weakly stratified regions alternated with thinner sheets having strong stratification (Cimatoribus & van Haren Reference Cimatoribus and van Haren2015; van Haren, Cimatoribus & Gostiaux Reference van Haren, Cimatoribus and Gostiaux2015). During the cooling phase, particularly strong stratification is often observed at the top of a cold temperature layer, more weakly stratified, extending to the bottom of the mooring (5 m above the sea floor). See figure 3 for an example of the temperature data recorded at the mooring, in particular (a), showing the ubiquitous presence of layers.
(a) Five vertical temperature profiles, randomly chosen in the data set, reordered to remove any inversion. (b) A contour plot with an example of the temperature data recorded by the mooring during three days, from 10 to 13 June 2013.
3.2 Thorpe scale flux estimation method
At the time of deployment and recovery, two shipborne conductivity–temperature–depth (CTD) surveys were performed, for comparison with the results from the thermistors and to check the temperature–density relationship. As discussed in more detail in Cimatoribus et al. (Reference Cimatoribus, van Haren and Gostiaux2014), the low noise level of the temperature data, in combination with the tight and linear temperature–density relationship observed in the CTD data, supports the use of temperature as a proxy for buoyancy. This enables the estimation of overturning scales, providing in turn an estimate of turbulent diffusivity.
The Thorpe scale (
$l_{T}$
) is the root mean square (r.m.s.) of the displacements necessary to adiabatically reorder a buoyancy profile, to remove all buoyancy inversions (Thorpe Reference Thorpe1977). It is a measure of the overturning caused by turbulent eddies, and can be related to the Ozmidov scale
$l_{O}$
, the scale at which buoyancy and inertial forces balance in a stratified flow (Thorpe Reference Thorpe1977; Dillon Reference Dillon1982). Using the empirical relation of Dillon (Reference Dillon1982), the Thorpe scale can provide an estimate of vertical diffusivity and flux:
$$\begin{eqnarray}K_{{\it\theta}}^{T}=0.128l_{T}^{2}N,\end{eqnarray}$$
where
${\it\theta}$
is the potential temperature,
$N$
is the buoyancy frequency and
$0.128$
is an empirical constant. The superscript is a reminder of the method used for obtaining the estimate. Multiplication of the diffusivity from (3.1) with the background temperature gradient provides a Thorpe-scale-based estimate of the flux, which will be indicated by
${\it\phi}_{{\it\theta}}^{T}$
. The flux estimate
${\it\phi}_{{\it\theta}}^{T}$
is the diffusive one, rather than the total, which includes the advective one (see, e.g., Scotti Reference Scotti2015). While experiments show that the advective flux may dominate the total (Martin & Rehmann Reference Martin and Rehmann2006), the irreversible flux is the component actually changing the available potential energy of the system irreversibly.
Since
$K_{{\it\theta}}^{T}$
would tend to the molecular value for vanishing
$l_{T}$
, molecular diffusivity may produce an N-shaped flux–gradient relation in highly stable conditions. However, this is not observed here, since the estimates of the flux obtained from (3.1) are at least three orders of magnitude higher than the molecular values. Equation (3.1) is based on the assumption of a constant mixing efficiency, whose validity is a subject of active discussion (see Mater & Venayagamoorthy Reference Mater and Venayagamoorthy2014, and references therein), and which is used here as it is the simplest option and the most common one in the oceanographic context.
Furthermore, (3.1) can be derived only for the case of fully developed turbulence outside boundary layers. Turbulence is indeed strong in this area, with internal waves constantly breaking and producing overturns throughout the depth range spanned by the mooring. We would expect the presence of the sea floor to bias the Thorpe-scale-based estimates of diffusivity and flux, due to the introduction of an external length scale (the distance from the boundary), as suggested by Thorpe (Reference Thorpe1977). This is checked by comparing estimates of
${\it\phi}_{{\it\theta}}^{T}$
computed excluding the deeper parts of the mooring in § 4.
We apply the Thorpe reordering approach in combination with (3.1) to our data, averaged in 200 s long independent time windows, obtaining a time series of vertically averaged estimates of diffusivity. No isolated overturns can be identified in this (highly turbulent) data set, and thus the vertical averaging is performed over the whole mooring depth (in some cases excluding the deepest 25 m of the mooring, see § 4).
3.3 Winters and D’Asaro flux estimation method
The Thorpe-scale-based method for estimating the temperature flux cannot resolve the vertical details of the turbulent fluxes, since it is based on the r.m.s. displacement in the whole profile. In this context, the method suggested by Winters & D’Asaro (Reference Winters and D’Asaro1996) provides a framework to study the vertical dependence of turbulent fluxes. We outline it here. The advection–diffusion equation for a scalar (in this case the potential temperature
${\it\theta}$
) in a fluid flow is
$$\begin{eqnarray}\frac{\partial {\it\theta}}{\partial t}+\boldsymbol{u}\boldsymbol{\cdot }\boldsymbol{{\rm\nabla}}{\it\theta}={\it\kappa}{\rm\nabla}^{2}{\it\theta},\end{eqnarray}$$
where
$t$
is time,
$\boldsymbol{u}$
is the velocity vector and
${\it\kappa}$
is the coefficient of molecular diffusivity. Winters & D’Asaro (Reference Winters and D’Asaro1996) considered the purely diffusive flux
${\it\kappa}\boldsymbol{{\rm\nabla}}{\it\theta}$
, i.e. the irreversible flux, averaged over an isothermal surface
$S$
:
$$\begin{eqnarray}{\it\phi}_{{\it\theta}}^{d}=\frac{1}{A}\int _{S}{\it\kappa}\boldsymbol{{\rm\nabla}}{\it\theta}\boldsymbol{\cdot }\hat{\boldsymbol{n}}\,\text{d}S.\end{eqnarray}$$
In (3.3),
$\hat{\boldsymbol{n}}$
is the unit vector normal to the isothermal surface
$S$
(pointing in the direction of increasing
${\it\theta}$
). The flux is normalised by
$A$
, the cross-sectional area of the surface we are considering, rather than by
$A_{S}$
, the surface area of
$S$
(
$A_{S}\geqslant A$
). Winters & D’Asaro (Reference Winters and D’Asaro1996) showed that
${\it\phi}_{{\it\theta}}^{d}$
can be expressed as
$$\begin{eqnarray}{\it\phi}_{{\it\theta}}^{d}={\it\kappa}\frac{\text{d}z^{\ast }}{\text{d}{\it\theta}}\langle |\boldsymbol{{\rm\nabla}}{\it\theta}|^{2}\rangle ,\end{eqnarray}$$
with brackets indicating isothermal averaging. The coordinate
$z^{\ast }$
is defined by
$$\begin{eqnarray}z^{\ast }(\boldsymbol{x},t)=\frac{1}{A}\int _{V}{\mathcal{H}}({\it\theta}(\boldsymbol{x},t)-{\it\theta}(\boldsymbol{x}^{\prime },t))\,\text{d}^{3}\boldsymbol{x}^{\prime },\end{eqnarray}$$
with
$V$
the volume containing the fluid,
$\boldsymbol{x}$
the position vector and
${\mathcal{H}}$
the Heaviside function. Equation (3.5) defines a reordering procedure, mapping the three-dimensional volume occupied by the fluid onto the one-dimensional coordinate
$z^{\ast }$
. The coordinate
$z^{\ast }$
has a unique value for all points on an isothermal surface, with
$z^{\ast }({\it\theta}_{2})>z^{\ast }({\it\theta}_{1})$
for
${\it\theta}_{2}>{\it\theta}_{1}$
.
It should be noted that the Osborn–Cox method provides an estimate of the average advective flux rather than the diffusive flux (3.3); these two fluxes are not necessarily equal (see Winters & D’Asaro Reference Winters and D’Asaro1996). While the flux estimates according to the Osborn–Cox method have been computed, their correlation coefficient with the Thorpe-scale-based estimates is lower (0.6 versus 0.8). This is not surprising, since both the Winters–D’Asaro and Thorpe methods provide an estimate of the diffusive flux, while the Osborn–Cox method provides an estimate of the advective flux.
Equation (3.4) is derived without any assumption on the flow, and is thus correct for scalars satisfying (3.2) in incompressible flows. However, its use requires the complete knowledge of the three-dimensional temperature field, and it is thus not directly applicable in an oceanographic context, because observations are usually one-dimensional profiles. Winters & D’Asaro (Reference Winters and D’Asaro1996) suggested a method to approximate (3.4) if a series of one-dimensional profiles is available. We adapt the method of Winters & D’Asaro (Reference Winters and D’Asaro1996) to our particular data set in the following way.
Each profile provides a series of potential temperature measurements
${\it\theta}_{i}(z_{j})$
at equally spaced vertical positions
$z_{j}$
, with
$j=1,2,\ldots ,M$
(
$z_{n}>z_{m}$
for
$n>m$
), and at equally spaced times
$i$
. At
$z_{j}$
, the temperature gradient
$|\boldsymbol{{\rm\nabla}}{\it\theta}|_{i}(z_{j})$
is estimated, using the Taylor hypothesis for obtaining its horizontal component (see Appendix). Here,
$P$
consecutive profiles in time are combined to provide an estimate of the flux. The procedure is as follows (Winters & D’Asaro Reference Winters and D’Asaro1996).
-
(i) Each profile, of both
${\it\theta}$
and its gradient, is sorted in ascending order of potential temperature. The original profiles
${\it\theta}_{i}(z_{j})$
and
$|\boldsymbol{{\rm\nabla}}{\it\theta}|_{i}(z_{j})$
thus become
${\it\theta}_{i}(z_{j}^{T})$
and
$|\boldsymbol{{\rm\nabla}}{\it\theta}|_{i}(z_{j}^{T})$
, where
$z_{j}^{T}$
is the equally spaced depth coordinate obtained from the ‘Thorpe reordering’ method of Thorpe (Reference Thorpe1977). -
(ii) Each profile
${\it\theta}_{i}(z_{j}^{T})$
and
$|\boldsymbol{{\rm\nabla}}{\it\theta}|_{i}(z_{j}^{T})$
, having constant sampling rate, is interpolated onto a set of equally spaced values of potential temperature
${\it\theta}_{j}$
, obtaining the new profiles
${\it\theta}_{i}({\it\theta}_{j})$
and
$|\boldsymbol{{\rm\nabla}}{\it\theta}|_{i}({\it\theta}_{j})$
. The sampling rate of these new profiles is
$1\times 10^{-3}~\text{K}$
, a safe upper limit of the instrument precision. -
(iii) Isothermal averages are formed:
The hypothesis is that the(3.6a )
$$\begin{eqnarray}\displaystyle & \displaystyle \langle |\boldsymbol{{\rm\nabla}}{\it\theta}|^{2}\rangle ^{P}({\it\theta}_{j})=\frac{1}{P}\mathop{\sum }_{i=1}^{P}|\boldsymbol{{\rm\nabla}}{\it\theta}|_{i}^{2}({\it\theta}_{j}), & \displaystyle\end{eqnarray}$$
(3.6b )
$$\begin{eqnarray}\displaystyle & \displaystyle \langle z^{T}\rangle ^{P}({\it\theta}_{j})=\frac{1}{P}\mathop{\sum }_{i=1}^{P}z_{i}^{T}({\it\theta}_{j}). & \displaystyle\end{eqnarray}$$
$P$
-averages of (3.6) provide a good estimate of the local isoscalar average of the temperature gradient, and of the
$z^{\ast }$
coordinate respectively.
-
(iv) Using the profile
$\langle z^{T}\rangle ^{P}({\it\theta}_{j})$
, an estimate
$(\text{d}z^{T}/\text{d}{\it\theta})^{P}({\it\theta}_{j})$
of the background gradient is computed using centred discrete derivatives. -
(v) The estimate of the flux is finally provided by rewriting (3.4) as
(3.7)The flux is positive downwards.
$$\begin{eqnarray}{\it\phi}_{{\it\theta}}^{d}({\it\theta}_{j})={\it\kappa}\left(\frac{\text{d}z^{T}}{\text{d}{\it\theta}}\right)^{P}({\it\theta}_{j})\langle |\boldsymbol{{\rm\nabla}}{\it\theta}|^{2}\rangle ^{P}({\it\theta}_{j}).\end{eqnarray}$$
-
(vi) The estimate of the flux of (3.7) is interpolated back on the depth coordinate
$z_{j}$
.
The results presented here use
$P=200$
, and thus the quantities are averages over 200 s of data. This time scale is chosen for being the longest one still shorter than all fast gravity wave periods (
$2{\rm\pi}/N$
). In other words, during 200 s, the measured temperature changes due to turbulence and horizontal advection, rather than due to internal waves. Since the mean velocity is of the order of
$0.1~\text{m}~\text{s}^{-1}$
, the horizontal displacement during this time is of the order of 10 m. The choice
$P=200$
thus avoids averaging regions with widely different dynamics, while reducing the uncertainty of the estimates (a turbulent overturn in the data can easily span tens of metres in the vertical). The use of different values of
$P$
(
$\pm 100$
) does not change the results significantly, except for an obvious increase of noise for lower values of
$P$
and a decrease of the dynamic range for higher values of
$P$
.
An essential element in the method is the interpolation from the depth coordinate to the temperature coordinate, which provides a consistent way to estimate the background gradient, on an isotherm rather than on an isobath, and thus provides a consistent definition of the flux at each depth (after interpolating back onto the depth coordinate).
It is important to note that the estimate of the temperature gradient
$\boldsymbol{{\rm\nabla}}{\it\theta}$
, obtained from the measurements, is biased low. This is a well-known issue, already noted by Thorpe (Reference Thorpe1977), related to the fact that the instruments do not resolve the length scales of molecular diffusion (see Cimatoribus & van Haren Reference Cimatoribus and van Haren2015, for further details). For a more general discussion of the wavenumber spectrum of scalar gradients, see Warhaft (Reference Warhaft2000). To estimate the impact of this underestimation, we compare the values of
${\it\phi}_{{\it\theta}}^{d}$
from (3.7) with those obtained using the classical Thorpe method, considered for the moment the ground truth. In the next section, we will show that the two estimates are correlated; a linear fit in the log–log space thus provides an empirical way to compensate for the underestimation of the gradient. While a completely independent fully resolved measurement would ideally be the only term of comparison, we argue that this approach is a plausible, albeit imperfect, one.
Comparison of the two flux estimates using the Thorpe scale method (
${\it\phi}_{{\it\theta}}^{T}$
, § 3.2) and the Winters–D’Asaro method (
${\it\phi}_{{\it\theta}}^{d}$
, § 3.3). (a) The comparison of the flux estimates vertically averaged over the entire mooring. (b) The comparison of the flux estimates averaged excluding the lower 25 m of the mooring. The red line in each panel is a linear fit in the log–log space; the result of the fit is shown.
4 Results
4.1 Comparison between the two flux estimates
Since the Thorpe-scale-based estimate is a vertical average, we consider the vertical average of
${\it\phi}_{{\it\theta}}^{d}$
as well. As already discussed,
${\it\phi}_{{\it\theta}}^{T}$
is expected to provide a good estimate of the average flux only far from the boundary. Close to the sea floor, the Thorpe scale is limited by the distance from the sea floor itself rather than by the buoyancy gradient, as hypothesised in the derivation of (3.1). On the other hand, the estimate (3.7) does not rely on any such hypothesis, and is thus expected to be valid everywhere, within the limits due to the under-resolution of temperature gradients discussed above.
Figure 4(a) shows the comparison of the two estimates in log–log space, vertically averaged over the whole mooring. Interpolation in time is performed to align the two time series, with a time step of 200 s. A linear fit of the two estimates in the log–log space gives
$\log _{10}({\it\phi}_{{\it\theta}}^{d})=-5.6+0.5\,\log _{10}({\it\phi}_{{\it\theta}}^{T})$
. The coefficient of determination of this fit is
$R^{2}=0.58$
. The two estimates are thus correlated, with
${\it\phi}_{{\it\theta}}^{d}$
approximately proportional to
$({\it\phi}_{{\it\theta}}^{T})^{1/2}$
, even if there is significant spread in the data. The spread may be due to the intrinsic limitations of the estimation methods (in particular for
${\it\phi}_{{\it\theta}}^{T}$
), or may be connected to the fact that we do not resolve the entire three-dimensional turbulent field with our measurements. The large underestimation of
${\it\phi}_{{\it\theta}}^{T}$
by
${\it\phi}_{{\it\theta}}^{d}$
is the effect of limited resolution, as discussed above.
With the aim of assessing whether the closeness to the boundary has a significant effect on the flux estimates, we show in figure 4(b) the two flux estimates vertically averaged excluding the bottom 25 m of the mooring. In the case of the Thorpe estimate, not defined at each depth, we compute a value of the diffusivity only using the Thorpe displacements in the top 75 m of the mooring. The results are very similar to those using the whole depth range, and the increase in the coefficient of determination is marginal (
$R^{2}=0.61$
). We thus suggest that, at least in this data set, the effect of the solid boundary on the Thorpe scale flux estimates is small. This result will be discussed in more detail in § 5. We must note, however, that our measurements do not extend below an HAB of approximately 5 m.
Assuming that the Thorpe scale flux estimate is correct (on average), we correct the underestimated
${\it\phi}_{{\it\theta}}^{d}$
using the results of the fit in figure 4(a). This approach can be compared with the common practice of fitting the Batchelor wavenumber spectrum on under-resolved temperature measurements to obtain an estimate of the temperature variance dissipation rates (Ruddick, Anis & Thompson Reference Ruddick, Anis and Thompson2000). By necessity, we are also assuming that the correction can be applied in the same way for each estimate
${\it\phi}_{{\it\theta}}^{d}$
(at all times and depths). This ‘corrected’ value of the flux will be denoted by
${\it\phi}_{{\it\theta}}^{\ast }$
.
4.2 Dependence of
${\it\phi}_{{\it\theta}}^{\ast }$
on depth and stratification
As discussed in § 3.3, the method proposed by Winters & D’Asaro (Reference Winters and D’Asaro1996) for estimating the temperature (or other scalar) flux enables a consistent definition of the flux and the background gradient at each depth. In this section, we fully exploit this advantage, studying the dependence of the flux on the local stratification and on the HAB.
The probability density functions (p.d.f.s) of the flux estimates
${\it\phi}_{{\it\theta}}^{\ast }$
are shown in figure 5, computed separately for the cooling and warming tidal phases ((a) and (b) respectively, see § 3.1 for the definition of the phases). During each tidal phase, two p.d.f.s are computed, one for data from the upper half of the mooring (empty symbols) and one for data from the lower half (filled symbols). The insets in the figure show the ratio between the empirical distribution and a normal distribution having the mean and the standard deviation of the log-transformed data. The p.d.f.s all show an approximately log-normal behaviour close to the mean, similar to what is observed for dissipation, both in the laboratory and in the ocean (see, e.g., Baker & Gibson Reference Baker and Gibson1987; Warhaft Reference Warhaft2000). Deviations from log-normality are strong only in the tails, where undersampling leads to shorter tails in the empirical distribution in comparison to the theoretical distribution. Deviations from log-normality are even smaller for the p.d.f. of the entire data set (not shown). The effect of noise, leading to fatter left tails, is barely visible only in the cooling phase close to the sea floor, testifying to the extremely low noise level of the instruments. Values of the flux higher by almost one order of magnitude are observed during the cooling phase; higher values of the flux are also observed in the lower half of the mooring.
The p.d.f.s of the logarithm of the flux (
$\log _{10}({\it\phi}_{{\it\theta}}^{\ast })$
), compensated for underestimation of the temperature gradient (see text). (a) The results during the cooling phase of the tide and (b) those during the warming phase. The insets above each panel show the ratio of the p.d.f. with a normal distribution, demonstrating that the flux has an approximately log-normal distribution. Filled symbols refer to data from the lower half of the mooring (continuous line in the inset) and empty symbols to data from the upper half of the mooring (dashed line in the inset). The horizontal scale is logarithmic.
Temperature flux
${\it\phi}_{{\it\theta}}^{\ast }$
as a function of background stratification
$\text{d}{\it\theta}/\text{d}z^{T}$
(a–c) and as a function of HAB (d). The error bars, visible only for large values of
$\text{d}{\it\theta}/\text{d}z^{T}$
, are 95 % confidence intervals, computed using Cox’s method as described and evaluated in Zhou & Gao (Reference Zhou and Gao1997). (a) Averages of all the flux estimates; (b) the estimates averaged separately during the cooling (blue triangles) and warming (red circles) tidal phases; (c) the estimates averaged separately in the lower half (red filled symbols) and upper half (blue empty symbols) of the mooring; (d) the depth dependence of the flux, estimated separately for the cooling (triangles) and warming (circles) tidal phases. The dashed lines in (a–c) are guides for the eye
$\propto (\text{d}{\it\theta}/\text{d}z^{T})^{1/2}$
(see text).
In figure 6(a), estimates of
${\it\phi}_{{\it\theta}}^{\ast }$
are shown, averaged in 32 bins of the background temperature gradient
$\text{d}{\it\theta}/\text{d}z^{T}$
(equally spaced in log space). Since
${\it\phi}_{{\it\theta}}^{\ast }$
is log-normally distributed, the mean is estimated from
$\exp ({\it\mu}+{\it\sigma}^{2}/2)$
, with
${\it\mu}$
and
${\it\sigma}^{2}$
the mean and variance of the log-transformed data respectively, as motivated, e.g., in Baker & Gibson (Reference Baker and Gibson1987). While this procedure is not essential, it improves the results in particular at the highest temperature gradients, for which fewer data points are available. The average flux has a tight nonlinear dependence on stratification, with a local maximum for
$\text{d}{\it\theta}/\text{d}z^{T}\approx 10^{-3}~\text{K}~\text{m}^{-1}$
, a local minimum at
$\text{d}{\it\theta}/\text{d}z^{T}\approx 10^{-2}~\text{K}~\text{m}^{-1}$
and a new increase for higher values of the stratification. It is impressive to see how closely this curve resembles the flux–gradient relationship hypothesised by Phillips (Reference Phillips1972), Posmentier (Reference Posmentier1977) and B98, and predicted by our model (§ 2, figure 1). As expected, the flux also grows as the square root of the gradient for high values of the gradient (dashed line). We stress that such a smooth curve can be obtained only after sufficient averaging (i.e. using a large data set such as the one here). We would expect such a relation to hold locally as well, had we sufficient knowledge of the flow to extract it reliably.
Figure 6(a) could suggest that the diffusive flux in a turbulent flow can be universally described by a function of stratification alone. However, this is obviously not the case, as is clear from a more detailed analysis of the data set. Figure 6(b) shows the mean flux averaged separately during the cooling and warming phases. Stronger fluxes take place during the cooling phase (blue triangles), with the positions of the stationary points virtually the same as the ones of the warming phase (red circles).
Figure 6(c) shows the mean flux averaged separately in the lower half (filled red squares) and the upper half of the mooring (empty blue squares). The curve from the data closer to the bottom has a local maximum at a higher value of stratification, and a less pronounced local minimum. The comparison of (a–c) suggests that the flux is a function of stratification alone for very strong stratification, at least within the large error bars of the estimates therein. For strong stratification (
$\text{d}{\it\theta}/\text{d}z^{T}>10^{-1}~\text{K}~\text{m}^{-1}$
), the flux increases as
$(\text{d}{\it\theta}/\text{d}z^{T})^{1/2}$
, as predicted by the theory in § 2. A small bump is visible in the flux–gradient relation at
$\text{d}{\it\theta}/\text{d}z^{T}\approx 5\times 10^{-2}~\text{K}~\text{m}^{-1}$
in (b,c) (red symbols, lower half of the mooring). This barely visible feature is robust to changes in the number of profiles used for averaging (
$P$
of the Winters–D’Asaroestimates).
Dependence of
${\it\phi}_{{\it\theta}}^{\ast }$
on stratification and HAB, shown in grey contours with a logarithmic scale. Red (blue) dots mark the positions of the local maxima (minima) of the flux. The red and blue dashed lines show the positions of the stationary points predicted by the model (see text).
Figure 6(d) shows the dependence of the flux on the HAB, as obtained by averaging the data in 24 equally spaced depth bins separately during the two tidal phases (triangles for cooling phase, circles for warming phase). The essential new information here is that the flux has a maximum at approximately 30–40 m above the bottom, which must be linked to horizontal advection if the profile is stationary. During each tidal phase, the maximum averaged flux value is approximately a factor 2 of the minimum.
The dependence of
${\it\phi}_{{\it\theta}}^{\ast }$
on the gradient and HAB is summarised as a contour plot in figure 7, showing also the estimated positions of the local maxima and minima of the flux. This figure will be compared with figure 2 in the next section.
5 Discussion and conclusions
We have presented observations of the temperature flux in the deep ocean, above a sloping sea floor. The observations provide detail on the dependence of the flux on the local vertical temperature gradient, and on the height above the sea floor. They support the predictions of the model for the flux–gradient–HAB relation discussed in § 2. This is surprising, in particular considering that we solved a highly simplified version of the model, assuming steady state and negligible diffusion of turbulent kinetic energy. It is, however, important to stress once more that averaging is essential here: only using a high-resolution long data set is it possible to extract the smooth curves of figure 6, which are otherwise obscured by the huge variability of the turbulent quantities. Furthermore, the use of the method outlined by Winters & D’Asaro (Reference Winters and D’Asaro1996) is essential for correctly disentangling the flux–gradient–HAB relation as well. While a more spatially complete observation of the turbulent flow is expected to reduce the need for averaging, and to bring different estimation methods to convergence, we are very far from this situation, both in the ocean and, possibly to a somewhat lesser extent, also in the laboratory.
The agreement between model and observations supports the mixing-length hypothesis, in particular when considering the similar patterns observed in the model and in the observations as a function of the HAB (compare figures 2 and 7). The agreement is also supported by the fact that the model and observations span a similar range of scales (the model results being presented in non-dimensional form) in terms of gradient, HAB and flux. Furthermore, the results strongly suggest that the solid boundary limits the mixing length. Thorpe (Reference Thorpe1977) argued that the Thorpe scale could not be used as a proxy for turbulence dissipation close to a boundary, since the distance from the boundary, rather than stratification, limits the overturn scale therein. While our results support the latter argument, they also suggest that the mixing length is limited by the wall in a way similar to the Thorpe scale, in particular considering the model and the comparison between
${\it\phi}_{{\it\theta}}^{d}$
and
${\it\phi}_{{\it\theta}}^{T}$
. This in turn supports the use of the Thorpe scale close to the boundary for obtaining an order-of-magnitude estimate of turbulence parameters, at least as long as stratification is present and averaging is performed.
It is interesting to consider the positions of the stationary flux points in the model and the observations. This is done in figure 7, reporting as red (blue) dots the positions of the maxima (minima) of the flux–gradient relation estimated from the observations, and with similarly coloured dashed lines those predicted by the model. The stationary points of the model, according to (2.9), are loci where the gradient is proportional to
$\log ^{2}({\it\eta}z)/z^{2}$
, using
${\it\eta}=0.5~\text{m}^{-1}$
. The model and the observations roughly agree.
It should be noted that the value of
${\it\eta}$
used in figure 7 is much smaller than the one computed from the values in table 1. The classical theory, which we followed, states that
${\it\eta}=u_{\ast }/{\it\nu}$
, with
$u_{\ast }$
of the same order of magnitude as the velocity fluctuations. If we use a larger value of
${\it\eta}$
, e.g.
$10^{3}~\text{m}^{-1}$
as implied by the values in table 1, a much worse agreement with the observations is obtained. This is confirmed by fitting a function
$\text{d}{\it\theta}/\text{d}z^{T}\propto \log ^{2}({\it\eta}z)/z^{2}$
to the estimated positions of the stationary points in the
$(\text{d}{\it\theta}/\text{d}z^{T},z)$
plane. We suggest that this mismatch results from a reduced vertical shear stress in the stratified case considered here in comparison with the uniform-density fluid considered in the classical theory for a log layer; buoyancy forces ‘isolate’ the stacked shear layers from each other.
No attempt has been made to fit the model to the observations and, in order to obtain a better quantitative agreement, the analysis of § 2 should be refined in the future. In particular, a more realistic functional dependence of the mixing length on the HAB is needed, something that could be obtained only through a dedicated investigation, most likely in a controlled environment. The geometric average is used here only for simplicity, due to the lack of a more physically motivated functional form.
A more realistic representation of the mixing length may also improve the dependence of the magnitude of the flux on the HAB. If we compare figure 1 with figure 6(c), we see that the model predicts a marked decrease in the value of the maximum of the flux for higher HAB, while the observations indicate that this decrease is marginal at best.
Furthermore, the transport of turbulent kinetic energy should not be neglected, for which the full differential equation system needs to be considered. In this respect, it is interesting to note that, in the steady state, (2.1a ) implies that the buoyancy flux is linear with depth. Observations indicate that this is not a realistic result (figure 6 d), suggesting that a dissipation or forcing term in the buoyancy gradient equation (2.1a ) may be missing. The derivation of the model from the Reynolds decomposed primitive equations would obviously provide a firmer ground for discussion.
Further progress also implies considering, in the forcing term
${\mathcal{P}}$
of (2.1), the effect of direct wave breaking, e.g. due to wave energy concentration after reflection of a wave train. This mechanism, whose energetic implications are discussed in some detail by Scotti (Reference Scotti2015), is not included in our forcing term, which rather assumes a classical wall-bounded shear layer. The forcing from breaking waves is probably related to the slope and horizontal structure of the sea floor as well, another element that is ignored in our model, and that may contribute to the discrepancies with the observations noted above.
As a last remark, we briefly consider the small bump visible in the flux–gradient relation at
$\text{d}{\it\theta}/\text{d}z^{T}\approx 5\times 10^{-2}~\text{K}~\text{m}^{-1}$
, in particular in figure 6(b,c) (filled symbols alone). This feature, not predicted by the model, may be linked to non-stationarities or, possibly, to a second instability of the flux–gradient relation, whose origin is at the moment unclear. A preliminary investigation suggests that it may be more prominent in the profiles that do not have a layered structure.
Further progress in the analysis of this data set implies a careful characterisation of the scales of the layers in the observations, and, possibly, their initial development and further evolution. While these issues go beyond the scope of this work, we would expect that their correct description by the model will prove to be essential for its further development. A closer description of the observations by the model may also provide the basis for applications, e.g. as a turbulence closure for a larger-scale ocean model.
Acknowledgement
We would like to thank the crew of the RV Pelagia for enabling the deployment and recovery of the instruments, and all of the technicians, in particular M. Laan, indispensable for the success of the measurements. We would like to thank Professor L. Maas for providing comments and suggestions for the improvement of the model. The authors are grateful to the referees and the editor for providing insightful and detailed comments on the first version of this paper.
Appendix
The velocity data recorded by the ADCP are used for transforming the measured temperature time series at each depth into a spatial series, assuming the validity of the ‘frozen turbulence’ (Taylor) hypothesis. Since the velocity is not constant, the approach of Pinton & Labbé (Reference Pinton and Labbé1994) is followed.
The velocity field
$v$
is decomposed into a mean component
$\langle v\rangle$
and a fluctuating one
$v^{\prime }$
using a low-pass Butterworth filter with stop band
$1/{\it\sigma}=3000~\text{s}$
. On average, this can be considered to be the forcing frequency, based on the frequency spectra of temperature and velocity fluctuations (not shown). The mean velocity is higher in the upslope tidal phase, while the velocity fluctuations have a similar distribution during both tidal phases.
The mean velocity obtained by filtering the ADCP data is vertically averaged in four bins, and these averages are used to transform the temperature time series of each thermistor to the spatial domain (using the average velocity from the appropriate vertical bin). The spatial series thus obtained has a non-constant sampling rate because the velocity is not constant in time. By linear interpolation to a constant time step of 0.2 m, a constant-step spatial series is obtained, reducing the size of the data set by a factor of approximately 2. The vertical resolution is fixed to 0.7 m by the thermistor spacing.
The turbulence intensity (
$\langle v^{\prime 2}\rangle ^{1/2}/\langle v\rangle$
) estimated from the ADCP data is much smaller than
$1$
most of the time, supporting the use of the Taylor frozen turbulence hypothesis. However, since the ADCP resolution is insufficient to resolve the complete turbulence inertial range, this value of the turbulence intensity is probably underestimated. On the other hand, the presence of noise in the ADCP data artificially increases the estimated turbulence intensity.
Only the first increments of the Taylor-transformed data are used in this work, for computing the horizontal component of the temperature gradient discussed in § 3.3.
