## Abstract

The interaction of shear, stratification, and turbulence in boundary layers on sloping topography is investigated with the help of an idealized theoretical model, assuming uniform bottom slope and homogeneity in the upslope direction. It is shown theoretically that the irreversible vertical buoyancy flux generated in the boundary layer is directly proportional to the molecular destruction rate of small-scale buoyancy variance, which can be inferred, for example, from microstructure observations. Dimensional analysis of the equations shows that, for harmonic boundary layer forcing and no rotation, the problem is governed by three nondimensional parameters (slope angle, roughness number, and ratio of forcing and buoyancy frequencies). Solution of the equations with a second-moment closure model for the turbulent fluxes reveals the periodic generation of gravitationally unstable boundary layers during upslope flow, consistent with available observations. Investigation of the nondimensional parameter space with the help of this model illustrates a systematic increase of the bulk mixing efficiencies for (i) steep slopes and (ii) low-frequency forcing. Except for very steep slopes, mixing efficiencies are substantially smaller than the classical value of Γ = 0.2.

## 1. Introduction

Since Munk (1966) suggested that processes near sloping topography may significantly contribute to basin-scale mixing, numerous investigations in the ocean, in lakes, and in the laboratory have provided support for this idea. Early observations of detached boundary layers advecting into the ocean’s interior (Armi 1978, 1979a) showed indirect evidence for the effect of boundary mixing that was later supported by turbulence microstructure measurements, revealing strongly enhanced dissipation levels in the bottom boundary layer (BBL) for various sites in the ocean (e.g., Nash et al. 2007; Rudnick et al. 2003; Lueck and Mudge 1997) and in lakes (Lemckert et al. 2004; Wüest and Lorke 2003; Gloor et al. 2000). The probably most direct demonstration of the basin-scale importance of boundary mixing, however, came from tracer experiments conducted in stratified basins of different scales, suggesting a dramatic increase of the vertical spreading rates of the tracer after it had touched the boundaries (Goudsmit et al. 1997; Ledwell and Bratkovich 1995; Ledwell and Hickey 1995).

Breaking and shoaling internal waves and bottom friction caused by near-bottom currents are believed to be the major energy sources for BBL mixing in stratified basins (Thorpe 2005; Wüest and Lorke 2003). Near sloping topography, an interesting phenomenon related to the oscillating shear caused by internal waves or internal seiching motions has been revealed by a number of recent measurements in lakes (Lorke et al. 2005, 2008). These data suggest that differential advection due to the near-bottom vertical shear acting on the upslope density gradient during upslope flow may create unstable stratification, with considerable implications for turbulence and mixing in the BBL. Similar occurrences of gravitationally unstable BBLs were found in the idealized numerical simulations of Slinn and Levine (1997), who studied mixing caused by cross-slope tidal currents on sloping topography at low latitudes. Although these studies have clearly shown that convection provides an important additional energy source for turbulence with a yet unknown impact on basin-scale mixing, an investigation of this effect and its implications for mixing over a broader parameter range is still lacking.

Idealized theoretical models of boundary mixing, assuming constant diffusivities and homogeneity in the upslope direction (Phillips 1970; Wunsch 1970; Thorpe 1987), and laboratory experiments (Phillips et al. 1986) have illustrated that mixing near a sloping boundary inevitably leads to a tilting of near-bottom isopycnals, thus setting up an internal pressure gradient that drives a residual or “secondary” circulation. Garrett (1990, 1991) recognized that this secondary circulation could have a crucial feed back on mixing in the BBL in at least two respects. First, he pointed out that the secondary circulation has a tendency to create stable stratification in the BBL, thus providing a mechanism for restratification without requiring an exchange of fluid between the BBL and the interior. This effect is likely to result in a higher mixing efficiency and thus in a larger buoyancy flux. Second, however, as shown by Garrett (1990, 1991), the vertical advective buoyancy flux associated with the secondary circulation is expected to be countergradient, thus reducing the net vertical buoyancy flux. This puzzling situation is further complicated by the fact that mixing parameters are tightly coupled to changes in shear and stratification.

Our aim here is to derive a framework for estimating the total vertical buoyancy flux in a sloping BBL from microstructure observations and numerical simulations that includes the shear-induced generation and destruction of stratification and its effect on net mixing. Following previous theoretical work (Phillips 1970; Thorpe 1987; Garrett 1990) and idealized modeling studies (Ramsden 1995a,b), the theory will be based on an idealized geometry with an infinite slope, ignoring flow variations in the upslope direction. The latter assumption eliminates any “tertiary” flow (e.g., intrusions), resulting from the upslope divergence of the secondary flow, as one potentially important mechanism for boundary layer restratification (Garrett et al. 1993; McDougall 1989; Phillips et al. 1986). Nevertheless, the investigations mentioned above have demonstrated that important insights into the basic mechanisms of boundary mixing can be gained in particular with such simplified models.

From the equations governing the problem, we will identify the key nondimensional parameters and explore the parameter space with the help of a second-moment turbulence model, similarly to previous boundary mixing studies by Weatherly and Martin (1978) and Middleton and Ramsden (1996) that have focused on the effect of rotation. Because here we focus on the derivation of the general theoretical framework, the presentation of the basic mixing mechanisms, and the role of gravitationally unstable boundary layers for net mixing, we leave the additional complications due to rotational effects to future work.

## 2. Geometry and governing equations

### a. Geometry

We investigate a Boussinesq fluid periodically moving up and down a uniform slope with slope angle *α* in an infinitely deep basin with variable stratification,

where *ẑ* is the vertical coordinate and

