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
 119         do {
 120                 fread(Rand_state, sizeof(Rand_state[0]), 4, fp);
 121         } while ((Rand_state[0] | Rand_state[1] | Rand_state[2] | Rand_state[3]) == 0);
 122
 123         fclose(fp);
 124
 125 #elif defined(WINDOWS)
 126
 127         HCRYPTPROV hProvider;
 128
 129         CryptAcquireContext(&hProvider, NULL, NULL, PROV_RSA_FULL, 0);
 130
 131         do {
 132                 CryptGenRandom(hProvider, sizeof(Rand_state[0]) * 4, (BYTE*)Rand_state);
 133         } while ((Rand_state[0] | Rand_state[1] | Rand_state[2] | Rand_state[3]) == 0);
 134
 135         CryptReleaseContext(hProvider, 0);
 136
 137 #else
 138
 139         /* Basic seed */
 140         u32b seed = (time(NULL));
 141 #ifdef SET_UID
 142         /* Mutate the seed on Unix machines */
 143         seed = ((seed >> 3) * (getpid() << 1));
 144 #endif
 145         /* Seed the RNG */
 146         Rand_state_set(seed);
 147
 148 #endif
 149 }
 150
 151 /*
 152  * Backup the RNG state
 153  */
 154 void Rand_state_backup(u32b* backup_state)
 155 {
 156         int i;
 157
 158         for (i = 0; i < 4; ++ i) {
 159                 backup_state[i] = Rand_state[i];
 160         }
 161 }
 162
 163 /*
 164  * Restore the RNG state
 165  */
 166 void Rand_state_restore(u32b* backup_state)
 167 {
 168         int i;
 169
 170         for (i = 0; i < 4; ++ i) {
 171                 Rand_state[i] = backup_state[i];
 172         }
 173 }
 174
 175
 176 /*
 177  * Extract a "random" number from 0 to m-1, via "division"
 178  */
 179 static s32b Rand_div_impl(s32b m, u32b* state)
 180 {
 181         u32b scaling;
 182         u32b past;
 183         u32b ret;
 184
 185         /* Hack -- simple case */
 186         if (m <= 1) return (0);
 187
 188         scaling = Rand_Xorshift_max / m;
 189         past = scaling * m;
 190
 191         do {
 192                 ret = Rand_Xorshift(state);
 193         } while (ret >= past);
 194
 195         return ret / scaling;
 196 }
 197
 198 s32b Rand_div(s32b m)
 199 {
 200         return Rand_div_impl(m, Rand_state);
 201 }
 202
 203
 204
 205
 206 /*
 207  * The number of entries in the "randnor_table"
 208  */
 209 #define RANDNOR_NUM     256
 210
 211 /*
 212  * The standard deviation of the "randnor_table"
 213  */
 214 #define RANDNOR_STD     64
 215
 216 /*
 217  * The normal distribution table for the "randnor()" function (below)
 218  */
 219 static s16b randnor_table[RANDNOR_NUM] =
 220 {
 221         206,     613,    1022,    1430,         1838,    2245,    2652,    3058,
 222         3463,    3867,    4271,    4673,        5075,    5475,    5874,    6271,
 223         6667,    7061,    7454,    7845,        8234,    8621,    9006,    9389,
 224         9770,   10148,   10524,   10898,   11269,       11638,   12004,   12367,
 225         12727,   13085,   13440,   13792,   14140,      14486,   14828,   15168,
 226         15504,   15836,   16166,   16492,   16814,      17133,   17449,   17761,
 227         18069,   18374,   18675,   18972,   19266,      19556,   19842,   20124,
 228         20403,   20678,   20949,   21216,   21479,      21738,   21994,   22245,
 229
 230         22493,   22737,   22977,   23213,   23446,      23674,   23899,   24120,
 231         24336,   24550,   24759,   24965,   25166,      25365,   25559,   25750,
 232         25937,   26120,   26300,   26476,   26649,      26818,   26983,   27146,
 233         27304,   27460,   27612,   27760,   27906,      28048,   28187,   28323,
 234         28455,   28585,   28711,   28835,   28955,      29073,   29188,   29299,
 235         29409,   29515,   29619,   29720,   29818,      29914,   30007,   30098,
 236         30186,   30272,   30356,   30437,   30516,      30593,   30668,   30740,
 237         30810,   30879,   30945,   31010,   31072,      31133,   31192,   31249,
 238
 239         31304,   31358,   31410,   31460,   31509,      31556,   31601,   31646,
 240         31688,   31730,   31770,   31808,   31846,      31882,   31917,   31950,
 241         31983,   32014,   32044,   32074,   32102,      32129,   32155,   32180,
 242         32205,   32228,   32251,   32273,   32294,      32314,   32333,   32352,
 243         32370,   32387,   32404,   32420,   32435,      32450,   32464,   32477,
 244         32490,   32503,   32515,   32526,   32537,      32548,   32558,   32568,
 245         32577,   32586,   32595,   32603,   32611,      32618,   32625,   32632,
 246         32639,   32645,   32651,   32657,   32662,      32667,   32672,   32677,
 247
 248         32682,   32686,   32690,   32694,   32698,      32702,   32705,   32708,
 249         32711,   32714,   32717,   32720,   32722,      32725,   32727,   32729,
 250         32731,   32733,   32735,   32737,   32739,      32740,   32742,   32743,
 251         32745,   32746,   32747,   32748,   32749,      32750,   32751,   32752,
 252         32753,   32754,   32755,   32756,   32757,      32757,   32758,   32758,
 253         32759,   32760,   32760,   32761,   32761,      32761,   32762,   32762,
 254         32763,   32763,   32763,   32764,   32764,      32764,   32764,   32765,
 255         32765,   32765,   32765,   32766,   32766,      32766,   32766,   32767,
 256 };
 257
 258
 259
 260 /*
 261  * Generate a random integer number of NORMAL distribution
 262  *
 263  * The table above is used to generate a pseudo-normal distribution,
 264  * in a manner which is much faster than calling a transcendental
 265  * function to calculate a true normal distribution.
 266  *
 267  * Basically, entry 64*N in the table above represents the number of
 268  * times out of 32767 that a random variable with normal distribution
 269  * will fall within N standard deviations of the mean.  That is, about
 270  * 68 percent of the time for N=1 and 95 percent of the time for N=2.
 271  *
 272  * The table above contains a "faked" final entry which allows us to
 273  * pretend that all values in a normal distribution are strictly less
 274  * than four standard deviations away from the mean.  This results in
 275  * "conservative" distribution of approximately 1/32768 values.
 276  *
 277  * Note that the binary search takes up to 16 quick iterations.
 278  */
 279 s16b randnor(int mean, int stand)
 280 {
 281         s16b tmp;
 282         s16b offset;
 283
 284         s16b low = 0;
 285         s16b high = RANDNOR_NUM;
 286
 287         /* Paranoia */
 288         if (stand < 1) return (mean);
 289
 290         /* Roll for probability */
 291         tmp = (s16b)randint0(32768);
 292
 293         /* Binary Search */
 294         while (low < high)
 295         {
 296                 int mid = (low + high) >> 1;
 297
 298                 /* Move right if forced */
 299                 if (randnor_table[mid] < tmp)
 300                 {
 301                         low = mid + 1;
 302                 }
 303
 304                 /* Move left otherwise */
 305                 else
 306                 {
 307                         high = mid;
 308                 }
 309         }
 310
 311         /* Convert the index into an offset */
 312         offset = (long)stand * (long)low / RANDNOR_STD;
 313
 314         /* One half should be negative */
 315         if (randint0(100) < 50) return (mean - offset);
 316
 317         /* One half should be positive */
 318         return (mean + offset);
 319 }
 320
 321
 322
 323 /*
 324  * Generates damage for "2d6" style dice rolls
 325  */
 326 s16b damroll(int num, int sides)
 327 {
 328         int i, sum = 0;
 329         for (i = 0; i < num; i++) sum += randint1(sides);
 330         return (sum);
 331 }
 332
 333
 334 /*
 335  * Same as above, but always maximal
 336  */
 337 s16b maxroll(int num, int sides)
 338 {
 339         return (num * sides);
 340 }
 341
 342
 343 /*
 344  * Given a numerator and a denominator, supply a properly rounded result,
 345  * using the RNG to smooth out remainders.  -LM-
 346  */
 347 s32b div_round(s32b n, s32b d)
 348 {
 349         s32b tmp;
 350
 351         /* Refuse to divide by zero */
 352         if (!d) return (n);
 353
 354         /* Division */
 355         tmp = n / d;
 356
 357         /* Rounding */
 358         if ((ABS(n) % ABS(d)) > randint0(ABS(d)))
 359         {
 360                 /* Increase the absolute value */
 361                 if (n * d > 0L) tmp += 1L;
 362                 else            tmp -= 1L;
 363         }
 364
 365         /* Return */
 366         return (tmp);
 367 }
 368
 369
 370
 371
 372 /*
 373  * Extract a "random" number from 0 to m-1, using the RNG.
 374  *
 375  * This function should be used when generating random numbers in
 376  * "external" program parts like the main-*.c files.  It preserves
 377  * the current RNG state to prevent influences on game-play.
 378  *
 379  * Could also use rand() from <stdlib.h> directly. XXX XXX XXX
 380  */
 381 s32b Rand_external(s32b m)
 382 {
 383         static bool initialized = FALSE;
 384         static u32b Rand_state_external[4];
 385
 386         if (!initialized)
 387         {
 388                 /* Initialize with new seed */
 389                 u32b seed = time(NULL);
 390                 Rand_Xorshift_seed(seed, Rand_state_external);
 391                 initialized = TRUE;
 392         }
 393
 394         return Rand_div_impl(m, Rand_state_external);
 395 }