A Study of Millimeter Variability in FUor Objects

John Wendeborn; Catherine C. Espaillat; Enrique Macías; Orsolya Fehér; Á. Kóspál; Lee Hartmann; Zhaohuan Zhu; Michael M. Dunham; Marina Kounkel

doi:10.3847/1538-4357/ab9304

1. Introduction

FU Ori objects were originally identified as systems in star-forming regions that exhibit large outbursts at optical wavelengths (Ambartsumyan 1971; Herbig 1977; Hartmann & Kenyon 1996; Audard et al. 2014). These outbursts have been explained by the onset of rapid mass accretion, rising from quiescent states of 10⁻⁸–10⁻⁷ M_⊙ yr⁻¹ to as much as 10⁻⁴ M_⊙ yr⁻¹ (Hartmann & Kenyon 1996). In their high states, FUors can also exhibit low-amplitude photometric/spectroscopic variability on ∼weekly timescales (Siwak et al. 2018; Takagi et al. 2018). Such high accretion rates, which can last for decades, may explain how stars accrete much of their mass, up to 10⁻² M_⊙ for a single outburst (Hartmann & Kenyon 1996).

The causes of these outbursts are not clearly understood. Several explanations, including perturbations due to nearby companions, thermal instabilities, and gravitational instabilities, have all been proposed (Bonnell & Bastien 1992; Hartmann & Kenyon 1996; Clarke et al. 2005; Reipurth et al. 2007). Companion perturbations, while an attractive solution for particular systems with specific orbital parameters, fail to explain many FUor events for isolated stars. Thermal and/or gravitational instabilities are more attractive alternatives for the underlying mechanism(s) behind FUor outbursts.

Thermal instabilities can occur when high disk opacities trap heat, leading to a runaway temperature increase within the disk. This heightened temperature then increases the effective disk viscosity, which in turn leads to high accretion rates (Bell & Lin 1994). This is considered necessary to explain the fast rise times of FU Ori and V1057 Cyg during their initial outbursts (Audard et al. 2014). However, thermal instabilities depend significantly on the disk viscosity and temperature, which are generally difficult to determine (Vorobyov & Basu 2005).

Gravitational instabilities are thought to occur when mass infalling from the surrounding envelope builds up in the disk, becoming gravitationally unstable. These instabilities can manifest in several ways, including magnetorotational instabilities (e.g., Armitage et al. 2001; Zhu et al. 2010) and/or clumps of material that then accrete onto the host star (e.g., Vorobyov & Basu 2005). If the envelope continues to dump mass into the disk over time, outbursts can be repetitive, something thought to be generally true of FUor objects (Hartmann & Kenyon 1996).

Despite the above caveats, it is possible that both thermal and gravitational instabilities work in concert to produce FUor events. However, both explanations require disk masses large enough to sustain high mass accretion rates (10⁻⁴ M_⊙ yr⁻¹) for decades, even hundreds of years, and for the gravitational instability to be triggered (M_disk/M_* ≳ 0.1; see Hartmann & Kenyon 1996; Liu et al. 2016; Cieza et al. 2018). Some measurements of FUor disk masses have called into question whether disks are, in general, massive enough for these instabilities to fully explain the observed outburst. For example, Dunham et al. (2012) and Kóspál et al. (2016) found upper limits on the disk mass of the FUor object HBC 722 to be 0.02 M_⊙ and 0.01 M_⊙, respectively, likely too small for gravitational instabilities to explain its outburst (Audard et al. 2014). Fehér et al. (2017, hereafter F17) also found that the disk masses of several FUor objects are quite low (e.g., 0.04 M_⊙ and <0.05 M_⊙ for V1515 Cyg and V710 Cas, respectively).

In order to gauge the viability of gravitational and/or thermal instabilities to explain FUor outbursts, accurate estimates of disk masses need to be made. This is not simple given that recent observations of FUor objects have shown that there may be up to 25%–60% variability at 1.3 mm (Liu et al. 2018). Liu et al. (2018) pointed to two possible explanations for the variability, including the perturbations of the thermal or density structure in the disk, as well as increased irradiation from the host star.

Coordinated millimeter and optical/NIR observations can provide insight into the underlying mechanism behind the observed millimeter variability. Most of the millimeter emission from FUor objects is thought to originate from the outermost regions of the disk, and traces cool, optically thin, millimeter-sized dust grains. The shorter optical/NIR emission originates from the innermost, hottest regions of the disk. A simultaneous change at millimeter and optical/NIR wavelengths would indicate that such variation is likely due to temperature changes in the disk, and hence thermal instabilities. This would further imply that future millimeter observations of FUor objects could, in general, benefit from coordinated optical/IR observations so as to better constrain the disk's properties (Reipurth et al. 2007). On the other hand, observing solely a change in millimeter flux density would indicate that the millimeter emission is disconnected from perturbations in the thermal structure of the disk, perhaps pointing toward gravitational instabilities or other alternative mechanisms.