is the buoyancy, with *ρ* denoting the density, *ρ*_{0} denoting a constant reference density, and *g* denoting the acceleration of gravity (Fig. 1a). Outside the BBL, we assume constant vertical stratification (*N* = *N*_{∞}) and horizontal isopycnals. Because rotation is ignored, the flow remains in a plane spanned by the geopotential coordinates *x̂* and *ẑ* (Fig. 1a).

The geometrical key simplification results from the assumption that all gradients in the upslope direction are negligible, except the upslope buoyancy gradient, which is assumed to be constant (Garrett et al. 1993). Defining a rotated coordinate system with the *x* axis pointing upslope and the *z* axis directed normal to the slope (Fig. 1), we find that the upslope buoyancy gradient is related to the interior stratification according to

The buoyancy gradient normal to the slope,

converges to the constant value

outside the BBL, where *b*_{∞} denotes the buoyancy that would be observed for perfectly horizontal isopycnals throughout the domain. As indicated in Fig. 1, in general the isolines of *b* and *b*_{∞} coincide only outside the BBL.

For the given geometry, it is straightforward to derive a relation between the different forms of the buoyancy frequency,

Because in a perfectly well-mixed BBL the buoyancy gradient normal to the slope vanishes by definition (*Ň* ^{2} ≡ 0), we conclude from (6) that *N* ^{2} = sin^{2}*α* in the BBL, illustrating that the terms well mixed and unstratified do not coincide (contrary to the situation in a basin with flat bottom). Outside the BBL, we find *Ň* ^{2} = cos*α*, which may be inserted into (6) to recover the trivial case *N* ^{2} = .

### b. Governing equations

Under the assumptions mentioned above, the Reynolds-averaged Boussinesq equations in rotated coordinates adopt the following simple form:

where *u* is the upslope velocity and *p* is the pressure minus the pressure in a reference state of rest (Garrett et al. 1993). The turbulent flux of momentum, *τ _{x}* =

*u*′

*w*′ (

*w*is the velocity normal to the slope, the prime denotes a fluctuating quantity, and the bar denotes the ensemble average), and the turbulent buoyancy flux,

*G*=

*w*′

*b*′, are computed from a second-moment turbulence model as described in detail below.

The upslope pressure gradient in (7) can be eliminated by taking the *x* derivative of the second relation in (7), inserting ∂*b*/∂*x* = tan*α* ∂*b*_{∞}/∂*z* obtained from combining (3) and (5), and integrating with respect to *z*. This yields

where in the buoyancy equation we have used (3) to express the constant upslope buoyancy gradient (Garrett et al. 1993), and assume

at the lower boundary. The new integration constant ∂*u*_{∞}/∂*t* appearing in (8) is a function of time, and it plays the role of a prescribed interior pressure gradient governing the evolution of the flow outside the BBL. Identifying this term with a periodic function of time will serve below as a simple representation of the oscillating pressure gradient due to internal waves. This approach is completely analogous to that used by Mellor (2002) to describe BBL turbulence driven by surface waves; to be valid, waves are required to be long enough to be compatible with our assumption of upslope homogeneity.

The first term on the right-hand side of the momentum budget in (8) represents the internal pressure gradient,

arising from the tendency of tilted isopycnals to relax back toward their equilibrium positions (indicated by the gray lines in Fig. 1).

The evolution of the equilibrium buoyancy *b*_{∞} appearing in (8) is computed from the special case of zero mixing,

## 3. Description of turbulence and mixing

### a. Turbulent fluxes

The turbulent fluxes of momentum and buoyancy appearing in (8) are computed from the second-moment turbulence closure model by Canuto et al. (2001), extended here by two prognostic equations for the turbulent kinetic energy (TKE) *k* and the dissipation rate ɛ, as described in Umlauf and Burchard (2005). This model is particularly suited for our study because (i) it does not contain dimensional parameters, which is a crucial requirement for the derivation of the nondimensional expressions discussed in section 4, and (ii) it does not contain any assumptions about the vertical structure of the bottom boundary layer. For these reasons, a number of frequently used parameterizations for ocean mixing are not applicable here (e.g., Durski et al. 2004; Large et al. 1994; Pacanowsci and Philander 1981).

From previous comparison with LES and laboratory data (Umlauf and Burchard 2005; Umlauf 2009) and turbulence measurements in entraining boundary layer flows (e.g., Arneborg et al. 2007; Umlauf et al. 2010), we expect that the model used here reproduces the basic interaction mechanisms between shear, stratification, and turbulence and thus yields credible estimates and correct trends for the mixing parameters.

It should be emphasized that the derivation of the second-moment equations for geophysical boundary layers typically involves two assumptions that conflict with the applicability of the resulting models for arbitrary slopes: (i) horizontal gradients of mean quantities in the second-moment equations are ignored, and (ii) the mean shear and gravity vectors are assumed to be aligned and vertical (Mellor and Yamada 1982; Canuto et al. 2001; Umlauf et al. 2005). Clearly, for the given geometry, these approximations are justified only for *α* ≪ 1. Therefore, we restrict our following investigations to *α* ≤ 0.1, which does, however, not impose a serious constraint for the typical slopes observed in the ocean and in lakes.

Under these conditions, the turbulent fluxes are computed from

where *S* = ∂*u*/∂*z* is the shear and *ν _{t}* and

*v*are the turbulent diffusivities computed from the second-moment turbulence model. We distinguish between vertical (

_{t}^{b}*B*) and slope-normal (

*G*) buoyancy fluxes, which is essential for the correct timing of the transition from weakly stable stratification to convection. Note that, by construction, the boundary layer turbulence model used here does not provide an expression for the upslope diffusivity, which we therefore ignore in the following. This is justified because (i) shear dispersion as an important upslope transport mechanism is taken into account and (ii) the vertical component of the upslope turbulent flux is likely to be negligible for the small slopes considered here (

*α*≪ 1) unless stratification is extremely weak or the upslope diffusivity is extremely large (Garrett et al. 1993).

