A&A 436, 1009-1020 (2005)
DOI: 10.1051/0004-6361:20052800

Evolution of V838 Monocerotis during and after the 2002 eruption[*]

R. Tylenda

Department for Astrophysics, N. Copernicus Astronomical Centre, Rabianska 8, 87-100 Torun, Poland

Received 1 February 2005 / Accepted 20 February 2005

By fitting the available photometric data on V838 Mon with standard supergiant spectra we have derived principal stellar parameters, i.e. effective temperature, radius and luminosity, and followed the evolution of the object since its discovery in early January 2002. Our analysis shows that the 2002 outburst of V838 Mon consisted of two major phases: pre-eruption, which was observed in January 2002 and a major outburst, the eruption, which started in the beginning of February 2002. During pre-eruption the object seemed to be relaxing after an initial event which had presumably taken place in last days of December 2001. The eruption phase, which lasted until mid-April 2002, resulted from a very strong energy burst, which presumably took place in last days of January at the base of the stellar envelope inflated in pre-eruption. The burst produced an energy wave, which was observed as a strong luminosity flash in the beginning of February, followed by a strong mass outflow in the form of two shells, which was observed as an expanding photosphere in later epochs. In mid-April, when the outflow became optically transparent and most of its energy radiated away, the object entered the decline phase during which V838 Mon was evolving along the Hayashi track. This we interpret as evidence that the main energy source during decline was gravitational contraction of the object envelope inflated in eruption. Late in 2002, dust formation started in the expanding shells which gave rise to a strong infrared excess observed in 2003.

Key words: stars: variables: general - stars: fundamental parameters - stars: winds, outflows - stars: individual: V838 Mon

1 Introduction

The eruption of V838 Mon was discovered in the beginning of January 2002 (Brown 2002). It appeared to be rather a pre-eruption brightening as the main eruption started at the beginning of February 2002 and lasted for about two months (e.g. Munari et al. 2002c). The rise to the main maximum in February 2002 caused to a light echo (Henden et al. 2002) whose spectacular images were obtained with HST (Bond et al. 2003). As shown in Tylenda (2004) an analysis of the echo images allows us to study the details of the dust distribution near V838 Mon, as well as to constrain the distance to the object.

Early classifications of the eruption of V838 Mon as nova-like quickly appeared to be inappriopriate (e.g. Munari et al. 2002c; Soker & Tylenda 2003). The light curve of V838 Mon was atypical for novae. Also, the observed mass loss seemed to be not enough violent for a nova. But the strongest argument against nova-like character came from the spectral evolution of V838 Mon. Classical novae usually reach a minimum effective temperature near their optical maximum brightness and then tend to evolve toward higher effective temperature, reaching values well above 105 K. V838 Mon, on the contrary, reached the maximum effective temperature, characteristic of an A-F spectral type, at the optical maximum in the beginning of February. The object then showed a general tendency to evolve to lower effective temperatures. In mid-April 2002 it almost disappeared from optical images but remained very bright in infrared, becoming one of the coolest M-type supergiants yet observed in astrophysics. Detailed descriptions of the spectral and photometric evolution of V838 Mon can be found in a number of papers including Munari et al. (2002c), Kimeswenger et al. (2002), Kolev et al. (2002), Osiwa\la et al. (2003), Wisniewski et al. (2003), Crause et al. (2003, 2005), Kipper et al. (2004).

The global fading of V838 Mon in optical after eruption enabled us to discover a faint hot continuum in short wavelengths (Desidera & Munari 2002; Wagner & Starrfield 2002) later classified as a normal B3 V star (Munari et al. 2002a). A detailed analysis of observational data available for the progenitor and the enviroment of V838 Mon has been carried out in Tylenda et al. (2005b). The main conclusion of this study is that V838 Mon is most likely a young binary system consisting of two intermediate mass stars. One component is the B-type companion discovered by Munari et al. (2002a). The other component that erupted in the 2002 event has a mass between 5-10 $M_{\odot}$ and before eruption was either a main sequence star of similar mass and spectral type as the B-type companion or a slightly less massive pre-main-sequence star. (However there is no direct evidence that both components form a real binary.) The system is embedded in the interstellar medium, which is most likely to be related to molecular clouds seen around the position of V838 Mon, rather than to be expelled by the object itself as argued by van Loon et al. (2004).

The aim of this paper is to study the evolution of the effective temperature and luminosity of V838 Mon during and after its eruption. This sort of study is crucial for testing possible eruption mechanisms and can be an important step towards a final understanding of the nature of V838 Mon. Our analysis is based on photometric results published in the literature. This is the only kind of observational data for V838 Mon where the available results are sufficiently rich to undertake a uniform and meaningful analysis over a time span since the discovery of the eruption untill the present (December 2004) state of the object.

2 Observational data and analysis

The analysis done in the present study is based on a comparison of the observed photometric magnitudes to a set of reference spectra. The least squares method has been used to get the best fit between the model spectrum and the observations. Since V838 Mon is considerably reddened (see Sect. 3) the model spectra have been reddened before fitting them to the observations. The standard extinction curve, $R \equiv A_{\rm V}/E_{B-V} = 3.1$, has been adopted. From the best fit one obtains an effective temperature, $T_{\rm eff}$, and an angular stellar radius, $\theta$. If the distance to the object is known, a linear radius, R, and a stellar luminosity, L, can then be calculated.

It is clear that the larger the wavelength range that is covered by the photometric measurements, the more reliable are the results obtained. Therefore we have primarily based our analysis on ${\it UBVR_{\rm c}I_{\rm c}JHKL}$ magnitudes available from Munari et al. (2002c) and Crause et al. (2003, 2005). For certain time periods no ${\it JHKL}$ measurements are available and we had to rely on ${\it UBVR_{\rm c}I_{\rm c}}$ data from the above references, supplemented by the ${\it BVRI}_{\rm c}$ measurements of Kimeswenger et al. (2002). In later epochs when the object became very cool, the infrared data were particularly important. Therefore we have also used photometric data from Henden & Munari (2002), Watson & Costero (2002), Lynch et al. (2003, 2004) and Tapia & Persi (2003).

Given the luminosity range obtained in the present study ( $10^4 {-} 10^6~L_{\odot}$) as well as the spectral appearance of V838 Mon during eruption (e.g. Kipper et al. 2004) we have used standard intrinsic colours for supergiants (luminosty class I) when constructing the set of reference spectra. The intrinsic photometric colours have been taken from Schmidt-Kaler (1982 - ${\it UBV}$), Johnson (1966 - ${\it VRI}$ for FGK types), Lee (1970 - ${\it VRI}$ for M types) and Koornneef (1983 - ${\it VJHKL}$). A calibration of the effective temperature and bolometric correction against the spectral type has been taken from Schmidt-Kaler (1982).

In this way we obtained a grid of 21 reference spectra ranging from F0 I to M5 I ( $T_{\rm eff}$ ranging between 7700 and 2800 K). Between the grid points a linear interpolation was done either in the colour - log  $T_{\rm eff}$ plane (e.g. B-V or V-K versus log  $T_{\rm eff}$) to get a reference spectrum for any given value of $T_{\rm eff}$ or in the colour - log $\lambda$ plane if the effective wavelength of a photometric band (e.g. $R_{\rm c}$, $I_{\rm c}$) differed from the Johonson standards.

After mid-April 2002 V838 Mon became cooler than the M5 I standard. In order to have estimates of the object parameters in later epochs we proceeded in two ways. First, we extrapolated our grid of standard spectra down to M 10. The extrapolation is in the log colour - log  $T_{\rm eff}$ plane (including $\log {\rm BC}- \log {\rm T}_{\rm eff}$). In the range of M types we explore the Wien part of the spectrum and consequently the relations are roughly linear in the above plane. Second, we fitted a blackbody spectrum to the observations. For the effective wavelengths of particular photomteric bands, the fluxes obtained from the blackbody law have been converted to magnitudes using the standard calibration of Vega (Schmidt-Kaler 1982; Tokunaga 2000). After having reddened them according to the extinction curve, the magnitudes were compared to the observations.

