Bryan M. Johnson

Nonlinear Outcome of Gravitational Instability in Disks with Realistic Cooling

We consider the nonlinear outcome of gravitational instability in optically thick disks with a realistic cooling function. We use a numerical model that is local, razor thin, and unmagnetized. External illumination is ignored. Cooling is calculated from a one-zone model using analytic fits to low-temperature Rosseland mean opacities. The model has two parameters: the initial surface density Σ0 and the rotation frequency Ω. We survey the parameter space and find the following. (1) The disk fragments when τcΩ ~ 1, where τc is an effective cooling time defined as the average internal energy of the model divided by the average cooling rate. This is consistent with earlier results that used a simplified cooling function. (2) The initial cooling time τc0 for a uniform disk with Q = 1 can differ by orders of magnitude from τc in the nonlinear outcome. The difference is caused by sharp variations in the opacity with temperature. The condition τc0Ω ~ 1 therefore does not necessarily indicate where fragmentation will occur. (3) The largest difference between τc and τc0 is near the opacity gap, where dust is absent and hydrogen is largely molecular. (4) In the limit of strong illumination the disk is isothermal; we find that an isothermal version of our model fragments for Q 1.4. Finally, we discuss some physical processes not included in our model and find that most are likely to make disks more susceptible to fragmentation. We conclude that disks with τcΩ 1 do not exist.

LOCALITY OF MHD TURBULENCE IN ISOTHERMAL DISKS

We numerically evolve turbulence driven by the magnetorotational instability (MRI) in a three-dimensional, unstratified shearing box and study its structure using two-point correlation functions. We confirm Fromang and Papaloizous 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 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.

Vortices in Thin, Compressible, Unmagnetized Disks

We consider the formation and evolution of vortices in a hydrodynamic shearing-sheet model. The evolution is done numerically using a version of the ZEUS code. Consistent with earlier results, an injected vorticity field evolves into a set of long-lived vortices, each of which has a radial extent comparable to the local scale height. But we also find that the resulting velocity field has a positive shear stress, Σδvrδv. This effect appears only at high resolution. The transport, which decays with time as t-1/2, arises primarily because the vortices drive compressive motions. This result suggests a possible mechanism for angular momentum transport in low-ionization disks, with two important caveats: a mechanism must be found to inject vorticity into the disk, and the vortices must not decay rapidly due to three-dimensional instabilities.

Orbital Advection by Interpolation: A Fast and Accurate Numerical Scheme for Super-Fast MHD Flows

In numerical models of thin astrophysical disks that use a Eulerian scheme, gas orbits supersonically through a fixed grid. As a result the time step is sharply limited by the Courant condition. Also, because the mean flow speed with respect to the grid varies with position, the truncation error varies systematically with position. For hydrodynamic (unmagnetized) disks an algorithm called FARGO has been developed that advects the gas along its mean orbit using a separate interpolation substep. This relaxes the constraint imposed by the Courant condition, which now depends only on the peculiar velocity of the gas, and results in a truncation error that is more nearly independent of position. This paper describes a FARGO-like algorithm suitable for evolving magnetized disks. Our method is second-order-accurate on a smooth flow and preserves ∇B = 0 to machine precision. The main restriction is that B must be discretized on a staggered mesh. We give a detailed description of an implementation of the code and demonstrate that it produces the expected results on linear and nonlinear problems. We also point out how the scheme might be generalized to make the integration of other supersonic/superfast flows more efficient. Although our scheme reduces the variation of truncation error with position, it does not eliminate it. We show that the residual position dependence leads to characteristic radial variations in the density over long integrations.

Analytical shock solutions at large and small Prandtl number

Exact one-dimensional solutions to the equations of fluid dynamics are derived in the large-Pr and small-Pr limits (where Pr is the Prandtl number). The solutions are analogous to the Pr = 3/4 solution discovered by Becker and analytically capture the profile of shock fronts in ideal gases. The large-Pr solution is very similar to Beckers solution, differing only by a scale factor. The small-Pr solution is qualitatively different, with an embedded isothermal shock occurring above a critical Mach number. Solutions are derived for constant viscosity and conductivity as well as for the case in which conduction is provided by a radiation field. For a completely general density- and temperature-dependent viscosity and conductivity, the system of equations in all three limits can be reduced to quadrature. The maximum error in the analytical solutions when compared to a numerical integration of the finite-Pr equations is O(1/Pr) for large Pr and O(Pr) for small Pr.

Closed-Form Shock Solutions

It is shown here that a subset of the implicit analytical shock solutions discovered by Becker and by Johnson can be inverted, yielding several exact closed-form solutions of the one-dimensional compressible Navier-Stokes equations for an ideal gas. For a constant dynamic viscosity and thermal conductivity, and at particular values of the shock Mach number, the velocity can be expressed in terms of a polynomial root. For a constant kinematic viscosity, independent of Mach number, the velocity can be expressed in terms of a hyperbolic tangent function. The remaining fluid variables are related to the velocity through simple algebraic expressions. The solutions derived here make excellent verification tests for numerical algorithms, since no source terms in the evolution equations are approximated, and the closed-form expressions are straightforward to implement. The solutions are also of some academic interest as they may provide insight into the non-linear character of the Navier-Stokes equations and may stimulate further analytical developments.

