Dust unveils the formation of a mini-Neptune planet in a protoplanetary ring

TL;DR

This study uses ALMA observations and hydrodynamical modeling to demonstrate that a low-mass migrating planet can create the complex ring structures observed in the protoplanetary disk of HD169142, indicating mini-Neptune formation.

Contribution

The paper provides evidence that a single low-mass planet can explain complex ring structures in protoplanetary disks, linking observations to planet formation theories.

Findings

ALMA images reveal three narrow, asymmetric rings in HD169142.

A migrating 10 Earth-mass planet explains the ring architecture.

Outer rings suggest formation of mini-Neptune-sized planets.

Abstract

Rings and radial gaps are ubiquitous in protoplanetary disks, yet their possible connection to planet formation is currently subject to intense debates. In principle, giant planet formation leads to wide gaps which separate the gas and dust mass reservoir in the outer disk, while lower mass planets lead to shallow gaps which are manifested mainly on the dust component. We used the Atacama Large Millimeter/submillimeter Array (ALMA) to observe the star HD169142, host to a prominent disk with deep wide gaps that sever the disk into inner and outer regions. The new ALMA high resolution images allow for the outer ring to be resolved as three narrow rings. The HD169142 disk thus hosts both the wide gaps trait of transition disks and a narrow ring system similar to those observed in the TW Hya and HL Tau systems. The mass reservoir beyond a deep gap can thus host ring systems. The observed rings are narrow in radial extent (width/radius of 1.5/57.3, 1.8/64.2 and 3.4/76.0, in au) and have asymmetric mutual separations: the first and middle ring are separated by 7 au while the middle and outermost ring are distanced by ~12 au. Using hydrodynamical modeling we found that a simple explanation, involving a single migrating low mass planet (10 M$_\oplus$), entirely accounts for such an apparently complex phenomenon. Inward migration of the planet naturally explains the ring's asymmetric mutual separation. The isolation of HD169142's outer rings thus allows a proof of concept to interpret the detailed architecture of the outer region of protoplanetary disks with low mass planet formation of mini-Neptune's size, i.e. as in the protosolar nebula.

Equations

$$\frac{r_{\mathrm{OG}}-r_{\mathrm{IG}}}{r_{\mathrm{p}}}\approx 2.9{\left(\frac{\gamma +1}{12/5}\frac{M_{\mathrm{p}}}{M_{\mathrm{th}}}\right)}^{-2/5}\left(\frac{h}{r}\right),$$

$$h_{i,\mathrm{d}}=h\times \sqrt{\frac{D_{\mathrm{z}}}{D_{\mathrm{z}}+{\mathrm{St}}_i}},$$

Full text

protoplanetary disks — planets and satellites: formation — planet-disk interactions — submillimeter: planetary systems

††facilities: ALMA Observatory††software: GPUVMEM (Cárcamo et al., 2018), CASA (McMullin et al., 2007), FARGO-ADSG (Baruteau & Masset, 2008a, b), RADMC-3D (Dullemond et al., 2012).

1 Introduction

Planet formation theories describe two main pathways for the formation of protoplanet embryos: gravitational instability followed by fragmentation, which forms gas giant planets, or (sub-)stellar mass companions (e.g., Boss, 1997; Kratter et al., 2010); and bottom-up growth by core accretion, which forms terrestrial and icy planets, or the cores of giant planets (e.g., Pollack et al., 1996). However, the timescale of core accretion increases dramatically with stellocentric distance due to the longer dynamical timescale and low surface densities at large radii (e.g., Goldreich et al., 2004), which poses a problem for the formation of planets before the protoplanetary disk dissipates. Mechanisms such as pebble accretion (Ormel & Klahr, 2010) are thought to alleviate the tension required to form the cores of Jupiter, Saturn, and the icy giants in the Solar System (all within 30 au) and the lifetime of the protosolar disk (Lambrechts & Johansen, 2012).

The formation and evolution of pebbles and planetesimals, and hence planetary cores, relies upon radial pressure bumps that curb the catastrophic drift of solids towards the star due to aerodynamic drag (e.g., Pinilla et al., 2012). These radial dust traps have been associated with the ring systems observed in disks in the thermal emission from cold dust grains at radio wavelengths. Forming planets are thought to open gaps in the dust distribution, and the masses of such putative protoplanets are commonly inferred from the gap’s depth and width (Rosotti et al., 2016; Dong & Fung, 2017). It is frequent to associate one planet per gap (Dipierro et al., 2015; Mentiplay et al., 2018; Clarke et al., 2018), multiple planets to one (wide) gap (Dong et al., 2015), and one planet to multiple (narrow) gaps (e.g., a 0.1$M_{\rm Jup}$ planet can explain the multiple gaps in AS 209, Zhang et al., 2018).