3 Distance and interstellar extinction

Initial estimates of the distance to V838 Mon done in Munari et al. (2002c) and Kimeswenger et al. 2002) gave values of 0.6-0.8 kpc. They were however based on a wrong interpretation of the observed expansion of the light echo. From a more realistic study of the echo structure Bond et al. (2003) concluded that the distance is at least 6 kpc. A detailed study of the light echo evolution done in Tylenda (2004) and Tylenda et al. (2005b) shows that the distance cannot be unambiguously determined from the observed expansion of the outer echo edge. Only a lower limit of $\sim $5 kpc can be established. However, from an analysis of the inner structure of the echo Tylenda (2004) estimated that V838 Mon is at a distance of $8.0\pm2.0$ kpc. More recently, Crause et al. (2005) made an analysis of the echo expansion observed at the SAAO between May 2002 and December 2004 and conclude that the distance is $\sim $9 kpc.

Wisniewski et al. (2003) have analysed radial velocity structures in the interstellar NaI D lines observed in the spectrum of V838 Mon. From a comparison with a velocity contour map in the Galaxy they conclude that the lower limit on the distance is $\sim $2.5 kpc. From observations of the same lines Kipper et al. (2004) determined a reddening distance >3 kpc and a kinematic distance >4 kpc.

The identification of the B-type companion in the spectrum of V838 Mon in Munari et al. (2002a) enables the application of the spectroscopic parallax method. Adopting the photomatric data, spectral classification and EB-V from Munari et al. (2002a) as well as the absolute magnitude $M_{\rm V} = -1.6$ for the spectral class B3 V from Schmidt-Kaler (1982) one gets a distance to the object of 9.4 kpc. (Munari et al. 2002a, give d = 10.5 kpc, presumably adopting a somewhat different calibration of $M_{\rm V}$.) A lower limit to the distance from the above data can be obtained adopting the ZAMS calibration of Schmidt-Kaler (1982) which gives $M_{\rm V} = -1.1$ for the B3 V type. The distance then decreases to 7.5 kpc. From the photometric data of Crause et al. (2005), Tylenda et al. (2005b) derived parameters for the B-type companion slightly different from those of Munari et al. (2002a), namely $V = 16.26 \pm 0.02$, EB-V = 0.71 and a B4 V spectral type. In this case one obtains a distance of 12.4 kpc and 9.4 kpc using the main sequence and ZAMS calibrations of Schmidt-Kaler (1982), respectively.

The conclusion of Tylenda et al. (2005b) that V838 Mon is very likely to be a young binary system allows us to assume that the system follows more or less the Galactic disc rotation and to estimate the distance from the radial velocity of the system. Kolev et al. (2002) from the fairly symmetric emission H$\alpha$ profile in January 2002 have derived a heliocentric wavelength of the line of 6564.155 Å which results in a heliocentric radial velocity of +62  $\mbox{km}\ \mbox{s}^{-1}$. Kipper et al. (2004) obtained $+59\pm6\ \mbox{km}\ \mbox{s}^{-1}$ from emission line profiles in February 2002. The line profiles were rather complex so this estimate is uncertain. T. Tomov (2004, private communication) derived a radial velocity of $\sim $ $+64\ \mbox{km}\ \mbox{s}^{-1}$ from the H$\beta$ absorption line in the spectrum of the B-type component. The above values agree very well with the radial velocity of a more redshifted component in the interstellar NaI D lines, i.e. $+64\pm1\ \mbox{km}\ \mbox{s}^{-1}$ (Zwitter & Munari 2002; Kolev et al. 2002; Kipper et al. 2004). This supports the conclusion of Tylenda et al. (2005b) that V838 Mon is embedded in the interstellar medium. Thus we can adopt that the heliocentric radial velocity of V838 Mon is +64  $\mbox{km}\ \mbox{s}^{-1}$. Using the results of Dehnen & Binney (1998) this can be transformed to $V_{\rm LSR} = +53\ \mbox{km}\ \mbox{s}^{-1}$ which, adopting the Galactic rotation curve of Brand & Blitz (1993), gives a distance of $\sim $6.7 kpc.

Munari et al. (2005) from an analysis of interstellar extinction, galactic kinematics and properties of the B-type companion have concluded that the distance to V838 Mon is $\sim $10 kpc.

Summarizing the above distance determinations we can conclude that observations of V838 Mon itself and of its B-type companion lead to consistent results and that the system is most likely at a distance between 5 and 12 kpc. For the purpose of the present study we assume that the distance to V838 Mon is 8.0 kpc.

From interstellar components of K I and Na I lines in the spectrum of V838 Mon Zwitter & Munari (2002) estimated an interstellar extinction $E_{B-V} = 0.80 \pm 0.05$ and concluded that this value would imply a distance $\ga$3 kpc. Munari et al. (2002c) have, however, argued that this value of extinction was too large as their distance estimate gave 790 pc and they adopted EB-V = 0.50. This distance estimate was however wrong, as discussed above, and V838 Mon is certainly more distant than 3 kpc. Therefore the EB-V estimate in Zwitter & Munari (2002) can be considered as reasonable.

Kimeswenger et al. (2002) have compared effective temperatures estimated from spectroscopy with photometric colours and found that $0.5 \la E_{B-V} \la 0.9$. For their final discussion Kimeswenger et al. adopted $E_{B-V} \simeq 0.7$.

Kipper et al. (2004) have applied different methods for estimating the reddening. From equivalent widths of the Na I interstellar lines they derived EB-V = 0.49. Using equivalent widths of diffuse interstellar bands they obtained EB-V between 0.41 and 2.1. From dust distribution maps these authors estimated an upper limit of $E_{B-V} \simeq 0.89$ in the direction of V838 Mon.

Probably the most reliable estimate of the reddening done in Kipper et al. is that resulting from their Table 3. These authors have classified their spectra of V838 Mon obtained between 4 February and 2 April 2002 using comparison spectra of supergiants. The (B-V) values expected from the classification are then compared with the observed values giving EB-V. From our compilations we have supplemented the missing data in their Table 3 for 26 March and 2 April and then from all the results we have obtained $E_{B-V} = 1.09 \pm 0.17$.

As discussed in Tylenda et al. (2005b) the available photometric data on the B-type companion give EB-V between 0.7 and 1.0.

In a recent paper Munari et al. (2005) have determined the interstellar extinction from interstellar NaI and KI lines, distribution of the extinction along the line of sight and properties of the B-type companion. They conclude that $E_{B-V} \simeq 0.87$.

In the present study we assume EB-V = 0.9. This value is a good compromise between the various results obtained for V838 Mon itself and its B-type companion, especially if estimates as low as 0.5 are not considered as being too low for a distance >5 kpc near the Galactic plane. We also note that the best fits between the photometric results analysed in the present study and the reference supergiant spectra are in most cases obtained when $E_{B-V} \simeq 0.9$ is adopted.

4 Results

Throughout the present paper the epoch of observations is given in days counted since 1 January 2002. We devide the time span discussed into three phases (for the light curve of V838 Mon see e.g. Munari et al. 2002b, 2002c or Crause et al. 2003, 2005). The first one, called pre-eruption, lasted untill 1 February 2002, i.e. it is includes epochs $\le 32$ days. The object then slowly evolved between V = 10 - 11. The second one, called eruption, lasted untill mid-April 2002 covering the epochs between 32 and $\sim $100 days. The object was then brighter than 9 mag, twice reaching $V \simeq 7$. In mid-April V838 Mon started a fast decline in V, the period with epochs $\ga$ 100 days is called decline.