Addendum: “Orbital Advection by Interpolation: A Fast and Accurate Numerical Scheme for Super-Fast MHD Flows” (ApJS, 177, 373 [2008])

The descriptions of some of the numerical tests in our original paper are incomplete, making reproduction of the results difficult. We provide the missing details here. The relevant tests are described in x 4 of the original paper (Figs. 8Y11). We use the analytical solutions outlined by B.M. Johnson (ApJ, 660, 1375 [2007]) as the initial conditions for the linear tests in B.M. Johnson et al. (ApJS, 177, 373, [2008]). The incompressive solution is given by the real parts of expressions (80)Y(82) of that paper. For imaginary ! and ! and a Keplerian rotation profile, these are

On the interaction between turbulence and a planar rarefaction

The modeling of turbulence, whether it be numerical or analytical, is a difficult challenge. Turbulence is amenable to analysis with linear theory if it is subject to rapid distortions, i.e., motions occurring on a timescale that is short compared to the timescale for nonlinear interactions. Such an approach (referred to as rapid distortion theory) could prove useful for understanding aspects of astrophysical turbulence, which is often subject to rapid distortions, such as supernova explosions or the free-fall associated with gravitational instability. As a proof of principle, a particularly simple problem is considered here: the evolution of vorticity due to a planar rarefaction in an ideal gas. Analytical solutions are obtained for incompressive modes having a wave vector perpendicular to the distortion; as in the case of gradient-driven instabilities, these are the modes that couple most strongly to the mean flow. Vorticity can either grow or decay in the wake of a rarefaction front, and there are two competing effects that determine which outcome occurs: entropy fluctuations couple to the mean pressure gradient to produce vorticity via baroclinic effects, whereas vorticity is damped due to the conservation of angular momentum as the fluid expands. Whether vorticity grows or decays depends upon the ratio of entropic to vortical fluctuations at the location of the front; growth occurs if this ratio is of order unity or larger. In the limit of purely entropic fluctuations in the ambient fluid, a strong rarefaction generates vorticity with a turbulent Mach number on the order of the rms of the ambient entropy fluctuations. The analytical results are shown to compare well with results from two- and three-dimensional numerical simulations. Analytical solutions are also derived in the linear regime of Reynolds-averaged turbulence models. This highlights an inconsistency in standard turbulence models that prevents them from accurately capturing the physics of rarefaction-turbulence interaction. In addition to providing physical insight, the solutions derived here can be used to verify algorithms of both the Reynolds-averaged and direct numerical simulation variety. Finally, dimensional analysis of the equations indicates that rapid distortion of turbulence can give rise to two distinct regimes in the turbulent spectrum: a distortion range at large scales where linear distortion effects dominate, and an inertial range at small scales where nonlinear effects dominate.

Three-temperature plasma shock solutions with gray radiation diffusion

The effects of radiation on the structure of shocks in a fully ionized plasma are investigated by solving the steady-state fluid equations for ions, electrons, and radiation. The electrons and ions are assumed to have the same bulk velocity but separate temperatures, and the radiation is modeled with the gray diffusion approximation. Both electron and ion conduction are included, as well as ion viscosity. When the material is optically thin, three-temperature behavior occurs. When the diffusive flux of radiation is important but radiation pressure is not, two-temperature behavior occurs, with the electrons strongly coupled to the radiation. Since the radiation heats the electrons on length scales that are much longer than the electron–ion Coulomb coupling length scale, these solutions resemble radiative shock solutions rather than plasma shock solutions that neglect radiation. When radiation pressure is important, all three components are strongly coupled. Results with constant values for the transport and coupling coefficients are compared to a full numerical simulation with a good match between the two, demonstrating that steady shock solutions constitute a straightforward and comprehensive verification test methodology for multi-physics numerical algorithms.

SIMPLE WAVES IN IDEAL RADIATION HYDRODYNAMICS

In the dynamic diffusion limit of radiation hydrodynamics, advection dominates diffusion; the latter primarily affects small scales and has negligible impact on the large-scale flow. The radiation can thus be accurately regarded as an ideal fluid, i.e., radiative diffusion can be neglected along with other forms of dissipation. This viewpoint is applied here to an analysis of simple waves in an ideal radiating fluid. It is shown that much of the hydrodynamic analysis carries over by simply replacing the material sound speed, pressure, and adiabatic index with the values appropriate for a radiating fluid. A complete analysis is performed for a centered rarefaction wave, and expressions are provided for the Riemann invariants and characteristic curves of the one-dimensional system of equations. The analytical solution is checked for consistency against a finite difference numerical integration, and the validity of neglecting the diffusion operator is demonstrated. An interesting physical result is that for a material component with a large number of internal degrees of freedom and an internal energy greater than that of the radiation, the sound speed increases as the fluid is rarefied. These solutions are an excellent test for radiation hydrodynamic codes operating in the dynamic diffusion regime. The general approach may be useful in the development of Godunov numerical schemes for radiation hydrodynamics.