In the following, we present millimeter observations of six FUor objects and coordinated millimeter and optical observations for a subset of three of these objects in order to determine the manner in which they vary over time. Details of our observations, as well as our data reduction techniques can be found in Section 2. Section 3 describes our analysis of these observations and the results. We discuss our findings in Section 4.

2. Observations and Data Reduction

2.1. Sample

Our entire sample consists of six known FUor objects (V1735 Cyg, V2494 Cyg, V2495 Cyg, V1057 Cyg, V1515 Cyg, and V733 Cep). All six targets had one observing run at 2.7 mm in 2017 May–June. In this paper, we focus primarily on V1735 Cyg, V2494 Cyg, and V2495 Cyg, which were selected for two follow-up observing runs including both millimeter and optical observations in 2018 June and August. V2494 Cyg and V2495 Cyg were chosen for follow-up observations because Liu et al. (2018) found tentative evidence for millimeter variability in these objects at 1.3 mm. V1735 Cyg was chosen because our 2017 data displayed variability relative to 2014 observations taken by F17. Our 2017 observations of V1057 Cyg, V1515 Cyg, and V733 Cep were consistent with 2014 flux densities reported by F17, so we elected not to carry out follow-up observations of these targets.

2.2. Millimeter Observations

Our millimeter observations were carried out using the NOrthern Extended Millimeter Array (NOEMA) in Plateau de Bure, France in 2017 May and/or June for all six objects in our sample. Three objects in our sample (V1735 Cyg, V2494 Cyg, and V2495 Cyg) were observed twice more in 2018 June and August. This gave us baselines of one year and two months. We also used the 2014 April observations from F17 of V1735 Cyg, which gives us a longer, three-year baseline for that object. Observing conditions were generally good throughout each observation, with stable precipitable water vapor typically below 10 mm. Conditions during our 2017 observations were worse, with precipitable water vapor levels sometimes above 10 mm. More details of all observations can be found in Table 1.

Table 1. NOEMA Observation Log

Object	R.A./Decl.	Date	Array Config.	Central Frequency	Distance^a
	(J2000)		(No. of Antennas)	(GHz)	(pc)
V1735 Cyg	21:47:20.66^a	2014 Apr 5	6Cq (6)	109.918	616^b
	+47:32:03.6	2017 Jun 11/12/14	8D-E10 (7)/8D (8)/8D-W12N13 (6)	108.502
		2018 Jun 9/10	8ant-Special (8)	106.744
		2018 Aug 14/15	8D-W05 (7)	106.744

V2494 Cyg	20:58:21.09^a	2017 May 31/Jun 11	8D-N09 (7)/8D-E10 (7)	108.502	800^b
	+52:29:27.7	2018 Jun 10	8ant-Special (8)	106.744
		2018 Aug 14/15	8D-W05 (7)	106.744

V2495 Cyg	21:00:25.38^c	2017 Jun 10	8D (8)	108.502	800^d
	+52:30:15.5	2018 Jun 10	8ant-Special (8)	106.744
		2018 Aug 14/15	8D-W05 (7)	106.744

V1057 Cyg	20:58:53.73^a	2017 May 31/Jun 2	8D-N09 (7)	108.502	898^c
	+44:15:28.38

V1515 Cyg	20:23:48.02^a	2017 Jun 3/4/7	8D-N09 (7)	108.502	981^b
	+42:12:25.78

V733 Cep	22:53:33.26^a	2017 Jun 5/8	8D-N09 (7)/8D (8)	108.502	661^b
	+62:32:23.63

Notes.

^aGaia Collaboration (2018). ^bBailer-Jones et al. (2018). ^cCutri (2012). ^dMagakian et al. (2013).

Download table as: ASCII Typeset image

All continuum observations reported here were centered near 108 GHz (∼2.8 mm), while the 2014 observations of V1735 Cyg from F17 were centered near 110 GHz (∼2.7 mm). The 2017 observations used NOEMA's WideX correlator tuned from 106.743 to 110.261 GHz, which covered the ¹³CO(1–0) (110.201 GHz) and C¹⁸O(1–0) (109.782 GHz) lines (see Table 2) at 78.118 kHz resolution. Our 2018 observations used the new PolyFix correlator tuned from 87.399 to 95.117 GHz (lower sideband) and 102.886 to 110.603 GHz (upper sideband). The two frequency ranges allowed us to measure the continuum emission at two different wavelengths, 2.81 and 3.29 mm. These frequencies covered the lines mentioned above, as well as the ¹³CS(2–1) (92.494 GHz) and HCO⁺(1–0) (89.189 GHz) lines, at a spectral resolution of 62.495 kHz.

Table 2. Summary of Observed Molecular Lines

Molecule	Line Frequency	Transition
¹³CO	110.201 GHz	J = 1-0
C¹⁸O	109.782 GHz	J = 1-0
¹³CS	92.494 GHz	J = 2-1
HCO⁺	89.189 GHz	J = 1-0

Download table as: ASCII Typeset image

