In 1985 Parisi and Rapuano reported that the susceptibility of the three-dimensional Ising model on a 24^3 lattice, measured with their generator, disagreed by four standard deviations with the result obtained on a special-purpose computer, and attributed the discrepancy to the random numbers; they proposed an exclusive-or (EOR) output mixing, un=(xn⊕xn−61)/232u_n=(x_n \oplus x_{n-61})/2^{32} on the lagged-Fibonacci sequence xn=xn−24+xn−55x_n=x_{n-24}+x_{n-55}, as a remedy against the order-three correlations of feedback-register generators. The later literature resolved much of the original question: two-tap feedback-register generators were shown to carry explicit three- and four-point correlations, to bias observables under cluster and single-spin update algorithms, and the exact Parisi-Rapuano generator fails eight Crush tests of the TestU01 battery while remaining in production use on Janus-class FPGA machines. What remains open is the gap documented by the subsequent literature and stated verbatim by the most recent authoritative review (James and Moneta, 2020): 'There is no general method to determine how good a RNG must be for a particular MC application' — namely, a general, predictive, quantitative diagnostic that anticipates, for a given pseudo-random number generator and a given Monte Carlo implementation, whether and by how much the extracted observables will be systematically biased, without having to run the production simulation. This problem asks for such a criterion or certifying protocol, demonstrated against the pinned historical benchmark incidents, and for its application to the standing case of the Parisi-Rapuano EOR-mixed generator across update schemes and lattice sizes.

Background

Monte Carlo (MC) simulation is the standard tool of statistical physics: expectation values such as the susceptibility of the Ising model (the standard model of spins on a lattice with nearest-neighbor ferromagnetic couplings, where susceptibility measures the response of the magnetization to a small field) are estimated by stochastic sampling. Near a continuous phase transition, finite-size scaling of such estimates yields universal quantities (critical exponents and amplitude ratios), and high-precision determinations of these are a core activity of computational statistical physics and lattice field theory. All stochastic sampling consumes pseudo-random numbers produced by deterministic algorithms: linear congruential generators (LCG, of the form xn+1=axn+c mod mx_{n+1}=ax_n+c \bmod m), feedback-register or lagged-Fibonacci generators in which a new word is combined from earlier ones (for example the very fast two-tap exclusive-or shift-register sequence generator of Kirkpatrick and Stoll, 1981), and composite schemes. These sequences are exactly reproducible but not truly independent, and their defects matter precisely in the delicate regime of long, correlated simulations.

A two-term recurrence such as xn=xn−24+xn−55x_n=x_{n-24}+x_{n-55} is an exact algebraic relation among three outputs, so the raw words of such a generator have non-uniform joint statistics at the third order (order-three or third-order correlations: the average of a product of three outputs at fixed lags differs from the value for ideal independent uniform random numbers). Parisi and Rapuano (Phys. Lett. B 157 (1985) 301) reported the empirical consequence: a four-standard-deviation discrepancy in the 24^3 three-dimensional Ising susceptibility between simulations that differed, in their assessment, only in the random number generator, the comparison target being a result from a special-purpose computer. They proposed the EOR (exclusive-or, XOR) output mixing Re(i)=EOR⁡(R(i),R(i−61))R_e(i)=\operatorname{EOR}(R(i),R(i-61)) as a cheap repair for feedback-register generators; in modern notation the generator that bears their name is xn=xn−24+xn−55 mod 232x_n=x_{n-24}+x_{n-55} \bmod 2^{32}, un=(xn⊕xn−61)/232u_n=(x_n \oplus x_{n-61})/2^{32}. A companion study the same year (Filk, Marcu and Fredenhagen, Phys. Lett. B 165 (1985) 125) analyzed the long-range correlations of congruential generators and their influence on lattice simulations near criticality.

