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

HAL
open science

Reduction of dust radial drift by turbulence in
protoplanetary disks

Fabiola Antonietta Gerosa, Jérémie Bec, Héloïse Méheut, Anand Utsav
Kapoor

► To cite this version:

Fabiola Antonietta Gerosa, Jérémie Bec, Héloïse Méheut, Anand Utsav Kapoor. Reduction of dust radial drift by turbulence in protoplanetary disks. Astronomy and Astrophysics - A&A, 2024, 685, pp.L4. 10.1051/0004-6361/202449660. hal-04570253

HAL Id: hal-04570253

https://hal.science/hal-04570253

Submitted on 7 May 2024

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

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

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

Distributed under a Creative Commons Attribution 4.0 International License
# Reduction of dust radial drift by turbulence in protoplanetary disks

Fabiola Antonietta Gerosa, J Bec, Héloïse Méheut, Anand Utsav Kapoor

► To cite this version:

Fabiola Antonietta Gerosa, J Bec, Héloïse Méheut, Anand Utsav Kapoor. Reduction of dust radial drift by turbulence in protoplanetary disks. Astronomy & Astrophysics - A&A, 2024, 685, pp.L4. 10.1051/0004-6361/202449660. [hal-04570253](https://hal.science/hal-04570253v1)

HAL Id: hal-04570253

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

Submitted on 7 May 2024

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

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

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

Distributed under a Creative Commons Attribution 4.0 International License
LETTER TO THE EDITOR

# Reduction of dust radial drift by turbulence in protoplanetary disks

Fabiola Antonietta Gerosa<fnref n="1" />, Jérémie Bec<fnref n="2" />,<fnref n="3" />, Héloïse Méheut<fnref n="1" />, and Anand Utsav Kapoor<fnref n="1" />

Received 19 February 2024 / Accepted 16 April 2024

## ABSTRACT

*Context.* Dust particles in protoplanetary disks rotate at velocities exceeding those of the surrounding gas due to a lack of pressure support. Consequently, they experience a headwind from the gas that drives them toward the central star. Radial drift occurs on timescales much shorter than those inferred from disk observations or those required for dust to aggregate and form planets. Additionally, turbulence is often assumed to amplify the radial drift of dust in planet-forming disks when modeled through an effective viscous transport. However, the local interactions between turbulent eddies and particles are known to be significantly more intricate than in a viscous fluid.

*Aims.* Our objective is to elucidate and characterize the dynamic effects of Keplerian turbulence on the mean radial and azimuthal velocities of dust particles.

*Methods.* We employed 2D shearing-box incompressible simulations of the gas, which is maintained in a developed turbulent state while rotating at a sub-Keplerian speed. Dust is modeled as Lagrangian particles set at a Keplerian velocity, therefore experiencing a radial force toward the star through drag.

*Results.* Turbulent eddies are found to reduce the radial drift, while simultaneously enhancing the azimuthal velocities of small particles. This dynamic behavior arises from the modification of dust trajectories due to turbulent eddies.

**Key words.** turbulence – planets and satellites: formation – protoplanetary disks

## 1. Introduction

Radial drift stands out as a pivotal process in the evolution of dust within planet-forming disks, playing a crucial role in planet formation. Dust drift arises from the sub-Keplerian rotation of the gas, exerting drag forces that slow down the inherently Keplerian motion of dust. The resulting depletion of dust angular momentum triggers its radial displacement toward the central star on a timescale significantly shorter than the disk's lifetime. This process is seemingly in contradiction with the observed presence of dust in much older disks. Moreover, planet formation is expected to occur over timescales on the order of a million years, and needs a large amount of dust left in the disk to form the observed massive objects. However, this timescale is significantly longer than the radial drift timescale, leading to an inconsistency in planet formation theories, known as the radial drift barrier (Weidenschilling 1977; Nakagawa et al. 1986). Consequently, planetesimal formation is believed to be localized in specific regions where radial drift is arrested, and dust concentrates.

Recent observational evidence suggests that disk turbulence exhibits lower amplitude than previously believed (Flaherty et al. 2015, 2018; Villenave et al. 2022). This is consistent with the saturated nonlinear state obtained from various hydrodynamical instabilities, such as the vertical shear instability (Arlt & Urpin 2004; Nelson et al. 2013) or the convective overstability (Klahr & Hubbard 2014; Lyra 2014), rather than situations displaying strong turbulence (e.g., through the magnet

torotational instability, Balbus & Hawley 1991). A weak turbulence modifies our understanding of turbulent particle concentration. In particular, rotation and shear become nonnegligible compared to the strength of turbulence, resulting in a reshape of the turbulent eddies. In this context, turbulence cannot be treated as homogeneous and isotropic. Our previous work (Gerosa et al. 2023) demonstrates that rotation, through the Coriolis force, can counterbalance the ejection of particles from turbulent eddies. This phenomenon fosters the clustering of dust and potential planetesimal formation within anticlonic eddies, in contrast with earlier studies that focused on lower rotation rates (Cuzzi et al. 2001; Pan et al. 2011). Furthermore, it shows that Keplerian turbulence cannot be simplistically approximated as a purely diffusive process acting on particles.

Turbulence is frequently invoked for its ability to transport gas angular momentum, often modeled as a viscosity (Shakura & Sunyaev 1973; Lesur 2021). This paradigm leads to gas radial motion, consequently amplifying dust radial drift through drag in an accreting disk (Takeuchi & Lin 2002). However, such a mean-field approach neglects the local modifications that turbulent eddies can impose on particle trajectories. The nuanced understanding of non-isotropic turbulence and its interplay with dust dynamics, as elucidated in Gerosa et al. (2023), prompts an investigation into whether Keplerian turbulence additionally alters dust drift. Turbulent gas could potentially slow down or even halt the radial drift of dust, thus overcoming the need for an azimuthally extended pressure bump, which

<footnote n="1">¹ Université Côte d'Azur, Observatoire de la Côte d'Azur, CNRS, Laboratoire Lagrange, Nice, France
e-mail: fabiola.gerosa@oca.eu</footnote>

<footnote n="2">² Université Côte d'Azur, CNRS, Institut de Physique de Nice, Nice, France</footnote>

<footnote n="3">³ Université Côte d'Azur, Inria, CNRS, Calisto team, Sophia Antipolis, France</footnote>
LETTER TO THE EDITOR

# Reduction of dust radial drift by turbulence in protoplanetary disks

Fabiola Antonietta Gerosa<fnref n="1" />, Jérémie Bec<fnref n="2" />,<fnref n="3" />, Héloïse Méheut<fnref n="1" />, and Anand Utsav Kapoor<fnref n="1" />

Received 19 February 2024 / Accepted 16 April 2024

## ABSTRACT

*Context.* Dust particles in protoplanetary disks rotate at velocities exceeding those of the surrounding gas due to a lack of pressure support. Consequently, they experience a headwind from the gas that drives them toward the central star. Radial drift occurs on timescales much shorter than those inferred from disk observations or those required for dust to aggregate and form planets. Additionally, turbulence is often assumed to amplify the radial drift of dust in planet-forming disks when modeled through an effective viscous transport. However, the local interactions between turbulent eddies and particles are known to be significantly more intricate than in a viscous fluid.

*Aims.* Our objective is to elucidate and characterize the dynamic effects of Keplerian turbulence on the mean radial and azimuthal velocities of dust particles.

*Methods.* We employed 2D shearing-box incompressible simulations of the gas, which is maintained in a developed turbulent state while rotating at a sub-Keplerian speed. Dust is modeled as Lagrangian particles set at a Keplerian velocity, therefore experiencing a radial force toward the star through drag.

*Results.* Turbulent eddies are found to reduce the radial drift, while simultaneously enhancing the azimuthal velocities of small particles. This dynamic behavior arises from the modification of dust trajectories due to turbulent eddies.

**Key words.** turbulence – planets and satellites: formation – protoplanetary disks

## 1. Introduction

Radial drift stands out as a pivotal process in the evolution of dust within planet-forming disks, playing a crucial role in planet formation. Dust drift arises from the sub-Keplerian rotation of the gas, exerting drag forces that slow down the inherently Keplerian motion of dust. The resulting depletion of dust angular momentum triggers its radial displacement toward the central star on a timescale significantly shorter than the disk's lifetime. This process is seemingly in contradiction with the observed presence of dust in much older disks. Moreover, planet formation is expected to occur over timescales on the order of a million years, and needs a large amount of dust left in the disk to form the observed massive objects. However, this timescale is significantly longer than the radial drift timescale, leading to an inconsistency in planet formation theories, known as the radial drift barrier (Weidenschilling 1977; Nakagawa et al. 1986). Consequently, planetesimal formation is believed to be localized in specific regions where radial drift is arrested, and dust concentrates.

Recent observational evidence suggests that disk turbulence exhibits lower amplitude than previously believed (Flaherty et al. 2015, 2018; Villenave et al. 2022). This is consistent with the saturated nonlinear state obtained from various hydrodynamical instabilities, such as the vertical shear instability (Arlt & Urpin 2004; Nelson et al. 2013) or the convective overstability (Klahr & Hubbard 2014; Lyra 2014), rather than situations displaying strong turbulence (e.g., through the magnet

torotational instability, Balbus & Hawley 1991). A weak turbulence modifies our understanding of turbulent particle concentration. In particular, rotation and shear become nonnegligible compared to the strength of turbulence, resulting in a reshape of the turbulent eddies. In this context, turbulence cannot be treated as homogeneous and isotropic. Our previous work (Gerosa et al. 2023) demonstrates that rotation, through the Coriolis force, can counterbalance the ejection of particles from turbulent eddies. This phenomenon fosters the clustering of dust and potential planetesimal formation within anticlonic eddies, in contrast with earlier studies that focused on lower rotation rates (Cuzzi et al. 2001; Pan et al. 2011). Furthermore, it shows that Keplerian turbulence cannot be simplistically approximated as a purely diffusive process acting on particles.

Turbulence is frequently invoked for its ability to transport gas angular momentum, often modeled as a viscosity (Shakura & Sunyaev 1973; Lesur 2021). This paradigm leads to gas radial motion, consequently amplifying dust radial drift through drag in an accreting disk (Takeuchi & Lin 2002). However, such a mean-field approach neglects the local modifications that turbulent eddies can impose on particle trajectories. The nuanced understanding of non-isotropic turbulence and its interplay with dust dynamics, as elucidated in Gerosa et al. (2023), prompts an investigation into whether Keplerian turbulence additionally alters dust drift. Turbulent gas could potentially slow down or even halt the radial drift of dust, thus overcoming the need for an azimuthally extended pressure bump, which

<footnote n="1">¹ Université Côte d'Azur, Observatoire de la Côte d'Azur, CNRS, Laboratoire Lagrange, Nice, France
e-mail: fabiola.gerosa@oca.eu</footnote>

<footnote n="2">² Université Côte d'Azur, CNRS, Institut de Physique de Nice, Nice, France</footnote>

<footnote n="3">³ Université Côte d'Azur, Inria, CNRS, Calisto team, Sophia Antipolis, France</footnote>
is the conventional solution to this barrier (Whipple et al. 1972; Pinilla et al. 2012).

The aim of this study is to investigate the impact of gas turbulence on the radial and azimuthal drift of dust. We employ 2D numerical simulations of forced Keplerian turbulence within a shearing box and track the dynamics of drifting Lagrangian particles in this flow, characterized by the absence of mean radial velocity. Our findings reveal how turbulence plays a crucial role in diminishing the radial drift of particles, while concurrently enhancing their azimuthal drift.

## 2. Model and numerical methods

### 2.1. Gas

We assumed that the gas is incompressible; thus, its dynamics are governed by a divergence-free velocity field $v$ that satisfies the Navier-Stokes equation, with the inclusion of the external gravitational potential from the star. Additionally, we considered the gas flow to be two-dimensional. We considered a small region located at distance $r_0$ from the star, rotating at a rate $\Omega$ along with the gas, and we employed the periodic shearing box approximation. If $(r, \theta)$ are the polar coordinates centered onto the star, we defined the local coordinates $x = r - r_0$ and $y = r_0\theta$. The total gas velocity solves as

$$ \partial_t v + v \cdot \nabla v = -\frac{1}{\rho_g} \nabla p + v \nabla^2 v - 2\Omega \times v + 3\Omega^2 x e_x + f - D, \quad (1) $$

where $\rho_g$ and $v$ represent the gas mass density and kinematic viscosity, respectively. The computational domain has size $L_x \times L_y$, with an aspect ratio of $L_y/L_x = 4$ to limit spurious geometrical effects at large $\Omega$. The fluctuating velocity of the gas $u = v + \frac{3}{2}\Omega \times e_y$ is assumed periodic in the frame distorted by the mean shear. We conducted numerical simulations with 256 $\times$ 1024 collocation points using the open-source pseudo-spectral solver SNOOPY (see Lesur & Longaretti 2005).

The turbulent flow was maintained in a statistically steady state by introducing both a forcing term, $f$, and a large-scale dissipation $D = \gamma(v + \frac{3}{2}\Omega \times e_y)$. The prescribed force is random, Gaussian, homogeneous, and isotropic in the sheared frame, with a zero mean, white noise in time, and with spatial correlations concentrated at large scales. This implementation allowed us to investigate a generalized form of 2D turbulence. Although the resulting flow exhibited some similarities with gas subjected to the subcritical baroclinic instability (Lesur & Papaloizou 2010), investigating the transition to a turbulent state through instabilities is beyond the scope of this paper. The forcing has a fixed amplitude, so that the mean injection rates of kinetic energy $\epsilon_l$ and entrophy $\eta_l$ are prescribed. This specifies a forcing length-scale $\ell_f = (\epsilon_l/\eta_l)^{1/2}$ and timescale $\tau_f = (\eta_l)^{-1/3}$, which remain constant across all runs. We can also define the small eddy turnover time $\tau_\omega = \langle\omega^2\rangle^{-1/2}$, where $\omega = \partial_x u_y - \partial_y u_x$ is the vorticity of the gas turbulent fluctuations. Using this timescale, we nondimensionalize the rotation rate $\Omega$, defining the Rossby number $Ro = 1/(\tau_\omega \Omega)$. It serves as a measure of turbulence strength relative to rotation. For our study, we selected values of $Ro$ in the range of [0.1, 10] (see Appendix A for details on the simulation parameters). We note that at small $Ro$, rotation becomes predominant over turbulence. The effects of such weak turbulence on particle dynamics have been understudied until now.

### 2.2. Dust

The position and velocity ($R_p, U_p$) of a particle in the rotating and sheared frame solve as

$$ \frac{dR_p}{dt} = U_p + \frac{3}{2} \Omega R_{p,x} e_y, \quad (2) $$

$$ \frac{dU_p}{dt} = -\frac{1}{\tau_p} [U_p - u(R_p, t)] - 2\Omega \times U_p + \frac{3}{2} \Omega U_{p,x} e_y - \epsilon e_x. \quad (3) $$

The equation governing the motion of dust particles is derived in Appendix B. The first term on the right-hand side of Eq. (3) comes from the drag with the gas. It involves the particle stopping time $\tau_p$ defined in Epstein regime as $\tau_p = (8/\pi)^{1/2}\rho_p a/(\rho_g c_s)$, with $\rho_p$ denoting the particle material mass density, $a$ its radius, and $c_s$ the speed of sound. We nondimensionalize $\tau_p$ by defining the Stokes number $St = \tau_p \Omega$. The subsequent term appearing in Eq. (3) accounts for the Coriolis force. The last term represents the adjustment for the difference in azimuthal velocity between gas and dust in the shearing box framework. It involves the parameter $\epsilon = (\Omega_K^2 - \Omega^2)r_0$, where $\Omega_K$ denotes the particle Keplerian rotation rate. When scaled relative to the two other fixed parameters quantities of the simulation, namely $\tau_f$ and $\ell_f$, which prescribe the level of gas turbulence, we can define the dimensionless drift parameter $\tilde{\epsilon} = \epsilon\tau_f^2/\ell_f$. We adopt the fiducial value $\tilde{\epsilon} = 0.1$ to examine the effects of weak but nonnegligible turbulence on dust drift. Nonetheless, we extensively explore the influence of this parameter on our results in the discussion, considering values in the range $[10^{-2}, 10]$. We neglect both particle-particle interactions and their feedback onto the gas.

## 3. Results

An overview of the results is depicted in Fig. 1. The strength of rotation with respect to turbulence is increased along the vertical axis, while the horizontal axis represents the particle Stokes number, varying from small to large grains. The orange background zone highlights the region where dust forms point clusters within anticyclones, whereas in the gray background area dust fills the whole space. In the other regions of the parameter space, dust particles tend to concentrate on filaments. For a more detailed discussion of distinct dust concentration regimes, we refer to Gerosa et al. (2023) or to the highlights given in Appendix C. In the lower region, indicative of strong turbulence and/or slow rotation, dust particles tend to be expelled from eddies. Conversely, in the upper region, dust predominantly resides within anticyclones. The transition between these two behaviors can be quantitatively characterized (see Appendix D), and is identified in Fig. 1 by the dashed line. To determine whether turbulence diminishes the radial drift of dust, we compare its mean speed in our turbulent runs to that in the laminar case. In the absence of turbulence, the particle radial drift velocity is given by $U_{p,x}^{NT} = -\tau_p \epsilon / (1 + \tau_p^2 \Omega^2)$ (as detailed in Appendix E). We therefore define the relative velocity $\Delta U_{p,x}/U_{p,x}^{NT} = ((U_{p,x}) - U_{p,x}^{NT})/U_{p,x}^{NT}$, where $\langle\cdot\rangle$ denotes the average over all particles and times. This quantity measures the variation in dust radial drift velocity attributable to turbulence. A negative (positive) value indicates a decrease (increase) in the radial drift of dust compared to the laminar case. It is observed in Fig. 1 that the radial drift velocity is significantly reduced at large $Ro^{-1}$ and small Stokes numbers. This indicates that, for these parameters, turbulence substantially slows down dust radial drift. It is noteworthy that, for $Ro^{-1} = 3.29$, the plot is the result of an average
is the conventional solution to this barrier (Whipple et al. 1972; Pinilla et al. 2012).

The aim of this study is to investigate the impact of gas turbulence on the radial and azimuthal drift of dust. We employ 2D numerical simulations of forced Keplerian turbulence within a shearing box and track the dynamics of drifting Lagrangian particles in this flow, characterized by the absence of mean radial velocity. Our findings reveal how turbulence plays a crucial role in diminishing the radial drift of particles, while concurrently enhancing their azimuthal drift.

## 2. Model and numerical methods

### 2.1. Gas

We assumed that the gas is incompressible; thus, its dynamics are governed by a divergence-free velocity field $v$ that satisfies the Navier-Stokes equation, with the inclusion of the external gravitational potential from the star. Additionally, we considered the gas flow to be two-dimensional. We considered a small region located at distance $r_0$ from the star, rotating at a rate $\Omega$ along with the gas, and we employed the periodic shearing box approximation. If $(r, \theta)$ are the polar coordinates centered onto the star, we defined the local coordinates $x = r - r_0$ and $y = r_0\theta$. The total gas velocity solves as

$$ \partial_t v + v \cdot \nabla v = -\frac{1}{\rho_g} \nabla p + v \nabla^2 v - 2\Omega \times v + 3\Omega^2 x e_x + f - D, \quad (1) $$

where $\rho_g$ and $v$ represent the gas mass density and kinematic viscosity, respectively. The computational domain has size $L_x \times L_y$, with an aspect ratio of $L_y/L_x = 4$ to limit spurious geometrical effects at large $\Omega$. The fluctuating velocity of the gas $u = v + \frac{3}{2}\Omega \times e_y$ is assumed periodic in the frame distorted by the mean shear. We conducted numerical simulations with 256 $\times$ 1024 collocation points using the open-source pseudo-spectral solver SNOOPY (see Lesur & Longaretti 2005).

The turbulent flow was maintained in a statistically steady state by introducing both a forcing term, $f$, and a large-scale dissipation $D = \gamma(v + \frac{3}{2}\Omega \times e_y)$. The prescribed force is random, Gaussian, homogeneous, and isotropic in the sheared frame, with a zero mean, white noise in time, and with spatial correlations concentrated at large scales. This implementation allowed us to investigate a generalized form of 2D turbulence. Although the resulting flow exhibited some similarities with gas subjected to the subcritical baroclinic instability (Lesur & Papaloizou 2010), investigating the transition to a turbulent state through instabilities is beyond the scope of this paper. The forcing has a fixed amplitude, so that the mean injection rates of kinetic energy $\epsilon_l$ and entrophy $\eta_l$ are prescribed. This specifies a forcing length-scale $\ell_f = (\epsilon_l/\eta_l)^{1/2}$ and timescale $\tau_f = (\eta_l)^{-1/3}$, which remain constant across all runs. We can also define the small eddy turnover time $\tau_\omega = \langle\omega^2\rangle^{-1/2}$, where $\omega = \partial_x u_y - \partial_y u_x$ is the vorticity of the gas turbulent fluctuations. Using this timescale, we nondimensionalize the rotation rate $\Omega$, defining the Rossby number $Ro = 1/(\tau_\omega\Omega)$. It serves as a measure of turbulence strength relative to rotation. For our study, we selected values of $Ro$ in the range of [0.1, 10] (see Appendix A for details on the simulation parameters). We note that at small $Ro$, rotation becomes predominant over turbulence. The effects of such weak turbulence on particle dynamics have been understudied until now.

### 2.2. Dust

The position and velocity ($R_p, U_p$) of a particle in the rotating and sheared frame solve as

$$ \frac{dR_p}{dt} = U_p + \frac{3}{2}\Omega R_{p,x} e_y, \quad (2) $$

$$ \frac{dU_p}{dt} = -\frac{1}{\tau_p} [U_p - u(R_p, t)] - 2\Omega \times U_p + \frac{3}{2} \Omega U_{p,x} e_y - \epsilon e_x. \quad (3) $$

The equation governing the motion of dust particles is derived in Appendix B. The first term on the right-hand side of Eq. (3) comes from the drag with the gas. It involves the particle stopping time $\tau_p$ defined in Epstein regime as $\tau_p = (8/\pi)^{1/2}\rho_p a/(\rho_g c_s)$, with $\rho_p$ denoting the particle material mass density, $a$ its radius, and $c_s$ the speed of sound. We nondimensionalize $\tau_p$ by defining the Stokes number $St = \tau_p \Omega$. The subsequent term appearing in Eq. (3) accounts for the Coriolis force. The last term represents the adjustment for the difference in azimuthal velocity between gas and dust in the shearing box framework. It involves the parameter $\epsilon = (\Omega_K^2 - \Omega^2)r_0$, where $\Omega_K$ denotes the particle Keplerian rotation rate. When scaled relative to the two other fixed parameters quantities of the simulation, namely $\tau_f$ and $\ell_f$, which prescribe the level of gas turbulence, we can define the dimensionless drift parameter $\tilde{\epsilon} = \epsilon\tau_f^2/\ell_f$. We adopt the fiducial value $\tilde{\epsilon} = 0.1$ to examine the effects of weak but nonnegligible turbulence on dust drift. Nonetheless, we extensively explore the influence of this parameter on our results in the discussion, considering values in the range $[10^{-2}, 10]$. We neglect both particle-particle interactions and their feedback onto the gas.

## 3. Results

An overview of the results is depicted in Fig. 1. The strength of rotation with respect to turbulence is increased along the vertical axis, while the horizontal axis represents the particle Stokes number, varying from small to large grains. The orange background zone highlights the region where dust forms point clusters within anticyclones, whereas in the gray background area dust fills the whole space. In the other regions of the parameter space, dust particles tend to concentrate on filaments. For a more detailed discussion of distinct dust concentration regimes, we refer to Gerosa et al. (2023) or to the highlights given in Appendix C. In the lower region, indicative of strong turbulence and/or slow rotation, dust particles tend to be expelled from eddies. Conversely, in the upper region, dust predominantly resides within anticyclones. The transition between these two behaviors can be quantitatively characterized (see Appendix D), and is identified in Fig. 1 by the dashed line. To determine whether turbulence diminishes the radial drift of dust, we compare its mean speed in our turbulent runs to that in the laminar case. In the absence of turbulence, the particle radial drift velocity is given by $U_{p,x}^{NT} = -\tau_p\epsilon/(1 + \tau_p^2\Omega^2)$ (as detailed in Appendix E). We therefore define the relative velocity $\Delta U_{p,x}/U_{p,x}^{NT} = ((U_{p,x}) - U_{p,x}^{NT})/U_{p,x}^{NT}$, where $\langle\cdot\rangle$ denotes the average over all particles and times. This quantity measures the variation in dust radial drift velocity attributable to turbulence. A negative (positive) value indicates a decrease (increase) in the radial drift of dust compared to the laminar case. It is observed in Fig. 1 that the radial drift velocity is significantly reduced at large $Ro^{-1}$ and small Stokes numbers. This indicates that, for these parameters, turbulence substantially slows down dust radial drift. It is noteworthy that, for $Ro^{-1} = 3.29$, the plot is the result of an average
<figure><img src="image_3.png" /><figcaption>**Fig. 1.** Variation in dust radial drift velocity due to turbulence for $\bar{\epsilon} = 0.1$ as a function of the inverse Rossby number and the Stokes number.</figcaption></figure>

over 12 runs. Because of particles clustering to a point in phase space, many runs are needed here to achieve sufficiently high statistics.

### 3.1. Radial drift

The mean radial drift velocity of particles is presented in Fig. 2. The dashed line corresponds to the nonturbulent case, $U_{p,x}^{NT}$. For small Stokes numbers, the radial drift velocity exhibits a substantial reduction for almost all the Rossby numbers, diminishing with faster rotation rates (larger $Ro^{-1}$). This can be clearly identified on this figure as the absolute value of the mean radial velocities are smaller than in the laminar case. For $Ro^{-1} > 0.7$, dust lies longer in anticyclones, where it gets locked to the gas radial velocity<fnref n="1" />. The turbulent decrease in dust drift reaches up to an order of magnitude in some cases (i.e., circles for which $\Delta U_{p,x}/U_{p,x}^{NT} = -90\%$ in Fig. 1). The yellow circles in Fig. 1, on the other hand, correspond to a complete cessation of radial drift due to turbulence, with particles clustered to a point and trapped inside an anticyclone. For $Ro^{-1} = 0.49$ particles mostly reside between eddies, yet their radial drift is still slowed down by turbulence. For $Ro^{-1} = 0.10$, instead, small particles drift faster toward the star, as usually assumed for strong turbulence. We delve into the interpretation of these last two phenomena in Sect. 4. Finally, for all Rossby numbers, $\langle U_{p,x} \rangle$ converges at large Stokes numbers to $U_{p,x}^{NT}$, indicating large solids weakly coupled to the turbulent flow.

### 3.2. Azimuthal drift

Figure 3 shows the mean azimuthal drift velocity $\langle U_{p,y} \rangle$ of dust particles. The laminar drift velocity in the azimuthal direction is given by $U_{p,y}^{NT} = \tau_p^2 \Omega \epsilon / (2(1 + \tau_p^2 \Omega^2))$ (see Appendix E), plotted as a dashed line. Turbulence enhances the azimuthal drift of small dust particles, while for large St, it remains unaltered. The values of $\langle U_{p,y} \rangle$ clearly appear larger at slower rotation rates

<figure><img src="image_4.png" /><figcaption>**Fig. 2.** Mean radial drift velocity as a function of the Stokes number, for various values of the Rossby number and for $\bar{\epsilon} = 0.1$. The dotted line represents the radial drift velocity particles would have in a laminar gas flow.</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>**Fig. 3.** Mean azimuthal drift velocity as a function of the Stokes number, for various values of the Rossby number and for $\bar{\epsilon} = 0.1$. The dotted line represents the azimuthal drift velocity particles would have in a laminar gas flow.</figcaption></figure>

(small $Ro^{-1}$). It is important to note that the magnitude of the strongest azimuthal drift velocity represents only a small fraction of the Keplerian velocity ($\langle U_{p,y} \rangle \approx 10^{-4} r_0 \Omega$ for $Ro^{-1} = 0.10$ and $St = 0.1$). Therefore, it would not currently be detectable from observational data. Conversely, the azimuthal drift is reduced at very high rotation rates and intermediate Stokes numbers. This is again consistent with the mechanism of dust trapping inside anticyclones, resulting in particle cluster velocities governed by the mean flow.

<footnote n="1">¹ By construction, there is no mean gas radial velocity in the shearing box.</footnote>
<figure><img src="image_2.png" /><figcaption>**Fig. 1.** Variation in dust radial drift velocity due to turbulence for $\bar{\epsilon} = 0.1$ as a function of the inverse Rossby number and the Stokes number.</figcaption></figure>

over 12 runs. Because of particles clustering to a point in phase space, many runs are needed here to achieve sufficiently high statistics.

## 3.1. Radial drift

The mean radial drift velocity of particles is presented in Fig. 2. The dashed line corresponds to the nonturbulent case, $U_{p,x}^{NT}$. For small Stokes numbers, the radial drift velocity exhibits a substantial reduction for almost all the Rossby numbers, diminishing with faster rotation rates (larger $Ro^{-1}$). This can be clearly identified on this figure as the absolute value of the mean radial velocities are smaller than in the laminar case. For $Ro^{-1} > 0.7$, dust lies longer in anticyclones, where it gets locked to the gas radial velocity<fnref n="1" />. The turbulent decrease in dust drift reaches up to an order of magnitude in some cases (i.e., circles for which $\Delta U_{p,x}/U_{p,x}^{NT} = -90\%$ in Fig. 1). The yellow circles in Fig. 1, on the other hand, correspond to a complete cessation of radial drift due to turbulence, with particles clustered to a point and trapped inside an anticyclone. For $Ro^{-1} = 0.49$ particles mostly reside between eddies, yet their radial drift is still slowed down by turbulence. For $Ro^{-1} = 0.10$, instead, small particles drift faster toward the star, as usually assumed for strong turbulence. We delve into the interpretation of these last two phenomena in Sect. 4. Finally, for all Rossby numbers, $\langle U_{p,x} \rangle$ converges at large Stokes numbers to $U_{p,x}^{NT}$, indicating large solids weakly coupled to the turbulent flow.

## 3.2. Azimuthal drift

Figure 3 shows the mean azimuthal drift velocity $\langle U_{p,y} \rangle$ of dust particles. The laminar drift velocity in the azimuthal direction is given by $U_{p,y}^{NT} = \tau_p^2 \Omega \epsilon / (2(1 + \tau_p^2 \Omega^2))$ (see Appendix E), plotted as a dashed line. Turbulence enhances the azimuthal drift of small dust particles, while for large St, it remains unaltered. The values of $\langle U_{p,y} \rangle$ clearly appear larger at slower rotation rates

<figure><img src="image_3.png" /><figcaption>**Fig. 2.** Mean radial drift velocity as a function of the Stokes number, for various values of the Rossby number and for $\bar{\epsilon} = 0.1$. The dotted line represents the radial drift velocity particles would have in a laminar gas flow.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>**Fig. 3.** Mean azimuthal drift velocity as a function of the Stokes number, for various values of the Rossby number and for $\bar{\epsilon} = 0.1$. The dotted line represents the azimuthal drift velocity particles would have in a laminar gas flow.</figcaption></figure>

(small $Ro^{-1}$). It is important to note that the magnitude of the strongest azimuthal drift velocity represents only a small fraction of the Keplerian velocity ($\langle U_{p,y} \rangle \approx 10^{-4} r_0 \Omega$ for $Ro^{-1} = 0.10$ and $St = 0.1$). Therefore, it would not currently be detectable from observational data. Conversely, the azimuthal drift is reduced at very high rotation rates and intermediate Stokes numbers. This is again consistent with the mechanism of dust trapping inside anticyclones, resulting in particle cluster velocities governed by the mean flow.

<footnote n="1">¹ By construction, there is no mean gas radial velocity in the shearing box.</footnote>
# 4. Discussion

## 4.1. Preferential sweeping

Preferential sweeping is a well-known process that induces variations in particle speed, due to turbulence in the gas phase. This phenomenon has been extensively studied and demonstrated in the context of droplets or aerosols settling in the atmosphere (Maxey 1987; see also Bec et al. 2024 for a recent review). In such scenarios the interaction between gravity and turbulence leads to an acceleration of particle settling as particles tend to preferentially sample down-welling regions of the turbulent flow. We anticipate a similar mechanism for dust particles moving within turbulent protoplanetary disks, where gravity drives radial drift of the dust toward the central star. Figure 2 shows, particularly for $Ro^{-1} = 0.10$, that the radial drift is increased by turbulence for intermediate Stokes numbers. A simple analysis using cellular eddies, as depicted in Fig. 4a, further elucidates this point. Cyclonic eddies (rotating in the direction of the disk) are shown in red, while anticyclonic eddies are shown in blue. It is evident that dust particles exclusively sample left-welling regions, attracted by the star’s gravity, resulting in increased radial drift velocities. When inspecting snapshots from our numerical simulations (Figs. 4b and d), the dynamics of dust appears to be more intricate. Here, particles predominantly sample left-welling regions, while also undergoing azimuthal drift. During this process, they are once again subject to preferential sweeping. Consequently, they experience an additional acceleration in the azimuthal direction, as seen for small enough Stokes numbers in Fig. 3.

## 4.2. Reduction of radial dust drift

Figure 2 shows that, for $Ro^{-1} > 0.10$, the radial velocity of dust particles is decreased by turbulence. Two different regimes should be considered here. For $Ro^{-1} > 0.7$, particles mostly remain trapped inside anticyclones. Therefore, the reduction in their drift can be readily explained: particles evolve with the radial velocity of eddies, which is slower than the dust velocity. On the other hand, at $Ro^{-1} = 0.49$, where particles are instead ejected from eddies, counterintuitively in regard to the effect of preferential sweeping on dust speed, our findings are in contrast with the typical enhancement of dust speed caused by standard preferential sweeping. The reason for this radial drift reduction lies in the Coriolis force, which induces an interchange of azimuthal and radial velocities:

$$f_C = -2\Omega \times U_p = 2\Omega U_{p,y} e_x - 2\Omega U_{p,x} e_y. \quad (4)$$

The decrease in radial drift velocities by turbulence can therefore be easily explained. While preferential sweeping would normally enhance the radial velocity of dust, it also increases $U_{p,y}$ that, when entering the Coriolis term, can counteract the acceleration of radial drift. This effect can even prevail, particularly at high rotation rates, resulting in an effective decrease in dust radial drift velocities. Moreover, as the rotation rate increases, the shear renders the turbulence increasingly non-isotropic (Gerosa et al. 2023). In particular, eddies become strongly stretched azimuthally, as can be seen when comparing Fig. 4b ($Ro^{-1} = 0.10$) and Fig. 4d ($Ro^{-1} = 0.49$). This, in turn, affects the dynamics of dust. To better understand this mechanism, in Figs. 4a and c, we illustrate dust dynamics in a model flow consisting of cellular eddies with two different aspect ratios to represent vortex stretching due to shear. In both cases, particles follow the preferential sweeping scenario, solely sampling

<figure><img src="image_6.png" /><figcaption>Fig. 4. Vorticity $\omega$, normalized by $\tau_\omega$, and dust particle position $R_p$ (black dots) for $St \approx 0.1$. The snapshots in a) and c) were obtained for a simple cellular flow, while from our direct numerical simulations at $\tilde{\epsilon} = 0.1$ and $Ro^{-1} = 0.10$ for b) and $Ro^{-1} = 0.49$ for d).</figcaption></figure>

