Vertical Shear, Diapycnal Shear and the Gradient Richardson Number

Abstract

In Cartesian coordinates $\left(x,y,z\right),$ the gradient Richardson number $Ri$ is the ratio between the square of the buoyancy frequency N and the square of the vertical shear S, $Ri=N^2/S^2$, where $N^2=-\left(g/\rho \right)\partial \rho /\partial z$ and $S^2={\left(\partial u/\partial z\right)}^2+{\left(\partial v/\partial z\right)}^2$, with ρ potential density, $\left(u,v\right)$ the horizontal velocity components and g gravity acceleration. In isopycnic coordinates $\left(x,y,\rho \right)$, $Ri$ is expressed as the ratio between ${M^2\equiv N}^{-2}$ and the squared diapycnal shear ${S_{\rho }}^2={\left(\rho /g\right)}^2\left[{\left(\partial u/\partial \rho \right)}^2+{\left(\partial v/\partial \rho \right)}^2\right]$, $Ri=M^2/{S_{\rho }}^2$. This could suggest that a decrease (increase) in stratification brings a decrease (increase) in dynamic stability in Cartesian coordinates, but a stability increase (decrease) in isopycnic coordinates. The apparently different role of stratification arises because S and $S_{\rho }$ are related through the stratification itself, $S_{\rho }=S/N^2$. In terms of characteristic times, this is equivalent to $\tau \equiv S_{\rho }={t_o}^2/t_d$, which is interpreted as a critical dynamic time $\tau$ that equals the buoyancy period $t_o\equiv N^{-1}$ normalized by the ratio $t_d/t_o$, where $t_d=S^{-1}$ is the deformation time. Here we follow simple arguments and use field data from three different regions (island shelf break, Gulf Stream and Mediterranean outflow) to endorse the usefulness of the isopycnal approach. In particular, we define the reduced squared diapycnal shear ${{\sigma }_{\rho }}^2={S_{\rho }}^2-M^2$ and compare it with the reduced squared vertical ${\sigma }^2=S^2-N^2$, both being positive (negative) for unstable (stable) conditions. While both Ri and ${\sigma }^2$ remain highly variable for all stratification conditions, the mean ${{\sigma }_{\rho }}^2$ values approach ${S_{\rho }}^2$ with increasing stratification. Further, the field data follow the relation ${{\sigma }_{\rho }}^2=\left(1-Ri\right)/\left(N^{2}Ri\right)$, with a subcritical $Ri=0.22$ for both the island shelf break and the Mediterranean outflow. We propose ${{\sigma }_{\rho }}^2$ and ${S_{\rho }}^2$ to be good indexes for the occurrence of effective mixing under highly stratified conditions.

1. Introduction

The analysis of small flow perturbations, for the case of uniformly stratified and sheared flow, goes back to Geoffrey Taylor’s 1915 essay for the Adams Prize [1], later followed by many other researchers [2,3,4,5]. Other authors employed energetic arguments to assess under which conditions the energy of fluctuations will increase [6,7,8], for a brief historical review see [9]. The conclusion from these works is that the growth of flow perturbations may be expressed in terms of the (gradient) Richardson number, Ri, with perturbations remaining stable when the Richardson number is above some critical value ${Ri}_c$ (supercritical, $Ri>{Ri}_c$) but becoming unstable when it goes below (subcritical, $Ri<{Ri}_c$). Further, the nonlinear stability analysis of three-dimensional stratified shear flow leads to a critical value of one [5], in agreement with energy arguments [8]. Despite this, there are still some theoretical and observational uncertainties on the actual Ri threshold value for the development of turbulence [10,11].

The most common representation for the Richardson number comes in the Cartesian $\left(x,y,z\right)$ coordinate system,

$$Ri={\displaystyle \frac{-\left(\frac{g}{\rho }\right)\frac{\partial \rho }{\partial z}}{{\left(\frac{\partial u}{\partial z}\right)}^2+{\left(\frac{\partial v}{\partial z}\right)}^2}}={\displaystyle \frac{N^2}{S^2}}$$

(1)

with ρ the potential density and $\left(u,v\right)$ the horizontal velocity components in the horizontal $\left(x,y\right)$ directions. Both the vertical-oscillation or buoyancy frequency $N={\left[-\left(g/\rho \right)\partial \rho /\partial z\right]}^{1/2}$ and vertical shear $S={\left[{\left(\partial u/\partial z\right)}^2+{\left(\partial v/\partial z\right)}^2\right]}^{1/2}$ have frequency units.

An alternative, possibly more natural representation of the Richardson number is in isopycnic coordinates. We apply the chain rule $\partial u/\partial z=\left(\partial u/\partial \rho \right)\left(\partial \rho /\partial z\right)$ as well as the inverse function rule $\partial \rho /\partial z={\left(\partial z/\partial \rho \right)}^{-1}$ and introduce the $\left(\rho /g\right)$ factor to obtain

