DISCOVERY OF MOLECULAR GAS AROUND HD 131835 IN AN APEX MOLECULAR LINE SURVEY OF BRIGHT DEBRIS DISKS

In our previous survey we searched for gas in 20 young (<50 Myr) debris disks with high fractional luminosity (${L}_{\mathrm{IR}}/{L}_{\mathrm{bol}}={f}_{d}\gt 5\times {10}^{-4}$) around A- and F-type stars (Moór et al. 2011a). The current project focuses on 20 additional debris disks. We included seven disks where the selection criteria were similar to those applied in the first survey, except that we allowed lower fractional luminosities down to f_d = 10⁻⁴. We added four targets (HD 10939, HD 17848, HD 161868, and HD 182681) that are similar to the gaseous 49 Cet and HD 21997 systems in that they harbor cold, very extended debris disks. By including six debris disks around G–K-type stars, we also extended the sample to Sun-like stars. Recently, Montgomery & Welsh (2012) and Welsh & Montgomery (2013) reported nightly variability of the Ca ii K absorption line in the spectra of some A-type debris disk host stars, indicating the presence of circumstellar gas. Note that among the known gaseous disks, β Pic, 49 Ceti, and HD 172555 also show similar variable Ca absorption (Ferlet et al. 1987; Montgomery & Welsh 2012; Kiefer et al. 2014). Thus, we selected three nearby (<75 pc) targets, HD 110411, HD 182919, and HD 183324, from this sample. Basic properties of the selected targets are summarized in Table 1.

Table 1.  Stellar Properties and CO Observations for Our Targets

ID	SpT	Dist.	Age	Memb.	${L}_{{\rm{IR}}}/{L}_{{\rm{bol}}}$	Prog.	${S}_{\mathrm{CO}(2-1)}$	${S}_{\mathrm{CO}(3-2)}$	${S}_{\mathrm{CO}(4-3)}$	MCO
		(pc)	(Myr)				(Jy km s−1)	(Jy km s−1)	(Jy km s−1)	(M⊕)
(1)	(2)	(3)	(4)	(5)	(6)	(7)	(8)	(9)	(10)	(11)
Young Debris Disks with ${f}_{d}\gt {10}^{-4}$
HD 16743	F0/F2III	58.9	30 (2)	⋯	5.1 × 10−4 (3)	2	⋯	<1.24	⋯	<4.6 × 10−5
HD 38206a	A0V	75.1	30	COL	1.4 × 10−4 (1)	2	⋯	<2.40	⋯	<1.4 × 10−4
HD 95086b	A8III	90.4	17	LCC	1.7 × 10−3 (3)	2	⋯	<1.95	⋯	<1.7 × 10−4
HD 121191	A5IV/V	130.0	17	LCC	4.9 × 10−3 (2)	4	⋯	<1.47	⋯	<2.6 × 10−4
HD 131488	A1V	150.0	16	UCL	6.0 × 10−3 (2)	4	⋯	<1.39	⋯	<3.3 × 10−4
HD 131835a	A2IVc	122.7	16	UCL	3.0 × 10−3 (5)	1, 2, 3	1.60 ± 0.78	2.74 ± 0.55	4.46 ± 2.99	>5.2 × 10−4d
η Tela	A0V	48.2	23	BPMG	2.4 × 10−4 (1)	4	⋯	<1.38	<5.36	<3.4 × 10−5
Older A-type Stars with Extended, Cold Debris Disks
HD 10939	A1V	62.0	346 (4)	⋯	1.2 × 10−4 (3)	2	⋯	<1.37	⋯	<5.6 × 10−5
HD 17848a	A2V	50.5	372 (4)	⋯	7.6 × 10−5 (3)	2	⋯	<1.42	⋯	<3.8 × 10−5
HD 161868	A0V	31.5	450 (3)	⋯	1.3 × 10−4 (3)	2	⋯	<1.82	⋯	<1.9 × 10−5
HD 182681	B8/B9V	69.9	144 (4)	⋯	2.8 × 10−4 (3)	2	⋯	<2.15	⋯	<1.1 × 10−4
Debris Disks around G–K-type Sun-like Stars
HD 105a	G0V	39.4	30	THA	1.9 × 10−4 (1)	3	⋯	<1.42	<6.08	<2.4 × 10−5
HD 377a	G2V	39.1	183 (5)e	⋯	2.6 × 10−4(1)	3	⋯	<1.03	<2.41	<1.7 × 10−5
HD 61005a	G3/G5V	35.3	40	ARG	1.6 × 10−3 (1)	3	⋯	<1.27	<4.37	<1.7 × 10−5
HD 92945	K1V	21.4	211 (5)e	⋯	4.9 × 10−4(1)	3	⋯	<1.33	<3.86	<0.7 × 10−5
HD 160305	F8/G0V	72.5	23	BPMG	8.6 × 10−5 (4)	3	⋯	<1.61	<5.54	<9.0 × 10−5
HD 202917	G5V	43.0	30	THA	2.8 × 10−4 (1)	3	⋯	<1.46	<7.22	<2.9 × 10−5
A-type Stars with Debris Disks and Variable Ca Absorption Line
HD 110411	A0V	36.3	90 (4)	⋯	5.9 × 10−5 (1)	5	<0.20	⋯	⋯	<0.6 × 10−5
HD 182919	A0V	72.9	198 (1)	⋯	3.4 × 10−5 (1)	5	<0.16	⋯	⋯	<2.0 × 10−5
HD 183324	A0V	61.2	140 (1)	⋯	2.2 × 10−5 (1)	5	<0.16	⋯	⋯	<1.4 × 10−5

Notes. Column (1): identification. Column (2): spectral type. Column (3): distance. Based on Hipparcos trigonometric parallax when it is available. For HD 121191 and HD 131488 distances were taken from Melis et al. (2013). Column (4): stellar age. For members of young moving groups and associations we adopted the age of the corresponding group. For the other objects the age information was taken from the literature (references are in brackets): (1) Chen et al. (2014); (2) Moór et al. (2011b); (3) Moór et al. (2015); (4) Nielsen et al. (2013); (5) Vican & Schneider (2014). Column (5): Membership status of the star. ARG: the Argus moving group; BPMG: the β Pic moving group; COL: the Columba moving group; LCC: the Lower Centaurus Crux association; THA: the Tucana-Horologium association; UCL: the Upper Centaurus Lupus association. Column (6): Fractional dust luminosity. Apart from HD 131835, the data were taken from the literature: (1) Chen et al. (2014); (2) Melis et al. (2013); (3) Moór et al. (2015); (4) Patel et al. (2014); (5) this work. Column (7): Observation projects. (1) E-083.C-0303; (2) M-087.F-0001; (3) M-092.F-0012; (4) M-093.F-0010; (5) IRAM-172-13. Columns (8)–(10): integrated line fluxes for CO J = 2–1, J = 3–2, J = 4–3, respectively (see Section 2.2). Column (11): estimated mass of the CO gas.

^aHD 17848, HD 105, HD 377, HD 61005, and η Tel were also observed at the CO (3–2) transition by Hales et al. (2014) using the APEX and ASTE radio telescopes; for HD 377 observation at the 2–1 transition is also available in Pascucci et al. (2006). HD 38206 and HD 131835 were searched for gas at the 2–1 transition of CO (Kastner et al. 2010; Zuckerman & Song 2012). All of these previous measurements produced nondetections for these targets. ^bCO observational data regarding this object are already published in Moór et al. (2013). ^cThe effective temperature of this star estimated by us (see Section 3.1) is rather more consistent with a spectral type of A4. ^dDerived from an LTE model (see Section 3.3). ^eVican & Schneider (2014) provided two different age estimates for these targets. Here we quoted the average of the two values.