The second-moment equations lead to the following expressions for the turbulent diffusivities (Canuto et al. 2001; Umlauf and Burchard 2005):

where, using the “quasi equilibrium” assumption of Galperin et al. (1988), the stability functions, *c _{μ}* and can be shown to depend only on the dimensionless buoyancy number

*Nk*/ɛ.

Different from Canuto et al. (2001), we do not assume that the TKE equation is in equilibrium, and we solve the full transport equation,

where *P* = *ν _{t}S*

^{2}denotes the shear production. The dissipation rate is computed from a transport equation of the form

with all model parameters summarized in Table 1.

Using the boundary conditions ∂*k*/∂*z* = 0 and ɛ = *u _{f}*

^{3}/(

*κz*

_{0}) at

*z*= 0, where

*u*

_{f}^{2}= |

*τ*|

_{x}_{z=0}denotes the square of the friction velocity,

*z*

_{0}denotes the bottom roughness length, and

*κ*denotes the von Kármán constant (see Table 1), it can be shown that the velocity follows the expected logarithmic wall profile,

Note that all boundary conditions are imposed at the “virtual origin” (*z* = 0) of the logarithmic profile in (16), where *u* = 0 consistent with (9). Details and numerical implementation are described, for example, in Burchard et al. (2005).

### b. Quantification of mixing

To derive the framework for a precise definition of mixing, we start from the equation for the variance of turbulent buoyancy fluctuations *b*′^{2} (Tennekes and Lumley 1972; Umlauf and Burchard 2005). With the geometrical constraints described in section 2 (but without any turbulence modeling assumptions), this equation can be written as

where *F _{b}* denotes the sum of the turbulent and viscous transport terms [the exact form of

*F*is given, e.g., in Eq. (A.3) of Umlauf and Burchard 2005]. For our idealized geometry, the production term is of the form

_{b}and the destruction of small-scale buoyancy variance by molecular mixing is given by

where *ν ^{b}* is the molecular diffusivity of buoyancy (corresponding to that of either heat or salt), and we sum over

*i*= 1, 2, 3 with the

*x*corresponding to some arbitrary orthogonal coordinate system. The destruction rate

_{i}*χ*in (19) is positive by definition and describes the molecular, irreversible mixing of fluid with different densities. This quantity is of central importance for the following analysis because, as discussed in detail below, it is closely related to the irreversible vertical buoyancy flux in the BBL.

_{b}Finally, it is important to note that the turbulence model by Canuto et al. (2001) used for our numerical simulations employs the equilibrium version of (17): *P _{b}* =

*χ*. Consistent with the boundary layer approximation discussed in the previous section, the production of buoyancy variance in (18) is dominated by vertical stratification: that is,

_{b}*P*= −2

_{b}*BN*

^{2}(see Canuto et al. 2001), which is a good approximation for both the total and slope-normal production as long as

*α*≪ 1. With these assumptions, the equilibrium form of (17) reduces to the familiar Osborn–Cox model that is frequently used for the interpretation of oceanic microstructure data (Osborn and Cox 1972; Gregg 1987). It is worth pointing out that the integral of the flux term

*F*in (17) vanishes because (i) there is no mixing outside the BBL and (ii) it can be shown that

_{b}*F*= 0 at

_{b}*z*= 0 because there is no buoyancy flux through the lower boundary (Umlauf and Burchard 2005). Much of the following analysis is based on time-averaged, vertically integrated equations for periodic flows (for which the rate terms vanish as well), such that the equilibrium assumption does not imply any restrictions.

## 4. Nondimensional description

We start our nondimensional analysis from (8), rewritten in the form

As a simple representation of internal-wave forcing, we assume that the outer flow is purely harmonic,

where *U* denotes a constant velocity scale and *ω* denotes the forcing frequency. This relationship introduces the characteristic time scale *ω*^{−1} to the problem, which we use to nondimensionalize the time according to

where here and in the following the asterisk denotes a nondimensional variable. From (20), we identify the slope angle *α* and the background stratification *N*_{∞} as independent parameters. The parameters *U* and *ω* are introduced by the forcing term in (21), and the bottom roughness *z*_{0} enters the problem via the boundary condition for the dissipation rate discussed in the context of (15).

With these parameters, we nondimensionalize the vertical coordinate as

and introduce the nondimensional variables

The nondimensional form of (20) is thus

where we have defined the nondimensional frequency,

Using the nondimensional friction velocity = *u _{f}*/

*U*, the boundary condition for ɛ becomes

which introduces the last nondimensional key quantity: the roughness parameter,

For harmonic forcing, the parameter space is thus spanned be three independent quantities, *α*, *R*, and Ω. It is straightforward to show that the equations for the turbulent quantities in (13), (14), and (15) can be made dimensionless in a similar way without introducing any new nondimensional parameters. Note that the existence of three independent nondimensional parameters also follows directly from combining the five dimensional parameters governing the problem (*α*, *U*, *ω*, *N*_{∞}, and *z*_{0}) into nondimensional products.

## 5. Properties

### a. Shear-induced stratification

The effect of the shear, *S* = ∂*u*/∂*z*, on vertical stratification inside the BBL is easiest understood from the transport equation of the buoyancy frequency,

which is derived from the time derivative of (6), inserting the *z* derivative of the buoyancy equation in (8).

