On the Diversity of Asymmetries in Gapped Protoplanetary Disks

The ALMA continuum images of the samples are presented in Figure 1. All images have a high S/N; the median S/N of these images is 55, and the lowest S/N is 40. The azimuthal and radial profiles are extracted for each image using the position angle and inclination of the outer disk. As most asymmetric disks have a moderate to face-on inclination, optical depth effects are not considered to be affecting the results significantly. The only exception is IRS 48 with a 50° inclination, which implies that the continuum contrast might be overestimated by a factor of a few. The radial profiles are taken by averaging the ±30° on either side of the angle of the peak of the asymmetry, if present, or around the major axis position angle for axisymmetric disks. The radial profiles provide the radial locations of the dust rings, in combination with more detailed analysis from the literature (after correction for the new Gaia distances). These radii are also listed in Table 2 under the R_dust column. At each radial location, the azimuthal profile is extracted using a radial width of half the beam size. Figure 2 presents both the radial and azimuthal profiles after normalization, where the latter are split into asymmetric and nonasymmetric structures. In the asymmetric curves, the profiles are normalized to the flux at the opposite side of the asymmetric peak; in case of a nondetection on that side, a 3σ upper limit is assumed for the normalization. Note that disks with multiple rings appear multiple times, with one curve for each ring. In the radial plots, the deprojected beam profile is overplotted at the location of the dust ring to show how well the ring is resolved radially.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Contrasts in the asymmetric rings between the peak and the opposite side range from ∼3 to 395. The disks of SR 21 and CQ Tau contain two asymmetric features along the same ring. The disk of SR 21 is not well resolved radially, and higher-resolution ALMA data show that these asymmetries are in fact more pronounced (Muro-Arena et al. 2020, T. Muto et al., in preparation). Also, for the CQ Tau disk, the asymmetries are moderate, and higher-resolution images (M. Benisty et al., in preparation) confirm that these asymmetric features are real. Dust continuum asymmetries are marginally optically thick (τ ∼ 0.5, using the brightness temperature and the expected temperature at the dust ring location, which is computed using Equation (4) in Francis & van der Marel 2020), so the intensity contrasts do not directly correspond to the density contrasts. The symmetric rings show moderate azimuthal variations with values between 1 and 2 that may be optical depth effects in combination with the orientation of the disks or minor asymmetric dust traps.

Both MWC 758 and AB Aur show a clear indication of an eccentric dust cavity, as the inner dust disk detected in the ALMA continuum image is offset from the center of the disk.

The azimuthal contrast in the images cannot be compared reliably across the different images, as its value depends on the beam size and S/N along the ring. Instead, the asymmetries can be described reasonably well by 2D Gaussian profiles in the radial and azimuthal direction, where the FWHM is a measure of the asymmetric nature of the ring. Each asymmetric disk is fit in the visibility plane with a parameterized model I(r, ϕ), including such Gaussians. The visibility curves, parameterization, and best fits are given in Appendix A. For HD 135344B, a detailed visibility analysis was already performed by Cazzoletti et al. (2018). Using the σ_ϕ values, we compute the FWHM of each asymmetric feature. For the axisymmetric disks, the FWHM is set to 360°. The FWHM is listed in Table 2.

The CO isotopologue images of transition disks reveal that the gas cavity radii are well within the dust cavity radii. Gas surface density profiles have been derived from these CO images in order to quantify the depth and width of these gas gaps that can be used to derive information about possible embedded companions. Unlike the dust ring, which generally has a sharp inner edge due to the trapping (Pinilla et al. 2018b), the gas gap edge is not expected to be sharp, but it has been shown to have a gradual drop in density consistent with clearing by a companion, with the minimum at the location of the companion and the dust trapped at the outer edge. We have defined R_dust as a maximum in the intensity profile of dust emission (see previous section) and the location of the gas gap R_gap as the minimum in the gas surface density profile.

In model fitting of CO images, the gap in the surface density profile is usually parameterized (e.g., van der Marel et al. 2016b) in order to limit the number of free parameters. Unfortunately, the parameterization is not the same in different studies, and a comparison across the sample is challenging. Furthermore, the parameterization in early studies usually contained unphysical sharp edges at the gas cavity radii. Therefore, we reevaluate the gas surface density profiles and gap radii R_gap by analyzing the normalized azimuthal averaged intensity profiles of ¹³CO for each target (Figure 3 and Appendix B). The properties and origin of each CO image are summarized in Table 3.

Zoom In Zoom Out Reset image size

None of the ¹³CO intensity profiles reveal a gap structure as expected from simple planet–disk interaction models; the minimum is located at the center of the disk. However, the spatial resolution of these images is limited (usually on the order of the size of the gap, ∼02–03 or ∼30–50 au), so the inner gas disk emission, if detected, is either unresolved (see also Figure 9 in van der Marel et al. 2018a) or lower than predicted by these models due to either a decrease in temperature, insufficient knowledge of processes inside the planet orbit leading to a further depletion of the gas, or both. Massive planets on eccentric orbits will significantly deplete the inner gas disk compared to regular planet–disk interaction models, where planets are held fixed on circular orbits (Muley et al. 2019).

The deep, high-resolution ¹²CO images of PDS 70 do reveal a clear gas gap (Keppler et al. 2019), and high-resolution dust continuum images reveal that inner dust disks are common in transition disks (Francis & van der Marel 2020), suggesting that in fact, many of these disks indeed harbor gaps rather than cavities. Also, the high accretion rates in transition disks (comparable to those of full disks) suggest a higher gas surface density close to the star (Manara et al. 2014; Francis & van der Marel 2020). Bosman et al. (2019) found evidence that the CO temperature in Herbig disks must be significantly lower than physical–chemical models predict to explain the ratios between different rovibrational lines. Such a decrease in temperature may also be a reason that ¹³CO remains undetectable in the inner parts of the disk; thus, we assume in this study that all disks in fact harbor gas gaps.

