A new code for astrophysical magnetohydrodynamics (MHD) is described. The code has been designed to be easily extensible for use with static and adaptive mesh refinement. It combines higher-order Godunov methods with the constrained transport (CT) technique to enforce the divergence-free constraint on the magnetic field. Discretization is based on cell-centered volume-averages for mass, momentum, and energy, and face-centered area-averages for the magnetic field. Novel features of the algorithm include (1) a consistent framework for computing the time- and edge-averaged electric fields used by CT to evolve the magnetic field from the time- and area-averaged Godunov fluxes, (2) the extension to MHD of spatial reconstruction schemes that involve a dimensionally-split time advance, and (3) the extension to MHD of two different dimensionally-unsplit integration methods. Implementation of the algorithm in both C and Fortran95 is detailed, including strategies for parallelization using domain decomposition. Results from a test suite which includes problems in one-, two-, and three-dimensions for both hydrodynamics and MHD are given, not only to demonstrate the fidelity of the algorithms, but also to enable comparisons to other methods. The source code is freely available for download on the web.Comment: 61 pages, 36 figures. accepted by ApJ

We study the formation of planetesimals in protoplanetary disks from the gravitational collapse of solid over-densities generated via the streaming instability. To carry out these studies, we implement and test a particle-mesh self-gravity module for the Athena code that enables the simulation of aerodynamically coupled systems of gas and collisionless self-gravitating solid particles. Upon employment of our algorithm to planetesimal formation simulations, we find that (when a direct comparison is possible) the Athena simulations yield predicted planetesimal properties that agree well with those found in prior work using different numerical techniques. In particular, the gravitational collapse of streaming-initiated clumps leads to an initial planetesimal mass function that is well-represented by a power-law, dN/dM p ∝ M −p p , with p ≃ 1.6 ± 0.1, which equates to a differential size distribution dN/dR p ∝ R −q p , with q ≃ 2.8 ± 0.1. We find no significant trends with resolution from a convergence study of up to 512 3 grid zones and N par ≈ 1.5 × 10 8 particles. Likewise, the power-law slope appears indifferent to changes in the relative strength of self-gravity and tidal shear, and to the time when (for reasons of numerical economy) self-gravity is turned on, though the strength of these claims is limited by small number statistics. For a typically assumed radial distribution of minimum mass solar nebula solids (assumed here to have dimensionless stopping time τ = 0.3), our results support the hypothesis that bodies on the scale of large asteroids or Kuiper Belt Objects could have formed as the high-mass tail of a primordial planetesimal population.

Turbulence can transport angular momentum in protoplanetary disks and influence the growth and evolution of planets. With spatially and spectrally resolved molecular emission line measurements provided by (sub)millimeter interferometric observations, it is possible to directly measure non-thermal motions in the disk gas that can be attributed to this turbulence. We report a new constraint on the turbulence in the disk around HD 163296, a nearby young A star, determined from ALMA Science Verification observations of four CO emission lines (the CO(3-2), CO(2-1), 13 CO(2-1), and C 18 O(2-1) transitions). The different optical depths for these lines permit probes of non-thermal line-widths at a range of physical conditions (temperature and density) and depths into the disk interior. We derive stringent limits on the non-thermal motions in the upper layers of the outer disk such that any contribution to the line-widths from turbulence is <3% of the local sound speed. These limits are approximately an order of magnitude lower than theoretical predictions for full-blown MHD turbulence driven by the magneto-rotational instability, potentially suggesting that this mechanism is less efficient in the outer (R 30AU) disk than has been previously considered.

Gas kinematics are an important part of the planet formation process. Turbulence influences planetesimal growth and migration from the scale of sub-micron dust grains through gas-giant planets. Radio observations of resolved molecular line emission can directly measure this non-thermal motion and, taking advantage of the layered chemical structure of disks, different molecular lines can be combined to map the turbulence throughout the vertical extent of a protoplanetary disk. Here we present ALMA observations of three molecules (DCO + (3-2), C 18 O(2-1) and CO(2-1)) from the disk around HD 163296. We are able to place stringent upper limits (v turb <0.06c s , <0.05c s and <0.04c s for CO(2-1), C 18 O(2-1) and DCO + (3-2) respectively), corresponding to α 3×10 −3 , similar to our prior limit derived from CO(3-2). This indicates that there is little turbulence throughout the vertical extent of the disk, contrary to theoretical predictions based on the magneto-rotational instability and gravito-turbulence. In modeling the DCO + emission we also find that it is confined to three concentric rings at 65.7±0.9 au, 149.9 +0.5 −0.7 au and 259±1 au, indicative of a complex chemical environment.

