Modeling the arc and ring structures in the HD 143006 disk

To date, there are more than 4000 exo-planets discovered using various techniques, such as radial velocity, transit surveys, and a recently developed method based on disk kinematics [1, 2, 3, 4]. How exactly these planets form remains one of the major questions in modern astrophysics. One direct way to solve the problem is to image the (proto)planets that are still in the making within their natal disks, as they can provide clear insights on where, when, and how planets are born. Unfortunately, such observations are very challenging because the embedded planets are much fainter in the optical/infrared regime compared to the host central star, and even to the emission knots, if present, from the disk material in the vicinity of the planet. Consequently, only a few young stars have been claimed to harbour (proto)planets so far [5, 6], and most of them still require further confirmation [7, 8, 9].

Nowadays, infrared and millimeter images of disks, taken by the Very Large Telescope Spectro-Polarimetric High-contrast Exoplanet REsearch (VLT/SPHERE) and Atacama Large Millimeter/submillimeter Array (ALMA), are able to achieve angular resolutions on the order of ${\sim}\,0.05^{\prime\prime}$ [10], corresponding to spatial scales of ${\sim}\,7\,\rm{AU}$ at the distance of nearby star formation regions (${\sim}\,140\,\rm{pc}$, for instance Taurus). The accumulating SPHERE and ALMA surveys have revealed a variety of disk substructures, such as rings [11, 12, 13], gaps [14, 15, 16], spirals [17, 18, 19, 20], and arcs [21, 22, 23, 24]. The interaction between the disk and embedded planets is one of the proposed explanations for the formation of these features [25]. The properties of substructures together with the disk morphology can be used to constrain the parameters of the embedded planets [26, 27]. Alternative mechanisms, which can produce substructures without planets, include zonal flows generated by the magnetorotational instability turbulence [28, 29], fast pebble growth near condensation fronts of volatiles [30], operation of a secular gravitational instability [31], and the pile-up of drifting pebbles [32, 33]. In these scenarios, disk substructures may play a role of dust trapping, which on the one hand prevents dust particles from rapid inwards drift [34], and on the other hand promotes grain growth and concentration to form planetesimals [35].

It is of particular interest to investigate the properties and stability of disk substructures, because such efforts help to understand how the disk evolves and for what timescales the substructures will impact the process of planet formation. Arcs, or their analogs, e.g., crescents and horseshoes, are often interpreted as vortices arising from the Rossby wave instability (RWI) at the edge of a gap created by a young planet [21, 36]. Some numerical simulations show that, in weakly turbulent disks, vortex-driven arcs and rings might survive for more than ${\sim}10^{3}$ planetary orbits, and even for disk evolutionary timescales (${\sim}\,\rm{Myr}$) [37, 38, 33]. However, observational constraints on the strength of turbulence exist for only a few disks [39, 40].

Rings and arcs in some cases appear extraordinarily far from the central star [22, 41]. In such cold and dense regions, gravitational instability may start to operate, which in turn complicates the determination for the fate of substructures. HD 143006 stands out from this disk category because of its complex geometry and very bright arc outside of the millimeter dust disk [41, 24]. This is a G7 young star located in the Upper Scorpius association at a distance of 166 pc [42, 43]. Its effective temperature and luminosity are ${\sim}\,5623\,\rm{K}$ and ${\sim}\,3.8\,L_{\odot}$, respectively [44, 13]. Comparing the location in the Hertzsprung-Russell diagram with models of pre-main-sequence evolutionary tracks returns a stellar mass of ${\sim}\,1.77\,M_{\odot}$, and an age of $4\,{-}\,12\,\rm{Myr}$ [13]. In this work, we focus on this intriguing object, and conduct three-dimensional radiative transfer modeling of the spectral energy distribution (SED), and high resolution infrared and millimeter images, with the goal of constraining the dust density and temperature distributions. On the basis of the modeling result, we investigate the properties and gravitational stability of the substructures.

