# The ALMA Survey of Gas Evolution of PROtoplanetary Disks (AGE-PRO). VI. Comparison of Dust Evolution Models to AGE-PRO Observations

Nicolás Kurtovic, Matías Gárate, Paola Pinilla, Ke Zhang, Giovanni Rosotti, Rossella Anania, Ilaria Pascucci, Benoît Tabone, Leon Trapman, Dingshan Deng, et al.

► To cite this version:

Nicolás Kurtovic, Matías Gárate, Paola Pinilla, Ke Zhang, Giovanni Rosotti, et al. The ALMA Survey of Gas Evolution of PROtoplanetary Disks (AGE-PRO). VI. Comparison of Dust Evolution Models to AGE-PRO Observations. The Astrophysical Journal, 2025, 989 (1), pp.6. 10.3847/1538-4357/add1d0. hal-05411943

HAL Id: hal-05411943

https://hal.science/hal-05411943v1

Submitted on 17 Dec 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

# The ALMA Survey of Gas Evolution of PROtoplanetary Disks (AGE-PRO). VI. Comparison of Dust Evolution Models to AGE-PRO Observations

Nicolás T. Kurtovic<sup>1,2</sup>, Matías Gárate<sup>2</sup>, Paola Pinilla<sup>2,3</sup>, Ke Zhang<sup>4</sup>, Giovanni P. Rosotti<sup>5</sup>, Rossella Anania<sup>5</sup>, Ilaria Pascucci<sup>6</sup>, Benoît Tabone<sup>7</sup>, Leon Trapman<sup>4</sup>, Dingshan Deng<sup>6</sup>, Miguel Vioque<sup>8,9</sup>, John Carpenter<sup>10</sup>, Lucas A. Cieza<sup>11</sup>, Laura M. Pérez<sup>12</sup>, Carolina Agurto-Gangas<sup>12</sup>, Anibal Sierra<sup>3,12</sup>, Dary A. Ruiz-Rodriguez<sup>13</sup>, James Miley<sup>14,15,16</sup>, Camilo González-Ruilova<sup>11,15,16</sup>, Estephani Torres-Villanueva<sup>4</sup>, and Aleksandra Kuznetsova<sup>17</sup>

<sup>1</sup>Max Planck Institute for Extraterrestrial Physics, Giessenbachstrasse 1, D-85748 Garching, Germany; kurtovic@mpe.mpg.de

<sup>2</sup>Max-Planck-Institut für Astronomie (MPIA), Königstuhl 17, 69117 Heidelberg, Germany

<sup>3</sup>Mullard Space Science Laboratory, University College London, Holmbury St Mary, Dorking, Surrey RH5 6NT, UK

<sup>4</sup>Department of Astronomy, University of Wisconsin-Madison, 475 N Charter Street, Madison, WI 53706, USA

<sup>5</sup>Dipartimento di Fisica, Università degli Studi di Milano, Via Celoria 16, I-20133 Milano, Italy

<sup>6</sup>Lunar and Planetary Laboratory, The University of Arizona, Tucson, AZ 85721, USA

<sup>7</sup>Université Paris-Saclay, CNRS, Institut d'Astrophysique Spatiale, Orsay, France

<sup>8</sup>European Southern Observatory, Karl-Schwarzschild-Str. 2, 85748 Garching bei München, Germany

<sup>9</sup>Joint ALMA Observatory, Alonso de Córdova 3107, Vitacura, Santiago 763-0355, Chile

<sup>10</sup>Joint ALMA Observatory, Avenida Alonso de Córdova 3107, Vitacura, Santiago, Chile

<sup>11</sup>Instituto de Estudios Astrofísicos, Universidad Diego Portales, Av. Ejercito 441, Santiago, Chile

<sup>12</sup>Departamento de Astronomía, Universidad de Chile, Camino El Observatorio 1515, Las Condes, Santiago, Chile

<sup>13</sup>National Radio Astronomy Observatory, 520 Edgemont Road, Charlottesville, VA 22903, USA

<sup>14</sup>Departamento de Física, Universidad de Santiago de Chile, Av. Victor Jara 3659, Santiago, Chile

<sup>15</sup>Millennium Nucleus on Young Exoplanets and their Moons (YEMS), Chile

<sup>16</sup>Center for Interdisciplinary Research in Astrophysics and Space Exploration (CIRAS), Universidad de Santiago de Chile, Chile

<sup>17</sup>Center for Computational Astrophysics, Flatiron Institute, 162 Fifth Avenue, New York, NY 10025, USA

Received 2024 September 6; revised 2025 March 17; accepted 2025 March 23; published 2025 July 31

## Abstract

The potential for planet formation of a circumstellar disk depends on the dust and gas reservoirs, which evolve as a function of the disk age. The Atacama Large Millimeter/submillimeter Array AGE-PRO Large Program has measured several disk properties across three star-forming regions of different ages, and in this study, we compare the observational results to dust evolution simulations. Using *DustPy* for the dust evolution, and *RADMC-3D* for the radiative transfer, we ran a large grid of models spanning stellar masses of 0.25, 0.50, 0.75, and 1.0 $M_\odot$, with different initial conditions, including: disk sizes, disk gas masses, and dust-to-gas ratio, and viscosity. Our models are performed assuming smooth, weakly, or strongly substructured disks, aiming to investigate if any observational trend can favor or exclude the presence of dust traps. The observed gas masses in the disks of the AGE-PRO sample are not reproducible with our models, which only consider viscous evolution with constant $\alpha$, suggesting that additional physical mechanisms play a role in the evolution of the gas mass of disks. When comparing the dust continuum emission fluxes and sizes at 1.3 mm, we find that most of the disks in the AGE-PRO sample are consistent with simulations that have either weak or strong dust traps. The evolution of spectral index in the AGE-PRO sample is also suggestive of an unresolved population of dust traps. Future observations at high angular resolution are still needed to test several hypotheses that result from comparing the observations to our simulations, including that more massive disks in gas mass have the potential to form dust traps at larger disk radii.

*Unified Astronomy Thesaurus concepts: Protoplanetary disks (1300); Stellar accretion disks (1579); Hydrodynamical simulations (767)*

## 1. Introduction

Planets are formed from the gas and dust in young circumstellar disks. The diverse properties observed in the exoplanet population, and the architecture of their planetary systems, will depend on the properties and evolution of their natal disks. Gas evolution can lead to (magneto)hydrodynamical instabilities that could create pressure bumps in the disks, halting the dust radial drift and allowing dust particles to efficiently grow from micron-sized particles to kilometer-sized bodies (see J. Bae et al. 2023; J. Drążkowska et al. 2023;

G. Lesur et al. 2023, for recent reviews). An important consideration is that the radial distribution of dust and gas in disks evolves as a function of time. The evolution of the dust particles strongly depends on the distribution and evolution of the gas (T. Birnstiel 2024; N. van der Marel & P. Pinilla 2023). The youngest disks (<1 Myr) have a larger mass reservoir, and also a higher frequency rate of mass replenishment through material infall (M. Kuffmeier et al. 2023). In middle-aged disks (~1–3 Myr), most of their disk material is concentrated in the midplane where dust particles are typically observed in ring-like substructures (e.g., F. Long et al. 2018; S. M. Andrews 2020; M. Villenave et al. 2020). Toward the end of a disk lifetime (3~10 Myr), the gas disk mass reservoir is expected to decrease to undetectable levels, as well as the dust mass, although the amount of dust mass over time strongly depends on the presence

<figure><img src="image_3.png" /></figure>

