Open Access
Issue
A&A
Volume 671, March 2023
Article Number A125
Number of page(s) 14
Section Interstellar and circumstellar matter
DOI https://doi.org/10.1051/0004-6361/202244869
Published online 17 March 2023

© The Authors 2023

Licence Creative CommonsOpen Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

This article is published in open access under the Subscribe to Open model.

Open access funding provided by Max Planck Society.

1 Introduction

Protoplanetary disks around young stars are the nursery of planetary systems. These disks are almost entirely composed of gas and a small mass fraction of dust. Small dust grains (~1 µm) are well coupled to the gas structure, extending, thus, vertically over several scale heights, following the gas distribution. Dust settling in turbulent disks was extensively studied (e.g., Dubrulle et al. 1995; Schräpler & Henning 2004; Fromang & Nelson 2009; Zsom et al. 2011; Woitke et al. 2016). On the contrary, large millimeter grains settle much faster towards the disk midplane, where the growth process is particularly efficient because of the high densities. By constraining the distribution of both large and small dust grains, we can extrapolate the grain’s vertical size distribution, probing fundamental disk properties such as the turbulent structure and the dust-to-gas ratio.

Disks tend to have a flared structure (Kenyon & Hartmann 1987; Bell et al. 1997; Chiang & Goldreich 1997). Flaring occurs when the midplane temperature falls off more slowly than r−1. The disk surface then flares outward with increasing radius as a consequence of vertical hydrostatic equilibrium. With the development of our observational techniques, we characterized several examples of disks with a flared structure (e.g., Pinte et al. 2016; Avenhaus et al. 2018; Villenave et al. 2020; Law et al. 2022). For instance, Law et al. (2022) recently observed a correlation between the CO emitting height and the disk size in a sample of disks observed with the Atacama Large Millimeter/submillimeter Array (ALMA). Using the dust emission, Avenhaus et al. (2018) measured the flaring of V4046 Sgr, RXJ 1615, and IM Lup with the Very Large Telescope (VLT) Spectro-Polarimetric High-contrast Exoplanet REsearch (SPHERE) instrument. These near-infrared polarized light images show unprecedented details of the disk’s surface layer by tracing the scattering of small, micron-sized dust grains. At these wavelengths, the disk is expected to be optically thick and the observed light is dominated by light scattered by the dust grains on the disk surface. Among these disks, IM Lup observations are arguably among the most impressive images of a protoplanetary disk, as part of the Disks ARound TTauri Stars with SPHERE (DARTTS-S) survey by Avenhaus et al. (2018), targeting eight TTauri stars, at 1.25 µm and 1.65 µm wavelengths. With its intermediate inclination (~56°), we have a three-dimensional (3D) perspective on the disk around IM Lup, which can be used to constrain the vertical disk structure and the small grain distribution in a model-independent way.

The large dust grain properties can be similarly constrained through observations. At (sub)millimeter wavelengths, ALMA data have provided a new understanding of the properties and distribution of large, millimeter-sized dust grains, found in the disk midplane (e.g., Andrews et al. 2016; Huang et al. 2017; Pinte et al. 2018). IM Lup was also part of the Disk Substructures at High Angular Resolution Project (DSHARP) survey (Andrews et al. 2018 and following papers), conducted with ALMA. The survey led to the characterization of substructures for 20 nearby protoplanetary disks using observations at 1.25 mm (240 GHz).

By combining ALMA and SPHERE images, it is possible to further characterize the structure of disks and examine how they affect the planet formation processes. The SPHERE and ALMA data show two very different geometries of the disk. In the near-infrared, the disk shows strong flaring and a prominent vertical structure. This implies that small grains are: well coupled to the gas structure, located in the upper layers of the disk, and prevented from settling (Avenhaus et al. 2018). The millimeter image, on the other hand, shows a flat disk, which suggests that the larger dust grains are well settled to the midplane. The disk emission differs also in its radial extent: the ALMA data extend up to 290 au, while the SPHERE data go as far as 400 au. These two images can be used to model the distribution of the largest and smallest dust grains, which we can then use to constrain the full grain size distribution.

Pinte et al. (2008) took a similar approach and combined infrared spectra and scattered light emission with millimeter thermal emission to constrain the parameter range of their IM Lup model. The data consist of scattered light detected with the WFPC2 instrument on board the Hubble Space Telescope (HST), along with Spitzer near and mid-infrared spectroscopy and SubMillimeter Array (SMA) millimeter emission. These authors combined a manual exploration of the model parameter space with a more rigorous fitting of the parameters that cannot be directly inferred from the observations. They found IM Lup to be a massive disk, with a mass of about ~0.1 M, assuming a 1% dust-to-gas ratio, extending up to ~400 au, but with a low Hα emission hinting at a low accretion rate. However, a low accretion rate would be a surprising result for a massive disk and these authors suggested that may be due to a period of low or moderate accretion. More recently, the accretion rate was measured by Alcalá et al. (2017, 2019), using the spectrum observed with the VLT/X-shooter spectrograph and the far-utraviolet continuum excess emission observed with the HST Cosmic Origins Spectrograph (HST-COS) and HST Space Telescope Imaging Spectrograph (HST-STIS). These studies reported a high mass-accretion rate of 10−8 M. Pinte et al. (2008) found that a dust population perfectly mixed with the gas cannot explain both the dust infrared spectral features and the millimeter continuum. Millimeter-sized grains must be present in the disk, but their emission cannot come from the same region as the infrared emission. The lower modeling power and quality of the data available at the time of this study, however, did not allow us to uniquely constrain the stratified dust disk structure; and so, other disk models apart from their best-fit model may match the applied data.

In this paper, we introduce a comprehensive modeling of protoplanetary disk dust distribution, focusing, in particular, on the vertical distribution of small grains. We then test this model by fitting it to the high-resolution millimeter ALMA data of IM Lup and near-infrared polarized radiation SPHERE data using a Markov chain Monte Carlo (MCMC) algorithm. This allows us to constrain the complete dust grain size distribution of the disk. Finally, we discuss our results and how they affect our understanding of the dust structure of protoplanetary disks.

2 Disk model

The gas distribution for our disk physical model is the axisymmmetric, self-similar Lynden-Bell & Pringle profile (Lynden-Bell & Pringle 1974):

(1)

where rc is the characteristic radius and γ is the gas surface density exponent. To recover the vertical gas distribution, we need to compute the vertical scale height of the disk. This scale height is expressed by the ratio of the sound speed, cs, to the Keplerian angular velocity Ω:

(2)

where k is the Boltzmann constant, Tmid(r) the midplane gas temperature distribution, µ = 2.3 the reduced mass of the H2 molecule, mH the hydrogen atomic mass, G the gravitational constant, and M the stellar mass. This equation comes from balancing the vertical pressure force of the gas with the gravitational force towards the midplane. Using this scale height, we obtain the gas density at a radius r and height z:

(3)

The two terms Hp(r) and Tmid(r) are coupled with each other: the amount of starlight penetrating to the disk midplane is determined by Hp. This radiation increases the midplane temperature Tmid, which results in a higher pressure scale height. Therefore, Hp and Tmid need to be consistent with each other, and the disk structure must be computed iteratively. To estimate the amount of heating from stellar irradiation, we need to compute the optical depth as a function of the location in the disk and the flaring angle φ(r), the angle at which the stellar radiation hits the disk surface. The flaring angle is, according to Chiang & Goldreich (1997):

(4)

where Hs is the surface height of the disk. We note that this is a different quantity than the pressure scale height, Hp, (which is usually lower) and is defined by the height above the midplane where the optical depth to the stellar radiation is unity.

Given the irradiation angle, φ(r), which is assumed to not depend on z, the optical depth to the stellar radiation is:

(5)

where κP, ⋆ is the Planck mean opacity at stellar wavelengths, ε is the assumed dust-to-gas ratio, and ρg is the gas density.

