Effects of small-scale dynamo and compressibility on the $\Lambda$ effect

The $\Lambda$ effect describes a rotation-induced non-diffusive contribution to the Reynolds stress. It is commonly held responsible for maintaining the observed differential rotation of the Sun and other late-type stars. Here the sensitivity of the $\Lambda$ effect to small-scale magnetic fields and compressibility is studied by means of forced turbulence simulations either with anisotropic forcing in fully periodic cubes or in density-stratified domains with isotropic forcing. Effects of small-scale magnetic fields are studied in cases where the magnetic fields are self-consistently generated by a small-scale dynamo. The results show that small-scale magnetic fields lead to a quenching of the $\Lambda$ effect which is milder than in cases where also a large-scale field is present. The effect of compressibility on the $\Lambda$ effect is negligible in the range of Mach numbers from 0.015 to 0.8. Density stratification induces a marked anisotropy in the turbulence and a vertical $\Lambda$ effect if the forcing scale is roughly two times larger than the density scale height.


Introduction
Solar and stellar differential rotation is thought to arise due to the interaction of density-stratified convective turbulence and global rotation of the star (e.g. Rüdiger 1980Rüdiger , 1989Rüdiger et al. 2013). While turbulence is often associated with only enhanced diffusion, several non-diffusive effects have also been discovered. Arguably the most well-known of these in the astrophysical context is the α effect which leads to the generation of large-scale magnetic fields in helical turbulence (Steenbeck et al. 1966). In the hydrodynamic (HD) context, a non-diffusive contribution to the Reynolds stress, also known as the Λ effect, is thought to be crucial for the maintenance of stellar differential rotation (e.g. Kichatinov & Rüdiger 1993;Kitchatinov & Rüdiger 1995;Kitchatinov & Rüdiger 2005;Rüdiger 1980). Considerable observational (e.g. Rüdiger et al. 2014, and references therein) and numerical (e.g. Käpylä 2019; Käpylä & Brandenburg 2008;Pulkkinen et al. 1993;Rüdiger et al. 2005) evidence support the existence of the Λ effect in flows akin to those in stellar convection zones. The Λ effect occurs in rotating anisotropic turbulence which means that angular momentum transport in accretion disks is also likely affected by it (e.g. Käpylä et al. 2010;Snellman et al. 2009). Other non-diffusive HD effects include the anisotropic kinetic alpha (AKA) effect (e.g. Frisch et al. 1987;Käpylä et al. 2018) and the inhomogeneous helicity effect (Yokoi & Brandenburg 2016), but their role in the maintenance of stellar differential rotation is likely to be subdominant to the Λ effect.
Corresponding author: pkaepyl@uni-goettingen.de Numerical simulations of magnetohydrodynamic (MHD) convection in spherical coordinates have reached sufficient spatial resolution that allow the excitation of small-scale dynamo action (e.g. Hotta et al. 2014;Käpylä et al. 2017;Nelson et al. 2013). The study of Käpylä et al. (2017) showed that differential rotation in simulations is strongly quenched at the highest magnetic Reynolds numbers where an efficient small-scale dynamo is excited. Furthermore, the turbulent Reynolds and Maxwell stresses were found to have similar spatial distributions and magnitudes but opposite signs. These findings can be interpreted as magnetic quenching of the Λ effect.
In a subsequent study (Käpylä 2019), the effect of largescale magnetic fields on the Λ effect was studied. These results show that the Λ effect is significantly quenched when the large-scale magnetic field reaches a substantial fraction of the equipartition strength. However, imposing a largescale field will also induce small-scale fields due to tangling by the turbulent motions and it is not possible to disentangle the two contributions. Here this caveat is avoided by self-consistent generation of small-scale magnetic fields by a small-scale dynamo in a setup where no simultaneous large-scale dynamo is present.
While the mean-field theory of the Λ effect is derived under the assumption of incompressibility, the numerical simulations used to compute the coefficients are most often fully compressible (e.g. Käpylä & Brandenburg 2008;Käpylä et al. 2004;Pulkkinen et al. 1993). Although the Mach numbers in these studies are still clearly subsonic, the effects of compressibility have not been studied in detail. Here such an effort is undertaken with a controlled set of simulations where the minimal ingredients (rotation and Copyright line will be provided by the publisher arXiv:1903.04363v2 [astro-ph.SR] 17 Feb 2020 anisotropic turbulence) for the Λ effect are included while the Mach number is varied.
Another aspect that has not received much attention is the contribution of density stratification to the anisotropy of turbulence (see, however Brandenburg et al. 2012) and the resulting Λ effect in rotating cases. Isolating these effects in convection is not possible because the forcing due to the convective instability is in itself highly anisotropic. Here this aspect is studied with isothermal but density stratified setups where the turbulence is driven by isotropic forcing and where the convective instability is absent.

