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