$$Ri={\displaystyle \frac{-\left(\frac{\rho }{g}\right)\frac{\partial z}{\partial \rho }}{{\left(\frac{\rho }{g}\frac{\partial u}{\partial \rho }\right)}^2+{\left(\frac{\rho }{g}\frac{\partial v}{\partial \rho }\right)}^2}}={\displaystyle \frac{M^2}{{S_{\rho }}^2}}$$

(2)

where ${M\equiv N}^{-1}={\left[-\left(\rho /g\right)\partial z/\partial \rho \right]}^{1/2}$ is the vertical-oscillation or buoyancy period and $S_{\rho }=\left(\rho /g\right){\left[{\left(\partial u/\partial \rho \right)}^2+{\left(\partial v/\partial \rho \right)}^2\right]}^{1/2}$ is the diapycnal shear, both with units of time. We introduce the “diapycnal shear” terminology to emphasize that the derivatives are with respect to density; when multiplied by $\left(\rho /g\right)$, the result has time units.

As a result of the above change of coordinates, a paradox arises. From the Cartesian perspective, Ri is directly proportional to stratification, while from the isopycnic viewpoint, Ri is inversely proportional to stratification [12,13,14]. The answer lies in the fact that the vertical and diapycnal shear are related through stratification, $S=N^2{S}_{\rho }$. For constant vertical shear, an increase in stratification enhances dynamic stability. In contrast, for constant diapycnal shear, the same increase in N actually reduces dynamic stability as it implies a quadratic increase in vertical shear. We are undoubtedly most used to the vertical perspective, so perhaps we are missing some lessons that may be learned from the isopycnic view. In this work, we will pursue such an isopycnic inspection, by means of both simple conceptual arguments and through the analysis of real data for different dynamic situations.

The appearance and extensive use of instruments that sense small-scale turbulence, hence providing data that can be translated into turbulent dissipation rates and vertical diffusivities, has progressively reduced the early emphasis on the Richardson number as an index for vertical instability. However, there are still many observational and numerical studies that, in the absence of microstructure measurements, resort to hydrographic-velocity data and the standard parameterizations of vertical diffusivity in terms of the Richardson number [15,16,17,18,19]. Hence, it seems timely to follow other authors [18,20,21,22,23,24] and further explore Ri-related alternative indexes that could possibly be used to assess the presence of vertical instability.

In the next section, we explore some relevant differences between the vertical and isopycnic perspectives, and introduce the concepts of natural oscillation time, reduced squared vertical shear and reduced squared diapycnal shear. We present the non-dimensional variables in Section 3, examining the interdependence between background stratification, vertical shear and diapycnal shear. The three datasets used in our study are briefly described in Section 4, and in Section 5 we look at these data in two different ways: first as cloud points plotted in terms of two variables (among background stratification, vertical shear and diapycnal shear) and next examining the dependence of both the reduced squared vertical shear and the reduced squared diapycnal shear in terms of the background stratification. Finally, in Section 6, we summarize the main conclusions.

2. Isopycnic Versus Vertical Perspectives

2.1. Basic Considerations

In vertical coordinates, the critical Richardson number criterion is commonly interpreted as if sufficiently high stratification (high $N$) leads to stable conditions ($Ri>1$), provided the vertical shear remains moderate, while low stratification is prone to unstable conditions ($Ri<1$) (Equation (1)). In a layered ocean, with near-constant density layers bounded by relatively thin steps of rapid density changes, this vertical-coordinate perspective would usually associate the steps with stable regions and the layers with unstable zones (Figure 1). This is true for a static ocean, but dynamic stability depends on the vertical shear. A layer will have most properties well mixed so that both N and S tend to zero; weak vertical shear will lead to instabilities, but these will only redistribute near-homogeneous water. In a step, on the other hand, N is large and instabilities can effectively redistribute density (and other water properties), but these will only develop for large enough vertical shear S.

Equation (2) in isopycnic coordinates, on the other hand, suggests that a highly stratified region (a step, $M=N^{-1}\to 0$) may be prone to mixing as long as the diapycnal shear $S_{\rho }$, which depends on the velocity of water parcels in adjacent layers, remains finite. Consider for example a 2D symmetrically stable front with the density field $\rho \left(x,z\right)$, and further consider that the along-front flow is in geostrophic balance. The thermal wind equation in Cartesian and isopycnic coordinates are $\rho f\partial v/\partial z=g\partial \rho /\partial x$ and $\rho f\partial v/\partial \rho =g\partial z/\partial x$, respectively, and the corresponding Richardson numbers are $Ri=C\left(\partial \rho /\partial z\right)/{\left(\partial \rho /\partial x\right)}^2$ and $Ri=C\left(\partial z/\partial \rho \right)/{\left(\partial z/\partial x\right)}^2$, where $C=-\left(f^2\rho /g\right)$. In a frontogenetic situation where the isopycnals get steeper with time (both $\partial \rho /\partial x$ and $\partial z/\partial x$ increase) and vertical stratification does not change ($\partial \rho /\partial z$ remains constant), Ri will decrease in both reference frames. However, in those locations where frontogenesis also leads to increased vertical stratification, the Cartesian expression suggests that stability may increase while the isopycnic perspective points to further instability.

