© P. Pinilla et al. 2020

1 Introduction

Censuses performed on the thousands of exoplanets that have been discovered thus far have revealed a large diversity of planetary architectures, along with a variety of trends that may be the result of the common physical processes ruling the formation and evolution of planets in disks around young stars. Exoplanet statistics show that giant planet occurrence increases with stellar metallicity and mass (e.g., Fischer & Valenti 2005; Johnson et al. 2010), while sub-Neptune type planets are more common around M-dwarfs than around Sun-like stars (e.g., Mulders 2018). These observed trends stand as open challenges for current theories of planet formation, examining stellar and disk properties; in particular, questioning how the effects of the ratio between stellar mass and disk mass affect the final outcome of planet formation.

Current observations of protoplanetary disks with the Atacama Large Millimeter/submillimete Array (ALMA) provide demographic surveys of disks (Ansdell et al. 2016, 2017, 2018; Barenfeld et al. 2016, 2017; Pascucci et al. 2016; Cazzoletti et al. 2019; Cieza et al. 2019; Long et al. 2019; van Terwisga et al. 2019) that aid in the understanding of the missing keys for current planet formation models. Using continuum observations of disks, we can currently access to different disk properties, including a large range of disk dust masses – an important property to understand the amount of material available in young disks to form different types of planets. The most common method to determine disk mass is to assume that the dust emission is optically thin over most of the disk volume and, therefore, that the detected flux is proportional to the dust disk mass (e.g., Beckwith et al. 1990). This dust mass is usually used as a tracer of the total disk mass by assuming a constant dust-to-gas ratio, conventionally the interstellar medium ratio, that is, 1/100.

Obtaining information about the amount of gas available in disks is challenging because faint optically thin molecular lines, such as ¹³CO and C¹⁸O, must be observed and their intensity can be also affected by isotope selective processes (Bruderer et al. 2012; Miotello et al. 2014, 2016) and the unknown carbon abundance relative to H₂ (e.g., Kama et al. 2016; Schwarz et al. 2016). In addition, even ¹³CO and C¹⁸O may be optically thick inside the CO ice line at the disk midplane, resulting in an underestimation of the disk mass (Zhang et al. 2017), especially for warm disks around Herbig stars. In these cases, observations of rare CO isotopologues, such as ¹³C¹⁷O are needed (although challenging) for a more robust estimation of gas disk masses (Booth et al. 2019) and, therefore, the best knowledge that we have so far about the disk material still comes from the dust emission.

Information obtained about the dusty material in disks from observations is, however, prone to uncertainties as well since it is based on several assumptions, such as the disk temperature (which is typically assumed to be 20 K but could vary depending on disk size; Hendler et al. 2017), and the dust opacity. The dust opacity is the major source of uncertainty because it depends on the grain size and composition, which are unknown in protoplanetary disks (e.g., Birnstiel et al. 2018). Nevertheless, when comparing millimeter-observations of protoplanetary disks with models of disk evolution, the same assumptions can be taken in order to facilitate the comparison and understand the effect of different crucial physical processes in the formation of planets, such as grain growth and drift.

Due to the fast radial drift that millimeter- and centimeter-sized particles undergo in the outer parts (>20 au) of protoplanetary disks, models of dust evolution at million-year timescales contradict millimeter observations of disks (Birnstiel et al. 2010a; Pinilla et al. 2012). This is because the models predict a millimeter spectral index (an indicator of grain growth) that is higher than that which is observed for most protoplanetary disks in different star forming regions (e.g., Ricci et al. 2010; Testi et al. 2014). Reducing or completing the suppression of the radial drift of dust particles, for instance, with the presence of particle traps or pressure bumps, is necessary in order to keep millimeter grains in protoplanetary disks for million-year timescales (Pinilla et al. 2012).

Observational surveys that measure the disk dust mass found a power-law relation between M_dust and M_⋆ that steepens with time (e.g., Pascucci et al. 2016). This relation was recovered by dust evolution models that include the growth, fragmentation, and drift of particles from Krijt et al. (2016). Drift is crucial for reproducing the steepness of the M_dust−M_⋆ relation because it is expected to reduce the dust mass with time and because drift is more effective around low-mass stars (Pinilla et al. 2013). However, in these models, the disks are short in solids compared with observations taken after ≲1 Myr of evolution. When considering disks with large inner cavities resolved at millimeter-wavelengths (the so-called transition disks; TDs), Pinilla et al. (2018) find that this relation is much flatter, suggesting that the flatness is due to an effective reduction of radial drift in pressure bumps.

In this paper, we aim to study the effect that dust evolution has on the M_dust−M_⋆ relation. We are particularly interested in investigating the effect of particle traps in the disk on this relation. Motivated by the fact that exoplanet statistics show a dependence on stellar mass and metallicity, we focus also on the effect of having different stellar and disk masses in the dust evolution models, assuming that the disk dust mass is a proxy of stellar metallicity. Furthermore, we investigate the effect of having pressure bumps of different amplitudes. In the context of planets being responsible for those pressure bumps, it is an open question of whether the M_dust−M_⋆ relations are, perhaps, tracing the presence of giant planets (hence stronger pressure bumps) around more massive stars, while sub-Neptune planets (hence weaker pressure bumps) around low-mass stars, in agreement with exoplanet statistics. In this picture, planets have already formed at the very first million years and the millimeter disk mass is the leftover mass, which does not drift inward and can be still observed by ALMA.

This paper is organized as follows. In Sect. 2, we explain the dust evolution models and the assumptions made when comparing them with observations. In Sect. 3, we present the results of these dust evolution models and the comparison with current millimeter observations of planet-forming disks, in particular in the context of the M_dust − M_⋆ relation. In Sect. 4, we discuss our results and assumptions, and the main conclusions are given in Sect. 5.

2 Methods

2.1 Dust evolution models

To study theeffect of pressure traps in the observed M_dust−M_⋆ relation, we performed dust evolution models that include the growth and dynamics of particles, based on the models from Birnstiel et al. (2010b).

In thesemodels, all the particles are initially micron-sized grains distributed as the gas with a dust-to-gas ratio of 1/100. We assumed a maximum grain size of 2 m, and the grain size grid was logarithmically spaced in 180 cells. For the growth process, we took into account the coagulation, fragmentation, and erosion of particles; while for the dynamics, particles are influenced by the gas drag and, hence, the radial and azimuthal drift, the turbulent mixing, vertical settling, and Brownian motion. Particles stick and grow when their relative velocities before collisions are below a threshold, which is set in the simulations to be 10 m s⁻¹. This value is based on numerical and laboratory experiments of icy particles, that show that above these velocities particles are expected to fragment (e.g., Blum & Wurm 2000; Blum 2018; Paszun & Dominik 2009; Gundlach & Blum 2015), although recent laboratory experiments show that ice particles may be as weak as silicate particles (e.g., Gundlach et al. 2018), with fragmentation velocities of ~1 m s⁻¹. The dust diffusion is assumed to be as the turbulent gas viscosity (Youdin & Lithwick 2007) and the turbulent velocities are proportional to $\sqrt{\alpha }$, which is the effective viscosity parameter for the disk evolution (Shakura & Sunyaev 1973). In the models, this parameter α has an influence, therefore, on the relative velocity of particles (and their fragmentation) due to turbulence and settling, as well as on the diffusion or mixing of grains.

