1.1 --- /dev/null	Thu Jan 01 00:00:00 1970 +0000
     1.2 +++ b/src/rng.cc	Sun Nov 01 00:09:12 2015 +0200
     1.3 @@ -0,0 +1,217 @@
     1.4 +#include <stdlib.h>
     1.5 +#include <string.h>
     1.6 +#include "rng.h"
     1.7 +
     1.8 +//#define USE_STD_RAND
     1.9 +
    1.10 +// defined at the end
    1.11 +/**
    1.12 + * tinymt32 internal state vector and parameters
    1.13 + */
    1.14 +struct tinymt32_t {
    1.15 +	uint32_t status[4];
    1.16 +	uint32_t mat1;
    1.17 +	uint32_t mat2;
    1.18 +	uint32_t tmat;
    1.19 +};
    1.20 +
    1.21 +#ifndef USE_STD_RAND
    1.22 +static void tinymt32_init(tinymt32_t * random, uint32_t seed);
    1.23 +inline static uint32_t tinymt32_generate_uint32(tinymt32_t * random);
    1.24 +inline static float tinymt32_generate_float(tinymt32_t * random);
    1.25 +#endif
    1.26 +
    1.27 +static RandGen defrng;
    1.28 +
    1.29 +struct RandGenState {
    1.30 +	tinymt32_t st;
    1.31 +};
    1.32 +
    1.33 +RandGen::RandGen()
    1.34 +{
    1.35 +	state = new RandGenState;
    1.36 +	memset(state, 0, sizeof *state);
    1.37 +#ifndef USE_STD_RAND
    1.38 +	tinymt32_init(&state->st, 0);
    1.39 +#endif
    1.40 +}
    1.41 +
    1.42 +RandGen::~RandGen()
    1.43 +{
    1.44 +	delete state;
    1.45 +}
    1.46 +
    1.47 +void RandGen::seed(uint32_t s)
    1.48 +{
    1.49 +#ifndef USE_STD_RAND
    1.50 +	tinymt32_init(&state->st, s);
    1.51 +#else
    1.52 +	srand(s);
    1.53 +#endif
    1.54 +}
    1.55 +
    1.56 +uint32_t RandGen::generate()
    1.57 +{
    1.58 +#ifndef USE_STD_RAND
    1.59 +	return tinymt32_generate_uint32(&state->st);
    1.60 +#else
    1.61 +	return rand();
    1.62 +#endif
    1.63 +}
    1.64 +
    1.65 +float RandGen::generate_float()
    1.66 +{
    1.67 +#ifndef USE_STD_RAND
    1.68 +	return tinymt32_generate_float(&state->st);
    1.69 +#else
    1.70 +	return (float)rand() / (float)RAND_MAX;
    1.71 +#endif
    1.72 +}
    1.73 +
    1.74 +void rng_srand(uint32_t s)
    1.75 +{
    1.76 +	defrng.seed(s);
    1.77 +}
    1.78 +
    1.79 +uint32_t rng_rand()
    1.80 +{
    1.81 +	return defrng.generate();
    1.82 +}
    1.83 +
    1.84 +float rng_frand()
    1.85 +{
    1.86 +	return defrng.generate_float();
    1.87 +}
    1.88 +
    1.89 +
    1.90 +#ifndef USE_STD_RAND
    1.91 +/**
    1.92 + * @file tinymt32.h
    1.93 + *
    1.94 + * @brief Tiny Mersenne Twister only 127 bit internal state
    1.95 + *
    1.96 + * @author Mutsuo Saito (Hiroshima University)
    1.97 + * @author Makoto Matsumoto (University of Tokyo)
    1.98 + *
    1.99 + * Copyright (C) 2011 Mutsuo Saito, Makoto Matsumoto,
   1.100 + * Hiroshima University and The University of Tokyo.
   1.101 + * All rights reserved.
   1.102 + *
   1.103 + * The 3-clause BSD License is applied to this software, see
   1.104 + * LICENSE.txt
   1.105 + */
   1.106 +#define TINYMT32_MEXP 127
   1.107 +#define TINYMT32_SH0 1
   1.108 +#define TINYMT32_SH1 10
   1.109 +#define TINYMT32_SH8 8
   1.110 +#define TINYMT32_MASK UINT32_C(0x7fffffff)
   1.111 +#define TINYMT32_MUL (1.0f / 4294967296.0f)
   1.112 +
   1.113 +/**
   1.114 + * This function changes internal state of tinymt32.
   1.115 + * Users should not call this function directly.
   1.116 + * @param random tinymt internal status
   1.117 + */
   1.118 +inline static void tinymt32_next_state(tinymt32_t * random) {
   1.119 +	uint32_t x;
   1.120 +	uint32_t y;
   1.121 +
   1.122 +	y = random->status[3];
   1.123 +	x = (random->status[0] & TINYMT32_MASK)
   1.124 +		^ random->status[1]
   1.125 +		^ random->status[2];
   1.126 +	x ^= (x << TINYMT32_SH0);
   1.127 +	y ^= (y >> TINYMT32_SH0) ^ x;
   1.128 +	random->status[0] = random->status[1];
   1.129 +	random->status[1] = random->status[2];
   1.130 +	random->status[2] = x ^ (y << TINYMT32_SH1);
   1.131 +	random->status[3] = y;
   1.132 +	random->status[1] ^= -((int32_t)(y & 1)) & random->mat1;
   1.133 +	random->status[2] ^= -((int32_t)(y & 1)) & random->mat2;
   1.134 +}
   1.135 +
   1.136 +/**
   1.137 + * This function outputs 32-bit unsigned integer from internal state.
   1.138 + * Users should not call this function directly.
   1.139 + * @param random tinymt internal status
   1.140 + * @return 32-bit unsigned pseudorandom number
   1.141 + */
   1.142 +inline static uint32_t tinymt32_temper(tinymt32_t * random) {
   1.143 +	uint32_t t0, t1;
   1.144 +	t0 = random->status[3];
   1.145 +#if defined(LINEARITY_CHECK)
   1.146 +	t1 = random->status[0]
   1.147 +		^ (random->status[2] >> TINYMT32_SH8);
   1.148 +#else
   1.149 +	t1 = random->status[0]
   1.150 +		+ (random->status[2] >> TINYMT32_SH8);
   1.151 +#endif
   1.152 +	t0 ^= t1;
   1.153 +	t0 ^= -((int32_t)(t1 & 1)) & random->tmat;
   1.154 +	return t0;
   1.155 +}
   1.156 +
   1.157 +/**
   1.158 + * This function outputs 32-bit unsigned integer from internal state.
   1.159 + * @param random tinymt internal status
   1.160 + * @return 32-bit unsigned integer r (0 <= r < 2^32)
   1.161 + */
   1.162 +inline static uint32_t tinymt32_generate_uint32(tinymt32_t * random) {
   1.163 +	tinymt32_next_state(random);
   1.164 +	return tinymt32_temper(random);
   1.165 +}
   1.166 +
   1.167 +/**
   1.168 + * This function outputs floating point number from internal state.
   1.169 + * This function is implemented using multiplying by 1 / 2^32.
   1.170 + * floating point multiplication is faster than using union trick in
   1.171 + * my Intel CPU.
   1.172 + * @param random tinymt internal status
   1.173 + * @return floating point number r (0.0 <= r < 1.0)
   1.174 + */
   1.175 +inline static float tinymt32_generate_float(tinymt32_t * random) {
   1.176 +	tinymt32_next_state(random);
   1.177 +	return tinymt32_temper(random) * TINYMT32_MUL;
   1.178 +}
   1.179 +
   1.180 +#define MIN_LOOP 8
   1.181 +#define PRE_LOOP 8
   1.182 +
   1.183 +/**
   1.184 + * This function certificate the period of 2^127-1.
   1.185 + * @param random tinymt state vector.
   1.186 + */
   1.187 +static void period_certification(tinymt32_t * random) {
   1.188 +	if ((random->status[0] & TINYMT32_MASK) == 0 &&
   1.189 +			random->status[1] == 0 &&
   1.190 +			random->status[2] == 0 &&
   1.191 +			random->status[3] == 0) {
   1.192 +		random->status[0] = 'T';
   1.193 +		random->status[1] = 'I';
   1.194 +		random->status[2] = 'N';
   1.195 +		random->status[3] = 'Y';
   1.196 +	}
   1.197 +}
   1.198 +
   1.199 +/**
   1.200 + * This function initializes the internal state array with a 32-bit
   1.201 + * unsigned integer seed.
   1.202 + * @param random tinymt state vector.
   1.203 + * @param seed a 32-bit unsigned integer used as a seed.
   1.204 + */
   1.205 +static void tinymt32_init(tinymt32_t * random, uint32_t seed) {
   1.206 +	random->status[0] = seed;
   1.207 +	random->status[1] = random->mat1;
   1.208 +	random->status[2] = random->mat2;
   1.209 +	random->status[3] = random->tmat;
   1.210 +	for (int i = 1; i < MIN_LOOP; i++) {
   1.211 +		random->status[i & 3] ^= i + UINT32_C(1812433253)
   1.212 +			* (random->status[(i - 1) & 3]
   1.213 +					^ (random->status[(i - 1) & 3] >> 30));
   1.214 +	}
   1.215 +	period_certification(random);
   1.216 +	for (int i = 0; i < PRE_LOOP; i++) {
   1.217 +		tinymt32_next_state(random);
   1.218 +	}
   1.219 +}
   1.220 +#endif	// !defined(USE_STD_RAND)