A&A 425, 725-728 (2004)
DOI: 10.1051/0004-6361:20040531

Vorticity and helicity of the solar supergranulation flow-field[*]

P. Egorov - G. Rüdiger - U. Ziegler

Astrophysikalisches Institut Potsdam, An der Sternwarte 16, 14482 Potsdam, Germany

Received 26 March 2004 / Accepted 9 June 2004

Box simulations of rotating convection are presented with the aim to study the divergence-vorticity correlation ${\cal C}$ computed from the horizontal velocity components of the turbulent flow field. Our investigation is motivated by new observations of the solar velocity field in the upper convection zone using the techniques of so-called time-distance helioseismology (Duvall et al. 1993). From such observations ${\cal C}$ can be derived directly and because of the close relation to the helicity it is an important quantity for solar dynamo theory. The divergence-vorticity correlation obtained from the simulations with Taylor number 103 is found to be in agreement with the observations of supergranulation. In particular, ${\cal C}$ vanishes at the equator and is negative (positive) in the northern (southern) hemisphere. For larger Taylor numbers the amplitude of ${\cal C}$ increases, which is consistent with the increase in strength of the Coriolis force acting on the turbulent flow. We also demonstrate that the correlation ${\cal C}$ and the kinetic helicity of the flow field have very similar spatial profiles. The observed values of ${\cal C}$, therefore, yield a clear indication of the existence of negative (positive) kinetic helicity in the surface layers of the northern (southern) hemisphere.

Key words: turbulence - Sun: granulation - Sun: helioseismology

1 Introduction

There is no doubt that the magnetic activity of the cool stars is driven by the turbulence in the outer convection zones which is subject to (density) stratification and the global rotation. This finding is a common result of analytical considerations and numerical simulations. The existence of kinetic helicity and current helicity which are both antisymmetric with respect to the equator is central in those computations.

It is therefore natural to look for observational confirmations of this concept. A big step in this direction was the statistics of the current helicity

 \begin{displaymath}{\cal H}_{\rm curr}= \langle \mbox{\boldmath$J$ }'\cdot \mbox{\boldmath$B$ }'\rangle
\end{displaymath} (1)

by Seehafer (1990) which has the same kind of equatorial (anti-) symmetry as the dynamo-$\alpha$. For homogeneous global magnetic fields the latter is related to the turbulent EMF according to $
{\mbox{\boldmath$\cal E$ }} = \langle \mbox{\boldmath$u$ }'\times \mbox{\boldmath$B$ }'\rangle \simeq \alpha
\circ \bar{\mbox{\boldmath$B$ }}
$. The antiphase relation

 \begin{displaymath}\alpha_{\phi\phi} {\cal H}_{\rm curr} < 0
\end{displaymath} (2)

by Keinigs (1983) between the $\alpha$-effect and the observed current helicity has suggested the existence of an $\alpha$-effect which is positive (negative) in the northern (southern) hemisphere.

The kinetic helicity

 \begin{displaymath}{\cal H}_{\rm kin}= \langle \mbox{\boldmath$u$ }' \cdot \nabla \times \mbox{\boldmath$u$ }'\rangle
\end{displaymath} (3)

which much earlier has been discussed as the physical background of the $\alpha$-effect seemed to withstand direct observation. Nevertheless, the rotational influence on the turbulence field should be observable on the solar surface. However, the effect is expected to be small and it is rather complicated to prove as data only exist for the horizontal components u and v of the flow.

In particular, the vertical component of the vorticity,

 \begin{displaymath}{\rm CURL} = \frac{\partial v}{\partial x} - \frac{\partial u}{\partial y},
\end{displaymath} (4)

and the (two-dimensional) divergence,

 \begin{displaymath}{\rm DIV} = \frac{\partial u}{\partial x} + \frac{\partial v}{\partial y},
\end{displaymath} (5)

