Radiative hydrodynamics simulations of red supergiant stars

A. Chiavassa; B. Plez; E. Josselin; B. Freytag

doi:10.1051/0004-6361/200911780

Home

All issues

Volume 506 / No 3 (November II 2009)

A&A, 506 3 (2009) 1351-1365

Full HTML

Highlight

Free Access

Issue		A&A Volume 506, Number 3, November II 2009


Page(s)		1351 - 1365
Section		Stellar atmospheres
DOI		https://doi.org/10.1051/0004-6361/200911780
Published online		27 August 2009

A&A 506, 1351-1365 (2009)

I. interpretation of interferometric observations

A. Chiavassa^1,2 - B. Plez^2,4 - E. Josselin² - B. Freytag^3,4

1 - Max-Planck-Institut für Astrophysik, Karl-Schwarzschild-Str. 1, Postfach 1317, 85741 Garching b. München, Germany
2 - GRAAL, Université de Montpellier II - IPM, CNRS, Place Eugéne Bataillon, 34095 Montpellier Cedex 05, France
3 - Centre de Recherche Astrophysique de Lyon, UMR 5574: CNRS, Université de Lyon, École Normale Supérieure de Lyon, 46 allée d'Italie, 69364 Lyon Cedex 07, France
4 - Department of Physics and Astronomy, Division of Astronomy and Space Physics, Uppsala University, Box 515, 75120 Uppsala, Sweden

Received 3 February 2009 / Accepted 14 July 2009

Abstract
Context. It has been proposed that convection in red supergiant (RSG) stars produces large-scale granules causing observable surface inhomogeneities. This convection is also extremely vigorous and is suspected to be one of the main causes of mass-loss in RSGs. It should thus be understood in detail. Evidence has accumulated of asymmetries in the photospheres of RSGs, but detailed studies of granulation are still lacking. Interferometric observations provide an innovative way of addressing this question, but they are still often interpreted using smooth symmetrical limb-darkened intensity distributions, or simple, spotted, ad hoc models.
Aims. We explore the impact of the granulation on visibility curves and closure phases using the radiative transfer code OPTIM3D. We simultaneously assess how 3D simulations of convection in RSG with CO⁵BOLD can be tested by comparing with these observations.
Methods. We use 3D radiative hydrodynamical (RHD) simulations of convection to compute intensity maps at various wavelengths and time, from which we derive interferometric visibility amplitudes and phases. We study their behaviour with time, position angle, and wavelength, and compare them to observations of the RSG $\alpha$ Ori.
Results. We provide average limb-darkening coefficients for RSGs. We describe the prospects for the detection and characterization of granulation (i.e., contrast, size) on RSGs. We demonstrate that our RHD simulations provide an excellent fit to existing interferometric observations of $\alpha$ Ori, in contrast to limb darkened disks. This confirms the existence of large convective cells on the surface of Betelgeuse.

Key words: stars: supergiants - stars: atmospheres - hydrodynamics - radiative transfer - techniques: interferometric

1 Introduction

Massive stars with masses between roughly 10 and 25 $M_{\odot}$ spend some time as red supergiant (RSG) stars, which represent the largest stars in the universe. They have effective temperatures, $T_{\rm eff}$ , ranging from 3450 to 4100 K, luminosities of between 20 000 and 300 000 $L_\odot$ and radii up to 1500 $R_\odot$ (Levesque et al. 2005). Their luminosities imply that they are among the brightest stars, visible to very large distances. However, a number of open issues remain. They shed large amounts of mass back to the interstellar medium, but their mass-loss mechanism is unidentified, although Alfvén and acoustic waves have been proposed (Pijpers & Hearn 1989; Cuntz 1997; Hartmann & Avrett 1984), as well as acoustic waves and radiation pressure on molecules (Josselin & Plez 2007). Their chemical composition is largely unknown, despite the work of e.g. Carr et al. (2000), and Cunha et al. (2007), because of difficulties in analysing their spectra, which contain broad, asymmetric lines that according to a convection pattern of large granules and (super-)sonic velocities (Gray 2008; Josselin & Plez 2007). Their $T_{\rm eff}$ -scale has been revised both at solar and Magellanic Cloud metallicities using 1D hydrostatic models (Massey et al. 2007; Levesque et al. 2005,2006,2007). Although these MARCS models (Gustafsson et al. 2008) provide an accurate fit to the optical spectra allowing the derivation of $T_{\rm eff}$ and reddening, problems remain. There is a blue-UV excess in many of the observed spectra, which may be related to either scattering by circumstellar dust or to an insufficiency in the models. There is also some disagreement in the IR colours, which may be related to atmospheric temperature inhomogeneities that are characteristic of convection (Levesque et al. 2006).

Hydrodynamical modeling of convection in RSGs has lagged behind that of solar type stars because of the need to include the entire star in the simulation box. Freytag et al. (2002) managed to develop numerical simulations of a typical RSG. We thus attempted to improve our understanding and description of RSGs using detailed numerical simulations and a large set of observational material.

This paper is the first in a series exploring the granulation pattern of RSGs and its impact on interferometric observations.

2 3D radiative transfer in a radiative hydrodynamical simulation

2.1 3D hydrodynamical simulations with CO $^{\sf\it 5}$ BOLD

The numerical simulations employed in this work were developed using CO⁵BOLD (Freytag & Höfner 2008; Freytag et al. 2002) in the star-in-a-box configuration: the computational domain is a cube, and the grid is equidistant in all directions. All six faces of the cube use the same open boundary conditions for material flows and emergent radiation. In addition, there is an ``inner boundary condition'': in a small spherical region in the centre of the cube, a source term to the internal energy provides the stellar luminosity and a drag force brakes dipolar flows through it. Otherwise, the hydrodynamics and the radiative transfer scheme are unaffected by the core and an integration can be completed without problem. Radiation transport is strictly in LTE. The grey Rosseland mean opacity is a function of both gas pressure and temperature. The necessary values are derived by interpolating in a table derived for temperatures around 12 000 K from high-temperature OPAL data (Iglesias et al. 1992) and low-temperature PHOENIX data (Hauschildt et al. 1997) by Hans-Günter Ludwig. Some more technical information can be found in Freytag & Höfner (2008), the CO⁵BOLD Online User Manual, and in a forthcoming paper by Freytag (2009).

The 12 $M_{\odot}$ model that we use in this paper, called st35gm03n07, was produced by of intensive calculations corresponding to about 7.5 years of simulated stellar time. It has a numerical resolution of 8.6 $R_\odot$ within a cube of 235³ grid points. The model parameters are a luminosity of $L=93~000 \pm 1300~{L}_{\odot}$ , an effective temperature of $T_{\rm {eff}}=3490\pm 13$ K, a radius of $R=832\pm0.7~{R}_{\odot}$ , and correspondingly a surface gravity of ${\rm log}~g=-0.337\pm0.001$ . These values are averages over both spherical shells and time (over the past year), and the errors are one standard deviation fluctuations with respect to the average over time. We define the stellar radius, R, and the effective temperature, $T_{\rm eff}$ , as follows. First, we compute the average temperature and luminosity over spherical shells, T(r), and L(r). We then search the radius R for which ${L(R)}/(4\pi R^2)=\sigma T^4(R)$ , where $\sigma$ is the Stefan-Boltzmann constant. The effective temperature is then $T_{\rm {eff}}=T(R)$ . Figure 1 shows the value of the radius, temperature, and luminosity over the past 3.5 years. The radius drifts by about $-0.5\%$ per year and seems to have stabilised to $R=832~{R}_{\odot}$ in the last year. Over the entire sequence $T_{\rm eff}$ fluctuates by $\pm 1\%$ , and has a constant average value. The luminosity fluctuations are of the order of $\pm 4\%$ , reflecting the temperature variations, with a decrease of about $1\%$ per year in the first few years, reflecting the radius decrease. These drifts indicate that the simulation has not completely converged in the first few years. In this work, we consider the entire 3.5 year sequence, despite the small radius drift, to obtain better statistics. The preceding 4 years of the simulation are not considered here, since they show larger drifts. The interferometric observables derived in this work are insensitive to the drift in the parameters.

This is our most successful RHD simulation so far because it has stellar parameters closest to those of real RSGs (e.g., 3650 K for $\alpha$ Ori, Levesque et al. 2005). We are developing additional simulations with different stellar parameters and will present the analysis of these simulations in a forthcoming paper.

$\begin{figure} \par\includegraphics[width=9cm,clip]{178fig1.eps}\vspace{2mm} \i... ...ig2.eps}\vspace{2mm} \includegraphics[width=9cm,clip]{178fig3.eps} \end{figure}$	Figure 1: Radius (Panel A), luminosity (Panel B) and temperature (Panel C) as a function of time for the simulation used in this work. The radius is fitted with the law: $R(t)=936.9 - 4.3\times t$ for $t\leq 23.8$ yrs and $R=832~{R}_{\odot}$ for t>23.8 yrs.
Open with DEXTER

2.2 Radiative transfer code: OPTIM3D

