EDP Sciences
Free Access
Issue
A&A
Volume 603, July 2017
Article Number A67
Number of page(s) 26
Section Interstellar and circumstellar matter
DOI https://doi.org/10.1051/0004-6361/201730385
Published online 07 July 2017

© ESO, 2017

1. Introduction

The yet unknown physical mechanisms responsible for the onset of asphericity and polar acceleration in planetary nebulae (PNe) are already active in early stages of the evolution beyond the asymptotic giant branch (AGB). Therefore, pre-PNe (pPNe) and young PNe (yPNe) hold the key to understanding the complex and fast (~1000 yr) nebular evolution from the AGB toward the PN phase. Some studies of pPNe support the idea that the multiple lobes and high-velocities observed are produced by the impact of collimated fast winds (CFWs or jets) on the spherical and slowly expanding circumstellar envelopes (CSEs) formed in the previous AGB phase (see, e.g., Balick & Frank 2002, for a review). However, this jet+“AGB CSE” two-wind interaction scenario remains unconfirmed by direct characterization of the post-AGB jets themselves and of the central nebular regions from which these jets would be launched (within ~few×100 au). Studying these central regions is difficult because of their small (sub-arcsec) angular sizes and, most importantly, because they are usually heavily obscured by optically thick circumstellar dust shells or disks. Moreover, post-AGB ejections, believed to be the main PN-shaping agents, are short-duration phenomena (with lifetimes of 100 yr) and they are already inactive (or, at least, much less energetic) in evolved PNe. Therefore, progress requires sensitive observations at long wavelengths of the crucial, short-lived pPN/yPN phases.

Table 1

Parameters of the sources observed in our pilot survey.

Central stars of pPNe start ionizing their surroundings at the mid-stages of their post-AGB evolution toward the central star of PN (CSPN) phase, typically when they reach a B-type spectral classification. Optical spectroscopic observations of pPNe and yPNe (e.g., Sánchez Contreras et al. 2002, 2008; Arrieta & Torres-Peimbert 2003, and references therein) have revealed the widespread presence of broad (FWZI ≈ 100–1000 km s-1) Hα emission, often with blue-shifted absorption features (P-cygni profiles) that are produced at their nuclei. In some pPNe, the Hα emission suggests the presence of active fast shocks in the stellar vicinity that could be linked, for example, to ongoing jets sculpting the inner regions of slower and older mass ejecta. Owing to the large amount of dust in most pPNe/yPNe, which is mainly concentrated in the equatorial regions, at optical wavelengths the light from the ionized nucleus can only be observed indirectly after it is reflected by the dust in the lobes. This situation not only impedes line detections but also complicates the interpretation of the observed scattered-Hα profiles enormously (Sánchez Contreras & Sahai 2001; Arrieta et al. 2005).

The just emerging central ionized cores of pPNe/yPNe, however, can be traced by radio recombination line (RRL) emission with the important advantage that dust extinction effects are minimal. The theory of RRL emission is well understood (e.g., Brown et al. 1978), however, RRL observations are challenging, especially in the cm-wavelength range where these lines are extremely weak in comparison to the continuum level (typically 1%). However, in the mm-wavelength range larger line-to-continuum flux ratios are expected, i.e., above ~80% and ~10% at 1 and 3 mm, respectively. To date, some mm-wavelength, but mainly cm-wavelength, RRL studies, have been carried out for some of the best-known and most luminous evolved PNe (e.g., Roelfsema et al. 1991; Bachiller et al. 1992; Vázquez et al. 1999). However, this field remains largely unexplored for pPNe/yPNe, where central ionized regions are in an early stage of development.

To our knowledge, RRLs have only been detected in the mm-wavelength range (H41α, H35α, and H30α) toward one pPN, the C-rich CRL 618 (Martín-Pintado et al. 1988). As shown by these authors, mm-wavelength RRLs (mm-RRLs) are excellent tools to study the structure, physical conditions (electron temperature, Te, and density, ne), and kinematics of the dense, central ionized regions of pPNe/yPNe. Indeed, even though single-dish observations do not enable us to spatially resolve the structure of the compact ionized cores, to some extent it is possible to constrain the properties of the current (post-AGB) mass-loss process from detailed radiative transfer modeling of RRL observations. Then, RRLs allow us to study post-AGB mass-loss, which remains very poorly known even when it is key to understanding the physics mediating the AGB-to-PN evolution.

In this paper, we report the results from our pilot study of the emerging central ionized regions of a sample of pPNe/yPNe candidates by means of single-dish observations of RRLs at mm-wavelengths. In Sect. 2, we describe our sample and provide a brief, yet somewhat detailed, introduction to the three sources with mm-RRL detections. The observations and observational results are described in Sects. 3 and 4, respectively. The analysis of the data, which includes radiative transfer model of the line and continuum emission, and the derived results are reported in Sect. 5. Results are interpreted and discussed in Sects. 6 and 7 and our main conclusions are summarized in Sect. 8.

2. Sample

A small sample of eight targets (Table 1) was selected as follows. First, we chose pPNe/yPNe candidates with compact Hα emission signalling active ionized winds at their cores (Sánchez Contreras et al. 2008; Arrieta & Torres-Peimbert 2003; van de Steene et al. 2000) observable with the IRAM 30 m dish. In contrast to evolved PNe, radio continuum flux measurements are lacking for the vast majority of pPNe/yPNe, however, we carefully dug into the literature and publicly available data archives to build the spectral energy distributions (SEDs) of our targets. We preferentially selected objects with indications of positive spectral indexes in the mm-to-cm range (Sννα, with α> 0), indicative of partially optically thick continuum emission. A positive spectral index is expected in the case of an isothermal spherical wind with a power-law radial density gradient (e.g., Panagia & Felli 1975; Reynolds 1986). In this case, the peak-line intensities of RRLs steeply increase with frequency, optimizing the chances of detection at mm-wavelengths. We also made a final sieve by choosing objects with the strongest Hα line and continuum fluxes. All targets in our sample have relatively hot central stars with effective temperatures in the range Teff~ 20 000–40 000 K, which correspond to early-B to late-O spectral types (references are given in Table 1).

In Fig. 1, we show the SEDs of the targets in our sample. In the mm-wavelength region, the continuum emission is a mixture of thermal dust and free-free emission. In the case of the well studied pPN CRL 618, with a strong variability of the radio continuum flux reported over the past three decades (e.g., Kwok & Feldman 1981; Martín-Pintado et al. 1988; Sánchez Contreras et al. 2004; Tafoya et al. 2013, and Sect. 2.1), we have used coeval and almost coeval mm-to-cm data as compiled by Tafoya et al. (2013) and Planck Collaboration Int. XVIII (2015). For the rest of the sources, ancillary data are not necessarily coeval and a certain degree of variability of the radio continuum cannot be ruled out; this could partially explain the scatter of the cm-continuum data points in some targets. Indeed, this is most likely the case of He 3-1475 as noticed by Cerrigone et al. (2011), who found a periodic pattern in the 3.6 cm continuum light curve of this object.

thumbnail Fig. 1

Spectral energy distributions of our targets showing literature photometry data (black circles and triangles indicate detections and upper limits, respectively) and the 1.3 and 2.7 mm continuum fluxes measured by us (star-like symbols; Table 1). The error bars represent a conservative 30% flux uncertainty assigned to all continuum data points. The dotted line indicates a modified blackbody fit to the far-IR thermal dust emission, which is subtracted from the total (dust+free-free) emission observed at mm-wavelengths. No attempt has been made to fit the warm dust emission blueward of the SED peak. The dashed line represents a fit of the type Sffνα to the mm-to-cm data points. The blue solid line represents the addition of these two fits.

Open with DEXTER

In the following subsections, we provide a broad introduction to the three sources with mm-RRL detections.

2.1. CRL 618

CRL 618 (aka RAFGL 618 = IRAS 04395+3601 = Westbrook Nebula) is a well-studied C-rich pPN currently going through the very early stages of PN development (e.g., Westbrook et al. 1975; Wynn-Williams 1977; Kwok & Feldman 1981). The optical nebula consists of multiple elongated (jet/finger-like) lobes, roughly oriented in the E-W direction, composed of shock-excited gas rapidly outflowing from the central star with velocities of up to Vexp~ 180 km s-1 (Sánchez Contreras et al. 2002; Trammell & Goodrich 2002; Riera et al. 2011; Balick et al. 2013). The lobes seem to be “cavities” excavated in the AGB CSE via an interaction with a spray of clumps/bullets that resulted from brief, episodic, and asymmetric ejection events along multiple axes about a century ago (Balick et al. 2013; Velázquez et al. 2014; Huang et al. 2016). In the optical, the remnant of the AGB CSE is seen as a large-scale, roughly round reflection halo around the lobes (Trammell et al. 1993; Sánchez Contreras et al. 2002).

Most of the nebular material in CRL 618 is in the form of molecular gas (~0.3 M), which has been very well studied based on high-angular resolution submm/mm-wavelength emission maps of CO and other molecules (Sánchez Contreras & Sahai 2004; Sánchez Contreras et al. 2004; Lee et al. 2013a,b; Nakashima et al. 2007). The molecular envelope includes fast compact outflows that are aligned with the optical lobes with expansion velocities increasing linearly with the distance to the nebula center, up to Vexp~ 340 km s-1 (at the lobe tips). In addition, there are, at least, three distinct structures expanding at low velocity: a large equatorial torus with a dense compact core (Vexp 12 km s-1), a thin-walled, bipolar shell encompassing the shock-excited optical lobes (Vexp~ 22 km s-1), and an extended tenuous halo (Vexp~ 16 km s-1) surrounding all other components.

The central star, which is hidden from direct view at optical wavelengths by a dusty equatorial torus, has been classified as ~B0 with an effective temperature of 30 000 Teff 40 000 K (Westbrook et al. 1975; Schmidt & Cohen 1981; Trammell et al. 1993; Balick et al. 2014).

Around the central star, there is a compact H ii region that produces free-free continuum emission (e.g., Wynn-Williams 1977; Kwok & Feldman 1981; Kwok & Bignell 1984; Martín-Pintado et al. 1988). At cm-wavelengths, the ionized region shows an overall elliptical brightness distribution elongated in the E-W direction (Martín-Pintado et al. 1993). The continuum flux and size of the ionized core at cm-wavelengths have been increasing monotonically in the past three decades. In their discovery paper, Kwok & Feldman (1981) interpreted this result as the advance of the ionization front within the molecular envelope as a consequence of the increasing stellar Teff. The rapid angular expansion of the nuclear H ii region, which had reached an angular size of ~0.̋5 × 0.̋2 (~450 × 180 au at d = 900 pc) at 22 GHz in 2008, has recently been studied in detail by Tafoya et al. (2013). These authors use high-angular resolution maps of the ~1.5–43 GHz continuum emission at multiple epochs spanning over 25 yr. Tafoya et al. (2013) demonstrate that the advance of the ionization front cannot explain the current expansion of the optically thick continuum-emitting region. These authors conclude that it is the optical depth of the free-free continuum at a given distance that increases with time, which is attributed to the rise of the electron density with time. This supports the idea that CRL 618 is expelling gas in the form of an ionized wind whose mass-loss rate has been increasing, from ~ 4 × 10-6 to 6 × 10-6M yr-1, during the last century. Tafoya et al. (2013) establish the onset of the ionization around 1971.

The free-free continuum at mm-wavelengths, which is optically thin and, thus, probes denser and deeper (150 au) regions than the continuum at cm-wavelengths (optically thick), displays much abrupt and unpredictable flux changes that cannot be simply attributed to the smooth growth of the optically thick layers of the H ii region but most likely denote alterations in the activity and/or physical conditions of the post-AGB wind at its core (Sánchez Contreras & Sahai 2004; Sánchez Contreras et al. 2004, see also Fig. 2).

CRL 618 is the only pPN in which mm-RRL emission has been reported prior to this work: the H41α, H35α, and H30α transitions (at 3.5, 2, and 1.3 mm, respectively) have been detected by Martín-Pintado et al. (1988). From the analysis of their single-dish spectra, which included LTE radiative transfer modeling of the line and continuum emission, these authors conclude that mm-RRLs trace a Te~ 13 000 K stellar wind with a terminal expansion velocity of Vexp 20 km s-1 and a mass-loss rate of ~ 7.6 × 10-6M yr-1.

