ocean.math.random.engines.KISS source code

1 /*******************************************************************************
2 
3         Copyright:
4             Copyright (c) 2008. Fawzi Mohamed
5             Some parts copyright (c) 2009-2016 dunnhumby Germany GmbH.
6             All rights reserved.
7 
8         License:
9             Tango Dual License: 3-Clause BSD License / Academic Free License v3.0.
10             See LICENSE_TANGO.txt for details.
11 
12         Version: Initial release: July 2008
13 
14         Authors: Fawzi Mohamed
15 
16 *******************************************************************************/
17 module ocean.math.random.engines.KISS;
18 import ocean.transition;
19 import Integer = ocean.text.convert.Integer_tango;
20 import ocean.core.Verify;
21 
22 /+ Kiss99 random number generator, by Marisaglia
23 + a simple RNG that passes all statistical tests
24 + This is the engine, *never* use it directly, always use it though a RandomG class
25 +/
26 struct Kiss99{
27     private uint kiss_x = 123456789;
28     private uint kiss_y = 362436000;
29     private uint kiss_z = 521288629;
30     private uint kiss_c = 7654321;
31     private uint nBytes = 0;
32     private uint restB  = 0;
33 
34     enum int canCheckpoint=true;
35     enum int canSeed=true;
36 
37     void skip(uint n){
38         for (int i=n;i!=n;--i){
39             next;
40         }
41     }
42     ubyte nextB(){
43         if (nBytes>0) {
44             ubyte res=cast(ubyte)(restB & 0xFF);
45             restB >>= 8;
46             --nBytes;
47             return res;
48         } else {
49             restB=next;
50             ubyte res=cast(ubyte)(restB & 0xFF);
51             restB >>= 8;
52             nBytes=3;
53             return res;
54         }
55     }
56     uint next(){
57         enum ulong a = 698769069UL;
58         ulong t;
59         kiss_x = 69069*kiss_x+12345;
60         kiss_y ^= (kiss_y<<13); kiss_y ^= (kiss_y>>17); kiss_y ^= (kiss_y<<5);
61         t = a*kiss_z+kiss_c; kiss_c = cast(uint)(t>>32);
62         kiss_z=cast(uint)t;
63         return kiss_x+kiss_y+kiss_z;
64     }
65     ulong nextL(){
66         return ((cast(ulong)next)<<32)+cast(ulong)next;
67     }
68 
69     void seed(scope uint delegate() r){
70         kiss_x = r();
71         for (int i=0;i<100;++i){
72             kiss_y=r();
73             if (kiss_y!=0) break;
74         }
75         if (kiss_y==0) kiss_y=362436000;
76         kiss_z=r();
77         /* Don’t really need to seed c as well (is reset after a next),
78            but doing it allows to completely restore a given internal state */
79         kiss_c = r() % 698769069; /* Should be less than 698769069 */
80         nBytes = 0;
81         restB=0;
82     }
83     /// writes the current status in a string
84     istring toString(){
85         char[] res=new char[6+6*9];
86         int i=0;
87         res[i..i+6]="KISS99";
88         i+=6;
89         foreach (val;[kiss_x,kiss_y,kiss_z,kiss_c,nBytes,restB]){
90             res[i]='_';
91             ++i;
92             Integer.format(res[i..i+8],val,"x8");
93             i+=8;
94         }
95         verify(i==res.length,"unexpected size");
96         return assumeUnique(res);
97     }
98     /// reads the current status from a string (that should have been trimmed)
99     /// returns the number of chars read
100     size_t fromString(cstring s){
101         size_t i=0;
102         verify(s[i..i+4]=="KISS","unexpected kind, expected KISS");
103         verify(s[i+4..i+7]=="99_","unexpected version, expected 99");
104         i+=6;
105         foreach (val;[&kiss_x,&kiss_y,&kiss_z,&kiss_c,&nBytes,&restB]){
106             verify(s[i]=='_',"no separator _ found");
107             ++i;
108             uint ate;
109             *val=cast(uint)Integer.convert(s[i..i+8],16,&ate);
110             verify(ate==8,"unexpected read size");
111             i+=8;
112         }
113         return i;
114     }
115 }