The model
The model is the same as that used in Käpylä & Brandenburg (2008) and Käpylä (2019), except in cases where gravity and density stratification are included.

Basic equations
Compressible HD or MHD turbulent flow in a fully periodic cube is modeled. An isothermal equation of state with p = ρc 2 s , where p is the pressure and c s is the constant speed of sound, is assumed. The following set of MHD equations is solved: where A is the magnetic vector potential, U is the fluid velocity, B = ∇ × A is the magnetic field, J = µ −1 0 ∇ × B is the current density, η is the magnetic diffusivity, µ 0 is the permeability of vacuum, ρ is the density, D/Dt = ∂/∂t − U · ∇ is the advective time derivative, g = (0, 0, −g) is the acceleration due to gravity, Ω is the rotation vector, and F visc and F force describe the viscous force and external forcing, respectively.
The viscous force is given by where ν is the kinematic viscosity and is the traceless rate of strain tensor. The forcing term on the rhs of Eq.
(2) is given by where x = (x, y, z), N = fc s (kc s /δt) 1/2 is a normalization factor, f contains the dimensionless amplitudes of the forcing, k = |k|, δt is the length of the time step, and −π < φ(t) < π is a random delta-correlated phase. The vector f k describes nonhelical transversal waves, with whereê is an arbitrary unit vector and where the wavenumber k is randomly chosen at every time step. The Pencil Code 1 was used to perform the simulations.

Units, system parameters, and boundary conditions
The units of length, time, density, and magnetic field are where k 1 is the wavenumber corresponding to the scale of the domain and ρ 0 is the initially uniform value of density.
The forcing amplitude f i j is given by where f 0 and f 1 are the amplitudes of the isotropic and anisotropic parts, respectively, δ i j is the Kronecker delta, and Θ k is the angle betweenê z and k. The forcing wavenumber is chosen from a narrow range around a predefined wavenumber k f . The Mach number of the flow is varied by adjusting the sound speed c s and the forcing amplitudes f 0 and f 1 . The rotation vector is given by Ω = Ω 0 (− sin θ, 0, cos θ) T , where θ is the angle with the vertical (z) direction. Viscosity and rotation can be combined into the Taylor number where L d = 2π/k 1 corresponds to the size of the computational domain. Furthermore, in the MHD cases the magnetic Prandtl number is an additional system parameter. The density scale height H ρ = c 2 s /g describes the stratification in cases with g 0. In cases with g = 0 the system is fully periodic. When g 0, impenetrable stress-free boundary conditions corresponding to: are enforced at the vertical (z) boundaries.

Diagnostics quantities
The following quantities are outcomes of the simulations that can only be determined a posteriori. The fluid and magnetic Reynolds numbers are given by The rotational influence on the flow is quantified by the Coriolis number based on the forcing scale Table 1 Summary of the SSD simulations. Each set consists of ten runs where θ is varied in steps of 10 • . All runs havek f = k f /k 1 = 10, Ta = 6.2 · 10 7 , ∆ρ = 1, c s = 3, Ω = 0.9, Ma = 0.05, and A V = −0.5.
The magnetic field strength is given in terms of the equipartition value Finally, the parameter characterises vertical anisotropy of turbulence. The density stratification is quantified by the ratio of the densities at the top and bottom of the domain, ∆ρ = ρ z bot /ρ z top where z bot k 1 = −π and z top k 1 = π.