The distance to CRL 618 adopted in this work is d = 900 pc, as determined using independent techniques by Goodrich (1991), from a total luminosity criteria, and by Sánchez Contreras & Sahai (2004), from the analysis of proper motions of the fast expansion of the molecular gas in the lobes. The inclination of the nebular axis with respect to the plane of the sky is i ~ 30° (Sánchez Contreras et al. 2002; Sánchez Contreras & Sahai 2004).

2.2. M 2-9

M 2-9, also known as the “Butterfly” or the “Twin Jet” nebula, is a bright bipolar nebula discovered by Minkowski (1947) and considered by many authors to be in the earliest stages of becoming a PN. In the optical, it shows a bright compact core from which two bilobed nested structures emerge oriented nearly north-south (Schwarz et al. 1997; Clyne et al. 2015). The brightness and morphology of the optical nebula have been changing notably on a timescale of a few years or less since its discovery (e.g., Allen & Swings 1972; van den Bergh 1974; Kohoutek & Surdej 1980; Doyle et al. 2000). In particular, the main bright condensations (knots) in the lobes progressively shift in the west-east direction, which has been most recently interpreted as produced by a rotating collimated spray of high velocity particles (jet) from the central binary system, which excites the walls of the inner cavity of M 2-9 (Corradi et al. 2011). From the rotation period of the beam particles, these authors infer an orbital period of P ~ 90 yr.

thumbnail Fig. 2

Light curve of the 1 mm-continuum of CRL 618 using ancillary data compiled by Sánchez Contreras et al. (2004; black circles); new measurements are added from Reuter et al. (1997), Nakashima et al. (2007), Planck Collaboration Int. XVIII (2015), and this work (red circles). Some reference dates are indicated. Together with significant abrupt changes in short timescales of less than 0.1 yr, a periodic pattern with a period of ~17 yr is hinted (dotted line) on top of an overall, long-term smooth increase.

Open with DEXTER

Optical spectroscopic observations point to an expansive kinematics in the lobes, as for most pPNe and young PNe, with expansion velocities ranging from ~23 km s-1 in the inner nebula up to ~164 km s-1 in the outer lobe tips (Schwarz et al. 1997; Torres-Peimbert et al. 2010; Clyne et al. 2015). The central core is dominated by Hα emission with ~1600 km s-1-wide wings. The broad Hα profile could indicate ongoing ultra-fast winds, although Hα wing emission could be partially (or totally) broadened by other mechanisms, such as Raman scattering (Torres-Peimbert et al. 2010; Arrieta & Torres-Peimbert 2003).

The central exciting star has been estimated to be of spectral type B1 and possibly late O, Teff~ 35 000 K (Calvet & Cohen 1978; Swings & Andrillat 1979). As in the case of CRL 618, the large extinction and high-infrared excess observed are attributed to a dusty torus in our line of sight that prevents us from seeing the central star directly but allows it to illuminate the lobes. The binarity of the nucleus of M 2-9 has been observationally established by Castro-Carrizo et al. (2012), hereafter CC12. Interferometric 12CO (J = 2 − 1) emission maps with subarcsecond-resolution by these authors show two coaxial, slightly off-centered rings lying in the equatorial plane of the nebula. These rings have been interpreted as the result of two short mass ejections produced at different positions in the binary orbit. The rings lie at a plane that is inclined with respect to the line of sight by ~17 ± 1° (Castro-Carrizo et al. 2017).

M 2-9 is known to be a source of thermal radio continuum emission. According to Kwok et al. (1985), two components are needed to properly fit the cm-wavelength data: a compact ionized wind at the center that produces free-free emission following a Sνν~ [ 0.6 − 0.7 ] law, and high-density condensations located in the extended ionized lobes that produce a nearly flat continuum emission distribution. At mm-wavelengths, the contribution from the extended lobes to the observed continuum is small and the free-free emission from the compact ionized core together with thermal emission by dust dominate the observed continuum flux (Sánchez Contreras et al. 1998, and Fig. 1). At cm-wavelengths, the continuum brightness distribution is elongated along the N-S (lobe’s) direction (Kwok et al. 1985; Lim & Kwok 2000, 2003). At mm-wavelengths, the ionized core is expected to be very compact, 0.̋1, but its geometry is unknown. Based on VLA cm-continuum emission maps obtained in 1982, Kwok et al. (1985) deduce a major-to-minor axis ratio that decreases monotonically as the frequency increases, implying a nearly spherical brightness distribution at ~80 GHz and beyond. However, based on a different set of VLA cm-continuum maps obtained in 1999, Lim & Kwok (2003) report a frequency-independent major-to-minor axis ratio of ~2.5 (in the 1.3-to-6 cm wavelength range) that, if extrapolated to mm-wavelengths, would imply an elongated geometry at these high frequencies.

The distance to M 2-9 is very uncertain. Values spanning almost two orders of magnitude, from 50 pc to 3 kpc, are found in the literature. The two most accepted values derived from recent works are d ~ 650 pc, from the analysis of the proper motions of the expanding molecular rings at the nebula equator (CC12), and d = 1.3 ± 0.2 kpc, from the analysis of the proper motions of the dusty blobs in the lobes (Corradi et al. 2011).

At d ~ 1.3 kpc, the total luminosity of M 2-9 is ~3000 L (Lykou et al. 2011). At d ~ 650 pc, the luminosity is a factor four smaller and, therefore, lower than the low limit of ~1000 L expected for the least massive AGB and post-AGB stars (Bloecker 1995, see also Fig. 8).

2.3. MWC 922

The evolutionary status of this, yet barely studied, dust/gas-enshrouded B[e] star is unclear. An evolved nature has been suggested based on similarities of its SED and nebular morphology with another well-known pPN, the Red Rectangle (e.g., Cohen et al. 1975; Waters et al. 1998; Bujarrabal et al. 2016). Based on the appearance of its surrounding nebulosity in the near-infrared, MWC 922 has been nicknamed the Red Square Nebula. Tuthill & Lloyd (2007) report H-band adaptive optics imaging showing an X-shaped nebular structure, which extends for about 5′′, with an equatorial dark band (running along PA ~ 46°) at the center. A similar square-like structure is observed in the mid-infrared (Lagadec et al. 2011). This nebular appearance may result from the projection of two twin opposing biconical lobes excavated in a dense dust-scattered halo, although illumination effects could be significant; for example, if the light from the central source is blocked by an inner dusty torus or disk at the nucleus.

The optical and infrared spectra of MWC 922 show a great amount of emission lines, including many recombination and forbidden Fe ii lines but also transitions from neutral elements, e.g., O i, H2, etc. (Allen & Swings 1976; Rudy et al. 1992; Pereira et al. 2003; Wehres et al. 2014). A B 3-B 6 spectral class for the exciting star has been suggested by Rudy et al. (1992) based on the nebular emission line spectrum. These authors acknowledge that this classification is very uncertain (in fact, it would imply Teff~ 12 000–17 000 K, which is not sufficient to produce significant photoionization) and propose an upper limit for the stellar temperature of Teff< 31 000 K based on the strength of the He iλ10830 Å line. A symbiotic nature is not favored by these authors given the lack of molecular absorption bands (very prominent in late-type giants) in the spectrum of MWC 922.

There is a bright compact radio continuum source associated with MWC 922, which has been studied by Rodríguez et al. (2012) using a 0.̋95 × 0.̋66 angular resolution 3.6 cm radio continuum VLA map. In these maps, the ionized core is angularly resolved with deconvolved dimensions ~0.̋18 × 0.̋20. Continuum observations at 6 cm from the CORNISH project exist (Purcell et al. 2013); the corresponding flux is included in the SED shown in Fig. 1. The 6 cm radio continuum source is unresolved with ~1.̋5-angular resolution.

The distance to MWC 922 is unknown. As mentioned by Tuthill & Lloyd (2007), if it lies within the Ser OB1 association then the distance is d = 1.7 kpc. If this (speculative) association is real then it may be consistent with a pre-main sequence nature for this object, which would still be linked to other young massive OB stars (but see discussion in Sect. 6.3). We estimate the kinematic distance to MWC 922 from its radial velocity, VLSR = +32.5 km s-1 (deduced from the centroids of the RRLs; see Sect. 4), and its galactic coordinates (l = 17.9275°, b = + 00.6336°) by assuming a simple galactic rotation law and adopting a value for the A Oort constant of 14.4 km s-1 kpc-1 and a galactocentric radius of 8.5 kpc (Kerr & Lynden-Bell 1986). We deduce values of ~3 and ~13 kpc, for the near and far kinematic distance, respectively.

The bolometric flux of MWC 922 computed by integrating its SED from ~1 μm to cm-wavelengths is ~5800 L kpc-2. At the far kinematic distance, this would imply a total luminosity on the order of a million L, approaching that of the most luminous and massive (100 M) stars known in the Galaxy. At d = 1.7–3 kpc, the total luminosity is ~[1.7–5.4] × 104L, which is at the high end of the range for post-AGB objects (Bloecker 1995). The foreground-ISM extinction toward the 7× 7 region around MWC 922 is AV ~ 1.7–2.7 mag at a distance of d ~ 2–3 kpc (Zasowski et al. 2015). After applying the corresponding ISM extinction correction, the total luminosity deduced at d = 1.7–3.0 kpc is ~[1.8–5.9] × 104L.

3. Observations and data reduction

Observations have been carried out with the IRAM 30 m radiotelescope (Pico Veleta, Granada, Spain) using the new generation heterodyne Eight MIxer Receiver (EMIR; Carter et al. 2012). Spectra were taken in two observational campaigns in June/July and September 2015 (project 050-15).

We simultaneously used the E090 (at 3 mm) and E230 (at 1 mm) bands operated in dual sideband (2SB) mode. For all targets, we observed the ~8 GHz-wide upper side band in dual (H+V) polarization. This enabled us to cover the ~105–113 and ~226–234 GHz frequency ranges simultaneously and, thus, to observe H39α, H48β, H49β, H55γ, and 13CO (J = 1 − 0) at 3 mm together with H30α and 12CO (J = 2 − 1) at 1 mm. For CRL 618, MWC 922, and M 2-9, we also used a different receiver setup to survey the lower and upper side bands of each receiver at the same time in single (either V or H) polarization. This enabled us to expand the frequency coverage to the ~89–97 and ~210–218 GHz ranges, which include H41α and H51β at 3 mm and H31α at 1 mm.

The image sideband was rejected with an average sideband rejection of ~14 dB. Each receiver was connected to the Fast Fourier Transform Spectrometer (FTS) in its 195 kHz spectral resolution mode (ΔÆ~ 0.25 and 0.5 km s-1 at 1 and 3 mm, respectively). Observations were performed in wobbler switching mode with a wobbler throw of 120′′ and frequencies of 0.5 Hz. Calibration scans on the standard two load system were taken every ~15 min. Pointing and focus were checked regularly on nearby continuum sources. After pointing corrections, the typical pointing accuracy was ~2′′–3′′. Average system temperatures were 130–200 K at 3 mm and 400–1000 K at 1.3 mm. On-source integration times are in the range ~ 2–6.5 h per target and receiver setup.

The beam parameters of the IRAM 30 m telescope used in this work, that is, the half-power beam width (HPBW), the main beam efficiency (ηeff), and the point source sensitivity (S/, i.e.,  the K-to-Jy conversion factor), are described to a good accuracy as a function of the frequency (ν, in GHz) by (1)according to measurement updates performed in August 20131.

We reduced the data using CLASS2 following the standard procedure, which includes killing bad channels, subtracting baseline, and averaging individual good-quality scans to produce the final spectra. We obtained and presented the spectra in antenna-temperature () scale, which can be converted to main-beam temperature (TMB) applying TMB = /ηeff or to a flux scale using the K-to-Jy conversion factor given above.

The uncertainty in the relative calibration of our observations was estimated by comparing the spectra in the same spectral ranges, including the 12CO (J = 2 − 1) and 13CO (J = 1 − 0)  transitions, every day. We estimate that the total line calibration uncertainty is 20% and 15% at 1 mm and 3 mm, respectively. The profiles of the bonus 12CO (J = 2 − 1) and 13CO (J = 1 − 0) transitions observed in this project are presented in Figs. A.1 and A.2.

We measured the continuum emission level (in a scale) within the different frequency ranges covered in these observations by fitting a low-order polynomial function to the spectral baseline in our data, which were obtained using the wobbler-switching method described above. For the sources that were observed in two epochs (June/July and September 2015), namely, CRL 618, MWC 922, M 2-9, and M 1-92, we obtained continuum measurements separately for both epochs and we find that the continuum levels in the two epochs are in agreement within 30% at 1 mm. The continuum flux uncertainty at 3 mm is probably lower than 20%.

