Gas-phase Molecules in Protoplanetary Nebulae with the 21 μm Emission Feature

It has been more than 30 years since the enigmatic 21 μm emission feature was first discovered in protoplanetary nebulae (PPNs). Although dozens of different dust carrier candidates have been proposed, there is as yet no widely accepted one. We present the results of molecular observations toward 21 μm objects using the 10 m Submillimeter Telescope of Arizona Radio Observatory in the 1.3 mm band and the 13.7 m telescope of Purple Mountain Observatory in the 3 mm band, aiming to investigate whether the gas-phase environments of these unusual sources have some peculiarities compared to normal PPNs. We detect 31 emission lines belonging to seven different molecular species, most of which are the first detection in 21 μm PPNs. The observations provide clues to the identification of the 21 μm feature. We report a correlation study between the fractional abundance of gas-phase molecules and the strengths of the 21 μm emission. Our study shows that, given the small sample size, the 21 μm feature has weak or no correlations with the gas-phase molecules. Future radio observations of high spatial and spectral resolution toward a large sample are desirable to elucidate the 21 μm emission phenomena.


INTRODUCTION
Low-and intermediate-mass (0.8-8.0 M ⊙ ) stars eject material into the interstellar medium (ISM) during their evolution from the asymptotic giant branch (AGB) to the white-dwarf phase, creating circumstellar envelopes (CSEs, see, e.g., Höfner & Olofsson 2018).CSEs can serve as a natural chemical laboratory and play a significant role on cosmic gas/dust cycle.However, we are still far from fully understanding circumstellar chemistry (see, e.g., Ziurys 2006;Cernicharo et al. 2011).The chemical processes in the protoplanetary nebula (PPN) phase that is a transition phase from the AGB to white-dwarf stages with a timescale of about ∼10 3 yr (Herwig 2005) are particularly unclear.A PPN consists of a post-AGB star with a effective temperature of less than 3 × 10 4 K and a separated shell.The central star is not hot enough to ionize the shell, but the growing ultraviolet (UV) radiation field could process the dust and gas in PPNs and stimulate photochemistry.
Protoplanetary nebule exhibit various infrared (IR) unidentified features, among which the socalled 21 µm emission feature is one of the most intriguing enigmas (Kwok 2021).This feature was first discovered in four PPNs by Kwok et al. (1989) from the Inf rared Astronomical Satellite (IRAS) Low-Resolution Spectrometer (LRS) database.Subsequent studies of Volk et al. (1999) expanded the sample of 21 µm feature objects and show that this feature has a consistent line profile from source to source and actually peaks at 20.1 µm.The objects with the 21 µm feature are rare, suggesting that this feature may arise from a transient molecule (Kwok et al. 1989).So far, this feature has been unambiguously detected only in 18 and 9 PPNs in the Galaxy and the Large/Small Magellanic Clouds, respectively (Cerrigone et al. 2011;Volk et al. 2011), suggesting that it must be a transient phenomenon.Accurate identification of its carrier is crucial for understanding circumstellar chemistry, but is a subject of considerable debate (see Volk et al. 2020 for a recent review).
Proposed carriers of the 21 µm feature include SiS 2 , Fe m O n , TiC, thiourea, hydrogenated fullerenes, polycyclic aromatic hydrocarbons, hydrogenated amorphous carbon, nano-diamonds, etc. (see, e.g., Volk et al. 2020, and references therein).No definite consensus has been reached.Zhang et al. (2009) conducted a quantitative investigation for nine kinds of inorganic candidates, including TiC nanoclusters, fullerenes coordinated with Ti atoms, SiS 2 , doped-SiC, SiO 2 -coated SiC dust, carbon and silicon mixtures, nano FeO, Fe 2 O 3 , and Fe 3 O 4 dust, and concluded that the cold FeO dust, which was first suggested by Posch et al. (2004), is a feasible candidate for the 21 µm emission feature carrier.However, the laboratory studies of Koike et al. (2017) suggested that the IR spectra of FeO vary depending on defects and disorders and the synthesis of FeO is sensitive to surrounding conditions, disfavoring FeO as the carrier of the 21 µm feature.Fullerenes have recently been discovered in circumstellar environments (Cami et al. 2010), which has renewed interest in the idea that hydrogenated fullerenes might be the 21 µm feature carrier (Zhang 2020a,b;Zhang et al. 2020).These fullerene derivatives may partly survive, be transferred into more complex compounds during the dispersion of CSEs, and finally contribute to other unidentified spectral features such as anomalous microwave emission (Iglesias-Groth 2006) and diffuse interstellar bands (Omont 2016).According to Mishra et al. (2016), there is a strong correlation between the mass-loss rates in the AGB phase and the flux of the 21 µm feature.They proposed that its carrier has been formed in the AGB phase but is unable to be excited by sufficient UV photons until the PPN phase.
It is established that the gas-phase chemistry and the dust processes depend on the physical conditions and the evolutionary stages of CSEs.Different astronomical environments are characterized by different gas-phase molecular patterns (e.g., Cernicharo et al. 2011;Aladro et al. 2015;Jørgensen et al. 2020;Taniguchi et al. 2023).Therefore, it is instructive to statistically investigate the gas-phase molecules in 21 µm sources, which could presumably provide some clues on the origin of the 21 µm feature.Such studies have been rather scarce because the molecular lines in 21 µm sources are typically very faint.The majority of earlier investigation were focused on CO only (e.g., Woodsworth et al. 1990;Hrivnak & Bieging 2005;Nakashima et al. 2012).Zhang (2020a) reported a line survey of the 21 µm source IRAS 22272+5435 at 3 and 1.3 mm bands, which resulted in the detection of ten molecules and three isotopologues.A large sample is still required to obtain statistically meaningful results.
In this paper, we present observations of gas-phase molecules in six 21 µm sources, aiming to investigate the relation between the 21 µm feature and the gas-phase molecules.The paper is organized as follows.In Sect.2, we present the sample, observations, and data reduction.In Sect.3, we present the observation results as well as the calculations of column densities, fractional abundances, and isotopic ratios.In Sect.4, we discuss the correlations between the strengths of the 21 µm emission feature and the measurements of the molecular lines in our sample.