Data analysis
The coefficients pertaining to the Λ effect were extracted by fitting the latitudinal profiles of the off-diagonal Reynolds stresses with the same procedure as in Käpylä (2019). The Reynolds stress is given by Q i j = u i u j where the overline denotes horizontal averaging and where u = U − U is the fluctuating velocity. In the homogeneous cases an additional z-averaging is performed. The fitting procedure assumes that the off-diagonal Reynolds stresses are solely generated by the Λ effect and that they can be represented as where ν t = 2 15 u rms is an estimate of the turbulent viscosity, and The expansions can in principle contain an arbitrary number of higher powers of sin 2 θ but here the simulations were made in such a regime that adding higher order contributions to the coefficients does not yield a significantly improved fit (see Käpylä 2019). Error estimates were computed by dividing the time series in three parts and averaging over each part. The greatest deviation of these from the average over the full data set was taken to represent the error.

Results
Three sets of simulations were made to study the effects of small-scale magnetic fields (Set SSD, Table 1), Mach number (Set MA, Table 2), and density stratification (Set STR, Table 3). In Sets SSD and MA the system is homogeneous and anisotropically forced while in Set STR the forcing is isotropic and a strong density stratification is present. Visualizations of typical flow patterns for three representative runs are shown in Figure 1. A somewhat surprising result is that the presence of strong density stratification is not clearly visible from the flow patterns, compare the left and right panels of Figure 1.   case the forcing is no longer solenoidal. The spectrum is clearly steeper than the Kolmogorov (1941) (hereafter K41) −5/3 prediction. This is most likely due to the insufficient scale separation between the forcing and viscous scales which does not allow the formation of an inertial range. Indeed, even the highest resolution simulations up to date with 8192 3 grid resolution are able recover only a rather modest well-defined inertial range (Iyer et al. 2017). The viscous scale is now well resolved due to the relatively modest Reynolds number (Re ≈ 15) in the current simulations. It is also evident that the turbulence is anisotropic at all scales all the way down to the grid scale, see the red, blue and yellow curves in Figure 2. This is quantified by a spectral analogy of the anisotropy parameter:

Anisotropy of turbulence
where E K (k) is the power spectrum of the total velocity and E i (k) are the power spectra of the individual velocity com-

Fig. 3
Anisotropy parameter A V (z) from the runs in Set STR with varyingk f . Figure 2 for Run MA7. A V (k) has a minimum at k = k f which is due to the fact that the forcing mainly puts energy in the z component of the velocity in this case.

ponents. A representative result for A V (k) is shown in the inset of
While A V (k) is mostly negative for large scales, it gradually increases and peaks around 0.5 fork ≈ 100. The fact that the anisotropy survives to the smallest resolved scales is in apparent disagreement with one of the cornerstones of the K41 theory which assumes that the turbulence is fully isotropic at small enough scales. However, the current simulations operate in a very modest Reynolds number regime which cannot be directly compared with the K41 theory which formally applies to fully developed turbulence at very high Re.
In addition to using explicitly anisotropic forcing, setups where anisotropy arises naturally due to gravity and density stratification are studied in Set STR. By virtue of the isothermal equation of state this setup is characterised by a constant density scale height H ρ = c 2 s /g = 9/10 such that simulation domain contains seven scale heights and a density contrast ∆ρ of more than a thousand. Although a modest resolution of 288 grid points was used, each scale height is covered by more than 40 grid points. This differs from the case of convection where the density (and pressure) scale height varies strongly as a function of depth and imposes much more restrictive constraints on the grid size (see, e.g. Käpylä et al. 2016). Similar setups, albeit with somewhat lower stratification, were used in an earlier study by Brandenburg et al. (2012). They showed that turbulence anisotropy remains small when isotropic forcing is used unless the forcing scale is larger than the density scale height. This is confirmed by the current simulations where the scale separation ratio, quantified by the ratio of the forcing and system scales k f = k f /k 1 , is varied between 1.5 and 10, see Figure 3 and Table 3 where a z-dependent variant of Equation (17) has been used. The density stratification-induced anisotropy is almost non-existent in the bulk of the domain in the case of the largest scale separationk f = 10 or /H ρ = 0.7. The stress-free and impenetrable boundary conditions enforce U z = 0 and lead to A V = 1 at the vertical boundaries. In the cases with poorer scale separation or larger /H ρ , A V tends to become more positive in the deep parts and obtains negative values near the surface. For the poorest scale separation (Run STR1, /H ρ = 4.5) the magnitude of the anisotropy is comparable to typical values achieved with anisotropic forcing in Sets SSD and MA.

