SEISMIC AND DYNAMICAL SOLAR MODELS. I. THE IMPACT OF THE SOLAR ROTATION HISTORY ON NEUTRINOS AND SEISMIC INDICATORS

The solar radiative zone represents 98% of the total mass of the Sun and the adopted microphysics plays a basic role in this region. The equilibrium between gravitational energy, nuclear energy production, and the energy escaping by photon interaction needs to be followed in the radiative interior over long timescales. This region is now permanently probed by helioseismology, which is a key to validate the various ingredients used in the construction of the standard solar model (SSM). One success of the SoHO space observatory has been obtained by measuring the Doppler velocity shifts down to the low frequency range of the acoustic spectrum with two instruments: GOLF, Global Oscillations at Low Frequency described by Gabriel et al. (1995) and MDI, a Michelson Doppler Imager built by Scherrer et al. (1995). This part of the spectrum, contrary to the high frequency range, is not sensitive to the variability of the subsurface layers along the solar cycle (Turck-Chièze et al. 2001b; Couvidat et al. 2003b). The corresponding modes have longer lifetimes but smaller intensities (Bertello et al. 2000; García et al. 2001). From these observations, a very clean sound speed profile has been established down to 0.06 R_☉ together with a reasonable density profile. It has been compared to the theoretical models for different solar compositions. For the latest sets of solar abundances with lowered oxygen abundance (Holweger 2001; Asplund et al. 2009; Caffau et al. 2008), the theoretical profiles show clear differences compared to the observed sound speed, that still need to be understood (Turck-Chièze et al. 2001b, 2004a).

2.1.1. The Nuclear Processes

2.1.2. The Impact of the Solar Composition on the Radiative Zone Properties

One important question has emerged from the new spectroscopic estimates of the carbon, nitrogen, and oxygen atmospheric abundances, leading to a reduction by 20%–30% of the related abundances: "Are we so sure that we understand properly the inner composition of the Sun in the radiative zone"? This question is an important issue because the solar abundance is chosen as a reference for the stellar evolution models. For instance, up until recently, the abundance of these nuclides was larger than that of any other heavy element by about a factor 10 in Population II stars, so the change in the solar CNO abundance has a strong impact on the evaluation of their chemical pattern (Turck-Chièze et al. 2004a). The situation is now clarified because the three independent analyses (Holweger 2001; Asplund et al. 2009; Caffau et al. 2008) converge to about the same reduction. A reduction by a similar factor (30%) of the iron content also appeared in the early nineties, but with rather small effects on the whole radiative zone (Turck-Chièze & Lopes 1993).

As far as CNO elements are concerned, this substantial reduction has a real impact on the radiative sound speed profile mainly due to the subsequent change in the oxygen opacity (see Figure 1) and modifies the depth of the convective envelope (Turck-Chièze et al. 2004a; Bahcall et al. 2005). The confrontation between helioseismic indicators and results of SSMs was examined extensively (Turck-Chièze et al. 2004a; Bahcall et al. 2005, Guzik et al. 2005). It is frustrating to note that the sound speed and density profiles of the model including the updated solar chemical composition are very similar to those corresponding to the case where we omitted to account for the slow gravitational settling of the nuclides (Thoul et al. 1994; Brun et al. 1998) and that this last update destroys the apparent agreement with the observed one. We may thus quantify the effect of the change of composition since the inclusion of the microscopic processes corresponds to a reduction by 10%–12% of all the heavy elements. Consequently a 20% variation of CNO abundances is almost equivalent to the 7%–9 % variation of all the heavy elements due to gravitational settling, or to an equivalent change of opacity coefficients if all the elements were concerned. See Turck-Chièze et al. (1997) where the impact of a lot of changes is considered and Table 1 for the consequence on neutrino predictions.

Zoom In Zoom Out Reset image size

Effectively, a change in the CNO abundances also affects the mean total opacity of the stellar plasma, and thus the evolution of the Sun up to its present status as well as its present structure. The knowledge of the opacity coefficients used in a solar model is up to now purely theoretical. The conditions of temperature and density in the radiative zone ensure that the plasma is only fully ionized for its main constituents, hydrogen and helium, but heavier species such as iron down to oxygen become partially ionized in the radiative zone as summarized in Figure 1. It is why the bound–bound and the bound–free interactions of photons with matter limit the evacuation of the energy produced within the inner 25% of the radius, corresponding to half the solar mass (Iglesias & Rogers 1996). The small amount of iron (only 2.8 × 10⁻⁵ of the hydrogen contribution in fraction number) contributes to about one-fifth of the opacity cross section in central conditions (see Figure 1 and Turck-Chièze et al. 1993 for more details). The oxygen is the second contributor and plays a crucial role in determining the base of the convective zone.

Up until now the details of the ion interactions have been verified only indirectly through acoustic pulsation eigenmodes. These modes probe plasma properties throughout the Sun from 6% R_☉ up to the limit of the convective zone. However, they depend on both the detailed composition and opacity. To disentangle both contributions would require a high radial resolution of the seismic data allowing us to capture some specific bound–bound effect through bumps in the radial sound speed profile in the radiative zone (Turck-Chieze 1998) that is still out of reach at present.