We developed a 3D pure LTE radiative transfer code, OPTIM3D, to generate synthetic spectra and intensity maps from snapshots of the 3D hydrodynamical simulations, taking into account the Doppler shifts caused by the convective motions. The radiation transfer is calculated in detail using pre-tabulated extinction coefficients generated with the MARCS code (Gustafsson et al. 2008). These tables are functions of temperature, density, and wavelength, and were computed by assuming the solar composition of Asplund et al. (2006). The tables include the same extensive atomic and molecular data as the MARCS models. They were constructed with no micro-turbulence broadening and the temperature and density distribution was optimized to cover the values encountered in the outer layers of the RHD simulations. The wavelength resolution is $R=\lambda/\Delta\lambda=500~000$ and we checked that this resolution is sufficient to ensure an accurate calculation of broadened line profiles of RSGs even after interpolation of the opacity at the Doppler shifted wavelengths.

The monochromatic intensity emerging towards the observer at a given position on the simulation can be computed by integrating the source function along a ray perpendicular to a face of the cube, at that position. In LTE it is given by:

$\begin{displaymath}{I}_{\lambda}\left(0\right)=\int^{\tau_\lambda}_0 S_{\lambda}\left({t_\lambda}\right){\rm e}^{-{t_\lambda}}{\rm d}{t_\lambda}, \end{displaymath}$

(1)

where $I_{\lambda}$ is the intensity, $t_{\lambda}$ is the optical depth along the ray increasing inwards, $\tau_\lambda$ is the maximum optical depth reached along the line-of-sight, and $S_{\lambda}=B_{\lambda}\left(T\right)$ , the Planck function at the temperature T, is the source function. A Gauss-Laguerre quadrature of order n can be performed to evaluate the integral, Eq. (1), when $\tau_\lambda\rightarrow\infty$ . This method is much faster than a detailed integration along the discretized ray, because it uses only the value of the source function at n depth points weighted with n predetermined weights. This method is reliable as long as the source function is sufficiently smooth along the optical depth scale, and is well known at the quadrature points. This is not always the case in our simulations, where the optical depth scale may jump by large amounts between 2 successive cells, e.g., from $\tau=1$ to $\tau=300$ for extreme cases. The source function must then be interpolated at intermediate optical depths, and the result is largely dependent on the way this interpolation is performed. We note however that this is also the case for a detailed integration, where too large jumps in the source function or optical depth scale will cause uncertainties in the resulting intensity. We checked for differences between a Gauss-Laguerre quadrature and a detailed summation of the contributions from all cells, with different kinds of interpolations, and found differences in intensities emerging from a single ray of less than $10\%$ on average, although there were some extreme cases that reached more than $100\%$ because of a particularly significant illconditioning of the source function. The average differences were within an acceptable range, we therefore rely on the Gauss-Laguerre quadrature, with a linear interpolation of the source function on the logarithmic $\tau$ -scale. The quadrature points and weights that we use are listed in Table 1 (Abramowitz & Stegun 1972). For the rays where the optical depth does not reach sufficiently high values, we complete a detailed summation of the contribution from all cells along the ray.

In practice, once the input simulation is read, OPTIM3D interpolates the opacity tables in temperature and logarithmic density for all the simulation grid points using a bi-linear interpolation. The interpolation coefficients are computed only once, and stored. Bi-linear interpolation was chosen instead of spline interpolation because: (i) spline interpolation is significantly more time consuming, and (ii) comparisons with other codes do not indicate that great improvements can be achieved using splines (see below). The logarithmic extinction coefficient is then linearly interpolated at each Doppler-shifted wavelength in each cell along the ray, and the optical depth scale along the ray is calculated. Equation (1) is then integrated, inferring the intensity emerging towards the observer at that wavelength and position. This calculation is performed for every line-of-sight that is perpendicular to the face of the computational box, for all the required wavelengths.

Table 1: Gauss-Laguerre quadrature weights for n=10.

Comparisons with existing codes were carried out. The spectral synthesis code Turbospectrum (Plez et al. 1993; Alvarez & Plez 1998, and further improvements by Plez) was used with one-dimensional MARCS models, where the source function was very well sampled on the $\tau$ -scale. OPTIM3D computations completed with bi-linear interpolation deviate by less than $5\%$ , and the deviation decreases to 0.2 $\%$ with spline interpolation. We also checked OPTIM3D against Linfor3D (Cayrel et al. 2007 for the Non-LTE version, and for the LTE version) using 3D CO⁵BOLD local models. We compared synthetic spectra computed for three artificial iron lines (of increasing strength) centered on a laboratory wavelength of 5500 Å using the same abundances. The discrepancy between the results of the codes was less than 3 $\%$ and was even less than 0.2 $\%$ when a spline interpolation of the opacity tables was used in OPTIM3D (with a significant increase in the CPU time). Finally, comparisons were made with the spectral line formation code used by, e.g., Asplund (2000) for 3D local convection simulations carried out with the code by Stein & Nordlund (1998) for giant stars (Collet et al. 2007). The tests were carried out on the [OI] line at 6300.3 Å and various Fe I and Fe II lines around 5000 Å. The discrepancies between the resulting synthetic spectra were less than 2 $\%$ , and became even less than 0.6 $\%$ when a spline interpolation of the opacity tables is used in OPTIM3D. Thus, the interpolation is the main source of error. In conclusion, if only a few lines are computed for, e.g., accurate abundance determinations, Linfor3D or the Asplund code are superior because they mostly avoid interpolations into opacity tables. The code used by Asplund performs bi-cubic interpolations of both the continuum opacity and the individual number densities, whereas we interpolate the total opacity from all lines contributing at a given wavelength, which is in principle less accurate. Therefore, when a large wavelength range must be calculated by taking into account many molecular and atomic lines simultaneously, OPTIM3D is superior, faster choice, which still provides results with an accuracy of a few percent.

3 Simulated images in the H and K bands: giant convective cells

We analyze the properties of the simulations in the H and K bands, where many interferometric observations were completed, and existing interferometers, e.g., VLTI/AMBER, operate routinely. We calculated intensity maps for a series of snapshots about 23 days apart covering 3.5 years of the model described above. We used the transmission curve of the four K band filters mounted on FLUOR (Fiber Linked Unit for Optical Recombination; Coude Du Foresto et al. 1998), and the H band filter mounted on IONIC (Integrated Optics Near-infrared Interferometric Camera; Berger et al. 2003) at the IOTA interferometer (Traub et al. 2003). The K band filters (Fig. 2) are: K203 (with a central wavelength of 2.03 $\mu$ m), K215 (2.15 $\mu$ m), K222 (2.2 $\mu$ m), and K239 (2.39 $\mu$ m). The H band filter has a central wavelength of 1.64 $\mu$ m (Fig. 3). The resulting intensities reported in this work are normalized to the filter transmission as $\frac{\int I_{\lambda} T\left(\lambda\right){\rm d}\lambda}{\int T\left(\lambda\right){\rm d}\lambda}$ where $I_{\lambda}$ is the intensity and $T\left(\lambda\right)$ is the transmission curve of the filter at a certain wavelength. The intensity maps are showed after applying a median [ $3\times3$ ] smoothing (see Sect. 5).

$\begin{figure} \par\includegraphics[width=9cm,clip]{178fig4.ps} \end{figure}$

Figure 2:

The transmission curves of the 4 narrow band filters mounted on the FLUOR instrument at IOTA together with the K band synthetic spectrum of a snapshot of the simulation and the corresponding continuum (bottom black curve). From the top, the spectra computed with only H₂O (red), only CO (green), and only CN (blue) are shown with an offset of, respectively, $0.9\times 10^{33}$ , $0.6\times 10^{33}$ , and $0.3\times 10^{33}$ erg cm^-2 s^-1 Å^-1.

Open with DEXTER

$\begin{figure} \par\includegraphics[width=9cm,clip]{178fig5.ps} \end{figure}$	Figure 3: The transmission curve of the filter mounted on IONIC at IOTA together with the H band synthetic spectra computed as in Fig. 2. From the top, the offset of the spectra is $3\times 10^{33}$ , $2\times 10^{33}$ , and $1\times 10^{33}$ erg cm^-2 s^-1 Å^-1.
Open with DEXTER