We now wish to find the value of z, at a given r, at which τ=1. Dullemond et al. (2001) found that this is equivalent to solving:

(6)

where erf is the error function:

(7)

This estimate for Hs can now be used to derive the flaring index from Eq. (4). Since the flaring angle determines Tmid, which determines Hp, these equations can be solved iteratively to derive Hp, Hs, φ, and Tmid consistently.

2.1 Vertical structure

To calculate the irradiation of the disk by the star in a 2D model, we use a ray-tracing technique, a 1+1D approach (Dullemond et al. 2002). This allows us to approximate a full 2D solution by calculating the transport of the stellar radiation through the disk, even though there is no radial transfer of thermal energy in the gas. This is a good approximation of the disk structure achievable within a reasonable computation time. Computational time is a strong constraint on our model: each new set of parameters we test against the observational data requires the calculation of a new disk structure and new radiative transfer calculations to produce the images to compare to the observations. A faster disk model allows us to test more parameters with the same computational resources and a better sampling of our posterior parameter distribution.

Ray-tracing consists of a simple integration on the lines ζ = z/r = const. The stellar flux can then rewritten as:

(8)

and the optical depth as:

(9)

From this radiation field, we can compute the temperature as:

(10)

where σSB is the Stephan-Boltzmann constant, κP⋆ the mean Planck opacity at stellar wavelength, and κP(T) is the mean Planck opacity of the dust at temperature, T. Values Jdiff and J are the mean intensity from dust re-emission and the intensity from stellar radiation, which can be derived from Eq. (8) and the dust opacity (discussed in Sect. 2.3). Once we have the temperature profile, we can solve for the vertical density profile from Eq. (3). This will be a different profile from the one used to solve Eq. (10), so we iterate these solutions until convergence.

2.2 Dust

We discuss above the fact that ALMA and SPHERE observations indicate that the disk around IM Lup is composed of at least two dust populations: one of small grains, coupled to the gas structure and vertically extended, and one of large grains settled to the midplane. These are just two parts of a continuous size distribution of grains, as predicted by dust coagulation models (e.g., Dullemond & Dominik 2005; Brauer et al. 2008; Birnstiel et al. 2010, 2012). Since a detailed dust coagulation and fragmentation model is computationally expensive, parametric functions to model the dust grain distribution are common in the literature.

We assume that the largest grain size at each radial location is given by a power-law distribution:

(11)

where a0 is the maximum grain size at radius, rc. From this maximum grain size at any radial location we can reconstruct the grain size distribution using:

(12)

where amin is 10−5 cm, amax is given by Eq. (11) and n(a) is normalized to the total dust volume density:

(13)

with m(a) the mass of a grain of radius a. The value of the index q in Eq. (12), one of our model parameters, has already been investigated in previous studies. Measurements of the interstellar extinction (e.g., Mathis et al. 1977) found q ≈ 3.5, which is also consistent with more recent submillimeter observations of debris disks (Ricci et al. 2015). However, the physics of gas-rich young protoplanetary disk is very different from the one of gas-depleted debris disks. Birnstiel et al. (2011) used dust coagulation and fragmentation models to fit an analytic multi-power law to dust evolution models.

The total dust surface density is derived from the gas density using a dust-to-gas ratio of:

(14)

where ε0 is the dust-to-gas ratio at the normalization radius, r0.

From this estimate of the dust density, we compute the vertical dust distribution by balancing turbulent diffusion and vertical settling. Assuming steady state, this is equivalent to solving (Dubrulle et al. 1995; Schräpler & Henning 2004; Fromang & Nelson 2009):

(15)

where τs is the grains stopping time, D is the turbulent diffusivity, and Ω the Keplerian angular velocity. A grain stopping time is the typical time it takes for a grain initially at rest to reach the local gas velocity. This depends on the grain bulk density, ρs, and its size, a:

(16)

where cs is the local sound speed. Different grain populations will have different vertical distributions, as larger grains quickly settle toward the midplane, while small grains remain coupled to the vertical structure of the gas. In our models, we divide the grain distribution from Eq. (12) in 30 populations of different grain sizes, logarithmically spaced between amin and amax, and solve their vertical distribution in Eq. (15).

The turbulent diffusivity, D, is a parametrization of how the grains are removed from the midplane by turbulent velocity fluctuations of the gas. The simplest case is to assume that it is constant:

(17)

where Sc is the Schmidt number, which has been measured to be of order one in zero net MHD turbulence (Johansen & Klahr 2005; Johansen et al. 2006). The diffusivity is intimately linked to the turbulence properties of the gas flows, therefore assuming it to be constant in space is a reasonable assumption when the turbulence is homogeneous. The most likely deviation from homogeneous turbulence is MHD effects, which are vertically and radially stratified (Dzyurkevich et al. 2013). As a consequence, the diffusivity coefficient, D, will change with the disk height, z, at any given radius. In our models, we do not explore MHD effects on the disk structure, and we use a constant diffusivity as in Eq. (17).

Table 1

DSHARP dust composition from Birnstiel et al. (2018)

2.3 Dust opacity

To compare our models to the observations, we produced synthetic emission maps of the dust structure. It requires us to calculate the opacity of the dust distribution derived in Sect. 2.2. We follow the approach of Birnstiel et al. (2018). This work analyzed the millimeter emission of DSHARP disks (including IM Lup) to study the dust properties of these sources.

The grain opacity is strongly affected by the grain composition. Here, we adopt the same composition used in the DSHARP project (Birnstiel et al. 2018), a mixture of water ice, astronomical silicates, troilite, and refractory organic material (see Table 1). Birnstiel et al. (2018) had assumed particles without porosity, however grains are expected to have larger porosities, at least in the initial stage of collisional growth (e.g., Krijt et al. 2015); therefore, we assumed grains with 30% porosity.

To calculate the dust mass absorption and scattering coefficients, and , as a function of frequency, v, we need the optical constants of the composing materials, n(λ) and k(λ), which depend on the observed wavelength λ. In Table 1 we report the dust composition assumed and the laboratory experiments providing the optical constants of the composing materials. Deriving the optical constants of a medium composed of multiple materials is a challenging task, as it depends on how we assume these materials are distributed in the dust grains. For general cases, some computationally expensive numerical models need to be used, but in some limit cases, analytic solutions are possible. For a homogeneous mix of the grain components, the solution is given by the Bruggemann rule (see Bohren & Huffman 1998 for more details). The resulting optical constants are shown in Fig. 1.

To calculate the size dependent opacities, and , from the optical constants in Fig. 1, we use the dshaгp_opac1 code, introduced by Birnstiel et al. (2018), based on the original Mie calculation code by Bohren & Huffman (1998).

To derive the total absorption opacity, we created a grid of 30 grain sizes, logarithmically spaced between 10−5 and 10 cm. The opacity (Eq. (12)) at frequency, v, is then averaged over the grain distribution n(α) from Eq. (12):

(18)

where αmin = 10−5 cm and amax is taken from Eq. (11). Similarly, we calculate the total scattering opacity, . Finally, to predict the map of the dust continuum and scattered light emission from the model described above using these opacities, we simulated the 3D radiative transfer with the Monte Carlo radiative transfer code RADMC–3D2 (Dullemond et al. 2012).

thumbnail Fig. 1

Medium optical constant for the dust composition in Table 1 (Birnstiel et al. 2018).

3 Observations

To constrain the vertical dust grain height in a disk, we need a disk meeting three requirements. First, the disk must be moderately inclined (~30–70°), so that its 3D structure can be reconstructed from its 2D projection on the sky. Secondly, it must show some ring structures in the scattered light emission to estimate the vertical extension of small grains as a function of radius. Finally, we required millimeter continuum observations to constrain the large grain distribution in the disk midplane. The faint, smooth, and extended emission observed in the scattered light in the disk around IM Lup makes this disk the ideal candidate for a test study. The methods developed here can then be applied in a future work to other axisymmetric disks and extended to include the effects of radial structures on the radial grain distribution due to mixing and drift, using methods similar to the ones presented in this work for the vertical grain distribution.