Original content from this work may be used under the terms of the [Creative Commons Attribution 4.0 licence](https://creativecommons.org/licenses/by/4.0/). Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.
or absence of substructures in disks (e.g., P. Pinilla et al. 2020). Thus, by contrasting gas and dust observations of planet-forming disks at different stages of evolution to gas and dust evolution models, we can identify evolutionary patterns and better understand the potential for planet formation of a disk.

Large observational samples of planet-forming disks have revealed correlations between the properties of disks and their stellar hosts. For example, stars of larger masses typically have disks of larger continuum luminosity ($M_{\star}-L_{\text{mm}}$ relation; M. Ansdell et al. 2016; S. A. Barenfeld et al. 2016; I. Pascucci et al. 2016; L. Testi et al. 2022). Similarly, more-massive stars typically have higher mass accretion rates, which also correlates with the disk (dust) mass (C. F. Manara et al. 2016, 2023). Correlations between the properties of disks have also been found, such as the size-luminosity relation, where the disk size in millimeter continuum emission positively correlates with the millimeter luminosity ($R_{\text{mm}}-L_{\text{mm}}$; A. Tripathi et al. 2017; N. Hendler et al. 2020). When studying these relations as a function of disk age, we can begin testing the mechanisms driving disk and planet-formation evolution, such as viscous evolution models, MHD-disk driven evolution models, and internal and external photoevaporation (e.g., M. Gárate et al. 2021, 2023; A. Somigliana et al. 2022, 2023; B. Tabone et al. 2022; C. Toci et al. 2023).

An important caveat to these observational trends is the calculation of disk masses and sizes from observations at millimeter wavelengths. The gas properties are often inferred from the dust properties. When measuring the dust mass content of a disk, the results rely on assuming optically thin emission, dust opacity, and dust temperature, all of them carrying significant uncertainties. Then, obtaining the gas mass requires an additional assumption on the gas-to-dust ratio in the disk. Therefore, observationally constraining the gas mass and its distribution has remained challenging, which also limits the comparison and interpretation of planet-formation models. This caveat also limits the comparison with dust evolution models. As the evolution of the dust strongly depends on its coupling with the gas, key processes of dust evolution (such as drift, settling, and diffusion) depend on the gas mass of the disk, which is still unknown for the large majority of systems.

The observations of the Atacama Large Millimeter/ submillimeter Array (ALMA) Large Program AGE-PRO have allowed for the measurement of the gas mass and size in a sample of 30 disks three different star-forming regions (SFRs), which are representative of young, middle-aged, and evolved disks (K. Zhang et al. 2025).

These measurements provide the first opportunity to compare global disk evolution models to observations of disk gas masses and sizes, for a large sample of disks. These comparisons are relevant not only to understanding current observational trends but also to exploring which ones we might be missing, and contribute to guiding future efforts toward better characterizing planet-forming disks.

In the present work, we investigate the role of dust traps and viscous evolution in influencing the evolution of the dust reservoir, and we compare them with observational trends across the AGE-PRO sample. In Section 2, we describe the models and the explored parameter space for the dust evolution and radiative transfer models. Section 3 summarizes the results of our models and the comparison with the AGE-PRO observations.

In Sections 4 and 5, we present the discussion and conclusions of this work, respectively.

## 2. Dust Model and Synthetic Observations

We simulate the temporal and radial evolution of the dust size distribution in planet-forming disks for a wide parameter space, along with the relevant observable quantities, such as the flux in the millimeter continuum, the corresponding gas and dust disk size, the dust and gas masses, the dust-to-gas ratio, pebble flux into the inner disk, and the accretion rates. For the gas and dust evolution, we used the 1D dust evolution code Dustpy<fnref n="18" /> (S. M. Stammler & T. Birnstiel 2022), which is based on the disk model of T. Birnstiel et al. (2010), and includes the effects of coagulation and fragmentation of multiple dust species (by dust species, we refer to dust grain of different sizes), their radial advection and diffusion, as well as the viscous evolution of the gas.

Then, we used the radiative transfer code RADMC-3D<fnref n="19" /> (C. P. Dullemond et al. 2012) to post-process the simulations at selected snapshots, obtaining the corresponding emission in the millimeter continuum. We compare the properties of these radiative transfer images directly to those measured by AGE-PRO observations, as the AGE-PRO disk sizes were obtained using visibility modeling (M. Vioque et al. 2025), which should alleviate the effect of beam convolution. In the Appendix, we also describe the comparison between the true properties of the radiative transfer images, and those that are measured from images convolved with a Gaussian beam and with thermal noise, to simulate the conditions of the CLEAN images. In this section, we provide an overview of our workflow, and describe the relevant parameters of each code. For a detailed description of Dustpy and RADMC-3D, we refer to S. M. Stammler & T. Birnstiel (2022) and C. P. Dullemond et al. (2012).

### 2.1. Dust Evolution

The advection, diffusion, and growth or fragmentation of dust particles are strongly determined by their Stokes number, St, which can be understood as the dimensionless coupling time to the gas motion or, more intuitively, as the dynamical grain size. The Stokes number is related to the physical grain size $a$ through the following equation:

$$ \text{St} = \frac{\pi a \rho_s}{2 \Sigma_g} \cdot \begin{cases} \frac{1}{4} & \lambda_{\text{mfp}}/a \ge 4/9 \\ \frac{a}{9} & \lambda_{\text{mfp}}/a < 4/9, \end{cases} \quad (1) $$

where $\rho_s$ is the grain material density, $\Sigma_g$ is the gas surface density, and $\lambda_{\text{mfp}}$ is the gas mean free path, where the latter is used to distinguish between the Epstein and Stokes I drag regimes (T. Birnstiel et al. 2010).

Then, the radial velocity of dust particles with a given size particle is defined through:

$$ v_d = \frac{1}{1 + \text{St}^2} v_\nu - \frac{2\text{St}}{1 + \text{St}^2} \eta v_k, \quad (2) $$

<footnote n="18">¹⁸ DustPy is available on github.com/stammler/DustPy. Version 0.5.6 was used for this work.</footnote>

<footnote n="19">¹⁹ www.ita.uni-heidelberg.de/~dullemond/software/radmc-3d/</footnote>
where the first term represents the coupling to the gas radial viscous velocity (J. E. Pringle 1981), and:

$$v_{\nu} = - \frac{3}{\Sigma_g \sqrt{r}} \frac{\partial}{\partial r} (\nu \Sigma_g \sqrt{r}), \quad (3)$$

with $\nu = \alpha c_s h_g$ being the kinematic viscosity (N. I. Shakura & R. A. Sunyaev 1973), $c_s$ is the sound speed, $h_g$ is the pressure scale height, $r$ is the radial coordinate, and $\alpha$ is the dimensionless turbulence parameter. The second term corresponds to the radial drift induced by the exchange of angular momentum between gas and dust (Y. Nakagawa et al. 1986; T. Takeuchi & D. N. C. Lin 2002), which depends on the difference between the gas azimuthal velocity and the local Keplerian speed $v_K$, and can be expressed in terms of the pressure gradient through $\eta = -(1/2) (h_g/r)^2 d \ln P/d \ln r$, where $P$ is the pressure.

In principle, the drift term leads to pebble-sized dust grains moving inward in a smooth disk, which will lead to very small dust disk sizes over time, very low millimeter fluxes, and high spectral indices (P. Pinilla et al. 2012, 2020). In contrast, when disks have pressure bumps, particles can be trapped in pressure maxima, which we commonly refer to as dust traps (S. J. Weidenschilling 1977; P. Pinilla et al. 2012), while dust diffusivity counterbalances the trapping efficiency, preventing the rings from becoming infinitely narrow.

In our models, the dust grains of different sizes interact between each other through coagulation and fragmentation following the Smoluchowski equation (see T. Birnstiel et al. 2010; S. M. Stammler & T. Birnstiel 2022, for details), which leads to two distinct regimes in which the local maximum grain size may be either limited by either by fragmentation or drift (see review by T. Birnstiel et al. 2016). We do not include the effect of back-reaction of dust onto the gas, which can affect disk evolution when there are high concentrations of dust in the disk (e.g., M. Gárate et al. 2020).

## 2.2. Initial Conditions and Substructure Model

We initialized the disk surface density profile using a variation of the D. Lynden-Bell & J. E. Pringle (1974) self-similar solution,

$$\Sigma_g(r) = \frac{M_{\text{disk}}}{2\pi r_c^2} \left(\frac{r}{r_c}\right)^{-1} \exp(-r/r_c) \frac{\alpha_0}{\alpha(r)}, \quad (4)$$

where $M_{\text{disk}}$ and $r_c$ correspond to the initial disk mass and characteristic radius. We included bump-like substructures $\alpha(r)$ in our turbulence profile that lead to self-sustained gaps in the gas surface density profile (as in J. Stadler et al. 2022). The radial turbulence profile is described by the following expression:

$$\alpha(r) = \alpha_0 \times \left( 1 + \sum_i A_{\text{gap}} \exp \left( - \frac{(r - r_{\text{gap},i})^2}{2w_{\text{gap}}^2} \right) \right), \quad (5)$$

which accounts for multiple Gaussian substructures of amplitude $A_{\text{gap}}$ and location $r_{\text{gap},i}$, imposed over a base turbulence parameter of $\alpha_0$. The width of the traps $w_{\text{gap}}$ is defined by setting the Gaussian FWHM to the local scale height value at the gap location $h_g(r_{\text{gap}})$. As the disk evolves, the viscous evolution will increase the disk characteristic

radius, but keep the gap-like substructures that act as dust traps. The work by P. Pinilla et al. (2021b) compared the effect of including static pressure bumps from the beginning of the simulations versus pressure bumps as assumed in our models. As pressure bumps introduced using Equation (5) can take more than 1 Myr to reach the final amplitude, they trap less dust than static traps that introduce early in the simulations. This affects the final observables, such as millimeter fluxes and spectral indices, which are discussed in Section 4.2.

This approach has the limitation that it does not capture the 3D structure of the dust dynamics, which, for example, can result in more efficient diffusion across the gap when the azimuthal dimension is considered (P. Weber et al. 2018). Similarly, our approach does not include mechanisms such as photoevaporation or MHD winds, which can contribute to the disk evolution (e.g., C. J. Clarke et al. 2001; D. Fedele et al. 2010; R. Alexander et al. 2014; Á. Ribas et al. 2015; B. Tabone et al. 2022; I. Pascucci et al. 2023).

We assumed that the disk temperature is set by stellar irradiation and neglected the contribution of other alternative sources, such as viscous heating or external irradiation. Therefore, the temperature profile is given by the following power-law profile:

$$T(r) = \left( \frac{L_*}{L_{\odot}} \right)^{1/4} \left( \frac{r}{1 \text{ au}} \right)^{-1/2} 263 \text{ K}, \quad (6)$$

which depends on the stellar luminosity $L_*$.

Finally, the initial dust size distribution follows the interstellar medium distribution (J. S. Mathis et al. 1977), with a radially uniform dust-to-gas ratio of $\epsilon_0$, which then evolves through the simulations. The initial dust distribution goes from 0.5–1 µm, and the material density of the dust grains is $\rho_s = 1.6 \text{ g cm}^{-3}$, with a fragmentation velocity of $v_{\text{frag}} = 10 \text{ m s}^{-1}$ (B. Gundlach et al. 2011; K. Wada et al. 2011; B. Gundlach & J. Blum 2015), though laboratory studies suggest that this value might be lower (B. Gundlach et al. 2018; G. Musiolik & G. Wurm 2019; T. Steinpilz et al. 2019). Dust grains could still grow to larger sizes even at lower fragmentation velocities, as long as the turbulence near the midplane is also reduced (P. Pinilla et al. 2021b).

## 2.3. Parameter Space

Our simulations explore the parameter space in the stellar mass $M_*$, the disk initial mass-to-star ratio ($M_{\text{disk}}/M_*$), the disk characteristic size $r_c$, and the different gap amplitudes $A_{\text{gap}}$. The ranges covered by these parameters are wide enough to represent the observed properties of our systems. For simplicity, we also refer to the trap amplitudes of Equation (5) as “None” ($A_{\text{gap}} = 0$), “Weak” ($A_{\text{gap}} = 1$), or “Strong” ($A_{\text{gap}} = 4$). We also explore the parameter space in the $\alpha_0$ turbulence, and initial dust-to-gas ratio $\epsilon_0$. The values explored are shown in Table 1. It is important to note that these simulations were performed during the time that the AGE-PRO program was receiving data, and for this reason, the initial gas disk mass ($M_{\text{gas}}$) was not tailored to the particular values that are inferred from these observations (as presented in L. Trapman et al. 2025b). Simulations that are tailored to the properties of each disks will be part of a future publication.

The stellar masses range from 0.25 $M_{\odot}$ to 1.0 $M_{\odot}$ in steps of 0.25, which covers the stellar masses for the AGE-PRO sources. For each stellar mass, we assign the corresponding
<table><caption>**Table 1**
Disk Model Parameters</caption><thead><tr><th>Symbol</th><th>Description [Units]</th><th>Value</th></tr></thead><tbody><tr><td>M<sub>*</sub></td><td>Stellar mass [M<sub>☉</sub>]</td><td>0.25, 0.5, 0.75, 1.0</td></tr><tr><td>L<sub>*</sub></td><td>Stellar Luminosity [L<sub>☉</sub>]</td><td>0.15, 0.25, 0.5, 1.0</td></tr><tr><td>M<sub>disk</sub></td><td>Initial disk mass [M<sub>*</sub>]</td><td>0.01, 0.05, 0.1</td></tr><tr><td>r<sub>c</sub></td><td>Characteristic radius [au]</td><td>15, 30, 60, 120</td></tr><tr><td>A<sub>gap</sub></td><td>Trap amplitude</td><td>0, 1, 4</td></tr><tr><td>r<sub>gap</sub></td><td>Gap locations [au]</td><td>10, 40, 70</td></tr><tr><td>α<sub>0</sub></td><td>Viscosity parameter</td><td>10<sup>-4</sup>, 10<sup>-3</sup></td></tr><tr><td>ε<sub>0</sub></td><td>Initial dust-to-gas ratio</td><td>0.01, 0.05</td></tr></tbody></table>

stellar luminosity (matched in the order shown in Table 1), which is used to determine the disk temperature in our dust evolution model. In addition, for each combination of disk mass, we initiate simulations with four different $r_c$, which therefore includes simulations of low-mass stars with very compact disks ($M_* = 0.25 M_☉$ and $r_c = 15$ au), as well as solar-type stars with very extended disks ($M_* = 1.0 M_☉$ and $r_c = 120$ au). All of the other possible combinations, such as low-mass stars with extended disks, and solar-type stars with compact disks, are also considered.

The locations of the gaps (when $A_{gap} > 0$) are fixed at 10, 40, and 70 au, provided that $r_{gap} \le 2r_c$. That means that the disks with $r_c = 15$ au only have one gap at $r_{gap} = 10$ au, the disks with $r_c = 30$ au have two gaps, while the disks with $r_c = 60, 120$ au have all three gaps at the same time. We include gaps only inside $2r_c$, since 83% of the total disk mass is contained within this radius in a disk that follows the same prescription as in D. Lynden-Bell & J. E. Pringle (1974), and therefore traps farther outside would only affect a minor fraction of the inward-drifting solids.

The turbulence value associated with the dust dynamics (settling, diffusion, and fragmentation) is typically the same as the base value of the gas turbulence $\alpha_0$. However, for the particular case with low gas turbulence $\alpha_0 = 10^{-4}$, we chose to decouple the turbulence value associated with the dust fragmentation, in order to avoid an overly efficient grain growth regime, which is computationally expensive to trace. For this reason, we impose that the fragmentation turbulence must be at least $\alpha_{frag} = 5 \times 10^{-4}$, while turbulence values associated with radial diffusivity and settling are always set to the value of the gas $\delta_{rad} = \delta_{vert} = \alpha_0$.

## 2.4. Numerical Grid

The standard radial grid of the simulations goes from $r_{in} = 5$ au to $r_{out} = 500$ au with $n_r = 250$ logarithmically spaced grid cells. This grid is used for most of our simulations, except those with a small characteristic radius $r_c = 15$ au, where the grid goes from 2.5 au to $r_{out} = 250$ au, to prevent losing a significant fraction of the mass in the inner regions.

The standard mass grid for the grain size distribution ranges from $m_{min} = 10^{-12}$ g to $m_{max} = 3 \times 10^6$ g, which, in terms of the grain sizes, corresponds to approximately to $a_{min} \approx 0.5$ µm and $a_{max} \approx 75$ cm, with $n_m = 134$ logarithmically spaced grid cells. For the simulations with low turbulence $\alpha_0 = 10^{-4}$, the upper limit of the mass grid is extended to $m_{max} = 10^8$ g (and $n_m = 141$ grid cells).

The simulations ran from time $t = 0$ and up to $t = 10$ Myr. The key output of the DustPy simulations is time evolution of the dust surface density, for every dust species in the grain size distribution. We show an example of the dust size

<figure><img src="image_4.png" /><figcaption>**Figure 1.** Dust size distribution for a disk with (top) and without (bottom) substructures at 1 Myr. The simulations correspond to a disk with $M_* = 1 M_☉$, $M_{disk} = 0.1 M_*$, $R_c = 60$ au, $\alpha_0 = 10^{-3}$, and $\epsilon_0 = 0.01$.</figcaption></figure>

distribution in a disk with gaps and without gaps in Figure 1, where dust particles are trapped at the outer edges of the gaps.

Considering the different combinations of $\alpha_0$ and $\epsilon_0$ values, we have two groups of simulations with $\alpha_0 = 10^{-3}$ and different $\epsilon_0$, and one with $\alpha_0 = 10^{-4}$ and $\epsilon_0 = 0.01$. Each of these three groups has 144 disks, which result from choosing the different $M_*$, $M_{disk}$, $R_c$, and $A_{gap}$. The state of each disk is saved at nine different ages (0.1, 0.5, 0.8, 1, 2, 3, 5, 7, and 10 Myr), thus resulting in 3888 different snapshots of the disks evolution.

## 2.5. Radiative Transfer

The millimeter flux of the dust continuum emission is obtained from the simulated dust size distribution with the radiative transfer code RADMC-3D (C. P. Dullemond et al. 2012). The main goal of obtaining these fluxes is to compare with the observations of AGE-PRO, and thus we chose to generate images at wavelengths of $\lambda = 1.05$ mm and $\lambda = 1.3$ mm, matching the two spectral setups of AGE-PRO. Given the large number of simulations, we select only five ages per disk to run the radiative transfer, at 0.5, 1, 2, 5, and 10 Myr, covering all of the relevant ages of the SFRs in the AGE-PRO sample (Ophiuchus, Lupus, and Upper Sco; K. Zhang et al. 2025). This implies a total of 19,440 radiative transfer models ($5 \times 3888$).

To obtain the dust opacities, we use the code OpTool$^{20}$ (C. Dominik et al. 2021), and assume that the grain composition follows the model of L. Ricci et al. (2010),

<footnote>$^{20}$github.com/cdominik/optool</footnote>
where the dust grains are composed of 10% Silicate (B. T. Draine 2003), 20% Carbon (V. G. Zubko et al. 1996), 30% water ice (S. G. Warren & R. E. Brandt 2008), and 40% vacuum.

Rather than calculating the opacities for all 134 (or 141) dust grain species, we chose to re-sample our grain size grid (and the corresponding dust surface densities) down to only 30 dust grain species, which still cover all of the key features of the opacity profile from L. Ricci et al. (2010) at our target wavelength. The size grid for the radiative transfer model is defined with three intervals of logarithmically spaced grid cells between 0.5 µm, 50 µm, 0.1 cm, and 20 cm, with 5, 15, and 10 grid cells, respectively. Dust sizes that fall outside the grid boundaries are ignored in the re-sampling, with no significant loss of mass.

The 1D dust surface density profiles are converted into 2D volumetric density profiles, in a spherical grid where the radial coordinates are the same as the Dustpy radial grid, and the colatitude coordinate goes from $\theta_{min} = 0$ to $\theta_{max} = \pi$, with $n_{\theta} = 180$ grid cells. The disk is assumed to be axisymmetric, and we also assume that the dust is in vertical hydrostatic equilibrium distributed according to the following Gaussian profile (S. Fromang & R. P. Nelson 2009):

$$ \rho_{d,i}(r, z) = \frac{\Sigma_{d,i}(r)}{\sqrt{2\pi} h_{d,i}(r)} \exp\left(-\frac{z^2}{2h_{d,i}^2(r)}\right), \quad (7) $$

where $\rho_{d,i}$ is the volumetric dust density for the dust grain species i, $\Sigma_{d,i}$ is the dust surface density, $h_{d,i}$ is the dust scale height from B. Dubrulle et al. (1995), and z is the vertical distance from the midplane. We illustrate the spatial distribution of small and large grains for the case of a disk with substructures in Figure 2, where we observe how larger grains have migrated inward and are settled closer to the midplane. The central star is assumed to emit as a perfect blackbody, where the total luminosity is the corresponding value of $L_{*}$, and the stellar temperature is obtained from the stellar tracks (A. Dotter et al. 2008).

The radiative transfer calculation was performed in two steps. First, we ran a Monte Carlo simulation to obtain the dust temperature for every grain species, using $n_{phot} = 10^7$ photon packages. In a second step, we obtained the simulated images at the target wavelengths, with $n_{phot} = 2 \cdot 10^5$ photon packages to account for the contribution of isotropic scattering. We ignore more complex scattering processes due to the high computational cost for covering our complete parameter space, and we verified the flux variations were well below 10% of the total flux.

For additional testing, we also further processed the radiative transfer images into synthetic observations to simulate the sensitivity and angular resolution of the AGE-PRO observational setup. In the main text, all of the comparisons are done between the radiative transfer models and the properties of the AGE-PRO disks as recovered from the observations, and additional discussion on the influence of the beam convolution when comparing disk radii is provided in the Appendix.

### 3. Results

In the following Section, the comparison between the AGE-PRO sources and our simulations is done by groups of different ages. The Ophiuchus sources are compared to our simulations at

<figure><img src="image_5.png" /><figcaption>Figure 2. Dust spatial distribution of small (top, 1 µm) and large (bottom 1mm) grains, for a disk with substructures at 1 Myr from the same model presented in the top panel of Figure 1.</figcaption></figure>

0.5 Myr. The Lupus sources are compared to the simulations at 1 Myr and 2 Myr. The Upper Sco sources are divided in two groups: For stellar hosts younger than 3 Myr, we compared them to our simulations at 2 Myr, while for the stellar hosts older than 3 Myr, we compared them to the simulations at 5 Myr. The decision to split at these ages comes from the recalculation of stellar properties for the sources within AGE-PRO, which resulted in two age clusters in the Upper Sco sample, one of which is more consistent with the ages of the Lupus disks, while another is consistent with an older age. The full discussion about the estimation of $M_{*}$ and age is given in K. Zhang et al. (2025), and we refer the reader to that work for a detailed description. We compared the observational results of AGE-PRO to the simulations with standard values of $\alpha = 10^{-3}$ and dust-to-gas ratio of 0.01. While the specific properties of the disks can change slightly for $\alpha = 10^{-4}$ or initial dust-to-gas ratio of 0.05, the qualitative evolution and spread of properties is the same as with the standard values.

The following subsections focus on comparing the results of the simulations to those of the observed disks in AGE-PRO. For each comparison, we present the results of the simulations in groups of disks with the same initial condition, such as their initial disk mass as a fraction of stellar mass ($M_{disk}/M_{*}$), or by their initial disk characteristic size $r_c$.

#### 3.1. Gas Mass Evolution

In our simulations, the evolution of $M_{gas}$ is only driven by viscous evolution. The presence of dust traps does not produce
<figure><img src="image_6.png" /><figcaption>**Figure 3.** Evolution of disk gas mass ($M_{gas}$) as a function of stellar mass ($M_{*}$) for different times of evolution, as shown at the top of each panel. The shaded regions enclose the disks with the same initial disk mass relative to their host star, and the markers show the maximum, median, and minimum value. The spread in $M_{gas}$ for a fixed stellar mass comes from the different initial disk characteristic size $R_{c}$ within the same group of initial disk mass, and its time evolution. The plot features only the simulations without traps ($A_{gap} = 0$), standard turbulence $\alpha_{0} = 10^{-3}$, and interstellar medium dust-to-gas ratio of $\epsilon_{0} = 0.01$. The disks with traps, and those from the alternative simulation setups show a comparable behavior. For comparison, the gray square markers show the estimated values for $M_{gas}$ for the disks in Oph (for the 0.5 Myr snapshot), Lupus (for the 1 and 2 Myr snapshots), and Upper Sco (for the 2 and 5 Myr snapshots) as obtained in L. Trapman et al. (2025b).</figcaption></figure>

a significant difference in the evolution of $M_{gas}$, and thus, we only compare the disks masses measured from AGE-PRO to those of our disks without traps ($A_{gap} = 0$), as shown in Figure 3. When comparing against $M_{gas}$ derived from L. Trapman et al. (2025b), we find that the latter have a larger dispersion in their values compared to those found in our simulations. In fact, the measured $M_{gas}$ from observations spans ~3 orders of magnitude in Ophiuchus, while our initial disk conditions ($t = 0$ Myr) only span 1 order of magnitude for each stellar mass ($M_{gas}$ from 0.01 $M_{*}$ to 0.1 $M_{*}$). The lowest disk masses in Ophiuchus have 2 orders of magnitude less gas mass than those of comparable $M_{*}$ in our simulations, and only about half of the targets in Lupus and Upper Sco are in agreement to the measured $M_{gas}$ of the models. This mismatch between the trends of $M_{gas}$ evolution prevents further comparisons between $M_{gas}$ and other disk or stellar parameters, such as the trends analyzed in L. Trapman et al. (2025b). This topic is further discussed in detail in Section 4, as well as in B. Tabone et al. (2025) and R. Anania et al. (2025).

## 3.2. Dust Mass Comparison

Despite the difference in distribution between modeled and observed gas masses, we find consistent values for the 1.3 mm dust continuum brightness among the AGE-PRO targets and our models, which were estimated from the radiative transfer models. In observations, the dust mass content in pebbles can be roughly approximated from the disk dust continuum emission using (R. H. Hildebrand 1983):

$$M_{\text{dust}} = \frac{d^2 F_{\nu}}{\kappa_{\nu} B_{\nu}(T)}, \qquad (8)$$

where it is assumed that the emission is optically thin, where d is the distance to each star, measured as the inverse of the parallax from Gaia DR3 (Gaia Collaboration et al. 2021), $\nu$ is the frequency of the AGE-PRO Band 6 observation (1.3 mm, 230 GHz), $B_{\nu}$ is the Planck function, and $\kappa_{\nu}$ is the mass absorption coefficient, which is typically assumed to be $\kappa_{\nu} = 2.3 \cdot (\nu/230 \text{ GHz})^{0.4} \text{ cm}^2 \text{ g}^{-1}$, as in S. M. Andrews et al. (2013). For the temperature, a constant $T = 20$ K is used for consistency with the literature (e.g., M. Ansdell et al. 2016). The dust content of the AGE-PRO disks was calculated directly from the measured millimeter fluxes, and they are

presented in K. Zhang et al. (2025), D. A. Ruiz-Rodriguez et al. (2025), D. Deng et al. (2025), and C. Agurto-Gangas et al. (2025). The pebble mass in our simulations ($M_{pebble}$) is measured by adding the mass content in particles of sizes ranging from 0.1 mm to 10 cm, which are the grains that should contribute more efficiently to the emission in the millimeter continuum in Bands 6 and 7 (L. Ricci et al. 2010). These results are presented in Figure 4.

The simulations with weak or strong traps ($A_{gap} = 1-4$) tend to match the observed dust masses closely over time, although there is a wide spread in $M_{pebble}$ from the simulations, given by the different values of initial disk mass, disk size, and amplitude of dust trapping. For the older sources in Upper Sco, the disks are more consistent with the simulations with traps and larger initial disk mass, of at least 0.05 $M_{*}$ (see panels showing the 5 Myr outputs). Overall, disk with strong dust traps and high initial disk masses are able to retain their dust mass over their lifetime evolution, remaining massive in dust even at the 10 Myr of age.

In contrast to the strong traps, $M_{pebble}$ in the simulations tends to drop quickly if no traps are present, reaching values that are several orders of magnitude below $0.1 M_{\oplus}$ even by 2 Myr for the lowest initial disk masses. This trend also holds at higher initial dust-to-gas ratios of $\epsilon_{0} = 0.05$, and at a lower alpha value of $\alpha = 10^{-4}$. Thus, the absence of traps results in the loss of dust mass for a wide range of initial conditions and disk viscosity, and these models are unable to explain observations.

An interesting trend from the models is that the $M_{pebble}$ decrease over time is more significant for lower $M_{*}$ in the absence of traps, as previously demonstrated by P. Pinilla et al. (2013) using a smaller parameter space. However, this trend is not seen in the AGE-PRO observations, potentially due to the presence of unresolved dust traps and the limited sample size.

## 3.3. Comparison with Observations: Dust Continuum Emission

The assumption of optically thin emission and constant dust temperature might oversimplify the conditions of the dust continuum emission in the AGE-PRO sample, especially considering that unresolved dust substructures could become more optically thick, the temperature might not be constant,
<figure><img src="image_7.png" /><figcaption>**Figure 4.** Pebble mass vs. stellar mass at different ages and trap amplitudes (strong, weak, and none from top to bottom, respectively) from simulations, considering pebbles as grain sizes between 0.1 mm and 10 cm. The shaded regions enclose the simulated disks of different initial disk masses. For comparison, the star markers show the estimated values for $M_{\text{dust}}$ for the disks in Ophiuchus (compared to 0.5 Myr), Lupus (compared to 1 and 2 Myr), and Upper Sco (compared to 3 or 5 Myr, depending on the age of the sources). The dashed horizontal line shows the approximated $3\sigma$ sensitivity of AGE-PRO.</figcaption></figure>

and the opacity law might be more complex than our simple assumption. Thus, comparing the observed flux of the AGE-PRO sources to the radiative transfer model of the simulations is a more robust approach than comparing $M_{\text{pebble}}$ to the estimated $M_{\text{dust}}$.

In Figure 5, the total flux of our simulations at 1.3 mm is compared to the measurements of AGE-PRO at the same wavelength. The overall conclusion from comparing dust masses remains the same when comparing the fluxes: disks with strong and weak traps are more representative of the evolutionary trend observed across the three regions. However, we also find that the total brightness of some disks could be explained by disks without any dust trapping. The key difference with the dust mass comparison from Section 3.2 comes from the contribution to the total brightness from small grains ($<0.1$ mm) in the opacities from L. Ricci et al. (2010). See Section 4 for a more detailed discussion on this topic.

The 1.3 mm flux of the three faintest sources in Lupus sample (Lupus 5, 7, 9) is barely consistent with the strong trapping dust models at 2 Myr of age, and they are not consistent with the 1 Myr strong trap disks. On the opposite side, the brightest sources in Lupus (Lupus 2, 10) cannot be explained in a model without dust traps at either 1 or 2 Myr of age. Similarly, the brightest sources of the Upper Sco sample (Upper Sco 1, 8, 9, 10) cannot be explained by models without dust traps if they are older than 5 Myr.

As in the case of the comparison with the dust mass in Section 3.2, we find that the fluxes observed in the AGE-PRO sample are more consistent with the simulations with strong and weak dust traps for ages later than 1 Myr. Differently from the comparison with pebble mass, where the dust trapping amplitude could be distinguished even at the early stages at 0.5 Myr, their dust continuum brightness is consistent with any of the dust trapping scenarios. Therefore, according to our models, disentangling the potential trapping efficiency in planet-forming disks purely from continuum flux measurements is only possible when looking at the evolution over time, as young disks ($<1$ Myr) have similar fluxes in the models independent of the presence and amplitude of pressure bumps.

### 3.4. Disk Size in Pebble Mass and Dust Continuum Emission

We estimate the size of the disks in dust content by measuring the radius at which 90% of the $M_{\text{pebble}}$ is contained, and we refer to this value as $R_{90\%_{\text{pebble}}}$, as shown in Figure 6. Because it is derived from $M_{\text{pebble}}$, this $R_{90\%_{\text{pebble}}}$ only considers pebbles larger than 0.1 mm. In the disks with strong and weak traps, the location of the radius containing the $R_{90\%_{\text{pebble}}}$ is completely determined by the location of the outer dust trap, even at early ages of 1 Myr, and also for the late ages of 10 Myr. In the case of disks with no dust trap, the $R_{90\%_{\text{pebble}}}$ decreases as a function of age as the larger pebbles drift unimpeded toward the center of the
<figure><img src="image_8.png" /><figcaption>Figure 5. Same as Figure 4, but with dust continuum fluxes at 1.3 mm (230 GHz) vs. the stellar mass.</figcaption></figure>

disk. We further discuss this result in comparison with the literature in Section 4.

From an observational perspective, calculating the R$_{90\%_{,pebble}}$ is challenging due to the poorly constrained distribution of optical depth, opacity, and temperature. Thus, an alternative measure for the disk size is the radius containing a certain fraction of the brightness, typically chosen to be 90% (e.g., G. P. Rosotti et al. 2019; N. Hendler et al. 2020; N. T. Kurtovic et al. 2021; E. Sanchis et al. 2021; F. Long et al. 2022). Further details about these measurements in the AGE-PRO sample are described in M. Vioque et al. (2025) and L. Trapman et al. (2025a). We compared these radii containing the 90% of the emission (R$_{90\%_{,mm}}$) from the observations of AGE-PRO to that of our radiative transfer images at 1.3 mm, as shown in Figure 7 as a function of M$_{\star}$ and age.

At the early age of 0.5 Myr, the distribution of R$_{90\%_{,mm}}$ is similar regardless of the strength of trapping, with the R$_{90\%_{,mm}}$ only slightly evolving from the initial disk conditions. In the case of the disks with traps, at an age of 1 Myr, the material in the dust traps represents the major contribution to the disk continuum brightness. Thus, R$_{90\%_{,mm}}$ also traces the location of the outermost trap, both for weak and strong traps.

The role of dust traps in determining R$_{90\%_{,mm}}$ changes as the radial drift starts depleting the disk of pebbles. At 2 Myr for the disks with weak traps, and at 5 Myr for the disks with strong traps, the most compact disks (with a starting critical radius of 15 au) show an increase in their R$_{90\%_{,mm}}$, as the relative contribution to the disk brightness reverses from the dust in the traps to the micron-sized grains in the outer regions

of the disk. By the age of 10 Myr, disks with very different initial conditions show similar R$_{90\%_{,mm}}$, even for groups of disks with the same trap strength. However, these disks can still be distinguished based on their millimeter brightness.

In the case of disks with no dust traps, the R$_{90\%_{,mm}}$ tends to increase as a function of age, as the micron-sized grains in the outer disk become the dominant contributing source to the millimeter emission when the larger grains are depleted by migrating inward in the disk. This is the opposite as the evolution of R$_{90\%_{,pebble}}$. Similarly to the simulations with dust traps, the disks with smaller initial r$_{c}$ show a faster increase as a function of age, as solids are depleted at earlier stages. The disks in Upper Sco that are older than 3 Myr can only be reproduced by the simulations where the outermost dust trap was located closer to the star (<20 au), and they are more consistent with strong trap models.

### 3.5. Spectral Index

The dual-band setup of the AGE-PRO observations (Band 6 and Band 7 at 1.3–1.05 mm) allows for comparison of the spatially integrated spectral index α$_{mm}$ with our radiative transfer images. We measure the spectral index from the total flux of the disks, assuming Rayleigh–Jeans limit and using the relation:

$$ \alpha_{mm} = \frac{\log_{10} F_{B7} - \log_{10} F_{B6}}{\log_{10} \nu_{B7} - \log_{10} \nu_{B6}}, \quad (9) $$

where F is the flux in each ALMA band, and ν is the representative frequency for that flux. Values of α$_{mm}$ ≈ 2 are
<figure><img src="image_9.png" /><figcaption>**Figure 6.** Radius enclosing 90% of the pebble mass vs. stellar mass. The shaded regions enclose the simulated disks of different initial disk critical radius, as shown by their different symbols. The dashed lines show the location of the dust traps in the simulations.</figcaption></figure>

expected for disks that are optically thick at the observed wavelengths, while higher values are typically associated with more optically thin emission (see M. Tazzari et al. 2021). The results are shown in Figure 8.

Interestingly, a few disks in the Lupus region have $\alpha_{mm} \approx 2$ or below, which are not matched by any of our simulations. These very low values are common in the Lupus SFR (M. Tazzari et al. 2021), and a more detailed description of $\alpha_{mm}$ with a wider wavelength coverage might be needed to understand the mechanisms driving these low values. Among the possible scenarios, free-free emission coming from the disk center could contribute to lower $\alpha$ when the total flux is considered (e.g., A. Garufi et al. 2025), especially in compact disks, and also self-scattering of millimeter emission in optically thick disks (H. B. Liu 2019). Thus, multiwavelength observations at higher angular resolution are needed to recover $\alpha_{mm}$ as a function of radii (e.g., J. Huang et al. 2018; F. Long et al. 2020; P. Pinilla et al. 2021a; K. Doi & A. Kataoka 2023; S. Das et al. 2024; A. Sierra et al. 2024). The remaining disks in Lupus seem consistent with either weak or strong traps, as their spectral index cannot be reproduced by disks with no traps.

A large dispersion in the values of spectral indices is also observed in Upper Sco, ranging from 1.9–3.2, which is probably tracing the diversity of locations and strength of dust traps in these disks. Almost all of the values for spectral index in this region are consistent with strong and weak dust traps, but a few disks with the highest spectral index could be consistent with weak traps. Even though some Upper Sco disks

have a matching spectral index with simulations with no traps, those disks are not consistent in the other properties (flux or radii) to those simulations, which contributes to rule-out that scenario. Our results are also consistent with those from J. Stadler et al. (2022) and L. Delussu et al. (2024), who also demonstrated that dust traps are needed to explain the spectral index of large disk populations.

The observed spectral indices of Lupus and Upper Sco are also shown in Figure 9 in comparison with multiple properties of those systems. These comparisons show no clear trend with disk gas or dust mass ($M_{gas}$, $M_{dust}$), nor with disk size ($R_{gas}$, $R_{dust}$). A tentative trend, however, is seen between the spectral index and the gas-to-dust mass ratio of the disks (panel (g) in Figure 9). In the Upper Sco SFR, there is a potential relation between larger $R_{gas}$ and higher spectral index, but a larger disk sample is needed to test it. Similarly, higher angular resolution observations would allow us to explore the relation between spectral index and dust emission morphology, including the location and shape of dust traps that could be driving the large range of total flux measurements.

# 4. Discussion

## 4.1. The Evolution of $M_{gas}$

The $M_{gas}$ of the AGE-PRO sample from L. Trapman et al. (2025b) are distributed over a wider range of masses compared to our models, which only include the pure viscous evolution scenario, as presented in Figure 3. The dispersion of 1 order of magnitude in the initial mass of our disks, which were set to be
<figure><img src="image_10.png" /><figcaption>**Figure 7.** Radius enclosing 90% of the millimeter flux continuum as a function of the stellar mass. The symbols and colors regions are the same as in Figure 6. The star markers show the measured R<sub>90%</sub> of the AGE-PRO sources, with the exception of V1094 Sco (Lupus 10). The dashed lines show the location of the dust traps in the simulations.</figcaption></figure>

initially 1%, 5%, and 10% of the stellar host mass $M_\star$, is not enough to cover the 3 orders of magnitude of difference that is observed in Ophiuchus (see Figure 3), which later increases to be ~4 orders of magnitude for Lupus and Upper Sco. Within the uncertainties, less than half of the disks are consistent with those of our simulations, and the disks with the lowest $M_{\text{gas}}$ do not even coincide with the simulations with initial $M_{\text{gas}} = 0.01 M_\star$.

The slow decline of $M_{\text{gas}}$ in our simulations is a result of the long viscous timescales in our disks, which are simulated with a constant $\alpha$ of $10^{-3}$ or $10^{-4}$ depending on the setup. While a value greater than $10^{-3}$ for $\alpha$ would accelerate the disk dispersal by viscous accretion, various indirect observational results suggest that, at a disk population level, the disk turbulence typically has low values of $\alpha \lesssim 3 \cdot 10^{-3}$ (see review by G. P. Rosotti 2023). This is suggested by the absence of evidence for viscous spreading (e.g., F. Long et al. 2022), the large spread in the $\dot{M}/M_\star$ relation (e.g., G. Lodato et al. 2017), or the low scale height at millimeter wavelengths measured from edge-on disks (M. Villenave et al. 2020, 2022).

Another caveat when considering viscous timescales in our simulations is that our models have a constant $\alpha$ not only in radii but also in time. A time evolution in the value of $\alpha$ is expected as the disks and their environment evolve, thus changing the distribution of temperature, density, and ionization (e.g., T. N. Delage et al. 2022, 2023). Very short viscous timescales for ages <1.0 Myr could contribute to removing disk masses, with $\alpha$ values decreasing for later disk stages,

explaining the small variation in disk masses between Lupus and Upper Sco, the low $M_{\text{gas}}$ for some Ophiuchus disks, and also the large variation between Ophiuchus and Lupus.

While the star formation process can potentially result in a large diversity of initial conditions, population synthesis models such as M. R. Bate (2018) have shown that the initial disk mass is typically in the range between $M_{\text{disk}}/M_\star = 0.1$ to 1 for ages <$10^5$ yr. Since it seems unrealistic to expect the initial disk-star-mass ratio to be <1% (especially considering that the sample selection of AGE-PRO required a previous $^{12}$CO detection), the mismatch between our simulations and the AGE-PRO observations suggests that viscous evolution with globally low-viscosity scenarios ($\alpha = 10^{-3} \sim 10^{-4}$) is not solely responsible for transporting angular momentum and/or removing gas from the disk over its lifetime. Instead, higher-viscosity scenarios, or other processes in addition to viscous evolution, must also be contributing, such as MHD winds or photoevaporation. These mechanisms must begin removing mass from the disks at a very early stage (<0.5 Myr), as even for a region such as Ophiuchus, the masses of our viscous evolution models are larger than some of those observed. This hypothesis is also supported by B. Tabone et al. (2025) and R. Anania et al. (2025), where they explore the role of MHD and (internal/external) photoevaporation in removing gas mass from the disk.

Alternatively, external perturbers in the form of companions or fly-bys could also contribute to dispersing disk masses through tidal interactions (see S. S. R. Offner et al. 2023;
<figure><img src="image_11.png" /><figcaption>**Figure 8.** Spectral index between 230 and 283 GHz (1.3 and 1.06 mm) vs. stellar mass, representative of the wavelengths covered by the AGE-PRO ALMA Band 6 and Band 7 observations. The colored regions enclose the minimum and maximum spectral indexes for each initial disk mass group. The squares show the spectral index as measured from the AGE-PRO observations for the Lupus and Upper Sco SFR.</figcaption></figure>

A. Zurlo et al. 2023, for reviews), especially at early ages when SFRs have a higher number density of stars. Several systems from the literature have been reported as having gas morphologies consistent with simulations of interactions over the last few hundred years (e.g., N. T. Kurtovic et al. 2018; J. E. Rodriguez et al. 2018; F. Ménard et al. 2020; L. A. Zapata et al. 2020); systems that seem isolated today could have suffered from a fly-by tens of thousands of years ago or earlier (e.g., A. J. Winter et al. 2018). Thus, the large range of $M_{\text{gas}}$ observed in AGE-PRO might originate from a combination of internal and external mechanisms.

Among the external mechanisms that could modify the mass of a disk, recent population synthesis simulations have shown that late infall of material could contribute to increasing the mass reservoir of a disk by more than 50% (M. Kuffmeier et al. 2023) within a disk lifetime. These infalls have already been observed toward young objects (e.g., H.-W. Yen et al. 2017; S. L. Grant et al. 2021; L. Cacciapuoti et al. 2024). These positive mass contributions put a higher tension between the AGE-PRO observations and explaining mass loss through only viscous evolution models.

Despite the low values of $M_{\text{gas}}$ that are not reproduced by our models, the comparison between $M_{\text{pebble}}$ (or $F_{\text{mm}}$) shows a better agreement, suggesting that early dispersal mechanisms at <1 Myr are affecting the gas more than the dust component. Tentatively, dust growth and evolution through radial drift could be protecting the dust reservoir at more compact radii in the midplane of disks compared to gas, which extends much

farther from the star. Since $M_{\text{gas}}$ is higher in our simulations, it implies that there are high dust-to-gas ratios in the observations that we do not reproduce in our models. The impact of this discrepancy on the trends derived from dust observables remains to be investigated.

## 4.2. Trends in AGE-PRO that Are Reproduced by Dust Evolution Models

We find a better agreement when we compare the properties of dust continuum emission of the AGE-PRO observations to the dust evolution simulations with either strong or weak traps. The combination of disk size in continuum emission, as well as the continuum brightness, allows us to exclude disks with no pressure traps as the dominant population in the observations of AGE-PRO, in line with the previous results from A. Zormpas et al. (2022) and L. Delussu et al. (2024).

At the very early disk stages, traced by the sources in Ophiuchus, the observed values of $F_{\text{mm}}$ or $R_{90\%\text{,mm}}$ are consistent with any of the trapping scenarios, suggesting there has not been enough time for dust radial drift to deplete the disks of solid material. Thus, the solid content at early stages (<1 Myr) might be a good tracer of the initial dust content of disks if the young disks have had a constant low viscosity in their midplanes. By the age of Lupus and Upper Sco, the disks are more consistent with the simulations with either weak or strong dust trapping, as those with no traps have already lost a significant amount of its solid reservoir. It is relevant to mention that our disks have static traps over time and location,
<figure><img src="image_12.png" /><figcaption>**Figure 9.** Observed spectral index between 230 and 283 GHz (1.3 and 1.06 mm), compared to other observational properties. The dashed lines in panel (i) show the median spectral index for each SFR.</figcaption></figure>

which contributes to trapping dust more efficiently. Variability in the dust trapping strength or location would result in more leaky traps, which would make the results more similar to those of simulations without traps (P. Pinilla et al. 2021b).

When comparing observations to radiative transfer models, it is relevant to consider that the contribution of small grains (<100 µm) with the opacity from L. Ricci et al. (2010) is higher than that from other opacities, such as the DSHARP opacities (T. Birnstiel et al. 2018; see, for example, J. Stadler et al. 2022), thus the radiative transfer image of simulations with very little pebble mass (<0.1 $M_⊕$) can have detectable millimeter brightness, considering the AGE-PRO sensitivity. Even though the millimeter flux of simulations with no traps at 10 Myr of age is comparable to fluxes of some Upper Sco disks, the $R_{90\%}$ of the simulations does not coincide with the observed radii. Further studies considering multimillimeter wavelength observations will be able to set stronger constraints on the grain size distribution (e.g., A. Sierra et al. 2024), thus allowing for better direct comparisons between radiative transfer models and the AGE-PRO observations.

Previous studies have shown that there is a decreasing trend in dust continuum millimeter brightness as a function of age (N. Hendler et al. 2020; C. F. Manara et al. 2023), which is also observed within the AGE-PRO sample (see Figure 5). This decreasing trend with age is also reproduced in our simulations in $M_{\text{pebble}}$ or $F_{1.3\text{mm}}$, but it can be significantly slowed down for disks with strong traps, thus increasing the range of observed millimeter fluxes for a given SFR. The disks

with high dust content and older than 3 Myr favor the scenario of traps, when comparing their dust masses to those of the dust evolution simulations.

When considering the evolution of size in disks, we find that the connection between $R_{90\%_{\text{pebble}}}$ and $R_{90\%_{\text{flux}}}$ has a high dependency on the choice of opacities, as also discussed by J. Stadler et al. (2022) and L. Delussu et al. (2024). Observationally, the sizes of disks are typically measured from millimeter continuum emission, and they are used as a tracer for dust mass distribution. However, these values can have large discrepancies if most of the dust mass is contained in small grains (<100 µm). For example, the $R_{90\%_{\text{pebble}}}$ (which only considers grains >100 µm) in disks without traps shows a monotonically decreasing trend as a function of age, which is a natural result of dust radial drift. However, when considering the $R_{90\%}$ of the disks in millimeter brightness, we see these disks becoming larger as a function of time, as the contributions from small grains in the outer part of the disks become the dominant millimeter flux source.

The distinction of pebbles used in this work is relevant to compare with previous works analyzing the evolution of disk sizes. For example, G. P. Rosotti et al. (2019) found that the radius containing the total dust mass of the disk grows as a function of time in the absence of traps, as a result of micron-sized grains being coupled to the expanding gas in viscous evolution. This is in contrast to the behavior of our $R_{90\%_{\text{pebble}}}$ in the absence of traps, which does not consider the contribution of micron-sized grains to the pebble-mass
content. Depending on the millimeter opacity of these small
grains (<100 µm-sized), the radiative transfer models of disks
with no traps can show shrinking radii as a function of time (as
in G. P. Rosotti et al. 2019), or increasing radii as a function of
time (this work; see Figure 7); thus, dust opacity is one of the
most relevant sources of uncertainty when comparing models
to observations (see also L. Delussu et al. 2024).

In the presence of traps, the connection between R<sub>90%,pebble</sub> and R<sub>90%,mm</sub> becomes more direct. The value of R<sub>90%,pebble</sub> is completely determined by the location of the outermost dust trap, which is also true for R<sub>90%,mm</sub> as long as the millimeter flux contribution of the material in the trap is larger than that of the micron-sized grains in the outer disk. Thus, the size of disks in millimeter emission might not be a good tracer of disk evolution, consistent with the findings of G. D. Mulders et al. (2017). In L. Trapman et al. (2025b), a positive correlation between M<sub>gas</sub> and R<sub>90%,mm</sub> is detected. If R<sub>90%,mm</sub> is indeed tracing the location of the outermost trap, this relation could be evidence of more-massive disks having a higher potential for generating dust traps at larger radii. Additional observations at higher angular resolution across the whole AGE-PRO sample are needed to test this hypothesis.

Finally, we observe a large spread in the spectral indices both in Lupus and in Upper Sco, as presented in Figure 9. From our simulations, an increase in the spectral index is expected as a function of age. In the absence of traps, this increase happens very early in a disk age (about 2 Myr), and it is only in the presence of dust traps that the spectral index distribution can have a distribution ranging from 2–4 for later ages, as observed by AGE-PRO. Even though the role of self-scattering, optical depth, and free–free emission needs to be further investigated with additional multiwavelength observations, our measurements of spectral indices offer additional support for the presence of either weak or strong traps in the AGE-PRO sample.

**5. Conclusions**

In this work, we compared dust evolution models to the AGE-PRO observations. Our models considered similar stellar masses and luminosities as the host stars of the AGE-PRO sample and included three different initial gas disk masses and four initial disk characteristic radii. Additionally, we assumed models that either do not have any dust traps or models that have weak or strong dust traps. The results of the dust evolution are taken to perform radiative transfer calculations, considering five different evolution times $t = [0.5, 1, 2, 5, 10]$ Myr. These outputs are compared to observations of disks in Ophiuchus (compared to output at 0.5 Myr), Lupus (1 and 2 Myr), and Upper Sco (2 and 5 Myr). The main conclusions of our work are as follows:

1. The pebble masses and millimeter fluxes obtained at 0.5 Myr of evolution are representative of the initial dust reservoir for planet formation in disks, if the disks have had a low viscosity from the beginning. For this reason, when comparing our models to observations of the young disks in Ophiuchus, the data does not strongly favor models with or without substructures, as little dust loss has occurred due to drift at these early stages.

2. When comparing our models at more evolved stages with the Lupus and Upper Sco data, the observations (in

particular, the observed millimeter fluxes and spectral
indices) favor the presence of dust traps in disks, as the
simulations with no traps are unable to follow the
observational trends. The preference between weak and
strong traps depends on the time of evolution when we
compare the data, and also on each specific disk. Our
models suggest that the dust content of the disks in
Upper Sco have mostly survived until their current age
due to the presence of traps. When traps are present in
disks, moderate drift can still happen, decreasing the
millimeter flux over time, but allowing Rdust to remain
stable over time, as found in the AGE-PRO observations.

3. The dust disk size is determined by the location of the outermost dust trap, even at early evolution times (<1 Myr). The AGE-PRO observations show a positive trend between the R<sub>90%,mm</sub> and M<sub>gas</sub>. In the context of the models of this work, this would imply that more-massive disks have a higher potential to form traps farther out. In the absence of traps, the expected shrink of R<sub>90%,mm</sub> over time is not obtained in our synthetic observations: when drift removes millimeter-sized particles in the outer parts of the disks, the micron-sized particles start to dominate the millimeter emission, which could contribute to hiding the shrinking of disk size in high sensitivity observations, as it is the case for AGE-PRO.

4. The spatially integrated spectral index obtained for Lupus and Upper Sco sources in AGE-PRO shows a large spread, ranging from disks that are consistent with being completely optically thick, to others where the emission is more consistent with being optically thin at the wavelength of the AGE-PRO observations. In simulations with stronger traps, the spectral index can remain closer to $\alpha_{mm} \approx 2$ until later ages, compared to weaker traps. Thus, the spectral index comparison also favors the presence of traps.

**Acknowledgments**

N.T.K. acknowledges the support of Ewine van Dishoeck while working on this project. The authors also thank the feedback provided by the referee, which contributed to improve this manuscript. N.T.K., M.G., and P.P. acknowledge the 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. P.P. and A.S. acknowledge the support from the UK Research and Innovation (UKRI) under the UK government's Horizon Europe funding guarantee from ERC (under grant agreement No. 101076489). K.Z. and L.T. acknowledge the support of the NSF AAG grant No. 2205617. G.R. acknowledges funding from the Fondazione Cariplo, grant No. 2022-1217, and the European Research Council (ERC) under the European Union's Horizon Europe Research & Innovation Program under grant agreement No. 101039651 (DiscEvol). Views and opinions expressed are, however, those of the author(s) only, and do not necessarily reflect those of the European Union or the European Research Council Executive Agency. Neither the European Union nor the granting authority can be held responsible for them. B.T. acknowledges support from the Programme National "Physique et Chimie du Milieu Interstellaire" (PCMI) of CNRS/INSU with INC/INP and co-funded by CNES. I.P. and D.D. acknowledge support
from Collaborative NSF Astronomy & Astrophysics Research grant (ID: 2205870). L.A.C and C.G.R. acknowledge support from the Millennium Nucleus on Young Exoplanets and their Moons (YEMS), ANID-Center Code NCN2024_001 and L.A.C. also acknowledges support from the FONDECYT grant No. 1241056. A.S. also acknowledges support from FONDECYT de Postdoctorado 2022 No. 3220495. L.P. also acknowledges support from ANID BASAL project FB210003. C.A.G. and L.P. acknowledge support from FONDECYT de Postdoctorado 2021 No. 3210520.

This paper makes use of the following ALMA data: ADS/JAO.ALMA#2021.1.00128.L. ALMA is a partnership of ESO (representing its member states), NSF (USA), and NINS (Japan), together with NRC (Canada), MOST and ASIAA (Taiwan), and KASI (Republic of Korea), in cooperation with the Republic of Chile. The Joint ALMA Observatory is operated by ESO, AUI/NRAO, and NAOJ. The National Radio Astronomy Observatory is a facility of the National Science Foundation operated under cooperative agreement by Associated Universities, Inc.

# Appendix

## Synthetic Observations: Comparison with CLEAN Images

There are three main products to consider when generating an image from interferometric data with the CLEAN algorithm: (a) the model image, which is generated from the visibilities; (b) the residual image, which represents the Fourier transform of the residual visibilities, whose amplitude should be comparable with the thermal noise background; and (c) the point-spread function (PSF) image, which shows the response of the visibility sampling to a point-source emission. In order to estimate the angular resolution of an observation, a Gaussian function is fitted to the PSF image, thus quantifying the smallest possible angular scale that can be recovered in the image. Spatial variations in the model image that are smaller than the angular resolution could be the result of the reconstruction process, and do not necessarily represent an existing structure detected by the observation. Thus, the model image is convolved with the Gaussian representative of the PSF shape, erasing brightness structures of higher spatial frequency than those recoverable with the visibility coverage. As a final step, the convolved model is added with the residual image, and the result is the reconstructed image of CLEAN.

In Section 3, properties such as total flux and disk size from radiative transfer models of the simulations have been compared to the AGE-PRO sources. In the case of radiative transfer models, those properties are easily obtained due to the absence of beam convolution or noise. However, for interferometric observations, they are more challenging to recover due to the limited angular resolution and sensitivity. Even though measurements with visibility modeling of the dust continuum emission contribute to increase the accuracy of results such as total flux and disk size, and partially overcoming the issues related to image reconstruction (see M. Vioque et al. 2025), they are still limited to the physical constraints of the data set, such as the signal-to-noise ratio at the longest baselines.

For additional testing of the influence of the beam convolution and the thermal noise in the measured radii of the disk, we generated synthetic observations directly in the image plane, using the radiative transfer image as the CLEAN model (J. A. Högbom 1974; B. G. Clark 1980), and convolving

it with a Gaussian kernel representative of the angular resolution. We added a thermal noise image with a comparable standard deviation to that of the AGE-PRO observations. In comparison with the CLEAN algorithm, the main difference is in the method used to generate the underlying model. With CLEAN, the model is constructed as a collection of point sources and Gaussians, while in our simulated observations, the model comes directly as an output of the radiative transfer.

The AGE-PRO observations have an average angular resolution of 0.25 and 0.8 in the 1.3 mm and 1.05 mm setups; thus, we use those values as the FWHM of the Gaussian kernel to convolve the radiative transfer images. In the CLEAN algorithm, the sensitivity is set by the standard deviation of the intensity distribution in the residual map. In the idealized scenario where all of the real emission is extracted in the CLEAN model, then the intensity of the CLEAN residuals should distribute as thermal noise. Thus, we represent the limited sensitivity in our synthetic observations by adding a thermal noise image, generated by assigning fluxes to each pixel following a Gaussian distribution centered at 0 and a standard deviation of $\mu = 1$. An example of a thermal noise image is shown in panel (a) of Figure 10. This thermal noise image is then convolved with the same beam as the radiative transfer image, and normalized to match a sensitivity of 25 $\mu$Jy beam$^{-1}$, comparable to AGE-PRO (see panel (b) of Figure 10). The convolved thermal noise image is then added to the convolved model image, thus obtaining a synthetic observation, as shown in panel (d) of Figure 10. From these final images, we obtain the synthetic disk fluxes, sizes, and infer spectral indices.

We calculated the $R_{90\%\_cont}$ from the synthetic observations similarly to the AGE-PRO observations (L. Trapman et al. 2025a). The disk size is measured by calculating an azimuthally symmetric radial profile ($I_\nu := I_\nu(r)$), and finding the radii that includes 90% of the total flux. However, the measurement of the total flux is dependent on the size of the elliptical aperture where flux is integrated. The emission from the disks would only contribute positively to an intensity profile, but the contribution from thermal noise will oscillate between positive and negative values. Thus, the radius containing 100% of a disk flux ($R_{100\%}$) is estimated from the radius at which the azimuthally averaged radial profile transitions from a positive value of $I_\nu$ to a negative value. The $R_{100\%}$ can also be interpreted as the radius at which the radial profile goes from being dominated by disk emission to a noise-dominated regime. An equivalent method is to find the $R_{100\%}$ by looking at the cumulative flux radial profile of a disk, and finding the radius at which the first derivative of $I_\nu$ becomes 0.

When measuring $R_{100\%}$ from the brightness profiles, an important caveat should be considered for disks with large depleted gaps or cavities, where the emission could become noise dominated before reaching the next emission substructure. In the synthetic images of our simulated disks, however, this is never the case, as there is no cavity or depleted gap larger than the beam size. Due to the relatively large beam compared to the disk sizes in the AGE-PRO observations, none of the cavities are resolved as empty either.

The synthetic images generated from the radiative transfer simulations do not have an inclination or rotation applied, and thus, they simulate disks that have been observed in a face-on configuration. Since the beam convolution is done with a circular Gaussian kernel, there is no distortion of the flux in a
<figure><img src="image_13.png" /><figcaption>**Figure 10.** Example of the different components of a synthetic observation. Panel (a) shows an image of the pixel noise following a Gaussian distribution. Panel (b) shows the same noise image, but convolved with a Gaussian of 0.25 of FWHM, and normalized to a standard deviation of 40μJy beam⁻¹. Panel (c) shows the radiative transfer model at 1.3 mm of one of the simulations. Panel (d) shows the model convolved with the 0.25 Gaussian, after adding the convolved noise, mimicking the output image of the CLEAN algorithm.</figcaption></figure>

<figure><img src="image_14.png" /><figcaption>**Figure 11.** Comparison between the R<sub>90</sub> in the radiative transfer images (R<sub>90,mod</sub>) and in the synthetic observations (R<sub>90,obs</sub>). Each dot represents a different disk, and images from all ages are included. The dots are color-coded by the total flux of the model.</figcaption></figure>

preferred axis, as is the case with interferometric observations, where there is typically one beam axis larger than the other. Thus, our measurements are not affected by either disk geometry or beam asymmetries. The measurement of disk sizes using the above-described method proves robust for bright disks that are larger than approximately 3 times the angular resolution of the image (the FWHM of the beam in our synthetic observations is 0.250), as shown in Figure 11.

## ORCID iDs

Nicolás T. Kurtovic https://orcid.org/0000-0002-2358-4796

Matías Gárate https://orcid.org/0000-0001-6802-834X

Paola Pinilla https://orcid.org/0000-0001-8764-1780

Ke Zhang https://orcid.org/0000-0002-0661-7517

Giovanni P. Rosotti https://orcid.org/0000-0003-4853-5736

Rossella Anania https://orcid.org/0009-0004-8091-5055

Ilaria Pascucci https://orcid.org/0000-0001-7962-1683

Benoît Tabone https://orcid.org/0000-0002-1103-3225

Leon Trapman https://orcid.org/0000-0002-8623-9703

Dingshan Deng https://orcid.org/0000-0003-0777-7392

Miguel Vioque https://orcid.org/0000-0002-4147-3846

John Carpenter https://orcid.org/0000-0003-2251-0602

Lucas A. Cieza https://orcid.org/0000-0002-2828-1153

Laura M. Pérez https://orcid.org/0000-0002-1199-9564

Carolina Agurto-Gangas https://orcid.org/0000-0002-7238-2306

Anibal Sierra https://orcid.org/0000-0002-5991-8073

Dary A. Ruiz-Rodriguez https://orcid.org/0000-0003-3573-8163

James Miley https://orcid.org/0000-0002-1575-680X

Camilo González-Ruilova https://orcid.org/0000-0003-4907-189X

Estephani Torres-Villanueva https://orcid.org/0000-0001-9961-8203

Aleksandra Kuznetsova https://orcid.org/0000-0002-6946-6787

## References

Agurto-Gangas, C., Pérez, L. M., Sierra, A., et al. 2025, ApJ, 989, 4

Alexander, R., Pascucci, I., Andrews, S., Armitage, P., & Cieza, L. 2014, in Protostars and Planets VI, ed. H. Beuther et al. (Tucson, AZ: Univ. Arizona Press), 475

Anania, R., Rosotti, G., Gárate, M., et al. 2025, ApJ, 989, 8

Andrews, S. M. 2020, ARA&A, 58, 483

Andrews, S. M., Rosenfeld, K. A., Kraus, A. L., & Wilner, D. J. 2013, ApJ, 771, 129

Ansdell, M., Williams, J. P., van der Marel, N., et al. 2016, ApJ, 828, 46

Bae, J., Isella, A., Zhu, Z., et al. 2023, in ASP Conf. Ser. 534, Protostars and Planets VII, ed. S. Inutsuka et al. (San Francisco, CA: ASP), 423

Barenfeld, S. A., Carpenter, J. M., Ricci, L., & Isella, A. 2016, ApJ, 827, 142

Bate, M. R. 2018, MNRAS, 475, 5618

Birnstiel, T. 2024, ARA&A, 62, 157

Birnstiel, T., Dullemond, C. P., & Brauer, F. 2010, A&A, 513, A79

Birnstiel, T., Dullemond, C. P., Zhu, Z., et al. 2018, ApJL, 869, L45

Birnstiel, T., Fang, M., & Johansen, A. 2016, SSRv, 205, 41

Cacciapuoti, L., Macias, E., Gupta, A., et al. 2024, A&A, 682, A61

Clark, B. G. 1980, A&A, 89, 377

Clarke, C. J., Gendrin, A., & Sotomayor, M. 2001, MNRAS, 328, 485

Das, S., Kurtovic, N. T., & Flock, M. 2024, A&A, 689, A104

Delage, T. N., Gárate, M., Okuzumi, S., et al. 2023, A&A, 674, A190

Delage, T. N., Okuzumi, S., Flock, M., Pinilla, P., & Dzyurkevich, N. 2022, A&A, 658, A97

Delussu, L., Birnstiel, T., Miotello, A., et al. 2024, A&A, 688, A81

Deng, D., Vioque, M., Pascucci, I., et al. 2025, ApJ, 989, 3

Doi, K., & Kataoka, A. 2023, ApJ, 957, 11

Dominik, C., Min, M., & Tazaki, R. 2021, OpTool: Command-line Driven Tool for Creating Complex Dust Opacities, Astrophysics Source Code Library, ascl:2104.010

Dotter, A., Chaboyer, B., Jevremović, D., et al. 2008, ApJS, 178, 89

Draine, B. T. 2003, ARA&A, 41, 241
Drążkowska, J., Bitsch, B., Lambrechts, M., et al. 2023, in ASP Conf. Ser. 534, Protostars and Planets VII, ed. S. Inutsuka et al. (San Francisco, CA: ASP), 717

Dubrulle, B., Morfill, G., & Sterzik, M. 1995, *Icar*, **114**, 237

Dullemond, C. P., Juhasz, A., Pohl, A., et al. 2012, RADMC-3D: A Multi-purpose Radiative Transfer Tool, Astrophysics Source Code Library, ascl:1202.015

Fedele, D., van den Ancker, M. E., Henning, T., Jayawardhana, R., & Oliveira, J. M. 2010, A&A, **510**, A72

Fromang, S., & Nelson, R. P. 2009, A&A, **496**, 597

Gaia Collaboration, Brown, A. G. A., Vallenari, A., et al. 2021, A&A, **649**, A1

Gárate, M., Birnstiel, T., Drążkowska, J., & Stammler, S. M. 2020, A&A, **635**, A149

Gárate, M., Birnstiel, T., Pinilla, P., et al. 2023, A&A, **679**, A15

Gárate, M., Delage, T. N., Stadler, J., et al. 2021, A&A, **655**, A18

Garufi, A., Carrasco-Gonzalez, C., Macias, E., et al. 2025, A&A, **694**, A290

Grant, S. L., Espaillat, C. C., Wendeborn, J., et al. 2021, ApJ, **913**, 123

Gundlach, B., & Blum, J. 2015, ApJ, **798**, 34

Gundlach, B., Kilias, S., Beitz, E., & Blum, J. 2011, *Icar*, **214**, 717

Gundlach, B., Schmidt, K. P., Kreuzig, C., et al. 2018, MNRAS, **479**, 1273

Hendler, N., Pascucci, I., Pinilla, P., et al. 2020, ApJ, **895**, 126

Hildebrand, R. H. 1983, QJRAS, **24**, 267

Högbom, J. A. 1974, A&AS, **15**, 417

Huang, J., Andrews, S. M., Cleeves, L. I., et al. 2018, ApJ, **852**, 122

Kuffmeier, M., Jensen, S. S., & Haugbølle, T. 2023, EPJP, **138**, 272

Kurtovic, N. T., Pérez, L. M., Benisty, M., et al. 2018, ApJL, **869**, L44

Kurtovic, N. T., Pinilla, P., Long, F., et al. 2021, A&A, **645**, A139

Lesur, G., Flock, M., Ercolano, B., et al. 2023, in ASP Conf. Ser. 534, Protostars and Planets VII, ed. S. Inutsuka et al. (San Francisco, CA: ASP), 465

Liu, H. B. 2019, ApJL, **877**, L22

Lodato, G., Scardoni, C. E., Manara, C. F., & Testi, L. 2017, MNRAS, **472**, 4700

Long, F., Andrews, S. M., Rosotti, G., et al. 2022, ApJ, **931**, 6

Long, F., Pinilla, P., Herczeg, G. J., et al. 2018, ApJ, **869**, 17

Long, F., Pinilla, P., Herczeg, G. J., et al. 2020, ApJ, **898**, 36

Lynden-Bell, D., & Pringle, J. E. 1974, MNRAS, **168**, 603

Manara, C. F., Ansdell, M., Rosotti, G. P., et al. 2023, in ASP Conf. Ser. 534, Protostars and Planets VII, ed. S. Inutsuka et al. (San Francisco, CA: ASP), 539

Manara, C. F., Rosotti, G., Testi, L., et al. 2016, A&A, **591**, L3

Mathis, J. S., Rumpl, W., & Nordsieck, K. H. 1977, ApJ, **217**, 425

Ménard, F., Cuello, N., Ginski, C., et al. 2020, A&A, **639**, L1

Mulders, G. D., Pascucci, I., Manara, C. F., et al. 2017, ApJ, **847**, 31

Musiolik, G., & Wurm, G. 2019, ApJ, **873**, 58

Nakagawa, Y., Sekiya, M., & Hayashi, C. 1986, *Icar*, **67**, 375

Offner, S. S. R., Moe, M., Kratter, K. M., et al. 2023, in ASP Conf. Ser. 534, Protostars and Planets VII, ed. S. Inutsuka et al. (San Francisco, CA: ASP), 275

Pascucci, I., Cabrit, S., Edwards, S., et al. 2023, in ASP Conf. Ser. 534, Protostars and Planets VII, ed. S. Inutsuka et al. (San Francisco, CA: ASP), 567

Pascucci, I., Testi, L., Herczeg, G. J., et al. 2016, ApJ, **831**, 125

Pinilla, P., Benisty, M., & Birnstiel, T. 2012, A&A, **545**, A81

Pinilla, P., Birnstiel, T., Benisty, M., et al. 2013, A&A, **554**, A95

Pinilla, P., Kurtovic, N. T., Benisty, M., et al. 2021a, A&A, **649**, A122

Pinilla, P., Lenz, C. T., & Stammler, S. M. 2021b, A&A, **645**, A70

Pinilla, P., Pascucci, I., & Marino, S. 2020, A&A, **635**, A105

Pringle, J. E. 1981, ARA&A, **19**, 137

Ribas, Á., Bouy, H., & Merín, B. 2015, A&A, **576**, A52

Ricci, L., Testi, L., Natta, A., et al. 2010, A&A, **512**, A15

Rodriguez, J. E., Loomis, R., Cabrit, S., et al. 2018, ApJ, **859**, 150

Rosotti, G. P. 2023, NewAR, **96**, 101674

Rosotti, G. P., Tazzari, M., Booth, R. A., et al. 2019, MNRAS, **486**, 4829

Ruiz-Rodriguez, D. A., González-Ruilova, C., Cieza, L. A., et al. 2025, ApJ, **989**, 2

Sanchis, E., Testi, L., Natta, A., et al. 2021, A&A, **649**, A19

Shakura, N. I., & Sunyaev, R. A. 1973, A&A, **24**, 337

Sierra, A., Pérez, L. M., Sotomayor, B., et al. 2024, ApJ, **974**, 306

Somigliana, A., Testi, L., Rosotti, G., et al. 2023, ApJL, **954**, L13

Somigliana, A., Toci, C., Rosotti, G., et al. 2022, MNRAS, **514**, 5927

Stadler, J., Gárate, M., Pinilla, P., et al. 2022, A&A, **668**, A104

Stammler, S. M., & Birnstiel, T. 2022, ApJ, **935**, 35

Steinpilz, T., Teiser, J., & Wurm, G. 2019, ApJ, **874**, 60

Tabone, B., Rosotti, G. P., Lodato, G., et al. 2022, MNRAS, **512**, L74

Tabone, B., Rosotti, G. P., Trapman, L., et al. 2025, ApJ, **989**, 7

Takeuchi, T., & Lin, D. N. C. 2002, ApJ, **581**, 1344

Tazzari, M., Testi, L., Natta, A., et al. 2021, MNRAS, **506**, 5117

Testi, L., Natta, A., Manara, C. F., et al. 2022, A&A, **663**, A98

Toci, C., Lodato, G., Livio, F. G., Rosotti, G., & Trapman, L. 2023, MNRAS, **518**, L69

Trapman, L., Vioque, M., Kurtovic, N., et al. 2025a, ApJ, **989**, 10

Trapman, L., Zhang, K., Rosotti, G., et al. 2025b, ApJ, **989**, 5

Tripathi, A., Andrews, S. M., Birnstiel, T., & Wilner, D. J. 2017, ApJ, **845**, 44 van der Marel, N., & Pinilla, P. 2023, arXiv:2310.09077

Villenave, M., Ménard, F., Dent, W. R. F., et al. 2020, A&A, **642**, A164

Villenave, M., Stapelfeldt, K. R., Duchêne, G., et al. 2022, ApJ, **930**, 11

Vioque, M., Kurtovic, N. T., Trapman, L., et al. 2025, ApJ, **989**, 9

Wada, K., Tanaka, H., Suyama, T., Kimura, H., & Yamamoto, T. 2011, ApJ, **737**, 36

Warren, S. G., & Brandt, R. E. 2008, JGRD, **113**, D14220

Weber, P., Benítez-Llambay, P., Gressel, O., Krapp, L., & Pessah, M. E. 2018, ApJ, **854**, 153

Weidenschilling, S. J. 1977, MNRAS, **180**, 57

Winter, A. J., Booth, R. A., & Clarke, C. J. 2018, MNRAS, **479**, 5522

Yen, H.-W., Takakuwa, S., Chu, Y.-H., et al. 2017, A&A, **608**, A134

Zapata, L. A., Rodríguez, L. F., Fernández-López, M., et al. 2020, ApJ, **896**, 132

Zhang, K., Pérez, L. M., Pascucci, I., et al. 2025, ApJ, **989**, 1

Zormpas, A., Birnstiel, T., Rosotti, G. P., & Andrews, S. M. 2022, A&A, **661**, A66

Zubko, V. G., Mennella, V., Colangeli, L., & Bussoletti, E. 1996, MNRAS, **282**, 1321

Zurlo, A., Gratton, R., Pérez, S., & Cieza, L. 2023, EPJP, **138**, 411