In this work, we derive the gap location directly from the ¹³CO profile across the sample. Since R_gap cannot be directly constrained for all sources, we rely on two additional quantities, R_CO and R_peak, to infer it. In model simulations, the gap edge R_CO is defined (Rosotti et al. 2016; Facchini et al. 2018b) as the radius R where the normalized intensity $\bar{I}(R)\lt 1-0.66(1-{\bar{I}}_{\min })$. The low spatial resolution of our ¹³CO observations does not allow us to measure this parameter directly, as the gap remains unresolved, and it is unclear whether potential inner disk gas emission is confused or highly decreased. As an alternative, we measure the location of R_peak, the peak in the integrated ¹³CO emission. Inspection of the results in Facchini et al. (2018b) shows that R_peak is approximately 1.4 times larger than R_CO. Figures 11 and 12 in Facchini et al. (2018b) provide the relations between R_CO and R_gap and, for planet masses 5–15 M_Jup, the ratio R_CO/R_gap = 1.3. Using these relations, we derive the location of R_peak, R_CO, and R_gap, as listed in Table 3. The derived gap radii are well inside the dust cavity, typically at 10–20 au radius. If the relations between R_peak and R_gap are invalid because of eccentricity, R_gap is likely even further in. Unfortunately, no grid of models exists for planet–disk interaction models including eccentricity with predictions for the CO emission, so we have to make the assumption here that the radial gap shape remains similar. For SR 21, no gap was resolved in ¹³CO, and the inner radius derived from the rovibrational emission (Pontoppidan et al. 2008) is assumed. For PDS 70, no ¹³CO data are available, but the ¹²CO profile directly reveals the gas gap in the image (Keppler et al. 2019). We notice that the ratios above still recover the gap radius at almost the exact same location for the R_peak of ¹²CO in this case. For HD 169142, DM Tau, HD 163296, and TW Hya, no gas gaps were resolved, and the gap locations are estimated to be located in between the dust ring radii. The drop in emission in the center of HD 163296 is caused by continuum oversubtraction (Isella et al. 2016).

Second, we present the first moment maps (velocity maps) of ¹²CO data in Appendix C for each of our targets. A twist pattern (deviation from Keplerian rotation) in the inner part of the disk points toward a misalignment between the inner and outer disk, or a warp. Such a misalignment can be explained by the presence of a massive companion (>1 M_Jup) that breaks the disk, leading to a different precession of the inner and outer disk (e.g., Facchini et al. 2018a; Zhu 2019). An even stronger misalignment can be induced as a result of a secular resonance between the companion and the disk (Owen & Lai 2017). Also, radial flows of gaseous material from the outer to the inner disk have been proposed to explain the twist pattern (e.g., Price et al. 2018), but it is almost impossible to distinguish observationally from a misalignment (Rosenfeld et al. 2014) and not unique for substellar companions such as found in HD 142527, as lower-mass companions can result in radial flows as well (Calcino et al. 2020). Misalignment has been discovered independently in several targets through shadows (e.g., Marino et al. 2015), dippers (e.g., Ansdell et al. 2016), and direct measurements of the inner dust disk orientation (e.g., Francis & van der Marel 2020).

For the targets in this study, a warp was confirmed for four targets that were previously found in the literature: IRS 48, HD 142527, MWC 758, and J1604–2130. Also, AB Aur appears to show non-Keplerian motion on larger scales, but this is most likely due to the strong contributions from the spiral arms detected in ¹²CO (Tang et al. 2017). For the other disks, no warp was detected, but we present a comprehensive overview of the properties of the observations suggesting that warps are impossible to detect with the available spatial/spectral resolution (Table 6). The table also provides references for other studies suggesting misalignment based on other data. Overall, all targets possibly have a misalignment between the inner and outer disk, pointing toward the presence of a massive companion. Derivation of the mass of the companion requires detailed knowledge of the viscosity, scale height, and precession time (Figure 12 in Zhu 2019), so no quantitative information can be derived from these maps.

Many transition disks have been studied in direct imaging searches for embedded companions. Figure 4 presents the best known limits for companions in the disks in our sample, in comparison with the locations of the dust rings and gas gaps. The right y-axis provides the mass ratio with respect to the stellar mass in percentage. References for high-contrast direct imaging searches are provided in Table 2, and additional data are discussed below. The mass upper-limit curves are computed in these works by comparing the brightness limit with evolutionary models such as BT-SETTL (Allard 2014) or COND (Baraffe et al. 2003), assuming the age of the system and hot-start models. When values for both models were computed, we adopted the BT-SETTL results. Companion detections are marked with red (unconfirmed) and blue (confirmed) circles. The evolutionary models for upper limits and most candidates do not include contributions from a circumplanetary disk (CPD), which could dominate the brightness and lower the mass limit estimates by a factor of 10 (Zhu 2015). Also, extinction, age estimate, and choice of evolutionary model may affect the derived mass estimates and limits.

Zoom In Zoom Out Reset image size

For PDS 70 and LkCa 15, Hα and multiwavelength data (spectral energy distribution, SED) provide constraints on the CPD, which means that their estimated companion mass is much lower than the typical contrasts in other disks. Although the SED of HD 142527B might be explained with a planetary companion with a CPD as well (Brittain et al. 2020), the companion mass is very likely substellar based on a proper-motion study of the primary star (Claudi et al. 2019) and previous fits of the SED and SINFONI H + K spectrum to BT-SETTL models (Lacour et al. 2016; Christiaens et al. 2018).