The protoplanetary disk around IM Lupus has been extensively studied in infrared scattered light (Pinte et al. 2008; Stapelfeldt et al. 2014; Avenhaus et al. 2018), the millimeter continuum (Pinte et al. 2008, 2018; Cleeves et al. 2016; Andrews et al. 2018), and gas line emission (Cleeves et al. 2016; Pinte et al. 2016). In particular, Zhang et al. (2021) used high-resolution ALMA observations to study CO isotopologue emission and retrieved the gas mass and distribution of several disks, including IM Lup. In our models, we assume the gas distribution derived in this work. This disk is the best candidate to test the simultaneous modeling of large and small grains, with the bonus of having extensive models of its gas structure to check the consistency of our results. Indeed, the efficient coupling of small grains to the gas implies similarities in their spatial distribution. To constrain our models, we use two sets of observations. One is the SPHERE measurements of the 1.65 µm near-infrared polarized emission, tracing the scattering of stellar radiation by small grains on the disk surface (Avenhaus et al. 2018). The second is the 1.25 mm continuum emission from the large grains in the disk midplane, detected by ALMA as part of the DSHARP disk survey (Andrews et al. 2018). A summary of the parameters of IM Lup can be found in Table 2.

There are two possible approaches to compare our model images to the observations: either we directly compare the model images to the observed images or we compare their radial emission profiles. We chose the latter to ensure we give the same weight to both the millimeter and near-infrared data when constraining the models. Indeed, a direct comparison of disk models to data from different instruments is a challenging task, especially in this case where the image resolution element is defined differently (beams for ALMA data and pixels for SPHERE data). On the other hand, using azimuthally averaged brightness profiles allows us to avoid the subtleties of working directly on the imaging data. Since either approach is valid and neither one offers a greater advantage over the other, we chose to work with the azimuthally averaged profiles to simplify our analysis.

In addition, we constrained our models using only the data outside of the 1 arcsec (158 au) radius. We are primarily interested in characterizing the difference between the small and large grain distribution, so we focus on the outer regions of the disk where this difference is maximized. Secondly, the SPHERE coronograph is obscuring the inner 0.1 arcsec (~ 16 au) of the disk , making the observed inner disk always dimmer than in the models. Finally, IM Lup observations show evidences both in the large grain (Huang et al. 2018a) and small grain emission (Avenhaus et al. 2018) up to 1 arcsec. Here, we do not address the effect these structures have on the dust distribution, which will be addressed in a follow-up study, and we focus on the region outside the 1 arcsec radius.

Table 2

IM Lupi stellar and disk properties.

Table 3

IM Lup small grain ring parameters, as in Avenhaus et al. (2018).

3.1 Near-infrared scattered light: Height of small grain distribution

The disk around IM Lup has been observed in its scattered light emission by HST (Pinte et al. 2008) and SPHERE (Avenhaus et al. 2018), as shown in Fig. 2. Scattered light images of an inclined disk such as IM Lup are expected to show azimuthal asymmetries due to the anisotropic nature of scattering by dust grains. These anisotropies are caused by both grain properties and the flared disk geometry. Grain scattering favors scattering in the forward direction, and the disk region closer to the observer is the brightest (e.g., Mulders et al. 2013). Moreover, the fraction of light that gets polarized by grain scattering depends on the grain structure, composition, and scattering angle.

The disk geometry also affects how much light is scattered by the disk surface. These images trace stellar radiation scattered by the small grains near the disk surface, because the disk is optically thick at infrared wavelengths (e.g., Grady et al. 2000, 2009; Hashimoto et al. 2011; Monnier et al. 2017) and a flared disk geometry can cause an excess in its infrared emission.

IM Lup scattered light emission has its outer edge at ≈340 au, and a faint halo extending up to ≈ 700 au (Avenhaus et al. 2018). The SPHERE data reveal several ring structures. If we assume these rings to be circular and centered on the star, we can fit ellipses to the peak flux of the rings to estimate radius, inclination, position angle, and the vertical displacement from the midplane due to the height of the small grains. These measurements are carried out in Avenhaus et al. (2018) and reported in Table 3 as well for convenience.

They find that the dust scattering surface is described by a power law in radius:

(19)

where H0 is the dust scale height at the fiducial radius r0, and β the flaring index. The height of the scattering surface is typically about three to four times the pressure scale height.

We use these results to extrapolate the radial scattered light emission profile of IM Lup. The flaring index gives us the information necessary to reconstruct the disk 3D structure from its projected 2D image on the sky plane, the SPHERE data. This can be done using the GoFish package (Teague 2019). To preserve the information about the non-isotropic nature of the polarized scattered light emission, we extract the emission profiles in four different directions. As shown in Fig. 2, we chose four cones 10º wide in PA in the forward and backward scattering direction, and the two side directions 90º away. We appliued this same procedure to the simulated disks, then we compared the models to the observations by fitting these four profiles.

thumbnail Fig. 2

IM Lup polarized emission at 1.65 µm wavelength observed by SPHERE (Avenhaus et al. 2018). The cones show the regions over which we extract the radial profiles we use to compare our models.

3.2 Millimeter continuum emission: Large grains in the disk midplane

The second set of observational constraints is the high resolution (~5 au) millimeter continuum emission from large grains in the disk midplane, at λobs = 1.25 mm. These data come from the ALMA DSHARP program (Andrews et al. 2018). Since both our model and the continuum emission are azimuthally axisymmetric, we compared our models to the azimuthal average of the observations. As with the scattered light emission, the millimeter emission shows radial structures in the inner disk, a spiral pattern extending from 25 au to 110 au (Huang et al. 2018b). Outside of this region, there is a gap in the emission profile. For this reason, we constrained our models with the observed emission profile outside the 1 arcsec radius, as with the scattered light data. The 1.25 mm emission has a sharp truncation at 295 au, also found in previous studies (Pinte et al. 2018), and extends up to about 400 au.

Tazzari et al. (2021) analyzed the 0.9 mm, 1.3 mm, and 3.1 mm emission for the 26 brightest disks in the Lupus region, including IM Lup. One of their most interesting results is that they measure very little variation in disk size across different wavelengths. This behavior is common for highly structured disks, such as IM Lup and the other DSHARP sources. For instance, it was also described by Pinilla et al. (2015, 2021), who studied the effect of disk radial structures on the dust grain size distribution. The most obvious explanation would be a physical drop in the dust density, but an emission drop can also be explained by an opacity feature. For example, Rosotti et al. (2019) and Powell et al. (2017) argued that the outer edge of the emission profile at a given millimeter wavelength is also the location where the maximum grain size corresponds to the opacity resonance at that wavelength (αmax ~ λ/2π). When the maximum grain size, decreasing with radius, falls below this size, the grain opacity drops quickly, and we observe a sharp truncation in the emission profile. This is in accord with the predictions of dust growth models and radial drift (e.g., Birnstiel et al. 2016; Powell et al. 2019; Franceschi et al. 2022). However, if this truncation were caused by radial drift, at different wavelengths we would observe truncation at different radii since radial drift affects different grain sizes with different efficiency. This is not what was measured by Tazzari et al. (2021), suggesting that the outer edge of the continuum emission is caused by a drop in the dust density – and not by an opacity feature.

While Tazzari et al. (2021) conclude that observations of millimeter continuum emission are not enough to distinguish which scenario is most likely in IM Lup, Avenhaus et al. (2018) are more specific and identify a outer edge at 340 au in the dust scattered light emission. This structure could trap grains at its location, and the grains could follow a different distribution inside and outside the ring. In our model, we recover the small grain distribution from the large grain distribution using the power law distribution in Eq. (12). However, if the grain distribution changes at the two sides of the ring and there is no large grain emission in the outer disk, we cannot recover the small grain abundance using this approach. Therefore, in our fitting routine we do not include the region outside the ring location, at 340 au.