\includegraphics[width=6.5cm]{2800f1b.eps} \end{figure} Figure 1: Examples of the fits of reference spectra to the observed magnitudes of V838 Mon. Left panel: pre-eruption and eruption phases. Right panel: decline phase. Symbols - observed magnitudes. Full curves - supergiant spectra. Dotted curves - blackbody distributions. All reference spectra have been reddened with EB-V = 0.9. The curves are labelled with the epoch of observations (days since 1 Jan. 2002). Note that the curves and the data corresponding to days 40.8, 47.6 and 72.8 ( left panel) have been shifted by -2, -4 and -6 mag, respectively, while those corresponding to days 301, 503, 861 and 1075 (right panel) have been shifted by 2, 5, 8 and 10 mag, respectively. Parameters of the fits can be found in Table 1.
Open with DEXTER

Figure 1 presents examples of the reference spectra (full curves: supergiants - dotted curves: blackbodies) fitted to the photometric data (symbols). The individual spectra have been labelled with the epoch of observations.

Table 1 shows the results from fitting the photometric data. Column (1) shows the epoch of observations (days since 1 January 2002). Column (2) gives the spectral range used in the fitting procedure (e.g. BI stands for ${\it BVRI}$ while BL means ${\it BVRIJHKL}$) and the references to the data (explained in the bottom of Table 1). The spectral type, the corresponding effective temperature and the angular photospheric radius ($\theta$ is in radians) resulting from fitting the standard supergiant spectra to the data are presented in Cols. (3)-(5), respectively. The (linear) photospheric radius and the luminosity, given in Cols. (6) and (7), have been calculated assuming a distance of 8 kpc. The second line, labelled BB in Col. (3), for epochs greater than 110 days gives the results of the fits using blackbodies as reference spectra.

Table 1: Results from fitting the supergiant and blackbody spectra (the latter only for day >110) to the photometric data. See text for explanations of the columns.

The left panel of Fig. 1 displays the spectra during the pre-eruption and eruption phases. The U magnitude has not been taken into account in the fitting procedure. This band is dominated by the Balmer continuum which is very sensitive to nonstandard phenomena such as departures from hydrostatic equilibrium, non-LTE and winds. As can be seen from Fig. 1 (day 12.2), in pre-eruption the U magnitude was below that expected for standard supergiants. However, during the sharp maximum at the begining of February the Balmer continuum was well above that in supergiants (see day 35.8 in Fig. 1). This period was also marked by the appearance of strong Balmer emission lines (see e.g. Osiwa\la 2003). Also the quality of the obtained spectral fit was then less good in the sense that the observed spectrum was somewhat flatter than the standard one. In later epochs, after the sharp maximum, the overall spectral fit became better, including the U magnitude (see days 40.8 and 47.6 in Fig. 1). Thus the object then resembled a standard supergiant quite closely.

As can be seen from Table 1, for a number of epochs in pre-eruption and eruption the fitting has been done using only the ${\it BVRI}$ magnitudes, as no ${\it JHKL}$ measurements were available for these dates. The results from these fits are expected to be less reliable. However, as discussed in Sect. 5, there is no significant systematic difference between them and those based on the ${\it BVRIJHKL}$ data.

The right panel of Fig. 1 shows examples of the fits obtained in the decline phase. As explained in Sect. 2, the fits for this period were obtained using both the extrapolated standard supergiant spectra and the blackbody distributions. In all fits in the decline period a contribution from the B-type companion, the same as in Tylenda et al. (2005b), has been taken into account. This contribution is important only in the ${\it UBV}$ bands. As can be seen from the right panel of Fig. 1, the observed spectra change their slope below the I band, i.e. for shorter wavelength bands the spectrum is steeper than for the longer ones. This effect, due to strong molecular bands significantly lowering the flux at shorter wavelengths, is roughly reproduced in the supergiant spectra. The blackbody spectra do not take this into account. Therefore, when fitting blackbodies we excluded the ${\it BVR}$ magnitudes from the fitting procedure. For the epochs later than 300 days the L magnitude was also omitted from the fitting (both using supergiant as well as blackbody spectra). This band was likely to be affected by the infrared excess discussed in Sect. 6. Inclusion or exclusion of the L magnitude in the fitting procedure does not significantly change the fit results.

As can be seen from Fig. 1, given the crudeness of the reference spectra used in the decline phase (extrapolated supergiants or blackbodies), the obtained fits to the observations are reasonable (although not as good as those in pre-eruption and eruption). The values of the fitting parameters ( $T_{\rm eff}$ and $\theta$ in Table 1) obtained using the supergiant spectra are not significantly different from those using the blackbodies in most cases. Only for two spectra obtained near day 120 does fitting a blackbody give a $T_{\rm eff}$ larger by almost 600 K and $\theta$ smaller by a factor of $\sim $1.6 compared to the results of fitting a supergiant. As discussed in Sect. 5.2, the object was probably then surrounded by an extended shell in a transition from optically thick to optically thin. In the case of the last spectrum (day 1075) the blackbody fitting also gave $T_{\rm eff}$ about 350 K lower and, consequently, a somewhat larger effective radius than the supergiant fitting.

Note that for the reasons discussed above we do not consider fits based only on ${\it BVRI}$ in the decline phase.

Lane et al. (2005) attempted to determine the angular size of V838 Mon from interferomatric observations done in the K band in November-December 2004. Their result of an angular diameter of 1.8 mas is slightly larger than our result obtained on day 1075 ($\sim $1.1 mas). One possible reason for this difference is that the uniform disc used in Lane et al. to deconvolve the data might be not a good model for the photosphere of V838 Mon. Another possible reason is that although the K band, according to our analysis, is dominated by the central object, a certain contribution from an extended infrared emission, discussed in Sect. 6, may also be important. This might explain why an elliptical model better fits the data of Lane et al. than a circular one. It thus would be desirable to make similar interferometric measurements but in shorter wavelengths, say in I or J, where the central object certainly dominates the flux, as well as in longer wavelengths, say in M or N, where the expanding matter is dominant.

5 Evolution of the object

...2b.eps}\vspace*{3mm} %
\includegraphics[width=7.1cm]{2800f2c.eps} %
\end{figure} Figure 2: Evolution of the effective temperature ( upper), luminosity ( middle) and radius ( bottom) of V838 Mon. The luminosity and radius are in solar units. The abscissa gives the logarithm of time in days counted since 1 Jan. 2002. See text for explanation of the symbols and curves (in the middle and bottom panels).
Open with DEXTER

Figure 2 shows the evolution of the principal parameters of V838 Mon with time. The effective temperature, luminosity and effective radius from Table 1 are displayed in the upper, middle and bottom panel, respectively. On the abscissa in each panel the epoch of observations (days since 1 Jan. 2002), in the logarithmic scale is given. Star-like symbols show the results obtained from fitting the supergiant spectra to the ${\it BVRIJHKL}$ (or ${\it BVRIJHK}$ for epochs > 300 days) magnitudes, i.e. those labelled with BL or BK in Col. (2) of Table 1. The same but derived from the ${\it BVRI}$ magnitudes (labelled BI in Col. (2) of Table 1) is represented with full dots. Open symbols display the results obtained from fitting the blackbody spectra, i.e. those labelled BB in Col. (3) of Table 1. Dotted and dashed curves are discussed Sects. 5.2 and 5.3.