For the models, we assumed different stellar parameters and disk masses and fixed other disk parameters. This approach is aimed at achieving a better understanding of the effect of stellar and disk mass on the radial drift and trapping of particles and, therefore, on the final M_dust that would be observed. Because radial drift is more efficient around low-mass stars (Pinilla et al. 2013; Zhu et al. 2018), we expect to obtain dust masses that decrease with decreasing stellar mass. To visualize this, we consider four different type of stars: (1) a very low-mass star with a mass of M_⋆ = 0.1 M_⊙, a luminosity of L_⋆ = 8 × 10⁻² L_⊙, and an effective temperature of 3000 K. These properties are motivated by the CIDA 1, a low-mass star whose disk has a large inner cavity resolved at millimeter wavelengths (Pinilla et al. 2018); (2) an intermediate-mass star with M_⋆ = 0.3 M_⊙, with a luminosity of L_⋆ = 0.2 L_⊙, and an effective temperature of 3700 K; (3) a solar-type star, that is, with M_⋆ = 1 M_⊙, L_⋆ = 2 L_⊙, and an effective temperature of 4700 K. And finally, (4) a Herbig star with M_⋆ = 2.0 M_⊙, a luminosity of L_⋆ = 20 L_⊙, and an effective temperature of 9000 K. In the first set of simulations, we keep the disk mass the same independently of the stellar mass, such that M_disk = 0.05 M_⊙. These models aimed to investigate the effect of radial drift around different stellar masses only. In the second set of simulations, we assumed M_disk = 0.05 M_⋆. According to protostellar disks formed in radiation hydrodynamical simulations of star-cluster formation, the disk mass is expected to be around 0.01–0.08 M_⋆ for ages higher than 10 000 yr (Bate 2018)

The temperature of the dust is assumed to be a power law that depends on the stellar luminosity, such that, $$T_{\text{d}}(r)=T_{10}{\left(\frac{r}{10\text{\hspace{0.17em}au}}\right)}^{-1/2}{\left(\frac{L_{\star }}{L_{\odot }}\right)}^{1/4},$$(1)

where T₁₀ is taken to be 30 K (Andrews et al. 2013; Tripathi et al. 2017) and L_⋆ is the stellarluminosity. Equation (1) is expected if the dust disk size is independent of stellar mass (Hendler et al. 2017).

Since the disk temperature is a function of stellar luminosity, we expect significant differences in some physical quantities across the range of stellar masses that we study. First, at a given distance, the sound speed is higher around more massive stars, making the disk scale height also higher (h(r) = c_sΩ⁻¹, with c_s being the sound speed and Ω the Keplerian frequency). Second, in dust evolution models, the relative velocities between particles are dominated by turbulence and radial drift. A change in temperature affects the maximum value of turbulent velocities, which are directly proportional to the sound speed, $v_{\text{turb}}\propto \sqrt{{\alpha }_{\text{turb}}} c_{\text{s}}$ (Ormel & Cuzzi 2007).

The disk surface density distribution is parametrized by an exponentially tapered power-law function, given by $${\Sigma }_{\text{gas}}(r)={\Sigma }_0{\left(\frac{r}{R_{\text{c}}}\right)}^{-\gamma }\exp \left[-{\left(\frac{r}{R_{\text{c}}}\right)}^{2-\gamma }\right],$$(2)

where γ = 1, R_c = 80 au, and Σ₀ is taken such that the disk mass is either 0.05 M_⊙ or 0.05 M_⋆. The radial grid is taken from 1 to 300 au and it is logarithmically spaced in 300 cells.

Particle traps

To include the pressure bumps in the disk, we assume Gaussian perturbations to the disk surface density given by $$B(r)=A\exp \left(-\frac{{(r-r_{\text{p}})}^2}{2w^2}\right),$$(3)

where A is the amplitude, r_p the center, and w the width of the Gaussian perturbation. We assume three pressure bumps in the disk motivated by the average amount of substructures observed in protoplanetary disks at high angular resolution (Huang et al. 2018; Long et al. 2018) and claimed to be particle traps (Dullemond et al. 2018). Hence, the gas surface density of the perturbed disk is given by $${\Sigma }_{\text{gas}}^{\prime }(r)={\Sigma }_{\text{gas}}(r)\times [1+B_1+B_2+B_3],$$(4)

where B₁, B₂, and B₃ are given by Eq. (3). We consider two different types of perturbations. One type we refer to as strong pressure bumps, with A = 4; and another that we call weak pressure bumps, with A = 1. This simple parametric form aims to mimic the presence of multiple giant planets or (sub-) Neptune planets creating pressure bumps exterior to their orbits. We note that, in reality, A depends notonly on the planet mass, but also on several disk parameters, such as viscosity and local scale height (e.g., Crida et al. 2006; Fung et al. 2014). When comparing the assumed gas surface density with the hydrodynamical simulations from Zhang et al. (2018), who studied three different disk scale heights and three different values of α for 4 different planet masses, the amplitude of our perturbations resemble planet masses of 33 M_Earth and 0.3 M_Jup in the case of α = 10⁻³ and h∕r = 0.05 when A = 4 and A = 1, respectively; or 0.3 M_Jup and 1 M_Jup when h = 0.07 − 0.1 (for A = 4 and A = 1, respectively). For this α value and planet masses, the gap width does not change significantly with the planet mass, so our assumption of keeping the same width independently of A is fair for this range of planet masses and disk scale heights.

In both cases, the Gaussian perturbations are centered at 20, 40, and 80 au; and the width is the local scale height h(r). Because we assumed different stellar masses and luminosities (and, hence, different disk temperatures), the local scale height is different for each case. However, we set the width of the perturbation as the disk scale height that corresponds to the case of a very low mass star; this means that for all cases, w has values of ~2, ~5, and, ~11 au at 20, 40, and 80 au, respectively. The width of the perturbation is, therefore, larger than the local scale height for the other three stellar masses, becoming almost a factor of two for the Herbig disk. Although the Herbig disk is warmer and, hence, it is expected that the disk scale height is higher at a given location, the dependence on the stellar mass makes that the disk scale height increases forthe low mass stars. Hence, the widths are larger or equal to the pressure scale height for all cases to ensure the stability of the pressure bumps (e.g., Pinilla et al. 2012; Dullemond et al. 2018). The width of the perturbation should not be confused with the width of the dust rings (see e.g., DSHARP) as the latter can be much smaller than the width of the pressure bump and of the local scale height.