Subsequent work made the phenomenon quantitative and famous. Ferrenberg, Landau and Wong (Phys. Rev. Lett. 69 (1992) 3382) showed that the Wolff single-cluster algorithm (an update scheme that flips connected clusters of same-sign spins and consumes a variable number of random numbers per step) produced systematically wrong two-dimensional Ising results with a shift-register generator then considered among the best, while older, known-defective generators performed acceptably. Vattulainen, Ala-Nissila and Kankaala (Phys. Rev. Lett. 73 (1994) 2513) introduced physical tests of randomness (Ising autocorrelation times, random walks, block statistics), traced these errors to short-range correlations of generalized feedback shift-register (GFSR) sequences, and measured the correlation length; their longer analysis (Phys. Rev. E 52 (1995) 3205) systematized the approach. Coddington (Int. J. Mod. Phys. C 5 (1994) 547) ran an independence-testing program on a benchmark: he compared many widely used generators by high-precision two-dimensional Ising simulations, for which exact finite-size results are known, under the Metropolis, Swendsen-Wang, and Wolff update algorithms, and found that many generators passing standard statistical tests fail these Monte Carlo tests. Ziff (Comput. Phys. 12 (1998) 385) computed the explicit three- and four-point correlations of two-tap GFSR rules and showed that virtually all two-tap generators fail a sensitive random-walk test, while four-tap rules with large tap offsets pass, which explains the order-three correlation problem of the whole family and supplies a structural fix. Shchur and Blöte (Phys. Rev. E 55 (1997) R4905) measured scaling relations for how the systematic deviations produced by binary shift-register generators in Wolff-cluster simulations of the two-dimensional Ising model depend on lattice size, shift-register length, and the number of bits correlated by the production rule — a family-specific, a-posteriori quantitative link. On the diagnostic side, Ballesteros and Martín-Mayor (Phys. Rev. E 58 (1998) 6787) showed that Schwinger-Dyson equations (exact integration-by-parts identities satisfied by Ising expectation values) provide a near-costless, sensitive test of random-number-induced bias that can be appended to any two- or three-dimensional Ising simulation. On the generator side, Lüscher's RANLUX (Comput. Phys. Commun. 79 (1994) 100, Fortran implementation by James, ibid. 111) derived quality from the mixing properties of an underlying chaotic dynamical system, and the TestU01 battery of L'Ecuyer and Simard (ACM Trans. Math. Softw. 33 (2007) 22) became the standard statistical yardstick (SmallCrush, Crush, BigCrush test suites).

Two facts define the current frontier. First, Manssen, Weigel and Hartmann (Eur. Phys. J. Spec. Top. 210 (2012) 53) tested the exact Parisi-Rapuano generator under TestU01 and found two suspicious tests in SmallCrush and eight failed tests in Crush, so the EOR mixing demonstrably does not remove all correlations; nevertheless the same generator is the production random number generator of the Janus and Janus II FPGA (field-programmable gate array)-based special-purpose machines used for large spin-glass simulations, spin glasses being disordered magnetic systems with frustrated interactions (Baity-Jesi et al., Comput. Phys. Commun. 185 (2014) 550), where roughly one thousand update engines per chip draw from it. Second, the field's answers to 'how do I know my generator is safe for my simulation' remain either a posteriori (run physical tests or cross-check with an independent generator, as still done in applied toolkits, e.g. the 2020 study of generator selection in the GEANT4-based GATE Monte Carlo toolkit, where GEANT4 is a widely used particle-transport simulation toolkit and GATE a Monte Carlo application built on it) or generator-level and a priori (choose generators with theoretical mixing guarantees); the modern review by James and Moneta (Comput. Softw. Big Sci. 4 (2020) 2) states that empirical testing of random number generators was shown by the 1992 incident not to be good enough, recommends theoretically grounded generators, and concludes flatly that 'There is no general method to determine how good a RNG must be for a particular MC application'; it gives no criterion that couples a generator's measurable correlation defects to the bias of a specific estimator in a specific Monte Carlo implementation. That coupling is what is still missing: a predictive, quantitative diagnostic rather than an after-the-fact detection.

Problem Statement

