Again, random numbers fall mainly in the planes: xorshift128+ generators
AAgain, random numbers fall mainly in the planes:xorshift128+ generators
HIROSHI HARAMOTO,
Faculty of Education, Ehime University
MAKOTO MATSUMOTO,
Graduate School of Science, Hiroshima UniversityXorshift128+ are pseudo random number generators with eight sets of parameters. Some of them are standardgenerators in many platforms, such as JavaScript V8 Engine. We show that in the 3D plots generated by thismethod, points concentrate on planes, ruining the randomness.CCS Concepts: •
Mathematics of computing → Mathematical software ; •
Theory of computation → Pseudorandomness and derandomization .Additional Key Words and Phrases: pseudo random number generators, xorshift128+, discrete mathematics
ACM Reference Format:
Hiroshi Haramoto and Makoto Matsumoto. 2019. Again, random numbers fall mainly in the planes: xorshift128+generators.
ACM Trans. Math. Softw.
0, 0, Article 0 (January 2019), 7 pages. https://doi.org/0
Let x be a 64-bit (unsigned) integer. Let F = { , } denote the two-element field, and x ∈ F isconsidered to be a 64-dimensional row vector. A linear transform x (cid:55)→ x L is defined as the leftshift by 1 bit, and x (cid:55)→ x R is the right shift by 1 bit. The matrix (cid:169)(cid:173)(cid:173)(cid:173)(cid:173)(cid:173)(cid:173)(cid:173)(cid:171)
01 01 . . .. . . . . . (cid:170)(cid:174)(cid:174)(cid:174)(cid:174)(cid:174)(cid:174)(cid:174)(cid:172) is L and (cid:169)(cid:173)(cid:173)(cid:173)(cid:173)(cid:173)(cid:173)(cid:173)(cid:171) . . . . . .. . . (cid:170)(cid:174)(cid:174)(cid:174)(cid:174)(cid:174)(cid:174)(cid:174)(cid:172) is R . Let I denote the identity matrix.The xorshift128+ [3] pseudo random number generator (PRNG) has 128-bit state space. A stateconsists of two 64-bit words ( s i , s i + ) , and the next state is ( s i + , s i + ) , where s i + = s i ( I ⊕ L a )( I ⊕ R b ) ⊕ s i + ( I ⊕ R c ) . (1) Authors’ addresses: Hiroshi Haramoto, [email protected], Faculty of Education, Ehime University, 3 Bunkyocho, Mat-suyama, Ehime, 790-8577; Makoto Matsumoto, [email protected], Graduate School of Science, HiroshimaUniversity, 1-3-1 Kagamiyama, Higashi-Hiroshima, Hiroshima, 739-8526.Permission to make digital or hard copies of all or part of this work for personal or classroom use is granted without feeprovided that copies are not made or distributed for profit or commercial advantage and that copies bear this notice andthe full citation on the first page. Copyrights for components of this work owned by others than ACM must be honored.Abstracting with credit is permitted. To copy otherwise, or republish, to post on servers or to redistribute to lists, requiresprior specific permission and/or a fee. Request permissions from [email protected].© 2019 Association for Computing Machinery.0098-3500/2019/1-ART0 $15.00https://doi.org/0 ACM Trans. Math. Softw., Vol. 0, No. 0, Article 0. Publication date: January 2019. a r X i v : . [ c s . CR ] A ug :2 Haramoto and Matsumoto Here the notation ⊕ is used for bitwise xor operation, or equivalently addition of vectors in F and that of 64 ×
64 matrices I and L a , etc. The output o i at the i -th state ( s i , s i + ) is given by o i = s i + s i + mod 2 , where + denotes addition of 64-bit integers. The 128-bit state ( s , s ) is the initial state.The consecutive three outputs from ( s i , s i + ) is x = s i + s i + mod 2 , y = s i + + s i + mod 2 and z = s i + + s i + mod 2 , and we shall show some relations among x , y and z . Let x , y be n -bit unsigned integers. We consider mainly n = x and y . Our claim is that x ⊕ y is with non negligible probability well-approximated byone of x + y , x − y or y − x , as analyzed below.We consider x = ( x , . . . , x n ) ∈ F n as an n -dimensional vector, which is also considered as an n -bitunsigned integer denoted by x = n (cid:205) i = x i n − i . In this situation, we write x = x . Let y = ( y , . . . , y n ) beanother n -dimensional vector. Then x ⊕ y is the addition of F vectors. We discuss when x ⊕ y = x + y holds, where the both sides are regarded as n -bit integers. Note that the operation x ⊕ y is similarto x + y , except that no over flow is reflected. Then x ⊕ y ≤ x + y holds and the equality holds ifand only if no overflow occurs. Equivalently, if and only if ( x i , y i ) ∈ {( , ) , ( , ) , ( , )} holds for i =
1, 2, . . . , n . Thus, among 4 n possibilities of pairs x , y , exactly 3 n pairs satisfy x + y = x ⊕ y . Ifplus is taken module 2 n , there are more cases with equality, e.g., x = ( , , . . . , ) = y is the case.This observation is summarized as follows.Theorem 2.1 (xor eqals sum). Let x , y ∈ F n be n -bit integers x , y , where x = ( x , . . . , x n ) , y = ( y , . . . , y n ) , x = n (cid:205) i = x i n − i and y = n (cid:205) i = y i n − i . Then x ⊕ y ≤ x + y holds, and the equality holds ifand only if ( x i , y i ) (cid:44) ( , ) holds for i = , , . . . , n . Among n pairs ( x , y ) , n pairs satisfy the equality.More pairs satisfy x ⊕ y = x + y mod 2 n . Proof. Proof follows from the previous observation. □ Another observation is about when the equality x ⊕ y = x − y occurs. Again, we do not take modulo 2 n for the right hand side. If we compute subtraction x − y in binary without borrows, we obtain x ⊕ y . There may be borrows, so we have inequality x ⊕ y ≥ x − y , with equality holds when no borrow occurs for each digit, or equivalently, the pair of bits ( x i , y i ) lies in {( , ) , ( , ) , ( , )} for each i = , , . . . , n . There are 3 n such pair ( x , y ) .Theorem 2.2 (xor eqals subtraction). Let x , y ∈ F n be as in Theorem 2.1. We have inequality x ⊕ y ≥ x − y , and the equality holds if and only if ( x i , y i ) (cid:44) ( , ) for each i = , , . . . , n . There are n such pairs. Theorem 2.3.
Let X be the set of pairs {( x , y ) | x , y ∈ F n } , and put A : = {( x , y ) ∈ X | x ⊕ y = x + y } B : = {( x , y ) ∈ X | x ⊕ y = x − y } C : = {( x , y ) ∈ X | x ⊕ y = y − x } . ACM Trans. Math. Softw., Vol. 0, No. 0, Article 0. Publication date: January 2019. gain, random numbers fall mainly in the planes: xorshift128+ generators 0:3
Then, X = n , A = B = C = n , ( A ∩ B ) = ( B ∩ C ) = ( C ∩ A ) = n , and ( A ∩ B ∩ C ) = hold. In particular, ( A ∪ B ∪ C ) = · n − · n + holds. Proof. We have A = {( x , y ) ∈ X | ∀ i ( x i , y i ) ∈ {( , ) , ( , ) , ( , )}} B = {( x , y ) ∈ X | ∀ i ( x i , y i ) ∈ {( , ) , ( , ) , ( , )}} C = {( x , y ) ∈ X | ∀ i ( x i , y i ) ∈ {( , ) , ( , ) , ( , )}} , and the second series of equalities hold. Since A ∩ B = {( x , y ) ∈ X | ∀ i ( x i , y i ) ∈ {( , ) , ( , )}} , ( A ∩ B ) = n follows, and ( B ∩ C ) = ( C ∩ A ) = n is similarly proved. We have A ∩ B ∩ C = {( x , y ) ∈ X | ∀ i ( x i , y i ) ∈ {( , )}} , and see that ( A ∩ B ∩ C ) =
1. The last equality follows fromthe standard inclusion-exclusion principle. □ Example 2.4.
Suppose n =
3, i.e., we consider three-bit precision. Then, among 4 =
64 pairs ( x , y ) , we showed that at least one of x ⊕ y = x + y , x − y , y − x occurs for 3 · − · + = − = n =
4, among 256 pairs, 196pairs satisfy one of the three relations. + IN XORSHIFT128+
We consider the consecutive three outputs x = s i + s i + mod 2 y = s i + + s i + mod 2 z = s i + + s i + mod 2 , where s i + and s i + are determined from ( s i , s i + ) by the recursion (1). We analyze z = s i + + s i + mod 2 = ( s i + ( I ⊕ R c ) ⊕ s i ( I ⊕ L a )( I ⊕ R b )) (2) + ( s i + ( I ⊕ R c ) ⊕ s i + ( I ⊕ L a )( I ⊕ R b )) mod 2 (3)We could not give an exact analysis, but give an intuitional approximation. Numbers b and c arelarger than 10, and the most significant b bits of x ( I + R b ) is identical with those of x . Thus, as faras we concentrate on the most significant several bits (we use Theorems mainly for n = I ⊕ R b and I ⊕ R c to be the identity matrix I . Thus, we have an approximation z ≈ ( s i + ⊕ s i ( I ⊕ L a )) (4) + ( s i + ⊕ s i + ( I ⊕ L a )) . (5)We denote by ≈ when the both sides coincide up to the most significant min { b , c } bits or somespecified n -bits, except that with small probability the matrices R b and/or R c , through the carry of + between (4) and (5), may affect on the MSBs.Now in (4), we have s i ( I ⊕ L a ) = s i ⊕ s i L a = s i ⊕ ( a s i mod 2 ) (6) ACM Trans. Math. Softw., Vol. 0, No. 0, Article 0. Publication date: January 2019. :4 Haramoto and Matsumoto which is in Theorem 2.3, according to the cases A , B , C , s i ⊕ ( a s i mod 2 ) ≈ ( + a ) s i mod 2 Case A i ( − a ) s i mod 2 Case B i ( a − ) s i mod 2 Case C i unknown otherwisewith respect to the most significant n bits. From now on we consider the most significant n bits,with mainly n =
3. Again in (4), s i + ⊕ s i ( I ⊕ L a ) ≈ s i + + s i ( I ⊕ L a ) Case A ′ i s i + − s i ( I ⊕ L a ) Case B ′ i − s i + + s i ( I ⊕ L a ) Case C ′ i unknown otherwiseThe same kind of case divisions are straight froward for i replaced with i +
1, denoted by A i + etc.We consider the following cases: Case + : Case A ′ i and Case A ′ i + occur, Case − : Case B ′ i and Case B ′ i + occur, Case t − : Case C ′ i and Case C ′ i + occur,and orthogonally the cases: Case + a : Case A i and Case A i + occur, Case − a : Case B i and Case B i + occur, Case a − C i and Case C i + occur.We assume that both one of the cases + , − , t − and one of the cases 1 + a , 1 − a , 2 a − n = (cid:32)(cid:18) (cid:19) × (cid:33) ≈ . · · · . For example, assume that the case − and 1 − a occur. Then, we have z ≈ ( s i + − s i ( − a )) + ( s i + − s i + ( − a )) = ( s i + + s i + ) + ( a − )( s i + s i + ) = y + ( a − ) x , where we omit modulo 2 . Table 1. z from x and y by case division Case\Case 1 + a − a a − + ( + a ) x + y ( − a ) x + y ( a − ) x + y − −( + a ) x + y ( a − ) x + y ( − a ) x + y t − ( + a ) x − y ( − a ) x − y ( a − ) x − y A straight forward computation by case division (Table 1) gives that with non negligible proba-bility one of z ≈ ±( + a ) x ± y , z ≈ ±( a − ) x ± y . ACM Trans. Math. Softw., Vol. 0, No. 0, Article 0. Publication date: January 2019. gain, random numbers fall mainly in the planes: xorshift128+ generators 0:5 y x 0z 0 0.2 0.4 0.6 0.8 1 0.20.40.60.810.2 0.4 0.60.801
Fig. 1. 3D plots by xorshift128+: x -axis magnified by a factor of hold. This shows that the consecutive three outputs ( x , y , z ) by xorshift128+ tend to lie on eightplanes, which give an explanation on Figure 1 [1]. We compare these planes with the outputs ofxorshift128+.Figure 2 describes four planes z = ±( + ) x ± y mod 1with restriction 0 ≤ x ≤ / , 0 ≤ y ≤
1. The x -axis is magnified with the factor 2 . The otherfour planes with coefficient 2 − +
1, and so omitted. Each planehas two connected components in this region. Figure 3 shows the union of these four planes.Figure 4 shows the outputs of xorshift128+ with parameter ( a , b , c ) = ( , , ) . Let ( x , y , z ) bethe consecutive outputs in [ , ) . We only pick up those with x ≤ / , and plot ( x , y , z ) . Werepeat this until we obtain 10000 points.Figure 5 contains both the four planes (Figure 3) and the outputs of xorshift128+ (Figure 4). Thiscoincidence justifies the approximated analysis done in this section. G. Marsaglia said “random numbers fall mainly in the planes” [2]. Experiments show that Vigna’sxorshift128+ have this property. An analysis of this phenomenon based on approximation of xor byarithmetic addition and subtraction is discussed.
ACKNOWLEDGMENTS
This work was supported by JSPS KAKENHI Grant Numbers 26310211, JP17K14234 and JP18K03213.
REFERENCES [1] Hiroshi Haramoto, Makoto Matsumoto, and Mutsuo Saito. 2019. Pseudo random number generators: attention for anewly proposed generator. arXiv:arXiv:1907.03251[2] George Marsaglia. 1968. RANDOM NUMBERS FALL MAINLY IN THE PLANES.
Proceedings of the National Academy ofSciences
61, 1 (1968), 25–28. https://doi.org/10.1073/pnas.61.1.25[3] Sebastiano Vigna. 2017. Further scramblings of MarsagliaâĂŹs xorshift generators.
J. Comput. Appl. Math.
315 (2017),175 – 181. https://doi.org/10.1016/j.cam.2016.11.006ACM Trans. Math. Softw., Vol. 0, No. 0, Article 0. Publication date: January 2019. :6 Haramoto and Matsumoto (a) (b)(c) (d)Fig. 2. Pictures of four planes: (a): z = ( + ) x + y mod 1 , (b): z = ( + ) x − y mod 1 , (c): z = −( + ) x + y mod 1 , (d): z = −( + ) x − y mod 1 Fig. 3. The union of four planes
ACM Trans. Math. Softw., Vol. 0, No. 0, Article 0. Publication date: January 2019. gain, random numbers fall mainly in the planes: xorshift128+ generators 0:7