Turbulent diffusion of large solids in a protoplanetary disc

We study the turbulent diffusion of solids in a protoplanetary disc, in order to discriminate between two existing analytical models of the turbulent diffusion process. These two models predict the same radial turbulent diffusion coefficient D-p,D-x for small particles (tau(s) << 1), but differ in the value of D-p,D-x for large particles (tau(s) >> 1, where tau(s) is the dimensionless particle stopping time, closely related to particle radius). The model given by Youdin & Lithwick (YL) takes into account orbital oscillations of the solids, while the other model given by Cuzzi, Dobrovolskis & Champney (CDC) does not. The CDC model predicts D-p,D-x similar to tau(-1)(s) for tau(s) >> 1, but the YL model gives D-p,D-x similar to tau(-2)(s). To investigate, we perform 3D, magnetohydrodynamic (MHD) numerical simulations. Turbulence in the disc is generated by the magnetorotational instability. The ATHENA code is used to solve the equations of ideal MHD in the shearing-box approximation, which allows us to model a local region of the disc with the relevant orbital dynamics. Solids are represented by Lagrangian particles that interact with the gas through drag, and are also subject to orbital forces. The aerodynamic coupling of particles to the gas is parametrized by tau(s). In one set of simulations, particle displacements along the radial direction are measured in a shearing box without vertical stratification of the gas density. In another simulation, the vertical component of stellar gravity is included, with a Gaussian gas density vertical profile, but the particle motion is restricted to fixed planes of constant height z. In both cases, the radial diffusion coefficient as a function of stopping time tau(s) is in very good agreement with the YL model. To study particle vertical diffusion, we use the unstratified shearing box, in which we allow the effects of vertical gravity and turbulence on the particles to balance out, resulting in particle layers whose scaleheight varies approximately as tau(-1/2)(s). Based on this result and YL, we calculate a vertical diffusion coefficient D-p,D-z that, in the limit tau(s) >> 1, varies as tau(-2)(s), similarly to radial diffusivity.