The SED has been widely used to constrain the grain size distribution and overall structure of the disk [45]. We compiled the photometry of HD 143006 from the GAIA [43], 2MASS [46], WISE [47], and IRAS catalogs [48]. In addition, individual flux measurements from optical to millimeter regimes reported in the literature are also collected [49, 50, 51, 52, 44, 13]. As shown in fig. 1, the SED shows a clear dip at $\lambda\,{\sim}\,10\,\mu{\rm m}$, suggesting that substructures (e.g., gaps) exist in the distribution of micron-sized dust grains. A linear fit to the data points longer than $0.8\,\rm{mm}$ results in a millimeter spectral index of $\alpha\,{\sim}\,2.9$, which yields a dust opacity slope of 0.9 under an assumption of optically thin emission in the Rayleigh-Jeans limit [53, 54]. This value is lower than that of the interstellar medium dust (i.e., ${\sim}\,1.7$), indicating that dust grains in the disk have already grown to millimeter sizes [55, 56].

Using the VLT/SPHERE instrument, Benisty et al. [41] observed HD 143006 at an angular resolution of ${\sim}\,0.037^{\prime\prime}$. fig. 2 shows the J-band ($\lambda_{0}\,{\sim}\,1.25\,\mu{\rm m}$) polarized intensity, which traces the scattered light from micron-sized dust grains populating the disk’s surface layer [41]. The disk reveals a strong east/west brightness asymmetries, and two rings. The outer Ring #1 extends from $R\,{\sim}\,0.3^{\prime\prime}$ to $0.5^{\prime\prime}$, and appears obscure in the western side. The inner Ring #2 has a central position of ${\sim}\,0.15^{\prime\prime}$. Along the azimuthal direction, its brightness is not continuously distributed, but shows two shadows at $\rm{PA}\,{\sim}\,190^{\circ}$ and $340^{\circ}$. These prominent features are indicative of a complex disk configuration with a misaligned inner disk [57]. The clear signature of flux depletion between rings (i.e., in the gap) and also in the central cavity is consistent with the dip shown in the mid-infrared SED.

HD 143006 was selected as one of the 20 disks observed in the Disk Substructures at High Angular Resolution Project (DSHARP) [13]. With an angular resolution of ${\sim}\,0.05^{\prime\prime}$, the 1.25 mm continuum image obtained with ALMA is resolved into three rings (B8, B40, and B64) at radial locations of $R\,{\sim}\,0.05^{\prime\prime}$, $0.24^{\prime\prime}$ and $0.39^{\prime\prime}$ respectively, see the left panel of fig. 3. Note that the nomenclature of these rings is directly taken from Pérez et al. [24]. Unlike the J-band scattered light, millimeter dust rings are broadly symmetric, except that the innermost B8 ring features two peaks along a PA of ${\sim}\,165^{\circ}$. Ring #2 identified in the SPHERE image lies closer to the central star than the B40 ring, whereas the position of the SPHERE Ring #1 matches well with the B64 ring. The most distinct feature in the ALMA image is the bright arc in the southeast, located just exterior to the B64 ring, and extends from ${\sim}\,118^{\circ}$ to ${\sim}\,165^{\circ}$ in the azimuthal direction. This arc has a counterpart in the SPHERE image, but it shows a broader radial and azimuthal extent.

The origin of the rich substructures of HD 143006 is yet not fully understood. Hydrodynamic simulations demonstrate that an inclined equal mass binary causes disc breaking into two distinct annuli, with the resulting inner disk misaligned with the outer region [58, 41]. This scenario can produce azimuthal asymmetries shown in the scattered light. Ballabio et al. [59] found that besides an inclined binary, a planetary companion located at $R\,{=}\,32\,\rm{AU}$ is required to simultaneously explain the morphological features in the scattered light and the millimeter dust rings. We do not attempt to re-investigate the same story. Instead, we are interested to constrain the temperature and density distribution of the disk with detailed radiative transfer modeling, and to infer the properties of the substructures.