Develop a general, predictive, quantitative diagnostic for random-number-generator-induced systematic error in Monte Carlo simulations, and settle the residual case of the Parisi-Rapuano generator with it. Concretely, answer the following, in a form usable before committing production computing resources:

  1. Formulate a criterion, formula, or certifying protocol that takes as inputs (i) measurable or computable correlation properties of a pseudo-random number generator (for example its lag-space three- and higher-order correlation functions, as computable for two-tap feedback-register rules, or its TestU01-family failure profile) and (ii) the random-number consumption pattern of a Monte Carlo implementation (update algorithm such as Metropolis single-spin, Swendsen-Wang or Wolff cluster moves; lattice size; how many random numbers are consumed per update and in what structure), and predicts whether, and at what magnitude relative to statistical error, the implementation's target observables (for example energies, susceptibilities, autocorrelation times, critical exponents or amplitude ratios obtained by finite-size scaling) will be biased.

  2. Demonstrate the criterion's predictive power against the historical benchmark record: the two-dimensional Ising model with exact finite-size results, reproducing the known incidents — the Ferrenberg-Landau-Wong failure (the Wolff cluster algorithm on the two-dimensional Ising model producing systematically wrong estimates with a two-tap shift-register generator then considered among the best available), and Coddington's multi-generator, multi-algorithm results — while correctly classifying generators that were found safe in the same studies. Its predictive content must additionally be exhibited out of sample: for at least one generator-algorithm-lattice combination not used in formulating or calibrating the criterion, the criterion's prediction of the sign and magnitude of the bias must be checked against the historical record or by new simulation.

  3. Apply it to the standing open case: whether the residual, battery-detected correlations of the exact Parisi-Rapuano generator xn=xn−24+xn−55 mod 232x_n=x_{n-24}+x_{n-55} \bmod 2^{32}, un=(xn⊕xn−61)/232u_n=(x_n \oplus x_{n-61})/2^{32} (known to fail eight Crush tests and in production use on Janus-class FPGA spin-glass machines) bias standard observables across update schemes (Metropolis, Swendsen-Wang, Wolff) and at lattice sizes beyond those already tested in the cited literature and relevant to current production use, and quantify any such bias against statistical errors.

The common objective uniting the three items is to replace after-the-fact detection of random-number-induced bias (physical tests, Schwinger-Dyson checks, statistical batteries, cross-run comparisons) with prediction or certification for a stated implementation. Items 2 and 3 are independent demonstrations of this single objective rather than illustrations of item 1.