The above simple arguments suggest that diapycnal shear may be a convenient variable to analyse whether sheared oceanic flows are prone to mixing. As an illustrative example, consider the two-layer flow between parallel planes such as that set up in the laboratory to simulate interfacial Kelvin–Helmholtz instabilities, e.g., [25] and references therein. The upper and lower layers have densities $\rho ={\rho }_0\left(1-∆\rho \right)$ and $\rho ={\rho }_0\left(1+∆\rho \right)$, respectively, and the corresponding velocities are $u=∆\rho gtsin\theta$ and $-∆\rho gtsin\theta$, where θ is the tilting of the table in the experimental setup. The vertical shear between both layers is not uniquely defined but, in contrast, the isopycnic derivative of the horizontal velocity is easily calculated from $\rho ={\rho }_0\left(1-u/\left(gtsin\theta \right)\right)$, as $\partial u/\partial \rho =-\left(gtsin\theta \right)/{\rho }_0$; notice this derivative is equal to the finite-difference fraction, $\delta u/\delta \rho =\left(2gt∆\rho sin\theta \right)/\left(2{\rho }_0∆\rho \right)=\left(gtsin\theta \right)/{\rho }_0$. The resulting expression for the diapycnal shear is simply $S_{\rho }=tsin\theta$, with time units, implying that instabilities will develop independently of stratification if sufficient time is allowed.

2.2. Characteristic Times

As mentioned in the Introduction, the Richardson number criterion arises from either flow-stability or energy arguments. A complementary interpretation comes in terms of characteristic times associated with the inverse of the buoyancy frequency and the vertical shear, as well as with the diapycnal shear itself. The inverse of the buoyancy frequency is the time characteristic for the vertical oscillation of a water parcel in a stratified medium, hereafter the vertical-oscillation or buoyancy period, $t_o\equiv N^{-1}$. The inverse of the vertical shear gives a characteristic time for the horizontal deformation of a water parcel (or temporal rate of strain), henceforth the vertical-deformation time, $t_d\equiv S^{-1}$.

The diapycnal shear is proportional to the derivative of the horizontal velocity with respect to density, e.g., $\partial u/\partial \rho$; when we multiply this derivative by $\rho /g$, we get the diapycnal shear with units of time, so we set $\tau \equiv S_{\rho }$ to simply remind us of its time units. Because the vertical and diapycnal shears are related through the stratification, $S_{\rho }=S/N^2$, the diapycnal shear is related to the buoyancy and deformation times as $\tau ={t_o}^2/t_d$. We envision the diapycnal shear as reflecting adjacent layers that maintain their velocity, so it naturally incorporates both stratification (through $t_o)$ and shear deformation (through $t_d$). We may hence interpret τ as a critical dynamic time for the system, equal to the vertical-oscillation period $t_0$ normalized by the ratio $\left(t_d/t_o\right)$.

Equations (1) and (2) may be written in terms of these characteristic times as

$$Ri={\displaystyle \frac{N^2}{S^2}}={\left({\displaystyle \frac{t_d}{t_o}}\right)}^2$$

(3)

$$Ri={\displaystyle \frac{1}{{N^2{S}_{\rho }}^2}}={\left({\displaystyle \frac{t_0}{\tau }}\right)}^2$$

(4)

or alternatively as

$$Ri={\displaystyle \frac{1}{S{S}_{\rho }}}={\displaystyle \frac{t_d}{\tau }}$$

(5)

We may now see that Equations (3) and (5) are not really that different. According to Equation (3), the Richardson number is subcritical (less than one) when $t_d<t_o$, i.e., when the vertical excursion of the water parcel takes long enough for it to become greatly distorted. Similarly, Equation (5) tells us that subcritical conditions require a short deformation time $t_d$ as compared with the critical dynamic time τ, $t_d<\tau$. The main difference lies on the linear and quadratic dependences: in vertical coordinates the dependence is quadratic ${\left(t_d/t_0\right)}^2$ while in isopycnal coordinates it is linear $\left(t_d/\tau \right)$. Equations (3) and (4) further tell us that the instability condition implies $t_d<t_0<\tau$.

2.3. Reduced Squared Shears

