A&A 368, 497-526 (2001)
DOI: 10.1051/0004-6361:20000554

Structure and physical properties of the rapidly evolving dusty envelope of IRC$~\!$+10216 reconstructed by detailed two-dimensional radiative transfer modeling

A. B. Men'shchikov1,3 - Y. Balega2 - T. Blöcker3 - R. Osterbart3 - G. Weigelt3


1 - Stockholm Observatory, 13336 Saltsjöbaden, Sweden
2 - Special Astrophysical Observatory, Nizhnij Arkhyz, 357147 Karachaevo-Cherkesia, Russia
3 - Max-Planck-Institut für Radioastronomie, Auf dem Hügel 69, 53121 Bonn, Germany

Received 13 July 2000 / Accepted 12 December 2000

Abstract
We present the first detailed, two-dimensional radiative transfer model of the dusty envelope around the carbon star IRC$~\!$+10216. Our goal was to find a self-consistent model of the star and its envelope which takes into account as many observational constraints as possible. The model reproduces very well the entire beam-matched spectral energy distribution of IRC$~\!$+10216 from optical to centimeter wavelengths (at several phases of stellar luminosity), observed intensity profiles of the object at 1.25, 2.2, 10.5, 50, 100${\mu }$m, and 1.3mm, a 10.5${\mu }$m lunar occultation intensity profile, our high-resolution J, H, K, and $H\!-\!K$ bispectrum speckle-interferometry images, and visibilities in J, H, K, L, M, and N bands. For the adopted distance of 130pc, the model of IRC$~\!$+10216 implies that the object changes its luminosity between 13000 and 5200$L_{\odot }$, its effective temperature between 2800 and 2500K, and its radius between 500 and 390$R_{\odot}$. There is a dense non-spherical dust shell around the star, with outflow cavities at position angle PA ${\,\approx \,}$20$^\circ $. The southern cavity with a full opening angle of 36$^\circ $ is tilted toward us by 40$^\circ $ from the plane of sky, causing the observed bipolar appearance of the object on a subarcsecond scale. If the envelope's outflow velocity of 15kms-1 applies to the material making up the dense core, then just ${\sim\,}$15 years ago the star was losing mass at a rate of 9 10-5$M_{\odot }$yr-1. Dust exists in the envelope of IRC$~\!$+10216 everywhere from the stellar photosphere up to a distance of 3pc from the star. The total mass of the envelope lost by the central star is 3$M_{\odot }$ and the dust-to-gas mass ratio is 0.004. The total optical depth $\tau_{V}$ toward the star in the visual is 40, in the polar cavities it is 10. The innermost parts of the envelope are optically thick even at 10.7${\mu }$m due to a strong resonance absorption of silicon carbide grains at that wavelength. In addition to SiC dust, the model contains inhomogeneous grains made of a mixture of SiC and incompletely amorphous carbon with thin [Mg0.5Fe0.5]S mantles. This is the simplest dust mixture required to fit all observations of IRC$~\!$+10216 and to correctly interpret the well-known 11.3${\mu }$m and 27${\mu }$m emission bands. The dust model found in this study can also be successfully applied to many other carbon stars exhibiting broad emission features in the 10.3-12.6${\mu }$m and 25-37${\mu }$m wavelength regions. An important and firm result of our modeling is that the brightest compact peak observed in IRC$~\!$+10216 is not the direct light from the underlying central star. In contrast to previous suggestions, the brightest southern component, labeled A in our high-resolution near-infrared images (Weigelt et al. 1998a, Weigelt et al. 1998b; Osterbart et al. 2000), is only the radiation emitted and scattered in the optically thinner southern cavity of the bipolar dense shell moving away from the central star. The carbon star is at the position of the fainter component B in our H and K images, which is 0 $.\!\!^{\prime\prime}$21 away from A along the symmetry axis. Direct stellar light (component B) is not seen at all in the Hubble Space Telescope 0.8${\mu }$m and 1.1${\mu }$m images, being absorbed by the dense dusty material. The even fainter components C and D in the Hand K images are probably due to smaller deviations of the dense shell from the spherical shape. IRC$~\!$+10216 seems to have entered a phase immediately before moving off the asymptotic giant branch and started developing asymmetries in its envelope.


   
1 Introduction