of the horizontal velocity seem to be suitable quantities to study the rotation-influenced surface flow field. At the top of the convection zone a cellular network is formed by the divergence of upflows and convergence of downflows. In the compressible medium the rising matter in the northern hemisphere expands and spins up under the basic rotation in order to conserve its angular momentum. In a righthand coordinate system with z pointing outwards these regions have ${\rm DIV} > 0$ and rotate clockwise so that ${\rm CURL} < 0$. The downflows converge ( ${\rm DIV} < 0$), compress and rotate in opposite direction ( ${\rm CURL} > 0$) under the influence of the Coriolis force. Therefore, lefthand spiral motions dominate in the northern hemisphere. The divergence-vorticity correlation

 \begin{displaymath}{\cal C} = \langle{\rm DIV}\!\cdot\!{\rm CURL}\rangle
\end{displaymath} (6)

is thus expected to be negative (positive) in the northern (southern) hemisphere. A quasilinear calculation of this (nontrivial) correlation for anisotropic turbulence with the $\tau$-approximation leads to the result

 \begin{displaymath}{\cal C} \simeq - \frac{{\it\Omega}}{\tau_{\rm corr}} \cos\theta
\end{displaymath} (7)

(Rüdiger et al. 1999) which is of order 10-11 s-2 for correlation times of 1 day.

The candidates to study the rotational influence on the turbulence at the solar surface are mesogranules and supergranules due to their nonmagnetic character. The average cell sizes of supergranules are in the range 15-30 Mm with predominantly horizontal velocities. Recent observations (e.g., Hathaway et al. 2002) have shown that the typical rms horizontal flow speed is $\sim$360 m/s while the rms vertical flow speed is $\sim$30 m/s.

The rotational effect is estimated by the Coriolis number

 \begin{displaymath}\Omega^*= 2\tau_{\rm corr} {\it\Omega},
\end{displaymath} (8)

where $\tau_{\rm corr}$ is a typical correlation (life-) time and ${\it\Omega}$ is the angular velocity of the basic rotation. The larger the scale of a velocity field the more sensitive is the flow to the action of the Coriolis force. For supergranulation one finds $\Omega^*\simeq 0.1$. For mesogranulation $\Omega^*$ is of the order 10-2.

We did not find observations for the magnitude of ${\cal
C}_{\rm meso}$, but only maximum values of the divergence and vorticity were given (Brandt et al. 1988; Simon et al. 1988; Simon et al. 1994). Wang et al. (1995) report indications of a negative divergence-vorticity correlation ${\cal C}$ for the northern hemisphere.

A clear detection of the influence of the Coriolis force is due to Duvall & Gizon (2000) and Gizon & Duvall (2003) using the method of time-distance helioseismology (Duvall et al. 1993). This technique provides both maps of the horizontal divergence of the flows and information about two individual horizontal components of the velocity. We shall use their results for comparison with our numerical simulations. The magnitude of the vorticity was found for the latitudinal interval of $\pm 45^{\circ}$ with zero at the equator and a value $\vert{\rm CURL}\vert\simeq 2\times10^{-6}$ s-1 at the interval boundaries (see Fig. 7b in Gizon & Duvall 2003). A vorticity of 10-6 s-1 corresponds to an angular velocity of $2.5^{\circ}\ {\rm day}^{-1}$ or a typical circular velocity of 10 m/s. Gizon & Duvall (2003) detected a "significant correlation of a few percent between the vertical vorticity and the horizontal divergence'' (their Fig. 7a). Positive (negative) divergence in the northern hemisphere is correlated with clockwise (counterclockwise) vorticity i.e. ${\cal C}<0$. ${\cal C}$ changes its sign in the southern hemisphere and vanishes at the equator. Both the sign and the latitudinal variation of ${\cal C}(\theta)$characterize the effect of the Coriolis force on the flow. The magnitude of the correlation ${\cal C}$ is obtained for the range of latitudes $[-45^{\circ}, 45^{\circ}]$ where it varies from $2\times10^{-10}$ s-2to $-2\times10^{-10}$ s-2.

Similar kinds of simulation, both the "local-box'' (e.g., Hathaway 1984; Brummell et al. 1996, 1998) and the full spherical shell (e.g., Miesch et al. 2000) have been used to study the effect of rotation on laminar/turbulent, incompressible/compressible convection. A correlation of counter-clockwise (cyclonic) flows with downdrafts and clockwise (anticyclonic) with updrafts as well as a switch in the direction of circulations in the bottom part of convection zone has been found for a wide range of Taylor number varying from $\sim$102 to $\sim$107. The second result in its turn yields a sign reverse of the kinetic helicity along the radial direction; this has also been confirmed by simulations of Brummell et al. (1998) and Miesch et al. (2000).

