Skip to content Skip to sidebar Skip to footer

Picking A, C And M For Linear Congruential Generator

I am looking to implement a simple pseudorandom number generator (PRNG) that has a specified period and guaranteed no collisions for the duration of that period. After doing some r

Solution 1:

From Wikipedia:

Provided that c is nonzero, the LCG will have a full period for all seed values if and only if:

  1. c and m are relatively prime,
  2. a-1 is divisible by all prime factors of m,
  3. a-1 is a multiple of 4 if m is a multiple of 4.

You said you want a period of 48-1, so you must choose m≥48-1. Let's try choosing m=48-1 and see where that takes us. The conditions from the Wikipedia article prohibit you from choosing c=0 if you want the period to be m.

Note that 11, 47, 541, and 911 are the prime factors of 48-1, since they're all prime and 11*47*541*911 = 48-1.

Let's go through each of those conditions:

  1. For c and m to be relatively prime, c and m must have no common prime factors. So, pick any prime numbers other than 11, 47, 541, and 911, then multiply them together to choose your c.
  2. You'll need to choose a such that a-1 is divisible by all the prime factors of m, i.e., a = x*11*47*541*911 + 1 for any x of your choosing.
  3. Your m is not a multiple of 4, so you can ignore the third condition.

In summary:

  • m = 48-1,
  • c = any product of primes other than 11, 47, 541, and 911 (also, c must be less than m),
  • a = x*11*47*541*911 + 1, for any nonnegative x of your choice (also, a must be less than m).

Here's a smaller test case (in Python) using a period of 48-1 (which has prime factors 7 and 47):

deflcg(state):
    x = 1
    a = x*7*47 + 1
    c = 100
    m = 48**2 - 1return (a * state + c) % m

expected_period = 48**2 - 1
seeds = [5]
for i inrange(expected_period):
    seeds.append(lcg(seeds[-1]))
print(len(set(seeds)) == expected_period)

It outputs True, as it should. (If you have any trouble reading Python, let me know and I can translate it to JavaScript.)

Solution 2:

Based on Snowball's answer and the comments I've created a complete example. You can use the set == list comparison for smaller numbers. I could not fit 48^5-1 into memory.

To circumvent the a < m problem, I'm incrementing the target a few times to find a number where a is able to be < m (where m has duplicated prime factors). Surprisingly +2 is enough for a lot of numbers. The few extra numbers are later skipped while iterating.

import random

def__prime_factors(n):
    """
    https://stackoverflow.com/a/412942/6078370
    Returns all the prime factors of a positive integer
    """
    factors = []
    d = 2while n > 1:
        while n % d == 0:
            factors.append(d)
            n //= d
        d += 1if d * d > n:
            if n > 1: factors.append(n)
            breakreturn factors

def__multiply_numbers(numbers):
    """multiply all numbers in array"""
    result = 1for n in numbers:
        result *= n
    return result

def__next_good_number(start):
    """
    https://en.wikipedia.org/wiki/Linear_congruential_generator#c%E2%89%A00
    some conditions apply for good/easy rotation
    """
    number = start
    factors = __prime_factors(number)
    whilelen(set(factors)) == len(factors) or number % 4 == 0:
        number += 1
        factors = __prime_factors(number)
    return number, set(factors)

# primes < 100 for coprime calculation. add more if your target is large
PRIMES = set([2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41,
        43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97])

defcreate_new_seed(target):
    """be aware, m might become > target"""
    m, factors = __next_good_number(target)
    a = __multiply_numbers(factors) + 1# https://en.wikipedia.org/wiki/Coprime_integers
    otherPrimes = [p for p in PRIMES if p notin factors]
    # the actual random part to get differnt results
    random.shuffle(otherPrimes)
    # I just used arbitary 3 of the other primes
    c = __multiply_numbers(otherPrimes[:3])
    # first number
    state = random.randint(0, target-1)
    return state, m, a, c

defnext_number(state, m, a ,c, limit):
    newState = (a * state + c) % m
    # skip out of range (__next_good_number increases original target)while newState >= limit:
        newState = (a * newState + c) % m

    return newState

if __name__ == "__main__":
    target = 48**5-1
    state, m, a, c = create_new_seed(target)
    print(state, m, a, c, 'target', target)

    # list and set can't fit into 16GB of memory
    checkSum = sum(range(target))
    randomSum = 0for i inrange(target):
        state = newState = next_number(state, m, a ,c, target)
        randomSum += newState
    print(checkSum == randomSum) # true

LCG is quite fascinating and usable in things like games. You can iterate a giant list of things in a deterministic random order. Shuffeling and saving the whole list is not required:

defriter(alist):
    """ iterate list using LCG """
    target = len(alist)
    state, m, a, c = create_new_seed(target)
    for _ inrange(target):
        yield alist[state]
        state = next_number(state, m, a ,c, target)

It is easy to save the state in between iteration steps:

savedState = '18:19:25:6:12047269:20'print('loading:', savedState)
i, state, m, a, c, target = (int(i) for i in savedState.split(':'))
state = next_number(state, m, a, c, target)
i += 1print('i:', i, 'is random number:', state, 'list done:', i+1 == target)
print('saving:', '{}:{}:{}:{}:{}:{}'.format(i, state, m, a, c, target))

Post a Comment for "Picking A, C And M For Linear Congruential Generator"