Some proposed solutions to reconcile the standard model with revised CNO abundances and helioseismology (Guzik et al. 2005; Lin et al. 2007; Basu & Antia 2008) have not been retained. Some others seem to favor the old composition (Basu et al. 2007; Serenelli et al. 2009). We explore, here, the idea that the SSM framework is questionable. Considering that any improvement has always been considered like a new informative fact and that a better determination of CNO has been awaited since the nineties (Turck-Chièze et al. 1993), the present facts suggest three reasonable solutions which could all exist simultaneously:

• 1.  
  The gravitational settling is not well constrained today for CNO and heavy elements. The photospheric helium content deduced from helioseismology (Vorontsov et al. 1991) of Y = 0.25 ± 0.03 (and confirmed by Basu 1995) for OPAL equation of state (EOS) has confirmed the need to take into account this slow atomic diffusion introduced first by Cox et al. (1989). This process leads today to a reduction of about 10%–15% of the He mass fraction at the solar surface (Christensen-Dalsgaard et al. 1993; Thoul et al. 1994; Brun et al. 1998; Michaud 2004) but we cannot check yet if the order of magnitude is correct also for CNO and iron.
• 2.  
  The opacity coefficients used in the models have never been verified by laboratory measurements under the physical conditions of the stellar interiors and may be a source of uncertainty. The high energy lasers developed in France (LIL, LMJ) and in the United States (NIF) will offer conditions corresponding to the solar radiative zone and might help in constraining these quantities. We are preparing experiments to measure the absorption energy spectrum of isolated elements and mixtures in plasma conditions which are more and more useful for the stellar community (Bailey et al. 2007; Loisel et al. 2009; Turck-Chièze et al. 2009). Laboratory opacity and EOS experiments might benefit in the coming decade from an equivalent effort to the one dedicated to nuclear reaction rate measurements. A detailed verification of the interaction between photons and plasma will contribute to disentangle the different processes in the deep interior of stars.
• 3.  
  The radiative zone of the Sun is not well reproduced by the hypotheses of the SSM, and we need to improve the models by the addition of dynamical processes. It is one of the reasons to develop the DSM which must introduce the dynamics of the tachocline, the gravity waves generated at the base of the convective envelope and propagating in the radiative interior, the presence of a potential magnetic field, and the different transport processes induced by the rotation. Guzik & Mussak (2010) examine the idea of mass loss and accretion, and in a coming work, we question the hypothesis of the energetic balance (Turck-Chièze et al. 2010).

Such a complete model of the Sun is not ready to be compared to seismic and neutrino probes; it is why we have developed an intermediate model, the seismic model described in the following section. A first step toward the DSM is then to confront all the probes with solar models including rotation, which is the purpose of the present paper. We shall produce models in the following sections which can be considered as the DSM1 step.

The quality of the seismic observations is such that one can build a seismic model. This model is designed to reproduce the measured solar sound speed using the classical stellar evolution equations and slightly changing some recommended specific absolute values of physical inputs within their own uncertainties (Turck-Chièze et al. 2001b, 2004a; Couvidat et al. 2003a). It can be established for the new solar composition and will converge to the same predictions by definition.

The interest of such a model comes from the idea that the standard framework is obviously too crude to reproduce all the existing observed quantities. This model is useful to improve the prediction of the neutrino fluxes (Tables 1 and 2) or the gravity mode frequencies because it implicitly includes the known sound speed profile (Turck-Chièze et al. 2004a; Couvidat et al. 2003a; Mathur et al. 2007).

Table 1 illustrates the time evolution of the predicted ⁸B neutrino flux which depends strongly on the central temperature and consequently, on the details of the plasma properties. The seismic model prediction (Turck-Chièze et al. 2001b; Couvidat et al. 2003a) agrees remarkably well with the measured value obtained with the SNO detector (filled with heavy water), which is sensitive to all the neutrino flavors (Ahmed et al. 2004). For the gallium or water detector predictions, one needs to inject the energy dependent reduction factor due to the fact that the electronic neutrinos are partially transformed into muon or tau neutrinos. This property of neutrinos was confirmed last year by the Borexino results (Arpesella et al. 2008). Doing so, the agreement between the predictions of the seismic model and all the detectors is extremely good (Turck-Chièze et al. 2005). Table 2 confirms that it is also true for the beryllium neutrinos.

Another interesting property of the SSeM is the fundamental periodicity P₀ which characterizes the asymptotic behavior of the gravity modes nearly equally spaced in period for frequencies below ∼100 μHz. $P_0=2{\pi }^2 \left(\int _{0}^{r_c} \frac{N}{r} dr \right)^{-1}$ where N is the Brunt-Väisälä frequency and r_c is the radius at the convective bottom. Before the launch of SoHO, there was a great dispersion in the theoretical predicted values of P₀. Following Hill et al. (1991), its value was varying from 29 minutes to 63 minutes depending on the models. Today, standard and seismic models agree within one minute, which helps to predict the general properties of the gravity modes (Mathur et al. 2007) because models have been improved by the understanding of the property of the modes and their sensitivity. Nevertheless the present standard model, including the new CNO composition, largely deviates from the seismic observations; consequently, the value of P₀ might be better determined by the seismic model or any coming model in agreement with seismic observations. Of course, the detection of individual gravity modes might improved the density profile in the core and consequently will put more pressure on the seismic model in the core.

