# IEEE754 doubles and uniform random distribution in C

IEEE754 floating point numbers have 3 parts: a sign (1bit), an exponent (11bits) and a mantissa (52bits).

To translate into C:

``````    union {
struct {
unsigned long long m : 52;
unsigned long long e : 11;
unsigned long long s : 1;
};

double d;
} X;
``````

To make a uniformly distributed number between -1 and 1, we need:

• A random mantissa `m`
• An exponent `e` starts at 1022 (which is actually -1) and has 50% chance of being 1022, 25% chance to be 1021, 12.5% chance to be 1020...
• a random sign `s`

``````#include <stdint.h>

// Please don't use C rand()
extern uint32_t decent_rand();

double uniform_rand_double_between_minus_one_and_one() {
union {
struct {
unsigned long long m : 52;
unsigned long long e : 11;
unsigned long long s : 1;
};

double d;
} X;

uint32_t bit_divider = decent_rand();
// Random mantissa
X.m = decent_rand() + ((uint64_t)bit_divider<<32);
// Base exponent
X.e = 1022;

// Recover unused entropy
bit_divider >>= 20;
// Infinite loop, for speed
int i = 12;
for(;;) {
for(; i; --i) {
// Consume entropy one bit at a time
if(bit_divider & 1) {
// 50% chance to reduce
X.e -= 1;
bit_divider >>= 1;
// Handle underflow of the exponent
if(X.e == 0) goto very_unlikely;
} else {
// 50% chance to be done
bit_divider >>= 1;
// Set the sign
X.s = bit_divider & 1;
return X.d;
}
}
// We consumed all the entropy, we start again
bit_divider = decent_rand();
i = 32;
}
very_unlikely:
X.s = bit_divider & 1;
return X.d;
}
``````

Then we can shift that space into another slice of the space, but we the sign is a hurdle here. Let's have the function without the sign:

``````double uniform_rand_double_between_zero_and_one_opti() {
union {
struct {
unsigned long long m : 52;
unsigned long long e : 11;
unsigned long long s : 1;
};

double d;
} X;

uint32_t bit_divider = decent_rand();
// Random mantissa
X.m = decent_rand() + ((uint64_t)bit_divider<<32);
// Base exponent
X.e = 1022;

// Recover unused entropy
bit_divider >>= 20;
// Infinite loop, for speed
int i = 12;
for(;;) {
for(; i; --i) {
// Consume entropy one bit at a time
// Handle underflow of the exponent
if(bit_divider & 1 && X.e != 0) [[unlikely]] {
// 50% chance to reduce
X.e -= 1;
bit_divider >>= 1;
} else {
// 50% chance to be done
goto out;
}
}
// We consumed all the entropy, we start again
bit_divider = decent_rand();
i=32;
}
out:
X.s = 0;
return X.d;
}
``````
``````    double uniform_rand_between(double a, double b) {
return uniform_rand_double_between_zero_and_one()*(b-a) + a
}
``````

There you have it, a uniform random number generator for real numbers stored as IEEE754 doubles. #### Ludovic Lagouardette

Founder of NekoIT, Developer and Cloud architect

IEEE754 doubles and uniform random distribution in C