diff --git a/src/bin/pgbench/pgbench.c b/src/bin/pgbench/pgbench.c
index b5c75ce1c6..06f1083d5a 100644
--- a/src/bin/pgbench/pgbench.c
+++ b/src/bin/pgbench/pgbench.c
@@ -185,7 +185,7 @@ int64		latency_limit = 0;
 char	   *tablespace = NULL;
 char	   *index_tablespace = NULL;
 
-/* random seed used when calling srandom() */
+/* random seed used to initialize the PRNG */
 int64		random_seed = -1;
 
 /*
@@ -279,14 +279,6 @@ typedef struct StatsData
 	SimpleStats lag;
 } StatsData;
 
-/*
- * Struct to keep random state.
- */
-typedef struct RandomState
-{
-	unsigned short xseed[3];
-} RandomState;
-
 /*
  * Connection state machine states.
  */
@@ -383,7 +375,7 @@ typedef struct
 	 * Separate randomness for each client. This is used for random functions
 	 * PGBENCH_RANDOM_* during the execution of the script.
 	 */
-	RandomState cs_func_rs;
+	pg_prng_state cs_func_rs;
 
 	int			use_file;		/* index in sql_script for this client */
 	int			command;		/* command number in script */
@@ -450,9 +442,9 @@ typedef struct
 	 * random state to make all of them independent of each other and therefore
 	 * deterministic at the thread level.
 	 */
-	RandomState ts_choose_rs;	/* random state for selecting a script */
-	RandomState ts_throttle_rs;	/* random state for transaction throttling */
-	RandomState ts_sample_rs;	/* random state for log sampling */
+	pg_prng_state ts_choose_rs;	/* random state for selecting a script */
+	pg_prng_state ts_throttle_rs;	/* random state for transaction throttling */
+	pg_prng_state ts_sample_rs;	/* random state for log sampling */
 
 	int64		throttle_trigger;	/* previous/next throttling (us) */
 	FILE	   *logfile;		/* where to log, or NULL */
@@ -832,30 +824,34 @@ strtodouble(const char *str, bool errorOK, double *dv)
 }
 
 /*
- * Initialize a random state struct.
+ * Initialize a pseudo-random state struct from main PRNG.
+ * This done at the beginning of process allows to have deterministic
+ * pgbench runs.
  */
 static void
-initRandomState(RandomState *random_state)
+pg_prng_state_init(pg_prng_state *rs)
 {
-	random_state->xseed[0] = random();
-	random_state->xseed[1] = random();
-	random_state->xseed[2] = random();
+	pg_prng_setstate_r(rs, pg_prng_u64());
 }
 
 /* random number generator: uniform distribution from min to max inclusive */
 static int64
-getrand(RandomState *random_state, int64 min, int64 max)
+getrand(pg_prng_state *random_state, int64 min, int64 max)
 {
 	/*
 	 * Odd coding is so that min and max have approximately the same chance of
 	 * being selected as do numbers between them.
 	 *
-	 * pg_erand48() is thread-safe and concurrent, which is why we use it
+	 * pg_prng_double_r() is thread-safe and concurrent, which is why we use it
 	 * rather than random(), which in glibc is non-reentrant, and therefore
 	 * protected by a mutex, and therefore a bottleneck on machines with many
 	 * CPUs.
+	 *
+	 * Note: pg_prng_double returns 52 bits (and pg_erand48 only 48 bits),
+	 * so this does not really work if more is needed between min & max.
+	 * Should rather use a modulo? Would be biased a little close to 64 bits.
 	 */
-	return min + (int64) ((max - min + 1) * pg_erand48(random_state->xseed));
+	return min + (int64) ((max - min + 1) * pg_prng_double_r(random_state));
 }
 
 /*
@@ -864,7 +860,7 @@ getrand(RandomState *random_state, int64 min, int64 max)
  * value is exp(-parameter).
  */
 static int64
