Accelerating equilibrium isotope effect calculations. II. Stochastic implementation of direct estimators

Path integral calculations of equilibrium isotope effects and isotopic fractionation are expensive due to the presence of path integral discretization errors, statistical errors, and thermodynamic integration errors. Whereas the discretization errors can be reduced by high-order factorization of the path integral and statistical errors by using centroid virial estimators, two recent papers proposed alternative ways to completely remove the thermodynamic integration errors: Cheng and Ceriotti [J. Chem. Phys. 141, 244112 (2015)] employed a variant of free-energy perturbation called "direct estimators," while Karandashev and Vaníček [J. Chem. Phys. 143, 194104 (2017)] combined the thermodynamic integration with a stochastic change of mass and piecewise-linear umbrella biasing potential. Here, we combine the former approach with the stochastic change in mass in order to decrease its statistical errors when applied to larger isotope effects and perform a thorough comparison of different methods by computing isotope effects first on a harmonic model and then on methane and methanium, where we evaluate all isotope effects of the form CH_4-xD_x/CH₄ and CH_5-xD_x ⁺/CH₅ ⁺, respectively. We discuss the reasons for a surprising behavior of the original method of direct estimators, which performed well for a much larger range of isotope effects than what had been expected previously, as well as some implications of our work for the more general problem of free energy difference calculations.