Turbulent Mixing in Solar Models

• Turbulent mixing in solar models is the process by which turbulence transfers heat, momentum, and chemical species between regions, affecting energy transport and rotation.
• It involves convective and shear instabilities, with models using turbulent diffusion coefficients and non-local convection theories to match helioseismic and spectroscopic data.
• Quantitative analyses show that incorporating turbulent mixing improves surface abundance predictions while posing challenges in reconciling core neutrino fluxes and sound-speed profiles.

Turbulent mixing in solar models describes the processes by which turbulence transports heat, angular momentum, chemical species, magnetic flux, and even radiation within the solar interior and envelope. These processes control energy transport, chemical stratification, dynamo action, evolution of the solar rotation profile, and match between solar models and helioseismic or spectroscopic measurements. Turbulent transport mechanisms operate in both convection zones and stably stratified radiative regions, with distinct physical regimes governed by rotation, magnetic fields, radiative transfer, and chemical gradients.

1. Theoretical Foundations: Turbulent Transport in the Solar Interior

Turbulent mixing arises from convective and shear instabilities. In the solar convection zone, turbulent transport is dominated by strong, nonlinear motions driven by superadiabatic temperature gradients. The effectiveness of turbulent mixing is characterized using turbulent transport coefficients, such as the turbulent thermal diffusivity (η), turbulent viscosity, and diffusivities for chemical species.

In the context of mixing-length theory (MLT), turbulent diffusivity is estimated as

$\eta \simeq \frac{1}{3} \alpha_p H_p v,$

where α_p is the mixing length parameter, H_p the local pressure scale height, and v the convective velocity. First-order turbulence models extend this by relating η to the turbulent kinetic energy and eddy turnover time (τ): $\eta = \frac{\tau}{3} \langle v^2 \rangle.$ For stably stratified radiative interiors, turbulence is highly anisotropic. Turbulent intensities satisfy for strong stratification: $$\frac{\langle u_\phi^2 \rangle}{\langle u_r^2 \rangle} \approx \frac{\tau^2 N^4}{\Omega^2},$$ where N is the buoyancy frequency and Ω the rotation rate. This anisotropy leads to much stronger horizontal mixing than radial mixing (Kitchatinov et al., 2012, Pipin et al., 2013). Angular momentum transport is dominated not by eddy viscosity but by non-diffusive Λ-effects arising from the correlation of radial and azimuthal turbulent motions under rotation.

Overshooting at the base of convection zones is modeled as turbulent diffusion with distance-dependent diffusivity, often parametrized as exponentially decaying from the convective boundary (Zhang, 2012, Zhang et al., 2022). The core overshoot diffusion coefficient, e.g.,

$D(r) = C D_0\, (P/P_{cz})^\theta,$

with parameters C and θ set by helioseismic and neutrino flux constraints.

Horizontal shear-driven instabilities at low Prandtl number can also generate turbulent mixing below the convection zone, as shown in DNS modeling (Garaud, 2020); the vertical eddy size and mixing coefficient in this regime scale nontrivially with stratification (see Section 4 below).

2. Turbulent Mixing, Element Transport, and Surface Abundances

Turbulent mixing competes with atomic diffusion (gravitational settling and radiative acceleration) to determine the evolving chemical composition in the solar interior and atmosphere (Moedas et al., 2022). Atomic diffusion alone predicts unphysically strong abundance gradients in solar-like stars. MESA and other models incorporate a parametrized turbulent diffusion coefficient to simulate efficient mixing: $$D_{\text{t}} = \omega\, D_{\text{He},0} \left(\frac{\rho_0}{\rho}\right)^n,$$ with ω and n as calibration constants and D_{\text{He},0} a reference helium diffusion coefficient at a certain depth.

In detailed solar evolutionary models, turbulent mixing below the base of the convection zone is essential to reconcile predicted and observed solar Li and Be abundances and to suppress excessive gravitational settling of helium. Solar models with only diffusion predict Li abundances ~30σ away from observations, while those including turbulent mixing reproduce both Li and Be surface abundances within ~0.6σ and 0.5σ, respectively, and improve surface helium agreement to ~0.3σ (Kunitomo et al., 1 Sep 2025). The diffusion is typically parameterized as

