THERMOHALINE MIXING: DOES IT REALLY GOVERN THE ATMOSPHERIC CHEMICAL COMPOSITION OF LOW-MASS RED GIANTS?

Red giants with M ≲ 1.5 M_☉ are known to experience extra-mixing in their convectively stable radiative zones that separate the hydrogen burning shell (HBS) from the bottom of the convective envelope (BCE). It operates during their first ascent along the red giant branch (RGB) above the bump luminosity, after the HBS has crossed a chemical composition discontinuity left behind by the BCE at the end of the first dredge-up. This causes the atmospheric ¹²C/¹³C ratio and carbon abundance to resume their declines with an increase in luminosity.

The most promising physical mechanism of RGB extra-mixing is thermohaline convection (Charbonnel & Zahn 2007a). It is driven by a double-diffusive instability that occurs when a scalar component that stabilizes a density stratification, e.g., the temperature T, diffuses away faster than a destabilizing component, e.g., the mean molecular weight μ. Eggleton et al. (2006) brought our attention to the fact that the reaction ³He(³He, 2p)⁴He produces a local depression of μ, Δμ ∼ −10⁻⁴, in the HBS tail above the bump luminosity. As a result, μ is increasing with the radius, which turns the μ profile into the destabilizing component. Moreover, the radiative thermal diffusivity K exceeds the molecular diffusivity ν_mol by many orders of magnitude (the inverse Lewis number τ = ν_mol/K ∼ 10⁻⁶ in the region of μ depression), which should lead to the double-diffusive instability. In the ocean, a similar instability develops where warm salty water overlies cold fresh water. The oceanic thermohaline mixing usually takes the form of vertically elongated fluid parcels of rising fresh and sinking salty water adjacent to each other, called "salt fingers" (Stern 1960). It has extensively been studied both experimentally and theoretically (e.g., Kunze 2003).

The efficiency of thermohaline mixing in stellar cases can be estimated only theoretically, e.g., via a linear stability analysis of its governing equations (Ulrich 1972) or via direct numerical simulations. Denissenkov (2010, hereafter Paper I) has recently conducted a comparative study of thermohaline convection in the oceanic and RGB cases using both a linear stability analysis and two-dimensional (2D) numerical simulations. Unfortunately, the linear theory yields a thermohaline diffusion coefficient that is proportional to the square of the (unknown) aspect ratio, a = l/d, of a finger, where l and d are the finger's length and diameter. Applications of this diffusion coefficient to model the evolutionary decline of carbon abundance in low-mass red giants demand that a>7 (Paper I), which is in agreement with the original result of Charbonnel & Zahn (2007a). However, our 2D numerical simulations of ³He-driven thermohaline convection in the vicinity of a μ depression yield a mixing rate equivalent to that approximated by the linear-theory diffusion coefficient with a < 1, i.e., it is a factor of ∼1/50th of its observationally constrained value. The difference in the finger aspect ratios between the oceanic (a>1) and RGB (a < 1) cases is most likely determined by the very low viscosity ν, or the Prandtl number Pr = ν/K ∼ 10⁻⁶, in the RGB case. This facilitates the development of secondary shear instabilities that strongly limit the growth of salt fingers (Radko 2010 and Paper I). On the other hand, the diffusion coefficient for thermohaline convection with a>7 has been shown to reproduce many of the RGB extra-mixing observational patterns at all metallicities so well that we had to postpone further discussion of these findings until we were able to confirm or refute the results of our 2D numerical simulations using three-dimensional (3D) computations.

In this Letter, we present first results of our 3D numerical simulations of thermohaline convection for the oceanic and RGB cases. It turns out that they do not differ much from the 2D results reported in Paper I. We also consider a possibility, mentioned in Paper I, that RGB thermohaline mixing may be enhanced by large-scale mean-field instabilities, such as the collective, γ, and intrusive ones, that are known to occur in the ocean. Finally, we speculate on possible alternative solutions of the problem of RGB extra-mixing.

The basic equations and numerical techniques that are employed in our 3D simulations of thermohaline convection for the oceanic and RGB cases are described in detail by Gargett et al. (2003). They are very similar to those used in Paper I. The main goal of Gargett et al. was to study differential scalar diffusion in 3D stratified turbulence. Therefore, they considered a density stratification stable against the double-diffusive instability and used a specified turbulent velocity field as an initial condition. To apply their method and computer code for our purposes, we have made the slower diffusing component, which represents salinity S in the ocean and μ in the current investigation, destabilizing by changing the sign of the vertical velocity term in their Equation (16). We have also reduced the amplitude of the initial velocity field, so that it serves as a mild perturbation to initiate the double-diffusive instability. Our available computational resources have allowed us to run the code within reasonable time intervals with a spatial resolution of 320³. Simple estimates of the relevant turbulent kinetic energy dissipation rate and the Kolmogorov and Batchelor length scales show that this resolution is close to marginal for our runs, however, not fully resolving the slower-diffusing component has been shown not to have a large effect on the computed fluxes (Traxler et al. 2010). The 3D simulations have been carried out for the same number of fastest growing fingers in the computational domain and parameter sets from Table 1 in Paper I that were used in our 2D computations.