For CRL 618, which is a strong (>2 Jy) millimeter continuum emitter, we also measured the continuum from cross-scans in azimuth and elevation direction performed for pointing. After discarding a few obvious bad/noisy scans, we obtained final, weighted average scans3 for the four frequency ranges observed, centered at 93, 109, 214, and 230 GHz. Gaussian fits were applied to the average scans to derive the antenna temperature at the peak, which was converted to flux density using the frequency-dependent conversion in Eq. (1).

In CRL 618, the line emission flux contribution to the total flux in the continuum backend (8 GHz-wide) used for the pointing cross-scans is estimated ~10% at ~230 GHz, where the strong 12CO (J = 2 − 1) line lies, and less than 4% at the lower frequency ranges observed, and has been considered to derive the (line-free) continuum flux level. The value of the continuum flux obtained from the two methods (baseline fitting and cross-scans) agree to a 20% level.

thumbnail Fig. 3

Recombination lines observed toward the pPN CRL 618 (histogram) and synthetic line profiles (pink) from our model in Table 4 (Sect. 5.1). The bottom and top X-axis represent VLSR (km s-1) and frequency (MHz). The frequency of He and C α-transitions, some of them detected, are also indicated. The many additional features observed in the spectrum, with emission and absorption profile components, are molecular transitions produced in the C-rich molecular envelope of CRL 618. Some of the mm-RRLs observed are partially blended with molecular lines.

Open with DEXTER

thumbnail Fig. 4

Recombination lines observed toward M 2-9 (histogram) and synthetic line profiles (pink) from our model in Table 4 (Sect. 5.2).

Open with DEXTER

thumbnail Fig. 5

Recombination lines observed toward MWC 922 (histogram) and synthetic line profiles (pink) from our model in Table 4 (Sect. 5.3).

Open with DEXTER

thumbnail Fig. 6

Line parameters of the RRLs detected in CRL 618, MWC 922, and M 2-9 derived from the observed profiles (Table 2; filled circles) and from the models in Table 4 including non-detections (empty circles). Dotted lines are fits to the observed line parameters of the Hnα transitions (red filled circles). For CRL 618, we show the mean VLSR of the Hnα centroids (dotted line) and a linear fit of VLSR as a function of frequency (dashed line).

Open with DEXTER

4. Observational results

4.1. Continuum

Our continuum flux measurements at frequencies adjacent to the H39α and H30α lines, that is, at 2.7 and 1.3 mm, respectively, are given in Table 1 and are shown in Fig. 1 (red symbols). In all cases, we find continuum fluxes that are in good agreement with previous measurements at similar wavelengths (whenever available) within the absolute flux uncertainties. A fit to the mm-to-cm free-free continuum fluxes, after removal of the dust thermal emission contribution, was performed and is shown in Fig. 1 (dashed line). In the case of M 2-9, two free-free continuum components, one nearly flat and another one as Sνν0.7, were used for the fit based on previous results at cm-wavelengths (Kwok et al. 1985, and Sect. 2.2). To our knowledge, our IRAM 30 m data represent the first detection of the continuum flux at mm-wavelengths for MWC 922, M 1-91 (tentatively detected at ~1 mm by Sánchez Contreras et al. 1998), IRAS 20462+3416, and M 2-56. For these sources, we find that the mm-to-cm free-free continuum is consistent with Sννα, with values of α in the range ~0.7 and 1.3.

The mm- and cm-wavelength continuum of CRL 618 has been observed many times during the last decades and is known to be variable, alternating periods of brightening and weakening at mm-wavelengths (Sect. 2.1). In Fig. 2 we show an updated light curve of the continuum near 1mm including our recent measurement from this work.

The free-free continuum of CRL 618 is nearly flat at mm-wavelengths, which is consistent with optically thin emission and follows Sνν1.6at longer wavelengths (Fig. 1). Our data (and model presented in Sect. 5) are consistent with a turnover frequency νt 80–90 GHz in 2015, in agreement with the value of νt obtained in 1998 (Wyrowski et al. 2003). At earlier epochs, the turnover frequency of the free-free continuum was larger, in particular, νt ~ 100–200 GHz in 1987 (Martín-Pintado et al. 1988). The decrease of νt with time could imply a decrease of the emission measure from EM ≈ (1.6 ± 0.9) × 1011 to 5.5 × 1010 cm-6pc (adopting Te = 13 000 K). However, this simplistic interpretation is only correct if the physical conditions in, and geometry of, the ionized core have remained unaltered, which seems unlikely given the rapid unforeseen changes of the mm-continuum emission (Fig. 2 and Sect. 2.1); for example, a variation of the electron temperature from Te~ 10 000 to 18 000 K would result in a shift of νt from ~150 to ~100 GHz for the same value of EM ≈ 1.1 × 1011  cm-6 pc.

Table 2

Line parameters from Gaussian profile fitting for RRL detections.

In the case of the pPN M 2-9, our continuum data (both at 1 and 3 mm) appear to be systematically lower than all previous published measurements at these wavelengths. This, together with the scatter of the continuum (ancillary) data points in the mm- and cm-wavelength domain could reflect time variability of the free-free emission in M 2-9.

The pPN Hen 3-1475 is another object in our sample with flux-density variations. A possibly periodic pattern and an inversion of its spectral index on a timescale of a few years have been suggested from previous multifrequency high-angular resolution radio observations (Cerrigone et al. 2011).

In the case of M 1-92, M 1-91, and IRAS 20462+3416, the continuum flux at 3–4 cm clearly exceeds that expected by extrapolating the trends observed at lower wavelengths, including the mm-range. This could be due to contaminating foreground/background emission sources, strong continuum variability, or additional superimposed free-free emission components with different contributions at cm- and mm-wavelengths (e.g., analog to the faint, extended emission from the lobes of M 2-9, Sect. 2.2).

4.2. Recombination lines

We have detected radio recombination line emission at mm-wavelengths for three (out of a total of eight) sources: CRL 618, M 2-9, and MWC 922 (Figs. 35). The main line parameters directly derived from the observed profiles are given in Table 2. In Fig. 6, we plot the flux, full width at half maximum (FWHM), and centroids of the hydrogen RRLs detected as a function of their rest frequencies. The rms noise of the spectra for RRL non-detections are given in Table 3. These undetected sources are rather weak mm-continuum emitters, which anticipate low mm-RRL intensities below our detections limits (see Sect. 5.4).

The strongest RRLs are observed toward CRL 618, which also shows the strongest mm-wavelength continuum. We detect a total of eight hydrogen recombination lines (α, β, and γ) together with helium α-transitions (Fig. 3). The RRLs appear to be partially blended with some of the many molecular transitions with composite emission and absorption profiles that crowd the submm/mm spectrum of this object (see, e.g., Pardo et al. 2007). The profile of the weak H55γ line is the most adversely affected by deep adjacent absorptions.

We find significant changes in the spectral profiles and total fluxes of the H41α and H30α lines of CRL 618 relative to the observations performed in 1987 by Martín-Pintado et al. (1988). In particular, these RRLs are now a factor ~2–3 more intense and ~20–60% broader than 29 yr ago. This is not totally unexpected given the obvious variability of the free-free continuum emission observed over many years (Sect. 2.1 and Fig. 2), and suggests changes with time of the physical conditions at the inner regions around the central source. In M 2-9, weak emission from the H30α and H39α lines and, tentatively, from the H31α and H41α transitions is detected (Fig. 4).

In MWC 922, we detect the four Hnα transitions observed in this work plus the weaker H48β and H49β lines (Fig. 5). In contrast to CRL 618 and M 2-9, in this source there is an evident drastic change from single-peaked to double-peaked profiles in the Hnα lines at 3 and 1 mm. A similar transition from single-peaked profiles of RRLs near 3 mm (and at cm-wavelengths) to double-peaked profiles of the RRLs near 1 mm was observed for the first time toward the ultra-compact H ii region of MWC 349A by Martín-Pintado et al. (1989a,b). These authors showed that, in this case, the two-horned line profiles result from maser amplification of the emission from a rotating ionized disk. As we see in Sect. 5.3, we reach a similar conclusion for MWC 922 that, thus, adds to the (yet short) list of sources where mm-RRL masers have been discovered to date; to our knowledge, these include MWC 349A, η Carinae, Cepheus A HW2, and MonR2-IRS2 ( see Martín-Pintado et al. 1989a; Cox et al. 1995; Jiménez-Serra et al. 2011, 2013; Abraham et al. 2014).

There are a few general observational results that can be readily seen in Fig. 6 and are summarized here:

  • The line strength of the RRLs with the sameΔn increases with the frequency with MWC 922 showing the steepest slope. In particular, for the α-transitions, the slope in this target clearly deviates from the ν1.1 dependency expected in case of optically thin LTE emission (e.g., Rodríguez et al. 2009). This is in agreement with the presence of maser amplification of the high-frequency lines. Non-LTE emission conditions are also inferred in the case of CRL 618 from the α-to-β line ratios, which depart from the theoretical predictions under LTE (e.g., Peters et al. 2012).

  • The line widths of the α-transitions range between FWHM~ 30 and 50 km s-1 and, thus, they are larger than expected from thermal motions alone (e.g., the sound speed is cs ~ 13 km s-1 for ionized gas with an average electron temperature of Te 10 000 K, typical of H ii regions). As discussed in Sect. 5, this indicates macroscopic ordered motions of the bulk of the gas at moderate speeds of 10–30 km s-1. We cannot rule out the presence of broad (>±30 km s-1) emission wings below our detection limits, i.e., at <20% of the line peak.

    Table 3

    rms noise (in units) for RRL non-detections.

  • The Hnα line width increases with increasing quantum number n in CRL 618 and (tentatively) M 2-9, but the opposite behavior is found for MWC 922. The high-n transitions trace more tenuous regions than lower-n lines (at higher frequencies and, thus, with continuum emission optically thinner). Therefore, assuming that the density in the ionized region falls off with the distance to the center (r), the trend observed in CRL 618 and M 2-9 suggests that the gas velocity increases with r, which is consistent with an expanding wind. In the case of MWC 922, the spectra seem more consistent with a rotating disk/wind, where one expects the rotation speed to decrease with r. In CRL 618, the widths of the β-lines are systematically higher than the α-lines at comparable frequencies. This is partially due to line blending and poor baseline subtraction related uncertainties, which happen to be on average more severe in the β-transitions, but also due to larger pressure broadening for the β-lines compared with the α-lines predicted by theory (e.g., Brocklehurst & Seaton 1972; Walmsley 1990; Strelnitski et al. 1996; Peters et al. 2012, and Sect. 5).

  • The centroids of RRLs provide us with an estimate of the line-of-sight velocity of the bulk of the ionized gas (and, plausibly, of the ionizing star that lies at the center). Using the Hnα profiles, with the highest S/N, we find VLSR = − 23.5 ± 0.2 km s-1 for CRL 618, VLSR = + 75 ± 2 km s-1 for M 2-9, and VLSR = + 32.5 ± 0.4 km s-1 for MWC 922. In the case of CRL 618, our data are consistent with a progressive velocity shift toward the blue at low frequencies (Table 2 and Fig. 6). Although tentative, the small velocity shift could be caused by the gradual decrease of the optical depth of the lines with the frequency (e.g.,  Martín-Pintado et al. 1988). In an outflowing gas, the optically thicker, lower frequency RRLs mainly trace a major part of the gas from the front (blue-shifted) of the nebula, since part of the emission from the back side (red-shifted) is absorbed through the line+continuum emitting region. If this is the case, the ~1 mm lines are expected to be centered closer to the line-of-sight velocity of the ionizing star of CRL 618 (see Sect. 5.1).

5. Radiative transfer model

In order to constrain the structure, physical conditions, and kinematics of the emerging ionized regions of CRL 618, MWC 922, and M 2-9, we modeled their free-free continuum and mm-RRL emission using the non-LTE radiative transfer code MORELI (MOdel for REcombination LInes), which is described in detail in Báez-Rubio et al. (2013).

