a fast

linear congruential generator (or specifically the

park-miller variation) is described below.

basic algorithm:

Code:

@init
rnd=1;
@sample
rnd=rnd*a%m;

but the main problem is the range in our case, since with some of the common lcg parameters we go out of range:

Code:

@init
rnd=1;
@sample
rnd=rnd*16807%2147483647;
// call 0
// 1*16807%2147483647=16807
// call 1
// 16807*16807%2147483647=2824275249
// call 2
// 2824275249*16807%2147483647=1 -> out of range prevention

therefore a larger modulus cannot be used, such as the 2^31-1 -> 'marsenne prime'

so here are some parameters i've calculated based on the "zx-spectrum" multiplier,(77) but with much larger modulus.

Code:

//all with initial seed 1
rnd0*=77; //$x4D
rnd0%=16777219; //$x1000003=2^24+3
rnd1*=76; //$x4C
rnd1%=16777221; //$x1000005=2^24+5
rnd2*=77; //$x4D
rnd2%=8388609; //$x800001=2^23+1

so basically all that is needed is a scaling factor in the end:

Code:

**//full code**
@init
rnd=1;
@sample
rnd*=77;
rnd%=16777219;
spl0=rnd*0.0000000593; //scale down to 0-1 range

comparison with the

marsenne twister implementation in rand().

100 executions of rand() = 61% cpu

100 executions of one recursion lane lcg = 26% cpu

@gfx plots (ermm...30sec freerun)

two rand() x,y:

[img]http://img195.**************/img195/7922/mtwister.png[/img]

two recursions lanes lcg for x,y:

[img]http://img195.**************/img195/5992/gauss.png[/img]

the more expensive gaussian noise version as a bonus:

Code:

**//gaussian version**
@init
rnd0=rnd1=rnd2=1;
@sample
rnd0*=77;
rnd0%=16777219;
rnd1*=76;
rnd1%=16777221;
rnd2*=77;
rnd2%=8388609;
spl0=(rnd0+rnd1+rnd2)*0.00000002403333333333333; //.... *k
//k=(1/3*0.0000000721)

----

lubomir