From (29), negative near-bottom shear (downslope flow) is seen to have a tendency to enhance stratification. The impact of this effect on mixing is, however, not obvious because there are two competing feedback mechanisms: (i) the suppression of turbulence due to the generation of stable stratification and (ii) the enhancement of mixing due to an increase of the mixing efficiency (e.g., Shih et al. 2005). Conversely, for positive near-bottom shear (upslope flow), stratification is destroyed, which may ultimately lead to unstable stratification, thus providing an additional energy source in the TKE budget (14). The relative importance of these effects will be studied in detail below.

### b. Boundary layer resonance

An equation for the boundary layer velocity relative to the interior flow, *ũ* = *u* − *u*_{∞}, can be obtained from the time derivative of the momentum equation in (20), inserting the buoyancy equation in (20) on the right-hand side. After some rearrangement, this yields

where the left-hand side is recognized as that of a harmonic oscillator with resonance at the critical frequency *ω _{c}* =

*N*

_{∞}sin

*α*. This generalizes a similar result derived by Thorpe (1987) for piecewise constant diffusivities, leading to the interesting conclusion that boundary layer resonance is reached for the frequency at which internal waves are critically reflected from the slope (e.g., Thorpe and Umlauf 2002).

### c. Residual transport

In the following, expressions for the residual volume and buoyancy fluxes are obtained, either from periodic averages for harmonic forcing or by integrating over an arbitrary “mixing event” of finite duration, assuming that, before the event (*t* < *t*_{1}) and after (*t* > *t*_{2}), the boundary layer has fully restratified and returned to a state of rest with *u* = 0 and horizontal isopycnals (*b* = *b*_{∞}).

Assuming negligible mixing outside the BBL, integration of the buoyancy equation in (8) from the bottom to an arbitrary reference level, *z* = *z _{r}*, outside the BBL yields

where we have used the facts that *G* = 0 at the bottom and no mixing occurs outside the BBL. Integrating (31) from *t*_{1} to *t*_{2} over a mixing event, the left-hand side vanishes, and thus there is no net transport, irrespective of the nature of the event and the value of *z _{r}*. Similarly, for monochromatic forcing, (21), and fully periodic flow we find