Despite the remarkable agreement between the two probes of the central region of the Sun (neutrinos and helioseismology), SSeM is not a physical model and SSM predictions with the new solar composition encounter a series of difficulties. In the following section, we present a first step to transform SSM or SSeM into a more realistic DSM1, by introducing rotation and the associated transport of matter and angular momentum in the computations.

The present Sun is a slowly rotating star (≃2.2 km s⁻¹ at the equator, 1.7 km s⁻¹ near the poles). It is thus a reasonable first-order approximation to treat rotation with a perturbative approach and to adopt the hypothesis of spherical symmetry in the standard model framework. We may nevertheless recall that the solar superficial deformations can be measured and that the quadrupolar moment is of the order of 1.84 × 10⁻⁷ (Lydon & Sofia 1995) and will be improved by the microsatellite PICARD (Thuillier 2005).

Even if departures from spherical symmetry can be considered small in the solar case, they are sufficient for rotation to induce large-scale circulations in the radiative and convective interior, that will simultaneously advect angular momentum, and chemical species (Busse 1982; Zahn 1992). Moreover, in case of radial differential rotation, which could be the case in the solar deepest interior, different hydrodynamical instabilities may develop that will generate hydrodynamical turbulence in the radiative regions. Laboratory Couette–Taylor experiments conducted for Reynolds numbers Re ≃ 10⁶ by Richard & Zahn (1999) indicate the development of turbulence in the flows. This gives a good hint that in stellar interiors, where the typical Reynolds number is of about 10¹² (solar case), the hydrodynamical instabilities associated with differential rotation eventually become turbulent.

Following Zahn (1992), Maeder & Zahn (1998), and Mathis & Zahn (2004), the transport of angular momentum and chemical species by meridional circulation and shear-induced turbulence has been added to the stellar structure equations in order to estimate the effects of such a transport on the evolution and the structure of the Sun up to the present time.

Under the assumption of shellular rotation ensured by a strong anisotropy of the turbulent diffusivities D_h ≫ D_v, all the relevant variables and vectorial fields may be described as the sum of the mean value over an isobar and of a second (or fourth) order perturbation (Zahn 1992; Mathis & Zahn 2004). This means that the variables are projected on a basis of Legendre polynomials, an approach very similar to that used in helioseismology. This allows us to separate angular and radial parts and thus to account for the rotational transport of angular momentum and chemicals over secular timescales in one-dimensional stellar evolution codes, while two-dimensional or three-dimensional secular stellar evolution is still far from maturity.

Using this formalism, the momentum equation

$$\begin{equation} \rho \left[\frac{\partial \vec{V}}{\partial t} + (\vec{V}. \vec{\nabla }) \vec{V} \right] = - \vec{\nabla }P - \vec{\nabla }\phi + \vec{\nabla }. || \tau || \end{equation} \tag{ 1 }$$

becomes, when averaging over the isobar and using an azimuthal projection:

$$\begin{equation} \rho \frac{d}{dt}{ (r^2 \overline{\Omega })} = \frac{1}{5 r^2}\frac{\partial }{\partial r} (\rho r^4 \overline{\Omega }U_2) + \frac{1}{r^2}\frac{\partial }{\partial r} \left(\rho \nu _v r^4 \frac{\partial \overline{\Omega }}{\partial r} \right), \end{equation} \tag{ 2 }$$

where $\overline{\Omega } = \Omega (r)$ is the mean angular velocity on the isobar of radius r and ν_v is the vertical turbulent viscosity associated to the shear instability (see Equation (8)).

U₂ is the vertical meridional velocity component:

$$\begin{equation} U_2 = \frac{P}{\rho g C_p T [\nabla _{\rm ad} - \nabla +(\phi /\delta) \nabla _\mu ]} \left(\frac{L}{M_\star }(E_\Omega + E_\mu)\right), \end{equation} \tag{ 3 }$$

with P is the pressure, C_p is the specific heat, E_Ω and E_μ are terms depending on the Ω- and μ-distributions respectively, up to the third-order derivatives and on various thermodynamic quantities (see details in Maeder & Zahn 1998).

Equation (2) is split into a system of four first-order equations in Ω complemented by an equation describing the evolution of the horizontal mean molecular weight fluctuations due to the vertical μ gradient for a given horizontal turbulence and a vertical velocity:

$$\begin{equation} \frac{d \Lambda }{dt} + U_2(r) \frac{\partial \ln \overline{\mu }}{\partial r}= -6 D_h \Lambda, \end{equation} \tag{ 4 }$$

where $\Lambda = \frac{\tilde{\mu }}{\mu }$ and D_h is the horizontal turbulent diffusivity.