We study how the structure and variability of magnetohydrodynamic (MHD) turbulence in accretion discs converge with domain size. Our results are based on a series of vertically stratified local simulations, computed using the athena MHD code, that have fixed spatial resolution, but varying radial and azimuthal extent (from ΔR= 0.5H to 16H, where H is the vertical scale height). We show that elementary local diagnostics of the turbulence, including the Shakura–Sunyaev α parameter, the ratio of Maxwell stress to magnetic energy and the ratio of magnetic to fluid stresses, converge to within the precision of our measurements for spatial domains of radial size Lx≥ 2H. We obtain α≃ 0.02–0.03, consistent with other recent determinations. Very small domains (Lx= 0.5H) return anomalous results, independent of spatial resolution. This convergence with domain size, however, is only valid for a limited set of diagnostics: larger spatial domains admit the emergence of dynamically important mesoscale structures. In our largest simulations, the Maxwell stress shows a significant large‐scale non‐local component, while the density develops long‐lived axisymmetric perturbations (‘zonal flows’) at the 20 per cent level. Most strikingly, the variability of the disc in fixed‐sized patches decreases strongly as the simulation volume increases, while variability in the magnetically dominated corona remains constant. Comparing our largest local simulations to global simulations with comparable spatial resolution, we find generally good agreement. There is no direct evidence that the presence of curvature terms or radial gradients in global calculations materially affect the turbulence, except to perhaps introduce an outer radial scale for mesoscale structures. The demonstrated importance of mean magnetic fields – seen in both large local and global simulations – implies, however, that the growth and saturation of these fields is likely of critical importance for the evolution of accretion discs.

We investigate the effect of shear viscosity, ν, and Ohmic resistivity, η on the magnetorotational instability (MRI) in vertically stratified accretion disks through a series of local simulations with the Athena code. First, we use a series of unstratified simulations to calibrate physical dissipation as a function of resolution and background field strength; the effect of the magnetic Prandtl number, P m = ν/η, on the turbulence is captured by ∼ 32 grid zones per disk scale height, H. In agreement with previous results, our stratified disk calculations are characterized by a subthermal, predominately toroidal magnetic field that produces MRI-driven turbulence for |z| 2H. Above |z| ∼ 2H, magnetic pressure dominates and the field is buoyantly unstable. Large scale radial and toroidal fields are also generated near the mid-plane and subsequently rise through the disk. The polarity of this mean field switches on a roughly 10 orbit period in a process that is well-modeled by an α-Ω dynamo. Turbulent stress increases with P m but with a shallower dependence compared to unstratified simulations. For sufficiently large resistivity, η ∼ c s H/1000 where c s is the sound speed, MRI turbulence within 2H of the mid-plane undergoes periods of resistive decay followed by regrowth. This regrowth is caused by amplification of toroidal field via the dynamo. This process results in large amplitude variability in the stress on 10 to 100 orbital timescales, which may have relevance for partially ionized disks that are observed to have high and low accretion states.

Turbulence is a fundamental parameter in models of grain growth during the early stages of planet formation. As such, observational constraints on its magnitude are crucial. Here we self-consistently analyze ALMA CO(2-1), SMA CO(3-2), and SMA CO(6-5) observations of the disk around TW Hya and find an upper limit on the turbulent broadening of <0.08c s (α <0.007 for α defined only within 2-3 pressure scale heights above the midplane), lower than the tentative detection previously found from an analysis of the CO(2-1) data. We examine in detail the challenges of image plane fitting vs directly fitting the visibilities, while also considering the role of the vertical temperature gradient, systematic uncertainty in the amplitude calibration, and assumptions about the CO abundance, as potential sources of the discrepancy in the turbulence measurements. These tests result in variations of the turbulence limit between <0.04c s and <0.13c s , consistently lower than the 0.2-0.4c s found previously. Having ruled out numerous factors, we restrict the source of the discrepancy to our assumed coupling between temperature and density through hydrostatic equilibrium in the presence of a vertical temperature gradient and/or the confinement of CO to a thin molecular layer above the midplane, although further work is needed to quantify the influence of these prescriptions. Assumptions about hydrostatic equilibrium and the CO distribution are physically motivated, and may have a small influence on measuring the kinematics of the gas, but they become important when constraining small effects such as the strength of the turbulence within a protoplanetary disk.

We numerically evolve turbulence driven by the magnetorotational instability (MRI) in a 3D, unstratified shearing box and study its structure using two-point correlation functions. We confirm Fromang & Papaloizou's result that shearing box models with zero net magnetic flux are not converged; the dimensionless shear stress α is proportional to the grid scale. We find that the two-point correlation of B shows that it is composed of narrow filaments that are swept back by differential rotation into a trailing spiral. The correlation lengths along each of the correlation function principal axes decrease monotonically with the grid scale. For mean azimuthal field models, which we argue are more relevant to astrophysical disks than the zero net field models, we find that: α increases weakly with increasing resolution at fixed box size; α increases slightly as the box size is increased; α increases linearly with net field strength, confirming earlier results; the two-point correlation function of the magnetic field is resolved and converged, and is composed of narrow filaments swept back by the shear; the major axis of the two-point increases slightly as the box size is increased; these results are code independent, based on a comparison of ATHENA and ZEUS runs. The velocity, density, and magnetic fields decorrelate over scales larger than ∼ H, as do the dynamical terms in the magnetic energy evolution equations. We conclude that MHD turbulence in disks is localized, subject to the limitations imposed by the absence of vertical stratification, the use of an isothermal equation of state, finite box size, finite run time, and finite resolution.