The MORELI code considers a given geometry (with spherical or cylindrical symmetry) for the ionized region. The electron density is assumed to decrease with the distance to the central star following a power law of the type ne(r) ∝ rαn. Also, MORELI takes the electron temperature, Te(r), and a given velocity field as inputs. For simplicity, and since the exact gas kinematics cannot be precisely determined from single-dish observations, in our models the velocity field is described by radial expansion (at constant velocity or increasing, linearly or asymptotically, with r) and/or Keplerian rotation (see Sects. 5.1 to 5.3 for more details). In addition to the gas macroscopic motions, thermal broadening and electron impact (pressure) broadening from high-density gas are included in the line profile function. Turbulence velocity dispersion is not considered as an additional form of line broadening since it is expected to be very small (σturb< 2–3 km s-1) in the envelopes of AGB and post-AGB objects (e.g., Schöier et al. 2004; Bujarrabal et al. 2005; Decin et al. 2010) and, thus, negligible compared to other broadening terms.

The MORELI code predicts the free-free RRLs and continuum emission by solving the equation of radiative transfer and integrating the intensity along the line of sight. The ionized region starts at minimum radius, Rin, which in principle could be the stellar photospheric radius; this region extends up to a distance that is a few times larger than the radius that encircles a major percentage (90%) of the observed fluxes, which is referred to as the effective outer radius Rout4. Thomson electron scattering is commonly neglected, to a first approximation, to solve the radiative transfer problem for recombination lines at millimeter wavelengths (e.g., Peters et al. 2012). The effect of this process is, in general, weak and high electron concentrations of ne 108–1012  cm-3in electron-scattering regions of sizes of a few × 10 au are needed to cause a measurable effect (Arrieta & Torres-Peimbert 2003; Sekeráš & Skopal 2012).

The hydrogen non-LTE level populations are computed in MORELI using departure coefficients, bn, which are defined as the ratio between the actual level population for the principal quantum number, n, and its value in LTE (bnNn/). Computations of the bn coefficients are available in the literature for n as low as 20 and for a broad range of electron densities and temperatures (typically, ne ~ 10–108.5 cm-3and Te~ 5000–15 000 K; Salem & Brocklehurst 1979, Walmsley 1990, Storey & Hummer 1995, and references therein). The bn factors computed by different authors for the same values of Te and ne do not fully agree. The uncertainties in the bn coefficients have a major impact on the accuracy to which we can reproduce the intensity of RRLs where stimulated emission appears to be present. It is particularly critical in the case of the H30α and H31α maser lines in MWC 922 (Sect. 4.2). The bn coefficients from Storey & Hummer (1995) were found to provide the best fit of the maser transitions in MWC 349A (Báez-Rubio et al. 2013) and we systematically used these for our targets.

The relative intensities of the RRLs and the line-to-continuum flux ratios observed in M 2-9 are consistent with LTE conditions (or, at most, with small non-LTE effects) and, therefore, the line profiles were modeled with MORELI in LTE mode in this case. In MWC 922 and CRL 618, however, there was a need for more general non-LTE calculations and we performed these calculations to reproduce simultaneously the α- and β-profiles observed at different frequencies.

Table 4

Parameters of the emerging H ii regions deduced from line and continuum radiative transfer modeling (Sect. 5).

In all cases, we explored a broad range of physical conditions until we obtained a reasonable match between the model predictions and mm-wavelength line and continuum observations. Our model also provides reasonable predictions for radio continuum data at longer wavelengths (whenever available), such as the total cm-to-mm free-free continuum flux and the angular size frequency dependence. In general, we tried to keep the number of model input parameters as low as possible given the data available.

As shown by, for example, Lee et al. (2007) and Báez-Rubio et al. (2013), the mass-loss rate is almost exclusively constrained by the flux of the free-free continuum, which has a moderate observational uncertainty of 30% (Sect. 3). The other two main sources of uncertainty in the mass-loss rate are the expansion velocity and distance to the source. Except for MWC 922 (see Sect. 5.3), the expansion velocity of the bulk of the ionized wind is rather well constrained, with a 20% accuracy, from the widths of the RRLs. For the distances adopted, we find that the range of possible mass-loss rates that can reproduce the observed continuum, and are also consistent with the line profiles, are always within a factor 2 of uncertainty. The spatio-kinematics and electron temperature stratification across the ionized mm-wavelength emitting layers critically determine the shape and intensity of the RRLs, but have a rather limited effect on the continuum flux. The expansion velocities and temperatures are then mainly constrained by the fit of the single-dish profiles and should be taken as representative average values within the ionized regions under study.

The modeling results and their main uncertainties are described for each of the modeled sources in the next subsections.

5.1. CRL 618

Based on previous works (Sect. 2.1), as a simple model of the ionized core of CRL 618, we adopted an isothermal stellar wind represented by a cylinder of ionized gas with a central cavity around the star, constant electron temperature, and an inverse-square density profile. This is similar to the model used by Martín-Pintado et al. (1988) to fit their H41α, H35α, and H30α spectra and continuum emission observed in 1987. The presence of a central cavity, with radius Rin, was originally proposed by these authors to explain with this simple model the nearly flat (optically thin) free-free continuum at mm-wavelengths and, at the same time, the strong continuum emission and large emission measure deduced (EMfew × 1010–1011 cm-6 pc).

The outer radius of the cylinder (along PA ~ 0°) was set to Rout~95 au, which is the size of the minor semi-axis of the H ii region expected in 2015 adopting the growth function θmin = 0.03 × (t − 1971)0.5 deduced by Tafoya et al. (2013). As shown by these authors, the measured minor axis of the cm-continuum, which remains basically constant at all frequencies, pinpoints the locus of the ionization front along these dense equatorial regions. In the direction of the lobes (PA ~ 90°), the ionization front was already at distances 600 au in 1998 and is probably near 800 au at present. In our model, the cylinder extends out to infinity along its major axis. It is not necessary to truncate the cylinder at 600–800 au because the ionized wind layers farther out would not contribute appreciably (<0.001%) to the mm-RRLs/continuum emission.

thumbnail Fig. 7

Same as in Fig. 1 with the free-free continuum emission predicted by our model (red line). Model parameters are given in Table 4 for CRL 618, MWC 922, and M 2-9, and in Table 5 for the rest of the targets.

Open with DEXTER

The input parameters of the best-fit model for CRL 618 are given in Table 4. The synthetic spectral Hnα profiles are overimposed to the observed spectra in Fig. 3. The line peak intensities and their FWHM are well reproduced within a  15%. Our model recreates the progressive blue-shift of the line centroids as the frequency decreases, as observed (see Sect. 4.2). This is due to the moderate line optical depth effects in this source, which produce an overall velocity shift toward the blue relative to the line-of-sight velocity of the ionizing star, adopted to be Vsys(LSR) = 19 km s-1 in our model. The free-free continuum emission in the mm-to-cm regime predicted by our model is also in very good agreement with the observational data (Fig. 7).

We find a characteristic electron temperature of Te~ 16 000–17 000 K, which is somewhat larger than the value obtained assuming LTE conditions by Martín-Pintado et al. (1988). The bulk of the mm-RRLs/continuum emission is produced inside the adopted cylindrical structure within a distance along its major axis (or semi-length) of Rout(90%)~110 au. The electron densities range between ne~ 2 × 107 and 4 × 106 cm-3in the mm-RRLs/continuum emitting region. This is in very good agreement with the values deduced by Tafoya et al. (2013) extrapolating their density law, which applies to the ~150–650 au layers, to the inner wind regions. We find that an inverse-square density profile yields satisfactory model predictions, however, both the RRL profiles and the mm-to-cm SED can also be reproduced with slightly different values of the power-law index, in particular, with the steeper gradient, ner-2.4 , proposed by Tafoya et al. (2013). Finally, since we are able to explain our single-dish data suitably with a latitude-independent density law, we have not attempted any models with equatorially enhanced density profiles. We do not conclude that there is no equatorial overdensity present in the ionized wind, but increasing the model parameter space is undesirable with the data available since it would not provide a more accurate characterization of the density structure of the inner wind regions under study.

We chose to model the wind kinematics using an asymptotic expansion velocity law similar to what is commonly used to model the dust/wind acceleration region of mass-losing stars (Table 4). This law reproduces the broader profiles of the ~3 mm Hnα lines compared to the ~1 mm transitions. We emphasize that, in the absence of spatially resolved maps of the emission of mm-RRLs, it is not possible to accurately determine the exact velocity field in the ionized core. From the observed profiles, and taking into account that the 1 mm and 3 mm lines arise in neighboring layers, probably only 10 au apart from each other, we deduce that there must be a steep rise of the expansion velocity at these small spatial scales. In particular, the data are consistent with a minimum expansion velocity of Vexp~ 7–12 km s-1 at the innermost regions of the wind (probably near 40–60 au) and Vexp~ 20–25 km s-1 near the outer boundary of the modeled region (~70–90 au). The observed profiles can also be reproduced with different velocity fields (e.g.,  with Vexpr or Vexpr2) as long as the full range of expansion velocities observed (10–20 km s-1) is preserved. The expansive kinematics inferred imply a kinematical age for the inner wind layers of only tk ~ 15–20 yr. Our model also reproduces very well the broader profiles of the β-transitions, which confirms that, as expected from theory, these lines undergo pressure broadening.

The radius of the central cavity required in the model adopted is rather uncertain and, if it exists, it can only be well determined with spatially resolved observations at mm-wavelengths. Although the model presented in Table 4 provides the best fit to the data, larger values of Rin are also acceptable (up to Rin~ 70 au). The larger the size of the cavity, the larger the mass-loss rate that is required to keep the mm-continuum emission at the observed flux level. For example, if Rin~ 70 au, a mass-loss rate of ~1.1 × 10-5M yr-1 is necessary. In the case of a bigger cavity the velocity gradient across the 1 mm- and 3 mm-emitting layers should also be larger.

Our final model also reproduces the frequency-dependence of the major and minor axes sizes of the cm-continuum deduced by Tafoya et al. (2013) based on their detailed analysis of VLA 1.4–23 GHz continuum maps observed in the period 1982–2007. In principle, a good fit to the cm-continuum maps observed in the past should not be a strict requirement for our models because our mm-wavelength observations probe the optically thin (110 au) part of the ionized wind that was ejected in the last ~15–20 yr. But the cm-wavelength observations trace the outer regions, ~150–630 au, of an older wind that was ejected at earlier times during the last century. We chose to roughly reproduce simultaneously ancillary cm-continuum data to avoid extremely abrupt changes in the inner and outer wind properties.

5.2. M 2-9

As explained in Sect. 2.2, the free-free continuum emission of M 2-9 at millimeter wavelengths is expected to be dominated by the compact ionized wind at the core of the nebula, with a positive spectral index, α ~ 0.7. The mm-RRLs observed are also most likely produced in this ionized wind. The spectral distribution of the continuum, with no apparent changes in the spectral slope over ~10–230 GHz range, suggests that the free-free emission is optically thick at mm-wavelengths. In fact, this emission probably becomes optically thin for frequencies 230 GHz alone since no transition to a flat continuum spectrum is observed at lower frequencies. However, without angularly resolved maps of the continuum at millimeter wavelengths, a separate, accurate determination of the spectral index of the core and extended lobes is not possible.

As in the case of CRL 618, we modeled the compact ionized core of M 2-9 as an isothermal cylindrical wind that expands radially, which is the simplest structure that is consistent with the elongation of the radio core observed at cm-wavelengths. The input parameters of our best-fit model are given in Table 4 and the predicted RRL profiles and free-free continuum emission from the compact ionized wind are shown in Figs. 4 and 7.

The spectral index of the compact ionized wind continuum is well reproduced with an electron density gradient with a power index of αn = 2.1 and electron density of ne~ 6 × 105  cm-3at 45 au, resulting in a current mass-loss rate of pAGB 3.9 × 10-7M yr-1. We also derive an upper limit for the inner radius of the H ii region of Rin 2 au, which is needed to explain the turnover frequency at νt> 230 GHz, as previously stated (Fig. 1). The bulk of the RRL and continuum emission at millimeter wavelengths arises within the ~5–50 au inner layers of the ionized wind.

We find that an LTE electron temperature of Te~ 7500 K reproduces well the velocity-integrated line to free-free continuum flux ratios (LTCR) observed for all α-transitions. This implies a value for the sound speed of cs ~ 11 km s-1. The isothermal assumption reproduces the LTCR frequency dependence observed, LTCR ν1.1, which is expected under LTE conditions (see, e.g., Mezger & Hoglund 1967).