In the so-called ‘gapped disks’, sometimes called transitional disks, dust is evacuated from central cavities or deep radial gaps, resulting in one or two bright rings. Accreting giant planets of masses comparable to or larger than that of Jupiter are able to open deep gaps in the gas and accumulate dust into an outer ring (Crida et al., 2006; Pinilla et al., 2012; Dipierro et al., 2016). By contrast, in disks without a deep gap, dust accumulation in mild pressure maxima produces sequences of several shallow rings (e.g., the TW Hya disk, Andrews et al., 2016). Clear signs of low-mass planet formation, as hinted by these fine shallow rings (Dong et al., 2017, 2018), have been absent in the outer regions of ‘gapped disks’.

The 1.7 $M_{\odot}$ star HD 169142 (estimated age 6${}^{+6}_{-3}$ Myr, Grady et al., 2007) hosts a prominent disk with deep wide gaps imaged in the near-infrared (Quanz et al., 2013; Momose et al., 2015; Pohl et al., 2017; Monnier et al., 2017; Ligi et al., 2018; Bertrang et al., 2018), mid-infrared (Honda et al., 2012), 1.3 millimeter (Fedele et al., 2017), and centimetre (Osorio et al., 2014; Macías et al., 2018) wavelengths. Circum-planetary features related to giant protoplanets have tentatively been reported in HD 169142 from radio/IR imaging (Reggiani et al., 2014; Biller et al., 2014; Osorio et al., 2014; Gratton et al., 2019), but the nature of these signals remains to be confirmed (Ligi et al., 2018).

In this paper, we present new 1.3 mm observations of the dust structures around HD 169142 at $\sim$2 au resolution (Sec. 2). In these new observations, the outer disk shows an intricate system of fine rings and narrow gaps (described in Sec. 2.2). Using hydrodynamical modeling (Sec. 3), a simple connection between the fine ring structure and a single, migrating, mini-Neptune protoplanet can be made, which we present and discuss in Sec. 3.3. Implications are discussed in Sec. 4.

2 Observations and data reduction

We obtained 1.3 millimeter observations by combining ALMA 12-m array extended (C40-9) and more compact (C40-6) configurations, resulting in baselines ranging from 19 meters to up to 13.9 kilometers and a total of 40-46 antennas. The combined observations are sensitive to spatial scales of up 2$\farcs$0, with a spatial resolution of $\sim$20 mas. The long baseline observations were acquired between Sept. and Nov. 2017 (Cycle 4) in four different blocks of $\sim$40 min each. Precipitable water vapor ranged between 0.7 and 1.8 mm. Observations of a phase calibrator (J1826-2924) were alternated with the science to calibrate the time-dependent variation of the complex gains. The cycling time for phase calibration was set to 8 minutes and 54 seconds for the compact and extended configurations, respectively. The ALMA correlator was configured in Frequency Division Mode (FDM). Two spectral windows with 1.875 GHz bandwidth were set up for detecting the dust continuum, centred at 232.0 GHz and 218.0 GHz, respectively.

All data were calibrated by the ALMA staff using the ALMA Pipeline version 40896 in the CASA package (McMullin et al., 2007), including offline Water Vapor Radiometer (WVR) calibration, system temperature correction, as well as bandpass, phase, and amplitude calibrations. The short baseline and long baseline datasets were calibrated independently. Self-calibration of the data was performed to improve coherence (a single iteration of phase-only), after which they were combined using the CASA task concat. A positional offset of (-20,-40) mas between the short and long baseline resulting images was measured and corrected prior to combining the datasets. This offset is consistent with the (-2.8,-38.0) mas/yr proper motion listed for this star (SIMBAD) and with the expected astrometric accuracy of ALMA (typically a 10th of the synthesized beam).

2.1 ALMA continuum imaging

Image reconstruction was initially performed using the CLEAN algorithm (CASA version 5.2, task tclean). As the source is relatively bright, we super-resolve the continuum data using our Maximum Entropy Method (MEM) package uvmem, which is part of the family of algorithms based on maximum-entropy regularization. Here we used the publicly-available GPU adaptation gpuvmem111The gpuvmem code is publicly available and can be found here: https://github.com/miguelcarcamov/gpuvmem. (Cárcamo et al., 2018), and regularized by minimizing the Laplacian of the model image with an objective function. The MEM reconstruction provides the highest resolution without compromising on sensitivity. We adopted the MEM image for our analysis. The final image has a peak of 67 mJy beam*-1* and an rms ($1\sigma$) of 14.7$\mu$Jy beam*-1* for a synthesized beam of 27$\times$20 mas. The total flux density over the whole image is 178$\pm$18 mJy (10% flux calibration uncertainty). The flux density of the central source at 1.3 millimeter is 194$\pm$20 $\mu$Jy. The bright ring at $\sim$0$\farcs$22 has an integrated flux of 59$\pm$6 mJy, and the flux density of the outer region is 126$\pm$13 mJy. A side-by-side comparison of the MEM versus CLEAN image reconstructions is shown in Fig. 5 (Appendix A).

2.2 A triple ring system in the disk’s outer region