Data reduction was done in GILDAS¹² using the NOEMA pipeline in the CLIC package. Minimal flagging was required to remove spurious data. The calibrated data were then transferred to CASA 5.3.0 (McMullin et al. 2007) for cleaning and further analysis. The flux calibrator MWC349 was used for all observations. J2037+511 was the phase calibrator for V2494 Cyg and V2495 Cyg, while J2201+508 was used for V1735 Cyg. The flux calibration source MWC 349 is time variable at 3 mm by <10%. Taking this and other factors into account, the nominal absolute flux uncertainty of NOEMA at 2.7 mm is estimated to be 10%.

We note that we cannot exclude the possibility that unaccounted instrumentation errors also affected the flux calibration (e.g., antenna pointing errors). We expect any other possible effect to have a lower impact than the nominal 10% absolute flux calibration uncertainty. However, quantifying their impact is extremely uncertain, and could lead to an underestimate of all of our cited flux calibration uncertainties.

Continuum data of V1735 Cyg, V2494 Cyg, and V2495 Cyg were first restricted to a uv range of 30–65 kλ in order to mitigate the effects of different uv-coverages and ensure that the amount of missing flux is the same between epochs. This also has the effect of removing extended emission, thereby ensuring that the measured flux densities are solely of the central, compact source. These continuum data were imaged using the clean task with natural weighting and then convolved to a 4'' × 4'' beam. After cleaning, but before uv-restriction/beam convolution, our angular resolution was about 3''–4''. After restricting the uv range to 30–65 kλ, our resolutions improved to about 2''–3''. And after beam convolution, the resolution was constant at 4''.

The resulting continuum images show compact emission for all the sources. The morphology of this continuum emission remained unchanged. No source was resolved, either fully or marginally, at 2.7 mm, before or after uv restriction or beam convolution. Continuum data of V1057 Cyg, V1515 Cyg, and V733 Cep were cleaned in the same manner, but were not uv-restricted nor convolved to a 4'' × 4'' beam. The line data were also processed using the clean task and natural weighting, but were corrected using imcontsub to remove continuum emission. The spectral resolution for all data cubes was about 0.25 km s⁻¹. Continuum rms values were obtained from emission-free regions, whereas line rms values were determined using emission-free channels. Note that the continuum rms values therefore include the effects of thermal and phase atmospheric noise as well as some contribution from the imaging reconstruction process due to the limited uv-coverage, whereas the line rms values do not include the latter effect.

2.3. Optical Observations

Optical observations of V1735 Cyg, V2494 Cyg, and V2495 Cyg were carried out using the Lowell Discovery Telescope's (LDT) Large Monolithic Imager (LMI) in Happy Jack, AZ (Bida et al. 2014), in 2018 June and August. These observations were coordinated with our 2018 NOEMA observations. Our June millimeter and optical observations were separated by about 24 hr, whereas our August observations were separated by about 36 hr. Details of these observations can be found in Table 6.

All observations were performed in photometric conditions using three optical filters: Johnson V, R, and I, with central wavelengths of 551 nm, 658 nm, and 806 nm respectively (Johnson et al. 1966). The FOV of the LMI is 12 farcm 5 × 12 farcm 5 (0 farcs 12 per unbinned pixel). We used 2 × 2 pixel binning for these observations (0 farcs 24 per pixel). The bias, flat-field, and overscan calibration of the CCD images were performed using a custom Python routine utilizing the Astropy package (Astropy Collaboration et al. 2018). The photometric calibration of all images was carried out interactively also using a custom Python routine utilizing the Astropy package. Exposure times varied among the targets, dependent on precise seeing conditions and object brightness. V2495 Cyg was too faint to detect in any of our observations. We were therefore unable to extract optical magnitudes for this object.

Using the standard stars SA92 253 and GD 246 (Landolt 2009), we obtained V, R, and I magnitudes for our targets in 2018 August. To obtain magnitudes for our 2018 June observations, we utilized differential photometry. This involved normalizing the flux density of our target to that of the total brightness of several background stars, which fluctuated by less than 3% throughout the observations. Then we directly compared the 2018 June flux densities to those of 2018 August. Finally, we used that ratio to obtain an optical magnitude for the 2018 June epoch. Magnitudes of our sample are listed in Table 6.

3. Results and Analysis

In the following, we present millimeter flux densities for all six FUors in our sample. We then search for variability in the 2.7 mm flux densities; ¹³CO, C¹⁸O, ¹³CS, and HCO⁺ molecular line fluxes; and optical magnitudes of V1735 Cyg, V2494 Cyg, and V2495 Cyg. We then describe the variability that is seen in V1735 Cyg for the 2.7 mm continuum.

3.1. Millimeter Continuum

Continuum flux densities were measured by 2D Gaussian fits within a 5 farcs 66 × 5 farcs 66 circular region, twice the convolved beam area. These results can be found in Table 3. To account for the frequency discrepancy, we adjusted our 2014 and 2017 measured flux densities using a spectral index of 2.5. The uncertainties stated in Table 3 were obtained via the root sum square of the 10% absolute flux calibration uncertainty of NOEMA (15% in the case of 2017 June observations) and the Gaussian fit uncertainties (typically 1%–5%).