4 Inferred disk model parameters

Starting from the geometrical structure of IM Lup, inferred directly from the observations as discussed in Sect. 3, and previous modeling efforts of its gas structure (Zhang et al. 2021), we can put quantitative constraints on the disk physical parameters using the Markov chain Monte Carlo sampler emcee3. The physical distributions explored by the fit are the grain size distribution, the maximum grain size radial distribution, the dust-to-gas ratio radial distribution, and the turbulent viscosity parameter αturb, according to the standard Shakura & Sunyaev (1973) viscosity recipe.

The grain size distribution is given by Eq. (12):

(20)

with the exponent of the distribution q as a fitting parameter. The maximum grain size radial distribution is given by Eq. (11):

(21)

with rc = 300 au and a0, p as fitting parameters. The dust-to-gas ratio radial distribution is given by Eq. (14):

(22)

with rin = 160 au and ε0 and γ as our fitting parameters. We assume the dust-to-gas ratio. and the maximum grain size distribution to have the same rcutoff. Finally, we have the turbulent viscosity parameter, αturb. This parameter sets the gas viscosity, v, as shown in Eq. (17). The turbulent diffusivity is used to compute the settling velocity and the turbulent diffusion coefficient of the dust species. These are used in turn to compute their settling mixing equilibrium, as discussed in Sect. 2.2.

This results in a model with six independent parameters4. The fitting is based on the definition of a reduced ;χ2. For the ALMA 1.25 mm emission, we measured the reduced;χ2 between the observed emission profile and the azimuthal average of the synthetic millimeter images. Concerning the scattered light, as described in Sect. 3.1, we measured the average brightness profiles within four segments 10° wide in PA, along the major and minor axes (as in Fig. 2). To give the same weight to the ALMA and SPHERE profiles, we considered the average ;χ2 of the scattered light brightness profiles:

(23)

where are the average brightness profiles in the four directions of the SPHERE images.

The best-fit parameters are reported in Table 4 and the comparison between the synthetic and observed profiles is shown in Fig. 3. In Fig. 4, we compare the simulated images to the observations. The best-fit model aptly reproduces both the ALMA and SPHERE images. The simulated emission profiles at both wavelengths (and the four directions for the scattered light images) match the observations. Models and observations are less in accord in the inner disk (< 1 arcsec). As discussed in Sect. 3, the inner disk shows evidence of structures that we have not included in our models, and the coronograph used in SPHERE observations masks the emission in the very inner disk of the scattered light image. Our 1+1D disk model (as discussed in Sect. 2) does not take into account radial transfer of energy and the differences between model and the data in the inner disk structure are within the model uncertainty of using a 1+1D model, rather than a full 2D model. Therefore, we do not include the profiles within 1 arc-second in our ;χ2 calculation, as they do not significantly affect our outer disk model.

thumbnail Fig. 3

Comparison between the observation and best fit radial profiles for the millimeter continuum emission (left) and near-infrared scattered light emission (right). The shaded region shows the observation noise level. In the right plot, the four colors indicate the profiles in the four different directions (same as in Fig. 2).

Table 4

IM Lup best fit parameters of Eqs. (20)(22).

4.1 Dust mass

A key result of our modeling is the constraints on the dust mass of the disk. In the optically thin case, the millimeter flux is proportional to the disk dust mass times the dust opacity. The opacity at millimeter wavelengths is the highest for millimeter-sized grains, but it is only a factor of a few higher than the millimeter opacity of submillimeter-sized grains (e.g., Woitke et al. 2016; Birnstiel et al. 2018, also Fig. 1). If millimeter grains are absent, we can still observe the millimeter emission from submillimeter-sized grains. However, dust mass estimates from the continuum millimeter emission are usually based on the assumption that the emission is coming from millimeter-sized grains. This assumption, however, cannot be justified by looking only at the millimeter emission. Our study takes into account both the emission from small and large grains, without assuming the size of the emitting grains. However, the reader should keep in mind that our grain distribution is constrained by the observations only in the region between 160 au and 340 au, as explained in Sects. 3.1–3.2, and our dust mass estimate may be off because the dust distribution in the inner disk is extrapolated from the best-fit dust distribution in the 160-340 au region.

The gas mass is fixed at 0.2 M, derived from thermochemical modeling of ALMA high resolution CO emission maps (Zhang et al. 2021). While it could be possible to find multiple grain size distributions reproducing the millimeter and near-infrared emission, the corner plot in Fig. 5 shows us that this is not the case, as the posterior probability distribution converge to a unique solution. As a further proof that the presence of millimeter grains is necessary to simultaneously match the millimeter and near-infrared emission, we also run our fitting routine using a model where the grain size is not allowed to be larger than 10−3 cm. In this case, the disk model predicts a high abundance of small grains to match the millimeter emission, while still poorly reproducing the data. The large abundance of small grains causes the disk to be more flared than it appears the observations, suggesting that millimeter grains are required to reproduce the observed geometry.

As a further test, we check how well our model constrains the maximum grain size distribution by changing the a0 parameter from Eq. (11) in our best-fit model. By halving a0, the millimeter image appears to be much less extended than the observations, since now millimeter-sized grains are less extended in the radial direction. The disk appears more flared in the scattered light, causing the disk surface in the forward direction to be more aligned to the observed line of sight. The scattered light polarization, in this direction, is at its minimum and therefore the disk appears much dimmer in this region compared to the observations. If we instead double a0, the millimeter emitting region is now more extended than in the observations, but also dimmer, since we have lost some millimeter-sized grains to produce larger grain sizes. In the scattered light, the disk is now less flared than the best-fit model. This causes the light in the forward scattering direction to be more polarized and here we observe a brightness excess compared to the observations. These effects of the maximum grain size distribution on both the millimeter and scattered-light profiles prove that this distribution is well constrained by our model and that the presence of millimeter-sized grains is necessary to simultaneously explain both the millimeter continuum and the near-infrared scattered light observations.

From the dust-to-gas ratio of our best fit model, we find a dust mass of 1.5 × 10−4M in the fitted region. When we extrapolate our model results to the whole disk structure, we get a total dust mass of 4 × 10−3 M. This dust mass estimate is about half the previous estimate of 0.01 M by Pinte et al. (2008), based on 1.3 mm continuum emission obtained with SMA (Panic et al. 2009). This difference comes from the different assumptions made on the grain size distribution, as they assumed well-mixed dust populations with a fixed maximum grain size. They found that a maximum size of 3 mm and a disk dust mass of 0.01 M reproduces the 1.3 mm emission, but not the silicate emission features from micron-sized grains observed in the mid-IR spectra. They suggest that to account for all observables a spatial dependence of the dust grain size distribution is necessary, with larger grains closer to the disk midplane and small grains on the disk surface. Indeed, our model features a vertical distribution of grain sizes, with a radial dependence of the maximum grain size. This results in a lower abundance of large grains in the outer disk. Since these grains carry most of the dust mass (see also the discussion in Sect. 4.2), our model results in a lower dust mass than the one found in Pinte et al. (2008). Moreover, Pinte et al. (2008) used an older estimate of the disk distance of 190 au, based on HIPPARCOS parallax measurement (Wichmann et al. 1998), which is higher than the value used in this work, namely, 158 au (Gaia Collaboration 2016). In their model, the emission of large grains (which carry most of the dust mass) is then brighter by a factor of about 1.4, leading to an overestimation of the dust mass of the same factor.