The equation for the transport of nuclides by meridional circulation and shear-induced turbulence can be written as a diffusion equation (see, e.g., Chaboyer & Zahn 1992):

$$\begin{eqnarray} \left(\frac{{d} X_i}{{d}t}\right)_{M_r}&=& \frac{\partial }{\partial M_r}\left[(4\pi r^2\rho)^2\left(D_v+D_{\rm eff}\right) \frac{\partial X_i}{\partial M_r}\right]\nonumber\\ && + \left(\frac{{d}X_i}{{d}t}\right)_{\rm nucl} + \left(\frac{{d}X_i}{{d}t}\right)_{\rm micro}, \end{eqnarray} \tag{ 5 }$$

where X_i is the mass fraction of the ith nuclide, the second and third terms are, respectively, the nuclear (nucl) and gravitational settling (micro) terms. D_v ≡ ν_v is the vertical component of the turbulent diffusivity. D_eff is the diffusion coefficient associated to the action of the meridional circulation, and it is given by

$$\begin{eqnarray} &&D_{\rm eff} = \frac{\left(r U_2\right) ^2}{30D_h}. \end{eqnarray} \tag{ 6 }$$

These transport equations have to be solved at each evolutionary time step, adding a set of five coupled nonlinear equations to the usual set of structure equations. This explains the small number of stellar evolution codes including this treatment: the STAREVOL code (Palacios et al. 2006) and the Geneva stellar evolution code (Eggenberger et al. 2008). It has been recently included in the CESAM stellar evolution code.

Following the effects of rotation according to the above described formalism is even more difficult in the solar case due to the small perturbation they cause on the solar structure. In order to test the quality of our results, we have decided to confront the results obtained with two different stellar evolution codes using distinct numerical approaches. We have computed rotating solar models with both STAREVOL and CESAM codes. In the following, we present the codes, the numerical implementation of the Maeder & Zahn formalism, and the main characteristics of the solar models that will be discussed in the forthcoming sections.

3.2.1. CESAM

3.2.2. STAREVOL

The Sun is the only star where we can perform quantitative comparison of the impact of the above prescriptions on the sound speed and the rotation profiles which can now be deduced from helioseismic data. So, it is interesting to discuss in details the role of the different terms and the choice for the prescriptions. In the following section, we shall present two types of models: one model with an extremely slow initial rotation which does not justify external magnetic braking (A or slow models), and two models with greater initial rotation and magnetic braking at the surface (models B and C, or intermediate and strong rotation models). The main characteristics of these models are given in Table 3. These cases allow us to discuss the order of magnitude of the meridional circulation, the diffusion coefficients, and the gradient of the internal angular velocity profile, together with their different roles in the building of the present solar rotation profile.

For these models, we have adopted the Mathis et al. (2004) prescription for the horizontal component of the turbulent diffusivity D_h:

$$\begin{equation} D_h =\sqrt{ \left(\frac{\beta }{10} \right)(r^2\overline{\Omega })\left[r | 2V_2 - \alpha U_2 |\right]}. \end{equation} \tag{ 7 }$$

The vertical component of the turbulent diffusivity D_v is from Talon & Zahn (1997):

$$\begin{equation} D_v = \frac{Ri_c}{N^2_T\big/(K_T + D_h) + N^2_\mu \big/D_h}(r \partial _r \overline{\Omega })^2 \end{equation} \tag{ 8 }$$

with Ri_c = 1/6 the critical Richardson number, K_T is the thermal diffusivity, N_T and N_μ are the chemical and thermal parts of the Brunt-Väisälä frequency N² = N²_T + N²_μ.

Referring back to Talon & Zahn (1997), let us recall here that for the shear instability to develop, not only the Richardson criterion must be satisfied, but the sheared flow shall be turbulent. This is controlled by the Reynolds criterion. In our models, we used a critical Reynolds number Re_c = 10ν_mol. When ν_v < Re_c, D_v is replaced by ν_mol in the transport equations for angular momentum and nuclides.

The choice of slow and mild rotators on the zero age main sequence (hereafter ZAMS) refers to the two identified populations of young stars as pinpointed by Bouvier (2009): the slow rotators for which the angular momentum evolves like Ω³ and the fast rotators for which the angular momentum varies like Ω.

The A models explore the extreme case of slow rotators on the pre-main sequence (hereafter PMS) and on the ZAMS, with a surface equatorial velocity of about 2.2 km s⁻¹ remaining almost constant from the ZAMS to the present age of the Sun. No braking is applied at the surface. These models do not actually reproduce observed stars, and can be considered more as academic cases used to show the respective impact of the meridional circulation and shear turbulence on the shape of the angular velocity profile even with such small velocities. They are also used to estimate interesting quantities relevant for the present observations.

The B and C models are mild rotators at their arrival on the ZAMS (approximately 20 km s⁻¹ and 50 km s⁻¹, respectively). These models undergo magnetic braking at the arrival on the main sequence, similar to the solar-mass models by Talon & Charbonnel (2005).