Description of the sample
The sample of 21 µm sources is listed in Table 1.They were chosen for observation because of their comparatively strong CO emission.The properties of the sample are described below.
IRAS 06530-0213: Its central star is an F5 I spectral type supergiant with an effective temperature of 6 900 K (Hrivnak & Reddy 2003).The star is strongly carbon-enriched with a C/O abundance ratio of 2.8 and has a sub-solar metallicity with [Fe/H] = −0.5 (Reyniers et al. 2004).The significantly enhanced s-process elements indicate that the star has experienced a strong third dredge-up (Reyniers et al. 2004).The optical image taken by the Hubble Space Telescope (HST) shows an elliptical nebula (∼2 ′′ × 1 ′′ ) with limb-brightened bipolar lobes (Ueta et al. 2000).
IRAS 07134+1005: This PPN, also named HD 56126 or SAO 96709, has an elliptical morphology with an extended circular reflection nebula (Ueta et al. 2000).It is the biggest 21 µm source, with an optical size of 4 ′′ , according to the HST image.A two-peaked, limb-brightened dust torus is plainly seen in the mid-IR (MIR) image (Meixner et al. 1997).Its central supergiant star has a variability period of 12.1 days with an effective temperature of 7 250 K and a spectral class of F5 Iab (Bakker et al. 1997;Ueta et al. 2000).The [Fe/H] and C/O ratio are −1.0 and 1.2, respectively (Reddy et al. 2002;Reyniers et al. 2004).It is one of the four sources in which the 21 µm emission feature was first discovered (Kwok et al. 1989).
IRAS 19500+1709: This PPN, also named HD 187885 or SAO 163075, has a central star with an effective temperature of 8 000 K and an optical spectral type of F2-6 Ia, a metallicity of [Fe/H] = −0.7 to −0.6, and a C/O ratio of 1.1 (Hrivnak et al. 1989;Van Winckel & Reyniers 2000;Reddy et al. 2002;Reyniers et al. 2004).It is an intrinsically polarized bipolar reflection nebula with a size of 3. ′′ 6 × 2. ′′ 2 at the near-IR (NIR) J−band, consisting of two lobes divided by a dark lane that is perpendicular to the major axis (Gledhill et al. 2001).The 21 µm emission feature in IRAS 19500+1709 was first detected by Justtanont et al. (1996) using the United Kingdom Infrared Telescope (UKIRT).

SMT 10 m Observations
Five 21 µm sources were observed in the 1.3 mm window with the Arizona Radio Observatory 10 m Submillimeter Telescope1 (SMT) at Mt. Graham on 2014 April 21-26 and 2015 April 15-21.The study is complemented by the 1.3 mm spectra of IRAS 22272+5435 reported by Zhang (2020a).We used the 2048-channel acousto-optical spectrometer and the 1024-channel Forbes Filterbanks simultaneously with a frequency channel spacing of 500 kHz and 1 MHz, respectively, to cover the frequency range from 216-273 GHz with several gaps.The selections of these frequency scans are based on the relatively strong lines discovered in IRC+10216 and CRL 2688 (He et al. 2008;Zhang et al. 2013).As the 21 µm sources are much fainter, it is reasonable to suppose that no line in the frequency gaps is detectable with our instrumental sensitivity.The half-power beam widths (HPBW) are in the range of 27. ′′ 7-34.′′ 8, which are larger than the source sizes derived from their mid-IR images.The standard wobbler switching mode was carried out with an azimuth beam throw of ±2 ′ .The on-source times of each band are about 1-3 hours.The dual-channel superconductor-insulator-superconductor (SIS) receivers in single-sideband dual polarization mode were applied.The typical system temperatures are about 200-500 K. Pointing and focus were checked every 2 hours with nearby planets.