The ALMA image (Fig. 1a) reveals a prominent bright ring (B1, see labelling in Fig. 1b) followed by substructure in the outer region in the form of a system of three fine rings. In previous observations the three rings (B2, B3, and B4) were seen as a single structure (we refer to this structure as the ‘outer region’), with evidence of ripples at the last scattering surface of micron-sized grains (Pohl et al., 2017; Gratton et al., 2019). In these new observations, we see that the outer region is a packed system of three rings, with mutual separations of $\sim$7 au between B2 and B3, and $\sim$12 au between B3 and B4, for a distance to HD 169142 of 113.9$\pm$0.8 pc (Gaia Collaboration et al., 2018).

Rings B2 and B3 are detected at a 5 to 20 $\sigma$ level (1$\sigma$ = 15$\mu$Jy beam*-1*), while emission from B4 reaches 4-5$\sigma$ over the image. For a quantitative characterization of the rings, we deproject the image by the disk inclination angle $i$ = 12.5^{\circ}$$\pm0.5∘ and transform to polar coordinates centred on the compact millimeter emission detected at the location of the star (see Appendix B). The polar deprojections zoomed on each ring are shown in Fig. 2. The disk azimuthally-averaged surface brightness profile is shown in Fig. 1b.

The locations and widths of the fine rings are determined from Gaussian fits to the average radial profile. Three Gaussian components were fit simultaneously plus a common low-order polynomial to account for low level background emission. Thus, the observed locations of B2, B3 and B4 are 0$\farcs$503$\pm$0$\farcs$005 (57.3 au), 0$\farcs$563$\pm$0$\farcs$005 (64.2 au), and 0$\farcs$667$\pm$0$\farcs$009 (76.0 au), respectively. Uncertainties only represent the error in fitting the Gaussian centroids. Since the rings are not circular (see below), the radii of each ring depend on azimuth, varying within $\pm$0$\farcs$01 around the mean locations quoted here.

The deconvolved width of each ring is determined by subtracting the beam width $\sigma_{\rm beam}={\rm FWHM}/\sqrt{8\log{2}}$ from the best fit Gaussian width in quadrature. Rings B2 and B3 have deconvolved widths of only 1.5 and 1.8 au, respectively, while B4 is more radially extended with a deconvolved width of 3.4 au. Rings B2 and B3 in HD 169142 are likely some of the narrowest structures in protoplanetary disks reported so far (see DSHARP ring characterization in Dullemond et al., 2018, where the narrowest ring is $\sim$3.4 au in deconvolved width).

Fig. 2 shows evidence for azimuthal dips along ring B2 spanning azimuthal angles $\sim$50∘ to 100∘ and $\sim$180∘ to 270∘, measured counter clockwise from the disk PA. Interestingly, all three outer rings have a finite eccentricity (0.09$\pm$0.03) and share a common focus. This is demonstrated in Fig. 2 which shows ellipses fit to the rings’ radii as a function of azimuthal angle (see Appendix B for details).

2.3 A structured ring and two deep wide gaps in the disk’s inner region

The brightest feature at 1.3 millimeter is the previously known inner ring B1, whose remarkable radial and azimuthal structure is further emphasized here at finer angular resolutions (Fig. 2). Previous hydrodynamical modelling have suggested that this structure is consistent with dynamical interactions with forming giant planets located inside the gaps D1 and D2 (Bertrang et al., 2018). The narrowness of ring B1 is explained by efficient radial trapping of dust in a pressure maximum, which puts an abrupt halt to the inward radial drift of dust particles $\geq$1 mm (Pinilla et al., 2012). Interestingly, as seen in Figs. 1a and 2, the structure in B1 may well be interpreted as a multiple ring. B1 is not split in the same way as the outer ring, rather it shows a very irregular and faint inner ring.

In addition, compact continuum emission is detected at the star’s location. We interpret this emission as coming from thermal emission from an inner disk of $\leq$1 au in radius, compatible with its near-infrared excess (Chen et al., 2018) as free-free emission is not expected to be significant in HD 169142 (Osorio et al., 2014).

The massive protoplanets required to explain the perturbed morphology of B1 are unlikely to split the outer dust disk into the observed triple ring structure (Bae et al., 2018). Instead, their dynamical interaction with the disk should lead to the clearing of deep and wide gaps such as D1 and D2, which are both $\sim$20 au in width and heavily depleted of dust and gas (Pohl et al., 2017; Fedele et al., 2017; Ligi et al., 2018; Bertrang et al., 2018). Yet, with the fine angular resolution provided by these new ALMA observations, we find that the outer region of HD 169142 is reminiscent of the ring systems in TW Hya (Andrews et al., 2016) and HL Tau (ALMA Partnership et al., 2015), both of which lack the presence of a deep wide gap in the dust continuum.

3 Hydrodynamical simulations and radiative transfer calculations

3.1 Theoretical background

Several models have been proposed to explain disks with rings. Enhanced dust growth at condensation fronts can produce fine rings (Banzatti et al., 2015; Zhang et al., 2015; Okuzumi et al., 2016), however, the location of rings and gaps are uncorrelated with the expected locations of snowlines (Long et al., 2018; Huang et al., 2018; van der Marel et al., 2019), and these rings are unlikely to be eccentric. Sharp variations in the gas viscosity associated with dead-zones (Flock et al., 2015; Miranda et al., 2017) and radially variable magnetic disk winds (Suriano et al., 2018) also yield ring-like structures, but they are driven by MHD instabilities that occur mostly within the central $\sim$20 au.

