ref: 3206bcd5977bc732622e5a5f408d64f23ebee097
parent: 64307b360c9b55982e9dca4984f8ca9e1007e915
author: Sigrid Solveig Haflínudóttir <[email protected]>
date: Fri Jan 26 10:16:59 EST 2024
PF_random: use mersenne twister
--- /dev/null
+++ b/m_random.c
@@ -1,0 +1,40 @@
+#include "quakedef.h"
+#include "mt19937-64.h"
+
+static mt19937_64 ctx;
+
+void
+m_random_init(u64int seed)
+{+ init_genrand64(&ctx, seed);
+}
+
+u64int
+m_random_64(void)
+{+ return genrand64_int64(&ctx);
+}
+
+s64int
+m_random_63(void)
+{+ return genrand64_int63(&ctx);
+}
+
+double
+m_random_i0i1(void)
+{+ return genrand64_real1(&ctx);
+}
+
+double
+m_random_i0e1(void)
+{+ return genrand64_real2(&ctx);
+}
+
+double
+m_random_e0e1(void)
+{+ return genrand64_real3(&ctx);
+}
--- /dev/null
+++ b/m_random.h
@@ -1,0 +1,6 @@
+void m_random_init(u64int seed);
+u64int m_random_64(void);
+s64int m_random_63(void);
+double m_random_i0i1(void);
+double m_random_i0e1(void);
+double m_random_e0e1(void);
--- a/meson.build
+++ b/meson.build
@@ -51,6 +51,7 @@
'i_tga.c',
'i_wad.c',
'keys.c',
+ 'm_random.c',
'mathlib.c',
'menu.c',
'model.c',
@@ -60,6 +61,7 @@
'model_bsp2.c',
'model_bsp30.c',
'model_sprite.c',
+ 'mt19937-64.c',
'net_loop.c',
'net_main.c',
'pal.c',
--- a/mkfile
+++ b/mkfile
@@ -39,7 +39,8 @@
in_plan9.$O\
isnanf.$O\
keys.$O\
- m_dotproduct`{test -f span_$objtype.s && echo -n _$objtype}.$O\+ m_dotproduct`{test -f m_dotproduct_$objtype.s && echo -n _$objtype}.$O\+ m_random.$O\
mathlib.$O\
menu.$O\
model.$O\
@@ -49,6 +50,7 @@
model_bsp2.$O\
model_bsp30.$O\
model_sprite.$O\
+ mt19937-64.$O\
nanosec.$O\
net_dgrm_plan9.$O\
net_loop.$O\
--- /dev/null
+++ b/mt19937-64.c
@@ -1,0 +1,182 @@
+/*
+ A C-program for MT19937-64 (2004/9/29 version).
+ Coded by Takuji Nishimura and Makoto Matsumoto.
+
+ This is a 64-bit version of Mersenne Twister pseudorandom number
+ generator.
+
+ Before using, initialize the state by using init_genrand64(seed)
+ or init_by_array64(init_key, key_length).
+
+ Copyright (C) 2004, Makoto Matsumoto and Takuji Nishimura,
+ All rights reserved.
+
+ Redistribution and use in source and binary forms, with or without
+ modification, are permitted provided that the following conditions
+ are met:
+
+ 1. Redistributions of source code must retain the above copyright
+ notice, this list of conditions and the following disclaimer.
+
+ 2. Redistributions in binary form must reproduce the above copyright
+ notice, this list of conditions and the following disclaimer in the
+ documentation and/or other materials provided with the distribution.
+
+ 3. The names of its contributors may not be used to endorse or promote
+ products derived from this software without specific prior written
+ permission.
+
+ THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
+ "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
+ LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
+ A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR
+ CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
+ EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
+ PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
+ PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
+ LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
+ NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
+ SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+
+ References:
+ T. Nishimura, ``Tables of 64-bit Mersenne Twisters''
+ ACM Transactions on Modeling and
+ Computer Simulation 10. (2000) 348--357.
+ M. Matsumoto and T. Nishimura,
+ ``Mersenne Twister: a 623-dimensionally equidistributed
+ uniform pseudorandom number generator''
+ ACM Transactions on Modeling and
+ Computer Simulation 8. (Jan. 1998) 3--30.
+
+ Any feedback is very welcome.
+ http://www.math.hiroshima-u.ac.jp/~m-mat/MT/emt.html
+ email: m-mat @ math.sci.hiroshima-u.ac.jp (remove spaces)
+*/
+
+#include "mt19937-64.h"
+
+#define NN 312
+#define MM 156
+#define MATRIX_A 0xB5026F5AA96619E9ULL
+#define UM 0xFFFFFFFF80000000ULL /* Most significant 33 bits */
+#define LM 0x7FFFFFFFULL /* Least significant 31 bits */
+
+/* initializes mt[NN] with a seed */
+void
+init_genrand64(mt19937_64 *context, unsigned long long seed)
+{+ context->mt[0] = seed;
+ for(context->mti=1; context->mti<NN; context->mti++)
+ context->mt[context->mti] = (6364136223846793005ULL * (context->mt[context->mti-1] ^ (context->mt[context->mti-1] >> 62)) + context->mti);
+}
+
+/* initialize by an array with array-length */
+/* init_key is the array for initializing keys */
+/* key_length is its length */
+void
+init_by_array64(mt19937_64 *context, unsigned long long init_key[], unsigned long long key_length)
+{+ unsigned long long i, j, k, m;
+
+ init_genrand64(context, 19650218ULL);
+ i = 1;
+ j = 0;
+ k = NN > key_length ? NN : key_length;
+ for(; k; k--){+ m = context->mt[i-1];
+ m = (m ^ (m >> 62)) * 3935559000370003845ULL;
+ context->mt[i] = context->mt[i] ^ m + init_key[j] + j; /* non linear */
+ i++;
+ j++;
+ if(i >= NN){+ context->mt[0] = context->mt[NN-1];
+ i = 1;
+ }
+ if(j >= key_length)
+ j = 0;
+ }
+ for(k = NN-1; k; k--){+ m = context->mt[i-1];
+ m = (m ^ (m >> 62)) * 2862933555777941757ULL;
+ context->mt[i] = context->mt[i] ^ m - i; /* non linear */
+ i++;
+ if(i >= NN){+ context->mt[0] = context->mt[NN-1];
+ i = 1;
+ }
+ }
+
+ context->mt[0] = 1ULL << 63; /* MSB is 1; assuring non-zero initial array */
+}
+
+/* generates a random number on [0, 2^64-1]-interval */
+unsigned long long
+genrand64_int64(mt19937_64 *context)
+{+ /* This is the altered Cocoa with Love implementation. */
+ int i, j;
+ unsigned long long result;
+
+ if(context->mti >= NN){ /* generate NN words at one time */+ int mid = NN / 2;
+ unsigned long long stateMid = context->mt[mid];
+ unsigned long long x;
+ unsigned long long y;
+
+ /* NOTE: this "untwist" code is modified from the original to improve
+ * performance, as described here:
+ * http://www.cocoawithlove.com/blog/2016/05/19/random-numbers.html
+ * These modifications are offered for use under the original icense at
+ * the top of this file.
+ */
+ for(i = 0, j = mid; i != mid - 1; i++, j++){+ x = (context->mt[i] & UM) | (context->mt[i + 1] & LM);
+ context->mt[i] = context->mt[i + mid] ^ (x >> 1) ^ ((context->mt[i + 1] & 1) * MATRIX_A);
+ y = (context->mt[j] & UM) | (context->mt[j + 1] & LM);
+ context->mt[j] = context->mt[j - mid] ^ (y >> 1) ^ ((context->mt[j + 1] & 1) * MATRIX_A);
+ }
+ x = (context->mt[mid - 1] & UM) | (stateMid & LM);
+ context->mt[mid - 1] = context->mt[NN - 1] ^ (x >> 1) ^ ((stateMid & 1) * MATRIX_A);
+ y = (context->mt[NN - 1] & UM) | (context->mt[0] & LM);
+ context->mt[NN - 1] = context->mt[mid - 1] ^ (y >> 1) ^ ((context->mt[0] & 1) * MATRIX_A);
+ context->mti = 0;
+ }
+
+ result = context->mt[context->mti];
+ context->mti = context->mti + 1;
+
+ result ^= (result >> 29) & 0x5555555555555555ULL;
+ result ^= (result << 17) & 0x71D67FFFEDA60000ULL;
+ result ^= (result << 37) & 0xFFF7EEE000000000ULL;
+ result ^= (result >> 43);
+
+ return result;
+}
+
+/* generates a random number on [0, 2^63-1]-interval */
+long long
+genrand64_int63(mt19937_64 *context)
+{+ return (long long)(genrand64_int64(context) >> 1);
+}
+
+/* generates a random number on [0,1]-real-interval */
+double
+genrand64_real1(mt19937_64 *context)
+{+ return (genrand64_int64(context) >> 11) * (1.0/9007199254740991.0);
+}
+
+/* generates a random number on [0,1)-real-interval */
+double
+genrand64_real2(mt19937_64 *context)
+{+ return (genrand64_int64(context) >> 11) * (1.0/9007199254740992.0);
+}
+
+/* generates a random number on (0,1)-real-interval */
+double
+genrand64_real3(mt19937_64 *context)
+{+ return ((genrand64_int64(context) >> 12) + 0.5) * (1.0/4503599627370496.0);
+}
--- /dev/null
+++ b/mt19937-64.h
@@ -1,0 +1,14 @@
+typedef struct mt19937_64 mt19937_64;
+
+struct mt19937_64 {+ unsigned long long mt[312];
+ int mti;
+};
+
+void init_genrand64(mt19937_64 *context, unsigned long long seed);
+void init_by_array64(mt19937_64 *context, unsigned long long init_key[], unsigned long long key_length);
+unsigned long long genrand64_int64(mt19937_64 *context);
+long long genrand64_int63(mt19937_64 *context);
+double genrand64_real1(mt19937_64 *context);
+double genrand64_real2(mt19937_64 *context);
+double genrand64_real3(mt19937_64 *context);
--- a/pr_cmds.c
+++ b/pr_cmds.c
@@ -457,20 +457,7 @@
static void
PF_random(pr_t *pr)
{- static double xmax;
- static long xand;
-
- if(xand == 0){- if((RAND_MAX & ((unsigned)RAND_MAX+1)) == 0){- xmax = (unsigned)RAND_MAX+1;
- xand = RAND_MAX;
- }else{- xmax = 0x8000;
- xand = 0x7fff;
- }
- }
-
- G_FLOAT(pr, OFS_RETURN) = ((double)(rand() & xand) + 0.5) / xmax;
+ G_FLOAT(pr, OFS_RETURN) = m_random_i0e1();
}
/*
--- a/quakedef.h
+++ b/quakedef.h
@@ -113,6 +113,7 @@
#include "zone.h"
#include "dat.h"
#include "mathlib.h"
+#include "m_random.h"
#include "fns.h"
#include "bspfile.h"
#include "vid.h"
--- a/sys_plan9.c
+++ b/sys_plan9.c
@@ -176,7 +176,9 @@
break;
default: usage();
}ARGEND
- srand(getpid());
+
+ m_random_init(time(nil));
+ srand(time(nil));
/* ignore fp exceptions: rendering shit assumes they are */
setfcr(getfcr() & ~(FPOVFL|FPUNFL|FPINVAL|FPZDIV));
notify(croak);
--- a/sys_unix.c
+++ b/sys_unix.c
@@ -57,12 +57,6 @@
return ts;
}
-int
-nrand(int n)
-{- return random() % n;
-}
-
void
fatal(char *fmt, ...)
{@@ -187,7 +181,8 @@
if(snailenabled)
initsnail();
- srand(nanosec() + time(nil));
+ m_random_init(time(nil));
+ srand(time(nil));
paths[1] = strdup(va("%s/.quake", getenv("HOME")));Host_Init(nargs, argv, paths);
--- a/unix/platform.h
+++ b/unix/platform.h
@@ -61,6 +61,5 @@
#define werrstr(fmt...) do{snprint(lasterr, sizeof(lasterr), fmt); }while(0)char *seprint(char *, char *, char *, ...);
-int nrand(int);
#define DotProduct(x,y) DotProduct_((x),(y))