Fixed-flux salt-finger convection in the small diffusivity ratio limit

Doubly diffusive systems, where two components with different diffusivities control buoyancy, are responsible for mixing processes that are important in geophysics and astrophysics.^8,26,35,37 For example, in the oceanic setting, heat and salinity both contribute to density stratification, while in stellar astrophysics, chemical concentration replaces salinity. We are interested in statically stable regimes in which one component plays a stabilizing role, while the other is destabilizing. Such configurations may become diffusively unstable once diffusive effects are incorporated.

In this paper, we focus on salt-finger convection, where the contribution to the stratification from the fast diffuser (temperature) is stabilizing, while that from the slow diffuser (salinity) is destabilizing. Stern²⁹ was the first to study the linear instability of the resulting stratification and established the existence of the salt-finger instability in which denser salt-rich fingers descend with less dense fresher fluid rising in between. However, the nonlinear evolution of this instability and, in particular, its nonlinear saturation is still imperfectly understood. Stern³⁰ suggested a collective instability mechanism where parallel fingers generate large scale gravity waves that disrupt the salt-finger field. This mechanism is a special case of a secondary instability identified by Holyer.¹¹ Based on this formulation, a saturation criterion, namely, that the growth rate of the secondary instability be comparable to that of the primary salt-finger instability, was proposed.^1,22 Other saturation mechanisms include finger collision,²⁷ finger clustering,¹⁸ and finger-mean shear interaction.^9,40,42

Reduced models of geophysical and astrophysical flows (cf. Ref. 10) are systematic reductions of the primitive equations that are valid in appropriate asymptotic parameter regimes. Such models are not only able to reach the extreme parameter regimes of interest in such flows but also do so at a reduced numerical cost. In doubly diffusive systems characterized by a small diffusivity ratio τ ≪ 1, the ratio of diffusivities of the competing salinity and temperature fields, such models are obtained by systematically exploring the limit τ → 0.^19,20,24,40 In this paper, we focus on the inertia-free salt convection model derived by Xie et al.⁴⁰ in the limit τ → 0 and large density ratio. Although the model filters out gravity waves, it captures both linear and secondary instabilities, as well as the strongly nonlinear dynamics that lead to the saturation of the salt-finger instability.

In the ocean, the characteristic spatial scale of salt-finger convection is in the range of several centimeters, and its characteristic time scale is the time scale for salinity diffusion, ∼10 days. These small spatial and slow temporal scales are very hard to resolve in ocean models, and salt-finger convection must therefore be parameterized. The parameterization of the impact of such small-scale structures on large-scale motions requires a relation between fluxes and mean gradients of temperature and salinity. To obtain such relations, most authors fix the mean gradients and study the resulting mean fluxes in a periodic domain, e.g., Refs. 20 and 34, or fix the temperature and salinity difference across a vertically bounded domain.^18,27,38,45

The approach of the present paper is different. We consider the conjugate situation in which the total (large-scale) fluxes are fixed. These should be understood as the sum of the mean fluxes generated by fluctuations and the impact of mean gradients. In a fixed-flux setup, the mean gradients are permitted to evolve according to the nonlinear dynamics of the fluctuations, and we focus on their final saturated mean values. From an energy perspective, the fixed-flux condition bounds the total energy of the system and therefore avoids blowup. We mention that fixed-flux conditions are relevant in many other situations, including wall-bounded shear flow³⁹ and the saturation of the magnetorotational instability.⁴¹ 

We focus on periodic domains, as appropriate for oceanic and astrophysical applications with distant boundaries.^21,28,33,34 Such domains allow for the existence of exact elevator-mode solutions, while ignoring complications arising from the presence of boundary layers. Thus, our work differs from bounded domain studies in which such modes are absent.^{12,18,23,27,38,43–45} Moreover, we are not only concerned with ocean-relevant large-domain properties but also study the dynamics in confined small domains, a setup applicable to rock fracture and crystallization in cooling magma chambers⁵ and other industrial and laboratory processes.³⁶ Our main conclusions can be illustrated using a two-dimensional formulation of the problem, and this is therefore the approach adopted here; the corresponding three-dimensional results will be reported elsewhere.