For IRS 48, upper limits on companion brightness were measured by Ratzka et al. (2005) with speckle imaging at 015 and 05 (20 and 67 au), converted to mass limits of 100 and 50 M_Jup by Wright et al. (2015), which are interpolated in our contrast curve. In addition, Simon et al. (1995) measured a K-band contrast in the regime 002–1'' (∼3–135 au) using the lunar occultation method, which was converted to a mass limit of 150 M_Jup (Wright et al. 2015). The latter is marked with an orange line in Figure 4. Also, for DoAr 44, we add limits from Ratzka et al. (2005) converted by Wright et al. (2015). For AB Aur, no contrast curves have been derived to our knowledge. For SR 21, Sallum et al. (2019) performed a detailed analysis using sparse aperture masking detecting features around 7 au, but modeling showed that these were more consistent with an inner dust ring rather than a companion. The features require a warped inner disk or spiral features. For LkCa 15, the companion candidates c and d identified by Sallum et al. (2015) have been suggested to originate from scattered light by inner disk material in follow-up studies (Thalmann et al. 2016; Currie et al. 2019), but we include LkCa 15b as a companion candidate due to its detection at Hα and the lack of polarized emission at its location.

In additions to the limits presented in Figure 4, we show additional constraints from a brown dwarf survey through sparse aperture masking for SR 21, DoAr 44, J1604–2130, LkCa 15, DM Tau, and TW Hya (Kraus et al. 2008, 2011; Evans et al. 2012; Cheetham et al. 2015) for the inner 015 of the disk and additional limits for SR 21 for the inner 5 au (S. Sallum, private communication) using the data from Sallum et al. (2019). Unlike coronagraphy, sparse aperture masking allows the detection of companions at angular separations well within the diffraction limit down to a few au at typical disk distances (e.g., Sallum & Skemer 2019). Mass limits were derived using a range of evolutionary models using the procedure described in Kraus & Hillenbrand (2007).

A handful of indirect estimates of companion candidates are known from the literature, using a range of techniques.

• 1.  
  Baines et al. (2006) claimed that AB Aur has an accreting binary companion at 82–489 au separation using Hα spectro-astrometry, which could be either inside or outside the cavity and with unknown mass.
• 2.  
  Boccaletti et al. (2020) deduced a planet of 4–13 M_Jup at 30 au in AB Aur based on a spiral arm twist.
• 3.  
  Gratton et al. (2019) found a tentative detection of a 3 M_Jup companion at 38 au separation in HD 169142.
• 4.  
  Willson et al. (2016) found a tentative detection of a companion at 6 au separation in DM Tau using sparse aperture masking without significant sky rotation.
• 5.  
  Pinte et al. (2018) and Teague et al. (2018) found evidence for ∼2 M_Jup planets at 83, 137, and 260 au in HD 163296 using deviations of Keplerian motion in ¹²CO channel maps.
• 6.  
  Calcino et al. (2019) and Poblete et al. (2020) claimed substellar companions with mass ratios of ∼0.2 at 10 and 30 au in IRS 48 and AB Aur, respectively, to explain the kinematics in the system.

Due to their more speculative nature, these candidates are not marked in the contrast curves in Figure 4. Follow-up observations are required to confirm their existence.

In Figure 5, we present the estimated gas surface density profile, described as a Gaussian centered on R_dust and an inner width consistent with the R_gap location. Here R_dust is chosen as the outer edge of the gas gap, as the pressure maximum is thought to be located there. The R_peak is inward of R_dust, but this can be understood by the radial temperature dependence; whereas the gas surface density is decreasing, the ¹³CO emission remains partially optically thick and peaks further inward. For DM Tau and SR 21, no R_peak could be measured from the CO profile, and it was estimated to be located at ∼75% of the dust ring radius. The profile is scaled to the derived gas surface density profile from the literature based on combined ¹³CO and C¹⁸O data (see references in the last column of Table 3) to match the surface density at peak, gap, and outer disk locations. The gap radii in the literature profiles were rescaled to the Gaia distances. The literature surface density profile is overplotted as a blue dotted line, and Table 3 lists the gas cavity radii R_gascav from these profiles from the literature. For PDS 70 and J1604.3–2130, the CO data were fit with a gap-like profile (material inside the companion orbit). For HD 142527 and MWC 758, the gas surface density profile was described as a Gaussian in the literature analysis (Boehler et al. 2017, 2018); R_gascav is taken as the radius where the density drops below 10⁻² g cm⁻², where ¹²CO becomes optically thin. For AB Aur, no detailed physical–chemical model was used, and the CO abundance was taken to be constant throughout the disk, so the derived profile remains highly uncertain. Because TW Hya and HD 163296 do not have a resolved inner gas gap, no values are provided here. Though HD 163296 does show gas gaps in the outer disk (Isella et al. 2016), the depth and width remain highly uncertain, and we refrain from including them in the plot (van der Marel et al. 2018b).

Zoom In Zoom Out Reset image size

The inner part of the gas disk remains largely unconstrained by our analysis of the CO images due to spatial resolution. We indicate the unconstrained regime with a gray area between zero and the beam radius in Figure 5. This representation reveals that, except for PDS 70, the distinction between a gas gap and gas cavity cannot be made for any of the disks. However, inner dust disks have been detected in about half of the sample (Francis & van der Marel 2020), and all disks show signs of significant gas accretion onto the star, which makes it most likely that gas is still present in the inner disk as well and the gas "cavities" are in fact gas gaps, consistent with clearing by a companion. The derived inner disk dust profiles from Francis & van der Marel (2020) and the expected gas surface density based on the accretion are included in the plot to reflect this. The expected gas surface density close to the star can be estimated indirectly assuming a viscous disk model from the stellar accretion rate by Manara et al. (2014),

$$\begin{eqnarray}&&{{\rm{\Sigma }}}_{g}(r)=\displaystyle \frac{\dot{M}2\,{m}_{p}}{3\pi \alpha {k}_{B}T(r)}\sqrt{\displaystyle \frac{{{GM}}_{* }}{{r}^{3}}},\end{eqnarray} \tag{ 2 }$$

