Make improvements

- Make rand64() thread safe
- Introduce lemur64 lcg prng
- Improve strace on New Technology
- Improve msync() on New Technology

libc/rand/getrandom.c

@@ -43,7 +43,7 @@ static bool have_getrandom;
 /**
  * Returns cryptographic random data.
  *
- * This random number seed generator blends information from:
+ * This random number seed generator obtains information from:
  *
  * - getrandom() on Linux
  * - RtlGenRandom() on Windows
@@ -61,6 +61,9 @@ static bool have_getrandom;
  * This function is safe to use with fork() and vfork(). It will also
  * close any file descriptor it ends up needing before it returns.
  *
+ * @note this function could block a nontrivial time on old computers
+ * @note this function is indeed intended for cryptography
+ * @note this function takes around 900 cycles
  * @asyncsignalsafe
  * @restartable
  * @vforksafe

libc/rand/lcg.internal.h

@@ -7,7 +7,7 @@ forceinline uint64_t KnuthLinearCongruentialGenerator(uint64_t prev[1]) {
      Seminumerical Algorithms, Third Edition, Addison-Wesley, 1998,
      p. 106 (line 26) & p. 108 */
   prev[0] = prev[0] * 6364136223846793005 + 1442695040888963407;
-  return prev[0];
+  return prev[0]; /* be sure to shift! */
 }
 
 #endif /* !(__ASSEMBLER__ + __LINKER__ + 0) */

libc/rand/lemur64.c

@@ -0,0 +1,45 @@
+/*-*- mode:c;indent-tabs-mode:nil;c-basic-offset:2;tab-width:8;coding:utf-8 -*-│
+│vi: set net ft=c ts=2 sts=2 sw=2 fenc=utf-8                                :vi│
+╞══════════════════════════════════════════════════════════════════════════════╡
+│ Copyright 2022 Justine Alexandra Roberts Tunney                              │
+│                                                                              │
+│ Permission to use, copy, modify, and/or distribute this software for         │
+│ any purpose with or without fee is hereby granted, provided that the         │
+│ above copyright notice and this permission notice appear in all copies.      │
+│                                                                              │
+│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL                │
+│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED                │
+│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE             │
+│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL         │
+│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR        │
+│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER               │
+│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR             │
+│ PERFORMANCE OF THIS SOFTWARE.                                                │
+╚─────────────────────────────────────────────────────────────────────────────*/
+#include "libc/rand/rand.h"
+
+/**
+ * Returns linear congruential deterministic pseudorandom data, e.g.
+ *
+ *     uint64_t x = lemur64();
+ *
+ * You can generate different types of numbers as follows:
+ *
+ *     int64_t x = lemur64() >> 1;    // make positive signed integer
+ *     double x = _real1(lemur64());  // make float on [0,1]-interval
+ *
+ * If you want a fast pseudorandom number generator that seeds itself
+ * automatically on startup and fork() then consider rand64(). If you
+ * want true random data then consider rdseed, rdrand, and getrandom.
+ *
+ * @return 64 bits of pseudorandom data
+ * @note this is Lemire's Lehmer generator
+ * @note this function takes at minimum 1 cycle
+ * @note this function passes bigcrush and practrand
+ * @note this function is not intended for cryptography
+ * @see rand64(), rngset(), _real1(), _real2(), _real3()
+ */
+uint64_t lemur64(void) {
+  static uint128_t s = 2131259787901769494;
+  return (s *= 15750249268501108917ull) >> 64;
+}

libc/rand/rand.c