The structure of this paper is as follows. We formulate the fixed-flux system in Sec. II. In Secs. III and IV, we study the single-mode theory and a two-mode truncation both theoretically and numerically. This mode-truncation model differs from the other mode-truncation models in doubly diffusive convection^{6,7,13–15,17,25,32,38} not only in the reduced model (2) from which our truncated model is derived but also in the presence of elevator modes whose growth is controlled by the fixed-flux condition (5a). We study the statistically steady large-domain dynamics in Sec. V. Section VI summarizes our results.

Motivated by oceanographic applications with large Schmidt numbers, i.e., large ratio of viscous to salt dissipation, Xie et al.⁴⁰ derived the following inertia-free salt convection (IFSC) model:

in order to describe two-dimensional nonlinear salt-finger convection in the limit of a small solute to heat diffusivity ratio τ ≪ 1 and a large density ratio R_ρ ≫ 1. Here, the density ratio R_ρ is the ratio between the contributions to the density from the background temperature and salinity, and R_ρ > 1 implies a stable stratification; S is a scaled salinity perturbation on top of a linear salinity background, while ψ is the stream function. The salinity in the background state increases upwards and is stabilized by a corresponding linear thermal stratification; the assumption τ ≪ 1 allows one to replace thermal fluctuations by corresponding fluctuations in ψ leading to a closed system for S and ψ at leading order in τ (for details, see Ref. 40).

All the quantities are nondimensionalized by the salinity diffusive temporal and spatial scales. The only external parameter in this model is the Rayleigh ratio Ra, defined as the ratio between (destabilizing) salinity and (stabilizing) temperature Rayleigh numbers, Ra ≡ Ra_S/Ra_T = 1/(R_ρτ), and is required to be O(1). Here, the Rayleigh numbers are defined in the usual fashion, e.g., Ra_T = gα_Tβ_TH⁴/(νκ_T), where g is the gravitational acceleration, α_T is the thermal expansion coefficient, β_T is the background temperature gradient, H is the domain height, ν is the kinematic viscosity, and κ_T is the thermal diffusivity. In the model, the temperature is slaved to the stream function, and the mean temperature gradient is unchanged by the dynamics. This reduced model is perhaps the simplest model that captures the key features of salt-finger convection including both linear and secondary instabilities and the various stages leading toward a statistically saturated state.

Xie et al.⁴⁰ used this model to study the dependence of the salinity flux on the Rayleigh ratio in the presence of a fixed mean salinity gradient, i.e., at a fixed Rayleigh ratio. We consider here the possibility that in the ocean, the mean gradients are instead controlled by large-scale mean fluxes and propose a conjugate scenario where the mean salinity flux is prescribed and the Rayleigh ratio is allowed to evolve, leading to the modified model,

Here, F is the prescribed salinity flux and $⋅¯$ and $⋅$ denote, respectively, the horizontal and vertical averages. The model (2) has the same prognostic equation (2a) and diagnostic equation (2b) as the IFSC model, but the mean gradient, the Rayleigh ratio Ra, is no longer fixed, and it evolves subject to the prescribed constant total flux F through relation (2c). In such a model, we expect to find a fixed point, i.e., a Rayleigh ratio that generates the right salinity flux, but this state is no longer guaranteed to be unique.

The definition (2c) of Ra stems from the horizontally averaged Eq. (1a),

For a (statistically) steady state, $S¯t=0$⁠, and therefore,

Thus, we can consistently define the salinity flux F and the Rayleigh ratio Ra as

leading to Eq. (2c). Since F is fixed, Ra must now evolve in response to changes in the fluctuation-induced transport term,