-getExponentialRand(RandomState *random_state, int64 min, int64 max,
+getExponentialRand(pg_prng_state *random_state, int64 min, int64 max,
 				   double parameter)
 {
 	double		cut,
@@ -875,7 +871,7 @@ getExponentialRand(RandomState *random_state, int64 min, int64 max,
 	Assert(parameter > 0.0);
 	cut = exp(-parameter);
 	/* erand in [0, 1), uniform in (0, 1] */
-	uniform = 1.0 - pg_erand48(random_state->xseed);
+	uniform = 1.0 - pg_prng_double_r(random_state);
 
 	/*
 	 * inner expression in (cut, 1] (if parameter > 0), rand in [0, 1)
@@ -888,7 +884,7 @@ getExponentialRand(RandomState *random_state, int64 min, int64 max,
 
 /* random number generator: gaussian distribution from min to max inclusive */
 static int64
-getGaussianRand(RandomState *random_state, int64 min, int64 max,
+getGaussianRand(pg_prng_state *random_state, int64 min, int64 max,
 				double parameter)
 {
 	double		stdev;
@@ -912,13 +908,13 @@ getGaussianRand(RandomState *random_state, int64 min, int64 max,
 	do
 	{
 		/*
-		 * pg_erand48 generates [0,1), but for the basic version of the
+		 * pseudo-random output in [0,1), but for the basic version of the
 		 * Box-Muller transform the two uniformly distributed random numbers
 		 * are expected in (0, 1] (see
 		 * https://en.wikipedia.org/wiki/Box-Muller_transform)
 		 */
-		double		rand1 = 1.0 - pg_erand48(random_state->xseed);
-		double		rand2 = 1.0 - pg_erand48(random_state->xseed);
+		double		rand1 = 1.0 - pg_prng_double_r(random_state);
+		double		rand2 = 1.0 - pg_prng_double_r(random_state);
 
 		/* Box-Muller basic form transform */
 		double		var_sqrt = sqrt(-2.0 * log(rand1));
@@ -948,7 +944,7 @@ getGaussianRand(RandomState *random_state, int64 min, int64 max,
  * not be one.
  */
 static int64
-getPoissonRand(RandomState *random_state, double center)
+getPoissonRand(pg_prng_state *random_state, double center)
 {
 	/*
 	 * Use inverse transform sampling to generate a value > 0, such that the
@@ -957,7 +953,7 @@ getPoissonRand(RandomState *random_state, double center)
 	double		uniform;
 
 	/* erand in [0, 1), uniform in (0, 1] */
-	uniform = 1.0 - pg_erand48(random_state->xseed);
+	uniform = 1.0 - pg_prng_double_r(random_state);
 
 	return (int64) (-log(uniform) * center + 0.5);
 }
@@ -1035,7 +1031,7 @@ zipfFindOrCreateCacheCell(ZipfCache *cache, int64 n, double s)
  * Luc Devroye, p. 550-551, Springer 1986.
  */
 static int64
-computeIterativeZipfian(RandomState *random_state, int64 n, double s)
+computeIterativeZipfian(pg_prng_state *random_state, int64 n, double s)
 {
 	double		b = pow(2.0, s - 1.0);
 	double		x,
@@ -1046,8 +1042,8 @@ computeIterativeZipfian(RandomState *random_state, int64 n, double s)
 	while (true)
 	{
 		/* random variates */
-		u = pg_erand48(random_state->xseed);
-		v = pg_erand48(random_state->xseed);
+		u = pg_prng_double_r(random_state);
+		v = pg_prng_double_r(random_state);
 
 		x = floor(pow(u, -1.0 / (s - 1.0)));
 
@@ -1065,11 +1061,11 @@ computeIterativeZipfian(RandomState *random_state, int64 n, double s)
  * Jim Gray et al, SIGMOD 1994
  */
 static int64
-computeHarmonicZipfian(ZipfCache *zipf_cache, RandomState *random_state,
+computeHarmonicZipfian(ZipfCache *zipf_cache, pg_prng_state *random_state,
 					   int64 n, double s)
 {
 	ZipfCell   *cell = zipfFindOrCreateCacheCell(zipf_cache, n, s);
-	double		uniform = pg_erand48(random_state->xseed);
+	double		uniform = pg_prng_double_r(random_state);
 	double		uz = uniform * cell->harmonicn;
 
 	if (uz < 1.0)
@@ -1081,7 +1077,7 @@ computeHarmonicZipfian(ZipfCache *zipf_cache, RandomState *random_state,
 
 /* random number generator: zipfian distribution from min to max inclusive */
 static int64
-getZipfianRand(ZipfCache *zipf_cache, RandomState *random_state, int64 min,
+getZipfianRand(ZipfCache *zipf_cache, pg_prng_state *random_state, int64 min,
 			   int64 max, double s)
 {
 	int64		n = max - min + 1;
@@ -3565,7 +3561,7 @@ doLog(TState *thread, CState *st,
 	 * to the random sample.
 	 */
 	if (sample_rate != 0.0 &&
-		pg_erand48(thread->ts_sample_rs.xseed) > sample_rate)
+		pg_prng_double_r(&thread->ts_sample_rs) > sample_rate)
 		return;
 
 	/* should we aggregate the results or not? */
@@ -5126,12 +5122,11 @@ printResults(TState *threads, StatsData *total, instr_time total_time,
 	}
 }
 
-/* call srandom based on some seed. NULL triggers the default behavior. */
+/* init prng based on some seed. NULL triggers the default behavior. */
 static bool
 set_random_seed(const char *seed)
 {
-	/* srandom expects an unsigned int */
-	unsigned int iseed;
+	uint64 iseed;
 
 	if (seed == NULL || strcmp(seed, "time") == 0)
 	{
@@ -5139,7 +5134,7 @@ set_random_seed(const char *seed)
 		instr_time	now;
 
 		INSTR_TIME_SET_CURRENT(now);
-		iseed = (unsigned int) INSTR_TIME_GET_MICROSEC(now);
+		iseed = INSTR_TIME_GET_MICROSEC(now);
 	}
 	else if (strcmp(seed, "rand") == 0)
 	{
@@ -5155,7 +5150,7 @@ set_random_seed(const char *seed)
 		/* parse seed unsigned int value */
 		char		garbage;
 
-		if (sscanf(seed, "%u%c", &iseed, &garbage) != 1)
+		if (sscanf(seed, UINT64_FORMAT "%c", &iseed, &garbage) != 1)
 		{
 			fprintf(stderr,
 					"unrecognized random seed option \"%s\": expecting an unsigned integer, \"time\" or \"rand\"\n",
@@ -5165,9 +5160,12 @@ set_random_seed(const char *seed)
 	}
 
 	if (seed != NULL)
-		fprintf(stderr, "setting random seed to %u\n", iseed);
-	srandom(iseed);
-	/* no precision loss: 32 bit unsigned int cast to 64 bit int */
+		fprintf(stderr, "setting random seed to " UINT64_FORMAT "\n", iseed);
+
+	/* initialize main PRNG */
+	pg_prng_setstate(iseed);
+
+	/* keep value, casted as signed */
 	random_seed = iseed;
 	return true;
 }
@@ -5775,7 +5773,7 @@ main(int argc, char **argv)
 	for (i = 0; i < nclients; i++)
 	{
 		state[i].cstack = conditional_stack_create();
-		initRandomState(&state[i].cs_func_rs);
+		pg_prng_state_init(&state[i].cs_func_rs);
 	}
 
 	if (debug)
@@ -5862,10 +5860,7 @@ main(int argc, char **argv)
 	/* set default seed for hash functions */
 	if (lookupVariable(&state[0], "default_seed") == NULL)
 	{
-		uint64		seed = ((uint64) (random() & 0xFFFF) << 48) |
-		((uint64) (random() & 0xFFFF) << 32) |
-		((uint64) (random() & 0xFFFF) << 16) |
-		(uint64) (random() & 0xFFFF);
+		uint64		seed = pg_prng_u64();
 
 		for (i = 0; i < nclients; i++)
 			if (!putVariableInt(&state[i], "startup", "default_seed", (int64) seed))
@@ -5909,9 +5904,9 @@ main(int argc, char **argv)
 		thread->state = &state[nclients_dealt];
 		thread->nstate =
 			(nclients - nclients_dealt + nthreads - i - 1) / (nthreads - i);
-		initRandomState(&thread->ts_choose_rs);
-		initRandomState(&thread->ts_throttle_rs);
-		initRandomState(&thread->ts_sample_rs);
+		pg_prng_state_init(&thread->ts_choose_rs);
+		pg_prng_state_init(&thread->ts_throttle_rs);
+		pg_prng_state_init(&thread->ts_sample_rs);
 		thread->logfile = NULL; /* filled in later */
 		thread->latency_late = 0;
 		thread->zipf_cache.nb_cells = 0;
diff --git a/src/bin/pgbench/t/001_pgbench_with_server.pl b/src/bin/pgbench/t/001_pgbench_with_server.pl
index 17323e6ca5..6e5e2cf997 100644
--- a/src/bin/pgbench/t/001_pgbench_with_server.pl
+++ b/src/bin/pgbench/t/001_pgbench_with_server.pl
@@ -260,10 +260,10 @@ pgbench(
 		qr{setting random seed to 5432\b},
 
 		# After explicit seeding, the four * random checks (1-3,20) should be
-		# deterministic, but not necessarily portable.
-		qr{command=1.: int 1\d\b},        # uniform random: 12 on linux
-		qr{command=2.: int 1\d\d\b},      # exponential random: 106 on linux
-		qr{command=3.: int 1\d\d\d\b},    # gaussian random: 1462 on linux
+		# deterministic on all platforms
+		qr{command=1.: int 18\b},     # uniform random
+		qr{command=2.: int 101\b},    # exponential random
+		qr{command=3.: int 1509\b},   # gaussian random
 		qr{command=4.: int 4\b},
 		qr{command=5.: int 5\b},
 		qr{command=6.: int 6\b},
@@ -276,7 +276,7 @@ pgbench(
 		qr{command=15.: double 15\b},
 		qr{command=16.: double 16\b},
 		qr{command=17.: double 17\b},
-		qr{command=20.: int \d\b},    # zipfian random: 1 on linux
+		qr{command=20.: int 3\b},     # zipfian random
 		qr{command=21.: double -27\b},
 		qr{command=22.: double 1024\b},
 		qr{command=23.: double 1\b},
diff --git a/src/include/port.h b/src/include/port.h
index a55c473262..7597f1a78b 100644
--- a/src/include/port.h
+++ b/src/include/port.h
@@ -348,6 +348,57 @@ extern long pg_lrand48(void);
 extern long pg_jrand48(unsigned short xseed[3]);
 extern void pg_srand48(long seed);
 
+/*----
+ * Postgres internal pseudo-random number generator (PRNG)
+ *
+ * 0: rand48 (pg historical linear congruential generator)
+ * 1: Knuth's MMIX LCG
+ * 2: splitmix64
+ * 3: xoshiro128+
+ * 4: xoshiro256** (default)
+ */
+#ifndef PG_PRNG_ALGORITHM
+#define PG_PRNG_ALGORITHM 0
+#endif
+
+/* set PRNG state byte size depending on algorithm */
+#if PG_PRNG_ALGORITHM == 0
+#define PG_PRNG_STATE_SIZE 6
+#elif PG_PRNG_ALGORITHM == 1 || PG_PRNG_ALGORITHM == 2
+#define PG_PRNG_STATE_SIZE 8
+#elif PG_PRNG_ALGORITHM == 3
+#define PG_PRNG_STATE_SIZE 16
+#elif PG_PRNG_ALGORITHM == 4
+#define PG_PRNG_STATE_SIZE 32
+#else
+#error "unexpected PG_PRNG_ALGORITHM value, must be in 0..4"
+#endif
+
+/* closed type seen as a byte array */
+typedef struct {
+	unsigned char data[PG_PRNG_STATE_SIZE];
+} pg_prng_state;
+
+/* interface with a reference to an external state */
+extern double pg_prng_double_r(pg_prng_state *state);
+extern int32 pg_prng_u31_r(pg_prng_state *state);
+extern int32 pg_prng_s32_r(pg_prng_state *state);
+extern uint32 pg_prng_u32_r(pg_prng_state *state);
+extern int64 pg_prng_s64_r(pg_prng_state *state);
+extern uint64 pg_prng_u64_r(pg_prng_state *state);
+extern void pg_prng_setstate_r(pg_prng_state *state, uint64 u);
+extern bool pg_prng_strong_setstate_r(pg_prng_state *state);
+
+/* interface with a shared internal state */
+extern double pg_prng_double(void);
+extern int32 pg_prng_u31(void);
+extern int32 pg_prng_s32(void);
+extern uint32 pg_prng_u32(void);
+extern int64 pg_prng_s64(void);
+extern uint64 pg_prng_u64(void);
+extern void pg_prng_setstate(uint64 u);
+extern bool pg_prng_strong_setstate(void);
+
 #ifndef HAVE_FLS
 extern int	fls(int mask);
 #endif
diff --git a/src/port/Makefile b/src/port/Makefile
index 9cfc0f9279..afa4161dde 100644
--- a/src/port/Makefile
+++ b/src/port/Makefile
@@ -37,7 +37,7 @@ LIBS += $(PTHREAD_LIBS)
 
 OBJS = $(LIBOBJS) $(PG_CRC32C_OBJS) chklocale.o erand48.o inet_net_ntop.o \
 	noblock.o path.o pgcheckdir.o pgmkdirp.o pgsleep.o \
-	pg_strong_random.o pgstrcasecmp.o pgstrsignal.o pqsignal.o \
+	pg_strong_random.o pgstrcasecmp.o pgstrsignal.o pqsignal.o prng.o \
 	qsort.o qsort_arg.o quotes.o snprintf.o sprompt.o strerror.o \
 	tar.o thread.o
 
diff --git a/src/port/prng.c b/src/port/prng.c
new file mode 100644
index 0000000000..4b8002270e
--- /dev/null
+++ b/src/port/prng.c
@@ -0,0 +1,438 @@
+/*-------------------------------------------------------------------------
+ *
+ * prng.c
+ *
+ * Postgres internal PRNG (pseudo-random number generator).
+ *
+ * This file supplies several implementations of a PRNG: POSIX-standard
+ * rand48, splitmix64, xoroshiro128+ and xoshiro256**, all wrapped into
+ * the same abstract interface so that changing the function only requires
+ * recompiling, without further changes. The variants with a reference to
+ * the state (*_r) are thread-safe if used with a thread-specific state.
+ *
+ * None of these functions should apply to anything which needs cryptographic
+ * security, even if seeded from strong random data: beware that the internal
+ * state can probably be recovered from a few outputs. If the state is too
+ * small or initialized from too small (random) data, they may be subject to a
+ * state enumeration attack. Most of these functions pass statistical tests,
+ * although they fail under some chosen transformation.
+ *
+ * Portions Copyright (c) 1996-2018, PostgreSQL Global Development Group
+ *
+ * Portions Copyright (c) 1993 Martin Birgmeier
+ * All rights reserved.
+ *
+ * You may redistribute unmodified or modified versions of this source
+ * code provided that the above copyright notice and this and the
+ * following conditions are retained.
+ *
+ * This software is provided ``as is'', and comes with no warranties
+ * of any kind. I shall in no event be liable for anything that happens
+ * to anyone/anything when using this software.
+ *
+ * IDENTIFICATION
+ *	  src/port/prng.c
+ *
+ *-------------------------------------------------------------------------
+ */
+
+#include "c.h"
+#include "port.h"
+
+#include <math.h>
+
+/*
+ * Left rotation, decent compilers should generate hardware instructions
+ */
+#define ROTL(x, n) ((x << n) | (x >> (64 - n)))
+
+/*
+ * PG PRNG ALGORITHM
+ */
+const char * pg_prng_algorithm =
+#if PG_PRNG_ALGORITHM == 0
+	/* POSIX standard */
+	"rand48"
+#elif PG_PRNG_ALGORITHM == 1
+	/* Donald Knuth's MMIX LCG */
+	"knuth-mmix-lcg"
+#elif PG_PRNG_ALGORITHM == 2
+	/*----
+	 * Basic 8 bytes state PRNG from Java 8,
+	 * see http://xoshiro.di.unimi.it/splitmix64.c
+	 */
+	"splitmix64"
+#elif PG_PRNG_ALGORITHM == 3
+	/*----
+	 * Algorithm by David Blackman and Sebastiano Vigna,
+	 * see http://xoshiro.di.unimi.it/xoroshiro128plus.c
+	 */
+	"xoroshiro128+"
+#elif PG_PRNG_ALGORITHM == 4
+	/*----
+	 * Algorithm by David Blackman and Sebastiano Vigna,
+	 * see http://xoshiro.di.unimi.it/xoshiro256starstar.c
+	 */
+	"xoshiro256**"
+#endif
+	;
+
+#if PG_PRNG_ALGORITHM >= 2
+/* also used to initialize large states */
+static uint64 _splitmix64(uint64 *state)
+{
+	uint64	v;
+	*state += UINT64CONST(0x9e3779b97f4a7c15);
+	v = *state;
+	v = (v ^ (v >> 30)) * UINT64CONST(0xbf58476d1ce4e5b9);
+	v = (v ^ (v >> 27)) * UINT64CONST(0x94d049bb133111eb);
+	return v ^ (v >> 31);
+}
+#endif
+
+/*
+ * INTERNAL STATE TYPE
+ */
+typedef struct {
+#if PG_PRNG_STATE_SIZE == 6
+	uint16	data[3];
+#else
+	uint64 data[PG_PRNG_STATE_SIZE / 8];
+#endif
+} _pg_prng_state;
+
+/*
+ * INTERNAL STATE DEFAULT
+ *
+ * pi in hexa is 3.243F6A8885 A308D31319 8A2E037073 44A4093822
+ * 299F31D008 2EFA98EC4E 6C89452821 E638D01377...
+ */
+static _pg_prng_state _prng_state = {
+	{
+
+#if PG_PRNG_STATE_SIZE == 6
+
+/* POSIX specifies 0x330e's use in srand48, but the other bits are arbitrary */
+#define RAND48_SEED_0	(0x330E)
+#define RAND48_SEED_1	(0x243F)
+#define RAND48_SEED_2	(0x6A88)
+
+		RAND48_SEED_0,
+		RAND48_SEED_1,
+		RAND48_SEED_2
+
+#elif PG_PRNG_STATE_SIZE == 8
+
+		UINT64CONST(0x243F6A8885A308D3)
+
+#elif PG_PRNG_STATE_SIZE == 16
+
+		UINT64CONST(0x243F6A8885A308D3),
+		UINT64CONST(0x13198A2E03707344)
+
+#else /* PG_PRNG_STATE_SIZE == 32 */
+
+		UINT64CONST(0x243F6A8885A308D3),
+		UINT64CONST(0x13198A2E03707344),
+		UINT64CONST(0xA4093822299F31D0),
+		UINT64CONST(0x082EFA98EC4E6C89)
+
+#endif
+
+	}
+};
+
+/*
+ * NEXT: advance internal state and return 64 bit value
+ */
+static uint64
+_next(_pg_prng_state *state)
+{
+#if PG_PRNG_ALGORITHM == 0
+
+/* These values are specified by POSIX */
+#define RAND48_MUL	UINT64CONST(0x0005deece66d)
+#define RAND48_ADD	UINT64CONST(0x000b)
+
+	/*
+	 * We do the arithmetic in uint64; any type wider than 48 bits would work.
+	 */
+	uint64		in;
+	uint64		out;
+
+	in = (uint64) state->data[2] << 32 |
+		(uint64) state->data[1] << 16 |
+		(uint64) state->data[0];
+
+	/* this is really 64 bits */
+	out = in * RAND48_MUL + RAND48_ADD;
+
+	/* truncated back to 48 bits */
+	state->data[0] = out & 0xFFFF;
+	state->data[1] = (out >> 16) & 0xFFFF;
+	state->data[2] = (out >> 32) & 0xFFFF;
+
+	return out;
+
+#elif PG_PRNG_ALGORITHM == 1
+
+	state->data[0] = UINT64CONST(6364136223846793005) * state->data[0] +
+					 UINT64CONST(1442695040888963407);
+
+	return state->data[0];
+
+#elif PG_PRNG_ALGORITHM == 2
+
+	return _splitmix64(state->data);
+
+#elif PG_PRNG_ALGORITHM == 3
+
+	const uint64 s0 = state->data[0];
+	uint64 s1 = state->data[1];
+	const uint64 next = s0 + s1;
+
+	s1 ^= s0;
+
+	state->data[0] = ROTL(s0, 24) ^ s1 ^ (s1 << 16);
+	state->data[1] = ROTL(s1, 37);
+
+	return next;
+
+#else /* PG_PRNG_ALGORITHM == 4 */
+
+	const uint64 v = state->data[1] * 5;
+	const uint64 next = ROTL(v, 7) * 9;
+	const uint64 t = state->data[1] << 17;
+
+	state->data[2] ^= state->data[0];
+	state->data[3] ^= state->data[1];
+	state->data[1] ^= state->data[2];
+	state->data[0] ^= state->data[3];
+
+	state->data[2] ^= t;
+	state->data[3] = ROTL(state->data[3], 45);
+
+	return next;
+
+#endif
+}
+
+/*----
+ * Initialize the internal state from a 64-bit seed.
+ *
+ * This makes sense for determinism, eg needed under debug, where the user
+ * can supply an integer for this purpose. If this is not needed, consider
+ * using the strong variant.
+ */
+void
+pg_prng_setstate_r(pg_prng_state *state, uint64 u)
+{
+	_pg_prng_state * s = (_pg_prng_state *) state;
+
+	Assert(sizeof(pg_prng_state) == sizeof(_pg_prng_state));
+
+#if PG_PRNG_STATE_SIZE == 6
+
+	/*
+	 * Per POSIX, this uses only 32 bits from "seed" even if "long" is wider.
+	 * Hence, the set of possible seed values is smaller than it could be.
+	 * Better practice is to use caller-supplied state and initialize it with
+	 * random bits obtained from a high-quality source of random bits.
+	 *
+	 * Note: POSIX specifies a function seed48() that allows all 48 bits
+	 * of the internal state to be set, but we don't currently support that.
+	 */
+	s->data[0] = RAND48_SEED_0;
+	s->data[1] = (uint16) u;
+	s->data[2] = (uint16) (u >> 16);
+
+#elif PG_PRNG_STATE_SIZE == 8
+
+	s->data[0] = u;
+
+#elif PG_PRNG_STATE_SIZE == 16
+
+	s->data[0] = u;
+	s->data[1] = _splitmix64(&u);
+
+#else /* PG_PRNG_STATE_SIZE == 32 */
+
+	s->data[0] = u;
+	s->data[1] = _splitmix64(&u);
+	s->data[2] = _splitmix64(&u);
+	s->data[3] = _splitmix64(&u);
+
+#endif
+}
+
+/*------------------------------------------------------------------------------
+ * Common pseudo-random data extraction functions
+ */
+
+/*----
+ * Generate a pseudo-random floating-point value using caller-supplied state.
+ * Values are uniformly distributed over the interval [0.0, 1.0).
+ */
+double
+pg_prng_double_r(pg_prng_state *state)
+{
+	uint64		x = _next((_pg_prng_state *) state);
+
+#if PG_PRNG_STATE_SIZE == 6
+	/* not enough bits to fill the mantissa?
+	 * in fact we have but that would depart from erand48
+	 */
+	return ldexp((double) (x & UINT64CONST(0x0000FFFFFFFFFFFF)), -48);
+#else /* sizes 8, 16 and 32 are ok, prefer higher order bits though */
+	return ldexp((double) (x >> 12), -52);
+#endif
+}
+
+/* Generate a pseudo-random double in [0.0, 1.0) from internal state. */
+double
+pg_prng_double(void)
+{
+	return pg_prng_double_r((pg_prng_state *) &_prng_state);
+}
+
+/*----
+ * Generate a pseudo-random non-negative integral value using caller-supplied
+ * state. Values are uniformly distributed over the interval [0, 2^31).
+ */
+int32
+pg_prng_u31_r(pg_prng_state *state)
+{
+	uint64		x = _next((_pg_prng_state *) state);
+
+#if PG_PRNG_STATE_SIZE == 6
+	/* LCG must avoid low-weight bits */
+	return (x >> 17) & UINT64CONST(0x7FFFFFFF);
+#else /* combine available 64 bits */
+	return (x + (x >> 32)) & UINT64CONST(0x7FFFFFFF);
+#endif
+}
+
+/* Generate pseudo-random non-negative integral from internal state. */
+int32
+pg_prng_u31(void)
+{
+	return pg_prng_u31_r((pg_prng_state *) &_prng_state);
+}
+
+/*---
+ * Generate a pseudo-random signed integral value using caller-supplied state.
+ * Values are uniformly distributed over the interval [-2^31, 2^31).
+ */
+int32
+pg_prng_s32_r(pg_prng_state *state)
+{
+	uint64		x = _next((_pg_prng_state *) state);
+
+#if PG_PRNG_STATE_SIZE == 6
+	/* LCG must avoid low-weight bits */
+	return (int32) ((x >> 16) & UINT64CONST(0xFFFFFFFF));
+#else
+	/* combine available 64 bits */
+	return (int32) ((x + (x >> 32)) & UINT64CONST(0xFFFFFFFF));
+#endif
+}
+
+/* Generate pseudo-random signed integral from internal state. */
+int32
+pg_prng_s32(void)
+{
+	return pg_prng_s32_r((pg_prng_state *) &_prng_state);
+}
+
+/*---
+ * Generate a pseudo-random unsigned integral value using caller-supplied state.
+ * Values are uniformly distributed over the interval [0, 2^32).
+ */
+uint32
+pg_prng_u32_r(pg_prng_state *state)
+{
+	uint64		x = _next((_pg_prng_state *) state);
+
+#if PG_PRNG_STATE_SIZE == 6
+	/* LCG must avoid low-weight bits */
+	return (x >> 16) & UINT64CONST(0xFFFFFFFF);
+#else
+	/* combine available 64 bits */
+	return ((x + (x >> 32)) & UINT64CONST(0xFFFFFFFF));
+#endif
+}
+
+/* Generate pseudo-random signed integral from internal state. */
+uint32
+pg_prng_u32(void)
+{
+	return pg_prng_u32_r((pg_prng_state *) &_prng_state);
+}
+
+/*
+ * Generate a pseudo-random non-negative integral value using caller-supplied
+ * state. Values are uniformly distributed over the interval [-2^63, 2^63).
+ * (note: not really for low-end LCG)
+ */
+int64
+pg_prng_s64_r(pg_prng_state *state)
+{
+#if PG_PRNG_STATE_SIZE == 6 /* rand48 */
+	/* not enough bits, call _next twice and drop some low bits */
+	return (int64) (_next((_pg_prng_state *) state) << 16  ^
+					_next((_pg_prng_state *) state) >> 20);
+#else
+	return (int64) _next((_pg_prng_state *) state);
+#endif
+}
+
+/* Generate a pseudo-random int64 from internal state */
+int64
+pg_prng_s64(void)
+{
+	return pg_prng_s64_r((pg_prng_state *) &_prng_state);
+}
+
+/*
+ * Generate a pseudo-random non-negative integral value using caller-supplied
+ * state. Values are uniformly distributed over the interval [0, 2^64).
+ */
+uint64
+pg_prng_u64_r(pg_prng_state *state)
+{
+#if PG_PRNG_STATE_SIZE == 6 /* rand48 */
+	/* not enough bits, call _next twice and drop some low bits */
+	return _next((_pg_prng_state *) state) << 16  ^
+		   _next((_pg_prng_state *) state) >> 20;
+#else
+	return _next((_pg_prng_state *) state);
+#endif
+}
+
+/* Generate a pseudo-random uint64 from internal state */
+uint64
+pg_prng_u64(void)
+{
+	return pg_prng_u64_r((pg_prng_state *) &_prng_state);
+}
+
+/* seed internal state from 64 bit integer */
+void
+pg_prng_setstate(uint64 u)
+{
+	pg_prng_setstate_r((pg_prng_state *) &_prng_state, u);
+}
+
+/* try to initialize prng state from a strong random source */
+bool
+pg_prng_strong_setstate_r(pg_prng_state *state)
+{
+	return pg_strong_random(&state->data, PG_PRNG_STATE_SIZE);
+}
+
+/* try to initialize the internal prng state from a strong random source */
+bool
+pg_prng_strong_setstate(void)
+{
+	return pg_prng_strong_setstate_r((pg_prng_state *) &_prng_state);
+}