In this paper the influence of the Coriolis force on the turbulence is simulated for the solar surface flow pattern of supergranulation. We shall show that the results of the simulations are close to both the quasilinear theory and the observations. In the last section the relation between the DIV-CURL correlation and the helicity is discussed with the result that the observations indicate a negative (positive) kinetic helicity in the northern (southern) hemisphere.

2 Basic equations and model

In this section we give a short description of the numerical model. The 3D numerical simulations of compressible, thermal convection under the influence of rotation are performed in a small rectangular box (Cartesian grid) which can be considered as a small piece of a spherical shell. The domain is placed tangentially at different colatitudes from $\theta =0^\circ $ (at the north pole) to $180^{\circ}$(at the south pole) with a step of $30^\circ $. A sandwich configuration which contains a convectively unstable layer in between stably stratified layers with overshooting convection has been constructed. The height d of the convection zone is chosen as our length unit.

The box rotates around the polar axis from west to east. The computational domain is $(x,y,z)\in[0,8d]\times[0,8d]\times[-2d,0]$ discretized by $100\times 100\times 128$ grid points which are uniformly distributed in each coordinate direction.

The governing equations describing thermal convection in a stratified medium involving the effect of the Coriolis force are

$\displaystyle \frac{\partial\rho}{\partial t} = -\nabla\cdot(\rho\mbox{\boldmath$u$ }),$     (9)
$\displaystyle \frac{\partial\rho\mbox{\boldmath$u$ }}{\partial t} = -\nabla\cdo...
...\boldmath$g$ } - 2\rho\mbox{\boldmath${\it\Omega}$ }\times\mbox{\boldmath$u$ },$     (10)
$\displaystyle \frac{\partial e}{\partial t} =
-\nabla\cdot(e\mbox{\boldmath$u$ })-P\nabla\cdot\mbox{\boldmath$u$ }+Q_{\rm vis}+
\nabla\cdot(\kappa\nabla T).$     (11)

The notation for the physical variables is the following: $\rho$ density, P pressure, e thermal energy density, T temperature, $\mbox{\boldmath$u$ }$ velocity and $\mbox{\boldmath$g$ }=-g\mbox{\boldmath$\hat z$ }, g>0$ the gravitational field. The remaining quantities in the equations are the viscous stress tensor $\tau$ given by

 \begin{displaymath}\tau_{ij}=\nu\rho\left(\frac{\partial u_i}{\partial x_j} +
...frac{2}{3}\frac{\partial u_k}{\partial x_k}\delta_{ij}\right),
\end{displaymath} (12)

where $\nu$ denotes the kinematic viscosity coefficient, the viscous heating term

\begin{displaymath}Q_{\rm vis} = \sum_{ij}\tau_{ij}\frac{\partial u_i}{\partial x_j}
\end{displaymath} (13)

and the thermal conductivity coefficient $\kappa$.

The set (9)-(11) is closed through the ideal gas equation of state

 \begin{displaymath}P=\frac{k}{\bar{\mu}m_{\rm p}}\rho T = (\gamma - 1)e,
\end{displaymath} (14)

where k is Boltzmann's constant, $m_{\rm p}$ is the atomic mass unit, $\bar{\mu}$ is the mean molecular weight and $\gamma=C_{\rm P}/C_{\rm V}=5/3$is the ratio of specific heats $C_{\rm P}$ and $C_{\rm V}$. The chemical composition of the gas is ignored, hence $\bar{\mu}=1$.

The initial distribution of the physical quantities $T,\rho,P$ corresponds to a 3-layer polytropic stratification. All quantities are periodic in the horizontal direction. At the bottom (z=-2) and the top (z=0) of the box impermeable conditions are imposed for the vertical velocity while the horizontal velocities satisfy a stress-free boundary condition. The temperature and density are fixed at the top of the domain, $T_3=1,\ \rho_3=1$, and a constant heat flux is injected at the bottom.