To better understand the best-fit dust distribution we look at the vertical distribution of the dust properties, and how they change with the radial location. In Fig. 6, we show the Stokes number, the dust vertical profile, ρd(z)/ρd(0), and the vertical dust-to-gas ratio, ρd(z)/ρg(z). The Stokes number is defined as the grain stopping timescale times the Keplerian frequency: S t = τstop ΩK, and it describes how fast the grains drift towards the inner disk. At smaller radii, we have a high Stokes number close to the midplane. This is due to large grains that are quickly drifting towards the inner disk, while at larger radii the midplane has already been depleted from the larger grains by drift. Above the midplane, we have a low Stokes number at all radii. This is caused to the vertical sedimentation of larger grains, suggesting that the vertical sedimentation happens on a much faster rate than the radial drift.

At different radii, we find significant differences in the dust-to-gas ratio vertical profiles. At larger radii the profiles get flatter, with a overall lower ratio, while at smaller radii the dust-to-gas ratio is strongly peaked at the mid-plane, where it gets as high as 0.3. This is an important result, as at such high dust-to-gas ratios streaming instabilities can start to develop, possibly setting up the right conditions for planetesimal formation (Youdin & Goodman 2005; Johansen et al. 2009; Bai & Stone 2010). Moreover, while not included in our disk model, observations show evidence of structures in the disk midplane inside the 160 au radius (as discussed in Sect. 3.2). If there are indeed active planet formation processes in the inner disk, the growing planets could then imprint non-axisymmetrical structures (e.g., Lovelace et al. 1999; Li et al. 2000). Our model, independently from this observational evidence, suggests that at this radius we find the right condition for the developing of streaming instabilities, possibly triggering planet-formation processes which, in turn, could be the origin of the structures observed in the inner disk.

thumbnail Fig. 4

Comparison between the simulated 1.25 mm continuum and 1.65 µm polarized emission, and the ALMA and SPHERE data. The black masks show the region inside 1 arcsec radius, not included in our fitting, where radial structures are observed. The blackcircle is the 2.5 arcsec radius, the outer edge of the disk.

thumbnail Fig. 5

Posterior distribution of the fitting parameters described in Sect. 2. The parameters follow a Gaussian distribution, showing that the parameters are well constrained and not correlated.

4.2 Distribution of the dust populations

The dust mass distribution of different grain populations is an important topic for the understanding of dust evolution and planet formation. Our best-fit model is in agreement with a vertical and radial stratification of dust grains. We find micron-sized grains on the disk surface (z/r ≈ 0.2) and in the outer disk. They are however missing in the midplane at the smaller radii, where they can efficiently grow to form larger grains, as shown in Fig. 7.

Submillimeter grains (Fig. 8) are more abundant than micrometric grains and are likewise depleted at smaller radii by the grain growth process. They also extend vertically in the disk, though they do not reach the disk surface (z/r ≈ 0.1). As these grains are a product of dust evolution, they are not found in the outer disk, where they are removed by radial drift.

Millimeter-sized grains (Fig. 9) are instead found exclusively in the midplane. In contrast to smaller grains, their density increases at smaller radii, as they are a product of dust evolution, which happens on shorter timescales in this context. They are more affected by radial drift, so they are less extended in radius than smaller grains.

The total dust mass, and the mass fraction in small, intermediate, and large grains, is tightly constrained by our fit. In Fig. 10, we show the mass ratio of each dust population to the total dust mass and the total dust mass of the ten highest likelihood models. The ratios and the total dust mass are in accord within all the models, proving that they are well constrained by observing both large and small grain emission.

In Sect. 4.1, we discuss how the millimeter emission does not necessary trace the emission of millimeter-sized grains, as is usually assumed. It is then worthwhile to check whether, according to our model, this is a good assumption for the IM Lup disk. The best-fit grain size distribution, in Fig. 11, shows that according to our model, the millimeter-sized distribution and the observed millimeter continuum emission overlap, and therefore the millimeter grains dominate the emission.

thumbnail Fig. 6

Stokes number, dust vertical profile, and dust-to-gas ratio at 160, 250, and 340 au radii (from top to bottom).

thumbnail Fig. 7

Density distribution of micrometric grains (agrain < 10−3 cm). These grains are more abundant in the regions where large grains are removed by radial drift (outer disk) and vertical settling (disk surface), and depleted in the midplane at smaller radii by dust growth processes.

thumbnail Fig. 8

Density distribution of submillimeter grains (10−3 cm < agrain < 10−1 cm). These grains have a moderate vertical extension and are depleted at both small and large radii, where they are removed by dust evolution processes.

thumbnail Fig. 9

Density distribution of millimeter grains (agrain>10−1 cm). These grains have a moderate vertical extension and are depleted at both small and large radii, where they are removed by dust evolution processes.

4.3 Grain evolution timescales

Since our dust structure is the outcome of a parametric model (and not of a dust evolution simulation), the grain distribution shown in Fig. 11 may not be physically consistent with the disk gas structure. Moreover, our model suggests a steep maximum grain size distribution, as shown both in Figs. 5 and 11. The parametrization assumed for the maximum grain size distribution, Eq. (11), does not hold for all grain sizes. Once grains reach centimeter sizes, they break each other apart by colliding, limiting the maximum grain size to 1 mm to 1 cm, depending on the turbulence strength and gas density (e.g., Birnstiel et al. 2011; Zsom et al. 2011). This prevents grains from growing to unphysically large sizes in the region inside 1 arcsecond, however, we did not include this region in our fit due to the presence of radial structures, as discussed in Sect. 3.2. However, it is worthwhile to investigate if the resulting grain distribution is physically consistent.

In Fig. 12, we compare the collisional, settling, and drift timescales to the disk age, averaged over the grain size distribution. With the exception of the disk midplane, where grains can efficiently grow, the collisional timescale is much longer than the disk age. Here, grains do not have enough time to grow through collision processes and they cannot produce grains larger than micron sizes. However, our model predicts the presence of submillimeter grains even outside of the midplane, which then must be produced by other processes. A possible explanation is coagulation driven by sedimentation (Zsom et al. 2011). To understand this process we must first discuss the grain settling timescales. Figure 12 shows that small grains have a settling timescale of about the disk age on the disk surface, and did not have enough time to start settling to the midplane. Deeper in the disk, the settling timescale gets shorter and grains start to settle towards the midplane. The difference in the settling velocity of grains drives their coagulation, as explored by several authors (e.g., Dullemond & Dominik 2004, 2005; Schräpler & Henning 2004). Zsom et al. (2011) showed that coagulation driven by sedimentation can be a very efficient process for grain growth, and can produce millimeter-sized grains as early as after 104 orbital timescales (104 yr at 1 au). Since the age estimates of IM Lup range from 0.5 Myr (Andrews et al. 2018) to 1 Myr (Andrews et al. 2018), this provides enough time for the formation of the large grains predicted by our model. The distribution of submillimeter grains shown in Fig. 8 is right inside the region where the settling timescale is shorter than the age of the disk, while the millimeter grain distribution corresponds to the region where the settling timescale is ten times shorter than the age of the disk. This anti-correlation between the grain size and the low level of turbulence suggests that in the disk grain growth may be driven by sedimentation instead of turbulent collisions. For a more thorough justification of our model grain distribution, however, we would need to use 2D dust evolution models, which we leave to a future work.

Finally, the drift happens on a much longer timescale than the disk age everywhere but in the midplane, while only being efficient at radii smaller than about 250 au. It is interesting to point out that these regions match the ones where submillimeter and millimeter grains are found (see Figs. 89). This is also in accord with the discussion in Sec. 4.1, where we show how the Stokes number gets higher at smaller radii.

thumbnail Fig. 10

Mass ratios of large (agrain > 10−2 cm), intermediate (10−3 < agrain < 10−2 cm), and small grains (agrain < 10−3 cm) of the ten best-fit models, from upper to lower panels.

thumbnail Fig. 11