The trend of decreasing line width with increasing frequency guessed from the Hα profiles of M 2-9 (Fig. 6) suggests radial expansion with a positive velocity gradient, as for CRL 618. As a first approximation we assumed constant radial expansion because in this case the trend is only tentative and the inferred velocity difference between the inner and outer layers of the ionized wind (where the H30α and H39α emission is produced, respectively) is small (~4 km s-1). Under this assumption, the RRL profiles are reasonably well fitted assuming Vexp 22 km s-1 and a systemic velocity of VLSR = +75 km s-1. Unlike for CRL 618, the optical depth of the mm-RRLs is very small (1) over the bulk of the emitting region of M 2-9 and, therefore, the line centroids are predicted to be well centered on the adopted systemic velocity.

Finally, our simple cylindrical model is able to reproduce the RRL profiles and free-free continuum at mm-wavelengths produced at the compact core of M 2-9, but our model (moderately) underestimates the continuum flux density from the core at wavelengths 3 cm. This is because a non-negligible fraction of the cm-continuum flux arises in the outer (>160 au) layers of the ionized wind, which are being truncated, and thus excluded, in the cylindrical model along the equatorial direction. We used a cylindrical geometry based on the elongated shape of the cm-continuum brightness distribution observed by Kwok et al. (1985) and Lim & Kwok (2003). However, the model presented here mainly aims to constrain the properties of the inner regions of the ionized wind (~5–50 au), where the free-free (RRL and continuum) emission at mm-wavelengths reported in this work arises and which are satisfactorily reproduced. Detailed modeling of the properties of the cm-continuum emission produced in the outermost regions of the ionized wind is beyond the scope of this paper.

We stress again that, as for the other targets, the simple model presented here is most likely not unique and that more complex geometries, kinematics, density, and temperature profiles (e.g., latitude-dependent) cannot be ruled out. In the particular case of M 2-9, we checked that a double-cone geometry, similar to that used for MWC 922 (see next section), could also reproduce the mm-wavelength data. In this case, to fit the widths of the mm-RRLs, the expansion velocity of the wind would have to be larger and the mass-loss rate moderately smaller (for example, ~ 1.4 × 10-7M yr-1 for a semi-opening angle of ~20° and Vexp = 25 km s-1). High-angular resolution maps of the mm-continuum emission are needed to accurately determine the spatio-kinematic structure and physical conditions of the inner ionized wind layers of our targets.

5.3. MWC 922

The ionized core of MWC 922 has been represented by a double-cone geometry like that used for modeling the ultra-compact H ii region MWC 349 A in Martín-Pintado et al. (2011) and, more recently, Báez-Rubio et al. (2013). This is based on the resemblance of the maser and non-maser RRL profiles of both sources and it is also consistent with the X-like morphology of the optical/NIR nebula surrounding MWC 922 (Sect. 2.3). The ionized emission is thus assumed to arise in two opposing conical structures. Each of these structures is formed by two nested cones co-axial with, and inscribed in, an extended neutral rotating disk. A sketch of this geometry is given by Báez-Rubio et al. (2013) in their Fig. 3. The outermost surface of the bicone (in contact with the neutral disk) has a semi-opening angle θa measured from its revolution axis to the midplane. The outer (shell-like) bicone has an angular width θd. The regions that are inside of (surrounded by) the innermost surface of the outer bicone, at colatitude θ<θaθd, form the inner bicone. Regions of colatitude θ>θa reside in the neutral disk, which is adopted to be edge-on.

The inner and outer ionized bicones can have different kinematics. The outer bicone is assumed to rotate in a Keplerian fashion, presumably sharing the same kinematics with the neutral disk in which this bicone is hypothetically inscribed. The outer bicone would then simply represent the thin conical surface of the rotating disk that is illuminated and photoionized by the central source. If the density above the disk is lower than that in the equatorial plane, the central star irradiates/ionizes the upper thin surfaces of the disk. We note that the disk material at low latitudes could remain largely neutral if the ionizing radiation is additionally blocked along the midplane by dust, for example, at the inner edge of a compact circumstellar disk. The inner bicone is modeled as a rotating and radially expanding wind; the same Keplerian rotation field adopted in the disk is vectorially added to a radial expansion field. The inner bicone could represent the stellar wind and/or gas photoevaporating from the ionized surface of the rotating disk. We refer to the inner and the surrounding outer bicone as the “ionized outflow” and the “ionized disk component”, respectively.

The gas in the rotating disk outside a certain radius rd is considered to be part of the outflow, that is, it is expanding and rotating. The value of rd is always chosen to be smaller than or equal to the disk gravitational radius, which is defined as the point where the Keplerian velocity (~escape velocity) and the sound speed are equal; the disk region outside the gravitational radius is no longer gravitationally bound to the central mass and it flows away freely (e.g., Hollenbach et al. 1994; Owen et al. 2012).

The electron density is assumed to vary throughout the disk and outflow as a function of the radial distance, following a power-law distribution and also as a function of the latitude, decreasing exponentially from the midplane toward the poles. In principle, our model enables the electron temperature to be different in the disk and outflow (but see below).

In Table 4 we give the input parameters of a model that successfully reproduces the mm-RRL profiles and the free-free continuum emission of MWC 922 (Figs. 5 and 7). This model predicts the observed transition from single-peak profiles at ~3 mm to double-peak profiles at ~1 mm, including the velocity peak-to-peak separation of the latter, and the steep increase of the Hnα-line FWHM and intensity as a function of frequency; the latter is indicative of non-LTE conditions and, in particular, of H30α maser emission in this object. It also recreates the spectral continuum emission distribution from the mm- to the cm-wavelength range and the angular diameter (~0.̋2 = 340 au at 1700 pc) of the 3.6 cm-continuum emitting region (Sect. 2.3).

Before discussing this model in detail, we stress that in the case of MWC 922, the geometry, kinematics, and electron temperature deduced are particularly uncertain. This is because of the extremely high sensitivity of the maser line intensity not only to small changes in the local physical conditions but also to the uncertain bn departure coefficients (Sect. 5). Our final model presented in Table 4 uses the coefficients by Storey & Hummer (1995), however, we also explored a large range of model parameters using the bn coefficients Walmsley (1990). Since the bn coefficients computed by these two groups differ significantly, the physical model needed to reproduce the observed line intensities are different when using one set of bn or the other. Irrespective of the bn set adopted, the biggest challenge is always to reproduce simultaneously the maser and non-maser RRL profiles observed.

Most of the mm-wavelength free-free emission observed in MWC 922 arises in a rather compact region within r ≲ 150 au. In particular, 90% of the RRL and continuum emission at 1 and 3 mm arise in the inner (~10–80 au) and outer (~20–150 au) parts of the central ionized core, respectively. The full extent of the H ii region is probably larger as suggested by the 3.6 cm-continuum maps. The inner radius of the ionized core of MWC 922 is poorly determined from these data since the deepest ionized layers at 10–15 au are very optically thick at mm-wavelengths and thus are not probed by our mm-RRL profiles. We chose a nominal inner radius in our model of Rin = 0.05 au, which is equivalent to Rin~ 1 R, adopting the values of Teff and luminosity given in Table 1. As for M 2-9, we derive an upper limit to the radius of an hypothetical central cavity, which could be devoid of ionized gas, of Rin 35 au; this is the largest cavity radius that would still reproduce the uniform spectral index of the free-free continuum as observed, which indicates optically thick emission from cm wavelengths up to ~230 GHz (at least).

The electron temperature distribution in the ionized core can only be poorly constrained from single-dish data. Initially, the temperature Te has been enabled to be different in the disk and in the outflow, however, for the geometry and opening angle adopted in our final model, the temperature in the disk and outflow have to be comparable to reproduce simultaneously the observed maser and non-maser profiles, Te~ 6500–7000 K.

The electron density rapidly declines from the ionized surface/layer of the rotating disk toward the poles. Along θ = θa + θd/ 2, i.e., within the thin ionized layer of the disk, the density ranges from ne 107 to 105 cm-3, following a r-2.4 power law, across the r ~ 15–150 au layers that emit the bulk of the continuum and emission of RRLs at mm-wavelengths. The density falls off by a factor of ~20–30 across the full width of the ionized disk, θd, and continues to decrease rather steeply across the ionized outflow, dropping by 5 dex at the revolution axis.

Table 5

Results from modeling of the free-free continuum for mm-RRL non-detections (Sect. 5.4).

As discussed in Sect. 4, the gradual increase of the width of the Hnα profiles as a function of frequency is consistent with Keplerian rotation. To explain the velocity peak-to-peak separation observed in the H30α maser line, we considered a central mass of ~8 M, which is equivalent to the rotation velocity field given in Table 4; we find that a range of ~5–10 M gives also acceptable predictions. The maser spike emission arises mainly in the Keplerian disk, at radial distances ~35–55 au, because this region has the geometry (a high coherence length) and the densities, of ne 106–107  cm-3, required to produce strong maser amplification of the H30α line. The double-peaked profile of the H31α line, which is also indicative of stimulated non-LTE amplification, is well reproduced by our model. However, the intensity of the H31α line relative to the rest of the RRLs detected is always notably lower than observed for all models regardless of the physical conditions and kinematics adopted. For this reason, and given the extremely high uncertainties of the bn coefficients, in our final model, we used the same departure b30-coefficient for the H31α and H30α transition, which results in a much better model-data agreement.

In the ionized outflow, we adopted a constant expansion velocity of Vexp~ 5 km s-1, which is consistent with the profile of the non-maser H39α and H41α lines, and, in particular, with the lack of emission wings. Given the dimensions of the outermost regions of the ionized wind probed by our observations, ~140 au, we derive a kinematical age for the outflow of tk ~ 150 yr. This is a lower limit, since the outer boundaries of the ionized core of MWC 922 are not observationally constrained and could reach larger distances. The innermost layers of the outflow, for example, at Rin = 15 au, would have been ejected only 10 yr ago. The age of the rotating disk cannot be estimated from these data but it, in principle, it could be a much older long-lived structure.

Finally, we already mentioned the large uncertainties of the model input parameters in the case of MWC 922. From our multiple model runs, we conclude that it is also possible to reproduce the data used in this work adopting a larger opening angle for the ionized disk+outflow structure. For example, using θa ~ 53°, we are able to find a set of physical parameters that satisfactorily explain the observations. One important difference with respect to the model in Table 4 is that, for a larger opening angle, the temperature in the ionized disk must be notably larger than in the outflow. In our opinion, the presence of a hot ionized thin disk sandwiched between two cooler structures (namely, the outflow and the neutral disk) is improbable5.

As for the other targets, we estimate an equivalent isotropic mass-loss rate for the outflow. For a constant expansion velocity of Vexp = 5 km s-1, and considering that we are adopting a non-isotropic (bipolar) outflow with a latitude-dependent density, we find ~ 1.8 × 10-5×(1cos(25°))~2 × 10-6M yr-1 for the value of θa = 25° used in our model. As expected, the equivalent isotropic mass-loss rate does not notably change when using a different opening-angle, since it is mainly constrained by the free-free continuum flux. In particular, for θa ~ 53°, we find ~ 5.2 × 10-6× (1cos(53°)) ~ 2.0 × 10-6M yr-1. However, the nature of the outflow in MWC 922 is unclear (stellar wind and/or photoevaporating gas from the ionized disk) and, therefore, the meaning/interpretation of the rate derived is not straightforward.

5.4. Remaining targets: mm-RRL non-detections

We modeled the free-free continuum emission of the sample targets with no mm-RRL detections. This mainly aims at obtaining an estimate of the mass-loss rate for these objects as well. In the simplified model used in this case, we assume that the free-free emission arises in a spherical isothermal wind expanding at constant velocity. In particular, we adopt Vexp = 15 km s-1, Te = 10 000 K, and an electron density power-law of the type ne(r) rαn. In this case, the free-free continuum is expected to vary with frequency as and, therefore, the αn parameter can be constrained from the spectral index of the free-free continuum (see, e.g., Rodríguez et al. 2009).

The main input parameters and synthetic free-free continuum emission of the best models are given in Table 5 and in Fig. 7. In all cases, we checked that MORELI predicts mm-RRLs with intensities below our detection limits. As we can see, most of the free-free continuum emission at mm-wavelengths arises at the very inner layers of the wind around the star, typically within Rout(99%)  130 au and, in some cases such as M 2-56, in regions as compact as Rout(99%)  30 au. The pPNe/yPNe with the most extended H ii region is He 3-1475 (with Rout(99%) ~ 800 au). In all but one case (He 3-1475), these inner layers represent very recent mass ejections that happened less than ~15 × yr ago.