The dimensionless parameters ${\rm Ra},\ {\rm Pr}$ and ${\rm Ta}$ are used to control the simulations. In our calculations ${\rm Ra}=3\times10^5$ and ${\rm Pr}=0.1$ are fixed whereas ${\rm Ta}=\{10^3,\ 10^5,\ 10^6\}$ characterize the magnitude ${\it\Omega}$ of the angular velocity vector.

Equations (9)-(11) are solved with the finite-difference, fractional-step code NIRVANA (Ziegler 1998, 1999).

3 Results and discussion

Despite the simplifications in the numerical model (local box approximation, gradients of temperature and density small compared to the real values in the solar convection zone) the results discussed in the following are quite well in agreement with the observations.

3.1 Anisotropy of the velocity field

The numerical simulations show an anisotropy between the vertical and horizontal turbulent velocities. Vertical motions dominate in the bulk of convection zone due to rising and descending flows. At the top the flow becomes mostly horizontal due to divergent and convergent motions. Figure 1 presents the computed anisotropy for ${\rm Ta}=10^6$and for different colatitudes in the northern hemisphere: in the upper part of the convection zone the intensity of the horizontal velocities $\langle u_{\phi}'^2\rangle$ dominates over the vertical ones. One can see there how the basic rotation modifies the intensity ratio with colatitude. For slower rotation rates the dependence on colatitude clearly becomes weaker. For the purpose of comparing our numerical results with the observations we consider a region near the top of the convection zone, for example at z=-0.3 ($5\%$ in depth) where the horizontal cellular network is still visible as in the observations. The ratio $\langle u_r'^2\rangle/\langle u_{\phi}'^2\rangle$ at this point varies as $\simeq $30-45% in all simulations.
\psfig{figure=0531fig1.ps,width=7cm,clip=} }
\end{figure} Figure 1: The intensity ratio $\langle u_r'^2\rangle/\langle u_{\phi}'^2\rangle$ in the box simulations for different colatitudes. The colatitudes are $\theta =0^\circ $ (poles, solid), $30^\circ $ (dashed) and $60^\circ $ (dot-dashed), ${\rm Ta}=10^6$.
Open with DEXTER

3.2 Divergence-vorticity correlation

To allow a direct comparison with the observed ${\cal C}$ a scaling factor $(D/c_{\rm s})^2$ is introduced where D is the depth of the convection zone and $c_{\rm s}=(\gamma p_3/\rho_3)^{1/2}$ is the sound speed at the top of the domain (z=0). Note that $c_{\rm s}$ is the same in all simulations. For solar parameters i.e. $D=0.3~R_{\odot}$ and $T_{\odot}=5700$ K, we find

 \begin{displaymath}\left(\frac{D}{c_{\rm s}}\right)^2=\frac{0.09~R_{\odot}^2 m_{\rm p}}{\gamma kT_{\odot}}\cdot
\end{displaymath} (15)

Figure 2 shows the computed correlation $\langle {\rm
DIV}\cdot{\rm CURL}\rangle$ versus depth for ${\rm Ta}=10^6$. At the top of the convection zone (z<-0.25) ${\cal C}$ is negative (positive) for the northern (southern) hemisphere. As expected, ${\cal C}$ vanishes at the equator and its magnitude increases toward the poles. At the bottom of the convection zone and in the lower overshoot region ${\cal C}$ is relatively small and it reverses its sign twice at $z\approx -0.75$ and $z\approx -1.15$. The results for ${\rm Ta}=10^5$ qualitatively look similar except that there is no sign reversal. The amplitude is reduced roughly by a factor of 2 compared with the ${\rm Ta}=10^6$ case.

\psfig{figure=0531fig2.ps,width=7cm,clip=} }
\end{figure} Figure 2: The correlation ${\cal C}$ in the box simulations for various latitudes after both horizontal and temporal averaging vs. depth for ${\rm Ta}=10^6$.
Open with DEXTER

Figure 3 shows horizontal and time averages of ${\cal C}$ scaled by Eq. (15) vs. latitude plotted in comparison with observational results indicated by the shaded area. In Fig. 3 triangles ( ${\rm Ta}=10^3$), asterisks ( ${\rm Ta}=10^5$) and diamonds ( ${\rm Ta}=10^6$) denote the computed values of ${\cal C}$at a depth z=-0.3 (see Fig. 2) and for different latitudes $90-\theta$.