Dust surface density of the best fit model, as a function of grain size and radial position.

thumbnail Fig. 12

Comparison of the collisional, settling, and drift timescales to the disk age, from top to bottom.

thumbnail Fig. 13

Comparison of the best-fit model millimeter emission profile but with low (red) and high turbulence (black) to the observed profile (blue).

4.4 Viscous turbulence parameter

Previous studies have constrained the value of αturb to be about αturb = 10−3, based on the observation of turbulent motions in the disks (Hughes et al. 2011; Flaherty et al. 2017). However, a more precise value is more complex to obtain as the relation between the turbulent motions and αturb depends on the nature of the turbulence (Cuzzi et al. 2001). Our best-fit model is in agreement with these findings, with αturb = 3 × 10−3. This is lower than the previous measurement by Pinte et al. (2018, αturb = 9 × 10−3). In this work, they fit a parametrical disk model to ALMA band 6 dust continuum emission and CO,13 CO, and C18 O line emission. The αturb is then estimated from the vertical dust settling in their best-fit model.

To check how the turbulence parameter affects the grain distribution, we ran two models with the same parameters as our best-fit model, but with a low (αturb = 10−4) and high (αturb = 10−2) level of turbulence. In Fig. 13 we show how the millimeter emission profiles of the new models compare to the profile of ALMA data. The millimeter profiles are unaffected by the change in αturb, and this is a further indication that the millimeter emission is coming from the large grains in the midplane, which is not impacted by turbulence.

The near-infrared emission, on the other hand, changes significantly with the turbulence strength. A highly turbulent viscosity results in more grains in the disk atmosphere, as it prevents the dust from settling into the midplane. When αturb = 10−2, the disk appears much more flared in the near-infrared, whereas when αturb = 10−4, the near-infrared emission originates very close to the disk midplane, as the turbulence is not strong enough to sustain micron-sized grains.

This behaviour is also descrived by Rich et al. (2021). While these authors did not measure the value of αturb (they assumed a nominal value of αturb = 10−3 ), they used an analytical model to manually explore the effect of αturb on the 12CO emission height and small dust grain scattering surface. By decreasing the value of αturb, they found a lower scattering surface for small grains, in accordance with our model prediction.

5 Conclusions

We present in this paper a model for the dust structure in IM Lup reproducing both ALMA midplane millimeter continuum emission and SPHERE near-infrared scattered light emission from the disk surface. While these two wavelength are particularly sensitive to large (ALMA) and small (SPHERE) particles, we used these data to build a model for the full grain size distribution. The geometrical disk parameters are taken from the SPHERE DARTTS-S survey (Avenhaus et al. 2018) and the gas mass from modeling CO isotopologues emission from ALMA data (Zhang et al. 2021). The dust radial and vertical grain size distribution, as well as the dust-to-gas ratio, successfully reproduce ALMA millimeter continuum emission from the disk midplane (Andrews et al. 2018) and SPHERE near-infrared scattered light emission from the disk surface (Avenhaus et al. 2018).

Our main results can be summarized as follows:

  1. The posterior probability distribution of the model parameters shows a single solution for the dust distribution, which can reproduce the observations. We find that ~99% of the dust mass is carried by millimeter or larger grains, and a 10−4 mass fraction of micron-sized grains is sufficient to reproduce the scattered light emission surface observed by SPHERE. Assuming a gas mass of 0.2 M inferred from CO emission (Zhang et al. 2021), we estimated a dust mass of 4 × 10−3 M.

  2. The turbulent parameter αturb determines which particle sizes are sedimented and which are well mixed to the gas. We constrained the turbulence parameter to 3 × 10−3 , which is slightly higher than the usual assumed value of 10−3 . The turbulence strength does not affect the model millimeter emission, dominated by large grains unaffected by turbulence. On the other hand, the disk flaring observed in the near-infrared is strongly affected by the turbulence. At low values of αturb, the disk can appear flat also at small wavelengths, since small grains are free to settle to the midplane.

  3. Large millimeter-sized grains in the midplane are necessary to reproduce the near-infrared observations. It requires a vertical stratification of the dust populations, rather than being well mixed to the gas, with millimeter-sized grains in the disk midplane and micron-sized grains on the disk surface, in accordance with dust coagulation models which include fragmentation (e.g., Dullemond & Dominik 2005; Birnstiel et al. 2010).

  4. The dust distribution shows a strong radial gradient in the particle size, a feature often associated with a pressure trap. It causes the disk to have a steep maximum grain size distribution, as the structures slow down the radial drift of the larger grains, which gather together on the inner side of the ring (e.g., Pinilla et al. 2015, 2021).

  5. We find a vertical dust distribution with most of the mass carried by large grains in the disk midplane. The existence of large particles and the low turbulence in such a young disk point towards sedimentation-driven coagulation.

  6. As a consequence of the strong radial gradient of the particle size distribution, we find millimeter and larger grains and a high dust-to-gas ratio in the disk midplane at smaller radii. It could be linked to the increase in brightness and the spiral structures observed in the inner disk, as in this environment dust drift is reduced and streaming instability can start to develop, possibly triggering planetesimal formation.

With these results, we demonstrate that the dust properties of a protoplanetary disk can be inferred from the properties of large and small grains. This limits the amount of observational diagnosis needed to characterize the dust structure of a protoplanetary disk.

Acknowledgements

R.F. and T.H. acknowledge support from the European Research Council under the Horizon 2020 Framework Program via the ERC Advanced Grant Origins 83 24 28. We thank Myriam Benisty for the helpful discussions and for providing the SPHERE observational data.

Appendix A Fitting the outer disk

In Sect. 3.2, we discuss how the observations justify our assumption that the dust distribution changes the 340 au annulus, where a outer edge in the dust emission was observed (Avenhaus et al. 2018). In our model, the abundance of small grain follows the one of the larger grains. However, if the dust distributions changes at the ring location, and there are no large grains in the outer disk, we cannot derive the abundance of small grains as the ring affects the large and the small grains differently. Another possible explanation for the sharp truncation in the millimeter profile could be an opacity feature. When the maximum grain size, decreasing with radius, drops below the one which emits most efficiently (a = λ/2π according to the Mie scattering theory), we observe a sharp truncation of the emission profile. In our tests, we explored each of these two cases to check if they can reproduce the observations until the outer edge of the scattered light emission, at 400 au, instead of constraining our fit to the region inside the location of the ring, at 340 au.

To test if the observations can be explained by an opacity feature, we follow the same approach discussed in Sect.2, but we add to the fit the region outside the ring. This new model cannot reproduce both the continuum and the scattered light profiles. The model profiles which best fit the continuum emission are shown in Fig.A.1. The grain distribution reproducing both the brightness and the size of the continuum emitting region does not produce a small grain population that would scatter light with enough efficiency to reproduce the scattered light emission.

The model that best matches the scattered light data is shown in Fig. A.2. This small grain distribution is coupled to an abundance of large grains that is too low to reproduce the continuum emission. Since this approach cannot reproduce, at the same time, the small and large grain emission, we conclude that the sharp drop in the continuum emission is not caused by an opacity effect, and a physical truncation of the dust distribution is necessary to explain the data.

To test whether a cutoff in the dust density can reproduce the observations, we modified Eq. (11) and Eq. (14) by adding an exponential cutoff:

(A.1)

(A.2)