PMO 13.7 m Observations
We also observed the the J = 1 → 0 lines from CO, 13 CO, and C 18 O toward the six 21 µm sources by using the Purple Mountain Observatory (PMO) 13.7 m millimeter telescope 2 at Delingha in China.The observations were performed during May 11-20 2022.A supplementary observation toward IRAS 06530-0213 was carried out on November 9 2023.The on-source time of each band is 0.5-1 hour.The receiver used is a two-sideband SIS mixer.The spectrometer backends for each sideband utilize a Fast Fourier Transform Spectrometer (FFTS) with a bandwidth of about 1 GHz including 16384 channels and a frequency resolution of 61 kHz.The corresponding velocity resolution is 0.16 km s −1 for the CO J = 1 → 0 line at the rest frequency of 115 GHz.The two-sideband SIS receiver was designed to cover the CO J = 1 → 0 line at the upper sideband and the J = 1 → 0 transitions of 13 CO and C 18 O at the lower sideband simultaneously with a specific local oscillator (LO) frequency (Shan et al. 2012).The observations were made using a beam switch method.Regular pointing checks were made with a precision better than 5 ′′ .The HPBWs of this telescope are about 55 ′′ and 52 ′′ at 110 and 115 GHz, respectively.All sources are thus unresolved.

Data reduction
The spectral data were reduced by using the Grenoble Image and Line Data Analysis Software (GILDAS) of the Continuum and Line Analysis Single-dish Software (CLASS) package3 (Pety 2005).Bad scans caused by bad atmospheric conditions and receiver instabilities were discarded.The artificial features at the edge of each band were removed.A low-order polynomial baseline subtraction was carried out for each scan.The scans were co-added by using the inverse noise as the weight.The antenna temperature (T * A ) is converted to the main-beam temperature (T mb ) following the expression of T mb = T * A /B eff , where B eff is the main beam efficiency of 0.7 and 0.6 for the SMT 10 m and PMO 13.7 m telescopes, respectively.To further increase the signal-to-noise ratio, we smoothed the spectra over several channels using the 'boxcar' method.The resultant root mean square (RMS) levels and velocity resolutions of the transitions are listed in Tables 2-7.

RESULTS
The spectra for the whole frequency range covered by our observations are shown in Figs. 1 and  2, where for a better view of the detected lines against the noise we have binned several adjacent channels.Thirty-one emission lines assigned to seven molecular species and isotopologues, including CO, 13 CO, HCN, H 13 CN, HNC, CN, and HC 3 N, are detected.Figs.3-8 display the profile of each line separately.SiS and SiO are marginally detected in a few sources.The density of the detected lines is about one line per GHz.Except for IRAS 22272+5435, none of the 21 µm sources display SiC 2 , a compound that is frequently detected in C-rich CSE.In previous work, CO, 13 CO, and HCN were the only molecules detected in these sources.The velocity-integrated intensities of strong lines are obtained through the stellar-shell model fitting using the CLASS package (Pety 2005).The fitting function of the stellar-shell model is where A is the integrated intensity in K MHz, ν 0 and ∆ν are the middle frequency and full width at zero level in MHz.H dictates the Horn/Center parameter, which represents the edge-to-center intensity ratio, with the expression of f (∆ν/2)/f (0) = 1 + H.As a result, we can estimate the expanding velocity by V exp = c∆ν/(2ν 0 ) with c the light speed.The faint lines and the blended CN N = 2 → 1 fine structures cannot be fitted by the stellar-shell model.Tables 2-7 present the fitting results, where the rest frequency, molecular species, upper and lower energy levels, expanding velocity (V exp ), Horn/Center parameter (H), intensity at the line center (T c ), velocity-integrated intensity ( T mb dV ), RMS of the baseline, line center velocity (V LSR ), velocity resolution (D v ), HPBW of the telescope at the corresponding frequency are listed.In the following, the observational results of each source are reported.

IRAS 06530-0213
The CO J = 2 → 1 line in this object was previously reported by Hrivnak & Bieging (2005) using the 10 m diameter single dish Heinrich Hertz Submillimeter Telescope (HHT), which exhibits a slightly asymmetric double-peak profile.In this work, we detect five emission lines from gas-phase molecules, including the CO J = 1 → 0 and 13 CO J = 2 → 1 lines, the J = 3 → 2 transitions of HCN and HNC, and the CN N = 2 → 1 line, as shown in Fig. 3.A narrow feature with a width of 2 km s −1 is visible at the blue side of the CO J = 1 → 0 line.Based on the Herschel archive data and CO map of the Milk Way (Dame et al. 2001), it is likely to be a contamination from a nearby interstellar cloud lying within the telescope beam.The broad CO J = 1 → 0 feature and the HCN J = 3 → 2 line show a parabolic profile, suggesting that they are optically thick.The CO J = 1 → 0 and 13 CO J = 2 → 1 lines exhibit a prominent absorption at about 45 km s −1 , which has also been seen in the CO J = 2 → 1 line (Hrivnak & Bieging 2005).The 13 CO J = 2 → 1 line shows an emission feature at 55-70 km s −1 , which is not detected in the other transitions.The mean line width of CO J = 1 → 0, 13 CO J = 2 → 1, and the J = 3 → 2 transitions of HCN and HNC indicate an expansion velocity of 9.7 ± 0.5 km s −1 .The HNC J = 3 → 2 line shows an asymmetric double-peak profile with a stronger red peak.

