pythonnumpydecimalprecisionmpmath

How to generate seeded random floats of at least 20 decimal places in python?


I would like to generate seeded random floats of at least 20 decimal places in the range of -1 to 1.

import numpy as np
np.random.seed(2000)
np.random.uniform(low=-1, high=1, size=2)

This is what I would do in the past, but my DOE research requires more precision.

I've looked in to mpmath and decimal modules in python

import mpmath
mpmath.mp.dps = 22
xis = (mpmath.rand() - mpmath.mpf(1), mpmath.rand())

This does not work because I could not find a way to seed mpmath.mp.rand().

import numpy as np
from decimal import Decimal, getcontext

getcontext().prec = 20
xis = np.asarray([Decimal(f"{x}") x for x in np.random.uniform(low=-1, high=1, size=2)], dtype=np.float64)

Here, my floats are dependent on numpy, and has max to 16-7 decimal places.

How can I create random floats that have 20+ decimal places and are seedable?


Solution

  • First off, you will need ~70 bits of precision to correctly store a number with 20 decimal digits. A long double (float96 or float128) in IEEE754 format has either 80 or 106 bits of precision, both of which are adequate for the task.

    Unfortunately, numpy.random.Generator.random only supports float32 and float64.

    Based on the current python implementation of mpmath.rand, the source of randomness is just the standard library random module, specifically random.getrandbits. You should be able to get reproducible results by calling the usual random.seed. Unfortunately, the state is global: there does not appear to be a mechanism to use individual random.Random objects, e.g., for different threads.

    Keep in mind that you need to set the working precision, which can be done in a context manager:

    random.seed(my_seed)
    with mp.workdps(22):
        data = [mp.rand() for _ in range(my_size)]
    

    If you have a collection of scripts that are being imported as modules, a single call to random.seed(my_seed) in your main driver will suffice to initialize the entire program.

    If you are calling scripts externally, e.g. via os or subprocess, you will want to put the same line in each script. It's safer to use an import guard at that point, so you can use the script both ways without thinking about it:

    if __name__ == '__main__':
        random.seed(my_seed)