The reduced vertical shear was introduced as a complementary way to understand the relevance of subcritical conditions for mixing [20,21]. By considering a critical value $Ri=1/4$, the subcritical condition would correspond to $S^2-4N^2=\left(S+2N\right)\left(S-2N\right)>0$ or, alternatively, to a reduced shear greater than zero, $S-2N>0$. A formal advantage of the reduced shear is that it is approximately proportional to the growth rate of Kelvin–Helmholtz instabilities [20,26]. A practical advantage is that reduced shear helps assess when effective mixing occurs. As mentioned above, growing instabilities are effective blenders in well-stratified conditions; in poorly stratified conditions, however, they can only redistribute water which is already well mixed. The Richardson number Ri is not capable of distinguishing cases of low stratification and low shear from those of high stratification and high shear, but the reduced shear does, i.e., for two different situations (high and low stratification) with equal Ri subcritical values, it turns out that $S-2N=N\left({Ri}^{-1/2}-2\right)$ is larger for the well-stratified condition.

Let us pursue these ideas and define a reduced squared vertical shear ${\sigma }^2$ and a reduced squared diapycnal shear ${{\sigma }_{\rho }}^2$, but using the critical Richardson value of one [5]:

$${\sigma }^2\equiv S^2-N^2=S^2\left(1-Ri\right)=N^2\left(1-Ri\right)/Ri$$

(6)

$${{\sigma }_{\rho }}^2\equiv {S_{\rho }}^2-M^2={S_{\rho }}^2\left(1-Ri\right)=\left(1-Ri\right)/\left(N^{2}Ri\right)$$

(7)

The definitions in (6) and (7) show that both squared shears are positive for unstable conditions ($Ri<1$) and negative for stable conditions ($Ri>1$). For constant shears, ${\sigma }^2$ decreases with stratification from $S^2$ to increasingly negative values, and ${{\sigma }_{\rho }}^2$ increases with stratification up to a maximum value of ${S_{\rho }}^2$. Equations (6) and (7) illustrate how, for constant Ri values, the reduced squared shears change with either stratification or vertical/diapycnal shear. Both reduced squared shears are equal to the corresponding squared shear but with a $\left(1-Ri\right)$ reduction factor, i.e., the reduced squared shears increase with decreasing Ri up to a maximum value equal to the squared shear, with ${\sigma }^2$ approaching $S^2$ and ${{\sigma }_{\rho }}^2$ approaching ${S_{\rho }}^2$. When expressed as a function of stratification $N^2$, both variables increase as Ri decreases but with ${\sigma }^2$ and ${{\sigma }_{\rho }}^2$ directly and inversely proportional to stratification, respectively.

3. Data Analysis

To see how the above ideas fit in the real world, we consider several different cases of stratification and vertical/diapycnal shear. In order to compare these cases, it is convenient to consider the nondimensional form of the dependent variables in Equations (3)–(7). For this purpose we use a background stratification, $N_b$, to set the following relations: $N=N_b\left(N/N_b\right)\equiv N_b{N}^{\prime }$, $S=N_b\left(S/N_b\right)\equiv N_b{S}^{\prime }$, $M=\left(1/N_b\right)\left(N_{b}M\right)\equiv \left(1/N_b\right)M^{\prime }$ and $S_{\rho }=\left(1/N_b\right)\left(N_b{S}_{\rho }\right)\equiv \left(1/N_b\right){S_{\rho }}^{\prime }$; notice there are only two independent variables as $M^{\prime }=1/N^{\prime }$ and ${S_{\rho }}^{\prime }=S^{\prime }/{N^{\prime }}^2$.

Hereafter, we drop primes for all dependent variables and will always refer (except where indicated) to their nondimensional forms. In particular, notice that $N^{-1}$, $S^{-1}$ and $S_{\rho }$, respectively, are the nondimensional forms of the vertical-oscillation period, the deformation time and the critical dynamic time. The Richardson number may hence be expressed in terms of the nondimensional variables exactly as in Equations (3)–(5). Analogously, the nondimensional reduced squared shears have the same form as in Equations (6) and (7).

The subcritical condition (Ri < 1) translates into the following instability conditions for the reduced squared shears:

$${\sigma }^2=\left(S^2-N^2\right)>0$$

(8)

$${{\sigma }_{\rho }}^2=\left({S_{\rho }}^2-{\displaystyle \frac{1}{N^2}}\right)>0$$

(9)

To illustrate the relations among the nondimensional dependent variables, we have produced several sets of plots, where we explore the relations between two variables as a function of the third variable and the Richardson number (Figure 2). Specifically, we look at the following dependences:

-

in the $\left(S,N\right)$ domain, $S=S_{\rho }{N}^2=\left(1/{Ri}^{1/2}\right)N$,

-

in the $\left(S_{\rho },N\right)$ domain, $S_{\rho }=S/N^2=\left(1/{Ri}^{1/2}\right)\left(1/N\right)$, and

-

in the $\left(S_{\rho },S\right)$ domain, $S_{\rho }=S/N^2=\left(1/Ri\right)\left(1/S\right)$.

