Generalized inversive congruential pseudorandom numbers

An approach to nonlinear congruential methods of generating uniform pseudorandom numbers in the interval [0,1) is the Inversive congruential generator with prime modulus. A generalization for arbitrary composite moduli [math]\displaystyle{ m=p_1,\dots p_r }[/math] with arbitrary distinct primes [math]\displaystyle{ p_1,\dots ,p_r \ge 5 }[/math] will be present here. Let [math]\displaystyle{ \mathbb{Z}_{m} = \{0,1,...,m-1\} }[/math]. For integers [math]\displaystyle{ a,b \in \mathbb{Z}_{m} }[/math] with gcd (a,m) = 1 a generalized inversive congruential sequence [math]\displaystyle{ (y_{n})_{n \geqslant 0} }[/math] of elements of [math]\displaystyle{ \mathbb{Z}_{m} }[/math] is defined by

[math]\displaystyle{ y_{0} = {\rm seed} }[/math]

[math]\displaystyle{ y_{n+1}\equiv a y_{n}^{\varphi(m)-1} + b \pmod m \text{, }n \geqslant 0 }[/math]

where [math]\displaystyle{ \varphi(m)=(p_{1}-1)\dots (p_{r}-1) }[/math] denotes the number of positive integers less than m which are relatively prime to m.

Example

Let take m = 15 = [math]\displaystyle{ 3 \times 5\, a=2 , b=3 }[/math] and [math]\displaystyle{ y_0= 1 }[/math]. Hence [math]\displaystyle{ \varphi (m)= 2 \times 4=8 \, }[/math] and the sequence [math]\displaystyle{ (y_{n})_{n \geqslant 0}=(1,5,13,2,4,7,1,\dots ) }[/math] is not maximum.

The result below shows that these sequences are closely related to the following inversive congruential sequence with prime moduli.

For [math]\displaystyle{ 1\le i \le r }[/math] let [math]\displaystyle{ \mathbb{Z}_{p_{i}}= \{0,1,\dots ,p_{i}-1\}, m_{i}= m / p_{i} }[/math] and [math]\displaystyle{ a_{i} ,b_{i} \in \mathbb{Z}_{p_{i}} }[/math] be integers with

[math]\displaystyle{ a\equiv m_{i}^{2} a_{i}\pmod {p_{i}} \;\text{and }\; b\equiv m_{i} b_{i}\pmod {p_{i}} \text{. } }[/math]

Let [math]\displaystyle{ (y_{n})_{n \geqslant 0} }[/math] be a sequence of elements of [math]\displaystyle{ \mathbb{Z}_{p_{i}} }[/math], given by

[math]\displaystyle{ y_{n+1}^{(i)}\equiv a_{i} (y_{n}^{(i)})^{p_{i}-2} + b_{i} \pmod {p_{i}}\; \text{, }n \geqslant 0 \;\text {where} \;y_{0}\equiv m_{i} (y_{0}^{(i)})\pmod {p_{i}}\;\text{is assumed. } }[/math]

Theorem 1

Let [math]\displaystyle{ (y_{n}^{(i)})_{n \geqslant 0} }[/math] for [math]\displaystyle{ 1\le i \le r }[/math] be defined as above. Then

[math]\displaystyle{ y_{n}\equiv m_{1}y_{n}^{(1)} + m_{2}y_{n}^{(2)} + \dots + m_{r}y_{n}^{(r)} \pmod m. }[/math]

This theorem shows that an implementation of Generalized Inversive Congruential Generator is possible, where exact integer computations have to be performed only in [math]\displaystyle{ \mathbb{Z}_{p_{1}},\dots , \mathbb{Z}_{p_{r}} }[/math] but not in [math]\displaystyle{ \mathbb{Z}_{m}. }[/math]

Proof:

First, observe that [math]\displaystyle{ m_{i}\equiv 0\pmod {p_{j}}, \; \text {for} \; i\ne j, }[/math] and hence [math]\displaystyle{ y_{n}\equiv m_{1}y_{n}^{(1)} + m_{2}y_{n}^{(2)} + \dots + m_{r}y_{n}^{(r)} \pmod m }[/math] if and only if [math]\displaystyle{ y_{n}\equiv m_{i} (y_{n}^{(i)})\pmod {p_{i}} }[/math], for [math]\displaystyle{ 1\le i \le r }[/math] which will be shown on induction on [math]\displaystyle{ n \geqslant 0 }[/math].

Recall that [math]\displaystyle{ y_{0}\equiv m_{i} (y_{0}^{(i)})\pmod {p_{i}} }[/math] is assumed for [math]\displaystyle{ 1\le i \le r }[/math]. Now, suppose that [math]\displaystyle{ 1\le i \le r }[/math] and [math]\displaystyle{ y_{n}\equiv m_{i} (y_{n}^{(i)})\pmod {p_{i}} }[/math] for some integer [math]\displaystyle{ n \geqslant 0 }[/math]. Then straightforward calculations and Fermat's Theorem yield