Download table as:  ASCIITypeset image

Observations with APEX. The majority of our observations were carried out with the 12 m Atacama Pathfinder Experiment radio telescope (APEX; Güsten et al. 2006) in the framework of three different programs. All of our targets were observed in the J = 3–2 transition of ¹²CO at a rest frequency of ν = 345.796 GHz. In program M-087.F-0001-2011 the observations were performed with the Swedish Heterodyne Facility Instrument/APEX2 (SHeFI; Vassilev et al. 2008) receiver, while in M-092.F-0012-2013 and M-093.F-0010-2014 the First Light APEX Submillimeter Heterodyne receiver (FLASH+; Klein et al. 2014) was used. The latter is a dual-frequency receiver that operates simultaneously in the 345 and 461 GHz atmospheric windows; therefore, we obtained ¹²CO J = 4–3 line (ν = 461.041 GHz) observations as well. For HD 131835, CO (2–1) line observations were also conducted using the SHeFI/APEX1 instrument.

For SHeFI observations acquired before 2011 June we utilized the Fast Fourier Transform Spectrometer (FFTS) backend with a spectral resolution of 488 kHz (0.42 km s⁻¹ at J = 3–2 transition). For later SHeFI measurements the eXtended bandwidth Fast Fourier Transform Spectrometer (XFFTS; Klein et al. 2012), with 32,768 channels providing a spectral resolution of 77 kHz (0.066 km s⁻¹ at J = 3–2 transition), was connected to the receiver. FLASH+ was always connected to the XFFTS backend. Beam sizes of APEX are ∼27'', ∼18'', and ∼14'' at 230, 345, and 460 GHz, respectively. All of our APEX observations were performed in on-off observing mode.

The CO (3–2) spectral line observations of HD 131835 were further supplemented by yet unpublished measurements obtained in our previous APEX program E-083.C-0303 (for details see Moór et al. 2011a).

Observations with the IRAM 30 m telescope. The three targets selected because of variable Ca ii K absorption line were observed with the IRAM 30 m telescope using the multiband heterodyne Eight MIxer Receiver (EMIR; Carter et al. 2012) as part of IRAM program no. 172-13. We searched for gas in our targets at 230.538 GHz, the J = 2–1 transition of the ¹²CO. The observations were conducted in wobbler-switching on-off mode with a wobbler throw of 60''. For the backend, we used the new Fast Fourier Transform Spectrometer (FTS) with a frequency resolution of 200 kHz, providing a velocity resolution of ∼0.25 km s⁻¹ in the J = 2–1 transition. Table 2 summarizes the main characteristics of the observing programs.

Data reduction. Both the APEX and IRAM spectra have been processed using the GILDAS/CLASS package.¹³ For the final average spectrum, we discarded noisy scans, and a baseline was subtracted from each individual scan. The baseline was typically linear, except for a few cases where we used second-order polynomials. The final spectrum was derived as an average of the individual spectra weighted by the inverse square of their rms noise. As a final step, the obtained antenna temperatures were converted to line flux densities. For the APEX data we used Kelvin-to-Jansky conversion factors of 39, 41, and 48 Jy K⁻¹ for CO (2–1), CO (3–2), and CO (4–3) transitions, respectively. For the IRAM CO (2–1) observations the corresponding number was 7.8 Jy K⁻¹. These conversion factors were taken from the relevant APEX¹⁴ and IRAM¹⁵ home pages.

Outcome of the survey. CO emission was evident only for HD 131835, where the integrated line flux of the 3–2 transition was successfully detected at a 5.0σ level. Figure 1 shows the baseline-subtracted CO (3–2) profile of this target. The measured line has a width of roughly 12 km s⁻¹ and is centered at a velocity of 7.2 km s⁻¹ with respect to the local standard of rest (LSR), in very good agreement with the systemic LSR velocity of the star (6.1 ± 1.1 km s⁻¹, computed from the radial velocity derived in Section 3.1). The measured line tentatively shows a double-peak profile in accordance with what we expect for a gas disk in Keplerian rotation around a star. The peak flux is ∼370 mJy, and the integrated line flux—obtained by integrating the line over a 12 km s⁻¹ interval—is 2.74 ± 0.55 Jy km s⁻¹. By inspecting the same velocity interval in the J = 2–1 and 4–3 spectra (Figure 1), we found no statistically significant evidence for lines; the obtained integrated line fluxes are 1.60 ± 0.78 and 4.46 ± 2.99 Jy km s⁻¹, respectively. We note that HD 131835 was previously observed in the CO 2–1 transition by Kastner et al. (2010), resulting in a nondetection. This observation was about two times less sensitive than ours. Assuming a line width of 5 km s⁻¹, Zuckerman & Song (2012) deduced a 5σ line flux upper limit of 5.1 Jy km s⁻¹ from this measurement, which is compatible with our result.

Zoom In Zoom Out Reset image size

None of our other targets were detected at any of the CO transitions. For these sources, upper limits were computed as ${S}_{\mathrm{rms}}{\rm{\Delta }}v\sqrt{N}$, where S_rms is the 1 σ measured noise, Δv is the velocity channel width, and N is the number of velocity channels over an interval of 10 km s⁻¹. The obtained line fluxes and upper limits are listed in Table 1.

Because of the detection of CO (Section 2.2), for HD 131835 we performed several additional observations. In order to better characterize its SED at far-infrared/submillimeter wavelengths and to search for O i and C ii emission from the disk, we obtained photometric and spectroscopic measurements using the Photodetector Array Camera and Spectrometer (PACS; Poglitsch et al. 2010) and the Spectral and Photometric Imaging Receiver (SPIRE; Griffin et al. 2010) on board the Herschel Space Observatory (program id OT2_amoor_3). These data were complemented by mid-infrared spectra obtained with the Spitzer Space Telescope (Werner et al. 2004). To investigate the stellar properties, a high-resolution ground-based optical spectrum was also taken.

2.3.1. Herschel Observations

