The instability of diffusive convection and its implication for the thermohaline staircases in the deep Arctic Ocean

In the present study, the classical description of diffusive convection is updated to interpret the instability of diffusive interfaces and the dynamical evolution of the bottom layer in the deep Arctic Ocean. In the new consideration of convective instability, both the background salinity stratification and rotation are involved. The critical Rayleigh number of diffusive convection is found to vary from 10 3 to 1011 in the deep Arctic Ocean as well as in other oceans and lakes. In such a wide range of conditions, the interfaceinduced thermal Rayleigh number is shown to be consistent with the critical Rayleigh number of diffusive convection. In most regions, background salinity stratification is found to be the main hindrance to the occurrence of convecting layers. With the new parameterization, it is predicted that the maximum thickness of the bottom layer is 1051 m in the deep Arctic Ocean, which is close to the observed value of 929 m. The evolution time of the bottom layer is predicted to be ∼ 100 yr, which is on the same order as that based on 14C isolation age estimation.


Introduction
Double diffusion is one of the most important nonmechanically driven mixing processes.It occurs in a fluid in which there are gradients of two (or more) properties with different molecular diffusivities and of opposing effects on the vertical density distribution.This phenomenon is of great interest to many disciplines in physical sciences and engineering (Turner, 1973), but the most active research area is exploration of thermocline staircases in oceans and lakes (Schmitt, 1994;Kelley et al., 2003).Basically, there exist two modes: salt finger (SF) and diffusive convection (DC).When cold and fresh water lies on top of warm and saline water, vertical mixing is triggered by different molecular diffusion rates of heat and salt, and DC forms (Kelley et al., 2003).DC is characterized by a series of thermohaline staircases: a stack of homogenous mixed layers of nearly constant temperature and salinity, separated by strongly stratified thin interfaces.
As sketched in Fig. 1, each convecting layer includes a well-mixed layer and is bounded by two adjacent interfaces.In each convecting layer of DC, the fluid properties and flow dynamics can be described by four dimensionless parameters.The first one is the thermal Rayleigh number Ra T = αgL 3 T /νκ T , where g is the gravitational acceleration, and α, ν, κ T , T , and L being, respectively, the thermal expansion coefficient, the kinematic viscosity, the thermal diffusivity, the temperature difference and the typical length scale of the convecting layer.The second one is the salinity Rayleigh number, Ra S = βgL 3 S/νκ T , where S is the salinity difference and β is the haline contraction coefficient.The other two parameters are the Prandtl number, P r = ν/κ T , and the Lewis number, Le = κ S /κ T , where κ S is the salinity diffusivity.As the interface is the boundary of two adjacent convecting layers, the thermal Rayleigh number of the interface, Ra TI , is proposed to be on the same order as the critical Rayleigh number of convection, which is on the order of 1000 (Turner, 1968(Turner, , 1973)).This argument has been found to work well in DC staircases in Lake Banyoles (Sánchez and Roget, 2007).
Recently, DC staircases have been observed in the deep Arctic Ocean (Timmermans et al., 2003;Timmermans and Garrett, 2006;Björk and Winsor, 2006;Carmack et al., 2012).These thermohaline staircases exhibit several unique Published by Copernicus Publications on behalf of the European Geosciences Union.characteristics.One of them is the thick diffusive interface.It is about 5-8 m, which is much larger than those observed in laboratory experiments (Turner, 1968;Huppert and Linden, 1979;Fernando, 1987), lakes (Sánchez and Roget, 2007;Schmid et al., 2010) and other ocean regions (Voorhis and Dorson, 1975;Larson and Gregg, 1983;Padman and Dillon, 1987;Anschutz and Blanc, 1996).Associated with the temperature difference across the interface, θ ∼ 1.3×10 −3 • C, and other fluid properties, as listed in Table 1 below, the thermal Rayleigh number of the interface, Ra TI , is on the order of 10 9 , which is much larger than the typical value reported in the literature (Turner, 1968(Turner, , 1973;;Sánchez and Roget, 2007).Therefore, additional understanding is needed to elucidate the difference between the result in the deep Arctic Ocean and those in previous studies.
Another prominent characteristic is that there exists a thick bottom mixed layer (Timmermans et al., 2003;Björk and Winsor, 2006;Timmermans and Garrett, 2006;Carmack et al., 2012).In the Canada Basin, this isothermal and isohaline bottom layer reaches to approximately 1000 m thick and extends about 1000 km across the basin (Timmermans et al., 2003).Björk and Winsor (2006) proposed a simple one-dimensional diffusive-convection model to understand the dynamics of the bottom layer.However, their proposed evolution process strongly depends on turbulent eddy diffusivity, which is unavailable at the present stage, and can only explain parts of the observation.Therefore, it is necessary to seek other mechanisms that may dominate in the deep-water evolution process.
In this paper, we explore the instability of the diffusive interface using hydrographic data measured by McLane Moored Profilers (MMP) at a fixed location in the Canada Basin.The classical description proposed by Turner (1968Turner ( , 1973) ) is extended here by considering the influences of background salinity stratification and rotation.Then the

Data
Hydrographic data were mainly obtained from the Beaufort Gyre Exploration Project (BGEP) (Ostrom et al., 2004;Proshutinsky et al., 2009).We focus on the variances of temperature and salinity in the deep ocean at a fixed location in the Beaufort Sea of the Canada Basin.As marked in Fig. 2, Mooring A was deployed at 75  (Timmermans et al., 2003;Björk and Winsor, 2006;Timmermans and Garrett, 2006).These structures are less pronounced in the salinity profiles, as shown in Fig. 3b, due to the instrument resolution.As analyzed in our previous work (Zhou and Lu, 2013), the interface properties could be determined with an averaging technique.In Fig. 3a, four DC steps can be identified from the potential temperature, θ, profile, which are referred to as the 1st, 2nd, 3rd and 4th steps from the bottom to the top, respectively.Each step includes a mixed layer and its overlying interface.Note that the mixed layer of the first step is the bottom layer.Usually, two parameters are used to characterize the susceptibility of water columns to DC staircases.One is the density ratio R ρ (Turner, 1965), which is expressed as R ρ is the ratio of the stabilizing force due to the salinity gradient and the destabilizing force due to the temperature gradient.The other one is the buoyancy frequency, which has the form N shows the stability of the stratification of water columns.
It has been proposed that the shape of DC structures and the vertical heat transfer strongly depend on these parameters (see the review paper of Kelley et al., 2003).In the case of well-developed DC steps, where the mixed layer is homogeneous in temperature and salinity, ∂S/∂z and ∂θ/∂z are reduced to the salinity and temperature gradients across the interface.According to the definitions of Ra T and Ra S , R ρ can be rewritten as R ρ = Ra S /Ra T .For detailed analysis of the determination of the diffusive interface properties, readers are referred to Zhou and Lu (2013).The mean properties of the four interfaces are listed in Table 1.

The classical description of DC instability
The instability of DC has been extensively studied based on the results in laboratory experiments (Veronis, 1965;Turner, 1973;Caldwell, 1974;Huppert and Moore, 1976; Pearlstein, 1981).Basically, two modes may occur at the onset of convection; one is steady convection and the other is oscillatory convection.In the case of the ocean (Le ∼ 0.01), it has been found that the onset Rayleigh number of steady convection appears to be less than the value for oscillatory convection (Huppert and Moore, 1976).This implies that the steady convection mode is responsible for convective instability in DC staircases.In terms of linear stability analysis (Turner, 1973), the convection occurs when the thermal Rayleigh number, Ra T , exceeds a critical value, Ra c .Ra c is written as In the case of homogeneous fluid, Ra S = 0 and Le = 0, and the equation is reduced to the onset Rayleigh number of Rayleigh-Bénard convection, Ra c = 27π 4 4 .In the deep Arctic Ocean, the DC staircases are thick and the corresponding Coriolis frequency, f , is large; one may expect that the rotation influences the occurrence of convection.The rotation is characterized by the Taylor number, T a = 4( f L 2 ν ) 2 .In linear stability analysis, it has been found that the rotation inhibits the onset of convective instability (Pearlstein, 1981).When the rotation strongly affects DC, the onset Rayleigh number, Ra c , would be in the form of Equation ( 4) implies that the convection occurs when the heat-induced buoyancy force gradient overcomes the resistance produced by the salt stratification and rotation.In the ocean, Le ∼ 0.01 and P r ∼ 4 − 13; Eqs. ( 3) and ( 4) are then simplified to for T a > 10 5 .

Results and discussions
There is not much work to address the instability of DC in the ocean because the involved convective flow state is typically far beyond the onset convection regime (Carpenter et al., 2012).However, findings about the instability mechanism of DC in previous laboratory experiments (Howard, 1964;Turner, 1968;Huppert and Linden, 1979) may actually help to understand the characteristics of diffusive interface and the bottom layer in the deep Arctic Ocean, both of which are discussed next.

Diffusive interface
For the diffusive interface in the deep Arctic Ocean, the corresponding thermal Rayleigh number, Ra TI , is obtained from its definition, with the typical scale L and temperature difference T being identified by the interface thickness h and temperature difference θ.As h and θ vary with time (Zhou and Lu, 2013), the deduced Ra TI is a function of time too.Its temporal distributions (shown in Fig. 4) indicate that the Ra TI of each interface is distributed approximately lognormal, which suggests that Ra TI is strongly intermittent.Similar distributions have been found in vertical heat flux, eddy diffusivity and other properties (Zhou and Lu, 2013).These results imply that the deep Arctic Ocean exhibits certain turbulent behaviors (Frisch, 1995).The mean Ra TI of each interface, as listed in Table 1, is on the order of ∼ 10 9 .
Although the diffusive interface is the internal boundary of two adjacent convecting layers, how to relate the interface to the boundary layer is arguable (Turner, 1973;Linden and Shirtcliffe, 1978;Padman and Dillon, 1989).Linden and Shirtcliffe (1978) proposed a boundary layer model.In this model, the boundary layer, δ, is defined as δ = ( h T − h S )/2, where h T and h S are the thicknesses of interface based on temperature and salinity.This model has been further employed in some studies, e.g., Worster (2004), andCarpenter et al. (2012).
On the other hand, many studies follow the classical boundary layer definition in diffusive convection (Turner, 1968(Turner, , 1973;;Padman and Dillon, 1989;Sánchez and Roget, 2007;Zaussinger and Spruit, 2013), which is similar to that in Rayleigh-Bénard convection.In Rayleigh-Bénard convection, two thermal boundary layers exist near the top and bottom ends.At the boundary, the temperature gradient is the largest.Away from the boundary, the temperature gradient becomes smaller and smaller till it reaches zero in the mixed layer.Such a temperature distribution is the result of convecting rolls.The thermal boundary layer thickness is defined as a small region between the boundary and the mixed layer (Lui and Xia, 1998).In DC, the temperature gradient is the largest around the middle of the interface, and it gets smaller and smaller as it approaches the edges of the interface.By analogy with Rayleigh-Bénard convection, each interface indeed consists of two boundary layers, as shown in Fig. 1, which are the top boundary layer of the lower convecting layer and the bottom one of the upper convecting layer.Here, we still use this definition of the boundary layer thickness.To the first-order approximation, the boundary layer thickness, δ, and the corresponding temperature difference, θ δ , can be taken as δ ∼ h/2 and θ δ θ/2.Consequently, the boundary layer thermal Rayleigh number, Ra Tδ , is Ra Tδ ∼ Ra TI /16.The same consideration is also applied to the salinity Rayleigh number, Ra Sδ , and the Taylor number, T a δ , of the boundary layer.Typical values of Ra Tδ and T a δ in the deep Arctic Ocean are listed in Table 2.
In Rayleigh-Bénard convection, it is assumed that the Rayleigh number based on the boundary layer thickness, Ra δ , is on the order of the critical value, namely, Ra Tδ ∼ Ra c (Howard, 1964;Siggia, 1994).This assumption has been confirmed to be marginally correct in experiments (e.g., Castaing et al., 1989;Lui and Xia, 1998).Analogously, the same argument has been employed to study the diffusive interface instability by Turner (1968Turner ( , 1973)).This argument has been found to work well in the DC staircases of the Arctic Ocean (Padman and Dillon, 1989) and Lake Banyoles (Sánchez and Roget, 2007).However, in the deep Arctic Ocean, Ra TI of interfaces is on the order of ∼ 10 9 ; then the deduced Ra Tδ is 10 8 , which is much larger than the observed value (∼ 10 3 ) Table 2. Data used in the analysis.The average values are listed when more than one data set are included in the data source.There are the temperature difference, θ , salinity difference, S, thickness, h, and Prandtl number, P r, of the interface; saline buoyancy frequency, N S , based on the background salinity profile; Taylor number, T a δ , thermal Rayleigh number, Ra Tδ , and salinity Rayleigh number, Ra Sδ , of the boundary layer.Deep Arctic -data of the present work; Upper Arctic -average data measured by the Ice-Tethered Profiler (ITP) 2 in the Beaufort Gyre observation system (BGOS); Atlantis II 1 -Fig.4 of Anschutz and Blanc (1996); Atlantis II 2 -Fig.4 of Swift et al. (2012); Lake Kivu -Fig. 1 and Table 2 of Schmid et al. (2010); Lake Banyoles -Figs. 2 and 8 of Sánchez and Roget (2007); and the Bahamas -Fig. 2 of Larson and Gregg (1983) in Rayleigh-Bénard convection (Castaing et al., 1989) and in DC (Sánchez and Roget, 2007).One may expect that the influences of salt stratification or rotation must be involved, that is, Eqs. ( 5) and ( 6) can be used to explain the results here.However, as listed in Table 1, to each interface, the salinity Rayleigh number Ra SI is always larger than Ra TI .Thus, when Ra Sδ (∼ Ra SI /16) is employed in Eqs. ( 5) and ( 6), the obtained Ra c would be larger than Ra Tδ (∼ Ra TI /16), which may mean that the assumption of Ra Tδ ∼ Ra c does not work here.Carpenter et al. (2012) suggested that Eqs. ( 5) and ( 6) cannot apply to the instability of diffusive convection of large density ratio R ρ (R ρ > 1.15).
An alternative consideration can be found in a previous laboratory study by Turner (1968).In this study, Eq. ( 5) was used to interpret the onset of the convecting layer in DC.When the breakdown of the unstable boundary layer is described, it is further argued that the original salinity gradient remains unchanged because the salt diffuses much more slowly than heat.In other words, the salinity Rayleigh number, Ra S , is in the form of In Eq. ( 7), N S = gβ dS dz is the buoyancy frequency based on the initial salinity gradient before being heated.This argument has been supported by the measurement of the advancing boundary layer in another laboratory experiment (Huppert and Linden, 1979).In the laboratory, however, the calculated Ra S is found to be very small due to the thin boundary layer, and it is reasonable to neglect it.Consequently, Ra c is reduced to the ordinary critical Rayleigh number based on heat alone (Turner, 1968).
It is expected that the same idea can be applied to the case of interest here.As it is impossible to obtain the initial salinity gradient of the deep Arctic Ocean, the background salinity gradient is alternatively used in Eq. ( 7).As shown in Fig. 3b, the background salinity gradient is obtained by the linear fitting to all DC staircases in salinity profiles, resulting in ds dz = 1.3×10 −5 m −1 and the corresponding N S = 3.1 × 10 −4 s −1 .Meanwhile, more data in the oceans and lakes have been collected to check the applicability of the argument.These data include the field observations in the upper layer of the Canada Basin (BGOS, 2013), deep Red Sea (Anschutz and Blanc, 1996;Swift et al., 2012), Lake Kivu (Schmid et al., 2010), Lake Banyoles (Sánchez and Roget, 2007), and the Bahamas (Larson and Gregg, 1983).Typical interface properties of these data sources are listed in Table 2.
In terms of Eqs. ( 6) and ( 7), the thermal critical Rayleigh number, Ra c , is obtained here from the δ-based salinity Rayleigh number, Ra Sδ , and the δ-based Taylor number, T a δ .Note that Eq. ( 5) is used when T a < 1000.Within the collected data, the thermal critical Rayleigh number, Ra c , is found to vary in a wide range from 10 3 in the Lake Banyoles to 10 11 in the deep Red Sea, as shown in Fig. 5a. Figure 5a also shows that the calculated thermal critical Rayleigh number, Ra c , is very close to the δ-based thermal Rayleigh number, Ra Tδ , in almost all the collected data.As some data were captured from the published figures (Sánchez and Roget, 2007;Larson and Gregg, 1983), the limited accuracy would be the most probable reason responsible for the data scatters in Fig. 5a.When the set of data is large enough to have good statistics, e.g., those in the Upper Arctic (BGOS, 2013), Ra c is in good agreement with Ra Tδ .As far as we know, this is the first time that the classical description of DC has been applied to interpret instability of diffusive interfaces in the oceans.
In addition, we examine the influences of background salinity stratification and the rotation on the onset of convection.The comparison between the first and the second terms in Eq. ( 6) is plotted in Fig. 5b.In the deep Arctic Ocean, the influence of rotation is supposed to be the largest because of the thickest diffusive interface and the largest Coriolis www.ocean-sci.net/10/127/2014/Ocean Sci., 10, 127-134, 2014 frequency in all the data.Even so, the contribution of rotation is found to be about 1000 times smaller than that of the salinity stratification.In other regions, the contribution of rotation is even smaller.Therefore, background salinity stratification would be the main hindrance to the occurrence of convective flow in the ocean.The influence of the salinity stratification, however, is not always important.In Lake Banyoles, as shown in Fig. 5b, it is found that Ra Sδ ∼ 55; thus, the inherent mechanism of thermal convection plays a dominating role there, which results in Ra c ∼ 1000 (Sánchez and Roget, 2007).

Bottom homogenous layer
As introduced in Sect. 1, the thick homogenous bottom layer is a unique feature of the deep Arctic Ocean (Timmermans et al., 2003;Björk and Winsor, 2006;Timmermans and Garrett, 2006;Carmack et al., 2012;Zhou and Lu, 2013).
As the model proposed by Björk and Winsor (2006) cannot explain the dynamics of the bottom layer well, we attempt to find other mechanisms instead.Recently, more evidence has been collected that geothermal heating plays an important role in the hydrographic configuration of deep water in the Arctic Ocean (Timmermans et al., 2003;Björk and Winsor, 2006;Timmermans and Garrett, 2006;Carmack et al., 2012;Zhou and Lu, 2013).The situation that the weakly stratified deep Arctic Ocean is being heated by geothermal heating is similar to those studied in the laboratory experiments and numerical simulations, where DC occurred in salt-stratified fluid heated from below (Turner, 1968;Huppert and Linden, 1979;Fernando, 1987).Thus, it is interesting to examine whether the results obtained from the laboratory experiments can be applied to the deep ocean and whether they can shed more light on the deep water dynamic process.
In the laboratory experiments where a salt-stratified fluid was heated from below (Turner, 1968(Turner, , 1973;;Huppert and Linden, 1979;Fernando, 1987), it has been found that the bottom layer is much thicker than the overlying staircases, which is similar to the observation in the deep Arctic Ocean.In Turner's (1968) theoretical work, he suggested that the homogeneous bottom layer reaches a critical thickness before a second convecting layer appears on top of it.This maximum thickness was proposed to be in the form of Here q 0 is the buoyancy flux, which is derived as q 0 = αgF /ρ 0 c p , and F is the heat flux supplied from the bottom boundary.Eq. ( 8) has been verified in Turner's laboratory, except that the fitted Ra c was about ∼ 2.4 × 10 4 , which is larger than the ordinary value (∼ 1000).As discussed in the last subsection, the instability of DC follows the arguments based on Eqs. ( 5), ( 6) and ( 7).Then we can evaluate the bottom layer thickness in the Arctic Ocean by using Eq. ( 8).
In the deep Arctic Ocean, geothermal heating is taken as F ∼ 50 m Wm −2 (Langseth et al., 1990).At Mooring A, the physical properties of sea water are respectively α With the related properties of the first interface, as listed in Table 1, h c1 is calculated to be 526 m, which is smaller than the observed value of 929 m (Zhou and Lu, 2013).However, as found in the laboratory experiment (Turner, 1968), the fitted Ra c was about 20 times the ordinary value, and so was close to the interface thermal Rayleigh number Ra TI .If Ra c is replaced by Ra TI in Eq. ( 8), h c1 is found to be 1051 m, which is much closer to the observed value.
In another laboratory study, Fernando (1987) found that the bottom mixed layer still grows upwards even after the second mixed layer forms.Based on the assumption that the kinetic and potential energies of turbulent eddies are Ocean Sci., 10, 127-134, 2014 www.ocean-sci.net/10/127/2014/balanced when the mixed layer grows to a critical height, an alternative expression for the layer thickness has been proposed, With the properties of the interface, the thickness h c2 is found to be 29 m, which is much smaller than the observation.This suggests that the assumption of energy balance cannot explain the formation of the bottom layer.
In a study of 14 C isolation age, the deep water of the Canada Basin is estimated to be about 500 yr old (Macdonald et al., 1993).According to the laboratory study of Turner (1968), the evolution time of the bottom layer can be predicted by τ = (N S h) 2 2q 0 (10) before the formation of the secondary convecting layers.Supposing that the bottom layer thickness, h, is about 1000 m, the time τ is calculated to be τ = 109 yr.In observations, however, there are more than four thermohaline steps overlying the bottom layer, which means that the evolution time of the bottom layer must be longer than τ .In some sense, τ can be taken as the lower bound of evolution time of the bottom layer.Based on the above results, it is suggested that the main evolution mechanism of the deep water in the deep Arctic Ocean is similar to that in the laboratory experiments (Turner, 1968(Turner, , 1973;;Huppert and Linden, 1979;Fernando, 1987).The homogenous bottom layer is the result of geothermal heating under salinity stratification.Its maximum thickness can be described by Eq. ( 8), where Ra TI is alternatively used.The time, τ , is inferred to be on the same order as that estimated from the 14 C isolation age detection.Both results imply that the evolution of the bottom layer is mainly controlled by the convective instability mechanism.

Conclusions
In summary, the classical description of DC has been updated to interpret the instability of diffusive interfaces and the dynamical characteristics of the bottom layer in the deep Arctic Ocean.When both background salinity stratification and rotation are considered, the critical Rayleigh number, Ra c , of DC is found to vary in a wide range from 10 3 to 10 11 within the collected data, and the thermal Rayleigh number of the boundary layer, Ra Tδ , is consistent with the critical Rayleigh number, Ra c , of DC.It is expected that this new parameterization of Eqs. ( 5), ( 6) and ( 7) can be extensively applied to DC in other oceans and lakes.In most cases, the salinity stratification is found to be the main hindrance to the occurrence of convective flow, except for regions where the diffusive interface is thin, e.g., in Lake Banyoles.When the interface thermal Rayleigh number, Ra TI , is alternatively employed in the old parameterization of Eq. ( 8), the predicted maximum thickness of the bottom layer is 1051 m, which is close to the observed value of 929 m.The evolution time of the bottom layer, τ , is predicted to be ∼ 100 yr, which is on the same order as the 14 C isolation age estimate.As several layers overlie the bottom layer, the evolution time of water in the deep Arctic Ocean must be longer than τ , and τ can be taken as the lower bound of the residence time.
In the formation of deep Arctic water, one perhaps cannot exclude the effects of other instabilities, e.g., thermobaric convection (Carmack et al., 2012) or topographic Rossby waves formed as a result of instability of a strong current (Timmermans et al., 2010).Nonetheless, according to the results in this study, the main evolution mechanism of the deep water would be similar to that in the laboratory experiments where salt-stratified fluid was heated from below (Turner, 1968(Turner, , 1973;;Huppert and Linden, 1979;Fernando, 1987).Therefore, when the classical description of DC is updated, it can be applied to interpret the flow dynamics of the deep Arctic Ocean.
Fig. 1.(a) A sketch of temperature distribution of a diffusive interface in DC staircases; δ is the boundary layer thickness of a convecting layer.(b) A cartoon of flow pattern near the interface.

Fig. 2 .
Fig. 2. Map of the Canada Basin in the Arctic Ocean.Isobaths are plotted using the ETOP data.The location of Mooring A is marked by a star.

Fig. 3 .
Fig. 3. Profiles of (a) potential temperature (θ ) and (b) salinity (S) measured by conductivity-temperature-depth (CTD) at Mooring A on 28 Septermber 2010.The dashed line is the linear fit of the DC staircases in the salinity profile.

Fig. 4 .
Fig. 4. Probability density functions (PDFs) of the thermal Rayleigh number, Ra TI , of the four diffusive interfaces.Data of the first, second, third and fourth interfaces are plotted using the dashed lines with solid circles, open circles, solid down triangles and open triangles, respectively.The log-normal fitting of the PDF of the first interface is shown as a gray curve.

Table 1 .
Mean properties of four diffusive interfaces at Mooring A between 2 October 2009 and 9 August 2010.They are potential temperature difference, θ , salinity difference, S, thickness, h, density ratio, R ρ , buoyancy frequency, N , thermal Rayleigh number, Ra TI , salinity Rayleigh number, Ra SI , and Taylor number, T a I .