Patent Application
20140337399
The present invention relates to methods of obtaining uniform and independent random numbers on computers. More specifically, the methods relate to the multiplicative congruential generation that includes:
1. a positive integer d called modulus,
2. a positive integer z called multiplier that is coprime with d, and
3. a positive integer n called initial value or seed coprime with d, emits a sequence of integers {r1, r2, r3, . . . } from recursive congruence relations
r
1
≡n,r
j+1
≡zr
j mod(d),0<rj<d,j=1,2,3, . . . , Eq. 1
and gives an output sequence {v1, v2, v3, . . . } in the interval (0, 1) by realizing the division
v
j
:=r
j
/d,j=1,2,3, . . . Eq. 2
A multiplicative congruential generator for uniform and independent random numbers with the modulus d, the multiplier z and the initial value n will also be noted as (d, z, n) symbolically.
A Random Number Generator (RNG) is a computational or physical device designed to generate a sequence of numbers or symbols that lack any pattern, i.e. appear random.
The many applications of randomness have led to the development of several different methods for generating random data. Many of these have existed since ancient times, including dice, coin flipping, the shuffling of playing cards, the use of yarrow stalks, and many other techniques. Because of the mechanical nature of these techniques, generating large numbers of sufficiently random numbers (important in statistics) required a lot of work and/or time. Thus, results would sometimes be collected and distributed as random number tables. Nowadays, after the advent of computational RNGs, a growing number of government-run lotteries, and lottery games, are using RNGs instead of more traditional drawing methods. The use of RNGs arises now in gambling, statistical sampling, computer simulation, cryptography, completely randomized design, and other areas where the required statistical precision, namely how precisely the generated random numbers should realize the assumed statistical distribution, has a large variety. The present invention is related to the rapid increase of applications of large-scale scientific and mathematical simulations that call more than ever for higher precision of the statistical nature of generated random numbers. In one extreme end, numerical simulations of stochastic differential equations (SDEs) arise in various scientific, technological or even economical applications. Schemes for SDEs with higher precision realize faster computation, but it is only on the basis of the precise realization of random number statistics. The precision of the statistics of random numbers is only remotely tied to the nature of the unpredictability which is fundamental, typically in cryptography. The now popular Mersenne twister is at the one end of the highest position for this unpredictability with its gigantic length of the period and the nominal equidistribution property. Yet, the precision can be realized solely and indispensably by statistical tests valid on properties over the whole period of generated random number sequences, and the tremendous period length of Mersenne twister impedes any systematic statistical tests of generated random numbers. Note also that the equidistribution property is insufficient in ensuring the statistical precision: There are a vast number of ways to choose primitive polynomials that realize the same equidistribution property, and an arbitrary choice can hardly hit on statistically excellent outputs. Mersenne twister still lacks any theoretical or empirical basis to selects one of such primitive polynomials as its own. In this regard a fundamental technological point is that any RNGs on computers should deal exclusively with integer sequences. This is necessitated by the requirement to reproduce the identical sequence of random numbers with arbitrarily large length on demand of users, which should be realized without the hindrance of truncation and round-off errors. Another reason for this necessity is the requirement of transportability, that one and the same random number sequence should be generated on any computers or in any computer languages. Though seemingly obstructive, these necessities in fact provide a great mathematical convenience by the fact that any finite sequence of integers {x1, x2, . . . , xT}, that is to give a sequence of uniform and independent random numbers as {x1/z, x2/z, . . . , xT/z} with a large integer z, may be approximated by the so-called multiplicative congruential RNGs. And any multiplicative congruential RNGs are admitted readily to spectral tests, one of the most powerful of tests on the precision of the statistics of random number sequences. The first of above noted points, significant by its mathematical and conceptual contents, should here be accounted for fully. Any sequence of integers {x1, x2, . . . , xT}, each of which are in the range from 0 to z−1 with a large integer z, may be regarded as corresponding to the base z expansion of a rational number x=0. {dot over (x)}1 x2 . . . {dot over (x)}T with period T. Multiplication of x with zT gives that x has the quotient expression x=n/d with d:=zT−1, n:=x1zT−1+x2zT−2+ . . . +xT z0. Integers {x1, x2, . . . , xT} are thus obtained by the arithmetic process of division of n by d, and there holds the equation zrj=xj+1d+rj+1 with 0<rj<d for j=1, 2, . . . implying that the division of the j-th remainder rj multiplied by the base z with d is to give the next quotient xj+1 and the next remainder rj+1. Another expression of this equation is rj+1≡zrj mod(d). Division of the original equation by zd gives |rj/d−xj+1/z|=rj+1/zd<1/z, the key estimate. Imagine the usual circumstance for random numbers with d and z large than 232. The term rj/d is precisely the j-th output of a multiplicative congruential RNG with the modulus d and the multiplier z, and the inequality proves that any sequence {x1/z, x2/z, . . . , xT/z} of uniform and independent random numbers, as it should be realized on computers, is approximated by a multiplicative congruential one {r1/d, r2/d, . . . , rT/d}. Since the sequence {x1/z, x2/z, . . . , xT/z} of uniform and independent random numbers is arbitrary, the finding of this structure dismisses some of our concerns on the representation of uniform and independent random numbers; if a sequence of such random numbers should exist at all, it should be represented by a multiplicatic congruential sequence. Efforts may be safely concentrate on good multiplicative congruential RNG for uniform and independent random numbers with a long period required.
The present invention provides a method to obtain uniform and independent random numbers with long periods and certified statistical quality, which may be realized in operating systems on computers with the smallest computational cost.
The invention claims a method of obtaining uniform and independent random numbers comprising:
1. two odd distinct primes p1 and p2, and their integer exponents i1≧1 and i2≧1, which are restricted in the ways (1a)-(1c) noted below,
(1a) the integer q1:=(p1−1)/2 is odd and the integer q2:=(p2−1)/2 is even,
(1b) integers q1{p1̂(i1−1)} and q2 {p2̂(i2−1)} are coprime,
(1c) the case i1=i2=1 is excluded,
2. primitive roots z1 and z2 of prime-power moduluses d1:=p1̂i1 and d2:=p2̂i2, respectively,
3. the modulus d=d1d2
4. the multiplier z coprime with d=d1d2 and determined by the system of congruence relations
z≡z
1 mod(d1=p1̂i1),z≡z2 mod(d2=p2̂i2) Eq. 3
or by the system of congruence relations
z≡−z
1 mod(d1=p1̂i1),z≡z2 mod(d2=p2̂i2), Eq. 4
both system of which specify z uniquely in the sense of modulo d,
5. an arbitrary integer n coprime with d.
The method constructs a sequence of integers {r1, r2, . . . } recursively by
r
1
≡n,r
j+1
≡zr
j mod(d),0<rj<d,j=1,2, . . . Eq. 5
and gives the output sequence {v1, v2, . . . } for uniform and independent random numbers by realizing the arithmetic,
v
j
:=r
j
/d,0<vj<1,j=1,2, . . . Eq. 6
Structures of the multiplicative congruential generator (d, z, n) as summarized in the preceding paragraph ensure that the emitted random number sequence has the usable period
T′=2q1{p1̂(i1−1)}×q2{p2̂(i2−1)} Eq. 7
and realize the optimal efficiency measure τ:=(the usable period/the modulus)=T′/d≈½;. The generator (d, z, n) with a modulus formed by two odd-prime-powers has another great merit originating from the inherent shuffling mechanism. For, the sequence {n, nz, nz2, . . . } satisfies
Nz
j−1
≡n
k(zk)j−1 mod(dk),j=1,2, . . . ,nk:≡n mod(dk),k=1,2, Eq. 8
and the following expression holds true by Sun Tzu's theorem,
nz
j−1
≡U
1
n
1(z1)j−1+U2n2(z2)j−1 mod(d),j=1,2, . . . , Eq. 9
with integers U1, U2 determined by component moduluses d1 and d2. The generator (d, z, n) may generate a sequence of random integers by shuffling sequences from component generators (d1, z1, n1) and (d2, z2, n2). In the event that component generators (d1, z1, n1) and (d2, z2, n2) are chosen by their respective spectral tests, the likelihood of hitting a good generator (d, z, n) as the synthesis is expected technologically to be well increased. Shuffling may be effective when d1, d2 are similar (but not too close) in magnitude, so that the adequate choice will be d1≈d2≈d1/2, with d1/2<<d in cases of practical importance. This setting facilitates decisively the computation of spectral tests of component multipliers, and of the combined multiplier as well.
Understanding of the present invention will be facilitated by consideration of the following detailed description of the preferred embodiments of the present invention taken in conjunction with the accompanying drawings, in which like numerals refer to like parts:
It is to be understood that the figures and descriptions of the present invention have been simplified to illustrate elements that are relevant for a clear understanding of the present invention, while eliminating, for the purpose of clarity, many other elements found in computing systems and number generators. Those of ordinary skill in the art may recognize that other elements and/or steps are desirable and/or required in implementing the present invention. However, because such elements and steps are well known in the art, and because they do not facilitate a better understanding of the present invention, a discussion of such elements and steps is not provided herein. The disclosure herein is directed to all such variations and modifications to such elements and methods known to those skilled in the art.
In order to keep the clarity of notations the power pi will be denoted as p̂i when confusion is feared. The present invention has its concern in the use of a specific form of the modulus d=(p1̂i1)×(p2̂i2), where p1 and p2 are odd primes and exponents i1≧1 and i2≧1 are integers. The integer q1:=(p1−1)/2 is assumed odd and the integer q2:=(p2−1)/2 is assumed even. The multiplier z is defined by the system of congruence relations:
z≡±z
1 mod(p1̂i1),z≡z2 mod(p2̂i2), Eq. 10
where z1 is a primitive root (or a generating element) modulo the odd-prime-power modulus p1̂i1 and z2 is a primitive root modulo the odd-prime-power modulus p2̂i2. Congruence relations in this system determine z uniquely modulo d by the theorem of Sun Tzu (or by the Chinese remainder theorem). The setting i1=i2=1 and the cases of the specification
z≡±z
1 mod(p1),z≡z2 mod(p2), Eq. 11
was submitted for patents to USPTO already as Ser. Nos. 12/379,964 and 13/105,351. The present invention generalizes exponents i1 and i2 to all cases of i1>1 and i2>1, with the exclusion of the case i1=i2=1.
Computers invariably have finite amounts of memory, and any generators of random numbers are destined to give only finite periods, if they are to emit random numbers in any reproducible way. In order for a generating method of random numbers to be of technological value, its period T should be large enough to sustain simulations. In this regard, consider a multiplicative congruential generator (d, z, n) with an odd prime modulus d=p and its primitive root multiplier z. The period T of the sequence is given as T=φ(p)=p−1 by Fermat's little theorem with the Euler's function φ. This stipulates that zT/2≡−1 mod(p) should be true, that the last half of the sequence (as viewed with the change of sign) is identical with the first half, and that only the first half is usable for independent random numbers. Any multiplicative congruential sequence generated by (d, z, n) has its usable period T′ for the independent random number sequence as T′:=T/2, if the cyclic sequence {1, z, z2, . . . } modulo d has an even period T and gives zT/2≡−1 mod(d). An even T does not always imply zT/2≡−1 mod(d). The factor τ:=T′/d characterizes the efficiency of a multiplicative congruential generator (d, z, n), because the computation of consecutive random numbers becomes heavy in proportion to d and also by statistical reasons to be described later. There exists a general upper bound for this efficiency in the form that z cannot exceed ½. (d, z, n) generators aimed in the present invention all realize the largest τ≈½.
A method of generating uniform and independent random numbers is given by comprising two distinct odd primes that give an odd integer and an even integer, together with by taking an integer exponent and an integer exponent, by forming the composite modulus by taking a primitive root modulo and a primitive modulo and giving the multiplier modulo by either the system of congruence relations, any of which determines the multiplier modulo uniquely, by taking an initial value coprime. The method generates the sequence of integers by recursive congruence relations and gives an output of uniform and independent random numbers.
Let p1 and p2 be distinct odd primes, exponents i1≧1 and i2≧1 be integers, the integer z1 be a primitive root modulo d1:=pîi1 and the integer z2 be a primitive root modulo d2:=p2̂i2. Define integers q1:=(p1−1)/2 and q2:=(p2−1)/2, and assume that q1 is odd and q2 is even. Assume that {p1̂(i1−1)}q1 and {p2̂(i2−1)}q2 are taken mutually coprime. Define finally d:=d1d2=(p1̂i1)×(p2̂i2). Then the following assertions hold true.
First, the integer z specified modulo d uniquely by the system of congruence relations
z≡z
1 mod(d1=p1̂i1),z≡z2 mod(d2=p2̂i2), Eq. 12
gives, together with an arbitrary initial value n coprime to d, the multiplicative congruential generator (d, z, n) that realizes the period T=2{p1̂(i1−1)}q1×{p2̂(i2−1)}q2 with the largest efficiency measure τ≈½.
Second, the integer z specified modulo d uniquely by the system of congruence relations
z≡−z
1 mod(d1=p1̂i1),z≡z2 mod(d2=p2̂i2), Eq. 13
gives, together with an arbitrary initial value n coprime to d, the multiplicative congruential generator (d, z, n) that realizes the period T=2{p1̂(i1−1)}q1×{p2̂(i2−1)}q2 with the largest efficiency measure τ≈½, though −z1 is not a primitive root modulo d1=p1̂i1.
Since component moduluses d1=p1̂i1 and d2=p2̂i2 are coprime mutually, Sun Tzu's theorem gives that z is determined modulo d=d1d2 uniquely by the noted system of congruence relations, and that an explicit form of z is given by
z≡U
1
z
1
+U
2
z
2 mod(d),Uj=δjk mod(dk), Eq. 14
with integers U1 and U2 depending only on component moduluses d1=p1̂i1 and d2=p2̂i2. The cyclic sequence {1, z, z2, . . . } generated by z modulo d is seen at once to be a shuffling of component cyclic sequences and gives,
z
j−1
≡U
1
z
1
j−1
+U
2
z
2
j−1 mod(d),j=1,2, . . . , Eq. 15
with the least common multiple period T that is even,
T=LCM(T1,T2)=T1T2/2,T1:=2q1{p1̂(i1−1)},T2:=2q2{p2(i2−1)}. Eq. 16
Therefore, zt≡−1≡d−1 mod(d) can arise, if at all, only at at t=T/2=T1T2/4. The assumption that q2 is even gives an integer for T2/4, and we have
z
T/2
≡z
1̂(T1T2/4)≡{z1̂(T1)}̂(T2/4)≡1 mod(d1=p1̂i1), Eq. 17
which proves that zT/2≠1 mod(d). Hence the whole length T of the sequence may be used for independent random numbers, and we have
τ=T/d=(T1T2/2)/{(p1̂i1)×(p2̂i2)}≈½. Eq. 18
Since z1 is a primitive root modulo p1̂i1, the cyclic sequence {1, z1, z12, . . . } has the following even period T1 in the same modulus,
T
1:=φ(p1̂i1)=p1̂(i1−1)φ(p1)=2q1{p1̂(i1−1)},q1:=(p1−1)/2, Eq. 19
and z1̂(T1/2)≡−1 mod(p1̂i1) holds true. However, the assumption that q1 is odd implies T1/2=q1{p1̂(i1−1)} is odd, and
(−z1)̂(T1/2)≡−(−1)≡1 mod(p1̂i1). Eq. 20
Therefore, for the case i1=1, the cyclic sequence {1, −z1, (−z1)2, . . . } modulo p1̂i1 realizes the period T1/2 and −z1 is not a primitive root modulo p1̂i1. For all these complications the multiplier z characterized by
z≡−z
1 mod(p1̂i1),z≡z2 mod(p2̂i2), Eq. 21
has its cyclic sequence modulo d with the least common multiple period T,
T:=LCM(T1/2,T2)=2q1{p1̂(i1−1)}×q2{p2̂(i2−1)}, Eq. 22
which is a multiple of 4 by the assumption that q2 is even. Therefore, T/2 is even, and zT/2 modulo p1̂i1 is equivalent to an even power of (−z1)̂T1. This proves zT/2≡1 mod(p1̂i1) together with zT/2≠−1 mod(d). Thus, the whole length T of the cyclic sequence generated by z may be used for independent random numbers, implying τ=T/d≈½.
The largest efficiency τ≈½ of a multiplicative congruential generator (d, z, n) is desirable but not sufficient for the generated random numbers to have excellent statistics. The multiplicative congruential generator (d, z, n) with a Mersenne prime modulus d=p=231−1 may be tested exhaustively for all of its primitive roots for the multiplier z and was found by Fishman and Moore (1986) that 7.669×10−5% of primitive roots pass the criterion. Significant implications of this small number will be that some tests over the whole length of the random number sequence are indispensable for us to grasp at excellent random number sequences out of the mediocre or bad multitudes, and that the circumstance will be universal to any sort of RNGs. Yet, merits of RNGs of multiplicative congruential type are that these may be run on computers with the smallest computational load, may readily be tested spectrally, and may be designed to have large efficiency τ≈½, not to mention the stated theoretical significance that they represent any sequences from RNGs on computers.
Let (d, z, n) be a multiplicative congruential generator with an arbitrary modulus d and an arbitrary multiplier z coprime with d. Denote {r1, r2, . . . } for the output sequence of integers defined at this place simply by
r
1
≡n mod(d),rk+1≡zrk mod(d),k=1,2, . . . Eq. 23
Denote the L-tuple starting from rm as Pm:≡(rm, rm+1, . . . , rm+L−1) and regard it as a point in the L-dimensional Euclidean space EL; Pm will also be noted loosely as a position vector. These points {P1, P2, . . . }, together with all points that have coordinates congruent to them modulo d, occupy a portion of lattice points in EL. Define L vectors in EL by the integer coordinates,
The point or the position vector Pm may be denoted Pm≡nzm−1e1 mod (d). Additions of vectors e2, e3, . . . , eL realize geometrically d-translations along coordinate axes 2, 3, . . . , L in EL. As regards the d-translation along the first axis the addition of the vector e1′ formed by an integral linear combination
e
1
′:=de
1
−ze
2
−z
2
e
3
− . . . −z
L−1
e
L=(d,0,0, . . . ,0,0) Eq. 29
plays the role. Thus, points with position vectors {Pm|m=1, 2, . . . } and the d-translations along axes are all expressed as integral linear combinations of vectors {e1, e2, . . . , eL} which are manifestly linearly independent (i.e. with positive determinant dL−1). These are all in the lattice spanned by these basis vectors. This lattice may be denoted as G(e1, e2, . . . , eL); we employ the simplification GL(d, z), because basis vectors are determined by d and z.
Denote CL(d) for an L-dimensional hypercube in El with sides 0≦xj<d along the j-th coordinate axes. The number of lattice points of GL(d, z) in the hypercube CL(d) is d. Take any lattice vector x=c1e1+c2e2+ . . . +CLeL with integers c1, c2, . . . , CL. The first coordinate of x is c1, and there exist d values of c1 that satisfy 0≦c1<d. Fix c1 at a value in this range. The second coordinate of x is then x2=c1z+c2d, and the requirement 0≦x2<d determines c2 uniquely by c1 and d. Integers c3, c4 . . . are likewise determined uniquely. Hence precisely d lattice vectors x or lattice points exist in the cube CL(d).
Now let (d, z, n) be any multiplicative congruential generator. For any dimension L the consecutive L-tuples from (d, z, n) in its usable period can occupy only less than half of points of GL(d, z) lattice, because of the efficiency measure z=(the usable period/d)=T′/d<½.
Spectral tests give valuations on the geometrical configuration of the lattice GL(d, z) introduced above, usually for 2≦L≦6. Inspections of possible lattice configurations for cases with L=2 in the Euclidean plane will give the intuitive comprehension of the content of spectral tests. If d lattice points in the hypercube CL(d) is distributed uniformly and isotropically to all directions, then the valuation of spectral test is high. A complication is that this valuation is not so easy with basis vectors or bases {e1, e2, . . . , eL} of the lattice GL(d, z). Rather, a different set of integer basis vectors, {f1, f2, . . . , fL} that form the reciprocal lattice, is more convenient. Basis vectors or bases of the reciprocal lattice are defined by inner products with basis vectors:
(ek,fl)=dδkl,1≦k,l≦L. Eq. 30
This implies that the matrix formed by {f1, f2, . . . , fL} transposed is the inverse matrix divided by d of the matrix formed by {e1, e2, . . . , eL}. Explicit forms of reciprocal bases are:
The reciprocal lattice consists of vectors that are integral linear combinations of these bases, typically v=c1f1+c2f2+ . . . +CLfL with any set of integers c1, c2, . . . , CL. Since all coordinates of v are integers, there exists the shortest vector among non-zero vectors with length greater or equal to 1. Denote bL(d, z) for this positive shortest length. Some geometrical restrictions intervene in the Euclidean space EL that give the theoretical upper bound BL(d) for this shortest value bL(d, z). This upper bound is listed in Table I below.
Since BL(d) refers to the geometrically ideal reciprocal lattice which (along with the geometrically ideal original lattice) requires irrational Cartesian components for their description, there holds bL(d, z)<BL(d) invariably.
The valuation of spectral tests may now be defined. The ratio is
μL(d,z):=BL(d)/bL(d,z)>1,2≦L≦6. Eq. 36
If this ratio μL(d, z) is closer to 1, the generator (d, z, n) has fewer reasons to be denied of the statistical uniformity and independence of random numbers they generate. Implications of values in L=2 cases,
μ2(d,z)≈1.05,1.10,1.15, . . . ,1.50, Eq. 37
have been depicted by Nakazawa and Nakazawa (2011), in “Spectral tests of primitive roots for primes up to 101027 and of multiplicative congruential random number sequences generated by pairs of primitive roots.” A glance over them will give clear geometrical image behind the above noted valuations. Fishman and Moore initiated the use of the criterion
μL(d,z)<1.25,2≦L≦6, Eq. 38
for the excellence of a primitive root multiplier z with respect to the Mersenne prime modulus d=231−1=2147483647, which is highly versatile and useful beyond this special case. Examinations of multipliers under this or other criterion require the search of the shortest reciprocal lattice vectors. Expressions of reciprocal lattice vectors with bases {f1, f2, . . . , fL} are yet not so convenient to this end; the expression of the reciprocal lattice vector v=(j1, j2, . . . , jL) with integer Cartesian coordinates should be employed by its simple relation to ∥v∥=(j12+j22+ . . . +jL2)1/2, the vector length. It is not that every such integer vector v belongs to the reciprocal lattice, but the discrimination may be performed readily. The congruence relation
j
1
+zj
2
+z
2
j
3
+ . . . +z
L−1
j
L≡0 mod(d), Eq. 39
is necessary and sufficient for an integer vector v=(j1, j2, . . . , jL) to be in the reciprocal lattice of dimension L≧2 for the modulus d and the multiplier z.
For, the reciprocal lattice vector has the expression v=c1f1+c2f2+ . . . +cLfL by a set of integers {c1, c2, . . . , cL}. Its Cartesian coordinates are
v=(j1,j2,j3, . . . ,jL),j1=dc1−zc2−z2c3− . . . −zL−1cL, Eq. 40
with j2=c2, j3=c3, . . . , jL=CL. This implies
j
1
+zj
2
+z
2
j
3
+ . . . +z
L−1
j
L
=dc
1≡0 mod(d). Eq. 41
Therefore, the congruence relation Eq. 40 is a necessary condition for v to be in the reciprocal lattice. Conversely, if Cartesian coordinates of v=(j1, j2, . . . , jL) fulfill this congruence relation, some integer j′ exists and gives
j
1
+zj
2
+ . . . +z
L−1
j
L
=j′d,j
1
=j′d−zj
2
−z
2
j
3
− . . . −z
L−1
j
L, Eq. 42
v=(j1,j2, . . . ,jL)=j′f1+j2f2+ . . . +jLfL. Eq. 43
Hence v is a reciprocal lattice vector, and the sufficiency follows.
In the usual circumstance for a multiplicative congruential generator (d, z, n) the modulus d is a very large integer. Therefore, d1/L<<d holds true for any L=2, 3, . . . . The geometrical restriction, that BL(d) z O(d1/L)>bL(d, z) holds for the shortest length bL(d, z) of the L-dimensional reciprocal lattice vector v=(j1, j2, . . . , jL) in Cartesian coordinates, allows us to search for bL(d, z) only over integers j1, j2, . . . , jL satisfying Eq. 40 within the range
−BL(d)<jk<BL(d),1≦k≦L,2≦L≦6. Eq. 44
In search for the shortest reciprocal lattice vector v=(j1, j2, . . . , jL) in Cartesian coordinates, therefore, the condition j1+zj2+z2j3+ . . . +zL−1jL≡0 mod(d) may be solved for j1 as j1≡−zj2−z2j3− . . . −zL−1jL mod(d). Noting that v and −v have the same length, a restriction of 0≦j1≦BL(d)≈d1/L may be performed. The search of such shortest vector should be on all combinations of L−1 integers j2, j3, . . . , jL sweeping over approximately 2BL(d) values in the interval (−BL(d), BL(d)), namely over the total number NL(d):={2BL(d)}L−1 of integer sets. They are listed below in Table II.
Consideration of the number of candidate multipliers to be tested should be given. We first simplify the problem to the case of Fishman and Moore, the odd prime modulus d=p and its primitive roots. There may be some symmetries. A primitive root z and its inverse modulo p give the identical performance in spectral tests. And a primitive root z gives −z that has the identical valuation in spectral tests, even though the latter is not a primitive root in the case that the integer q=(p−1)/2 is odd; note that both of +z give the same performance measure τ=(usable period)/p≈½. The spectral tests should sweep over all primitive roots of the modulus d=p. The expression of all primitive roots by powers of one of them readily shows that an odd prime d=p has φ(φ(p))=φ(p−1) primitive roots. A general evaluation of this Euler's function φ(p−1) is not easy. Since p−1 is even, however, the estimate φ(p−1)≦(p−1)/2≈d/2 is certain. There may of course be an odd prime with many prime factors for p−1, and φ(p−1) is far smaller. However, the choice of such a prime modulus may imply too few primitive roots to have abundant excellent multipliers. A rough estimate φ(p−1)≈p/2 provides Mersenne prime p=231−1 with φ(p−1)≈0.249p. Table III below tabulates the number OL(d=p)=φ(p−1)NL(d=p) of integer sets to be examined in the L-th degree spectral tests.
Take now a modulus d=d1d2, dk=(pk)̂ik, k=1, 2. The strategy is as follows.
1. Design the magnitude of d z 2T′ for the usable period T′ so as to be sufficient for needs of planned simulations.
2. Choose odd primes p1, p2 and their exponents i1, i2 so as for d1≈d2≈d1/2 to hold. Values of d1 and d2 too near should be avoided. Also, primes p1, p2 should be chosen in such a way that q1=(p1−1)/2 is odd and q2=(p2−1)/2 is even.
3. Perform exhaustive spectral tests of degrees L=2, 3, . . . , 6 over every primitive root z1 of d1, select top N excellent ones, and save them in memory. N≈210 might suffice for practical use. Perform likewise exhaustive spectral tests for every primitive root z2 of d2, select, and store N excellent ones.
4. Take consecutively N2 pairs (z1, z2) of selected excellent primitive roots, and solve the system of congruence relations
z≡±z
1 mod(d1),z≡z2 mod(d2) Eq. 45
to obtain the multiplier z for the modulus d=d1d2. Perform second stage spectral tests of the generator for the modulus d and the multiplier z, starting from the 2nd degree and recording passers, then moving to third degree tests over remainders, and so forth.
In the first stage exhaustive spectral tests of this strategy, a primitive root multiplier z1 under the modulus d1=(p1)̂(i1)≈d1/2 involves NL(d1)≈NL(d1/2) integer sets to be examined. And the total number of primitive roots is
φ(φ(d1))=qφ(φ{(p1)̂(i1−1)(p1−1)})={(p1)̂(i1−2)}(p1−1)φ(p1−1)≈(p1)̂(i1−1)φ(p1−1)<d1/2. Eq. 46
Since the initial value is irrelevant in spectral tests, we denote briefly (d, z) for the (d, z, n) generator. Approximations dk≈d1/2 and φ(φ(dk))≈dk/2 Z d1/2/2 for k=1, 2 give Table IV below for the total number PL(d) of integer sets to be examined in the first stage L-th degree spectral tests.
Finally, consider the 2nd stage spectral test. Assume a choice of the best-N primitive roots of the modulus d1=(p1)̂(i1)≈d1/2 and best-N primitive roots of the modulus d2=(p2)̂(i2)≈d1/2. N≈210 will presumably be sufficient. Taking N2 pairs (z1, z2) of primitive roots to construct the multiplier z under the modulus d=d1d2 and proceed to L-th degree spectral test for L=2, 3, . . . , 6. Tabulation of the number of sets of integers to be swept is now easy; we need to multiply Table II with N2 to obtain Table V shown below.
The ratio FL(d):=PL(d)/OL(d) gives the reduction rate of the total number of integer sets to be swept over in the pair of first stage L-th degree spectral tests for the two odd-prime power modulus design versus the number of integer sets to be treated in the single odd prime modulus scheme. The ratio SL(d):=QL(d)/OL(d) similarly gives the rate of reduction of numbers of integer sets to be swept over in the second stage L-th degree spectral tests of N2 composite multipliers versus the single odd prime modulus case. These ratios are tabulated below in Table VI.
Recent supercomputers realize 1.6×1016=253.15 flops. If a computing program requires one random number in 23 floating point operations, such a computer will consume 250.15, 258.96 or 263.87 random numbers in a second, in a day or in a month, respectively. If these random numbers are to be supplied from a (d, z, n) generator with the efficiency τ≈½, the modulus d should be as large as d=265 at present. In the future, the speed of computers will increase, and the modulus d of a (d, z, n) generator should be increased in proportion. The matter may be looked at inversely. The number of cases of integers to be treated in spectral tests may be increased only in proportion to the computer speed. Hence the form of the dependence on d of OL(d) in Table III indicates computational difficulties arising with exhaustive spectral tests in single odd prime modulus schemes. In contrast, Tables IV and V reveal the recovery of computability by two odd-prime-power modulus designs. For all such achievements, however, all tables show that spectral tests become harder for larger L. Therefore, the preparation of two-odd-prime-power generators by spectral tests should be started from the degree L=2 and be discarded of non-passable multipliers preventing time-consumptive higher degree computations. Thus, passers of the desired high performance in L=2 tests should be tabulated in a file. Then the file should be read into subfiles, candidate multipliers in different subfiles should be processed in parallel in the test of the next degree L, and procedures should be repeated with L increased. Many other minute optimizations will simultaneously be needed: Multiplication of numbers should be replaced with addition if possible; use of function mod had better be avoided if simple subtraction suffices; and so forth. Spectral tests for larger moduluses are still a challenging subject. A simpler strategy to tighten the criterion for passers can be useful since more retirees will ease the burden of spectral tests. But this might be risky since a too tight criterion may extinguish passers.
The memory device 104 may be or include a device such as a Dynamic Random Access Memory (D-RAM), Static RAM (S-RAM), or other RAM or a flash memory. The data storage device 108 may be or include a hard disk, a magneto-optical medium, an optical medium such as a CD-ROM, a digital versatile disk (DVDs), or Blu-Ray disc (BD), or other type of device for electronic data storage.
The communication interface 106 may be, for example, a communications port, a wired transceiver, a wireless transceiver, and/or a network card. The communication interface 106 may be capable of communicating using technologies such as Ethernet, fiber optics, microwave, xDSL (Digital Subscriber Line), Wireless Local Area Network (WLAN) technology, wireless cellular technology, and/or any other appropriate technology.
The peripheral device interface 112 is configured to communicate with one or more peripheral devices. The peripheral device interface 112 operates using a technology such as Universal Serial Bus (USB), PS/2, Bluetooth, infrared, serial port, parallel port, and/or other appropriate technology. The peripheral device interface 112 may, for example, receive input data from an input device such as a keyboard, a mouse, a trackball, a touch screen, a touch pad, a stylus pad, and/or other device. Alternatively or additionally, the peripheral device interface 112 may communicate output data to a printer that is attached to the computing device 100 via the peripheral device interface 112.
The display device interface 110 may be an interface configured to communicate data to display device 110. The display device 110 may be, for example, a monitor or television display, a plasma display, a liquid crystal display (LCD), and/or a display based on a technology such as front or rear projection, light emitting diodes (LEDs), organic light-emitting diodes (OLEDs), or Digital Light Processing (DLP). The display device interface 110 may operate using technology such as Video Graphics Array (VGA), Super VGA (S-VGA), Digital Visual Interface (DVI), High-Definition Multimedia Interface (HDMI), or other appropriate technology. The display device interface 110 may communicate display data from the processor 102 to the display device 110 for display by the display device 110. As shown in
An instance of the computing device 100 of
As used herein, the term “processor” broadly refers to and is not limited to a single- or multi-core processor, a special purpose processor, a conventional processor, a Graphics Processing Unit (GPU), a digital signal processor (DSP), a plurality of microprocessors, one or more microprocessors in association with a DSP core, a controller, a microcontroller, one or more Application Specific Integrated Circuits (ASICs), one or more Field Programmable Gate Array (FPGA) circuits, any other type of integrated circuit (IC), a system-on-a-chip (SOC), and/or a state machine.
As used to herein, the term “computer-readable medium” broadly refers to and is not limited to a register, a cache memory, a ROM, a semiconductor memory device (such as a D-RAM, S-RAM, or other RAM), a magnetic medium such as a flash memory, a hard disk, a magneto-optical medium, an optical medium such as a CD-ROM, a DVDs, or BD, or other type of device for electronic data storage.
Although the methods and features are described above, the methods and features described above may be performed, mutatis mutandis, using any appropriate architecture and/or computing environment. Although features and elements are described above in particular combinations, each feature or element can be used alone or in any combination with or without the other features and elements. For example, each feature or element as described may be used alone without the other features and elements or in various combinations with or without other features and elements. Sub-elements and/or sub-steps of the methods described above may be performed in any arbitrary order (including concurrently), in any combination or sub-combination.
A system of generating uniform and independent random numbers for use within a computer system is disclosed. The system includes a processor and a communication device. The processor that operates to: take a positive integer d called modulus; take a positive integer z called multiplier coprime with d; take a positive integer n called initial value or seed coprime with d; and generates a sequence {r1, r2, . . . } by realizing congruence relations
r
1
≡n mod(d),rj+1≡zrj mod(d),0<rj<d,j=1,2, . . . Eq. 47
The communication device outputs the random number sequence {v1, v2, . . . } by realizing the arithmetic vj=rj/d for j=1, 2, . . . . The modulus d has the form of a product d=(p1)̂(i1)×(p2)̂(i2) of powers of distinct odd primes p1, p2 with exponents i1 and i2 that may take arbitrary integral values i1≧1 and i2≧1 excluding the case i1=i2=1, the odd prime p1 gives an odd integer q1=(p1−1)/2, the odd prime p2 gives an even integer q2=(p2−1)/2, the integers p1, q1, i1, p2, q2, i2 give mutually coprime integer q1 (p1)̂(i1−1) and integer q2(p2)̂(i2−1), the multiplier z is determined modulo d with a primitive root z1 modulo (p1)̂(i1) and with a primitive root z2 modulo (p2)̂(i2) either by congruence relations
z≡z
1 mod {(p1)̂(i1)},z≡z2 mod {(p2)̂(i2)} Eq. 48
or by congruence relations
z≡−z
1 mod {(p1)̂(i1)},z≡z2 mod {(p2)̂(i2)}, Eq. 49
and the modulus d and the multiplier z pass the 2nd and 3rd spectral test criterion of Fishman and Moore, which is the same as that the two dimensional vector v=(j1, j2) with integer coordinates fulfilling
j
1
+zj
2≡0 mod(d) Eq. 50
and having its length ∥v∥={(j1)2+(j2)2}1/2 has the smallest positive length b2(d, z) that satisfies b2(d, z)>(2d)1/2/(31/4μ) for μ=1.25, as well as that the three dimensional vector v=(j1, j2, j3), with integer coordinates fulfilling
j
1
+zj
2
+z
2
j
3≡0 mod(d) Eq. 51
and with its length ∥v∥:={(j1)2+(j2)2+(j3)2}1/2, has the smallest positive length b3(d, z) that satisfies b3(d, z)>21/6d1/3/μ.
Although the invention has been described and pictured in an exemplary form with a certain degree of particularity, it is understood that the present disclosure of the exemplary form has been made by way of example, and that numerous changes in the details of construction and combination and arrangement of parts and steps may be made without departing from the spirit and scope of the invention as set forth in the claims hereinafter.