[math]\displaystyle{ y_{n+1}\equiv a y_{n}^{\varphi(m)-1} + b \equiv m_{i}(a_{i}m_{i}^{\varphi(m)}(y_{n}^{(i)})^{\varphi(m)-1}+ b_{i})\equiv m_{i}(a_{i}(y_{n}^{(i)})^{p_{i}-2} + b_{i})\equiv m_{i}(y_{n+1}^{(i)}) \pmod {p_{i}} }[/math],

which implies the desired result.

Generalized Inversive Congruential Pseudorandom Numbers are well equidistributed in one dimension. A reliable theoretical approach for assessing their statistical independence properties is based on the discrepancy of s-tuples of pseudorandom numbers.

Discrepancy bounds of the GIC Generator

We use the notation [math]\displaystyle{ D_m ^{s}=D_m (x_0,\dots , x_m-1) }[/math] where [math]\displaystyle{ x_n=(x_n, x_n+1 ,\dots , x_n+s-1) }[/math] ∈ [math]\displaystyle{ [0,1)^s }[/math] of Generalized Inversive Congruential Pseudorandom Numbers for [math]\displaystyle{ s\ge 2 }[/math].

Higher bound

Let [math]\displaystyle{ s\ge 2 }[/math]

Then the discrepancy [math]\displaystyle{ D_m ^s }[/math] satisfies

[math]\displaystyle{ D_m }[/math] [math]\displaystyle{ ^s }[/math] < [math]\displaystyle{ m^{-1/2} }[/math] × [math]\displaystyle{ ( \frac{2}{\pi} }[/math] × [math]\displaystyle{ \log m + \frac{7}{5})^s }[/math] × [math]\displaystyle{ \textstyle \prod_{i=1}^r (2s-2 + s(p_i)^{-1/2})+ s_{m}^{-1} }[/math] for any Generalized Inversive Congruential operator.

Lower bound:

There exist Generalized Inversive Congruential Generators with

[math]\displaystyle{ D_m }[/math] [math]\displaystyle{ ^s }[/math] ≥ [math]\displaystyle{ \frac{1}{2(\pi+2)} }[/math] × [math]\displaystyle{ m^{-1/2} }[/math] : × [math]\displaystyle{ \textstyle \prod_{i=1}^r (\frac {p_{i} - 3}{p_{i} - 1})^{1/2} }[/math] for all dimension s  :≥ 2.

For a fixed number r of prime factors of m, Theorem 2 shows that [math]\displaystyle{ D_m^{(s)} = O (m^{-1/2} (\log m )^s) }[/math] for any Generalized Inversive Congruential Sequence. In this case Theorem 3 implies that there exist Generalized Inversive Congruential Generators having a discrepancy [math]\displaystyle{ D_m^{(s)} }[/math] which is at least of the order of magnitude [math]\displaystyle{ m^{-1/2} }[/math] for all dimension [math]\displaystyle{ s\ge 2 }[/math]. However, if m is composed only of small primes, then r can be of an order of magnitude [math]\displaystyle{ (\log m)/\log\log m }[/math] and hence [math]\displaystyle{ \textstyle \prod_{i=1}^r (2s-2 + s(p_i)^{-1/2})= O{(m^\epsilon )} }[/math] for every [math]\displaystyle{ \epsilon\gt 0 }[/math]. Therefore, one obtains in the general case [math]\displaystyle{ D_m^{s}=O(m^{-1/2+\epsilon}) }[/math] for every [math]\displaystyle{ \epsilon\gt 0 }[/math].

Since [math]\displaystyle{ \textstyle \prod_{i=1}^r ((p_{i} - 3)/(p_{i} - 1))^{1/2} \geqslant 2^{-r/2} }[/math], similar arguments imply that in the general case the lower bound in Theorem 3 is at least of the order of magnitude [math]\displaystyle{ m^{-1/2 - \epsilon} }[/math] for every [math]\displaystyle{ \epsilon\gt 0 }[/math]. It is this range of magnitudes where one also finds the discrepancy of m independent and uniformly distributed random points which almost always has the order of magnitude [math]\displaystyle{ m^{-1/2} (\log\log m)^{1/2} }[/math] according to the law of the iterated logarithm for discrepancies. In this sense, Generalized Inversive Congruential Pseudo-random Numbers model true random numbers very closely.

See also

• Pseudorandom number generator
• List of random number generators
• Linear congruential generator
• Inversive congruential generator
• Naor-Reingold Pseudorandom Function

Notes

• Eichenauer-Herrmann, Jürgen (1994), "On Generalized Inversive Congruential Pseudorandom Numbers", Mathematics of Computation (American Mathematical Society) 63 (207): 293–299, doi:10.1090/S0025-5718-1994-1242056-8