Global 3D radiation-hydrodynamics models of AGB stars. Effects of convection and radial pulsations on atmospheric structures
AAstronomy & Astrophysics manuscript no. aaagb3dfirstgrid c (cid:13)
ESO 2018October 8, 2018
Global 3D radiation-hydrodynamics models of AGB stars.
Effects of convection and radial pulsations on atmospheric structures
B. Freytag, S. Liljegren, and S. Höfner
Division of Astronomy and Space Physics, Department of Physics and Astronomy, Uppsala University, Box 516, SE-751 20 Uppsala,Sweden (e-mail:
[email protected] )Received ...; accepted ...
ABSTRACT
Context.
Observations of asymptotic giant branch (AGB) stars with increasing spatial resolution reveal new layers of complexity ofatmospheric processes on a variety of scales.
Aims.
To analyze the physical mechanisms that cause asymmetries and surface structures in observed images, we use detailed 3Ddynamical simulations of AGB stars; these simulations self-consistently describe convection and pulsations.
Methods.
We used the CO5BOLD radiation-hydrodynamics code to produce an exploratory grid of global "star-in-a-box" modelsof the outer convective envelope and the inner atmosphere of AGB stars to study convection, pulsations, and shock waves and theirdependence on stellar and numerical parameters.
Results.
The model dynamics are governed by the interaction of long-lasting giant convection cells, short-lived surface granules, andstrong, radial, fundamental-mode pulsations. Radial pulsations and shorter wavelength, traveling, acoustic waves induce shocks onvarious scales in the atmosphere. Convection, waves, and shocks all contribute to the dynamical pressure and, thus, to an increaseof the stellar radius and to a levitation of material into layers where dust can form. Consequently, the resulting relation of pulsationperiod and stellar radius is shifted toward larger radii compared to that of non-linear 1D models. The dependence of pulsation period onluminosity agrees well with observed relations. The interaction of the pulsation mode with the non-stationary convective flow causesoccasional amplitude changes and phase shifts. The regularity of the pulsations decreases with decreasing gravity as the relative sizeof convection cells increases. The model stars do not have a well-defined surface. Instead, the light is emitted from a very extendedinhomogeneous atmosphere with a complex dynamic pattern of high-contrast features.
Conclusions.
Our models self-consistently describe convection, convectively generated acoustic noise, fundamental-mode radialpulsations, and atmospheric shocks of various scales, which give rise to complex changing structures in the atmospheres of AGBstars.
Key words. convection – shock waves – methods: numerical – stars: AGB and post-AGB – stars: atmospheres – stars: oscillations(including pulsations)
1. Introduction
Variability with typical periods of 100 – 1000 days is a character-istic feature of stars on the asymptotic giant branch (AGB). Thepronounced changes of the stellar luminosity are generally at-tributed to low-order large-amplitude pulsations. Various corre-lations between pulsation properties and stellar parameters havebeen derived from observations, for example, period-luminosity(P-L) relations. Evolved AGB stars, known as Mira variables,follow a linear P-L relation, similar to that of Cepheids, where alarger luminosity results in a longer period (see, e.g., Whitelocket al. 2008, 2009). Less evolved AGB stars fall on a series of P-L sequences, parallel to the sequence of Mira stars (Wood et al.1999; Wood 2015). These sequences are interpreted as represent-ing stars pulsating in various modes, with Mira variables pul-sating in the fundamental mode while low-amplitude variablespulsate in overtones.Pulsation and convection seem to play a decisive role in theheavy mass loss experienced by AGB stars. According to thestandard picture, the observed outflows – with typical veloci-ties of 5 – 20 km s − and mass loss rates in the range of 10 − – 10 − M (cid:12) yr − – are accelerated by radiative pressure on dustgrains, which are formed at about 2 – 3 stellar radii. Pulsations and non-stationary convective flows trigger strong atmosphericshock waves, which lift gas out to distances where temperaturesare low enough and gas densities are su ffi ciently high to allowfor the condensation of dust particles (for a recent review on thismass loss scenario, see, e.g., Höfner 2015).Classically, variable stars have been analyzed via their lightcurves, which have been obtained for a large sample of objectsand are usually regular for evolved AGB stars. However, high-resolution observations of nearby stars reveal complex irregularstructures and dynamical phenomena on various spatial and tem-poral scales.In the case of Mira (o Cet), Karovska et al. (1991) usedspeckle interferometry with various telescopes to detect asym-metries in the extended atmosphere that changed over time.Later, Karovska et al. (1997) speculated that the asymmetriesdetected with the Faint Object Camera on the Hubble SpaceTelescope (HST) in the UV and optical might be due to spotsor non-radial pulsations. Lopez et al. (1997) found a best fit toobservations in the IR taken with the Infrared Spatial Interferom-eter (ISI) by models with inhomogeneities or clumps. Chandleret al. (2007) presented a number of explanations for their ISIobservations, among them shock waves and non-uniform dustformation. Recently, Stewart et al. (2016b) produced images of Article number, page 1 of 14 a r X i v : . [ a s t r o - ph . S R ] F e b & A proofs: manuscript no. aaagb3dfirstgrid
Mira by analyzing occultations by the rings of Saturn observedwith the Cassini spacecraft showing layers and asymmetries inthe stellar atmosphere.Nevertheless, some phenomena may not be intrinsic to thestar itself, since Mira is part of a binary system and the complexlarge-scale structure of the envelope may be attributed to wind-wind interaction (Ramstedt et al. 2014).However, similar clumpy and non-spherical structures canbe found around other AGB stars. For instance, very recently,Ohnaka et al. (2016) observed the semi-regular (type a) starW Hya with VLT / SPHERE-ZIMPOL and VLTI / AMBER andfind supporting evidence for clumpy dust clouds caused by pul-sations and large convective cells, as predicted by 3D simula-tions for AGB stars (Freytag & Höfner 2008).A large number of observations exist of the carbon starIRC + ff & Buscher (1998) attribute asymmetries indi ff raction-limited interferometric images to envelope clearingalong a bipolar axis. Weigelt et al. (1998) used speckle-maskinginterferometry with the SAO 6 m telescope to find a clumpy dustshell, that was considered as direct evidence for an inhomoge-neous mass-loss process caused large-scale surface convection-cells, whose presence was suggested by Schwarzschild (1975).Recently, Stewart et al. (2016a) investigated the dynamical evo-lution of dust clouds in images reconstructed from aperture-masking interferometric observations using the Keck and VLTand from occultation measurements by Cassini, resulting ina complete change in the asymmetry and distribution of theclumps compared to the observations about 20 years earlier byKastner & Weintraub (1994); Hani ff & Buscher (1998); Weigeltet al. (1998).The reason for these observed asymmetries could be the un-derlying irregular global convective flows in AGB stars, althoughthe star itself or the stellar surface are not directly visible.A realistic modeling of the underlying physical processes isa notoriously di ffi cult problem and theoretical studies found inthe literature are usually restricted to 1D simulations. Still, dy-namical 1D atmosphere and wind models have been used suc-cessfully to explore the basic dust-driven mass-loss process, re-lying on a parameterized description of sub-photospheric ve-locities due to radial pulsations (so-called piston models; see,e.g., Bowen 1988; Fleischer et al. 1992; Winters et al. 2000;Höfner et al. 2003; Jeong et al. 2003; Höfner 2008). Mass lossrates and wind velocities along with spectral energy distribu-tions, variations of photometric colors with pulsation phase, andother observable properties, resulting from the latest generationof such time-dependent atmosphere and wind models computedwith the DARWIN code, show good agreement with observa-tions (Nowotny et al. 2010; Eriksson et al. 2014; Bladh et al.2015; Höfner et al. 2016).One-dimensional CODEX models of gas dynamics in thestellar interior and atmosphere, simulating (radial) Mira-typepulsations as self-excited oscillations and following the propaga-tion of the resulting shock waves through the stellar atmosphere,have been presented by Ireland et al. (2008, 2011), for instance.Despite their success in reproducing various observationalresults, 1D dynamical models are not su ffi cient to give a com-prehensive picture of the physical processes leading to mass losson the AGB. These 1D models require a number of free param-eters that have to be carefully adjusted, for example, to describethe variable inner boundary ("piston") in dynamical atmosphereand wind models or the time-dependent extension of the mixing-length theory (MLT) in non-linear pulsation models. Further, such models cannot simulate intrinsically 3D phenomena suchas stellar convection and they therefore cannot describe the gi-ant convection cells that are presumably giving rise to the non-spherical and clumpy morphology of the atmosphere.Turbulent flows, in general, and stellar convection, in par-ticular, are known to produce acoustic waves if the Mach num-ber is su ffi ciently large (Lighthill 1952). This excitation processhas been studied with local 3D radiation-hydrodynamics simula-tions for the case of the Sun, for example, by Nordlund & Stein(2001) and Stein & Nordlund (2001). The steepening of thesewaves while traveling upward into the chromosphere was inves-tigated by Wedemeyer et al. (2004). Such local 3D radiation-hydrodynamics simulations have been used for decades to modelsmall patches on the surface of (more or less) solar-type stars.Now a number of grids produced by di ff erent groups with vari-ous codes are available (Ludwig et al. 2009; Magic et al. 2013;Trampedach et al. 2013; Beeck et al. 2013), which also extendto the regime of white dwarfs (Tremblay et al. 2015).In contrast to these compact types of stars, giants might becovered by only a small number of huge convective cells, as sug-gested by Stothers & Leung (1971) as an explanation for irregu-larities in the light curves. Schwarzschild (1975) argued that, ifthe size of surface convection cells is governed by some charac-teristic length scale, such as the density scale height, the counter-part of small solar granules should be huge convective structureson giants. Consequently, these stars should produce large-scaleacoustic waves. Surface convection and associated waves, pulsa-tion, and shocks have been investigated with global 3D radiation-hydrodynamics simulations with the CO5BOLD code (Freytaget al. 2012), both for the case of red supergiants, for instance, byFreytag et al. (2002) and Chiavassa et al. (2009), and with an ex-ploratory study for an AGB star by Freytag & Höfner (2008).These studies confirm the existence of large-scale convectioncells and acoustic waves.In this paper, we present the first grid of AGB star modelsproduced with a new, improved version of CO5BOLD. We an-alyze the dependence of the properties of convection and pul-sations on stellar parameters and also look at the influence ofrotation.
2. Setup of global AGB star models with CO5BOLD
Following our earlier work (Freytag & Höfner 2008), we presentnew radiation-hydrodynamics simulations of AGB stars withCO5BOLD (Freytag et al. 2012). Improvements regarding theaccuracy and stability of the hydrodynamics solver are outlinedin Freytag (2013).The CO5BOLD code solves the coupled non-linear equa-tions of compressible hydrodynamics (with an approximate Roesolver) and non-local radiative energy transfer (for global mod-els with a short-characteristics scheme) in the presence of a fixedexternal gravitational field. The numerical grid is Cartesian. Inall models presented here, the computational domain and all in-dividual grid cells are cubical. The tabulated equation of statetakes the ionization of hydrogen and helium and the formation ofH molecules into account. It is assumed that solar abundancesare appropriate for M-type AGB stars. The tabulated gray opac-ities (very similar to those used in Freytag & Höfner 2008) aremerged from Phoenix (Hauschildt et al. 1997) and OPAL (Igle-sias et al. 1992) data at around 12 000 K, with a slight reductionin the very cool layers (that are hardly reached in the currentmodel grid) to remove any influence of dust onto the opacities.There are no source terms or dedicated opacities for dust: no dust Article number, page 2 of 14. Freytag et al.: Global 3D radiation-hydrodynamics models of AGB stars.
Table 1.
Basic model parameters model M (cid:63) M env L (cid:63) n x × n y × n z x box P rot t avg R (cid:63) T e ff log g P puls σ puls M (cid:12) M (cid:12) L (cid:12) R (cid:12) yr yr R (cid:12) K (cgs) yr yrst28gm06n02 1.0 0.196 7079 127 ∞ ∞ ∞ ∞ ∞ ∞ ∞ ∞ ∞ ∞ ∞ ∞ ∞ ∞ ∞ ∞ ∞ ∞ ff ective temperature and surface gravity and of a running number;the mass M (cid:63) , used for the external potential; the envelope mass M env , derived from integrating the mass density of all grid cellswithin the computational box; the average emitted luminosity L (cid:63) , close but not identical to the inserted luminosity of either 5000,7000, or 10000 L (cid:12) in the core; the model dimensions n x × n y × n z ; the edge length of the cubical computational box x box ; the rotationalperiod P rot ; the time t avg , used for the averaging of the rest of the quantities in this table and for the further analysis; the averageapproximate stellar radius R (cid:63) ; the average approximate e ff ective temperature T e ff ; the logarithm of the average approximate surfacegravity log g ; the pulsation period P puls ; and the half width of the distribution of the pulsation frequencies σ puls . The first (“old”)group of simulations comprises models used in Freytag & Höfner (2008). In the second (“test”) group, numerical parameters (e.g.,the box size or the number of grid points) or the rotation period were varied. The last (“grid”) group comprises models with slightlyvarying stellar parameters ( M env , L (cid:63) ).is included in any of the current models in contrast to the two oldmodels st28gm06n05 and n06.The gravitational potential is spherically symmetric, (seeEq. (41) in Freytag et al. 2012) corresponding to a 1 M (cid:12) corein the outer layers and smoothed in the center at r (cid:46) R = R (cid:12) (see Fig. 4 in Freytag & Höfner 2008). The tiny central nuclear-reaction region cannot possibly be resolved with grid cells ofconstant size. Instead, in the smoothed-core region, a source termfeeds in energy corresponding to 5000, 7000, or 10000 L (cid:12) . Adrag force is only active in this core to prevent dipolar flowstraversing the entire star. All outer boundaries are open for theflow of matter and for radiation.The start models for the first AGB simulations presented inFreytag & Höfner (2008) were based on hydrostatic 1D strati-fications (see Fig. 1 of Chiavassa et al. 2011a, showing imagesfrom the initial phase of a simulation of a red supergiant startedthis way). All later AGB models were derived from a predeces-sor 3D model by either refining or coarsening the numerical grid,adding or removing grid layers to change the size of the compu-tational box, increasing or decreasing the core luminosity in timein steps of 1000 L (cid:12) , or by modifying the density in each gridcells in each time step by a tiny amount to change the initial tothe final envelope mass. While the initial adjustment time lastedtypically a few years, the typical total simulation time span is10 s (just above 30 yr). We monitored the time evolution of sev-eral quantities, for instance, the total emerging luminosity and the mass density in the very center, to ensure that the models arewell relaxed when we compute averaged quantities for furtheranalysis (see below).Table 1 and Fig. 1 give an overview of the simulations splitinto three groups: old models from Freytag & Höfner (2008),new test models used for parameter studies (including rotationrate), and a grid of models with di ff erent stellar parameters.While, for instance, the mass M (cid:63) (controlling the gravita-tional potential), the resolution and the extent of the numericalgrid, and the rotation rate are pre-chosen fixed parameters (sec-ond group of rows in Table 1), other model properties are deter-mined after a simulation is finished (third group in Table 1). Thestated stellar luminosity is a time average of the luminosity foreach “fine” snapshot without the full 3D information, but withlots of preprocessed – and thus compressed – data, saved every250 000 s ≈ M env , however, is calculatedfrom the integrated density for each fine snapshot averaged overtime. Nevertheless, the radius is more di ffi cult to determine andless well defined. The radius is chosen as that point R (cid:63) wherethe spherically (abbreviated as (cid:104) . (cid:105) Ω ) and temporally (denoted as (cid:104) . (cid:105) t ) averaged temperature and the average luminosity (cid:104) L (cid:105) Ω , t ful-fill (cid:104) L (cid:105) Ω , t = πσ R (cid:63) (cid:104) T (cid:105) Ω , t . Then, e ff ective temperature and sur-face gravity follow.To investigate purely radial motions we take averages overspherical shells for each fine snapshot for the radial mass flux (cid:104) ρ v radial (cid:105) Ω ( r , t ) and the mass density (cid:104) ρ (cid:105) Ω ( r , t ) and take the ratio Article number, page 3 of 14 & A proofs: manuscript no. aaagb3dfirstgrid eff [K]4000500060007000800090001000011000 L [ L ⊙ ] Old modelsTest modelsGrid models
Fig. 1.
Luminosity and e ff ective temperature of all models. The colorsof the markers represent the temperature of the models while the sizeof the markers represent the luminosity. The shape of the markers in-dicates which group in Table 1 the models belong to; the “old” modelsfrom (Freytag & Höfner 2008) are circles, the “test” models with variednumerical parameters are diamonds, and the “grid” models with di ff er-ent stellar parameters are squares. as the radial velocity, which is now a function (cid:104) v radial (cid:105) ( r , t ) of ra-dial distance and time. The derivation of the pulsation quantities P puls and σ puls of Table 1 is described in Sect. 3.3.
3. Results
The time evolution of various quantities can be followed in thesnapshot sequences in Figs. 2 and 3; the plots of radial velocityas function of radius and time are shown in Fig. 4 (along a partic-ular radius vector: v radial ( r , t )) and Fig. 5 (averaged over spheres: (cid:104) v radial (cid:105) ( r , t )). The convection zone is essentially indicated by the bright (highentropy) irregular inner part of the entropy snapshots in Fig. 2with a radius around 300 R (cid:12) , which is smaller than the inferredstellar radius of R (cid:63) = R (cid:12) given in Table 1; i.e., light is emittedfrom layers further out. The drop in entropy at the top of theconvective layers is accompanied by a drop in temperature andeven a thin density-inversion layer.Huge convective cells can span 90 degrees or more inthe cross sections, which is in line with extrapolations bySchwarzschild (1975) from solar granulation. The cells are out-lined by non-stationary downdrafts reaching from the surface ofthe convection zone to the center of the model star. While theflow travel time from the surface to the center is around halfa year, the convective cells can have a lifetime of many years,causing long intervals of one preferred flow direction in the con-vection zone in Fig. 4. These downdrafts even have a tendencyto traverse the artificial stellar core and to create a global dipo-lar flow field. However, in the (non-rotating) models presentedhere, the drag force, which is somewhat arbitrarily applied in thestellar core, prevents these flows.The surface of the convection zone appears corrugated; thisis caused by many smaller short-lived convection cells close tothe surface. Their number increases with numerical resolution, as a comparison of Fig. 2 with Fig. 1 of Freytag & Höfner (2008)reveals.Simulations of red supergiants (RSGs) with CO5BOLD(Freytag et al. 2002; Chiavassa et al. 2009; Arroyo-Torres et al.2015) show large-scale convection cells as well. But for RSGsthe ratio of surface pressure scale height, and therefore the sizeof the convective structures and the local radius fluctuations, toradius is smaller than for AGB stars, so that RSGs appear morespherical. The density snapshots in Fig. 2 show irregular structures withconvection cells in the interior and a network of shocks in theatmosphere. To investigate purely radial motions we considerthe density-weighted spherical average of the radial velocity (cid:104) v radial (cid:105) ( r , t ) as described in Sect. 2 and plotted in Fig. 5 for thestandard model st28gm06n26.The behavior of the inner part of the model di ff ers from thatof the atmosphere, which is particularly evident in Fig. 5. Belowabout 400 R (cid:12) (the nominal radius is R (cid:63) = R (cid:12) ), the pulsationis rather regular and coherent over all layers, close to a standingwave. The fundamental mode dominates as there are no nodes inthe velocity map. In the outer layers, however, the slopes in thevelocity map indicate propagating shock waves (see Sect. 3.4).In the right panel in Fig. 5, a close-up of part of the velocityfield is shown, with movements of mass shells overlaid, analo-gous to plots for the 1D models, for example, in Höfner et al.(2003) or Nowotny et al. (2010). The amplitude of the 3D mass-shell oscillations is smaller than for corresponding 1D modelsby a factor of two or more. This is at least partly due to the aver-aging over the full 3D model, which smoothes the amplitude asthe shock waves are not exactly spherical.To quantify the periodic behavior a Fourier analysis is per-formed using the averaged radial velocity (cid:104) v radial (cid:105) ( r , t ). An exam-ple for r = R (cid:12) is shown in the left panel of Fig. 6, again formodel st28gm06n26 with the corresponding power spectrum inthe right panel. Frequencies around 0.6 yr − dominate.To explore the behavior all through the star the power spectraof the averaged radial velocity (cid:104) v radial (cid:105) ( r , t ) at all radii are plottedin the middle panel of Fig. 7. Again the two di ff erent behaviorsin the inner and outer layers of the model are evident. There isa dominant frequency in the stellar interior, as suspected fromFig. 5. However, the outer layers beyond r ∼ R (cid:12) no longerpulsate with the same period as the interior of the star. Lower fre-quency signals become more prominent and beyond r ∼ R (cid:12) (far out in the shock-dominated atmosphere), and the frequencywith the largest amplitude is significantly smaller than that of theinterior of the star.In the power spectra of models with di ff erent stellar parame-ters in Fig. 7, the dominant mode generally becomes more dif-fuse with increasing radius of the model. For more compactmodels, for instance model st29gm04n001 with R (cid:63) = R (cid:12) inthe left panel in Fig. 7, there is a very clear dominant mode inthe interior. For the standard model (middle panel in Fig. 7 withradius R (cid:63) = R (cid:12) ) there is still a dominant frequency, but thespread around this frequency is larger. For the largest model(model st28gm07n001 with R (cid:63) = R (cid:12) in the right panel inFig. 7), the power spectrum seems to be equally distributed overa large frequency range, lacking the clear dominant frequencythat is present in the two other models. Article number, page 4 of 14. Freytag et al.: Global 3D radiation-hydrodynamics models of AGB stars. -5000500 y [ R O • ] -5000500 y [ R O • ] -5000500 y [ R O • ] -5000500 y [ R O • ] -5000500 y [ R O • ] -500 0 500x [R O • ]-5000500 y [ R O • ] -500 0 500x [R O • ] -500 0 500x [R O • ] -500 0 500x [R O • ] t= 0.0 d t= 92.6 d t= 185.2 d t= 277.8 d t= 370.3 d -500 0 500x [R O • ] t= 463.0 d Fig. 2.
Time sequences (for the extended model st28gm06n25) of temperature slices, density slices, entropy slices, pseudo-streamlines (integratedover about 3 months), and bolometric intensity maps, where black corresponds to zero intensity. The snapshots are nearly 3 months apart (seethe counter in the right-most panels), so that the entire sequence covers about one pulsational cycle. The color scales are kept the same for allsnapshots. This figure can be directly compared to Fig. 1 of Freytag & Höfner (2008). Article number, page 5 of 14 & A proofs: manuscript no. aaagb3dfirstgrid -600-400-2000200400600 y [ R O • ] t= 0.0 d t= 28.9 d t= 57.9 d t= 86.8 d t= 112.9 d -600-400-200 0 200 400 600x [R O • ]-600-400-2000200400600 y [ R O • ] t= 141.8 d -600-400-200 0 200 400 600x [R O • ] t= 170.7 d -600-400-200 0 200 400 600x [R O • ] t= 199.7 d -600-400-200 0 200 400 600x [R O • ] t= 228.6 d -600-400-200 0 200 400 600x [R O • ] t= 257.5 d Fig. 3.
Time sequence (for the standard model st28gm06n26) of bolometric intensity maps, where black corresponds to zero intensity. Thesnapshots are about 1 month apart to demonstrate the relative short timescale for changes of the small surface features. r [ R ⊙ ] −30−24−18−12−60612182430 v r a d i a l [ k m s − ] Fig. 4.
Plot (for the standard model st28gm06n26) of radial velocity v radial ( r , t ) for all grid points in one column from the center to the side of thecomputational box for the entire simulation time. Blue indicates outward and red inward flow. r [ R ⊙ ] v r a d i a l [ k m s − ] Fig. 5.
Examples for the standard model st28gm06n26:
Left : The spherically averaged radial velocity (cid:104) v radial (cid:105) ( r , t ) for the full run time and radialdistance. The di ff erent colors show the average vertical velocity at that time and radial distance. Right:
Part of the velocity field from the rightimage, indicated with the rectangle with mass-shell movements plotted as iso-mass contour lines.Article number, page 6 of 14. Freytag et al.: Global 3D radiation-hydrodynamics models of AGB stars. v r a d i a l [ k m s − ] −1 ]0.00.51.01.52.02.5 A r e a - n o r m a li s e d a m p li t u d e [ y r ] Fig. 6.
Examples for the standard model st28gm06n26:
Left:
The vertical velocity at a constant r = R (cid:12) , over 30 years. Right:
The powerspectrum of the vertical velocity, shown to the left, showing the clearly dominant frequency. At the bottom of the panel, the density plot for thispower spectrum is shown. This type of data is used in Fig. 7 for all radial points. r [ R ⊙ ] −1 ] 0.0 0.5 1.0 1.5 2.0 2.5 Fig. 7.
Power spectra derived from the velocity fields of three di ff erent models, mapped over frequency and radial distance. For easier comparisonall the di ff erent power spectra have been area normalized. Left:
Model st29gm04n001, the hottest and smallest model, is shown.
Middle:
Thestandard model st28gm06n26 is shown.
Right:
Model st28gm07n001, the coolest and most luminous model with the largest radius, is shown. -600 -400 -200 0 200 400 600x [R O • ]-600-400-2000200400600 y [ R O • ] -600 -400 -200 0 200 400 600x [R O • ] -600 -400 -200 0 200 400 600x [R O • ] Fig. 8.
Density slices for the three models used in Fig. 7. The range used for all color tables is -16 ≤ log ρ ≤ -6.7. Article number, page 7 of 14 & A proofs: manuscript no. aaagb3dfirstgrid
To find the dominant frequency and therefore the pulsationperiod and to investigate the spread in frequencies, the area-normalized power spectra of the radial velocities for radial dis-tances r = R (cid:63) were added. For instance, this correspondsto radial points in the range r = R (cid:12) for the standardmodel. A Gaussian distribution was fitted in the frequency do-main containing the strongest signal. The central value for the fitfor each model is taken as the dominant frequency f dom , whilethe spread in frequencies is represented by the standard devia-tion σ puls . The resulting periods P puls = / f dom and spreads σ puls are listed in Table 1 and plotted in Fig. 9, where the colors of thesquares represent the temperature with lighter colors indicatinghigher temperature and the size represents the luminosity of themodel, as seen in Fig. 1. The irregular spread of the mode frequencies in Fig. 7 is likelydue to interactions between the pulsations and large-scale con-vective motions causing occasional amplitude changes and phaseshifts. With larger radii, the convective cells increase further inrelative size resulting in stronger disturbances of the pulsationmode. Within a luminosity group the frequency spread growswith decreasing radius, which is likely due to the increase ofconvective velocities with increasing e ff ective temperature; cf.the bottom right panel in Fig. 10.The power spectra in Fig. 7 can be compared to the powerspectrum derived from a local solar model in Fig. 3 of Stein &Nordlund (2001). The similarity is an indication of a commonexcitation mechanism.Analyzing light curves, Christensen-Dalsgaard et al. (2001)attributed oscillations in semi-regular variables to stochastic ex-citation by convection. Bedding et al. (2005) distinguished sev-eral cases. These authors attributed large phase fluctuations inthe light curve of the semi-regular star W Cyg to stochastic exci-tation, whereas the very stable phase of the true Mira star X Camwas found to be consistent with the excitation by the κ mecha-nism. L Pup was classified as an intermediate case in whichboth mechanisms play a role.Analyzing the work integral in 1D pulsation models, Lat-tanzio & Wood (2004) concluded that the driving mechanism,at least for Mira variables, is likely a κ mechanism acting in thepartial hydrogen and helium I ionization zone.Our grid of 3D models, which only covers a small part of theentire range of the AGB star parameters, already shows a rangeof di ff erent behaviors of the oscillations (see Fig. 7 and the dis-cussion in the previous section), where the trend clearly pointsto the role of convection and the size of the convective cellsfor the excitation of, or at least interaction with, the pulsation.The non-stationary transonic convective flows with Mach num-bers, which in the downdrafts often exceed unity, produce acous-tic noise as described by Lighthill (1952) for turbulent flows.This excitation mechanism has been investigated for instance byNordlund & Stein (2001) and Stein & Nordlund (2001) withlocal 3D radiation-hydrodynamics simulations of solar granu-lation. The Mach numbers in near-surface convective flows ofAGB stars can be even larger causing more e ffi cient wave exci-tation (Lighthill 1952). As the relative sizes of the exciting con-vective structures are very large, the generated waves have muchlower wave numbers than on the Sun.To check that the radial pulsations are not just (long-lasting)transient phenomena introduced by the initial conditions, weadded a strong (purely artificial) drag force in the entire modelvolume. This drag force reduced the amplitude of the pulsa- tions (and the convective flows) but did not lead to an exponen-tial decay of the radial mode, indicating that an e ffi cient mode-excitation mechanism must be at work.In the models of Freytag & Höfner (2008), these pulsationsexisted and their amplitude was extracted as a description of thepiston boundary in 1D models. However, the global pulsationswere harder to distinguish from the local shock network becauseof the lower numerical resolution; this led to larger sizes of con-vective and wave structures and shorter time sequences, whichmade a Fourier analysis less reliable than what is possible forthe current model grid. As has been pointed out by Fox & Wood (1982) and Wood(1990), AGB stars do not seem to follow the simple period-mean-density relation, P puls × ( ¯ ρ/ ¯ ρ (cid:12) ) / = Q , where Q is the pul-sation constant. This is not very surprising because the derivationof the period-mean-density relation relies on the assumptionsthat displacements are adiabatic and non-linear e ff ects are small,both of which are probably incorrect for AGB stars. Instead,these works, using 1D pulsation models, find that P puls ∝ R α(cid:63) M − β(cid:63) with α ∼ β ∼ M (cid:63) = M (cid:12) is compared to that ofFox & Wood (1982) and Wood (1990) in the left panel of Fig. 9.We findlog( P puls ) = .
39 log( R (cid:63) ) − . , (1)which gives generally a larger radius for a given period. Theremight be several reasons for this di ff erence in addition to uncer-tainties in the 1D models.There is a contribution to the extension of the atmosphereof the 3D models (see Sect. 4.2) due to the convectively gen-erated, small-scale shocks (see Sect. 3.4) that do not exist in the1D models. This would a ff ect the photospheric radius but not thepulsation period if the shocks occur above the top of the acousticcavity. This contribution might even lower the slope because thelargest models have the most extended atmospheres. The convec-tive envelope is not in hydrostatic equilibrium but is a ff ected bythe convective dynamical pressure as well. In addition, the treat-ment of the artificial core in the 3D models might play a role.However, the e ff ect should be small because the sound-crossingtime for the core, and therefore the contribution to the period, isrelatively small. Finally, there is some uncertainty in the deter-mination of radius and period.In the left panel of Fig. 9, observations of the radii for di ff er-ent AGB stars are plotted against the periods with C stars shownas crosses and M stars as circles. The diameter observations arefrom Richichi et al. (2005), periods from Samus et al. (2009),and parallaxes from van Leeuwen (2007). The 3D models agreefairly well with observations, especially the lower luminositymodels. The period-radius relation from the 1D models mightproduce periods that are too long for the range of radii exploredby the 3D models when compared to observations. However, 1Dmodels predict that altering the mass leads to a change of theperiod-radius relation. While only the relations for 1 M (cid:12) are plot-ted here for direct comparison to the results from the 3D models,1D models from Fox & Wood (1982) can reach the shorter peri-ods, if the mass of the models is increased.It is di ffi cult to draw any final conclusion as the errors in de-termining fundamental parameters for field AGB stars are verylarge. The uncertainties in observed absolute magnitudes origi-nate mainly from uncertainties in the parallaxes, which are dif-ficult to determine for AGB stars as the photo-centers of these Article number, page 8 of 14. Freytag et al.: Global 3D radiation-hydrodynamics models of AGB stars.
Fig. 9.
Left:
Logarithm of the period of the models plotted against the logarithm of the resulting radius with three di ff erent period-radius relationsplotted. The green line indicates period-radius relations from this work, the orange dashed line indicates Wood (1990), and the purple dotted lineindicates Fox & Wood (1982) for M = M (cid:12) . The error bars of the pulsation period of the models refer to the width of the peak in the powerspectrum, that is, larger than the uncertainty of the period. The statistical error in the model radius is small. The typical error in the observed radiiis about 31 % due to uncertainties in the parallax. The crosses and circles are observations of C stars and M stars, respectively, with radius from theCHARM2 catalog (Richichi et al. 2005) and periods from GCVS catalog (Samus et al. 2009). Only stars with measured parallaxes were pickedso that the distance determination is independent of the measured period (Ramstedt & Olofsson 2014). Right : The absolute magnitude against thelogarithm of the period for all models. The error in the absolute magnitude is only due to the finite simulation time and is tiny (a few mmag). Theline is taken from Whitelock et al. (2009) and is the P-L relation for AGB stars in the LMC, where the gray area is the 1 σ error of the fit to theobservations. stars are variable (for a statistical analysis of photocentric vari-ability, see Ludwig 2006; Chiavassa et al. 2011b). Also, the un-certainty of observations of the radii of AGB stars is fairly large,sometimes on the order of the stellar radius. In addition, the ra-dius varies significantly during a pulsation cycle and the phasesare not always well known. Here, we use the mean radius as ref-erence. However, the radius varies by around 20 % during onepulsation period for our models.It is also likely that the stellar masses a ff ect the period-radiusrelation. However, unless there are well-observed companionstars, it is usually not possible to determine masses from obser-vations. A theoretical prediction of the e ff ect of di ff erent massesis not yet possible, as all the models in the current grid have amass of 1 M (cid:12) .A property that is better constrained by observation is the P-Lrelation, which has been extensively studied. A comparison be-tween the 3D models and observation by Whitelock et al. (2009)is shown in the right panel of Fig. 9. The models from the gridfollow a trend of brighter absolute magnitude with larger period,which is qualitatively similar to that of the observations. Thereis, however, a spread in the periods because of the di ff erent ef-fective temperatures and radii for a given luminosity, giving con-straints for our – so far a bit arbitrary – choice of the combinationof the main control parameters M (cid:63) , M env , L (cid:63) , and metallicity ofthe 3D models. The steepening of acoustic waves in the solar chromosphere andthe transformation into a network of shocks was modeled byWedemeyer et al. (2004) and Muthsam et al. (2007). AGB starshave a larger temperature drop from the convection zone to theatmosphere (cf. Figs. 5 and 6 in Freytag & Chiavassa 2013, for amodel sequence from the Sun with log g = M (cid:12) star with log g = -0.43) that is accompanied by a larger change inpressure or density scale height. This leads to a stronger com-pression and amplification of the waves in the cool atmosphere,so that the waves very early turn into shocks, leaving no room foran essentially undisturbed photosphere. Sound-speed variations,particularly at the rugged surface of the convection zone, andtransonic convective flow speeds shape the waves and contributeto the generation of small-scale shock structures (at r > R (cid:12) inFig. 4), which give rise to ballistic gas motions with peak heightsof a few ten to a few hundred solar radii. The strongest shockscan even leave the computational box.The plot of entropy versus radius in Fig. 10 shows that part ofthe atmosphere is a zone of convective instability – with negativeentropy gradient – separate from the normal interior convectionzone. The image sequences in Fig. 2 and the plot of v radial ( r , t )in Fig. 4 show a very di ff erent behavior of the two zones: onlythe inner zone is governed by the normal overturning motions,whereas propagating shocks dominate in the outer zone (as dis-cussed in Sect. 3.3.1). But still, the convective instability mightdestabilize the shocks favoring a tendency toward smaller struc-tures as seen in the intensity snapshots in Figs. 2 and 3. For mod-els with L = L (cid:12) , this instability zone lies two times furtherout relative to the stellar radius, compared to the higher lumi-nosity models. Woitke (2006) showed that in 2D simulations ofthe atmosphere of an AGB star, where the radial stellar pulsationand inhomogeneities generated by convection are ignored andshocks are generated by an external κ mechanism, instabilitieswithin the shock fronts can cause small-scale structures to form.A combination of relatively small-scale acoustic noise andglobal, radial pulsations generates the network of shocks with asize spectrum regulated by several complex processes; a waveemitted from a small region close to the surface might turn intoan expanding shock wave, which fills a large part of the atmo-sphere if the surrounding flow permits it. The overtaking and Article number, page 9 of 14 & A proofs: manuscript no. aaagb3dfirstgrid O • ]-14-12-10-8 l og ( ρ /[ g c m - ] ) O • ]3.54.04.5 l og ( T /[ K ] ) O • ]1.61.82.02.22.42.62.8 s [ e r g K - g - ] L=10028 L O • , R=531 R O • , P rot = 0 yr : st28gm07n001L= 6986 L O • , R=436 R O • , P rot = 0 yr : st26gm07n002L= 6952 L O • , R=400 R O • , P rot = 0 yr : st26gm07n001L= 6954 L O • , R=370 R O • , P rot = 0 yr : st28gm06n26L= 6948 L O • , R=348 R O • , P rot = 0 yr : st29gm06n001L= 4982 L O • , R=345 R O • , P rot = 0 yr : st27gm06n001L= 4977 L O • , R=312 R O • , P rot = 0 yr : st28gm05n002L= 4989 L O • , R=300 R O • , P rot = 0 yr : st28gm05n001L= 4982 L O • , R=294 R O • , P rot = 0 yr : st29gm04n001L= 6955 L O • , R=384 R O • , P rot =20 yr : st28gm06n29L= 6950 L O • , R=394 R O • , P rot =10 yr : st28gm06n30 O • ]12345 v r ad i a l [ k m s - ] O • ]-1.0-0.50.00.51.0 l og ( P d y n / P ) O • ]24681012 v r ad i a l , r m s [ k m s - ] Fig. 10.
Plots of selected quantities, averaged over spherical shells and time vs. radii for all grid models. The plus signs roughly indicate thedepth with τ Ross = Top left:
Logarithm of the mass density ρ . Top right : Logarithm of the temperature T . Center left : Entropy per mass unit s . Center right : Temporal rms value of the spherically averaged radial velocity (includes mostly contributions from the radial pulsations).
Bottomleft : Logarithm of the ratio of dynamical pressure P dyn and gas pressure P . Bottom right : Temporal and spatial rms value of the radial velocity(includes contributions both from convection and radial pulsations). merging of shocks increases the typical size. On the other hand,large-scale waves in the interior are shaped by the convectivebackground flow and sound-speed distribution and can becomeas rugged as the surface of the convection zone itself. And theconvective instability in the atmosphere also favors small scales.The shock fronts show up prominently in the density slices inFig. 2 and are only occasionally visible in the temperature slices(e.g., in the upper left quadrant in the fourth row) because theradiative relaxation is fast enough to cool down the heated post-shock material to the local equilibrium value, even with the grayradiative transfer used. The radiative relaxation time decreasesfurther when non-gray radiative transfer is employed (in a first simulation), essentially wiping out any temperature signature ofthe shocks. That might change again with a drastic increase innumerical resolution, which is currently not achievable.While the temperature stratifications of the various models inFig. 10 show many di ff erences in the convection zone and inneratmosphere, they converge into one of three profiles in the outeratmosphere, depending only on stellar luminosity. This and therelatively small temperature fluctuations in the outer part of thetemperature slices in Fig. 2 indicate that the thermal structure ofthe outer atmosphere is dominated by radiative processes and isnot far from radiative equilibrium. The shocks show up in thevelocities and the density but hardly in the temperature. Article number, page 10 of 14. Freytag et al.: Global 3D radiation-hydrodynamics models of AGB stars.
Expectedly, the size spectrum of the shocks is extended to-ward smaller scales with increasing resolution as a comparisonof Fig. 2 of the current paper with Fig. 1 of Freytag & Höfner(2008) demonstrates.
100 200 300 400 500 600r [R O • ]-600-400-2000200400600 100 200 300 400 500 600r [R O • ]-600-400-2000200400600 Fig. 11.
Azimuthally and temporally averaged velocity fields for modelsst28gm06n29 ( P rot =
20 yr, left) and st28gm06n30 ( P rot =
10 yr, right).The color plot shows the angular momentum in the corotating frame;bright blue means rotation faster than the mean given by P rot , and darkblue means rotation slower than the mean. Rotation plays an important role in stellar activity and pos-sibly also in the dynamics of stellar winds. The loss of angu-lar momentum from magnetic coupling with the surroundingsor from a wind likely slows down evolved stars considerably.Dorfi & Höfner (1996) performed 1D stationary-wind models ofa 1 M (cid:12) / / L (cid:12) AGB star and found that a rotation pe-riod of 40 yr modifies the isotropic mass loss marginally, while aperiod of 10 yr results in a drastic increase of the mass loss rateand causes a significant axial asymmetry of the wind.In hydrodynamical simulations of convection in the interiorof rotating stars, the Mach numbers are so low that the flow isessentially incompressible and the optical thickness of individ-ual grid cells is so large that radiation transport is adequatelydescribed by the di ff usion approximation.However, these conditions are not met in the convective en-velopes and surrounding atmospheres of AGB stars, where thenumerical treatment in CO5BOLD (e.g., of non-local radiationtransport and compressible hydrodynamics) is required, so wemodified the CO5BOLD setup to account for rotation. While it ispossible to just add a rotational velocity field to the initial model,we chose to perform the simulation in a corotating frame. A cen-trifugal potential is added to the gravitational potential. Beforeeach hydrodynamics step the Coriolis force is applied, whichjust rotates the velocity vectors, but with twice the amount sim-ply suggested by period and time step. The artificial drag force inthe model core was chosen to act only radially so as not to a ff ectthe angular momentum. We computed a first exploratory pair ofmodels with rotation periods of 20 and 10 years (see Table 1).Longer periods would require even longer integration times. The P rot =
20 yr simulation was started from a snapshot from a non-rotating run, and the period was decreased (by increasing cen-trifugal potential and Coriolis force) gradually over a few stellaryears to avoid transient oscillations due to a too rapid change ofthe potential.As expected for slow rotators, where the rotational period islonger than typical convective turnover timescales, angular mo-mentum is advected inward into a small region close to the coreof the model (see Fig. 11). In spite of the drag force in the core,a global dipole flow develops with typical velocities in the atmo-sphere of v meridional = / s (for st28gm06n29 with P rot =
20 yr)and v meridional = / s (for st28gm06n30 with P rot =
10 yr).While the core generally rotates very rapidly, the part of theconvection zone that is not close to the axis rotates slower thanthe nominal rate. Part of the material close to the core of themodel with P rot =
20 yr even shows a slow retrograde rotation.However, all azimuthally averaged velocity components showlarge fluctuations and the averaging time or relaxation time forthis model might not yet be su ffi cient.The mean atmospheric stratification is a ff ected by the rotat-ing star: while the temperature stratification shows hardly any ef-fect, the average density in the atmosphere increases with shorterrotation period (see Fig. 10).Our rotating models have a number of shortcomings. Theapproximation of the smoothed stellar core plays a larger rolethan for purely convective (and pulsating) flows that are not ro-tating. Would the angular momentum in a real star be advectedeven further in and leave only a very slowly rotating convectiveenvelope and atmosphere behind? What role would magneticfields play in coupling the interior to the convective envelope?The outer boundaries might also influence the results becausethe slowly rotating atmosphere moves with respect to the com-putational box and might exchange angular momentum throughthe boundaries. Clearly, improved models need a larger compu-tational domain and update of the treatment of the stellar core.Still, the presented models support the results of Dorfi &Höfner (1996) that rotation in AGB stars – if the period is onthe order of 20 yr or shorter – can influence the atmospheric ve-locity field and the wind, and might be responsible for the shapesof some planetary nebulae.
4. Discussion
The current models have an improved resolution over those pre-sented in Freytag & Höfner (2008) because of an increased num-ber of grid points and a higher order reconstruction scheme of thehydrodynamics solver (Freytag 2013). While sub-surface con-vection cells and above-surface shocks show finer structures,there are no qualitative changes in the general results. Still, thenewer models allow a better separation of the radial pulsationmode from the convective noise and also present a longer timebase for the Fourier analysis.The averaged density stratifications for selected models inFig. 12 give an idea about the size of e ff ects due to changes in nu-merical parameters. The two smallest and oldest models, n02 andn06 (referring to the last three letters of the model names), whichwere already used in Freytag & Höfner (2008), have the sameresolution but di ff erent extensions of the computational box andagree well with each other. The same is true for the pair n24 / n25from the current “test” group of models, indicating that the boxsize does not have a major impact. The curves for n24 and ourstandard model n26 are indistinguishable. The di ff erence is a Article number, page 11 of 14 & A proofs: manuscript no. aaagb3dfirstgrid O • ]-14-12-10-8-6 l og ( ρ /[ g c m - ] ) T eff =2531K, logg=-0.85, R=437 R O • , x box =1243.6 R O • , 127 pts : st28gm06n02T eff =2538K, logg=-0.83, R=429 R O • , x box =1674.4 R O • , 171 pts : st28gm06n06T eff =2687K, logg=-0.73, R=383 R O • , x box =1380.6 R O • , 281 pts : st28gm06n13T eff =2616K, logg=-0.76, R=394 R O • , x box =1380.6 R O • , 401 pts : st28gm06n16T eff =2733K, logg=-0.71, R=371 R O • , x box =1380.6 R O • , 281 pts : st28gm06n24T eff =2727K, logg=-0.71, R=371 R O • , x box =1970.2 R O • , 401 pts : st28gm06n25T eff =2737K, logg=-0.70, R=370 R O • , x box =1380.6 R O • , 281 pts : st28gm06n26 Fig. 12.
Logarithm of the gas density averaged over spherical shells andtime for six models. The ordering in the legend essentially correspondsto the order of the curves. The plus signs roughly indicate the depth with τ Ross =
1. In the legend, the e ff ective temperature, surface gravity, edgelength of the computational box, number of grid points, and the nameare given for each model. change in the version of the code, which mainly entails a modifi-cation of the WENO scheme used for the high-order reconstruc-tion in the hydrodynamics solver. We implemented a reductionof the order of the reconstruction, and therefore an increase indissipation, in case of large Mach numbers or large gradients inpressure or entropy. The numerical resolution of model n16 wasincreased by a factor of 1.43 in each direction, compared to n13with only a slight e ff ect on the mean-density stratification.The noticeable decrease of the density by about one order ofmagnitude in the outer layers from n13 / n16 to n24 / n25 / n26 iscaused by a change in the outer boundary conditions. Assuminga less steep density decrease while filling the ghost cells at theboundaries of the computational box induces a stronger infall ofmaterial in the phase between outward moving shocks and some-what impedes propagation of the shocks, leading to lower aver-aged densities in the outer layers. An additional reduction mightbe due to a di ff erence in the treatment of the velocity damping inthe core region of the models, which prevents large-scale dipo-lar flows that could develop and dominate the entire convectiveenvelope. Between the old models n02 / n06 and the new models,there are di ff erences in the envelope mass M env and too manychanges in the numerics to allow a disentanglement of the influ-ence of individual settings. The bottom left panel in Fig. 10 shows the ratio of the ra-dial component of the dynamical pressure and the gas pres-sure, both averaged over spheres and time, (cid:104) ρ v (cid:105) Ω , t / (cid:104) P (cid:105) Ω , t .From the total dynamical pressure, we derive a density-weighted rms value of the radial velocities according to( (cid:104) ρ v (cid:105) Ω , t / (cid:104) ρ (cid:105) Ω , t ) / , and the (smaller) contribution by the ra-dial pulsations ( (cid:104) ( (cid:104) ρ v radial (cid:105) Ω / (cid:104) ρ (cid:105) Ω ) (cid:105) t / (cid:104) ρ (cid:105) Ω , t ) / , plotted in the bot-tom right and the middle right panels, respectively.In the convection zone, where the pressure ratio is below one,it already stays above 20 % for a number of models for a largefraction of the radius. Thus, the dynamical pressure is not negli-gibly small and influences the stratification at least of the outer convective envelope. The peak of the radial velocities near thesurface of the convection zone is accompanied by a peak in thepressure ratio and the dynamical pressure becomes even largerthan the gas pressure.Further out in the atmosphere, the dynamical pressure dom-inates over the gas pressure (radiation pressure is not includedin the current models) by factor of 5 to 10. It significantly in-creases the density and pressure scale heights compared to thecase with only gas pressure, causing a levitation of dense mate-rial into cool layers. This allows molecules to form and possiblycreates a highly irregular, non-spherical, and dynamic versionof a MOLsphere (an extended sphere around the star that is op-tically thick in molecular lines, for instance, due to water; seeTsuji 2000). In addition, the conditions become even su ffi cientlycool and dense for dust to form (see, e.g., Freytag & Höfner2008). The contributions of purely radially symmetric motions(middle right panel in Fig. 10) and spatially fluctuating flows tothe total radial velocities (bottom right panel in Fig. 10) are bothsignificant, but only the radially symmetric motions can be ac-counted for by dynamical 1D models.While we get very extended atmospheres in our models ofAGB stars, Arroyo-Torres et al. (2015) concluded that for red su-pergiant stars (with much larger masses in the range 5 – 25 M (cid:12) )the dynamical pressure in pulsating 1D models and the evenlarger dynamical pressure in 3D CO5BOLD models is not su ffi -cient to enlarge the photosphere to the observed sizes. We presented some observations (e.g., by VLT / SPHERE, HST,VLTI, and various other interferometers, or with the Cassinispacecraft) of asymmetries and clumps in the dust envelopes ofnear-by AGB stars in the Introduction.The bolometric-intensity maps in Figs. 2 and 3 derived fromthe 3D models show that the smallest scale patterns change ontimescales of less than a month, while intensity changes of largerareas occur on timescales of about a year.The surface of the normal stellar convection zone (seeSect. 3.2) sits too deep to directly a ff ect the emergent intensity:surface granules and larger convection cells themselves are notobservable. Instead, the visible structures are caused by shockson various scales. However, as described in Sect. 3.4, the shocksare shaped by the underlying convective structures. A dimmingand brightening of a large area (see Fig. 3) might well indirectlyreflect the dynamics of the convection.A detailed comparison of results from the 3D models withobservations has to await simulations with non-gray opacitiesand a detailed treatment of dust as well as time series of high-angular-resolution observations. A comparison of the pulsation period with the various typicaltimescales of convective structures (see Sect. 3.2) gives t granule (cid:28) t downdraft ∼ t turnover (cid:46) P puls (cid:28) t giant − cell , (2)i.e., that the granules adjust to the pulsation whereas the giantcells are more or less frozen in. However, a restructuring of thegiant cells can lead to a variation of the detailed pulsation be-havior on rather long timescales of several years. The strongestinteraction is to be expected between pulsations and downdrafts.There is no such thing as “the” convective timescale. Article number, page 12 of 14. Freytag et al.: Global 3D radiation-hydrodynamics models of AGB stars.
5. Conclusions
We presented a first exploratory grid of 3D radiation-hydrodynamics models of AGB stars computed with a new ver-sion of CO5BOLD, which is improved in terms of accuracy, sta-bility, and boundary treatment compared to the version used forthe simulations presented in Freytag & Höfner (2008). The in-creased e ff ective resolution leads to additional finer structuresin the convective flow (surface granules and deep turbulent ed-dies) and in near-surface shocks. However, there is no significantchange in the dynamical behavior of the models and only a smallone in quantities that are averaged spatially (over spheres) andtemporally.Several interacting processes govern the dynamics of the at-mosphere of the model AGB stars: non-stationary convectionmanifests as giant convective cells, which change topology onvery long timescales, and short-lived small surface granules. Thegiant convections cells are outlined by turbulent downdrafts thatreach from the surface of the convection zone down into the verycore of the model. Convective motions emit acoustic noise (i.e.,waves, that get trapped inside the star, giving rise to standingacoustic modes, comparable to solar p modes) and they shape,i.e., distort, wave fronts. In addition, there are large-amplitude,radial, fundamental-mode pulsations. The small-scale acousticwaves steepen when they reach the thin cool atmosphere andturn into shocks. The shocks interact and merge, so that the scaleof the atmospheric shocks increases with radial distance, from asmall-scale shock network close to the surface of the convec-tion zone to distorted but almost global, more or less radiallyexpanding shock fronts in the outer layers. The cycles of out-ward moving shocks and material falling back toward the starhave a longer period than the pulsations themselves.The radial pulsations have realistic properties in spite of thecrude treatment of the stellar core in the models. The modelsreproduce the correct period for a given luminosity compared toobservation, if we chose an appropriate ratio of envelope mass tototal stellar mass. The radius of the 3D models is however largerfor a given period compared to previously found period-radiusrelations from 1D pulsation models. The reason for this is notentirely clear, as the 3D models might appear larger owing to amore extended atmosphere inflated by the dynamical pressure ofsmall-scale shocks, or the di ff erence in the representation of theinterior (description of convective energy flux, dynamical pres-sure in the interior, treatment of the stellar core, etc.) might playa role. Higher gravity models have a clearly defined pulsationperiod, whereas lower gravity objects show a much more irreg-ular behavior, depending on the relative size of the convectioncells and the typical convective flow speed.The convective cells themselves do not reach out into visi-ble layers. However, the network of shocks propagating into the(partially convectively unstable) atmosphere gives rise to short-lived spatial inhomogeneities across the stellar disk, which mightbe the cause for observed dynamical features.In the convection zone, the dynamical pressure alreadyreaches significant values ( P dyn / P gas ∼ Acknowledgements.
This work has been supported by the Swedish ResearchCouncil (Vetenskapsrådet). The computations were performed on resources(“milou”) provided by SNIC through Uppsala Multidisciplinary Center for Ad-vanced Computational Science (UPPMAX) under Project p2013234. We thankKjell Eriksson for his helpful comments on the manuscript.
References
Arroyo-Torres, B., Wittkowski, M., Chiavassa, A., et al. 2015, A&A, 575, A50Bedding, T. R., Kiss, L. L., Kjeldsen, H., et al. 2005, MNRAS, 361, 1375Beeck, B., Cameron, R. H., Reiners, A., & Schüssler, M. 2013, A&A, 558, A48Bladh, S., Höfner, S., Aringer, B., & Eriksson, K. 2015, A&A, 575, A105Bowen, G. H. 1988, ApJ, 329, 299Chandler, A. A., Tatebe, K., Wishnow, E. H., Hale, D. D. S., & Townes, C. H.2007, ApJ, 670, 1347Chiavassa, A., Freytag, B., Masseron, T., & Plez, B. 2011a, A&A, 535, A22Chiavassa, A., Pasquato, E., Jorissen, A., et al. 2011b, A&A, 528, A120 + Chiavassa, A., Plez, B., Josselin, E., & Freytag, B. 2009, A&A, 506, 1351Christensen-Dalsgaard, J., Kjeldsen, H., & Mattei, J. A. 2001, ApJ, 562, L141Dorfi, E. A. & Höfner, S. 1996, A&A, 313, 605Eriksson, K., Nowotny, W., Höfner, S., Aringer, B., & Wachter, A. 2014, A&A,566, A95Fleischer, A. J., Gauger, A., & Sedlmayr, E. 1992, A&A, 266, 321Fox, M. W. & Wood, P. R. 1982, ApJ, 259, 198Freytag, B. 2013, Mem. Soc. Astron. Italiana Suppl., 24, 26Freytag, B. & Chiavassa, A. 2013, in EAS Publications Series, Vol. 60, EASPublications Series, ed. P. Kervella, T. Le Bertre, & G. Perrin, 137–144Freytag, B. & Höfner, S. 2008, A&A, 483, 571Freytag, B., Ste ff en, M., & Dorch, B. 2002, Astronomische Nachrichten, 323,213Freytag, B., Ste ff en, M., Ludwig, H.-G., et al. 2012, J.Comp.Phys., 231, 919Hani ff , C. A. & Buscher, D. F. 1998, A&A, 334, L5Hauschildt, P. H., Baron, E., & Allard, F. 1997, ApJ, 483, 390Höfner, S. 2008, A&A, 491, L1Höfner, S. 2015, in Astronomical Society of the Pacific Conference Series, Vol.497, Why Galaxies Care about AGB Stars III: A Closer Look in Space andTime, ed. F. Kerschbaum, R. F. Wing, & J. Hron, 333Höfner, S., Bladh, S., Aringer, B., & Ahuja, R. 2016, A&A, 594, A108Höfner, S., Gautschy-Loidl, R., Aringer, B., & Jørgensen, U. G. 2003, A&A,399, 589Iglesias, C. A., Rogers, F. J., & Wilson, B. G. 1992, ApJ, 397, 717Ireland, M. J., Scholz, M., & Wood, P. R. 2008, MNRAS, 391, 1994Ireland, M. J., Scholz, M., & Wood, P. R. 2011, MNRAS, 418, 114Jeong, K. S., Winters, J. M., Le Bertre, T., & Sedlmayr, E. 2003, A&A, 407, 191Karovska, M., Hack, W., Raymond, J., & Guinan, E. 1997, ApJ, 482, L175Karovska, M., Nisenson, P., Papaliolios, C., & Boyle, R. P. 1991, ApJ, 374, L51Kastner, J. H. & Weintraub, D. A. 1994, ApJ, 434, 719Lattanzio, J. C. & Wood, P. 2004, in “Asymptotic Giant Branch Stars”, eds.Habing H.J., Olofsson H., Springer, 23–104Lighthill, M. J. 1952, Royal Society of London Proceedings Series A, 211, 564Lopez, B., Danchi, W. C., Bester, M., et al. 1997, ApJ, 488, 807Ludwig, H.-G. 2006, A&A, 445, 661Ludwig, H.-G., Ca ff au, E., Ste ff en, M., et al. 2009, Mem. Soc. Astron. Italiana,80, 711Magic, Z., Collet, R., Asplund, M., et al. 2013, A&A, 557, A26Muthsam, H. J., Löw-Baselli, B., Obertscheider, C., et al. 2007, MNRAS, 380,1335Nordlund, Å. & Stein, R. F. 2001, ApJ, 546, 576Nowotny, W., Höfner, S., & Aringer, B. 2010, A&A, 514, A35Ohnaka, K., Weigelt, G., & Hofmann, K.-H. 2016, A&A, 589, A91Ramstedt, S., Mohamed, S., Vlemmings, W. H. T., et al. 2014, A&A, 570, L14Ramstedt, S. & Olofsson, H. 2014, A&A, 566, A145Richichi, A., Percheron, I., & Khristoforova, M. 2005, A&A, 431, 773Samus, N. N., Durlevich, O. V., & et al. 2009, VizieR Online Data Catalog, 1Schwarzschild, M. 1975, ApJ, 195, 137Stein, R. F. & Nordlund, Å. 2001, ApJ, 546, 585Stewart, P. N., Tuthill, P. G., Monnier, J. D., et al. 2016a, MNRAS, 455, 3102Stewart, P. N., Tuthill, P. G., Nicholson, P. D., & Hedman, M. M. 2016b, MN-RAS, 457, 1410Stothers, R. & Leung, K.-C. 1971, A&A, 10, 290 Article number, page 13 of 14 & A proofs: manuscript no. aaagb3dfirstgrid
Trampedach, R., Asplund, M., Collet, R., Nordlund, Å., & Stein, R. F. 2013,ApJ, 769, 18Tremblay, P.-E., Ludwig, H.-G., Freytag, B., et al. 2015, ApJ, 799, 142Tsuji, T. 2000, ApJ, 540, L99van Leeuwen, F. 2007, A&A, 474, 653Wedemeyer, S., Freytag, B., Ste ff en, M., Ludwig, H.-G., & Holweger, H. 2004,A&A, 414, 1121Weigelt, G., Balega, Y., Blöcker, T., et al. 1998, A&A, 333, L51Whitelock, P. A., Feast, M. W., & van Leeuwen, F. 2008, MNRAS, 386, 313Whitelock, P. A., Menzies, J. W., Feast, M. W., et al. 2009, MNRAS, 394, 795Winters, J. M., Le Bertre, T., Jeong, K. S., Helling, C., & Sedlmayr, E. 2000,A&A, 361, 641Woitke, P. 2006, A&A, 452, 537Wood, P. R. 1990, in From Miras to Planetary Nebulae: Which Path for StellarEvolution?, ed. M. O. Mennessier & A. Omont, 67–84Wood, P. R. 2015, MNRAS, 448, 3829Wood, P. R., Alcock, C., Allsman, R. A., et al. 1999, in IAU Symposium,Vol. 191, Asymptotic Giant Branch Stars, ed. T. Le Bertre, A. Lebre, &C. Waelkens, 151en, M., Ludwig, H.-G., & Holweger, H. 2004,A&A, 414, 1121Weigelt, G., Balega, Y., Blöcker, T., et al. 1998, A&A, 333, L51Whitelock, P. A., Feast, M. W., & van Leeuwen, F. 2008, MNRAS, 386, 313Whitelock, P. A., Menzies, J. W., Feast, M. W., et al. 2009, MNRAS, 394, 795Winters, J. M., Le Bertre, T., Jeong, K. S., Helling, C., & Sedlmayr, E. 2000,A&A, 361, 641Woitke, P. 2006, A&A, 452, 537Wood, P. R. 1990, in From Miras to Planetary Nebulae: Which Path for StellarEvolution?, ed. M. O. Mennessier & A. Omont, 67–84Wood, P. R. 2015, MNRAS, 448, 3829Wood, P. R., Alcock, C., Allsman, R. A., et al. 1999, in IAU Symposium,Vol. 191, Asymptotic Giant Branch Stars, ed. T. Le Bertre, A. Lebre, &C. Waelkens, 151