It can be seen from our simulations (see Fig. 4) that the surface of the stellar model is covered by a few large convective cells of size between about 400 and 500 $R_\odot$ that evolve on a timescale of years. These cells have strong downdrafts that can penetrate to the stellar core (Freytag et al. 2002, 2009, in preparation). Close to the surface, there are short-lived (a few months to one year) small-scale (50 to 100 $R_\odot$ ) granules (bottom panels of Fig. 4). Freytag et al. (1997) found a relation between the mean horizontal size of convective granules $x_{\rm {gran}}$ and the atmospheric pressure scale-height defined as $H_{\rm p0}=\frac{kT_{\rm {eff}}}{g\mu m_{\rm {H}}}$ for GK dwarfs and subgiants. It is unclear if this relation can be extrapolated to 3D simulations of RSGs. Using it we find that $x_{\rm {gran}}/R_{\star}=10\times H_{\rm p0}/R_{\star}=0.1$ , for parameters appropriate to a RSG atmosphere dominated by gas pressure. Obviously, this leads to a size much smaller than can be seen in Fig. 4. Freytag et al. (1997) found that a value of $x_{\rm {gran}}/H_{\rm p0} = 10$ would fit 2D simulations for GK dwarfs and subgiants, but they also showed that A-type and F-type stars lie above the curve indicating that they have larger granules. These stars have high turbulent pressure that may dominate over the gas pressure in turbulent convective layers. Following Gustafsson et al. (2008), we write $P_{\rm {turb}}=\beta \rho v_{\rm {turb}}^2$ , where $v_{\rm {turb}}$ is the turbulent velocity, $\rho$ is the gas density, and $\beta$ is a parameter close to one, whose value depends on the anisotropy of the velocity field. A clearer way to express $H_{\rm p0}$ is thus $H_{\rm p0}= \frac{kT_{\rm {eff}}}{g\mu m_{\rm {H}}}\left(1+\beta\gamma(\frac{v_{\rm {turb}}}{c_{\rm s}})^2\right)$ , where $\gamma$ is the adiabatic exponent, and $c_{\rm s}$ the sound speed. If $v_{\rm {turb}}$ is only a factor of 2 larger than $c_{\rm s}$ , $H_{\rm p0}$ is increased by a factor of about 5. This is the case for our RSG simulation, where $P_{\rm {turb}}/P_{\rm {gas}}\sim2$ at the surface, R_*, as determined in Sect. 2.1. This implies then that $x_{\rm {gran}}/R_{\star}=0.5$ , after extrapolating Freytag et al. (1997) formula. This is more consistent with the properties of the large granules visible on intensity maps in Fig. 4. Additional mechanisms might influence the size of the granules: (i) in RSGs, most of the downdrafts will not grow fast enough to reach any significant depth before they are swept into the existing deep and strong downdrafts enhancing the strength of neighboring downdrafts; (ii) radiative effects and smoothing of small fluctuations can cause an enhancement in the growth time for small downdrafts, while the granule crossing time is short because of the large horizontal velocities; (iii) sphericity effects, see for example (Steffen & Freytag 2007; Freytag & Ludwig 2007); (iv) Freytag et al. (1997) use both the effective temperature and the pressure scale height at the bottom of the photosphere as a reference, although layers below the photosphere also matter; and (v) numerical resolution (or lack of it) could also have an effect.

$\begin{figure} \par\includegraphics[width=18.5cm,clip]{178fig6.eps}\par\includegraphics[width=18.5cm,clip]{178fig7.eps} \end{figure}$	Figure 4: Top 6 panels: maps of the intensity in the IONIC filter (linear scale with a range of [0; $2.5\times 10^5$ ] erg cm^-2 s^-1 Å^-1). The different panels correspond to snapshots separated by 230 days ( $\sim$ 3.5 years covered). Bottom 6 panels: successive snapshots separated by 23 days ( $\sim$ 140 days covered).
Open with DEXTER

4 Intensity profiles

The simulated RSG atmospheres appear very irregular both in their structure and dynamics. The surface inhomogeneities and their temporal evolution induce strong variations in the emerging spectra, and intensity profiles. In this section, we analyse the average centre-to-limb intensity profiles, and their time variations.

4.1 Surface inhomogeneities and temporal evolution

The top left panel of Fig. 5 shows a three-dimensional image representation of the intensity emerging from one face of a snapshot of the simulation in the H band. The K band appearance is similar. This image shows very sharp intensity peaks, two to three pixels wide. This is also evident in the top right panel of the figure as small bright (up to $40\%$ brighter than the surrounding points) patches. These patches reflect the ill-conditioning of the source function, because of a lack of spatial resolution around $\tau_\lambda=1$ along some lines of sight, where the source function may have a significant discontinuity (see Sect. 2.2). Attempts were made to solve this problem by interpolating the source function and opacity inside CO⁵BOLD, although this led to numerical instabilities. The only possible solution is to increase the number of grid points, which would require larger and faster computers.

Radial intensity profiles within a given snapshot exhibit large variations with position angle in their radial extension of about $10\%$ (see bottom left panel of Fig. 5). The variation with time in the intensity profiles are of the same order of magnitude (10 $\%$ , see bottom right panel of figure).

$\begin{figure} \par\hspace*{1cm}\includegraphics[width=7cm,clip]{178fig8.ps}\hsp... ...ig10.ps}\hspace*{3mm} \includegraphics[width=8cm,clip]{178fig11.ps} \end{figure}$

Figure 5:

Top left panel: three-dimensional image of a snapshot from Fig. 4. Top right panel: intensity map of the same snapshot represented using the histogram equalization algorithm to underline the thin bright patches because of the undersampling of the $\tau$ scale. Bottom left panel: intensity profiles for three position angles of the same snapshot. The numerical box edge is at impact parameter $r/R_{\star }\sim 1.3$ . The intensity is normalized to the intensity at disk center. Bottom right panel: radially averaged intensity profiles for all the snapshots of the simulation (grey); one snapshot of the simulation is emphasized with a solid black line. The intensity is normalized to the area subtended by the curve, area $=\int_0^{1.3} I\left(r/R_{\star}\right) {\rm d}r/R_{\star}$ .

Open with DEXTER

4.2 The limb darkening law

Despite the large azimuthal variations in the intensity profiles, and their temporal variations, it is interesting to derive radially averaged intensity profiles for each snapshot. These may be be used, e.g., as a first approximation when interpreting interferometric observations, instead of limb-darkening (LD) laws computed from hydrostatic models (Claret 2000). The bottom right panel of Fig. 5 shows all the radially averaged intensity profiles obtained from the simulation.

We use a LD law of the form

$\begin{displaymath} \frac{I(\mu)}{I(1)}=\sum_{k=0}^3 a_k\left(1-\mu\right)^k, \end{displaymath}$

(2)

where $I(\mu)$ is the intensity, a_k are the LD coefficients, and $\mu=\cos\theta$ , where $\theta$ is the angle between the line of sight and the radial direction. $\mu$ is related to the impact parameter $r/R_{\star}$ by means of $r/R_{\star}=\sqrt{1-\mu^2}$ , where $R_{\star}$ is the stellar radius determined as Sect. 2.1. The average intensity profiles were constructed by using rings regularly spaced in $\mu$ for $\mu\le1$ (i.e., $r/R_{\star}\le1$ ), and adding a few points for $r>R_{\star}$ up to the numerical box limit. The standard deviation of the average intensity, $\sigma_{I\left(\mu\right)}$ , was computed within each ring. There is a small tail at $r>R_{\star}$ , which provides a minor contribution to the total flux (less than 1 $\%$ , see bottom right panel of Fig. 5), and cannot be fitted with Eq. (2). We fitted the radially averaged profiles of all the snapshots of the simulation (57 profiles, 23 days apart covering 3.5 years). The fit was weighted by $1/\sigma_{I\left(\mu\right)}$ to decrease the importance of central points with poor statistics. The fit was first made to $I(\mu)/I_{\rm {norm}}$ , where $I_{\rm {norm}}$ the intensity normalized to the area subtended by the curve: $I_{\rm {norm}}=\frac{I\left(r/R_\star\right)}{\int_0^{1.3} I\left(r/R_\star\right) {\rm d}r/R_\star}$ . This was done to diminish the impact of intensity fluctuations between snapshots on the fitting coefficients.

Table 2: Time-averaged limb-darkening coefficients for the RHD simulation.

In Table 2, we present the values of the four LD coefficients averaged over all 3.5 years, renormalized to the disk center, for both the IONIC H-band filter, and the K222 filter (because the sensitivity of the FLUOR instrument is always superior in the continuum than in the molecular bands (Perrin et al. 2004b), and it samples the maximum transmission region of the K band). Figure 6 shows an example of a LD fit. The intensity profiles for different position angles in the same snapshot are very different (see Fig. 5, bottom left panel), the fitting coefficients having large scatter. The time averaged LD fits provides, however, an indication of the shape of the intensity profile in the H and K bands (we note that they are very similar). They differ of course very much from simple first order LD laws. They also differ from LD laws calculated by Claret (2000) for parameters appropriate to RSGs (see Fig. 7). When comparing with observations of RSGs, we recommend the use of our fits. Ideally, one should use single snapshots, as we do below in our analysis of Betelgeuse, because they may deviate from the average LD fit by large amounts (see Table 2).

$\begin{figure} \par\includegraphics[width=8.5cm,clip]{178fig12.ps} \end{figure}$	Figure 6: Example of a LD fit (dashed line) using the LD law described in the text for the radially averaged intensity profile (solid line) emphasized in Fig. 5 (bottom right panel). The intensity is normalized to the area subtended by the curve. This best-fit model has a $\chi ^2=0.02$ .
Open with DEXTER

$\begin{figure} \par\includegraphics[width=8.5cm,clip]{178fig13.ps} \end{figure}$	Figure 7: The time-averaged H-band radial intensity profile of our simulation (solid line), and the fit of Table 2 (dashed line). A fully LD (dash-dotted line), a partially LD (triple dot-dashed line), and a LD fit from Claret (2000, dotted line) for comparable RSG parameters are plotted for comparison.
Open with DEXTER

5 Visibility curves and phases

The granulation pattern has a significant impact on the interferometric visibility curves and phases. We attempted here to derive their characteristic signature.

5.1 Computation