Table 3. Millimeter Continuum Flux Densities

		${I}_{\nu ,\mathrm{Peak}}$	${S}_{\nu ,\mathrm{Gauss}}$	Rms
Object	Date
		(mJy/beam)	(mJy)	(mJy/beam)
V1735 Cyg	2014 Apr	1.63 ± 0.16^a	1.33 ± 0.16^a	0.10
	2017 Jun	2.27 ± 0.34^a	2.44 ± 0.39^a	0.07
	2018 Jun	1.94 ± 0..19	1.84 ± 0.20	0.08
	2018 Aug	1.68 ± 0.17	1.60 ± 0.19	0.09

V2494 Cyg	2017 May/Jun	16.7 ± 2.5^a	16.6 ± 2.5^a	0.16
	2018 Jun	16.0 ± 1.6	15.8 ± 1.6	0.07
	2018 Aug	13.8 ± 1.4	12.5 ± 1.4	0.40

V2495 Cyg	2017 Jun	14.6 ± 2.2^a	13.9 ± 2.2^a	0.17
	2018 Jun	14.0 ± 1.4	13.8 ± 1.4	0.07
	2018 Aug	12.9 ± 1.3	12.2 ± 1.2	0.24

V1057 Cyg	2017 May/Jun	3.97 ± 0.60	5.40 ± 0.81	0.08

V1515 Cyg	2017 Jun	0.66 ± 0.10	0.84 ± 0.13	0.04

V733 Cep	2017 Jun	0.54 ± 0.08	1.50 ± 0.23	0.05

Notes. ${I}_{\nu ,\mathrm{Peak}}$ and ${S}_{\nu ,\mathrm{Gauss}}$ are obtained from Gaussian fitting and the rms is obtained from an emission-free region.

^aCorrected for frequency discrepancy, using a spectral index of 2.5.

Download table as: ASCII Typeset image

Figure 1 shows the continuum maps of all six targets. While flux densities did vary somewhat, intensity distributions remained unchanged, thus we only show one epoch. The maps of V1735 Cyg, V2494 Cyg, and V2495 Cyg were made using a restricted uv coverage, whereas the rest were not. We note that despite steps taken to mitigate the differing uv coverages, a portion of the variability reported may still be due to remaining differences.

Also of note, our measured flux density for V1057 Cyg (5.4 ± 0.1 mJy) agrees with that of F17 (4.9 ± 0.2 mJy). V1515 Cyg and V733 Cep were weakly detected in F17, where F17 estimates peak intensities of 0.18 ± 0.03 and 0.38 ± 0.10 mJy/beam, respectively. F17 could not, however, estimate flux densities. The F17 peak intensities are somewhat lower than what we report here for V1515 Cyg and V733 Cep, 0.66 ± 0.10 and 0.54 ± 0.08 mJy/beam, respectively.

For V1735 Cyg, V2494 Cyg, and V2495 Cyg, we show the continuum flux densities and optical magnitudes for the two epochs in 2018 in Figures 2, 3, and 4, respectively. The only object to display any millimeter variability in our sample, V1735 Cyg, exhibited an ∼80% increase in flux density from 2014 to 2017. This flux density increase falls outside our stated uncertainties, thus, we conclude that the observed variability is intrinsic to the source.

**Figure 2.** Top: V1735 Cyg millimeter continuum flux density vs. time. Error bars are the uncertainties and are listed in Table 3. Bottom: V1735 Cyg VRI magnitudes vs. time. Error bars are roughly the size of the points.
Download figure:
Standard image High-resolution image

**Figure 3.** Top: V2494 Cyg millimeter continuum flux density vs. time. Error bars are the uncertainties and are listed in Table 3. Bottom: V2494 Cyg VRI magnitudes vs. time. Error bars are roughly the size of the points.
Download figure:
Standard image High-resolution image

**Figure 4.** V2495 Cyg millimeter continuum flux density vs. time. Error bars are the uncertainties and are listed in Table 3.
Download figure:
Standard image High-resolution image

We also see that following its rise from 2014 to 2017, V1735 Cyg possibly dimmed in 2018 June and then again in 2018 August. This may be a sign that it is returning to some quiescent state, from some burstlike, heightened state. However, this downward trend from 2017 to 2018 is within or close to the flux uncertainties of NOEMA, and is also seen in both V2494 Cyg and V2495 Cyg. More observations are needed to confirm if V1735 Cyg's flux density at 3 mm has truly decreased since 2017 June. No other objects in our sample, over any time period, show signs of millimeter variability outside of the measurement and flux calibration uncertainties.

3.2. Molecular Lines