As can be seen from Fig. 2, in pre-eruption and eruption (log time $\la$ 2.0), there is no significant systematic difference between the results obtained from fitting only the ${\it BVRI}$ magnitudes and those derived form the ${\it BVRIJHKL}$ measurements. The full dots consistently follow the same evolutionary trend as the asterisks. In the decline phase, given the uncertainties due to the reference spectra, the results obtained from fitting the blackbodies (open circles) are fairly consistent with those based on the supergiant spectra (asterisks).

In various estimates done later in this paper, requiring parameters of the central star, we adopt, following Tylenda et al. (2005b), standard values for a B3 V star, i.e. $M_\ast = 8~M_{\odot}$ and $R_\ast = 5~R_{\odot}$.

5.1 Pre-eruption: Epochs earlier than 32 days

In the pre-eruption phase (log time $\la$1.5 in Fig. 2) V838 Mon seemed to have reached an equilibrium or relaxation configuration. This phase followed an initial event which had presumably taken place a few days before 1 January 2002 (e.g. Munari et al. 2002c). The photospheric radius between days 10-32 remained at an almost constant value of $\sim $ $350~R_{\odot}$, suggesting that the photosphere was in a quasi-hydrostatic equilibrium. The effective temperature reached a maximum at days 12-15 and subsequently slowly declined. The luminosity showed a similar trend and at the end of the pre-eruption phase the objects was half as luminous as at the beginning.

This behaviour suggests that the initial event, which inflated the stellar envelope to $\sim $ $350~R_{\odot}$, had properties of a short-lived energy burst, after which the energy sources dropped significantly, provoking a cooling of the photospheric regions and a decline in the apparent luminosity. In principle, cooling and a luminosity decline should be followed by a contraction of the photosphere. However, for an $8~M_{\odot}$ star and the observed photospheric radius the dynamical time is $\sim $30 days. Thus the inflated stellar envelope did not have enough time to react during the pre-eruption phase.

Within the above hypothesis one can estimate the mass of the inflated envelope assuming that during the pre-eruption phase it radiated away a significant part of its internal energy. The latter can be estimated from Eq. (A.16). Integrating the observed luminosity over the pre-eruption phase and assuming (unperturbed) stellar parameters, $M_\ast = 8~M_{\odot}$ and $R_\ast = 5~R_{\odot}$, one gets a mass of the inflated envelope, $M_{\rm e} \simeq 4 \times 10^{-3}~M_{\odot}$. This value should be considered with caution. On the one hand, the observed luminosity and temperature evolved only slightly during pre-eruption, so only a part of the envelope energy radiated away. On the other hand, it is likely the energy sources that had presumably initiated pre-eruption did not extinguish completely and were contributing to the observed luminosity in the course of the pre-eruption phase.

Spectroscopic observations obtained in the pre-eruption phase showed numerous lines with P-Cygni profiles indicating significant mass loss. The typical outflow velocity was of 200-300  $\mbox{km}\ \mbox{s}^{-1}$ with a terminal one reaching 500  $\mbox{km}\ \mbox{s}^{-1}$ (Munari et al. 2002c; Osiwa\la et al. 2003). This mass loss, observed during the whole pre-eruption phase, is not surprising, as the stellar luminosity was then of a quarter of the Eddington limit for an 8 $M_{\odot}$ star.

5.2 Eruption: Epochs between 33 and $\sim $100 days

The gentle relaxation in the pre-eruption phase was suddenly interupted on 2 February 2002 (day 33). The event was very dramatic. For almost 3 days V838 Mon brightened by a factor of $\sim $25 in luminosity while the effective temperature increased from $\sim $4600 K to $\sim $7200 K. The event was also followed by a fast rise of the effective radius. A major qualitative difference between the evolution of V838 Mon in the whole eruption phase compared to that in pre-eruption is a more or less steady increase of the effective radius in eruption versus a roughly constant radius in pre-eruption. At the end of the eruption phase the effective radius reached values approaching $3000~R_{\odot}$, compared to $\sim $ $350~R_{\odot}$ in pre-eruption.

However, as can be seen from Fig. 2 (bottom panel), the rise of the effective radius was not monotonic. One can easily distinguish three periods during which the radius was steadily increasing, interupted by two short phases of contraction. A simplistic interpretation of this behaviour might be that the effective radius evolution reflects the behaviour of the stellar envelope, namely that each phase of the envelope expansion is followed by a period during which the envelope is contracting. However, the expansion velocities as estimated below are much higher than the local escape velocities. Thus there is no way to stop the expanding matter at a certain time and to force it to contract. A more likely hypothesis is that during the contraction phase the hitherto photospheric layers become, due to expansion, optically thin and the photosphere starts to form in deeper denser regions. In this case the observed evolution of the effective radius during eruption would suggest ejection of three consecutive shells.

To futher investigate the evolution of the effective radius we have plotted three dotted lines in the bottom panel of Fig. 2 showing radii increasing with constant velocity. The parameters of the lines, i.e. velocity and epoch of the zero radius, are 800  $\mbox{km}\ \mbox{s}^{-1}$ and 28.4 days, 400  $\mbox{km}\ \mbox{s}^{-1}$ and 30.5 days, 270  $\mbox{km}\ \mbox{s}^{-1}$ and 28.0 days for lines (a), (b) and (c), respectively. The fact that the three lines start at almost (within 2.5 days) the same time moment strongly suggests that the whole eruption phase resulted from a single outburst event which took place in the last days of January 2002, deep in the stellar envelope already inflated by the pre-outburst event.

The above suggestion that during the eruption phase V838 Mon ejected three consecutive shells poses some problems when spectroscopic data are considered. The velocity of the first shell (line (a) in the bottom panel of Fig. 2) of 800  $\mbox{km}\ \mbox{s}^{-1}$ is much higher than the observed values. In general, the observed velocities during eruption were similar or even lower than during pre-eruption. The terminal velocity derived by Kipper et al. (2004) on 4/5 Feb. 2002 (days 35-36) was between 240 and 390  $\mbox{km}\ \mbox{s}^{-1}$. Even if the systematic velocity of the object, $\sim $ $65\ \mbox{km}\ \mbox{s}^{-1}$ as discussed in Sect. 3, is taken into account, this is much lower than the velocity of our hypothetical first shell. Typical velocities of maximum absorption in the P-Cygni profiles in eruption were between 60 and 300  $\mbox{km}\ \mbox{s}^{-1}$ with a general tendecy of decreasing (or the appearance of lower velocity components) with time (Kolev et al. 2002; Crause et al. 2003; Wisniewski et al. 2003; Kipper et al. 2004). If the systematic velocity of the object is added, these values are consistent with the velocities of the shells (b) and (c) in Fig. 2.

A varying photospheric radius need not reflect a mouvement of the matter in the photospheric layers. As discussed above, the shrinkage of the photosphere, observed on days $\sim $38 and $\sim $60, cannot be considered as due to a collapse of the envelope. Similarily an increase of the photospheric radius does not necessarily require expansion of a dense, optically thick shell. In the case of an optically thick wind, e.g. as discussed in Bath & Shaviv (1976) for classical novae, the photosphere is formed in the wind and its position is mainly determined by the mass loss rate, $\dot{M}$. If $\dot{M}$ is increasing the photosphere moves outward, while a decrease in $\dot{M}$ leads to photospheric contraction. However, even if $\dot{M}$ remains constant, but the temperature of the wind matter varies, the position of the photosphere would also vary because of the varying opacity. We argue that this is what happened in the luminosity peak observed during the first $\sim $10 days of the eruption phase.

Let us assume a steady-state spherically symmetric wind with a mass loss rate, $\dot{M}$, and expanding with a velocity, $v_{\rm w}$. Then the density, $\rho$, of the wind matter is distributed with a radius, r, as

 \begin{displaymath}\rho = \frac{\dot{M}}{4 \pi~r^2~v_{\rm w}}\cdot