We computed visibilities and phases using the IDL data visualization and analysis platform. For each image, a discrete Fourier transform was calculated. To reduce the problems caused by the finite size of the object and avoid edge effects, the resolution in the Fourier plane was increased by resampling the input $235\times 235$ pixel image to a grid size of $2048\times 2048$ pixels. The visibility V was defined as the modulus |z| of the complex Fourier transform, $z = x + {\rm i}y$ , normalized to the modulus at the origin of the frequency plane, |z₀|, with the phase $\tan\theta = \Im(z)/\Re(z)$ . When dealing with observations, the natural spatial frequency ( $\nu$ ) unit is arcsec^-1. Since we study theoretical models, we instead use $R_\odot^{-1}$ units. The conversion factor between these is

$\begin{displaymath} \nu~[{\rm arcsec}^{-1}]=\nu~[R_\odot^{-1}]\cdot d~[{\rm pc}]\cdot214.9, \end{displaymath}$

(3)

where 214.9 is the astronomical unit expressed in solar radii, and d is the distance of the observed star. The relation between the baseline, B (in m), of an interferometer, and the spatial frequency $\nu$ (in arcsec^-1) observed at a wavelength $\lambda$ (in $\mu{\rm m}$ ) is $\nu=B/\lambda/0.206265$ . Since our calculated images are affected by source function discontinuities, we investigated how the visibility curves are affected by the resulting bright spikes (Fig. 5). In Fig. 8, we compare the visibility curves computed for one projected baseline from the raw image, and after applying a median $[3\times 3]$ smoothing that effectively erases the artifacts. The visibility curves are affected by these artifacts mostly for frequencies higher than 0.03 $R^{-1}_\odot$ (corresponding to 33 $R_\odot$ , i.e., $\sim$ 4 pixels). We can therefore apply this cosmetic median filter, because it will not affect the visibilities at lower frequencies, which are the only ones observable in practice with modern interferometers.

$\begin{figure} \par\includegraphics[width=8.5cm,clip]{178fig14.eps} \end{figure}$	Figure 8: The solid line is the visibility curve for the IONIC filter intensity map of Fig. 5 (top right). The dotted line is computed for the same map after applying a $[3\times 3]$ median smoothing.
Open with DEXTER

We now study the first few lobes of the visibility curves of our simulations, and how they are affected by asymmetries and surface structure.

5.2 The first lobe

The first lobe of the visibility curve is mostly sensitive to the radial extension of the observed source. Figure 9 (bottom panel) shows the visibility curves computed for 36 different angles from the intensity map of one snapshot in the IONIC filter (top panel). A dispersion of the visibility curves (thin grey lines in Fig. 9) is noticeable. This behavior is similar for all the snapshots. These synthetic visibilities have been compared with a uniform disk (UD) model (solid line in Fig. 9), and with limb-darkened (LD) models. We use both a fully limb-darkened disk ( $I_{\mu}/I_1=\mu$ , dotted-dashed line in the figure), and a partially limb-darkened model with a₁=-0.5 ( $I_{\mu}/I_1=0.5+0.5\cdot\mu$ , dashed line in the figure). The radius determined by fitting a UD disk model to the computed visibilities ranges from 794 to 845 $R_\odot$ for the 36 angles, and is up to 5 $\%$ smaller than $R_{\star}=836.5~R_\odot$ , the radius of the simulation determined as described in Sect. 2.1. The partially-, and fully-darkened models are respectively $\sim$ 2%, and only $\sim$ 1% smaller than $R_{\star}$ . In Fig. 9, we also show the visibility amplitude resulting from our average LD fit of Table 2. The resulting diameter is then 842 $R_\odot$ , very close to the simulation radius. We note, that Nardetto et al. (2006) also found that the UD radius is about 4 to 5 $\%$ smaller than the photospheric radius in their simulation of Cepheids, and that the LD radius is much closer to the radius of their simulation. Stellar diameters determined with UD or LD fits to the observed first visibility lobe of RSGs will be affected by these systematic errors. As shown below, observations of higher spatial frequencies will greatly improve the knowledge of both limb-darkening and asymmetries, thus helping the radius to be more tightly constrained.

$\begin{figure} \par\includegraphics[width=8.5cm,clip]{178fig15.eps}\vspace*{2mm} \par\includegraphics[width=8.5cm,clip]{178fig16.eps} \end{figure}$

Figure 9:

Top panel: intensity map in the IONIC filter (the range is [0; $2.5\times 10^5$ ] erg cm^-2 s^-1 Å^-1). Bottom panel: visibility curves from the above snapshot computed for 36 different angles, 5 $^\circ$ apart (thin grey lines). Note the logarithm visibility scale. The solid black curve is a UD model, with a radius of 810 $R_\odot$ . The dashed black line is a partially LD disk with a radius of 822 $R_\odot$ . The dot-dashed line is a fully LD disk with a radius of 830 $R_\odot$ . the triple-dot-dashed line is our average LD law (cf. Table 2) for a radius of 842 $R_\odot$ . The stellar parameters of this snapshot are: $L = 98~400~{L}_\odot$ , $R_* = 836.5~{R}_\odot$ , $T_{\rm eff} = 3534$ K, and $\log g = -0.34$ .

Open with DEXTER

It is interesting to compare the angular and temporal visibility fluctuations at one standard deviation, defined as $F=\sigma_{\rm {vis}}/\rm {vis}$ , where vis is the visibility: (i) the temporal evolution, fixing one angle and following the RHD simulation for 3.5 years with a time-step of $\sim$ 23 days; (ii) the angular evolution, considering a single snapshot and computing the visibilities for 36 different angles 5 $^\circ$ apart. In the first lobe, Fig. 10 shows that temporal and angular fluctuations have the same order of magnitude. The fluctuations are less than 1% at frequency $\sim$ 0.00040 $R^{-1}_\odot$ (at this frequency, the visibility is greater than $50\%$ ), they are $\sim$ 3% at frequency 0.00057 $R^{-1}_\odot$ , and are close to $\sim$ 10% at 0.00069 $R^{-1}_\odot$ .

$\begin{figure} \par\includegraphics[width=8.5cm,clip]{178fig17.eps} \end{figure}$	Figure 10: Standard deviation of the visibility normalized to the visibility in the first lobe. The solid line indicates the temporal fluctuations for one fixed angle over 3.5 years. The trend is similar for all other angles. The dashed curve corresponds to the angular fluctuations of the snapshot in Fig. 9.
Open with DEXTER

$\begin{figure} \par\includegraphics[width=8.5cm,clip]{178fig18.eps} \end{figure}$	Figure 11: Same as Fig. 9 for the second, third, and fourth lobes. In addition, the dotted line is the visibility curve for the particular azimuth parallel to the x-axis of the IONIC intensity map of Fig. 9.
Open with DEXTER

5.3 The second, third, and fourth lobes: signature of convection

As in Sect. 5.2, we analyze the angular and temporal visibility fluctuations at one sigma with respect to the average value in the H band (IONIC filter). Figure 11 shows an enlargement of the second, third, and fourth lobes of the visibility curves computed for different position angles. The dispersion clearly increases with spatial frequency, and visibilities deviate significantly from the UD or LD cases, because of the small scale structure in the model stellar disk. The same is true for temporal fluctuations of the visibility at a given position angle. Figure 12 shows the temporal fluctuations of the visibilities for one fixed position angle, as well as the angular fluctuations for the snapshot of Fig. 9. As for the first lobe, there is no clear distinction between angular and temporal fluctuations. Relative fluctuations are of course large around the minima of visibility, where visibilities are also more difficult to measure. However, with the precision of current interferometers (e.g., 1% for visibilities of $\sim$ 5-10% for VLTI-AMBER), it should be possible to characterize the granulation pattern of RSGs. This requires observing the third and the fourth lobes and not limiting the observation at the first and second lobes, which only provide information about the radius and LD. The signal to be expected in these lobes is higher than the UD or LD predictions (see dashed line in Fig. 11). Efforts should therefore be directed toward observing at these frequencies.

$\begin{figure} \par\includegraphics[width=8.5cm,clip]{178fig19.eps} \end{figure}$	Figure 12: Same as Fig. 10 for the second, third, and fourth lobe.
Open with DEXTER

It may however be that approximations in our model (e.g., limited spatial resolution, grey radiative transfer) significantly affect the intensity contrast of the granulation. The radiation transfer in our RHD models indeed uses a frequency-independent grey treatment to speed up the calculations. This approximation leads to errors in the mean temperature structure within the optically thin layers that are difficult to quantify. The implementation of non-grey opacities can decrease the temperature fluctuations compared to the grey case (e.g., Ludwig et al. 1994, for local RHD models). As a consequence, the intensity contrast will decrease, reducing the visibility fluctuations. To investigate its impact on visibilities, we artificially decrease the intensity contrast in one of our images, and use the snapshot with its nominal intensities as reference. We first fit a LD law (as in Sect. 4) to the radially average intensity profile. After subtracting this average profile from the intensity map, we obtain the fluctuations caused by granulation, and measure the contrast $C_{\rm {ref}}$ = $\frac{I_{\rm {max}}-I_{\rm {min}}}{I_{\rm {max}}+I_{\rm {min}}}$ . It is then easy to scale this contrast before again adding the LD profile, to produce a reduced contrast image. An example of the resulting intensities is shown in Fig. 13 (top row). With a contrast of only 1% with respect to the nominal value, small surface structures are hardly visible. As previously, we determine $\sigma_{\rm {vis}}/\rm {vis}$ for all the images with various contrasts, around the top of the second ( $\sim$ 0.0010 $R^{-1}_\odot$ ), third ( $\sim$ 0.0016 $R^{-1}_\odot$ ), and fourth lobes ( $\sim$ 0.0022 $R^{-1}_\odot$ ). The bottom left panel of Fig. 13 shows that when the contrast is reduced, and the surface structures fade, the resulting visibility fluctuations similarly decrease in all the lobes (almost in proportion to the intensity contrast decrease). Reducing the contrast of course brings the visibility curves towards the visibility of the LD profile (Fig. 13, bottom right panel). This proportionality can be used to determine the granulation contrast from observations of the visibility fluctuations with time or position angle.