Thus, the fixed-flux condition captures the feedback from the fluctuation–fluctuation interaction on Ra, which is missing in fixed-gradient setups.

Multiplying the salinity Eq. (2a) by S and taking temporal and spatial averages, we obtain an equation for the energy balance. In a steady or a statistically steady state, this equation reads

where ϵ > 0 is the dissipation rate. In this energy relation, we can see that the fixed-flux and fixed-gradient conditions are conjugate to one another. However, the one-to-one correspondence between Ra and F no longer holds. Moreover, for a fixed energy dissipation rate ϵ, there are two possible values of Ra that correspond to a given F. Finding the relation between F and Ra is one of the primary aims of this paper.

From the energy equation (7), the positiveness of the dissipation rate implies that the mean salinity gradient is positive, indicating that salinity is indeed destabilizing, and so generates a negative perturbation flux, reflecting the conversion of the saline available potential energy into kinetic energy. Quantitatively, 0 < Ra < −F, and so

According to the Poincaré inequality, the saline available potential energy

is therefore bounded,

where Λ ∼ L⁻² is the smallest eigenvalue of the Laplacian on the domain. Thus, the fixed-flux condition also avoids blowup. Here, equality is achieved by an exact elevator-mode solution whose wavelength equals the domain width. We study such exact solutions in Sec. III.

Linearizing Eqs. (2) with respect to the zero static state, we obtain the linear growth rate of the fixed-flux system,

where k and m are the horizontal and vertical wavenumbers, respectively, and $K=k2+m2$⁠. Thus, linear instability requires F < −1, a fact consistent with the instability condition Ra > 1 in a fixed-gradient system.⁴⁰ This correspondence is mathematically trivial because the linearization of (2c) is Ra = −F.

We can also write the equation for the saline available potential energy in spectral space,

where

$⋅^$ denotes a Fourier transform, and k denotes the wavenumber. Evidently, there are two mechanisms that are potentially responsible for reaching a (statistically) steady state. One is the modification of the mean gradient Ra to a marginally stable state through relation (2c), while the other is the transfer of energy from linearly growing modes to linearly damped modes through the nonlinear term in Eq. (12).

In this section, we present an analytical saturated single-mode solution subject to the fixed-flux requirement. We emphasize that the saturation mechanism is strongly nonlinear owing to the finite-amplitude modification of Ra through the perturbation flux [cf. Eq. (2c)]. For this purpose, we focus on elevator-type solutions since these single-mode solutions are known to play a potentially important role in the statistically stationary state even when the flow is turbulent. For example, elevator modes are responsible for abrupt but episodic energy growth in Rayleigh–Bénard convection in a three-dimensional periodic domain.⁴ 

We assume an inclined single-mode solution in the form

From the diagnostic relation (2b), we see that this single-mode solution generates a negative-definite salinity flux,

The fixed-flux condition (2c) thus indicates the presence of a modified mean salinity gradient,

Substituting (16) into the prognostic equation (2a), we obtain

where

Note that nonlinear advection vanishes identically for single-mode solutions and that the fixed-flux condition is responsible for the nonlinear term in Eq. (17).

We need a positive linear growth rate λ to excite nonzero solutions of (17), i.e., F < −1. Since D is positive definite, λ > 0 indicates a stable saturate state with

On substituting (19) into (16), we find that this saturated state corresponds to the saturated Rayleigh ratio,

i.e., to the marginal instability of the fixed-gradient state. It is interesting to note that in this case, Ra is flux-independent but scale-dependent. This is not true in the statistically steady state discussed farther below.

Elevator modes are exact solutions of (2), but owing to the nonlinear term, it is possible that an elevator mode will not be selected in the presence of perturbations. Hence, in this section, we consider a mode-truncated model that combines one elevator mode, one inclined mode, and one mean component,

where