The blue and red solid curves in the upper left panel of Figure 1 show transitions of the ratio D_S/k_T to its equilibrium values in our 2D and 3D simulations of the oceanic thermohaline convection. Here, D_S and k_T are the turbulent salt and microscopic thermal diffusivities. The 2D simulations were performed with the resolution 1024² (Paper I). Our estimated 3D equilibrium value of D_S/k_T exceeds its corresponding 2D value merely by a factor of 1.5 (compare the blue and red dashed lines). It is in good agreement with the result obtained using a better resolution by Traxler et al. (2010) that we have interpolated for the value of density ratio R_ρ = αΔT/βΔS = 1.6 (α and β are the coefficients of thermal expansion and haline contraction) from their Table 1 (star symbol). The upper right panel makes a similar comparison for the RGB thermohaline convection. In this case, we do not find any difference between the 2D and 3D equilibrium values of the ratio D_μ/K (the blue and red dashed lines coincide). Therefore, we conclude that our results reported in Paper I seem to be correct. In particular, RGB thermohaline convection driven by ³He burning yields a turbulent mixing rate that is by a factor of 50 lower than the value required by observations.

Zoom In Zoom Out Reset image size

Observations, laboratory experiments, and direct numerical simulations show that small-scale salt-fingering convection in the ocean can lead to instabilities and the formation of dynamical structures on much larger scales. These include the so-called collective instability leading to the generation of internal gravity waves, the intrusive instability developing when fluid is stratified both vertically and horizontally (in stars, such a situation may result from the rotational distortion of level surfaces), and the so-called γ-instability that produces thermohaline staircases (e.g., see Traxler et al. 2010, and references therein). The latter structures consist of thick well-mixed convective layers separated by thin salt-fingering interfaces. The staircases are of special interest to us because they have been observed to enhance vertical mixing in the ocean by up to an order of magnitude (Schmitt et al. 2005). Moreover, our preliminary estimates in Paper I indicate that the turbulent flux ratio γ = F_T/F_μ is a decreasing function of R_ρ in the RGB case, which is a sufficient condition for the γ-instability, according to Radko (2003).

The large-scale instabilities are studied by averaging the Boussinesq and advection–diffusion equations for T and S (μ) over length and timescales of many salt fingers. As a result, a system of mean-field equations is obtained in which the small-scale turbulent velocity is incorporated into the mean heat and salinity fluxes F_T and F_S (Traxler et al. 2010). The fluxes are related to one another and to the large-scale temperature gradient via the flux ratio γ = F_T/F_S and the Nusselt number Nu, which are assumed to depend only on the local value of the density ratio R_ρ; specification of these dependences closes the mean-field equations. Large-scale changes in R_ρ thus result in modulations of F_T and F_S via the dependences of γ and Nu on R_ρ that can either amplify or damp initial perturbations.

A linear stability analysis of the mean-field equations yields a third-order dispersion relationship for the growth rates of all the three large-scale instabilities mentioned above (Equation (2.13) of Traxler et al. 2010). Coefficients of the cubic equation depend on the Prandtl and inverse Lewis numbers as well as on Nu, γ, and their derivatives with respect to R_ρ. The coefficients are also functions of the horizontal and vertical wave numbers, l and k. To calculate them, we first need to estimate the dependences of Nu and γ on R_ρ. Traxler et al. (2010) estimated those dependences via a body of 3D numerical simulations of the oceanic salt-fingering convection. We have used data from their Table 2 to reproduce, for a comparison, their "flower plot" depicting the growth rates of fastest growing instabilities at R_ρ = 1.5 in the absence of lateral gradients (the lower left panel of Figure 1). The bulb of the flower corresponds to the salt-fingering mode that turns out to be a solution of the general dispersion equation as well. The leaves of the flower represent the gravity waves of the collective instability, whereas the γ-instability modes are the ones dominating in the region beneath the leaves. To obtain estimates of γ(R_ρ) and Nu(R_ρ) for the RGB case (Figure 2), we have performed additional 2D numerical simulations of the ³He-driven thermohaline convection for the same parameter set that was used in Paper I but for a number of different values of R_ρ. The resulting flower plot is shown in the lower right panel of Figure 1. We see that in the RGB case only the salt-fingering mode is unstable in spite of the flux ratio γ being a decreasing function of R_ρ, as in the oceanic case. This difference is explained by the fact that the Nusselt number is extremely small, Nu ∼ 10⁻⁶, in the RGB case. Because of this, it is necessary to include horizontal diffusive fluxes in the mean-field stability analysis—something which was not done by Radko (2003) but was treated by Traxler et al. (2010). In the RGB case, the inclusion of reasonable lateral gradients proportional to their corresponding vertical gradients multiplied by the ratio of centrifugal to gravitational acceleration only produces an asymmetry of the bulb with respect to the vertical axis.

Zoom In Zoom Out Reset image size