IRAS 07134+1005
The first detection of CO in this object was made by Zuckerman et al. (1986) using the Five College Radio Astronomy Observatory (FCRAO) 14 m telescope through the J = 1 → 0 transition.The CO J = 1 → 0, 2 → 1, and 4 → 3, 13 CO J = 2 → 1, and HCN J = 1 → 0 transitions were previously detected in this object by Caltech Submillimeter Observatory (CSO) 10.4 m telescope, Institut de Radioastronomie Millimétrique (IRAM) 30 m telescope, and HHT (Bujarrabal et al. 1992;Omont et al. 1993;Knapp et al. 1998Knapp et al. , 2000;;Hrivnak & Bieging 2005).The HCO + J = 1 → 0 line was possibly detected by Bujarrabal et al. (1992) using the IRAM 30 m telescope.As shown in Fig. 4, our observations clearly show the presence of CO J = 1 → 0, 13 CO J = 2 → 1, HCN J = 3 → 2, and CN N = 2 → 1 lines.The profiles of the CO J = 1 → 0 and 13 CO J = 2 → 1 transitions show a smoothed parabolic and a double-peaked profile, corresponding to optically thick and thin cases, respectively.Bujarrabal et al. (1992) found that the HCN J = 1 → 0 line, in spite of the poor signal-to-noise ratio, seems to be shifted to the red side, which can be ascribed to self-absorption.Our higher signal-to-noise ratio observations show that the HCN J = 3 → 2 line clearly has a stronger emission at the red side, supporting the result of Bujarrabal et al. (1992).Based on an interferometric observation of the Berkeley-Illinois-Maryland-Association (BIMA) millimeter array, Meixner et al. (2004) found that IRAS 07134+1005 is a low-mass clumpy PPN with a size extending to about 7 ′′ in radius.Nakashima et al. (2009) constructed a morphokinematic model based on the Submillimeter Array (SMA) interferometric observations of the CO J = 3 → 2 line, which indicates the presence of a central torus with an inner and outer radius of 1. ′′ 2 and 3. ′′ 0, respectively.The mean width of the CO J = 1 → 0, 13 CO J = 2 → 1, and HCN J = 3 → 2 lines in our spectra suggests an expansion velocity of 9.9 ± 0.2 km s −1 .

IRAS 19500-1709
The first detection of CO species in IRAS 19500-1709 was made by Likkel et al. (1987) in a survey of the CO J = 1 → 0 line toward the coldest IRAS objects performed with the IRAM 30 m telescope.They found that IRAS 19500-1709 has the lowest circumstellar expansion velocity (∼11 km s −1 ) among their sample.Subsequent observations of higher-J CO transitions indicate the existence of a speedy outflow with a velocity of ∼30 km s −1 (Knapp et al. 1989;Likkel et al. 1991;Bujarrabal et al. 1992;Omont et al. 1993;Hrivnak & Bieging 2005).The 13 CO and HCN species have also been detected in previous observations (Bujarrabal et al. 1992).The CO J = 1 → 0, 13 CO J = 2 → 1, HCN J = 3 → 2, and CN N = 2 → 1 transitions, are visible in our spectra (Fig. 5).The CO J = 1 → 0 and 13 CO J = 2 → 1 lines show broad and narrow components.Their peak velocities show a small offset (5.2 km s −1 for CO and 4.2 km s −1 for 13 CO).This is primarily due to the rather sensitive dependence of the fitting on the faint broad wing, and cannot be considered as real.Nevertheless, the broad component suggests a outflow velocity of 29-33 km s −1 , which is consistent with previous results.The two components are not distinct for the HCN J = 3 → 2 line due to low signal-to-noise ratio.An inspection of the Fig. 5 clearly indicates that the stellar-shell fitting for the HCN line is unreliable.We thus mark it as an uncertain measurement, and assume a 50% uncertainty for the analysis.Bujarrabal et al. (1992) discovered a CO J = 2 → 1 absorption at 13 km s −1 , and ascribed this to a slow envelope.Such an absorption feature is also detected in our 13 CO J = 2 → 1 spectrum.Differing from Bujarrabal et al. (1992), we discover a 13 CO emission component at ∼10 km s −1 .The mean line widths suggest the expansion velocities of the broad and narrow components to be 9.2 ± 0.5 and 30.0 ± 1.8 km s −1 , respectively.