This choice is motivated by the observation of inclined fingers in two-dimensional numerical simulations^9,42 and in ocean observations.¹⁶ The diagnostic equation (2b) then results in the stream functions

and the amplitudes of these modes evolve according to the truncated equations,

where

The system (24) has four steady solutions,

and

hereafter referred to as I–IV, respectively.

As exact solutions, the solutions I (26) and II (27) have the form of single-mode theory (19), and their existence requires that the following conditions hold:

and

respectively; these solutions may be stable. In contrast, solution III is always unstable (see the  Appendix) and so is the conduction state, solution IV. Note that solution IV is a saddle point because the mean mode C represents a stable eigendirection.

Within the mode-truncation model (24), the saturated state I (26) is always linearly stable. On the other hand, solution II (27) may be stable or unstable, depending on parameters, and its linear growth rate λ_II is determined by

Since $K2+βF+32β2BII2>0$ whenever solution II exists, the condition for a stable solution II is that both solutions of equation

have negative real parts.

In Fig. 1, we show the stability of the solutions in the (k, m) plane for two different values of F: F = −1.05 and F = −2.7. Both solutions are supercritical. Solution I is stable to the left of the vertical line, corresponding to the condition for the existence of solution I, while solution II can either be stable or unstable. This solution exists below the curve: when −F is small, e.g., panel (a) with F = −1.05, solution II is stable everywhere in its region of existence (region II_s), while for larger −F, an unstable region (region II_u) emerges, e.g., panel (b) with F = −2.7.

We check the stability condition in the case F = −2.7 and m = 0.4 by simulating Eqs. (24) numerically. The results are shown in Fig. 2, where we separate the two possible states using the perturbative saline available potential energy,

Since the saturated states are exact solutions of the original system (2), we can observe them in direct numerical simulations of Eqs. (2). Figure 3 shows a contour plot of the saturated salinity field S when F = −2.7 computed in a domain of size l × l, where l = 2π/0.4. This parameter regime corresponds to that used in Fig. 2 and permits stable solutions of both types I and II. Figure 4 shows an unstable solution II in a larger domain with l = 2π/0.1. This failure of the stability analysis in Fig. 1(b) is due to the presence of a sheared inclined Rayleigh–Taylor instability that takes place at interfaces where the salinity gradient is positive (Fig. 4).

We need to mention that solutions I and II only exist in the fixed Ra setup if the wavenumber is carefully selected as the marginal wavenumber. In the fixed-flux setup, the evolution of these solutions is less constrained and the failure to capture important mechanisms, such as the sheared inclined Rayleigh–Taylor instability, limits the validity of the truncated model. In addition, the three-mode model does not capture other important dynamics such as the collective instability³⁰ or finger collisions²⁷ that can also lead to saturation. Consequently, in Sec. V, we turn to a study of the partial differential equation system (2) in a large domain.

In the fixed-gradient setup, Xie et al.⁴⁰ found that the saturated state arising from an initially small amplitude random field within the reduced model (1) is maintained by the combined effect of secondary instabilities and finger collision. In addition, they found that the most unstable elevator mode (or optimal mode) determines the characteristic length scale of the saturated state. In our fixed-flux setup, exact single-mode solutions are easily reached in small domains, but in large domains, these single-mode solutions are subject to the same secondary instabilities as in the fixed-gradient setup, and we therefore expect that the above-mentioned saturation mechanism persists and hence that the relation between the salinity flux and Rayleigh ratio remains the same.

To obtain information on the statistically steady state, we numerically solve Eq. (2). We use a Fourier pseudospectral method with 2/3 dealiasing in space with resolution 512 × 512 in a square domain with a side length 32 times the optimal wavelength and a fourth-order explicit Runge–Kutta scheme in time in which the nonlinear terms are treated explicitly and linear terms implicitly using an integrating factor method. All calculations are initialized with small amplitude random Gaussian fields.

