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