Theory of disk-planet interactions predicts that a low-mass planet can produce multiple dust rings by the propagation and dissipation of spiral density waves (Goodman & Rafikov, 2001; Dong et al., 2017, 2018). The dissipation of waves progressively moves dust away from the planet’s orbit, with the consequence that two narrow dust gaps form, one on each side of the planet’s orbit. The expelled dust thus forms two rings exterior to the gaps. Also, dust can remain at the planet’s (co-)orbital radius, forming a third ring between the two dust gaps (Dong et al., 2017). Previous hydrodynamic simulations have explored this scenario and found that, in a low-viscosity disk ($\alpha$ viscosity of 10*-4* or lower), a single low-mass planet indeed produces multiple rings in the dust distribution separated by annular gaps (Dong et al., 2017, 2018; Bae et al., 2018), observable at ALMA resolutions. The resulting ring system depends mainly on disk viscosity, disk temperature (aspect ratio), planet mass, and the elapsed time, which will dictate the number of additional spiral arms and hence the potential number of gaps (Bae et al., 2018).

Given the sharpness and small widths of the rings in HD 169142 (the deconvolved width of ring B2 is only $\sim$1.5 au), the $\alpha$ turbulent viscosity should be fairly small allowing a single low-mass planet to produce multiple rings and gaps. Assuming standard radial profiles for the disk structure, a relation that connects the distance between the two dust gaps, relative to the orbit of the planet, and the planet’s mass and disk scale-height can be found. Informed by hydrodynamic calculations (assuming $\alpha\lesssim 10^{-4}$, Dong et al., 2018), the relation reads

[TABLE]

where $r_{\rm IG}$ and $r_{\rm OG}$ are the locations of the inner and outer gaps around the planet, respectively, and $r_{\rm p}$ is the planet’s orbit. $M_{\rm th}=M_{\star}(h/r)^{3}$ is the disk thermal mass, $\gamma$ is the polytropic index (equal to 1 for isothermal gas), and $h/r$ is the disk aspect ratio. Fig. 1b shows that the narrow gaps D3 and D4 are separated by $\sim$9 au. In the outer region, the observed flux density of the rings translates into a brightness temperature of 8-12 K, consistent with a disk aspect ratio of 3%. Thus, for HD 169142’s stellar mass ($M_{\star}\!\approx\!1.7\,M_{\odot}$), disk aspect ratio ($h/r\!\approx\!0.03$), and approximating $r_{\rm p}$ by B3’s location at $\sim$64 au, a planet of $\lesssim$10 $M_{\oplus}$ could produce gaps separated by $r_{\rm D4}-r_{\rm D3}\approx\,9$ au.

3.2 Physical model and numerical setup

We perform hydrodynamic simulations of an outer disk composed of gas and a distribution of dust particles, to test whether the three bright rings (B2, B3, B4) unveiled by our ALMA observations can be explained by a single forming planet. The main observables are the distance between the gaps, the mutual separation between the rings, and the observed fluxes. In this study, we do not attempt to model the previously known inner features (D1 and B1) but rather focus on the new substructure of the outer region.

3.2.1 Hydrodynamical simulations

Our 2D hydrodynamical simulation of the gas and dust in the HD 169142 disk has been carried out with the code Dusty FARGO-ADSG. It is an extended version of the publicly available code FARGO-ADSG, which solves the gas hydrodynamics equations on a polar mesh (Masset, 2000; Baruteau & Masset, 2008a, b), and which models dust as Lagrangian test particles222The FARGO-ADSG code is publicly available at the address http://fargo.in2p3.fr/-FARGO-ADSG-. The version including a Lagrangian treatment of the dust particles, Dusty FARGO-ADSG, can be made available upon request to co-author Clement Baruteau. (Baruteau & Zhu, 2016; Fuente et al., 2017). The simulation has been designed to test whether the three bright rings (B2, B3, B4) unveiled by our ALMA observations can be explained by a single forming planet. In our simulation the star’s mass is assumed to be 2 $M_{\odot}$. We perform a set of simulations sampling the main parameter space {planet mass, disk scale-height, disk alpha viscosity}, although no exhaustive exploration is intended as finding a perfect fit is not in the scope of this work.

The gas momentum and continuity equations are solved on a polar mesh with 900 cells logarithmically spaced between 35 and 105 au in radius, and 1,200 cells evenly spaced between 0 and 2$\pi$ in azimuth. Non-reflecting boundary conditions are used to avoid reflections of the planet wakes at the radial edges of the computational grid.