$\begin{figure} \par {\hspace*{1cm}\includegraphics[width=7cm,clip]{178fig20.eps}... ...g22.ps}\hspace*{3mm} \includegraphics[width=8cm,clip]{178fig23.eps} \end{figure}$

Figure 13:

Top left panel: three-dimensional image of a snapshot with nominal intensities. Top right panel: same snapshot with a feature contrast reduced to 1%. Bottom left panel: standard deviation of the visibility curves at 36 angles 5 $^\circ$ ) apart, for decreasing feature contrast. Solid line is for the top of the second lobe ( ${\sim}0.0010~{R}^{-1}_\odot$ ), dashed line is for the top of the third lobe ( ${\sim}0.0016~{R}^{-1}_\odot$ ), and dotted line is for the top of the fourth lobe ( ${\sim}0.0022~{R}^{-1}_\odot$ ). Bottom right panel: visibility in the second, third, and fourth lobes for one particular position angle. The dot-dashed line shows the original simulation contrast. The dashed line, and the dotted line show the visibility with a feature contrast reduced to 50 $\%$ , and 1% respectively. The solid line is the fitted LD profile computed for this snapshot.

Open with DEXTER

5.4 Importance of spectral resolution in interferometry: the H and K bands

Interferometric observations completed with a broad-band filter combine information about the lines and continuum. Higher spectral resolution allows us to recover for more information, about both visibility moduli and phases. The VLTI-AMBER interferometer provides spectral resolutions of R=35, 1500, and 12 000. To illustrates the differences between these resolutions, we compute intensity maps around the CO first overtone line at 23041.75 Å ( $\log(gf)=-5.527$ and $\chi_{\rm ex}=0.180$ eV) for the three resolutions (Fig. 14). The resulting images are shown in the central row of the figure, and the spectrum in the top row. The contrast, defined as in the previous section, is similar for the low and medium resolution images but is $\sim$ 30% lower in the high resolution image, at the CO line wavelength. Large fluctuations are visible in all lobes, but they are smaller than those seen in the H band IONIC filter (see Fig. 12). The visibility fluctuations in the high spectral resolution image of the CO line are larger (second, third, and fourth lobes; see dotted line in bottom panel), in spite of a lower intensity contrast, presumably because of the darkening of large patches of the simulated stellar surface.

$\begin{figure} \par\includegraphics[width=18cm,clip]{178fig24.eps}\vspace*{3mm} ... ...4cm}\includegraphics[width=10cm,clip]{178fig26.eps}\hspace*{4cm}}\end{figure}$

Figure 14:

Top row: synthetic spectrum centered on the CO first overtone band head. The left panel shows the range of wavelengths spanned by one resolution element at the VLTI-AMBER low spectral resolution of 35. The right panel shows the same for the VLTI-AMBER medium spectral resolution of 1500, and the high spectral resolution of 12 000 (thick mark). Central row: intensity maps for those three spectral resolution elements. The intensity range is [0;10⁵] erg cm^-2 s^-1 Å^-1. Bottom row: standard deviation of the visibility in the second, third, and fourth lobes. Solid and dashed lines correspond to low and medium resolution, respectively, and the dotted line to high resolution.

Open with DEXTER

We also computed wavelength-dependent visibility curves in the H band for the high and medium VLTI-AMBER resolutions. Figure 15 displays a three-dimensional view of the visibility curves with a resolution of 12 000 and 1500 (top panels). The simulated star was scaled to an apparent diameter of $\sim$ 43.6 mas (the observed diameter of $\alpha$ Ori, Perrin et al. 2004a). The displacement of the zeropoints with wavelength is easily seen, as well as the amplitude variations in the higher frequency lobes.

To simulate differential observations with an interferometer at medium and high spectral resolution, in Fig. 15, we also show the variation in the visibility modulus with wavelength for a fixed baseline (15 m, i.e., in the second lobe, at $\nu=45~{\rm arcsec}^{-1}$ ). The visibility shows variations related to the flux spectrum, which decreases in absorption lines. At these wavelengths, wiggles and dark spots appear on intensity maps (Fig. 14, central right panel) increasing the visibility signal at frequencies higher than the second lobe. The visibility variations are attenuated significantly at lower spectral resolution. Observations at wavelengths of a spectral line and the nearby continuum probe different atmospheric depths, and thus layers at different temperatures. They thus provide important information about the wavelength dependence of limb darkening. Since the horizontal temperature and density fluctuations depend on the atmospheric depth, differential observations with relative phase determination provide unique constraints on the granulation pattern. The visibility variations in Fig. 15, such as the steep visibility jump from 0.123 to 0.107 between 1.5975 and 1.5980 $\mu$ m, could be measured in differential interferometric mode at high spectral resolution with the current precision of VLTI-AMBER (1% for visibilities of $\sim$ 5-10%), with optimal sky conditions.

$\begin{figure} \par\includegraphics[width=8.5cm,clip]{178fig27.eps}\hspace*{3mm}... ...m}\includegraphics[width=8.5cm,clip]{178fig29.ps}\hspace*{4.1cm}} \end{figure}$

Figure 15:

Top left panel: three-dimensional view of the visibility curves as a function of wavelength for a particular position angle. The spectral resolution is 12 000. Top right panel: same as in top left panel at a spectral resolution of 1 500. Bottom panel: visibility as a function of wavelength for a baseline of $\sim$ 15 m, i.e., the top of the second lobe, at one particular position angle. The simulation has been scaled to an apparent diameter of $\sim$ 43.6 mas. The synthetic spectrum convolved to a resolution of 12 000 is overplotted (thin solid line). The small black dots correspond to the highest resolution with AMBER (12 000), while the large red dots correspond to medium resolution (1500). The crosses show the uniform disk of 43.6 mas. When changing the position angle, the expected standard deviation is about 10 $\%$ of the visibility (see Fig. 12).

Open with DEXTER

5.5 Closure phase: departure from circular symmetry

Since terrestrial atmospheric turbulence affects the phases of the complex visibilities with random errors, it is impossible to derive them for individual pairs of telescopes. One instead uses closure phase between three telescopes, because the sum of all phase differences removes the atmospheric contribution,but not the phase information of the object visibility (see e.g., Monnier 2007). The closure phase is thus an important complementary piece of information, which can detect asymmetries in the RSG atmospheres. Figure 16 shows the scatter plot of the closure phase of one snapshot of the RHD model computed in the IONIC filter (the scatter is similar for the K222 filter). The behavior is similar for all snapshots. We used 500 random baseline triangles with a maximum linear extension of 40 m, and plot the closure phase as a function of the triangle maximum baseline. The closure phases deviate from zero or $\pm\pi$ already at $\sim$ 10 m (0.0008 $R^{-1}_\odot$ , if we scale the model to an apparent diameter of 43.6 mas at a distance of 174.3 pc). At higher baselines, it clearly differs from zero or $\pm\pi$ , values that are indicating of a point symmetric brightness distribution. This is a clear signature of surface inhomogenities. The characteristic size distribution on the stellar surface can also be derived from the closure phase: the contribution of small-scale convection-related surface structures increases with frequency. The first deviation at $\sim$ 0.0008 $R^{-1}_\odot$ (just beyond the first zero, see Fig. 11) corresponds to the deviation of the stellar disk from circular symmetry. It may be very efficient to constrain the level of asymmetry of RSG atmospheres by accumulating statistics on closure phase at short and long baselines, since they are easily measured to high precision. A small departure from zero immediately inters a departure from symmetry.

$\begin{figure} \par\includegraphics[width=8.5cm,clip]{178fig30.eps} \end{figure}$

Figure 16:

Scatter plot of closure phases in the IONIC filter centered on 1.64 $\mu$ m of 500 random baseline triangles with a maximum linear extension of 40 m. Closure phases are plotted against the longest baseline of the triangle. The upper x-axis corresponds to synthetic observations of the simulation at an apparent diameter of 43.6 mas (which corresponds to $\alpha$ Ori at a distance of 174.3 pc). The axisymmetric case is represented by the grey lines.

Open with DEXTER

We also computed the closure phase for the different K band VLTI-AMBER spectral resolution intensity maps of Fig. 14. The large deviations from circular symmetry are already noticeable at low spectral resolution (Fig. 17, left panels) and the closure phase scatter does not differ much from that at high spectral resolution (right panel). This provides the possibility of detecting asymmetries caused by granulation without being compelled to observe at high spectral resolution.