with Σ_g(r) the gas surface density, $\dot{M}$ the accretion rate, m_p the proton mass, α the viscosity, k_B the Boltzmann constant, G the gravitational constant, M_* the stellar mass, and T(r) the temperature profile, for which we use (Dullemond et al. 2001)

$$\begin{eqnarray}&&T(r)={\left(\displaystyle \frac{\phi {L}_{* }}{8\pi {\sigma }_{B}{r}^{2}}\right)}^{1/4}=\sqrt[4]{\displaystyle \frac{\phi {L}_{* }}{8\pi {\sigma }_{B}}}\displaystyle \frac{1}{\sqrt{r}},\end{eqnarray} \tag{ 3 }$$

with L_* the stellar luminosity, ϕ the flaring angle (taken as 0.02), and σ_B the Stefan–Boltzmann constant. The stellar properties are taken from Francis & van der Marel (2020). We compute the expected local gas surface density Σ_acc at 1 au derived from the accretion rate assuming α = 10⁻³ (see Table 3) and overplot this value on the derived profiles (Figure 5).

In order to estimate the companion mass assuming a single planet, the millimeter-dust radius R_dust is compared with the derived gas gap radius R_gap. Figure 6 presents this relation. Best-fit relations between the ratio and the planet mass were derived by Facchini et al. (2018b) for planet masses between 1 and 15 M_Jup for an average between α = 10⁻³ and 10⁻⁴, but as no models were run for higher-mass companions, this relation cannot be used to estimate accurate masses for these ratios. Also, Facchini et al. did not consider eccentric orbits, which significantly increase the separation between the dust ring and the gas gap (Muley et al. 2019). The majority of our disks lie in the regime >15 M_Jup, suggesting that they contain planets above this threshold, in the brown dwarf regime. The ratio between R_dust − R_gap and R_gap for our sources is typically between 1.3 and 2.5. Outliers are CQ Tau and SR 21 with even higher ratios (3–4), suggesting very massive companions, potentially (sub)stellar.

Zoom In Zoom Out Reset image size

A similar comparison was made between R_dust and the inner edge of the scattered-light gap for a sample of transition disks (Villenave et al. 2019) using the planet mass relations derived by de Juan Ovelar et al. (2013) for α = 10⁻³. For the overlapping targets, their estimates of companion mass are consistent with ours, with several substellar mass companions. Exceptions are IRS 48 and Sz 91, for which they used nonscattered-light observations for which the planet relation does not hold, and LkCa 15, for which the inner edge is more challenging to determine in scattered light due to its high inclination. Also, several other disks in their sample (not in our work) that are claimed to be in the planetary regime might suffer from high inclination.

With the azimuthal profiles fit in Appendix A and the gas surface density profiles described in Section 4.1, we construct a plot of the azimuthal FWHM as a function of the Stokes number of the traced dust grains for the radial location of the dust. The Stokes number St as defined in Equation (1) cannot be used directly, as dust continuum emission is originating from a large range of grain sizes (and Stokes numbers). Therefore, we use a simplification with the assumption that particles with size a_grain = λ_obs/2π (Draine 2006) are the primary contributors at observing wavelength λ_obs and introduce the observational Stokes number,

$$\begin{eqnarray}&&{{St}}_{\mathrm{obs}}=\displaystyle \frac{{\lambda }_{\mathrm{obs}}{\rho }_{s}}{4{{\rm{\Sigma }}}_{\mathrm{gas}}({R}_{\mathrm{dust}})},\end{eqnarray} \tag{ 4 }$$

with the gas surface density Σ_gas at the location of the dust ring R_dust. The result is shown in Figure 7. For the uncertainties, we assume an uncertainty of a factor of 3 on Σ_gas(r), based on the typical uncertainty on the gas surface density based on CO isotopologue data as derived by Woitke et al. (2019). Although other grain sizes than a_grain may contribute to the emission, which might add additional uncertainty, this is not considered an issue, as all Stokes numbers are computed in the same way. As it is reasonable that the grain size distributions are similar across the sample (under the assumptions that the disks have similar ages and are evolving in similar physical environments), it would thus shift all data points in the same direction, and the trend would remain the same.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Figure 7 shows that axisymmetric disks have low values of the observational Stokes number, but the asymmetric features are all located at St_obs > 10⁻². The derivation of gas surface density from CO isotopologue data remains uncertain, in particular due to problems with our knowledge of the carbon budget in disks (Kama et al. 2016b; Miotello et al. 2019). In addition, we show the dependence of the azimuthal contrast from the images on the observational Stokes number, which also shows a distinction between symmetric and asymmetric disks. As this trend might be affected by imaging artifacts, it is not further discussed.

To explore the physical implications of these results, we employ a model presented in Birnstiel et al. (2013) that analytically solves for the equilibrium of azimuthal drift and mixing of dust particles. The azimuthal contrast of this model was shown to be in good agreement with the 2D calculations of Lyra & Lin (2013). Our model consists of three steps: (1) constructing a particle size distribution, (2) calculating the azimuthal equilibrium density distribution of each particle size according to Birnstiel et al. (2013), and (3) calculating the dust intensity profiles and thus the azimuthal dust intensity contrast. For the first step, the particle size distribution, we tried two different choices: first, a truncated power-law size distribution with an MRN exponent (Mathis et al. 1977) up to a maximum particle size a_max, and second, the steady-state distributions of Birnstiel et al. (2011). For the second step, the azimuthal density distribution, we employed Equation (8) of Birnstiel et al. (2013) and parameterized the azimuthal gas density as a constant density plus a Gaussian overdensity,