The disk has an initial surface density profile decreasing with radius $r$ as $r^{-1}$ with an exponential cutoff beyond 100 au. The initial surface density is 2.9 g cm*-2* at 100 au, very close to the radiative transfer model which fits the CO istotopologue and continuum observations at 1.3 mm (Fedele et al., 2017). This large surface density means that Type I migration cannot be discarded, especially for the low planet masses we are exploring. A locally isothermal equation of state is used with the disk temperature decreasing with radius as $r^{-1/2}$, and equal to 8 K at 70 au (motivated by the observed brightness temperature), slightly lower than previous models. This low temperature, which is required to match the small mutual separations between the bright rings, translates into a disk aspect ratio profile in $r^{1/4}$ and equal to 0.033 at 69 au. Gas self-gravity is included since the disk’s Toomre $Q$-parameter varies from $\sim$6 to $\sim$9 throughout the computational grid. The $\alpha$ turbulent viscosity is set to $10^{-5}$. This low viscosity is motivated by the sharpness of the observed rings. Accretion onto the planet is not included.

The code solves the equations of motion for 200,000 dust particles with radii between 10 $\mu$m and 3 mm. Dust particles feel the gravity of the star, the planet, the disk gas (since gas self-gravity is included) and gas drag. Dust turbulent diffusion is also included as stochastic kicks on the particles position vector (Charnoz et al., 2011). However, dust self-gravity, dust drag, growth and fragmentation are not taken into account. Dust particles (assumed to be compact spheres of 2 g cm*-3* internal density) are inserted at the beginning of the simulation with a number density profile decreasing as $r^{-2}$ between 59 and 90 au. This corresponds to an initial dust surface density decreasing as $r^{-4}$. Although quite steep, the surface density evolves into a profile which increases with $r$ between the rings after 150 orbits. See the perturbed gas density as well as the distribution of dust particles shown in Fig. 3. Dust feedback onto the gas is discarded as the dust’s surface density along the B2, B3, and B4 rings remains comfortably smaller than the gas surface density in our model (see Fig. 3), for the gas-to-dust mass ratio assumed in the radiative transfer calculation (see Section 3.2.2).

A low-mass planet of a mini-Neptune mass (motivated by Sec. 3.1) is inserted at 69 au at the beginning of the simulation. The planet-to-star mass ratio is set to 1.7$\times$10*-5*, which for a 2 $M_{\odot}$ star translates into 11 $M_{\oplus}$333The stellar mass of HD 169142 is closer to $M_{\star}=1.7~{}M_{\odot}$, which yields a planet mass of 9.6 $M_{\oplus}$ in our simulation. Hence we round the mini-Neptune’s mass off to $\sim$10 $M_{\oplus}$ throughout the paper.. The planet migrates due to disk-planet interactions and reaches $\sim$64 au after 150 orbits, which is very close to the location of the B3 ring. Higher mass planets (here we considered planet-to-star mass ratio of 6$\times$10*-5*, i.e., a $\sim$33 $M_{\oplus}$ planet) produce significantly wider gaps even at higher viscosities (a simulation with $\alpha=10^{-4}$ was performed, not shown here).

The time that best reproduces our observations of the outer rings is approximately at 150 orbits of the mini Neptune, i.e., $\sim$0.05 Myr. We note that the rings are present since orbit $\sim$100 onwards (see Fig. 3). The model is stable for at least 1000 planet orbits ($\sim$0.4 Myr). After this, a global $m$=$1$ mode grows in the disk, leading both the planet and the disk to develop some eccentricities. After this, the co-rotating ring gradually diffuses with time as the gas surface density decreases.

3.2.2 Radiative transfer calculation

We compute the continuum emission from our dust simulations using the public radiative transfer (RT) code radmc3d (version 0.41). Twenty logarithmically-spaced bins are used to allocate the 200,000 dust particles from the simulation. The dust sizes ($s$) follow a power-law distribution $n(s)\propto s^{-3.5}$, with minimum and maximum sizes of $1\,\mu{\rm m}$ and $1{\rm mm}$, respectively. Assuming a gas-to-dust ratio of 33 yields a total dust mass in the outer region of $\sim$100$M_{\oplus}$.

Each 2D distribution of dust particles is turned into a surface density and expanded in the vertical direction assuming hydrostatic equilibrium. The dust’s scale height $h_{i,{\rm d}}$ is size-dependent and follows

[TABLE]

where $h$ is the initial gas pressure scale-height, ${\rm St}_{i}$ is the averaged Stokes number for the $i^{\rm th}$ dust size bin, and $D_{\rm z}$ is a turbulent diffusion coefficient in the vertical direction which depends on the level of turbulent activity across the vertical extent of the disk (Yang et al., 2018). We assume $D_{z}$ to be proportional to the $\alpha$ turbulent viscosity in our 2D hydrodynamic simulation (see Baruteau et al., 2019, for more details). A $D_{z}$ value equal to 10 $\alpha$ yields dust temperatures which match the observed brightness temperature of rings B2 and B3. In the RT calculation, we use 24 cells logarithmically-spaced to sample 3 scale-heights in colatitude (with finer sampling towards the midplane).