We find spectral indexes roughly consistent with inverse-square density laws except for M 2-56 and IRAS 20462+3416. In these two cases, a steeper radial density fall, as ner[ − 3.1: − 3.5 ], is inferred, which may suggest a gradual increase of the mass-loss rate in recent times.

We derive characteristic mass-loss rates in the range pAGB~ [0.3–4] × 10-6M yr-1. We do not observe a clear trend between pAGB and the bolometric luminosity; there are rather luminous objects, namely, M 1-91 and M 2-56, with relatively low mass-loss rates of pAGB~ [2–3] × 10-7M yr-1, as well as dim sources, namely, IRAS 20462+3416, with relatively high mass-loss rates of pAGB~ [1–2] × 10-6M yr-1.

6. Mass-loss history and evolutionary status

The emission from the mm-RRLs reported and modeled here allowed us to characterize the relatively dense, ne~ 106−8 cm-3, inner (10–100 au) ionized regions of three of our program objects: CRL 618, M 2-9, and MWC 922. Considering their expansion velocities and small extent (Table 4 and Sect. 5), these inner regions probe very recent mass ejections events that happened less than ~15–20 yr. This is probably also true for our sources with no mm-RRL detections, since there is no reason to expect remarkably different bulk gas motion speeds in these regions.

The mass-loss history in more ancient epochs (back to few thousand years ago) is relatively well known for CRL 618 and M 2-9 from previous works but, unfortunately, not for MWC 922 (Sect. 1). In the next subsections, we discuss our findings in the context of these former mass-loss episodes, the nature and evolutionary status of these objects.

6.1. CRL 618

In the case of CRL 618, we found a significant difference between the H30α and H41α line profiles (widths and peak intensity) and those observed in 1987 by Martín-Pintado et al. (1988). These changes are not totally surprising since our observations do not trace the same wind material as observed in the 1980s. This is because the expanding wind layers probed by the mm-RRLs observed ~30 yr ago by Martín-Pintado et al. (1988) have traveled ~125 au outward during this time (at Vexp 20 km s-1) starting from their position in 1987 (at Rin~ 70 au). Therefore, they are now beyond ~195 au and outside the mm-RRL emitting regions. The changes in the line profiles observed are then indicative of wind variability, i.e., evolution with time of the wind properties6.

The kinematic age of the layers of the post-AGB wind of CRL 618 observed by us is ~10–20 yr (Table 4). The mass-loss rate deduced from our mm-RRL observations, pAGB~ 8.4 × 10-6M yr-1, is that of a post-AGB wind that was ejected in the ~10–20 yr previous to 2015. Similarly, the mass-loss rate deduced by Martín-Pintado et al. (1988), pAGB~ 7.6 × 10-6M yr-1, probably corresponds to a wind ejected in the ~10–20 yr previous to 1987. The slightly larger mass-loss rate estimated by us in more recent epochs is then in agreement with the trend deduced by Tafoya et al. (2013), who proposed that the mass-loss rate had increased (from pAGB~ 4 × 10-6 to ~6 × 10-6M yr-1) in the ~130 yr previous to 1998. In this context, our result is consistent with a mass-loss rate that has continued to rise and, indeed, at nearly the same rate contemplated by Tafoya et al. (2013).

Martín-Pintado et al. (1988) considered a central cavity with radius Rin~ 70 au in their LTE model of the ionized core of CRL 618. As explained in Sect. 5, the size of this cavity, which could result from a short period of diminished mass-loss rate, is rather uncertain because of the degeneracy of the turnover frequency with Rin and pAGB. However, we would like to mention that the presence of a cavity with a radius of Rin~ 70 au in 1987 is not inconsistent with our data and model. This is because, as explained in the first paragraph of this section, the cavity (as the rest of the layers of the ionized wind) would have moved forward by ~125 au in the last ~30 yr. Therefore, if the cavity hypothesized by Martín-Pintado et al. (1988) was real and its size in 1987 roughly correct, 30 yr later, it should constitute a low-density shell at ~195–265 au, i.e., beyond the ionized wind layers traced by our mm-wavelength observations. The fact that a cavity is also needed to explain the observations in 2015 could indicate another (periodically recurrent?) recent short intermission of the mass loss. In any case, we caution that the cavity is needed in the context of the simple model used here (Table 4). The cavity may not be a requirement if one adopts a more sophisticated model with a sufficiently complex geometry, density, temperature, and velocity stratification (for example, an inhomogeous wind with dense clumps).

Of course, the comparison of the mass-loss rates (and other parameters) derived from different works has to be carried out with extreme caution since these are affected by hardly quantifiable uncertainties that are inherent to the different analysis methods and assumptions adopted. Therefore, our conclusions about the overall increase of pAGB and intermittent reduced mass-loss rate intervals should be considered tentative.

The mass-loss history of CRL 618 in earlier epochs is written in its nebular architecture (Sect. 2.1). From detailed observations and modeling of its molecular envelope, two major large-scale, mass-loss episodes in the form of slow winds (Vexp 17 km s-1) have been identified. The first took place 2500 yr ago at a rate of  10-5M yr-1 and generated a tenuous extended halo. The second mass-loss event started ~400 yr ago at very high rate of ~ few × 10-4–10-3M yr-1 resulting in the formation of a dense central core (Sánchez Contreras et al. 2004; Lee et al. 2013a; Soria-Ruiz et al. 2013). More recently, ~200 yr ago, the interaction between fast bullet-like ejections and the pre-existing AGB circumstellar envelope probably shaped the optical lobes (Balick et al. 2013; Velázquez et al. 2014). The compact, fast (100 km s-1) molecular outflows of CRL 618, with a kinematical age of ~30–80 yr, could have been partially shaped by the mentioned bullet-like ejections or by a distinct, more recent (80 yr) collimated fast wind (Sánchez Contreras et al. 2004; Lee et al. 2013b; Huang et al. 2016).

The mass-loss rate at which the post-AGB wind was ejected over the last ~150 yr in CRL 618 is pAGB~ [4–8] × 10-6M yr-1. Considering the whole range of values deduced including our work (see above), this mass-loss rate is then only moderately lower than the mass-loss rates of the slow AGB winds that led to the molecular extended halo and dense core. Given their low expansion velocity and high mass-loss rates, the halo and core were most likely ejected when the central star was ascending the AGB toward its tip and immediately before leaving this phase, respectively; this is when the most intense winds are believed to be blown out (e.g., Bloecker 1995; Vassiliadis & Wood 1993). Therefore, in the case of CRL 618, the drop of the mass-loss rate, which is believed to mark the end of the AGB phase (after the so-called “superwind”) and the beginning of the AGB-to-PN evolution, is rather modest, in particular, only about one order of magnitude lower than that of the core-wind, and comparable to that of the halo-wind.

6.2. M 2-9

Our mm-RRL observations trace the inner layers of a dense ionized wind at the core of M 2-9 that has been ejected at an average rate of pAGB~ 3.5 × 10-7M yr-1 and expansion velocity Vexp~ 22 km s-1 over the last 15 yr (Table 4). The mass-loss rate estimated by us is consistent (within uncertainties) with that obtained by Kwok et al. (1985) for the outer layers of the ionized core from their analysis of the cm-continuum emission maps (1.4–22 GHz, see also Sect. 2.2). These authors find pAGB = 3.3 × 10-5(Vexp/1600 km s-1)(d/kpc)3/2M yr-1, which is equivalent to pAGB = 2.4 × 10-7M yr-1 after scaling it to the values of the distance and the expansion velocity in Table 4.

In 1982, the outer boundary of the ionized core of M 2-9 along the nebula axis was probably near 0.̋3 (Rout~ 200 au); beyond this radius most of the material remains neutral, as deduced from the cm-continuum maps presented by Kwok et al. (1985). Adopting the same expansion velocity measured by us, Vexp~ 22 km s-1, the wind layers traced by the cm-continuum observed by these authors must have been ejected in the 40–50 yr previous to 1982 (the observing date), that is, sometime after 1932. Therefore, the current wind observed by us could have been ongoing with a constant or, at least, not markedly variable mass-loss rate of ~3 × 10-7M yr-1 since it began at least ~50 yr ago.

Prior to this, two other major mass-loss episodes have been identified in this object leading to two ring-shaped, eccentric molecular structures at the nebular equator mapped in CO emission by CC12; see also Sect. 2.2. Recent 0.̋18-resolution CO emission maps obtained with ALMA confirm the spatio-kinematics of the two rings (Castro-Carrizo et al. 2017). These authors conclude that these two expanding rings were formed during two short mass-loss episodes (of duration ~40 yr) produced when the mass-losing star was at different positions in the orbit of the central binary. This scenario is also supported by the different systemic velocities found for the two rings, VLSR = 80.0 and 80.6 km s-1 (±0.1 km s-1) for the outer and inner ring, respectively, which is readily explained if the mass-losing star changed its velocity due to its orbital motion within a binary (or multiple) system. CC12 proposed that the outer ring, with Rout~ 2340 au and Vexp = 7.8 ± 0.1 km s-1, was ejected at a rate of ~9 × 10-5M yr-1 about 1400 yr ago, whereas the inner ring, with Rout~ 780 au and Vexp = 3.9 ± 0.1 km s-1, would have occurred ~900 yr ago at a slightly smaller rate. It is unknown whether the mass-losing central star of M 2-9 was in the late AGB or early post-AGB phase when the CO-rings were ejected, but the relatively high mass-loss rates (and low expansion velocities) observed are most consistent with a late-AGB status. The mass-loss rate of the younger ionized core-wind studied here (Table 4) is about two orders of magnitude smaller than that of the preceding mass outburst that led to the inner CO ring, and thus most consistent with a post-AGB evolutionary stage at present.

The central binary. We found a significant difference between the centroids of the mm-RRLs, VLSR = 75 ± 2 km s-1 (Sect. 4), and those of the CO lines. This VLSR difference is consistent with the scenario proposed by CC12 in which the mass-losing star is moving along its orbit; the centroid of the mm-RRLs would then correspond to the line-of-sight velocity of the star averaged during the last 10–20 yr (when the inner layers of the ionized wind were ejected). The kinematical age of the inner wind layers probed by our mm-RRL observations only represents a small fraction (1/5) of the orbital period. We can obtain a lower limit to the orbital velocity of the mass-losing (primary) star of  (80.6 – 75)/(2 cos(17°))  2.9 km s-1, which would imply an average distance to the center of mass of a1 8.7 au for an orbital period of P ~ 90 yr. With this limited knowledge of the orbital elements for only one of the stars (m1, the mass-losing star), we cannot derive the mass (m2) and orbit of the secondary with certainty, but we can obtain some constraints by applying the well-known relationship between the binary mass function, f(m1,m2), and the observables a1 and P. Assuming a circular orbit, the mass function, which is derived from the third Kepler’s law, is given by (2)where m2 is the mass of the companion, q is the dimensionless secondary-to-primary mass ratio q = m2/m1; units are given in M, au, and years for the masses (including f(m1,m2)), semi-major axis, and period, respectively.

Using the lower limit to a1 deduced above, we can also obtain a lower limit to the left hand of this equation, f(m1,m2) 0.08 M and, thus, a lower limit to the mass of the companion relative to the primary. The mass of the primary star (m1), which is the stellar remnant after most of the envelope was ejected in the earlier red giant and AGB phases, is expected to be within ~0.5–1 M range (Schönberner 1983; Bloecker 1995; Vassiliadis & Wood 1994). In particular, the low luminosity of M 2-9 (~700–3000 L at d = 650–1700 pc, Table 1) is most consistent with a low-mass progenitor (with an initial mass of 1 M) and, thus, with a low-mass stellar post-AGB core of about 0.5 M (see Sect. 7 and Fig. 8). Adopting m1 = 0.5M, the limit to f(m1,m2) implies that the companion mass must be m2 0.4 M. The uncertainty in the systemic velocity of the RRLs, Vsys = 75 ± 2 km s-1, translate into a relative error of 50% in the mass of the companion, which could then be as low as m2 ~ 0.2M as derived from these data.

We now analyze the possible evolutionary status of the secondary for different values of the mass ratio. A main sequence companion would be consistent with a mass in the range m2 0.2 and 1 M. Larger masses of the secondary (assuming a dwarf star) would be in conflict with its slower evolution relative to the low-mass (1 M) primary. In this case, the relative orbital separation would be a = a1 + a2~ 18–23 au. If the companion is a white dwarf (WD), the upper limit to its mass is m2 ~ 1.4M, which would imply a slightly larger orbital separation of a~ 25 au. If this is the case, then M 2-9 could have been a classic symbiotic system (AGB+WD) until recently, when the primary mass-losing star would have left the AGB; at present, the whole system could be dying a second death as as post-symbiotic (“post-AGB”+WD) system.