where we have introduced the symbol 〈· · ·〉 to denote the average over one period, *T* = 2*π*/*ω*. This generalizes the result derived by Phillips (1970) for stationary flows with zero mixing outside the BBL. Note that in (32) and all following results, the averaging procedure removes any dependency on the reference level *z _{r}*.

For periodic flow, the vertical structure of the residual currents follows from averaging the buoyancy equation in (8),

Because 〈*G*〉 = 0 at the bottom, it is clear from (33) that, for 〈*G*〉 < 0 near the bottom, the residual near-bottom circulation will be upslope, consistent with the available analytical solutions for constant or piecewise constant diffusivities summarized by Garrett (1990, 1991) and Garrett et al. (1993). Conversely, for 〈*G*〉 > 0 (i.e., convection dominates), we expect a reversal of the near-bottom circulation.

To obtain an expression for the residual upslope buoyancy flux, we start from the evolution equation for the square of the buoyancy *b*^{2}, which is easily derived by multiplying the buoyancy equation in (8) with *b*. The result is integrated across the BBL and added to the integral of the equation for the buoyancy variance, (17). Using (18) for the production term and the fact that *F _{b}* vanishes at the integration limits (see above), this finally yields

where we have introduced

for the total (i.e., advective plus turbulent) upslope buoyancy flux and the integrated mixing rate, respectively.

It is evident from (34) that the instantaneous upslope buoyancy flux *Q _{b}* results from mixing (last term on the right-hand side), as well as from temporal variations of the variance terms on the left-hand side. The latter terms, however, vanish if integrated over an arbitrary mixing event as described above, and (34) thus yields

relating the irreversible buoyancy flux to the mixing rate. Similarly, for periodic forcing and fully periodic flow, averaging (34) over one period yields

Note that the local buoyancy *b* is only known up to an arbitrary constant *b*_{0}, which implies that the averaged advective buoyancy flux,

contains an arbitrary component that is proportional to the averaged speed 〈*u*〉. For fully periodic conditions, however, 〈*Q _{b}*〉 becomes independent of the reference buoyancy because the vertical integral of 〈

*u*〉 vanishes according to (32).

It should be noted that both relations, (36) and (37), imply that the irreversible buoyancy flux is strictly downslope and thus associated with an irreversible increase of potential energy due to mixing. To quantify this important conclusion, we exploit the fact that 〈*Q _{b}*〉, the total upslope buoyancy flux integrated over a line

*normal*to the slope coincides with the total vertical buoyancy flux integrated over a

*horizontal*line. To see this, we note that, analogous to the derivation for stationary flows (Garrett 1991; Garrett et al. 1993), for periodic flows it can be shown that

where we have arbitrarily chosen a position on the slope that coincides with the origin of the coordinate systems shown in Fig. 1. The terms on the right-hand side of (39) are recognized as (i) the vertical advective flux, (ii) the vertical component of the slope-normal turbulent flux, and (iii) the vertical component of the upslope turbulent flux. The sum of the latter two terms thus corresponds to the total vertical turbulent flux, referred to as 〈*Q _{G}*〉 in the following.

With the interpretation of 〈*Q _{b}*〉 as the total vertical buoyancy flux, (36) and (37) form one of the key results of this paper, linking the molecular fluxes at the smallest (Batchelor) scales of the flow to the large-scale net buoyancy flux and thus to the change of background potential energy (Winters et al. 1995). In this context, it is worth pointing out that, based on an analysis of a similar problem in isopycnal coordinates, Garrett (2001) arrived at the related conclusion that only the turbulent buoyancy flux

*across*mean isopycnals (i.e., the mean diapycnal flux) has an effect on the net vertical buoyancy flux.

## 6. Boundary layer evolution

To gain some first insight into the properties of the system described by (8), (9), and (11), we study the evolution of a BBL driven by harmonic forcing, (21), with *M*_{2} tidal period (*T* = 12.4 h). The turbulent fluxes are computed as described in section 3. We consider a moderate slope (*α* = 10^{−2}) and set the stratification and flow parameters to the reference values listed in Table 2, representing “typical” values observed in the deeper parts of lakes, fjords, and the ocean, where the BBL is forced by internal waves or internal seiching motions (recall that some restrictions are implied by that fact that we ignore rotational effects). A larger parameter space will be investigated in section 7.

### a. Periodic evolution

Our simulations start from a state of rest with uniform stratification *N*_{∞}, where we chose the reference buoyancy in (2) arbitrarily such that *b* = 0 at the bottom, initially. After the system has reached a fully periodic state, the cyclic evolution of the flow field corresponds to the situation displayed in Fig. 2. At this stage, the boundary layer has reached a thickness of approximately 6 m, where we distinguish between BBL and interior using a threshold of ɛ = 10^{−14} W kg^{−1} (this is the value assigned to nonturbulent regions in the turbulence model). A more standard criterion based on stratification turned out to be impractical because of the strongly variable density structure in the BBL.

Figure 2a illustrates a small lead of the near-bottom velocities with respect to the oscillating interior flow, which increases toward the bottom and mirrors the transition from inertially to frictionally dominated BBL dynamics. Buoyancy is periodically advected up and down the slope, where near-bottom mixing leads to a strong reduction of the buoyancy gradient normal to the slope (Fig. 2b). The combined effects of mixing and restratification lead to a periodic downward bending of near-bottom isopycnals and thus, according to (10), to a persistent internal pressure gradient (Fig. 2c). Although varying periodically, on the average this pressure gradient is of opposite sign in the upper and lower part of the BBL and has a tendency to create a residual flow that is directed upslope near the bottom (Garrett 1991).

Stratification in the BBL is periodically generated and destroyed because of the competing effects of mixing and restratification, with exception of the strong and permanent pycnocline near the upper edge of the BBL (Fig. 2d). Consistent with the discussion after (29), the lower part of the BBL is seen to become statically unstable during upslope flow, when dense fluid is advected on top of lighter fluid because of the strong near-bottom shear (Figs. 2a,d). Note that a similar effect occurring in the strongly sheared transition region between the interior and the BBL also leads to unstable stratification in confined regions, similar to the analytical results by Thorpe (1987). A second, lower pycnocline periodically evolves at the top of the (vigorously turbulent) unstable region and periodically merges with the permanent pycnocline (Fig. 2d). The thickness *d _{u}* of the unstable near-bottom layer may include more than half of the total BBL thickness

*d*

_{BBL}(Fig. 2e).

Along with the inversion of the boundary layer stratification, also the turbulent buoyancy flux periodically changes its sign (Fig. 3a). This alignment of boundary layer stratification and buoyancy flux is a simple consequence of the down-gradient formulation in (12), which is supported for both stable and unstable stratification by the detailed turbulence measurements in the oscillating BBL of a lake by Lorke et al. (2008). The dissipation rate (Fig. 3b) oscillates at twice the frequency of the imposed interior flow and exhibits a strong increase toward the bottom with near-bottom values consistent with the law of the wall (not shown). The values in the lower part of the BBL show a remarkable symmetry between stable and unstable periods, suggesting that the additional production of TKE during unstable periods does not have a dominant effect on the near-bottom energy dissipation. Based on their in situ turbulence measurements, Lorke et al. (2008) have drawn a similar conclusion.

As remarked in the context of (36) and (37) above, the molecular mixing rate *χ _{b}* is closely linked to the net buoyancy flux in the BBL, and it is therefore interesting to study its variability. Figure 3c illustrates that strong mixing occurs mainly during downslope flow, where stable stratification leads to efficient mixing. During upslope flow, strong mixing is confined to the periodically generated pycnocline at the top of the unstable layer; whereas the contribution of the unstable layer itself is negligible (fluid in this layer is already well mixed). This is corroborated by Fig. 3d, demonstrating that the by far strongest contribution to the integral of

*χ*[i.e., to defined in (35)] occurs during downslope flow.

_{b}Figure 4 compares profiles of selected parameters during periods of strongest stable and unstable stratification, separated in time by *T*/2. Figures 4a,c illustrate that the generation of stable stratification during downslope flow reduces the turbulent diffusivity, in spite of the stronger shear compared to the situation with upslope flow. This effect is compensated by an increase in the mixing efficiency such that substantially more net mixing occurs compared to upslope flow (Figs. 3c,d).

### b. Residual flow and mixing efficiency

The residual flow (Fig. 5a) for fully periodic conditions is approximately one order of magnitude smaller than the instantaneous velocities, and it exhibits a vertical structure with a near-bottom upslope transport and a return current of equal magnitude in the upper part of the BBL that is in qualitative agreement with the predictions made by Garrett (1991). The structure of the residual slope-normal buoyancy flux 〈*G*〉 (Fig. 5b) is consistent with structure of 〈*u*〉 implied by (33) (e.g., the minimum in 〈*G*〉 coincides with 〈*u*〉 = 0).

Figure 6a reveals a damped oscillation of the residual volume flux 〈*Q*〉, quickly decaying to zero as the flow becomes fully periodic in accordance with (32). The period of the damped oscillation is close to *T _{c}* = 2