For the two types of pressure bumps, the perturbed gas surface density is normalized such that the disk mass remains as the unperturbed cases, that is, 0.05 M_⊙ or 0.05 M_⋆. Figure 1 shows the gas surface density profile for the case without any traps (unperturbed density) and the cases with strong and weak pressure traps when M_disk = 0.05 M_⊙. Table 1 summarizes the model parameters considered in this work.

Fig. 1 Assumed profiles for the gas surface density distribution when the disk mass is assumed to be 0.05 M⊙ for all cases.

2.2 Synthetic dust mass

In order to compare the results of the models with observations, we took the dust density distribution at different times of evolution and calculated the vertical optical depth (τ_ν) at 0.87 mm, which corresponds to ALMA Band 7, the band adopted for most of the medium-resolution surveys of star-forming regions.

To calculate τ_ν, we took the vertically integrated dust density distribution σ(r, a) at a given time of evolution and calculate τ_ν = σ(r, a)κ_ν∕cosi, where κ_ν is the opacity for each grain size and at a given frequency or wavelength (κ_ν). To obtain κ_ν, we use Mie theory and we assume the same volume fractions and optical properties for the dust particles as in Ricci et al. (2010). Furthermore, we consider face-on disks only, that is, i = 0. This assumption has no effects for optically thin discs, but it can lead to overall higher disc fluxes and thus higher predicted masses if discs are close to being optically thick. With the optical depth, the intensity profile at a given wavelength is given by $$I(r)=B_{\nu }(T(r))\left[1-e^{-{\tau }_{\nu }(r)}\right].$$(5)

We obtain the total flux as $$F_{\nu }=\frac{2\pi \cos i}{d^2}{\displaystyle {\int }_{r_{\text{in}}}^{r_{\text{out}}}I}(r)r\text{d}r,$$(6)

where d is the distance to the source, which we set to 140 pc, and B_λ(T(r)) is the Planck function, for which we take the temperature profile of Eq. (1).

To obtain the dust mass as is usually calculated from the millimeter fluxes, we assume optically thin emission (Hildebrand 1983) and, hence, $$M_{\text{dust}}\simeq \frac{d^2{F}_{\nu }}{{\kappa }_{\nu }{B}_{\nu }(T(r))}.$$(7)

For κ_ν, we assume what is usually taken in disk surveys, which is a frequency-dependent relation given by κ_ν = 2.3 cm² g⁻¹ × (ν/230 GHz)^0.4 (Andrews et al. 2013), and for B_ν we assume a blackbody surface brightness corresponding to 20 K.

Fig. 2 Dust density distribution after 1 Myr of evolution for no traps (left panel), weak traps (middle panel), and strong traps (right panel), for the case of M⋆ = 0.1M⊙ and α = 10−3. Top panels: the disk mass is Mdisk = 0.05 M⊙. Bottom panels: the disk mass is Mdisk = 0.05 M⋆.

3 Results

3.1 Stellar mass effect

In this section, we first describe the results when M_disk = 0.05 M_⊙ for all cases and focus on the effect of having different stellar mass. The top panels of Fig. 2 show the dust density distribution after 1 Myr of evolution and α = 10⁻³ for the cases of no traps (left panel), weak traps (middle panel), and strong traps (right panel) for a stellar mass of 0.1 M_⊙. This figure shows the effect of particle trapping in pressure bumps and the effect of their amplitude.

Unperturbed density

In the absence of particle traps, grains grow from micron-sized particles to the maximum grain size that is seen in the case of 0.1 M_⊙ set by fragmentation in most of the disk. When particles reach centimeter or millimeter sizes, they drift towards the star, depleting the disk of pebbles in million-year timescales. This is shown in the left panel of Fig. 3, which illustrates the effect of stellar mass on the evolution of the total mass in particles larger than 0.1 mm (hereafter, called mass in large dust). In this case, the maximum mass in large dust is quickly reached after few thousand of years of evolution and decreases with time until few large dust particles remain in the disk. The maximum of the total mass in large dust is around 50 M_⊕ (30% of the initial mass in solids) for the Herbig star and it is around 25 M_⊕ for the low-mass star, even when in all cases, the disk (and dust) mass was initially the same. This difference is a result of the growth and more efficient drift around low-mass stars. However, after ~1 Myr of evolution, the mass in large dust is almost the same for the four stellar masses (~ 20 M_⊕). At 5 Myr of evolution, the dust mass has decreased to less than 1 M_⊕.

Fig. 3 Effect of stellar mass on the evolution of dust mass of particles larger than 0.1 mm from the models with no traps (left panel), weak traps (middle panel), and strong traps (right panel). In all cases α = 10−3 and Mdisk = 0.05 M⊙.

Particle traps

When particle traps are included in the disk, particles can grow to millimeter or centimeter sizes and remain trapped in pressure maxima (middle and right panels of Figs. 2 and 3). In the case of weak traps (with A = 1.0 in Eq. (3)), the amplitude of the pressure bumps is enough tokeep a high amount of large dust in the disk for the cases of a Solar-type star and a Herbig star (2 M_⊙). The mass in large dust is around 130 M_⊕ and 110 M_⊕ for the 2 M_⊙ stellar mass at 1 and 5 Myr of evolution, respectively (Fig. 3). For the Solar type star the mass in large dust is around 85 M_⊕ and 70 M_⊕ at 1 and 5 Myr of evolution, respectively.For the 0.1 M_⊙ and 0.3 M_⊙ stars, these pressure bumps are not enough to retain millimeter- or centimeter-sized particles and after around 0.5 Myr of evolutionthe mass in large dust decreases with time reaching values below 2 M_⊕ at 5 Myr of evolution. This is because of the more efficient radial drift around low-mass stars (Pinilla et al. 2013). The radial drift of particles is a consequence of the difference between the gas azimuthal velocity and the Keplerian speed, which is higher for disks around low-mass stars. Therefore, the weak traps cannot efficiently trap particles for the two low-mass stars that we consider (middle panel of Fig. 3).

When strong pressure bumps are assumed in the dust evolution models, their amplitude is high enough to retain most of the large dust independent of the stellar mass at million-year timescales (right panels of Figs. 2 and 3). From 1 Myr of evolution, the mass in large dust range between 150 to 120 M_⊕ for the most massive star (2 M_⊙) and the lowest mass star (0.1 M_⊙), respectively. When trapping is effective, increasing the amplitude does not change significantly the mass in large dust with time. This can be seen when comparing the results from weak and strong traps, where the mass in large dust in the cases of disks around 1 and 2 M_⊙ is very similar independently of the amplitude of the pressure bumps (middle and right panel of Fig. 3). Therefore, with an amplitude of A = 4, we reach a saturation in trapping for all the stellar masses (see also Fig. 6 in Pinilla et al. 2012).