Values of q larger than those discussed above are, in principle, possible based on our observables and interpretation (a1 8.7 au and P ~ 90 yr), but it would imply that the companion is a dark, compact stellar remnant such as a neutron star (NS, with a typical mass in the range m2 ~ 1.1–3.2 M) or a black hole (BH, if m2> 3M). If this is the case, the orbital separation of the system would be larger than for the dwarf or WD companion scenario (for example, a ~ 44 au for m2 ~ 10M). If the mass-losing star is accompanied by such a NS/BH then, given the relatively large orbital separation, it would be most likely in a “silent” mode; that is, this star would not be undergoing significant mass-transfer that, if present, should be accompanied by high-energetic phenomena. This is indeed consistent with the non-detection of X-ray emission toward M 2-9 (Ruiz 2014).

The mass of the companion deduced here (m2 ≳ 0.2M) is somewhat larger than that estimated by CC12, who proposed a low-mass companion of m2 ≲ 0.1–0.2 M orbiting around a m1 ~ 1M mass-losing star when the CO rings were ejected. The low secondary-to-primary mass ratio proposed by these authors results from the low orbital velocity of the primary used in their calculations, ~ 1 km s-1, which is smaller than the value inferred from the separation between the centroids of the mm-RRLs and the CO line ( 2.9±1 km s-1). We note that CC12 interpret the difference between the CO line centroids (~0.6 km s-1) and velocity gradients observed in the two rings as the result of two short mass-ejection events that took place when the mass-losing star was at two specific positions in a circular orbit. In particular, they proposed that the large CO ring was ejected when the system was near conjunction (with the stars moving perpendicularly to the line of sight) and the small ring was blown out ~500 yr later, when the stars were at quadrature (and, thus moving in a direction very close to the line of sight); see red and green marks in their Fig. 4.

The low radial velocity shifts deduced from the CO emission from the two rings are not inconsistent with a moderately larger orbital velocity, as deduced here, if the rings were expelled at positions slightly different from those adopted by CC12 or, most importantly, if the orbit is not circular. For elliptical orbits, the radial velocity curve exhibit a marked skew symmetry with different maximum and minimum absolute radial velocities. In fact, for eccentric orbits the observed radial velocity can be much smaller than the semi-amplitude of the curve during most of the orbital phase. If this is the case, the relatively low radial velocities derived from the spatio-kinematics of the CO rings would not be incompatible with a larger radial velocity measured at a third epoch from our RRLs, probably nearer to periastron.

Finally, one can also consider the possibility that the ionized gas traced by our mm-RRLs is not centered around the primary mass-losing star but around the companion. In this case, the difference between the centroids of the RRLs and CO transitions would represent the orbital velocity of the less massive secondary (moving faster than the primary about the mass center). If the companion is a very low-mass star, m2 0.2, an external source of ionization, other than the companion itself, would be needed. This is because a main sequence star of this mass or lower would have a spectral type later than ~M 4, i.e., Teff 3200 K, and would not be able to ionize the gas in its vicinity. High-speed shocks in a putative accretion disk around and/or within fast winds launched by the secondary (presumably, after wind mass transfer from the primary) may constitute an alternative ionization agent. If this is correct, at present, the companion must be undergoing active accretion from the wind of the primary (given the short lifetime, 10–20 yr, of the ionized core) and broader RRL profiles from the hot shocked emitting gas would be expected. In the case of high-speed shocks, one would also expect partial ionization of the slow dense wind around the mass-losing star by the UV radiation generated by the shocks.

With the observables available to date and, in particular, in the absence of direct empirical determination of its orbital or stellar parameters, an accurate characterization of the system at the core of M 2-9 is not possible.

6.3. MWC 922

The B [e] star MWC 922 and its surrounding nebulosity remain poorly characterized to date and, therefore, the long-term history of the mass loss and nebular shaping cannot be reconstructed as for CRL 618 and M 2-9. Our mm-RRL observations are consistent with the presence of an ionized outflow and a Keplerian rotating disk within a radius of ~150 au from the center. The origin of these two structures is unknown and cannot be determined from our data.

The rotating kinematics of the disk, if it is Keplerian as adopted in our model, provides an estimate of the mass of the central object. This parameter is poorly constrained, mainly owing to the large uncertainties of the departure coefficients bn on which the maser line profiles critically depend, but values of ~5–10 M are consistent with our observations. This large mass could represent the total mass of a binary or, more generally, multiple system.

The location of MWC 922 in the HR diagram is uncertain, mainly due to the unknown distance. The temperature of the central star is also not well established, although it probably lies within the range Teff~20 000–30 000 K (Sect. 2.3). As discussed by Tuthill & Lloyd (2007), the presence of gas and dust around MWC 922 suggests either a pre-main sequence star still partially embedded in its natal cloud or a post-main sequence star surrounded by matter ejected in its late evolutionary stages.

The range of Teff and luminosities given in Table 1 are consistent, in principle, with MWC 922 being a post-AGB object evolved from a massive progenitor (Fig. 8). Considering the post-AGB evolutionary tracks by Bloecker (1995), which include stars with initial masses of up to 7 M, MWC 922 could be the evolutionary product of a ~5 M main sequence star (adopting d = 1.7 kpc and log (L/L) = 4.25) or even more massive, >7 M (if d = 3.0 kpc and log (L/L) = 4.77). The mass of the remnant post-AGB core of MWC 922 would be around 1 M at present. Since the total mass at the center of the rotating disk of MWC 922 deduced from our analysis is significantly larger, ~5–10 M, either a massive companion (m2 ~ 4–9 M) or multiple lower mass companions must exist. As explained earlier, the mass of the central object in uncertain. If we take the lowest value to be consistent with the observations, which is 5 M, then the secondary could be a ~4 M star still on the main sequence [with log (L/L) ~ 2.2 and spectral type ~B 6-B7], while the more massive primary (with an initial mass 5 M) would have evolved off the main sequence already (in less than 10–100 Myr).

If this scenario is correct and MWC 922 is an evolved/post-main sequence star, its age would still be sufficiently low to be associated with other young less massive stars born in the same natal cloud, which could also remain not fully dispersed. This is not inconsistent with the strong ISM contamination of the CO spectra toward this target (see Appendix A), which could be completely masking out the emission from the presumptive mass-ejecta from the star. In any case, given its location close to the Galactic plane, strong ISM contamination is expected toward MWC 922.

If the total mass of the system at the core of MWC 922 is 8–10 M, and m1 ~ 1M, then the companion (or companions) should account for the remaining stellar mass of 7–9 M. In this case, and if there is only one companion, a pre-main sequence nature would be probably favored. This is the case unless, of course, the initial mass of MWC 922 was even larger, i.e., 10 M, and it is now in a pre-oxygen-neon WD or pre-supernova stage; this stage is expected to be extremely short, which makes this scenario less probable.

7. Post-AGB mass-loss rates and post-AGB evolution

The post-AGB phase is, beyond question, one of the least understood phases of the evolution of low-to-intermediate mass stars. Mass loss is the dominant mechanism impelling these stars across the HR diagram on their way to the CSPN phase, before fading along the white dwarf cooling sequence. The total AGB-to-PN transition time depends on the heating rate of the central star, which is, indeed, critically dictated by the mass-loss rate, especially in the early post-AGB stages: the more intense the post-AGB wind, the faster the evolution to the final (hottest) point of the CSPN phase. In spite of its decisive role in stellar evolution, empirical data on post-AGB mass-loss rates are notably lacking and, therefore, evolutionary models are bound to adopt completely unconstrained mass-loss prescriptions, for example, interpolating between those theoretically or semi-empirically determined for the preceding AGB and the following CSPN phases (e.g., Schönberner 1983; Bloecker 1995; Vassiliadis & Wood 1993; Miller Bertolami 2016).

thumbnail Fig. 8

Adapted from Figs. 12 and 5 of Bloecker (1995). Top: Post-AGB evolutionary tracks with observational data for RRLs and non-RRL detections (filled and empty symbols, respectively – Table 1). For MWC 922 and M 2-9, the uncertainty of the luminosity due to the uncertainty of the distance is indicated by the arrows. The stellar mass of the post-AGB remnant core of each track is indicated to its right. Vertical dotted lines delimit the ionization onset region (Teff~ 20 000–40 000 K) where our targets lie. Bottom: Post-AGB mass-loss vs. effective temperatures adopted by post-AGB models of 0.605, 0.696, 0.836, and 0.894 M (Bloecker 1995) and of 0.546 M and 0.565 M (from Schönberner 1983). As in the top panel, the symbols indicate the location of our targets using the values of pAGB deduced in this work (Tables 4 and 5; vertical error bars correspond to an uncertainty factor of 2). The thick dash- and dot-dashed lines indicate the mass-loss rate typically assumed by evolutionary models for our three objects with RRL detections.

Open with DEXTER

Commonly, the post-AGB mass-loss is described by a sudden abrupt drop immediately after the departure from the AGB phase (where very high rates of 10-3–10-4M yr-1 are applied) followed by a less steep, linear decline, and final flattening in the CSPN regime (see Fig. 8, bottom panel, and references above). The mass-loss rates derived for CRL 618, MWC 922, and M 2-9 in this work (Table 4) are well above the typical values adopted in post-AGB evolutionary models, which vary between pAGB~ 1 × 10-9 and 2 × 10-7M yr-1 for the range of Teff considered here.

The largest discrepancy, by three orders of magnitude, is found for the well-known pPN CRL 618. The location of this object in the HR diagram (Fig. 8, upper panel) is consistent with a remnant core of ~0.625 M, which is the descendant of a star with an initial mass of ~3 M on the main sequence, for which pAGB~ 1 × 10-8M yr-1 would be assumed by evolutionary models currently in use (Fig. 8, bottom panel).

The mass-loss rate deduced for M 2-9 is also much higher than the values adopted by models, although the progenitor mass in this case is less certain because the distance to this object is not well known. Adopting the largest distance proposed in previous works, d = 1.7 kpc, the position of M 2-9 in the HR diagram would be consistent with a remnant mass of 0.565 M, which is the descendant of a 1 M star, for which pAGB~ 10-9M yr-1 or less is normally assumed (Fig. 8). The lowest luminosity point in the HR diagram (log (L/L) ~ 2.8, adopting d = 650 pc) falls below the evolutionary tracks of the least massive stars (~0.8 M), which are expected to evolve off the main sequence in less than a Hubble time.

In the case of MWC 922, the mass-loss rate is also above the model prescriptions assuming that it is the evolutionary product of a 5 M star (with a 0.836 M remnant core). However, not only the post-AGB nature of this object is questionable (Sect. 6.3); the origin of the expanding outflow (+rotating disk), which is consistent with the observed mm-RRL line profiles, is not clear since it could represent a stellar wind or gas photoevaporating from the ionized surface of a circumstellar disk, or a mixture of both.

Although the mass-loss rates inferred for the five objects in our sample with no mm-RRL detections are uncertain (given that the spatio-kinematics of their ionized cores is unconstrained), for all of these objects we obtain relatively large values of pAGB 10-6–10-7M yr-1 (Table 5) that clearly exceed the rates adopted by models. For example, the most luminous pPNe/yPNe in our sample with no mm-RRL detections, namely, He3-1475, M 1-92, and M 2-56, would have similar ~2–3 M progenitors to CRL 618 for which models would use pAGB10-8M yr-1. IRAS 20462+3416 and M 2-56 have luminosities that are comparable to M 2-9 and, therefore, they also probably descend from a low-mass (1 M) progenitor. For such a low-mass progenitor, the post-AGB mass-loss rates used by models are 3 orders of magnitude below the rates inferred by us.

As mentioned in the beginning of this section, the mass-loss rates applied by models during the transition between the AGB and hot CSPN stage are expected to influence the HR crossing time dramatically. Some models, for example,  those by Schönberner (1983) for ~0.8–1 M stars, assume that most of the envelope mass is removed at very high rates of ~10-3–10-4M yr-1 in a short “superwind” phase that happens only in the vicinity of the AGB, i.e. for Teff~ 4000–6000 K. After that point, which is arbitrarily defined by most theorists as the beginning of the post-AGB phase, the adopted mass-loss rates fall well below the growth rate of the core due to hydrogen burning, and therefore, they are considered to influence only moderately the transition of the star to the CSPN region. Other evolutionary models, such as those by Bloecker (1995) for stars with initial masses between 1 and 7 M, however, indicate that the transition time from the AGB to the CSPN region depends strongly on the treatment of mass loss beyond the AGB phase, although rather low post-AGB mass-loss rates (inferred by radiation-driven wind theories) are still applied in these computations.