We expand the hydrodynamic simulation’s grid to include the inner regions and then artificially add the inner disk and the bright ring B1. Their parameterisation follows previous RT modelling which fits the 1.3 millimeter data (Fedele et al., 2017). The presence of these structures that block stellar radiation are needed to recover the low temperatures of the outer regions. A 7320 K star with a radius of 1.5 $R_{\odot}$ (Gaia Collaboration et al., 2018) is placed at the centre of the grid and $10^{8}$ photon packages are propagated to compute the temperature of each dust bin size via the Monte Carlo radmc3d mctherm task. The temperatures range between $\sim$13 and $\sim$20 K in the outer rings, consistent with previous multi-wavelength modelling (Fedele et al., 2017). These values are 2$\times$ higher than the locally isothermal prescription used in the hydrodynamical simulation, which is representative of the gas scale-height, rather than a physical temperature of dust grains. Ray-tracing is performed to solve the radiative transfer equation for continuum emission assuming the dust is a mix of 30% amorphous carbons and 70% silicates. The resulting synthetic image is convolved with the MEM resolution beam and it is shown in Fig. 3. The high fidelity of the ALMA image allows direct comparison with the data in the image plane. The hydrodynamical model reproduces the locations of B2, B3, and B4 as well as the level of flux along rings B2 and B3, and to a lesser extend ring B4’s. This can be best appreciated in Fig. 8 (Appendix C).

3.3 Results

A hydrodynamical model including a 10 $M_{\oplus}$ planet embedded in a disk with an aspect ratio $h=0.033$ reproduces the observations (Fig. 4). The initial condition is such that the planet is introduced at 69 au. After about fifty thousand years, the planet migrates inwards to a $\sim$64 au orbit, shepherding ring B3 and “pushing” rings B2 and B4, thus explaining the asymmetry in the rings mutual separations. This can be appreciated in Fig. 3. A simulation with a planet on a fixed orbit (not included) produces rings that are too equidistant for the disk and planet parameters assumed here. If the 10 $M_{\oplus}$ core formed closer to the star and then migrated outwards, it would not be consistent with the rings locations: the observed asymmetry reveals the direction of migration. This migration could have distinct observational signatures in multiwavelength observations, for example, by probing the spatial segregation of larger grains with longer wavelength observations. Such predictions have recently been reported by Meru et al. (2019) for a planet migrating in a high viscosity disk. Despite the little knowledge we have about the initial conditions of the dust density distribution, the model reproduces the observed fluxes along B2 and B3, albeit with a slightly brighter ring B4. The initial amount of large particles near B4’s location is likely overestimated in our simulations.

As the planet continues to interact dynamically with the dust reservoir, it also starts clearing an azimuthal opening in its vicinity. This produces a dip in azimuth along the middle ring at the planet’s location. Whether this dip is observable at 1.3 millimeters depends on the planet’s mass, the maximum size in the dust distribution, and the elapsed time. The more massive the planet, the faster and deeper is the opening of the azimuthal dip. On the other hand, a large maximum particle size produces a more pronounced dip. Dust growth and fragmentation, which are not taken into account in our simulations, probably also have an effect on the observability of this azimuthal opening. The lack of a clear azimuthal dip along ring B3 indicates that the bulk of particles we observe have sizes $<$1 mm, the planet is rather small ($\lesssim$10 $M_{\oplus}$), or that it may have formed less than a few hundred orbits ago ($\sim$50 kyr), or most likely a combination of these options. Note that the signal-to-noise ratio at the rings location might not be sufficient to probe the azimuthal structure along B3. Longer wavelength observations might be able to trace the location of the mini-Neptune as the larger grains pile-up behind the planet (note mm particles in Fig. 3).

3.4 Discussion

3.4.1 Could rings B2, B3 and B4 be shaped by two low-mass planets?

A scenario where two low mass planets (one for each narrow gap) results in the triple ringed structure is unlikely for several reasons. The small separation between the gaps (only 9 au) requires scale-heights at the location of the planets which would be too small for a realistic disk temperature, these translate into $<$5 K at 64 au (following eq. 16 in Dong & Fung, 2017). In other words, if you put two planets into the two gaps, they will likely open one common gap instead of two. Moreover, the similarity between D3 and D4 suggests that if they are produced by two planets, these ought to have formed at similar times and with similar masses. The dynamical interactions between the closely orbiting planets might also make them unstable.

3.4.2 How does an inner giant planet impact the structure of rings B2, B3 and B4?