where p1, p2, and rCutoff are new fitting parameter, for a total of eight parameters. The physical drop in the continuum emission at a specific wavelength can be explained by either a cutoff in the maximum grain size distribution or in the dust-to-gas ratio. Since we do not have a theoretical argument to confirm the distribution for which the cutoff is more physically justifiable, we allow for the possibility of a cutoff in both distributions, and allowed the fit to find the one that fits the observed small and large grain emission. The profiles of best model produced by this approach is shown in Fig. A.3. This model places a cutoff of both the maximum grain size and the dust to gas ratio at 300 au, the location of the outer edge of the continuum emission. This model well reproduces both the continuum and scattered light data, with the exception of the outer part of the disk. Here, the truncation removes the small grains emitting in the scattered light and in the model, we observe a sharp drop in the emission not present in the observations. While a truncation of the large grain distribution is necessary to reproduce the continuum observations, this would cause a similar truncation in the small grain distribution. This is not, however, what we observe in the data, because the ring in the dust structure affects the distribution of large grains, but not the one of small grains. Since there are no large grains in the outer part of the disk, we cannot use our model to constrain the small grain distribution from the one of the large grains; thus we constrained our fit to the region where the continuum emission indicates the presence of large grains.

thumbnail Fig. A.1

Model without a truncation in the dust distribution best reproducing the continuum data. This model does not include enough small grains to reproduce the scattered light observations.

thumbnail Fig. A.2

Model without a truncation in the dust distribution best reproducing the scattered light data. This model does not have enough millimeter grains to match the observed continuum emission.

thumbnail Fig. A.3

Best-fit model with a truncation of the maximum grain size and dust-to-gas ratio distribution. This model well reproduces both the continuum and scattered light data in the region inside the truncation, but adds a sharp drop in the scattered light profiles not seen in the observations.

References

  1. Alcalá, J.M., Manara, C.F., Natta, A., et al. 2017, A&A, 600, A20 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  2. Alcalá, J.M., Manara, C.F., France, K., et al. 2019, A&A, 629, A108 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  3. Andrews, S.M., Wilner, D.J., Zhu, Z., et al. 2016, ApJ, 820, L40 [NASA ADS] [CrossRef] [Google Scholar]
  4. Andrews, S.M., Huang, J., Pérez, L.M., et al. 2018, ApJ, 869, L41 [NASA ADS] [CrossRef] [Google Scholar]
  5. Avenhaus, H., Quanz, S.P., Garufi, A., et al. 2018, ApJ, 863, 44 [NASA ADS] [CrossRef] [Google Scholar]
  6. Bai, X.-N., & Stone, J.M. 2010, ApJ, 722, 1437 [NASA ADS] [CrossRef] [Google Scholar]
  7. Bell, K.R., Cassen, P.M., Klahr, H.H., & Henning, T. 1997, ApJ, 486, 372 [NASA ADS] [CrossRef] [Google Scholar]
  8. Birnstiel, T., Dullemond, C.P., & Brauer, F. 2010, A&A, 513, A79 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  9. Birnstiel, T., Ormel, C.W., & Dullemond, C.P. 2011, A&A, 525, A11 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  10. Birnstiel, T., Klahr, H., & Ercolano, B. 2012, A&A, 539, A148 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  11. Birnstiel, T., Fang, M., & Johansen, A. 2016, Space Sci. Rev., 205, 41 [Google Scholar]
  12. Birnstiel, T., Dullemond, C.P., Zhu, Z., et al. 2018, ApJ, 869, L45 [CrossRef] [Google Scholar]
  13. Bohren, C.F., & Huffman, D.R. 1998, Absorption and Scattering of Light by Small Particles, Wiley Science Paperback Series (New Jersey: John Wiley & Sons) [Google Scholar]
  14. Brauer, F., Dullemond, C.P., & Henning, T. 2008, A&A, 480, 859 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  15. Chiang, E.I., & Goldreich, P. 1997, ApJ, 490, 368 [Google Scholar]
  16. Cleeves, L.I., Öberg, K.I., Wilner, D.J., et al. 2016, ApJ, 832, 110 [NASA ADS] [CrossRef] [Google Scholar]
  17. Cuzzi, J.N., Hogan, R.C., Paque, J.M., & Dobrovolskis, A.R. 2001, ApJ, 546, 496 [NASA ADS] [CrossRef] [Google Scholar]
  18. Draine, B.T. 2003, ARA&A, 41, 241 [NASA ADS] [CrossRef] [Google Scholar]
  19. Dubrulle, B., Morfill, G., & Sterzik, M. 1995, Icarus, 114, 237 [Google Scholar]
  20. Dullemond, C.P., & Dominik, C. 2004, A&A, 421, 1075 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  21. Dullemond, C.P., & Dominik, C. 2005, A&A, 434, 971 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  22. Dullemond, C.P., Dominik, C., & Natta, A. 2001, ApJ, 560, 957 [NASA ADS] [CrossRef] [Google Scholar]
  23. Dullemond, C.P., van Zadelhoff, G.J., & Natta, A. 2002, A&A, 389, 464 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  24. Dullemond, C.P., Juhasz, A., Pohl, A., et al. 2012, Astrophysics Source Code Library [record ascl:1202.015] [Google Scholar]
  25. Dzyurkevich, N., Turner, N.J., Henning, T., & Kley, W. 2013, ApJ, 765, 114 [NASA ADS] [CrossRef] [Google Scholar]
  26. Flaherty, K.M., Hughes, A.M., Rose, S.C., et al. 2017, ApJ, 843, 150 [NASA ADS] [CrossRef] [Google Scholar]
  27. Franceschi, R., Birnstiel, T., Henning, T., et al. 2022, A&A, 657, A74 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  28. Fromang, S., & Nelson, R.P. 2009, A&A, 496, 597 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  29. Gaia Collaboration (Prusti, T., et al.) 2016, A&A, 595, A1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  30. Grady, C.A., Devine, D., Woodgate, B., et al. 2000, ApJ, 544, 895 [NASA ADS] [CrossRef] [Google Scholar]
  31. Grady, C.A., Schneider, G., Sitko, M.L., et al. 2009, ApJ, 699, 1822 [NASA ADS] [CrossRef] [Google Scholar]
  32. Hashimoto, J., Tamura, M., Muto, T., et al. 2011, ApJ, 729, L17 [Google Scholar]
  33. Henning, T., & Stognienko, R. 1996, A&A, 311, 291 [NASA ADS] [Google Scholar]
  34. Huang, J., Öberg, K.I., Qi, C., et al. 2017, ApJ, 835, 231 [NASA ADS] [CrossRef] [Google Scholar]
  35. Huang, J., Andrews, S.M., Dullemond, C.P., et al. 2018a, ApJ, 869, L42 [NASA ADS] [CrossRef] [Google Scholar]
  36. Huang, J., Andrews, S.M., Pérez, L.M., et al. 2018b, ApJ, 869, L43 [NASA ADS] [CrossRef] [Google Scholar]
  37. Hughes, A.M., Wilner, D.J., Andrews, S.M., Qi, C., & Hogerheijde, M.R. 2011, ApJ, 727, 85 [NASA ADS] [CrossRef] [Google Scholar]
  38. Johansen, A., & Klahr, H. 2005, ApJ, 634, 1353 [NASA ADS] [CrossRef] [Google Scholar]
  39. Johansen, A., Klahr, H., & Mee, A.J. 2006, MNRAS, 370, L71 [NASA ADS] [CrossRef] [Google Scholar]
  40. Johansen, A., Youdin, A., & Mac Low, M.-M. 2009, ApJ, 704, L75 [NASA ADS] [CrossRef] [Google Scholar]
  41. Kenyon, S.J., & Hartmann, L. 1987, ApJ, 323, 714 [NASA ADS] [CrossRef] [Google Scholar]
  42. Krijt, S., Ormel, C.W., Dominik, C., & Tielens, A.G.G.M. 2015, A&A, 574, A83 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  43. Law, C.J., Crystian, S., Teague, R., et al. 2022, ApJ, 932, 114 [NASA ADS] [CrossRef] [Google Scholar]
  44. Li, H., Finn, J.M., Lovelace, R.V.E., & Colgate, S.A. 2000, ApJ, 533, 1023 [NASA ADS] [CrossRef] [Google Scholar]
  45. Lovelace, R.V.E., Li, H., Colgate, S.A., & Nelson, A.F. 1999, ApJ, 513, 805 [Google Scholar]
  46. Lynden-Bell, D., & Pringle, J.E. 1974, MNRAS, 168, 603 [Google Scholar]
  47. Mathis, J.S., Rumpl, W., & Nordsieck, K.H. 1977, ApJ, 217, 425 [NASA ADS] [CrossRef] [Google Scholar]
  48. Monnier, J.D., Harries, T.J., Aarnio, A., et al. 2017, ApJ, 838, 20 [NASA ADS] [CrossRef] [Google Scholar]
  49. Mulders, G.D., Min, M., Dominik, C., Debes, J.H., & Schneider, G. 2013, A&A, 549, A112 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  50. Panic, O., Hogerheijde, M.R., Wilner, D., & Qi, C. 2009, A&A, 501, 269 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  51. Pinilla, P., van der Marel, N., Pérez, L.M., et al. 2015, A&A, 584, A16 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  52. Pinilla, P., Lenz, C.T., & Stammler, S.M. 2021, A&A, 645, A70 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  53. Pinte, C., Padgett, D.L., Ménard, F., et al. 2008, A&A, 489, 633 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  54. Pinte, C., Dent, W.R.F., Ménard, F., et al. 2016, ApJ, 816, 25 [Google Scholar]
  55. Pinte, C., Ménard, F., Duchêne, G., et al. 2018, A&A, 609, A47 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  56. Powell, D., Murray-Clay, R., & Schlichting, H.E. 2017, ApJ, 840, 93 [NASA ADS] [CrossRef] [Google Scholar]
  57. Powell, D., Murray-Clay, R., Pérez, L.M., Schlichting, H.E., & Rosenthal, M. 2019, ApJ, 878, 116 [NASA ADS] [CrossRef] [Google Scholar]
  58. Ricci, L., Maddison, S.T., Wilner, D., et al. 2015, ApJ, 813, 138 [NASA ADS] [CrossRef] [Google Scholar]
  59. Rich, E.A., Teague, R., Monnier, J.D., et al. 2021, ApJ, 913, 138 [NASA ADS] [CrossRef] [Google Scholar]
  60. Rosotti, G.P., Booth, R.A., Tazzari, M., et al. 2019, MNRAS, 486, L63 [CrossRef] [Google Scholar]
  61. Schräpler, R., & Henning, T. 2004, ApJ, 614, 960 [CrossRef] [Google Scholar]
  62. Shakura, N.I., & Sunyaev, R.A. 1973, A&A, 500, 33 [NASA ADS] [Google Scholar]
  63. Stapelfeldt, K.R., Duchêne, G., Perrin, M., et al. 2014, in Exploring the Formation and Evolution of Planetary Systems, eds. M. Booth, B.C. Matthews, & J.R. Graham, 299, 99 [NASA ADS] [Google Scholar]
  64. Tazzari, M., Clarke, C.J., Testi, L., et al. 2021, MNRAS, 506, 2804 [NASA ADS] [CrossRef] [Google Scholar]
  65. Teague, R. 2019, J. Open Source Softw., 4, 1632 [NASA ADS] [CrossRef] [Google Scholar]
  66. Villenave, M., Ménard, F., Dent, W.R.F., et al. 2020, A&A, 642, A164 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  67. Warren, S.G., & Brandt, R.E. 2008, J. Geophys. Res. (Atmos.), 113, D14220 [NASA ADS] [CrossRef] [Google Scholar]
  68. Wichmann, R., Bastian, U., Krautter, J., Jankovics, I., & Rucinski, S.M. 1998, MNRAS, 301, L39 [NASA ADS] [CrossRef] [Google Scholar]
  69. Woitke, P., Min, M., Pinte, C., et al. 2016, A&A, 586, A103 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  70. Youdin, A.N., & Goodman, J. 2005, ApJ, 620, 459 [NASA ADS] [CrossRef] [Google Scholar]
  71. Zhang, K., Booth, A.S., Law, C.J., et al. 2021, ApJS, 257, 5 [NASA ADS] [CrossRef] [Google Scholar]
  72. Zsom, A., Ormel, C.W., Dullemond, C.P., & Henning, T. 2011, A&A, 534, A73 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