Figure 5 tests the above conjecture by plotting the dependence of the saline potential energy E_S and the perturbation flux F_S on the supercriticality Ra − 1. We find that this dependence continues to exhibit two regimes just as in the fixed-gradient setup and that these can still be interpreted within the theory developed by Xie et al.,⁴⁰ indicating the equivalence of the statistically steady states in the two conjugate systems. The expressions (black curves in Fig. 5) for small Ra are

while for large Ra,

where C is a constant and the optimal wavenumber is given by

The fixed-flux condition provides us with a new perspective for distinguishing between these two statistically steady regimes. For fixed flux F and dissipation ϵ, the energy balance, Eq. (7), leads to two possible values of Ra,

one of which must be selected by the statistically steady dynamics. Figure 6 compares the simulation results with the two solution branches in (38) as a function of the prescribed flux F. We see that regimes I and II defined by Xie et al.⁴⁰ for the fixed-gradient setup correspond to Ra₊ and Ra₋, respectively. We conclude that for a fixed flux and a fixed dissipation rate, the saturated statistically steady state of the system (2) tends to minimize Ra. In other words, since F is fixed, a statistically steady state tends to maximize the perturbation flux F_S, as supposed by Stern³¹ in his determination of the saturated state of salt-finger convection. Thus, even though the solution Ra₋ always provides a smaller Ra compared with Ra₊, it is not realized in regime I owing to the constraint Ra > 1. In addition, it is easy to check that the single-mode theory follows branch Ra₋.

We showed above that the relation between the mean salinity flux and Rayleigh ratio remains the same with fixed-gradient and fixed-flux conditions, suggesting that the feedback from perturbations on Ra through the fixed-flux condition is unimportant for the statistically steady states. However, despite this, the two cases are quite different, as revealed in different statistical measures. To exhibit these distinctions, we need to identify two statistically steady states with the same Ra and salinity flux F, satisfying fixed-gradient and fixed-flux conditions, respectively. We perform simulations in both cases, with the same domain size, time step and resolution, and the same initial conditions: a small amplitude random salinity field and zero velocity. With the fixed-flux condition F = −1000, we find a saturated Rayleigh ratio Ra = 3. This is shown in the left panel of Fig. 7. Next, we ran a numerical simulation with the fixed-gradient condition Ra = 3 until it reached a statistically steady state. In the middle panel of Fig. 7, we confirm that the saturated salinity flux matches the prescribed value F = −1000 in the fixed-flux simulation. We thereby confirm that these averaged quantities are the same in the two conjugate cases. It is therefore reasonable to compare different statistical quantities from these two simulations.

The third panel of Fig. 7 shows that the evolution of the salinity potential energy E_S in the two setups is quite different, even though they share the same mean value in the statistically steady state. There are in fact two fundamental differences in the evolution of E_S in the two cases. First, the time taken from the initial condition to the statistically steady state in the fixed-flux simulation is much shorter than that in the fixed-gradient simulation. The reason is that in the fixed-flux setup, the initial salinity field is weak and so the initial salinity flux is small. As a result, the fixed-flux condition (2c) implies a very large initial Rayleigh ratio [Ra = O(10³)], which leads to linear growth rates much larger than those with Ra = 3.

The second difference appears in the statistics of the statistically steady state. In Fig. 8, we show the pdf and the frequency spectra of E_S in the two cases. The left panel shows that even though E_S has the same mean value in both cases, the fluctuations in the fixed-gradient case are much larger than those in the fixed-flux case, indicating that the mechanism maintaining the mean value is much stronger in the presence of a fixed flux. This difference stems from the different mechanisms present in the two setups: in the fixed-gradient setup, E_S increases due to the growing modes, and to stop this growth, nonlinear fluctuation–fluctuation interactions are required to transfer energy from the growing modes to damped modes, and this nonlinear effect is slow. In contrast, in the fixed-flux setup, the growth of E_S leads to an increase in salinity flux, requiring a reduction of Ra [cf. Eq. (2c)]. This modification of the background gradient occurs through a fluctuation–mean interaction, which abruptly damps some of the growing modes. The different restoring mechanisms, i.e., the difference between the fluctuation–fluctuation interaction and the fluctuation–mean interaction, also explain why the skewness of E_S is negative in the fixed-gradient setup and is much larger than that in the fixed-flux setup. The stronger restoring mechanism in the fixed-flux setup also implies a shorter characteristic restoring time scale, as confirmed in the right panel of Fig. 8 showing that the peak of the frequency spectrum in the fixed-flux setup is at a higher frequency compared with that in the fixed-gradient setup.

