src/z-rand.c

   1 /* File: z-rand.c */
   2
   3 /*
   4  * Copyright (c) 1997 Ben Harrison, and others
   5  *
   6  * This software may be copied and distributed for educational, research,
   7  * and not for profit purposes provided that this copyright and statement
   8  * are included in all such copies.  Other copyrights may also apply.
   9  */
  10
  11
  12 /* Purpose: a simple random number generator -BEN- */
  13
  14 #if defined(WINDOWS)
  15 #include <Windows.h>
  16 #endif
  17
  18 #include "z-rand.h"
  19
  20
  21
  22
  23 /*
  24  * Angband 2.7.9 introduced a new (optimized) random number generator,
  25  * based loosely on the old "random.c" from Berkeley but with some major
  26  * optimizations and algorithm changes.  See below for more details.
  27  *
  28  * Code by myself (benh@phial.com) and Randy (randy@stat.tamu.edu).
  29  *
  30  * This code provides (1) a "decent" RNG, based on the "BSD-degree-63-RNG"
  31  * used in Angband 2.7.8, but rather optimized, and (2) a "simple" RNG,
  32  * based on the simple "LCRNG" currently used in Angband, but "corrected"
  33  * to give slightly better values.  Both of these are available in two
  34  * flavors, first, the simple "mod" flavor, which is fast, but slightly
  35  * biased at high values, and second, the simple "div" flavor, which is
  36  * less fast (and potentially non-terminating) but which is not biased
  37  * and is much less subject to low-bit-non-randomness problems.
  38  *
  39  * You can select your favorite flavor by proper definition of the
  40  * "randint0()" macro in the "defines.h" file.
  41  *
  42  * Note that, in Angband 2.8.0, the "state" table will be saved in the
  43  * savefile, so a special "initialization" phase will be necessary.
  44  *
  45  * Note the use of the "simple" RNG, first you activate it via
  46  * "Rand_quick = TRUE" and "Rand_value = seed" and then it is used
  47  * automatically used instead of the "complex" RNG, and when you are
  48  * done, you de-activate it via "Rand_quick = FALSE" or choose a new
  49  * seed via "Rand_value = seed".
  50  *
  51  *
  52  * RNG algorithm was fully rewritten. Upper comment is OLD.
  53  */
  54
  55
  56 /*
  57  * Currently unused
  58  */
  59 u16b Rand_place;
  60
  61 /*
  62  * Current "state" table for the RNG
  63  * Only index 0 to 3 are used
  64  */
  65 u32b Rand_state[RAND_DEG] = {
  66         123456789,
  67         362436069,
  68         521288629,
  69         88675123,
  70 };
  71
  72
  73 /*
  74  * Initialize Xorshift Algorithm state
  75  */
  76 static void Rand_Xorshift_seed(u32b seed, u32b* state)
  77 {
  78         int i;
  79
  80         /* Initialize Xorshift Algorithm RNG */
  81         for (i = 1; i <= 4; ++ i) {
  82                 seed = 1812433253UL * (seed ^ (seed >> 30)) + i;
  83                 state[i-1] = seed;
  84         }
  85 }
  86
  87 /*
  88  * Xorshift Algorithm
  89  */
  90 static u32b Rand_Xorshift(u32b* state)
  91 {
  92         u32b t = state[0] ^ (state[0] << 11);
  93
  94         state[0] = state[1];
  95         state[1] = state[2];
  96         state[2] = state[3];
  97
  98         state[3] = (state[3] ^ (state[3] >> 19)) ^ (t ^ (t >> 8));
  99
 100         return state[3];
 101 }
 102
 103 static const u32b Rand_Xorshift_max = 0xFFFFFFFF;
 104
 105 /*
 106  * Initialize the RNG using a new seed
 107  */
 108 void Rand_state_set(u32b seed)
 109 {
 110         Rand_Xorshift_seed(seed, Rand_state);
 111 }
 112
 113 void Rand_state_init(void)
 114 {
 115 #ifdef RNG_DEVICE
 116
 117         FILE *fp = fopen(RNG_DEVICE, "r");
 118         u32b buf[4];
 119         do {
 120                 fread(buf, sizeof(buf[0]), 4, fp);
 121         } while ((buf[0] | buf[1] | buf[2] | buf[3]) == 0);
 122         memcpy(Rand_state, buf, sizeof(buf));
 123         fclose(fp);
 124
 125 #elif defined(WINDOWS)
 126
 127         HCRYPTPROV hProvider;
 128         u32b buf[4];
 129
 130         CryptAcquireContext(&hProvider, NULL, NULL, PROV_RSA_FULL, 0);
 131
 132         do {
 133                 CryptGenRandom(hProvider, sizeof(buf[0]) * 4, (BYTE*)buf);
 134         } while ((buf[0] | buf[1] | buf[2] | buf[3]) == 0);
 135
 136         CryptReleaseContext(hProvider, 0);
 137
 138 #else
 139
 140         /* Basic seed */
 141         u32b seed = (time(NULL));
 142 #ifdef SET_UID
 143         /* Mutate the seed on Unix machines */
 144         seed = ((seed >> 3) * (getpid() << 1));
 145 #endif
 146         /* Seed the RNG */
 147         Rand_state_set(seed);
 148
 149 #endif
 150 }
 151
 152 /*
 153  * Backup the RNG state
 154  */
 155 void Rand_state_backup(u32b* backup_state)
 156 {
 157         int i;
 158
 159         for (i = 0; i < 4; ++ i) {
 160                 backup_state[i] = Rand_state[i];
 161         }
 162 }
 163
 164 /*
 165  * Restore the RNG state
 166  */
 167 void Rand_state_restore(u32b* backup_state)
 168 {
 169         int i;
 170
 171         for (i = 0; i < 4; ++ i) {
 172                 Rand_state[i] = backup_state[i];
 173         }
 174 }
 175
 176
 177 /*
 178  * Extract a "random" number from 0 to m-1, via "division"
 179  */
 180 static s32b Rand_div_impl(s32b m, u32b* state)
 181 {
 182         u32b scaling;
 183         u32b past;
 184         u32b ret;
 185
 186         /* Hack -- simple case */
 187         if (m <= 1) return (0);
 188
 189         scaling = Rand_Xorshift_max / m;
 190         past = scaling * m;
 191
 192         do {
 193                 ret = Rand_Xorshift(state);
 194         } while (ret >= past);
 195
 196         return ret / scaling;
 197 }
 198
 199 s32b Rand_div(s32b m)
 200 {
 201         return Rand_div_impl(m, Rand_state);
 202 }
 203
 204
 205
 206
 207 /*
 208  * The number of entries in the "randnor_table"
 209  */
 210 #define RANDNOR_NUM     256
 211
 212 /*
 213  * The standard deviation of the "randnor_table"
 214  */
 215 #define RANDNOR_STD     64
 216
 217 /*
 218  * The normal distribution table for the "randnor()" function (below)
 219  */
 220 static s16b randnor_table[RANDNOR_NUM] =
 221 {
 222         206,     613,    1022,    1430,         1838,    2245,    2652,    3058,
 223         3463,    3867,    4271,    4673,        5075,    5475,    5874,    6271,
 224         6667,    7061,    7454,    7845,        8234,    8621,    9006,    9389,
 225         9770,   10148,   10524,   10898,   11269,       11638,   12004,   12367,
 226         12727,   13085,   13440,   13792,   14140,      14486,   14828,   15168,
 227         15504,   15836,   16166,   16492,   16814,      17133,   17449,   17761,
 228         18069,   18374,   18675,   18972,   19266,      19556,   19842,   20124,
 229         20403,   20678,   20949,   21216,   21479,      21738,   21994,   22245,
 230
 231         22493,   22737,   22977,   23213,   23446,      23674,   23899,   24120,
 232         24336,   24550,   24759,   24965,   25166,      25365,   25559,   25750,
 233         25937,   26120,   26300,   26476,   26649,      26818,   26983,   27146,
 234         27304,   27460,   27612,   27760,   27906,      28048,   28187,   28323,
 235         28455,   28585,   28711,   28835,   28955,      29073,   29188,   29299,
 236         29409,   29515,   29619,   29720,   29818,      29914,   30007,   30098,
 237         30186,   30272,   30356,   30437,   30516,      30593,   30668,   30740,
 238         30810,   30879,   30945,   31010,   31072,      31133,   31192,   31249,
 239
 240         31304,   31358,   31410,   31460,   31509,      31556,   31601,   31646,
 241         31688,   31730,   31770,   31808,   31846,      31882,   31917,   31950,
 242         31983,   32014,   32044,   32074,   32102,      32129,   32155,   32180,
 243         32205,   32228,   32251,   32273,   32294,      32314,   32333,   32352,
 244         32370,   32387,   32404,   32420,   32435,      32450,   32464,   32477,
 245         32490,   32503,   32515,   32526,   32537,      32548,   32558,   32568,
 246         32577,   32586,   32595,   32603,   32611,      32618,   32625,   32632,
 247         32639,   32645,   32651,   32657,   32662,      32667,   32672,   32677,
 248
 249         32682,   32686,   32690,   32694,   32698,      32702,   32705,   32708,
 250         32711,   32714,   32717,   32720,   32722,      32725,   32727,   32729,
 251         32731,   32733,   32735,   32737,   32739,      32740,   32742,   32743,
 252         32745,   32746,   32747,   32748,   32749,      32750,   32751,   32752,
 253         32753,   32754,   32755,   32756,   32757,      32757,   32758,   32758,
 254         32759,   32760,   32760,   32761,   32761,      32761,   32762,   32762,
 255         32763,   32763,   32763,   32764,   32764,      32764,   32764,   32765,
 256         32765,   32765,   32765,   32766,   32766,      32766,   32766,   32767,
 257 };
 258
 259
 260
 261 /*
 262  * Generate a random integer number of NORMAL distribution
 263  *
 264  * The table above is used to generate a pseudo-normal distribution,
 265  * in a manner which is much faster than calling a transcendental
 266  * function to calculate a true normal distribution.
 267  *
 268  * Basically, entry 64*N in the table above represents the number of
 269  * times out of 32767 that a random variable with normal distribution
 270  * will fall within N standard deviations of the mean.  That is, about
 271  * 68 percent of the time for N=1 and 95 percent of the time for N=2.
 272  *
 273  * The table above contains a "faked" final entry which allows us to
 274  * pretend that all values in a normal distribution are strictly less
 275  * than four standard deviations away from the mean.  This results in
 276  * "conservative" distribution of approximately 1/32768 values.
 277  *
 278  * Note that the binary search takes up to 16 quick iterations.
 279  */
 280 s16b randnor(int mean, int stand)
 281 {
 282         s16b tmp;
 283         s16b offset;
 284
 285         s16b low = 0;
 286         s16b high = RANDNOR_NUM;
 287
 288         /* Paranoia */
 289         if (stand < 1) return (mean);
 290
 291         /* Roll for probability */
 292         tmp = (s16b)randint0(32768);
 293
 294         /* Binary Search */
 295         while (low < high)
 296         {
 297                 int mid = (low + high) >> 1;
 298
 299                 /* Move right if forced */
 300                 if (randnor_table[mid] < tmp)
 301                 {
 302                         low = mid + 1;
 303                 }
 304
 305                 /* Move left otherwise */
 306                 else
 307                 {
 308                         high = mid;
 309                 }
 310         }
 311
 312         /* Convert the index into an offset */
 313         offset = (long)stand * (long)low / RANDNOR_STD;
 314
 315         /* One half should be negative */
 316         if (randint0(100) < 50) return (mean - offset);
 317
 318         /* One half should be positive */
 319         return (mean + offset);
 320 }
 321
 322
 323
 324 /*
 325  * Generates damage for "2d6" style dice rolls
 326  */
 327 s16b damroll(int num, int sides)
 328 {
 329         int i, sum = 0;
 330         for (i = 0; i < num; i++) sum += randint1(sides);
 331         return (sum);
 332 }
 333
 334
 335 /*
 336  * Same as above, but always maximal
 337  */
 338 s16b maxroll(int num, int sides)
 339 {
 340         return (num * sides);
 341 }
 342
 343
 344 /*
 345  * Given a numerator and a denominator, supply a properly rounded result,
 346  * using the RNG to smooth out remainders.  -LM-
 347  */
 348 s32b div_round(s32b n, s32b d)
 349 {
 350         s32b tmp;
 351
 352         /* Refuse to divide by zero */
 353         if (!d) return (n);
 354
 355         /* Division */
 356         tmp = n / d;
 357
 358         /* Rounding */
 359         if ((ABS(n) % ABS(d)) > randint0(ABS(d)))
 360         {
 361                 /* Increase the absolute value */
 362                 if (n * d > 0L) tmp += 1L;
 363                 else            tmp -= 1L;
 364         }
 365
 366         /* Return */
 367         return (tmp);
 368 }
 369
 370
 371
 372
 373 /*
 374  * Extract a "random" number from 0 to m-1, using the RNG.
 375  *
 376  * This function should be used when generating random numbers in
 377  * "external" program parts like the main-*.c files.  It preserves
 378  * the current RNG state to prevent influences on game-play.
 379  *
 380  * Could also use rand() from <stdlib.h> directly. XXX XXX XXX
 381  */
 382 s32b Rand_external(s32b m)
 383 {
 384         static bool initialized = FALSE;
 385         static u32b Rand_state_external[4];
 386
 387         if (!initialized)
 388         {
 389                 /* Initialize with new seed */
 390                 u32b seed = time(NULL);
 391                 Rand_Xorshift_seed(seed, Rand_state_external);
 392                 initialized = TRUE;
 393         }
 394
 395         return Rand_div_impl(m, Rand_state_external);
 396 }