*π*/(

*N*

_{∞}sin

*α*) = 55.2 h, thus corresponding to the resonance frequency of the BBL predicted from (30). Similarly, when fully periodic conditions are reached, the residual buoyancy flux 〈

*Q*〉 becomes stationary and, in agreement with (37), approaches the irreversible limit, 〈〉 (Fig. 6b). Note that 〈〉 approaches the asymptotic value much faster than 〈

_{b}*Q*〉, indicating that irreversible mixing is not sensitive with respect to the realization of perfectly periodic conditions in the BBL. Also, for the larger parameter space investigated in the following section, we found that 〈〉 converges already after a few periods, suggesting that our mixing estimates are robust enough to be valid also for real systems that almost never reach perfectly periodic conditions.

_{b}Interestingly, the magnitude of the total vertical buoyancy flux is smaller than the turbulent vertical flux 〈*Q _{G}*〉 alone (Fig. 6b). From (39), it thus follows that the vertical advective contribution is countergradient, resulting in a substantial reduction of the total flux. For stationary boundary layers, Garrett (1991) has predicted this behavior as a consequence of the interaction between secondary circulation and the density structure in the BBL. Thus, from observations of the turbulent flux alone (e.g., Lorke et al. 2008), the total contribution of boundary mixing is likely to be overestimated.

Because the residual flow has a tendency to restratify the BBL, it is likely that the associated increase of the mixing efficiency at least partly compensates for the countergradient fluxes discussed above. To study such variations, the mixing efficiency is defined here as the ratio of the irreversible vertical buoyancy flux and the rate at which energy is dissipated into heat, both integrated across a horizontal line:

where we have again used *dz* = −*dx̂*sin*α* to convert between horizontal and slope-normal integration. Note that this is different from the notion of “mixing effectiveness” introduced by Garrett (1990) to compare the total vertical flux to the turbulent flux that would occur if *N* = *N*_{∞} everywhere (i.e., if mixing would not reduce stratification in the BBL). Nevertheless, Garrett (1990) points out that also the mixing effectiveness is expected to increase because of the restratifying tendency of the secondary circulation, unless restratification is strong enough to suppress turbulence.

From Fig. 6c, values around 0.5% are found for the mixing efficiency Γ, indicating that for this parameter set boundary mixing is not very efficient. We will see below that this conclusion does not hold in general.

## 7. Parameter space

Extending the results derived in the preceding section to a larger parameter space, we now study the effect of variations in the key nondimensional parameters *α*, *R*, and Ω identified in section 4. We start our investigations by varying the slope angle over the range 10^{−3} ≤ *α* ≤ 10^{−1}, covering typical slopes observed in the ocean and in lakes while keeping the roughness parameter *R* and the frequency ratio Ω fixed at the values given in Table 2 (corresponding to the values used in the previous section). We resolve this range of slopes with 200 simulations, each of which is analyzed after the flow has reached fully periodic conditions. To obtain an estimate for the effect of the roughness parameter, these experiments were complemented by two additional series of simulations over the same range of slope angles, however, with *R* one order of magnitude larger and smaller than the value mentioned in Table 2. For each simulation, we computed the irreversible buoyancy flux 〈〉 from (37). As theoretically shown in (37) and numerically verified in Fig. 6b, for periodic conditions 〈〉 converges to the total vertical buoyancy flux 〈*Q _{b}*〉. The contribution of the upslope turbulent flux in (39) turned out to be negligible in all our simulations. Note that we do not discuss model predictions in the immediate vicinity of the critical value for BBL resonance (gray-shaded areas in Figs. 7, 8) because there (i) fully periodic conditions are reached only after an unrealistically large number of periods and (ii) near-critical reflection of internal waves occurring for

*α*≈

*α*is hardly compatible with the geometric assumptions outlined in section 2.

_{c}As illustrated in Figs. 7a,b, the nondimensional total and turbulent buoyancy fluxes, 〈〉 = 〈*Q _{b}*〉/

*U*

^{3}and 〈〉 = 〈

*Q*〉/

_{G}*U*

^{3}, show a systematic increase for increasing roughness numbers

*R*over the whole spectrum of slopes. However, given that both

*R*and

*α*vary by two orders of magnitude, the overall variability of the nondimensional buoyancy fluxes is surprisingly small with typical values on the order of 〈〉 ≈ 10

^{−3}, except near the critical slope, where both fluxes collapse.

For *α* < *α _{c}*, we find that the magnitude of the total flux 〈〉 is generally smaller than the magnitude of the turbulent flux 〈〉. Similar to the observations made in the context of Fig. 6b, this implies that the vertical advective flux is countergradient and may lead to a substantial reduction of the total flux. Because the two-layer structure of the secondary circulation remains qualitatively similar to that shown in Fig. 5a for subcritical slopes, these results support the suggestion of Garrett (1990, 1991) that the advective contribution of the secondary circulation results in a reduction of the total vertical buoyancy flux. For supercritical slopes,

*α*<

*αc*, however, we observe a reversal of the secondary circulation (not shown), with the consequence that both the turbulent and the advective vertical fluxes change their signs (Figs. 7a,b).