left-welling regions. However, with elongated eddies, dust paths are predominantly azimuthally oriented, thus resulting in a significant reduction in particle mean radial motion.

## 4.3. Estimation of the turbulent parameter $\alpha$

Turbulent intensity is classically defined as the ratio of the amplitude of the turbulent fluctuations $u$ of the gas to its mean velocity. This quantity is the most relevant for estimating the effects of turbulence on particle dynamics. From its squared value, we can compute the turbulent $\alpha$ parameter as

$$\alpha = \frac{\langle |u|^2 \rangle}{(\Omega H)^2}. \quad (5)$$

The modifications of dust velocities, due to preferential sweeping, and of particle concentrations can be estimated through this value. Considering $H = 10 l_f$, the $\alpha$ parameter can be as large as $10^{-2}$ for $Ro^{-1} = 0.10$. At this rotation rate, turbulence is ejecting particles from eddies and accelerating their drift. Its effect on particle dynamics is therefore similar to what is already known in the literature for strong MRI-like turbulence (Yang et al. 2018). On the other hand, at $Ro^{-1} = 6.19$ we obtain $\alpha = 10^{-5}$. At these low turbulence levels, not extensively studied in the literature, turbulent anticyclones concentrate dust and slow down drift. In this case, the Coriolis force due to rotation is able to counteract the diffusive power of eddies. We therefore observe a substantial difference in dust dynamics in weak turbulence compared to the purely laminar case ($Ro^{-1} = \infty$), where particles would be randomly distributed and would drift radially at a considerably faster velocity.