Due to the accumulation of dust particles in these strong pressure maxima, the disk emission becomes optically thick in several disk locations. In addition, dust grows as large as the maximum grain size that we allow in the simulations (right panels of Fig. 2). Both effects result on an underestimation of the dust disk masses when they are calculated from sub-millimeter fluxes, as discussed in Sect. 3.3.

Effect of the stellar mass on the fragmentation barrier

The maximum size that particles can reach before they fragment, or the fragmentation barrier, is given by (Birnstiel et al. 2012): $$a_{\text{frag}}=\frac{2}{3\pi }\frac{\Sigma }{{\rho }_{\text{s}}{\alpha }_{\text{turb}}}\frac{v_{\text{frag}}^2}{c_{\text{s}}^2}.$$(8)

When the disk surface density and α−viscosity are the same in the simulations, the effect of changing the stellar mass and luminosity is reflected in the sound speed. For a more luminous star, that is, with a warmer disk, the sound speed is higher, increasing the maximum turbulent velocity and therefore decreasing a_frag. As a result, in a warmer disk, particles reach lower sizes. Figure 4 shows the case of strong pressure bumps in the case of a disk around a Solar type star. The maximum grain size that particles can reach, is lower than in the case of a 0.1 M_⊙ stellar mass (top right panel of Fig. 2), an effect that is more clear inside the pressure bumps. In the case of a Solar-type star, more millimeter or centimeter grains remain in the disk at million-year timescales because fragmentation keeps the particles smaller and closer to the sizes where they are more affected by the gas drag, which means that they are easier to trap in pressure bumps. This effect will influence the synthetic sub-millimeter fluxes that are derived from the models because when growth is very efficient (as in the case of strong bumps and low-mass stars), the very large grains (≳ 10 cm) that are formed inside the pressure maxima have very low opacities contributing very little to the millimeter fluxes.

Fig. 4 Dust density distribution after 1 Myr of evolution for the case of strong traps, M⋆ = 1.0 M⊙, Mdisk = 0.05 M⊙, and α = 10−3.

3.2 Disk mass effect

There are two main effects of changing the disk mass to be a fraction of the stellar mass (0.05 M_⋆) instead of assuming 0.05 M_⊙ for all the cases. First of all, the fragmentation barrier (a_frag, Eq. (8)) is directly proportional to the gas surface density (Σ), implying that a_frag is lower in low-mass disks even though the sounds speeds are lower. On the other hand, the radial drift of particles comes from the gas drag forces, which depend on dust properties such as their mass and cross-sectional area, but also on the gas disk via the gas density. The dimensionless Stokes number (St) measures the importance of drag force by comparing the stopping time of the particles to the Keplerian frequency at a given radius. Particles with the Stokes number equal unity are subject to the highest radial drift. At the midplane, the Stokes number is inversely proportional to the gas surface density and as a consequence, a millimeter-sized particle in the outer disk has a higher Stokes number in alower mass disk. Therefore, the drift barrier is also lower in lower mass disks, which means that particles grow to smaller sizes in a low-mass disk, but also a particle with a given size drifts faster in a low-mass disk.

Figure 2 compares the case of 0.1 M_⊙ stellar mass when M_disk = 0.05 M_⊙ vs. M_disk = 0.05 M_⋆ and shows how for the case when M_disk = 0.05 M_⋆, the maximum grain size is lower and due to more efficient drift of millimeter- or centimeter-sized particles, the disk is more depleted of these kind of grains. As result, in the cases with traps, there is less effective trapping and fewer millimeter grains remain inside pressure maxima.

3.3 Comparison with the observations: M_dust−M_⋆ relation

In order to compare the results from dust evolution models with current millimeter observations of protoplanetary disks, in particular in the context of the M_dust−M_⋆ relation, we used the values for the slope and intercept for different star-forming regions from Pascucci et al. (2016) as stellar masses were homogeneously derived in this study, except for σ-Orionis, which was not included; hence, for this parameter, we took the values from Ansdell et al. (2017). In all regions, T =20 K was assumed (see solid lines in different colors in Fig. 5). We also collected from the literature a sample of disks that have been observed with ALMA, which have either a large inner cavity or gap, defined as transition disks (TDs, empty circles), and disks with substructures, such as rings, gaps, or asymmetries, but still exhibiting substantial emission from the inner disk (full triangles). We selected the TDs that have been observed with ALMAand whose cavities are resolved, which implies that the TDs in our sample have large cavities (≳20 au). We gathered a total of 73 disks (43 TDs and 30 disks with substructures, see Table A.1). The sample of TDs double the number of targets presented previously in Pinilla et al. (2018).

We placed the TDs and disks with substructures in the M_dust−M_⋆ relation and fit these points with the following relation: log₁₀(M_dust∕M_⊕) = βlog₁₀(M_⋆∕M_⊙) + α. We obtained $\beta ={0.72}_{-0.04}^{+0.04}$ and $\alpha ={1.86}_{-0.01}^{+0.01}$, and a Pearson coefficient of 0.50. The fit takes into account the uncertainties of the data, dominated by the 10% of uncertainty from flux calibration. When only the TDs are considered for the fit, the Pearson coefficient is higher (0.67) and, therefore, this correlation is stronger when only TDs are included. This is because for this correlation, we mainly have data for TDs around low-mass stars. If only disks with substructures are considered, the correlation is weak with a Pearson coefficient of 0.28. However, more observations in the low stellar mass regime that include the resolution of substructures are needed to confirm if there is indeed a correlation.

The values for the slope and intercept change when including or excluding disks with substructures. When we only consider TDs, the values remain similar with ${\beta }_{\text{TD}}={0.94}_{-0.06}^{+0.06}$ and ${\alpha }_{\text{TD}}={1.80}_{-0.02}^{+0.02}$, but the relation is shallower when only disks with substructures are considered, with ${\beta }_{\text{SB}}={0.42}_{-0.1}^{+0.1}$ and ${\alpha }_{\text{SB}}={1.90}_{-0.02}^{+0.02}$. With the current data, there is not any substantial evidence that the population of TDs is different than the population of disks with substructures based on the dust disk mass. To compare the two samples, we calculate the Kolmogorov-Smirnov (KS) two-sided test, and we find a maximum deviation between the cumulative distributions of only 0.15 with a high probability (~80%) that the two samples are similar. Overall, the M_dust−M_⋆ correlation of the whole sample is dominated by TDs, and it is shallower than seen from all disks in nearby star-forming as found in Pinilla et al. (2018), to which some of these disks actually belong to.

Next, we overplotted the M_dust obtained from the models at 1 and 5 Myr of evolution with the observed correlations. We plotted M_dust that were directly obtained from the dust evolution models for large dust without any post-processing assumption (as e.g., the values in Fig. 3) and, after dust temperature and opacities were assumed, we calculated the dust mass expected from millimeter observations as explained in Sect. 2.2 (hereafter called synthetic dust mass). Figure 6 shows the results for all the models that assumed M_disk = 0.05 M_⊙ and Fig. 7 shows the results for the cases where M_disk = 0.05 M_⋆.

