Dusty gas with smoothed particle hydrodynamics - I. Algorithm and test suite

We present a new algorithm for simulating two-fluid gas and dust mixtures in smoothed particle hydrodynamics (SPH), systematically addressing a number of key issues including the generalized SPH density estimate in multifluid systems, the consistent treatment of variable smoothing-length terms, finite particle size, timestep stability, thermal coupling terms and the choice of kernel and smoothing length used in the drag operator. We find that using double-hump-shaped kernels improves the accuracy of the drag interpolation by a factor of several hundred compared to the use of standard SPH bell-shaped kernels, at no additional computational expense. In order to benchmark our algorithm, we have developed a comprehensive suite of standardized, simple test problems for gas and dust mixtures: DUSTYBOX, DUSTYWAVE, DUSTYSHOCK, DUSTYSEDOV and DUSTYDISC, the first three of which have known analytic solutions. We present the validation of our algorithm against all of these tests. In doing so, we show that the spatial resolution criterion Delta less than or similar to c(s)t(s) is a necessary condition in all gas+ dust codes that becomes critical at high drag (i.e. small stopping time t(s)) in order to correctly predict the dynamics. Implicit timestepping and the implementation of realistic astrophysical drag regimes are addressed in a companion paper.