When instead considering the angular momentum transport due to turbulence (Shakura & Sunyaev 1973), we can compute a viscous $\alpha$ using the Reynolds stress $\langle u_x u_y \rangle$ (in place of $\langle |u|^2 \rangle$). For this parameter we find values that are two orders of magnitude lower than those given above. This indicates a small amount of turbulent momentum transport in our simulations, in line with
# 4. Discussion

## 4.1. Preferential sweeping

Preferential sweeping is a well-known process that induces variations in particle speed, due to turbulence in the gas phase. This phenomenon has been extensively studied and demonstrated in the context of droplets or aerosols settling in the atmosphere (Maxey 1987; see also Bec et al. 2024 for a recent review). In such scenarios the interaction between gravity and turbulence leads to an acceleration of particle settling as particles tend to preferentially sample down-welling regions of the turbulent flow. We anticipate a similar mechanism for dust particles moving within turbulent protoplanetary disks, where gravity drives radial drift of the dust toward the central star. Figure 2 shows, particularly for $Ro^{-1} = 0.10$, that the radial drift is increased by turbulence for intermediate Stokes numbers. A simple analysis using cellular eddies, as depicted in Fig. 4a, further elucidates this point. Cyclonic eddies (rotating in the direction of the disk) are shown in red, while anticyclonic eddies are shown in blue. It is evident that dust particles exclusively sample left-welling regions, attracted by the star’s gravity, resulting in increased radial drift velocities. When inspecting snapshots from our numerical simulations (Figs. 4b and d), the dynamics of dust appears to be more intricate. Here, particles predominantly sample left-welling regions, while also undergoing azimuthal drift. During this process, they are once again subject to preferential sweeping. Consequently, they experience an additional acceleration in the azimuthal direction, as seen for small enough Stokes numbers in Fig. 3.