Fig. 5 Mdust−M⋆ relation in different star-forming regions (colors) and for the disks with resolved cavities or transition disks (TDs, open circles) and disks with structures (full triangles, see appendix). The values for the slope and intercept for the different star forming regions are taken from Pascucci et al. (2016), except for σ-Orionis, which values are from Ansdell et al. (2017). From fitting this relation for TDs and disk with structures, we used log10 (Mdust∕M⊕) = βlog10(M⋆∕M⊙) + α, and we obtained $\beta ={0.72}_{-0.04}^{+0.04}$ and $\alpha ={1.86}_{-0.01}^{+0.01}$ and a Pearson coefficient of 0.50. The fit takes into account the uncertainties of the data, which are dominated by the 10% of uncertainty from flux calibration.

Fig. 6 Comparison with observations, models with α = 10−3 and Mdisk = 0.05 M⊙ in all cases. Top panels: total mass of the large dust directly obtained from models. Bottom panels: synthetic dust mass after considering assumptions from observations.

3.3.1 Synthetic dust mass compared to mass in large dust from models

In general, the synthetic dust mass is lower than the mass in large dust that is obtained directly from the models (Figs. 6 and 7). This difference is because of two different reasons. First, some parts of the disk become optically thick. In the case of an unperturbed density profile, the optical depth of the disk (τ) reaches values of τ~ 0.1 in the case of low-mass stars (<M_⊙), and τ ~ 0.5 in the case of higher stellar masses at 1 Myr. The higher τ around more massive stars, even when the disk mass is assumed to be the same independently of the stellar mass, is due to the drift more efficient removal of large particles in the case of low-mass stars. At 5 Myr of evolution, the optical depth is much lower than unity in all the cases without particle traps, while the synthetic dust mass and the mass in large dust from the models are fairly in agreement.

In the case of particle traps, the optical depth can reach values of between 1 to 10 inside the bumps. For the weak pressure bumps, these high values for τ are reached only in the innermost bump, while in the case of strong bumps, these values are reached at the peak of all the bumps. In the cases of particle traps, a second reason plays an important role to understand why the synthetic dust mass is lower than the mass in large dust from models, and it is due to the growth of solids beyond 10 cm, hereafter referred to as boulders (inthe minimum mass solar nebula model, the (sub-)meter-sized objects that are in between being strongly influenced by the gas and particles that are completely decoupled and moving with Keplerian speed are called “boulders”).

Boulders have very low opacities and their contribution to the millimeter fluxes is insignificant, which means that their mass is hidden when calculated from the millimeter-fluxes. This effect is clear when comparing the masses in the case of strong bumps. For low-mass stars (0.1 or 0.3 M_⊙), boulders areformed inside the bumps as discussed in the previous section (Fig. 2), while in the case of disks around 1 or 2M_⊙, the fragmentation barrier keeps the maximum grain size around millimeter and centimeter sizes in all the pressure bumps (Fig. 4). The formation of boulders reduces the optical depth inside the pressure bumps, but simultaneously, a lot of the solid material becomes invisible when observed at millimeter-wavelengths. As a consequence, while the synthetic dust mass is lower by a factor of ~ 1.5−2 for the disks around stars M_⋆ ≥ M_⊙, the synthetic dust mass is a factor ~10−100 lower than the total mass in large dust for M_⋆ < M_⊙. Hence, in the case of traps, the difference between the synthetic dust mass and the mass from the models originates from a combination of these two effects: optical thickness and boulder formation. For the weak traps, these two effects appear only in the innermost bump at 1 Myr of evolution. While for strong traps, these two effects appear in all the pressure traps at ≥ 1 Myr.

In summary, both optically thick regions and boulder formation lead to a large difference between the synthetic dust mass and the mass in large dust obtained directly from the models. Our result show that boulder formation has a more significant influence. Recently, Stammler et al. (2019) showed that to explain the optical depth between 0.2 and 0.5 inside the rings observed in the DSHARP sample (Andrews et al. 2018a), boulder and planetesimal formation must occur inside the pressure bump to reduce the optical depth. In our models, we do not reduce the mass of the dust inside the pressure bumps by assuming that boulders were formed when high dust-to-gas density ratio is reached, but rather particles naturally grow to sizes beyond decimeter when drift and turbulence are inefficient to fragment or keep the particles with millimeter to centimeter size.

Fig. 7 As in Fig 6, but with Mdisk = 0.05 M⋆ in all cases.

3.3.2 Synthetic dust mass versus observations: M_disk = 0.05 M_⊙ case

In the less realistic case of disk masses independent of stellar masses, the observed slope in different star forming regions is not reproduced by the models, not when comparing directly with the models nor after post-processing to compare with the observations (see Fig. 6). Only in the case of weak traps, the slope of the sample of TDs and disks with substructures is reproduced when comparing the mass directly from the models at 1 Myr of evolution. When comparing it with the synthetic dust mass, the slope of the TDs and disks with substructures is similar to the models with and without traps, but values are lower than observed, meaning that the initial dust mass should be higher for all stellar masses to have an agreement with the relation of TDs and substructured disks, an assumption that may be unlikely since with the assumed disk mass of 0.05 M_⊙, the stellar-disk mass ratio is already very high forthe disks around 0.1 M_⊙ (stellar-disk mass ratio of 0.5) and 0.3 M_⊙ (stellar-disk mass ratio of 0.17).

Strong pressure bumps are needed to produce the high disk dust masses that are observed around stars with masses equal or higher than the Solar mass (when comparing with the synthetic dust mass). With weak or no pressure bumps, the millimeter fluxes and hence the dust disk masses are lower than observed already at 1 Myr of dust evolution. In the context of planets causing these pressure bumps, this implies that if the disk mass is initially the same for all disks around different stellar masses, only the strong pressure bumps (or giant planets) can explain the observed values of M_dust for stellar masses higher than a solar mass. For low-mass stars no conclusion can be given besides that in all scenarios the synthetic dust mass is overestimated at 1 Myr. At 5Myr, the synthetic dust mass values are in agreement with the values in Lupus, Taurus and ChaI, although the averaged age of these regions is lower than 5 Myr.

3.3.3 Synthetic dust mass versus observations: M_disk = 0.05 M_⋆ case

In the more realistic case of disk masses scaling with stellar masses, the unperturbed density can explain the observed M_dust−M_⋆ relation for different star forming regions only at 1 Myr when compared directly with the mass in large dust obtained from the models (see Fig. 7). When comparing with the synthetic dust mass predicted by our model, the slope at 1 Myr is slightly lower than in the actual observations. We find M_dust values that are too low for the case of 1 and 2 M_⊙ at 5 Myr (both dust mass: directly from models or the synthetic dust mass as inferred from the simulations).