\end{displaymath} (1)

The effective photosphere is formed at a radius, R0, satisfying the condition

 \begin{displaymath}\int_{R_0}^{\infty} \kappa~\rho~{\rm d}r = \frac{2}{3},
\end{displaymath} (2)

where $\kappa$ is the opacity. Assuming a constant mass loss rate, velocity and opacity above the photosphere, Eq. (2) can be integrated, using Eq. (1), leading to

 \begin{displaymath}\kappa~\rho_0~R_0 = \frac{2}{3},
\end{displaymath} (3)

where $\rho_0 = \dot{M}/(4\pi~R_0^2~v_{\rm w})$ is the density in the photospheric region. For a given chemical composition, the opacity, $\kappa$, is a function of temperature and density. Therefore, if values for the effective temperature and radius are given, the density, $\rho_0$, can be derived from Eq. (3). Next, using the wind velocity, $v_{\rm w}$, the mass loss rate $\dot{M}$ can be calculated from Eq. (1).

\par\resizebox{8.1cm}{!}{\includegraphics[clip]{2800f3.eps}} \end{figure} Figure 3: Mass loss rate derived assuming an optically thick wind model (see text for more details) - large symbols. Small symbols - effective radius (the same as in Fig. 2 but plotted as 3 log $~(R/R_{\odot}) - 10.5$).
Open with DEXTER

Figure 3 shows the mass loss rate derived from the values of the effective temperature and radius listed in Table 1, using the above steady-state wind model. As the opacity we have used the Rosseland means as given by Alexander & Ferguson (1994) for the chemical composition X =0.7 and Z = 0.02. A constant wind velocity, $v_{\rm w} = 300$  $\mbox{km}\ \mbox{s}^{-1}$, has been assumed. To make the discussion easier the effective radius has also been plotted in the figure with small symbols.

The values of $\dot{M}$ plotted in Fig. 3 are reasonable estimates of the mass loss rate only if the wind is optically thick (so the effective photosphere is formed in the wind) and builds up an $\sim $r-2 density distribution above the photosphere. The latter condition takes place if the wind does not vary significantly on a time scale $\sim $ $R_0/v_{\rm w}$. This time scale is $\sim $10 days in pre-eruption. Therefore the very fast drop in mass loss rate, by more than 2 orders of magnitude on a time scale of $\la$1 day, obtained on day 32 is not real. If the wind in pre-eruption was optically thick and the drop in mass loss rate was real, the density distribution above the photosphere would have significantly changed after $\sim $10 days. However, this was the beginning of the eruption phase and the rise in luminosity should be followed by an increase in mass loss rate rather than a drop. However what is striking in Fig. 3 is a constant mass loss rate obtained during the first luminosity peak in eruption, i.e. between days 32 and 42. During these 10 days all the stellar parameters, i.e. luminosity, effective temperature and radius, varied by significant factors. For instance, the luminosity first increased by a factor of 25 and then dropped by a factor of 2. Yet the mass loss rate was constant within 10-15%. We therefore suggest that this was the time period during which the wind was optically thick and our model of a steady-state wind gives a good interpretation of the evolution of the object. However, the infered mass loss rate charactarises the wind during the pre-eruption phase rather than in eruption. Thus we interepret the results obtained in pre-eruption and eruption as follows.

In the pre-eruption phase, as discussed in Sect. 5.1, the photosphere was formed in layers approximately in hydrostatic equilibrium. The object had a wind as indicated by the observed P-Cygni profiles. However the mass loss rate was much lower than that shown in Fig. 3 during this phase. Instead we argue that $\dot{M}$ was then close to the value derived from the analysis of the first 10 days in eruption, i.e. it was $\sim $ $5 \times 10^{-4}~M_{\odot}$/yr. Note that, as discussed in Sect. 5.1, the luminosity of the object in pre-eruption was of a quarter of the Eddington limit so a significant mass loss could be expected. At the above mass loss rate and a photospheric temperature of $\sim $5000 K, the wind, during the pre-eruption phase, was optically thin and seen only in the lines. This roughly steady-state wind built up an $\sim $r-2 density distribution above the photosphere. On day 32 a strong luminosity wave (or a series of strong shock waves), presumably generated at the base of the inflated envelope, reached the photosphere and quickly heated up the regions above it. This resulted in an abrupt increase of the opacity (by 2 orders of magnitude if the temperature increased from 5000 K to 7000 K - Alexander & Ferguson 1994) and the hitherto optically thin regions built by the wind became optically thick. This was observed as a very fast expansion of the photospheric radius leading to an apparent velocity of $\sim $ $800\ \mbox{km}\ \mbox{s}^{-1}$, much faster than the wind velocity observed in the lines. Subsequent cooling of the wind matter decreased the opacity, which resulted in a contraction of the effective photosphere observed between days 36 and 42.

The luminosity wave was followed by an intense mass loss (the observed luminosity between days 33 and $\sim $110 was above the Eddington limit for an $8~M_{\odot}$ star) which reached the shrinking photosphere near day 42, as indicated by a large increase of $\dot{M}$ in Fig. 3. Note, however, that the obtained values of $\dot{M}$ should again be regarded with caution. The assumption of an r-2 density distribution in the model wind is not satisfied if the mass loss rate is rapidly changing, as seen in Fig. 3 after day 42. In this case the obtained evolution of $\dot{M}$ in the later epochs, when combined with the evolution of the effective radius, seems to suggests that mass loss following the fast luminosity wave took place in two shells. These shells were presumably formed by the initial eruption event which took place at the end of January, as discussed above, and followed the initial strong luminosity wave. Interaction of the luminosity wave with the envelope inflated in pre-eruption also might have been important for the formation of the shells.

The first shell, somewhat faster than the second one, became visible near day 42, when it met the photosphere contracting in the pre-eruption wind. Subsequently, it was observed as an expanding and cooling photosphere radiating at a roughly constant luminosity of $\sim $ $6 \times 10^5~L_{\odot}$. Due to expansion and cooling this shell became partly transparent around day 58 and the second shell started to be visible. This resulted in a shrinkage of the effective photosphere observed between days 58 and 63. At the same time the photosphere became hotter while the luminosity started to increase and quickly reached a value of $\sim $ $1 \times 10^6~L_{\odot}$remaining more or less constant between days 66 and 100. Judging from the evolution of the effective radius and luminosity, the second shell became partly transparent between days 100-120. The discrepancy in the effective radius derived from fitting supergiant spectra and blackbodies on day 120 (see Table 1 and Fig. 2) can be explained by the presence of an extended, partly transparent envelope. At a temperature of 2500-3000 K the matter is opaque in the optical due to molecular absorption so the effective radius derived from the ${\it BVR}$ bands (important when fitting a supergiant spectrum) would correspond to the outer radius of the envelope. In the infrared the envelope is more transparent and it is possible to see deeper and hotter layers. Thus a higher effective temperature and a smaller effective radius is obtained from fitting a blackbody to the ${\it IJHKL}$ measurements.

From the parameters observed when a shell becomes transparent one can estimate the mass of the shell, as the optical thickness of the shell is then expected to be

 \begin{displaymath}\tau = \kappa~\rho~\Delta R_{\rm s} \simeq \frac{2}{3}
\end{displaymath} (4)

where $\Delta R_{\rm s}$ is the radial thickness of the shell. If $\Delta R_{\rm s}$, temperature, T, and $\kappa(T,\rho)$ are known, the density, $\rho$, can be derived from Eq. (4). Then, knowing the shell radius, $R_{\rm s}$, the mass of the shell, $M_{\rm s}$, can be obtained.