PACS maps. Far-infrared maps were obtained with the PACS detector on 2012 September 9 and 11 in mini-scan-map mode (PACS Observer's Manual v2.5)¹⁶ using a scan speed of 20'' s⁻¹ with 10 scan legs of 3' length separated by 4''. We made two scans with scan angles of 70° and 110° both at 70 μm (OBSID: 1342250784, 1342250785) and at 100 μm (OBSID: 1342250868, 1342250869). Since the PACS photometer observed in two bands simultaneously (at 160 μm in addition to 70 μm or 100 μm), these measurements provided four 160 μm scans as well.

Data processing was done with the Herschel Interactive Processing Environment (HIPE, Ott 2010) version 13 using PACS calibration tree No. 65 and the standard HIPE script optimized for the reduction of bright sources. Additionally, we used the recently developed "gyro" correction to reduce the pointing jitter. We selected only those data frames from the timeline where the actual scan speed of the spacecraft was between 15'' and 25'' s⁻¹. To eliminate the marked low-frequency (1/f) noise, we applied high-pass filtering with filter width parameters of 15, 20, and 35 for the 70, 100, and 160 μm data, respectively. In order to avoid flux loss, the 25'' radius vicinity of our targets was excluded from the filtering. For glitch removal we used the second-level deglitching algorithm. As a final step, in each band we applied the "photProject" task to combine all frames into a map using the default pixel fraction (1.0) and pixel sizes of 11, 14, and 21 at 70, 100, and 160 μm, respectively.

HD 131835 was clearly detected in all PACS bands. We performed aperture photometry on the source using a radius of 18'', while the sky background was estimated in an annulus between 40'' and 50''. The aperture was placed at the source's centroid position. The offsets between the derived centroids and the targets' optical position (corrected for the proper motion using the epochs of PACS observations) were <1'' in all bands. For sky noise determination, we distributed 16 apertures with radii of 18'' (identical to the the source aperture) randomly along the background annulus. By performing aperture photometry without background subtraction in each aperture, we computed the sky noise as the standard deviation of these background flux values. We applied aperture correction in each band, by using correction factors taken from the calibration file attached to the measurement. To derive the final uncertainty of our photometry, we added the measurement errors and the calibration uncertainty (7%; Balog et al. 2014) quadratically. The obtained flux density values and their uncertainties are listed in Table 3.

Table 3.  Measured and Predicted Photospheric IR Flux Densities for HD 131835

λ	Meas. fluxa	Instr.	Pred. Flux	Reference
(μm)	(mJy)		(mJy)
3.35	297.0 ± 10.3	WISE	302.9	Wright et al. (2010)
4.60	165.2 ± 5.5	WISE	168.1	Wright et al. (2010)
6.00	104.1 ± 10.2	IRS	101.6	This work
7.00	80.2 ± 7.0	IRS	75.5	This work
7.97	67.7 ± 4.5	IRS	58.7	This work
9.00	58.6 ± 3.6	IRS	46.4	This work
9.00	68.7 ± 7.0	IRC	46.4	Ishihara et al. (2010)
10.00	52.8 ± 3.4	IRS	37.8	This work
10.99	50.4 ± 2.6	IRS	31.3	This work
11.56	49.1 ± 2.3	WISE	28.4	Wright et al. (2010)
13.05	50.8 ± 2.8	IRS	22.3	This work
14.83	56.2 ± 4.9	IRS	17.3	This work
16.99	70.4 ± 6.5	IRS	13.3	This work
19.02	88.5 ± 7.7	IRS	10.6	This work
21.28	115.2 ± 15.6	IRS	8.5	This work
22.09	161.5 ± 9.7	WISE	7.9	Wright et al. (2010)
23.67	153.1 ± 3.1	MIPS	6.8	Chen et al. (2012)
24.48	160.8 ± 15.2	IRS	6.4	This work
25.00	186.0 ± 33.5	IRAS	6.1	Moshir (1990)
27.45	207.5 ± 19.7	IRS	5.1	This work
30.50	249.5 ± 19.1	IRS	4.1	This work
33.55	309.7 ± 23.4	IRS	3.4	This work
60.00	684.0 ± 61.6	IRAS	1.1	Moshir (1990)
70.00b	738.7 ± 52.5	PACS	0.8	This work
71.42	659.2 ± 33.4	MIPS	0.7	Chen et al. (2012)
100.00b	637.0 ± 45.5	PACS	0.37	This work
160.00	382.3 ± 27.9	PACS	0.14	This work
250.00	156.4 ± 11.5	SPIRE	0.06	This work
350.00	84.3 ± 8.3	SPIRE	0.03	This work
500.00	35.4 ± 8.9	SPIRE	0.01	This work
870.00	8.5 ± 4.4	LABOCA	0.005	Nilsson et al. (2010)

Notes.

^aThe quoted flux densities are not color corrected. ^bThe emission is marginally resolved at these wavelengths.

Download table as:  ASCIITypeset image

PACS spectroscopy. We carried out PACS spectroscopic observations centered on the O i 63 μm (OBSID: 1342248686) and the C ii 158 μm (OBSID: 1342248687) lines. These fine-structure lines, if arising from circumstellar disks, are unresolved in Herschel PACS observations. Observations were performed on 2012 July 29 using the Line Spectroscopy observing template that enables the coverage of a small wavelength interval around these lines. To efficiently eliminate telescope and sky background, the measurements were obtained in chop-nod mode. The PACS spectrometer array consists of 5 × 5 spectral pixels (spaxel), each having a size of 94. We made small 3 × 3 and 2 × 2 raster maps with raster point/line step sizes of 3'' and 2'' for the O i and C ii observations, respectively. The maps were centered on the target, and they were repeated two and four times in the case of O i and C ii observations, respectively.

Raster map observing mode has been selected to mitigate the risk of eventual pointing offsets of the telescope significantly larger than 2''. While the 47'' × 47'' spatial footprint of the integral field unit is sufficiently large to almost entirely cover the spectrometer spatial beam, a single spaxel of 94 × 94 size is strongly undersampling the beam, especially at the short wavelength of the O i line. As a consequence of undersampling, in the baseline pointed observing scheme the spectrum of a point source is extracted from the central spaxel of the 5 × 5 unit applying a correction for fluxes due to small pointing offsets from the center of the central spaxel and point jitter. However, if the source flux is lower than about a few janskys per spaxel, then correction factors estimated from the source continuum sampled by the eight neighbors of the central spaxel may become unreliable and the spectrum extracted from the central spaxel may have to be left uncorrected for mispointing.

To overcome this problem, we apply mapping observation of the same field of view with sub-spaxel raster step size. When data from raster positions are combined, the sampling of the beam is highly improved, even in the case of an eventual pointing offset of the raster central position. This technique directly provides a high-resolution image of the source where the spectrum can be extracted within a synthetic aperture centered on the measured peak position.

We apply the telescope normalization scheme for spectrophotometric flux calibration, where the Herschel telescope background calibrated on Neptune is used as an absolute radiometric reference. The method efficiently eliminates any drifts with time in the system response at the frequency of individual chopper cycles, resulting in the offset signals in nod A and nod B positions perfectly canceling each other. This method may improve signal-to-noise ratio—especially for faint sources in the sub-jansky per spaxel regime—comparing to the baseline flux calibration method, which relies on response estimates from calibration blocks executed at the beginning of the observation and propagated to the covered wavelength range through the relative spectral response function.

Flux-calibrated IFU cubes obtained at subsequent raster positions are combined into a single spectral cube using a 3D drizzling algorithm (Regibo 2012) available in the PACS interactive pipeline in HIPE. The strength of the drizzling algorithm is that the convolution with the detector footprint is minimized, especially when we apply the "gyro" improved high spatial resolution reconstruction of instantaneous pointing jitter. The projected pixel size is adjusted to Nyquist sample the beam at the line central wavelength.

In the last data processing step, spectra from the drizzled spectral cube are combined within an extraction aperture of 63 radius centered on the source position.

Figure 2 shows the obtained O i and C ii spectra. Though in the C ii spectrum there is a peak at 157.7 μm, its significance is only 2.3σ; thus, neither the O i nor C ii line were detected toward our target. In order to derive upper limits for the line fluxes, we used the continuum-subtracted spectra and determined the rms noise (σ_rms) in a 20-pixel-wide window centered on the expected line wavelength. The upper limits were computed as $3\times {\sigma }_{{\rm{rms}}}\times {\rm{\Delta }}\nu \times \sqrt{p},$ where ${\rm{\Delta }}\nu$ is the width of one pixel in Hz, while p is the width of an unresolved emission line in pixels. Values of p were calculated as the ratio of the FWHM of unresolved lines (taken from PACS Observer's Manual v2.5, 0.020 μm and 0.126 μm for O i and C ii spectra, respectively) to the appropriate pixel widths (in μm). The obtained upper limits for the O i and C ii lines are 1.5 × 10⁻¹⁸ Wm⁻² and 5.3 × 10⁻¹⁹ Wm⁻², respectively.

Zoom In Zoom Out Reset image size

SPIRE maps. We observed HD 131835 with the SPIRE photometer on 2012 August 9, in Small Scan Map mode (Spire Handbook v2.5¹⁷ ) at 250, 350, and 500 μm simultaneously. Data were processed with HIPE 13 following the standard pipeline processing steps with default values for all applied tasks. The final maps were produced using the "naiveMapper" task. The beam sizes are 176, 239, and 352 at 250, 350, and 500 μm, respectively, and the maps were resampled to pixel sizes of 6'', 10'', and 14'' at these wavelengths.

Our target is detected in all bands. We used point-spread function (PSF) photometry to extract the source fluxes; the corresponding SPIRE beam profiles were taken from the calibration context. This method yielded flux densities of 156.4 ± 7.6 mJy, 84.3 ± 6.9 mJy, and 35.4 ± 8.6 mJy at 250, 350, and 500 μm, respectively. In order to validate our flux values, the photometry was also performed using the SPIRE Timeline Fitter task in HIPE. Instead of using the final maps, this task fits two-dimensional elliptical or circular Gaussian functions to the baseline-subtracted timeline data of the photometer at the coordinates of the source (Bendo et al. 2013). The radius of the region that includes the peak of the source was set to 22'', 30'', and 42'' for 250, 350, and 500 μm maps, respectively. The background was estimated from an annulus between 300'' and 350'' centered on the source. At 250 and 350 μm we used elliptical Gaussian during the fitting; at 500 μm, because of the source's lower signal-to-noise ratio, we applied a circular Gaussian. Using the Timeline Fitter, we obtained flux densities of 165.4 ± 4.4 mJy, 77.8 ± 4.4 mJy, and 44.1 ± 8.3 mJy at 250, 350, and 500 μm, respectively, in good accordance with the ones derived from PSF photometry. In the following analysis, we will use the PSF photometry. The final uncertainties were derived as the quadratic sum of the measurement errors and the overall calibration uncertainty of the SPIRE photometer (5.5%; Bendo et al. 2013). The derived results are quoted in Table 3.

2.3.2. Spitzer/IRS Observation

2.3.3. Fiber-fed Extended Range Optical Spectrograph (FEROS) Observation

HD 131835 is an A2-type star (Houk 1982) located at ${122.7}_{-12.8}^{+16.2}$ pc away (van Leeuwen 2007). Based on astrometric data made with ESA's Hipparcos satellite, de Zeeuw et al. (1999) proposed the star to be a member of the ∼16 Myr old (Pecaut et al. 2012) UCL subgroup of the Sco–Cen association with a membership probability of 95%. According to our current knowledge, HD 131835 has no stellar companion (Kouwenhoven et al. 2005; Wahhaj et al. 2013).

We used spectroscopic and photometric data to model the stellar photosphere and to estimate the fundamental stellar properties of HD 131835. Moreover, supplementing the Hipparcos astrometric data by our ground-based radial velocity measurement, its membership in UCL was also reconsidered.

Analysis of the FEROS spectrum. In order to estimate the radial velocity and different stellar parameters of our target, the measured FEROS spectrum was compared with templates taken from the spectral library of Munari et al. (2005). In this spectral library the continuum normalization was performed by dividing the absolute flux spectrum by its calculated continuum. Because of the applied method, in some spectral regions the continuum of the model spectra can significantly differ from a continuum that would be obtained using a regular normalization scheme of observational spectroscopy. Therefore, first we used high-order polynomials to flatten the model spectra, and then the same method was applied in the continuum normalization of our observed spectrum as well. We used an iterative method to estimate the stellar parameters. We first applied the cross-correlation technique to derive the radial velocity of the star by convolving the measured spectrum with a template taken from the spectral library. We selected a template spectrum with [Fe/H] = 0.0, $\mathrm{log}g=4.0,$ and ${T}_{{\rm{eff}}}=8750\;{\rm{K}}$ (corresponding to the A2 spectral type of the target quoted by the SIMBAD). The cross-correlation function was calculated by the fxcor task in IRAF. After transforming the measured spectrum to the laboratory system, we compared it with a grid of Munari synthetic spectra, using the 4000–6200 Å wavelength range, excluding the Hβ region and Na D lines. The grid was compiled by varying the effective temperature, surface gravity, and projected rotational velocity of the model spectra, while the metallicity was fixed to [Fe/H] = 0.0. By finding the global minimum, we used the derived parameters to repeat the radial velocity determination and then recompute the stellar parameters. This process yielded the following best-fit parameters: ${T}_{{\rm{eff}}}=8250\pm 250\;{\rm{K}}$, $\mathrm{log}g=4.0\pm 0.5$, $v\mathrm{sin}i=100\pm 15\;\mathrm{km}\;{{\rm{s}}}^{-1}$. Based on dwarf spectral type versus effective temperature scale presented in Pecaut & Mamajek (2013), the derived T_eff rather corresponds to a spectral type of A4. The radial velocity was found to be 1.6 ± 1.4 km s⁻¹, consistent with our previous measurement of 3.3 ± 1.7 km s⁻¹ (Moór et al. 2006). Taking the weighted average of the two radial velocity values, we obtained v_r = 2.3 ± 1.1 km s⁻¹.

Analysis of photometric and near-IR data. To provide an additional independent estimate for fundamental stellar properties and model the stellar photosphere, we fitted an ATLAS9 atmosphere model (Castelli & Kurucz 2004) to the optical and near-IR observations of the target. Photometric data were taken from the Tycho 2 (Høg et al. 2000), Hipparcos (Perryman et al. 1997), and Two Micron All Sky Survey (2MASS; Cutri et al. 2003) catalogs. These data were further supplemented by Wide-field Infrared Survey Explorer (WISE) W1-band photometry at 3.4 μm from the WISE All-sky Database (Wright et al. 2010) and Johnson B, V photometry from the catalog of Slawson et al. (1992). With a distance of ∼123 pc, HD 131835 could be outside the Local Bubble; thus, its reddening might not be negligible. In order to provide a rough estimate of the reddening at the distance of HD 131835, we collected those stars in its vicinity (with separations <3°) that have both Hipparcos-based trigonometric parallax and measured Strömgren color indices and Hβ index in the photometric catalog of Hauck & Mermilliod (1998). For these stars we derived $E(B-V)$ color excesses (as $E(B-V)=E(b-y)/0.74$) from the Strömgren data by using the appropriate calibration processes (Crawford 1975, 1979; Olsen 1984). Figure 3 shows the derived E(B − V) values as a function of trigonometric parallaxes. At the distance of HD 131835 (marked by a vertical line) the reddening is between ∼0.015 and ∼0.075 mag. We note that, consistent with our results, the extinction map of Schlafly & Finkbeiner (2011) shows a reddening of ∼0.083 mag in the direction of our target. This value gives the total reddening within the Milky Way for the line of sight; thus, it can be considered an upper limit. Therefore, in the photosphere modeling we fitted both the effective temperature and the reddening (the latter was limited between 0.0 and 0.08 mag). By adopting solar metallicity and $\mathrm{log}g=4.25,$ our χ² minimization yielded T_eff = ${8250}_{-100}^{+250}\;{\rm{K}}$ and E(B − V) = ${0.025}_{-0.02}^{+0.05}.$ Considering these data and the ${\text{}}{\text{}}{Hipparcos}$-based trigonometric distance, we derived a luminosity of L_bol = 9.2 ± 2.6 L_⊙ for HD 131835.

Zoom In Zoom Out Reset image size

UCL membership. Using the new combined radial velocity and Hipparcos astrometric data of HD 131835, we computed a Galactic space motion of U = −5.8 ± 1.3 km s⁻¹, V = −18.4 ± 2.1 km s⁻¹, and W = −4.8 ± 0.8 km s⁻¹ with respect to the Sun. This space velocity is in perfect agreement with the charateristic space motion of UCL (U = −5.1 ± 0.6 km s⁻¹, V = −19.7 ± 0.4 km s⁻¹, W = −4.6 ± 0.3 km s⁻¹; Chen et al. 2011). The position of HD 131835 in the color–magnitude diagram (Figure 4) matches well the locus defined by UCL members, indicating that its isochrone age is consistent with that of the cluster. Our results confirm the previously proposed membership of HD 131835 in UCL; therefore, we adopted the age of the cluster, ∼16 Myr, for our star. With this age and our previous effective temperature and stellar luminosity estimates, we derived a stellar mass of 1.77 ± 0.08 M_⊙ by using solar-metallicity isochrones from Siess et al. (2000).

Zoom In Zoom Out Reset image size

SED modeling. In order to construct the SED of HD 131835, the previously mentioned PACS, SPIRE, and IRS data were further supplemented by infrared and submillimeter photometry from the literature (Table 3). The compiled SED is plotted in Figure 5. HD 131835 exhibits strong excess at far-IR wavelengths, but as the IRS spectrum shows, the SED of the source already starts to depart from photospheric emission at ∼7 μm, indicating the presence of warmer circumstellar dust as well. We used the ATLAS9 atmosphere model of HD 131835 (Section 3.1) to predict photospheric flux contribution at different mid- and far-IR wavelengths. The average accuracy of the predicted photospheric fluxes is around 5%.

Zoom In Zoom Out Reset image size

The excess emission of debris disks is generally well fitted by a single (modified) blackbody or a combination of two different temperature blackbody components (Chen et al. 2014; Kennedy & Wyatt 2014). This simple model can provide estimates of fundamental disk properties such as the characteristic dust temperature(s) and the fractional luminosity (f_d is the ratio of the luminosity of the dust emission to the bolometric luminosity of the host star).

We fitted the excess emission of HD 131835 by single- and two-temperature models. For the single-temperature fitting we used a modified blackbody:

$$\begin{eqnarray}&&{F}_{\nu ,\mathrm{exc}}={{AB}}_{\nu }({T}_{d}){X}_{\lambda }^{-1},\end{eqnarray} \tag{ 1 }$$

where ${F}_{\nu ,\mathrm{exc}}$ is the measured excess emission at ν, ${B}_{\nu }$ is the Planck function, T_d is the dust temperature, A is a scaling factor that is proportional to the solid angle of the emitting region, and X_λ = 1 if $\lambda \leqslant {\lambda }_{0}$ and ${X}_{\lambda }={\left(\frac{\lambda }{{\lambda }_{0}}\right)}^{\beta }$ if $\lambda \gt {\lambda }_{0}.$ In the two-temperature model we used a combination of a warmer simple blackbody (there is no sense in using a modified blackbody since the values of λ₀ and β cannot be constrained well from the SED) and a modified blackbody. We applied a Levenberg–Marquardt algorithm (Markwardt et al. 2009) to find the best-fitting model. An iterative way was used to compute and apply color corrections for the photometric data during the fitting process (see, e.g., Moór et al. 2006). We found that only the two-temperature model (plotted in Figure 5) could reasonably fit both the mid- and far-IR/submillimeter excess emission. This model provides dust temperatures of T_d,w = 176 ± 22 K and T_d,c = 71 ± 3 K for the warm and cold dust components, respectively. We obtained λ₀ = 136 ± 30 μm and β = 0.50 ± 0.16 for the cold belt. The warm component has a fractional luminosity of (8.2 ± 1.9) × 10⁻⁴, while the cold component has a fractional luminosity of (2.2 ± 0.2) × 10⁻³.

We estimated the dust mass in the disk based on the submillimeter brightness using the standard formula:

$$\begin{eqnarray*}&&{M}_{d}=\frac{{F}_{\nu ,\mathrm{excess}}{d}^{2}}{{B}_{\nu }({T}_{{\rm{d,c}}}){\kappa }_{\nu }},\end{eqnarray*}$$

where ${F}_{\nu ,{\rm{excess}}}$ is the measured excess at 500 μm (we used this data point because the longest-wavelength LABOCA measurement has very low signal-to-noise ratio), d is the distance to the source, ${\kappa }_{\nu }={\kappa }_{0}{\left(\frac{\nu }{{\nu }_{0}}\right)}^{\beta }$ is the mass absorption coefficient, and ${B}_{\nu }$ is the Planck function. By adopting ${\kappa }_{0}\;=$ 2 cm² g⁻¹ at ${\nu }_{0}=345\;{\rm{GHz}}$ (e.g., Nilsson et al. 2010) and taking β and T_d,c values from the previous SED modeling, we derived a dust mass of 0.47 ± 0.18 M_⊕ for the disk.

Spatial extent. In order to evaluate whether the disk is spatially extended, its profiles measured at 70, 100, and 160 μm were compared with appropriate PSFs of the PACS instrument. The PSF profiles were compiled using mini-scan-map observations of four stars (α Boo, α Tau, α Cet, β And) that do not exhibit infrared excesses and served as fiducial standards in the calibration of PACS photometry (Balog et al. 2014). For the data processing of these observations we used identical reduction steps as in the case of HD 131835, and finally the obtained PSFs were rotated to match the roll angle of the telescope at the time of observing our target. By fitting a two-dimensional Gaussian function to the PSF profiles, we derive averaged FWHMs of (570 ± 002) × (559 ± 003) and (686 ± 004) × (676 ± 003) at 70 and 100 μm, respectively, while for the target we obtained (649 ± 020) × (587 ± 019) and (774 ± 034) × (719 ± 032), implying that the disk is slightly elongated and marginally extended along its major axis at these wavelengths. At the longest wavelength the target's profile was consistent with that of the PSF measurements. In order to derive the characteristic size, inclination, and position angle of the disk, the 70 and 100 μm PACS images were fitted using a simple, nonphysical disk model grid in the same way as described in Moór et al. (2015). In this model we assumed that the dust emitting at these wavelengths is located in a narrow outer ring around the central star. The model has three free parameters, the average radius (R_avg), the position angle (PA), and the inclination (i) of the disk, while the width of the disk was fixed to 0.1R_avg following Booth et al. (2013). We used a Bayesian analysis in the selection of the best-fitting model. Our best solution has R_avg = 160 ± 20 AU, PA = 54° ± 8°, and $i={86}_{-10}^{+4}$° at 70 μm and ${R}_{{\rm{avg}}}=207\pm 40\;{\rm{AU}}$, PA = 39° ± 12°, and i = 66° ± 24° at 100 μm. The final disk parameters computed as a weighted average of those derived in the two bands are ${R}_{{\rm{avg}}}=169\pm 18\;{\rm{AU}}$, PA = 49 ± 7°, and $i={83}_{-9}^{+7}$°. It is worth noting that the measured profiles can also be fitted using a more extended disk model with a smaller inner radius and a larger outer radius. Such a model would provide very similar PA and inclination parameters as above, and the average disk radius would also not change significantly.

The disk around HD 131835 has also been successfully resolved at 11.7 and 18.3 μm using the Gemini South telescope (Hung et al. 2015), deriving a PA of ∼61° and an inclination of ∼74°. These parameters are broadly consistent with our results. By modeling the SED and the mid-IR images simultaneously, Hung et al. (2015) proposed a combination of an extended continuous power-law disk between 35 and 310 AU and two narrow rings at 105 ± 5 and 220 ± 40 AU stellocentric distances as the best-fitting solution. The three disk components are made of three different grain populations. All of the emitting grains are proposed to be hotter than blackbodies. The radius of the outer narrow ring (220 ± 40 AU), which is predominantly responsible for the disk emission at wavelengths >50 μm, is broadly consistent with the characteristic disk radius derived from the marginally resolved PACS images.

3.3.1. Basic Gas Mass Estimates

Observations of a single CO isotopologue (¹²CO in our case) cannot provide information on the optical depth of the radiation. Assuming that the measured CO emission of HD 131835 is optically thin, the mass of CO gas in the disk can be estimated as

$$\begin{eqnarray}&&{M}_{\mathrm{CO}}=\displaystyle \frac{4\pi {{md}}^{2}}{h{\nu }_{{ul}}{A}_{{ul}}}\displaystyle \frac{{S}_{{ul}}}{{x}_{u}},\end{eqnarray} \tag{ 2 }$$

where m is the mass of the CO molecule, d is the distance of the object, h is the Planck constant, ${\nu }_{{ul}}$ and A_ul are the rest frequency and the Einstein coefficient for the given transition between the u upper and l lower levels, S_ul is the observed integrated line flux, and x_u is the fractional population of the upper level. Apart from x_u, all of the parameters from the right side of the equation are known. In the case of local thermodynamic equilibrium (LTE), the level populations are thermalized and governed by the Boltzmann equation, allowing the computation of x_u if we know the gas temperature. In optically thin disks the gas kinetic temperature can significantly differ from the temperature of dust grains (see, e.g., Kamp & van Zadelhoff 2001); thus, the dust emission does not constrain the gas temperature. Since in the 2–1 and 4–3 transitions we have upper limits for the CO emission (Figure 1), from the line ratios we can provide only a poor constraint for the excitation temperature, yielding T_ex > 8 K. Taking T_ex = 20 K, the value measured in β Pic (Roberge et al. 2006), we derived a total CO gas mass of (4.4 ± 2.2) × 10⁻⁴ M_⊕ from our J = 3–2 line flux. In LTE T_ex and T_kin are equal. As Figure 6 demonstrates, the derived CO mass is not very sensitive to the excitation temperature. For a range between 8 and 250 K, M_CO varies between 3.7 × 10⁻⁴ and 2.3 × 10⁻³ M_⊕. In principle, using similar assumptions, we can derive upper limits for both O i and C ii masses; however, for these constituents—especially for O i—the estimates have stronger dependency on the excitation temperature (see Figure 6). Additionally, LTE excitation of the O i 63 μm line requires a dense medium, e.g., considering H₂ molecules as collisional partners, the critical density is typically ∼5 × 10⁵ cm⁻³ (Lequeux 2005). In a lower-density environment this line becomes subthermal with excitation temperatures far below what would be expected from LTE. Since these issues would make the mass estimate of O i unreliable, we focused on C ii, where the critical densities are lower (e.g., ∼3 × 10³ cm⁻³, assuming hydrogen molecules as collisional partners; Lequeux 2005). For T_ex > 8 K, we obtained ${M}_{{\rm{C}}\;{\rm{II}}}=1.4$ ${M}_{\oplus }$ as an upper mass limit for C ii.

Zoom In Zoom Out Reset image size

3.3.2. Simple Gas Disk Model

It is important to note, however, that our basic assumptions concerning the optically thin gas emission and the LTE could be questioned. For instance, in the disk around HD 21997 not only the ¹²CO but even the ¹³CO lines turned out to be optically thick (Kóspál et al. 2013), while in debris disks with low amounts of gas, radiative excitation can dominate over collisional processes, leading to subthermal non-LTE level populations (Matrà et al. 2015). To partly consider these caveats and to further constrain the gas disk fundamental properties, the measured line profile was modeled with a simple disk geometry using the LIME radiation transfer code (Brinch & Hogerheijde 2010).

We adopted a prescription for the structure of the disk that is frequently used to model passively irradiated protoplanetary disks. The radial surface density distribution was taken to be a power law. We assumed that the disk is vertically isothermal and in hydrostatic equilibrium, resulting in a Gaussian vertical density distribution, with radially dependent pressure scale height. Assuming Keplerian rotation around a 1.77 M_⊙ star and a radial temperature profile of ${r}^{-0.5},$ the ratio of the pressure scale height to the radius has a profile of ${r}^{-0.25}.$ Using these approximations, the disk can be described as

$$\begin{eqnarray}&&{n}_{{\rm{CO}}}(r,z)=\displaystyle \frac{{{\rm{\Sigma }}}_{{\rm{CO}}}}{\sqrt{2\pi }{H}_{{\rm{p}}}}\mathrm{exp}\left[-\displaystyle \frac{1}{2}{\left(\displaystyle \frac{z}{{H}_{{\rm{p}}}}\right)}^{2}\right],\end{eqnarray} \tag{ 3 }$$

$$\begin{eqnarray}&&\displaystyle \frac{{H}_{{\rm{p}}}}{r}=0.1{\left(\displaystyle \frac{r}{{R}_{{\rm{out}}}}\right)}^{-0.25},\end{eqnarray} \tag{ 4 }$$

$$\begin{eqnarray}&&{{\rm{\Sigma }}}_{{\rm{CO}}}(r)={{\rm{\Sigma }}}_{{\rm{CO,out}}}{\left(\displaystyle \frac{r}{{R}_{{\rm{out}}}}\right)}^{-1},\end{eqnarray} \tag{ 5 }$$

where R_out is the outer disk radius, ${{\rm{\Sigma }}}_{{\rm{CO,out}}}$ is the gas surface density at the outer disk radius in the midplane, and H_p is the scale height. Further parameters are the inner radius (R_in), the inclination (i) of the disk, and the gas temperature T_g(r). In the modeling we always assumed that the gas disk is coplanar with the dust disk and therefore its inclination was fixed to 74° (Section 3.2). The origin of gas in this system is not clarified yet; it could be either primordial or secondary. In the modeling, therefore, we used two different scenarios corresponding to these alternative hypotheses.

3.3.3. Primordial Gas Disk Model

3.3.4. Secondary Gas Disk Model

3.3.5. Search for Accretional Signatures

For the nondetected sources CO mass upper limits were estimated assuming optically thin emission and LTE. The excitation temperature was assumed to be 20 K, and the line flux upper limits were taken from Table 1. In those cases where both CO (3–2) and (4–3) line observations were available, we used the lower transition in the calculations. We note again that our basic assumptions might not be fulfilled for all of the studied systems. With the achieved sensitivity we cannot exclude that some of our targets harbor optically thick gas disks. In very tenuous gas disks the excitation could be subthermal, even leading to a very low excitation temperature of <10 K (e.g., Matrà et al. 2015). In both of these cases our upper limits would underestimate the possible total CO mass in the disk.

With the discovery of CO gas in HD 131835 the number of known gaseous debris disks around A-type stars is increased to seven. Tables 4 and 5 present the main stellar and disk properties. All of these systems are likely younger than 50 Myr;¹⁹  thus, they represent the very early phase of debris disk evolution. HD 131835 may be one of the youngest objects in this sample, and based on their nearly identical stellar properties, it can be considered as a comparably young sibling of β Pic.

Table 5.  Fundamental Disk Properties of Known Gaseous Debris Systems around A-type Stars

ID	${L}_{{\rm{IR}}}/$ ${L}_{{\rm{bol}}}$	M${}_{{\rm{dust}}}$	Gas Detection	MCO	Variable	Proposed
			CO/O i/C ii/C i		Abs. Lines	Origin
		(M⊕)		(M⊕)
49 Ceti	1.1 × 10−3 (14)	2.7 ± 0.4 × 10−1 (14)	Y/Y/Y/Y (22, 19, 20)	>2.5 × 10−4 (4)	Y (12)	Secondary (19)
HD 21997	5.7 × 10−4 (14)	1.6 ± 0.5 × 10−1 (14)	Y/N/N/N (13)	(4–8) × 10−2 (10)	N	Primordial (10)
β Pic	2.6 × 10−3 (15)	7.8 ± 0.8 × 10−2 (5)	Y/Y/Y/Y (5, 18, 2)	2.9 × 10−5 (5)	Y (1)	Secondary (7, 5)
HD 131835	3.0 × 10−3 (21)	4.7 ± 1.6 × 10−1 (21)	Y/N/N/N (21)	>5.2 × 10−4 (21)	N	Undecided
HD 32297	4.4 × 10−3 (3)	3.7 ± 1.1 × 10−1 (11)	N/N/Y/N (6)	<2.2 × 10−4 (13)	N	Secondary (6)
HD 172555	7.3 × 10−4 (16)	4.8 ± 0.06 × 10−4 (16)	N/N/N/Y (16)	<2.4 × 10−5 (13)	Y (9)	Secondary (8, 16)
η Tel	2.4 × 10−4 (3)	1.3 ± 0.08 × 10−2 (17)	N/N/Y/N (17)	<6.8 × 10−5 (21)	N	Secondary (17)

Note. Column (1): identification. Column (2): fractional luminosity of the disks. Column (3): dust mass. In the case of HD 32297 the mass was estimated from the 1.3 mm data published by Meeus et al. (2012), using the outline described in Section 3.2. Column (4): detection of CO, O i, C ii, and C i gas. Column (5): mass of CO gas. For HD 21997 and β Pic the mass estimates were taken from the literature or from Table 1. For the sake of homogeneity, in the case of 49 Ceti, HD 32297, and HD 172555 the previous mass estimates were recomputed using their CO (3–2) line flux data taken from the literature and assuming optically thin emission and T_ex = 20 K. Because of unknown ¹²CO optical depth, the mass estimates for 49 Ceti and HD 131835 are lower limits. Column (6): detection of variable optical absorption lines. Column (7): proposed origin of the observed gas. References for the literature data used in this table: (1) Beust et al. (1994); (2) Cataldi et al. (2014); (3) Chen et al. (2014); (4) Dent et al. (2005); (5) Dent et al. (2014); (6) Donaldson et al. (2013); (7) Fernández et al. (2006); (8) Johnson et al. (2012); (9) Kiefer et al. (2014); (10) Kóspál et al. (2013); (11) Meeus et al. (2012); (12) Montgomery & Welsh (2012); (13) Moór et al. (2011a); (14) Moór et al. (2015); (15) Rhee et al. (2007); (16) Riviere-Marichalar et al. (2012); (17) Riviere-Marichalar et al. (2014); (18) Roberge et al. (2006); (19) Roberge et al. (2013); (20) Roberge et al. (2014); (21) this work; (22) Zuckerman et al. (1995).

Download table as:  ASCIITypeset image

Figure 7 shows the integrated CO (3–2) fluxes and upper limits for these seven sources, normalized to 100 pc, and plotted against fractional luminosities. Partly based on our two surveys, additional debris disks around main-sequence stars located within 125 pc, as well as HD 141569, a nearby transitional disk that may represent a very final phase of protoplanetary disk evolution (Wyatt et al. 2015), are also displayed. Although there is no clear trend with the fractional luminosity, apart from η Tel all known gaseous debris disks have a fractional luminosity >5 × 10⁻⁴. Based on debris disk catalogs of Moór et al. (2006), Rhee et al. (2007), and Chen et al. (2014), we identified 11 A-type stars within 125 pc that harbor disks with L_IR/ L_bol above 5 × 10⁻⁴. All of these systems are younger than 50 Myr, and nearly all of them have already been observed in CO rotational transitions in different surveys (the two exceptions are HD 98363 and HD 143675). Among the observed nine disks, four harbor CO gas (49 Ceti, HD 21997, β Pic, and HD 131835), while in the case of HD 32297 and HD 172555 atomic gas tracers were detected. The gas components of HD 95086, HR 4796, and HD 110058 have remained undetected so far. This suggests that the most massive debris disks around young A-type stars quite commonly possess a detectable amount of gas. It is important to note that among these disks only HD 21997 and β Pic were observed with ALMA, and the sensitivity of most of the other observations did not allow the detection of even a β-Pic-like gas disk. From the four CO-bearing debris disks HD 131835 is the brightest; its normalized line flux is ∼17 times higher than that of the disk around β Pic, but falls below the normalized line flux of HD 141569.

Zoom In Zoom Out Reset image size

Within 125 pc there are many F-type stars with age of <50 Myr that possess debris disks with a fractional luminosity of >5 × 10⁻⁴ (Figure 7), i.e., having similar properties to that of the above-mentioned A-type sample. We observed 15 such systems in the framework of our two surveys, but none of them were detected; thus, based on the current data sets, the presence of CO molecules may be the characteristic of debris disks around young A-type stars. Rigliaco et al. (2015) detected H i lines in mid-infrared spectrum of eight young debris disks around F–K-type stars. The presence of these lines can be explained by either a low rate of gas accretion onto these stars or chromospheric activity. In the case of the former scenario these findings would indicate that some disks around F–K-type stars can also retain primordial gas for their debris phase.

Based on its infrared and submillimeter continuum data, the dust disk properties of HD 131835 are similar to those of the other young debris disks. Although the disk has a high total fractional luminosity of ∼3 × 10⁻³, this value is still less than the usually invoked threshold value of 0.01 that divides protoplanetary and debris disks. Recently, Wyatt et al. (2015) introduced a new classification scheme to differentiate between Herbig Ae disks and debris disks with A-type host stars. This method is based on the ratios of the observed to the stellar photospheric fluxes at 12 and 70 μm (R₁₂ and R₇₀). For debris disks R₁₂ and R₇₀ have to be lower than 3 and 2000, respectively. Based on color-corrected WISE 12 and PACS 70 μm data and predicted photospheric fluxes, we derived R₁₂ = 1.8 and R₇₀ = 955 for HD 131835, i.e., it can be classified as a debris disk. The majority of Herbig Ae systems show PAH emission. The absence of PAH features in the IRS spectrum of HD 131835 is also in accordance with its classification as a debris disk. These arguments suggest that the observed excess emission is due to dust grains that are probably second generation, produced via erosion of larger unseen bodies.

According to the current paradigm, the gas content of debris disks may also be second generation, derived from previously formed planetesimals. However, in young disks we cannot exclude the possibility that the evolutions of dust and gas were not parallel and the observed gas is predominantly composed of residual primordial material (of course, we never exclude that a part of the gas is produced from already-emerged icy bodies). Indeed, though most known gaseous debris disks are proposed to be rather secondary, one of the oldest among them, HD 21997, may likely harbor a hybrid disk with a significant amount of primordial gas (Kóspál et al. 2013). In the following, we investigate the origin of gas in the HD 131835 system. In this analysis we considered the β Pic and HD 21997 systems as references: the former may be the best established representantive of young debris disks with secondary gas (Fernández et al. 2006; Dent et al. 2014), while the latter is the sole known example for a hybrid disk. Furthermore, both systems have already been observed with ALMA (Kóspál et al. 2013; Dent et al. 2014), providing very detailed information on the spatial distribution of CO gas and dust, and their host stars—particularly β Pic—are quite similar to HD 131835.

In the absence of significant H₂ and dust content in a debris disk with secondary gas, only self-shielding can protect CO molecules againt the UV photons of stellar and interstellar radiation field. To judge the efficiency of self-shielding, we took our secondary gas disk model with a CO mass of 1.0 × 10⁻³ M_⊕ (i.e., the model with the minimum CO mass) and computed the shielding factors using the photodissociation model of Visser et al. (2009). The stellar UV flux was determined from the fitted ATLAS model (Section 3.1) and was found to be dominant at <45 AU. The contributions of interstellar radiation field were computed from Draine (1978). Our calculations yielded a CO lifetime of 40 yr at the directly illuminated inner edge of the disk and of 500 yr in the most efficiently self-shielded midplane regions. Assuming that the gas is continuously replenished and taking the CO mass of 1.0 × 10⁻³ M_⊕, we obtained a CO production rate of ≳2 × 10⁻⁶ M_⊕ yr⁻¹ or 1.2 × 10¹⁹ kg yr⁻¹. This is at least eight times higher than the gas production rate obtained for the β Pic disk from its ALMA observation (Dent et al. 2014). Assuming a CO mass abundance of 10% in planetesimals, it would require destruction at least of 2 × 10⁻⁵ M_⊕ yr⁻¹ of icy bodies. The CO gas can be released from icy grains/planetesimals in different processes. The temperature significantly exceeds the value needed for the sublimation of pure CO ice (∼20 K) everywhere in the disk. In such an environment CO ice is thought to be mainly present in deeper layers of icy planetesimals and as an admixture in the amorphous water-ice on the surface of grains and larger bodies. Thus, photodesorption or, in warmer regions (>110 K), sublimation of surface water-ices can produce CO gas as well. CO entrapped in a water-ice matrix can also be released via collisions between grains and planetesimals; moreover, fragmentation of larger bodies can lead to the excavation of CO ices persisting in deeper layers. Figure 8 shows the ratio of CO mass to the dust mass as a function of age for debris disks with a known gaseous component and for some selected disks (additional debris disks around A-type stars with ${f}_{{\rm{d}}}\gt 5\times {10}^{-4};$ Fomalhaut and HD 107146, two debris systems recently observed by ALMA; and HD 141569). Based on ALMA observations of Dent et al. (2014), the M_CO/M_dust ratio measured for β Pic disk is $3.6\times {10}^{-4}$ (note that they used an assumption of LTE, and if it is not fulfilled, the CO mass could be higher; see also Matrà et al. 2015). The spatially resolved images also revealed a clump at ∼85 AU from the star, in which an enhanced CO-to-dust mass ratio was measured. Sensitive ALMA observations are available only for a few debris disks. Upper limits for M_CO/M_dust at Fomalhaut (Matrà et al. 2015) and at HD 107146 (Ricci et al. 2015), however, hint that the ratio measured at β Pic is already quite high. Interestingly, the other three debris disks with CO gas show even higher M_CO/M_dust ratios. For HD 21997 we measured a ∼1000 times higher M_CO/M_dust ratio (note that here we used the optically thin ¹²C¹⁸O line to estimate the CO mass) than for β Pic. This object clearly differs from the others and likely harbors a significant fraction of primordial gas. The lower limits of CO-to-dust mass ratios for HD 131835 (in calculation of its M_CO/M_dust we used the minimum CO mass of 5.2 × 10⁻⁴ M_⊕ that we obtained in our models) and 49 Ceti are ≳2.5× higher than that of β Pic. In a secondary disk the M_CO/M_dust ratio depends both on the dust and gas production rates and on the removal time of the two elements. The observed higher CO-to-dust mass ratios thus can be explained, e.g., by a more effective self-shielding due to larger CO density, which would result in a higher CO lifetime, or by a higher fraction of volatile material in icy planetesimals/grains.

Zoom In Zoom Out Reset image size

Stellar and interstellar UV photons can permeate debris disks without any hindrance and ionize carbon atoms formed via photodissociation of CO molecules. The fraction of neutral and ionized carbon depends on the strength of the local radiation field. Owing to their longer lifetime in such environments, the amount of C i and C ii can significantly exceed that of CO molecules, e.g., in the disk of β Pic there is approximately 200 times more C ii gas than CO gas in terms of mass (Cataldi et al. 2014; Dent et al. 2014). Figure 9 displays the ${S}_{{\rm{C}}\;{\rm{II}}}/$ ${S}_{\mathrm{CO}(3-2)}$ line flux ratios for debris disks where C ii and/or CO gas has been detected. The ratio of β Pic is about 30 times higher than the upper limit we obtained for HD 131835. This may reflect a real difference in the mass ratio, or—as our modeling showed (Section 3.3)—may be a consequence of the low gas temperature in HD 131835. In the latter models the total gas mass and the gas-to-dust mass ratio would be rather high.

Zoom In Zoom Out Reset image size

In the secondary gas scenario, gas and dust are released at the same location from planetesimals; thus, the two components are expected to be co-located, just as we see at β Pic. Based on our current data and modeling (Section 3.3), we cannot determine the accurate location of gas; we just can say that the results do no exclude that the two components are co-located. It is worth noting, however, that our gas disk model coincides with that part of the dust disk model of Hung et al. (2015), which is composed by very hot grains. These grains give only a low fraction of the dust disk mass, and it is questionable whether they are accompanied by a large number of planetesimals that could be the source of observed CO gas.

All in all, in the secondary scenario, HD 131835 could be considered as of comparable age and—in terms of its disk—a more massive analog of the β Pic system. However, none of our current data exclude that HD 131835 in fact possesses a hybrid disk. Our modeling in Section 3.3 shows that the CO (3–2) emission may be optically thick; thus, we cannot exclude that the CO mass is significantly higher, reaching even the same level as in the case of HD 21997. The low ${S}_{{\rm{C}}\;{\rm{II}}}/$ ${S}_{\mathrm{CO}(3-2)}$ flux ratio can be explained in the secondary scenario, but as Figure 9 shows, among the currently known gaseous debris disks only HD 21997, the sole hybrid disk, exhibits similarly low ratio.

For a better clarification of the origin of gas, additional observations are necessary. By measuring the emission of less abundant ¹³CO and C¹⁸O isotopologues from the disk, a more reliable CO gas mass estimate could be obtained. Moreover, mapping and comparing the spatial distribution of gas and dust would allow us to investigate whether the two components are co-located, which is a prerequisite in a secondary gas scenario. These observations requires better sensitivity and spatial resolution than that of a single-dish radio telescope but could definitely be carried out using the Atacama Large Millimeter/submillimeter Array.