In Fig. 9, we show snapshots of the salinity field S at a maximum and a minimum of E_S(t) (see Fig. 7, third panel) in the fixed-gradient and fixed-flux cases. We see immediately that at maximum E_S in the fixed-gradient case, the finger structure (both descending salt fingers and rising fresh fingers) is much more prominent than at the corresponding peak state in the fixed-flux case, with substantial power on large scales. The situation is reversed at E_S minima, where the finger structure in the fixed-flux case is more pronounced than in the fixed-gradient case.

To provide a quantitative comparison of the snapshots in Fig. 9, we present in the left panel of Fig. 10 the normalized salinity correlation function in the horizontal direction. This correlation function is an average over several realizations of the flow as well as over both x and z in order to reduce noise. The panel shows that the horizontal salinity correlation oscillates with a period close to the optimal wavelength l_opt, indicating the importance of linear dynamics in the saturated statistically steady state in both cases, and, in particular, of the characteristic length scale of the finger state. Moreover, in the fixed-gradient setup, the salinity field at peak E_S has a larger positive correlation than in the fixed-flux case, particularly at large scales, indicating larger patchiness, reflecting the visual impression from Fig. 9. We surmise that this is a consequence of the weaker restoring force present in the fixed-gradient case. The opposite is typically the case at E_S minima. As a result, the correlations at the peak and trough in the fixed-flux case lie in between those in the fixed-gradient case. These tendencies are confirmed by the corresponding Fourier spectra (Fig. 10, right panel), all of which show structures up to scales of order 8l_opt and then increase in a power–law fashion as k decreases toward zero, suggesting that a characteristic patch size is something like 8l_opt with very substantial fluctuations on yet larger scales. All four curves have a similar bimodal shape: the spectra peak around k = 0 and k = k_opt ≈ 0.76 for Ra = 3. Dividing the spectra at the maxima by those at the minima reveals a difference between the two formulations (inset in the right panel of Fig. 10): the ratio shows a peak around the optimal wavenumber in the fixed-gradient formulation, indicating the prominence of finger-like structures at E_S maxima, while in the fixed-flux formulation, this peak is largely absent, showing that the salinity field at a maximum is not that different from that at the minimum. We have also checked that the spectra of positive anomaly regions (descending salt fingers with S > 0) are indistinguishable from the spectra of negative anomaly regions (rising fresh fingers with S < 0) and conclude that no spontaneous symmetry breaking between rising and sinking fingers, either in the finger strength or their scale, takes place. This is so both at E_S maxima and at E_S minima.

In this paper, we studied salt-finger convection in the limit of small diffusivity and large density ratios in a fixed-flux setup. The model studied here is a modification of the inertia-free salt convection model derived by Xie et al.⁴⁰ for a fixed salinity gradient. Our setup differs from the classical setup with a fixed-gradient condition (fixed Rayleigh ratio Ra in the reduced model) by allowing the evolution of the Rayleigh ratio due to feedback from the perturbation salinity flux. We studied the single-mode dynamics and the interaction between two modes in a mode-truncated model both theoretically and numerically. The existence and stability of possible saturated states were determined, and the solutions were confirmed by numerical simulations of the original system (2) in a small domain. The fixed-flux setup avoids blowup of the single-mode solutions in a broad parameter regime, a property that distinguishes this setup from the more common fixed-gradient setup.