Our 2D and 3D numerical simulations of thermohaline convection driven by ³He burning have demonstrated that it is unlikely to be the sole mechanism of RGB extra-mixing. Its corresponding mixing rate in the vicinity of a μ depression, D_μ ≈ 2 × 10⁻³ K (the upper right panel of Figure 1), is too low compared to the observationally constrained value of (D_μ)_obs ≈ 0.1 K (Paper I). Therefore, for its efficiency to be consistent with observations, the RGB thermohaline mixing has to be either assisted or modified by other processes. In Paper I, we have shown that a turbulent viscosity ν_t exceeding the microscopic one by a factor of 10⁴ could enhance the rate of ³He-driven thermohaline mixing up to the observed value, provided that it was not itself associated with strong turbulent mixing. The enhancement of thermohaline mixing in this case is caused by the fact that the higher viscosity suppresses the development of secondary shear instabilities that limit the growth of salt fingers. On the other hand, if the increase of viscosity is accompanied by a similar increase of the rate of turbulent mixing, the latter will reduce salt-finger buoyancy and vertical transports by smoothing out the chemical composition contrast between rising and sinking fluid parcels.

The additional source of turbulence in the radiative zones of RGB stars could arise from the differential rotation that is produced by the mass inflow that feeds the HBS (Paper I). We have employed the COMSOL Multiphysics software package to solve the angular momentum transport equation for the radiative zone of the low-metallicity bump-luminosity RGB model from Paper I in the presence of mass inflow and turbulent shear mixing. Their respective rate $\dot{r}$ and diffusion coefficient ν_v are given by Equation (23) from Paper I and Equation (3) with the uncertainty factor f_v from the paper of Denissenkov et al. (2006). The angular velocity at the BCE was assumed to have the same value, Ω_BCE = 6.3 × 10⁻⁷ rad s⁻¹, as in the M2 bump-luminosity RGB model of Palacios et al. (2006). Simple estimates show that for f_v ≳ 1 the transport of angular momentum by the rotation-induced turbulent diffusion dominates over its advection by meridional circulation, therefore, and for the sake of simplicity, the latter was neglected. Our computations show that the competition between the terms proportional to $\dot{r}$ and ν_v leads to a stationary solution with ν_v ≈ 3 × 10⁻⁴ K in the vicinity of a μ depression that very weakly depends on both f_v and Ω_BCE, although the final rotational shear is a decreasing function of f_v. This value exceeds the microscopic viscosity by a factor of 10², which is not enough to enhance D_μ to its observed value (Paper I). It is also difficult to understand why the rotation-induced small-scale turbulence should not produce a mixing rate of the same magnitude as ν_v, in which case D_μ will not be enhanced at all (Paper I). Therefore, we do not think that the proposed viscosity enhancement is a possible solution of the problem.

Our angular momentum transport computations confirm the conclusion that radiative zones of low-mass RGB stars above the bump luminosity should possess strong differential rotation. This may lead to the generation of high-amplitude toroidal magnetic fields, especially close to the HBS, provided that a weak poloidal field is present in the radiative zone. Under these circumstances, a possible alternative solution of the problem could be magneto-thermohaline mixing (Denissenkov et al. 2009). We intend to investigate this possibility using 2D numerical simulations in one of our forthcoming papers.

Toroidal magnetic fields in radiative zones may also filter out salt fingers with small diameters, so that the maximum growth rate is shifted toward thicker fingers. Since the finger growth rate decreases with an increase in diameter as σ ∝ d⁻² (e.g., Paper I), the square of the vertical velocity shear between the axes of neighboring rising and sinking fingers, ∼(σl/d)², is proportional to a²d⁻⁴. This parameter is likely to control the development of secondary shear instabilities (Kunze 1987); consequently, thicker fingers will probably reach higher aspect ratios before they are destroyed. Following Charbonnel & Zahn (2007b), we have complemented the dispersion relationship (13) from Paper I (note that the term νk_Tk_Sk⁶ was erroneously omitted) with terms that take into account the Lorentz force associated with a specified toroidal magnetic field B₀. Its solutions are shown in Figure 3 for the field strengths B₀ = 0, 10, and 100 G (the solid blue curves from left to right). In the case of B₀ = 100 G, the fastest growing fingers transverse to the field have a diameter that is nearly 40 times as large as that of fingers in the absence of a magnetic field. This diameter is also almost four times larger than the diameter of fastest growing fingers in the case of a viscosity that is amplified by the factor of 10⁴ (the red solid curve). However, the influence of such a field on the planform of the fastest growing fingers still needs to be studied. In our next study, we will perform 3D numerical simulations of the RGB thermohaline convection in the presence of a toroidal magnetic field to find out if the magnetic shift of the salt-finger growth-rate spectrum toward larger diameters can really enhance its associated mixing rate up to the empirically constrained value.

Zoom In Zoom Out Reset image size

P.A.D. is grateful to Falk Herwig and Don VandenBerg who have supported this work through their Discovery Grants from Natural Sciences and Engineering Research Council of Canada.