¹³CO, C¹⁸O, ¹³CS, and HCO⁺ line fluxes for V1735 Cyg, V2494 Cyg, and V2495 Cyg were extracted from velocity-integrated spectral cubes after continuum subtraction and cleaning. The velocity range used for integration varied per object and per species, but was always centered in the frames that contained emission. In all cases, the same 5 farcs 66 × 5 farcs 66 circular aperture (centered at the primary source of emission) was used to measure the flux. Given the extended and asymmetric morphology of the line emission, we chose not to use 2D Gaussian fits to obtain line fluxes. These results can be found in Tables 4 and 5. We note that for the line fluxes, due to the variable uv coverages, lower SNR, extension of the emission, and possible contamination from the surrounding envelope, we assume uncertainties of 25%. In one case (V2495 Cyg) ¹³CS emission was not detected on or near the target's location, but further from the target at about ∼15'' away. We measure 0.039 and 0.051 Jy km s⁻¹ in 2018 June and August, respectively.

Table 4. ¹³CO, C¹⁸O Fluxes

		${I}_{\nu ,13\mathrm{CO}}$	Rms_13CO	Δv_13CO	${I}_{\nu ,{\rm{C}}18{\rm{O}}}$	Rms_C18O	Δv_C18O
Object	Date
		(Jy km s⁻¹)	(Jy/beam km s⁻¹)	(km s⁻¹)	(Jy km s⁻¹)	(Jy/beam km s⁻¹)	(km s⁻¹)
V1735 Cyg	2014 Apr	0.7 ± 0.2	0.12	−4.43, +5.87	0.21 ± 0.05	0.03	−2.73, +1.52
	2017 Jun	0.7 ± 0.2	0.12	−4.50, +6.00	0.22 ± 0.06	0.03	−2.75, +1.50
	2018 Jun	1.0 ± 0.3	0.11	−4.50, +6.00	0.21 ± 0.05	0.04	−2.75, +1.50
	2018 Aug	1.3 ± 0.3	0.17	−4.50, +6.00	0.39 ± 0.10	0.06	−2.75, +1.50

V2494 Cyg	2017 May/Jun	1.2 ± 0.3	0.09	−2.75, +2.25	0.40 ± 0.10	0.05	−3.00, +1.50
	2018 Jun	1.2 ± 0.3	0.08	−2.75, +2.25	0.40 ± 0.10	0.04	−3.00, +1.50
	2018 Aug	1.2 ± 0.3	0.11	−2.75, +2.25	0.38 ± 0.10	0.06	−3.00, +1.50

V2495 Cyg	2017 Jun	0.19 ± 0.05	0.03	−2.00, +3.75	0.22 ± 0.06	0.03	−2.00, +2.50
	2018 Jun	0.23 ± 0.06	0.05	−2.00, +3.75	0.28 ± 0.07	0.04	−2.00, +2.50
	2018 Aug	0.31 ± 0.08	0.03	−2.00, +3.75	0.27 ± 0.07	0.05	−2.00, +2.50

V1057 Cyg	2017 May/Jun	2.2 ± 0.6	0.15	−4.32, +1.84	0.49 ± 0.12	0.03	−2.12, +1.08

V1515 Cyg^a	2017 Jun	⋯	⋯	⋯	⋯	⋯	⋯

V733 Cep^a	2017 Jun	⋯	⋯	⋯	⋯	⋯	⋯

Notes. I_ν is the flux obtained from a 5 farcs 66 × 5 farcs 66 aperture on the object location. Rms_13CO and rms_C18O are the background rms obtained from emission-free regions. Δv is the integrated velocity range.

^aNo data, due to poor quality.

Download table as: ASCII Typeset image

¹³CO and C¹⁸O fluxes (Table 4) were generally consistent for all objects across all epochs. The only exception may be the C¹⁸O emission of V1735 Cyg from 2018 June to August (see Figure 5). The flux appears to have risen by about 86%, though we emphasize that the line fluxes are highly uncertain because of the lack of short baselines to recover the extended emission. Additionally, the slightly different uv coverages between epochs also result in artificial differences in the morphology of the extended line emission. Regardless, these differences have a relatively small (yet hard to quantify) effect on our flux measurements since we focus only on the compact line emission at the position of each object. This is partly shown by the fact that the 2014 observations (which included IRAM 30 m observations to cover short uv spacings) display similar fluxes to our observations.

**Figure 5.** V1735 Cyg C¹⁸O flux vs. time. Error bars are the uncertainties and are listed in Table 4.
Download figure:
Standard image High-resolution image

¹³CS and HCO⁺ was observed in V1735 Cyg, V2494 Cyg, and V2495 Cyg in 2018 June and August (Table 5). These lines show no signs of variability, either in flux or spatial morphology. We did not detect ¹³CS or HCO⁺ in V2495 Cyg. The emission morphology for ¹³CS and HCO⁺ for all objects was consistent throughout 2018, thus we display only one epoch in Figures 6 and 7.

**Figure 6.** Moment-0 molecular line maps of V1735 Cyg. Top: ¹³CO. Middle: C¹⁸O. Bottom: ¹³CS (left) and HCO⁺ (right). Solid (dashed) gray contours denote positive (negative) 3, 6, 12, 24, 48, 96, and 192σ levels. σ for each epoch is equivalent to the rms of each image, which can be found in Table 4. The central "+" denotes the target's location (see Table 1). The "×" denotes the location of V1735 Cyg SM1 (Harvey et al. 2008). Morphological differences between 2014 April and other epochs is due to differing uv coverages. Note that we do not include the ¹³CS and HCO⁺ line maps from 2018 August since they are very similar to those of 2018 June shown here.
Download figure:
Standard image High-resolution image