As discussed above the first shell started to become partely transparent on day $\sim $58. We can thus assume the effective radius and temperature obtained at this epoch (see Table 1) as estimates of the shell radius and temperature, i.e. $R_{\rm s} = 1.3 \times 10^3~R_{\odot}$ and T = 4500 K. From a difference between curves (b) and (c) in the bottom panel of Fig. 2 at this epoch, we can estimate that $\Delta R_{\rm s} \simeq 0.2~R_{\rm s}$. Adopting the Rosseland mean opacity as tabulated in Alexander & Ferguson (1994) we finally obtained $M_{\rm s} \simeq 0.08~M_{\odot}$.

The second shell presumably became transparent between days 100-120. Assuming that Eq. (4) was satisfied on day 100 we have (see Table 1) $R_{\rm s} = 2.5~ \times~ 10^3~R_{\odot}$ and T = 3500 K, which results in $M_{\rm s} \simeq 0.45~M_{\odot}$. If we assume that the shell became partly transparent on day 120 and that the results obtained from fitting a supergiant spectrum (see Table 1) are estimates of the shell parameters, i.e. $R_{\rm s} = 3.0 \times 10^3~R_{\odot}$ and T = 2500 K, the result is $M_{\rm s} \simeq 0.65~M_{\odot}$.

The above estimates show that the second shell was $\sim $7 times more massive than the first one. However, the obtained absolute values probably overestimate the real shell masses. The reason is that the Rosseland mean, which very well approximates the true opacity in optically thick conditions, is likely to underestimate the real opacity in partly transparent media, particularly in conditions where the Planck mean is much higher than the Rosseland mean. This is the case for the temperature and density range considered above where, according to the results of Alexander & Ferguson (1994), the Planck mean is $\sim $2 orders of magnitude higher than the Rosseland mean. If the Planck mean is used instead of the Rosseland mean in Eq. (4) the results are $M_{\rm s} \simeq 8 \times 10^{-4}~M_{\odot}$ and $\simeq$ $ 4 \times 10^{-3}~M_{\odot}$for the masses of the first and second shell, respectively. These values can be considered as lower limits to the shell masses.

5.3 Decline: Epochs later than $\sim $100 days

This phase is first marked by a continuous decline in luminosity. This decline was initially fast. For about 50 days V838 Mon faded by an order of magnitude. Next the time scale became longer and after almost three years the object was $\sim $40 times less luminous than in eruption.

The initial decline, untill day $\sim $230, was also characterised by a decreasing effective temperature which dropped to $\la$2000 K. The results of our analysis during this phase are particularly uncertain. The spectral fitting was based on an uncertain extrapolation well beyond the range of the standard supergiant spectra. Also the blackbody spectra were not successful in fitting the optical and infrared magnitudes at the same time. The evolution of the effective radius was particularly uncertain, as can be seen from the bottom panel of Fig. 2.

In the later decline, i.e. days $\ga$250, the object shrunk in effective radius. Also the trend in the effective temperature was reversed. The object, although still very cool for a star, slowly became hotter. If the observed evolution of V838 Mon is plotted on the HR diagram the conclusion is clear: in the later decline the object declined along the Hayashi track. The only known stellar eruptive object that shows a similar evolution is V4332 Sgr (Tylenda et al. 2005a) which is usually considered to be of the same class as V838 Mon.

In the stellar evolution the only case of decline along the Hayashi track is the protostellar phase when the main source of energy is gravitational contraction. This strongly suggests that V838 Mon in the late decline also is powered by a gravitational collapse of its envelope inflated during eruption. As shown in Tylenda et al. (2005a) this hypothesis can satisfactorily explain the observed evolution of V4332 Sgr for almost 10 years after its eruption. Within this hypothesis one can estimate the mass of the inflated envelope.

Appendix A presents an approach to the problem of a gravitationally contracting envelope. The structure of the envelope is treated in a polytropic approximation which seems to be better than the r-5/2 density distribution assumed in Tylenda et al. (2005a). For the effective temperatures observed during decline it is reasonable to assume that the envelope is convective so n = 3/2. The dashed curves in the middle and bottom panels of Fig. 2 were derived from a numerical integration of Eq. (A.20) (using Eqs. (A.12) and (A.13)). We found that $T_{\rm eff} = 6120 {-} 1200~\log~ (R/R_{\odot})$is a good approximation to our results obtained from fitting the supergiant spectra to the observations for days $\ga$250. This relation was used in our modelling of the envelope contraction. $M_\ast = 8~M_{\odot}$ and $R_\ast = 5~R_{\odot}$ have been assumed as the parameters of the central star. The mass of the envelope was chosen to fit the observed radius and in the case of the dashed curves in Fig. 2 it is $M_{\rm e} = 0.22~M_{\odot}$.

6 Infrared excess

Infrared spectra of V838 Mon have been analysed by Lynch et al. (2004). The observations taken in January 2002 (i.e. in pre-eruption) they have been able to fit with a single blackbody of 2700 K. They, however, note a possible small excess near 10 $\mu $m. The data obtained a year later, between December 2002 and February 2003 (i.e. in decline), show a more complex nature. In particular, the spectrum shows strong absorption bands due to numerous molecules characteristic of an oxygen-rich ( ${\rm C/O} < 1$) matter. Lynch et al. (2004) conclude that the overall spectrum at short wavelengths can then be fitted with a 2000 K blackbody, whereas at longer wavelengths (8-13 $\mu $m) the continuum indicates 650-800 K.

Table 2: Infrared excess in V838 Mon: data and results. See text for explanations of the columns.

We have collected the available measurements in the range of $\lambda \ga 5~\mu$m. With the shorter wavelength photometry (discussed in Sect. 4) this allowed us to investigate the presence of an infrared excess for six epochs, as given in Col. (1) of Table 2. Two of them are in pre-eruption, one in eruption and three in decline. The data and results are summarized in Table 2 and Fig. 4.

\includegraphics[width=6.5cm]{2800f4b.eps} \end{figure} Figure 4: Infrared excess in pre-eruption and eruption ( left panel) as well as in decline ( right panel). Symbols - observed magnitudes. Full curves - fitted spectra. Dotted and dashed curves ( right panel) - stellar and dust spectral components, respectively. The curves are labelled with the epoch of observations (days since 1 Jan. 2002). A downward arrow for day 258 shows the upper observational limit at 18 $\mu $m. Note that the curve and the data on day 15.0 ( left panel) have been shifted by 2 mag, while those corresponding to days 395 and 650 ( right pannel) have been shifted by 4 and 8 mag, respectively. The data sources and the parameters of the fitted spectra can be found in Table 2.
Open with DEXTER

Standard supergiant colours are poorly known for the ${\it MN}$ bands. Therefore in order to determine the stellar contribution for $\lambda \ga 5~\mu$m we have used blackbody distributions fitted to shorter wavelengths. For the decline period these are simply the same blackbody spectra as those discussed in Sect. 4 and whose parameters are given in Table 1. They are shown with dotted curves in the right panel of Fig. 4. In the case of pre-eruption and eruption we have fitted blackbodies to the observed magnitudes shortward of 5 $\mu $m. They are shown with full curves in the left panel of Fig. 4. The values of the stellar blackbody temperatures, $T_{\rm BB,s}$, are given in Col. (3) of Table 2. Note that the blackbody temperatures obtained for the first three epochs in Table 2 are not significantly different from the values of $T_{\rm eff}$ obtained from fitting supergiant spectra (see Table 1).

As can be seen from Fig. 4 (left panel) the data in pre-eruption and eruption show a 0.5-1.0 mag excess at $\sim $$10~\mu$m. This is in agreement with the analysis of Lynch et al. (2004) of their January 2002 spectrum[*]. As noted in Lynch et al. the 10 $\mu $m excess is likely to be due to silicate emission.