$\begin{figure} \par\includegraphics[width=8.8cm,clip]{178fig31.eps}\includegraphics[width=8.8cm,clip]{178fig32.eps} \end{figure}$	Figure 17: Scatter plot of closure phases (cf. Fig. 16 for details) obtained from the VLTI-AMBER K band low, and high spectral resolution intensity maps (Fig. 14, central left, and right panels).
Open with DEXTER

6 Comparison with interferometric observations of $\alpha$ Ori

A preliminary comparison of our model predictions with true observations is possible for $\alpha$ Ori. We compare the synthetic visibilities derived from our RHD simulations in the continuum filter K222 (Fig. 2) with the observation of $\alpha$ Ori by Perrin et al. (2004a) that reach the third lobe in the K band. The absolute model dimensions were scaled to match the interferometric observation in the first lobe. This corresponds to an apparent diameter of 43.6 mas at a distance of 179 pc. These values agree with those of Perrin et al., who found a diameter of $43.64\pm0.10$ mas, and Harper et al. (2008), who reported a distance of $197\pm45$ pc.

We computed over 2000 visibility curves and we found that the data can be described by visibility fluctuations caused by the granulation of the simulation (Fig. 18), as already shown in Chiavassa et al. (2007). Among this large number of visibility curves, we find some that match all the observational datapoints with greater accuracy than the uniform disk (with a diameter of 43.33 mas; Perrin et al. 2004a), or limb-darkened disk model (linear-limb darkening law, $I\left(\mu\right)=1-a\left(1-\mu\right)$ , with a diameter of 43.64 mas and a=0.09 also in Perrin et al., see Fig. 18). The best-fit model solution corresponded to a reduced $\chi ^2=0.21$ , and all the visibility curves fall within a $\chi^2$ range of [0.21, 18.1]. Our RHD simulations represent a great improvement over parametric models (the UD model with reduced $\chi^2=19.9$ , and the LD model with $\chi ^2=22.3$ ) when interpreting these interferometric observations. The observation datapoints in the first, second, and third lobes can be reproduced by a single visibility curve, from the projection at a particular position angle of one of our snapshots (see Fig. 18). There is one observed point in the first lobe at 24.5 arcsec^-1 that is difficult to reproduce. Adjustments to the absolute model dimensions of the star required to reproduce this datapoint, would lead to mismatch of the other observations at higher frequencies. However, this may be a problem when calibrating the observation.

A more detailed comparison with $\alpha$ Ori data in the H band (Haubois et al. 2006) will be presented in a forthcoming paper (Chiavassa et al. 2009, in prep.)

$\begin{figure} \par {\hspace*{4.75cm}\includegraphics[width=7.5cm,clip]{178fig33... ...}\hspace*{3mm} \includegraphics[width=8.5cm,clip]{178fig37.eps} } \end{figure}$

Figure 18:

Comparison of the RHD simulation with the $\alpha$ Ori observations (red dots with error bars) by Perrin et al. (2004a). Top panel: intensity map in the K222-FLUOR filter of the snapshot that most closely matches the interferometric data. The range is [0; 10⁵] erg cm^-2 s^-1 Å^-1. The stellar parameters of this snapshot are: $L = 93~480~{L}_\odot$ , $R = 833~{R}_\odot$ , $T_{\rm eff} = 3497$ K, and $\log~g = -0.34$ . The simulation has been scaled to an apparent diameter of 43.6 mas at a distance of 179 pc. The white line indicates the position angle of the projected baselines. Other panels: synthetic visibilities from the simulation, compared with observations. Panel B covers the entire observational range, and panels C-E are close-ups of clusters of observations. The thick, solid line corresponds to the most closely matching visibility curve(reduced $\chi ^2=0.21$ ). The thin solid lines show the minimum and maximum extent of variations in the visibilities. The dot-dashed, and the dashed lines are the LD, and the UD models used by Perrin et al. (reduced $\chi ^2=22.3$ , and 19.9). Note the logarithmic visibility scale.

Open with DEXTER

7 Conclusions

Our radiation hydrodynamical simulations have confirmed that only a few large granules cover the surfaces of RSG stars. The granules in the simulation that we analyze here are 400-500 $R_\odot$ in diameter, and have lifetimes of years. Smaller scale structures develop and evolve within these large granules, on shorter timescales (of a month).

We have demonstrated that RHD simulations are essential to a proper quantitative analysis of interferometric observations of the surface of RSGs beyond the smooth, symmetrical, limb-darkened intensity profiles. We present new average limb darkening coefficients within the H and K bands, which differ significantly from those commonly used in UD or LD profiles. However, these LD coefficients fluctuate with time, and the average is only indicative. Our model surface granulation causes angular and temporal variations in both visibility amplitudes and phases. In the first lobe, which is sensitive to radius, fluctuations can be as high as 5 $\%$ , and radii determinations can be affected to this extent: the radius determined with a UD fit is 3-5% smaller than the radius of the simulation, while the radius determined with a fully LD model is 1% smaller.

The second, third, and fourth lobes, which carry a signature of both limb-darkening and smaller scale structure, are very different from the simple LD case. The visibility amplitudes can be greater than the UD or LD case, and closure phases differ significantly from 0 and $\pm\pi$ , because of a departure from circular symmetry. The visibilities also show fluctuations with both position angle and time, which are directly related to the granulation contrast. We also want to emphasize that data at high spectral resolution provides extremely valuable information. The stellar surface differs dramatically between images of absorption line and its nearby continuum, and differential observations should be easier to carry out with high precision. At lower resolution (e.g. R=1500), continuum and line information become mixed and there is a considerable loss of differential signal.

We conclude that the detection of the signature of granulation, as predicted by our simulation, is measurable with todays interferometers, if observations of both amplitude and closure phase are made at high spatial frequencies (second, third, and fourth lobes, or even higher if possible). These observations will provide direct information about the timescale of the variation, and the size and contrast of granulation.

A few RSGs are prime targets for interferometry, because of their large diameter, proximity, and high infrared luminosity. Only 4 or 5 are sufficiently close and bright for imaging to be possible, but a larger number (10-20) are within reach of less ambitious programs, such as closure phase measurements.

Three approaches can be combined to characterize the granulation pattern by

searching for angular visibility variations, observing with the same telescope configuration (covering high spatial frequencies) and using the Earth's rotation to study 6-7 different position angles in one night, should be sufficient, if measurement errors can be kept below 10%, for visibilities of the order of between 5 and 10%. One or two other telescope configurations would provide more frequency points, but the change of configuration would then have to be made within days, which is possible at VLTI;
looking for temporal visibility fluctuations by observing at two (preferably more) epochs $\sim$ 1 month apart with the same telescope configuration. This can easily be scheduled with existing interferometers, such as CHARA or the VLTI;
looking for visibilities as a function of wavelength, at high spectral resolution, in different spectral regions belonging to both spectral features and continuum. If measurement errors can be kept close to 1%, for visibility of the order of $\sim$ 10%, variations in the visibility correlated with the flux spectrum could be detected, indicating variations in the radius, the limb-darkening, or the granulation pattern. These relative measurements are more easily completed at the required precision than absolute visibility measurements.

These observations provide a wealth of information about both the stars and our RHD models. which we know are affected by the limitations. Comparisons with observations will help us decide which approximations must be relaxed. The simulations are primarily constrained by execution time, which requires several approximations, the most important ones being the limited spatial resolution and the complete lack of wavelength resolution, i.e., grey radiative transfer. This speeds up the simulations to managable execution times of between several months and one year of intensive calculations for about seven years of stellar time. A higher numerical resolution shows smaller scale structures appearing within the granules that are already present in lower resolution simulations (Chiavassa et al. 2006, Fig. 10). However, this should not affect the first few visibility lobes.

The approximation of grey radiative transfer is only justified in the stellar interior and is inappropriate in the optically thin layers. The implementation of non-grey opacities (e.g., five wavelength groups employed to describe the wavelength dependency of the radiation field within a multigroup radiative transfer scheme, see Nordlund 1982, for details) would provide an important improvement to the hydrodynamical simulations. For local RHD simulations, Ludwig et al. (1994) found that frequency-dependent radiative transfer causes an intensified heat exchange of a fluid element with its environment, which tends to reduce the temperature differences. Consequently, the temperature fluctuations in the non-grey local models are smaller than in the grey case. This is also expected for global RSG models, and the result of such a decrease in the temperature fluctuations, would be a decrease in both the intensity contrast and the visibility fluctuations.

Acknowledgements

This project was supported by the French Ministry of Higher Education through an ACI (PhD fellowship of Andrea Chiavassa, postdoc fellowship of Bernd Freytag, and computational resources). Present support is ensured by a grant from ANR (ANR-06-BLAN-0105). We are also grateful to the PNPS and CNRS for their financial support through the years. We thank the CINES for providing some of the computational resources necessary for this work. We thank Michel Belmas and Philippe Falandry for their help. Part of this work was made while BPz was on sabbatical at Uppsala Astronomical Observatory. We thank Bengt Gustafsson, Hans-Gunter Ludwig, John Monnier, Martin Asplund, Nik Piskunov, and Nils Ryde for enlightening discussions.