## 4.2. Reduction of radial dust drift

Figure 2 shows that, for $Ro^{-1} > 0.10$, the radial velocity of dust particles is decreased by turbulence. Two different regimes should be considered here. For $Ro^{-1} > 0.7$, particles mostly remain trapped inside anticyclones. Therefore, the reduction in their drift can be readily explained: particles evolve with the radial velocity of eddies, which is slower than the dust velocity. On the other hand, at $Ro^{-1} = 0.49$, where particles are instead ejected from eddies, counterintuitively in regard to the effect of preferential sweeping on dust speed, our findings are in contrast with the typical enhancement of dust speed caused by standard preferential sweeping. The reason for this radial drift reduction lies in the Coriolis force, which induces an interchange of azimuthal and radial velocities:

$$f_C = -2\Omega \times U_p = 2\Omega U_{p,y} e_x - 2\Omega U_{p,x} e_y. \quad (4)$$

The decrease in radial drift velocities by turbulence can therefore be easily explained. While preferential sweeping would normally enhance the radial velocity of dust, it also increases $U_{p,y}$ that, when entering the Coriolis term, can counteract the acceleration of radial drift. This effect can even prevail, particularly at high rotation rates, resulting in an effective decrease in dust radial drift velocities. Moreover, as the rotation rate increases, the shear renders the turbulence increasingly non-isotropic (Gerosa et al. 2023). In particular, eddies become strongly stretched azimuthally, as can be seen when comparing Fig. 4b ($Ro^{-1} = 0.10$) and Fig. 4d ($Ro^{-1} = 0.49$). This, in turn, affects the dynamics of dust. To better understand this mechanism, in Figs. 4a and c, we illustrate dust dynamics in a model flow consisting of cellular eddies with two different aspect ratios to represent vortex stretching due to shear. In both cases, particles follow the preferential sweeping scenario, solely sampling

