diff options
Diffstat (limited to 'src/common/pg_prng.c')
| -rw-r--r-- | src/common/pg_prng.c | 37 |
1 files changed, 36 insertions, 1 deletions
diff --git a/src/common/pg_prng.c b/src/common/pg_prng.c index e58b471cff0..c7bb92ede38 100644 --- a/src/common/pg_prng.c +++ b/src/common/pg_prng.c @@ -19,11 +19,17 @@ #include "c.h" -#include <math.h> /* for ldexp() */ +#include <math.h> #include "common/pg_prng.h" #include "port/pg_bitutils.h" +/* X/Open (XSI) requires <math.h> to provide M_PI, but core POSIX does not */ +#ifndef M_PI +#define M_PI 3.14159265358979323846 +#endif + + /* process-wide state vector */ pg_prng_state pg_global_prng_state; @@ -236,6 +242,35 @@ pg_prng_double(pg_prng_state *state) } /* + * Select a random double from the normal distribution with + * mean = 0.0 and stddev = 1.0. + * + * To get a result from a different normal distribution use + * STDDEV * pg_prng_double_normal() + MEAN + * + * Uses https://en.wikipedia.org/wiki/Box%E2%80%93Muller_transform + */ +double +pg_prng_double_normal(pg_prng_state *state) +{ + double u1, + u2, + z0; + + /* + * pg_prng_double generates [0, 1), but for the basic version of the + * Box-Muller transform the two uniformly distributed random numbers are + * expected to be in (0, 1]; in particular we'd better not compute log(0). + */ + u1 = 1.0 - pg_prng_double(state); + u2 = 1.0 - pg_prng_double(state); + + /* Apply Box-Muller transform to get one normal-valued output */ + z0 = sqrt(-2.0 * log(u1)) * sin(2.0 * M_PI * u2); + return z0; +} + +/* * Select a random boolean value. */ bool |