References

Abramowitz, M., & Stegun, I. A. 1972, Handbook of mathematical functions, 9th edn. (New York: Dover publications, INC.)
Alvarez, R., & Plez, B. 1998, A&A, 330, 1109 [NASA ADS]
Asplund, M. 2000, A&A, 359, 755 [NASA ADS]
Asplund, M., Grevesse, N., & Sauval, A. J. 2006, Commun. Asteroseismol., 147, 76 [NASA ADS] [CrossRef]
Berger, J.-P., Haguenauer, P., Kern, P. Y., et al. 2003, in Interferometry for Optical Astronomy II, ed. W. A. Traub. Proc. SPIE, 4838, 1099
Carr, J. S., Sellgren, K., & Balachandran, S. C. 2000, ApJ, 530, 307 [NASA ADS] [CrossRef]
Cayrel, R., Steffen, M., Chand, H., et al. 2007, A&A, 473, L37 [NASA ADS] [CrossRef] [EDP Sciences]
Chiavassa, A., Plez, B., Josselin, E., & Freytag, B. 2006, in EAS Publ. Ser. 18, ed. P. Stee, 177
Chiavassa, A., Plez, B., Josselin, E., & Freytag, B. 2007, in SF2A-2007: Proceedings of the Annual meeting of the French Society of Astronomy and Astrophysics held in Grenoble, France, ed. J. Bouvier, A. Chalabaev, & C. Charbonnel, 447
Claret, A. 2000, A&A, 363, 1081 [NASA ADS]
Collet, R., Asplund, M., & Trampedach, R. 2007, A&A, 469, 687 [NASA ADS] [CrossRef] [EDP Sciences]
Coude Du Foresto, V., Perrin, G., Ruilier, C., et al. 1998, in Astronomical Interferometry, ed. R. D. Reasenberg, Proc. SPIE, 350, 856
Cunha, K., Sellgren, K., Smith, V. V., et al. 2007, ApJ, 669, 1011 [NASA ADS] [CrossRef]
Cuntz, M. 1997, A&A, 325, 709 [NASA ADS]
Freytag, B., & Höfner, S. 2008, A&A, 483, 571 [NASA ADS] [CrossRef] [EDP Sciences]
Freytag, B., & Ludwig, H.-G. 2007, in SF2A-2007: Proceedings of the Annual meeting of the French Society of Astronomy and Astrophysics held in Grenoble, France, ed. J. Bouvier, A. Chalabaev, & C. Charbonnel, 481
Freytag, B., Holweger, H., Steffen, M., & Ludwig, H.-G. 1997, in Science with the VLT Interferometer, ed. F. Paresce, 316
Freytag, B., Steffen, M., & Dorch, B. 2002, Astron. Nachr., 323, 213 [NASA ADS] [CrossRef]
Gray, D. F. 2008, AJ, 135, 1450 [NASA ADS] [CrossRef]
Gustafsson, B., Edvardsson, B., Eriksson, K., et al. 2008, A&A, 486, 951 [NASA ADS] [CrossRef] [EDP Sciences]
Harper, G. M., Brown, A., & Guinan, E. F. 2008, AJ, 135, 1430 [NASA ADS] [CrossRef]
Hartmann, L., & Avrett, E. H. 1984, ApJ, 284, 238 [NASA ADS] [CrossRef]
Haubois, X., Perrin, G., Lacour, S., et al. 2006, in SF2A-2006: Semaine de l'Astrophysique Française, ed. D. Barret, F. Casoli, G. Lagache, A. Lecavelier, & L. Pagani, 471
Hauschildt, P. H., Baron, E., & Allard, F. 1997, ApJ, 483, 390 [NASA ADS] [CrossRef]
Iglesias, C. A., Rogers, F. J., & Wilson, B. G. 1992, ApJ, 397, 717 [NASA ADS] [CrossRef]
Josselin, E., & Plez, B. 2007, A&A, 469, 671 [NASA ADS] [CrossRef] [EDP Sciences]
Levesque, E. M., Massey, P., Olsen, K. A. G., et al. 2005, ApJ, 628, 973 [NASA ADS] [CrossRef]
Levesque, E. M., Massey, P., Olsen, K. A. G., et al. 2006, ApJ, 645, 1102 [NASA ADS] [CrossRef]
Levesque, E. M., Massey, P., Olsen, K. A. G., & Plez, B. 2007, ApJ, 667, 202 [NASA ADS] [CrossRef]
Ludwig, H.-G., Jordan, S., & Steffen, M. 1994, A&A, 284, 105 [NASA ADS]
Massey, P., Levesque, E. M., Olsen, K. A. G., Plez, B., & Skiff, B. A. 2007, ApJ, 660, 301 [NASA ADS] [CrossRef]
Monnier, J. D. 2007, New Astron. Rev., 51, 604 [NASA ADS] [CrossRef]
Nardetto, N., Fokin, A., Mourard, D., & Mathias, P. 2006, A&A, 454, 327 [NASA ADS] [CrossRef] [EDP Sciences]
Nordlund, A. 1982, A&A, 107, 1 [NASA ADS]
Perrin, G., Ridgway, S. T., Coudé du Foresto, V., et al. 2004a, A&A, 418, 675 [NASA ADS] [CrossRef] [EDP Sciences]
Perrin, G., Ridgway, S. T., Mennesson, B., et al. 2004b, A&A, 426, 279 [NASA ADS] [CrossRef] [EDP Sciences]
Pijpers, F. P., & Hearn, A. G. 1989, A&A, 209, 198 [NASA ADS]
Plez, B., Smith, V. V., & Lambert, D. L. 1993, ApJ, 418, 812 [NASA ADS] [CrossRef]
Steffen, M., & Freytag, B. 2007, Astron. Nachr., 328, 1054 [NASA ADS] [CrossRef]
Stein, R. F., & Nordlund, A. 1998, ApJ, 499, 914 [NASA ADS] [CrossRef]
Traub, W. A., Ahearn, A., Carleton, N. P., et al. 2003, in Interferometry for Optical Astronomy II, ed. W. A. Traub, Proc SPIE, 4838, 45

Footnotes

... Manual: www.astro.uu.se/ bf/co5bold_main.html
... and: http://www.aip.de/ mst/Linfor3D/linfor_3D_manual.pdf

All Tables

Table 1: Gauss-Laguerre quadrature weights for n=10.

Table 2: Time-averaged limb-darkening coefficients for the RHD simulation.

All Figures

$\begin{figure} \par\includegraphics[width=9cm,clip]{178fig1.eps}\vspace{2mm} \i... ...ig2.eps}\vspace{2mm} \includegraphics[width=9cm,clip]{178fig3.eps} \end{figure}$	Figure 1: Radius (Panel A), luminosity (Panel B) and temperature (Panel C) as a function of time for the simulation used in this work. The radius is fitted with the law: $R(t)=936.9 - 4.3\times t$ for $t\leq 23.8$ yrs and $R=832~{R}_{\odot}$ for t>23.8 yrs.
Open with DEXTER
In the text

$\begin{figure} \par\includegraphics[width=9cm,clip]{178fig4.ps} \end{figure}$	Figure 2: The transmission curves of the 4 narrow band filters mounted on the FLUOR instrument at IOTA together with the K band synthetic spectrum of a snapshot of the simulation and the corresponding continuum (bottom black curve). From the top, the spectra computed with only H₂O (red), only CO (green), and only CN (blue) are shown with an offset of, respectively, $0.9\times 10^{33}$ , $0.6\times 10^{33}$ , and $0.3\times 10^{33}$ erg cm^-2 s^-1 Å^-1.
Open with DEXTER
In the text

$\begin{figure} \par\includegraphics[width=9cm,clip]{178fig5.ps} \end{figure}$	Figure 3: The transmission curve of the filter mounted on IONIC at IOTA together with the H band synthetic spectra computed as in Fig. 2. From the top, the offset of the spectra is $3\times 10^{33}$ , $2\times 10^{33}$ , and $1\times 10^{33}$ erg cm^-2 s^-1 Å^-1.
Open with DEXTER
In the text

$\begin{figure} \par\includegraphics[width=18.5cm,clip]{178fig6.eps}\par\includegraphics[width=18.5cm,clip]{178fig7.eps} \end{figure}$	Figure 4: Top 6 panels: maps of the intensity in the IONIC filter (linear scale with a range of [0; $2.5\times 10^5$ ] erg cm^-2 s^-1 Å^-1). The different panels correspond to snapshots separated by 230 days ( $\sim$ 3.5 years covered). Bottom 6 panels: successive snapshots separated by 23 days ( $\sim$ 140 days covered).
Open with DEXTER
In the text

$\begin{figure} \par\hspace{1cm}\includegraphics[width=7cm,clip]{178fig8.ps}\hsp... ...ig10.ps}\hspace{3mm} \includegraphics[width=8cm,clip]{178fig11.ps} \end{figure}$	Figure 5: Top left panel: three-dimensional image of a snapshot from Fig. 4. Top right panel: intensity map of the same snapshot represented using the histogram equalization algorithm to underline the thin bright patches because of the undersampling of the $\tau$ scale. Bottom left panel: intensity profiles for three position angles of the same snapshot. The numerical box edge is at impact parameter $r/R_{\star }\sim 1.3$ . The intensity is normalized to the intensity at disk center. Bottom right panel: radially averaged intensity profiles for all the snapshots of the simulation (grey); one snapshot of the simulation is emphasized with a solid black line. The intensity is normalized to the area subtended by the curve, area $=\int_0^{1.3} I\left(r/R_{\star}\right) {\rm d}r/R_{\star}$ .
Open with DEXTER
In the text