In recent years, it has been shown that planet perturbations can extend far beyond the Hill radius if the disk viscosity is low (Bae et al., 2018). Since the model presented here assumes a very low viscosity, any inner giant planet (in D2, for example) could have an impact on the outer triple ring system. In a debris disk or a protoplanetary disk with low gas surface density, giant planets could potentially produce fine structures by mean motion resonances with the dust. However, HD 169142 is gas rich and the dust in the outer region is still coupled to the gas distribution. To further test the possibility that the fine structure is somehow related to inner giant planets, we included a giant planet in our simulation, of roughly a Jupiter mass located in the middle of D2 at 38 au (see Appendix C). The giant indeed opens a deep gap in the gas, but fails to generate any additional narrow rings in the outer disk. However, the spiral wakes excited by the giant and the vortices at the gap edge (produced by the Rossby wave instability) induce eccentricity on the outer triple rings associated with the mini-Neptune (Fig. 6). Recent scattered light imaging indeed shows faint spiral features at the disk surface which can be associated to putative giant planets in the wide gaps (Gratton et al., 2019). The dust particles around the mini-Neptune can acquire non-circular trajectories either (i) by direct gravitational interaction with the mini-Neptune, whose eccentricity varies on account of the inner giant’s wakes depositing energy and angular momentum in the planet’s horseshoe region (Fig. 7), and/or (ii) by being deflected by the (shock) wakes upon crossing them.

The addition of a giant planet adds significant complexities as the final configuration of rings would strongly depend on the initial conditions for the dust distribution and the precise timing of the formation of the mini Neptune and the giant planet. As protoplanetary disk observations grow in sensitivity and resolution, more and new ingredients need to be considered in the modeling efforts. This paper is meant as a proof of concept to show that a simple one-planet simulation can account for seemingly complex sub-structure in the outer ring of a transition disk. Finding a perfect fit to the data is not in the scope of this work. Indeed, some important aspect have been left aside such as: dust growth (e.g., Bae et al., 2018) and fragmentation, dust back-reaction onto the gas (e.g., Gonzalez et al., 2017; Dipierro et al., 2018), dust torque due to scattered pebbles (Benítez-Llambay & Pessah, 2018), 3-D and MHD effects (e.g., Flock et al., 2015; Miranda et al., 2017), and possibly many more. The inclusion of planet migration was needed here to explain the asymmetry in the rings’ mutual separations, a step forward towards understanding how low-mass planet migration can be studied from observations of dust radial structures. However, a proper account of dust dynamics, planet accretion and thermodynamics, in a self-consistent way, is necessary to build the migration history of these low-mass protoplanets (Benítez-Llambay & Pessah, 2018).

3.4.3 What can these observations tell us about the planet formation process?

The low mass of the putative protoplanet in HD 169142, in addition to the lack of clear evidence for structures associated with gravitational instabilities in (abundant) multi-wavelength observations of this source, suggests that a bottom-up process such as core accretion could be responsible for its formation. This implies that core accretion can potentially operate at $\sim$65 au within the age of the system ($\sim$6 Myr, Grady et al., 2007), possibly assisted by pebble accretion as the outer region banks more than 100 $M_{\oplus}$ worth of pebble-sized solids (subject to uncertainties of the dust opacities). The modeling presented here also shows that the mini-Neptune should have formed well after the giant planets carved gaps D1 and D2. The presence of a giant planet in D2 would have produced a pressure maxima in the outer region which could have enhanced core accretion via a dust trap (Pinilla et al., 2012).

4 Concluding remarks

The new ALMA observations presented here show that a transition disk with wide deep gaps can also host narrow-ring structures in its outer region, similar to those observed in the HL Tau and the TW Hya systems.

The HD 169142 observation allows to link the architecture of protoplanetary disks with low mass planets. The interpretation via hydrodynamics is a proof of concept that links the structure of closely packed double gaps and tripple rings with a single and migrating low-mass planet. Planetary migration naturally explains the distinct mutual separations between the narrow rings. The connection was made possible in HD 169142 thanks to the isolation of its outer region. In the absence of a clear gap that separates an outer ring, the superposition of multiple rings due to several planets hampers simple and clear explanations such as that found for HD 169142.

In HD 169142, we have thus found evidence that suggests that low-mass planet formation can occur in the outer regions of disks bearing evidence for giant protoplanets. The planet formation mechanism, likely core accretion or any bottom-up process, can thus produce planet embryos at $\sim$65 au and outside the orbit of inner giants, at least in certain disks.

We thank Ed Fomalont, Anya Yermakova and Philipp Weber for useful discussions, as well as our anonymous referee for their constructive comments. Financial support was provided by the government of Chile grants Millennium Scientific Initiative RC130007, CONICYT-Gemini 32130007, and CONICYT-FONDECYT grant numbers 1171624, 1171246 and 1191934. S.P acknowledges support from the Joint Committee of ESO and the Government of Chile. The data analysis and some of the simulations were carried out in the Brelka cluster, hosted at DAS/U. de Chile (Fondequip EQM140101). Numerical simulations were carried out on the CalMip machine of the Centre Interuniversitaire de Calcul de Toulouse, which is gratefully acknowledged. This work was performed in part at the Aspen Center for Physics, which is supported by National Science Foundation grant PHY-1607611. This paper makes use of the following ALMA data: ADS/JAO.ALMA#2016.1.00344.S. 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. This work has made use of data from the European Space Agency (ESA) mission Gaia, processed by the Gaia Data Processing and Analysis Consortium (DPAC). Funding for the DPAC has been provided by national institutions, in particular the institutions participating in the Gaia Multilateral Agreement.

Appendix A Comparison between CLEAN and MEM algorithm