In all these plots, we show an overall perspective that goes all the way to stratification values as large as 10 times the background stratification (Figure 2, upper panels), and display a close-up on that region with stratification values between well mixed and twice the background stratification (Figure 2, lower panels). As justified above, we are most interested in what leads to mixing in relatively well-stratified conditions, e.g., for N $\in \left[1,2\right]$.

In the left panels of Figure 2, the domain is halved in regions of supercritical and subcritical conditions, illustrating that subcritical conditions may only be reached if the vertical shear is sufficiently large, $S>N$. For low stratification values, subcritical conditions ($Ri\le 1$) can only be attained for relatively large diapycnal shears (Equation (1)), e.g., $S_{\rho }\ge 2$ when $N\le 0.5$; however, as stratification increases, we find that subcritical conditions are achieved with a progressively smaller diapycnal shear, e.g., when $N\ge 2$ it is sufficient to have $S_{\rho }\le 0.5$.

The middle panels of Figure 2 show a quite different partition between subcritical and supercritical conditions, which is the result of the quadratic inverse dependence of Ri on both stratification and diapycnal shear (Equation (2)). Only a small portion of the domain corresponds to supercritical conditions, tending to zero as either $N$ or $S_{\rho }$ become large. In particular, for well-stratified conditions, the domain of subcritical conditions reaches moderate and even low diapycnal shears, necessarily coincident with high values of vertical shear.

Figure 2. Relations among nondimensional dependent variables. (Left panels) Vertical shear as a function of stratification, (Middle panels) diapycnal shear as a function of stratification, (Right panels) diapycnal shear as a function of vertical shear, for different Ri values (red curves) and different values of the third variable (blue curves) as shown; the dashed, solid and dotted lines, respectively, correspond to values less than, equal to and greater than one. In the top panels, we present the results in the [0:10, 0:10] domain, while in the bottom panels we zoom into the [0:2, 0:2] domain; the colored arrows illustrate the direction for decreasing Ri (red curves in all panels), increasing $S_{\rho }$ (blue in left panels), increasing $S$ (blue in middle panels) and increasing $N$ (blue in right panel).

Figure 2. Relations among nondimensional dependent variables. (Left panels) Vertical shear as a function of stratification, (Middle panels) diapycnal shear as a function of stratification, (Right panels) diapycnal shear as a function of vertical shear, for different Ri values (red curves) and different values of the third variable (blue curves) as shown; the dashed, solid and dotted lines, respectively, correspond to values less than, equal to and greater than one. In the top panels, we present the results in the [0:10, 0:10] domain, while in the bottom panels we zoom into the [0:2, 0:2] domain; the colored arrows illustrate the direction for decreasing Ri (red curves in all panels), increasing $S_{\rho }$ (blue in left panels), increasing $S$ (blue in middle panels) and increasing $N$ (blue in right panel).

Finally, the right panels in Figure 2 show the inverse linear dependence of Ri on both vertical and diapycnal shear (Equation (3)), with subcritical conditions requiring sufficiently large values of either one or both variables, depending on the size of vertical stratification. In particular, subcritical flow requires only moderate diapycnal shear during high stratification conditions.

The principal idea arising from Figure 2 is that moderate diapycnal shear in well-stratified conditions is a guarantee for actual mixing of distinct waters. We may further explore this idea by examining the dependence of the nondimensional reduced squared shears as a function of the nondimensional stratification $N$ for different $S$, $S_{\rho }$ and Ri values:

-

in the $\left({\sigma }^2,N\right)$ domain, ${\sigma }^2=S^2-N^2=\left({S_{\rho }}^2{N}^2-1\right)N^2=N^2\left(1-Ri\right)/Ri$, and

-

in the $\left({{\sigma }_{\rho }}^2,N\right)$ domain, ${{\sigma }_{\rho }}^2={S_{\rho }}^2-1/N^2=\left(S^2-N^2\right)N^4=\left(1-Ri\right)/\left(N^{2}Ri\right)$.

The corresponding plots are shown in Figure 3, again up to $N=10$ (upper panels) and for $N\in \left[0,2\right]$ (lower panels). Recall that subcritical conditions correspond to positive values of either reduced squared shear.

The plots for the reduced squared vertical shear (Figure 3, left panels) illustrate that the same Richardson number may correspond to many different ${\sigma }^2$ values, each of them with a different combination of $N$, $S$ and $S_{\rho }$ values (except Ri = 1, which corresponds to ${\sigma }^2=0$). In particular, positive values of ${\sigma }^2$ are obtained under highly stratified conditions for moderate and even low diapycnal shear but require high values of vertical shear; for example, a value as low as $S_{\rho }=0.5$ guarantees the subcritical conditions for $N\ge 2$, implying $S\ge 2$.

