good random permutation?

Started by
14 comments, last by Lode 9 years, 8 months ago

Hello,

I'm looking for a way to map a 31-bit positive integer (leaving some room for Java!) to another 31-bit positive integer, for purpose of random world generation, with the following properties:

-if you go through the inputs in sequence, the resulting output sequence should be a good pseudo-random sequence

-it should not require state. So a linear congruence generator or similar will not work afaik. Instead, it should be able to give you the output for any of the inputs at any time, no requirement to go through it in order.

-it should be fast, so cryptographic hashes are probably not it

-it should be stable: always give the same output for the same input (so no true random generator smile.png)

-each result is unique, that is, every possible integer is output by one input integer (not the most important requirement, but would be nice, each output value should have equal probability). So some arbitrary bit twiddling and multiplications won't do it either

I tried to implement something similar to the following one:

http://preshing.com/20121224/how-to-generate-a-sequence-of-unique-random-integers/

but it's actually not that great, I want to use a seed value as well, and while this one works great for some seeds, it does not for many as well and to make it even more annoying, it gives too similar patterns for seeds that are near each other (so if you generate for coordinates x and y a random world temperature with seed 1000 and a random height with seed 1001, they will be too correlated. Or something like that, in any case the patterns are ugly with nearby seeds).

Does anyone have suggestions for such a function?

Thanks a lot! smile.png

Advertisement
I don't think a stateless PRNG makes any sense. Is there not any way to store the state of a PRNG with the terrain generator?