Small-scale magnetic fields
Testing the dependence of purely small-scale magnetic fields is possible with the current homogeneous setups in cases where the magnetic Reynolds number exceeds that of the critical value for the excitation of a small-scale dynamo. Due to the absence of inhomogeneities, large-scale shear or helicity, no large-scale magnetic fields are expected to develop. A limited range of magnetic field strengths has been studied in the Set SSD, see Table 1. Run SSD1 with Re M = 27 corresponds to a slightly supercritical case whereas Run SSD4 corresponds to the highest Re M (≈ 135) that can be resolved with the adopted grid resolution. The saturation levelb rms = B rms /B eq , where B rms is the rmsvalue of the magnetic field, increases from roughly 10 to 50 per cent of the equipartition value in this range. In practice the characteristics of the flow, that is the Reynolds and Coriolis numbers and the degree of anisotropy, were kept fixed in all runs by adjusting f 0 and f 1 while Re M was varied by changing Pm.
The results for the Λ coefficients corresponding to Equations (24) to (26) from Sets SSD1-4 are shown in Figure 4. All of the coefficients are quenched as the small-scale magnetic fields increase: the values forb rms ≈ 0.5 are typically roughly half of their hydrodynamic values. It is also evident that the quenching as a function of pure small-scale fields is weaker than that in the cases where an imposed large-scale vertical field is present, see the dashed lines in Figure 4 for the corresponding data from Käpylä (2019).

Dependence on Mach number
Although the Mach numbers in the foregoing studies (Käpylä 2019;Käpylä & Brandenburg 2008) were typically relatively low (O(0.1)), it cannot be ruled out that a contribution due to compressibility is present. Furthermore, results from low-Reynolds number shear flows (Rogachevskii et al. 2011) suggest that compressibility significantly affects turbulent pumping for Ma 0.1.
Results for the Mach number dependence from Set MA where Ω = 0.8, A V = −0.4, and Re = 15 are kept fixed are shown in Figure 5. The range of Mach numbers spans from 0.015 to 0.8. Only the vertical stress Q yz for Ma ≈ 0.8 is statistically significantly different from the values obtained for lower Ma and even there the change is only on the or- Fig. 5 Off-diagonal Reynolds stresses Q xy (black), Q xz (red), and Q yz (blue) as functions of Mach number from Set MA, see Table 2. der of ten per cent. Temporal fluctuations of Q xz increase, manifested by the drastically increased error estimates, as a function of Ma but the time-averaged values for all runs are still consistent with a Ma-independent value.