**Figure 7.** Moment-0 molecular line maps of V2494 Cyg (top) and V2495 Cyg (bottom). From left to right: ¹³CO, C¹⁸O, ¹³CS, and HCO⁺. Solid (dashed) gray contours denote positive (negative) 3, 6, 12, 24, 48, 96, and 192σ levels. σ for each epoch is equivalent to the rms of each image, which can be found in Table 4. The central "+" denotes the source's location (see Table 1). Note that here we only show ¹³CO and C¹⁸O line maps from 2018 June since the maps from 2017 June and 2018 August are very similar. Likewise, we only show line maps of ¹³CS and HCO⁺ from 2018 June since the maps from 2018 August are very similar.
Download figure:
Standard image High-resolution image

Table 5. ¹³CS, HCO⁺ Fluxes

		${I}_{\nu ,13\mathrm{CS}}$	Rms_13CS	Δv_13CS	${I}_{\nu ,{\mathrm{HCO}}^{+}}$	Rms ${}_{{\mathrm{HCO}}^{+}}$	Δv ${}_{\mathrm{HCO}+}$
Object	Date
		(Jy km s⁻¹)	(Jy/beam km s⁻¹)	(km s⁻¹)	(Jy km s⁻¹)	(Jy/beam km s⁻¹)	(km s⁻¹)
V1735 Cyg	2018 Jun	0.033 ± 0.009	0.008	−1.00, +0.50	1.14 ± 0.29	0.07	−1.00, +7.50
	2018 Aug	0.040 ± 0.010	0.014	−1.00, +0.50	1.35 ± 0.34	0.11	−1.00, +7.50

V2494 Cyg	2018 May/Jun	⋯	0.004	−1.50, −0.50	0.18 ± 0.05	0.01	−2.25, +0.25
	2018 Aug	⋯	0.007	−1.50, −0.50	0.17 ± 0.04	0.03	−2.25, +0.25

V2495 Cyg	2018 Jun	⋯	0.01	−1.25, +3.25	⋯	0.02	−2.00, +1.00
	2018 Aug	⋯	0.03	−1.25, +3.25	⋯	0.03	−2.00, +1.00

Note. Our 2017 correlator setup did not include these lines. Therefore, the 2017 observations are not included in this table. I_ν is the flux obtained from a 5 farcs 66 × 5 farcs 66 aperture on the object location. Rms_13CS and rms ${}_{\mathrm{HCO}+}$ are the background rms obtained from emission-free regions. Δv is the integrated velocity range. We do not detect ¹³CS emission from V2494 Cyg or V2495 Cyg. We do not detect HCO⁺ emission from V2495 Cyg.

Download table as: ASCII Typeset image

3.3. LDT-LMI Photometry

Optical photometry taken in 2018 June and August for V1735 Cyg and V2494 Cyg do not show variability (Table 6) and are consistent with previous observations. Our measurements of V1735 Cyg generally agree with those of Peneva et al. (2009) from 2003 to 2009. They measured V ∼ 18.9 and R ∼ 16.6, but found I ∼ 13.8, about half a magnitude brighter than reported here. Our measurements of V2494 Cyg agree with those of Magakian et al. (2013) from 2003 to 2010 in R and I. They found R ∼ 16.4 and I ∼ 14.7, but did not report V. We note that V2495 Cyg was observed, but not detected (see Section 2.3).

Table 6. LDT Photometry

Object	Date	Start Time (UT)	V	R	I
V1735 Cyg	2018 Jun 11	10:44	19.02 ± 0.09	16.61 ± 0.06	14.24 ± 0.04
	2018 Jun 12	10:45	19.07 ± 0.09	16.61 ± 0.06	14.24 ± 0.04
	2018 Aug 13	11:10	19.13 ± 0.09	16.62 ± 0.07	14.28 ± 0.05

V2494 Cyg	2018 Jun 11	10:51	18.47 ± 0.09	16.65 ± 0.07	14.79 ± 0.04
	2018 Jun 12	10:53	18.47 ± 0.09	16.63 ± 0.06	14.80 ± 0.04
	2018 Aug 13	11:30	18.38 ± 0.09	16.55 ± 0.07	14.73 ± 0.05

Download table as: ASCII Typeset image

4. Discussion

Prior to our observations, V2494 Cyg and V2495 Cyg were the only two FUor objects thought to be variable at millimeter wavelengths, displaying 1.3 mm flux density changes of ∼25%–60% on a timescale of about one year (Liu et al. 2018). Here we report that V1735 Cyg has also exhibited variability in the millimeter, but at 2.7 mm, and over a timescale of ∼3 yr, from 2014 to 2017. We discuss here possible underlying mechanisms for this variability.