<figure><img src="image_5.png" /><figcaption>Fig. 4. Vorticity $\omega$, normalized by $\tau_\omega$, and dust particle position $R_p$ (black dots) for $St \approx 0.1$. The snapshots in a) and c) were obtained for a simple cellular flow, while from our direct numerical simulations at $\tilde{\epsilon} = 0.1$ and $Ro^{-1} = 0.10$ for b) and $Ro^{-1} = 0.49$ for d).</figcaption></figure>

left-welling regions. However, with elongated eddies, dust paths are predominantly azimuthally oriented, thus resulting in a significant reduction in particle mean radial motion.

## 4.3. Estimation of the turbulent parameter $\alpha$

Turbulent intensity is classically defined as the ratio of the amplitude of the turbulent fluctuations $u$ of the gas to its mean velocity. This quantity is the most relevant for estimating the effects of turbulence on particle dynamics. From its squared value, we can compute the turbulent $\alpha$ parameter as

$$\alpha = \frac{\langle |u|^2 \rangle}{(\Omega H)^2}. \quad (5)$$

The modifications of dust velocities, due to preferential sweeping, and of particle concentrations can be estimated through this value. Considering $H = 10 l_f$, the $\alpha$ parameter can be as large as $10^{-2}$ for $Ro^{-1} = 0.10$. At this rotation rate, turbulence is ejecting particles from eddies and accelerating their drift. Its effect on particle dynamics is therefore similar to what is already known in the literature for strong MRI-like turbulence (Yang et al. 2018). On the other hand, at $Ro^{-1} = 6.19$ we obtain $\alpha = 10^{-5}$. At these low turbulence levels, not extensively studied in the literature, turbulent anticyclones concentrate dust and slow down drift. In this case, the Coriolis force due to rotation is able to counteract the diffusive power of eddies. We therefore observe a substantial difference in dust dynamics in weak turbulence compared to the purely laminar case ($Ro^{-1} = \infty$), where particles would be randomly distributed and would drift radially at a considerably faster velocity.

When instead considering the angular momentum transport due to turbulence (Shakura & Sunyaev 1973), we can compute a viscous $\alpha$ using the Reynolds stress $\langle u_x u_y \rangle$ (in place of $\langle |u|^2 \rangle$). For this parameter we find values that are two orders of magnitude lower than those given above. This indicates a small amount of turbulent momentum transport in our simulations, in line with
<figure><img src="image_7.png" /><figcaption>Fig. 5. Phase diagram of the drift parameter vs. the inverse Rossby number. The highlighted parameters present the results of the dust concentration and radial drift for St = 0.1.</figcaption></figure>

recent observations that tend to infer low turbulent viscosity in disks.

### 4.4. Dependence on the drift parameter $\tilde{\epsilon}$

The shearing box approximation holds true only when the box is located at a distance $r_0$ from the star much greater than the largest scales of the simulation. This requires $r_0 \gg l_f$. Consequently, we obtain the condition $\epsilon \gg (\Omega_K^2 - \Omega^2) l_f$. The acceptable values of $(\Omega_K^2 - \Omega^2)$ are determined from

$$v_y^2/v_K^2 = 1 - n c_s^2/v_K^2, \quad (6)$$

where n is the index of the power law describing pressure as a function of radius. In the literature the typical disk properties yield 0.996 < $v_y^2/v_K^2$ < 0.999. For our analysis we make the arbitrary choice $v_y^2/v_K^2 = 0.998$. Therefore, for slow rotation rates, $\tilde{\epsilon}$ can be as small as $10^{-4}$, while faster rotations necessitate considering a larger drift parameter (e.g., $\tilde{\epsilon}_{min} \approx 0.2$ for $Ro^{-1} = 10$). In Fig. 5 the black hatched region shows the inaccessible values of $\tilde{\epsilon}$ for the shearing box approximation as a function of the inverse Rossby number. Choosing a larger (smaller) value for $v_y^2/v_K^2$ shifts the inaccessible region to smaller (larger) values of the drift parameter $\tilde{\epsilon}$.

In our exploration of the parameter space, we investigated various values of $\tilde{\epsilon}$ within the accessible range. Figure 5 presents the results of dust concentration and radial velocity variation as a function of $\tilde{\epsilon}$ and $Ro^{-1}$ for $St = 0.1$. It demonstrates that, for $\tilde{\epsilon} \le 1$, dust is ejected from eddies, and its velocity is enhanced at small $Ro^{-1}$. Conversely, dust tends to cluster more readily and experiences a decrease in radial drift for $Ro^{-1} > 0.7$. The results at $\tilde{\epsilon} = 1$ further corroborate the findings presented in this paper, as they exhibit the same behavior as those at $\tilde{\epsilon} = 0.1$, while being much farther away from the inaccessible region (see Appendix F). However, at $\tilde{\epsilon} = 10$, for which turbulence becomes a secondary mechanism compared to drift, a higher rotation rate is required for the slowdown of dust and particle concentration to happen. It is important to note that the decrease in dust radial

drift is quantitatively smaller for the same Ro and St at larger $\tilde{\epsilon}$. Therefore, a small value of the drift parameter is essential for the reduction of radial drift due to turbulence to significantly impact the long-term evolution of dust in disks.

## 5. Conclusions

In this letter we addressed the issue of dust particle drift in a turbulent gas, usually considered as a problem in theories of planetesimal formation. We demonstrated that turbulent eddies can reduce the radial drift velocity of small dust particles through two distinct mechanisms: preferential sweeping of particles between elongated eddies at slow rotation rates and concentration within anticyclones for faster rotations. In the latter scenario, dust particles may even come to a complete halt in their radial drift. Concurrently, turbulence tends to enhance the azimuthal drift of dust in most cases, although azimuthal velocities of particles are reduced when they aggregate inside anticyclonic eddies.

The reduction of dust radial drift has significant implications, potentially reconciling the masses of dust disks inferred from observations with the predicted mass loss rates from theory, particularly when turbulence prolongs the timescale of dust drift by an order of magnitude. Furthermore, our findings suggest that a turbulent region in the disk could act as a traffic jam, potentially explaining the formation of certain observable substructures, such as those detected by telescopes like ALMA. Slowing down radial drift may also facilitate planetesimal formation by assisting dust particles to overcome the radial drift barrier and by reducing collision speeds, both crucial factors for a favorable outcome of their interactions.

In our simulations we made several assumptions. The shearing box approach is relevant for our local numerical simulations, focusing specifically on dust particles dynamics. The turbulence intensity and considered scales are also small enough for density perturbations to be weak in the box (Méheut et al. 2015), therefore supporting the incompressible approximation. The choice of two-dimensional simulations is justified for small Rossby numbers, since numerous studies have indicated a two-dimensionalization of the flow at high rotation rates (Yeung & Zhou 1998; Biferale et al. 2016). Additionally, for large Ro, the 2D approximation remains valid due to the highly stratified nature of protoplanetary disks (Cambon 2001). Finally, in this study we did not consider the back-reaction of dust on gas or self-interactions between particles. While these effects can be significant when dust particles are concentrated, they are often considered secondary aspects when diluted. However, dust feedback can have strong consequences on both particle clustering (Johansen & Youdin 2007) and dust radial drift (Dipierro et al. 2018), topics that are left for future work.

Acknowledgements. This work was supported by the “Programme National de Physique Stellaire” (PNPS) and “Programme national de planétologie” (PNP) of CNRS/INSU co-funded by CEA and CNES, and Observatoire de la Côte d’Azur. This work was supported by the French government, through the UCAJEDI Investments in the Future project managed by the National Research Agency (ANR) under reference number ANR-15-IDEX-01. The authors are grateful to the OPAL infrastructure and the Université Côte d’Azur’s Center for High-Performance Computing for providing resources and support.

