]> git.kaiwu.me - klib.git/commitdiff
removed MT19937-64.
authorHeng Li <lh3@me.com>
Fri, 31 Aug 2018 08:58:52 +0000 (22:58 -1000)
committerHeng Li <lh3@me.com>
Fri, 31 Aug 2018 08:58:52 +0000 (22:58 -1000)
kmath.c
kmath.h

diff --git a/kmath.c b/kmath.c
index e51b3a78bc07c8b20d3e32879d8397e0aa40bea9..9368b2c9dc7074ca8f4d9f753a547e758eef6551 100644 (file)
--- a/kmath.c
+++ b/kmath.c
@@ -3,105 +3,6 @@
 #include <math.h>
 #include "kmath.h"
 
-/**************************************
- *** Pseudo-random number generator ***
- **************************************/
-
-/* 
-   64-bit Mersenne Twister pseudorandom number generator. Adapted from:
-
-     http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/VERSIONS/C-LANG/mt19937-64.c
-
-   which was written by Takuji Nishimura and Makoto Matsumoto and released
-   under the 3-clause BSD license.
-*/
-
-#define KR_NN 312
-#define KR_MM 156
-#define KR_UM 0xFFFFFFFF80000000ULL /* Most significant 33 bits */
-#define KR_LM 0x7FFFFFFFULL /* Least significant 31 bits */
-
-struct _krand_t {
-       int mti, n_iset;
-       double n_gset;
-       krint64_t mt[KR_NN];
-};
-
-static void kr_srand0(krint64_t seed, krand_t *kr)
-{
-       kr->mt[0] = seed;
-       for (kr->mti = 1; kr->mti < KR_NN; ++kr->mti) 
-               kr->mt[kr->mti] = 6364136223846793005ULL * (kr->mt[kr->mti - 1] ^ (kr->mt[kr->mti - 1] >> 62)) + kr->mti;
-}
-
-krand_t *kr_srand(krint64_t seed)
-{
-       krand_t *kr;
-       kr = (krand_t*)calloc(1, sizeof(krand_t));
-       kr_srand0(seed, kr);
-       return kr;
-}
-
-krint64_t kr_rand(krand_t *kr)
-{
-       krint64_t x;
-       static const krint64_t mag01[2] = { 0, 0xB5026F5AA96619E9ULL };
-    if (kr->mti >= KR_NN) {
-               int i;
-               if (kr->mti == KR_NN + 1) kr_srand0(5489ULL, kr);
-        for (i = 0; i < KR_NN - KR_MM; ++i) {
-            x = (kr->mt[i] & KR_UM) | (kr->mt[i+1] & KR_LM);
-            kr->mt[i] = kr->mt[i + KR_MM] ^ (x>>1) ^ mag01[(int)(x&1)];
-        }
-        for (; i < KR_NN - 1; ++i) {
-            x = (kr->mt[i] & KR_UM) | (kr->mt[i+1] & KR_LM);
-            kr->mt[i] = kr->mt[i + (KR_MM - KR_NN)] ^ (x>>1) ^ mag01[(int)(x&1)];
-        }
-        x = (kr->mt[KR_NN - 1] & KR_UM) | (kr->mt[0] & KR_LM);
-        kr->mt[KR_NN - 1] = kr->mt[KR_MM - 1] ^ (x>>1) ^ mag01[(int)(x&1)];
-        kr->mti = 0;
-    }
-    x = kr->mt[kr->mti++];
-    x ^= (x >> 29) & 0x5555555555555555ULL;
-    x ^= (x << 17) & 0x71D67FFFEDA60000ULL;
-    x ^= (x << 37) & 0xFFF7EEE000000000ULL;
-    x ^= (x >> 43);
-    return x;
-}
-
-double kr_normal(krand_t *kr)
-{ 
-       if (kr->n_iset == 0) {
-               double fac, rsq, v1, v2; 
-               do { 
-                       v1 = 2.0 * kr_drand(kr) - 1.0;
-                       v2 = 2.0 * kr_drand(kr) - 1.0; 
-                       rsq = v1 * v1 + v2 * v2;
-               } while (rsq >= 1.0 || rsq == 0.0);
-               fac = sqrt(-2.0 * log(rsq) / rsq); 
-               kr->n_gset = v1 * fac; 
-               kr->n_iset = 1;
-               return v2 * fac;
-       } else {
-               kr->n_iset = 0;
-               return kr->n_gset;
-       }
-}
-
-#ifdef _KR_MAIN
-int main(int argc, char *argv[])
-{
-       long i, N = 200000000;
-       krand_t *kr;
-       if (argc > 1) N = atol(argv[1]);
-       kr = kr_srand(11);
-       for (i = 0; i < N; ++i) kr_rand(kr);
-//     for (i = 0; i < N; ++i) lrand48();
-       free(kr);
-       return 0;
-}
-#endif
-
 /******************************
  *** Non-linear programming ***
  ******************************/
diff --git a/kmath.h b/kmath.h
index c5f7933bb6234c0a6540d63f2630ef84ab7655aa..1815a14c8acac37c5008f3903cbc589dc6994ef0 100644 (file)
--- a/kmath.h
+++ b/kmath.h
@@ -7,22 +7,6 @@
 extern "C" {
 #endif
 
-       /**********************************
-        * Pseudo-random number generator *
-        **********************************/
-
-       typedef uint64_t krint64_t;
-
-       struct _krand_t;
-       typedef struct _krand_t krand_t;
-
-       #define kr_drand(_kr) ((kr_rand(_kr) >> 11) * (1.0/9007199254740992.0))
-       #define kr_sample(_kr, _k, _cnt) ((*(_cnt))++ < (_k)? *(_cnt) - 1 : kr_rand(_kr) % *(_cnt))
-
-       krand_t *kr_srand(krint64_t seed);
-       krint64_t kr_rand(krand_t *kr);
-       double kr_normal(krand_t *kr);
-
        /**************************
         * Non-linear programming *
         **************************/