10
typedef signed int sint;
11
typedef unsigned int uint;
12
typedef signed long slong;
13
typedef unsigned long ulong;
19
RNG(ulong z_, ulong w_, ulong jsr_, ulong jcong_ ) { };
22
void init(ulong z_, ulong w_, ulong jsr_, ulong jcong_ ) { }
23
long rand_int31() { return RNG_INT31(); }
24
double rand_halfopen01() // (0,1]
25
{ return RNG_UNIF01(); }
26
int binomial(double pp, int n) { return RNG_BINOM(n,pp); }
31
/* This was the original RNG, but now we use the igraph version */
33
// __________________________________________________________________________
34
// random.h - a Random Number Generator Class
35
// random.cpp - contains the non-inline class methods
37
// __________________________________________________________________________
38
// This C++ code uses the simple, very fast "KISS" (Keep It Simple
39
// Stupid) random number generator suggested by George Marsaglia in a
40
// Usenet posting from 1999. He describes it as "one of my favorite
41
// generators". It generates high-quality random numbers that
42
// apparently pass all commonly used tests for randomness. In fact, it
43
// generates random numbers by combining the results of three other good
44
// random number generators that have different periods and are
45
// constructed from completely different algorithms. It does not have
46
// the ultra-long period of some other generators - a "problem" that can
47
// be fixed fairly easily - but that seems to be its only potential
48
// problem. The period is about 2^123.
50
// The ziggurat method of Marsaglia is used to generate exponential and
51
// normal variates. The method as well as source code can be found in
52
// the article "The Ziggurat Method for Generating Random Variables" by
53
// Marsaglia and Tsang, Journal of Statistical Software 5, 2000.
55
// The method for generating gamma variables appears in "A Simple Method
56
// for Generating Gamma Variables" by Marsaglia and Tsang, ACM
57
// Transactions on Mathematical Software, Vol. 26, No 3, Sep 2000, pages
60
// The code for Poisson and Binomial random numbers comes from
61
// Numerical Recipes in C.
63
// Some of this code is unlikely to work correctly as is on 64 bit
69
// #include <process.h>
70
// #define getpid _getpid
72
// #include <unistd.h>
76
// static const double PI = 3.1415926535897932;
77
// static const double AD_l = 0.6931471805599453;
78
// static const double AD_a = 5.7133631526454228;
79
// static const double AD_b = 3.4142135623730950;
80
// static const double AD_c = -1.6734053240284925;
81
// static const double AD_p = 0.9802581434685472;
82
// static const double AD_A = 5.6005707569738080;
83
// static const double AD_B = 3.3468106480569850;
84
// static const double AD_H = 0.0026106723602095;
85
// static const double AD_D = 0.0857864376269050;
93
// ulong z, w, jsr, jcong; // Seeds
95
// ulong kn[128], ke[256];
96
// double wn[128],fn[128], we[256],fe[256];
100
// static const double PI = 3.1415926535897932;
101
// static const double AD_l = 0.6931471805599453;
102
// static const double AD_a = 5.7133631526454228;
103
// static const double AD_b = 3.4142135623730950;
104
// static const double AD_c = -1.6734053240284925;
105
// static const double AD_p = 0.9802581434685472;
106
// static const double AD_A = 5.6005707569738080;
107
// static const double AD_B = 3.3468106480569850;
108
// static const double AD_H = 0.0026106723602095;
109
// static const double AD_D = 0.0857864376269050;
114
// RNG() { init(); zigset(); }
115
// RNG(ulong z_, ulong w_, ulong jsr_, ulong jcong_ ) :
116
// z(z_), w(w_), jsr(jsr_), jcong(jcong_) { zigset(); }
120
// inline ulong znew()
121
// { return (z = 36969 * (z & 65535) + (z >> 16)); }
122
// inline ulong wnew()
123
// { return (w = 18000 * (w & 65535) + (w >> 16)); }
124
// inline ulong MWC()
125
// { return (((znew() & 65535) << 16) + wnew()); }
126
// inline ulong SHR3()
127
// { jsr ^= ((jsr & 32767) << 17); jsr ^= (jsr >> 13); return (jsr ^= ((jsr << 5) & 0xFFFFFFFF)); }
128
// inline ulong CONG()
129
// { return (jcong = (69069 * jcong + 1234567) & 0xFFFFFFFF); }
130
// inline double RNOR() {
131
// slong h = rand_int32();
132
// ulong i = h & 127;
133
// return (((ulong) abs((sint) h) < kn[i]) ? h * wn[i] : nfix(h, i));
135
// inline double REXP() {
136
// ulong j = rand_int32();
137
// ulong i = j & 255;
138
// return ((j < ke[i]) ? j * we[i] : efix(j, i));
141
// double nfix(slong h, ulong i);
142
// double efix(ulong j, ulong i);
145
// inline void init()
146
// { ulong yo = time(0) + getpid();
147
// z = w = jsr = jcong = yo; }
148
// inline void init(ulong z_, ulong w_, ulong jsr_, ulong jcong_ )
149
// { z = z_; w = w_; jsr = jsr_; jcong = jcong_; }
151
// inline ulong rand_int32() // [0,2^32-1]
152
// { return ((MWC() ^ CONG()) + SHR3()) & 0xFFFFFFFF; }
153
// inline long rand_int31() // [0,2^31-1]
154
// { return long(rand_int32() >> 1);}
155
// inline double rand_closed01() // [0,1]
156
// { return ((double) rand_int32() / 4294967295.0); }
157
// inline double rand_open01() // (0,1)
158
// { return (((double) rand_int32() + 0.5) / 4294967296.0); }
159
// inline double rand_halfclosed01() // [0,1)
160
// { return ((double) rand_int32() / 4294967296.0); }
161
// inline double rand_halfopen01() // (0,1]
162
// { return (((double) rand_int32() + 0.5) / 4294967295.5); }
164
// // Continuous Distributions
165
// inline double uniform(double x = 0.0, double y = 1.0)
166
// { return rand_closed01() * (y - x) + x; }
167
// inline double normal(double mu = 0.0, double sd = 1.0)
168
// { return RNOR() * sd + mu; }
169
// inline double exponential(double lambda = 1)
170
// { return REXP() / lambda; }
171
// double gamma(double shape = 1, double scale = 1);
172
// double chi_square(double df)
173
// { return gamma(df / 2.0, 0.5); }
174
// double beta(double a1, double a2)
175
// { double x1 = gamma(a1, 1); return (x1 / (x1 + gamma(a2, 1))); }
177
// // Discrete Distributions
178
// double poisson(double lambda);
179
// int binomial(double pp, int n);