## References

Arlt, R., & Urpin, V. 2004, A&A, 426, 755

Balbus, S., & Hawley, J. 1991, ApJ, 376, 214

Bec, J., Gustavsson, K., & Mehlig, B. 2024, Annu. Rev. Fluid Mech., 56, 189
<figure><img src="image_6.png" /><figcaption>Fig. 5. Phase diagram of the drift parameter vs. the inverse Rossby number. The highlighted parameters present the results of the dust concentration and radial drift for $St = 0.1$.</figcaption></figure>

recent observations that tend to infer low turbulent viscosity in disks.

### 4.4. Dependence on the drift parameter $\tilde{\epsilon}$

The shearing box approximation holds true only when the box is located at a distance $r_0$ from the star much greater than the largest scales of the simulation. This requires $r_0 \gg l_f$. Consequently, we obtain the condition $\epsilon \gg (\Omega_K^2 - \Omega^2) l_f$. The acceptable values of $(\Omega_K^2 - \Omega^2)$ are determined from

$$v_y^2/v_K^2 = 1 - n c_s^2/v_K^2, \quad (6)$$

where n is the index of the power law describing pressure as a function of radius. In the literature the typical disk properties yield $0.996 < v_y^2/v_K^2 < 0.999$. For our analysis we make the arbitrary choice $v_y^2/v_K^2 = 0.998$. Therefore, for slow rotation rates, $\tilde{\epsilon}$ can be as small as $10^{-4}$, while faster rotations necessitate considering a larger drift parameter (e.g., $\tilde{\epsilon}_{min} \approx 0.2$ for $Ro^{-1} = 10$). In Fig. 5 the black hatched region shows the inaccessible values of $\tilde{\epsilon}$ for the shearing box approximation as a function of the inverse Rossby number. Choosing a larger (smaller) value for $v_y^2/v_K^2$ shifts the inaccessible region to smaller (larger) values of the drift parameter $\tilde{\epsilon}$.

In our exploration of the parameter space, we investigated various values of $\tilde{\epsilon}$ within the accessible range. Figure 5 presents the results of dust concentration and radial velocity variation as a function of $\tilde{\epsilon}$ and $Ro^{-1}$ for $St = 0.1$. It demonstrates that, for $\tilde{\epsilon} \le 1$, dust is ejected from eddies, and its velocity is enhanced at small $Ro^{-1}$. Conversely, dust tends to cluster more readily and experiences a decrease in radial drift for $Ro^{-1} > 0.7$. The results at $\tilde{\epsilon} = 1$ further corroborate the findings presented in this paper, as they exhibit the same behavior as those at $\tilde{\epsilon} = 0.1$, while being much farther away from the inaccessible region (see Appendix F). However, at $\tilde{\epsilon} = 10$, for which turbulence becomes a secondary mechanism compared to drift, a higher rotation rate is required for the slowdown of dust and particle concentration to happen. It is important to note that the decrease in dust radial

drift is quantitatively smaller for the same $Ro$ and $St$ at larger $\tilde{\epsilon}$. Therefore, a small value of the drift parameter is essential for the reduction of radial drift due to turbulence to significantly impact the long-term evolution of dust in disks.

## 5. Conclusions

In this letter we addressed the issue of dust particle drift in a turbulent gas, usually considered as a problem in theories of planetesimal formation. We demonstrated that turbulent eddies can reduce the radial drift velocity of small dust particles through two distinct mechanisms: preferential sweeping of particles between elongated eddies at slow rotation rates and concentration within anticyclones for faster rotations. In the latter scenario, dust particles may even come to a complete halt in their radial drift. Concurrently, turbulence tends to enhance the azimuthal drift of dust in most cases, although azimuthal velocities of particles are reduced when they aggregate inside anticyclonic eddies.

The reduction of dust radial drift has significant implications, potentially reconciling the masses of dust disks inferred from observations with the predicted mass loss rates from theory, particularly when turbulence prolongs the timescale of dust drift by an order of magnitude. Furthermore, our findings suggest that a turbulent region in the disk could act as a traffic jam, potentially explaining the formation of certain observable substructures, such as those detected by telescopes like ALMA. Slowing down radial drift may also facilitate planetesimal formation by assisting dust particles to overcome the radial drift barrier and by reducing collision speeds, both crucial factors for a favorable outcome of their interactions.

In our simulations we made several assumptions. The shearing box approach is relevant for our local numerical simulations, focusing specifically on dust particles dynamics. The turbulence intensity and considered scales are also small enough for density perturbations to be weak in the box (Méheut et al. 2015), therefore supporting the incompressible approximation. The choice of two-dimensional simulations is justified for small Rossby numbers, since numerous studies have indicated a two-dimensionalization of the flow at high rotation rates (Yeung & Zhou 1998; Biferale et al. 2016). Additionally, for large $Ro$, the 2D approximation remains valid due to the highly stratified nature of protoplanetary disks (Cambon 2001). Finally, in this study we did not consider the back-reaction of dust on gas or self-interactions between particles. While these effects can be significant when dust particles are concentrated, they are often considered secondary aspects when diluted. However, dust feedback can have strong consequences on both particle clustering (Johansen & Youdin 2007) and dust radial drift (Dipierro et al. 2018), topics that are left for future work.

*Acknowledgements.* This work was supported by the “Programme National de Physique Stellaire” (PNPS) and “Programme national de planétologie” (PNP) of CNRS/INSU co-funded by CEA and CNES, and Observatoire de la Côte d’Azur. This work was supported by the French government, through the UCAJEDI Investments in the Future project managed by the National Research Agency (ANR) under reference number ANR-15-IDEX-01. The authors are grateful to the OPAL infrastructure and the Université Côte d’Azur’s Center for High-Performance Computing for providing resources and support.

## References

Arlt, R., & Urpin, V. 2004, A&A, 426, 755

Balbus, S., & Hawley, J. 1991, ApJ, 376, 214

Bec, J., Gustavsson, K., & Mehlig, B. 2024, Annu. Rev. Fluid Mech., 56, 189
Biferale, L., Bonaccorso, F., Mazzitelli, I., et al. 2016, Phys. Rev. X, 6, 041036

Cambon, C. 2001, Eur. J. Mech. B Fluids, 20, 489

Cuzzi, J. N., Hogan, R. C., Paque, J. M., & Dobrovolskis, A. R. 2001, ApJ, 546, 496

Dipierro, G., Laibe, G., Alexander, R., & Hutchison, M. 2018, MNRAS, 479, 4187

Flaherty, K. M., Hughes, A. M., Rosenfeld, K. A., et al. 2015, ApJ, 813, 99

Flaherty, K. M., Hughes, A. M., Teague, R., et al. 2018, ApJ, 856, 117

Gerosa, F. A., Méheut, H., & Bec, J. 2023, Eur. Phys. J. Plus, 138, 9

Johansen, A., & Youdin, A. 2007, ApJ, 662, 627

Klahr, H., & Hubbard, A. 2014, ApJ, 788, 21

Lesur, G. R. 2021, J. Plasma Phys., 87, 205870101

Lesur, G., & Longaretti, P.-Y. 2005, A&A, 444, 25

Lesur, G., & Papaloizou, J. C. 2010, A&A, 513, A60

Lyra, W. 2014, ApJ, 789, 77

Maxey, M. R. 1987, J. Fluid Mech., 174, 441

Méheut, H., Fromang, S., Lesur, G., Joos, M., & Longaretti, P.-Y. 2015, A&A, 579, A117

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

Nelson, R. P., Gressel, O., & Umurhan, O. M. 2013, MNRAS, 435, 2610

Pan, L., Padoan, P., Scalo, J., Kritsuk, A. G., & Norman, M. L. 2011, ApJ, 740, 6

Pinilla, P., Birnstiel, T., Ricci, L., et al. 2012, A&A, 538, A114

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

Takeuchi, T., & Lin, D. 2002, ApJ, 581, 1344

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

Weidenschilling, S. 1977, MNRAS, 180, 57

Whipple, F. L. 1972, in From Plasma to Planet, ed. A. Elvius, 211

Yang, C.-C., Mac Low, M.-M., & Johansen, A. 2018, ApJ, 868, 27

Yeung, P. K., & Zhou, Y. 1998, Phys. Fluids, 10, 2895
Biferale, L., Bonaccorso, F., Mazzitelli, I., et al. 2016, Phys. Rev. X, 6, 041036

Cambon, C. 2001, Eur. J. Mech. B Fluids, 20, 489

Cuzzi, J. N., Hogan, R. C., Paque, J. M., & Dobrovolskis, A. R. 2001, ApJ, 546, 496

Dipierro, G., Laibe, G., Alexander, R., & Hutchison, M. 2018, MNRAS, 479, 4187

Flaherty, K. M., Hughes, A. M., Rosenfeld, K. A., et al. 2015, ApJ, 813, 99

Flaherty, K. M., Hughes, A. M., Teague, R., et al. 2018, ApJ, 856, 117

Gerosa, F. A., Méheut, H., & Bec, J. 2023, Eur. Phys. J. Plus, 138, 9

Johansen, A., & Youdin, A. 2007, ApJ, 662, 627

Klahr, H., & Hubbard, A. 2014, ApJ, 788, 21

Lesur, G. R. 2021, J. Plasma Phys., 87, 205870101

Lesur, G., & Longaretti, P.-Y. 2005, A&A, 444, 25

Lesur, G., & Papaloizou, J. C. 2010, A&A, 513, A60

Lyra, W. 2014, ApJ, 789, 77

Maxey, M. R. 1987, J. Fluid Mech., 174, 441