Rotation is included from the PMS in all our models. We assume solid-body rotation as the initial state in the completely convective PMS star. When magnetic braking is applied at the arrival on the main sequence, we use the formalism developed by Kawaler (1988; see also Bouvier et al. 1997; Palacios et al. 2003) in order to obtain the solar surface equatorial velocity when the model reaches 4.6 Gyr:

$$\begin{equation} \left(\frac{dJ}{dt}\right) \, = \, -K \Omega ^3\left(\frac{R}{R_\odot }\right)^{1/2}\left(\frac{M}{M_\odot }\right)^{-1/2} (\Omega < \Omega _{{\rm sat}}), \end{equation} \tag{ 9 }$$

$$\begin{equation} \left(\frac{dJ}{dt}\right) \, = \, -K \Omega \Omega _{{\rm sat}}^2\left(\frac{R}{R_\odot }\right)^{1/2} \left(\frac{M}{M_\odot }\right)^{-1/2} (\Omega \ge \Omega _{{\rm sat}}). \end{equation} \tag{ 10 }$$

This formulation corresponds to a field geometry intermediate between a dipolar and a radial field (Kawaler 1988), and it is widely used in the literature (Chaboyer et al. 1995; Krishnamurthi et al. 1997; Bouvier et al. 1997; Sills et al. 2000). The parameter K in Kawaler's law is related to the magnitude of the magnetic field strength, and is adjusted accordingly with the adopted initial rotation velocity. Ω_sat expresses the fact that the magnetic field generation saturates beyond a certain evolutionary point. This saturation is actually required in order to retain a sufficient amount of fast rotators in young clusters, as originally suggested by Stauffer et al. (1987). Here, we adopt Ω_sat =  14Ω_☉ following Bouvier et al. (1997). Table 3 presents the main characteristics of models A, B, and C computed with STAREVOL and CESAM.

We consider first the extreme case where the initial angular momentum is small so that the model arrives on the ZAMS as a slow rotator (υ_ZAMS ≈ 2.2 km s⁻¹) without undergoing any magnetic braking.

Table 4 shows the evolution of several quantities, among them the surface velocity in the center and at the surface of model A_C. With an initial rotation of about 0.1 km s⁻¹, the surface velocity is seen to increase to reach 2.19 km s⁻¹ at the arrival on the ZAMS (t = 390 Myr), and then remains almost constant up to the present age of the Sun. In the central regions, the increase of the rotation velocity is maintained during the main sequence, so that at the present age, υ_C is about 4 times larger than υ_S. This velocity gradient clearly appears in Figure 2 and on the left panel of Figure 7, where the evolution of the angular velocity Ω is displayed for models A_C and A_S. We may note here the very good agreement between the profiles obtained with CESAM and STAREVOL for solar models calibrated within 1%, thus validating the numerics of both codes. The building of the Ω-gradient is due on the PMS to the rapid contraction of the central regions when radiation becomes the dominant energy transport process in the core, as well as to the combined effects of meridional circulation and shear-induced turbulence.

Zoom In Zoom Out Reset image size

Despite the very slow rotation of these models, a small temperature gradient is generated between the pole and the equator and a meridional circulation develops in the radiative interior. The velocity of the meridional currents, U₂, is of the order of 10⁻⁹ cm s⁻¹ in absolute value at 4.6 Gyr as can be seen in Figure 3, left. This flow is extremely small compared to the meridional circulation in the convective zone. The circulation consists in a single counterclockwise loop, e.g., U₂ negative. The flow is directed downward, peaking near the base of the convective zone, extracting angular momentum from the central region to the base of the convective envelope. The profiles obtained with both codes are of similar shape and amplitude. They present irregularities mainly due to the mean molecular weight variations $\Lambda = \frac{\tilde{\mu }}{\mu }$ that introduce nonlinearities in the system because some noise is coming from the mean molecular weight gradient. Let us here recall that the transport of angular momentum introduced in the codes is treated in a self-consistent but simplified manner. Such models give hints, order of magnitudes, and global shapes, but should not be expected to deliver exact profiles.

Zoom In Zoom Out Reset image size

The diffusion coefficient profiles are also similar in the two codes (see Figure 4, left). The effective diffusion coefficient D_eff is proportional to (U₂)²/D_h and consequently is very small (about 10⁻⁵–10⁻³ cm² s⁻¹) due to the slow meridional circulation velocity and the large horizontal component of the turbulent diffusivity ensured by the use of the Mathis et al. (2004) prescription. The Reynolds number of the sheared flow is lower than the critical Reynolds number Re_c ∝ 10ν_mol in the bulk of the radiative zone. The flow is not turbulent, so the vertical component of the turbulent diffusivity D_v is essentially set equal to the local molecular viscosity (see Section 3.3). The effective diffusion coefficient is actually much smaller than the microscopic diffusivity (see, e.g., Brun et al. 1999), which is of the order of 10 cm² s⁻¹ in solar models. Consequently, the microscopic diffusion is the most efficient process to transport chemical species below the convective envelope in the models with slow rotation. This shows up in the helium mass fraction profile presented in quadrant III of Figure 5. There, the step in the helium profile at the base of the convective zone reveals the gravitational settling of helium in this region, and this, regardless of the existence of meridional flows and shear. If the young Sun was a slow rotator, the chemical stratification probed by helioseismic data in the region below the convective envelope would be the same as if the Sun had not rotated (provided that no other dynamical process than rotation is included).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Of course, one cannot ensure that the diffusivity coefficients are properly estimated, so it is interesting to see the impact of a different prescription for them. Using the STAREVOL code, we compute model A'_S using the prescription of Zahn (1992) for D_h:

$$\begin{equation} D_h = r [2 V_2 - \alpha U_2]. \end{equation} \tag{ 11 }$$

Figure 5 shows the profiles of the angular velocity Ω, the meridional circulation velocity U₂, the helium mass fraction Y, and the diffusion coefficients D_h, D_v, and D_eff for models A_S and A'_S. As expected, the horizontal turbulent diffusion coefficient D_h is 3 orders of magnitudes smaller when using Zahn's (1992) prescription (see also Mathis et al. 2004). This translates into an increase by the same amount of the effective diffusivity D_eff. However, D_eff remains smaller than D_v and ν_mic in model A'_S, so that the helium and angular velocity profiles remain unchanged compared to model A_S. From this comparison, we may conclude that in the case of a Sun that has always been rotating slowly, the choice of the prescription for the horizontal component of the shear-induced turbulent diffusivity does not affect significantly the structure, rotation and initial chemical stratification.

We have also computed two models with a higher initial angular momentum, corresponding to a surface velocity of about 20 km s⁻¹ and 50 km s⁻¹ at the arrival on the ZAMS. These models undergo a strong angular momentum extraction on the early main sequence to reach a surface velocity of about 2 km s⁻¹ at the age of the present Sun. This phenomenon associated to the idea of a magnetic braking due to the decoupling from the disk environment is modeled using the Kawaler formula as detailed in the previous section.

Figure 6 compares the results of the two models B_S and C_S at 4.6 Gyr. The diffusion coefficients D_v and D_eff entering the diffusion equation for the chemical species are similar in models B_S and C_S, so they lead to the same profile for the helium mass fraction. The angular velocity is slightly larger in the central regions in model C_S. This model is submitted to a more efficient braking in order to reach a surface velocity of about 2–2.2 km s⁻¹ at the age of the Sun. In fact, when magnetic braking is initiated, the meridional circulation velocity U₂ is five times smaller in model C_S than that found in model B_S in absolute value, but it is positive (negative in model B_S) generating an efficient inward transport of angular momentum and leading to a faster and larger increase of the central velocity than in model B_S. The behavior found for model C_S is in all points similar to that reported by Eggenberger et al. (2005), and also resembles that described by Meynet & Maeder (2000) for more massive and fast rotating stars. The evolution further on the main sequence is similar for models B_S and C_S, with U₂ rapidly becoming negative in the later (see figure top right of Figure 6). The profiles displayed in Figure 6 show no significant difference for the diffusivities, meridional circulation, and helium profile at 4.6 Gyr, so that in the following we will focus on models B_S and model B_C.

Zoom In Zoom Out Reset image size

The meridional velocity profiles at 4.6 Gyr for the B models are shown on the right panel of Figure 3. The velocity of the meridional currents is 3 orders of magnitude larger compared to models A_S and A_C previously described. At the ZAMS, the meridional circulation is maximal below the convective envelope as a result of the strong magnetic torque applied at the surface. U₂ is then 2 orders of magnitude larger than it will be later on the main sequence, when the efficient spin down of the surface layers is accompanied by a global weakening of the meridional circulation. The form of the meridional circulation velocity profile U₂ is maintained over the main sequence evolution, and mainly consists in one counterclockwise cell peaking below the convective envelope and extracting angular momentum from the central regions, this time more efficiently than in the slowly rotating solar models. Similar to what was obtained by Decressin et al. (2009) for their 1.5 M_☉ model, the meridional circulation is mainly driven by the local losses of angular momentum due to the magnetic braking. Meridional circulation and shear turbulence do not ensure an efficient transport of angular momentum. One notes that U_2,max is larger in model B_S by 25% compared to model B_C.

The gradient of angular velocity is much larger in these models compared to models A_C and A_S, and Ω increases by more than 1 order of magnitude from the surface to the solar core. We will discuss in more detail the implications of such a steep profile when confronted with helioseismic data. The match between the two codes shown in Figure 7 is not as good as that found for slowly rotating models. This is due to the difference in the numerical treatment used in both codes and the way the braking is introduced. The rotation profile is even steeper for models computed with CESAM, in which braking occurs very efficiently on the ZAMS so that the Ω-profile does not evolve at all in the last 3.6 Gyr. In model B_S, angular momentum is extracted from the central parts as the star evolves and both the central and surface values of the angular velocity decrease with increasing time so that the gradient is less steep at 4.6 Gyr than it is at 1 Gyr.

Zoom In Zoom Out Reset image size

