function out = CMWC32(s1,s2,n) # 32 bit CMWC
a = 18782;
c = s1;
x = s2;
b = 2^32-1;
# initialize state vector
s = zeros(4097,1);
for i=1:4096
x = b-bitand(uint64(a*x+c),b);
c = bitshift(uint64(a*x+c),-32);
s(i) = x;
endfor
# Now compute new values for use
for i = 1:n
x = s(1);
s(4097) = b-bitand(uint64(a*x+c),b);
c = bitshift(uint64(a*x+c),-32);
disp(double(s(4097))/2^32)
s = s(2:4097);
endfor
--------------------
No doubt this could be streamlined and optimized considerably. For the very simple C code which describes this PRNG, see
here.
-Alasdair
On Tue, Jun 5, 2012 at 12:25 AM, Alasdair McAndrew
<address@hidden> wrote:
Here's a simple function which generates n pseudo-random numbers according to the 2006 Wichmann-Hill method, which is much stronger than the one linked to above, and which passes the "Big Crush" test, as well as having a period of 2^121:
whrandom.m
-----------------
# Implements the 2006 Wichmann-Hill PRNG
# Computes n pseudo-random numbers starting with a seed of 4 user-chosen
# values.
function out = whrandom(ix,iy,iz,iw,n)
out = zeros(n,1);
for i=1:n
ix = mod(11600*ix,p1);
iy = mod(47003*iy,p2);
iz = mod(23000*iz,p3);
iw = mod(30000*iw,p4);
out(i) = mod(ix/p1+iy/p2+iz/p3+iw/p4,1);
end
-----------------
-Alasdair
On Mon, Jun 4, 2012 at 10:26 PM, Alasdair McAndrew
<address@hidden> wrote:
Writing your own is a good idea. Cleve Moler, here, speaks of the linear congruential generator
x_{k+1} = ax_k (mod m)
with a,m = 16807, 2^31-1
which was in fact Matlab's internal random number generator until replaced by the Mersenne Twister. You'll find some good links here; there are some PRNGs developed by the late George Marsaglia which are both far simpler, and have a far longer period, than the Mersenne Twister.
For a Matlab implementation of the Wichmann-Hill PRNG (which works basically by adding up three floating point values and taking the fractional part as output), seehere. This is a simple program which should work equally well in Matlab and Octave.
-Alasdair
On Mon, Jun 4, 2012 at 12:03 AM, Mathieu Dubois
<address@hidden> wrote:
Hello,
I want to comapre the output of the Self-Organinzin Map algorithm under Matlab and Octave.
The algorithm is stochastic because, as in most training algorithm, the data are presented inrandom order.
In order to compare as precisely as possible the Matlab and Octave version it seems necessary to have the same sequence of random numbers.
Is there a way to do that?
I have tried the following:
rand('twister', 10);
s = rand(10, 1);
disp(s);
disp(mean(s));
disp(var(s));
The sequences are different (but the statistics are close).
Mathieu
_______________________________________________
Help-octave mailing list
address@hidden
https://mailman.cae.wisc.edu/listinfo/help-octave
--
Blog:
http://amca01.wordpress.comWeb:
http://sites.google.com/site/amca01/
Facebook:
http://www.facebook.com/alasdair.mcandrew
--
Blog:
http://amca01.wordpress.comWeb:
http://sites.google.com/site/amca01/
Facebook:
http://www.facebook.com/alasdair.mcandrew