The fitting results are presented in figs. 1, 3 and 7, whereas table 1 summarizes the corresponding parameter set. Our model well matches the ALMA data and the morphology of the infrared polarized intensity map, see for instance how high the quality of fit is for the arc and rings in the right panel of fig. 3.

The model SED is broadly consistent with the observation, except for the data points at ${\sim}\,3\,\mu{\rm m}$ and ${\sim}\,40\,\mu{\rm m}$. The geometry of the inner disk for HD 143006 might be more complicated than what we assumed (see fig. 4), for instance there might be an puffed-up inner rim [74]. Moreover, a variability in accretion or an inclined companion can lead to changes in the structure of the inner disk [75]. These factors can explain the discrepancies in the near-infrared flux. For simplicity, we assume that the flaring index $\beta$ and scale height $H_{100}$ are the same in the inner and outer disk. The mild overestimation of far-infrared flux can be reconciled by reducing $H_{100}$ (or $\beta$, or both) for the outer disk while maintaining the vertical structure of the inner disk. Nevertheless, such finely refinements are not expected to have a significant impact on the overall properties (i.e., temperature, surface density, and Toomre parameter) of the disk.

The left panel of fig. 5 shows the surface densities projected onto the plane of sky. A projection is performed to have a one-to-one comparison with the observation. The reconstructed surface densities show obvious fluctuations along both the radial and azimuthal directions. We checked the quality of the ALMA data, and found that in most disk regions (excepting for the gap around $R\,{\sim}\,22\,\rm{AU}$), the signal-to-noise ratio is higher than 3. In the arc and rings where we are most interested, the ratio is higher than ${\sim}\,10$. These mean that the derived surface densities for the arc and rings are reliable since our iterative fitting is driven by the data. Nevertheless, surface density fluctuations in the gap around $R\,{\sim}\,22\,\rm{AU}$ may reflect the noise of the ALMA observation. The right panel presents the optical depth at 1.25 mm and the mass-averaged temperature, derived from an azimuthal average between ${\rm PA}\,{=}\,118^{\circ}$ and $165^{\circ}$. As can be seen, all the three ALMA rings (B8, B40, B64) are marginally optically thick: $\tau_{1.25\,\rm{mm}}\,{\sim}\,1$. This is different from the result presented by Huang et al. [60], who obtained an optical depth of $\tau_{1.25\,\rm{mm}}\,{\lesssim}\,0.2$. When calculating $\tau_{1.25\,\rm{mm}}$, they assumed a temperature profile being the simplified expression for a passively heated, flared disk: $T\,{\propto}\,R^{-0.5}$ [76, 74]. It is clear that such an approximation is not sufficient for HD 143006, see the blue solid line in the figure. Instead, due to the misalignment, the temperature exhibits a discontinuity at the boundary ($R_{\rm b}\,{=}\,18\,\rm{AU}$) between the inner and outer disk. Along the radial direction towards the arc (${\rm PA}\,{=}\,141.5^{\circ}$), stellar photons can directly impinge upon the inner rim of the outer disk, resulting in a sharp increase in temperature. Beyond that radius, the temperature decreases by roughly following a power law with $T\,{\sim}\,20\,\rm{K}$ at the arc location. Moreover, in the analysis of Huang et al. [60], dust scattering is not taken into account. The difference in optical depth between our work and theirs can be understood in terms of these two reasons. The bright arc is optically thick with $\tau_{1.25\,\rm{mm}}\,{\sim}\,2.7$.