Diffusion coefficients are displayed in the right panel of Figure 4 at 4.6 Gyr. They have similar shapes and amplitude, which implies that the chemical stratification should also be similar in models B_S and B_C. Although U₂ is larger by 3 orders of magnitude compared to A models, the horizontal component of the turbulent diffusivity is also much larger leading to D_eff that is again smaller than 1 cm s⁻¹, which remains much smaller than both the local molecular diffusivity and the microscopic diffusivity. Meridional circulation is thus not contributing significantly to the transport of chemical species. On the other hand, in these mild rotating models, the Reynolds number associated with the sheared flow is larger than the critical Reynolds number Re_c, and the flow eventually becomes turbulent. The vertical turbulent diffusivity reaches values of the order of 100–1000 cm s⁻² below the convective envelope, which is at least 10 times larger than the microscopic diffusivity in this region. As a consequence, the abundance profiles are flattened in both models B_S and B_C (see, e.g., see solid line in quadrant III of Figure 6 for model B_S), and thus differ significantly from those derived for SSMs. If the young Sun was a fast rotator and experienced magnetic braking during the early main sequence without any inhibition of the action of the shear-induced turbulent mixing, the chemical stratification of the present Sun below the convective zone might bear the signature of this past rotational history.

After more than 10 years of observations with GOLF and MDI instruments located onboard SoHO, and with the facilities of the GONG ground-based network, we have derived very important constraints on the radiative zone rotation profile. The determination of the splittings of a large number of acoustic modes has definitely established that the rotation profile in the part of the radiative zone which is not influenced by the nuclear reaction rates, e.g., the region between 0.25 and 0.68 R_☉, is really flat with invisible latitudinal differential variation (Couvidat et al. 2003b; Eff-Darwich et al. 2008), in great contrast with the latitudinal differential profile of the convective zone. In the central region, below 0.25 R_☉, only gravity modes may inform the rotation profile. Some individual candidate modes at high frequency have been observed by the GOLF instrument (Turck-Chièze et al. 2004b; Mathur et al. 2007; García et al. 2008), and two solutions have been extracted from these first data depending on the interpretation of the pattern attributed to an ℓ = 2, n = 2 mode: a slightly reduced rotation rate or an increase by about a factor 3 in the central region (Mathur et al. 2008). Nevertheless, only the solution of a rapid rotating core is compatible with the other kind of detection using the asymptotic behavior at low frequency which is detected with a very high probability (García et al. 2007, 2008). It could presume a rather complex solar core rotation larger than the rest of the radiative zone, with a different axis, and maybe with some manifestation of the radiative zone magnetic field. The other solution cannot be totally excluded but is not the most probable. Again, the core profile needs to be confirmed with extensive data from SoHO and with an improved instrument like the GOLF-NG concept (Turck-Chièze et al. 2006, 2008) which will measure velocity displacements at eight heights in the solar atmosphere. Considering that the complete solar rotation profile might be accessible to observation, it is important to predict it properly down to the central region.

This paper is a new step toward a detailed understanding of how the different dynamical processes can influence the present rotation of the Sun. In this study we have followed three different cases because we can only have some indirect information on the solar internal rotation profile at the end of its contraction phase. We have already seen in the previous section that the initial angular momentum content and the rotation history, in particular the fact of undergoing magnetic braking on the ZAMS, shape the angular velocity profile, determining the absolute value as well as the gradient. In this section, we wish to focus on the evolution of the angular velocity profile. We present hereafter a detailed analysis of the construction of the predicted internal profile at the age of the Sun and confront these predictions with helioseismology.

Compared to previous studies by Eggenberger et al. (2005) and Talon & Charbonnel (2005) where rotation was only included from the ZAMS, with an initial flat profile, we have decided to follow rotation from the PMS when the star is actually completely convective. We also assume solid-body rotation as the initial state, but by the time the models reach the ZAMS, the convective region has retreated to the surface and the Ω-profile is not flat anymore.

Tables 4 and 5 list the evolution of the radius, the luminosity, the central and the superficial velocity for models A_C and B_C, and some other indicators of the contraction phase and of the involved nuclear reactions. The evolution of the angular velocity profile for these models is also presented in the left panels of Figure 7.

For model A_C, where the initial rotation rate is of about 0.08 km s⁻¹, the contraction of the star during the PMS rapidly generates a radial gradient in the Ω-profile. The very slow meridional currents and the small amount of shear generated in such slowly rotating models are not very efficient to transport the angular momentum. The differential rotation established during the early evolution by the contraction of the inner regions is maintained and slightly amplified during the main sequence by the advection term. The predicted overall contrast between core and surface velocities is of about a factor of 4, which is comparable to that obtained by Eggenberger et al. (2005) in their slowly rotating model.