A side-by-side comparison of HD 169142 image synthesis using CLEAN and MEM algorithms is shown in Fig. 5. The triple ring structure in the outer region is robust in all image reconstructions. The ‘uniform’ weighting CLEAN image has a comparable resolution to MEM but at the cost of decreasing sensitivity.

Appendix B Polar deprojection and ellipse fitting

To study the eccentricity of the rings, we extract the radial profile of each ring as a function of azimuth, $R(\theta)$, by fitting gaussian functions in radius (using the Python iminuit optimization library). Only signal above $5\sigma$ is included in the process. Ring B4 is not considered as it only has a small number of points with high signal-to-noise. B1 is also left aside as its perturbed morphology cannot be reproduced by a simple elliptical curve. The optimization yields the values for $R(\theta)$ and their error bars $\Delta R(\theta)$, for each ring. The ellipse equation is then fitted to $R(\theta)$ for a given ring using a direct least square procedure (Fitzgibbon et al., 1996). The fitting procedure yields the ellipse centre, its angle of rotation and its eccentricity ($e=\sqrt{1-(b/a)^{2}}$ where $a$ and $b$ are the semi-major and semi-minor axes, respectively). The best fit ellipse for each ring are plotted on the polar projection in Fig. 2. The uncertainties in the eccentricity values are calculated using a Monte-Carlo approach where data points are modified by random numbers of the order of $3\Delta\theta$ and then fitted following the aforementioned procedure. This is repeated a thousand times and the final errors are drawn from the standard deviation. Correlated noise is accounted for by multiplying the centroids’ error by a factor $\sqrt{N_{\rm pix}}$, where $N_{\rm pix}=5$ is the number of pixels along the synthesized beam’s major axis.

The uncertainty of the eccentricity is dominated by our imprecise knowledge of the inclination angle of the disk. In order to determine the significance of the fitted values, we explore a range of inclination angles ranging between 9∘ and 16∘, repeating the fitting procedure described above for each inclination value. The results are shown in Fig. 2b. Rings B2 and B3 show eccentricities which indeed vary with inclination angle, with a minimum eccentricity of e$$\approx$$0.085. Interestingly, the rings have different eccentricities for a given inclination, except at i$$\approx$$12.5^{\circ}, where both rings have the same eccentricity. This inclination angle is close to the standard value of i$$\sim$$13^{\circ} used in previous multi-wavelength modelling (e.g. Fedele et al., 2017; Bertrang et al., 2018). The intrinsic position angles of the ellipses also coincides at i$$\approx$$12.5^{\circ}. This fitting procedure thus suggests a common inclination for both B2 and B3 of i$$\approx$$12.5^{\circ}\pm 0.5\deg (1$\sigma$). At this inclination, the eccentricity of B2 and B3 share the common value of 0.09$\pm$0.02 (2$\sigma$ uncertainty, same as the shaded area in Fig. 2). A disk inclination angle different from $12.5^{\circ}$ yields different eccentricity values for B2 and B3, which would suggest the occurrence of warping in the outer regions.

Appendix C Addition of a giant planet associated with gap D2

The impact that an inner giant could have on the outer narrow rings is explored by including a 0.5 $M_{\rm Jup}$ planet fixed at 38 au, in the middle of the observed gap D2, in addition to the mini-Neptune. The general setup inherits the same parameters as the single mini-Neptune simulation, with some changes. The initial radius of the orbit of the mini-Neptune planet is decreased to 68 au, this is to account for a slightly lower migration rate in the presence of the inner giant. The grid resolution is lowered to 600$\times$900 cells in radius and azimuth, respectively. The grid’s inner edge is extended down to 20 au as to include gap D2 in the radial domain. The smaller inner edge translates into a significant computational cost. The number of particles is decreased by half (i.e only 100,000 particles are used) so as to make the computation less expensive. Dust particles are distributed over the same region as in the previous simulation. Dust around the giant planet is not included as this would require increasing the number of particles, making the simulation too expensive. As mentioned in the main text, fitting the already known features B1 and D1 is beyond the scope of this work. Modeling all features in HD 169142 system requires an extensive search of a large parameter space, but the main difficulty lies on the little knowledge available of the initial conditions and the relative age of each planet. In our simple two-planet setup, both companions are introduced simultaneously at the beginning of the simulation.

Fig. 6 shows that the presence of a giant planet at 38 au carves a wide gap in the gas at the location of the D2 dust gap. At the same time, the mini-Neptune shapes the outer region into the triple ring system, consistent with the single planet simulation. The presence of the giant planet makes the morphology of the outer rings appear eccentric. The asymmetry in the outer rings mutual separation is also reproduced. The migration of the mini-Neptune is slightly affected by the gas giant’s wakes (see Fig. 7). Vortices develop at the outer edge of the giant’s gap but these start decaying shortly after a couple of hundred orbits of the mini-Neptune (see right panel in Fig. 6). The lack of observational evidence for these vortices may suggest that the mini-Neptune formed long after the giant planet carved the gap D2, at a time where the vortices had already decayed.