Méheut, H., Fromang, S., Lesur, G., Joos, M., & Longaretti, P.-Y. 2015, A&A, 579, A117

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

Nelson, R. P., Gressel, O., & Umurhan, O. M. 2013, MNRAS, 435, 2610

Pan, L., Padoan, P., Scalo, J., Kritsuk, A. G., & Norman, M. L. 2011, ApJ, 740, 6

Pinilla, P., Birnstiel, T., Ricci, L., et al. 2012, A&A, 538, A114

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

Takeuchi, T., & Lin, D. 2002, ApJ, 581, 1344

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

Weidenschilling, S. 1977, MNRAS, 180, 57

Whipple, F. L. 1972, in From Plasma to Planet, ed. A. Elvius, 211

Yang, C.-C., Mac Low, M.-M., & Johansen, A. 2018, ApJ, 868, 27

Yeung, P. K., & Zhou, Y. 1998, Phys. Fluids, 10, 2895
# Appendix A: Run parameters

<table><caption>**Table A.1.** Simulation parameters (in simulation units).</caption><thead><tr><th>L<sub>x</sub> × L<sub>y</sub></th><th>N<sub>x</sub> × N<sub>y</sub></th><th>ε<sub>1</sub></th><th>η<sub>1</sub></th><th>ν</th><th>γ</th></tr></thead><tbody><tr><td>2π × 8π</td><td>256 × 1024</td><td>1.38</td><td>22.12</td><td>0.005</td><td>0.05</td></tr></tbody></table>

<table><caption>**Table A.2.** Variation in the eddy turnover time and Rossby number as a function of rotation (in simulation units).</caption><thead><tr><th>Ω</th><th>0.7</th><th>2.7</th><th>5.3</th><th>8.7</th><th>10.7</th><th>13.3</th><th>21.3</th></tr></thead><tbody><tr><th>τ<sub>ω</sub></th><td>0.16</td><td>0.18</td><td>0.20</td><td>0.22</td><td>0.23</td><td>0.24</td><td>0.29</td></tr><tr><th>Ro<sup>-1</sup></th><td>0.10</td><td>0.49</td><td>1.09</td><td>1.92</td><td>2.49</td><td>3.29</td><td>6.19</td></tr></tbody></table>

Table A.1 presents the parameters used in the direct numerical simulations with the spectral code Snoopy.

We conducted seven distinct runs, each time varying the rotation rate $\Omega$, which resulted in specific values for the eddy turnover time $\tau_{\omega}$ (as listed in Table A.2), and, consequently, for the Rossby number.

Finally, each of these runs is seeded with 48 sets of $10^4$ Lagrangian particles. These families are associated with 4 different values of $\epsilon$ (0.1, 1, 10, 100) and 12 values of $\tau_p$ (0.008, 0.014, 0.025, 0.044, 0.080, 0.144, 0.260, 0.470, 0.849, 1.533, 2.769, and 5.000).

# Appendix B: Shearing box and drift

We consider dust particles in Keplerian rotation, thus rotating with a rotation rate $\Omega_K$ faster than the gas, which rotates with the shearing box at a sub-Keplerian speed $\Omega$. Considering the box centered at a distance $r_0$ from the star, the equation for passive particles is

$$ \frac{dV_p}{dt} = -\frac{1}{\tau_p} [V_p - v(\mathbf{R}_p, t)] - 2\Omega \times V_p + (2\Omega_K^2 + \Omega^2) R_{p,x} e_x - (\Omega_K^2 - \Omega^2)r_0 e_x. \quad (B.1) $$

For numerical reasons, it is easier to compute the particle velocity in the frame where the background negative shear has been subtracted. Therefore, we perform the change of variables $U_{p,y} = V_{p,y} + (3/2)\Omega x_p$ giving

$$ \frac{dU_{p,x}}{dt} = -\frac{1}{\tau_p} [U_{p,x} - u_x(\mathbf{R}_p, t)] + 2\Omega U_{p,y} + (2\Omega_K^2 + \Omega^2) R_{p,x} - (\Omega_K^2 - \Omega^2)r_0 - 3\Omega^2 R_{p,x}, \quad (B.2) $$

$$ \frac{dU_{p,y}}{dt} = -\frac{1}{\tau_p} [U_{p,y} - u_y(\mathbf{R}_p, t)] - \frac{1}{2} \Omega U_{p,x}. $$

Finally, if $R_{p,x} \ll r_0$, we obtain $2(\Omega_K^2 - \Omega^2)x_p \ll (\Omega_K^2 - \Omega^2)r_0$, so that

$$ \frac{dU_p}{dt} = -\frac{1}{\tau_p} [U_p - u(\mathbf{R}_p, t)] - 2\Omega \times U_p + \frac{3}{2} \Omega U_{p,x} e_y - \epsilon e_x, \quad (B.3) $$

with $\epsilon = (\Omega_K^2 - \Omega^2)r_0$.

<figure><img src="image_8.png" /><figcaption>**Fig. C.1.** Vorticity $\omega$, normalized by $\tau_{\omega}$, and dust particle position $\mathbf{R}_p$ (black dots) for $\tilde{\epsilon} = 0.1$, $\text{St} = 1$, $\text{Ro}^{-1} = 3.29$ and a) $t=200\tau_{\omega}$, b) $t=2000\tau_{\omega}$, c) $t=20000\tau_{\omega}$. The clusters in b) and c) are highlighted with pink circles.</figcaption></figure>

# Appendix C: Dust clustering

The clustering properties of particles in a turbulent Keplerian flow can vary significantly depending on the flow and dust parameters. In a previous paper where no drift was considered (Gerosa et al. 2023), we identified three possible turbulent concentration paths:

1. Ejection from eddies and concentration on filamentary structures in between, for small rotation rates and dust sizes.

2. Distribution in spirals inside the anticyclones, but avoiding their core, for intermediate rotation rates.

3. Formation of point clusters in the cores of anticyclones, for large rotation rates and small dust sizes.

When the drift is considered, the picture remains qualitatively the same. As an example, we show in Fig. C.1 the concentration of dust particles over time, for $\tilde{\epsilon} = 0.1$, $\text{St} = 1$, and $\text{Ro}^{-1} = 3.29$. Due to the Coriolis force, dust initially gathers in anticyclonic eddies (Fig. C.1a). Over time, each of these clouds coalesces into a point cluster in the core of an anticyclone (Fig. C.1b). Clusters eventually merge to form a single cluster (Fig. C.1c). In the extreme case of $\tilde{\epsilon} = 10$, however, the process of dust clustering is halted by drift, which carries particles away from anticyclones on a timescale too fast for dust to concentrate efficiently within them.

# Appendix D: Okubo-Weiss parameter

The Okubo-Weiss parameter estimates the properties of a gas flow, helping to identify its local structures. It is defined as

$$ \mathrm{OW} = -(\partial_x u_x)^2 - \partial_x u_y \partial_y u_x. \quad (\mathrm{D.1}) $$

When positive, it describes rotating elliptic structures of the flow (e.g., eddies). For straining hyperbolic regions, the Okubo-Weiss parameter instead assumes negative values. This parameter has been measured along the particles path and its averaged value, as a function of $\text{Ro}^{-1}$ for various Stokes number, is shown in Fig. D.1. We note that $\langle \mathrm{OW}(\mathbf{R}_p) \rangle$ assumes negative values at small
# Appendix A: Run parameters

<table><caption>**Table A.1.** Simulation parameters (in simulation units).</caption><thead><tr><th>L<sub>x</sub> × L<sub>y</sub></th><th>N<sub>x</sub> × N<sub>y</sub></th><th>ε<sub>1</sub></th><th>η<sub>1</sub></th><th>ν</th><th>γ</th></tr></thead><tbody><tr><td>2π × 8π</td><td>256 × 1024</td><td>1.38</td><td>22.12</td><td>0.005</td><td>0.05</td></tr></tbody></table>

<table><caption>**Table A.2.** Variation in the eddy turnover time and Rossby number as a function of rotation (in simulation units).</caption><thead><tr><th>Ω</th><th>0.7</th><th>2.7</th><th>5.3</th><th>8.7</th><th>10.7</th><th>13.3</th><th>21.3</th></tr></thead><tbody><tr><th>τ<sub>ω</sub></th><td>0.16</td><td>0.18</td><td>0.20</td><td>0.22</td><td>0.23</td><td>0.24</td><td>0.29</td></tr><tr><th>Ro<sup>-1</sup></th><td>0.10</td><td>0.49</td><td>1.09</td><td>1.92</td><td>2.49</td><td>3.29</td><td>6.19</td></tr></tbody></table>

Table A.1 presents the parameters used in the direct numerical simulations with the spectral code Snoopy.

We conducted seven distinct runs, each time varying the rotation rate $\Omega$, which resulted in specific values for the eddy turnover time $\tau_{\omega}$ (as listed in Table A.2), and, consequently, for the Rossby number.

Finally, each of these runs is seeded with 48 sets of $10^4$ Lagrangian particles. These families are associated with 4 different values of $\epsilon$ (0.1, 1, 10, 100) and 12 values of $\tau_p$ (0.008, 0.014, 0.025, 0.044, 0.080, 0.144, 0.260, 0.470, 0.849, 1.533, 2.769, and 5.000).

# Appendix B: Shearing box and drift

We consider dust particles in Keplerian rotation, thus rotating with a rotation rate $\Omega_K$ faster than the gas, which rotates with the shearing box at a sub-Keplerian speed $\Omega$. Considering the box centered at a distance $r_0$ from the star, the equation for passive particles is

$$ \frac{dV_p}{dt} = -\frac{1}{\tau_p} [V_p - v(\mathbf{R}_p, t)] - 2\Omega \times V_p + (2\Omega_K^2 + \Omega^2) R_{p,x} e_x - (\Omega_K^2 - \Omega^2)r_0 e_x. \quad (B.1) $$

