Recent ALMA observations have revealed that a large fraction of protoplanetary discs contain bright rings at submillimetre wavelengths. Dust trapping induced by pressure maxima in the gas disc is a popular explanation for these rings. However it is unclear whether such pressure bumps can survive for evolutionary timescales of the disc. In this work we investigate an alternative scenario which involves only dustgas interactions in a smooth gas disc. We postulate that ALMA rings are a manifestation of a dense clumpy midplane that is actively forming planetesimals. The clumpy medium itself hardly experiences radial drift but clumps lose mass by disintegration and vertical transport and planetesimal formation. Starting from a seed ring we numerically solve the transport equations to investigate the rings survival. In general rings move outwards due to diffusion of the clump component. Without pressure support rings leak material at rates 40 rm Moplus  mathrmMyr1 and in order for rings to survive they must feed from an external mass reservoir of pebbles. In the case where the pebble size is constant in the disc a cycle between ring formation and dispersion emerges. Rings produce large quantities of planetesimals which could be material for planet formation and explain the massive budget inferred debris disc. Mock images of ALMA observations compare well to the rings of Elias 24 and AS 209 from DSHARPs sample.
['protoplanetary discs', 'circumstellar matter', 'planets and satellites: formation', 'submillimetre: planetary systems', 'Astrophysics - Earth and Planetary Astrophysics']
Survival of ALMA rings in the absence of pressure maxima
{'No Header': ', 1-18 (2021)', 'Haochang Jiang ( 蒋 昊昌 ) 1 , Chris W. Ormel 1 ★': '1 Department of Astronomy, Tsinghua University, Haidian DS 100084, Beĳing, China  \nAccepted XXX. Received YYY; in original form ZZZ', 'ABSTRACT': "Recent ALMA observations have revealed that a large fraction of protoplanetary discs contain bright rings at (sub)millimeter wavelengths. Dust trapping induced by pressure maxima in the gas disc is a popular explanation for these rings. However, it is unclear whether such pressure bumps can survive for evolutionary time-scales of the disc. In this work, we investigate an alternative scenario, which involves only dust-gas interactions in a smooth gas disc. We postulate that ALMA rings are a manifestation of a dense, clumpy mid-plane that is actively forming planetesimals. The clumpy medium itself hardly experiences radial drift, but clumps lose mass by disintegration and vertical transport and planetesimal formation. Starting from a seed ring, we numerically solve the transport equations to investigate the ring's survival. In general, rings move outward, due to diffusion of the clump component. Without pressure support, rings leak material at rates GLYPH<24> 40 𝑀 GLYPH<8> Myr GLYPH<0> 1 and in order for rings to survive, they must feed from an external mass reservoir of pebbles. In the case where the pebble size is constant in the disk, a cycle between ring formation and dispersion emerges. Rings produce large quantities of planetesimals, which could be material for planet formation and explain the massive budget inferred debris disc. Mock images of ALMA observations compare well to the rings of Elias 24 and AS 209 from DSHARP's sample.  \nKey words: protoplanetary discs - circumstellar matter - planets and satellites: formation - submillimetre: planetary systems", '1 INTRODUCTION': 'Thanks to the high sensitivity and spatial resolution, the Atacama Large Millimeter/submillimeter Array (ALMA) has provided us with unprecedented imagery starting with HL tau (ALMA Partnership et al. 2015). At (sub)millimeter wavelengths, these images reveal the disc mid-plane materials, which harbors most of the discs\' mass. At continuum bands, the emission originates from so-called pebblesized particles - particles of moderate size, significantly larger than micron-size ISM dust grains but small enough to provide an opacity at GLYPH<24> mm wavelengths. These are the precursors of planetary objects: planetesimals and the cores of giant planets.  \nIntriguingly, ALMA found that these planet-building blocks are not homogeneously distributed, contrasting decades-long beliefs that disks are \'smooth\' (Weidenschilling 1977a; Hayashi 1981). Instead, disks manifest structure in diverse forms, ranging from spirals, gaps, lopsided arcs to rings. Particularly abundant are the axisymmetric rings, which appear in at least 17 out of 18 samples of the DSHARP survey (Huang et al. 2018). Even though the DSHARP survey is biased towards bright sources, these annular substructures are hosted by stars across a range of luminosities, masses, and accretion rates (Zhang et al. 2016; Andrews et al. 2018; Long et al. 2018; van der Marel et al. 2019; Cieza et al. 2019; Francis & van der Marel 2020). The morphology of the substructure shows remarkable diversity as well. Rings can appears at almost any radius (5 GLYPH<0> 150 au) and their widths range from a few astronomical units to tens of astronomical units even in the disc. This observational diversity perhaps suggests that the origins of substructures may also be various (Andrews 2020).  \nThere are many hypotheses aiming to explain the observed substructures. Arcs may reveal that pebbles are concentrated by vortices, while spirals or other non-axisymmetric structure may hint at the presence of gravitationally unstable discs (Toomre 1964; Boss 1997; Kratter & Lodato 2016) or massive companions (Dong et al. 2015; Zhu et al. 2015b; Bae & Zhu 2018). However, nonaxisymmetric substructures are relatively rare (Andrews 2020). In this work, we focus on the most common substructure - rings. Explanation of substructure formation can be grouped in two camps: with and without pressure maxima.  \nThe iceline scenario is one possible mechanism where the pressure bump is not needed for ring formation. The iceline explanation exploits that the mechanical properties of particles depend on the composition of the surface material which makes up the dust grains. For example, ice mantles can evaporate off grains when particles are exposed to higher ambient temperatures (e.g., Zhang et al. 2015; Pinilla et al. 2017) or the particles become more brittle due to sintering (Okuzumi et al. 2016; Sirono & Ueno 2017). This results in a radial dependence of -say- the fragmentation threshold, which in turn changes the physical characteristics (the particle size), their concentration (smaller solids drift slower and pileup), and hence the appearance. Therefore, the particles size can change rapidly after they cross condensation fronts. Various hydrodynamic models caused by gas-particle coupling can also generate annular substructure, for example, self-induced pileup of particles by aerodynamical feedback (Drążkowska et al. 2016; Gonzalez et al. 2017). When gas drag slows down the solids\' self-gravitational collapse, secular gravitational instability (SGI) happens and breaks down the over-densities into narrow rings (e.g., Takahashi & Inutsuka 2014; Tominaga et al. 2020).  \nAmore direct explanation is particles trapping due to a reversal of the pressure gradient (Kretke & Lin 2007; Pinilla et al. 2012; Dullemond et al. 2018). In smooth disks, the pressure radially decreases from the hot and dense inner disk to the cool and more tenuous outer disk. This pressure gradient provides hydrostatic support, causing the disk to rotate at sub-Keplerian velocity and causing particles to drift inwards (Weidenschilling 1977b; Brauer et al. 2007). A local reversal of the pressure gradient, however, would result in the opposite behaviour: gas would rotate super-Keplerian and particles would drift outwards. Hence, particles pile up at pressure maxima. Applied to ALMA rings, the idea is that rings are associated with these pressure maxima. The boundary of magnetically dead and active zones is one possible site to set such a pressure bump for ring formation (e.g., Flock et al. 2015; Pinilla et al. 2016). MHD zonal flows, creating narrow enhancements in pressure, are expected to concentrate solids (e.g., Uribe et al. 2011; Bai & Stone 2014; Suriano et al. 2017, 2018). A more popular explanation are disc-planet gravitational interactions. As planets gravitational push aside the gas in the disc (Lin & Papaloizou 1986), they are naturally considered to be the \'unseen\' driving forces behind pressure bumps (e.g., Rice et al. 2006; Paardekooper & Mellema 2006; Zhu et al. 2012; Dipierro et al. 2016; Bae et al. 2017). In its vicinity, the planet changes the rotation profile in a non-axisymmetric fashion. When this kinematic signal is observed (Teague et al. 2018; Pinte et al. 2020) the case for the planet explanation is particularly compelling.  \nThere are several reasons why it is unlikely that all rings are associated with planets. First, planets need to be of appreciable mass to perturb the disk, unless the disk is extremely inviscid (Rosotti et al. 2016; Dong et al. 2017). However, from GPI exoplanet survey and SPHERE infrared survey for exoplanets (SHINE), it is known that wide-orbit ( ¡ 5 au) giant planets are rare (occurrence glyph[lessorsimilar] 5% for FGK stars, Nielsen et al. 2019; Vigan et al. 2020). Another problem is that rings, due to their ubiquity, need to last for evolutionary times. In other words, how do rings survive for million years. Specifically, a planet large enough to open a gap in the gas disk, is also in its migration "sweet spot" (between the Type I and Type II migration regimes) and as such also tends to migrate on time-scales shorter than the evolutionary times. In addition, the pileup of material in the ring region should provide the optimal conditions of ongoing planet formation - a high density of material in a pressure bump begging the question how this material just sits around for GLYPH<24> Myr. Finally, blocking all pebbles at large distance is clearly unfavorable for the ability to form planets in the interior disc, conflicting the high occurrence of close-in planets inferred from the exoplanet census (e.g., Zhu et al. 2018).  \nIn this work, we investigate whether the survival of rings can yet be reconciled by a model where the gas disk is none the less smooth. Moreprecisely, we focus not ontheformation of the annular structure, which may well be due to one of the afore mentioned ring formation processes, like a pressure bump, but could also be attributed to the star formation process (e.g., Bate et al. 2010), but rather on its longterm survival in a disk characterized by a "standard" smooth pressure gradient. In such setups rings will leak a great amount of material and, in order to survive, must be replenished by an external reservoir of pebbles. Our goal is to quantify the conditions that must be fulfilled in order to match the ALMA and DSHARP observations, and to put the implications of our ring model in the context of planet formation.  \nParticles clumping has been observed in many shearing-box simulations studying the streaming instability (e.g., Johansen & Youdin 2007; Bai & Stone 2010). Systematic studies on clumping in stratified discs reveal that axisymmetric pebble filaments form via streaming instabilities, before they collapse into planetesimals (Yang & Jo-  \nhansen 2014; Yang et al. 2017; Li et al. 2018). Motivated by the pebble clumping, we conceive a scenario where solids are divided into two components: freely-drifting "disc pebbles" and nearly-stationary (mid-plane) "clump pebbles". In our simulations, a seed ring is initialized at 10-100 au. This perturbation will grow by absorbing drifting pebbles at its leading edge, but also lose a great number of pebbles at its back as there is no pressure confinement. In the clumpy ring model (hereafter CRM), planetesimal formation is also accounted for and plays a crucial role during the rings evolution. By numerically solving the 1D mass transport equations, we assess the conditions for which a ring-like feature can be sustained over several million years.  \nThe plan of the paper is as follows. The dust transport models are described in Section 2. The numerical results are presented in Section 3. We apply this model to DSHARP samples and compare with observations in Section 4. Section 5 presents an assessment of the model and outlines implications for planet formation. We summarize the main results and conclusions in Section 6.', '2 MODEL': 'In the CRM, it is assumed that a single particle size dominates the mass budget of the pebbles. We do not consider a particle size distribution. Solids are divided among three components:  \n- (i) Disc pebbles. Dust that is vertically extended, homogeneously distributed, and which drifts inwards.\n- (ii) Clump pebbles. Dust that resides in a clumpy state in the midplane regions. Clumps hardly drift, but are subject to radial diffusion.\n- (iii) Planetesimals. Km-sized particles formed after clumps have experienced gravitational collapse. They do not interact with the other two components.  \nAcartoon of the CRM is shown in Fig. 1. The disc gas surface density profile has no substructure and follows a standard power-law. Disc pebbles and clump pebbles continually interchange material vertically by settling and diffusion. Clump pebbles collapse into planetesimals. For simplicity in this paper, planetesimals do not interact with the other components. We treat disc pebbles as vertically extended, at a scale height no larger than the gas scale height. As the dust-to-gas ratio of disc pebble increases, the settling of disc pebble creates a new sublayer in the protoplanet disc (Sekiya 1998; Youdin & Shu 2002). Within this layer the local dust-to-gas ratio exceeds unity; it rotates Keplerian and pebbles trapped in this layer hardly drift toward the central star. Clumps yet diffuse radially. Clumps either disperse to return pebbles back to the disc component or collapse to trigger planetesimal formation.', '2.1 Transport Model': 'We consider the time evolution of the surface density of disc dust and clump dust in 1D cylindrical geometry. The transport equation for the two pebble components (disc pebbles and clump pebbles) are  \n𝜕 Σ d 𝜕𝑡 , 1 𝑟 𝜕 𝜕𝑟 GLYPH<0> 𝑟 Σ d 𝑣 dr GLYPH<150> d GLYPH<1> GLYPH<0> 1 𝑟 𝜕 𝜕𝑟 GLYPH<20> 𝑟𝐷 d Σ g 𝜕 𝜕𝑟 GLYPH<18> Σ d Σ g GLYPH<19> GLYPH<21> = 𝑆 cp 𝜕 Σ c 𝜕𝑡 , 1 𝑟 𝜕 𝜕𝑟 GLYPH<0> 𝑟 Σ c 𝑣 dr GLYPH<150> c GLYPH<1> GLYPH<0> 1 𝑟 𝜕 𝜕𝑟 GLYPH<20> 𝑟𝐷 c Σ g 𝜕 𝜕𝑟 GLYPH<18> Σ c Σ g GLYPH<19> GLYPH<21> = GLYPH<0> 𝑆 cp GLYPH<0> 𝑆 plt (1)  \nwhere Σ g, Σ d , Σ c are, respectively, the gas, disc pebble, and clump pebble surface densities. In analogy to the Shakura & Sunyaev  \nFigure 1. Sketch of the model in this paper. The three components are separately shown in black (disc pebbles), purple (clump pebbles; the ring) and brown (planetesimals). A constant external mass flux / 𝑀 ext is imposed on the outer boundary 𝑟 out where disc pebbles drift in at a velocity 𝑣 dr GLYPH<150> p . Disc and clump pebbles continually interchange material by settling (followed by their incorporation into clumps) and clump dispersal (followed by vertical diffusion). Clump pebbles inside the ring hardly drift toward the central star ( 𝑣 dr GLYPH<150> c GLYPH<28> 𝑣 dr GLYPH<150> p ) and can collapse into planetesimals, which in the CRM do not interact with the other components. The clump layer has a thickness GLYPH<24> 𝜂𝑟 , and the clump particles are characterized by a low but non-zero radial diffusivity ( 𝐷 c) due to clump dispersal and vertical diffusion. A leaking mass flux / 𝑀 leak is released from the ring.  \n<!-- image -->  \n(1973)\'s parametrization for the shear stress, the turbulent particle diffusivity of the clump pebbles is defined as 𝐷 c = 𝛼 c Ω 𝐾 𝐻 2 g (Youdin & Lithwick 2007; Yang et al. 2018), where Ω 𝐾 \' 𝑟 \' = √︁ 𝐺𝑀 ★ GLYPH<157> 𝑟 3 is the local Keplerian frequency with host stellar mass 𝑀 ★ . Similarly, for disc pebble it is 𝐷 d = 𝛼 d Ω 𝐾 𝐻 2 g . As the radial transport of disc pebble is dominated by advection, we omit its radial diffusion in most of our simulations ( 𝛼 d = 0), but see Sect. 3.5 for the non-zero 𝛼 d cases.  \nThe radial drift velocity is (e.g. Nakagawa et al. 1986, hereafter referred to as the NSH solution)  \n𝑣 dr = St \' 1 , 𝑍 \' 2 , St 2 𝑐 2 𝑠 Ω 𝐾 𝑟 𝜕 log 𝑃 𝜕 log 𝑟 = GLYPH<0> 2St \' 1 , 𝑍 \' 2 , St 2 𝜂𝑣 𝐾 (2)  \nwhere 𝑣 𝐾 = Ω 𝐾 𝑟 is the Keplerian velocity, 𝑍 is the local dust-to-gas mass ratio of each component, and 𝜂 GLYPH<17> GLYPH<0> \' 𝑐 2 𝑠 GLYPH<157> 2 𝑣 2 𝐾 \'\' 𝜕 log 𝑃 GLYPH<157> 𝜕 log 𝑟 \' is a dimensionless measure of the radial pressure gradient. For disc pebble and clump pebble, 𝑍 takes on different values according to their surface densities and scale heights (see Sect. 2.2 and Sect. 2.3).  \nThedrift velocity of disk pebbles is always given by Eq. (2). However, for clump pebbles we consider several options for their drift motion:  \n- (i) Clumps drift at the NSH solution\n- (ii) Clumps drift faster than the NSH solution\n- (iii) Clumps do not drift at all, 𝑣 dr = 0  \nIn the first two cases the drift velocity is determined by the 𝑓 ca parameter  \n𝑣 dr GLYPH<150> c = 𝑓 ca 𝑣 dr \' 𝑍 \' , \' 1 GLYPH<0> 𝑓 ca \' 𝑣 idv (3)  \nthat is, a mixture between the individual velocity ( 𝑓 ca = 0), where 𝑣 idv = 𝑣 dr \' 𝑍 = 0 \' , and the NSH solution ( 𝑓 ca = 1). The value of 𝑓 ca that we used is an input parameter to the CRM. It is listed in Table 1. When it is absent, the clump advection velocity used in Eq. (1) is 0 (case iii).  \nA clump velocity exceeding the NSH solution is motivated by the simulations of Johansen & Youdin (2007), where it was found that radial drift can be faster than the NSH solution. However, these simulations were conducted without vertical gravity. On the other hand, Yang & Johansen (2014) and Yang et al. (2017), accounting  \nfor vertical gravity, found that the regions with above-average pebble concentrations (a.k.a. clumps), were characterized by a virtually-zero clump drift.  \nThe RHS of Equation (1) are source terms. The material exchange between mid-plane and upper layer disc is  \n𝑆 cp = 𝑓 𝜏 Σ c 𝑡 diff GLYPH<0> 𝑓 c Σ d 𝑡 sett GLYPH<149> (4)  \nThe first term on the right of Equation (4) is the material lost from clumps due to turbulent diffusion and the second term describes the gain of clumps from disc pebbles that settle to the mid-plane. There are two modulation factors, 𝑓 𝜏 and 𝑓 c, with values between 0 and 1. The factor 𝑓 𝜏 captures the suppression of the (net) diffusion due to clumping depth effects within the clump sublayer and 𝑓 c describes the likelihood that disc pebbles become "clumpy" when they reach the mid-plane. These factors will be further discussed in Section 2.3. The other source term represents loss by clumps due to planetesimal formation  \n𝑆 plt = Σ c 𝑡 plts = 𝜖 Σ c 𝑡 sett (5)  \nwhere the settling time-scale is  \n𝑡 sett = 1 St Ω 𝐾 (6)  \nand St is the Stokes number. In the clump pebble region, where the dust-to-gas ratio exceeds unity, planetesimal formation via gravitational collapse is expected (Johansen et al. 2006; Chiang & Youdin 2010; Simon et al. 2016; Schoonenberg et al. 2018). We assume that the clumpy mid-plane dust forms planetesimals on a time-scale 𝑡 plts = 𝑡 sett GLYPH<157> 𝜖 (Drążkowska et al. 2016; Schoonenberg et al. 2018; Stammler et al. 2019), where 𝑡 sett is the timescale on which filaments form via streaming instability in shearing box (Yang et al. 2017). The planetesimal formation efficiency 𝜖 = 0 GLYPH<149> 01, unless otherwise stated, is taken similar to the value argued by Drążkowska et al. (2016) for the typical St GLYPH<24> 0 GLYPH<149> 01 in our simulations. For smaller particles, the time-scales for settling and streaming instability increases (Yang et al. 2017), while the formation rate of planetesimals is lower. After their formation, planetesimals do not interact with the other components.', '2.2 Disc Model: Gas and Disc Pebbles': "We consider gas and dust discs around 1 𝑀 GLYPH<12> mass star. The gas surface density profile is set as  \nΣ g = Σ g GLYPH<150> 0 GLYPH<16> 𝑟 100 au GLYPH<17> GLYPH<0> 1 GLYPH<149> (7)  \nwhere Σ g GLYPH<150> 0 is taken 3 g cm GLYPH<0> 2 in all of our simulations for simplicity, comparable with the gas surface density required to reproduce the DSHARP observations (Zhang et al. 2018). The gas disc in our simulations is taken isothermal in the vertical direction and does not evolve. The aspect ratio of the disc is chosen as  \nℎ = 𝐻 g 𝑟 GLYPH<17> 𝑐 𝑠 Ω 𝐾 𝑟 = ℎ 0 GLYPH<16> 𝑟 100 au GLYPH<17> 0 GLYPH<149> 25 (8)  \nwhere 𝑐 𝑠 = √︁ 𝑘 𝐵 𝑇 GLYPH<157> 𝜇𝑚 𝑝 is the isothermal sound speed with 𝑚 𝑝 the proton mass and 𝜇 = 2 GLYPH<149> 3 the mean molecular weight in atomic units. The aspect ratio at 100 au, ℎ 0 , is taken 0 GLYPH<149> 075 unless otherwise stated. Equation (8) corresponds to a disc temperature profile of  \n𝑇 = 𝜇𝑚 𝑝 𝐺𝑀 ★ ℎ 2 𝑘 𝐵 𝑟 = 𝑇 0 GLYPH<16> 𝑟 100 au GLYPH<17> GLYPH<0> 0 GLYPH<149> 5 (9)  \nMNRAS 000 , 1-18 (2021)  \nwhere 𝑇 0 = 13 GLYPH<149> 9 K is the disc temperature at 100 au for a solar mass star.  \nAt the outer boundary disc pebbles flow in at imposed rates / 𝑀 ext. The initial disc pebble surface density is taken to be the steady-state value  \nΣ d = / 𝑀 ext 2 𝜋𝑟𝑣 dr (10)  \nThe external mass flux / 𝑀 ext is a key model parameter.  \nIn the Epstein drag regime, - particle radius less than mean free path; a condition valid for the disc and dust parameters of interest here - the normalized stopping time of particles with radius 𝑎 GLYPH<15> and internal density 𝜌 GLYPH<15> in the mid-plane of the disc is the Stokes number (Weidenschilling 1977b; Birnstiel et al. 2012)  \nSt = 𝜋 2 𝜌 GLYPH<15> 𝑎 GLYPH<15> Σ g (11)  \nWe consider two setups:  \n(i) runs characterized by a constant Stokes number (the particle size changes with orbital radius);  \n(ii) runs characterized by constant particle size (the Stokes number changes with radius).  \nThe former is consistent with a scenario in which particles attain their size by a velocity fragmentation threshold (Birnstiel et al. 2010); the latter by a bouncing threshold (Güttler et al. 2010). We summarize these parameters in Table 1.  \nUsing the 𝛼 -prescription for vertical transport of the disc pebbles with dimensionless diffusivity parameter 𝛿 d , we set 𝛿 d = 10 GLYPH<0> 3 in all of our simulations as suggested from the vertical extent of dust in edge-on observations of ALMA disks (e.g., Pinte et al. 2016; Villenave et al. 2020). The disk pebbles' vertical diffusivity can be different from their radial diffusivity defined earlier (see Section 5.1 for discussion). Assuming a vertical Gaussian distribution for both disc pebble and gas density, the pebble disc scale height  \n𝐻 d = √︄ 𝛿 d 𝛿 d , St 𝐻 g (12)  \nis obtained by balancing dust settling with vertical diffusion (Dubrulle et al. 1995; Birnstiel et al. 2010).", '2.3 Clump Model': 'Following Sekiya (1998) and Youdin & Shu (2002), we assume that clump pebbles reside in the mid-plane region with the scale height  \n𝐻 c = √︁ 𝑅𝑖 crit 𝜂𝑟 Ψ \' 𝜓 \' (13)  \nhereafter referred to as the clump layer. For the critical Richardson number 𝑅𝑖 crit , Chiang (2008) points out that when the vertically-integrated dust-to-gas ratio is between 1 and 5 times solar abundance, the Richardson number 𝑅𝑖 GLYPH<24> 0 GLYPH<149> 1. Thus we take 𝑅𝑖 crit = 0 GLYPH<149> 1. The dimensionless self-gravitational term Ψ \' 𝜓 \' GLYPH<17> √︁ 1 , 2 𝜓 GLYPH<0> 𝜓 ln »\' 1 , 𝜓 , √︁ 1 , 2 𝜓 \'GLYPH<157> 𝜓 … is an order of unity, where 𝜓 GLYPH<17> 4 𝜋𝐺𝜌 g GLYPH<157> Ω 2 𝐾 = 0 GLYPH<149> 23 \' 𝑟 GLYPH<157> 100 au \' 0 GLYPH<149> 75 is the self-gravity of gas. By equating the diffusion time-scale from clump pebbles to disc pebbles  \n𝑡 diff = 𝐻 2 c 𝛼 c 𝑐 𝑠 𝐻 g (14)  \nand the settling time-scale (Equation (6)), we can calculate the sublayer diffusivity with default values of parameters for clumps. Accordingly, we choose 𝛼 c = 10 GLYPH<0> 5 in the thin mid-plane layer in our default model (see Sect. 5.1 for a more detailed discussion).  \nSpontaneous local overdense clumping of dust by streaming instability was firstly found by Johansen et al. (2006) and Johansen & Youdin (2007). A transition from a homogeneous to a clumpy medium happens at volume dust-to-gas ratio equal to unity (Johansen &Youdin 2007). And strong radial concentration of solids suddenly operates, leading to roughly axisymmetric substructures with vertically stratified dust layer (Yang & Johansen 2014; Li et al. 2018). To qualitatively mimic the appearance of clumps, we introduce the efficiency factor  \n𝑓 c = 1 GLYPH<0> exp \'GLYPH<0> 𝑍 3 mid \' (15)  \nwhere the mid-plane dust-to-gas ratio is  \n𝑍 mid = 𝜌 d GLYPH<150> 0 𝜌 g GLYPH<150> 0 = 𝑍 mid GLYPH<150> dd , 𝑍 mid GLYPH<150> cd = Σ c Σ g 𝐻 g 𝐻 c , Σ d Σ g 𝐻 g 𝐻 d (16)  \nin which 𝜌 d GLYPH<150> 0 and 𝜌 g GLYPH<150> 0 are the mid-plane density of dust and gas respectively. When the dust-to-gas ratio is GLYPH<28> 1, Eq. (15) simply implies that a clumpy mid-plane layer will not form: 𝑓 c GLYPH<28> 1 and there is no disc-to-clump conversion. Otherwise, when 𝑍 mid glyph[greaterorsimilar] 1, a clump layer is present and disc pebbles transform into clump pebbles on a settling time-scale. In order for the ring to survive, the transition of 𝑓 c from 0 to 1 must be rapid; otherwise the dichotomy between a clumpy ring and free-drifting dust will not be established. Exponents larger than 2 in Eq. (15) will not change the conclusions of this work significantly.  \nApart from planetesimal formation, vertical diffusive transport from the clump layer to the disc pebble component also contributes to the loss of clump pebble. However, the rate at which the clump sublayer "erodes" depends on the drain of clumping in clump layer. Pebbles, liberated from dispersing clumps, may yet be prevented from escaping the clump layer, as they can collide with other particles in the clump sublayer (Krĳt & Ciesla 2016).  \nThis suppression of diffusion by the collisional "self-shielding" is modeledthrough a reduction in the diffusivity by a factor 𝑓 𝜏 Eq. (A3):  \n𝑓 𝜏 = 1 GLYPH<0> 𝑒 GLYPH<0> 𝜏 𝜏 GLYPH<149> (17)  \nwhere 𝜏 is the "optical depth" of clump layer. In the optical thin limit ( 𝜏 GLYPH<28> 1) 𝑓 𝜏 \' 1 and there is no reduction: pebbles in the clump layer can freely exchange with disc pebbles through diffusion. In the other limit ( 𝜏 GLYPH<29> 1) 𝑓 𝜏 \' 1 GLYPH<157> 𝜏 , meaning that only a fraction of the clump layer ( Σ c GLYPH<157> 𝜏 ) is able to diffuse out of the clump mid-plane without being impeded. That is, only the surface layers do exchange material. This is analogous to the escape of radiation from layers of optical depth concept below unity.  \nFor the clump depth of the mid-plane layer, one would naively write 𝜏 = 𝜅 GLYPH<15> Σ c where 𝜅 GLYPH<15> the "collision opacity" of the pebbles  \n𝜅 GLYPH<15> GLYPH<17> 4 𝜋𝑎 2 GLYPH<15> 𝑚 GLYPH<15> = 3 𝜌 GLYPH<15> 𝑎 GLYPH<15> GLYPH<149> (18)  \nA complication, however, is the clumpiness of the medium as the medium becomes more transparent when particles clump. To encapsulate these effect, we need to present a macroscopic clump model. In Appendix A we present a simple clump model, accounting for these effects. In this model, it is assumed that the physical properties of clumps (their size, filling factor, etc) are constant, such that their density increases with mass loading ( 𝑍 mid or Σ c). Then, at a critical loading the clumps saturate and no longer contribute to the effective clumping depth of the medium. Correspondingly, we write for the effective clumping depth  \n𝜏 eff = min \' 𝜅 GLYPH<15> Σ c GLYPH<150> 𝜏 max \' (19)  \nand use this expression in Eq. (A3) to evaluate the reduction in diffusivity factor 𝑓 𝜏 . The parameter 𝜏 max signifies at which point clumping starts to make the medium more transparent. When 𝜏 max GLYPH<159> 1 it means that the clumps fill a small fraction of the space and that the particles liberated from dispersing clumps would easily escape, which does not favour ring survival. On the other hand 𝜏 max GLYPH<29> 1 meansthat collisional reduction of the diffusion (Krĳt & Ciesla 2016) would promote the survival of clumping state, which is the desired situation. In this paper we adopt 𝜏 max = 15. Although an (unknown) free parameter, we have verified that the conclusions of this paper are unaffected as long as 𝜏 max GLYPH<29> 1.', '3 RESULTS': 'We conduct two classes of simulations, models with constant Stokes number and models with constant particle sizes. The parameters are summarized in Table 1. Column (1) to (8) are input parameters of the model. Column (9) to (13) are output and corresponding measured values. We list the section in which each model is mainly discussed in Column (14). The grid resolution is fixed at Δ 𝑥 = 0 GLYPH<149> 005 au. For our default model, the outer boundary is at 150 au and the inner boundary is at 50 au. The gas density profile does not evolve in time. The evolution for each dust component is determined by the transport equations described above. The dust outer boundary is set according to a fixed external mass flux / 𝑀 ext. At the inner boundary of the domain, a free-outflow condition is set. We first describe our default model with constant Stokes number in Sect. 3.1. We quantify the mass flux leaking from the ring downstream in Sect. 3.2 and discuss the motion of the ring in Sect. 3.3 and the morphology of the ring in Sect. 3.4. Then we conduct a parameter study about the survival of the ring in Sect. 3.5. In Sect. 3.6, we turn to the constant size model.', '3.1 Constant Stokes Runs': 'In our default model t1f1 , we fix the Stokes number at 0.01, and set the constant external mass flux 100 M GLYPH<8> Myr GLYPH<0> 1 (See Sect. 5.1 for assessment). The planetesimal formation efficiency is 𝜖 = 0 GLYPH<149> 01 such that 1% of the clump collapse into planetesimal in a settling time-scale. Dimensionless diffusivities are 𝛿 d = 10 GLYPH<0> 3 for the disc pebbles and 𝛼 c = 10 GLYPH<0> 5 for the clumps respectively. At 𝑡 = 0, the simulation is initialized by inserting a Gaussian-shaped "ring" of clump material on top of the smooth power-law dust disc. The disc pebble background plus the radial Gaussian distribution of the clump give the total dust density:  \nΣ tot \' 𝑟 \' = Σ d , Σ c = Σ d \' 𝑟 \' , 𝐴 Σ 0 exp GLYPH<18> GLYPH<0> 1 2 \' 𝑟 GLYPH<0> 𝑟 0 \' 2 𝑤 2 GLYPH<19> (20)  \nwhere the reference radius 𝑟 0 = 80 au, and Σ 0 = Σ d \' 𝑟 0 \' (Eq. (10)). Hence, the relative amplitude 𝐴 = 1 means that the peak of total dust density is twice that of the background value. The width of the initial perturbation is taken as 𝑤 = 2 au. The main parameters are listed in Table 1.  \nWe ran the simulation for 5 Myr. We divide the evolution of this ring into two stages: the ring growth phase and the later evolution phase. The upper two panels of Fig. 2 show the radial profile of the total dust surface density and the mass flux during the ring growth phase. Within GLYPH<25> 0.1 Myr, the initial perturbation transforms into a ring-like shape. Since the differential radial drift caused by concentration of clump around the initial perturbation location of 𝑟 = 80 au, a pileup of pebbles happens and make the ring narrow. During GLYPH<24> 0.1-1 Myr, disc pebbles are being "eaten" by the ring, and ring  \nTable 1. Model runs. (1) Stokes number; (2) particle size; (3) initial clump perturbation contrast to background disc pebble density; (4) advection velocity for clumps; (5) external mass flux; (6) planetesimal formation efficiency; (7) clump pebble turbulence 𝛼 -parameter; (8) aspect ratio at 100 au (9) result of rings; (10) predicted velocity of the ring; (11) measured velocity of the ring; (12) predicted leaking mass flux; (13) measured leaking mass flux; (14) section reference where the model is mainly discussed.  \n- 𝑎 "1" for runs where clump advection follows the NSH solution. Empty for runs where there is no clump advection, see Sect. 2.1.\n- 𝑏 constant Stokes number, particle size changes with location\n- 𝑐 fraction of NSH velocity in clump drift, see Eq. (3)\n- 𝑑 constant particle size, Stokes number is variable. 103 𝜇 mcorresponds to St = 0 GLYPH<149> 01 at 𝑟 = 100 au  \ngrows as clump pebbles accumulate locally. The height and width of the ring are still evolving but the location of the ring peak does not change significantly. As the concentration of dust increases inside the ring, the measured mass flux leaking from the ring downstream gradually decreases to approach 65 𝑀 GLYPH<8> Myr GLYPH<0> 1 , which is lower than the external mass flux because of consumption of the ring and planetesimal formation inside. At GLYPH<24> 1 Myr, the ring shape has reached equilibrium reflecting a balance between the consumption of disc pebbles by the ring and the loss due to planetesimal formation.  \nA sufficiently strong seed ring is required to get sustained clumping. Otherwise the initial perturbation will just fade away and all the dust will end up in the disc state. In order to assess this quantitatively, we decrease the initial amplitude 𝐴 and find that the bifurcation point lies around 𝐴 = 0 GLYPH<149> 7 for the parameter setup of t1f1 . Once the initial perturbation\'s amplitude exceeds this value, further increasing 𝐴 will not affect the ring characteristics. (see t1f1A8 , t1f1A6 in Table 1)  \nTo better understand the material exchange, we show in Figure 3  \nthe surface density and radial mass flux contributed by disc pebbles and clump pebbles, respectively, at 𝑡 = 2 Myr. From the upper panel, it is clear that clumps only exist in the ring, where they constitute the vast majority of the ring\'s mass. The high concentration of clumps within the ring region slows down their drift according to Eq. (2); and the faster drifting disc pebbles are "eaten" by the slowly drifting ring. Therefore, both the disc pebble density and pebble flux are lower in the ring than outside the ring. Figure 3 also shows the diffusion suppression factor 𝑓 𝜏 and the settling factor 𝑓 c, which are indicators for the presence of the ring. The settling factor 𝑓 c quantifies the ability of disc dust to transfer into clumpy mid-plane dust, while the diffusion suppression factor 𝑓 𝜏 quantifies the ability of clump material to diffuse back to the disk component. Both of these factors depend on the mid-plane dust-to-gas ratio (Eq. (16); dashed line). Because 𝑍 mid GLYPH<29> 1 withing the ring region, 𝑓 c GLYPH<25> 1 over a broad spatial scale, allowing the incoming disk dust to be incorporated into the clumpy mid-plane. In addition, 𝑓 𝜏 has reached  \nFigure 2. Dust radial profile and mass flux of default model t1f1 for the ring growth phase (upper two panels) and later evolution phase (lower two panels). The Stokes number is constant at St = 0 GLYPH<149> 01. Ring growth happens within GLYPH<25> 1 Myr. After that, the ring\'s shape stays self-similar, while moving outward at a velocity of GLYPH<25> 1 au Myr GLYPH<0> 1 . Other parameters are listed in Table 1  \n<!-- image -->  \nr  \n[au]  \n.  \nits minimum within the disk region, signifying that the clumpy midplane has become "optical thick" to collisions and that exchange to the disk component only happens at the surface, since the liberated pebbles bump into each other and cannot escape the clump layer (see Sect. 2.3).  \nAs a result of these "clumpy effects", pebbles in the disc turn into nearly-stationary clumpy state efficiently and are trapped in clump layer. After the ring\'s growing phase, the combined surface density of pebble and clump dust stops increasing and the morphology of the ring stabilizes. The ring starts to "migrate" outward at speed GLYPH<24> au GLYPH<157> Myr as shown in the lower two panels of Fig. 2. The reason for the outward movement is that clump material is being added to the front side of the ring and lost at the back, whereas the ring as a whole suffers little from inward drifting due to its high concentration. At the outer edge of the ring, diffusion of clump pebbles catches up with the upstream disc pebble flow, whereas flux leaks from the ring downstream. As a result, the ring\'s morphology stays preserved and we see the ring "moving" outward. A more detailed discussion of the ring motion is given in Sect. 3.3.  \nWe show the mass budget evolution of different dust components in Fig. 4. The red and green line show the accumulated inflow and outflow of mass respectively. The blue region is the mass of clump pebbles in the simulation at time 𝑡 , which correspond to the mass of ring. The orange region shows the accumulated mass of formed planetesimal. The inflow slope is constant at 100 𝑀 GLYPH<8> Myr GLYPH<0> 1 consistent  \nFigure 3. Snapshots of run t1f1 at 𝑡 = 2 Myr,  \n<!-- image -->  \nUp: disc pebble (red) and clump (blue) surface density distribution around the ring. Middle: mass flux contributed by the two components. Clumps exist only inside the ring and keep moving inward. The outward motion of the ring is not an overall migration of the dust, but a phase velocity effect. The solid black arrow marks the point where negative mass flux of clump at ring\'s outer boundary is observed due to the steep Σ -gradient. Down: diffusion suppression factor 𝑓 𝜏 (blue), settling factor 𝑓 c (red) and mid-plane dust-togas ratio 𝑍 mid (black). Low 𝑓 c indicates suppression from disc to clump pebble conversion and low 𝑓 𝜏 represent inefficient transport of the clump sublayer.Figure 4. Time evolution of the mass budget. The red and green line show the accumulated inflow and outflow of mass respectively. The slope of the line denoted "inflow" corresponds to the external mass flux. The blue block represents the mass of clump pebbles in the simulation at time 𝑡 and the orange block on top represents the total mass in planetesimals. After 𝑡 = 2 Myr the ring\'s mass (clump pebble mass) has converged to \' 35 𝑀 GLYPH<8> .  \n<!-- image -->  \nwith the external mass flux. Starting with a small amount of mass, the clump pebble mass first increase rapidly during the ring growth phase. In the mean time, the outflow mass and planetesimal mass is still evolving as well. At 𝑡 = 2 Myr, both the outflow and planetesimal formation rate converge to constant values. Together with the constant external mass flux, the mass budget of different components enters equilibrium - the mass of the ring (clump mass) saturates to \' 35 𝑀 GLYPH<8> and the planetesimal production approaches \' 30 𝑀 GLYPH<8> GLYPH<157> Myr. After 5 Myr, the amount of pebbles inside the ring is 41 𝑀 GLYPH<8> and the planetesimal belt is 140 𝑀 GLYPH<8> massive. As shown by the plot, the clumpy ring does not absorb all the incoming mass flux but continue to leak a significant mass to the inner disk.', '3.2 The Leaking Mass Flux': "Here, we use a simple model to quantify the value of the leaking mass flux / 𝑀 leak - the rate at which the ring loses material to the downstream regions.  \nWe consider the inner edge of the ring at the boundary between the clumpy ring and freely-drifting downstream. We assume that near this point, 𝑍 GLYPH<25> 𝑍 crit = 1, most of the clump material gets released. The corresponding clump density is  \n𝜌 c 𝜌 g = Σ c Σ g 𝐻 g 𝐻 c = 𝑍 crit = 1 GLYPH<149> (21)  \nFor simplicity, the mid-plane contributed by disc pebbles is neglected comparing to the clumps'. Further assuming that all of the clumps at the this point will be converted into radially-drifting disc pebbles, the surface densities of disc pebble should equal the released clump pebbles. Defining the leaking mass flux as / 𝑀 leak = 2 𝜋𝑟𝑣 dr GLYPH<150> d Σ d we obtain using the expressions in Eq. (21)  \n/ 𝑀 leak = 2 𝜋𝑟𝑣 dr GLYPH<150> d Σ g 𝐻 c 𝐻 g (22)  \nSince the time-scale of materials exchange is comparable to the drift time-scale pass the the ring, the exact location of the edge is uncertain in our simulations. We take a numerical factor 𝑓 edg to quantify the uncertainty introduced by choosing 𝑍 = 1 for the edge. With further substitutions for 𝑣 dr GLYPH<150> d , 𝐻 c and 𝐻 g, we obtain  \n/ 𝑀 leak ' 𝑓 edg St Σ g GLYPH<150> 0 ℎ 3 0 𝑀 0 GLYPH<149> 5 ★ 𝑟 0 GLYPH<149> 25 = 40 𝑀 GLYPH<8> Myr GLYPH<0> 1 𝑓 edg GLYPH<18> St 10 GLYPH<0> 2 GLYPH<19> GLYPH<18> Σ g GLYPH<150> 0 3 g cm GLYPH<0> 2 GLYPH<19> GLYPH<2> GLYPH<18> ℎ 0 0 GLYPH<149> 075 GLYPH<19> 3 GLYPH<18> 𝑀 ★ 𝑀 GLYPH<12> GLYPH<19> 0 GLYPH<149> 5 GLYPH<16> 𝑟 80 au GLYPH<17> 0 GLYPH<149> 25 (23)  \nIn Table 1, we compute the leaking mass flux predicted from Eq. (23), and measure the exact mass flux downstream for each run as / 𝑀 leak GLYPH<150> ms . Bycomparison Eq. (23) to the measured / 𝑀 leak GLYPH<150> ms we obtain that 𝑓 edg is within » 0 GLYPH<149> 6 GLYPH<150> 0 GLYPH<149> 8 … for runs where clumps have no advection term and within » 1 GLYPH<149> 5 GLYPH<150> 1 GLYPH<149> 6 … for runs with clump drift (see more comperision in Sect. 3.3). The dependence on Stokes number and gas surface density signify that larger particles are easier to leak from the ring. The expression shows that the leaking mass flux is rather sensitive to the aspect ratio (temperature) of the disc. Thicker/hotter discs result in less efficient pileups because smaller mid-plane concentration and faster pebble drifting. In order for rings to survive, the mass flow entering the ring region must exceed / 𝑀 leak . Therefore, when / 𝑀 ext GLYPH<159> / 𝑀 leak the ring will disperse. On the other hand, as the leaking mass flux is fixed, / 𝑀 ext ¡ / 𝑀 leak signifies that there is an excess mass flux, which in the CRM is converted into planetesimals. We see the expression Eq. (23) is independent of external mass flux. Therefore,  \nFigure 5. Run t1f2 . Same as Figure 2 but no clump advection. Ring profile moves faster than that in default model t1f1 . Parameters are listed in Table 1.  \n<!-- image -->  \nonce a ring is present, the mass flux downstream is set and will not be influenced by any fluctuations in / 𝑀 ext upstream.  \nThe analytical values calculated basing on Eq. (23) are consistent with the numerical results in magnitude ( 𝑓 edg ' 1), but the leaking mass flux in runs where clumps don't drift is only half of the cases where clumps with advection. The difference between these two sets tells that clump advection increases / 𝑀 leak by contributing to the leaking. Small 𝛼 c clump will not affect as we expect in Eq. (23), but in run t1f2a4 (see Sect. 3.4), we see the leaking mass flux is larger than t1f2 (see Sect. 3.3) although they are different only on 𝛼 c. This is because clump radial diffusion cannot be ignored with such a high 𝛼 c. Similar to clump advection, radial diffusion of clump promotes the leaking mass flux as well.", '3.3 Ring Migration': 'As can be seen from the mass flux panel of Fig. 3 (middle panel), the clump\'s mass flux is overall positive, indicating that clump pebbles are drifting towards the host star. However, rather than drifting inward, the ring profile actually moves outwards as shown in Fig. 2. This tells us that the motion of the ring is not an overall migration of the dust inside the ring, but a phase velocity of the morphology. In Fig. 3 we find that there is a small negative feature in the mass flux at the outer edge of the ring, which corresponds to the outer edge of the ring. At this point, the radial diffusion of clumps is stronger than its inward drift. These clumps move outward and will grow by feeding from the incoming dust flow. We anticipate the ring to move outward faster when the radial diffusion of the clump is stronger (higher 𝛼 c) or the radial drift of disc pebbles is lower. The former helps with the outward diffusion of clumps at the ring edge. The latter weakens the inward drifting to render the outward moving more significant.  \nIn order to investigate how the clump drift prescription influences the ring\'s movement, we take several runs with different options for the clump\'s advection (see Sect. 2.1). In runs t1f1 , t1f1c9 , t1f1c6 , t1f1c3 and t1f1c0 , we vary 𝑓 ca from 1 to 0 to increase the clump drift velocity, and measure the ring migration velocity if there is ring survive. When 𝑓 ca is lower than 0.3, the initial seed ring does not feed from the external mass flux and decays. For 𝑓 ca = 0 GLYPH<149> 6 the radial drift of the clumps has become GLYPH<24> 5 times that of the NSH solution. A ring forms, but due to the high advective velocity it moves inward rather than outward.  \nNext, we consider the case without any clump advection. Another reason to study this 𝑣 dr GLYPH<150> c = 0 case is that it allows us to study the outward ring motion in isolation. Accordingly, we keep all parameters  \nFigure 6. Snapshots of rings at 𝑡 = 0 GLYPH<149> 5 Myr, varying the clump diffusivity 𝛼 c. The Gray line is the initial pebble distribution for all of the three runs. Higher 𝛼 c broaden the rings and result in stronger outward motion, until at 𝛼 c = 10 GLYPH<0> 3 the ring no longer survives.  \n<!-- image -->  \nFigure 7. Ring survival depends on Stokes number and mass influx rate. From up to down, run t1f3 , t1f4 , t2f1 and t3f1 . In the upper two panels, the external mass flux is lower (50 𝑀 GLYPH<8> Myr GLYPH<0> 1 ) and higher (200 𝑀 GLYPH<8> Myr GLYPH<0> 1 ) compared with the the default model t1f1 , while in the lower two panels the Stokes number is larger (0 GLYPH<149> 02) and smaller (0 GLYPH<149> 005). The green bar\'s upper side indicates the external mass flux, and lower side of green bar indicates the leaking mass flux in the default model. In order to preserve the ring, the external mass flux must be larger than the leaking mass flux. Therefore, larger / 𝑀 ext and lower St aid the long-term survival of the ring.  \n<!-- image -->  \nthe same as in run t1f1 , but switch off the advection term. The results of this simulation ( t1f2 ) is shown in Fig. 5. Similar to run t1f1 , within 1 Myr stage the location of the ring peak did not change but the height grows. On the other hand, after 1 Myr a massive ring forms with the front edge of the ring moving much faster than in run t1f1 . In this case, the pileup of clumps at the outer edge happens more easily, because the clumps which accumulate at the ring region do not drift inward.  \nWe estimate the magnitude of the ring migration as follows. Because of the sharp gradient at the outer edge of the ring, clumps diffuse outward radially. After their diffusion, clumps will be transformed  \nFigure 8. Parameter study about the survival of the ring. Gray : Regions where the initial background dust-to-gas (without a ring) is larger than 0 GLYPH<149> 5 (indicated by the red line). Excluded from our analysis. Blue : Regions where rings form. The colour indicates the ring surface density-to-background ratio (contrast); White : Regions where the initial perturbation has decayed. The green dashed line corresponds to the leaking mass flux for each St. When St = 0 GLYPH<149> 04 and / 𝑀 ext = 400 M GLYPH<8> Myr GLYPH<0> 1 maximum contrast is obtained.  \n<!-- image -->  \nto disk pebbles on a vertical diffusion time-scale 𝑡 diff = 𝐻 2 d GLYPH<157> 𝐷 𝑧 = 𝑡 sett = 1 GLYPH<157> St Ω 𝐾 . At the same time, these clumps at the front edge of ring have spread radially over a region Δ 𝑥 = p 𝐷 c 𝑡 diff before they disperse. The upstream-scattered clumps intercept a fraction 𝑃 of the disc pebble mass flux, where 𝑃 = 𝑓 acc \' Δ 𝑥 GLYPH<157> 𝑣 dr GLYPH<150> d \'GLYPH<157> 𝑡 sett is the ratio of the time pebbles drift over the region Δ 𝑥 and the time it takes them to settle into the clumpy sublayer. Considering the "accretion efficiency", we leave an order of unity factor 𝑓 acc. Therefore, the accumulation of clumpy material at the front edge of the ring progresses at a rate:  \nd 𝑀 d 𝑡 GLYPH<12> GLYPH<12> GLYPH<12> GLYPH<12> front = 𝑃 / 𝑀 ext = 𝑓 acc √︁ St 𝐷 c Ω 𝐾 / 𝑀 ext 𝑣 dr GLYPH<150> d GLYPH<149> (24)  \nEqual the outwards mass flux of Eq. (24) to 2 𝜋𝑟 Σ c 𝑣 ring , we obtain the ring migration velocity  \n𝑣 ring = 𝑓 acc √︁ St 𝐷 c Ω 𝐾 / 𝑀 ext 2 𝜋𝑟𝑣 dr GLYPH<150> d Σ c (25)  \nThe clump surface density at the ring front Σ c can be estimated realizing that the clumping depth of clump 𝜏 c GLYPH<25> 1. This indicates the transition from clumpy ring region to background where there is no clumps. In the Epstein regime (Equation (11))  \n𝜏 c = 𝜅 Σ c GLYPH<24> Σ c 𝜌 GLYPH<15> 𝑎 GLYPH<15> GLYPH<24> Σ c St Σ g (26)  \nTherefore, Σ c GLYPH<24> St Σ g and with further substitutions for 𝐷 c and 𝑣 dr GLYPH<150> d we obtain  \n𝑣 ring \' 20 au Myr GLYPH<0> 1 𝑓 acc GLYPH<18> 𝛼 c 10 GLYPH<0> 5 GLYPH<19> 0 GLYPH<149> 5 GLYPH<18> St 10 GLYPH<0> 2 GLYPH<19> GLYPH<0> 1 GLYPH<149> 5 GLYPH<18> ℎ 0 0 GLYPH<149> 075 GLYPH<19> GLYPH<0> 1 GLYPH<2> GLYPH<18> Σ g GLYPH<150> 0 3g cm GLYPH<0> 2 GLYPH<19> GLYPH<0> 1 GLYPH<18> / 𝑀 ext 100 𝑀 GLYPH<8> Myr GLYPH<0> 1 GLYPH<19> GLYPH<16> 𝑟 100 au GLYPH<17> GLYPH<0> 0 GLYPH<149> 25 GLYPH<149> (27)  \nUsing this equation, we successfully retrieve the migration velocity of the ring in Fig. 5, which is 9 GLYPH<149> 5 au Myr GLYPH<0> 1 at 𝑟 = 90 au, with 𝑓 acc = 0 GLYPH<149> 5. We therefore take 𝑓 acc = 0 GLYPH<149> 5 for all of our simulation and calculate the estimated values of the ring velocity basing on Eq. (27),  \nwhich are listed in Table 1 as 𝑣 ring . In addition, Table 1 lists the measured velocity 𝑣 ring GLYPH<150> ms , which is obtained from the position of the leading edge of clump surface density.  \nComparing the movement of rings in Fig. 2 and Fig. 5, one notes that the two ring velocities are different ( t1f1 : 1 GLYPH<149> 3 au Myr GLYPH<0> 1 , t1f2 : 9 GLYPH<149> 5 au Myr GLYPH<0> 1 ; measured velocity 𝑣 ring GLYPH<150> ms , see Table 1). The ring in the non-clump-advection run t1f2 moves much faster than the ring in the standard model t1f1 . The reason is simply that clumps in t1f1 have a non-zero radial drift velocity 𝑣 dr GLYPH<150> c . Thus, the ring itself drifts inward as well. A more accurate prediction is therefore GLYPH<24> \' 𝑣 ring GLYPH<0> 𝑣 dr GLYPH<150> c \' where 𝑣 dr GLYPH<150> c is the clump drift velocity. Since the clump radial drift depends on dust concentration 𝑍 , which varies with disk radius, a precise analytical expression is hard to obtain, however. But crudely, with 𝑣 dr GLYPH<150> c \' 9 au Myr GLYPH<0> 1 (individual radial drift velocity when 𝑍 mid = 0 GLYPH<149> 5 in Eq. (2)), the net velocity of the ring amounts to \' 1 au GLYPH<157> Myr for run t1f1 which approximately correspond to the measured value. Given the fact that 𝑣 ring and 𝑣 dr GLYPH<150> c are of the same magnitude, the net ring velocity and sign become hard to predict. For example, in run t1f1h8 , with a higher aspect ratio and faster dust drifting, the clump advection dominates over the ring\'s outward movement. The ring moves inward in this case (see t1f1h6 , t1f1h9 , t1f2h6 , t1f2h8 , t1f2h9 for more comparison in Table 1). For the cases where clumps do not drift, 𝑣 dr GLYPH<150> c = 0, the predicted ring velocity 𝑣 ring is more consistent with the measured velocity. And for the runs where the clump drift is faster than the NSH solution ( 𝑓 ca GLYPH<159> 1), the ring\'s movement is dominated by the clump inward advective motion.  \nEquation (27) shows that the ring velocity increases with higher clump diffusivity 𝛼 c and higher external mass flux, but that the ring migrates faster with decreasing Stokes number. A higher 𝛼 c helps with the diffusion of clumps; clumps spread faster at the front edge of the ring and speed up its migration as in run t1f2a4 . However, further increasing the 𝛼 c is detrimental to the ring survival. In run t1f2a3 , 𝛼 c = 10 GLYPH<0> 3 and no ring can survive due to the fast spreading of dust. Besides, a higher external mass flux supplies more materials to feed the ring, which leads to more dust mass pileup and then faster migration as well. The Stokes number influences the migration speed significantly. Larger St means that the clump radial drift is faster, which is bad for dust pileup and also tells that the interaction time-scale available for material exchanging is shorter.', '3.4 Ring Morphology': 'The CRM produces rings that generally are asymmetric: the outer edge of the ring is usually sharper than the inner edge. Several parameters will change the ring morphology, but the most important one is the clump diffusivity 𝛼 c. Figure 6 shows snapshots of runs t1f2a6 , t1f2 , t1f2a4 and t1f2a3 , where 𝛼 c increases from 10 GLYPH<0> 6 (10x less than the default) to 10 GLYPH<0> 3 (100x higher than the default), at 𝑡 = 0 GLYPH<149> 5 Myr. As we discussed in Sect. 3.3, larger 𝛼 c speeds up the outward migration of rings. But it also flattens the ring and makes the overall morphology more Gaussian. In run t1f2a3 , the clump diffusivity is too high ( 𝛼 𝑐 = 10 GLYPH<0> 3 ) that the ring cannot survive. Clump advection and disc pebble diffusivity change the motion of ring but do not affect its shape.', "3.5 Parameter Study of Ring's Survival": "To investigate the survival and migration of the ring, we study how parameters influence its evolution. As discussed in the previous sections, the two most important parameters are the external mass flux  \nand the Stokes number. The external mass flux, which present the amount of material available to feed the ring, decides the disc pebble density according to Equation (10). The Stokes number controls the drifting and settling time-scale, which determines the mass loss from the ring by advection - the leaking mass - and the planetesimal formation rate.  \nIn Fig. 7, we vary the external mass flux. In model t1f3 , the constant external mass flux decreases from 100 M GLYPH<8> Myr GLYPH<0> 1 to 50 M GLYPH<8> Myr GLYPH<0> 1 and in model t1f4 , it increases to 200 M GLYPH<8> Myr GLYPH<0> 1 (see Fig. 7 upper two panels). Since the external mass flux only matters the supplemental mass to the ring, it does not affect the leaking mass. In t1f3 , the external mass flux is smaller than the estimated leaking mass flux and the ring disperses because of the lack of dust supply. In contrast, with more material feeding the ring, the ring grows more massive, becomes wider and migrates faster in t1f4 .  \nIn run t2f1 , St = 0 GLYPH<149> 02 and both drifting and settling is two times faster than in the default run. As shown in the lower two panels of Fig. 7, the ring disappears very quickly in t2f1 with larger St. Larger St make the initial pebble density lower according to Eq. (10) and also raises the leaking mass flux above the external mass flux. Thus there is no ring. Conversely, for the case St = 0 GLYPH<149> 005 (run t3f1 ), we get a more massive ring than t1f1 .  \nIn Fig. 8 we conduct a parameter study varying the external mass flux(10 GLYPH<20> / 𝑀 ext GLYPH<20> 1000 𝑀 GLYPH<8> ) and Stokes number (0 GLYPH<149> 001 GLYPH<20> St GLYPH<20> 0 GLYPH<149> 1). We identify three regions:  \n- (i) Regions where the background mid-plane dust-to-gas (without a ring) is larger than 0 GLYPH<149> 5 (gray). We exclude them from our analysis since the pileup conditions are triggered everywhere in the disc, which indicates an inconsistent initial profile. 1\n- (ii) Regions where rings form. The colour indicates the ring surface density-to-background ratio (contrast) (blue).\n- (iii) Regions where the initial perturbation has decayed within 1 Myr (white).  \nFrom Fig. 8, it is clear that there is a sharp transition between the ring-surviving blue area (2) and ring decaying white area (3). We draw a green dashed line, which is consistent with the leaking mass flux in Eq. (23) with 𝑓 edg = 1 GLYPH<149> 6. Below this criterion, the leaking mass flux dominates the evolution and rings disperse within 1 Myr. Abovethis criterion, the external mass flux is larger than the mass loss due to planetesimal formation and ring leaking at the initial position. The contrast between the ring and background is a key observable. In Fig. 8, the contrast approaches the maximum around St = 0 GLYPH<149> 04 and / 𝑀 ext = 400 M GLYPH<8> Myr GLYPH<0> 1 .  \nThe survival of the rings also depends on the planetesimal formation efficiency 𝜖 , since it influences the mass loss of the ring directly. We vary the 𝜖 from 10 GLYPH<0> 1 to 10 GLYPH<0> 4 in run t1f1 and t1f2 (run t1f1e1 , t1f1e3 , t1f1e4 , t1f2e1 , t1f2e3 , t1f2e4 in Table 1). The results are shown in the two sides of Fig. 9 separately. In the t1f1 group, clumps have an advection velocity and the leaking mass flux is higher than that in the t1f2 group. In other words, as the external mass flux is the same in these two cases, the allowed upper limit for mass loss by planetesimal formation is lower in the t1f1 group than in the t1f2 group. For example, the top panels show that run t1f1 is unable to keep the ring structure with 𝜖 = 0 GLYPH<149> 1 whereas  \n- 1 Reducing the initial pebble surface density would limit the pileup to the (arbitrary) location of the outer boundary, which is also inconsistent. For runs in the gray region, our simulation setup is inapplicable as clumps (ring formation) is not limited to the ring region. Physically, these regions also correspond to extremely high dust-to-gas surface density ratios ( glyph[greaterorsimilar] 0 GLYPH<149> 3).  \nFigure 9. Effects of planetesimal formation rate on the outcome of the ring. From top-down, planetesimal formation efficiency 𝜖 = 10 GLYPH<0> 1 , 10 GLYPH<0> 2 to 10 GLYPH<0> 3 and 10 GLYPH<0> 4 . left: runs with clump advection (group t1f1 ). Higher 𝜖 causes the ring to collapse into planetesimal more rapidly. right: runs without clump advection (group t1f2 ). Higher 𝜖 will significantly weaken the width of the ring.  \n<!-- image -->  \nFigure 10. Effects of disc pebble diffusivity on the outcome of the ring. From top-down, disc pebble radial diffusivity coefficient 𝛼 d increase from 0, 10 GLYPH<0> 5 to 10 GLYPH<0> 4 and 10 GLYPH<0> 3 . left: runs with clump advection (group t1f1 ). Higher 𝛼 d suppress the outward movement of the ring. Ring move inward when 𝛼 d = 10 GLYPH<0> 3 right: runs without clump advection (group t1f2 ). 𝛼 d hardly influence the evolution.  \n<!-- image -->  \nt1f2 can. In both groups, decreasing 𝜖 from 10 GLYPH<0> 2 to 10 GLYPH<0> 4 will not change the ring evolution significantly. This is because for low 𝜖 mass loss is dominated by pebble drift rather than planetesimal formation. To summarize, inefficient planetesimal formation does not guarantee ring survival, but too efficient planetesimal formation ( 𝜖 glyph[greaterorsimilar] 0 GLYPH<149> 1) will certainly destroy the ring.  \nFinally, we consider the effect of the disc pebble diffusivity on the ring velocity. Until now 𝛼 d has been put zero for simplicity, because the transport of disc pebble is dominated by advection. We set up two groups of simulations with 𝛼 d ≠ 0 as shown in Fig. 10. From top to bottom, the values of 𝛼 d increase from 0, 10 GLYPH<0> 5 , 10 GLYPH<0> 4 to 10 GLYPH<0> 3 . As the  \ndisc pebble diffusivity increases, the ring's outwards motion slows down. In the group where clumps can drift, the direction of migration changes sign at 𝛼 d = 10 GLYPH<0> 3 , which is the same as 𝛿 d but much higher than the clump diffusivity 𝛼 c = 10 GLYPH<0> 5 . Similarly, in the no-clump advection group (right panels) the rate of ring migration slows down. This reduction of the outward migration can be understood as follows. According to the CRM in Sect. 3.3, the migration of the ring is a result of disc pebbles interacting with outward-diffusing clumps. As shown in Figure 3, in the CRM, the disc pebble surface density is lower inside the ring than outside because pebbles are getting incorporated into clumps (and have difficulty to escape). Accounting for a disk  \n<!-- image -->  \nFigure 11. Same as Figure 2 but with fixed particle size. Parameters are listed in Table 1 as s1f1 . Clumps that diffuse outward drift back more quickly than t1f1 because of the increasing St at larger radius.  \n<!-- image -->  \ndiffusivity, a diffusive flux directed towards the minimum of Σ d -i.e., the ring location - arises. This effectively increases the radial motion of disc pebbles at the leading edge, and thus decreases the probability of pebbles to be absorbed into the clumps. Therefore, the clump's inward drift velocity 𝑣 dr GLYPH<150> c dominates over the outward 𝑣 ring , resulting in a net inward movement of the ring. Rings still move outward in the case without clump advection, albeit at a slower pace.", '3.6 Constant Size Runs: Cyclical Rings': "Another choice for the pebble properties in the simulations is to fix the size of the pebbles, which is consistent with particle growth limited by bouncing barrier (Güttler et al. 2010; Zsom et al. 2010).  \nIn simulation s1f1 , we keep all parameters identical to the default run t1f1 . The only difference is that this model characterized by a fixed particle radius of 𝑎 GLYPH<15> = 103 𝜇 m, which corresponds to a Stokes number of 0.01 at 80 au. The result of this s1f1 run is shown in Fig. 11. Unlike the constant St case, the front of the ring features no significant outwards movement. We attribute this effect due the increase of the Stokes number with radius: clumps diffuse outward drift back more quickly because of the increasing St at larger radius. In addition, lower St at smaller radius means that the leaking mass flux is also smaller. Both these effects help the stability of the ring.  \nAccording to Eq. (11), the Stokes number decreases as dust drifts from the outer disc towards the host star, which promotes the traffic jam of drifting pebbles at the inner disc (Youdin & Shu 2002; Youdin &Chiang 2004; Drążkowska et al. 2016). With the same parameters as in run s1f1 , we set up run s1f1A0 in which the surface density is a smooth power-law profile at 𝑡 = 0. Even without an initial perturbation ' 𝐴 = 0 ' , we find that a pileup of pebbles happens in the inner region. After their formation, the ring moves outward as shown in the space-time plot of Figure 12.  \nAs was discussed in the motion of ring part (see Section 3.3), the ring's migration is faster for larger St and higher 𝛼 c. With a higher external mass flux of / 𝑀 ext = 400 𝑀 GLYPH<8> Myr GLYPH<0> 1 , larger fixed particle size 𝑎 GLYPH<15> = 1 mm and higher clump diffusivity 𝛼 c = 10 GLYPH<0> 4 , we get a cyclical process run s2f2 with period GLYPH<24> 0 GLYPH<149> 3 Myr (see Fig. 13). As St increases with disc radius, the leaking mass flux also increases. The enhanced particle flux leaking from the primary ring towards the secondary ring can form a significant pileup of pebbles at low St in the inner disc, which supports this cycling. The cycling period get shorter with larger St, higher 𝛼 c and location closer to the star. The  \nFigure 12. Space-time of surface density, illustrating spontaneous ring formation and outward migration (run s1f1A0 ). Particle size is held fixed in this simulation and there is no initial perturbation. Pileup of ring happens at inner disc region around 10 au and the ring moves outward to 60 au.Figure 13. Cycling of multiple rings (run s2f2 ). The evolution of the ring is periodic. Ring is formed by pileup at inner disc and move outward.  \n<!-- image -->  \nphysical reason behind this cyclical process can be understood as follows. According to the definition of Stokes number in Eq. (11), St increases with disc radius along the ring moving outward. Thus the leaking mass flux also increases. An enhanced particle flux is released from the primary ring towards the secondary ring when the primary ring decays, which can assist the secondary ring's external mass flux as well. Consequently, a pileup of material naturally happens in inner disc again, which supports this cycling.", '4 COMPARISON WITH ALMA': "In order to investigate the spatial dust distribution responsible for the ring emission, we generate a synthetic intensity map corresponding to the CRM output. The procedure is the following:  \nStarting from the luminosity of the host star 𝐿 ★ , the dust temperature 𝑇 𝑑 can be estimated by the irradiated flaring disc model  \n𝑇 d ' 𝑟 ' = GLYPH<18> 𝜙𝐿 ★ 8 𝜋𝑟 2 𝜎 SB GLYPH<19> 0 GLYPH<149> 25 (28)  \nFigure 14. Application of the CRM to Elias 24. Panels: (a) Convolved image resulting from our simulation with constant particle size, used 2D Gaussian beam is plotted at bottom left. (b) Observed continuum emission at 239 GHz, beam sizes is shown in the lower left. (Andrews et al. 2018). (c) Total dust surface density (green) and scaled gas density profile (gray). (d) Azimuthally averaged intensity profile from the observation (blue) and intensity profile of the our model (orange).  \n<!-- image -->  \nwhere 𝜎 SB is the Stefan-Boltzmann constant and 𝜙 = 0 GLYPH<149> 02 is the flaring angle (Chiang & Goldreich 1997; D'Alessio et al. 2001). The aspect ratio of the disc is therefore (Eq. (8) and Eq. (9)):  \nℎ = √︄ 𝑘 𝐵 𝑇 𝑑 𝑟 𝜇𝑚 𝑝 𝐺𝑀 ★ (29)  \nwhich is a function of host stellar mass 𝑀 ★ , stellar luminosity 𝐿 ★ , and radius 𝑟 .  \nThen we apply the CRM with constant particle size with parameters given in Sect. 4.1 and Sect. 4.2 and obtain the surface density profile Σ d ' 𝑟 ' and Σ c ' 𝑟 ' , from which we obtain the vertical distribution assumption  \n𝜌 d ' 𝑧GLYPH<150> 𝑟 ' = Σ d p 2 𝜋𝐻 d exp GLYPH<0> 1 2 𝑧 2 𝐻 2 d ! 𝜌 c ' 𝑧GLYPH<150> 𝑟 ' = Σ c p 2 𝜋𝐻 c exp GLYPH<0> 1 2 𝑧 2 𝐻 2 c ! (30)  \nWe rotate and incline the disc basing on the position angle and inclination from Andrews et al. (2018) and calculate the intensity.  \nFor simplify, we assume that the opacity of all dust in our simulation is a constant 𝜅 𝜈 GLYPH<24> 1 GLYPH<149> 0 cm 2 g GLYPH<0> 1 , corresponding to a grain radius of 𝑎 GLYPH<24> 100 𝜇 m at ALMA band 6 ( 𝜆 = 0 GLYPH<149> 125cm) according to the DSHARP opacity model (Birnstiel et al. 2018). Without accounting for scattering, the intensity profile 𝐼 𝜈 ' 𝑟 ' by radiative transfer is  \nd 𝐼 𝜈 ' 𝑟 ' d 𝜏 = GLYPH<0> 𝐼 𝜈 , 𝐵 𝜈 ' 𝑇 d ' 𝑑𝜏 = 𝜅 ' 𝜌 d , 𝜌 c ' d 𝑧 (31)  \nwhere 𝐵 𝜈 ' 𝑇 𝑑 ' is the Planck function. Then, we convolve these intensity maps with 2D Gaussian FWHM 𝛿 GLYPH<24> 0 GLYPH<149> 00 035 ( GLYPH<24> 5 au), which correspond to a beam size and distance from the DSHARP data (Andrews et al. 2018). The used parameters of host star are listed in Table 2.  \nWe apply the CRM with constant particle size towards two protoplanet discs, Elias 24 and AS 209 (Fig. 14 and Fig. 15). The two  \n2  \nFigure 15. Application of the CRM to AS 209. Panels same as Fig. 14  \n<!-- image -->  \nTable 2. Stellar parameters and image setups used for simulations (Andrews et al. 2018). (1) distance in parsec; (2) (3) stellar luminosity and mass in units of the solar values; (4) beam size used for 2D Gaussian kernel convolution; (5) opacity of dust.  \nsample discs are selected from ALMA DSHARP survey (Andrews et al. 2018; Huang et al. 2018; Guzmán et al. 2018). Both of them have well-resolved annular substructures at 239 GHz dust continuum and they are one of the sharpest rings who hold substantial peak-to-valley contrast among all DSHARP sources (Dullemond et al. 2018). We briefly discuss how well the CRM reproduces these substructures.", '4.1 Single ring: Elias 24': "First, we show a case featuring a single ring. Elias 24 is a system with a narrow bright ring at 77 au and a deep gap at 57 au (Huang et al. 2018). We choose an external mass flux of / 𝑀 ext = 200 𝑀 GLYPH<8> Myr GLYPH<0> 1 , a clump diffusivity 𝛼 c = 3 GLYPH<2> 10 GLYPH<0> 5 and a particle size of 93 GLYPH<149> 4 𝜇 m. According to the temperature profile, the aspect ratio ℎ 0 = 0 GLYPH<149> 1 at 100 au. The gas surface density is taken to be Σ g = 3gcm GLYPH<0> 2 ' 𝑟 GLYPH<157> 100 au ' GLYPH<0> 1 , hence St = 8 GLYPH<149> 66 GLYPH<2> 10 GLYPH<0> 3 at the ring location ( 𝑟 = 77 au). We initialize a small perturbation at 𝑟 = 73 au with relative amplitude 𝐴 = 0 GLYPH<149> 5 (see Eq. (20)). By 0 GLYPH<149> 45 Myr the seed ring has transformed into the outward moving state. The ring peak has arrived at the observed location of the ring at submillimeter wavelengths, where we take a snapshot. In order to scale the intensity, the opacity is chosen to be 𝜅 = 1 GLYPH<149> 46 cm 2 g GLYPH<0> 1 .  \nFigure 14a shows the convolved continuum map and Figure 14b shows the observational dust continuum at ALMA band 6 (239 GHz) and the radial profile, respectively. Panel (c) shows the unconvolved surface density profile corresponding to the simulations and the gas surface density scaled by a factor 0 GLYPH<149> 1. In panel (d), the orange line shows the azimuthally averaged intensity from the convolved continuum map (a) in orange while the blue line shows the radial profile inferred from the observations (Huang et al. 2018). The continuum  \nemission in panel (b) shows a perfect annular structure and high contrast between the ring to the background. After convolution, the radial profile turns to be smoother and more Gaussian, which make the edges of the ring not as sharp as in the surface density profile. This simulation reproduces not only the positions and intensities but also the shape of the ring and gap well.", '4.2 Two rings: AS 209': "Most DSHARP disks exhibit multiple annular substructures (Huang et al. 2018). AS 209 is such a system hosting a group of ring (Guzmán et al. 2018). For simplicity, we focus the CRM towards the two most prominent, well-separated outer rings centered at 74 au and 120 au. We set up the simulation with parameters / 𝑀 ext = 300 𝑀 GLYPH<8> Myr GLYPH<0> 1 , clump diffusivity 𝛼 c = 10 GLYPH<0> 5 and particle size 180 𝜇𝑚 . The aspect ratio ℎ = 0 GLYPH<149> 08 at 100 au correspond to the stellar luminosity of 0 GLYPH<149> 83 𝐿 GLYPH<12> (Eq. (29)). The gas surface density in power-law profile is Σ g = 3 g cm GLYPH<0> 2 ' 𝑟 GLYPH<157> 100 au ' GLYPH<0> 1 , hence St = 3 GLYPH<149> 75 GLYPH<2> 10 GLYPH<0> 2 at the primary ring location 𝑟 = 120 au and St = 2 GLYPH<149> 3 GLYPH<2> 10 GLYPH<0> 2 at the secondary ring location 𝑟 = 74 au. We initiate two small perturbation with relative amplitude 𝐴 = 0 GLYPH<149> 2 at 𝑟 = 69 au and with relative amplitude 𝐴 = 1 at 𝑟 = 113 au. After 0 GLYPH<149> 8 Myr, we take a snapshot when both seed rings have grown into a steady morphology. A dust opacity 𝜅 = 1 GLYPH<149> 75 cm 2 g GLYPH<0> 1 is used to scale the intensity.  \nThe observed dust continuum and radial profile of AS 209 are shown in Fig. 15 panels (b) and the blue line in (d). Panels (a) and the orange line in (d) are the profiles from the simulation above. The positions, intensities and profile of the rings and gaps are reproduced. The intensity contrast between the two ring is reproduced as well.  \nThe moving ring in the our models show a non-symmetric profile in radial direction, that is, the inner back edge of the ring is always wider and smoother than its outer front edge. New pileup happens at the front rapidly and its motion is opposite to the mass flux direction. Due to the current resolution of observation and the influence of data convolution, the ring is not well enough resolved to estimate the detailed morphology of the ring. But we expect such a shape should be observed in future observation.  \nFinally, we briefly discuss how the image changes with the further evolution of the ring. In general, since the simulations have reached a quasi-steady state, the key features in the synthetic image will only slowly change with time. The width of the ring will become a bit wider as the outer edge moves faster than the inner edge of the ring (see Fig. 5). After 1 Myr, for example, the rings in AS 209 will migrate out by a further 10 au, but the intensity and width change by only GLYPH<24> 10%. However, it is likely that other model parameters, such as the external mass flux, will change as well on evolutionary timescales, in turn affecting the ring observational characteristics. We discuss how the non-constant external mass flux will influence the results in Section 5.1.", '5.1 Assessment': 'Our scenario to maintain rings on evolutionary time-scales ( GLYPH<24> Myr) without a pressure bump relies on the following assumptions: a high mass flux of pebbles ( GLYPH<24> 100 𝑀 GLYPH<8> Myr GLYPH<0> 1 ), a relatively low but nonzero clump diffusivity ( 𝛼 c GLYPH<24> 10 GLYPH<0> 5 ) and a low planetesimal formation efficiency ( 𝜖 glyph[lessorsimilar] 0 GLYPH<149> 01). We discuss these parameters and compare them with previous work.  \nFirst, the CRM requires a relatively high mass budget for pebbles.  \nThe value of the pebble flux in the CRM ( / 𝑀 ext glyph[greaterorsimilar] 100 𝑀 GLYPH<8> Myr GLYPH<0> 1 ) suggests that discs have reservoirs of at least GLYPH<24> 100 𝑀 GLYPH<8> in pebbles. These values may seem high compared to typical Class II discs (Andrews et al. 2013; Ansdell et al. 2016; Barenfeld et al. 2016; Ansdell et al. 2017), but may be more common for the DSHARP sample. In addition, the disc dust mass inferred from radio continuum emission could be underestimated because of dust scattering (Zhu et al. 2019). Also, it has been suggested that the disk gas mass and by implication the dust mass that "hides" in the outer disk is higher than what is inferred from CO observations (Powell et al. 2017, 2019). Planet formation models require mass reservoirs of GLYPH<24> 100 𝑀 GLYPH<8> , as the efficiency of pebble accretion, for example, is low (Liu & Ormel 2018; Ormel & Liu 2018). For example, for the formation of the cores of the solar system\'s giant planets (Lambrechts &Johansen 2014) and super-Earth (Lambrechts et al. 2019), typical pebble budgets amount to 130 and 190 Earth mass respectively. Thus, it is not unreasonable to anticipate such high mass budgets.  \nAnother simplifying assumption is that the external mass flux is constant. In reality, the mass flux may both fluctuate and and decay with time. The ring will disappear when the mass flux drops below the leaking mass flux / 𝑀 leak (see Sect. 3.2), as there is not enough material to support the ring any longer. Thus, the lifetime of the ring is determined by the duration over which the external mass flux exceeds / 𝑀 leak . However, the ring is insensitive to small fluctuations in / 𝑀 ext. Even if there was zero external mass flux, rings still takes a time GLYPH<24> 𝑀 ring GLYPH<157> / 𝑀 leak GLYPH<25> 0 GLYPH<149> 1-1 Myr to disappear. In addition, the sudden release of pebbles may locally enhance the mass flux downstream to create or supply inner rings, see discussion in Sect. 5.2.  \nA novel feature of the CRM is the hypothesis that the rings are clumpy in nature, and that these clumps are characterized by a nonzero diffusivity. In (magneto)hydrodynamical simulations of streaming instability a clumpy distribution of solids naturally happen whenever 𝜌 d GLYPH<157> 𝜌 g GLYPH<21> 1 (e.g., Johansen et al. 2006; Johansen & Youdin 2007; Bai & Stone 2010). Recently, Yang et al. (2017) and Li et al. (2018) have measured the diffusivity of the particles undergoing the streaming stability. These characteristics - a clumpy, diffusive medium motivate the non-zero 𝛼 c we used. Following the CRM about the clump scale height in Sect. 2.3, we express this diffusivity in term of an 𝛼 -viscosity ( 𝐷 c = 𝛼 c 𝐻 2 g Ω 𝐾 ). Equating the diffusion time-scale of clumps 𝑡 diff = 𝐻 2 c GLYPH<157> 𝐷 c with the settling time-scale Eq. (6), we obtain  \n𝛿 c = St 𝐻 2 c 𝐻 2 g \' 6 GLYPH<2> 10 GLYPH<0> 6 GLYPH<18> St 10 GLYPH<0> 2 GLYPH<19> GLYPH<18> ℎ 0 0 GLYPH<149> 075 GLYPH<19> 2 (32)  \n(cf. Youdin & Chiang 2004; Takeuchi et al. 2012) where we, for simplicity, have assumed that the diffusion is isotropic 𝛼 c = 𝛿 c. For the CRM to produce rings, however, the value of 𝛼 c cannot be too high: when it were to exceed the disc pebble vertical diffusivity 𝛿 d = 10 GLYPH<0> 3 , the ring would simply spread out, rendering the midplane dust-to-gas ratio to fall below unity. On the other hand, if the clump radial diffusivity is lower than 𝛼 c = 10 GLYPH<0> 5 , the outward ring migration proceeds slower, but it does not affect the survival of the ring.  \nPlanetesimal formation is another mechanism that influences the survival of the ring (Stammler et al. 2019). It is expected that streaming instability and planetesimal formation happen when the local volume dust-to-gas ratio ¡ 1 (Umurhan et al. 2020) after which we expect the clumps to exist for a settling time-scale (Shariff & Cuzzi 2015). Hence 𝜖 = 1 implies optimal conversion of (clumpy) pebbles to planetesimals, in which case there would be no ring. Generally, fast planetesimal formation ( 𝜖 GLYPH<24> 0 GLYPH<149> 1) is problematic for every planetesimal-forming ring model. For example, in Drążkowska et al.  \n(2016) 𝜖 \' 0 GLYPH<149> 001 (more precisely, 0 GLYPH<149> 01 GLYPH<157> 2 𝜋 ) is taken as a standard value for St = 0 GLYPH<149> 01 particles, smaller than our standard value of 𝜖 = 0 GLYPH<149> 01, which corresponds to a time-scale glyph[greaterorsimilar] 1 Myr for planetesimal formation from clumps. Nevertheless, high 𝜖 seems in line with hydrodynamical simulations, which typically feature rapid planetesimal formation (e.g., Abod et al. 2019). However, the ensuing collapse of a bound clump into much more compact planetesimals is under active study with some studies pointing out that turbulent diffusion may suppress the collapse of low-mass clumps (Klahr & Schreiber 2020).  \nOf course, the simplest way to make massive rings for the ALMA rings is to avoid the conditions for planetesimal formation (streaming instability), i.e., by having a low mid-plane dust-to-gas ratio. But this requires a more turbulent environment with rings likely to become vertically thick (see discussion in Sect. 5.2).', '5.2 Comparison to Pressure Bump scenario': 'In the pressure bump scenario, the width of the ring is determined by the radial diffusivity and the value of the pressure gradient reversal (Dullemond et al. 2018). Recently, Rosotti et al. (2020) have combined kinematic data (the rotation profile) with the measured ring widths, deriving an aerodynamic constraint of 𝛼 GLYPH<157> St GLYPH<24> 0 GLYPH<149> 1 for the ring particles, which suggests that the dust radial diffusivity is at least moderate. On the other hand, geometric reasoning suggests a low value of vertical 𝛿 GLYPH<159> 10 GLYPH<0> 3 as deduced from the pebble scale height GLYPH<24> 1 au at 𝑟 = 100 au (e.g., Pinte et al. 2016; Villenave et al. 2020). Several suggestions have been proposed to fix this disagreement. One solution mentioned by Rosotti et al. (2020) is that the particles are very small (St glyph[lessorsimilar] 10 GLYPH<0> 3 ). Another is to invoke that the diffusivity in the radial and vertical direction is different (Zhu et al. 2015a; Xu et al. 2017; Yang et al. 2018; Baehr & Zhu 2021). the CRMdoes not suffer from this tension, as clumps are located within a thin vertical scale height "by definition" (Sekiya 1998) and the radial width of rings is independent of 𝛼 .  \nWhile the CRM maintains a long-lived ring, it is unclear whether the mechanisms responsible for the pressure bump ring can be sustained on GLYPH<24> Myr time-scales. Taki et al. (2016) point out that the outward flow of gas created by the back-reaction from dust to gas will flatten the gas pressure bump. Therefore, a mechanism supporting the bump must be present. One such mechanism is the torques exerted on the disc by planets (Lin & Papaloizou 1986). Previous work by Zhang et al. (2018) suggests that giant planets in the disc can explain annular substructures in some of DSHARP discs (Andrews et al. 2018). However, in that case one expects a kinematic imprint - velocity kinks in the rotation curves (Teague et al. 2018; Pinte et al. 2020). Recent efforts to obtain kinematic imprints of embedded protoplanet only found convincing detections in 8 out of 18 DSHARP sources (Pinte et al. 2020). In particular, there is no obvious CO emission kink in AS 209, despite the high data quality for this source (Pinte et al. 2020). Therefore, not every ring is necessarily associated with a planet, justifying exploring alternative mechanisms.  \nA common feature in the DSHARP observations is the presence of multiple rings (e.g. AS 209, see Sect. 4.2). Different from the pressure bump scenario, which traps essentially all of the incoming pebbles, the inward pebble flux is only partially halted by the ring in the CRM. The leaking of inward-drifting pebbles can explain such a multi-ring feature naturally, which is a particularly favored outcome in the constant size model, because the mass flux required to support the ring increases with disc radius (see Sect. 3.2). In other words, although the upstream flux is partly reduced by the outer ring, the  \nremaining flux is still large enough to maintain a second ring or multiple rings downstream.  \nWediscuss several ways to distinguish our scenario from the pressure bump scenario observationally. First, because the CRM does not rely on a pressure maximum, we do not expect ALMA rings to be accompanied with concomitant substructures in gas tracers (such as COisotopologues). Some early studies to link dust substructure with gas observations are now being made, which is a promising avenue to differentiate between these two scenarios (e.g, Kim et al. 2020; Rosotti et al. 2021). However, interpretation of gas substructure is hard because of the relatively low spatial resolution and sensitivity of gas emission lines and may be further affected by radiation transfer effects. Second, the size distribution downstream of the ALMA ring is expected to be different from the distribution inside the ring. In the pressure-bump induced ring scenario, large grains are more efficiently trapped by the pressure bump, whereas the small grains filter through (Pinilla et al. 2015; Powell et al. 2019). In contrast, in the CRM any particle can leak from the ring regardless of size. Therefore, measuring the change in particle size between the ring location and the region downstream of it is another way to differentiate between the two scenarios. Third, the morphology of the ring in our simulation displays a typical asymmetric profile - i.e., the leading (outer) edge is steeper than the trailing (inner) edge (see Sect. 3.4). Higher spatial resolution observation are needed to resolve the ring morphology. For example, the steep "cliffs" edge in AS 209 (Fig. 15c) would show up much stronger when the spatial resolution would improve to GLYPH<24> 0 GLYPH<149> 00 01.  \nFinally, it is possible for the two ring mechanism to work in tandem. Apressure bump is certainly a compelling way to concentrate solids and to quickly build rings, particularly in the outer disk. However, pressure bumps may not last for evolutionary time-scales - planets, for example migrate. The CRM then addresses how rings survive in the long run.', '5.3 Implications for planet formation': 'The outward migration of rings seen in most model setups carries important implications for planet formation and debris discs. Strictly speaking, the outward movement of the ring we discuss in Sect. 3.3 is only the migration of the leading edge. This outward motion means that the leading edge may not be greatly influenced by (planet formation) processes happening downstream. Pebbles released by the "propagating" ring (as discussed in Sect. 3.2) can aid planet formation in the inner disk. For example, rings could start in the inner disk regions ( GLYPH<24> 10 au) but then move outwards ( GLYPH<24> 50 au) in run s1f1A0 (Fig. 12). Then in the early stage, the massive ring can help with the formation of planet embryos, while in the later phases it still leaks pebbles to feed the possible inner planet. Therfore, in our scenario outer rings are beneficial towards forming interior planets.  \nOur outward propagating pebble ring leaves behind a planetesimal belt up to GLYPH<24> 100 Earth mass at a surface density GLYPH<24> 1 g cm GLYPH<0> 2 . Although these are massive belts, the collision time-scale of planetesimals at 𝑟 GLYPH<25> 100 au, 𝑡 coll GLYPH<24> 𝜌 GLYPH<15> 𝑅 plts GLYPH<157> ΣΩ , is still long especially when planetesimals are "born big" ( 𝑅 plts GLYPH<24> 100 km; Morbidelli & Nesvorný 2020). In other words, we expect migrating rings in the outer disk ( glyph[greaterorsimilar] 50 au) to leave behind massive planetesimal belts. Upon dispersal of the disk these planetesimal belts may be stirred to appear as debris disks. Intriguingly, it has recently been suggested that certain debris disks could indeed be very massive. By fitting the observed luminosity to state-of-the-art collision evolution models, one deduces that the primordial planetesimal belts could be as massive as 1000 𝑀 GLYPH<8> (Krivov et al. 2018; Krivov & Wyatt 2021)! Although these numbers  \nare associated with stars more massive stars than the DSHARP sample, rendering a 1-1 comparison somewhat problematic, the question is raised which mechanism could have created such massive belts. In our view, ALMA rings could well be cradle for such massive planetesimal belts.  \nFinally, we briefly discuss the potential of planet formation by pebble accretion in the ring region. In the standard pressure supported disk without pebble clumping, the headwind velocity always contributes to the relative velocity between pebble flow and the planetesimal. Pebble accretion would be ineffective at distances GLYPH<24> 50-100 au, unless the bodies have grown to masses similar to Pluto (Visser & Ormel 2016). However, in the CRM the pebble drift is significantly decelerated because of the high dust-to-gas ratios. Therefore, the pebble accretion threshold will be much lower inside the clumpy ring. Thus, in the interval where clump pebbles and planetesimals co-exist, (small) planetesimals could efficiently feed from pebbles. On the other hand, mutual scattering of planetesimals would oppose growth (Liu et al. 2019; Schoonenberg et al. 2019). In a future work, we intend to investigate the potential for planet formation in these ALMA rings.', '6 CONCLUSIONS': 'In this work, we have developed a model to understand the axisymmetric substructure seen by ALMA in protoplanetary disks rings - without relying on a pressure maximum. This clumpy ring model (CRM) assumes that rings are the manifestation of a dense and clumpy medium, which is actively forming planetesimals. We have conducted a series of 1D transport simulations where pebbles are being characterized either by a fixed Stokes number or by a fixed particle size. In these, we have investigated how key parameters as the pebble size, the external mass flux, and planetesimal formation efficiency affect the survival of the ring. Our main findings are:  \n- (i) The CRM relies on an initial non-linear perturbation or order unity above the background level. Once it is triggered and the conditions for ring survival are met, a (quasi)-steady state emerges where the ring can exist on evolutionary time-scale.\n- (ii) Without containment (from a pressure maximum), rings leak material at a rate given by Eq. (23). Ring survival demands the existence of an external flux of disc pebbles, exceeding this number. In addition, the time-scale on which planetesimals form must be much longer than the settling time-scale. Once these conditions are met, a ring can be maintained for GLYPH<24> Myr. A higher external mass flux, smaller particle size, and lower aspect ratio facilitate the long-time survival of rings.\n- (iii) Aplanetesimal formation efficiency lower than the default of 𝜖 = 10 GLYPH<0> 2 does not affect the survival of the ring. On the other hand, efficiencies 𝜖 glyph[greaterorsimilar] 0 GLYPH<149> 1 either destroy or else significantly diminish the extent of the ring.\n- (iv) Due to the diffusive nature of clumpy material, rings move outwards at rates up to GLYPH<24> 10 au Myr GLYPH<0> 1 , while material leaks away from the ring downstream. Outward ring movement is favoured by a higher clump diffusivity (but no more than 𝑎 c = 10 GLYPH<0> 3 ), small clump advection velocity, a larger external mass flux and a lower disc pebble Stokes number.\n- (v) Compared to the pressure bump ring formation model, our rings stand out for their sharp edges - a distinction between the CRM and the pressure bump scenario, which can potentially be observed at increased spatial resolution.\n- (vi) Mock images of DSHARP sample Elias 24 and AS 209 are  \nmade, which reproduces not only the positions and intensities but also the shape of the ring and gap well.  \n(vii) The high leaking pebble mass flux implies that in the CRM ALMArings can aid planet formation interior to the observed rings.  \n(viii) When the particle size is constant throughout the disk, a pileup in the inner disk naturally occurs, which spontaneously triggers the CRM. A cyclical process emerges where rings appear in the inner disk ( GLYPH<24> 10 au), move outwards over several tens of au before dissolving entirely to re-appear in the inner disk. A constant pebble size also naturally gives rise to multiple rings.  \n(ix) Up to GLYPH<24> 100 𝑀 GLYPH<8> in planetesimals are expected to form over a period of GLYPH<24> Myr. In the outer disk, this could be a possible explanation for the high mass budget of planetesimals inferred in debris discs.  \nTofurther investigate the viability of the CRM to reproduce ALMA morphology, we plan to develop a more detailed fitting procedure, accounting for the dust size distribution, self-consistent temperature profile and more precise radiation transfer.', 'ACKNOWLEDGEMENTS': 'We thank the anonymous referee for their thoughtful and constructive comments. H.J. and C.W.O. would like to thank Ruobing Dong, AkimasaKataoka, Seongjoong Kim, Sebastiaan Krĳt and Chao-Chin Yang for constructive suggestions and helpful discussions during this work. The authors also appreciate feedback on an earlier version of the manuscript by Joanna Drążkowska, Beibei Liu, Satoshi Okuzumi and Zhaohuan Zhu. This work has used Astropy (Astropy Collaboration et al. 2013), Matplotlib (Hunter 2007), Numpy (Harris et al. 2020), Scipy (Virtanen et al. 2020) software packages.', 'DATA AVAILABILITY': 'The data underlying this article will be shared on reasonable requests to the corresponding author.', 'APPENDIX A: CLUMP MICROPHYSICAL MODEL': 'In analogy to radiation transport, the variation of vertical mass flux 𝐹 with height d 𝑧 can be simplified as two parts: the loss term GLYPH<0> 𝜅𝜌 c 𝐹 and source term 𝜌 c GLYPH<157> 𝑡 diff . The first term describes the flux reduction caused by self-shielding and the second term is the new contribution due to diffusion of the local layer d 𝑧 :  \nd 𝐹 d 𝑧 = GLYPH<0> 𝜅𝜌 c 𝐹 , 𝜌 c 𝑡 diff GLYPH<149> (A1)  \nIntegrating Eq. (A1) over the vertical dimension, we obtain  \nd Σ c d 𝑡 GLYPH<17> 𝐹 = 1 GLYPH<0> 𝑒 GLYPH<0> 𝜏 𝜅 1 𝑡 diff = 1 GLYPH<0> 𝑒 GLYPH<0> 𝜏 𝜏 Σ c 𝑡 diff (A2)  \nwhere  \n𝜏  \n=  \n∫  \n𝜅𝜌  \ncd  \n𝑧  \n=  \n𝜅  \nΣ  \nc. Accordingly, we define  \n𝑓 𝜏 = 1 GLYPH<0> 𝑒 GLYPH<0> 𝜏 𝜏 (A3)  \nas the diffusion suppression factor in Eq. (4). It describes how efficiently clumps can diffuse back to disc pebbles. When the "clumping depth" (optical depth for collisions) of the medium is low, replenishment of disc pebble from clumps is unopposed, this factor attains its maximum value 𝑓 𝜏 \' 1. On the other hand, for 𝜏 GLYPH<29> 1, 𝑓 𝜏 GLYPH<28> 1: particles in dispersing clumps do not leave the clump layer, but \'bump\' into other clumps.  \nIn the preceding we have not yet specified the opacity 𝜅 of the clumpy medium, which is an average over the inhomogeneous space and can therefore be lower than the geometrical opacity 𝜅 GLYPH<15> of the particles. For monodisperse-sized spherical particles with size 𝑎 GLYPH<15> , the clumping depth is defined as the collisional cross section 4 𝜋𝑎 2 GLYPH<15> over the single particle mass 𝑚 GLYPH<15> Eq. (18):  \n𝜅 GLYPH<15> GLYPH<17> 4 𝜋𝑎 2 GLYPH<15> 𝑚 GLYPH<15> = 3 𝜌 GLYPH<15> 𝑎 GLYPH<15> (A4)  \nIn terms of the Stokes number (Equation (11)) this reads  \n𝜅 GLYPH<15> = 3 𝜋 2 1 St Σ 𝑔 \' 157 cm 2 g GLYPH<0> 1 GLYPH<18> St 0 GLYPH<149> 01 GLYPH<19> GLYPH<0> 1 GLYPH<18> Σ g 3 g cm GLYPH<0> 2 GLYPH<19> GLYPH<0> 1 GLYPH<149> (A5)  \nLet clumps be characterized by a radius 𝑎 cl , volume density 𝜌 cl and let the filling factor 𝑓 fill be the fraction of clumps in the medium. Then  \n𝜌 cl = 1 𝑓 fill 𝜌 c = 𝑍 mid 𝑓 fill 𝜌 g (A6)  \nWe can calculate the clumping depth for a single clump 𝜏 cl = 𝜅 GLYPH<15> 𝜌 cl 𝑎 cl . With further substitution of 𝜅 GLYPH<15> and 𝜌 cl and assuming a vertical Gaussian distribution for gas density Σ g = p 2 𝜋𝜌 g 𝐻 g, we hence get that at a critical mid-plane dust-to-gas ratio  \n𝑍 GLYPH<3> mid = 2 3 √︂ 2 𝜋 St 𝑓 fill 𝐻 g 𝑎 cl (A7)  \nthe clumps reach 𝜏 cl = 1. Below this value, dust within an individual clump can easily escape and the transport of dust performs the same as single-particle disc pebbles. Thus we could keep the opacity of single dust particles as the medium opacity. On the other hand, for collisionally thick clumps, clumps act as the target particles. In that situation, denser clumps (higher 𝑍 ) no longer increase the collision rate. We therefore write for the net opacity:  \n𝜅 \' 𝑍 mid \' = 8 > > < > > 3 4 1 𝜌 GLYPH<15> 𝑎 GLYPH<15> = 𝜅 GLYPH<15> 𝜏 cl GLYPH<28> 1 3 4 1 𝜌 cl 𝑎 cl = 𝜅 GLYPH<15> GLYPH<16> 𝑍 GLYPH<3> mid 𝑍 mid GLYPH<17> 𝜏 cl GLYPH<29> 1 (A8)  \n:  \nThis implies that the medium clumping depth 𝜏 = 𝜅 Σ c will also hit a maximum:  \n𝜏 = min \' 𝜅 GLYPH<15> Σ c GLYPH<150> 𝜏 max \' (A9)  \nwith  \n𝜏 max = 𝜅 GLYPH<15> 𝑍 GLYPH<3> mid 𝑍 mid Σ c = p 2 𝜋 𝑓 fill 𝐻 c 𝑎 cl (A10)  \nwhere the max clumping depth 𝜏 is independent of the dense clump density. Physically, the upper limit of 𝜏 describes that only dust at the surface of clump layer ( 𝜏 glyph[lessorsimilar] 𝜏 max) can diffuse out. The value of 𝜏 max is determined by several uncertain clump properties. In general, bigger clumps that fill a smaller fraction of the space will result in a more transparent medium, and therefore result in a lower 𝜏 max. Study of microphysics on the properties of clumps would be necessary to model 𝜏 max, which is beyond the scope of our work. Here, we simply parametrize this uncertainty in terms of 𝜏 max, for which we take 𝜏 max = 15 in our simulations. A larger 𝜏 max means that the  \nclump sublayer will be even more clumpy for dust to escape. Particles are more trapped in the ring by collision with each other. A lower 𝜏 max glyph[lessorsimilar] 10 will inversely limit the the clump layer more transparent for clump diffusing out, which will hinder the ring\'s survival. The conclusions of this paper are unaffected as long as 𝜏 max GLYPH<29> 1.  \nThis paper has been typeset from a T E X/L A T E X file prepared by the author.'}