IRAS 22272+5435
The detection of CO in IRAS 22272+5435 was first reported by Woodsworth et al. (1990) through the J = 2 → 1 transition using the JCMT.Bujarrabal et al. (2001) observed the J = 1 → 0 and 2 → 1 lines of CO and 13 CO using the IRMA 30 m telescope.The CO J = 2 → 1 and 4 → 3 lines were detected by Hrivnak & Bieging (2005) using the HHT.We observe the J = 1 → 0 transitions of CO and 13 CO, which yield an expansion velocity of 7.5 ± 0.7 km s −1 .As shown in Fig. 6, their line profiles are smoothed parabolic and double-peaked, respectively, corresponding to optically thick and thin cases, respectively.Interferometric observations toward IRAS 22272+5435 have been carried out by Fong et al. (2006) using the BIMA array in the CO J = 1 → 0 line and by Nakashima et al. (2012) using the Combined Array for Research in Millimeter-wave Astronomy (CARMA) in the CO J = 2 → 1 line.The morphologies of the CO J = 2 → 1 and J = 1 → 0 transitions appear to be different.The CO J = 1 → 0 mapping shows a roughly spherical envelope with a size of about 20 ′′ and some unusual protrusions.The CO J = 2 → 1 image clearly resolves the central structure, which exhibits two intensity peaks with a distance of about 1. ′′ 4 and that their connecting line is almost in the vertical direction of the two 12.5 µm peaks of Ueta et al. (2001).Thus far, the detected molecules in this object include C 2 H, C 4 H, HC 3 N, SiC 2 , CS, C 34 S, CO, 13 CO, CN, CH 3 CN, HCN, H 13 CN, and HNC (see the details in Zhang 2020a).

IRAS 23304+6147
The first detection of CO in IRAS 23304+6147 was made by Woodsworth et al. (1990) using the JCMT through the J = 2 → 1 transition.Likkel et al. (1991) and Hrivnak & Bieging (2005) observed the CO J = 1 → 0 and 4 → 3 transitions.Omont et al. (1993) tentatively detected the CS J = 3 → 2 transition using the IRAM 30 m telescope.Our observations result in the detection of six lines, including the CO J = 1 → 0, 13 CO J = 2 → 1, J = 3 → 2 transitions of HCN, HCN, and H 13 CN, and CN N = 2 → 1, as shown in Fig. 7.The 13 CO J = 2 → 1 and HCN J = 3 → 2 lines are composed of narrow and broad components.Similar to the CO J = 1 → 0 line in IRAS 19500-1709, the two components of the HCN J = 3 → 2 line show slightly different peak velocities with an offset of 3.2 km s −1 .It ought to be regarded as a fitting error because the fitting of the broad component is substantially influenced by its faint wing.The CO J = 1 → 0 line shows an asymmetric shape with stronger emission at the blue side.A similar profile has been reported by Hrivnak & Bieging (2005) for the CO J = 2 → 1 line.We observe a narrow emission feature at the velocity of about 0 km s −1 for the CO J = 1 → 0 line, which may stem from the interstellar medium close to the Galactic plane.The mean line widths result in expansion velocities of 10.3 ± 0.9 and 22 ± 0.8 km s −1 for the narrow and broad components, respectively.

IRAS Z02229+6208
Prior to this work, CO was the only gas-phase molecule detected in this object, whose J = 2 → 1 and 3 → 2, 4 → 3 transitions were observed by Hrivnak & Kwok (1999) and Hrivnak & Bieging (2005).The J = 1 → 0 transitions of CO and 13 CO, the J = 3 → 2 transitions of HCN, H 13 CN, and HNC, CN N = 2 → 1, and HC 3 N J = 24 → 23, 25 → 24, and 28 → 27 lines are detected in our observations.HCN, H 13 CN, HNC, and HC 3 N are newly discovered in this object.As shown in Fig. 8, the J = 1 → 0 transitions of CO and 13 CO and the J = 3 → 2 transitions of HCN and HNC show both broad and narrow components.The broad component is likely to originate from fast outflows.There is an emission feature at the velocity of about -5 km s −1 for CO J = 1 → 0 and 13 CO J = 1 → 0 and J = 2 → 1 lines; we infer that it may stem from the interstellar medium.The mean line widths indicate expansion velocities of 10.3 ± 0.9 and 22 ± 0.8 km s −1 for the broad and narrow components, respectively.