The misalignment leads to a highly asymmetric temperature distribution, see the middle panel of fig. 5. The eastern/western sides of the outer disk’s rim (located at $R_{\rm b}\,{=}\,18\,\rm{AU}$, the green dashed circle in the plot) are directly exposed to the central star. The typical temperature there reaches $110\,\rm{K}$. Conversely, the northern/southern parts of the rim are cool with $T_{\rm dust}\,{\sim}\,60\,\rm{K}$. table 2 gives the mean temperature of each substructure. We defined the arc region as the area enclosed with the white dashed lines in the left panel of fig. 3. The radial range of the arc overlaps with that of the B64 ring, because we considered that both features are tightly connected. Adopting the condensation temperature reported in the literature [77, 78, 30], the B8, B40, and B64 rings are associated with the ice lines of ${\rm NH_{3}}$, ${\rm H_{2}S}$, and ${\rm CH_{4}}$, respectively. By integrating the surface density, we obtained the dust mass in the substructures. The total dust mass of the best model is $M_{\rm dust}\,{=}\,9.7\times 10^{-5}\,M_{\odot}$ (or $32\,M_{\oplus}$). The arc contains about 20% of the dust mass. This number should be regarded as a lower limit due to the optical thickness $\tau_{1.25\,\rm{mm}}\,{\sim}\,2.7$. The mass fraction is lower than that of the arc in HD 135344B (i.e., ${\sim}\,50\%$), which shares a similar disk morphology [23]. The ring masses range from 0.6 to $16\,M_{\oplus}$. They are systematically lower than those derived in the younger (${\sim}\,1\,\rm{Myr}$) HL Tau disk. In HL Tau, the dust masses in the rings are from 16 to ${\sim}\,120\,M_{\oplus}$ [79]. If the rings are created by planets, this may suggest that much of the dusty material in the HD 143006 disk have already been converted to the mass of planets (or planetary cores).

As can be seen from table 2, the rings, particularly the B64 ring, and the arc occupy a large amount of dust mass. In addition, the arc is located so distant from the central star that the temperature is as low as $23\,\rm{K}$. Gravitational instability may be in action at these cold and dense regions. We checked this possibility using the Toomre criterion [80]. The Toomre $Q$ parameter is defined as

where $c_{s}$ is the local sound speed, $\Omega_{K}(R)$ is the angular frequency, and $\Sigma_{\rm gas}$ stands for the gas surface density. In terms of $Q$, a disk/substructure is unstable to its own self-gravity if $Q\,{<}\,1$, and stable if $Q\,{>}\,1$. As there is no constraints on $\Sigma_{\rm gas}$, for simplicity, we scaled up the constructed $\Sigma_{\rm dust}$ by a constant gas-to-dust mass ratio of 30 that is found to be representative via millimeter dust and gas observations of a sample of protoplanetary disks [81, 82, 83]. The two-dimensional distribution of $Q$ is shown in the left panel of fig. 6. The value is below unity at the center of the arc, as indicated by the green contour. The right panel of fig. 6 displays the azimuthally averaged $Q$ as a function of radius. Two ranges of azimuthal angles were considered in order to highlight the difference at various locations. While most regions including the three rings are gravitationally stable, $Q$ approaches ${\sim}1.3$ at the radial position of the arc. Pohl et al. [84] showed that an embedded massive planet (planet-to-star mass ratio of 0.01) can induce gravitational collapses in massive disks when $Q$ is around unity. Nevertheless, analyzing infrared interferometric data excludes the presence of a companion with mass ratios of ${\gtrsim}\,0.1$ [41]. Future searches for less massive (proto)planets with high-contrast imaging are desired to test the hypothesis.

The gas-to-dust mass ratio reflects how quickly the dust and gas components dissipate, in particular relative to each other. The inherited interstellar medium gas-to-dust mass ratio is 100. ALMA observations of dust and gas lines reveal that the ratio can drop to 10 in some disks [82, 85]. Arcs are often interpreted as a result of dust trapping into vortices that preferentially concentrate dust grains. Numerical simulations of vortex formation and dust trapping have shown that the gas-to-dust mass ratio in the center of vortices can decrease down to 10 and even ${\sim}\,1$ [86, 87]. In such a situation, $Q$ will be above unity, and gravitational instability is suppressed accordingly. Future high resolution and high sensitivity observations of gas lines are required to constrain the gas surface density, and therefore the gas-to-dust mass ratio, particularly in the local arc region. Furthermore, hydrodynamic simulations dedicated on HD 143006, with the inclusion of the feedback of dust onto gas [88], will help to understand the fate of the arc.

