High performance wrapper for CSPRNG arc4random(3)

Question

This function is a wrapper for arc4random(3), a cryptographic pseudo-random number generator on macOS and BSDs, but you can also get it by installing libbsd-devel on most Linux distros.

It serves the same purpose as arc4random_uniform(3), the official wrapper for arc4random(3), i.e. to generate a nonnegative random integer less than range.

static inline unsigned random_uniform(unsigned range) {
 static uint32_t random_32b;
 static uint32_t full_32b;
 if (full_32b < range) {
 random_32b = arc4random();
 full_32b = UINT32_MAX;
 }
 unsigned result = random_32b % range;
 random_32b /= range;
 full_32b /= range;
 return result;
}

You can find the implementation of arc4random_uniform(3) on macOS at the bottom of this page. Basically, it makes use of the rejection method to produce uniform deviates and eliminate the modulo bias. On the other hand, my implementation tries to eliminate the modulo bias by keeping track of the maximum possible range of the remaining random bytes in full_32b. What's more, it saves the unused random bytes in a static variable for future use, as a call to arc4random(3) can be costly. It appears to be faster according to my benchmark, but I'm not sure if it's correct.

• esote: I recommend viewing OpenBSD's arc4random_uniform(3), which is far cleaner than the macOS one you link. I'd be interested in it's performance according to your benchmarking.

Answer 1

result = random % range

You can only do this if range is constrained to be much smaller than the maximum value of random. Imagine a 2-bit PRNG that returns a value from 0 to 3. When your function is invoked with range 3, return value 0 will be twice as common as 1 or 2.

random_32b /= range;
full_32b /= range;

This is going to exhibit obvious bias when the result has a small number of bits. Imagine a range of 3 and a 4-bit random pool. On exit, full_32b will equal 5 but random_32b will only equal 5 one time in 15.

---

There are three good ways to go about this:

1. • choose an acceptable margin of bias (expressed as a fraction);
   • divide range by that fraction;
   • select enough bits to make max(random) larger than the result of the division;
   • return a simple modulus.

For example, 1% bias is 0.01. If the range is 10, the random number needs a maximum larger than \10ドル/0.01 = 1000\$. That's ten bits.

2. • select enough bits to represent range
   • repeat until the result is smaller than range.

Example: for a range of 10, consume four bits and discard results larger than 9.

3. • select enough bits to represent range, plus a couple more, to represent a multiple of range. Call the multiplier \$i\$.
   • repeat until the result is smaller than \$i * \$range
   • divide result by \$i\$

This last approach combines the advanatages of first two: you get bias-free numbers and mostly predictable runtime.

Example: approach #2 with a range of 9 consumes four bits per attempt, and discards almost half of the results. With 5 bits, instead of 4, you keep results from 0-26 (i.e., range\$*3\$) and divide them by 3. Discard 27 or larger—that's only 5/32 attempts wasted. With 6 bits and range of 9, values 0-62 are usable (range\$*7\$), and only one possible value (63) will be discarded.

In all cases, consume bits by bit-shifting (or dividing by powers of two). Integer division with arbitrary denominators won't work.