Effects of small-scale dynamo and compressibility on the Λ effect

Powered by TCPDF (www.tcpdf.org) This material is protected by copyright and other intellectual property rights, and duplication or sale of all or part of any of the repository collections is not permitted, except that material may be duplicated by you for your research use or educational purposes in electronic or print form. You must obtain permission for any other use. Electronic or print copies may not be offered, whether for sale or otherwise to anyone who is not an authorised user. Käpylä, Petri J.


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 1980;Rüdiger 1989;Rüdiger et al. 2013). While turbulence is often associated with only enhanced diffusion, several nondiffusive 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 nondiffusive 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;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 nondiffusive HD effects include the anisotropic kinetic alpha 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.
Numerical simulations of magnetohydrodynamic (MHD) convection in spherical coordinates have reached sufficient spatial resolution that allows 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 , the effect of large-scale 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 large-scale 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ä et al. 2004;Käpylä & Brandenburg 2008;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 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 , 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 = 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 Equation (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 ij is given by where f 0 and f 1 are the amplitudes of the isotropic and anisotropic parts, respectively, ij is the Kronecker delta, and Θ k is the angle betweenê 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 = 2 s ∕ 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 where = L d k 1 /k f = 2 /k f . The Mach number is given by The magnetic field strength is given in terms of the equipartition value Finally, the parameter characterizes vertical anisotropy of turbulence. The density stratification is quantified by the ratio of the densities at the top and bottom of the domain, Δ = bot ∕ 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 . The Reynolds stress is given by = where the overline denotes horizontal averaging and where u = U − U is the fluctuating velocity. In the homogeneous cases Sets SSD and MA 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 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 .
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.

T A B L E 2 Summary of the MA simulations
Set  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. Figure 2 shows the power spectra of the velocity from Run MA7. The power peaks near the forcing wavenumber for the total and vertical velocity. The peaks at the overtones of the forcing wavenumber arise because in the anisotropic 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 components. A representative result for A V (k) is shown in the inset of 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. While A V (k) is mostly negative for large scales, it gradually increases and peaks around 0.5 for̃≈ 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 Note: Grid resolution 288 3 , c s = 3, f 1 = 0, H k 1 = 9/10, Δ = 10 3 , and Re = 15. The values of Ω ⋆ and A V indicate the extrema from the range |z| < 9 /10. The starred runs were repeated at the same colatitudes as the SSD sets. F I G U R E 2 Power spectra of the total velocity (black), and its x (blue), y (yellow), and z (red) components. The inset show the spectral anisotropy parameter A V (k) according to Equation (27) for Run MA7 stratification are studied in Set STR. By virtue of the isothermal equation of state this setup is characterized by a constant density scale height = 2 s ∕ = 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).

Ruñf
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̃f = f ∕ 1 , is varied between 1.5 and 10, see Figure 3 and Table 3 where a F I G U R E 3 Anisotropy parameter A V (z) from the runs in Set STR with varying̃f 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 separa-tioñ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. Owing 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 level̃r ms = rms ∕ eq , where B rms is the rms-value 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 for̃r ms ≈ 0.5 F I G U R E 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 are typically roughly half of their HD 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ä & Brandenburg 2008) were typically relatively low ((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 order of 10%. 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. Pipin & Kosovichev 2018). Here this scenario is tested with a density-stratified setup with isotropic forcing in the Set STR, see Table 3.

F I G U R E 6
Vertical Reynolds stress Q yz normalized by the squared rms-velocity 2 rms (z) in density-stratified runs with varying /H The simulations considered here have Ω ⋆ ≈ 0.8 which is close to the Coriolis number where the Reynolds stresses obtain a maximum in . 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 and 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 near-surface regions where A V < 0. The results indicate that the scale separation ratio plays 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 wheñf increases. The values of Q yz are of the order of a few per cent of 2 rms for /H = 4.5 whereas for /H = 0.7 the effect is no longer statistically significant at the equator. Figure 7a shows the normalized vertical stress̃y z from seven latitudes for the STR4 runs. The maximal values are generally on the order of 1-2% 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 7b. 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 can 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  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 small-scale 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., . 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 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 forcing 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.

ACKNOWLEDGMENTS
The computations were performed on the facilities hosted by CSC -IT Center for Science Ltd. in Espoo, Finland, who are administered by the Finnish Ministry of Education. This work was supported in part by the Deutsche Forschungsgemeinschaft Heisenberg programme (grant No. KA 4825/1-1) and by the Academy of Finland ReSoLVE Centre of Excellence (grant No. 307411).