The plots for the reduced squared diapycnal shear (Figure 3, right panels) serve to emphasize the idea that, under well-stratified conditions, vertical instability and mixing will only occur if the vertical shear is very high, but even a small diapycnal shear can cause such stirring. In particular, for relatively large stratifications, the approximation ${{\sigma }_{\rho }}^2={S_{\rho }}^2-1/N^2\cong {S_{\rho }}^2>0$ holds always, which suggests that the diapycnal shear is a good stability index. In these conditions, the critical dynamic time is large enough to guarantee that the flow deformation will create instability.

In the following sections, we will use three very different datasets to explore the relevance of the isopycnic perspective. Our approach will be empirical: we will display scatter plots of the sort shown in Figure 2 and Figure 3 and will explore the different Ri-regimes in each dynamical system.

4. Datasets

The three different datasets represent very diverse oceanographic settings and dynamic conditions: (1) the shelf break of Gran Canaria Island (GCI), a deep-ocean island where internal waves are commonly present over the slope and shelf break (Figure 4); (2) the Gulf Stream (GS), an intense western boundary baroclinic current (Figure 5); and (3) the Mediterranean outflow (MO), a density-driven intense jet (Figure 6). All three datasets share the characteristic that a significant fraction of the data corresponds to well-stratified waters (Figure 7) which are experiencing an intense flow. Further, in all cases the velocity field is actually sampled, therefore allowing a good assessment of vertical and diapycnal shears. Therefore, these datasets may be quite useful for examining the relationships between stratification, vertical shear and diapycnal shear and to explore how the Richardson number and the reduced squared (vertical and diapycnal) shears change with stratification.

The first set of data corresponds to repeated vertical sampling at a fixed location on the 100 m isobath (yo-yo type) in the shelf break south of GCI (Figure 4, left panels) in fall 1994 [27]. Sampling consisted of a conductivity–temperature–depth (CTD) probe which included an acoustic Doppler point current meter. Density and velocity profiles were obtained at 1 m vertical resolution every hour during nearly four days, two days of neap tides and two days of spring tides, for a total of about 9000 samples.

Gran Canaria Island has a steep slope (about 0.15) that goes down to depths over 2500 m (Figure 4, bottom left). The dominant semidiurnal tide propagates northwards into the slope and, as a result of the interaction with the topography, packets of internal waves with semidiurnal and lower periodicities are generated. The structure of the water column consisted of a 50–70 m surface mixed layer, a 20–30 m thick seasonal thermocline where the potential density increased from 1025.5 to 1026.0 kg m⁻³, and a relatively thin bottom mixed layer (Figure 4). The buoyancy frequency reflects the fairly well-mixed region in the upper half of the water column, on top of varyingly stratified waters in the bottom half, with a mean buoyancy frequency equal to 0.010 s⁻¹ (Figure 7, left panel). The cross-slope speed often exceeds 0.2 m s⁻¹ over the entire water column (Figure 4, bottom right), leading to subcritical conditions not only in the mixed layers but also in the thermocline [24].

The second dataset consists of 20 repeated hydrographic sections across the GS (Figure 5, left panel), carried out between September 1980 and May 1983. These sections include temperature data obtained with expendable bathythermographs (XBT) and velocity data obtained from free-falling Pegasus instruments [31], originally reported and analyzed by [28], see also [29] and references therein. Temperature values were converted into potential density by means of a recursive utilization of algorithms [32,33]. Each section contained between four and ten stations of velocity and potential density from the sea surface down to 2000 m at 25 m intervals, for a total of some 10,000 samples.

The permanent thermocline of the GS, with temperatures between about 6 and 16 °C, slopes sharply towards the continental slope, rising some 600–700 m in 100 km. As a result, the GS core velocities increase from about 0.25 to more than 1.5 m s⁻¹ over this same temperature range. On both sides of the GS core, the mean vertical stratification remains similar but the velocities are substantially smaller (Figure 5, middle panel). The relatively low vertical data resolution limits the sampling of high-shear events associated to internal waves, which would be responsible for a significant increase in the gradient Richardson values, and hence restricts our analysis to near-geostrophic flow [9]. Despite this, the background geostrophic velocity field is intense enough for the core of the GS to experience several instances of subcritical conditions (Ri < 1) which are responsible for vertical mixing [12,34]. The buoyancy frequency decreases largely in the top 300 m and remains approximately constant down to 1000 m, with low variability characteristics for the subtropical permanent thermocline. The mean buoyancy frequency for this second dataset is 0.004 s⁻¹ (Figure 7, middle panel).

The third set of data consists of 56 stations with full-depth CTD profiles and simultaneous velocity measurements with expendable current profilers (XCP) that sampled the eastern Gulf of Cádiz during fall 1988 [30,35,36], with about half of the stations over the MO (Figure 6, left panel). Velocity and density data were recorded every 2 m down to an average depth of about 500 m so that some 14,000 measurements were collected.