For numerical reasons, it is easier to compute the particle velocity in the frame where the background negative shear has been subtracted. Therefore, we perform the change of variables $U_{p,y} = V_{p,y} + (3/2)\Omega x_p$ giving

$$ \frac{dU_{p,x}}{dt} = -\frac{1}{\tau_p} [U_{p,x} - u_x(\mathbf{R}_p, t)] + 2\Omega U_{p,y} + (2\Omega_K^2 + \Omega^2) R_{p,x} - (\Omega_K^2 - \Omega^2)r_0 - 3\Omega^2 R_{p,x}, \quad (B.2) $$

$$ \frac{dU_{p,y}}{dt} = -\frac{1}{\tau_p} [U_{p,y} - u_y(\mathbf{R}_p, t)] - \frac{1}{2} \Omega U_{p,x}. $$

Finally, if $R_{p,x} \ll r_0$, we obtain $2(\Omega_K^2 - \Omega^2)x_p \ll (\Omega_K^2 - \Omega^2)r_0$, so that

$$ \frac{dU_p}{dt} = -\frac{1}{\tau_p} [U_p - u(\mathbf{R}_p, t)] - 2\Omega \times U_p + \frac{3}{2}\Omega U_{p,x} e_y - \epsilon e_x, \quad (B.3) $$

with $\epsilon = (\Omega_K^2 - \Omega^2)r_0$.

<figure><img src="image_7.png" /><figcaption>**Fig. C.1.** Vorticity $\omega$, normalized by $\tau_{\omega}$, and dust particle position $\mathbf{R}_p$ (black dots) for $\tilde{\epsilon} = 0.1$, St = 1, $Ro^{-1} = 3.29$ and a) $t=200\tau_{\omega}$, b) $t=2000\tau_{\omega}$, c) $t=20000\tau_{\omega}$. The clusters in b) and c) are highlighted with pink circles.</figcaption></figure>

# Appendix C: Dust clustering

The clustering properties of particles in a turbulent Keplerian flow can vary significantly depending on the flow and dust parameters. In a previous paper where no drift was considered (Gerosa et al. 2023), we identified three possible turbulent concentration paths:

1. Ejection from eddies and concentration on filamentary structures in between, for small rotation rates and dust sizes.

2. Distribution in spirals inside the anticyclones, but avoiding their core, for intermediate rotation rates.

3. Formation of point clusters in the cores of anticyclones, for large rotation rates and small dust sizes.

When the drift is considered, the picture remains qualitatively the same. As an example, we show in Fig. C.1 the concentration of dust particles over time, for $\tilde{\epsilon} = 0.1$, St = 1, and $Ro^{-1} = 3.29$. Due to the Coriolis force, dust initially gathers in anticyclonic eddies (Fig. C.1a). Over time, each of these clouds coalesces into a point cluster in the core of an anticyclone (Fig. C.1b). Clusters eventually merge to form a single cluster (Fig. C.1c). In the extreme case of $\tilde{\epsilon} = 10$, however, the process of dust clustering is halted by drift, which carries particles away from anticyclones on a timescale too fast for dust to concentrate efficiently within them.

# Appendix D: Okubo-Weiss parameter

The Okubo-Weiss parameter estimates the properties of a gas flow, helping to identify its local structures. It is defined as

$$ OW = -(\partial_x u_x)^2 - \partial_x u_y \partial_y u_x. \quad (D.1) $$

When positive, it describes rotating elliptic structures of the flow (e.g., eddies). For straining hyperbolic regions, the Okubo-Weiss parameter instead assumes negative values. This parameter has been measured along the particles path and its averaged value, as a function of $Ro^{-1}$ for various Stokes number, is shown in Fig. D.1. We note that $\langle OW(\mathbf{R}_p) \rangle$ assumes negative values at small
<figure><img src="image_9.png" /><figcaption>Fig. D.1. Okubo-Weiss parameter computed at the particle position as a function of the inverse Rossby number, for various Stokes numbers and $\tilde{\epsilon} = 0.1$</figcaption></figure>

$Ro^{-1}$, while it become positive at high rotation rates. This means that, independently of the Stokes number, particles in a slow rotating gas flow predominantly lie between eddies. Instead, for $Ro^{-1} \ge 0.7$, dust resides longer inside anticyclonic eddies.

## Appendix E: Mean drift velocity

An expression for the mean particle velocity can be obtained by averaging Eq. (B.3) assuming statistical stationarity of turbulent fluctuations along particle paths. This leads to the following expressions for the mean radial and azimuthal particle velocities:

$$ \begin{aligned} \langle U_{p,x} \rangle &= \frac{\langle u_x(\mathbf{R}_p, t) \rangle + 2\tau_p\Omega \langle u_y(\mathbf{R}_p, t) \rangle - \tau_p\epsilon}{1 + \tau_p^2\Omega^2}, \\ \langle U_{p,y} \rangle &= -\frac{\tau_p\Omega \langle u_x(\mathbf{R}_p, t) \rangle + \langle u_y(\mathbf{R}_p, t) \rangle - \tau_p^2\Omega\epsilon}{2(1 + \tau_p^2\Omega^2)}. \end{aligned} \quad (E.1) $$

In the absence of gas turbulence ($u = 0$), these equations yield the mean particle velocity $U_p^{NT}$ in a laminar flow.

It should be noted that without drift ($\epsilon = 0$), even if the gas flow is turbulent, the average gas velocity at particle position $\langle u(\mathbf{R}_p, t) \rangle$ would vanish, leading to $\langle U_p \rangle = 0$. On the other hand, when drift is present, as in our case, the mean gas velocity at the particle position does not vanish, highlighting the phenomenon of preferential sweeping. This mechanism indeed arises from a preferential sampling of specific zones and structures of the flow along particle paths, modifying the mean radial and azimuthal drifts of the particles.

## Appendix F: Radial drift at $\tilde{\epsilon} = 1$

Figure F.1 shows the relative increase in dust radial drift due to turbulence $\Delta U_{p,x}/U_{p,x}^{NT} = (\langle U_{p,x} \rangle - U_{p,x}^{NT})/U_{p,x}^{NT}$ for $\tilde{\epsilon} = 1$. We note that the results closely resemble those presented in Fig. 1 for $\tilde{\epsilon} = 0.1$. A subtle reduction of radial drift can be noted at $Ro^{-1} = 6.19$ and small St compared to $Ro^{-1} = 3.29$, indicating a reversal in trend. In this case, turbulence is too weak to efficiently counteract dust drift.

<figure><img src="image_10.png" /><figcaption>Fig. F.1. Variation in dust radial drift velocity due to turbulence for $\tilde{\epsilon} = 1$ as a function of the inverse Rossby number and the Stokes number.</figcaption></figure>
<figure><img src="image_8.png" /><figcaption>Fig. D.1. Okubo-Weiss parameter computed at the particle position as a function of the inverse Rossby number, for various Stokes numbers and $\tilde{\epsilon} = 0.1$</figcaption></figure>

$Ro^{-1}$, while it become positive at high rotation rates. This means that, independently of the Stokes number, particles in a slow rotating gas flow predominantly lie between eddies. Instead, for $Ro^{-1} \ge 0.7$, dust resides longer inside anticyclonic eddies.

## Appendix E: Mean drift velocity

An expression for the mean particle velocity can be obtained by averaging Eq. (B.3) assuming statistical stationarity of turbulent fluctuations along particle paths. This leads to the following expressions for the mean radial and azimuthal particle velocities:

$$ \begin{aligned} \langle U_{p,x} \rangle &= \frac{\langle u_x(\mathbf{R}_p, t) \rangle + 2\tau_p\Omega \langle u_y(\mathbf{R}_p, t) \rangle - \tau_p\epsilon}{1 + \tau_p^2\Omega^2}, \\ \langle U_{p,y} \rangle &= -\frac{\tau_p\Omega \langle u_x(\mathbf{R}_p, t) \rangle + \langle u_y(\mathbf{R}_p, t) \rangle - \tau_p^2\Omega\epsilon}{2(1 + \tau_p^2\Omega^2)}. \end{aligned} \quad (E.1) $$

In the absence of gas turbulence ($u = 0$), these equations yield the mean particle velocity $U_p^{NT}$ in a laminar flow.

It should be noted that without drift ($\epsilon = 0$), even if the gas flow is turbulent, the average gas velocity at particle position $\langle u(\mathbf{R}_p, t) \rangle$ would vanish, leading to $\langle U_p \rangle = 0$. On the other hand, when drift is present, as in our case, the mean gas velocity at the particle position does not vanish, highlighting the phenomenon of preferential sweeping. This mechanism indeed arises from a preferential sampling of specific zones and structures of the flow along particle paths, modifying the mean radial and azimuthal drifts of the particles.

## Appendix F: Radial drift at $\tilde{\epsilon} = 1$

Figure F.1 shows the relative increase in dust radial drift due to turbulence $\Delta U_{p,x}/U_{p,x}^{NT} = (\langle U_{p,x} \rangle - U_{p,x}^{NT})/U_{p,x}^{NT}$ for $\tilde{\epsilon} = 1$. We note that the results closely resemble those presented in Fig. 1 for $\tilde{\epsilon} = 0.1$. A subtle reduction of radial drift can be noted at $Ro^{-1} = 6.19$ and small St compared to $Ro^{-1} = 3.29$, indicating a reversal in trend. In this case, turbulence is too weak to efficiently counteract dust drift.

<figure><img src="image_9.png" /><figcaption>Fig. F.1. Variation in dust radial drift velocity due to turbulence for $\tilde{\epsilon} = 1$ as a function of the inverse Rossby number and the Stokes number.</figcaption></figure>