\psfig{figure=0531fig3.ps,width=7cm,clip=} }
\end{figure} Figure 3: Comparison of computed $\langle {\rm
DIV}\cdot{\rm CURL}\rangle$ with observational data for ${\rm Ta}=10^3$ (triple-dot-dashed), ${\rm Ta}=10^5$ (solid) and ${\rm Ta}=10^6$ (dashed); the numerical values are taken near the top of the convection zone (z=-0.3). The shaded area covers the range of the observed data.
Open with DEXTER

For slow rotation ( ${\rm Ta}=10^3$) the correlation ${\cal C}$ is very close to the observations. For higher Taylor numbers the amplitude of $\langle {\rm DIV}\!\cdot\!{\rm CURL}\rangle$ is about one order of magnitude greater than the observed value for supergranules. The properties of the turbulence are shown in Table 1 for different rotation rates and box inclinations. The rms value of $u_{\phi}$ normalized by the sound speed $c_{\rm s}$ is given as characteristic for the flow field near the top of the convection zone. With increasing rotation rate the rms velocities change from almost uniform i.e. independent of $\theta $ in the case ${\rm Ta}=10^3$ to a nonuniform distribution for ${\rm Ta}=10^6$.

The coincidence between the simulations and the observations illustrated in Fig. 3 is expected when comparing the Coriolis numbers in the simulations with those prevailing in the supergranulation. The slower the box rotates the smaller the influence of the Coriolis force on the velocity field and, as a consequence, ${\cal C}$ would become smaller. Although the Sun has ${\rm Ta}\sim 10^7$, the supergranules rotate more slowly than suggested by their Taylor number. Indeed, for ${\rm Ta}=10^3$ we have a $\Omega^*$ similar to the value estimated for the supergranulation flow field. This explains the good agreement of the computed $\langle {\rm DIV}\!\cdot\!{\rm CURL}\rangle$ with the observations. The obtained Coriolis numbers for the higher Taylor number simulations are more related to giant cells (also the typical cell size corresponds to giant cells), hence, the effect of rotation is at least one order of magnitude greater than the observed one. Extrapolation of our results to the mesogranulation pattern leads to the value of 10-11 s-2 for the maximum(poles) $\langle {\rm
DIV}\cdot{\rm CURL}\rangle$ correlation.

3.3 Kinetic helicity

For the ${\rm Ta}=10^6$ case Fig. 4 shows the vertical distribution of the kinetic helicity (3) obtained by averaging over the horizontal directions and over time. The comparison with Fig. 2 reveals a close relationship of ${\cal C}$ and ${\cal H}_{\rm kin}$. The signs of both quantities are clearly correlated. Note that the vanishing of ${\cal H}_{\rm kin}$ at zero depth is forced by the boundary condition at the outer edge of the upper overshooting region. A relation ${\cal C} \propto {\cal H}_{\rm kin}$ has already been formulated by Rüdiger et al. (1999), the presented arguments can now be considered as numerically confirmed.Obviously, the observation of the correlation ${\cal C}$ for supergranulation clearly indicates the existence of negative (positive) kinetic helicity in the northern (southern) hemisphere even close to the surface. The standard relation between the $\alpha$-effect and kinetic helicity for hydrodynamic (!) turbulence is

 \begin{displaymath}\alpha \simeq - \tau_{\rm corr}{\cal H}_{\rm kin}
\end{displaymath} (16)

(Krause & Rädler 1980). According to Eq. (16) the observational results by Gizon & Duvall (2003) would indicate the existence of a positive (negative) $\alpha$-effect beneath the solar surface for the northern (southern) hemisphere.
\psfig{figure=0531fig4.ps,width=7cm,clip=} }
\end{figure} Figure 4: Horizontally- and time-averaged kinetic helicity ${\cal H}_{\rm kin}$ vs. depth for different latitudes in the box simulations for ${\rm Ta}=10^6$. Note the close relation of this quantity to the simulated correlation ${\cal C}$ given in Fig. 2. The shape of the helicity curve is very close to results of Miesch et al. (2000, their Fig. 22).
Open with DEXTER



Online Material

Table 1: Properties of turbulence for different Taylor numbers and colatitudes $\theta $ at z=-0.3.

Copyright ESO 2004