Column densities and fractional abundances
Three HC 3 N transitions are clearly detected in IRAS Z02229+6208 (Fig. 7), allowing us to perform a rotation-diagram analysis to derive the excitation temperature (T ex ) and total column density (N) of HC 3 N in this PPN.Based on the assumptions of local thermodynamical equilibrium (LTE) and optical thinness, the level populations can be expressed in the function of ln where k and h represent the Boltzmann constant and Planck constant, respectively.The rest frequency ν (in Hz), upper-level energy E u /k (in K), partition function Q(T ex ), spontaneous emission coefficient A ul , and upper state degeneracy g u of each transition are taken from the CDMS catalog (Müller et al. 2001(Müller et al. , 2005)).Based on the assumption of both the source brightness distribution and the antenna beam being Gaussian profiles, the line intensity T mb in the scale of main beam temperature is converted to the source brightness temperature T S by the expression of The θ b value is equal to the HPBW of each transition, as listed in Table 7.The source size θ s is simply assumed to be 10.′′ 4, a value that has been used to estimate the mass loss rate of IRAS Z02229+6208 by Mishra et al. (2016).
T S dV is the velocity-integrated intensity of each line.We derive T ex = 88.4 ± 27.7 K and N = (1.37 ± 0.76) × 10 13 cm −2 for HC 3 N in IRAS Z02229+6208.The resultant rotation diagram is shown in Fig. 9. Assuming T ex = 80 K and taking the source sizes from Mishra et al. (2016) (see Table 1), we can estimate the column densities of other molecules through Eq. 2. The analysis results are presented in Table 8.The assumption of a constant T ex is only a rough approximation, which however does not cause significant errors under the condition of E u /k < 80 K.
To determine the fractional abundances with respect of H 2 (f X ), we utilize the equation proposed by Olofsson (1996), where A ul is in s −1 , V exp is the expansion velocity in km s −1 , θ b is the beam size in arcsecond, D is the distance in pc, ν is the rest frequency in GHz, T mb dV is the velocity-integrated intensity in K km s −1 , E l /k is the energy of the lower level in K, ṀH 2 is the H 2 mass loss rate in M ⊙ /yr, and x i and x e represent the ratios of the inner (R i ) and outer (R e ) radii over the beam size, that is, x e,i = R e,i /(θ b D) (see Table 1 for these values).If more than one transition from a molecule is detected, we adopt the average abundance.The obtained f X values are listed in Table 8.The indicated uncertainties come only from the measurement uncertainties, and may be significantly larger if considering the uncertainties associated with the distance estimation and the assumed gasto-dust ratio.Nevertheless, the relative abundances between different molecules are free from the uncertainties of the distances and gas-to-dust ratios.It should be noted that if the lines are optically thick, the f X values should be considered as the lower limits.Ramstedt & Olofsson (2014) investigated 19 C stars and found the fractional abundance of 13 CO to be (0.1-5.5) × 10 −4 with a mean value of 1.3 × 10 −4 , systematically higher than the f13 CO values of 21 µm sources of (1-2) × 10 −5 .The CSE expansion leads to a decreased optical depth, and hence the shielding of the CO against the ambient interstellar radiation field is declining.As a result, more 13 CO may be photo-dissociated.On the other hand, the more extensive penetration of UV photons, either from the center star (or a companion) or the external radiation field, may ionize more molecules and drive photochemistry during the evolution from the AGB to PPN phase.
Table 9 compares the molecular abundances of the 21 µm sources and those of the prototypical Cstar IRC+10216 and PPN CRL 2688.We find that despite having molecular abundances more than one order of magnitude lower than the two prototype objects, the 21 µm sources' abundance pattern more closely resembles CRL 2688 than IRC+10216.Therefore, it seems that molecular spectra could serve as valuable means for classifying CSEs.
3.8. 12C/ 13 C abundance ratios Isotopic abundance ratios can shed some insights on the stellar nucleosynthesis.As an intermediate reaction product of the CNO cycle, 13 C can be brought to the surface of a red giant through the first dredge-up process. 12C is mainly generated by the triple-α reaction in AGB phase (Herwig 2005;Karakas & Lattanzio 2014).The 12 C/ 13 C ratios are estimated through CO and HCN lines and their 13 C isotopic lines.As no C 18 O lines are discovered, we estimate the lower limits of the 16 O/ 18 O ratios based on the observational sensitivity.The results are presented in Table 10.Because the CO and HCN are optically thicker than their isotopic lines, the 12 C/ 13 C ratios listed in Table 10 should be treated as lower limits.For IRAS Z02229+6208, the 12 C/ 13 C ratios obtained from the broad and narrow line components are listed.Their divergent values are most likely the cause of the two components having different optical depths.The lower limits of the 12 C/ 13 C ratios in all 21 µm sources are much lower than the solar value (∼89.3,Lodders 2003), which cannot be solely attributed to the influence of optical thickness.Therefore, the 21 µm sources exhibit a 13 C enhancements probably as a result of nonstandard mixing processes (Karakas 2010).We come to the conclusion that the isotopic ratios of 21 µm sources are representative of evolved stars in general.We also note that the extremely low 12 C/ 13 C ratios (<15) is a characteristic of J-type carbon stars that might be associated with binary evolution.