As the MO exits the Strait of Gibraltar, it experiences three different dynamic regions over a distance of less than 50 km. In the first phase, located between Spartel Sill (360 m) and Western Spartel Sill (420 m, located at about 6°20′W, 35°47′N), the MO follows west along a relatively narrow (about 5 km wide) and gently sloping channel that deepens from 360 to 420 m in about 20 km. In the second phase, the MO continues west along a channel of variable amplitude (typically 5 to 10 km) and undergoes three abrupt topography-driven accelerations (the bathymetry deepens from 420 to 700 m over a distance of another 20 km) with core velocities in excess of 1 m s⁻¹, resulting in substantial dilution of the Mediterranean water properties. During these two phases, the MO approximately occupies the bottom third of the water column (Figure 6, right panels). In the third phase, the MO follows northwest, its diluted core slowing down and leaning against the Iberian continental slope [30,37,38,39]. Therefore, the local profile of the buoyancy frequency depends on whether the station is close to Gibraltar and whether it is found on the path or away from the path of the MO; for this reason, the variability is largest at depths between 200 and 750 m, and decreases in the bottom 100 m simply because of the few stations we have for deep waters. The mean buoyancy frequency for this third dataset is 0.004 s⁻¹ (Figure 7, right panel).

5. Results

We nondimensionalize the variables with a different constant background frequency for each dataset, as explained in Section 3. This background frequency, which is defined as the mean buoyancy frequency over each region (

Table 1. Dimensional mean buoyancy frequency N (s⁻¹) and best-fit parameters to the nondimensional mean ${{\sigma }_{\rho }}^2$ values for two different functions, with the $R^2$ correlation coefficients.

		Gran Canaria Island	Gulf Stream	Mediterranean Outflow
Number of samples		9000	10,000	14,000
Mean N (s−1)		0.010	0.004	0.004
best exponential adjustment ${{\sigma }_{\rho }}^2={{\sigma }_0}^{2}exp\left(-rN\right)$	${{\sigma }_0}^2$	174.2	−1.9	250.3
	r	2.90	1.05	4.03
	$R^2$	0.66	0.83	0.97
best adjustment as ${{\sigma }_{\rho }}^2=\left(1-Ri\right)/\left(N^{2}Ri\right)$	Ri	0.217	2.067	0.222
	$R^2$	0.70	0.86	0.81

), allows consideration of the dynamic changes within each single region. Additionally, the utilization of different background stratifications (one per region) facilitates a unified interpretation of all three, potentially quite different, cases.

Consider first the cloud points for the nondimensional vertical stratification, vertical shear and diapycnal shear (Figure 8), in particular the zoom for the high stratification and high shear domain (Figure 9). The overall shape of S against N changes greatly in all three cases, with no noticeable dependence of vertical shear on stratification. The common denominator in all three cases is that S takes a relatively high range of values for any N (Figure 8 and Figure 9, top panels). In contrast, all $S_{\rho }$ against N cloud points have similar shapes, with relatively high values and high variability in low stratification conditions and the opposite for well-stratified waters (N > 3) (Figure 8 and Figure 9, middle panels); these cloud points actually suggest a decay of the maximum diapycnal shears with increasing vertical stratification. Finally, the $S_{\rho }$ against S scatter plots also display very large differences between all three cases, with the GCI and GS cases displaying completely different distributions and the MO possibly as an intermediate situation (Figure 8 and Figure 9, bottom panels).

These results may be interpreted as if each flow has its own dynamics, setting a range of characteristic diapycnal shears. We may step forward and explore if each dynamical system is characterized by a controlling parameter. In Section 2, we argued that the reduced squared diapycnal shear ${{\sigma }_{\rho }}^2$ may be such parameter. This variable is an index for the stability of the system, in the same way as the Richardson number, but with the characteristic that its maximum value is the squared diapycnal shear ${S_{\rho }}^2$, to be attained during high stratification conditions. These ideas are sustained by the cloud points of ${\sigma }^2$ and ${{\sigma }_{\rho }}^2$ as a function of stratification (Figure 10 and Figure 11). The distribution of ${\sigma }^2$ is quite different in each case, with large scattering for all N values. In contrast, the distribution of ${{\sigma }_{\rho }}^2$ does hint at a dependence on N, with its maximum (positive) values decreasing with increasing vertical stratification, with a threshold value ${S_{\rho }}^2$ that is itself conditioned by N, i.e., ${{\sigma }_{\rho }}^2={S_{\rho }}^2-N^{-2}=S^2{N}^{-4}-N^{-2}$.