$$D_{\text{mix}} = D_T \left( \frac{\rho(r)}{\rho(R_{cz})} \right)^{-n}$$

with D_T ≈ 5000 cm²/s and n = 4.

Turbulent mixing, however, also affects core properties and neutrino yields: it lowers the central metallicity and temperature, resulting in tensions with observed ^8B and CNO neutrino fluxes (discrepancies up to 6–8σ in recent models). This underscores that while turbulent mixing resolves surface abundance anomalies, it can exacerbate problems matching interior-sensitive solar observations.

Parametrization of turbulent mixing also enables simulation of the effect of radiative acceleration on heavy elements such as iron without the computational expense of full diffusion+acceleration modeling. For instance, in FG stars, surface [Fe/H] variations are fit by adjusting the turbulent diffusivity to mimic radiative levitation effects, though this approximation fails for elements like Ca or O whose radiative support differs significantly from that of iron (Moedas et al., 2022).

3. Turbulent Mixing and Solar Magnetism: Dynamo Action and Anisotropy

Turbulent mixing, especially turbulent magnetic diffusivity, is a central ingredient in mean-field solar dynamo models. Mixing-length theory predicts η values orders of magnitude higher than those presumed in kinematic dynamo models (which typically employ double-step diffusivity profiles with low values in the bulk of the convection zone). However, when magnetic quenching is included—nonlinear suppression of turbulence by strong fields—the effective (spatiotemporally geometrized) diffusivity profile can be reconciled with empirical models used in dynamo codes (Muñoz-Jaramillo et al., 2010). The effective diffusivity after magnetic quenching is fitted by a double-step analytic profile, validating the operational use of such profiles in global solar dynamo simulations.

Turbulent pumping, an advective process distinct from diffusion, also plays a critical role in Babcock–Leighton–type dynamo models. Downward (radial) and equatorward (latitudinal) pumping can replace deep meridional circulation as the dominant flux transport mechanism, maintaining solar-like cycle periodicity and parity even with only shallow or no meridional flow (Hazra et al., 2016).

Anisotropy in turbulent mixing modifies large-scale magnetic field distributions. The parameter

$$a = \frac{\langle u_h^2 \rangle - 2\langle u_r^2 \rangle}{\langle u_r^2 \rangle}$$

quantifies the horizontal/radial intensity ratio of turbulent mixing. Higher values result in flux concentration at the boundaries and alter cycle periods and overlap in butterfly diagrams (Pipin et al., 2013). This is consistent with increased horizontal mixing in the convectively stable radiative zone, as also shown by τ-approximation calculations (Kitchatinov et al., 2012).

4. Turbulent Mixing at the Convective Boundary, Tachocline, and Radiative Interior

Mixing at the base of the convection zone, in the tachocline, and in radiative regions is strongly influenced by the nature of turbulence, rotation, and magnetic fields.

Direct numerical simulations of horizontal shear instabilities at low Prandtl number reveal distinct turbulent regimes defined by the product $BPe$ (with $B = N^2L^2/U^2$, $Pe = UL/\kappa_T$):

• For high BPe (solar tachocline regime), the vertical eddy scale and vertical turbulent diffusivity scale as

$$l_z \simeq 2.1 \left(\frac{N^2L^2}{U^2}\right)^{-1/3} L, \qquad D_{\text{turb}} \sim 2.7 \left(\frac{N^2L^2}{U^2}\right)^{-2/3} UL.$$

The predicted vertical diffusion coefficient can be orders of magnitude higher than assumed in classical models, challenging thin tachocline models that require extreme anisotropy and weak vertical mixing (Garaud, 2020). These results, however, do not include the effects of rotation (which reduces the Rossby number to ~0.1 on tachocline eddy scales) or magnetic fields, both of which are expected to further suppress or modify turbulence. The development of a comprehensive theory for the tachocline requires incorporating these effects.