4

The source code can be found at: https://github.com/rfranceschi/IMLup_SPHERE_fit

All Tables

Table 1

DSHARP dust composition from Birnstiel et al. (2018)

Table 2

IM Lupi stellar and disk properties.

Table 3

IM Lup small grain ring parameters, as in Avenhaus et al. (2018).

Table 4

IM Lup best fit parameters of Eqs. (20)(22).

All Figures

thumbnail Fig. 1

Medium optical constant for the dust composition in Table 1 (Birnstiel et al. 2018).

In the text
thumbnail Fig. 2

IM Lup polarized emission at 1.65 µm wavelength observed by SPHERE (Avenhaus et al. 2018). The cones show the regions over which we extract the radial profiles we use to compare our models.

In the text
thumbnail Fig. 3

Comparison between the observation and best fit radial profiles for the millimeter continuum emission (left) and near-infrared scattered light emission (right). The shaded region shows the observation noise level. In the right plot, the four colors indicate the profiles in the four different directions (same as in Fig. 2).

In the text
thumbnail Fig. 4

Comparison between the simulated 1.25 mm continuum and 1.65 µm polarized emission, and the ALMA and SPHERE data. The black masks show the region inside 1 arcsec radius, not included in our fitting, where radial structures are observed. The blackcircle is the 2.5 arcsec radius, the outer edge of the disk.

In the text
thumbnail Fig. 5

Posterior distribution of the fitting parameters described in Sect. 2. The parameters follow a Gaussian distribution, showing that the parameters are well constrained and not correlated.

In the text
thumbnail Fig. 6

Stokes number, dust vertical profile, and dust-to-gas ratio at 160, 250, and 340 au radii (from top to bottom).

In the text
thumbnail Fig. 7

Density distribution of micrometric grains (agrain < 10−3 cm). These grains are more abundant in the regions where large grains are removed by radial drift (outer disk) and vertical settling (disk surface), and depleted in the midplane at smaller radii by dust growth processes.

In the text
thumbnail Fig. 8

Density distribution of submillimeter grains (10−3 cm < agrain < 10−1 cm). These grains have a moderate vertical extension and are depleted at both small and large radii, where they are removed by dust evolution processes.

In the text
thumbnail Fig. 9

Density distribution of millimeter grains (agrain>10−1 cm). These grains have a moderate vertical extension and are depleted at both small and large radii, where they are removed by dust evolution processes.

In the text
thumbnail Fig. 10

Mass ratios of large (agrain > 10−2 cm), intermediate (10−3 < agrain < 10−2 cm), and small grains (agrain < 10−3 cm) of the ten best-fit models, from upper to lower panels.

In the text
thumbnail Fig. 11

Dust surface density of the best fit model, as a function of grain size and radial position.

In the text
thumbnail Fig. 12

Comparison of the collisional, settling, and drift timescales to the disk age, from top to bottom.

In the text
thumbnail Fig. 13

Comparison of the best-fit model millimeter emission profile but with low (red) and high turbulence (black) to the observed profile (blue).

In the text
thumbnail Fig. A.1

Model without a truncation in the dust distribution best reproducing the continuum data. This model does not include enough small grains to reproduce the scattered light observations.

In the text
thumbnail Fig. A.2

Model without a truncation in the dust distribution best reproducing the scattered light data. This model does not have enough millimeter grains to match the observed continuum emission.

In the text
thumbnail Fig. A.3

Best-fit model with a truncation of the maximum grain size and dust-to-gas ratio distribution. This model well reproduces both the continuum and scattered light data in the region inside the truncation, but adds a sharp drop in the scattered light profiles not seen in the observations.

In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.