$\begin{figure} \par\includegraphics[width=8.5cm,clip]{178fig12.ps} \end{figure}$	Figure 6: Example of a LD fit (dashed line) using the LD law described in the text for the radially averaged intensity profile (solid line) emphasized in Fig. 5 (bottom right panel). The intensity is normalized to the area subtended by the curve. This best-fit model has a $\chi ^2=0.02$ .
Open with DEXTER
In the text

$\begin{figure} \par\includegraphics[width=8.5cm,clip]{178fig13.ps} \end{figure}$	Figure 7: The time-averaged H-band radial intensity profile of our simulation (solid line), and the fit of Table 2 (dashed line). A fully LD (dash-dotted line), a partially LD (triple dot-dashed line), and a LD fit from Claret (2000, dotted line) for comparable RSG parameters are plotted for comparison.
Open with DEXTER
In the text

$\begin{figure} \par\includegraphics[width=8.5cm,clip]{178fig14.eps} \end{figure}$	Figure 8: The solid line is the visibility curve for the IONIC filter intensity map of Fig. 5 (top right). The dotted line is computed for the same map after applying a $[3\times 3]$ median smoothing.
Open with DEXTER
In the text

$\begin{figure} \par\includegraphics[width=8.5cm,clip]{178fig15.eps}\vspace*{2mm} \par\includegraphics[width=8.5cm,clip]{178fig16.eps} \end{figure}$	Figure 9: Top panel: intensity map in the IONIC filter (the range is [0; $2.5\times 10^5$ ] erg cm^-2 s^-1 Å^-1). Bottom panel: visibility curves from the above snapshot computed for 36 different angles, 5 $^\circ$ apart (thin grey lines). Note the logarithm visibility scale. The solid black curve is a UD model, with a radius of 810 $R_\odot$ . The dashed black line is a partially LD disk with a radius of 822 $R_\odot$ . The dot-dashed line is a fully LD disk with a radius of 830 $R_\odot$ . the triple-dot-dashed line is our average LD law (cf. Table 2) for a radius of 842 $R_\odot$ . The stellar parameters of this snapshot are: $L = 98~400~{L}_\odot$ , $R_* = 836.5~{R}_\odot$ , $T_{\rm eff} = 3534$ K, and $\log g = -0.34$ .
Open with DEXTER
In the text

$\begin{figure} \par\includegraphics[width=8.5cm,clip]{178fig17.eps} \end{figure}$	Figure 10: Standard deviation of the visibility normalized to the visibility in the first lobe. The solid line indicates the temporal fluctuations for one fixed angle over 3.5 years. The trend is similar for all other angles. The dashed curve corresponds to the angular fluctuations of the snapshot in Fig. 9.
Open with DEXTER
In the text

$\begin{figure} \par\includegraphics[width=8.5cm,clip]{178fig18.eps} \end{figure}$	Figure 11: Same as Fig. 9 for the second, third, and fourth lobes. In addition, the dotted line is the visibility curve for the particular azimuth parallel to the x-axis of the IONIC intensity map of Fig. 9.
Open with DEXTER
In the text

$\begin{figure} \par\includegraphics[width=8.5cm,clip]{178fig19.eps} \end{figure}$	Figure 12: Same as Fig. 10 for the second, third, and fourth lobe.
Open with DEXTER
In the text

$\begin{figure} \par {\hspace{1cm}\includegraphics[width=7cm,clip]{178fig20.eps}... ...g22.ps}\hspace{3mm} \includegraphics[width=8cm,clip]{178fig23.eps} \end{figure}$	Figure 13: Top left panel: three-dimensional image of a snapshot with nominal intensities. Top right panel: same snapshot with a feature contrast reduced to 1%. Bottom left panel: standard deviation of the visibility curves at 36 angles 5 $^\circ$ ) apart, for decreasing feature contrast. Solid line is for the top of the second lobe ( ${\sim}0.0010~{R}^{-1}_\odot$ ), dashed line is for the top of the third lobe ( ${\sim}0.0016~{R}^{-1}_\odot$ ), and dotted line is for the top of the fourth lobe ( ${\sim}0.0022~{R}^{-1}_\odot$ ). Bottom right panel: visibility in the second, third, and fourth lobes for one particular position angle. The dot-dashed line shows the original simulation contrast. The dashed line, and the dotted line show the visibility with a feature contrast reduced to 50 $\%$ , and 1% respectively. The solid line is the fitted LD profile computed for this snapshot.
Open with DEXTER
In the text

$\begin{figure} \par\includegraphics[width=18cm,clip]{178fig24.eps}\vspace{3mm} ... ...4cm}\includegraphics[width=10cm,clip]{178fig26.eps}\hspace{4cm}}\end{figure}$	Figure 14: Top row: synthetic spectrum centered on the CO first overtone band head. The left panel shows the range of wavelengths spanned by one resolution element at the VLTI-AMBER low spectral resolution of 35. The right panel shows the same for the VLTI-AMBER medium spectral resolution of 1500, and the high spectral resolution of 12 000 (thick mark). Central row: intensity maps for those three spectral resolution elements. The intensity range is [0;10⁵] erg cm^-2 s^-1 Å^-1. Bottom row: standard deviation of the visibility in the second, third, and fourth lobes. Solid and dashed lines correspond to low and medium resolution, respectively, and the dotted line to high resolution.
Open with DEXTER
In the text

$\begin{figure} \par\includegraphics[width=8.5cm,clip]{178fig27.eps}\hspace{3mm}... ...m}\includegraphics[width=8.5cm,clip]{178fig29.ps}\hspace{4.1cm}} \end{figure}$	Figure 15: Top left panel: three-dimensional view of the visibility curves as a function of wavelength for a particular position angle. The spectral resolution is 12 000. Top right panel: same as in top left panel at a spectral resolution of 1 500. Bottom panel: visibility as a function of wavelength for a baseline of $\sim$ 15 m, i.e., the top of the second lobe, at one particular position angle. The simulation has been scaled to an apparent diameter of $\sim$ 43.6 mas. The synthetic spectrum convolved to a resolution of 12 000 is overplotted (thin solid line). The small black dots correspond to the highest resolution with AMBER (12 000), while the large red dots correspond to medium resolution (1500). The crosses show the uniform disk of 43.6 mas. When changing the position angle, the expected standard deviation is about 10 $\%$ of the visibility (see Fig. 12).
Open with DEXTER
In the text

$\begin{figure} \par\includegraphics[width=8.5cm,clip]{178fig30.eps} \end{figure}$	Figure 16: Scatter plot of closure phases in the IONIC filter centered on 1.64 $\mu$ m of 500 random baseline triangles with a maximum linear extension of 40 m. Closure phases are plotted against the longest baseline of the triangle. The upper x-axis corresponds to synthetic observations of the simulation at an apparent diameter of 43.6 mas (which corresponds to $\alpha$ Ori at a distance of 174.3 pc). The axisymmetric case is represented by the grey lines.
Open with DEXTER
In the text

$\begin{figure} \par\includegraphics[width=8.8cm,clip]{178fig31.eps}\includegraphics[width=8.8cm,clip]{178fig32.eps} \end{figure}$	Figure 17: Scatter plot of closure phases (cf. Fig. 16 for details) obtained from the VLTI-AMBER K band low, and high spectral resolution intensity maps (Fig. 14, central left, and right panels).
Open with DEXTER
In the text

$\begin{figure} \par {\hspace{4.75cm}\includegraphics[width=7.5cm,clip]{178fig33... ...}\hspace{3mm} \includegraphics[width=8.5cm,clip]{178fig37.eps} } \end{figure}$	Figure 18: Comparison of the RHD simulation with the $\alpha$ Ori observations (red dots with error bars) by Perrin et al. (2004a). Top panel: intensity map in the K222-FLUOR filter of the snapshot that most closely matches the interferometric data. The range is [0; 10⁵] erg cm^-2 s^-1 Å^-1. The stellar parameters of this snapshot are: $L = 93~480~{L}_\odot$ , $R = 833~{R}_\odot$ , $T_{\rm eff} = 3497$ K, and $\log~g = -0.34$ . The simulation has been scaled to an apparent diameter of 43.6 mas at a distance of 179 pc. The white line indicates the position angle of the projected baselines. Other panels: synthetic visibilities from the simulation, compared with observations. Panel B covers the entire observational range, and panels C-E are close-ups of clusters of observations. The thick, solid line corresponds to the most closely matching visibility curve(reduced $\chi ^2=0.21$ ). The thin solid lines show the minimum and maximum extent of variations in the visibilities. The dot-dashed, and the dashed lines are the LD, and the UD models used by Perrin et al. (reduced $\chi ^2=22.3$ , and 19.9). Note the logarithmic visibility scale.
Open with DEXTER
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.