Spectrally resolved 3D simulations confirm that turbulent mixing properties and the dominant scale of motion change with depth, with convective scales growing below granulation and changing in character near the bottom of the hydrogen ionization zone; this alters the diffusivity structure and therefore the mixing efficiency (Kitiashvili et al., 3 Feb 2025).

5. Non-Local Turbulent Convection Models and Mean-Field Theories

Non-local turbulent convection models attempt to supersede local mixing-length theory by modeling the spatial (and temporal) propagation of turbulence, including turbulent kinetic energy and in some models, convective flux and entropy fluctuations as dynamic variables. The 1-equation Kuhfuss model, implemented in GARSTEC, evolves turbulent kinetic energy (ω) and assumes the convective flux π is proportional to the entropy gradient (a downgradient closure): $$\pi = -\alpha_s \Lambda \omega^{1/2} \frac{\partial s}{\partial r}.$$ This approach automatically produces convective boundary mixing (CBM, or overshooting) determined by nonlinear transport and advection of ω beyond the formal Schwarzschild boundary (Braun et al., 4 Jul 2024). However, in solar models, the 1-equation Kuhfuss model with downgradient closure enforces a nearly adiabatic stratification too deeply below the convection zone, placing the base of the adiabatic region at substantially lower radius (e.g., 0.6845 R☉ versus the helioseismic value of ≈0.713 R☉), and resulting in significant discrepancies in the sound-speed profile. The underlying issue is a tight, local coupling of the convective heat flux to the entropy gradient; more sophisticated (e.g., 3-equation) nonlocal models, which evolve the heat flux as an independent variable, can reproduce a smoother, non-monotonic transition in the boundary region and improve agreement with seismic constraints.

6. Turbulent Mixing, Radiative Transport, and Penetration of Radiation

Turbulent mixing not only advects and diffuses chemical species and momentum, but also modifies radiative transfer by increasing the effective penetration of radiation through fluctuating opacity. Mean-field radiative transfer theory for turbulent media yields an effective penetration length

$$L_{\text{eff}} = L_r \left[1 + \frac{2 \langle \kappa'^2 \rangle}{\kappa^2} \left(1 + L_r/\ell_0\right)^{-1}\right],$$

where L_r is the mean radiative penetration length, κ the mean absorption coefficient, and ℓ_0 the turbulent correlation length (Rogachevskii et al., 2021). In the upper solar convection zone, where compressibility enhances temperature and density fluctuations (and thus $\kappa'$), the effective penetration length increases by a factor up to 2.5 compared to the mean. This turbulence-induced “softening” of radiative absorption is most significant near the photosphere, but can be more dramatic in cooler stars with higher Mach numbers. These findings motivate incorporation of turbulent radiative transport corrections in solar and stellar envelope models.

7. Open Issues and Future Prospects

Despite significant progress, multiple challenges remain:

• Achieving simultaneous agreement between solar models and all observational diagnostics (neutrino flux, atmospheric abundances, sound-speed profile, depth of the convection zone) remains elusive when turbulent mixing is incorporated; improvements in input microphysics (opacities, electron screening, nuclear rates) and better understanding of accretion history, wind enrichment, and the tachocline are necessary (Kunitomo et al., 1 Sep 2025).
• Parametric forms of turbulent mixing and their calibration often depend on the method and observational probe (asteroseismology, helioseismology, spectroscopic abundances).
• The impact of rotation, magnetism, and time dependence on mixing efficiency—particularly in radiative interiors and boundary layers—requires further direct numerical simulations and improved analytic theory.

Reductions in the structural “surface effect” in oscillation frequencies have been achieved by appending ⟨3D⟩ radiative hydrodynamical envelope models with turbulent pressure to 1D stellar models on the fly (Jørgensen et al., 2019). However, comprehensive grids of 3D simulations and more advanced coupling algorithms are still being developed to fully exploit these methodologies across the HR diagram.

In conclusion, turbulent mixing is a multifaceted, crucial ingredient of solar modeling. It mediates the interior structure, composition, rotation, and observables of the Sun and solar-type stars, and will remain central to the evolution of stellar astrophysics as constraints and models become increasingly precise and physically complete.