Dependence on density stratification
In the foregoing analysis the Λ effect resulted from the forcing that was designed to be anisotropic. While this is the case also in natural convection, a background density stratification also leads to anisotropy and should hence support a Λ effect. This is indeed predicted by analytic theories (e.g. Kitchatinov & Rüdiger 2005;Pipin & Kosovichev 2018). Here this scenario is tested with a density-stratified setup with isotropic forcing in the Set STR, see Table 3.
The simulations considered here have Ω ≈ 0.8 which is close to the Coriolis number where the Reynolds stresses obtain a maximum in Käpylä (2019). Furthermore, the simulation domain is situated at the equator at θ = 90 • where only the vertical Reynolds stress Q yz is non-zero. In order to isolate the contribution relevant for the Λ effect, the horizontally averaged horizontal mean flows U x and U y are artificially removed from the solution similarly as in Rüdiger et al. (2019).
The results for the vertical stress Q yz are shown in Figure 6. The stress is positive everywhere even in the nearsurface regions where A V < 0. The results indicate that the scale separation ratio can play an important role for the density-induced anisotropy and consequently for the associated Λ effect: both are substantial in the case of large-scale forcing and tend to approach zero whenk f increases. The values of Q yz are of the order of a few per cent of u 2 rms for /H ρ = 4.5 whereas for /H ρ = 0.7 the effect is no longer statistically significant at the equator. Figure 7(a) shows the normalized vertical stressQ yz from seven latitudes for the STR4 runs. The maximal values are generally on the order of one to two per cent of the squared rms-velocity. In this set the maximum value is obtained at θ = 45 • and only very small values are obtained at the equator (θ = 90 • ), see Figure 7(b). This, however, depends on the scale separation ratio because for a corresponding set with /H ρ = 4.5, Q yz is consistent with a monotonic increase toward the equator. In all of these runs the sign of A V differs from the sign of Q yz . The mismatch of the signs of Q yz and A V also suggests that non-locality may play a significant role when the scale separation ratio is small.
Another contribution to the Λ effect due to a vertical gradient of the Coriolis number was discussed recently by Pipin & Kosovichev (2018). However, this effect is relevant only for large Coriolis numbers and thus not applicable here. Furthermore, the local Coriolis number Co = 2Ω 0 /u rms (z) varies by a factor between two and three in the current simulations (see the sixth column of Table 3) which is mild in comparison to the variation of four orders of magnitude in the solar convection zone as considered by Pipin & Kosovichev (2018).

Conclusions
The effects of small-scale magnetic fields, compressibility and background density stratification on the Λ effect were studied with numerical simulations of forced turbulence with an isothermal equation of state.
The small-scale magnetic fields generated by a smallscale dynamo lead to a significantly milder quenching of the Λ effect in comparison to cases where also a uniform large-scale field is imposed (see, e.g. Käpylä 2019). Thus is appears that the small-scale dynamo alone could not explain the severely quenched differential rotation in recent semi-global convection simulations (Käpylä et al. 2017). It is also conceivable that other MHD instabilities, such as the magnetorotational instability (e.g. Masada 2011), can be excited in the high-resolution convection simulations, leading to repercussions for differential rotation.
Another aspect that has hitherto received little attention is the Mach number dependence of the Λ effect although Fig. 7 (a) Normalized vertical Reynolds stress Q yz as a function of height from runs with /H ρ = 0.7 and Ω = 0.8 from different colatitudes as indicated by the legend. (b)Q yz from z k 1 = 2 as functions of θ from two sets with /H ρ = 0.7 (black) and /H ρ = 4.5 (blue). The curves show best fits according to Eq. (25) with V (0) = 0.011, V (1) = 0.027, V (2) = −0.036 for /H ρ = 0.7 and V (0) = 0.010 for /H ρ = 4.5 most numerical studies of the subject operate in a fully compressible regime (e.g. Käpylä 2019; Käpylä & Brandenburg 2008). The current results indicate that while the fluctuations of the Λ coefficients tend to increase with Ma, the mean values are consistent with a Ma-independent value at least until Ma ≈ 0.8. Thus the effects of compressibility have most likely not had a significant contribution to the results regarding the Λ effect in the previous numerical studies.
In the mean-field-theoretical treatment the Λ effect has distinct contributions from the anisotropy of turbulence and from background density stratification. The former has been modeled by an anisotropic forcing in homogeneous and fully periodic setups (e.g. Käpylä 2019) while the latter requires a mean density gradient and inevitably leads to inhomogeneity. The latter setup was studied with a set of strongly stratified simulations where turbulence was driven by isotropic forcing. Thus the anisotropy of the turbulence was induced by the density stratification. The current results indicate that the anisotropy is weak in cases where the forc-ing scale is smaller or comparable with the density scale height. The vertical velocities are suppressed (enhanced) over the horizontal components in the deep (near-surface) parts of the simulations. The vertical Reynolds stress and hence the Λ effect are, however, positive everywhere.
These results suggest an opposite sign of the vertical Λ effect due to the density gradient in comparison to contribution from the vertically dominated homogeneous turbulence at a comparable Coriolis number. Further studies involving more realistic flows (e.g. convection) are needed to study the relevance of these findings for solar and stellar differential rotation.