[Back] Blum Blum Shub (BBS) is used as a pseudo-random number generator (it is pseudo as it is not a truly random number, and where its randomisation depends on a random seed). It was created by Lenore Blum, Manuel Blum and Michael Shub in 1968. Blum Blum Shub uses the form of \(x_{n+1}=x_{n}^{2} \pmod M \), and where \(x_0\) is a random seed. The value of \(M\) is equal to \(pq\), and where \(p\) and \(q\) are prime numbers. These values of \(p\) and \(q\) are both congruent to 3 mod 4 (\(p=q=3 \pmod 4\)). The security of the method involves the difficulty in factorizing \(M\). It is slow, but is the strongest proven random number generator. For each step, we extract some of the information from \(x_{n+1}\), and which is typically the least significant bit. It would not be used within a cipher application, but could be used in key generation.

## Blum Blum Shub |

## BBS

The Blum Blum Shub (BBS) method is as pseudorandom number generator and was created by Lenore Blum, Manuel Blum and Michael Shub in 1968. It uses the form of:

\(x_{n+1}=x_{n}^{2} \pmod M \)

and where \(x_0\) is a random seed. The value of \(M\) is equal to \(pq\), and where \(p\) and \(q\) are prime numbers. These values of \(p\) and \(q\) are both congruent to 3 mod 4 (p=q=3 (mod4) ). What does that mean? Well when I take the values of p or q and divide them by 4, I will get a remainder of 3.

So, \(p=7\) is possible (as 7 divided by 4 is 1 remainder 3), and \(p=11\) is also possible (as 11 divided by 4 is 2 remainder 3). A value of 13 is not possible is at will be 3 remainder 1. Let's try a simple example in Python:

>>> p=7 >>> q=11 >>> M=p*q >>> x0=5 >>> x1=(x0**2)%M >>> x2=(x1**2)%M >>> x3=(x2**2)%M >>> x4=(x3**2)%M >>> print (x1,x2,x3,x4) 25 9 4 16

In this case we are using p=7 and q=11, and then a seed of \(x_0=5\), and the random sequence is 25, 9, 4 and 16. We normally just take the least significant bit, so the output would be '1100'. Note that the seed value \(x_0\) should not share a factor with either \(p\) or \(q\).

Some sample code is:

import sympy import random import re import sys x = 3*10**10 y = 4*10**10 seed = random.randint(1,1e10) def lcm(a, b): """Compute the lowest common multiple of a and b""" return a * b / gcd(a, b) def gcd(a,b): """Compute the greatest common divisor of a and b""" while b > 0: a, b = b, a % b return a def next_usable_prime(x): p = sympy.nextprime(x) while (p % 4 != 3): p = sympy.nextprime(p) return p p = next_usable_prime(x) q = next_usable_prime(y) M = p*q N = 1000 if (len(sys.argv)>1): N=int(sys.argv[1]) print "\np:\t",p print "q:\t",q print "M:\t",M print "Seed:\t",seed x = seed bit_output = "" for _ in range(N): x = x*x % M b = x % 2 bit_output += str(b) print bit_output print "\nNumber of zeros:\t",bit_output.count("0") print "Number of ones:\t\t",bit_output.count("1") xi='' print "\nPredicting 1st 10 with Euler's Theorem:" for i in range(10): xi += str((seed ** (2**i % lcm((p-1),(p-1))) % M) %2) print xi

A sample run is:

p: 30000000091 q: 40000000003 M: 1200000003730000000273 Seed: 4882516701 1000100001000010001100111000101110111110000111010100110001100101101110111111001111010001000010010100 Number of zeros: 52 Number of ones: 48 Predicting 1st 10 with Euler's Theorem: 1000100001

An interesting feature with the Blum-Blum-Shub method is that we can calculate the \(i\)th value using Euler's Theorem:

\(x_i = (x_0^{2^i \pmod {λ(M)}}) \pmod M\)

and where \(λ(M)\) is lcm(p-1,q-1).

This allows us to directly calculate any value in the sequence, and also run the method in parallel.