linear congruential generator
- To: mathgroup at smc.vnet.net
- Subject: [mg4770] linear congruential generator
- From: cabana at flinet.com (David Cabana)
- Date: Fri, 13 Sep 1996 13:54:52 -0400
- Organization: Florida Internet Corp
- Sender: owner-wri-mathgroup at wolfram.com
This is in response to a recent posting regarding (among other things), a
linear congrential generator. Here is what I got when I ran the proposed
code:
In[1]:=
x[0] = 1;
x[n+1] := Mod[(7^7*x[n]), 2^31 -1]//N
u[n] := x[n]/(2^31-1)
Table[x[n], {n,1,10}]
Out[4]=
{x[1], x[2], x[3], x[4], x[5], x[6], x[7], x[8], x[9], x[10]}
Here is roughly the same code, with some changes to make it work.
x[0] = 1;
x[n_] := Mod[(7^7 * x[n-1]), 2^31 - 1]
u[n_] := x[n]/(2^31-1) //N
Table[x[n], {n,1,10}]
Notice the changes in the definitions of x and u. The underscores are
critical. I also expressed the recursion a bit differently. The generator
must use exact integer arithmetic to be correct, so I moved the //N to the
definition of u.
There is another problem. As written, this program is very slow. The
reason is that to compute x[100] it must compute x[99], x[98], ..., x[1].
Then to compute x[101], it must compute x[100], x[99], x[98], ... x[1]. It
computes the same thing over and over again. There are various fixes to
this. One is to use memoization, say by writing
x[n_] := x[n] = Mod[(7^7 * x[n-1]), 2^31 - 1]
instead of
x[n_] := Mod[(7^7 * x[n-1]), 2^31 - 1]
Another is to take an entirely different approach:
generator[n_]:= Mod[7^7 * n, 2^31 - 1]
uniformVariates[n_, seed_] :=
NestList[generator, seed, n] / (2^31 - 1) //N
I did not run timing tests, but I suspect the 2nd approach is quicker.
--
David Cabana cabana at flinet.com
==== [MESSAGE SEPARATOR] ====