In the case of weak traps, we obtained two interesting results. First, when comparing it with the mass in large dust directly obtained from the models, the M_dust−M_⋆ relation of TDs and disks with substructures is reproduced at 1 Myr of evolution. The sample of TDs and disks with substructures (Table A.1) is a combination of disks at different stages of evolution, from 0.5 to 10 Myr, but most of them are older than 1 Myr (Garufi et al. 2018). Second, due to the reasons explained above (optical thickness and boulder formation), when comparing it with the synthetic dust mass, we find that the model becomes steeper and closer to the slope of different star forming regions, especially at 1 Myr. At 5 Myr, the synthetic dust mass is lower than in any star forming region, in particular for disks around 1 and 2 M_⊙ stellar masses.

The case of strong traps shows very interesting results but they must be interpreted with caution. First, the M_dust−M_⋆ relation observed in TDs and disks with substructures is fairly reproduced at any time of evolution when comparing with the mass in large dust obtained directly from the models. However, due to boulder formation, in particular for the 0.1 and 0.3 M_⊙, the synthetic dust mass is much lower reproducing again the steeper trends at any time of evolution. This implies that to explain the trend of TDs and disks with substructures, strong pressure bumps (perhaps due to giant planets) are needed, but in addition, boulder formation must be inhibited inside these pressure bumps, otherwise the observed fluxes and masses are too low, reproducing instead the steeper relations in different star-forming regions.

Therefore, in the context of planets causing these pressure bumps, we found similar results as in the case of M_disk = 0.05 M_⊙ for the disks around 1 and 2 M_⊙, meaning those for which only strong pressure bumps (or giant planets) can explain the observed values of dust mass at any time of evolution. Contrary to the case of M_disk = 0.05 M_⊙, any of the three scenarios (no traps, weak or strong traps) can reproduce the M_dust−M_⋆ relations of different star formation regions for 0.1 and 0.3 M_⊙ at 1 or 5 Myr of evolution. The agreement with the strong traps scenario is surprising, but it is due to boulder formation, while in the other two cases, it is due to the actual loss of particles attributed to the drift.

4 Discussion

4.1 Balance between trapping and boulder formation

Our findings suggest that strong pressure traps can reproduce the M_dust−M_⋆ relation observed in different star forming regions (bottom right panel of Fig. 7) because of the trapping of dust particles and boulder formation inside pressure bumps, the latter in particular for M_⋆ < 1 M_⊙. In contrast, the flatter relation observed for TDs and disks with substructures has not been reproduced by any model, specifically one based on post-processing for the purpose of comparison with the observations. In the framework of our models, one solution to reproduce the relation of TDs and disks with substructures is to consider that boulder formation is inhibited inside pressure bumps in the case of low-mass stars (M_⋆ < 1 M_⊙).

Our results seem to contradict those of Stammler et al. (2019), who suggests that dust must be transformed into planetesimals to obtain the optical depths values inferred inside the DSHARP rings (Andrews et al. 2018a), using dust evolution models that include trapping in a pressure bump. However, that work was done assuming a Herbig star, in particular to reproduce the observations of one of the brightest disk (HD 163296), while our conclusion about hindering boulder formation inside pressure bumps applies only to low-mass stars (M_⋆ < 1 M_⊙). In our case of a disk around a Herbig star, we do not have any natural production of boulders in pressure maxima and the small differences (a factor of ~1.5−2) between the mass in large dust directly from the models and the observed total dust mass originate from the high optical depth (~1–2) inside the pressure maxima. According to our results, inducing boulder formation in the pressure bumps in the Herbig case can help to reduce τ and be in agreement with the values obtained in the DSHARP sample. However, as found in our models, boulder formation reduces by a lager factor the synthetic dust mass (~10−100, comparison between right panels of Fig. 7). Assuming a higher initial disk mass may help to have a better agreement with the synthetic dust mass as in Stammler et al. (2019), where the initial condition for the disk mass is 0.4 M_⊙ (four times than what we assume in this work) and they do reproduce the intensity of one of the observed rings of HD 163296.

There are several possibilities that may hinder boulder formation in pressure maxima in disks around low-mass stars. The first possibility is that the velocity threshold for which particles fragment (fragmentation velocity) is lower than what we considered (10 ms⁻¹). Laboratory experiments and numerical simulations show that these values are possible (e.g., Paszun & Dominik 2009; Gundlach & Blum 2015; Musiolik et al. 2016) for water ice grains, although recent experiments suggest that the fragmentation velocities for water ice grains may be as low or lower as for silicates, that is, 1 ms⁻¹ (e.g., Musiolik & Wurm 2019). However, if this is the case, the fragmentation barrier would be two orders of magnitude lower (Eq. (8)), which will imply that the maximum grain size would be set by fragmentation with typical values of ~100 μm in the outer disk (>20 au). These grain sizes are in agreement with the values proposed by Liu (2019) and Zhu et al. (2019), who suggest that the measured optical depth in several DSHARP disks can be explained by dust scattering when the optically thick dust has an albedo of ~0.9 at 1.25 mm, which corresponds to a maximum grain size of ~0.1–1 mm. Nevertheless, if this is the case, dust trapping may be inefficient in reproducing observations. Investigating the required conditions for particle trapping, and an agreement with observations when the fragmentation velocities are 1 ms⁻¹ in the entire disk, is a subject of future research.

Other possibilities that may explain the inhibition of the formation of boulders include bouncing in the models or an increase in the disk α viscosity. Nevertheless, to include bouncing would keep the grains too small (Stokes numbers lower than 10⁻⁴, e.g., Zsom et al. 2010) such that they could not effectively be trapped in pressure maxima. The other alternative is to have a higher value of the α viscosity, which sets the fragmentation barrier and turbulence mixing strength; hence, determining how much grains can be diffused out from pressure bumps, which is discussed in the next subsection.

Fig. 8 Effect of the turbulence mixing strength on the evolution of the dust mass of particles larger than 0.1 mm from the models with no traps (left panel), weak traps (middle panel), and strong traps (right panel). In all the cases, M⋆ = 0.1 M⊙ and Mdisk = 0.05 M⊙.

4.2 Effect of turbulence mixing strength

As discussed in the previous subsection, increasing the turbulence mixing strength may help to hinder boulder formation inside pressure bumps. In general, changing the value of the α viscosity in our models has a strong effect in the results. This effect has been studied and discussed in detail in de Juan Ovelar et al. (2016) and Bae et al. (2018), who show that by just varying the values of α and keeping the rest of the model parameters the same, the synthetic images from models can have very different outputs, such as: multiple rings and gaps, a cavity with a single ring, or compact disks.

