Reconditioning your quantile function

Monte Carlo simulation is an important tool for modeling highly nonlinear systems (like particle colliders and cellular membranes), and random, floating-point numbers are their fuel. These random samples are frequently generated via the inversion method, which harnesses the mapping of the quantile function Q(u) (e.g. to generate proposal variates for rejection sampling). Yet the increasingly large sample size of these simulations makes them vulnerable to a flaw in the inversion method; Q(u) is ill-conditioned in a distribution's tails, stripping precision from its sample. This flaw stems from limitations in machine arithmetic which are often overlooked during implementation (e.g. in popular C++ and Python libraries). This paper introduces a robust inversion method, which reconditions Q(u) by carefully drawing and using uniform variates. pqRand, a free C++ and Python package, implements this novel method for a number of popular distributions (exponential, normal, gamma, and more).