In order to further examine these ideas, we look at the way the mean values and standard deviations of all three stability indices (Ri, ${\sigma }^2$ and ${{\sigma }_{\rho }}^2$) change with stratification (Figure 12; for this calculation the mean values and standard deviations are calculated over 0.2 intervals of N). The results show no clear relation of ${\sigma }^2$ and Ri with stratification, but suggest the existence of a univocal dependence of ${{\sigma }_{\rho }}^2$ with stratification for each particular flow, with ${{\sigma }_{\rho }}^2$ proportional to $N^{-1}$ in GCI and the MO, and ${{\sigma }_{\rho }}^2$ proportional to ${-N}^{-1}$ in the GS (Figure 13).

Finally, we explore if there are simple functions that can provide a good fit to the observed mean ${{\sigma }_{\rho }}^2$ dependence on stratification N (Figure 13; mean values are here calculated over 0.1N intervals). Our first attempt involves simple exponential decay with stratification, given by ${{\sigma }_{\rho }}^2={{\sigma }_0}^{2}exp\left(-rN\right)$. It turns out that the adjustment is fairly good in all cases, with the best-fit parameters (as calculated separately for each dataset) presented in Table 1; because of the excellent fit for low and high N values, the highest correlation with the data corresponds to the MO case. For GCI and the MO, the exponential coefficient is relatively large (2.90 and 4.03) as compared with the one for the GS (1.05); the main difference, however, is that for GCI and the MO the mean ${{\sigma }_{\rho }}^2$ values remain always positive, indicative of unstable conditions, while for the GS they are negative, indicative of stable conditions for shear instabilities. The negative values in the GS are likely the consequence of the low vertical resolution of these data, which cannot resolve the relevant turbulent scales.

An equally good ${{\sigma }_{\rho }}^2$ fit corresponds to the Ri-dependence given by Equation (7). We may again appreciate a similar behavior for GCI and the MO, quite different from what happens for the GS. Both former cases are best adjusted with a subcritical Ri = 0.22, while the GS reflects a supercritical flow Ri = 2.1 (Table 1). It is remarkable that this happens despite the mean Richardson values being almost always significantly greater than one (top panels of Figure 12). Our interpretation is that, despite all the variability reflected by the scattered points in Figure 11, each case responds to a different flow dynamic. This dynamic setting is characterized by a ${{\sigma }_{\rho }}^2\left(N\right)$ relation with a characteristic flow-dependent Ri value, in agreement with the idea that strongly stratified flows self-organize around some local critical state [40].

6. Conclusions

In this work, we have endorsed the advantages of isopycnic thinking to better understand the dynamic conditions leading to effective mixing in stratified sheared flows. A nondimensional analysis, in terms of the background stratification, proves to be useful to explore the interdependences between stratification (expressed in terms of the buoyancy frequency $N$), vertical shear $S$, diapycnal shear $S_{\rho }$ and the (gradient) Richardson number $Ri$. In particular, we have discussed vertical instability in terms of characteristic times, showing that instability occurs when the critical dynamic time $\tau$—which is equal to the diapycnal shear and hence includes both stratification and vertical shear, $\tau \equiv S_{\rho }={t_o}^2/t_d$—is longer than the vertical-oscillation period $t_0=N^{-1}$, which itself is longer than the vertical-deformation time $t_d=S^{-1}$.

The Richardson number Ri is the classical index to assess the occurrence of instability, but it cannot differentiate between instances when mixing effectively redistributes density (and other properties, which occurs when stratification is high) from cases when instabilities simply displace water that is already homogeneous. Here, we have extended previous work [20,21] by defining the reduced squared vertical shear ${\sigma }^2$ in vertical coordinates and the reduced squared diapycnal shear ${{\sigma }_{\rho }}^2$ in isopycnic coordinates. The three indices for instability (Ri, ${\sigma }^2$ and ${{\sigma }_{\rho }}^2$) are sensitive to stratification, but ${{\sigma }_{\rho }}^2$ has the special property of decreasing with increasing stratification towards a limiting value equal to the squared diapycnal shear ${S_{\rho }}^2$.

Finally, we have examined the distribution of the nondimensional dependent variables ($N$, $S$, $S_{\rho }$) and indices (Ri, ${\sigma }^2$, ${{\sigma }_{\rho }}^2$) for three quite different dynamic regimes: the Gran Canaria Island shelf break characterized by intense internal waves, the strongly baroclinic Gulf Stream and the gravity-current Mediterranean outflow. The most remarkable feature is the decay of both the mean and maximum ${S_{\rho }}^2$ and ${{\sigma }_{\rho }}^2$ values with increasing stratification, which changes for each different flow regime, endorsing the idea that ${{\sigma }_{\rho }}^2>0$ is a good index for assessing flow stability: ${{\sigma }_{\rho }}^2$ has to be large for low stratification values, but small positive values are sufficient for high stratification. The non-dimensional expression ${{\sigma }_{\rho }}^2=\left(1-Ri\right)/\left(N^{2}Ri\right)$ fits well to the mean ${{\sigma }_{\rho }}^2$ experimental values, leading to a family of instability functions with a characteristic Ri value for each flow regime.