Accelerated Monte Carlo Methods for Coulomb Collisions

As an alternative to binary-collision models for simulating Coulomb collisions in the Fokker-Planck limit, we present a new numerical higher-order-accurate time integration scheme for Langevin-equation-based collisions. A Taylor-series expansion of the stochastic differential equations is used to improve upon the standard Euler time integration. Additional Milstein terms arise in the time-discretization due to both the velocity dependence of the diffusion coefficients, and the aggregation of angular deflections. We introduce an accurate, easily computable direct sampling method for the multidimensional terms -- an approximation to the double integral over products of Gaussian random processes. Including these terms improves the strong convergence of the time integration of the particle trajectories from O($\Delta $t1/2) to O($\Delta $t). This is useful as a both a first step towards direct higher-order weak schemes (for computing average quantities), and as a key component in a ``multi-level'' scheme that returns a computationally efficient estimate of averaged quantities. The latter is maximally efficient, in the asymptotic sense, when used with Milstein terms, and therefore the optimal choice of multi-level scheme. We present results showing both the improved strong convergence of the new integration method, and the increased efficiency multi-level scheme.