@@ -21,11 +21,18 @@
 #include "libc/rand/rand.h"
 
 /**
- * Returns 31-bit random number using a linear congruential generator.
+ * Returns 31-bit linear congruential pseudorandom number, e.g.
  *
- * Please note that, unlike rand32(), the rand() function uses the same
- * seed at startup by default, unless srand() is called. This makes it
- * useful in cases where deterministic behavior is needed.
+ *     int x = rand();
+ *     assert(x >= 0);
+ *
+ * This function always returns a positive number. If srand() isn't
+ * called, then it'll return the same sequence each time your program
+ * runs. Faster and more modern alternatives exist to this function.
+ *
+ * @note this function does well on bigcrush and practrand
+ * @note this function is not intended for cryptography
+ * @see lemur64(), rand64(), rdrand()
  */
 int rand(void) {
   return KnuthLinearCongruentialGenerator(&g_rando) >> 33;

libc/rand/rand.h

@@ -23,12 +23,13 @@ char *setstate(char *);
 long random(void);
 void srandom(unsigned);
 
+uint64_t lemur64(void);
+uint64_t rand64(void);
 uint64_t vigna(void);
 uint64_t vigna_r(uint64_t[hasatleast 1]);
 void svigna(uint64_t);
 uint64_t rdrand(void);
 uint64_t rdseed(void);
-uint64_t rand64(void);
 void _smt19937(uint64_t);
 void _Smt19937(uint64_t[], size_t);
 uint64_t _mt19937(void);

libc/rand/rand64.c

@@ -16,59 +16,80 @@
 │ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR             │
 │ PERFORMANCE OF THIS SOFTWARE.                                                │
 ╚─────────────────────────────────────────────────────────────────────────────*/
+#include "libc/assert.h"
 #include "libc/bits/bits.h"
-#include "libc/bits/xadd.h"
+#include "libc/bits/lockcmpxchg16b.internal.h"
+#include "libc/bits/lockxadd.internal.h"
+#include "libc/bits/weaken.h"
 #include "libc/calls/calls.h"
+#include "libc/intrin/kprintf.h"
 #include "libc/nexgen32e/rdtsc.h"
-#include "libc/nexgen32e/x86feature.h"
+#include "libc/nt/thread.h"
 #include "libc/rand/rand.h"
 #include "libc/runtime/runtime.h"
 #include "libc/sysv/consts/auxv.h"
+#include "libc/thread/create.h"
 
-static uint64_t thepool;
+extern int __pid;
+static int thepid;
+static int thecount;
+static uint128_t thepool;
 
 /**
  * Returns nondeterministic random data.
  *
- * This function automatically seeds itself on startup and reseeds
- * itself after fork() and vfork(). It takes about nanosecond to run.
- * That makes it much slower than vigna() and rand() but much faster
- * than rdrand() and rdseed().
+ * This function is similar to lemur64() except it'thepool intended to
+ * be unpredictable. This PRNG automatically seeds itself on startup
+ * using a much stronger and faster random source than `srand(time(0))`.
+ * This function will automatically reseed itself when new processes and
+ * threads are spawned. This function is thread safe in the sense that a
+ * race condition can't happen where two threads return the same result.
  *
  * @see rdseed(), rdrand(), rand(), random(), rngset()
- * @note based on vigna's algorithm
+ * @note this function is not intended for cryptography
+ * @note this function passes bigcrush and practrand
+ * @note this function takes at minimum 30 cycles
  * @asyncsignalsafe
+ * @threadsafe
  * @vforksafe
  */
 uint64_t rand64(void) {
-  bool cf;
-  register uint64_t t;
-  if (X86_HAVE(RDSEED)) {
-    asm volatile(CFLAG_ASM("rdseed\t%1")
-                 : CFLAG_CONSTRAINT(cf), "=r"(t)
-                 : /* no inputs */
-                 : "cc");
-    if (cf) {
-      thepool ^= t;
-      return t;
+  void *d;
+  int c1, p1, p2;
+  uint128_t s1, s2;
+  do {
+    p1 = __pid;
+    p2 = thepid;
+    c1 = thecount;
+    asm volatile("" ::: "memory");
+    s1 = thepool;
+    if (p1 == p2) {
+      // fast path
+      s2 = s1;
+    } else {
+      // slow path
+      if (!p2) {
+        // first call so get some cheap entropy
+        if ((d = (void *)getauxval(AT_RANDOM))) {
+          memcpy(&s2, d, 16);  // kernel entropy
+        } else {
+          s2 = kStartTsc;  // rdtsc() @ _start()
+        }
+      } else {
+        // process contention so blend a timestamp
+        s2 = s1 ^ rdtsc();
+      }
+      // toss the new pid in there
+      s2 ^= p1;
+      // ordering for thepid probably doesn't matter
+      thepid = p1;
     }
-  }
-  t = _xadd(&thepool, 0x9e3779b97f4a7c15);
-  t ^= (getpid() * 0x1001111111110001ull + 0xdeaadead) >> 31;
-  t = (t ^ (t >> 30)) * 0xbf58476d1ce4e5b9;
-  t = (t ^ (t >> 27)) * 0x94d049bb133111eb;
-  return t ^ (t >> 31);
+    // lemur64 pseudorandom number generator
+    s2 *= 15750249268501108917ull;
+    // sadly 128-bit values aren't atomic on x86
+    _lockcmpxchg16b(&thepool, &s1, s2);
+    // do it again if there's thread contention
+  } while (_lockxadd(&thecount, 1) != c1);
+  // the most important step in the prng
+  return s2 >> 64;
 }
-
-static textstartup void rand64_init(int argc, char **argv, char **envp,
-                                    intptr_t *auxv) {
-  for (; auxv[0]; auxv += 2) {
-    if (auxv[0] == AT_RANDOM) {
-      thepool = READ64LE((const char *)auxv[1]);
-      return;
-    }
-  }
-  thepool = kStartTsc;
-}
-
-const void *const g_rand64_init[] initarray = {rand64_init};

libc/rand/rdrand.c

@@ -55,9 +55,12 @@ static dontinline uint64_t rdrand_failover(void) {
  * aren't available then we try /dev/urandom and if that fails, we try
  * getauxval(AT_RANDOM), and if not we finally use RDTSC and getpid().
  *
- * This function takes between 10 nanoseconds to several microseconds.
- *
+ * @note this function could block a nontrivial time on old computers
+ * @note this function is indeed intended for cryptography
+ * @note this function takes around 300 cycles
  * @see rngset(), rdseed(), rand64()
+ * @asyncsignalsafe
+ * @vforksafe
  */
 uint64_t rdrand(void) {
   int i;

libc/rand/rdseed.c

@@ -31,9 +31,12 @@
  * sysctl(KERN_ARND). If those aren't available then we try /dev/urandom
  * and if that fails, we use RDTSC and getpid().
  *
- * This function takes about 32 nanoseconds.
- *
+ * @note this function could block a nontrivial time on old computers
+ * @note this function is indeed intended for cryptography
+ * @note this function takes around 800 cycles
  * @see rngset(), rdrand(), rand64()
+ * @asyncsignalsafe
+ * @vforksafe
  */
 uint64_t rdseed(void) {
   int i;

libc/rand/real1.c

@@ -21,9 +21,9 @@
 /**
  * Generates number on [0,1]-real-interval, e.g.
  *
- *     double x = _real1(vigna())
+ *     double x = _real1(lemur64())
  *
- * @see vigna(), mt19937()
+ * @see lemur64(), mt19937()
  */
 double _real1(uint64_t x) {
   return 1. / 9007199254740991. * (x >> 11);

libc/rand/real2.c

@@ -21,9 +21,9 @@
 /**
  * Generates number on [0,1)-real-interval, e.g.
  *
- *     double x = _real2(vigna())
+ *     double x = _real2(lemur64())
  *
- * @see vigna(), mt19937()
+ * @see lemur64(), mt19937()
  */
 double _real2(uint64_t x) {
   return 1. / 9007199254740992. * (x >> 11);

libc/rand/real3.c

@@ -21,9 +21,9 @@
 /**
  * Generates number on (0,1)-real-interval, e.g.
  *
- *     double x = _real3(vigna())
+ *     double x = _real3(lemur64())
  *
- * @see vigna(), mt19937()
+ * @see lemur64(), mt19937()
  */
 double _real3(uint64_t x) {
   return 1. / 4503599627370496. * ((x >> 12) + .5);

libc/rand/vigna.c

@@ -16,7 +16,6 @@
 │ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR             │
 │ PERFORMANCE OF THIS SOFTWARE.                                                │
 ╚─────────────────────────────────────────────────────────────────────────────*/
-#include "libc/bits/xadd.h"
 #include "libc/rand/rand.h"
 
 static uint64_t g_vigna;
@@ -40,9 +39,8 @@ static uint64_t g_vigna;
  *     static uint64_t s = 0;
  *     uint64_t x = svigna_r(&s);
  *
- * This function is the fastest way to generate good scalar pseudorandom
- * numbers that aren't truncated. If you want to fill a buffer with data
- * then rngset() implements vigna's algorithm to do that extremely well:
+ * If you want to fill a buffer with data then rngset() implements
+ * vigna's algorithm to do that extremely well:
  *
  *     char buf[4096];
  *     rngset(buf, sizeof(buf), vigna, 0);
@@ -52,6 +50,9 @@ static uint64_t g_vigna;
  * want true random data then consider rdseed, rdrand, and getrandom.
  *
  * @return 64 bits of pseudorandom data
+ * @note this function is not intended for cryptography
+ * @note this function passes bigcrush and practrand
+ * @note this function takes at minimum 4 cycles
  */
 uint64_t vigna(void) {
   return vigna_r(&g_vigna);