A&A 488, 9-23 (2008)
P. J. Käpylä1,2 - A. Brandenburg2
1 - Observatory, PO Box 14, 00014 University of Helsinki, Finland
2 - NORDITA, AlbaNova University Center, Roslagstullsbacken 23, 10691 Stockholm, Sweden
Received 19 November 2007 / Accepted 22 May 2008
Aims. We determine the components of the -effect tensor that quantifies the contributions to the turbulent momentum transport even for uniform rotation.
Methods. Three-dimensional numerical simulations are used to study turbulent transport in triply periodic cubes under the influence of rotation and anisotropic forcing. Comparison is made with analytical results obtained via the so-called minimal tau-approximation.
Results. In the case where the turbulence intensity in the vertical direction dominates, the vertical stress is always negative. This situation is expected to occur in stellar convection zones. The horizontal component of the stress is weaker and exhibits a maximum at latitude - regardless of how rapid the rotation is. The minimal tau-approximation captures many of the qualitative features of the numerical results, provided the relaxation time tau is close to the turnover time, i.e. the Strouhal number is of order unity.
Key words: hydrodynamics - turbulence - Sun: rotation - stars: rotation
Differential rotation plays a crucial role in dynamo processes
that sustain large-scale magnetic activity in stars like the
Sun (e.g. Moffatt 1978; Krause & Rädler
1980). The internal rotation of the Sun is familiar
from helioseismology (e.g. Thompson et al. 2003), but
the processes sustaining the observed rotation profile
are not understood well. The angular momentum balance of a rotating star is determined by the
The meridional flow can also be directly affected by the Reynolds stresses (e.g. Rüdiger 1989), but it is more strongly determined by the baroclinic term that arises if the isocontours of density and pressure do not coincide. Such a configuration can appear because of latitude-dependent turbulent heat fluxes that arise naturally in rotating convection (e.g. Rüdiger 1982; Pulkkinen et al. 1993; Käpylä et al. 2004; Rüdiger et al. 2005a) or from a subadiabatic tachocline (Rempel 2005) which is likely to occur in the Sun (Rempel 2004; Käpylä et al. 2006b). The flows due to thermodynamic effects are likely to be needed to avoid the Taylor-Proudman balance in the Sun (e.g. Durney 1989; Brandenburg et al. 1992; Kitchatinov & Rüdiger 1995; Rempel 2005). The overall importance of the meridional flow in the angular momentum balance of the Sun is, however, still unclear since no definite observational information about it is available below roughly 20 Mm depth (e.g. Zhao & Kosovichev 2004).
Although not much more is known about the Reynolds stresses from
observations, already this limited knowledge can be used to gain
insight into the theory of turbulent transport. Solar surface
observations indicate that there is an equatorward flux of angular
momentum, as described by the Reynolds stress component
103 m2 s-2 in the latitude range where sunspots are
observable (e.g. Ward 1965; Pulkkinen & Tuominen
1998; Stix 2002). In mean-field theory the
simplest approximation that can be made concerning the Reynolds
stresses is to assume them proportional to the gradient of mean velocity
(the Boussinesq ansatz):
Much effort has been devoted to determining Reynolds stresses
from convection simulations (Hathaway & Somerville
1983; Pulkkinen et al. 1993;
Rieutord et al. 1994; Brummell et al. 1998; Chan 2001; Käpylä et al. 2004; Rüdiger et al. 2005b;
Hupfer et al. 2005, 2006; Giesecke
2007). These studies have confirmed the existence of
the -effect and revealed some surprising results that are at
odds with theoretical considerations (e.g. Kitchatinov & Rüdiger
1993, 2005) derived under the second-order
correlation approximation (SOCA). The discrepancies are most prominent
in the rapid rotation regime,
Often the Reynolds stress realized in the simulation is taken to solely represent the -effect. This approach seems like a reasonable starting point but in an inhomogeneous system large scale mean flows are generated when rotation becomes important. These flows affect the Reynolds stresses via the turbulent viscosity. Furthermore, in the presence of stratification, heat fluxes can also significantly affect the stresses (Kleeorin & Rogachevskii 2006). In the present study we simplify the situation as much as possible in order to disentangle the effect of the turbulent velocity field from other effects. Thus we neglect stratification and heat fluxes by adopting a periodic isothermal setup. Turbulence is driven by external forcing, which provides clear scale separation between the turbulent eddies and the system size, which is typically not achieved in convection simulations. Further insight is sought from comparison of simple analytical closure models with numerical data.
Preliminary results on the -effect are presented in Käpylä & Brandenburg (2007). In the present paper, numerical datasets covering a significantly larger part of the parameter space are analyzed, and a more thorough study of the validity and results of the minimal tau-approximation are presented.
The remainder of the paper is organized as follows. Section 2 summarizes the model and the methods of the study, and in Sects. 3 and 4 the results and the conclusions are given.
We model compressible hydrodynamic turbulence in a triply periodic
cube of size .
The gas obeys an isothermal equation of state
characterized by a constant speed of sound, .
continuity and Navier-Stokes equations read
To make the resulting turbulence anisotropic the forcing
amplitude depends on direction
At this point it is important to emphasize that our forcing function is designed to capture the effects that lead to a finite -effect. The implementation of anisotropy should therefore be a simple one. In stars, anisotropy is produced by stratification and convection. Our goal is clearly not to simulate properties of convection other than its tendency to make the turbulence anisotropic.
The numerical computations were made with the P ENCIL-C ODE, which uses sixth-order accurate finite differences in space, and a third-order accurate time-stepping scheme (Brandenburg & Dobler 2002; Brandenburg 2003). Resolutions up to 2563 grid points were used in the simulations.
In the following we use non-dimensional variables by setting
The simulated domain is thought to represent a small rectangular
portion of a spherical body of gas. We choose (x,y,z) to correspond
of spherical coordinates. With this choice the
rotation vector can be written as
Since the turbulence is homogeneous, volume averages are employed and denoted by overbars.
An additional time average over the statistically saturated state of
the calculation is also taken. We define the Coriolis number as
Errors are estimated by dividing the time series into three equally long parts and computing mean values for each part individually. The largest departure from the mean value computed for the whole time series is taken to represent the error.
|Figure 1: Ux at the periphery of the simulation domain from a slowly rotating run with , , and , resolution 2563.|
|Open with DEXTER|
To have some understanding of the numerical results, we compare with the simplistic tau-approximation (hereafter MTA) in real space (e.g. Blackman & Field 2002; Brandenburg et al. 2004). Unlike the usual first-order smoothing approximation where nonlinearities in the fluctuations are neglected, they are retained in an approximate manner in MTA.
In order to develop a theory for the Reynolds stress,
we derive an equation for its time derivative,
The rotational influence in the MTA-model is measured by the Coriolis
On account of previous results on similar systems of forced turbulence (e.g. Brandenburg & Subramanian 2005) is used as our reference model. To reproduce the numerical results, we introduce an empirical isotropization term in the MTA-model; see Sect. 3.1 and in particular Eqs. (29) and (30). The use of this term is regulated by another free parameter , which can obtain the values zero or unity.
The MTA-model results were obtained by advancing the time-dependent Eqs. (22) with the parameterizations (24) and (25) until a stationary solution was reached.
Consider turbulence where the intensity in the vertical direction is stronger than the intensities in the horizontal directions. This situation is encountered if the turbulence is caused by the convective instability. First we consider the case with the maximum anisotropy that was achieved with the present model. Table B.1 summarizes the results for seven values of the Coriolis numbers ranging from 0.06 to roughly 5.4. Calculations at seven latitudes from the pole ( ) to the equator ( ) with equidistant intervals of were made with each Coriolis number, in all runs. In this set of runs, f0 = 10-6 and f1 = 0.2 were used in Eq. (13).
|Figure 2: Volume-averaged Reynolds stress components Qxx ( top), Qyy ( middle), and Qzz ( bottom), normalized by the square of the rms-velocity, as functions of latitude and rotation from the turbulence simulations listed in Table B.1. Coriolis number, as defined in Eq. (17), varies as indicated in the legend in the middle panel.|
|Open with DEXTER|
|Figure 3: Diagonal Reynolds stress components from the MTA-model. Here the Coriolis number is defined via Eq. (27) with and .|
|Open with DEXTER|
Figure 2 shows the diagonal components of the Reynolds tensor as functions of latitude and Coriolis number from the numerical turbulence simulations. As rotation is increased, the magnitudes of the horizontal components Qxx and Qyy increase monotonically while Qzz decreases, which illustrates an isotropizing effect of rotation on the turbulence. In the following, we refer to this effect as the rotational isotropization of turbulence. This effect is a purely empirical and refers to the observation that increased rotation leads to stronger mixing which washes out anisotropies that are caused by other effects. Of course, rotation itself can cause the turbulence to become anisotropic, but this seems to play a role only at much higher rotation rates. This is indeed seen in the most rapidly rotating case where the behavior is more complex. The MTA-model, Fig. 3, on the other hand, fails to reproduce isotropization of turbulence when rotation is included. The behavior is most obvious at the pole where rotation contributes no net effect to linear order, see Appendix A, Eqs. (A.1) to (A.6). The same is true for the Qxx component at the equator.
Some persistent trends arise as a function of latitude in the numerical simulations; for instance, Qxx peaks at the pole and decreases toward the equator except again for the fastest rotation. An approximately opposite trend is seen for Qyy for slow and rapid rotation, whereas in the intermediate range a minimum appears at mid-latitudes. For slow rotation, Qzz behaves approximately in the opposite way to Qyy, having a maximum at the pole and a minimum at the equator, whereas for rapid rotation the trend is reversed. There is also a persistent maximum at mid-latitudes. The MTA-model, however, fails to reproduce most of the characteristics of the latitude distribution of the diagonal stresses. This indicates that nonlinear terms, which are manifestly not described adequately well by the MTA-assumption, are the deciding factor in determining the behavior of the diagonal stresses.
To capture the rotational isotropization with the MTA-model at least qualitatively, we experimented by adding a term
Figure 4 shows the results for the diagonal components with . Now the magnitudes of the turbulence intensities are more in line with the full numerical simulations, although the latitude distribution is still manifestly wrong. Although the off-diagonal components are much better represented by the linear terms appearing in the equation of the Reynolds stress (see Sect. 3.3), the rotational isotropization term helps for reducing their magnitudes closer to the levels seen in the direct simulations also in that case.
The functional form of in Eq. (29) indicates that for rapid rotation. This behaviour cannot be justified based on the present numerical data (see Fig. 2).
Turbulence anisotropies can be characterized by the quantities
Table B.1 shows that for slow rotation increases monotonically from the pole to the equator. For , exhibits a negative minimum at mid-latitudes and reaches a positive maximum at the equator. The maximum at the equator can be explained as resulting from the surviving (y-)component of the Coriolis force at low latitudes. The relation between the horizontal -effect and the anisotropy parameter is at best poorly confirmed by the numerical results. It must also be noted that is smaller by at least an order of magnitude in comparison to , which also enters the equation for Qxy but with a higher order in (see Eq. (A.7)).
The vertical anisotropy, , on the other hand, retains its sign in all models. For slow rotation the absolute value of decreases monotonically from the pole toward the equator. The latitude distribution for rapid rotation is approximately opposite, although an additional minimum ( ) or a maximum ( ) can occur at mid-latitudes. Here the correspondence between and holds at least for the sign for all calculations, bar the few cases in the intermediate rotation regime where Qyz is positive (see below).
In stellar convection zones, the off-diagonal Reynolds stresses contribute to the angular momentum transport and work to generate (-effect) or to smooth out (turbulent viscosity) differential rotation. In the present case only the former effect is in operation. Figure 5 summarizes the results for the runs listed in Table B.1.
The Qxy component of the stress corresponds to in spherical coordinates and is responsible for latitudinal angular momentum transport. In the simulations this component is always positive and the latitude distribution peaks at latitude (see the uppermost panel of Fig. 5). The sign is in accordance with solar observations (Ward 1965; Pulkkinen & Tuominen 1998) and analytical turbulence models (Kitchatinov & Rüdiger 1993, 2005). The latitude distribution in the rapid rotation regime is significantly different from convection simulations where Qxy is sharply concentrated near the equator (Chan 2001; Käpylä et al. 2004; Hupfer et al. 2005; Rüdiger et al. 2005a). The reason for this difference is still unclear but is evidently related to the physics that have been omitted in the present study. The MTA-model captures the latitude dependence rather well, but the magnitude of the stress is clearly too great (see Fig. 6). If the rotational isotropization term, Eq. (29), is taken into account (see the uppermost panel of Fig. 7), the agreement is also better for the magnitude. Keeping the forcing and rotation rate fixed, the best agreement with the 3D models is found if is used, see Fig. 8.
|Figure 5: Same as Fig. 2 but for the off-diagonal stress components Qxy ( top), Qxz ( middle), and Qyz ( bottom). Linestyles as in Fig. 2.|
|Open with DEXTER|
|Figure 6: Same as Fig. 3 but for the off-diagonal stress components Qxy ( top), Qxz ( middle), and Qyz ( bottom). , . Linestyles as in Fig. 3.|
|Open with DEXTER|
While the analogue of the Qxz component does not play a direct role in the angular momentum balance in stars, it can still contribute via generating meridional flows (e.g. Rüdiger 1989). In the numerical simulations we find that, for all cases except the most rapidly rotating one, Qxz is negative and peaks at (middle panel of Fig. 5). For , however, the sign changes near the equator, with positive values toward the equator and negative ones toward the pole. This is at odds with the analytical result of Rüdiger et al. (2005b), but does agree with the results of convection simulations (Pulkkinen et al. 1993; Käpylä et al. 2004).
The MTA-model gives qualitatively similar results, although the sign change occurs at significantly slower rotation; see the middle panel of Fig. 6. Rotational isotropization helps to correct the magnitude, but not the earlier occurrence of the sign change (Figs. 7 and 8).
Since the turbulent intensity of the vertical (z-)motions is greater than the horizontal ones, the expectation is that the stress component Qyz is negative, i.e. that (Biermann 1951). This is indeed seen in the simulations quite consistently (lowermost panel of Fig. 5), although at intermediate rotation low positive values can occur at high latitudes. The highest values of Qyz occur for as opposed to for Qxy. Rotational quenching of Qyz seems to be stronger and occur for lower than for Qxy(see also Fig. 10). A similar trend was seen in convection simulations by Käpylä et al. (2004). The consistently positive values of Qyz for rapid rotation seen in convection simulations (Käpylä et al. 2004; Chan 2007, private communication) do not occur in the present calculations.
The value of Qyz always reaches a maximum at the equator. This contradicts with analytical results derived under first-order smoothing (Kitchatinov & Rüdiger 1993, 2005) and the MTA-model that predicts a maximum around for intermediate and rapid rotation. The rotational isotropization term is again needed to reduce the magnitude of the stress. It seems that to reproduce Qyz correctly, one would need to apply somewhat stronger isotropization than what is presently used (see Figs. 7 and 8).
One conclusion that can be drawn from the MTA-model results is that, although the diagonal Reynolds stresses are quite poorly reproduced in comparison to the 3D simulations, the off-diagonals have most of the qualitative features correct. This seems to imply that, to model the off-diagonals, the exact parameterization of the nonlinearities in the equation of the Reynolds stress is not crucial. The empirical rotational isotropization term helps to capture some of the missing features for the diagonal components and reduces the magnitudes of the off-diagonals to the level that is also seen in the simulations.
The homogeneous setup used so far prevents any mean flows from being generated. This is good for the purpose of testing the sole effect of turbulent velocity field on the Reynolds stresses, but can be argued to be unphysical because astronomical objects where turbulence is important; e.g., the solar convection zone, have boundaries and cannot be considered homogeneous.
To test how much the assumption of homogeneity affects the results,
we made a set of simulations with a setup where the z-boundaries are
impenetrable. For the horizontal velocity components, stress-free boundary
conditions are used, i.e.
|Ux,z = Uy,z = Uz = 0.||(33)|
Mean flows are generated in the runs that, however, are small in comparison to the
rms-velocity of the turbulence. The exception to this trend is the
where a large mean shear flow develops, for which
The results for the off-diagonal stresses are compared in Fig. 9. For slow rotation, , the differences are minute at all other latitudes except at the equator. The tendency for Qyz to have a minimum at the equator is reminiscent of results from convection simulations (e.g. Chan 2001; Käpylä et al. 2004; Rüdiger et al. 2005b). For more rapid rotation, the trends for different components seem to diverge: Qxy is somewhat reduced whereas Qxz seems to increase somewhat and there is hardly any change for Qyz apart from the equatorial case. For the stresses at the equator are not saturated because the large scale flow is not fully developed.
To connect to earlier studies, the 3D simulation data and MTA-model results are compared to the analytical SOCA results of Kitchatinov & Rüdiger (2005, hereafter KR05; see also Kitchatinov & Rüdiger 1993). Since analytical results are only available for the components relevant to the -effect, only the Qxy and Qyz components are thus considered. Furthermore, we restrict the comparison to a subset of the models with .
|Figure 9: Same as Fig. 5 but for homogeneous ( left panel) and inhomogeneous ( right) simulations.|
|Open with DEXTER|
|Figure 10: Comparison of the numerical simulations (solid lines), two MTA-models with , (dotted), and , (dashed), with the SOCA results of KR05 (dash-dotted). Left ( right) hand panels show ( ) for Coriolis numbers 0.15 (uppermost panels), 1.5 ( middle), and 5.4 ( lower panels).|
|Open with DEXTER|
The conventional way of writing the -effect is (e.g. Rüdiger 1989)
The normalized stresses are known from the simulations and the MTA-model, whereas H and V can be computed analytically for the turbulence model of KR05 (see Appendix B). The results are shown in Fig. 10 for three Coriolis numbers. Our Coriolis number is smaller than that of KR05 by a factor of .
For the horizontal stress, the SOCA results and the rotationally quenched MTA-model seem to fare similarly well. The former underestimates the magnitude for small and overestimates it for large , whereas for the latter the trend is exactly the opposite. If rotational isotropization is not taken into account, the agreement is poor for all Coriolis numbers considered here.
For Qyz the standard MTA-model is almost spot on for , but overestimates the magnitudes by at least a factor of two for the other cases. As discussed in the previous section, the latitude distribution shows a mid-latitude maximum that is not present in the numerical simulations. When rotational isotropization is taken into account, at least the magnitude can be reconciled with the numerical results. The SOCA result does not fare very well in this case, predicting, in general, values that are too low and an incompatible latitude distribution with zero stress at the equator.
|Figure 11: Volume-averaged Reynolds stress components Qxy ( top), Qxz ( middle), and Qyz ( bottom) as functions of latitude and vertical turbulence anisotropy , see Eq. (32). Line styles as indicated by the legend in the lower-most panel. Coriolis number in each case is roughly 0.3.|
|Open with DEXTER|
Table B.2 and Fig. 11 show the results for four sets of calculations in which the Coriolis number is kept approximately constant at 0.3 whereas the turbulence anisotropy is varied between -0.07 and -0.47. This is done by choosing suitable values for f0 and f1 so that the rms-velocity stays approximately constant.
From the figure it is seen that the ratio of the stress to the amount of anisotropy decreases monotonically as a function of . Although part of the difference can be explained by the somewhat smaller Co (0.31 as opposed to 0.34 in the other cases) in the calculation with the largest , the trend still persists.
This trend can be understood as follows: from the approximate relation we obtain . The decreasing trend of seen in the results suggests that the correlation time changes when the turbulence anisotropy is varied, i.e. decreases when is increased. This is plausible since to change different values of the forcing amplitudes, f0 and f1 (see Eq. (13)), need to be used resulting in differences in the turbulence.
|Figure 12: Off-diagonal Reynolds stresses Qxy ( top), Qxz ( middle), and Qyz ( bottom) as functions of latitude and Reynolds number, see the legend in the lower-most panel.|
|Open with DEXTER|
Figure 12 shows the off-diagonal stresses as functions of latitude and Reynolds number for a constant (see also Table B.3). There is no a priori reason to expect that the stresses should depend on the Reynolds number if and if the turbulence anisotropy is kept constant. This is essentially what is borne out of the simulations, although there seems to be a weak decreasing trend as a function of Re for for Qxy and Qxz, although the results are within error bars. For Qyz, however, the largest Reynolds number case shows a distinct drop in comparison to the less turbulent cases.
The decrease in the stresses as a function of the Reynolds number is likely to have the same origin as the decrease seen when the turbulence anisotropy is increased (see the previous section). In order to obtain the same in the simulations different values, f0 and f1 are needed for different values of , the kinematic viscosity. More precisely, the less the viscosity, the more difficult it is to obtain large anisotropy. Thus, for smaller , higher ratio f1/f0 is required to achieve a given . Thus, if the interpretation that the higher the ratio f1/f0, the smaller becomes is correct, the decreasing trend seen as a function of might be an artefact due to the small differences in the forcing between the different runs.
The MTA-model reproduces the simulation results for the off-diagonal stresses reasonably well when is used. Now we turn to the numerical simulations to determine independently. Although we have not been able to provide any direct support of the basic MTA-assumption, , the MTA-model is still able to reproduce many of the features of the numerical simulations adequately. This implies that a value for could be extracted from the MTA-relations.
We consider two methods to determine the Strouhal number from the simulations: (i) MTA-relations derived for the Reynolds stresses; and (ii) similar relations for passive scalar transport under the influence of rotation.
Using the minimal tau-approximation and assuming a stationary state where , the off-diagonal Reynolds stresses can be derived from Eq. (8), yielding
|Figure 13: Strouhal numbers computed from Eq. (47) at colatitude as functions of the Coriolis number for the runs listed in Table B.1. A power law proportional to shown for reference.|
|Open with DEXTER|
As an independent check of the dependence of the Strouhal number on rotation, we expand the passive scalar transport case, which was studied by Brandenburg et al. (2004) to cases where rotation is included. The numerical model is the same as in the runs presented so far, except that isotropic forcing is used with f0=0.01-0.03 and f1 = 0 in Eq. (13).
The turbulent passive scalar flux is denoted by
where c is the
fluctuation of passive scalar density, i.e. the passive scalar
concentration per unit volume. Following the MTA-approach, we
solve first for the time derivative
In the passive scalar cases we consider isotropically forced turbulence for which even when rotation is included. Equations (53)-(55) yield
In the passive scalar case, the second and third order terms are
indeed correlated (Brandenburg et al. 2004)
according to the basic MTA-assumption, and a Strouhal number can be thus computed using
|Figure 14: Uppermost panel: St from triple correlation , middle panel: , corresponding to the tau in Eq. (58), lowermost panel: Strouhal numbers from Eq. (59) as functions of the Coriolis number for a constant Reynolds number of roughly nine.|
|Open with DEXTER|
For slow rotation, the Strouhal number is consistently between one and three when computed from the Reynolds stresses (Fig. 13). Similar values are obtained from the passive scalar transport (Fig. 14) when the Reynolds number is ten or larger. The Strouhal number computed from the triple correlations, , is more strongly dependent on the Reynolds number, but it seems to converge slowly towards a constant value near unity. These results are in line with the values required in the MTA-model and earlier studies in different contexts employing similar turbulence calculations (Brandenburg et al. 2004; Brandenburg & Subramanian 2005, 2007).
However, when the Coriolis number approaches or exceeds unity, the Strouhal numbers computed from the equations of the Reynolds stresses, Eqs. (42)-(44), decrease rapidly so that for it has dropped at least by an order of magnitude (see Fig. 13). Similar results are obtained for the passive scalar transport under the influence of rotation, see Eqs. (56)-(58) and the lower panels of Fig. 14. The trend is clearer for and , whereas for the decreasing trend is seen only for rapid rotation, i.e. when . For slow rotation, however, is almost constant and increases somewhat when the Coriolis number approaches unity. These results seem to confirm the trend seen earlier in convection simulations (Käpylä et al. 2005, 2006a).
The Strouhal number from the triple correlations follows a trend similar to , with increasing values up to after which there is a rapid decrease, see the uppermost panel of Fig. 14. For low Reynolds number St can become negative in the range , hence the gaps in the corresponding data in Fig. 14.
Turbulent momentum fluxes, which are described by the Reynolds stresses, were determined from numerical simulations of homogeneous rotating anisotropic turbulence. Since no large-scale shear is present, the generated Reynolds stresses correspond to contributions that are already present for uniform rotation. The resulting term is known as the -effect (Krause & Rüdiger 1974). The component responsible for the horizontal transport, Qxy, is positive and peaks around latitude regardless of the Coriolis number. The vertical component is predominantly negative and it always peaks at the equator.
Although the numerical results for the -effect broadly agree with analytical SOCA calculations (Kitchatinov & Rüdiger 1993, 2005), the MTA-model seems to reproduce certain features of the numerical results somewhat more closely. The present numerical results do not show the enigmatic results, such as the extreme latitude distribution of Qxy or a positive Qyz for rapid rotation, which have been reported from convection simulations (e.g. Chan 2001; Käpylä et al. 2004). The difference lies most likely in our neglecting stratification and heat fluxes. The exact manner in which they affect the Reynolds stresses is not within the scope of the present paper, but should be investigated more closely in the future.
By applying the minimal tau approximation closure relation to the Reynolds stress equation, qualitatively similar results are obtained, but the magnitude of the stresses is in general too large. The vertical flux in the MTA-model, however, has a maximum at mid-latitudes for intermediate and rapid rotation. Adding an empirical rotational isotropization term (motivated in Sect. 3.1) also brings the magnitude in line with the 3D simulations. Although adding this term with this particular form has no rigorous theoretical basis, we can see that phenomenological effects of isotropization of turbulence due to rotation are indeed seen in the simulations and that the term is thus justified.
Another drawback of the MTA-model is that the diagonal components of the Reynolds tensor are rather badly reproduced since the nonlinear effects of rotation manifest in the numerical simulations are not explicitly taken into account. The empirically added rotational isotropization term augments the magnitudes, but not the latitude distribution. Furthermore, no direct evidence of the validity of the MTA-assumption was found in the numerical simulations. Contrasting the behavior of the diagonal components to the fairly good correspondence between the numerical simulations and the MTA-model for the off-diagonal components leads us to conclude that, whereas the behavior of the diagonal components is dominated by the inadequately modeled nonlinear effects, the off-diagonals are fairly well presented by the linear terms.
A Strouhal number of order unity in the MTA-model gives best fits to the numerical results. Fitting the numerical results to expressions derived under the MTA, similar values of are found for slow rotation. For Coriolis numbers of order unity or larger, however, the Strouhal number obtained in this manner decreases rapidly. In the passive scalar case, the situation is somewhat more complex, although a similar decreasing trend of the Strouhal number is recovered for rapid rotation, see Fig. 14. These results are in accordance with earlier results from convection simulations (Käpylä et al. 2005, 2006a) using Reynolds stresses or correlation analysis of the velocity field.
A related aspect in turbulent transport that requires closer study is the turbulent viscosity (see preliminary results in Käpylä & Brandenburg 2007) and the possibility of a -effect due to the anisotropy induced by a large-scale shear flow (Leprovost & Kim 2007). These matters will be considered in more detail in a future publication.
The computations were performed on the facilities hosted by the Center of Scientific Computing in Espoo, Finland, who are financed by the Finnish ministry of education. P.J.K. acknowledges the financial support from Helsingin Sanomat foundation and the Academy of Finland grant No. 121431. P.J.K. acknowledges the hospitality of Nordita during the program ``Turbulence and Dynamos'' during which this work was finalized. The anonymous referee is acknowledged for the critical reading and helpful comments on the manuscript.
and using the MTA-closure
and parameterizing the contributions of the forcing as
Eq. (22) gives six equations for the six unknowns of
the symmetric tensor Qij, i.e.
In the present study the forcing is such that Qyy(0) - Qxx(0) = 0, so the equations reduce to
The fluxes of angular momentum have commonly been parameterized by Eqs. (35) and (36), and the normalized fluxes by Eqs. (37) and (38). Kitchatinov & Rüdiger (2005) computed these coefficients using SOCA
The functions Ii and Ji are given by
Table B.1: Summary of the turbulence anisotropies and normalized Reynolds stresses, , for the calculations. , , and grid resolution 643 was used in all runs.
Table B.2: Summary of the turbulence anisotropies and normalized Reynolds stresses, , for a set of runs with varying turbulence anisotropy. , , , and grid resolution 643 was used in all runs.
Table B.3: Turbulence anisotropies and Reynolds stresses for varying Reynolds numbers. and were used in all runs. .