In the case of planets, usually a higher mass planet is required to open a gap when the disk viscosity is higher. This balance between planet mass and viscosity has an effect in the potential trapping of particles (Fig. 1 in de Juan Ovelar et al. 2016). The amplitude of A = 4 of our perturbation in Eq. (3) corresponds to planet masses similar to 0.1–1 M_Jup for 10⁻³−10⁻² α−viscosity. We note that increasing the amplitude of our perturbation, which would resemble gaps formed by more massive planets or less viscous disks, does not have a strong effect in our models, since with A = 4 effective trapping is already reached and therefore increasing A would only saturate the amount of dust trapped in pressure maxima. In our assumptions of the particle traps, the parametric perturbation does not depend on α, and hence the shape of the pressure bumps is independent of gas viscosity. As a consequence, the only effect of varying α in our models is reflected in the dust evolution.

Varying α directly influences the fragmentation barrier a_frag ∝ 1∕α, which means that for higher α, grains reach lower sizes when the maximum grain size is set by fragmentation. In addition, we assume dust diffusion D as Youdin & Lithwick (2007), which means D ∝ α∕(1 + St²) and, thus, higher dust diffusion for higher values of α. Increasing α implies lower a_frag and higher D, which influence the potential trapping because if grains do not grow to values that attain a high Stokes number, they are more difficult to trap and they diffuse out of pressure traps easily. To see the effect of α on the mass in large dust, we simulated the case of M_⋆ = 0.1 M_⊙ and M_disk = 0.05 M_⊙ with three different values of α (10⁻⁴, 10⁻³, 10⁻²). Figure 8 shows the main results of this numerical experiment.

When the disk viscosity is high (α = 10⁻²) trapping is less efficient and it does not occur even for the case of strong traps (Fig. 8). In these cases, the peak of the mass in large dust is reached earlier in the simulations (~ 10⁵ yr). On the other hand, when disk viscosity is low (α = 10⁻⁴), the peak of the mass in large dust is reached later (~1 Myr), and in this case trapping is very efficient even for weak traps. In this case, the fragmentation barrier is one order of magnitude higher than our standard models, which implies more efficient formation of boulders inside pressure bumps, creating more invisible solid material when observed at millimeter-emission and a higher discrepancy with the observed dust mass.

4.3 Disk size-luminosity relation

Tripathi et al. (2017) and Andrews et al. (2018b) find a significant correlation between the disk size (R_eff, usually defined as the radii that encloses 68% of the millimeter emission) and the continuum luminosity of the disk. Hendler et al. (2020) performed a homogeneous analysis of five star forming regions, finding that this relation flattens for Upper Sco, which is the oldest region in the sample (left panel of Fig. 9). The observed relation could originate from the growth and radial drift of dust particles, where the dust disk extension becomes smaller with time due to drift. Depending on the initial conditions, dust evolution models can reproduce the observed relations (Tripathi et al. 2017; Rosotti et al. 2019). It is, therefore, possible that the R_eff−L_mm flattens for disks in which radial drift is reduced or ineffective.

We checked the R_eff−L_mm relation for our sample of TDs and disks with substructures and found a significan correlation (Pearson coefficient of 0.57) with an intercept (A) equal to 2.29 ± 0.02 and a slope (B) equal to 0.36 ± 0.03, assuming log₁₀(R_eff∕au) = B log₁₀(L_mm∕mJy) + A. For this calculation, the millimeter luminosity is scaled to a 140 pc distance for all the disks. The values for R_eff are the ones that enclose 90 or 95% of the millimeter flux in Band 6 or Band 7 with ALMA (the wavelength is between ~0.8 and ~1.3 mm). The relations found in Hendler et al. (2020), Tripathi et al. (2017) and Andrews et al. (2018b) used the radius that encloses the 68% of the flux, but Hendler et al. (2020) demonstrates that there is a linear relation close to 1–1 between the two values.

The values for R_eff of the sample are larger than the averaged values found in different star forming regions, with a mean value of 116 au (see Table 2 from Hendler et al. (2020) for comparison). This is because most of the ALMA observations of TDs and disks with substructures are biased towards disks with larger radii. The slope of R_eff−L_mm in our sample is smaller than the values found in the youngest star forming regions studied in Hendler et al. (2020), that is, Ophiucus, Lupus, and Taurus/Auriga with ages between 1 and 3 Myr. The slope of our sample lies in between the values found for Chamaeleon I (2–3 Myr, slope of 0.4) and Upper-Sco (5–1 Myr, slope of 0.22), while the ages of the TDs and disks with substructures in our sample range between 0.1 and 10 Myr.

Based on our dust evolution models, we calculate the radius that encloses 95% of the mass of the dust particles that are larger than 0.1 mm. The evolution of this radius is shown in the right panel of Fig. 9 for the models that assume strong pressure bumps, M_disk = 0.05 M_⋆, and α = 10⁻³. Independent of the stellar mass, after ~0.5 Myr of evolution the value of this radius converges to the location of the farther pressure maximum, that is, ~ 75 au (similar values are found when the disk radius is calculated from the synthetic intensity profile). Based on these models, we expect a flat relation between R_eff and L_mm at any evolutionary stage. It is possible that the observed R_eff − L_mm relation for TDs and disks with substructures provides information about how far pressure bumps may form and how the less luminous disks that are around stars of lower mass may have the pressure bumps closer in. Exoplanet statistics suggest that low-mass stars have more planets with radii of 1–4 R_⊕ at the same semi-major axis (e.g., Mulders et al. 2015) but to test this idea, we require more detailed models to investigate different locations and the morphology of pressure bumps, along with a greater number of high-resolution observations around low-mass stars.

Fig. 9 Left panel: fit of the disk size-luminosity relation from sub-millimeter observations of protoplanetary disks in differentstar forming regions (colors) as reported in Hendler et al. (2020) compare to the relation for TDs (open circles) and disks with substructures (full triangles). Right panel: evolution of the radius that encloses 95% of the mass of dust particles larger than 0.1 mm from the models that assume strong pressure bumps, Mdisk = 0.05 M⋆, and α = 10−3.

5 Conclusions

We investigated the effect of dust evolution and particle trapping in the observed M_dust−M_⋆ relation of different star-forming regions. We focused our attention on the effect of the stellar and disk mass. Both properties affect the mass in large dust at million-year timescales due to: (1) the more effective drift around low-mass stars; (2) the more efficient fragmentation for disks around stars of higher mass, and it is less efficient in more massive disks; and (3) the drift being more efficient in disks of lower mass. The combination of these effects influence the efficiency of trapping and the potential to form boulders in pressure bumps, which influence the final M_dust. Our conclusions are:

• 1.

  Independent of the assumption for the disk mass, either M_disk = 0.05 M_⊙ or M_disk = 0.05 M_⋆, strong traps are required to reproduce the observed M_dust of disks around M_⋆ ≥ 1 M_⊙ at different time of evolution. In the context of planets as the origin of pressure bumps, this result agrees with current exoplanet statistics about giant planets, which are found to be more common around more massive stars.