Usually the go-to choice for alternative RNGs (if not already provided by the language's libraries) is usually mersenne twister.

Sounds like you're looking for a perfect hash.

I don't think a stateless PRNG makes any sense. Is there not any way to store the state of a PRNG with the terrain generator?

Usually the go-to choice for alternative RNGs (if not already provided by the language's libraries) is usually mersenne twister.

He needs a random-access permutation, which is a different problem. But, yes, unless the seed is hardcoded, there is going to be state (by definition).

One possible candidate I can think of is the function:

\[f(x) = \begin{cases} x ~ ~ &\text{if} ~ ~ x = 0 ~ ~ \text{or} ~ ~ x = 2^{31} - 1 \\ g^x \bmod{2^{31} - 1} ~ ~ &\text{otherwise} \end{cases}\]

which has pretty good pseudorandom properties and is a permutation (there are a couple obvious fixed points by construction, which can be eliminated by e.g. multiplying the input with a suitable constant beforehand). Here g is supposed to be a generator of 2^31 - 1, the smallest one is 7 and about 25% of all possible bases are generators, so you won't run out of seeds (it is easy to test if you have a generator, which means you can efficiently find one by instantiating any PRNG with your seed and looking for a generator that way). It is also reasonably fast, because modular exponentiation can be implemented efficiently using modular square and multiply, and, furthermore, reducing an integer modulo 2^31 - 1 thankfully doesn't require a division. I might write up a full implementation this afternoon if anyone is interested (or just for fun)

“If I understand the standard right it is legal and safe to do this but the resulting value could be anything.”

There are a few fast operations you can do on an `unsigned' value that are permutations of Z/(2^32):
- add a constant
- XOR with a constant
- multiply by an odd constant
- rotate bits (e.g. `x = (x>>18) | (x<<14)')

By composing a few of these operations together I think you can make something that satisfies all of your requirements.

Since you want to do it with 31 bits instead of 32, you'll have to AND the result of each operation with 0x7fffffff to emulate 31-bit arithmetic. If you were to do it with 32 bits, it would be a bit faster.

Care to explain what you mean by "leaving some room for Java!"?

Thanks! Yeah I didn't really know the exact scientific name of what I was looking for which made googling hard. Perfect hash, and random access random permutation, sound indeed like it.

It needs to be stateless because the whole world will not be generated at once, it generates a part only when you reach it. E.g. if you reach a part with a coordinate x + 65536 * y = 5000, it needs to give you the value for 5000.

About Java: I'm not even using it but if I would the world should be the same. So, it doesn't have unsigned integers. I don't like that, at all, but it is what it is. 64-bit integers can of course store the 32-bit value, but truth is the upper value doesn't really matter a lot so going up to 2 billion seems fine.

Yes, it would be much easier if you had asked for a 32-bit permutation. I think Java does let you access all 32 bits (signed integers just use a bit for the sign bit, after all) but bitwise operations are apparently rather tricky. I always say signed arithmetic is evil, but Java does not let you choose here, so it's harder than it needs to be.

“If I understand the standard right it is legal and safe to do this but the resulting value could be anything.”

Here's my (heavily optimized) implementation of the algorithm I mentioned above, as a Java class with some test/benchmark code:

package permutation;

import java.util.*;

/**
 * Provides a deterministic pseudorandom permutation of [0, 2^31).
 *
 * WARNING: heavily tuned for the specific 2^31 - 1 case - do not edit
 *          the prime and expect it to work (you will at least have to
 *          update the reduce() method and the cofactor array and also
 *          make sure that the modPow() 2-ary algorithm still works)!
 */
public class Permutation {
    private static final int p = 2147483647;         /* Prime = 2^31 - 1  */
    private static final int s = 31;                 /* Mersenne exponent */
    private static final int[] cofactors = {         /* List of cofactors */
        (p - 1) / 2,
        (p - 1) / 3,
        (p - 1) / 7,
        (p - 1) / 11,
        (p - 1) / 31,
        (p - 1) / 151,
        (p - 1) / 331
    };

    /**
     * Computes the value of x mod p, where p = 2^s - 1 is prime.
     *
     * This is a division-free optimization of the usual x % p
     * method, with performance improvements of around 10-20%.
     *
     * Only works with p = 2147483647 (it can be implemented for
     * smaller Mersenne primes but it becomes less efficient, as
     * one needs to check more carefully for modular overflow).
     */
    private static long reduce(long x) {
        long r = (x & p) + (x >> s);
        return r >= p ? r - p : r;
    }

    /**
     * Computes the value of g^e mod p, using precomputed 2-ary
     * modular exponentiation. The generator array must contain
     * the 2-ary powers of the generator, beginning at g^1.
     */
    private static int modPow(int[] generator, int e) {
        long r = 1;
        int k = 0;
        
        while (e != 0) {
            if ((e & 1) != 0) {
                r = reduce(r * generator[k]);
            }

            e >>= 1;
            ++k;
        }

        return (int)r;
    }

    /**
     * Returns whether g (as a 2-ary power table) is a generator
     * modulo p. Does this efficiently via the factors of p - 1.
     */
    private static boolean isGenerator(int[] g) {
        if (g[0] <= 1) return false;

        for (int cofactor : cofactors) {
           if (modPow(g, cofactor) == 1)
               return false;
        }

        return true;
    }

    /**
     * Deterministically converts a seed to a generator, and
     * returns it as a 2-ary power table (starting at g^1).
     *
     * Note: a sizeable fraction of values in [2, p - 2] are
     *       generators of p so this will terminate quickly.
     */
    private static int[] seedToGenerator(int seed) {
        Random random = new Random(seed);
        int[] powers = new int[31];

        while (true) {
            long g = powers[0] = 2 + random.nextInt(p - 3);

            for (int k = 1; k < 31; ++k) {
                powers[k] = (int)(g = reduce(g * g));
            }

            if (isGenerator(powers))
                return powers;
        }
    }
    
    /**
     * The lower bound of the permutation (inclusive).
     */
    public static final int Lower = 0;

    /**
     * The upper bound of the permutation (inclusive).
     */
    public static final int Upper = p;

    /**
     * Precomputed values of g^(2^k) mod p, used for modular
     * exponentiation. The actual generator can be read from
     * generator[0] = g^1.
     */
    private int[] generator;

    /**
     * Instantiates a new permutation of [0, p] using
     * the given seed (same seed = same permutation).
     */
    public Permutation(int seed) {
        generator = seedToGenerator(seed);
    }

    /**
     * Permutes the input in [0, p] and returns the
     * corresponding output - will throw an illegal
     * argument exception if x is out of bounds.
     */
    public int permute(int x) {
        if ((x < Lower) || (x > Upper))
            throw new IllegalArgumentException();
        else {
            /* This is for randomizing the fixed points at 0 and p. Remember
             * that the composition of two permutations is a permutation. */
            x ^= generator[0];

            if ((x == 0) || (x == p))
                return x;
            else
                return modPow(generator, x);
        }
    }

    /**
     * Little test driver for the class.
     */
    public static void main(String[] args) {
        final int benchTrials = 5, benchSamples = 10000000;
        final int dupTrials   = 5, dupSamples   = 10000000;
        final int arbitrarySeed = 0xDEADBEEF;

        System.out.printf("Instantiating permutation on [%d, %d] (seed %d)\n",
                          Permutation.Lower, Permutation.Upper, arbitrarySeed);

        Permutation perm = new Permutation(arbitrarySeed);
        
        System.out.println("\nDisplaying first few values:");

        for (int x = 0; x < 20; ++x)
        {
            System.out.printf("%d => %d\n", x, perm.permute(x));
        }

        System.out.println("\nBenchmarking performance:");

        for (int t = 0; t < benchTrials; ++t) {
            long begin = System.nanoTime();

            for (int x = t * benchSamples; x < (t + 1) * benchSamples; ++x)
                perm.permute(x);

            double elapsed = (System.nanoTime() - begin) / 1000000000.0;
            System.out.printf(" [+] %.1f perms/s\n", benchSamples / elapsed);
        }

        System.out.println("\nChecking for duplicates:");

        for (int t = 0; t < dupTrials; ++t) {
            Set<Integer> dups = new HashSet<Integer>();
            int dupCount = 0;

            for (int x = t * dupSamples; x < (t + 1) * dupSamples; ++x) {
                int output = perm.permute(x);

                if (!dups.contains(output))
                    dups.add(output);
                else
                    ++dupCount;
            }

            if (dupCount > 0)
                System.out.printf(" [!] Found %d duplicates!\n", dupCount);
            else
                System.out.printf(" [+] No duplicates\n");
        }
    }
}

It meets all of your requirements:

- with the added xor operation the fixed points are gone, and the permutation can now be called pseudorandom (it is at least as good as an LCG, and probably much better - someone else can weigh in here)

- it can give you the output permuted integer for any input in the [0, 2^31 - 1] interval in constant time/memory

- it is rather fast (on the PC I am on I am getting between 15.5 and 17.5 million lookups per second across the entire input range, in Java)

- it is deterministic ("stable") and lets you specify a seed (note: not all seeds produce a unique permutation, but most do, i.e. you have around 1 billion seeds to play with)

- it is, in fact, a permutation (so each output corresponds to one, and only one, input)

You don't have to use it, but that's okay, I really enjoyed writing it anyway :)

“If I understand the standard right it is legal and safe to do this but the resulting value could be anything.”

I don't think you should have any trouble using 32-bit arithmetic in Java. For arithmetic operations (addition and multiplication), it makes no difference whether you are interpreting your bit patterns as unsigned or signed with two's complement representation. XOR only has one possible interpretation. When shifting, just make sure to use >>> instead of >> so the sign bit isn't propagated.

Disclaimer: I am not a Java programmer.

I'll throw in my minimalistic GLSL PRNG. I'm using this in my particle code, which requires very basic random inital values. I stole it from here and modified it to make use of the VAO's vertex ID to generate unique seeds for each particle. It uses two variables (fTime and iSeed) as components to calculate the initial seed when the shader is run. iSeed has to be set once (or it could be omitted), fTime is modified per-update anyway. Modifying the code to return a 32-bit integer should be trivial ;).


uniform float fTime;
uniform int iSeed;

uint seed = uint(fTime * 1000.0) + uint(gl_VertexID) + iSeed;

const float InverseMaxInt = 1.0 / 4294967295.0;

//returns a result in the range [0, 1]
float drnd()
{
    seed++;

    uint i=(seed^12345391u)*2654435769u;
    i^=(i<<6u)^(i>>26u);
    i*=2654435769u;
    i+=(i<<5u)^(i>>12u);
    return float(i) * InverseMaxInt;
}

EDIT: you mentioned you need to map one 32-bit int to another, in which case you can simple pass the previous value in as the seed. As for whether the result is always unique, I cannot say.

This topic is closed to new replies.

Advertisement