DISCUSSION AND CONCLUSION
The CN N = 2 → 1 transition consists of eleven hyperfine components, which are divided into two groups located at frequency ranges of 226.6-226.7 and 226.8-226.9GHz, respectively, as shown in Figs.3-8.Their intensity ratios may differ from intrinsic values due to the influence of optical depth.However, the observed intensity ratios between the 226.8 and the 226.6 GHz groups are lower than 1.2, which is significantly less than the intrinsic strength ratio (∼1.8).To replicate the huge deviation, unreasonably large optical depths are needed.The same circumstances have been discovered for other CSEs (Bachiller et al. 1997;Zhang 2020a).The pumping of CN lines by optical and/or near-IR radiation provides a more likely explanation for the strength abnormalities.
The interferometric observations of IRAS 07134+1005 show a geometrically thick expanding torus with an expanding velocity of 8 km s −1 (Nakashima et al. 2009).Based on a comparison of two sources, Zhang (2020a) finds that the 21 µm source IRAS 22272+5435 may have a relatively narrow line width, probably tracing a slowly expanding torus.The narrow line profile is also shown in our spectra of all the 21 µm sources.It is reasonable to speculate that the presence of a torus might provide a favorable condition for the formation of the 21 µm carrier.However, apart from the narrow component, the lines of IRAS 19500-1709, IRAS 23304+6147, and IRAS Z02229+6208 exhibit a broad wing with a width of 20-30 km s −1 , tracing velocity-enhanced components.The CO observations with the Atacama Large Millimetre/submillimetre Array (ALMA) reveal that IRAS 07134+1005 has a hollow shell with a slight velocity enhancement at the poles, while another 21 µm source, IRAS 16594-4646 (not part of the current sample), does indeed exhibit a dense torus (Ueta et al. 2018).Therefore, the 21 µm sources seem to have a range of morphologies.This could also be inferred from the mid-IR observations of Ueta et al. (2000).Follow-up interferometric observations of a larger sample are required to elucidate the details of the three-dimensional structures and dynamics of 21 µm sources.
As shown in Table 9, the relative molecular abundances exhibit a similar degree of dispersion as the strengths of the 21 µm feature.To establish a link between the 21 µm emission feature and the gas-phase molecular species, we investigate the correlations between the relative intensities of the 21 µm feature respective to total IR emission (F 21µm /F IR ) and molecular fractional abundance f X for the six 21 µm feature sources, as shown in Fig. 10.The F 21µm /F IR values equal to 4.0%, 8.1%, 1.0%, 2.3%, 4.7%, and 1.3% for IRAS 06530-0213, IRAS 07134+1005, IRAS 19500+1709, IRAS 23304+6147, IRAS 22272+5435, and IRAS Z02229+6208, respectively.We note that different authors have reported different intensities of the 21 µm feature because of different spectral decomposition approaches (Zhang et al. 2010;Mishra et al. 2015Mishra et al. , 2016)).For the purpose of consistency, the values and uncertainties presented in Mishra et al. (2016) are used in our work.It must bear in mind that the adopted uncertainties of the F 21µm /F IR ratio in the analysis are likely to be underestimated.We fit the data with a linearized power law of the form log(f , where f X is the fractional abundances of CO, 13 CO, HCN, H 13 CN, HNC, CN, SiS, and HC 3 N, as shown in Fig. 10.The fittings are performed using the LINMIX Bayesian regression package Kelly (2007), which takes into account the uncertainties in x-and y-directions and is able to handle the upper limits of the y values.The fitting results are presented in Table 11.The correlation coefficients are too low to have a strong statistical significance.However, some of the correlations would be enhanced if the points for the feeblest 21 µm source IRAS 19500+1709 are excluded.
In Fig. 10, we also examine the correlation between the 21 µm feature and the HNC/HCN and CN/HCN abundance ratios.[HNC]/[HCN] has been found to decrease with increasing gas kinetic temperature (Graninger et al. 2014, and the references therein) and thus has been taken as a thermometer of astronomical environments (e.g., Jin et al. 2015;Hacar et al. 2020;Long et al. 2021).We find a strong inverse correlation between [HNC]/[HCN] and F 21µm /F IR , suggesting that a warm environment may favor the formation of the 21 µm carrier.The modelling of Huggins & Glassgold (1982) suggests that CN in C-rich CSEs is the direct photodissociation product of HCN.This is supported by the observations that CN is enhanced with respect to HCN during the PPN-PN evolution (Bachiller et al. 1997).If one excludes the points for IRAS 07134+1005 and IRAS 19500+1709, which have the strongest and weakest relative intensities of the 21 µm feature in our sample, the [CN]/[HCN] values lie within a narrow range.If we assume that [CN]/[HCN] is an indication of a change in the UV radiation field, this might suggest that these sources have a similar radiation environment.The rough trend shown in Fig. 10 can be tentatively explained in terms of a scheme in which the 21 µm carrier is being formed/excited with increasing UV radiation and is destroyed once the radiation strength reaches a certain level.This seems to be in line with the notion of Kwok et al. (1989) that the 21 µm emission is a transient phenomenon.However, the interpretation is highly speculative and can be further complicated by other factors, such as the properties of the central star and the dust-to-gas ratio.
The intensity ratio between the 13 CO J = 2 → 1 line and the total IR emission can be taken as a proxy of the gas-to-dust ratio.If we assume a constant mass-loss ratio ( ṀH 2 ) and expansion velocity (V exp ), the ṀH 2 /V exp could be used to trace the average density of the CSE.A correlation study may shed new insights into the physical conditions of synthesizing the 21 µm carrier.However, we cannot find any strong correlation, as shown in Fig. 10.A large sample study is needed to draw meaningful conclusions.

Qiu et al.
Table 1.Stellar Parameters.Note-The values of F 21µm , F IR , expansion velocity (V exp ), and inner and outer radius of the circumstellar shell (r i and r e ) are taken from Mishra et al. (2015) and Mishra et al. (2016).The total dust mass (M dust ) and the duration age from AGB stage are taken from Mishra et al. (2015) and Mishra et al. (2016), which have been scaled to the Gaia DR3 distance.The average mass-loss rate Ṁ is deduced through the dust mass divided by the duration age under an assumption of the gas-to-dust mass ratio of 100.Notea Obtained by integrating the main-beam temperature over the velocity coverage from 20-50 km s −1 .

Source name
b Obtained by integrating the main-beam temperature over the velocity coverage from -20-110 km s −1 .
c Obtained by integrating the main-beam temperature over the velocity coverage from -30-60 km s −1 .

Symbol:
? means the values that cannot be derived from our observations because of insufficient sensitivity.
× means the values that cannot be obtained because of unsolved hyperfine structures.: means uncertain measurements.Notea Obtained by integrating the main-beam temperature over the velocity coverage from 60-85 km s −1 .
b Obtained by integrating the main-beam temperature over the velocity coverage from 20-150 km s −1 .
c Obtained by integrating the main-beam temperature over the velocity coverage from 30-1000 km s −1 .Notea Obtained by integrating the main-beam temperature over the velocity coverage from -15-60 km s −1 .
b Obtained by integrating the main-beam temperature over the velocity coverage from -30-110 km s −1 .
c Obtained by integrating the main-beam temperature over the velocity coverage from -30-50 km s −1 .Notea Obtained by integrating the main-beam temperature over the velocity coverage from 5-40 km s −1 .
a Obtained by integrating the main-beam temperature over the velocity coverage from -10-90 km s −1 .
a Obtained by integrating the main-beam temperature over the velocity coverage from -20-50 km s −1 .

Notea
The average values of fractional abundances of each species with a standard deviation.b The relative standard deviation (RSD) is the ratio between the standard deviation and the average value.c Taken from Qiu et al. (2022).
f Obtained by scaling the f (CO) value to the 12 C/ 13 C ratio of 35 from Agundez et al. (2012).g In M ⊙ yr −1 km s −1 .

Notea
The non-detection of the 13 CO J = 1 → 0 line suggests that its integrated intensity is lower than 3σ.b Taken from Zhang (2020a).
c The non-detection of the H 13 CN J = 3 → 2 line suggests that its integrated intensity is lower than 3σ.d Derived from the broad line component.
e Derived from the narrow line component.
Table 11.LINMIX fitting results between the fraction abundances of gas-phase molecules and the relative intensities of the 21 µm features.

Figure 1 .
Figure1.Spectra of the 21 µm sources obtained with the PMO 13.7 m telescope.The dotted vertical lines mark the strong lines in the spectrum of CRL 2688.To bring out the lines against the noise, 36 adjacent channels have been binned, resulting in a frequency resolutions of about 2 MHz.To avoid smearing we cut off part of the spectra of IRAS Z0229+6208, and the zoom-out spectra are shown in the insets.

Figure 2 .
Figure2.Spectra of the 21 µm sources obtained with the SMT 10 m telescope.The dotted vertical lines mark the strong lines in the spectrum of CRL 2688.To bring out the lines against the noise, three adjacent channels have been binned, resulting in a frequency resolutions of about 3 MHz.
Obtained by integrating the main-beam temperature over the velocity coverage from -80-80 km s −1 .bObtained by integrating the main-beam temperature over the velocity coverage from -70-10 km s −1 .c Obtained by integrating the main-beam temperature over the velocity coverage from -40-5 km s −1 .

Table 10 .
Lower limits of 12 C/ 13 C and 16 O/ 18 O isotopic abundance ratios.Objects CO/ 13 CO CO/ 13 CO HCN/H 13 CN CO/C 18 O Note-α is the linear correlation coefficient.