A similar conclusion can also be drawn from a decline spectrum of day 258 (mid-September 2002). For this epoch we have more detailed measurements near 10 $\mu $m from Watson & Costero (2002). As can be seen from Fig. 4 (right panel), the stellar blackbody fits the data well up to $\sim $$8~\mu$m. Above this wavelength the observed flux steeply increases and between 9-13 $\mu $m it stays $\sim $0.8 mag above the stellar continuum. This feature can indeed be interpreted as due to optically thin emission from silicate dust. To show this we have added a component proportional to the dust emission coefficient, i.e. to $Q_{{\rm abs},\nu}~B_\nu(T_{\rm d}$), where $Q_{{\rm abs},\nu}$ is a monochromatic absorption coefficient for "astronomical silicates'' as tabulated in Draine (1985). This component is shown as a dashed curve on day 258 in Fig. 4. The full curve presents the sum of the stellar blackbody and the dust contribution and, as can be seen from the figure, it fits the observations very well. The dust temperature, $T_{\rm d}$, has been assumed to be 800 K as this is the expected temperature of silicate dust condensation. This value is also close to the blackbody temperature estimated in Lynch et al. (2004) and in the present study from the observations in December 2002-January 2003 (day 395, see below). The value of $T_{\rm d}$ cannot be precisely constrained from the observational data on day 258, however it cannot be lower than $\sim $400 K as then the predicted brightness at 18 $\mu $m is above the observational upper limit (shown by a downward arrow in Fig. 4).

We can therefore conclude that in pre-eruption, eruption and early decline V838 Mon showed an infrared excess at $\sim $$10~\mu$m. This excess was not significant for the global energetics of the object and was likely to be due to silicate dust condensation in the wind. The dust was probably optically thin, suggesting a moderate mass loss rate (at the distance of dust condensation). Most likely, following the discussion in Sect. 5.2, this concerns the pre-eruption wind matter expanding above the two dense shells ejected at the onset of eruption.

In later epochs the situation changed considerably. On days 395 (data from observations in December 2002 and January 2003) and 650 (October 2003) the infrared excess increased and cannot be fitted by an optically thin silicate emission. In the former epoch Lynch et al (2004) note a silicate feature at $\sim $$10~\mu$m but it is superimposed on a blackbody of 650-800 K. In the latter epoch, as noted in Lynch et al. (2003), the feature became even less pronounced with a continuum emission dominating the observed spectrum.

We have been able to fit the whole spectrum of V838 Mon on days 395 and 650 by a sum of two blackbody distributions, as shown in Fig. 4 (right panel). The parameters of the fits can be found in Table 2. The cooler component (dashed curves in the figure) can be interpreted as due to an optically thick dusty region reprocessing the radiation from the central object. The fact that we see the central object signifies that dust does not significantly obscure it. In other words, the covering factor of the dusty region is <1, in accordance with the discussion below.

As can be seen from Table 2, the dust temperature, $T_{\rm BB,d}$, decreased while the radius, $R_{\rm d}$, increased between the two epochs. This can be easily interpreted as due to expansion of the dusty region.

The values of log  $\theta_{\rm d}$ and $R_{\rm d}$given in Table 2 are lower limits to the dust distance from the central object, as they have been derived assuming a radiating sphere. The values of $L_{\rm d}/L_{\rm s}$ given in the last column of Table 2 can be considered as an estimate of the covering factor of the dusty region. If this is the case, the dust distance should be increased by the inverse square root of the covering factor and $R_{\rm d}$becomes $1.17~ \times~ 10^4~R_{\odot}$ and $1.78 \times 10^4~R_{\odot}$ on days 395 and 650, respectively. These values of $R_{\rm d}$ and the derived temperature values closely follow the expected relation $T_{\rm BB,d}/T_{\rm BB,s} = (R_{\rm s}/R_{\rm d})^{1/2}$($R_{\rm s}$ taken from Table 1). Assuming that, as discussed in Sect. 5.2, the dominant mass loss was initiated on day $\sim $30 one can estimate a mean expansion velocity from the above values of $R_{\rm d}$. The result is 260  $\mbox{km}\ \mbox{s}^{-1}$ and 230  $\mbox{km}\ \mbox{s}^{-1}$ for days 395 and 650, respectively. We can therefore interpret the strong infrared excess observed in the later decline as due to dust formed in the second, more massive, shell ejected in eruption, as discussed in Sect. 5.2. From the condition that the shell should be optically thick we can estimate a lower limit to its mass. Following Alexander & Ferguson (1994) the opacity $\kappa \simeq 2.7$for $T \simeq 700$ K (the Rosseland mean, but the Planck mean is not significantly different). Taking the shell radius of $1.8 \times 10^4~R_{\odot}$, a radial thickness of 20% and the filling factor of 0.5 the condition $\tau \ga 1$ implies $M_{\rm s} \ga 2 \times 10^{-3}~M_{\odot}$. This value is consistent with the estimates obtained in Sect. 5.2 ( $4 \times 10^{-3} {-} 0.6~M_{\odot}$).

The high value of $T_{\rm BB,d}$ in Dec. 2002-Jan. 2003 (800 K) suggests that dust was being formed during this epoch. Probably the whole process of dust formation was completed in a somewhat later epoch. This might be the reason why the covering factor as measured by $L_{\rm d}/L_{\rm s}$ on day 650 was larger than that on day 395.

7 Summary

We present the following description of the eruption of V838 Mon observed in 2002.

The eruption consisted of two distinct phases, pre-eruption and eruption. Both phases are likely to have resulted from rather short-lived but strong energy bursts.

The first one, which produced the pre-eruption, occured before the discovery of the object so we have no data for it. Presumably it took place in the last days of December 2001 and its main effect was to produce an envelope inflated to $\sim $ $350~R_{\odot}$. In January 2002 the object was in a relatively quiet phase, radiating at $\sim $ $7{-}8 \times 10^4~L_{\odot}$. The slowly decreasing effective temperature and luminosity suggest relaxation after the initial event. The object had a wind that built up an optically thin envelope, extending to $\sim $ $1500~R_{\odot}$ in last days of January 2002.

Presumably, another eruption event occured on 28-30 January at the base of the inflated envelope. This second energy burst was much stronger than the first one and generated the main phase, eruption. As a result an energy wave, probably in the form of shocks and/or a radiation wave, was formed, followed by a massive outflow of matter. The energy wave reached the photosphere on 2 February and heated up the extended envelope built up by the pre-eruption wind. This increased the opacity and made the inner wind layers optically thick, which was observed as a rapid expansion of the effective photosphere. In 3 days the effective radius increased to $\sim $ $700~R_{\odot}$, while the luminosity reached $\sim $ $1 \times 10^6~L_{\odot}$. Subsequent cooling of the wind envelope resulted in a decrease in luminosity and a contraction of the photosphere. On 10-11 February the first shell of mass outflow met the contractiong photosphere and became visible. This first, less massive shell, likely to have been produced from the pre-eruption envelope compressed by the second energy burst, was observed as an expanding and cooling photosphere untill $\sim $27 February. Then it became partly transparent, so that deeper and hotter layers started to be visible. In this way the main mass outflow was observed between 4 March and $\sim $10 April as an initially hotter, but then expanding and cooling, photosphere radiating at a luminosity of $\sim $ $1 \times 10^6~L_{\odot}$. Thus we show that the evolution of the object in the eruption phase, i.e. in February-mid-April, was due to a single event rather than a series of outbursts as often suggested in the literature.

Because of no significant energy supply from internal sources the object declined in luminosity when the shell of the main outflow became transparent and radiated away most of its internal energy. A certain part of the matter involved in the eruption did not gain enough energy to overcome the gravitational field of the central star. It formed an inflated envelope whose gravitational contraction was the main source of the luminosity observed in the decline phase. We estimate that the mass of this envelope is $\sim $ $0.2~M_{\odot}$.