As shown in Figs. 7c,d, the nondimensional dissipation rate 〈*〉 = 〈〉/*U*^{3} and the bulk mixing efficiency Γ defined in (40) are weak functions of the roughness parameter *R* but exhibit strong, opposing trends for increasing slope angle *α*. The strong decrease of the dissipation rate over the range of slopes investigated here can be understood from the following argument: The dissipation rate integrated in the slope-normal direction across the BBL can be expressed as *C _{D}U*

^{3}, where

*C*is a (not necessarily constant) dimensionless drag coefficient. Thus, using (40), we find a simple nondimensional relationship of the form 〈*〉sin

_{D}*α*=

*C*, where the appearance of the trigonometric factor mirrors the conversion from horizontal to slope-normal integration. The drag coefficient exhibits little variability in our simulations, which can be explained from the fact that most of the dissipation rate occurs very close to the bottom in a region that is not strongly affected by stratification (see Fig. 3b). Hence, for constant

_{D}*C*, an inverse relationship between 〈*〉 and the slope angle is expected.

_{D}This is largely consistent with the behavior of 〈*〉 displayed in Fig. 7c, from which we also infer *C _{D}* =

*O*(10

^{−3}), slightly variable for different roughness numbers but well inside the range of accepted values. The mixing efficiency increases from values less than Γ = 10

^{−3}for mild slopes up to efficiencies larger than 1% for slope angles exceeding

*α*= (1–2) × 10

^{−2}. This trend is consistent with the observations of Gloor et al. (2000), who carefully analyzed the energy budget of a small lake and reported increasing mixing efficiencies for steeper slopes.

The role played by the secondary circulation in this process can be understood from averaging the evolution equation for *N* ^{2} in (29) over one period. Assuming periodic conditions, the rate term vanishes, and the shear and mixing terms on the right-hand side must balance. Because negative residual shear, 〈*S*〉 < 0, is seen to have a tendency to generate stable stratification (here in the core of the BBL; see Fig. 5a), a compensation by the second term on the right-hand side of (29) (i.e., by mixing) is required to insure that 〈*N* ^{2}〉 remains constant. However, because the shear term in (29) is proportional to the upslope buoyancy gradient, sin *α*, we expect that the generation of stable stratification occurs at a higher rate for steeper slopes. The magnitude of the secondary circulation, defined here as the (nondimensional) standard deviation of the residual flow = *σ _{u}*/

*U*, is on the order of a few percent of the interior flow and shows relatively little variability for varying slope angles (Fig. 7e). The rate at which stable stratification inside the BBL is generated thus depends to a large extent on the upslope density gradient, which increases for steeper slopes.

The net stabilizing effect of the secondary circulation is mirrored in an increased mixing efficiency (Fig. 7d) but, as shown by our results, also in a suppression of turbulence (Fig. 7c) for steeper slopes. Both effects nearly cancel, and the net impact on the total buoyancy flux (Fig. 7a) is comparatively small. This is one of the key results of this section.

Similarly to the above analysis, we isolate the effects of variable forcing frequency Ω for fixed *R* and *α* (see Table 2) with 200 simulations spanning a range of two orders of magnitude around the reference value for Ω given in Table 2. Values outside this range are either too close to the upper limit for internal waves (*ω* = *N*) to be compatible with the geometrical assumptions made in section 2 or describe wave frequencies that are too small to reach fully periodic conditions on realistic time scales.

The key information obtained from these simulations is that both the buoyancy flux (Fig. 8a) and the mixing efficiency (Fig. 8b) become vanishingly small in the high-frequency limit. This behavior can be understood by recalling from the discussion of Figs. 3c,d that mixing is mainly associated with phases of downslope flow, where shear-induced stable stratification leads to more efficient mixing. For high-frequency forcing, however, increasingly shorter periods of downslope flow inhibit the full development of this process.

The nondimensional thickness, 〈〉 = 〈*d _{u}*〉

*N*

_{∞}/

*U*, of the unstable layer is only a small fraction of the total layer thickness, 〈〉 = 〈

*d*

_{BBL}〉

*N*

_{∞}/

*U*, but convection is seen to occur for all frequencies (Fig. 8c). In fact, in all our simulations we observed periods with unstable boundary layers, suggesting that this phenomenon is an intrinsic component of the BBL dynamics, at least in the nonrotating case. The secondary circulation (Fig. 8e) remains in the range of 1%–3% of the interior flow

*U*over the whole frequency spectrum and is thus comparable to the values observed in Fig. 7e. For supercritical Ω, we find again that the advective buoyancy flux associated with the secondary circulation leads to a reduction of the magnitude of the total flux (Fig. 8a).

## 8. Discussion and conclusions

Consistent with recent observations of mixing on sloping topography in lakes (Lorke et al. 2005, 2008) and numerical simulations of internal tides near ocean ridges (Slinn and Levine 1997), our results demonstrate that the frictional shear caused by oscillating currents may periodically generate gravitationally unstable BBLs. In spite of the additional source of TKE due to convection, the contribution of unstable regions to net mixing in the BBL was found to be negligible because these regions are already rather well mixed. Conversely, strong mixing occurred during periods of downslope flow, when shear-induced periodic stabilization of the boundary layer leads to efficient mixing.

The molecular mixing rate *χ _{b}* of small-scale buoyancy variance was identified as the key parameter determining the net (i.e., advective plus turbulent) vertical buoyancy flux in the boundary layer. This result does not involve any turbulence modeling assumptions, is valid for both rotating and nonrotating flows, and is consistent with the isopycnal view of boundary mixing suggested by Garrett (2001) and the energetic interpretation of mixing by Winters et al. (1995). The parameter

*χ*has the advantage that it can be estimated directly and rather precisely from turbulence microstructure profiles (Lueck et al. 2002; Luketina and Imberger 2001; Gregg 1999) and that it is the only parameter that fully captures the combined effects of turbulence and stratification for mixing.