4.1. Variable Disk Heating

In FUors, the inner disk is significantly heated by viscous heating from the accretion process and produces strong optical/IR emission (Hartmann & Kenyon 1996) and possibly millimeter emission as well (Takami et al. 2019). This hot inner disk irradiates the outer disk. Therefore, changes in the temperature of the inner disk may lead to changes in the heating of the outer disk, which we can trace with millimeter emission.

If temperature changes in the inner disk were the cause of the millimeter variability we see in V1735 Cyg, we would also expect to see a corresponding increase in its optical and/or IR flux as well. To the best of our knowledge, there is no optical data of V1735 Cyg taken close in time to the 2014 NOEMA data and so we cannot test this using our 2018 optical data. However, archival Wide-field Infrared Survey Explorer (WISE) data of V1735 Cyg exists at 3.4 and 4.6 μm from 2014 through 2020 (Figure 8), with data from 2014 June (∼2 months after the 2014 April NOEMA observations) and 2017 June (taken within a week of the 2017 June NOEMA observations). The WISE photometry displays no significant variability, indicating that the disk irradiation may have remained relatively constant during that time, and therefore that the millimeter variability is not tied to disk temperature changes. We can also rule out a change in disk temperature being responsible for the millimeter variability given that the ∼doubling of the millimeter flux in V1735 Cyg would imply an equivalent ∼doubling in disk irradiation, which is unlikely.

**Figure 8.** WISE photometry of V1735 Cyg from 2014 to 2019. Red circles are Band 1 (3.4 μm). Blue squares are Band 2 (4.6 μm). Dashed gray bars indicate dates of NOEMA observations.
Download figure:
Standard image High-resolution image

4.2. Gas and Dust Buildup in the Disk

Because we see solely millimeter variability and no optical changes, one may speculate that this may be evidence that material is building up in the disk from the envelope. Using the equation

$\begin{eqnarray}&&{M}_{\mathrm{cont}}=\displaystyle \frac{{{gS}}_{\nu }{d}^{2}}{{\kappa }_{\nu }{B}_{\nu }(T)}\end{eqnarray} \tag{ 1 }$

(which assumes an optically thin disk), where M_cont is the continuum mass, g = 100 is the gas-to-dust ratio, S_ν is the measured flux density at 2.7 mm, d is the distance, κ_ν = 0.2 cm² g⁻¹ is the dust opacity coefficient at 2.7 mm, and B_ν(T) is the Planck function for a blackbody with a temperature of T = 30 K, we find that the disk mass of V1735 Cyg (at a distance d = 616 pc; Bailer-Jones et al. 2018) must have increased from 0.13 to 0.21 M_⊙ from 2014 to 2017. We note that this is consistent with a previous disk mass (0.20 M_⊙) estimated with SED modeling (Gramajo et al. 2014). Our measured disk mass change would correspond to a mass infall rate of 0.027 M_⊙ yr⁻¹ from the envelope, which is highly unlikely (Ohtani & Tsuribe 2013; White et al. 2019). Therefore, given the degree to which the continuum flux density changes, an unrealistic rate of mass infall would be necessary to account for the magnitude of the millimeter variability seen in V1735 Cyg. In addition, the (likely optically thin) C¹⁸O emission was relatively constant from 2014 to 2017, implying that no C¹⁸O has built up during that time. Thus, material buildup does not seem to be the source of the variability of V1735 Cyg at 2.7 mm.

4.3. Free–Free Emission

One other potential source of millimeter variability is changes in the free–free emission of the system. Free–free emission can be identified by analysis of the spectral index, α, of the millimeter emission. The more significant the free–free emission, the shallower the spectral index, down to −0.1 to 0.6 for purely free–free (Reynolds 1986). Scattering in an optically thick disk can act to lower the spectral index as well, though this effect is generally strongest at the innermost regions of the disk (Liu 2019; Zhu et al. 2019).

The change in α of V1735 Cyg can be measured using existing data from 2013 June, 2014 April, 2018 June, and 2018 August. Liu et al. (2018) weakly detected V1735 Cyg at 1.3 mm in 2013 June. Given their 3σ upper limits, and using flux density estimates from 2014 April (F17), Liu et al. (2018) determined an upper limit on α of 1.7–2.0. Using the upper (106.7 GHz) and lower (91.3 GHz) sidebands described in Section 2, we are able to determine spectral indices of our 2018 June and August observations. We find tentative evidence of shallower slopes than Liu et al. (2018), α = 1.4 ± 0.4 in 2018 June and α = 1.3 ± 0.7 in 2018 August. These slopes are somewhat lower than the expected spectral index of most circumstellar disks, where generally α = 2–3 (Beckwith & Sargent 1991; Ubach et al. 2012; Liu et al. 2018), and are consistent with free–free emission. We note that the spectral indices we measure with our NOEMA data in V2494 Cyg (α = 2.5–2.6) and V2495 Cyg (α = 2.3–2.5) are in line with those of most circumstellar disks, thus free–free emission was likely not a significant contributor during those observations.