$$\begin{eqnarray}&&{{\rm{\Sigma }}}_{g}(r,y)={\bar{{\rm{\Sigma }}}}_{g}(r)\displaystyle \frac{1+(A-1){e}^{-\tfrac{{y}^{2}}{2{r}^{2}{\sigma }^{2}}}}{1+(A-1)\tfrac{\sigma }{2\sqrt{2\pi }}\,\mathrm{Erf}\left(\tfrac{\sqrt{2}\pi }{\sigma }\right)},\end{eqnarray} \tag{ 5 }$$

where y is the azimuthal coordinate going from zero to 2π r and σ_y is the azimuthal extent of the bump. Here Σ_g(r, y) is normalized such that the azimuthal average gives Σ_g(r).

The emission profile was then calculated assuming absorption opacity and a face-on, vertically isothermal disk, such that the intensity becomes ${I}_{\nu }={B}_{\nu }({T}_{\mathrm{dust}})\left(1-{e}^{-\tau }\right)$, where τ = ∑_iΣ_i κ_abs,i is the optical depth. Here Σ_i and κ_abs,i are the surface density and absorption opacity of grain size i. Optical depth is thus taken into account for the comparison with the data. The final azimuthal FWHM is computed directly for the intensity profile of the model. Although this is a rather simple approach, it is sufficient for the purpose of reproducing the trend of azimuthal extent with regard to observational Stokes number.

We explore a number of parameters, including α, gas-to-dust ratio, gas contrast, and gas azimuthal extent, and find the best fit, accompanied by two values on either side to show the dependence of the curve on the parameter (Figure 8). The threshold of the observational Stokes number where disks become asymmetric (∼10⁻² for most of our targets) depends on the disk properties and may vary somewhat from disk to disk. Furthermore, we explore whether the FWHM is set by fragmentation (different fragmentation velocities) or a default grain size distribution with a maximum grain size. Fragmentation velocities in lab experiments range from 1 to 10 m s⁻¹ (e.g., Blum & Wurm 2008) and even higher outside the snowline (Wada et al. 2009), although the latter has been called into question by recent experiments (Steinpilz et al. 2019). The fiducial model (orange) shows a possible combination based on a manual fitting procedure: α = 10⁻³, gas-to-dust ratio = 10, gas contrast = 1.2, a_max = 1 mm, and azimuthal σ_y = 10°. In addition, we plot the azimuthal FWHM as a function of temperature (at the dust ring) in this figure.

A final aspect that is relevant for this discussion is the presence of spiral arms and the link with asymmetries. Table 2 indicates which disks show spiral arms in scattered light, with the references provided in the last column. About half of the sample (all asymmetric disks) shows spiral arms. For DM Tau, no scattered-light imaging data are available, and the spiral nature of HD 169142 remains uncertain, as the spiral arms are seen in total scattered light (Gratton et al. 2019) through angular differential imaging but not in polarized scattered light (Pohl et al. 2017; Bertrang et al. 2018). The former technique suffers, unlike polarized differential imaging, from possible biases deriving from the disk emission self-subtraction, in particular if the disk is seen face-on.

4.3.1. Link with Asymmetries

Garufi et al. (2018) noticed that all asymmetric disks in ALMA show spiral arms in scattered light, and the apparent origin of one of the spiral arms in HD 135344B in the dust asymmetry (van der Marel et al. 2016a) suggests that these phenomena are physically linked if a vortex triggers the spiral arm (Lovelace & Romanova 2014). However, recent simulations show that spiral arms triggered by a vortex are unlikely to be detectable in scattered light (Huang et al. 2019), and the spirals must have a different origin, such as a companion.

Figure 9 presents a number of trend plots between the azimuthal extent and parameters that have been linked to spiral arms. The pitch angle was derived from the deprojected scattered-light images from the literature (for references, see Table 2). When no spiral arms were detected, the pitch angle is set to 1°. No difference is shown between "single" and "double" spirals, as the secondary spiral can easily be hidden in part of the disk. Ages, near-IR (NIR) excess, and stellar masses are taken from Garufi et al. (2018). For each plot, we compute the correlation coefficient r_corr using the linear regression procedure by Kelly (2007), resulting in values of r_corr = 0.0 ± 0.6 for each plot consistent with a lack of correlation.

Zoom In Zoom Out Reset image size

Although asymmetries are only seen in disks with spiral arms, as already demonstrated by Garufi et al. (2018), there is no trend between azimuthal extent and pitch angle. The pitch angle itself depends primarily on the disk temperature (hence aspect ratio) and to a lesser degree on planet mass (Fung & Dong 2015; Zhu et al. 2015). The asymmetric nature of CQ Tau is debatable, but spiral arms have been found. Asymmetries and spiral arms are only seen in stars with stellar masses >1.5 M_⊙, but not exclusively; HD 163296 is an intermediate-mass star without spiral arms or azimuthal asymmetries. Disks with spirals and asymmetries exist for a range of ages within the sample, although no young disks with spiral arms or asymmetries are known. The sample is intrinsically biased, as it only contains disks with ages of ∼4–10 Myr, which represent a minority disk population that lives longer than the average disk lifetime of ∼3 Myr (Mamajek 2009). The presence of wide gaps and thus dust traps is thought to be the main reason for the longer lifetime, as radial drift is efficiently reduced (Pinilla et al. 2020). The NIR excess has been linked before to spiral arms (Garufi et al. 2018), possibly due to changes in the temperature structure in the shadows in the outer disk (Montesinos et al. 2016). Most disks with high NIR excess also show asymmetries, but not for the entire sample.

In the bottom right panel of Figure 9, we test whether the detection of spiral arms can be linked to the local gas surface density in the same way as the millimeter asymmetries through the Stokes number, following Veronesi et al. (2019). We thus aim to test whether the Stokes numbers of micron-sized dust grains are significantly different in disks with spiral arms and disks with rings in scattered light (see Table 2). In each disk, the radial range of the scattered-light features (rings and spiral arms) is estimated from the literature, and the Stokes number of a micron-sized dust grain is computed using the gas surface density profile in that range with Equation (4). The Stokes number of the micron-sized grains is similar throughout the sample, regardless of contrast or spiral arms. An observed trend would contradict the result of Veronesi et al. (2019), who predicted that rings are only visible at (much) higher Stokes numbers at millimeter wavelengths. The presence of spiral arms might thus be unrelated to the local gas surface density.