After the first study by Becklin et al. (1969) and a series of photometric, spectroscopic, and polarimetric measurements by others (Lockwood 1970; Miller 1970; Herbig & Zappala 1970; Shawl & Zellner 1970), it has become clear that IRC$~\!$+10216 (CW Leo, AFGL 1381) deserves special attention. An extended object ( $2{\hbox{$^{\prime\prime}$ }}{\times}4{\hbox{$^{\prime\prime}$ }}$) outside the galactic plane, IRC$~\!$+10216 was the brightest source on the sky in the near infrared (IR). Its almost featureless spectral energy distribution (SED) resembled that of a 650K blackbody (or 1200K at wavelengths below 1${\mu }$m); its brightness at 2.2${\mu }$m varied by 2$^{\rm m}$ with a period of ${\sim\,}$600 days; its small proper motion implied a distance $D{\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\displaystyl...
...\offinterlineskip\halign{\hfil$\scriptscriptstyle ...100pc. The 1${\mu }$m linear polarization was in excess of 20% at a position angle PA ${\,\approx \,}$120$^\circ $, perpendicular to the elongation of the optical image (PA ${\,\approx \,}$30$^\circ $), indicating scattering by solid particles in a non-spherical environment. IRC$~\!$+10216 was identified as a longest-period, pulsating carbon star surrounded by an opaque dust shell, perhaps in a state preceding that of a planetary nebula.

During the three decades of studies, dust continuum radiative transfer calculations have been repeatedly applied to this well-observed object (Crabtree & Martin 1979; Mitchell & Robinson 1980; Rowan-Robinson & Harris 1983; Le Bertre 1987; Martin & Rogers 1987; Orofino et al. 1990; Griffin 1990; Lorenz-Martins & Lefèvre 1993; Winters et al. 1994; Sloan & Egan 1995; Bagnulo et al. 1995; Ivezic & Elitzur 1996; Groenewegen 1997) in order to understand in more detail physical, chemical, and evolutionary properties of the star and its dusty molecular envelope. All the models used an assumption of spherical geometry which seemed to be justified for this object by the presence of the large-scale (several arcminutes) apparently symmetric envelope. The latter appears to have, however, not very smooth density distribution in the deep optical images by Crabtree et al. (1987) and Mauron & Huggins (1999), showing a series of incomplete concentric shells at radii of ${\sim }$15-60 $^{\prime \prime }$.

On the other hand, several speckle observations (McCarthy et al. 1980; Mariotti et al. 1983; Ridgway & Keady 1988; Dyck et al. 1991; Danchi et al. 1994) and especially recent highest-resolution observations (Osterbart et al. 1997; Weigelt et al. 1998a,b; Haniff & Buscher 1998; Osterbart et al. 2000; Tuthill et al. 2000) have demonstrated that the inner dust shell of IRC$~\!$+10216 is actually non-spherical and extremely clumpy. This well-established observational fact emphasizes the need of multidimensional radiative transfer calculations for proper interpretation of the observations.

A lot of circumstellar structures seen in the high-resolution images of the last decade demonstrate that the spherical symmetry is no more a realistic basic assumption for understanding the observations. Clumpy and bipolar structures (often intrinsically related to the jets and outflows) are being widely observed in protostars, young stellar objects (YSO), protoplanetary nebulae, in the winds of evolved post-asymptotic giant branch (post-AGB) stars, and in active galactic nuclei. An axially symmetric geometry of an optically thick toroidal distribution of matter around a central energy source emerges as a good basic approximation for the modeling of the physically different non-spherical objects.

From a practical point of view, axisymmetric geometries are much more time-consuming for modeling, because of two additional parameters, the density structure in the polar direction and the viewing angle. Furthermore, to find a self-consistent and realistic model, one has to derive also the properties of dust grains, which makes the whole problem enormously difficult and requiring many more observational data sets to constrain the modeling. In a series of recent papers (Men'shchikov & Henning 1997; Men'shchikov et al. 1998; Men'shchikov et al. 1999; White et al. 2000), we applied our two-dimensional radiative transfer code to detailed modeling of toroidal envelopes around YSOs (L1551 IRS5, HL Tau) and around a close binary post-AGB star (the Red Rectangle). The models demonstrate that a derivation of reliable physical parameters requires using all available (spatial) information in a detailed radiative transfer modeling. The less constraints are used, the easier is to find a model and the higher is the risk that the model is unrealistic (Men'shchikov & Henning 2000). The common practice of sketchy fitting of SEDs alone can only give source parameters that are highly uncertain.

In the present study, we utilized our approach and two-dimensional radiative transfer code to find a self-consistent model of the dusty environment of IRC$~\!$+10216, which would match quantitatively all available dust continuum observations. In Sect. 2, we briefly describe the most substantial results obtained for IRC$~\!$+10216 during the three decades of very active studies, that will enable the reader synthesize a large amount of information and follow easier our subsequent discussion. In Sect. 3, all details of our approach and radiative transfer model are described. In Sect. 4, we confront the model with available observations of IRC$~\!$+10216. In Sect. 5, we briefly touch upon its present evolutionary state. In Sect. 6, we summarize results obtained in this study.

   
2 IRC$~\!$+10216: Three decades of observations

In the years following the initial discovery, a great number of observations in a wider wavelength range and with higher sensitivity, spectral and angular resolution painted quite a detailed and complex (yet incomplete) picture of this carbon star.

   
2.1 Properties of the envelope

GEOMETRY. Although on scales of more than several arcseconds the envelope of IRC$~\!$+10216 seems to have a roughly spherical shape, large optical and near-IR polarization measured by Dyck et al. (1971) and Capps & Knacke (1976) suggested significant deviations from spherical symmetry. Near-IR speckle interferometry by Mariotti et al. (1983), Dyck et al. (1984), and Ridgway & Keady (1988) revealed asymmetries in the envelope on subarcsecond scales and two scattering lobes of a lower density. Polarimetric imaging by Tamura et al. (1988) detected an extended envelope also in highly polarized light. A weakly-bipolar reflection nebula (Kastner & Weintraub 1994) or a geometrically thick disk (Dyck et al. 1987) or a dusty torus (Skinner et al. 1998) must be significantly inclined toward an observer. A complex, clumpy subarcsecond structure in the near IR is changing on time scales of a few years, as has been demonstrated by interferometric observations of Weigelt et al. (1998a) and Haniff & Buscher (1998). Osterbart et al. (2000) concluded that it is the southern scattering lobe that dominates the optical Hubble Space Telescope (HST) and near-IR images and polarization pattern. See also Sects. 3.2, 4.5.1, 4.5.2, 4.6, and 6.

DENSITY PROFILE. Near-IR interferometry observations by McCarthy et al. (1980) indicated large departures from the expected smooth $\rho{\,\propto\,}r^{-2}$ density distribution on an arcsecond scale. Deviations indicating an enhancement by a factor of 3 were found also over much larger distances by Bieging et al. (1984) and a flatter density distribution $\rho{\,\propto\,}(1+r/25{\hbox{$^{\prime\prime}$ }})\,r^{-2}$ was deduced by Fazio et al. (1980). Having analysed far-IR scans of the envelope on large scales, Harvey et al. (1991) concluded that dust density distribution is inconsistent with the picture of a smooth, constant-velocity outflow. Keady & Ridgway (1993) presented evidence for by a factor of 2-3 higher densities ${\sim }$1000 years ago. Two periods of enhanced density ${\sim }$200 and 850 years ago were identified by Groenewegen et al. (1997). Direct evidence of periodic outbursts with higher densities on a time scale of 200-800 years was found in deep optical images by Groenewegen et al. (1987) and Mauron & Huggin (1999). See also Sects. 4.3, 4.4, 4.5.3, and 6.

MASS-LOSS RATE. From a radiative transfer modeling of the CO line emission, Groenewegen et al. (1998) derived a mass-loss rate of 1.5 10-4$M_{\odot }$yr-1. A similar estimate of 2 10-5$M_{\odot }$yr-1 was obtained by Keady et al. (1988) and Kastner (1992). Somewhat higher mass-loss rate of 3.25 10-5$M_{\odot }$yr-1 was derived by Crosas & Menten (1997) from their modeling of the CO radiation, whereas Truong-Bach et al. (1991) obtained 4 10-5$M_{\odot }$yr-1 for the outer envelope. Molecular spectroscopy and a spectral synthesis modeling allowed Keady & Ridgway (1993) to derive an upper limit of 4 10-5$M_{\odot }$yr-1. Sahai (1987) estimated 4.8 10-5$M_{\odot }$yr-1 in the outer envelope, whereas Knapp & Morris (1985) obtained 5.5 10-5$M_{\odot }$yr-1. Note that the above estimates are based on different distances and need to be rescaled (see Table 5). See also Sects. 4.4 and 6.

OUTFLOW VELOCITY. From molecular line observations, Morris et al. (1975) derived an expansion velocity of the envelope of 12kms-1. Betz et al. (1979) detected NH3 in the inner dense region, concluding that the gas flow was accelerated to 14kms-1 on subarcsecond scales. Olofsson et al. (1982) observed several molecules in the envelope with an outflow velocity of 14.4kms-1. Kuiper et al. (1976) resolved a large CO envleope expanding at 15kms-1 and Knapp & Morris (1985) obtained almost the same velocity of 15.2kms-1, whereas Knapp et al. (1982) measured an outflow velocity of 17kms-1. See also Sect. 3.2 and Appendix A.

   
2.2 Properties of dust grains

COMPOSITION. A broad emission feature at about 11${\mu }$m in the spectra of IRC$~\!$+10216 obtained by Treffers & Cohen (1974) indicated that silicon carbide dust must be present, as well as other carbonaceous particles with featureless emissivity. A broad far-IR emission band beginning at 24${\mu }$m was discovered by Forrest et al. (1981). It was shown to extend beyond 37${\mu }$m and peak at ${\sim }$30${\mu }$m by Herter et al. (1982). Magnesium sulfide dust grains were proposed by Goebel & Moseley (1985) as the carrier of the band. See also Sects. 3.3.1, 3.3.4, 3.3.5, 3.3.6, 4.2, and 6.

SIZES. Spectrophotometry of the 11${\mu }$m feature by Treffers & Cohen (1974) showed that there exist small grains with radii $a{\,\ll\,}$1${\mu }$m. Near-IR low-resolution spectra obtained by Witteborn et al. (1980) indicated a presence of grains with radii $a{\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\displaystyl...
...\offinterlineskip\halign{\hfil$\scriptscriptstyle ...0.25${\mu }$m. Jura (1983) argued that grains are rather small, having $a{\,\approx\,}$0.05${\mu }$m, and later Jura (1994) proposed that dust in the outer envelope has a wide distribution of sizes, including both large and small grains ( $a{\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\displaystyl...
...\offinterlineskip\halign{\hfil$\scriptscriptstyle ...0.015${\mu }$m, $a{\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\displaystyl...
...\offinterlineskip\halign{\hfil$\scriptscriptstyle ...0.5${\mu }$m). See also Sects. 3.3.1, 3.3.6, and 4.2.

CONDENSATION. Modeling of dust formation by McCabe (1982) demonstrated that small SiC grains can condense very close to the stellar photosphere ( $r{\,\approx\,}$1.5$R_{\star}$, $T{\,\approx\,}$1500K), whereas graphite grains form at lower temperatures further away from the star (5$R_{\star}$ ${\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\displaystyle...
...\offinterlineskip\halign{\hfil$\scriptscriptstyle ...20$R_{\star}$, 1000K ${\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\displaystyle...
...\offinterlineskip\halign{\hfil$\scriptscriptstyle ...600K). Danchi et al. (1990) interpreted interferometry data as indicating that dust grains form at higher temperatures (1000-1800K), and closer to the star (1.5-3$R_{\star}$) than usually assumed. See also Sects. 3.3.2, 4.3, and 6.

EMISSIVITY. Submillimeter photometry by Phillips et al. (1982) showed that the far-IR opacity of the grains is close to $\kappa{\,\propto\,}{\lambda}^{-1}$. Comparing the far-IR and ultraviolet extinction of dust in the envelope, Jura (1983) argued that $\kappa{\,\propto\,}{\lambda}^{-1.3}$. Photometry in the far IR by Rengarajan et al. (1985) and at submm wavelengths by Sopka et al. (1985) confirmed that the dust emissivity is rather flat at long wavelengths ( $\kappa{\,\propto\,}{\lambda}^{-1}$ to ${\lambda}^{-1.2}$). A similar result was also found by Le Bertre (1987) from his radiative transfer fitting of the observed SED. See also Sects. 3.3.3 and 4.2.

   
2.3 Uncomfortable reality

Huge amounts of detailed observational information for IRC$~\!$+10216 present a real challenge to theoretical models attempting to explain the data in a consistent way. In addition to the usual difficulties and uncertainties of modeling and interpretation, the star and its envelope have significantly varied in time due to both periodic pulsations and non-periodic, long-term changes. What we observe is a complex result of a very long evolution of the central star, with its history recorded in the structure and chemical composition of the parsec-sized envelope.

Many usual simplifying assumptions are invalid for IRC$~\!$+10216: (1) we clearly have a non-stationary, time-dependent problem, (2) spherical symmetry cannot be adopted for dense inner regions and even axial symmetry is questionable there, (3) the mass-loss rate and the wind velocity are not constant, (4) the envelope's density and chemical composition cannot be described by simple power laws, (5) dust grains are very likely to have inhomogeneous structure and composition, non-spherical shapes, and different size distributions in various parts of the envelope, (6) the grains are accelerated by the radiation pressure to different velocities, which leads to their collisions and possible growth by coagulation or to destruction in the dense inner regions, (7) the observational constraints used in models are not coeval.

The last point is potentially very troublesome, because of the observed non-periodic changes which imply that the data from different epochs may well be fundamentally incomparable. An important question is how far one can go making simplifying assumptions while modeling very complex phenomena. The finer details are observed, the more simplifications we have to abandon, and the more complex the models should become.

   
3 IRC$~\!$+10216: The radiative transfer model

   
3.1 General formulation

Having mentioned many difficulties which may complicate thorough and detailed analyses, we are going to make a number of assumptions and simplifications, which we believe are reasonable, and present a model for IRC$~\!$+10216, which takes into account all existing dust continuum observations. Our goal is to find a snapshot structure of the dusty envelope, consistent with the recent high-resolution optical and near-IR images of the object (Weigelt et al. 1998a,b; Haniff & Buscher 1998; Osterbart et al. 2000), as well as with most other measurements of the dust radiation. Although our model describes only the present structure of the object (and as such is time-independent), in Sect. 5 we consider its evolutionary implications. The complex molecular chemistry of the envelope is largely ignored here and the gas component is described by means of a spatially constant dust-to-gas mass ratio. The assumed constancy of ${\rho _{\rm d}/\rho }$ may not be realistic enough, thus contributing to the uncertainties of the derived properties of the gaseous envelope (e.g., the total mass M of the envelope, lost by the central star during its evolution).

We utilized our 2D radiative transfer code based on a ray-tracing method which provides an accurate solution to the frequency-dependent radiative transfer problem, including isotropic scattering (we refer to Men'shchikov & Henning 1997, for a detailed description). Large parameter space was explored in many hundreds of runs by changing model parameters and comparing their SEDs, images, and visibilities to available observational constraints (see Appendix B for details).

   
3.2 Geometry and assumptions

The model assumes axially-symmetric geometry shown schematically in Figs. 1 and 2. Departures from spherical symmetry appear only on a subarcsecond scale, in the inner dense core. In our recent images (Osterbart et al. 2000) the bipolar holes appear very small, comparable to the stellar diameter (see also Sect. 4.5). The perfect conical shape of the polar cavities should be considered as only a first approximation to reality. In this 2D radiative transfer modeling, we had to ignore the possibility that the immediate environment of IRC$~\!$+10216 may be three-dimensional. One needs to keep this in mind, however, when comparing our results with observations.

Although it may contradict intuitive perception, in our model the star is located at the position of the fainter component B, whereas the brightest peak A is produced by the hot dust emission from the far side of the outflow cavity (Fig. 1). The even fainter components C and D may be, for example, optically thinner zones due to density fluctuations along a ring surrounding the outflow cavity (Fig. 2). High-resolution observations, in particular a 1.1${\mu }$m polarization map and the components' proper motions presented in our previous paper (Osterbart et al. 2000), support this picture. We have shown that the most natural and symmetric velocity field can be derived if one assumes the star at the position of B. In fact, only in this case, both symmetrically located components C and D are moving away from B (within the plane of sky) with the same apparent velocity $v_{\rm C}{\,\approx\,}v_{\rm
D}{\,\approx\,}$5kms-1; the component A is moving away from the star with an apparent velocity of $v_{\rm
A}{\,\approx\,}$14kms-1 (Table 2 in Osterbart et al. 2000). A velocity being different by factor of 3 is most likely a pure projection effect implying that the radius-vectors $r_{\rm C}$ and $r_{\rm D}$ to C and D are more inclined toward the observer than is the vector $r_{\rm A}$ to component A (Fig. 2). The actual deprojected radial velocity $v_{\rm r}{\,\approx\,}$15kms-1 (see Appendix A) is the same for all the components and equal to the general expansion velocity of the envelope around IRC$~\!$+10216.

It is instructive to list here our general simplifying assumptions, which are relevant (and essential) for the model: (1) dust distribution is axially-symmetric close to the star, (2) outer parts of the envelope are spherically-symmetric, (3) dust density depends only on the radial coordinate, (4) dust population consists of spherical, compact solid grains, (5) size distribution of the grains can be described by a power law, (6) wherever a dust component exists, its composition, structure, and grain sizes are spatially invariant, (7) light scattering by the dust particles is isotropic,

  \begin{figure}
\par\resizebox{\hsize}{!}{\includegraphics{h232401.eps}}\end{figure} Figure 1: Geometry of the circumstellar envelope of IRC$~\!$+10216. Schematically shown are three regions of the model envelope - the innermost dense core with bipolar cavities (dark color), the less denser envelope where molecules are observed (medium color), and the outer extended envelope (light color). The model is basically spherically-symmetric, except for the inner axially-symmetric region with bipolar outflow cavities (see Fig. 2). The geometry is defined by the opening angle of the cavities, $\omega {\,=\,}\pi {\,-\,}\psi {\,\approx \,}$36$^\circ $ ( $\psi {\,\approx \,}$144$^\circ $) and the viewing angle, $\theta _{\rm v}{\,\approx \,}$40$^\circ $, between the equatorial plane and the line of sight. The bright components A and B labeled by Weigelt et al. (1998a) correspond to those shown in Fig. 2
Open with DEXTER


  \begin{figure}
\par\resizebox{\hsize}{!}{\includegraphics{couv-new.eps}}\end{figure} Figure 2: Three-dimensional representation of the innermost dense core of the circumstellar envelope of IRC$~\!$+10216 as it appears in the highest-resolution near-IR speckle images in projection onto the sky plane. Components A, B, C, D observed by Weigelt et al. (1998a) are shown schematically, as well as much fainter components E, F, G found by Osterbart et al. (2000). Contrary to intuitive expectations, our model puts the star at the position of the fainter component B. The model identifies the brightest peak A with the light escaping from the southern cavity. The components C and D (not modeled here) may appear on a ring around the cavity most likely due to dust density (optical depth) fluctuations along the ring. Weak emission from the obscured northern cavity and inhomogeneous dusty environment cause a somewhat irregular appearance of the components E, F, and G in the speckle images
Open with DEXTER

(8) dust-to-gas mass ratio of a dust component is spatially homogeneous, (9) dust is in radiative equilibrium with the radiation field, (10) there are no sources or sinks of radiative energy in the envelope, (11) the contribution of molecular line emission to the observed fluxes is small, (12) the star radiates as a blackbody at long wavelengths.

There are no good reasons to believe that the above statements are extremely good approximations. One cannot avoid using them, however, primarily because of the obvious lack of sufficiently detailed and reliable observational constraints. The influence of the assumptions on model results is quite uncertain apriori and, therefore, it must be carefully investigated.

   
3.3 Dust particles

   
3.3.1 Chemical composition

Spectrophotometric observations of IRC$~\!$+10216 show an almost featureless continuum with only a couple of clearly detected broad bands (Sect. 2.2). The latter are the so-called 11.3${\mu }$m emission attributed to SiC grains (Treffers & Cohen 1974) and a very broad emission band at ${\sim }$30${\mu }$m, which is usually thought to be produced by MgS grains (Goebel & Moseley 1985). Also clearly visible is a 3.1${\mu }$m absorption feature attributed to the absorption by molecules (C2H2, HCN) in stellar photospheres (Ridgway et al. 1978).

 

 
Table 1: Model estimates of the grain temperature T1 at the inner dust boundary R1, dust grain size a and power-law exponent p (from $Q_{\rm abs}{\,\propto \,}{\lambda }^{-p}$) adopted in previous models of IRC$~\!$+10216. Bracketed values are uncertainties of the parameters given in the referenced papers. Fifth column roughly classifies chemical composition of dust components used in the models: "Gr'' graphite, "amC'' amorphous carbon, "SiC'' silicon carbide, "MgS'' magnesium sulfide

T1(K)
R1($R_{\star}$) a(${\mu }$m) p Composition Reference

1700
    2 Gr Crabtree & Martin (1979)
1500 1.5     Gr+SiC McCabe (1982)
1500 1.9 $\le$0.05 1 amC+SiC Griffin (1990)
1500 3.4 $\le$0.25   amC+SiC+MgS Skinner et al. (1999)
1300 2.5-3       Danchi et al. (1990)
1300 2.5     amC Winters et al. (1994)
1100 2.2-3 0.008 1.3 amC+SiC Bagnulo et al. (1995)
1075(5%) 4.5(10%) 0.16(6%)   amC+SiC Groenewege (1997)
1000 4.0     amC+SiC Sloan & Egan (1995)
1000 4.5 0.05 1 amC Martin & Rogers (1987)
950(25%) 6.6   1.3   Le Bertre (1987)
850 6 0.05 1.2 amC+SiC Lorenz-Martins & Lefèvre (1993)
750 5.6 0.1 1 amC Rowan-Robinson & Harris (1983)
750(7%) 5.1-8.5 $\le$0.2   amC+SiC Ivezic & Elitzur (1996)
600 20 0.05 2 Gr+SiC Mitchell & Robinson (1980)
    0.01   amC+SiC Sloan & Egan (1995)

The continuum is undoubtedly produced by some kind of carbonaceous solid particles, although details of the chemical composition, internal structure of the dust material, and of the condensation process itself are still a matter of debate. Most frequently, the observations have been interpreted in terms of amorphous carbon grains (which are conglomerates of highly disordered planar graphitic structures), although it seems very unlikely that any single form of carbon may be able to explain all observations. Moreover, there are good reasons to believe that properties of circumstellar dust are extremely complex.

In the present modeling we assumed that dust is a mixture of compact spherical solid particles of radius a, composed of inhomogeneous materials which include different forms of carbon, silicon carbide, and magnesium-iron sulfides. Although polycyclic aromatic hydrocarbon (PAH) molecules may also exist in the carbon-rich envelope around IRC$~\!$+10216, the observations show no evidence of their emission features. The three components considered in this work constitute the simplest realistic set of materials; the addition of other solids or PAHs would not be justified.

The radiative transfer code employed in this study (Men'shchikov & Henning 1997) treats separately an arbitrary number of dust components which may differ by the composition, structure, shapes and sizes of grains. Following Kim et al. (1994) and Jura (1994), we adopted a size distribution in the form ${\rm
d}n/{\rm d}a{\,\propto\,}a^{\gamma}{\,\exp\,}(-a/a_{\rm exp})$, which is defined by the minimum grain radius and by the exponential cutoff radius for the largest particles ( $a_{\rm min}, a_{\rm exp}$). Although we explored a very large range of various grain parameters and mixtures in the modeling, unknown details of the dust properties in IRC$~\!$+10216 prevented us from using more complex dust models in final runs.

   
3.3.2 Condensation

Laboratory experiments on dust condensation (Frenklach et al. 1989) suggest that SiC may be the first material to nucleate from the gas phase at fairly high temperatures $T_{\rm g}{\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\dis...
...\offinterlineskip\halign{\hfil$\scriptscriptstyle ...2000K, i.e. very close to (or even inside of) the stellar photosphere. At lower temperatures of $T_{\rm g}{\,\sim\,}$1500K (distances r of only a few $R_{\star}$) these very small SiC particles may provide sites for their further heterogeneous growth by the deposition of amorphous carbon (diamond, for $T_{\rm g}{\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\dis...
...\offinterlineskip\halign{\hfil$\scriptscriptstyle ...1300K) layers. When temperatures in the stellar wind drop to $T_{\rm g}{\,\sim\,}$1000-900K, PAH molecules can form in the gas phase (Frenklach & Feigelson 1989). Their clustering may produce more nucleation sites, provided that chemical, thermal, and dynamic conditions in the wind are favorable for the formation of high enough number density of PAHs. The existing SiC-C core-mantle grains could also grow substantially by the deposition of PAHs on their surfaces. At even larger distances from the star ( $r{\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\displaystyl...
...\offinterlineskip\halign{\hfil$\scriptscriptstyle ...10$R_{\star}$), MgS could condense at temperatures $T_{\rm g}{\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\dis...
...\offinterlineskip\halign{\hfil$\scriptscriptstyle ...800K, as indicated by nucleation calculations (Gail & Sedlmayr 1986).

Although laboratory experiments and theory suggest that condensation temperatures as high as 2000K are possible for the abundant refractory C and SiC dust, radiative transfer models have been adopting much lower temperatures T1 at the inner dust boundary R1. It is clear from Table 1, however, that the observationally based model estimates of the dust temperature and the dust-free zone radius in IRC$~\!$+10216 also substantially differ from each other. Note a very large scatter, by factors of ${\sim }$3-4 in the derived values for the two key model parameters. It is reasonable to consider the factors as rough estimates of the absolute uncertainties, or "error bars'', of typical simplified models that fit only SEDs. Less abundant (or more transparent) refractory components may condense at higher temperatures, closer to the star than the models assume, yet their presence in the "dust-free'' zone is almost impossible to deduce from simple fitting of observed SEDs.

As usual, our model assumes that dust appears instantaneously at a distance R1 from the star, where the gas temperature $T_{\rm g}$equals to the expected nucleation temperature $T_{\rm nuc}$ of the adopted dust component. In a multicomponent approach utilized in the present study, different components of dust grains (having different abundances) form in a considerable range of temperatures and distances (factors of ${\sim }$3 and 10, correspondingly). The stellar surface oscillates in our model between 1.8 and 2.3AU with the pulsational period of the star. We put the inner boundary of the dusty envelope as close as possible to the star, at $R_{1}{\,\approx\,}$3.3AU ${\,\approx\,}1.6\,R_{\star}$, and separated it from the latter for purely numerical reasons. This particular choice of R1 does not influence the model results in any way, because the innermost regions of the envelope ( $r{\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\displaystyl...
...\offinterlineskip\halign{\hfil$\scriptscriptstyle ...20AU) populated by SiC dust are optically thin.

Condensation temperature of SiC grains, the most refractory dust component ( $T_{\rm SiC}{\,\approx\,}$2000K), is higher than their radiative equilibrium temperatures at the inner boundary ( $T{\,\approx\,}$1000K), suggesting that SiC dust can indeed nucleate even inside the stellar photosphere. At larger distances, $R_{\rm
[SiC,C]}\approx 30$AU ${\,\approx\,}15\,R_{\star}$, carbonaceous dust condenses in the stellar wind ( $T_{\rm [SiC,C]}{\,\approx\,}$900K), forming [SiC,C] conglomerates. The formation zone of magnesium-iron sulfides, the least refractory material ( $T_{\rm [Mg,Fe]S}
{\,\approx\,}$600 K), is located further downstream, at $R_{\rm [Mg,Fe]S}{\,\approx\,}$45AU ${\approx\,}23\,R_{\star}$. Pre-existing [SiC,C] grains serve as the condensation sites for [Mg,Fe]S mantles in our model.

   
3.3.3 Optical properties

Most previous radiative transfer calculations tried to fit the SED of IRC$~\!$+10216 by postulating small spherical graphite or amorphous carbon grains, with a power-law long-wavelength emissivity $\kappa{\,\propto\,}{\lambda}^{-p}$ (Table 1). Recent radiative transfer models that took into account near-IR visibilities (in addition to the SED), derived somewhat larger grain sizes, up to $a{\,\sim\,}$0.2${\mu }$m. Some of the studies used also pure SiC grains as a minor additional component necessary to explain the 11.3${\mu }$m emission feature. All of them adopted optical constants of SiC published by Pegourie (1988); see, however, Papoular et al. (1998) for a discussion of problems with derivation of those constants. Lack of published measurements for MgS was the main reason why this dust component has not been included in previous models, despite the proposed association of the broad 30${\mu }$m emission with this material (Sect. 2.2).

In this study, relevant materials were identified using known optical constants of different laboratory analogs of cosmic dust. Mie theory was applied to compute the size- and wavelength-dependent absorption and scattering efficiencies, $Q_{\rm abs}(a,{\lambda})$ and $Q_{\rm
sca}(a,{\lambda})$. We utilized the optical constants of graphite (Draine 1987), amorphous carbon AC1 (Rouleau & Martin 1991), ACH2 (Zubko et al. 1996), carbon samples cel400, cel600, cel800, cel1000 (Jaege et al. 1998), $\alpha $-SiC (Choyke & Palik 1985; Pegourie 1988), and [Mg,Fe] sulfides of various compositions ([MgxFe1-x]S with x=0.9, 0.75, 0.5, 0.1, 0; Begemann et al. 1994).

The optical data of the laboratory samples of [Mg,Fe] sulfides have been measured by Begemann et al. (1994) in the wavelengths range 10-500${\mu }$m. To enable calculations of the radiative equilibrium dust temperatures, we extrapolated the optical constants to the wavelengths 1.2-10${\mu }$m using the Lorentz oscillator model in the limit ${\omega}{\,\gg\,}{\omega}_0$, where ${\omega}_0$ is the resonance frequency (e.g., Bohren & Huffman 1983. The data were extended to the wavelength range 0.3-1.1${\mu }$m using the optical constants of troilite (FeS) measured by Egan & Hilgeman (1977). A similar extrapolation to the wavelength range 17-300${\mu }$m, in the limit ${\omega}{\,\ll\,}{\omega}_0$, was applied to the optical constants of $\alpha $-SiC compiled for the wavelengths 0.04-17${\mu }$m by Choyke & Palik (1985). In this way, we were able to compute the efficiencies $Q_{\rm abs}$ and $Q_{\rm sca}$ in the interval 0.3-300${\mu }$m from the optical constants.

At longer wavelengths, we extrapolated $Q_{\rm abs}$ by a power law ${\lambda}^{-p}$ with $p{\,=\,}$2 for SiC and $p{\,=\,}$1 for the other materials. Although $Q_{\rm abs}{\,\propto\,}{\lambda}^{-2}$ at ${\lambda}{\,\gg\,}a$ in the ideal case of compact homogeneous spherical grains, the latter simple model is unlikely to describe the real cosmic dust in full detail. In fact, we know that the real dust grains have various non-spherical shapes, rough surfaces, inhomogeneous chemical composition, and aggregate (porous) structure. A lower exponent approaching unity is more appropriate for coagulated particles with fluffy structure (Wright 1987; Ossenkopf 1993; Ossenkopf & Henning 1994; Henning & Stognienko 1996) or layered materials (Tielens & Allamandola 1987). Moreover, a lot of observational material has been accumulated over the last decades indicating that $p{\,\approx\,}$1 describes the optical properties of circumstellar dust much better than $p{\,=\,}$2 at long wavelengths (cf. Table 1 for IRC$~\!$+10216).

In addition to the homogeneous spherical grains, we considered more realistic core-mantle particles and studied shape effects using the model of continuous distribution of ellipsoids (CDE, Bohren & Huffman 1983). To account for the composite, non-uniform internal structure of the real grains, we utilized the effective medium calculator, written by Volker Ossenkopf, which simulates the effects of pollution, porosity, and aggregation on the optical constants of the input materials (Ossenkopf 1991). We explored inclusions of one material into the other (Maxwell-Garnett effective medium) and of unorganized conglomerates (Bruggeman effective medium). Simple spherical impurities and their shapes in the extreme cases of disks and needles, have been investigated. Exploring the parameter space of the mixtures, we tried to keep a balance between the obvious complexity of the reality and simplifications necessary to avoid too many free parameters.

   
3.3.4 Emission feature at 11.3${\mu }$m

The broad emission band seen in IRC$~\!$+10216 and many other carbon stars between 10${\mu }$m and 13${\mu }$m (Fig. 3) is usually called the 11.3${\mu }$m feature.

  \begin{figure}
\par\resizebox{\hsize}{!}{\includegraphics{h232403.eps}}\end{figure} Figure 3: Silicon carbide feature and the 11.3${\mu }$m emission band in IRC$~\!$+10216 (Yamamura et al. 1998). The upper set of curves shows the absorption efficiencies of small (a=0.01${\mu }$m) spherical $\alpha $-SiC grains computed using the optical constants from Choyke & Palik (1985) and Pegourie (1988). The shape effects on $Q_{\rm abs}$ are also illustrated in the CDE approximation. The lower set of curves (shifted down for clarity) shows efficiencies for the Bruggeman effective medium grains of the same size, composed of $\alpha $-SiC (Choyke & Palik 1985) and cel800 (Jaeger et al. 1998) and mixed using the volume fractions of 9:1, 5:5, and 15:85. The CDE model efficiencies are also plotted for comparison. One should be cautious in making direct comparions of the dust optical properties with the observed SEDs: self-consistent radiative transfer modeling is necessary for correct analysis and reliable conclusions
Open with DEXTER

Several detailed profiles of the feature have been measured in the past 25 years (Treffers & Cohen 1974; Merrill & Stein 1976; Olnon et al. 1986; Yamamura et al. 1998). Although it is well established that the band is due to a SiC dust component, there are still considerable debates on what causes the shapes of the broad emission in carbon stars.

Blanco et al. (1994) attempted to explain them by a mixture of two crystalline forms of silicon carbide, $\alpha $- and $\beta$-SiC. Kozasa et al. (1996) proposed core-mantle grains consisting of $\alpha $-SiC cores and carbon mantles (or $\alpha $-SiC inclusions in carbon particles) as the most plausible carrier of the feature. Papoular et al. (1998 )argued, however, that the shape effects alone can explain the emission bands, favoring CDE as an adequate approximation. Speck et al. (1999) concluded that it is actually $\beta$-SiC, not $\alpha $-SiC, that reproduces the observed shapes in carbon stars. None of the above explanations is convincing, since they are based on a highly simplified approach and a narrow wavelength range of just one dust feature. We addressed the problem anew, in the frame of our detailed radiative transfer calculations, aiming to construct a dust model that is consistent with both broad features observed in IRC$~\!$+10216 as well as with all other observations of this object.

Note that the 11.3${\mu }$m band in Fig. 3 is significantly asymmetric and that continuum levels at the blue and red sides of the feature (${\sim }$9.5${\mu }$m and 14${\mu }$m, respectively) differ by 30-50%. This characteristic shape, which is also observed in other carbon stars (see, e.g., V Cyg, RU Vir, IRC+40540 in Yamamura et al. 1998), should be reproduced by successful models.

Many mixtures of SiC with different samples of carbon (volume fractions 0.01-0.99) have been tested in our study. Some of the results are shown in Fig. 3, where we compare shapes of the absorption efficiencies for small grains with the Infrared Space Observatory (ISO) spectrum of IRC$~\!$+10216. It is quite clear that the optical constants of SiC derived by Pegourie (1988) cannot explain the observed 11.3${\mu }$m feature. Both the spherical grains and the CDE model predict much broader shapes of the band and higher level of its red wing, opposite to what is observed. Spherical $\alpha $-SiC grains with the optical constants published by Choyke & Palik (1985) can also be ruled out as they produce a very strong single resonance at 10.75${\mu }$m. The CDE model applied to their data predicts a flatter, much broader profile which looks somewhat more similar to the observed band.

Yet, shapes of $\alpha $-SiC grains in the CDE approximation cannot explain the observations, as SiC is clearly not the only grain material present in carbon stars. In fact, amounts of SiC must be small compared to other carbonaceous dust, because of the great strength of its resonance feature. If we added, however, any small amount of the CDE grains to the most abundant carbon materials (featureless in this wavelength region), there would be no asymmetry in the continuum levels at both sides of the 11.3${\mu }$m feature. This would be an obvious inconsistency with the observations, implying that the CDE model is incapable of explaining the 11.3${\mu }$m feature, contrary to the conclusions of Papoular et al. (1998). By the same reasons, any additional component of spherical SiC grains would also fail to reproduce the essential characteristics of the observed profile.

Results of our modeling plotted in Fig. 3 demonstrate that inhomogeneous grains are the carrier of the band. We confirm the general idea of Kozasa et al. (1996) that a composite material containing both carbon and SiC is indeed responsible for the observed profiles in carbon stars. The results disagree, however, with their suggestion that carbon mantles or Maxwell-Garnett-type $\alpha $-SiC impurities in amorphous carbon grains provide a good model for the 11.3${\mu }$m band. We have found that unorganized Bruggeman-type aggregates of $\alpha $-SiC (constants from Choyke & Palik 1985) and cel800 with volume fractions 15:85 resemble the observed band (Fig. 3). Shape effects in the CDE approximation for such grains do not work so well in explaining the 11.3${\mu }$m profile.

We conclude that the [SiC,C] conglomerates with 80-90% of volume in carbon, not pure SiC grains, provide the best match for the 11.3${\mu }$m band in IRC$~\!$+10216 and other carbon stars. Moreover, as we show in the next sections, the composite [SiC,C] grains even become consistent with both features observed in IRC$~\!$+10216, provided that they are covered by mantles of [Mg,Fe] sulfides.

   
3.3.5 Emission feature at 27${\mu }$m (30${\mu }$m)
Firstly, note that the so-called 30${\mu }$m emission feature in IRC$~\!$+10216 does not actually peak at 30${\mu }$m, but rather at 27${\mu }$m (when the fluxes are plotted as $F_{\nu}$, Fig. 4) or even at 26${\mu }$m (in the $\nu
F_{\nu}$-type plots used in Sect. 4.2). Both the old observations presented by Goebel & Moseley (1985) and the most recent ISO spectrophotometry by Yamamura et al. (1998) display the same peak position.

The originally proposed carrier of the feature, MgS (Nuth et al. 1985), shows a very broad band at 36${\mu }$m, almost 10${\mu }$m away from the observed peak, which is obviously unacceptable. A sample of [Mg0.9Fe0.1]S, the closest laboratory analog of MgS measured by Begemann et al. (1994), displays a sharper peak at 28.2${\mu }$m. When used in radiative transfer calculations, however, this material produces an emission feature centered at 30${\mu }$m, almost 4${\mu }$m away from the peak position in IRC$~\!$+10216. We have to conclude that the identification of the band carrier with MgS should be reconsidered.

  \begin{figure}
\par\resizebox{\hsize}{!}{\includegraphics{h232404.eps}}\end{figure} Figure 4: Magnesium-iron sulfide dust features and the 27${\mu }$m emission band in IRC$~\!$+10216 (Yamamura et al. 1998). The upper set of curves shows the absorption efficiencies of small (a=0.01${\mu }$m) spherical core-mantle grains composed of different carbon materials cel400, cel600, cel800, cel1000 (Jaeger et al. 1998), graphite (Draine 1987), and AC1 (Rouleau & Martin 1991) and covered by [Mg0.9Fe0.1]S (Begemann et al. 1994) mantles with a thickness of 15% above the core's surface. The lower set of curves (shifted down for clarity) shows the efficiencies for the same cores covered by [Mg0.5Fe0.5]S mantles. The double-peak features almost disappear for the carbon cores containing graphite or a large fraction of ordered graphitic structures. One should be cautious in making direct comparions of the dust optical properties with the observed SEDs: self-consistent radiative transfer modeling is necessary for correct analysis and reliable conclusions
Open with DEXTER

In order to resolve the problem, we tried many combinations of [MgxFe1-x]S with different samples of carbon (volume fractions 0.01-0.99). We explored effects of the core-mantle grains, the Maxwell-Garnett and Bruggeman effective media, and the shape effects in the CDE approximation. The latter failed to agree with the observations of IRC$~\!$+10216, in most cases producing much too broad the emission feature (25-45${\mu }$m) centered at ${\lambda}{\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\dis...
...\offinterlineskip\halign{\hfil$\scriptscriptstyle ...30${\mu }$m. Moreover, the optical constants of all the [MgxFe1-x]S sulfides, alone or in a combination with different types of amorphous carbon (Sect. 3.3.3), failed to produce any dust materials having an emission peak between 26 and 27${\mu }$m.

In some cases, when the amorphous carbon cores were covered by the [Mg0.9Fe0.1]S mantles, the outcome was a minimum at 30${\mu }$m, between two peaks at 24${\mu }$m and 37${\mu }$m. Such behavior of the core-mantle grains has already been noticed (e.g., Szczerba et al. 1997). Selected results of our effective medium calculations displayed in Fig. 4 show the absorption efficiencies of relatively thin (15% in radius) [MgxFe1-x]S mantles (x=0.9,0.5) on small carbon grains made of AC1, cel400, cel600, cel800, cel1000, and graphite. The sequence of four laboratory carbon samples, from cel400 to cel1000, is characterized by an increasing amount of ordered graphitic structures (Jaeger et al. 1998). To reproduce the observed band characteristics, the [Mg0.5Fe0.5]S sulfides should be mixed with (or added as a mantle to) either graphite or incompletely amorphous carbon material. Due to changes of the optical properties of the inhomogeneous materials, the 35${\mu }$m maximum of the double-peak feature becomes heavily suppressed (Fig. 4) with only one peak remaining, centered at 27${\mu }$m.

We conclude, therefore, that such carbon grains covered by [Mg0.5Fe0.5]S mantles explain the 27${\mu }$m band observed in IRC$~\!$+10216, provided that the carbon material contains a large fraction of graphitic structures. Note also (see Fig. 4) that the same set of solid materials can be used to explain a family of broad emission bands observed in carbon stars, with peaks between 25 and 37${\mu }$m.

   
3.3.6 The dust model

An interesting outcome of our modeling is that a single population of core-mantle grains can reproduce both features in IRC$~\!$+10216, provided that it is the [SiC,C] aggregates (Sect. 3.3.4) that are covered by the [Mg0.5Fe0.5]S mantles (Sect. 3.3.5); hereafter we refer to such composite grains as to the [SiC,C]+[Mg,Fe]S dust. Extinction efficiencies of such grains for different thicknesses of mantles are compared with the ISO data for IRC$~\!$+10216 in Fig. 5. The shapes produced by moderately thin mantles ( ${\Delta}a/a{\,\sim\,}$0.1) of [Mg0.5Fe0.5]S around the [SiC,C] cores (SiC:C ${\,\sim\,}$1:9 by volume) best resemble the observed flux distribution. Close similarities between their total extinction and the SED, including all essential characteristics of the two bands (peak positions and widths, as well as the continuum levels), suggests that our identification of the materials is indeed a good approximation to the reality. This is further demonstrated by the results of our self-consistent radiative transfer calculations in Sect. 4.

  \begin{figure}
\par\resizebox{\hsize}{!}{\includegraphics{h232405.eps}}\end{figure} Figure 5: Inhomogeneous dust grains made of silicon carbide, carbon, and magnesium-iron sulfides explain both the 11.3${\mu }$m and 27${\mu }$m emission bands in the spectrum of IRC$~\!$+10216. The lower set of curves shows extinction efficiencies of small (a=0.01${\mu }$m) spherical core-mantle grains composed of [SiC,C] (volume ratio SiC:C=15:85) cores and [Mg0.5Fe0.5]S mantles with relative thicknesses of 7, 15, and 26%. The upper curves (shifted up for clarity) show that a mixture of carbon, [SiC,C], and core-mantle [SiC,C]+[Mg,Fe]S grains can also match the observed SED. The amount of SiC in the core-mantle grains is diluted (volume ratio SiC:C=4:96) due to admixture of additional carbon. One should be cautious in making direct comparions of the dust optical properties with the observed SEDs: self-consistent radiative transfer modeling is necessary forcorrect analysis and reliable conclusions
Open with DEXTER


 

 
Table 2: Dust properties in our model of IRC$~\!$+10216. Listed are the dust-to-gas mass ratios ${\rho }_i/\rho $, material densities ${\rho }_{{\rm gr}\,i}$, limiting grain sizes $a_{{\rm min\,}i}$ and $a_{{\rm exp\,}i}$, and exponents ${\gamma }_i$ of the power-law size distribution for three dust components. Arrows indicate the same values as in a column to the left

r (AU)
3-30 30-45 45-100 100-6 105

SiC

$\rho_1/\rho$
0.00042 6  10-5 $\leftarrow$ $\leftarrow$
$\rho_{\rm gr\,1}$(gcm-3) 3.3 $\leftarrow$ $\leftarrow$ $\leftarrow$
$a_{\rm min\,1\,}$(${\mu }$m) 0.005 $\leftarrow$ $\leftarrow$ $\leftarrow$
$a_{\rm exp\,1\,}$(${\mu }$m) 1.0 $\leftarrow$ $\leftarrow$ $\leftarrow$
$\gamma_1$ -3 $\leftarrow$ $\leftarrow$ $\leftarrow$

[SiC,C]

$\rho_2/\rho$
- 0.00165 3  10-5 $\leftarrow$
$\rho_{\rm gr\,2\,}$(gcm-3) - 2.2 $\leftarrow$ $\leftarrow$
$a_{\rm min\,2\,}$(${\mu }$m) - 0.01 $\leftarrow$ $\leftarrow$
$a_{\rm exp\,2\,}$(${\mu }$m) - 1.0 $\leftarrow$ $\leftarrow$
$\gamma_2$ - -4 $\leftarrow$ $\leftarrow$

[SiC,C]+[Mg,Fe]S

$\rho_3/\rho$
- - 0.0039 $\leftarrow$
$\rho_{\rm gr\,3\,}$(gcm-3) - - 3.5 $\leftarrow$
$a_{\rm min\,3\,}$(${\mu }$m) - - 0.015 0.05
$a_{\rm exp\,3\,}$(${\mu }$m) - - 1.15 3.0
$\gamma_3$ - - -6 -3


As is shown by the upper set of curves in Fig. 5, a mixture of carbon, [SiC,C], and core-mantle [Mg0.5Fe0.5]S grains can match the observed shape as well as does the single population of composite grains. It is clear that all the three major dust components are present in the envelope of IRC$~\!$+10216. There are, however, no observational constraints on the relative radial distributions of the components' abundances. We had, therefore, no choice but to reduce the free parameter space by assuming that there are only three dust components:

Table 2 lists all physical parameters of our model as a function of the component number i and the radial distance r. Instead of assuming specific (arbitrary) dust-to-gas mass ratios ${\rho }_i/\rho $, we use constraints of detailed stellar evolution calculations to fix the total mass of the envelope at $M{\,\approx\,}$3.0$M_{\odot }$ (see Sect. 5). Then our radiative transfer modeling predicts both the total mass of dust $M_{\rm d}$ and $\rho_3/\rho$. The latter is very similar to other estimates of the dust-to-gas mass ratio in IRC$~\!$+10216 (cf. Table 5); our grain model is further used to define $\rho_1/\rho$ and $\rho_2/\rho$ for the other components. Adopted bulk densities $\rho_{{\rm gr}}$ of C, SiC, and [Mg0.5Fe0.5]S (2, 3.3, and 6gcm-3) were used to compute densities of the other (inhomogeneous) grain materials. Assumed extremes of the power-law size distributions $a_{\rm min}$, $a_{\rm exp}$ are basically the same for all components. Our modeling shows that exponents $\gamma$ in the innermost and outer parts of the envelope can be similar to the standard interstellar value (-3.5), whereas the model requires larger relative amounts of small grains in dust condensation zones.

The overall physical picture implied by the model is quite natural, as expected (cf. Sect. 3.3.2). The SiC grains first condense in the wind very close to stellar surface, at $T_{\rm SiC}{\,\approx\,}$2000K. The most abundant carbon materials begin to nucleate in the wind when it reaches greater distances and lower temperatures ( $T_{\rm C}{\,\approx\,}$1000K), leading eventually to the formation of [SiC,C] aggregates. When they enter the condensation zone of the [Mg,Fe] sulfides ( $T_{\rm [Mg,Fe]S}
{\,\approx\,}$600K), the latter form mantles on the existing inhomogeneous particles.

In Fig. 6, we plotted optical depths toward the star as a function of wavelength and opacities of all dust components for only the smallest grains.

  \begin{figure}
\par\resizebox{\hsize}{!}{\includegraphics{h232406.eps}}\end{figure} Figure 6: Properties of dust in our model of IRC$~\!$+10216. Opacities for smallest grains (a=0.01${\mu }$m) of all three dust components as well as the total dust opacity are shown. Wavelength-dependent optical depth toward the star is also plotted for both absorption and scattering, along with the total extinction through the whole envelope produced by the entire population of dust grains
Open with DEXTER

The envelope of IRC$~\!$+10216 is optically thick at 10.7${\mu }$m due to a very strong resonance asorption of SiC dust at these wavelengths. In contrast to what has been assumed before, the optical depth at this wavelength ( $\tau{\,\sim\,}$10) is as high as that in the near-IR H and K bands. On the other hand, SiC grains are transparent at optical and near-IR wavelengths. This is the reason for their radiative equilibrium temperature distribution to markedly differ from other dust components (Sect. 4.3).

Note that one should be extremely cautious in making direct comparisons of dust absorption efficiencies to observed SEDs for a wide wavelength range, where radiation transfer effects are clearly important. Any sensible quantitative analysis should be done only in the frame of accurate, self-consistent radiative transfer modeling. We would also like to emphasize here again (as did also Kozasa et al. 1996) dangers of the widely used method of drawing a continuum and subtracting it to analyze the observed shapes. The actual continuum levels are essentially unknown and the contribution of different materials may well not be additive. Moreover, a wavelength-dependent effect of finite apertures used in observations modifies shapes of the wide dust features and needs to be properly taken into account (see Sect. 4.2).

Our exploration of the dust grains' parameter space presented in this and two previous sections included many runs of the radiative transfer code to verify our conclusions drawn on the basis of $Q_{\rm abs}$ curves alone. It is essential that experience accumulated from model calculations directs the simplified analysis described above. Without such a guidance, the comparisons of $Q_{\rm abs}$ to SEDs may easily be misleading.

   
3.4 Other input parameters

In this section we discuss the remaining basic parameters of our model of IRC$~\!$+10216: distance D, bolometric luminosity $L_{\star}$, effective temperature $T_{\star }$, and location R2 of the outer boundary.

   
3.4.1 Distance and luminosity

Estimates of the distance to IRC$~\!$+10216 range from a lower limit of ${\sim }$100pc to 290pc. As is evident from Table 3, a number of recent detailed studies point toward the lower end of the range.

 

 
Table 3: Estimates of D for IRC$~\!$+10216

D(pc)
Reference

290
Herbig & Zappala (1970)
150 Kastner (1992)
150 Crosas & Menten (1997)
100-150 Zuckerman et al. (1986)
110-135 Groenewegen et al. (1998)
$\ge$ 100 Becklin et al. (1969)


In the present model, the value of D=130pc is adopted (see also Sect. 4.1).

Since deviations from spherical symmetry in IRC$~\!$+10216 are moderate, luminosity of the central star obtained from the observed SED will depend weakly on the (unknown a priori) optical depth and viewing direction. Our model yields 13000 and 5200$L_{\odot }$ for the phases of maximum and minimum luminosity, respectively, with uncertainties of about $\pm$20%. These luminosities are based upon the assumption of a sinusoidal periodic variation of fluxes from IRC$~\!$+10216 and on our model fits of the observed SEDs at phases $\phi $=0.24 and $\phi $=0.56 (see Sect. 4.2). Specifically, the following formula was used to describe the variations of stellar luminosity with phase:

\begin{displaymath}L_{\star}(\phi){\,=\,}[\,0.3035\,\cos\,(2\pi\phi)+0.6965\,]\,13280\,L_{\odot}.
\end{displaymath} (1)

For the period of stellar pulsations we adopted the value of P=649 days determined by Le Bertre (1992). A sinusoidal fit to a compilation of all photometric measurements (from 1965 to 1998) in the K band from Dyck et al. (1991), Le Bertre (1992), and from our own data, shows that the period and a 2$^{\rm m}$ amplitude of the light variations were almost constant during the last two decades. Estimates of the phases corresponding to all observations used in our study are based on the period given above and the reference maximum phase of the light curve, $\phi $=0, on November 17, 1988 (Julian date JD=2447483).

Significant changes in the object occurred during the first 15 years of observations, when an average K brightness declined by about 1$^{\rm m}$ to the present level. An extrapolation back in time (to the late 1970s) using our sinusoidal fit starts to give phases deviating by ${\Delta}\phi{\,\approx\,}$0.2 from the old papers' estimates.

   
3.4.2 Effective temperature

The intrinsic stellar energy distribution of IRC$~\!$+10216 is reprocessed by dust grains having wavelength-dependent absorption and scattering properties, and as such is uncertain. According to the spectral classification C9,5 by Cohen (1979), the opaque dusty envelope obscures a star with $T_{\star }$ ${\,\approx \,}$2200K. A number of other indirect estimates from model fits to the observed flux distribution are collected in Table 4.

 

 
Table 4: Estimates of $T_{\star }$ for IRC$~\!$+10216. Bracketed values are uncertainties given in the referenced papers

$T_{\star }$(K)
Reference

2460
Sloan & Egan (1995)
2330(15%) Ridgway & Keady (1988)
2330(15%) Griffin (1990)
2326 Ridgway & Keady (1993)
2300 Keady et al. (1988)
2230 Cohen (1979)
2200(7%) Ivezic & Elitzur (1996)
2100 Lorenz-Martins & Lefèvre (1993)
2010 Winters et al. (1994)
2010 Skinner et al. (1999)
2000(5%) Groenewegen (1997)
2000 Crabtree & Martin (1979)
2000 Mitchell & Robinson (1980)
2000 Rowan-Robinson & Harris (1983)
2000 Martin & Rogers (1987)
2000-2500 Lucy (1976)
1800 Phillips et al. (1982)
1800 Sopka et al. (1985)
1800 Sahai et al. (1989)
1500-2700 Danchi et al. (1994)
1600-2100 Witteborn et al. (1980)
1100-1400 Orofino et al. (1990)


A brief comparison of the numbers suggests again (cf. also Sect. 3.3.2) that a factor of ${\sim }$2 gives an idea of "error bars'' for the parameters derived from a typical radiative transfer modeling.

One major problem with the estimates is that they are all based on the assumption that the central star in IRC$~\!$+10216 radiates as a blackbody. The latter has been known, however, as a poor approximation for the cool, extended molecular atmospheres of AGB stars; thus, the accuracy some estimates appear to claim may be misleading. The absorption-line spectra of stellar atmospheres strongly deviate from a blackbody in the optical to mid-IR range, for the extremely low effective temperatures $T_{\star }$ ${\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\displaystyle...
...\offinterlineskip\halign{\hfil$\scriptscriptstyle ...3000K.

Realistic spectra of stellar atmospheres have not been incorporated in previous radiative transfer models of IRC$~\!$+10216. This may have been justified by the apparent smoothness of the observed SED, which sharply contrasts with the absorption spectra of the model atmospheres (see Fig. 7). Our modeling has shown that, indeed, a blackbody central source with $T_{\star }$ ${\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\displaystyle...
...\offinterlineskip\halign{\hfil$\scriptscriptstyle ...2000K could give a very good fit to the observed short-wavelength SED. Note, however, that such low effective temperatures are not in line with stellar evolution calculations predicting 2500-3000K as a typical range of $T_{\star }$ for AGB stars. We would like to stress that the opaque dusty envelope of IRC$~\!$+10216 makes it impossible to derive $T_{\star }$ with any acceptable accuracy from the dust continuum radiative transfer calculations alone. Stellar spectrum is reprocessed by dust grains with poorly known properties.

It has been suggested by studies of stellar pulsations (Bessell et al. 1996; Hofmann et al. 1998) that stellar radius varies in a relatively narrow range, by ${\sim }$20%, between the phases of maximum and minimum luminosity. We assumed in the final model that the effective temperature of IRC$~\!$+10216 varies between 2800K and 2500K at those phases. More specifically, we used the following formula to describe phase variations of effective temperature:

\begin{displaymath}T_{\star}(\phi){\,=\,}[\,0.0536\,\cos\,(2\pi\phi)+0.9464\,]\,2800\,{\rm K}.
\end{displaymath} (2)

In this study of IRC$~\!$+10216, hydrostatic stellar atmosphere models were utilized, kindly computed for us by Rita Loidl using the MARCS code (Gustafsson et al. 1975) in the version of Jorgensen et al. (1992). Several plane-parallel and spherical atmospheres (calculated at wavelengths 0.25-12${\mu }$m with $T_{\star }$ in the range1500-3000K) were used. In the final models, we adopted spherical atmospheres with $T_{\star }$=2660, 2510, and 2520K for $\phi $=0.24, 0.56, and 0.61, respectively. Carbon-to-oxygen ratio C/O=2, gravity $\log\,g{\,=-\,}0.5$, and solar metallicity, [M/H]=0, have been fixed in all models. Although a $\log\,g{\,=-\,}0.8$ would be more appropriate for the central star in IRC$~\!$+10216, our choice is sufficiently accurate for the modeling, given the uncertainties of model atmospheres and applicability of the latter to this particular object. High-resolution spectra of the atmospheres had been smoothed out to ${\lambda}/\Delta{\lambda}{\,\sim\,}$100 and extrapolated using their slopes at $\lambda{\,\approx\,}$12${\mu }$m to longer wavelengths until they intersect the Rayleigh-Jeans tail of a blackbody of the same temperature and luminosity. However, reality may be more complicated, since radio photospheres of evolved stars are known to be typically twice as large as the optical photospheres (see Skinner et al. 1997; Reid & Menten 1997, for $\alpha $ Ori and various AGB stars).

   
3.4.3 Outer boundary

The outer boundary radius of the envelope of IRC$~\!$+10216 is certainly larger than the maximum distance of 600 $^{\prime \prime }$ to which the far-IR emission of the object has been traced by Infrared Astronomical Satellite (IRAS) (Young et al. 1993). The observed size does not give any information about the physical boundary, as it is a lower limit set by IRAS sensitivity. Molecules that have been traced up to radial distances of only 240 $^{\prime \prime }$ (Groenewegen et al. 1998) also provide no useful constraints on the real extent of the envelope. They set a lower limit because of either a detector sensitivity, or too low a number density of molecules for collisional excitation, or a photodissociation by the interstellar radiation field.

 

 
Table 5: Estimates of $\dot{M}$(10-5$M_{\odot }$yr-1), R2(104AU), and ${\rho _{\rm d}/\rho }$ for IRC$~\!$+10216 ($\dot{M}$ and R2 scaled to $D{\,=\,}$130pc). Bracketed values are uncertainties of $\dot{M}$ given in the papers

$\dot{M}$
R2 ${\rho _{\rm d}/\rho }$ Reference

Dust continuum ( $\dot{M}{\,\propto\,}D$)

6.1
0.9   Winters et al. (1994)
3.8 3.1 0.0045 Skinner et al. (1999)
2.9   0.0036 Danchi et al. (1994)
2.7 3.5 0.0016 Ivezic & Elitzur (1996)
2.6   0.0030 Sopka et al. (1985)
2.1 3.2 0.0040 Bagnulo et al. (1995)
2.1(10%) 14 0.0050 Groenewegen (1997)
2.0 2.6 0.0030 Griffin (1990)
1.1 19   Le Bertre (1997)

Molecular lines ( $\dot{M}{\,\propto\,}D^2$)

2.4
    Crosas & Menten (1997)
1.7     Kwan & Linke (1982)
1.7     Truong-Bach et al. (1991)
1.5   0.0049 Kastner (1992)
1.4(20%)   0.0014 Groenewegen et al. (1998)
1.2     Knapp & Morris (1985)
0.9     Sahai (1987)
0.8   0.0058 Keady & Ridgway (1993)
0.8     Kwan & Hill (1977)


As we are constructing a quantitative model for the entire envelope, having a good estimate of R2 is important for accuracy and consistency of the model parameters. For an expected $\rho_{\rm d}{\,\propto\,}r^{-2}$ density distribution in the outermost parts of our model (containing 98.4% of the total envelope's mass), $M{\,\propto\,}R_{2}$.

In contrast to previous radiative transfer models that used more or less arbitrary R2, we derive it by adopting an observationally determined value for the mass-loss rate $\dot{M}$ in IRC$~\!$+10216 (Table 5). Considering the mass conservation equation for a stationary, spherically-symmetric outflow, ${\rm
d}M/{\rm d}r{\,=\,}4{\pi}{r^2}\rho$, and assuming that $\rho{\,\propto\,}r^{-2}$ and $v{\,=\,}$15kms-1, we readily obtain the location of the physical outer boundary:

\begin{displaymath}R_2{\,=\,}\frac{M}{\dot{M}}\,v.
\end{displaymath} (3)

As explained in Sects. 3.3.6 and 5, we utilized stellar evolution theory constraints to define $M{\,=\,}$3.0$M_{\odot }$ in our model. The only remaining quantity to be specified, $\dot{M}$, is determined from observations of both dust continuum and molecular line emission. Averaging different estimates of $\dot{M}$ from Table 5, we obtain 2.8  10-5$M_{\odot }$yr-1 from dust emission and 1.4  10-5$M_{\odot }$yr-1 from molecular lines. We believe that molecular gas observations are likely to give more reliable estimates of gas mass-loss rates than models of dust radiation. Therefore, we adopt a similar value of $\dot{M}{\,=\,}1.6{~}10^{-5}$$M_{\odot }$yr-1 in this study, assuming also that it applies equally well to the (unobserved) early mass-loss history of IRC$~\!$+10216.

The resulting outer boundary of IRC$~\!$+10216 is located at $R_{2}{\,=\,}6{~}10^5$AU, which corresponds to an angular radius of 4615 $^{\prime \prime }$ (=77$^\prime$=1.3 108$R_{\odot}$=3pc$\,\approx$3 105$R_{\star}$). The outermost material at R2 implies a kinematic age of 2 105 years, a reasonable time for an evolved AGB star with an initial mass of $M_{\star}$ ${\,\approx \,}$4$M_{\odot }$ that has experienced several thermal pulses (Bloecker 1999; Bloecker et al. 2000).

It is clear that there remains some uncertainty in our determination of R2, arising from possible changes of the outflow velocity v or the mass-loss time scale $M/\dot{M}$ in the past, during the long evolution of the central star. A similar uncertainty would be present, however, in any other model adopting the commonly used $\rho{\,\propto\,}r^{-2}$density distribution in the outer parts of the envelope.

   
4 IRC$~\!$+10216: The model results

   
4.1 Model parameters

The model parameters, their uncertainties, and the region of the model parameter space which we have explored are summarized in Table 6.

 

 
Table 6: Model parameters of IRC$~\!$+10216. The only parameters that have been varied extensively are density exponents

Parameter
Symbol Units Range Value Error Comment

Distance
D pc -- 130 -- assumed (see Sect. 3.4.1)
Maximum luminosity $L_{\star}$ $L_{\odot }$ 1...2 104 13000 $\pm$20% model (from SED at $\phi=0.24$, Sect. 4.2)
Minimum luminosity $L_{\star}$ $L_{\odot }$ 3...6 103 5200 $\pm$20% model (from SED at $\phi=0.56$, Sect. 4.2)
Eff. temperature (max) $T_{\star }$ K 1200...3000 2800 $\pm$10% model (see Sect. 3.4.2)
Eff. temperature (min) $T_{\star }$ K 1200...3000 2500 $\pm$10% model (see Sect. 3.4.2)
Stellar radius (max) $R_{\star}$ $R_{\odot}$ -- 500 -- derived ( $R_{\star}{\,=\,}{(4\pi\sigma)^{-1/2}}{L_{\star}^{1/2}}\,{T_{\star}^{-2}}$)
Stellar radius (min) $R_{\star}$ $R_{\odot}$ -- 390 -- derived ( $R_{\star}{\,=\,}{(4\pi\sigma)^{-1/2}}{L_{\star}^{1/2}}\,{T_{\star}^{-2}}$)
Inner boundary R1 AU -- 3.3 -- derived ( $R_{1}{\,\approx\,}$1.6$R_{\star}$, Sect. 3.3.2)
Outer boundary R2 AU -- 6 105 -- derived ( $R_2{\,=\,}(M/\dot{M})\,v$, Sect. 3.4.3)
Envelope's total mass M $M_{\odot }$ -- 3 -- assumed (see Sect. 4.3)
Early mass-loss rate $\dot{M}$ $M_{\odot }$yr-1 -- 1.6  10-5 -- assumed (see Sect. 3.4.3)
Outflow velocity v kms-1 -- 15 -- assumed (see Sect. 3.2, Appendix A)
Density exponent ${\alpha}_1$ -- -9...+9 +1 -- model (3.3-30AU, see Sect. 4.3)
Density exponent ${\alpha}_2$ -- -9...-1 -2 -- model (30-45AU, see Sect. 4.3)
Density exponent ${\alpha}_3$ -- -9...-1 -5 -- model (45-70AU, see Sect. 4.3)
Density exponent ${\alpha}_4$ -- -9...-1 -3 $\pm$20% model (70-176AU, see Sect. 4.3)
Density exponent ${\alpha}_5$ -- -2...-1 -1.3 $\pm$20% model (176-3000AU, see Sect. 4.3)
Density exponent ${\alpha}_6$ -- -6...-2 -4 -- model (3000-5700AU, see Sect. 4.3)
Density exponent ${\alpha}_7$ -- -2...-2 -2 -- model (5700- 6 105AU, see Sect. 4.3)
Cavity opening angle ${\omega}$ $^\circ $ 10-80 36 $\pm$10% model (see Sects. 4.2, 4.5.2)
Viewing angle $\theta _{\rm v}$ $^\circ $ 0-90 40 $\pm$10% model (see Sects. 4.2, 4.5.2)


Density distribution is assumed to be a broken power law $\rho{\,\propto\,}r^{{\alpha}_i}$, the exponents ${\alpha}_i$(i=1, 2,..., 7) being allowed to change in different radial zones. "Early mass-loss rate'' corresponds to the outer envelope ( $r{\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\displaystyl...
...\offinterlineskip\halign{\hfil$\scriptscriptstyle ...5000AU) having a $\rho{\,\propto\,}r^{-2}$ density distribution. The fourth column shows a typical range of the parameter variations in our modeling. The fifth column lists the values adopted in the final model. The sixth column contains estimates of the uncertainties associated with the final model parameters. They indicate very approximately the maximum changes of the parameters, such that the fit to the observations can still be restored by adjusting some other model parameters. The uncertainties should not be interpreted as absolute error bars. If some of our general assumptions turn out to be insufficiently realistic, this would possibly affect the results more than we estimated.

The strategy of our exploration of the model parameter space was the following. In the first step, we used a single dust component in the model to get an idea of the approximate density structure and parameters that are necessary to come close to the observed SEDs and intensity maps. Very soon we found the general geometry of the dense core and cavity used in our final model. In the second step, we explored the dust grain model based mainly on the approximate comparisons of $Q_{\rm abs}$ with the SED; many full radiative transfer calculations were also necessary to test different possibilities. In the third step, we refined the envelope's model by plugging the dust model into the radiative transfer code and comparing the model results with all observational constraints. This was the most time-consuming part of the modeling that required many hundreds of runs to converge to the final model.

We would like to emphasize that not all of the parameters in Tables 2 and 6 are independent free parameters and that all of those that are independent and free are constrained by the model fits to available observations. The most extensively varied parameters in our modeling of the envelope are ${\alpha}_i$ (i=1, 2,.., 5), $\rho_{j}/\rho$, $a_{{\rm
min}\,j}$, $a_{{\rm exp}\,j}$, and ${\gamma}_j$ (j=1, 2, 3). Less extensively varied parameters are ${\omega}$, $\theta _{\rm v}$, ${\alpha}_i$ (i=6, 7), and $T_{\star }$. The least extensively varied parameters are $L_{\star}$, dust model, and numerical resolution in terms of the distributions and numbers of intervals over radii, wavelengths, and grain sizes (176, 220, and 19, respectively, in the final model).

To give the reader an idea of how "unique'' the model is, we can say that it was very difficult to find the model that quantitatively fits all observational constraints. It was much easier to find a model that fits only a slightly smaller number of observations. The unprecedented number of constraints used in a self-consistent modeling does make a real difference with other much simpler models. This gives us a strong feeling of confidence in the results described below.

It may be useful to remind here that results of our radiative transfer modeling are fully invariant with respect to the distance D, except for several parameters which are scaled in the following simple way:

\begin{displaymath}L{\,\propto\,}D^2, \,\,M{\,\propto\,}D^2, \,\,\dot{M}{\,\propto\,}D,
\,\,R{\,\propto\,}D, \,\,\rho{\,\propto\,}D^{-1}.
\end{displaymath} (4)

These relations make it easy to scale the model results to another distance, if necessary.

   
4.2 Spectral energy distribution

Figures 7-9 compare model energy distributions with the observed SEDs of IRC$~\!$+10216 at two reference phases $\phi $=0.24 and $\phi $=0.56.

  \begin{figure}
\par\resizebox{\hsize}{!}{\includegraphics{h232407.eps}}\end{figure} Figure 7: Observed SED of IRC$~\!$+10216 compared to our model with $L_{\star }{\,=\,}9500$$L_{\odot }$ ($\phi $=0.24; for references, see Appendix B). Stellar continuum, total and beam-matched fluxes, and SED for an equivalent spherical envelope are plotted. The latter emphasizes effect of the bipolar outflow cavities. The observed broad-band fluxes are shown by different symbols to differentiate between those measured at around the maximum luminosity, intermediate phases, and minimum luminosity. If available, error bars for the fluxes are drawn, whenever they are larger than the symbols. Thick curves (dots) show all available spectrometry and spectrophotometry data; the upper and lower continua between ${\sim }$2${\mu }$m and 8.5${\mu }$m correspond to phases $\phi $=0.24 and $\phi $=0.56, respectively. The model assumes that we observe the envelope at $\theta _{\rm v}$=40$^\circ $ (relative to the midplane). Effect of beam sizes is shown by the vertical lines, clearly visible in the model SED at $\lambda {\,>\,}$50${\mu }$m. Whereas only the lower points of the vertical lines are relevant, we have connected them to the adjacent continuum by straight lines, to better visualize the effect. The insert compares the millimeter- and centimeter-wavelength portion of the model SED with available radio observations. Three arrows at the top axis indicate the central frequencies of the J, H, and K photometric bands
Open with DEXTER


  \begin{figure}
\par\resizebox{\hsize}{!}{\includegraphics{h232408.eps}}\end{figure} Figure 8: Observed SED of IRC$~\!$+10216 compared to our model with $L_{\star }{\,=\,}5500$$L_{\odot }$ ($\phi $=0.56). See also caption to Fig. 7
Open with DEXTER


  \begin{figure}
\par\resizebox{\hsize}{!}{\includegraphics{h232409.eps}}\end{figure} Figure 9: More detailed comparison of our model to the observed SED of IRC$~\!$+10216. Both total and beam-matched continua are shown here for phases $\phi $=0.24 and $\phi $=0.56. Total model fluxes at maximum luminosity and the SED of a stellar atmosphere at minimum phase are also displayed. The three arrows at the top abscissa indicate frequencies of the K, L, and M photometric bands
Open with DEXTER

The latter were chosen because good spectrophotometry data covering top portions of the SEDs have been obtained at those phases. To better visualize the observed shapes, we found it useful to scale to these phases several spectrophotometry segments (measured at different epochs) as described in Appendix B.

In Fig. 7, we compare the model SED at $\phi $=0.24 to the observed spectrophotometry and broad-band flux distribution of IRC$~\!$+10216 in the entire range from optical to radio wavelengths. For reference, the SED of an equivalent spherical envelope is also shown. The latter differs from the model only in that the polar outflow cavities are closed, i.e. $\psi$=180$^\circ $ (Fig. 1). The effect of beam sizes is shown in Fig. 7 at all wavelengths where fluxes were measured; it is visible mainly at the far-IR and millimeter wavelengths. The beam-matched fluxes are by a factor of 3-4 lower than the total model fluxes, emphasizing the very large extent of the envelope that is much larger than any photometric apertures. Note that care must be taken to compare our model with only those (few) measurements, which were obtained close to the luminosity phase considered. At $\lambda{\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\displ...
...\offinterlineskip\halign{\hfil$\scriptscriptstyle ...100${\mu }$m, the dust radiation is compact enough compared to the beam sizes of the observations.

The fit is almost perfect at all wavelengths, except for two fluxes (at 3.6cm and 6cm), where some other sources of radio emission may be necessary to invoke to explain the observed fluxes which are by a factor of ${\sim }$2 in excess of the model predictions. Sahai et al. (1989) found that the contribution from molecular lines to the radio fluxes is negligible. Knapp et al. (1995) estimated that the radio flux at 3.6cm is too low to be from an HII region, yet too high to originate in a stellar photosphere, concluding that an extended chromosphere was most likely source of the excess. Griffin (1990) also found in his radiative transfer model that there must be additional free-free radiation of the stellar chromosphere at these wavelengths.

The model fits the observations at phases of both high and low stellar luminosity. In Fig. 8 we present the same model with $L_{\star }{\,=\,}5500$$L_{\odot }$ ($\phi $=0.56), close to the minimum luminosity. Comparison of the model SED with the flux distribution of an equivalent spherical model emphasizes the influence of the bipolar geometry on subarcsecond scales (Fig. 1). The spherical model produces fluxes by orders of magnitude lower than those observed at $\lambda{\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\displ...
...\offinterlineskip\halign{\hfil$\scriptscriptstyle ...2${\mu }$m (Figs. 7,8). Ivezic & Elitzur (1996) proposed that larger dust grains ( $a{\,\sim\,}$0.2${\mu }$m) alone can explain excess short-wavelength fluxes in a spherical model. Although it is easy to construct an optically thinner spherical model that perfectly fits the observed shape, such a model is unable to explain the entire set of observational constraints considered in the present study. Our model demonstrates that it is the radiation escaping through the optically thinner bipolar cavities, which produces the excess fluxes over the spherical SED at $\lambda{\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\displ...
...\offinterlineskip\halign{\hfil$\scriptscriptstyle ...2${\mu }$m.

It is true, however, that a fairly wide distribution of grain sizes (from $a_{\rm min} \approx 0.01$ ${\mu }$m to $a_{\rm
max}{\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\d...
...\offinterlineskip\halign{\hfil$\scriptscriptstyle ...1${\mu }$m) is essential for successful fits of all the data in our model of IRC$~\!$+10216. Although one can reasonably expect plenty of small grains in the dust formation zone, being poor scatterers they cannot reproduce the short-wavelength excess of the SED. On the other hand, large grains alone scatter too much to fit both the shapes of the broad emission features and our near-IR speckle images (Sect. 4.5). We are left with a physically quite natural choice of a wide distribution of grain sizes, consistently with what has been deduced from observations on the basis of other considerations (Sect. 2.2).

In Fig. 9, we compare to observations only those portions of the model SEDs (at $\phi $=0.24 and $\phi $=0.56) which include two broad dust emission features at 11.3${\mu }$m and 27${\mu }$m, as well as the 3.1${\mu }$m absorption feature (Sect. 3.3). Being a superposition of many molecular absorption lines in stellar photosphere, the 3.1${\mu }$m band is present in the intrinsic stellar spectra at both phases. It is also clearly visible in the model SED (Fig. 9) at phases around the luminosity minimum, due to scattered light. Increasing luminosity enhances hot dust continuum emission, veiling the absorption feature. In fact, it almost disappears in both observed and model SEDs already at $\phi $=0.24.

 

 
Table 7: Wavelengh-dependent amplitudes of photometric variations of IRC$~\!$+10216 between $\phi $=0 and $\phi $=0.5

Wavelength $\lambda$ (${\mu }$m)
1.24 1.63 2.19 3.79 4.64 8.4 9.7 12.9 18.1 Comment

Observed amplitude (mag)
2.20 2.14 2.03 1.65 1.47 0.84 0.94 0.63 0.70 from Le Bertre (1992)
Model amplitude (mag) 1.06 1.26 2.00 1.63 1.42 0.91 0.86 0.76 0.70 $2.5{\,\log{(F_{\lambda\,0.0}/F_{\lambda\,0.5})}}$


Shapes of the broad dust emission features are also reproduced remarkably well in Fig. 9. Note that the emission bands are fitted by the beam-matched model fluxes (a 25 $^{\prime \prime }$ circular beam was used in the model). Moreover, both the continuum level and shape significantly change in the 20-50 ${\mu }$m range when the effect of finite beams is taken into account. Differences between the total and beam-matched fluxes emphasize the fact that the spatial distribution of intensity is strongly wavelength-dependent. Fits to such wide emission or absorption features should be carried out using accurate radiative transfer calculations (see also Sects. 3.3.4, 3.3.5).

  \begin{figure}
\par\resizebox{\hsize}{!}{\includegraphics{h232410.eps}}\end{figure} Figure 10: Temperatures and densities of dust grains in the envelope of IRC$~\!$+10216. The distributions are displayed for the smallest grains only (a=0.01${\mu }$m) of all dust components; temperatures correspond to the luminosity maximum ($\phi $=0.0). Also shown are the total (gas+dust) density profiles in the midplane and in bipolar outflow cavities along the symmetry axis. As discussed in Sect. 3.3.6, different dust components nucleate sequentially in spatially separate zones (at radii of ${\sim }$30 and 45AU) according to their respective condensation temperatures of ${\sim }$1000 and 600K. Total mass of the huge envelope is 3$M_{\odot }$, assuming the dust-to-gas mass ratios given in Table 2. Additional labeling above the lower abscissa shows expansion times assuming a constant outflow velocity of 15kms-1, facilitating approximate kinematic considerations. The two vertical lines drawn on the left side of the diagram (labeled with an asterisk) indicate stellar radii at the phases $\phi $=0.0 and $\phi $=0.5
Open with DEXTER

In Table 7, we compare wavelength-dependent amplitudes of photometric variations of IRC$~\!$+10216 obtained by Le Bertre (1992) to the predictions of our model. The latter is obviously in a very good agreement with the observations, except for only two near-IR bands, J and H. An inspection of Figs. 7 and 8 shows that the deviations are related to a "bump'' on the model SED between 1${\mu }$m and 2.2${\mu }$m, which is also clearly visible in the spectrum of the central star. The bump is the light from the star scattered by dust grains. The deviation is, therefore, a consequence of the adopted model atmosphere; it may suggest that the real atmosphere of IRC$~\!$+10216 differs from the model adopted in this study.

   
4.3 Densities and temperatures

In Figs. 10 and 11, we display density and temperature distributions in the model of IRC$~\!$+10216 for only the smallest grains with a radius a=0.01${\mu }$m. The density distribution is rather well constrained by available spatially-resolved observations, except for both the innermost and outermost regions. In the immediate vicinity of the star ( $\theta{\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\displa...
...\offinterlineskip\halign{\hfil$\scriptscriptstyle ...0.1 $^{\prime \prime }$), we assumed a $\rho_{\rm d}{\,\propto\,}r$ gradient of the dust density. As there are no observational constraints on the radial dependence of the dust-to-gas mass ratio, we have no other choice but to assume that ${\rho _{\rm d}/\rho }$ is constant there. This assumption implies that the gas density $\rho$ is also proportional to r (Fig. 10).

We should emphasize that the gas density distribution close to the star is neither observationally constrained nor predicted by the model. The situation in the innermost parts of the envelope is very complicated due to both stellar pulsations and unknown details of the temporal and radial dependences of the dust condensation process there. It is fairly clear that in reality, ${\rho _{\rm d}/\rho }$ is a function of r and it may change a lot due to condensation of new dust components out of the gas phase. This means that the gas density may also increase toward the stellar surface; both our model and observations do not "resolve'' the innermost parts of the envelope with $r{\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\displaystyl...
...\offinterlineskip\halign{\hfil$\scriptscriptstyle ...5$R_{\star}$. Our model requires, however, that this transition zone must be relatively transparent in the optical and near-IR wavelengths.

In the outer envelope ( $\theta{\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\displa...
...\offinterlineskip\halign{\hfil$\scriptscriptstyle ...40 $^{\prime \prime }$), we adopted the standard $\rho{\,\propto\,}r^{-2}$ profile implying a spherical outflow with a constant velocity and mass-loss rate. Gas density at R2 in our model is very low ( $n_{\rm H}{\,\approx\,}$0.5cm-3), of the order of the average density of interstellar medium. Both the SED and 1.3mm intensity profile (Sects. 4.2, 4.5.3) require that the dust density distribution is noticeably flatter ( $\rho{\,\propto\,}r^{-1.3}$) at distances 200AU ${\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\displaystyle...
...interlineskip\halign{\hfil$\scriptscriptstyle ...AU (1 $.\!\!^{\prime\prime}$5 ${\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\displaystyle...
...\offinterlineskip\halign{\hfil$\scriptscriptstyle ...23 $^{\prime \prime }$). Our model confirms the presence of a large density enhancement at $r{\,\approx\,}$3000AU (cf. Sect. 2.1). Assuming that the gas and dust density distributions are similar, such a density bump implies an episode of higher mass loss (by a factor of 4, see Fig. 11) approximately 1000 years ago.

Our model predicts that there is a pronounced density enhancement at $r{\,\approx\,}$45AU from the star, most likely caused by a temporal variation of the mass-loss rate during the last decades (see Sect. 4.4). The wave extending over the distances ${\sim }$10-150AU is best visible in the total density distribution (Fig. 10, middle) and the midplane mass-loss rate (Fig. 11). It is possible that dust acceleration has also played a role in the formation of the steep density gradient $\rho{\,\propto\,}r^{-5}$ in the innermost, subarcsecond environment of IRC$~\!$+10216. If we assume that the density wave is moving outward at the observed expansion velocity (${\sim }$15kms-1, see Sects. 2.1, 3.2), we would conclude that the mass-loss increase began ${\sim }$30 years before our observations. We believe that it is more than a mere coincidence that, at the same time (in the early 1970s), the mean K-band brightness of the object started to gradually decline by 1 $.\!\!^{\rm m}$2 to its minimum in 1977, slowly recovering by ${\sim }$0 $.\!\!^{\rm m}$2 afterwards (Dyck et al. 1991). Near-IR colors also became redder in the years of the minimum, corroborating this picture.

In our model, there are SiC grains all the way down to the stellar photoshere, where they are assumed to condense at temperatures ${\sim }$2000K. Although we had to separate the inner boundary of the dusty envelope from the stellar surface for numerical reasons, it may well be that in reality the most refractive dust components nucleate even inside the stellar photosphere (Sect. 3.3.2). Depending on the unknown actual properties and amounts of the grains, structure and effective temperature of such a dusty stellar photosphere may be different from what has been assumed in this work and from what the current model atmospheres predict.

The profiles in Figs. 10 and 11 reflect different temperatures of the components, according to the optical properties of the dust materials (Fig. 6). In particular, the SiC dust with its strong absorption feature at 10.75${\mu }$m is significantly cooler than the other components in the innermost parts of the envelope, producing a 30% jump between their temperature profiles at $r{\,=\,}$30AU. Note that the same SiC grains are by a factor of 3 warmer than the other components in the outermost regions of the envelope of IRC$~\!$+10216, due to a reprocessed (cooler) radiation field there.

Two vertical jumps in the density profiles of individual dust components at 30 and 45AU, appearing in Fig. 10 on top of the broad expanding density wave described above, are due to the condensation of C and [Mg0.5Fe0.5]S, correspondingly, from the gas phase. Our model assumes that almost all SiC grains are incorporated into the [SiC,C] aggregates at $r{\,\approx\,}$30AU, which in turn are almost all incorporated into the [SiC,C]+[Mg,Fe]S grains at $r{\,\approx\,}$45AU (see Table 2). The dust density in the polar outflow cavities is by a factor of 30 lower than that in other parts of the dense core. Such an axially-symmetric density distribution causes the well-known bipolar appearance of IRC$~\!$+10216 on subarcsecond scales (see Sect. 4.5).

The model implies that the dust condensation radii do not change with stellar luminosity (Figs. 10,11). Our high-resolution speckle images obtained at different phases over a period of 3 years do not show any structural changes that can be attributed to the period of stellar pulsations, either (Osterbart et al. 2000). In contrast, previous radiative transfer models and some observational works adopted the view that the dust formation radius (the inner boundary) moves back and forth (by ${\sim }$40%) over half of the stellar luminosity period (Sect. 2.2). Our model is very sensitive to the location of the dust formation zones. Changes of the optical depths associated with the oscillations of the condensation radii lead to significant changes of the model SEDs, images, and visibilities (Sects. 4.2, 4.5, 4.6). We could not find any oscillating model that would be consistent with all observational constraints at different photometric phases. On the other hand, our final model with fixed condensation radii at 30 and 45AU explains all the observations very well.

We believe that both our images and radiative transfer models strongly suggest that the dust formation radii do not oscillate during stellar pulsations (Winters et al. 1995).

  \begin{figure}
\par\resizebox{\hsize}{!}{\includegraphics{h232411.eps}}\end{figure} Figure 11: Temperatures of the smallest dust grains at minimum luminosity ($\phi $=0.5) and mass-loss history of IRC$~\!$+10216. The latter is derived from the total density distribution displayed in Fig. 10. The mass-loss rates along the model symmetry axis (outflow cavities) and in a perpedicular direction (midplane) are plotted, assuming that the outflow velocity v=15kms-1 is constant throughout the envelope
Open with DEXTER

The oscillating models implicitly assume that the radiative-equilibrium dust temperatures are equal to the gas temperature. However, thermodynamic properties of gas may well be decoupled from those of dust and even to some extent from the radiation field in the vicinity of a (pulsating) star. It seems that some necessary conditions for condensation of solids out of the gas phase do not depend directly on the luminosity phase; they may be determined by other processes, such as, e.g., shock waves. This may be one of the reasons, why the dust condensation zones do not move toward the star when its luminosity decreases to a minimum.

   
4.4 Mass-loss rate

A time-dependent quantity of general interest can be readily reconstructed from our stationary model: the mass-loss history that is recorded in the density structure of the parsec-sized envelope over ${\sim}10^{5}$ years. The mass-loss rate

\begin{displaymath}\dot{M}(r){\,=\,}{4\pi}{r^2}{\rho(r)\,}v
\end{displaymath} (5)

is plotted in Fig. 11 as a function of the radial distance and outflow time, assuming that velocity v=15kms-1 is constant throughout the envelope.

As it is clearly visible in Fig. 11, the model of IRC$~\!$+10216 reveals two epochs of very high mass loss (${\sim }$15 and 1000 years ago) up to 9  10-5$M_{\odot }$yr-1. The mass loss rate was by an order of magnitude lower just ${\sim }$100 years ago. Over the earlier AGB evolution of the star, where we have no observational constraints, the model assumes a constant $\dot{M}{\,\approx\,}1.6{~}10^{-5}$$M_{\odot }$yr-1 and a $\rho{\,\propto\,}r^{-2}$ density profile (Sect. 3.4.3). The episode of enhanced mass loss 1000 years ago is "recorded'' in a 1.3mm image (Sect. 4.5.3), if we assume that the millimeter radiation is not significantly contaminated by molecular line emission; the peak mass-loss rate would be reduced otherwise.

The recent episode of high mass loss (terminated 15 years ago) is reconstructed by the present modeling on the basis of our high-resolution images (Osterbart et al. 2000). The sharp decrease of $\dot{M}$ below 30AU (during the last decade) is partly the real variation of the mass loss, or the density wave (Fig. 10). In part, especially in the immediate environment of the star ( $r{\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\displaystyl...
...\offinterlineskip\halign{\hfil$\scriptscriptstyle ...5$R_{\star}$, or during the last 3 years), it may be artificial, a consequence of the unknown radial distribution of the dust-to-gas mass ratio ${\rho _{\rm d}/\rho }$, as discussed in Sect. 4.3. The innermost parts of the envelope are "unresolved'' by both models and observations. We found in the modeling, however, that this zone must be relatively transparent in the optical and near-IR wavelengths.

   
4.5 Model images

   
4.5.1 0.8${\mu }$m HST image and J-band speckle image

In Figs. 12 and 13, we compare our model of IRC$~\!$+10216 with the 0.8${\mu }$m HST image and with our J-band speckle image.

  \begin{figure}
\par\mbox{\includegraphics[width=6.6cm,clip]{h232412a.eps} \hspace{4mm}}
\includegraphics[width=6.6cm,clip]{h232412b.eps} \end{figure} Figure 12: Comparison of the 0.8${\mu }$m model image with a 100mas resolution ( left panel) with an archival HST image of IRC$~\!$+10216 at the same wavelength ( right panel). The brightest southern cavity dominates both images, whereas the scattered light from the opposite (northern) cavity is visible only in the observed image. The northern cavity, as well as the direct stellar light are obscured by a dense dusty core. The position of the central star is marked by an asterisk and the location of the northern cavity is indicated by a circle. The distance between them amounts to 210mas
Open with DEXTER


  \begin{figure}
\par\mbox{\includegraphics[width=6.6cm,clip]{h232413a.eps} \hspace{4mm}
\includegraphics[width=6.6cm,clip]{h232413b.eps} }
\end{figure} Figure 13: Comparison of the J-band model image with a 150mas resolution ( left panel) with our speckle image of IRC$~\!$+10216 ( right panel; Fig. 3a in Osterbart et al. 2000). As in the optical HST image in Fig. 12, light from both (bright) southern and (faint) northern cavities is visible in the speckle image. In the model image, the star (asterisk) and the northern cavity (circle) are obscured by the dense dusty core
Open with DEXTER

The model images were convolved at these wavelengths with Gaussian point-spread functions (PSF) of 100 and 150mas (full width at half-maximum, FWHM). Although we rebinned the convolved 0.8${\mu }$m image to match the HST resolution of 46mas per pixel, we did not try to reproduce in more detail the actual complex PSF of the WFPC2 instrument. Due to the dust scattering which dominates at short wavelengths ( $\lambda{\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\displ...
...\offinterlineskip\halign{\hfil$\scriptscriptstyle ...1${\mu }$m), theoretical images depend critically on the details of the actual density distribution, and dust properties inside and near the outflow cavities. As our model geometry and distribution of dust properties are rather approximate, only a qualitative comparison between the model and HST images can be done. The effect of asymmetries is reduced at longer wavelength, enabling us to make more quantitative comparisons in Sects. 4.5.2 and 4.5.3.

The star is completely obscured in Figs. 12 and 13 by the optically thick dusty core (Fig. 1). As it is expected, faint scattered light from the northern outflow cavity is visible much better in the HST image than in our J-band speckle image. The characteristic cometary shape of the light from the southern outflow cavity is reminiscent of other bipolar outflow sources. If we applied the real PSF of WFPC2, its wider wings would spread the bright light of the southern cavity to the dark surroundings, thus making the model image even more similar to the observed one. There are some significant density inhomogeneities on subarcsecond scales, seen in the HST image as distortions of the light distribution symmetry, mainly on the right side of both cavities (see also Fig. 16); they can only be reproduced by a full 3D model.

Note that the brightest peak A is shifted by 0 $.\!\!^{\prime\prime}$4 away from the star at 0.8${\mu }$m compared to the J-band image. Furthermore, the separation of the peaks A and B in both H and K images (Sect. 4.5.2) is smaller than that in J band by another 0 $.\!\!^{\prime\prime}$3. Such wavelength-dependent shifts are the well-known behavior observed in many optically thick bipolar outflow sources (e.g., HL Tau, L1551 IRS5; Men'shchikov et al. 1999; White et al. 2000). One can easily predict that at even longer wavelengths, where optical depths for both absorption and scattering are smaller, the star will become more visible, whereas the cavity and other "clumps'' will become fainter. The latter will completely vanish at some mid-IR wavelengths leaving only one bright peak at the position of the star - a mixture of the diluted direct stellar light with the radiation coming from the hot innermost parts of the dense core.

   
4.5.2 $\mathsfsl{H,\; K}$, and $\mathsfsl{H-K}$ speckle images

In Figs. 14 and 15, our H- and K-band speckle images of IRC$~\!$+10216 are compared to the model convolved with a circular Gaussian PSF of 95mas (FWHM).

  \begin{figure}
\par\mbox{\includegraphics[width=6.6cm,clip]{h232414a.eps} \hspace{4mm}
\includegraphics[width=6.6cm,clip]{h232414b.eps} }
\end{figure} Figure 14: Comparison of the H-band model image with a 95mas resolution ( left panel) with our reconstructed speckle image of IRC$~\!$+10216 ( right panel; Fig. 3b in Osterbart et al. 2000). Direct light from the central star (located at the position of a black asterisk) has appeared at this wavelength, whereas the northern outflow cavity (white circle) remains obscured by the opaque dense core. The fainter components C and D (cf. Fig. 2) are not modeled in this paper; they are most likely due to local density (optical depth) variations
Open with DEXTER


  \begin{figure}
\par\mbox{\includegraphics[width=6.6cm,clip]{h232415a.eps} \hspace{4mm}
\includegraphics[width=6.6cm,clip]{h232415b.eps} }
\end{figure} Figure 15: Comparison of the K-band model image with a 95mas resolution ( left panel) with our speckle image of IRC$~\!$+10216 ( right panel; Fig. 3c in Osterbart et al. 2000). Direct light from the central star (located at the position of a black asterisk) appears relatively brighter due to a lower optical depth at this wavelength, whereas the northern cavity (white circle) still remains obscured. The fainter components C and D (cf. Fig. 2) are not modeled in this paper; they are most likely due to local density (optical depth) variations
Open with DEXTER


  \begin{figure}
\par\mbox{\includegraphics[width=6.6cm,clip]{h232416a.eps} \hspace{4mm}
\includegraphics[width=6.6cm,clip]{h232416b.eps} }
\end{figure} Figure 16: Comparison of the model $H\!-\!K$ color image with a 95mas resolution ( left panel) with our $H\!-\!K$ speckle image of IRC$~\!$+10216 right panel; Fig. 5 in Osterbart et al. 2000). As in the previous model images, the asterisk marks the position of the star and the circle indicates the position of the northern outflow cavity. Except for the perfectly black background outside the sharp edge of the gray central area, the darker the gray scale is, the bluer is the $H\!-\!K$ color; the brighter the gray scale, the redder the $H\!-\!K$ color. The bluest colors at the bottom parts of the images are due to the bright southern outflow cavity visible here as an elliptic darker area just below the image center (see also Figs. 14, 15). The irregular shape of the outer edge of the observed image is caused by local inhomogeneities of the dusty envelope seen also in Fig. 12
Open with DEXTER

At these longer near-IR wavelengths, the star is becoming visible, along with the dust emission of the hottest and densest inner parts of the envelope. The opposite outflow cavity is completely obscured by the opaque circumstellar environment. Note how the direct stellar light becomes more visible as we move from J to H to K filters (Figs. 13, 14, 15). This is exactly the appearance one would expect for an optically thick bipolar outflow source (cf. Sect. 4.5.1).
  \begin{figure}
\par\resizebox{\hsize}{!}{\includegraphics{h232417.eps}}\end{figure} Figure 17: Normalized intensity profiles in the H band (through the peaks A and B) from the model image in Fig. 14 compared to the observed profiles from our speckle image of IRC$~\!$+10216. Also shown are a 95mas PSF and true (unconvolved) model intensity profiles
Open with DEXTER


  \begin{figure}
\par\resizebox{\hsize}{!}{\includegraphics{h232418.eps}}\end{figure} Figure 18: Normalized intensity profiles in the K band (through the peaks A and B) from the model image in Fig. 15 compared to the observed profiles from our speckle image of IRC$~\!$+10216. Also shown are a 95mas PSF and true (unconvolved) model intensity profiles
Open with DEXTER

The radiation of the hot, dense dusty core clearly appears in the model K image as a faint circular halo around the star (Fig. 15). Note that the same halo is also present in our speckle image, visible only on the east and north-east sides, although its intensity level in the model is somewhat higher (Fig. 18). The shape of the faint emission in the observed image appears incompletely circular because of a patchy foreground extinction in the envelope. Obvious asymmetries in the light distribution of the optical 0.8${\mu }$m HST image and our speckle J-band image (Figs. 12, 13), with respect to the symmetry axis of the H, K images (PA ${\,\approx \,}$20$^\circ $), are most likely caused by the inhomogeneities of the dust density distribution.

In addition to the asymmetries present in the HST and J images of IRC$~\!$+10216, two dark patches of extinction are most clearly visible in our high-resolution $H\!-\!K$ speckle image (Fig. 16). Consistently with the observed color image, the model is shown for only those points of the area, where intensities in both H and K images are higher than 2% of the peak level. The (irregular) sharp edge of the observed $H\!-\!K$ image emphasizes the darker areas north and west of the brightest peak; the same zones of enhanced extinction distort the observed short-wavelength images.

The $H\!-\!K$ color images in Fig. 16 nicely delineates and strongly supports the axially-symmetric geometry of the dense core adopted in our model (Fig. 2). The model displays essentially the same color distribution as in the observed image. The redder color is concentrated around the star (component B), whereas the bright southern outflow cavity appears much bluer, as expected in a bipolar nebula. If the star were at the position of the brightest peak A as, e.g., Haniff & Buscher (1998) assumed, the reddest color would only be located on one side, north of the star. That would be a highly unusual situation, very difficult to understand. Our model demonstrates that a bluer color of the brightest peak does not necessarily imply that we observe direct light of the star. In optically thick, non-spherical (bipolar) environments, this would more likely be the light coming from their optically thinner regions or cavities.

  \begin{figure}
\par\resizebox{\hsize}{!}{\includegraphics{h232419.eps}}\end{figure} Figure 19: Model intensity profiles in the K band (as in Fig. 18) at maximum and minimum luminosities illustrate the phase variations in the near-IR images predicted by our model
Open with DEXTER

In Figs. 17 and 18, our model is compared to the observations in a more quantitative way. The model intensity profiles are plotted in two directions, parallel (PA ${\,\approx \,}$20$^\circ $) and perpendicular (PA ${\,\approx \,}$110$^\circ $) to the projected symmetry axis, as well as the observed intensity profiles extracted from our H and K speckle images. Unconvolved intensity cuts through the star and the southern cavity are also plotted to show the "true'' intensity distribution reconstructed by our modeling with a numerical resolution of 4mas. One can easily recognize the (rectangular) profile of the star below the fainter peak and the true shapes of the intensity distributions which are otherwise buried under a much wider PSF.

There are several characteristics of the observed profiles that are important to note. First of all, the brighter peak is very compact, similar to that produced by the direct stellar light. This means that the cavities are rather narrow, comparable in size to stellar diameter. Full opening angle $\omega{\,=\,}$36$^\circ $ in our model produces somewhat wider profiles at PA ${\,\approx \,}$110$^\circ $, so that we feel that the real cavities in IRC$~\!$+10216 should be slightly narrower ( $\omega{\,\approx\,}$30$^\circ $). We adopted the wider cavity in the model because otherwise the brightest peak A would be not as bright as the observed one relative to the fainter peak B. This small discrepancy results most likely from a simplified density distribution in the compact dense core (Figs. 1, 2), not depending on the polar direction in the vicinity of the outflow cavities. Another feature of the observed profiles is that the separation of the peaks A and B is the same in both H and K images. This implies that optical depths at the wavelengths must be comparable.

The intensity profiles are very similar to those derived from our high-resolution speckle images of the Red Rectangle (Men'shchikov et al. 1998). Bipolar geometry of the latter is similar to our present model, except for significantly wider outflow cavities ( $\omega{\,=\,}$70$^\circ $ vs. $\omega{\,=\,}$36$^\circ $). The dense torus of the Red Rectangle is oriented almost edge-on toward the observer, with a viewing angle $\theta _{\rm v}{\,\approx \,}$7$^\circ $ (relative to the midplane). One may be tempted to extrapolate that model to IRC$~\!$+10216 and interpret the bright components A and B as the light from the southern and northern cavities, respectively. This would imply that the star is between A and B and that the viewing direction is also close to edge-on.

  \begin{figure}
\par\resizebox{\hsize}{!}{\includegraphics{h232420.eps}}\end{figure} Figure 20: Normalized intensity profiles from the model images of IRC$~\!$+10216 in J and K bands compared to the observations by Kastner & Weintraub (1994). Cuts through a 2 $^{\prime \prime }$ PSF and true (unconvolved) model intensity profiles are also shown
Open with DEXTER

Although this interpretation might seem plausible, our modeling has shown that the star cannot be located between the two peaks. The central star in IRC$~\!$+10216 would be too luminous for this model to work. A dusty core at a twice smaller distance from the star would be too dense, too hot, and emit too much radiation to reproduce the observed images. The location of the star between the peaks A and B can definitely be ruled out.

In Fig. 19, we illustrate phase variations predicted by our model in the near-IR images. Plotted are the model intensity profiles in K band at the phases of maximum and minimum of luminosity. The effect of higher luminosity is to increase dust temperature and, therefore, enhance dust emission from the inner parts of the envelope. The outflow cavity becomes relatively brighter than the star, whose direct light gets heavily diluted in the hot dust emission. Another effect is that the distance between the star and the cavity becomes by ${\sim }$10mas larger. This apparent periodic "motion'' is purely the effect of increased stellar luminosity in a bipolar dust shell.

  \begin{figure}
\par\resizebox{\hsize}{!}{\includegraphics{h232421.eps}}\end{figure} Figure 21: Normalized intensity profiles from the model images of IRC$~\!$+10216 at 50${\mu }$m and 100${\mu }$m compared to the slit-scan observations by Lester et al. (1986). The 11 $.\!\!^{\prime\prime}$5 and 22 $.\!\!^{\prime\prime}$9 cuts through PSFs and true (unconvolved) model intensity profiles are also shown
Open with DEXTER


  \begin{figure}
\par\resizebox{\hsize}{!}{\includegraphics{h232422.eps}}\end{figure} Figure 22: Normalized intensity profiles from the model images of IRC$~\!$+10216 at 10.7${\mu }$m and 1.3cm compared to the observations by Bloemhof et al. (1988) and Groenewegen et al. (1997). The observed profiles have been mirrored and averaged before they were plotted. The 0 $.\!\!^{\prime\prime}$68 and 11 $^{\prime \prime }$ PSFs, as well as true (unconvolved) model intensity profiles are also shown
Open with DEXTER

   
4.5.3 ${J,\; K,}$ 10${\mu }$m, 50${\mu }$m, 100${\mu }$m, and 1.3cm profiles

In Fig. 20, we show the model intensity profiles in J and K bands, which are compared to the profiles derived from images of IRC$~\!$+10216 by Kastner & Weintraub (1994). The model is convolved with a circular Gaussian PSF of 2 $^{\prime \prime }$ (FWHM). Consistently with what has been observed in the polarimetric imaging, our model shows a very extended scattering envelope on spatial scales of ${\sim }$2-20 $^{\prime \prime }$ from the star.

In Fig. 21, we show the model intensity profiles at 50${\mu }$m and 100${\mu }$m, which are compared to the slit scans obtained by Lester et al. (1986). The model is convolved with circular Gaussian PSF of 11 $.\!\!^{\prime\prime}$5 and 22 $.\!\!^{\prime\prime}$9, respectively, using the slit sizes reported by the observers ( $11{\hbox{$.\!\!^{\prime\prime}$ }}5{\times}46{\hbox{$^{\prime\prime}$ }}$ and $22{\hbox{$.\!\!^{\prime\prime}$ }}9{\times}46{\hbox{$^{\prime\prime}$ }}$). The profiles along the symmetry axis are shown (PA ${\,\approx \,}$20$^\circ $), although the images are essentially circular at long wavelengths.

The model intensity profile at 50${\mu }$m agrees quite well with the observations, whereas the 100${\mu }$m intensity distribution in the model seems to be more extended. Similar discrepancy has been noted by Lester etal. (1986) between their measurements of the almost unresolved source and the observations by Fazio et al. (1980) of an extended envelope at 60${\mu }$m. We do not have any definite explanation of the differences; most likely, they are related to different observational techniques. While the unconvolved intensity profile seems to be dominated by the unresolved peak from the inner dense parts of the envelope (Fig. 21), low-level emission from $\theta{\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\displa...
...\offinterlineskip\halign{\hfil$\scriptscriptstyle ...10 $^{\prime \prime }$ makes significant contribution to the convolved intensity profile.

In Fig. 22, the model intensity profiles at 10.7${\mu }$m and 1.3mm are compared to those obtained by Bloemhof et al. (1988) and Groenewegen et al. (1997). Constrained by the 1.3mm image, our model of IRC$~\!$+10216 implies a flatter density profile at 200AU ${\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\displaystyle...
...\offinterlineskip\halign{\hfil$\scriptscriptstyle ...3000AU than that in the outermost envelope (Fig. 10). The dust density (and probably $\dot{M}$, too) was enhanced by a factor of ${\sim }$3 approximately 1000 years ago. This result is consistent with the density profiles reconstructed by Fazio et al. (1980) and Groenewegen et al. (1997).

The 10.7${\mu }$m model intensity profile agrees very well with the observed profile when convolved with the empirical PSF of 0 $.\!\!^{\prime\prime}$68 given by Bloemhof et al. (1988). At the same time, the unconvolved (true) intensity profile shows a bright core of ${\sim }$300 mas in radius (Fig. 22), in contrast to the unconvolved intensity profile predicted by the model of Ivezic & Elitzur (1996), which is dominated by direct stellar light. The difference is not surprising if one recalls that their spherical model is optically thin at 10.5${\mu }$m, whereas our two-dimensional model demonstrates that IRC$~\!$+10216 is optically thick (Fig. 6).

In Fig. 23, our model is compared to the recent lunar occultation observations at 10.5${\mu }$m by Stecklum et al. (1999) with an estimated resolution of ${\sim }$50mas. The eastern and western branches of the observed (reconstructed) intensity profile were averaged and mirrored in the figure. Comparing the model with the lunar occultation measurements, one should bear in mind that the profile published by Stecklum et al. was a preliminary reconstruction. We used it in our modeling because the authors were confident that its main characteristics are real and because it supplies additional constraints to the model.

Intensity profiles deconvolved from lunar occultation observations are integrated in the direction along the lunar limb. Thus, while the resolution perpendicular to the limb may be as high as ${\sim }$10mas, all structural information along the limb is essentially lost. Therefore, one cannot expect lunar occultation observations to resolve as fine structural details of an extended object as by high-resolution observations with a true circular PSF of 10mas.

The concept of a one-dimensional "resolution'' for the integrated intensity profile is different from the standard concept of resolution and it may be quite misleading. Critically depending on the (unknown) real intensity distribution in the object under consideration, the corresponding true resolution of the lunar occultation data may be much lower. Almost all information about the immediate surroundings of the star is buried in the emission coming from more distant regions. This is the main reason why Stecklum et al. (1999) found no trace of the star with an expected angular diameter of 40mas in their lunar occultation data.

Another problem is that the occultation observations were carried out very close to the strong absorption resonance of SiC dust at 10.75${\mu }$m which is not very suitable wavelength for resolving the star and the structure of the innermost parts of its circumstellar envelope. Our model shows that optical depth of the envelope of IRC$~\!$+10216 containing the SiC grains at this wavelength is as high as in the near-IR J and H bands ( $\tau{\,\sim\,}10$; Fig. 6). Thus, the flat-topped intensity distribution is to be expected for the lunar occultation observations of an extended optically thick envelope around the carbon star at 10.5${\mu }$m.

  \begin{figure}
\par\resizebox{\hsize}{!}{\includegraphics{h232423.eps}}\end{figure} Figure 23: Normalized intensity profiles from the model images of IRC$~\!$+10216 at 10.7${\mu }$m compared to reconstructed radial intensity profile from lunar occultation observations by Stecklum et al. (1999). The eastern and western branches of the data were mirrored and averaged. The model was convolved with a 0 $.\!\!^{\prime\prime}$03 PSF using a long slit of 0 $.\!\!^{\prime\prime}$03${\times }$300 $^{\prime \prime }$; true (unconvolved) model intensity profiles are also shown
Open with DEXTER

In fact, the smooth true model intensity profile shown in Fig. 23 demonstrates that there is no chance to resolve either the star or any subarcsecond-size structure (except for maybe the bright cavity) by this technique.

   
4.6 Near-infrared visibilities

Although we have compared our model to observations of IRC$~\!$+10216 in terms of the SEDs, images, and intensity profiles, in this section we are going to confront the model and observations also in terms of visibilities. The two-dimensional Fourier transform of an intensity map gives a much more detailed and sensitive measure of the intensity distribution over all spatial frequencies than do convolved images and intensity profiles. This would enable a more precise, quantitative assessment of the accuracy of our model.

In Figs. 24-27, we compare the model with available near-IR visibilities obtained by Mariotti et al. (1983), Dyck et al. (1987), Ridgway & Keady (1988), Weigelt et al. (1998a), and Osterbart et al. (2000). The model visibility is plotted for only two orthogonal directions, parallel (PA ${\,\approx \,}$20$^\circ $) and perpendicular (PA ${\,\approx \,}$110$^\circ $) to the symmetry axis. To some extent, comparisons with observations are ambiguous due to the variable nature of IRC$~\!$+10216 that changes on a time scale of one year, due to both the stellar pulsations and to non-periodic dynamic processes in its environment.

To avoid the problem, our approach was to mainly fit our recent speckle visibility data (Weigelt et al. 1998a; Osterbart et al. 2000) which are coeval with the reconstructed images we modeled above (Sects. 4.5.1, 4.5.2). To give the reader a measure of differences between the observations obtained in different epochs, we also plotted older visibilities obtained by Dyck et al. (1987) and Ridgway & Keady (1988). Some visibility data show apparent problems in their low-frequency range, where they should never exceed unity. For completeness, we show all the data sets which are useful at higher spatial frequencies and for an assessment of the overall uncertainties of the observational material.

Overall, our model fits the near-IR visibilities reasonably well, whereas in the L and N bands it deviates noticeably from available measurements. The short-wavelength J-band visibility (Fig. 24) is more difficult to fit, because of the large role of scattering that strongly depends on how realistic the basic model geometry is and how well we know the spatial distribution and properties of dust grains.

  \begin{figure}
\par\resizebox{\hsize}{!}{\includegraphics{h232424.eps}}\end{figure} Figure 24: Model visibilities of IRC$~\!$+10216 in J and H bands are plotted for only two orthogonal directions, PA ${\,\approx \,}$20$^\circ $ and PA ${\,\approx \,}$110$^\circ $. Visibilities from Ridgway & Keady (1988) are also shown, for reference, in the lower panel
Open with DEXTER


  \begin{figure}
\par\includegraphics[width=8cm]{h232425.eps} %
\end{figure} Figure 25: Model visibilities of IRC$~\!$+10216 in $K^{\prime }$ and K bands are plotted for only two orthogonal directions, PA ${\,\approx \,}$20$^\circ $ and PA ${\,\approx \,}$110$^\circ $. In general, data at low spatial frequencies may be less reliable ( $K^{\prime }$ band, upper panel)
Open with DEXTER


  \begin{figure}
\par\resizebox{\hsize}{!}{\includegraphics{h232426.eps}}\end{figure} Figure 26: Model visibilities of IRC$~\!$+10216 in K and L bands are plotted for only two orthogonal directions, PA ${\,\approx \,}$20$^\circ $ and PA ${\,\approx \,}$110$^\circ $. Data at low spatial frequencies may be less reliable (L band, lower panel)
Open with DEXTER


  \begin{figure}
\par\resizebox{\hsize}{!}{\includegraphics{h232427.eps}}\end{figure} Figure 27: Model visibilities of IRC$~\!$+10216 in M and N bands are plotted for only two orthogonal directions, PA ${\,\approx \,}$20$^\circ $ and PA ${\,\approx \,}$110$^\circ $. Data at lower and upper extremes of spatial frequencies may be less reliable (M band, upper panel). The ${\sim }$10${\mu }$m visibilities which we denoted here as N-band data, were actually measured in the 9-11${\mu }$m atmospheric window ( $\delta \lambda {\,\sim }$ several microns) which includes a strong SiC absorption resonance; the envelope is opaque in this band
Open with DEXTER

The visibility shape is very similar to the observed one, although the absolute values at high frequencies are smaller, indicating that the relative contribution of the largest spatial scales is somewhat too high in our model.

In general, the model visibilities in H and K bands are in much better quantitative agreement with our speckle data. There are some deviations from the measured visibilities, which are directly related to the departures of the model intensity profiles from the observed brightness distribution (Figs. 17, 18). In particular, we were unable to better reproduce the small width of the brightest outflow cavity. The model intensity profiles are noticeably wider at PA ${\,\approx \,}$110$^\circ $, resulting in lower intensities at largest spatial frequencies (Figs. 24, 25). Some of the deviations of the model visibilities from the observed ones are caused by the presence of peaks C and D in our speckle images (Osterbart et al. 2000), which do not exist in the model.

Also in L and N bands (Figs. 26, 27), our model shows some departures from the visibilities measured by Mariotti et al. (1983) and Ridgway & Keady (1988). Note, however, that the measurements themselves differ from each other, too. We believe that the main reason for the discrepancies is that the visibility data for IRC$~\!$+10216 are too old to rely on. Given the large non-periodic variations of the object over the last years, it is not surprising that our recent speckle images are fundamentally incomparable to much older observations. Thus, our model which has been designed to reconstruct the present structure and properties of IRC$~\!$+10216, may not match the old data very well.

We would like to emphasize that the mixture of data from different epochs is a real problem that highly complicates any detailed modeling and interpretation of IRC$~\!$+10216 on subarcsecond scales. Stationary modeling requires coeval observational data; one cannot get meaningful constraints on the model if the time scale for the object's changes is smaller than the time interval over which the observations were obtained. One solution to the problem would be to explicitly include time into the models. Sufficiently detailed and realistic time-dependent modeling of individual objects is, however, out of question in the foreseeable future at the present state of the theory, observations, and numerical modeling.

   
5 IRC$~\!$+10216: Evolutionary considerations

The star is in a very advanced stage of its AGB evolution, as indicated by its long pulsational period (650 days), high mass-loss rate, and carbon-rich chemistry of its dusty envelope. Based on the isotopic ratios of C, N, and O in the envelope, the star's initial mass $M^0_{\star}$ can be estimated to be close to 4$M_{\odot }$ due to moderate hot bottom burning (Guelin et al. 1995). This is in line with a mass estimate based on stellar luminosity that gives $M^0_{\star}$ ${\,\approx \,}$4-4.5$M_{\odot }$ (for $D{\,=\,}$170pc, Weigelt et al. 1998a).

Recently, Kahane et al. (2000) presented observational results on a chlorine isotopic abundance ratio that is sensitive to the s-process nucleosynthesis. Combining these observations with evolutionary models and s-process calculations, the authors derived $M^0_{\star}$ ${\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\displaystyle...
...\offinterlineskip\halign{\hfil$\scriptscriptstyle ...3$M_{\odot }$. Current stellar evolution models do not predict hot bottom burning to take place for this mass range (although the lithium enrichment of Galactic carbon stars of rather low luminosities seems to indicate that additional mixing processes may be present, Abia & Isern 1997). Given the observed CNO isotopic abundance ratios, this suggests that IRC$~\!$+10216 had a low-mass progenitor ( $M^0_{\star}$ ${\,\mathrel{\mathchoice {\vcenter{\offinterlineskip\halign{\hfil
$\displaystyle...
...\offinterlineskip\halign{\hfil$\scriptscriptstyle ...2$M_{\odot }$). For instance, only then the 12C/13C ratio appears to be in agreement with C/O$\,>\,$1 see Kahane et al. 2000.

In the model of IRC$~\!$+10216 presented in this study, we adopted a total mass of the envelope of 3.0$M_{\odot }$. Assuming a present-day stellar core mass of 0.7-0.8$M_{\odot }$, the model implies an initial mass of $M^0_{\star}$ ${\,\approx \,}$3.7-3.8$M_{\odot }$. We would like to stress here again that the model predicts only the dust mass $M_{\rm d}$ of the envelope, not its total mass M. In contrast to other dust radiative transfer models adopting a specific dust-to-gas mass ratio ${\rho _{\rm d}/\rho }$ to convert dust mass to gas mass, we prefer to fix $M{\,=\,}$3.0$M_{\odot }$ based on stellar evolution models for IRC$~\!$+10216. Then our model predicts the values of ${\rho _{\rm d}/\rho }$ ${\,=\,}M_{\rm d}/M$ in both the outer and inner parts of the envelope (Table 2).

It is important to note that the model results presented in this paper are independent of the adopted values of $M^0_{\star}$ and M for IRC$~\!$+10216, provided that simple scaling relations are applied to the following parameters:

 \begin{displaymath}\dot{M}{\,\propto\,}M, \,\,\,\rho{\,\propto\,}M,
\,\,\,\rho_{\rm d}/\rho{\,\propto\,}M^{-1}.
\end{displaymath} (6)

If we assumed a low-mass progenitor and, for example, $M{\,=\,}$1.5$M_{\odot }$, the present core mass can be expected to be close to 0.6$M_{\odot }$ and accordingly $M^0_{\star}$${\,=\,}$2.1$M_{\odot }$. This would increase dust-to-gas mass ratios given in Table 2 by a factor of 2, whereas densities (Fig. 10) and mass-loss rates (Table 6) would be decreased by the same factor. No other changes would be required.

Independently of the actual initial masss of this object, the latter appears to be at a very advanced evolutionary phase. Given the very high mass-loss rate derived in our modeling, we can expect only a few thermal pulses to occur until the AGB evolution of IRC$~\!$+10216 is terminated.

   
6 Conclusions

Recent high-resolution imaging of IRC$~\!$+10216 with a resolution better than 100mas revealed its remarkably non-spherical, "clumpy'' appearance on subarcsecond scales (Weigelt et al. 1998a). In our previous paper, the observations were interpreted in terms of a bipolar, dynamically evolving dusty envelope around the carbon star being at the very end of its AGB evolution (Osterbart et al. 2000).

In this study, we constructed the first two-dimensional radiative transfer model for IRC$~\!$+10216, which describes in detail many aspects of its present-day structure. Reconstructing the density distribution throughout the huge, parsec-sized envelope, the model also gives insights into the past AGB evolution of the star. Quantitatively explaining all relevant manifestations of circumstellar dust grains in a self-consistent, extensive exploration of a very large parameter space (Sect. 4.1), the model gives a reliable, coherent picture of the object.

STAR'S LOCATION. An important and firm result of the modeling is that the brightest component A, visible in the HST and our speckle interferometry images of IRC$~\!$+10216 (Sects. 4.5.1, 4.5.2), does not contain direct light from the central star. In contrast to an intuitive picture adopted in all previous works, the star is actually located at the projected position of the fainter component B. The brightest peak is only the light emitted and scattered by a narrow bipolar outflow cavity which is tilted toward us by 40$^\circ $. Fainter components C and D appear most likely due to density deviations around the outflow cavities from axial symmetry on smaller scales (Sect. 3.2).

STELLAR PARAMETERS. For the adopted distance of 130pc, the stellar luminosity varies from 13000 to 5200$L_{\odot }$ between the maximum and minimum brightness of IRC$~\!$+10216 (Sect. 3.4.1). We assumed that the effective temperature of the central star is changing from 2800 to 2500K between the two phases (Sect. 3.4.2). Being somewhat higher than the temperature usually adopted in radiative transfer calculations, it agrees better with stellar evolution models. The stellar radius varies in a relatively narrow range, from 500 to 390$R_{\odot}$. During the last 1000 years, this star has experienced at least two episodes of very high mass loss with $\dot{M}{\,\approx\,}9{~}10^{-5}$$M_{\odot }$yr-1(Sect. 4.4).

DENSITY DISTRIBUTION. The density distribution of dust cannot be described everywhere by a simple power law $\rho{\,\propto\,}r^{-2}$ characteristic for a stationary spherical mass loss. The density profile found in our modeling (Sect. 4.3) results most likely from a complex interplay between stellar pulsations, changing mass-loss rate, and radiative acceleration of dust and gas in the innermost parts of the envelope. All these processes are changing the structure of the dense subarcsecond region of the envelope on a time scale of years. Polar outflow cavities in the axially-symmetric dense shell cause the observed "clumpy'' appearance of IRC$~\!$+10216 in recent highest-resolution images. The huge envelope of 3pc in radius has a total mass of 3$M_{\odot }$, which corresponds to a dust-to-gas mass ratio of 0.0039 (Table 2). The mass is consistent with an initial mass of the central star of ${\sim\,}$4$M_{\odot }$ and the present-day core mass of 0.7$M_{\odot }$ (Sect. 5).

DUST COMPOSITION. In agreement with theoretical expectations, but in contrast to all previous fits of the object's SED, our model predicts that both the radial distribution and composition of dust around IRC$~\!$+10216 are strongly inhomogeneous (Sect. 3.3.6). The inner dust boundary of the envelope, a subject of considerable debate over decades, is located very close to or even inside stellar photosphere, where the newly formed SiC grains have temperatures of $T_{\rm SiC}{\,\approx\,}$2000K. Most of the carbonaceous dust nucleates further away from the star (at 30AU), whereas the most volatile components condense at even larger distances (at 45AU). Abundances and size distribution of dust grains in our model depend on the radial distance from the star (Table 2).

BROAD EMISSION BANDS. First detailed exploration of a large parameter space of various grain models enabled us to conclude that it is inhomogeneous conglomerate particles that are responsible for the well-known, yet poorly explained 11.3${\mu }$m and 27${\mu }$m (30${\mu }$m) emission features. Our modeling strongly suggests that the 11.3${\mu }$m feature is produced by [SiC,C] grains, unorganized aggregates made of incompletely amorphous carbon with significant graphitic content and silicon carbide (Sect. 3.3.4). The 27${\mu }$m emission band, also known as a 30${\mu }$m feature, is perfectly reproduced by the core-mantle [SiC,C]+[Mg,Fe]S grains (Sect. 3.3.5). Although it is absolutely clear that the real situation in IRC$~\!$+10216 is more complex, our dust model (Sect. 3.3.6) is the simplest approximation that consistently and quantitatively explains dust emission observed in this object. The same set of solid materials can be used to successfully model a family of broad emission bands observed in many other carbon stars, with peaks in the 10.3-12.6${\mu }$m and 25-37${\mu }$m wavelength regions.

IRC+10216 AND CRL2688. We would like to draw several parallels between the object studied in this work and the proto-planetary nebula CRL2688, also known as the Egg Nebula (e.g., Sahai et al. 1998). Both objects have distinct dense, optically thick dusty circumstellar cores with biconical cavities presumably produced by high-velocity collimated outflows. The opening angle of the holes in CRL2688 is $\omega{\,\approx\,}$30$^\circ $, almost the same as in our model of IRC$~\!$+10216 (cf. Sect. 4.5.2). Appearances of the two objects noticeably differ, however, primarily because they are inclined at different angles toward the observer (almost edge-on for CRL2688 vs. $\theta _{\rm v}{\,\approx \,}$40$^\circ $ for IRC$~\!$+10216). In both objects, shell-like density enhancements are observed in scattered light as incomplete arcs, which correspond to periodic increases of similarly high mass-loss rates ( ${\sim\,}10^{-4}$$M_{\odot }$yr-1) by a factor of 2-3. One can possibly identify the recent episode of a very high mass loss in IRC$~\!$+10216 (Sect. 4.4) with the formation of one of such shells. CRL2688 is in a significantly more advanced stage of a rapid transition into the planetary nebulae phase than IRC$~\!$+10216. The above brief (and incomplete) comparison suggests that similar physical processes are operating within this evolutionary sequence.

SYMMETRY, CLOUDS, AND CAVITIES. We should stress that the widely used basic assumption of spherical symmetry appears to have a very limited applicability, as we "turn on'' higher and higher angular resolutions. Furthermore, our model demonstrated that one should be very cautious while interpreting asymmetric optically thick dusty environments on the basis of common sense, instead of accurate quantitative modeling. Extending present results to other non-spherical objects, we believe that interpretations of observations in terms of dust clouds around mass-losing stars may be misleading. We expect that in many cases such "clouds'' will turn out to actually be outflow cavities or other types of inhomogeneities in the density distributions of optically thick circumstellar shells.

CONCLUDING REMARKS. IRC$~\!$+10216 is the most difficult object we took for modeling, primarily because of its variability and enormous amounts of observational material that should be quantitatively explained by a realistic model. Having constructed the model that fits a wide variety of observational constraints, we are confident that the general structure and properties of IRC$~\!$+10216 obtained in this study are indeed good approximations to the extremely complex reality. We should make it absolutely clear, however, that no model can provide ultimate answers. No matter how sophisticated, any model would ignore many important aspects of nature. Having explained all the observations of dust continuum self-consistently, our detailed model can serve as just a first step to a more complete and reliable understanding of IRC$~\!$+10216 and evolved AGB stars in general.

Acknowledgements
We would like to thank Rita Loidl for computations of several model atmospheres used in our study and Volker Ossenkopf for discussions on properties of dust grains. The HST image of IRC$~\!$+10216 has been retrieved from the Hubble Data Archive operated at the STScI, Baltimore, USA. Extensive radiative transfer modeling described in this paper has benefited a lot from the excellent computer network of Stockholm Observatory, in particular of its Infrared Group. Financial support of the latter is highly appreciated. We would like to thank Pawel Artymowicz for permission to use his Hydra cluster of workstations. We are grateful to the referee, Claudine Kahane, whose comments stimulated our work on further improvements of clarity and content of several sections of this paper.

   
Appendix A: Deprojected velocities

Geometry of our model (Fig. 1) implies that we look deep into the bright cavity excavated in the dense subarcsecond core of IRC$~\!$+10216. It is the far side of the outflow cavity, which emits and scatters most of the short-wavelength radiation we receive in the brightest peak A. Our model does not incorporate the fainter components C and D, which we believe are optically thinner spots due to column density fluctuations inside the dense core. We have, therefore, no specific information on their relevant physical parameters. In this section, we consider their apparent velocities in the plane of sky ( $v_{\rm
A}{\,\approx\,}$14kms-1, $v_{\rm C}{\,\approx\,}v_{\rm
D}{\,\approx\,}$5kms-1), measured recently by 2000 over 5 epochs. Simple geometrical estimates below demonstrate that deprojected radial velocities of all the components A, C, and D are equal to the general expansion velocity $v_{\rm r}{\,\approx\,}$15kms-1 of the envelope.

From the detailed modeling of IRC$~\!$+10216 reported in this paper, we know the viewing angle ${\theta_{\rm v}}{\,=\,}$40$^\circ $. The radius-vector $r_{\rm A}$ from B to A makes an angle $\theta_{\rm v}-{\omega}/2$ with the plane of sky, where $\omega{\,=\,}{\pi}-{\psi}{\,=\,}$36$^\circ $ is the opening angle of the cavity (radiation of A is emitted on the far side of the cone). The deprojected radial velocity of the brightest component is readily available:

\begin{displaymath}v_{\rm r,A} = \frac{v_{\rm A}}
{\cos(\theta_{\rm v}-{\omega}/2)}
\approx 15.1\,{\rm km\,s}^{-1}.
\end{displaymath} (A.1)

The radius-vectors $r_{\rm C}$ and $r_{\rm D}$ make angles $\zeta_{\rm C}$ and $\zeta_{\rm D}$ with the plane of sky. Due to the symmetry with which C and D appear in the observed images with respect to the cavity axis, it is reasonable to assume that $r_{\rm
C}{\,=\,}r_{\rm D}$ and $\zeta_{\rm C}{\,=\,}\zeta_{\rm D}$; let us denote them as $r_{\Delta}$ and $\zeta_{\Delta}$. The only other assumption is that the radiation which we observe in the fainter peaks C and D is emitted close to the "surface'' of the dense core, i.e. $r_{\Delta}{\,\approx\,}r_{\rm core}$.

From our K-band speckle images (Figs. 1f-j, Osterbart et al. 2000), the projected distances $x_{\rm A}$, $x_{\rm C}$, and $x_{\rm D}$ between B and each of the three components can be easily measured. Since the distances to C and D are very similar, we can average them and denote as $x_{\Delta}$. Averaging also over the images, we find a ratio of $x_{\rm
A}/x_{\Delta}{\,\approx\,}$1.56.

Considering the line of sight passing through component A (Fig. 1), one can write

 \begin{displaymath}\frac{r_{\rm A}}{r_{\Delta}} =\frac{\cos(\theta_{\rm v}+{\omega}/2)}
{\cos(\theta_{\rm v}-{\omega}/2)}\cdot
\end{displaymath} (A.2)

Considering also the line of sight passing through C or D, one can write

 \begin{displaymath}\frac{x_{\rm A}}{x_{\Delta}} = \frac{r_{\rm A}}{r_{\Delta}}
\frac{\cos(\theta_{\rm v}-{\omega}/2)}
{\cos\zeta_{\Delta}}
\end{displaymath} (A.3)

and finally get

 \begin{displaymath}\cos\zeta_{\Delta} = \frac{x_{\Delta}}{x_{\rm A}}
\cos(\theta_{\rm v}+{\omega}/2).
\end{displaymath} (A.4)

Substituting the numbers, we derive $\zeta_{\Delta}{\,\approx\,}$71$^\circ $ and the (average) deprojected radial velocity of C and D

 \begin{displaymath}v_{{\rm r,}\Delta} = \frac{x_{\rm A}}{x_{\Delta}}
\frac{v_{\...
...s(\theta_{\rm v}+{\omega}/2)}
\approx 14.7\,{\rm km\,s}^{-1}.
\end{displaymath} (A.5)

   
Appendix B: Observational constraints

We used all published broad-band fluxes in the modeling, taking into account beam sizes of the corresponding instruments (Becklin et al. 1969; Strecker & Ney 1974; Campbell et al. 1976; Schwartz & Spencer 1977; Shivanandan et al. 1977; Elias et al. 1978; Fazio et al. 1980; Phillips et al. 1982; Spergel et al. 1983; Sopka et al. 1985; Rengarajan et al. 1985; Goebel & Moseley 1985; Le Bertre 1987,1988; Sahai et al. 1989; Woody et al. 1989; Drake et al. 1991; Harvey et al. 1991; Walmsley et al. 1991; Le Bertre 1992; Altenhoff et al. 1994; Knapp et al. 1995; van der Veen et al. 1995; Groenewegen et al. 1997; Osterbart et al. 2000).

Whenever possible, we preferred to rely on existing spectrometry or spectrophotometry data (Becklin et al. 1969; Miller 1970; Merrill & Stein 1976; Witteborn et al. 1980; Forrest et al. 1981; Herter et al. 1982; Goebel & Moseley 1985; Olnon et al. 1986; Le Bertre 1988; Speck et al. 1997; Yamamura et al. 1998), rather than on the broad-band fluxes, because the spectra are coeval and they consistently cover much wider wavelength intervals. Calibration of the IRAS low-resolution spectrum (see Olnon et al. 1986) was corrected as prescribed by Cohen et al. (1992).

Unfortunately, there are only a few near-IR fluxes measured at phases close to $\phi $=0.56. In order to fill the gap and better visualize the SED (Figs. 7-9), we chose to scale down intermediate-phase spectrophotometry in the intervals 1.8-8.4${\mu }$m ( $\phi{\,\approx\,}$0.38, Witteborn et al. 1980) and 16-36${\mu }$m ( $\phi{\,\approx\,}$0.12, Forrest et al. 1981). The scaling factors of 0.74 and 0.53, respectively, fitted these two segments to the only available spectrum of IRC$~\!$+10216 at $\phi $=0.56 published by Speck et al. (1997) and to several broad-band fluxes in the 2-12${\mu }$m wavelength range ($\phi $=0.5) measured by Le Bertre (1992). To better define the SED shape at $\phi $=0.24, we also scaled (by a factor of 0.575) spectra in the 1.5-2.4${\mu }$m range from Beckli et al. (1969) and Merrill & Stein (1976), which were obtained close to the phase of maximum luminosity.

Important constraints for our modeling of IRC$~\!$+10216 are given by high-resolution optical and near-IR images (Weigelt et al. 1998a; Haniff & Buscher 1998; Osterbart et al. 2000), as well as by longer-wavelengths scans and images (Lester et al. 1986; Bloemhof et al. 1988; Groenewegen et al. 1997; Stecklum et al. 1999). Valuable spatial information has been provided by near-IR visibilities (Mariotti et al. 1983; Dyck et al. 1987; Ridgway & Keady 1988; Weigelt et al. 1998a; Osterbart et al. 2000).

References

 


Copyright ESO 2001