Vortices are regions acting as local pressure maxima in the disk, and they can form at the edge of a planet-induced gap due to the excitation of RWI. An important prediction of dust trapping toward vortices is that larger dust grains are expected to be concentrated into a narrower region [89, 90], which has been observationally identified in several disks, e.g., Oph IRS 48 [91], HD 142527 [92], MWC 758 [93], and HD135344B [23]. Investigating the spatial distribution of dust grains with different sizes requires high resolution observations at multi-wavelengths. Due to the lack of such datasets, we distributed the SGP and LGP into the same azimuthal extent in the radiative transfer models. Future high resolution observations at different wavelengths, and follow-up search for (proto)planets are indispensable to clarify the nature of the arc in the HD 143006 disk.

We performed three dimensional radiative transfer modeling of the HD 143006 disk that shows three dust rings and a bright arc in the outermost region on the DSHARP ALMA map. The goal of our work is to constrain the density and temperature distributions, and to investigate the possibility of triggering a gravitational instability in the substructures. In order to account for the ALMA and SPHERE observations, our model features a complex geometry with a misalignment between the inner and outer disk.

All of the three rings (B8, B40, B64) are marginally optically thick ($\tau_{1.25\,\rm{mm}}\,{\sim}\,1$) at the observed ALMA wavelength, while the optical thickness of the arc is $\tau_{1.25\,\rm{mm}}\,{\sim}\,2.7$. Our self-consistent modeling gives higher optical depth than those derived from a simple analytic method in the literature. We found that the temperature as a function of radius resembles two connecting power laws with a sharp jump at the boundary of the misalignment. The temperature in the rings are from 74 to 29 K, and it drops down to 23 K in the arc. The total dust mass in the disk is about $32\,M_{\oplus}$, with the arc sharing about 20% of the amount. The dust masses in the rings range from $0.6\,M_{\oplus}$ to $16\,M_{\oplus}$, which are systematically lower than those found in the younger HL Tau disk.

While all the rings are gravitationally stable, the Toomre parameter $Q$ in the arc is around unity when assuming a gas-to-dust mass ratio of 30 that is constant in the entire disk. The arc is close to the threshold of gravitational collapse in this case. If gravitational instability is in action and the cooling is efficient enough, giant planets may quickly form [94]. However, if the arc indeed traces a vortex, millimeter dust grains preferentially concentrate in the arc, which will decrease the local gas-to-dust mass ratio. Given that the gas-to-dust mass ratio is an unknown assumption, adopting a lower value would increase our estimate of $Q$. In this situation, gravitational instability will not take place. Future gas and dust observations at high resolutions are required to better constrain the gas-to-dust mass ratio, and to clarify the origin of the arc in the HD 143006 disk.

We thank the anonymous referees for their constructive comments that highly improved the manuscript. We sincerely thank Myriam Benisty for sharing the reduced SPHERE data of HD 143006. Y.L. acknowledges the financial support by the Natural Science Foundation of China (Grant No. 11973090). M.F. acknowledge funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation program (grant agreement No. 757957). We acknowledge the science research grants from the China Manned Space Project with NO. CMS-CSST-2021-B06. Herschel is an ESA space observatory with science instruments provided by European-led Principal Investigator consortia and with important participation from NASA. 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. This work has made use of data from the European Space Agency (ESA) mission Gaia (https://www.cosmos.esa.int/gaia), processed by the Gaia Data Processing and Analysis Consortium (DPAC, https://www.cosmos.esa.int/web/gaia/dpac/consortium). Funding for the DPAC has been provided by national institutions, in particular the institutions participating in the Gaia Multilateral Agreement.

The authors declare that they have no conflict of interest.