4.3.2. Link with Gaps and Stellar Properties

We explore the link between spirals and stars further in a wider sample comparison in Figure 10, as more massive stars are generally associated with high luminosities as well. Targets are taken from the scattered-light demographics study by Garufi et al. (2018) for which high-resolution ALMA data are available in Francis & van der Marel (2020) and Andrews et al. (2018) to estimate the dust gap width. Here R_dustgapwidth is defined as R_cav for the transition disks (Francis & van der Marel 2020) and as the gap width of the inner gap in Zhang et al. (2018) for the ring disks. The assessment of the presence of spirals is primarily based on the classification by Garufi et al. (2018, Figure 1), where spirals and giants are marked as "spiral" in our sample, rings and rims as "no spiral," and faint, small, or inclined as "unconfirmed," since the detection of spirals in these disks is hindered by the observational sensitivity, angular resolution, and disk geometry, respectively. Four giants from Garufi et al. (2018) were marked by those authors as controversial due to their high inclination, which makes the detectability of spiral arms more challenging. These are marked as "unconfirmed" in our plot. All data are provided in Appendix D.

Zoom In Zoom Out Reset image size

Figure 10 shows that all spiral disks lie in the upper right corner of the diagram, with high luminosity and a large gap width, in contrast to nonspirals, with either low luminosity or a narrow gap width. This difference can be understood as the detectability of spiral arms increases with pitch angle; the pitch angle is correlated with the aspect ratio (disk temperature) and, to some extent, companion mass (Dong et al. 2015; Fung & Dong 2015; Zhu et al. 2015). As the disk temperature generally scales with stellar luminosity (Dullemond et al. 2001) and the gap width roughly with companion mass (Varnière et al. 2004), the pitch angle is expected to be larger for more luminous disks with the most massive companions. This link supports a view where the scarcity of spirals around T Tauri stars is due to their low luminosity rather than their young age (Garufi et al. 2020). A larger pitch angle is more easily resolved and thus more likely to be detectable in NIR observations at low inclination. High inclination angles make the detectability of spiral arms more challenging (Dong et al. 2016), and indeed, the three purple data points in the upper right corner of Figure 10 are all disks with inclination i ∼ 40°–75°. The two almost face-on disks with tightly wound spirals and a large empty cavity (HD 142527 and GG Tau) may remain a separate category, as these are the only confirmed binaries in this comparison. Another possible connection is that more luminous stars are generally also more massive, and higher-mass stars are more likely to have a binary companion (Raghavan et al. 2010), although binary companions have not been detected yet in the majority of these disks (see Section 3.3).

This result demonstrates that spiral arms are possibly present in all disks with gaps, assuming all gaps are opened by planets, but only detected when the planet is sufficiently massive, the star sufficiently luminous, and the inclination angle not too high. This explains why spirals are uniquely found in low-inclination disks with wide gaps and a high-luminosity star. The thresholds appears to be at L_* ≳ 1.5 L_⊙ and R_gapwidth ≳ 15 au, but more data are required to confirm this.

4.3.3. Link with Morphology

This connection between detectability, luminosity, and gap width thus explains the locations of the spiral disks in Figure 9 in mass and luminosity, but the link with azimuthal extent remains unclear. Figure 11 presents an overlay of scattered-light images and ALMA data. The spiral arm(s) and dust asymmetries appear to be spatially connected, suggesting that the spiral arm and dust asymmetry may be physically related. Such a comparison has been made before for IRS 48 (Follette et al. 2015, Figure 12), MWC 758 (Dong et al. 2018b, Figure 1(c)), HD 135344B (Cazzoletti et al. 2018, Figure 1(b)), SR 21 (Muro-Arena et al. 2020, Figure 3(d)), and CQ Tau (Uyama et al. 2020, Figure 4) for the disks in our sample, but also for, e.g., V1247 Ori (Kraus et al. 2017, Figure 1(b)) and HD 100453 (Rosotti et al. 2020, Figure 5). Proposed scenarios for these connections include, e.g., the launching of a spiral by the vortex (van der Marel et al. 2016a) and the detection of part of the spiral in millimeter emission (e.g., Rosotti et al. 2020).

Zoom In Zoom Out Reset image size

No physical connection is visible in HD 142527 (Avenhaus et al. 2017), AB Aur (Tang et al. 2017), GG Tau (Keppler et al. 2020), or HD 143006 (Pérez et al. 2018), but this could be related to the limited detectability of the spiral arms in these systems around the radius of the dust asymmetry.

In Section 4.2, we compared the azimuthal extents of the dust asymmetries with the local gas surface density through the observational Stokes number and found a steplike trend that can be matched to a simple dust evolution model with an azimuthal pressure bump in the gas. This result suggests that the diversity of asymmetries and nonasymmetries is not linked to the disk, the companion, or a limited lifetime but rather to the local gas surface density at the location of the pressure bump. This implies that minor azimuthal pressure bumps may be very common in disks, but they are only detected as dust asymmetries when the Stokes number is sufficiently high, i.e., when the local gas surface density is sufficiently low. This also leads to the prediction that dust observations at centimeter wavelengths, such as the ngVLA, will show a much larger number of asymmetric disks, as a higher observational Stokes number is traced at these wavelengths.