Models B and C are similar to those presented in the three main previous studies of the effect of rotation on the solar evolution by Chaboyer et al. (1995), Eggenberger et al. (2005), and Charbonnel & Talon (2005). The evolution of the angular velocity profile is in agreement with that obtained in the later study, where STAREVOL models were also used. Comparing the evolution displayed in Figure 7 with that of Figures 1 and 2 of Charbonnel & Talon (2005) and Eggenberger et al. (2005), respectively, one can measure the impact of the numerical approach. In all these cases, if the general forms and amplitudes of the angular velocity profile at the age of the Sun show encouraging similarities (see also Figure 8), the in between evolution can be quite different. Models B_S and B_C undergo a very strong braking at the arrival on the ZAMS leading to a quick spin down of the convective envelope and a contrast between the surface and central angular velocities of about a factor 15 for model B_C already at 1 Gyr, and of about 7 in model B_S at the same age. In the later model, the contrast increases during the main sequence to reach a factor of 10 between the core and the surface at 4.6 Gyr. In both models, the Ω-gradient steepens in the bulk of the radiative zone as the star evolves. The profiles are flatter in model B_S, yet not as flat as those presented in Figure 2 of Eggenberger et al. (2005).

Zoom In Zoom Out Reset image size

A detailed analysis of the relative importance of meridional circulation versus shear-turbulence in our models A_S, B_S, and C_S using the tools presented in Decressin et al. (2009), demonstrates that the behavior of these solar models is in all points similar to that obtained for the 1.5 M_☉ model shown in that paper. The angular momentum transport is ensured by meridional currents. In the case of models B and C, these currents, although slow, are primarily generated by the action of the applied torques resulting from the action of magnetic braking. In the case of the slow model A, no torque is applied and the angular momentum flux is smaller yet also dominated by the meridional circulation. The flux of angular momentum due to shear turbulence is more than 5 orders of magnitude smaller than that attributed to the meridional circulation in all three cases.

These detailed calculations are self-consistent and rely on the use of a set of parameters that have been extensively tested throughout the H-R diagram. They consider initial rotation generally smaller than used in the previous literature (Chaboyer et al. 1995; Krishnamurthi et al. 1997) and show slower differential rotations in the core than found in previous works (Pinsonneault et al. 1989) or the most recent of Charbonnel & Talon (2005). Figure 8 shows a first comparison between the different profiles and the presumed observed rotation profile. It is clear that the models with an initial reasonable rotation velocity + magnetic braking are not supported by the observations. As indicated by García et al. (2007), the helioseismic data reject an increase by a factor 10 for the core rotation. On the contrary, the slow initial rotation leads to a profile for the present Sun for which the contrast between the core and the surface angular velocity is very similar to the presumed observed one. The contrast is similar but of course none of the different models investigated here create a flat rotation between 0.25 and 0.7 solar radius. Within the framework of the present study, models A_S and A_C corresponding to an extremely slow rotating young Sun are more compatible with the central observations. However, this cannot be taken at face value because this model seems rather unrealistic because it did not contain any magnetic braking (see Figure 2 of Bouvier 2009). So, one can imagine a slightly greater initial rotation, typically 5 or 10 km s⁻¹ at the ZAMS, with a central rotation value slightly eroded by some other process. The discrepancy between helioseismology and the profiles presented here calls for the inclusion of additional dynamical processes that efficiently extract angular momentum from the radiative interior. This result confirms the most recent works. Probably the radiative magnetic field (Eggenberger et al. 2005; Yang & Bi 2006; Denissenkov & Pinsonneault 2007) or the internal gravity waves (Charbonnel & Talon 2005), or both, could flatten the profile, but their role could be smaller than thought previously.

Figure 9 compares the squared sound speed profile difference between seismic observations and the models computed with CESAM. The used acoustic modes are given in Couvidat et al. (2003a). We present the three models (standard, low rotation, and moderate rotation) computed with the version of CESAM used in the present study. These models slightly differ from the standard models already published due to the large reduction of the (¹⁴N, p) reaction rate mentioned in Section 2. The CNO combustion is rather small, but this strong decrease, previously studied in Turck-Chièze et al. (2001a), leads to a negative sound speed difference in the core partly compensated by a small increase in the rest of the radiative zone. These two differences increase the difference with observations. But what we would like to emphasize in this figure is the intercomparison between models with and without rotation.

Zoom In Zoom Out Reset image size

Even if these models must be more representative of the real Sun, the sound speed profile obtained for the two calibrated CESAM rotating models A and B is in less agreement with the observed profile than the standard model. In fact, the rotation-induced mixing has reduced the effect of the microscopic diffusion (see also Section 4) partly inhibited by the macroscopic motions. As discussed earlier, such an effect is practically not visible in models A, where microscopic diffusivity remains larger than both the vertical turbulent and the effective diffusivity, but it is clearly visible in models B and C rotating faster on the ZAMS.

As a natural consequence, Table 6 shows also that the central temperature is slightly reduced in rotating models in comparison with the SSM. For the CESAM models calibrated in radius and luminosity at 10⁻⁴ level, the predictions for the neutrinos are always worse than those of the seismic model (see details in Table 1). Moreover, in CESAM models where we have kept superficial Z/X constant, one may notice that the surface helium content also increases with the initial rotation due to turbulent flow just below the convective zone, as previously mentioned in Section 4 for the STAREVOL results but by a smaller amount (comparison of the two codes in Table 6). This will contribute to reconciling the observed helium content to the predicted one when using the most recent composition of Asplund et al. (2009), if the initial rotation is sufficient to allow turbulent flow at the base of the convective zone and in the tachocline.