• 2.

  When alldisks have the same initial mass independent of the stellar mass (M_disk = 0.05 M_⊙), the dust evolution models produce a flatter M_dust−M_⋆ relationwhen compared to the observed ones from different star-forming regions. Therefore, we conclude that the initial disk mass has to scale with the stellar mass to reproduce the slope of the M_dust−M_⋆ relation, otherwise it is flatter than what is observed. In the case of M_disk = 0.05 M_⊙ for all cases, the slope obtained from the models is similar to the one observed for TDs and disks with substructures, however the modeled M_dust is overall lower than the observed one.

• 3.

  When the disk mass is a fraction of the stellar mass (5%), our models reproduce the observed M_dust−M_⋆ relation from different star-forming regions. However, in the cases of an unperturbed density and weak traps, the trends are reproduced only at 1 Myr while at later times (5 Myr), M_dust is too low when compared to observations for M_⋆ ≥ 1 M_⊙.

• 4.

  Strong pressure traps can reproduced the observed M_dust−M_⋆ relations of different star forming regions in the case of M_disk = 0.05 M_⋆. This result arises from dust trapping and dust growth beyond centimeter-sizes inside pressure traps, an effect that decreases the millimeter fluxes and reduces the expected M_dust, in particular for M_⋆ < 1 M_⊙.

• 5.

  To explain the flatter relation of M_dust−M_⋆ for TDs and disks with substructures, strong pressure bumps (perhaps caused by giant planets) are needed. However, efficient growth must be inhibited inside these pressure bumps for low-mass stars, otherwise the observed fluxes and masses are too low, reproducing instead the steeper trends in different star-forming regions. Different possibilities can hinder very efficient growth inside pressure bumps, such as a lower fragmentation velocity or including the effect of bouncing.

• 6.

  Due to the efficient formation of decimenter bodies in pressure bumps for disks around low-mass stars M_⋆ < 1 M_⊙, we cannot give a definitive conclusion on what type of traps are present in these disks. Indeed, all three cases we tested (unperturbed, weak, and strong traps) can reproduced the observed values in the more realistic case of disk mass scaling with stellar mass.

• 7.

  The R_eff−L_mm relation for TDs and disks with substructures is flatter than observed in Ophiucus, Lupus, and Taurus/Auriga star-forming regions and the slope value of our sample lies in the range found for older regions, that is, Chamaeleon I and Upper-Sco. This relation may flatten due to inefficient radial drift when pressure bumps exist in the disks, in which case R_eff traces the location of the farther pressure bump at any time of evolution. Currently, high angular resolution observations of TDs and disks with substructures are biased toward the brightest and larger disks, so more observations of small and faint ones are required to test this idea.

Future multi-wavelength and high-resolution observations of disks around low-mass stars are needed to investigate if gaps or rings are present and how much grains have growth there in comparison to the rings and gaps in disks around more massive stars. Our results suggest that boulder formation must be inhibited inside pressure bumps in TDs and disks with substructures in order to explain the observed M_dust. Further theoretical research is needed to investigate what conditions are required to explain the observed properties of this group of disks, in particular if boulder formation is hindered when the fragmentation velocities are lower than assumed in this work – as has been suggested by recent laboratory experiments.

Acknowledgements

The authors are thankful to the referee for the constructive referee report, to F. Long for all the enthusiastic discussions about the results of this paper, and to N. van der Marel for sharing some data. P.P. acknowledges support provided by the Alexander von Humboldt Foundation in the framework of the Sofja Kovalevskaja Award endowed by the Federal Ministry of Education and Research. I.P. acknowledges support from an NSF Astronomy & Astrophysics Research Grant (ID: 1515392). This material is based on work supported by the National Aeronautics and Space Administration under Agreement No. NNX15AD94G for the program Earths in Other Solar Systems. The results reported herein benefited from collaborations and/ or information exchange within NASA’s Nexus for Exoplanet System Science (NExSS) research coordination network sponsored by NASA’s Science Mission Directorate.

Appendix A Target list of TDs and disks with substructures

References

All Tables

All Figures

	Fig. 1 Assumed profiles for the gas surface density distribution when the disk mass is assumed to be 0.05 M⊙ for all cases.
In the text

	Fig. 2 Dust density distribution after 1 Myr of evolution for no traps (left panel), weak traps (middle panel), and strong traps (right panel), for the case of M⋆ = 0.1M⊙ and α = 10−3. Top panels: the disk mass is Mdisk = 0.05 M⊙. Bottom panels: the disk mass is Mdisk = 0.05 M⋆.
In the text

	Fig. 3 Effect of stellar mass on the evolution of dust mass of particles larger than 0.1 mm from the models with no traps (left panel), weak traps (middle panel), and strong traps (right panel). In all cases α = 10−3 and Mdisk = 0.05 M⊙.
In the text

	Fig. 4 Dust density distribution after 1 Myr of evolution for the case of strong traps, M⋆ = 1.0 M⊙, Mdisk = 0.05 M⊙, and α = 10−3.
In the text

Fig. 5 Mdust−M⋆ relation in different star-forming regions (colors) and for the disks with resolved cavities or transition disks (TDs, open circles) and disks with structures (full triangles, see appendix). The values for the slope and intercept for the different star forming regions are taken from Pascucci et al. (2016), except for σ-Orionis, which values are from Ansdell et al. (2017). From fitting this relation for TDs and disk with structures, we used log10 (Mdust∕M⊕) = βlog10(M⋆∕M⊙) + α, and we obtained $\beta ={0.72}_{-0.04}^{+0.04}$ and $\alpha ={1.86}_{-0.01}^{+0.01}$ and a Pearson coefficient of 0.50. The fit takes into account the uncertainties of the data, which are dominated by the 10% of uncertainty from flux calibration.

In the text

	Fig. 6 Comparison with observations, models with α = 10−3 and Mdisk = 0.05 M⊙ in all cases. Top panels: total mass of the large dust directly obtained from models. Bottom panels: synthetic dust mass after considering assumptions from observations.
In the text

	Fig. 7 As in Fig 6, but with Mdisk = 0.05 M⋆ in all cases.
In the text

	Fig. 8 Effect of the turbulence mixing strength on the evolution of the dust mass of particles larger than 0.1 mm from the models with no traps (left panel), weak traps (middle panel), and strong traps (right panel). In all the cases, M⋆ = 0.1 M⊙ and Mdisk = 0.05 M⊙.
In the text

Fig. 9 Left panel: fit of the disk size-luminosity relation from sub-millimeter observations of protoplanetary disks in differentstar forming regions (colors) as reported in Hendler et al. (2020) compare to the relation for TDs (open circles) and disks with substructures (full triangles). Right panel: evolution of the radius that encloses 95% of the mass of dust particles larger than 0.1 mm from the models that assume strong pressure bumps, Mdisk = 0.05 M⋆, and α = 10−3.

In the text