This scenario also explains the existence of dust asymmetries in outer rings, such as seen in, e.g., HD 135344B (this study) but also V1247 Ori (Kraus et al. 2017) and HD 143006 (Andrews et al. 2018), which were not included in this study due to a lack of gas analysis. Furthermore, a dust feature identified in high-resolution data of TW Hya at 1.3 mm at 52 au was interpreted as either a CPD or a small azimuthal dust trap (Tsukagoshi et al. 2019). Our data have insufficient sensitivity to reveal this feature. The Stokes number (using our gas surface density profile) is ∼0.3 at the location of the feature, which follows the same trend as the other data points in Figure 7, and the feature could thus indeed be another azimuthal dust trap. As there are no clear correlations between azimuthal extent and typical spiral arm properties (Figure 9), asymmetries may be unrelated to the location or mass of the companion.

If this scenario is correct, asymmetries are not related to the lifetime of a vortex or gas horseshoe. Previous studies have suggested that vortices dissipate on relatively short timescales due to dust feedback (Fu et al. 2014; Miranda et al. 2017), although 3D simulations do not reproduce rapid vortex dissipation (Lyra et al. 2018). Hammer et al. (2017) showed that vortices induced by planets may have limited lifetimes when the planet mass growth is not sufficiently fast. Both scenarios have been used to argue that the occurrence rate of asymmetries is caused by the limited lifetime. Our work demonstrates that a timescale may be irrelevant for the occurrence rate. This implies that dissipation of vortices and/or horseshoes could happen on much longer timescales than the lifetime of the disk.

Figure 8 also shows typical values for the gas overdensity consistent with the observations. Several parameters are redundant with each other; a lower viscosity requires a lower gas overdensity and/or gas-to-dust ratio. As vortices are thought to survive only when α ≲ 10⁻⁴ (de Val-Borro et al. 2007), this implies that for vortices, the gas-to-dust ratio is likely to be closer to unity, whereas gas horseshoes also survive at high α and may have higher gas-to-dust ratios and/or gas overdensities. Furthermore, the extent dependence becomes shallower for a wider azimuthal width of the gas bump.

Another interesting aspect is the choice of the dust grain size distribution. Both the equilibrium model using a fragmentation velocity and a grain size distribution with a fixed maximum grain size can reproduce the curve. This means that we can neither rule out nor confirm that fragmentation of grains is the limiting effect for dust growth inside the dust trap. Future multiwavelength data might be able to probe the material properties of these dust particles. In the fragmentation limit, the maximum Stokes number should be inversely proportional to α · T (Birnstiel et al. 2011), but the measured asymmetric contrast does not show a clear dependence on the temperature (Figure 8); linear regression analysis results in an r_corr coefficient of −0.3 ± 0.4.

The scenario described in Section 5.1 where dust asymmetries are linked to the local gas surface densities leaves the question of the origin of the gas asymmetry open. Gas horseshoes only appear at the inner edge of a wide eccentric gap and are unable to trigger secondary asymmetries, so they can be excluded for asymmetries in disks with multiple rings with asymmetries and multiple asymmetries, such as seen in SR 21 and CQ Tau. The disks with single dust rings and asymmetric features (HD 142527, IRS 48, and AB Aur) could still be explained by either gas horseshoes or vortices as the result of RWI.

The main distinctions between the gas horseshoe and vortex scenarios are the disk viscosity and the mass of the companion; vortices require α ≲ 10⁻⁴ to survive for a sufficient amount of time (e.g., Godon & Livio 1999; Regály et al. 2012), whereas the gas horseshoe can exist at higher α, and an RWI occurs at the edge of a planet gap as long as the planet is massive enough to carve a deep gap (≳1 M_Jup), whereas the gas horseshoe requires a mass ratio q > 0.05, corresponding to ≳50–100 M_Jup for 1–2 M_⊙ stellar mass. The RWI also develops at the edge of an eccentric gap as long as α ≲ 10⁻⁴ (Ataiee et al. 2013). The companion in the HD 142527 system has been estimated to be ∼150–440 M_Jup or q = 0.07–0.21, consistent with the horseshoe scenario (Price et al. 2018), but in other disks, no such companion has been identified. As the viscosity remains very challenging to constrain observationally, the companion mass provides the best constraint on the origin of the single dust asymmetries. The possible companion masses are discussed in the next section.

For the disks with multiple rings and asymmetric features, a vortex remains a likely explanation under the assumption that disks have a viscosity α ≲ 10⁻⁴. Hydrodynamic simulations show that the RWI always develops in these conditions at the edges of gaps, and vortices should be very common in disks. Our results demonstrate that the lack of detections of these vortices could simply be due to the local gas surface density.

Another possibility is that the underlying gas asymmetry is in fact part of the spiral density wave, which is supported by the physical connection between spiral arms and the millimeter-dust features discussed in Section 4.3. The reason that only a small part of the spiral arm is visible in the millimeter is that millimeter grains are only present at the edge of the planet gap, where the spiral density wave can lead to a further concentration of the dust azimuthally, but again, only when the local radial gas surface density is low enough. Small mismatches between the curve of the millimeter-dust and scattered-light features such as seen in HD 135344B (van der Marel et al. 2016a) can be understood, as the emission originates from different heights in the disk (Rosotti et al. 2020).

Spiral waves launched by a planet are rotating with respect to the background gas flow; they run over the dust particles with little time for the particles to react. Dust particles thus cannot get trapped and get carried along in spiral density waves; the timescales for dust accumulation in dust traps are at least 100 times longer than the local orbital period (Birnstiel et al. 2013). However, spiral waves still lead to changes in the pressure scale height and vertical flows, which may also lead to different spatial distributions of different particle sizes that reproduce the observed morphologies. Whether asymmetric features in the millimeter continuum really can be part of a spiral density wave remains a question.

