Here is a comparison of the time it takes to generate 100,000,000 random
numbers with three different generators, compiled with various compilers
on Linux.

rand31: linear congruential generator (m = 0x7FFFFFFF, a = 742938285)
rand: ANSI C rand()
(note that it is much slower if not compiled with -static)
MT: Mersenne Twister (MT19937)

------------------------------------------------------------------------

GCC 4.1 beta (-O2 -march=pentium3 -static)

rand31: 1.703, rand: 4.559, MT: 2.501

GCC 4.1 beta (-O2 -march=pentium3 -fomit-frame-pointer -static)

rand31: 1.731, rand: 4.561, MT: 2.321

GCC 4.1 beta (-O2 -march=pentium3)

rand31: 1.702, rand: 7.004, MT: 2.500

GCC 4.0.0 (-O2 -march=pentium3 -static)

rand31: 1.702, rand: 4.561, MT: 2.405

GCC 4.0.0 (-O2 -march=pentium3 -fomit-frame-pointer -static)

rand31: 1.732, rand: 4.562, MT: 2.339

GCC 3.3.4 (-O2 -march=pentium3 -static)

rand31: 1.863, rand: 4.172, MT: 2.579

GCC 3.3.4 (-O2 -march=pentium3 -fomit-frame-pointer -static)

rand31: 1.730, rand: 4.172, MT: 2.385

ICC 9.0 (-O2 -xK -static)

rand31: 1.841, rand: 4.410, MT: 2.474

------------------------------------------------------------------------

#include <stdio.h>
#include <stdlib.h>
#include <stdint.h>
#include <sys/time.h>

#define CS_PURE __attribute__ ((__pure__))
#define CS_NOINLINE __attribute__ ((__noinline__))
#define PUBLIC

typedef struct CsoundRandMTState_ {
int mti;
uint32_t mt[624];
} CsoundRandMTState;

PUBLIC CS_PURE int csoundRand31(int seedVal)
{
uint64_t tmp1;
uint32_t tmp2;

/* x = (742938285 * x) % 0x7FFFFFFF */
tmp1 = (uint64_t) ((int32_t) seedVal * (int64_t) 742938285);
tmp2 = (uint32_t) tmp1 & (uint32_t) 0x7FFFFFFF;
tmp2 += (uint32_t) (tmp1 >> 31);
if ((int32_t) tmp2 < (int32_t) 0)
tmp2 = (tmp2 + (uint32_t) 1) & (uint32_t) 0x7FFFFFFF;
return (int) tmp2;
}

/* Period parameters */

#define N 624
#define M 397
#define MATRIX_A 0x9908B0DFU /* constant vector a */
#define UPPER_MASK 0x80000000U /* most significant w-r bits */
#define LOWER_MASK 0x7FFFFFFFU /* least significant r bits */

/* initializes mt[N] with a seed */

PUBLIC void csoundSeedRandMT(CsoundRandMTState *p, uint32_t seedVal)
{
int i;
uint32_t x;

p->mt[0] = x = seedVal;
for (i = 1; i < N; i++) {
/* See Knuth TAOCP Vol2. 3rd Ed. P.106 for multiplier. */
/* In the previous versions, MSBs of the seed affect */
/* only MSBs of the array mt[]. */
/* 2002/01/09 modified by Makoto Matsumoto */
x = ((uint32_t) 1812433253 * (x ^ (x >> 30)) + (uint32_t) i);
p->mt[i] = x;
}
p->mti = N;
}

static CS_NOINLINE void MT_update_state(CsoundRandMTState *p)
{
/* mag01[x] = x * MATRIX_A for x=0,1 */
const uint32_t mag01[2] = { (uint32_t) 0, (uint32_t) MATRIX_A };
int i;
uint32_t y;
for (i = 0; i < (N - M); i++) {
y = (p->mt[i] & UPPER_MASK) | (p->mt[i + 1] & LOWER_MASK);
p->mt[i] = p->mt[i + M] ^ (y >> 1) ^ mag01[y & (uint32_t) 1];
}
for ( ; i < (N - 1); i++) {
y = (p->mt[i] & UPPER_MASK) | (p->mt[i + 1] & LOWER_MASK);
p->mt[i] = p->mt[i + (M - N)] ^ (y >> 1) ^ mag01[y & (uint32_t) 1];
}
y = (p->mt[N - 1] & UPPER_MASK) | (p->mt[0] & LOWER_MASK);
p->mt[N - 1] = p->mt[M - 1] ^ (y >> 1) ^ mag01[y & (uint32_t) 1];
}

/* generates a random number on [0,0xffffffff]-interval */

PUBLIC uint32_t csoundRandMT(CsoundRandMTState *p)
{
int i = p->mti;
uint32_t y;

if (i >= N) { /* generate N words at one time */
MT_update_state(p);
i = 0;
}
y = p->mt[i];
p->mti = i + 1;
/* Tempering */
y ^= (y >> 11);
y ^= (y << 7) & (uint32_t) 0x9D2C5680U;
y ^= (y << 15) & (uint32_t) 0xEFC60000U;
y ^= (y >> 18);

return y;
}

volatile uint32_t seed_;

int main(void)
{
CsoundRandMTState st;
double t0, t1, t2, t3;
struct timeval tv;
int i;

gettimeofday(&tv, NULL);
t0 = (double) tv.tv_usec * 0.000001 + (double) tv.tv_sec;
seed_ = 12345;
for (i = 0; i < 100000000; i++)
seed_ = csoundRand31(seed_);
gettimeofday(&tv, NULL);
t1 = (double) tv.tv_usec * 0.000001 + (double) tv.tv_sec;
srand(314159);
for (i = 0; i < 100000000; i++)
seed_ = rand();
gettimeofday(&tv, NULL);
t2 = (double) tv.tv_usec * 0.000001 + (double) tv.tv_sec;
csoundSeedRandMT(&st, 5489);
for (i = 0; i < 100000000; i++) {
seed_ = csoundRandMT(&st);
/* if (i < 500) {
printf("%10u ", seed_);
if ((i % 5) == 4) printf("\n");
} */
}
gettimeofday(&tv, NULL);
t3 = (double) tv.tv_usec * 0.000001 + (double) tv.tv_sec;
printf("rand31: %.3f, rand: %.3f, MT: %.3f\n",
(t1 - t0), (t2 - t1), (t3 - t2));
return 0;
}
--
Send bugs reports to this list.
To unsubscribe, send email to csound-unsubscribe@...