In large domains, these single-mode states are no longer stable and the system reaches a statistically steady regime. We found, based on relations between domain-averaged quantities, that the fixed-flux and the fixed-gradient conditions are identical in the saturated state. This is a nontrivial finding since in the two cases, the potential saturated states are different: for example, the elevator modes, which are exact solutions, can blow up in a fixed Ra setup but are bounded when the flux is fixed. As in the study of Ref. 40, two statistically steady regimes are identified as a function of increasing Ra and shown to coincide with the corresponding regimes in the fixed-gradient setup despite the imposition of a fixed flux and dissipation rate. Moreover, we saw that the statistically steady state selects the minimum Rayleigh ratio Ra or, equivalently, maximizes the perturbation flux. In the minimization process, we needed to impose the constraint that Ra > 1, i.e., the condition that there exist modes that supply energy to the statistically steady state. The branch selection is qualitatively consistent with the analysis of single-mode saturation, which also follows the smaller Ra branch.

The fact that the final states in the fixed-gradient and fixed-flux formulations of the problem yield the same saturated state is not surprising. After all, both states obey the same balance relations. However, the fact that there are two possible states for the given values of the flux −F and dissipation ϵ is not a priori obvious. It would be of great interest to characterize the competing unstable states, but this cannot be done using direct numerical simulation. In fact, the two stationary states that can be reached by numerical simulation in the fixed-gradient and fixed-flux cases have substantially different statistical properties. As shown in this paper, this is a consequence of the fact that the fixed-flux formulation includes the fluctuation–mean interaction, which instantaneously modifies the background Rayleigh ratio Ra, thereby changing the linear stability properties of the background state. For this reason, the fixed-flux formulation leads to a stronger restoring mechanism that is reflected in the resulting statistically steady state. In contrast, the restoring mechanism in the fixed-gradient formulation is via fluctuation–fluctuation interactions that transfer energy from linearly growing modes to damped modes and is therefore weaker. This difference in the restoring mechanism further leads to distinctive signatures in the statistics of the time series of the potential energy E_S: the fixed-flux setup has a smaller variance, skewness, and characteristic time scale than the fixed-gradient formulation. We believe that these observations are particularly significant: given a specific measured time series of salinity and velocity, we can, in principle, use these distinctive features to determine which formulation is more appropriate for the specific scenario. We can thus determine whether observed salt-finger convection is flux driven or gradient driven.

We mention finally that it is often assumed that the fixed-gradient and fixed-flux formulations are equivalent. This is indeed the case in bounded systems in which the velocity at the upper and lower bounding plates must vanish. This is not the case, however, in doubly periodic (or indeed triply periodic) setups such as ours. We believe that the ability provided by the fixed-flux formulation for the flux to select the associated mean salinity gradient is closer to what happens in real applications where horizontal boundaries are absent and no mechanism exists that can maintain a fixed value of the mean gradient. Moreover, the present approach opens the subject of doubly diffusive convection to formulations that seek to identify flows that maximize fluxes or minimize dissipation.² 

In future work, we plan to extend the approach of this paper to three dimensions, as well as applying it to the primitive equations describing salt-finger convection. We mention that related computations have been performed for triply periodic Rayleigh–Bénard convection forced by a fixed destabilizing temperature gradient (i.e., for a pure fluid) by Calzavarini et al.³ 

We thank Benjamin Miquel for useful discussions. This work was supported, in part, by the National Science Foundation under Grant Nos. OCE-2023499 (KJ) and OCE 2023541 (EK).

The data that support the findings of this study are available from the corresponding author upon reasonable request.

To study the local behavior near solution III, we consider small perturbations of this state satisfying the linear equation

where

We find that the determinant of M is negative,

implying that solution III is unstable.