Furthermore, it is unclear whether the disks without a clear physical connection between the spiral arm and millimeter-dust features listed above are intrinsically different (e.g., because of a substellar rather than a planetary companion) or whether these large-scale asymmetries are perhaps connected to the spiral arm after all. Spiral arms traced in high-resolution ¹²CO observations (such as seen on larger scales in, e.g., HD 142527; Christiaens et al. 2014) might help to reveal this connection, as ¹²CO remains visible at lower surface densities than small dust grains.

Multi-epoch observations of asymmetric millimeter-dust features in disks are required to measure their rotational speed, in order to see if they move along with the spiral (with the orbital speed of the companion) or on their own Keplerian orbit. The latter would directly rule out trapping in a spiral density wave.

Figure 4 shows the known limits for companions in each of the disks. The three confirmed companions (HD 142527B, PDS 70 b, and PDS 70 c) and three companion candidates (MWC 758 b, LkCa 15 b, and HD 169142 b) are located at or around the gas gap radii well inside the dust ring radii. The HD 142527B companion was detected at a small separation of 12 au, but orbital fitting indicates that the orbit is highly eccentric, and the companion may reach a separation of at least 57 au at apoapsis (Claudi et al. 2019), very close to the derived gas gap radius. For the 10 systems where the limits are derived through a contrast curve, the limits are known at the gas gap location for all except one system (CQ Tau), and for SR 21, the limits are very marginal at the gap location. Only AB Aur has no limits on companions. If the gap radii are overestimated and the gap radius is even closer in to the star, the contrast curves only provide limits for about half of the sample.

The contrast curves, which have been derived using hot-start models, rule out mass ratios q > 0.05 (M_p > 50 M_Jup, the minimum threshold for the formation of gas horseshoes) in the targeted regimes. Table 4 provides the limits for the possible companions and expected structure at the gap edge. However, high mass ratios are still possible in the inner parts of the disk where no contrast could be measured, which is particularly relevant for the disks where the derived gas gap radius is not covered by the contrast curve, such as CQ Tau, SR 21, AB Aur, and perhaps IRS 48. Such a high mass ratio was recently suggested for IRS 48 and AB Aur for reproducing the CO kinematics and dust contrast through a circumbinary simulation (Calcino et al. 2019; Poblete et al. 2020).

The large derived ratios between the gas gap radius and dust ring radius from Figure 6 imply minimum companion masses >15 M_Jup, which are in the brown dwarf and stellar regime, for the assumption that the gaps are cleared by a single companion on a circular orbit. On the other hand, the contrast curves generally rule out companion masses >50 M_Jup at the gap location. This limits the companion masses (at R_gap) at the gap location to the brown dwarf regime. This would be consistent with the derived companion candidate masses in MWC 758 and HD 169142, although these masses remain highly uncertain due to the lack of available data for analysis of the contributions by a CPD, if present.

On the other hand, it is very likely that simple planet–disk interaction models with a fixed orbit, such as those used by Facchini et al. (2018b), are insufficient to derive planet masses from CO-versus-dust images. A disk gap is thought to become eccentric when the mass ratio q ≳ 0.003 due to eccentric Lindblad resonances (Kley & Dirksen 2006); the back-reaction the disk exerts on the companion has been shown to grow its orbital eccentricity (e.g., Papaloizou et al. 2001; Ragusa et al. 2018), which is efficient down to super-Jovian planet masses when accretion onto the planet is included (D'Angelo et al. 2006). Muley et al. (2019) demonstrated that a proper planet–disk interaction model with a single planet, including accretion onto the planet and migration, is able to reproduce the observed gas gap and dust ring in PDS 70 with a planet with a mass of ∼4 M_Jup after 4 Myr, with a natural eccentricity growth up to e ∼ 0.3. This simulation was run to explain PDS 70 with only PDS 70 b, as PDS 70 c was still unknown at the time. This scenario has recently been proposed to explain MWC 758 as well (Calcino et al. 2020).

We hypothesize that the wide transition disk cavities in our sample are also caused by eccentric, super-Jovian planets; these planets are in the 3–15 M_Jup regime, and their eccentric orbits have developed naturally, as shown in Muley et al. (2019). Due to the eccentric planet orbit, the disk may no longer appear eccentric. Eccentric disk cavities have been observed in MWC 758 (Dong et al. 2018b) and AB Aur (this work) but are generally hard to determine observationally. The main motivation for this scenario is thus the large separation between the dust cavity radius and the deduced gas gap, which is thought to be representative of the companion orbit. Also note that brown dwarfs are expected to carve eccentric gaps considering their high mass ratios. It is also possible that multiple companions (such as seen in PDS 70) are responsible for the wide gaps. The (sub)stellar companions required for gas horseshoes are ruled out in the majority of our disks, with the exception of IRS 48, HD 142527, AB Aur, and CQ Tau.

The occurrence rate of massive companions from direct imaging surveys in older systems argues against brown dwarf companions as a common explanation for transition disks. Super-Jovians (5–13 M_Jup) have an occurrence rate of 8.9% at 10–100 au for intermediate-mass (1.5–5 M_⊙) stars (Nielsen et al. 2019), which is the stellar mass range of the majority of our sample. Above that mass threshold, the occurrence rate of brown dwarfs (13–80 M_Jup) at wide orbits is much lower (∼1%, also known as the "brown dwarf desert"; Nielsen et al. 2019), but the occurrence of stellar companions (>80 M_Jup) or a binarity rate at 10–100 au is again increased, with a fraction of ∼15% in the 1–2 M_⊙ stellar mass range (Moe & Kratter 2019). This suggests that transition disks are more likely caused by super-Jovians or stellar companions than brown dwarfs.

The narrow gaps in HD 163296 and TW Hya are consistent with lower-mass planets, <5 M_Jup, for which eccentricity is unlikely to develop. This is possibly a distinction between so-called transition disks and ring disks (van der Marel et al. 2019); wide gaps only develop when the planet is sufficiently massive to develop an eccentricity, which requires q > 0.003. A ring disk may host multiple lower-mass companions.