We have shown that the mass-loss rate in the post-AGB phase, when the star has reached Teff 20 000 K and the ionization has begun, can be orders of magnitude larger than the values assumed by evolutionary models. If this is found to be a common property of most post-AGB stars (which needs to be investigated using larger samples), much faster evolutionary speeds to the PN region may be expected.

A similar result was obtained by Trams et al. (1989) from the analysis of Hα P-Cygni profiles in four candidate post-AGB stars with Teff~ 6500–10 000 K, which yielded order-of-magnitude estimates of the mass-loss rate of ~10-8 to 10-7M yr-1. These authors discussed the drastic consequences of high mass-loss rates for the formation of PNe, especially for objects with low-mass progenitors. For example, for a low-luminosity object like M 2-9, with a ~0.55 M core, the time needed to evolve from Teff~ 5000 K to its current temperature is ~100 000 yr when adopting the model mass-loss prescriptions, 2 × 10-9M yr-1, but less than ~5000 yr when adopting the observed rate (see Fig. 8). The latter transition time leads to a much better agreement with the dynamical lifetime of the nebulosity in this and other similar objects.

8. Conclusions

We have detected RRL emission at mm-wavelengths from the central ionized regions of CRL 618, M 2-9, and MWC 922 at a moment when these H ii region are in a very early stage of development after the ionization onset. The RRLs at mm-wavelengths trace the inner regions of these targets at scales 10–100 au, indicating bulk motion speeds of ~10–30 km s-1.

Our analysis of the RRL spectra and free-free continuum data included detailed modeling of CRL 618, M 2-9, and MWC 922 using the non-LTE radiative transfer code MORELI. This has enabled us to set constraints of the density and temperature profiles and morpho-kinematics in the ionized cores of these objects. Typical densities are found to range between ne 106 and 108  cm-3, following radial power-law distributions ne(r) r[ − 2.1: − 2.4 ]. The total mass of ionized gas at the core of these objects is 10-6–10-4M, which is currently being ejected at rates of 10-6–10-7M yr-1. Simple models of the free-free continuum emission have also been performed for the remaining five targets with no mm-RRL detections under the assumption that their ionized cores are spherical isothermal winds. The post-AGB mass-loss rates derived in these cases range between pAGB~ 3 × 10-7 and 4 × 10-6M yr-1. These rates are significantly higher than those adopted by current post-AGB evolutionary models. New evolutionary tracks should be made allowing post-AGB stars to lose mass at higher rates, to study their influence on the star heating rate, AGB-to-PN transition times, initial-final mass relation, etc.

Future prospects. Our models are necessarily simplified representations of the presumably more complex spatio-kinematical and physical structure of the dense stellar surroundings in these objects. We need RRL emission maps with 0.̋05-angular resolution to better constrain the properties of these very inner winds, in particular, to test our conjectures on nebular morphology and kinematics.

For example, these maps are needed to check the rotating “disk+wind” scenario in MWC 922; if confirmed, it will be possible to obtain from the disk kinematics a more reliable estimate of the mass of the central source(s) and, thus, to clarify its nature. Disentangling the composite nebular structure of MWC 922 and, in particular, isolating the wind from the disk is needed for an accurate interpretation of the mass loss in this case. In parallel, measuring the outflow rotation and studying the disk dynamics can provide fundamental information regarding the origin of the outflow and the role of the disk in the wind launching/collimation process. High-angular resolution maps of RRLs are also needed to confirm or not the expansive kinematics and the suspected velocity gradient in the 100 au inner regions of CRL 618 and M 2-9: Does it follow a Vexpr law as observed at much larger linear scales in these, and most, pPNe? Also, these maps are necessary to unveil the presence of jet-like features, of which indirect observational indicia exist but whose direct characterization remains elusive. High angular-resolution of the ionized core of M 2-9 is also critical for a truthful judgement of the binary system at its nucleus and to unravel the different ionized gas (and dust) components that contribute to the observed mm-continuum.

With its unprecedented capabilities, ALMA is the only facility in the world that is able to map with sufficient sensitivity and angular resolution the dense inner winds of young PNe at the precise moment when these come into existence and begin to be ionized.


2

CLASS is a world-wide software used to process, reduce and analyze heterodyne line observations maintained by the Institut de Radioastronomie Millimétrique (IRAM) and distributed with the GILDAS software, see http://www.iram.fr/IRAMFR/GILDAS

3

We used 1/σ2 weights, where σ is the rms of the baseline of the cross-scan.

4

The terms Rout(90%) and Rout(99%), which are also explicitly given in some cases (e.g., Table 5), refer to the radius that encircles ~90% and 99% of the mm-continuum flux, respectively.

5

For example, based on our theoretical and observational current understanding of photoevaporating disks; see, for example, Hollenbach et al. (1994), Yorke & Welz (1996), Báez-Rubio et al. (2013).

6

Incidentally, in the case of unperturbed wind expansion, the layers studied by Martín-Pintado et al. (1988) would have now entered the optically thick part of the ionized wind (encompassing ~150–650 au; see Fig. 5 in Tafoya et al. 2013) and, thus, should be best traced at present by RRLs and continuum observations at cm-wavelengths.

Acknowledgments

We thank the referee for his/her comments and very valuable suggestions. The data presented in this paper were reduced using the packages available in the GILDAS software (http://www.iram.fr/IRAMFR/GILDAS). This work has been partially supported by the Spanish MINECO through grants CSD2009-00038, AYA2009-07304, AYA2010-2169-C04-01, AYA2012-32032, AYA2016-78994-P, FIS2012-39162-C06-01, ESP2013-47809-C03-01 and ESP2015-65597-C4-1. A. Báez-Rubio acknowledges support from DGAPA postdoctoral grant (year 2015) to UNAM. This research has made use of the SIMBAD database and VizieR catalog access tool (CDS, Strasbourg, France), the NASA Astrophysics Data System, and Aladin.

References

Appendix A: Spectra of the bonus 12CO (J = 2–1) and 13CO (J = 1–0) lines

In this Appendix we report the 12CO (J = 2 − 1) and 13CO (J = 1 − 0) spectral profiles observed toward our targets in this work (see Sect. 3). Spectra are channel-smoothed to a final velocity resolution of ~2 km s-1 and are shown in units of (K) and VLSR (km s-1) in Figs. A.1 and A.2.

All our sample targets are CO emission detections except for IRAS 20462+3416, which is also non-detected in previous surveys (e.g., Sánchez Contreras & Sahai 2012; He et al. 2014), and MWC 922. This is the first detection of CO emission in the young PN M 1-91 after some unsuccessful previous searches (e.g., Josselin et al. 2000; Bujarrabal et al. 2001). To our knowledge, these are also the first published CO spectra toward MWC 922. The narrow absorption and fewer emission features observed toward this object are most likely produced by intervening ISM clouds along the line of sight. For CO detections, we find line profiles and line fluxes in good agreement with previous observations, whenever available within typical ~20–30% calibration uncertainties (e.g., Bujarrabal et al. 2001; Sánchez Contreras & Sahai 2012; He et al. 2014).

thumbnail Fig. A.1

Spectra of the 12CO (J = 2 − 1) transition observed toward our sample targets (Table 1).

Open with DEXTER

thumbnail Fig. A.2

Spectra of the 13CO (J = 1 − 0) transition observed toward our sample targets (Table 1).

Open with DEXTER

All Tables

Table 1

Parameters of the sources observed in our pilot survey.

Table 2

Line parameters from Gaussian profile fitting for RRL detections.

Table 3

rms noise (in units) for RRL non-detections.

Table 4

Parameters of the emerging H ii regions deduced from line and continuum radiative transfer modeling (Sect. 5).

Table 5

Results from modeling of the free-free continuum for mm-RRL non-detections (Sect. 5.4).

All Figures

thumbnail Fig. 1

Spectral energy distributions of our targets showing literature photometry data (black circles and triangles indicate detections and upper limits, respectively) and the 1.3 and 2.7 mm continuum fluxes measured by us (star-like symbols; Table 1). The error bars represent a conservative 30% flux uncertainty assigned to all continuum data points. The dotted line indicates a modified blackbody fit to the far-IR thermal dust emission, which is subtracted from the total (dust+free-free) emission observed at mm-wavelengths. No attempt has been made to fit the warm dust emission blueward of the SED peak. The dashed line represents a fit of the type Sffνα to the mm-to-cm data points. The blue solid line represents the addition of these two fits.

Open with DEXTER
In the text
thumbnail Fig. 2

Light curve of the 1 mm-continuum of CRL 618 using ancillary data compiled by Sánchez Contreras et al. (2004; black circles); new measurements are added from Reuter et al. (1997), Nakashima et al. (2007), Planck Collaboration Int. XVIII (2015), and this work (red circles). Some reference dates are indicated. Together with significant abrupt changes in short timescales of less than 0.1 yr, a periodic pattern with a period of ~17 yr is hinted (dotted line) on top of an overall, long-term smooth increase.

Open with DEXTER
In the text
thumbnail Fig. 3

Recombination lines observed toward the pPN CRL 618 (histogram) and synthetic line profiles (pink) from our model in Table 4 (Sect. 5.1). The bottom and top X-axis represent VLSR (km s-1) and frequency (MHz). The frequency of He and C α-transitions, some of them detected, are also indicated. The many additional features observed in the spectrum, with emission and absorption profile components, are molecular transitions produced in the C-rich molecular envelope of CRL 618. Some of the mm-RRLs observed are partially blended with molecular lines.

Open with DEXTER
In the text
thumbnail Fig. 4

Recombination lines observed toward M 2-9 (histogram) and synthetic line profiles (pink) from our model in Table 4 (Sect. 5.2).

Open with DEXTER
In the text
thumbnail Fig. 5

Recombination lines observed toward MWC 922 (histogram) and synthetic line profiles (pink) from our model in Table 4 (Sect. 5.3).

Open with DEXTER
In the text
thumbnail Fig. 6

Line parameters of the RRLs detected in CRL 618, MWC 922, and M 2-9 derived from the observed profiles (Table 2; filled circles) and from the models in Table 4 including non-detections (empty circles). Dotted lines are fits to the observed line parameters of the Hnα transitions (red filled circles). For CRL 618, we show the mean VLSR of the Hnα centroids (dotted line) and a linear fit of VLSR as a function of frequency (dashed line).

Open with DEXTER
In the text
thumbnail Fig. 7

Same as in Fig. 1 with the free-free continuum emission predicted by our model (red line). Model parameters are given in Table 4 for CRL 618, MWC 922, and M 2-9, and in Table 5 for the rest of the targets.

Open with DEXTER
In the text
thumbnail Fig. 8

Adapted from Figs. 12 and 5 of Bloecker (1995). Top: Post-AGB evolutionary tracks with observational data for RRLs and non-RRL detections (filled and empty symbols, respectively – Table 1). For MWC 922 and M 2-9, the uncertainty of the luminosity due to the uncertainty of the distance is indicated by the arrows. The stellar mass of the post-AGB remnant core of each track is indicated to its right. Vertical dotted lines delimit the ionization onset region (Teff~ 20 000–40 000 K) where our targets lie. Bottom: Post-AGB mass-loss vs. effective temperatures adopted by post-AGB models of 0.605, 0.696, 0.836, and 0.894 M (Bloecker 1995) and of 0.546 M and 0.565 M (from Schönberner 1983). As in the top panel, the symbols indicate the location of our targets using the values of pAGB deduced in this work (Tables 4 and 5; vertical error bars correspond to an uncertainty factor of 2). The thick dash- and dot-dashed lines indicate the mass-loss rate typically assumed by evolutionary models for our three objects with RRL detections.

Open with DEXTER
In the text
thumbnail Fig. A.1

Spectra of the 12CO (J = 2 − 1) transition observed toward our sample targets (Table 1).

Open with DEXTER
In the text
thumbnail Fig. A.2

Spectra of the 13CO (J = 1 − 0) transition observed toward our sample targets (Table 1).

Open with DEXTER
In the text

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

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

Initial download of the metrics may take a while.