These possibly shallower spectral indices we find are suggestive that the slope of the millimeter emission of V1735 Cyg decreased from 2014 to 2017 while we see an increase in millimeter emission, and may indicate that the free–free emission of V1735 Cyg has increased to become a significant contributor to the overall SED near 2.7 mm. Free–free emission has been tied to ionized jets/winds in objects with disks (e.g., Macías et al. 2016; Ubach et al. 2017; Espaillat et al. 2019), and these jets/winds are linked to accretion (Frank et al. 2014). One would then expect to see signatures of accretion variability in V1735 Cyg, which may be traced in the IR. However, the WISE photometry shows no significant variability between 2014 and 2017 (Figure 8). It may be the case that the IR emission is variable due to accretion, but was not detected with the cadence of WISE.

Liu et al. (2018) note that free–free emission is not thought to be significant in FUor objects based on previous observations (see Rodriguez et al. 1990; Liu et al. 2014, 2017; Dzib et al. 2015). This may indeed be the case for certain objects and/or during quiescent states without enhanced accretion, but free–free emission may become significant following an accretion event. As such, future observations of FUor objects would benefit not only from multiepoch observations, but also from multiwavelength millimeter and centimeter observations (Liu et al. 2017). This will help inform how significant, if at all, free–free emission is for a given object. If significant, free–free emission may lead to overestimated disk masses.

5. Summary

We observed six FUor objects (V1735 Cyg, V2494 Cyg, V2495 Cyg, V1057 Cyg, V1515 Cyg, and V733 Cep) in 2017 at 2.7 mm. Motivated by comparison to previously published works, we then followed up with coordinated 2.7 mm and optical (V, R, I) observations for three objects (V1735 Cyg, V2494 Cyg, and V2495 Cyg) to probe for flux variability. We did not see variability outside our stated uncertainties (∼10%–15%) from 2017 to 2018 in either our millimeter or optical observations. However, we do see a ∼80% increase in the 2.7 mm flux density of V1735 Cyg in our 2017 June data relative to archival 2014 April data from F17. Although we took steps to mitigate the effect of differing uv coverages for each observation, it should be noted that they may still have had effects on our measurements.

We can likely rule out thermal changes in the disk as the source of millimeter variability in V1735 Cyg since 3.4 and 4.6 μm WISE photometry from 2014 to 2017 displayed no signs of corresponding variability, indicating that the millimeter variability is not related to temperature changes in the inner disk. Gas and dust buildup in the disk is also unlikely to be the sole mechanism behind the observed millimeter variability given that the mass transfer rate from the envelope to the disk necessary to account for the continuum flux density changes we see would be unreasonably large (∼0.027 M_⊙ yr⁻¹).

We find that the spectral slope of V1735 Cyg is shallower than expected for pure thermal dust emission at 3 mm, which may indicate a significant contribution from free–free emission. We also find that the 3 mm spectral index may have decreased since 2014, indicating a significant increase in the free–free emission. We hypothesize that V1735 Cyg may have experienced a small accretion event, leading to the ionization of ejected material, increasing the free–free emission and leading to the observed millimeter variability. If confirmed, this could imply that previously reported disk masses of FUor objects measured during enhanced accretion activity may be overestimated. Future study of FUor objects will benefit from both multiepoch and multiwavelength observations to disentangle the free–free component from that of thermal dust emission and allow for more accurate disk mass estimates, which will help constrain what role thermal/gravitational instabilities have in triggering FUor outbursts.

We thank the anonymous referee for a careful review and suggestions that greatly improved this paper. J.W., C.C.E., and E.M. acknowledge support from the National Science Foundation under CAREER grant AST-1455042. Á.K. acknowledges funding from the European Research Council under the European Union's Horizon 2020 research and innovation program under grant agreement 716155 (SACCRED). This work is based on observations carried out under project number SX17AG and S18AX with the IRAM NOEMA interferometer. IRAM is supported by INSU/CNRS (France), MPG (Germany), and IGN (Spain). These results made use of the Lowell Discovery Telescope (formerly Discovery Channel Telescope) at Lowell Observatory. Lowell is a private, nonprofit institution dedicated to astrophysical research and public appreciation of astronomy and operates the LDT in partnership with Boston University, the University of Maryland, the University of Toledo, Northern Arizona University, and Yale University. The Large Monolithic Imager was built by Lowell Observatory using funds provided by the National Science Foundation (AST-1005313). This publication makes use of data products from the Near-Earth Object Wide-field Infrared Survey Explorer (NEOWISE), a project of the Jet Propulsion Laboratory/California Institute of Technology. NEOWISE is funded by the National Aeronautics and Space Administration. This research made use of APLpy, an open-source plotting package for Python (Robitaille & Bressert 2012).

Facilities: IRAM:NOEMA - , LDT. -

Software: CASA, Python, Astropy, APLpy.

A Study of Millimeter Variability in FUor Objects

Article metrics

MathJax

Permissions

Author affiliations

ORCID iDs

Dates

Abstract

1. Introduction