At the same time the matter ejected during eruption expanded and cooled, so that dust formation was possible. As a result a strong infrared excess at $\lambda \ga 3~\mu$m was observed since October-November 2002. We estimate that the total mass lost from V838 Mon is $\la$ $ 0.6~M_{\odot}$. This is rather a conservative upper limit. The true value is likely to be significantly lower.

We postpone the discussion of the present results in terms of possible eruption mechanisms to another paper.

The author thanks to Noam Soker for discussions and comments on the text of the paper. This work was partly supported from grant No. 2.P03D.002.25 financed by the Polish State Committee for Scientific Research.



Online Material

Appendix A: Polytropic envelope

Suppose that a star of a mass, $M_\ast$, and a radius, $R_\ast$, has an envelope of mass, $M_{\rm e}$, inflated to a radius, $R_{\rm e}$. The envelope is in hydrostatic equilibrium. Suppose also that $M_{\rm e} \ll M_\ast$, so the gravitational potential within the envelope, V(r), is

 \begin{displaymath}V(r) = -\frac{G~M_\ast}{r}\cdot
\end{displaymath} (A.1)

Let us assume a polytropic equation of state, i.e.

 \begin{displaymath}P = K~\rho^{1+(1/n)},
\end{displaymath} (A.2)

where P is pressure, $\rho$ - density, K - polytropic constant, n - polytropic index. (A detailed theory of polytropes can be found in e.g. Chandrasekhar 1957.) Then the standard equation of hydrostatic equilibrium can be integrated giving

 \begin{displaymath}(n+1)\frac{P}{\rho} = \frac{G~M_\ast}{r} - \frac{G~M_\ast}{R_{\rm e}}\cdot
\end{displaymath} (A.3)

Equation (A.3) can be used with Eq. (A.2) to obtain the density distribution, i.e.

 \begin{displaymath}\rho = \rho_0 \left( \frac{x_0}{1-x_0} \right)^n \left( \frac{1-x}{x} \right)^n,
\end{displaymath} (A.4)

where $x = r/R_{\rm e}$, $x_0 = R_\ast/R_{\rm e}$, and $\rho_0 = \rho(x_0)$is the density at the base of the envelope. The mass of the envelope, $M_{\rm e}$, is

 \begin{displaymath}M_{\rm e}\! =\! \int_{R_\ast}^{R_{\rm e}}\! 4\pi~r^2\rho~{\rm...
...o_0 \left(\! \frac{x_0}{1-x_0} \!\right)^n
\!I_{\rm m}(n,x_0)
\end{displaymath} (A.5)


 \begin{displaymath}I_{\rm m}(n,x_0) = \int_{x_0}^{1} (1-x)^n~x^{2-n}~{\rm d} x .
\end{displaymath} (A.6)

The potential energy of the envelope, $\Omega$, is
$\displaystyle \Omega = \int_{R_\ast}^{R_{\rm e}} 4\pi~r^2~V(r)~\rho~
{\rm d} r
...~M_\ast~R_{\rm e}^2~\rho_0
\left( \frac{x_0}{1-x_0} \right)^n I_{\rm e}(n,x_0),$     (A.7)


 \begin{displaymath}I_{\rm e}(n,x_0) = \int_{x_0}^{1} (1-x)^n~x^{1-n}~ {\rm d} x .
\end{displaymath} (A.8)

Using Eq. (A.5), Eq. (A.7) can be rewritten as

 \begin{displaymath}\Omega = - \frac{G~M_\ast~M_{\rm e}}{R_{\rm e}}
\frac{I_{\rm e}(n,x_0)}{I_{\rm m}(n,x_0)}\cdot
\end{displaymath} (A.9)

According to the virial theorem the internal energy of the envelope, U, is

 \begin{displaymath}U = - \frac{1}{2} \Omega.
\end{displaymath} (A.10)

Thus the total energy of the envelope, E, is

 \begin{displaymath}E = \Omega + U = - \frac{G~M_\ast~M_{\rm e}}{2~R_{\rm e}}
\frac{I_{\rm e}(n,x_0)}{I_{\rm m}(n,x_0)}\cdot
\end{displaymath} (A.11)

In astrophysical applications two values of the polytropic index, i.e. n = 3/2 and n = 3, are of particular interest. The integrals in Eqs. (A.6) and (A.8) can be analytically evaluated in these cases giving
$\displaystyle { I_{\rm m}(3/2,x_0) = \frac{\pi}{16} - \frac{1}{8}~\arcsin~(x_0^{1/2})
- \frac{1}{8}~x_0^{1/2}~(1-x_0)^{1/2} {} }$
    $\displaystyle {} \qquad - \frac{1}{12}~x_0^{1/2}~(1-x_0)^{3/2}
+ \frac{1}{3}~x_0^{1/2}~(1-x_0)^{5/2},$ (A.12)

$\displaystyle { I_{\rm e}(3/2,x_0) = \frac{3}{8}\pi - \frac{3}{4}~\arcsin~(x_0^{1/2})
- \frac{3}{4}~x_0^{1/2}~(1-x_0)^{1/2} {} }$
    $\displaystyle {} \qquad \qquad - \frac{1}{2}~x_0^{1/2}~(1-x_0)^{3/2},$ (A.13)

 \begin{displaymath}I_{\rm m}(3,x_0) = - \frac{11}{6} - \ln~x_0 + 3 x_0 - \frac{3}{2} x_0^2
+ \frac{1}{3} x_0^3,
\end{displaymath} (A.14)

 \begin{displaymath}I_{\rm e}(3,x_0) = \frac{1}{x_0} + \frac{3}{2} + 3~\ln~x_0 - 3 x_0
+ \frac{1}{2} x_0^2.
\end{displaymath} (A.15)

Note that in the case of extended envelopes ($x_0 \ll 1$) the envelope energy becomes

 \begin{displaymath}E \simeq - \frac{3~G~M_\ast~M_{\rm e}}{R_{\rm e}}
\qquad \mbox{for} \qquad n = 3/2
\end{displaymath} (A.16)


 \begin{displaymath}E \simeq - \frac{G~M_\ast~M_{\rm e}}{2~R_\ast~\ln(R_{\rm e}/R_\ast)}
\qquad \mbox{for} \qquad n = 3.
\end{displaymath} (A.17)

Assuming that the luminosity of the envelope, L, is due to gravitational contraction of the envelope and neglecting mass loss from the envelope, we can write

$\displaystyle { L = - \frac{{\rm d}E}{{\rm d}t} {} }$
    $\displaystyle {} = \frac{G~M_\ast~M_{\rm e}}{2~R_{\rm e}}
\left[ \frac{{\rm d}}...
...(n,x_0)}{I_{\rm m}(n,x_0)}~
\frac{{\rm d}}{{\rm d}t}
~ (\ln R_{\rm e}) \right].$  

On the other hand, we have

 \begin{displaymath}L = 4\pi~R_{\rm e}^2~\sigma~T_{\rm eff}^4,
\end{displaymath} (A.19)

where $T_{\rm eff}$ is the effective temperature. Combining Eqs. (A.18) and (A.19) we obtain
$\displaystyle \frac{{\rm d}}{{\rm d}t} (\ln R_{\rm e})
= \frac{{\rm d}}{{\rm d}...
...t~M_{\rm e}}
\left( \frac{I_{\rm e}(n,x_0)}{I_{\rm m}(n,x_0)} \right)^{-1}\cdot$     (A.20)

If the behaviour of $T_{\rm eff}$ is known, Eq. (A.20) can be numerically integrated to obtain the evolution of the envelope radius.

Copyright ESO 2005