_{b}From dimensional analysis of the equations of motion we showed that only three nondimensional parameters (slope angle, roughness parameter, and nondimensional forcing frequency) govern the problem for harmonically forced, nonrotating boundary layers. Numerical simulations suggest that the solutions most strongly depend on the slope angle *α* for the physically interesting range 10^{−3} ≤ *α* ≤ 10^{−1}. For this range, the bulk mixing efficiency Γ varies by approximately two orders of magnitude, whereas the nondimensional buoyancy flux 〈〉 ≈ 10^{−3} only exhibits a small variability. The latter result is rather surprising, because it suggests that the net dimensional buoyancy flux is roughly proportional to *U*^{3} and only weakly sensitive to other parameters like the slope angle or the value of the background stratification. This forms one of the main conclusions of the paper.

Using these results and following previous suggestions by Armi (1979a) and Garrett (1979, 1990), it is possible to quantify the basin-scale effect of boundary layer mixing with an “effective” diffusivity *ν*_{eff}. This quantity is computed from the assumption that, at any given depth, the apparent vertical buoyancy flux across a horizontal surface with area *A* corresponds to the vertical buoyancy flux generated inside the BBL on the surrounding topography. If the length of the surrounding depth contour is denoted by *L*, this approach amounts to identifying the apparent vertical flux *ν*_{eff }*A* with the total vertical flux in the boundary layer 〈〉*L*. Expressing the total vertical buoyancy flux in the BBL in terms of the integrated mixing rate as in (37), our model predicts the effective diffusivity as

where all effects associated with the variability of the mixing efficiency, the dissipation rate, and the secondary circulation are implicitly taken into account. This model is applicable also for rotating flows, because the momentum equation is not involved in its derivation. Moreover, following the discussion around (36), it is legitimate to interpret 〈〉 appearing in (41) more generally as the mixing rate averaged over a series of arbitrary mixing events, which are estimated, for example, from microstructure measurements at different locations along a given depth contour.

If such measurements are not available, a similar approach can be used to express *ν*_{eff} in terms of the nondimensional buoyancy flux, 〈〉 = 〈*Q _{b}*〉

*U*

^{−3}, which yields

Using 〈〉 ≈ 10^{−3} (see above), order-of-magnitude estimates for the basin-scale diffusivity can be obtained from (42). Note, however, that the exchange of fluid between BBL and interior is not taken into account in our simple one-dimensional modeling approach. This clearly conflicts with the frequently reported evidence for intrusions, advecting mixed BBL fluid toward the interior near sloping topography (Armi 1978, 1979a; Phillips et al. 1986; Wain and Rehmann 2010). This additional mechanism for boundary layer restratification is likely to lead to more efficient mixing in the BBL (Wain and Rehmann 2010; Gloor et al. 2000; Armi 1979b), and we expect that (42) only provides a lower bound for the basin-scale diffusivity.

It is worth pointing out that Garrett (1990), summarizing different approaches to compute the basin-scale mixing rate, arrived at an expression similar to (42). Assuming that a constant fraction Γ of the vertically integrated dissipation rate, modeled as *C _{D}U*

^{3}, reappears as an increase in basin-scale potential energy, Garrett (1990) suggests

where *A _{b}* is the horizontal area covered by a BBL of thickness

*d*

_{BBL}. In the second step in (43), we have used the geometric relations

*A*=

_{b}*WL*(

*W*is the lateral width of the BBL) and

*d*

_{BBL}=

*W*sin

*α*. For constant

*C*, the expressions in (42) and (43) are consistent only if variations in the bottom slope are compensated by changes in the mixing efficiency Γ. Evidently, this condition cannot be satisfied for any constant Γ. One of the main conclusions is thus that the mixing efficiency increases with the slope angle and is much smaller than the canonical value Γ = 0.2 for interior mixing.

_{D}For typical values of the deep-water stratification ( = 10^{−6} to 10^{−5} s^{−2}), typical internal seiching velocities in lakes (*U* = 5 × 10^{−2} m s^{−1}), assuming a circular topography (*L*/*A* = 2/*r*) with radius *r* = 10 km, (42) predicts *ν*_{eff} = 0.25–2.5 × 10^{−5} m^{2} s^{−1}, using the order-of-magnitude estimate 〈〉 ≈ 10^{−3}. These values include the range of measured basin-scale diffusivities in small to medium-sized lakes (Wüest and Lorke 2003).

Using a similar computation for large-scale ocean basins, however, it is easy to show that the processes discussed here cannot explain the canonical value *ν*_{eff} = 10^{−4} m^{2} s^{−1}. Previous investigations suggest that boundary mixing in the abyssal ocean is more likely to be related to the near-critical reflection of internal waves in regions extending far beyond the well-mixed BBL, where stratification is strong and mixing is efficient (e.g., Klymak et al. 2008; Nash et al. 2007; Lueck and Mudge 1997; Toole et al. 1997). As a simple consequence of the one-dimensional nature of our model, any processes related to internal-wave motions are not taken into account. Finally, the effect of rotation may significantly change the simple pictured developed here by adding the complications of an along-slope flow component and associated upslope or downslope Ekman transports (Middleton and Ramsden 1996; Garrett et al. 1993). Future work is required to clarify the role of these additional processes for the mixing in stratified bottom boundary layers.

## Acknowledgments

This work was funded by the German Research Foundation (DFG) as part of the ShIC project (UM79/3-1). We are grateful for the helpful and constructive comments by Chris Garrett and two anonymous reviewers. Computations have been performed with a modified version of the turbulence modeling tool box GOTM (available online at http://www.gotm.net).

## REFERENCES

**,**(

**,**

**,**(

**,**

**,**

**,**

**,**

**,**(

**,**(

**,**

**,**

**,**

**,**(

**,**(

**,**(

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**(

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

## Footnotes

*Corresponding author address:* Lars Umlauf, Leibniz-Institute for Baltic Sea Research, Seestrasse 15, 18119 Warnemünde, Germany. Email: lars.umlauf@io-warnemuende.de