References

  1. G. Parisi and F. Rapuano, 'Effects of the random number generator on computer simulations', Physics Letters B 157 (1985) 301-302. DOI: 10.1016/0370-2693(85)90670-7. https://doi.org/10.1016/0370-2693(85)90670-7
  2. T. Filk, M. Marcu and K. Fredenhagen, 'Long range correlations in random number generators and their influence on Monte Carlo simulations', Physics Letters B 165 (1985) 125-130. DOI: 10.1016/0370-2693(85)90705-1. https://doi.org/10.1016/0370-2693(85)90705-1
  3. S. Kirkpatrick and E. P. Stoll, 'A very fast shift-register sequence random number generator', Journal of Computational Physics 40 (1981) 517-526. DOI: 10.1016/0021-9991(81)90227-8. https://doi.org/10.1016/0021-9991(81)90227-8
  4. A. M. Ferrenberg, D. P. Landau and Y. J. Wong, 'Monte Carlo simulations: Hidden errors from good random number generators', Physical Review Letters 69 (1992) 3382-3384. DOI: 10.1103/PhysRevLett.69.3382. https://doi.org/10.1103/PhysRevLett.69.3382
  5. I. Vattulainen, T. Ala-Nissila and K. Kankaala, 'Physical tests for random numbers in simulations', Physical Review Letters 73 (1994) 2513-2516. DOI: 10.1103/PhysRevLett.73.2513. https://doi.org/10.1103/PhysRevLett.73.2513
  6. P. D. Coddington, 'Analysis of random number generators using Monte Carlo simulation', International Journal of Modern Physics C 5 (1994) 547-560. DOI: 10.1142/S0129183194000726. https://doi.org/10.1142/S0129183194000726
  7. I. Vattulainen, T. Ala-Nissila and K. Kankaala, 'Physical models as tests of randomness', Physical Review E 52 (1995) 3205-3214. DOI: 10.1103/PhysRevE.52.3205. https://doi.org/10.1103/PhysRevE.52.3205
  8. H. G. Ballesteros and V. Martín-Mayor, 'Test for random number generators: Schwinger-Dyson equations for the Ising model', Physical Review E 58 (1998) 6787-6791. DOI: 10.1103/PhysRevE.58.6787. https://doi.org/10.1103/PhysRevE.58.6787
  9. R. M. Ziff, 'Four-tap shift-register-sequence random-number generators', Computers in Physics 12 (1998) 385-392. DOI: 10.1063/1.168692. https://doi.org/10.1063/1.168692
  10. M. Lüscher, 'A portable high-quality random number generator for lattice field theory simulations', Computer Physics Communications 79 (1994) 100-110. DOI: 10.1016/0010-4655(94)90232-1. https://doi.org/10.1016/0010-4655(94)90232-1
  11. F. James, 'RANLUX: A Fortran implementation of the high-quality pseudorandom number generator of Lüscher', Computer Physics Communications 79 (1994) 111-114. DOI: 10.1016/0010-4655(94)90233-X. https://doi.org/10.1016/0010-4655(94)90233-X
  12. F. James, 'A review of pseudorandom number generators', Computer Physics Communications 60 (1990) 329-344. DOI: 10.1016/0010-4655(90)90032-V. https://doi.org/10.1016/0010-4655(90)90032-V
  13. P. L'Ecuyer and R. Simard, 'TestU01: A C library for empirical testing of random number generators', ACM Transactions on Mathematical Software 33 (2007) article 22, 1-40. DOI: 10.1145/1268776.1268777. https://doi.org/10.1145/1268776.1268777
  14. M. Manssen, M. Weigel and A. K. Hartmann, 'Random number generators for massively parallel simulations on GPU', The European Physical Journal Special Topics 210 (2012) 53-71. DOI: 10.1140/epjst/e2012-01637-8. arXiv:1204.6193. https://doi.org/10.1140/epjst/e2012-01637-8
  15. M. Baity-Jesi, R. A. Baños, A. Cruz, L. A. Fernandez, J. M. Gil-Narvion, A. Gordillo-Guerrero, D. Iñiguez, A. Maiorano, F. Mantovani, E. Marinari, V. Martín-Mayor, J. Monforte-Garcia, A. Muñoz Sudupe, D. Navarro, G. Parisi, S. Perez-Gaviro, M. Pivanti, F. Ricci-Tersenghi, J. J. Ruiz-Lorenzo, S. F. Schifano, B. Seoane, A. Tarancon, R. Tripiccione and D. Yllanes, 'Janus II: a new generation application-driven computer for spin-system simulations', Computer Physics Communications 185 (2014) 550-559. DOI: 10.1016/j.cpc.2013.10.019. arXiv:1310.1032. https://doi.org/10.1016/j.cpc.2013.10.019
  16. P. D. Coddington and S.-H. Ko, 'Techniques for empirical testing of parallel random number generators', Proceedings of the 12th International Conference on Supercomputing (ICS'98), ACM (1998) 282-288. DOI: 10.1145/277830.277895. https://doi.org/10.1145/277830.277895
  17. K. Hongo, R. Maezono and K. Miura, 'Random number generators tested on quantum Monte Carlo simulations', Journal of Computational Chemistry 31 (2010) 2186-2194. DOI: 10.1002/jcc.21509. https://doi.org/10.1002/jcc.21509
  18. Y. Lin, F. Wang and B. Liu, 'Random number generators for large-scale parallel Monte Carlo simulations on FPGA', Journal of Computational Physics 360 (2018) 93-103. DOI: 10.1016/j.jcp.2018.01.029. https://doi.org/10.1016/j.jcp.2018.01.029
  19. F. Sepehri, M. Hajivaliei and H. Rajabi, 'Selection of random number generators in GATE Monte Carlo toolkit', Nuclear Instruments and Methods in Physics Research A 973 (2020) 164172. DOI: 10.1016/j.nima.2020.164172. https://doi.org/10.1016/j.nima.2020.164172
  20. F. James and L. Moneta, 'Review of High-Quality Random Number Generators', Computing and Software for Big Science 4 (2020) 2. DOI: 10.1007/s41781-019-0034-3. https://doi.org/10.1007/s41781-019-0034-3
  21. L. N. Shchur and H. W. J. Blöte, 'Cluster Monte Carlo: Scaling of Systematic Errors in the 2D Ising Model', Physical Review E 55 (1997) R4905. DOI: 10.1103/PhysRevE.55.R4905. arXiv:cond-mat/9703050. https://doi.org/10.1103/PhysRevE.55.R4905