• Home
  • Line#
  • Scopes#
  • Navigate#
  • Raw
  • Download
1 /* Random objects */
2 
3 /* ------------------------------------------------------------------
4    The code in this module was based on a download from:
5       http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/MT2002/emt19937ar.html
6 
7    It was modified in 2002 by Raymond Hettinger as follows:
8 
9     * the principal computational lines untouched.
10 
11     * renamed genrand_res53() to random_random() and wrapped
12       in python calling/return code.
13 
14     * genrand_int32() and the helper functions, init_genrand()
15       and init_by_array(), were declared static, wrapped in
16       Python calling/return code.  also, their global data
17       references were replaced with structure references.
18 
19     * unused functions from the original were deleted.
20       new, original C python code was added to implement the
21       Random() interface.
22 
23    The following are the verbatim comments from the original code:
24 
25    A C-program for MT19937, with initialization improved 2002/1/26.
26    Coded by Takuji Nishimura and Makoto Matsumoto.
27 
28    Before using, initialize the state by using init_genrand(seed)
29    or init_by_array(init_key, key_length).
30 
31    Copyright (C) 1997 - 2002, Makoto Matsumoto and Takuji Nishimura,
32    All rights reserved.
33 
34    Redistribution and use in source and binary forms, with or without
35    modification, are permitted provided that the following conditions
36    are met:
37 
38      1. Redistributions of source code must retain the above copyright
39     notice, this list of conditions and the following disclaimer.
40 
41      2. Redistributions in binary form must reproduce the above copyright
42     notice, this list of conditions and the following disclaimer in the
43     documentation and/or other materials provided with the distribution.
44 
45      3. The names of its contributors may not be used to endorse or promote
46     products derived from this software without specific prior written
47     permission.
48 
49    THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
50    "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
51    LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
52    A PARTICULAR PURPOSE ARE DISCLAIMED.  IN NO EVENT SHALL THE COPYRIGHT OWNER OR
53    CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
54    EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
55    PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
56    PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
57    LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
58    NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
59    SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
60 
61 
62    Any feedback is very welcome.
63    http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/emt.html
64    email: m-mat @ math.sci.hiroshima-u.ac.jp (remove space)
65 */
66 
67 /* ---------------------------------------------------------------*/
68 
69 #include "Python.h"
70 #include <time.h>               /* for seeding to current time */
71 #ifdef HAVE_PROCESS_H
72 #  include <process.h>          /* needed for getpid() */
73 #endif
74 
75 /* Period parameters -- These are all magic.  Don't change. */
76 #define N 624
77 #define M 397
78 #define MATRIX_A 0x9908b0dfU    /* constant vector a */
79 #define UPPER_MASK 0x80000000U  /* most significant w-r bits */
80 #define LOWER_MASK 0x7fffffffU  /* least significant r bits */
81 
82 typedef struct {
83     PyObject_HEAD
84     int index;
85     uint32_t state[N];
86 } RandomObject;
87 
88 static PyTypeObject Random_Type;
89 
90 #define RandomObject_Check(v)      (Py_TYPE(v) == &Random_Type)
91 
92 #include "clinic/_randommodule.c.h"
93 
94 /*[clinic input]
95 module _random
96 class _random.Random "RandomObject *" "&Random_Type"
97 [clinic start generated code]*/
98 /*[clinic end generated code: output=da39a3ee5e6b4b0d input=f79898ae7847c321]*/
99 
100 /* Random methods */
101 
102 
103 /* generates a random number on [0,0xffffffff]-interval */
104 static uint32_t
genrand_int32(RandomObject * self)105 genrand_int32(RandomObject *self)
106 {
107     uint32_t y;
108     static const uint32_t mag01[2] = {0x0U, MATRIX_A};
109     /* mag01[x] = x * MATRIX_A  for x=0,1 */
110     uint32_t *mt;
111 
112     mt = self->state;
113     if (self->index >= N) { /* generate N words at one time */
114         int kk;
115 
116         for (kk=0;kk<N-M;kk++) {
117             y = (mt[kk]&UPPER_MASK)|(mt[kk+1]&LOWER_MASK);
118             mt[kk] = mt[kk+M] ^ (y >> 1) ^ mag01[y & 0x1U];
119         }
120         for (;kk<N-1;kk++) {
121             y = (mt[kk]&UPPER_MASK)|(mt[kk+1]&LOWER_MASK);
122             mt[kk] = mt[kk+(M-N)] ^ (y >> 1) ^ mag01[y & 0x1U];
123         }
124         y = (mt[N-1]&UPPER_MASK)|(mt[0]&LOWER_MASK);
125         mt[N-1] = mt[M-1] ^ (y >> 1) ^ mag01[y & 0x1U];
126 
127         self->index = 0;
128     }
129 
130     y = mt[self->index++];
131     y ^= (y >> 11);
132     y ^= (y << 7) & 0x9d2c5680U;
133     y ^= (y << 15) & 0xefc60000U;
134     y ^= (y >> 18);
135     return y;
136 }
137 
138 /* random_random is the function named genrand_res53 in the original code;
139  * generates a random number on [0,1) with 53-bit resolution; note that
140  * 9007199254740992 == 2**53; I assume they're spelling "/2**53" as
141  * multiply-by-reciprocal in the (likely vain) hope that the compiler will
142  * optimize the division away at compile-time.  67108864 is 2**26.  In
143  * effect, a contains 27 random bits shifted left 26, and b fills in the
144  * lower 26 bits of the 53-bit numerator.
145  * The original code credited Isaku Wada for this algorithm, 2002/01/09.
146  */
147 
148 /*[clinic input]
149 _random.Random.random
150 
151   self: self(type="RandomObject *")
152 
153 random() -> x in the interval [0, 1).
154 [clinic start generated code]*/
155 
156 static PyObject *
_random_Random_random_impl(RandomObject * self)157 _random_Random_random_impl(RandomObject *self)
158 /*[clinic end generated code: output=117ff99ee53d755c input=afb2a59cbbb00349]*/
159 {
160     uint32_t a=genrand_int32(self)>>5, b=genrand_int32(self)>>6;
161     return PyFloat_FromDouble((a*67108864.0+b)*(1.0/9007199254740992.0));
162 }
163 
164 /* initializes mt[N] with a seed */
165 static void
init_genrand(RandomObject * self,uint32_t s)166 init_genrand(RandomObject *self, uint32_t s)
167 {
168     int mti;
169     uint32_t *mt;
170 
171     mt = self->state;
172     mt[0]= s;
173     for (mti=1; mti<N; mti++) {
174         mt[mti] =
175         (1812433253U * (mt[mti-1] ^ (mt[mti-1] >> 30)) + mti);
176         /* See Knuth TAOCP Vol2. 3rd Ed. P.106 for multiplier. */
177         /* In the previous versions, MSBs of the seed affect   */
178         /* only MSBs of the array mt[].                                */
179         /* 2002/01/09 modified by Makoto Matsumoto                     */
180     }
181     self->index = mti;
182     return;
183 }
184 
185 /* initialize by an array with array-length */
186 /* init_key is the array for initializing keys */
187 /* key_length is its length */
188 static void
init_by_array(RandomObject * self,uint32_t init_key[],size_t key_length)189 init_by_array(RandomObject *self, uint32_t init_key[], size_t key_length)
190 {
191     size_t i, j, k;       /* was signed in the original code. RDH 12/16/2002 */
192     uint32_t *mt;
193 
194     mt = self->state;
195     init_genrand(self, 19650218U);
196     i=1; j=0;
197     k = (N>key_length ? N : key_length);
198     for (; k; k--) {
199         mt[i] = (mt[i] ^ ((mt[i-1] ^ (mt[i-1] >> 30)) * 1664525U))
200                  + init_key[j] + (uint32_t)j; /* non linear */
201         i++; j++;
202         if (i>=N) { mt[0] = mt[N-1]; i=1; }
203         if (j>=key_length) j=0;
204     }
205     for (k=N-1; k; k--) {
206         mt[i] = (mt[i] ^ ((mt[i-1] ^ (mt[i-1] >> 30)) * 1566083941U))
207                  - (uint32_t)i; /* non linear */
208         i++;
209         if (i>=N) { mt[0] = mt[N-1]; i=1; }
210     }
211 
212     mt[0] = 0x80000000U; /* MSB is 1; assuring non-zero initial array */
213 }
214 
215 /*
216  * The rest is Python-specific code, neither part of, nor derived from, the
217  * Twister download.
218  */
219 
220 static int
random_seed_urandom(RandomObject * self)221 random_seed_urandom(RandomObject *self)
222 {
223     PY_UINT32_T key[N];
224 
225     if (_PyOS_URandomNonblock(key, sizeof(key)) < 0) {
226         return -1;
227     }
228     init_by_array(self, key, Py_ARRAY_LENGTH(key));
229     return 0;
230 }
231 
232 static void
random_seed_time_pid(RandomObject * self)233 random_seed_time_pid(RandomObject *self)
234 {
235     _PyTime_t now;
236     uint32_t key[5];
237 
238     now = _PyTime_GetSystemClock();
239     key[0] = (PY_UINT32_T)(now & 0xffffffffU);
240     key[1] = (PY_UINT32_T)(now >> 32);
241 
242     key[2] = (PY_UINT32_T)getpid();
243 
244     now = _PyTime_GetMonotonicClock();
245     key[3] = (PY_UINT32_T)(now & 0xffffffffU);
246     key[4] = (PY_UINT32_T)(now >> 32);
247 
248     init_by_array(self, key, Py_ARRAY_LENGTH(key));
249 }
250 
251 static PyObject *
random_seed(RandomObject * self,PyObject * arg)252 random_seed(RandomObject *self, PyObject *arg)
253 {
254     PyObject *result = NULL;            /* guilty until proved innocent */
255     PyObject *n = NULL;
256     uint32_t *key = NULL;
257     size_t bits, keyused;
258     int res;
259 
260     if (arg == NULL || arg == Py_None) {
261        if (random_seed_urandom(self) < 0) {
262             PyErr_Clear();
263 
264             /* Reading system entropy failed, fall back on the worst entropy:
265                use the current time and process identifier. */
266             random_seed_time_pid(self);
267         }
268         Py_RETURN_NONE;
269     }
270 
271     /* This algorithm relies on the number being unsigned.
272      * So: if the arg is a PyLong, use its absolute value.
273      * Otherwise use its hash value, cast to unsigned.
274      */
275     if (PyLong_Check(arg)) {
276         /* Calling int.__abs__() prevents calling arg.__abs__(), which might
277            return an invalid value. See issue #31478. */
278         n = PyLong_Type.tp_as_number->nb_absolute(arg);
279     }
280     else {
281         Py_hash_t hash = PyObject_Hash(arg);
282         if (hash == -1)
283             goto Done;
284         n = PyLong_FromSize_t((size_t)hash);
285     }
286     if (n == NULL)
287         goto Done;
288 
289     /* Now split n into 32-bit chunks, from the right. */
290     bits = _PyLong_NumBits(n);
291     if (bits == (size_t)-1 && PyErr_Occurred())
292         goto Done;
293 
294     /* Figure out how many 32-bit chunks this gives us. */
295     keyused = bits == 0 ? 1 : (bits - 1) / 32 + 1;
296 
297     /* Convert seed to byte sequence. */
298     key = (uint32_t *)PyMem_Malloc((size_t)4 * keyused);
299     if (key == NULL) {
300         PyErr_NoMemory();
301         goto Done;
302     }
303     res = _PyLong_AsByteArray((PyLongObject *)n,
304                               (unsigned char *)key, keyused * 4,
305                               PY_LITTLE_ENDIAN,
306                               0); /* unsigned */
307     if (res == -1) {
308         goto Done;
309     }
310 
311 #if PY_BIG_ENDIAN
312     {
313         size_t i, j;
314         /* Reverse an array. */
315         for (i = 0, j = keyused - 1; i < j; i++, j--) {
316             uint32_t tmp = key[i];
317             key[i] = key[j];
318             key[j] = tmp;
319         }
320     }
321 #endif
322     init_by_array(self, key, keyused);
323 
324     Py_INCREF(Py_None);
325     result = Py_None;
326 
327 Done:
328     Py_XDECREF(n);
329     PyMem_Free(key);
330     return result;
331 }
332 
333 /*[clinic input]
334 _random.Random.seed
335 
336   self: self(type="RandomObject *")
337   n: object = None
338   /
339 
340 seed([n]) -> None.
341 
342 Defaults to use urandom and falls back to a combination
343 of the current time and the process identifier.
344 [clinic start generated code]*/
345 
346 static PyObject *
_random_Random_seed_impl(RandomObject * self,PyObject * n)347 _random_Random_seed_impl(RandomObject *self, PyObject *n)
348 /*[clinic end generated code: output=0fad1e16ba883681 input=78d6ef0d52532a54]*/
349 {
350     return random_seed(self, n);
351 }
352 
353 /*[clinic input]
354 _random.Random.getstate
355 
356   self: self(type="RandomObject *")
357 
358 getstate() -> tuple containing the current state.
359 [clinic start generated code]*/
360 
361 static PyObject *
_random_Random_getstate_impl(RandomObject * self)362 _random_Random_getstate_impl(RandomObject *self)
363 /*[clinic end generated code: output=bf6cef0c092c7180 input=b937a487928c0e89]*/
364 {
365     PyObject *state;
366     PyObject *element;
367     int i;
368 
369     state = PyTuple_New(N+1);
370     if (state == NULL)
371         return NULL;
372     for (i=0; i<N ; i++) {
373         element = PyLong_FromUnsignedLong(self->state[i]);
374         if (element == NULL)
375             goto Fail;
376         PyTuple_SET_ITEM(state, i, element);
377     }
378     element = PyLong_FromLong((long)(self->index));
379     if (element == NULL)
380         goto Fail;
381     PyTuple_SET_ITEM(state, i, element);
382     return state;
383 
384 Fail:
385     Py_DECREF(state);
386     return NULL;
387 }
388 
389 
390 /*[clinic input]
391 _random.Random.setstate
392 
393   self: self(type="RandomObject *")
394   state: object
395   /
396 
397 setstate(state) -> None.  Restores generator state.
398 [clinic start generated code]*/
399 
400 static PyObject *
_random_Random_setstate(RandomObject * self,PyObject * state)401 _random_Random_setstate(RandomObject *self, PyObject *state)
402 /*[clinic end generated code: output=fd1c3cd0037b6681 input=b3b4efbb1bc66af8]*/
403 {
404     int i;
405     unsigned long element;
406     long index;
407     uint32_t new_state[N];
408 
409     if (!PyTuple_Check(state)) {
410         PyErr_SetString(PyExc_TypeError,
411             "state vector must be a tuple");
412         return NULL;
413     }
414     if (PyTuple_Size(state) != N+1) {
415         PyErr_SetString(PyExc_ValueError,
416             "state vector is the wrong size");
417         return NULL;
418     }
419 
420     for (i=0; i<N ; i++) {
421         element = PyLong_AsUnsignedLong(PyTuple_GET_ITEM(state, i));
422         if (element == (unsigned long)-1 && PyErr_Occurred())
423             return NULL;
424         new_state[i] = (uint32_t)element;
425     }
426 
427     index = PyLong_AsLong(PyTuple_GET_ITEM(state, i));
428     if (index == -1 && PyErr_Occurred())
429         return NULL;
430     if (index < 0 || index > N) {
431         PyErr_SetString(PyExc_ValueError, "invalid state");
432         return NULL;
433     }
434     self->index = (int)index;
435     for (i = 0; i < N; i++)
436         self->state[i] = new_state[i];
437 
438     Py_RETURN_NONE;
439 }
440 
441 /*[clinic input]
442 
443 _random.Random.getrandbits
444 
445   self: self(type="RandomObject *")
446   k: int
447   /
448 
449 getrandbits(k) -> x.  Generates an int with k random bits.
450 [clinic start generated code]*/
451 
452 static PyObject *
_random_Random_getrandbits_impl(RandomObject * self,int k)453 _random_Random_getrandbits_impl(RandomObject *self, int k)
454 /*[clinic end generated code: output=b402f82a2158887f input=8c0e6396dd176fc0]*/
455 {
456     int i, words;
457     uint32_t r;
458     uint32_t *wordarray;
459     PyObject *result;
460 
461     if (k <= 0) {
462         PyErr_SetString(PyExc_ValueError,
463                         "number of bits must be greater than zero");
464         return NULL;
465     }
466 
467     if (k <= 32)  /* Fast path */
468         return PyLong_FromUnsignedLong(genrand_int32(self) >> (32 - k));
469 
470     words = (k - 1) / 32 + 1;
471     wordarray = (uint32_t *)PyMem_Malloc(words * 4);
472     if (wordarray == NULL) {
473         PyErr_NoMemory();
474         return NULL;
475     }
476 
477     /* Fill-out bits of long integer, by 32-bit words, from least significant
478        to most significant. */
479 #if PY_LITTLE_ENDIAN
480     for (i = 0; i < words; i++, k -= 32)
481 #else
482     for (i = words - 1; i >= 0; i--, k -= 32)
483 #endif
484     {
485         r = genrand_int32(self);
486         if (k < 32)
487             r >>= (32 - k);  /* Drop least significant bits */
488         wordarray[i] = r;
489     }
490 
491     result = _PyLong_FromByteArray((unsigned char *)wordarray, words * 4,
492                                    PY_LITTLE_ENDIAN, 0 /* unsigned */);
493     PyMem_Free(wordarray);
494     return result;
495 }
496 
497 static PyObject *
random_new(PyTypeObject * type,PyObject * args,PyObject * kwds)498 random_new(PyTypeObject *type, PyObject *args, PyObject *kwds)
499 {
500     RandomObject *self;
501     PyObject *tmp;
502 
503     if (type == &Random_Type && !_PyArg_NoKeywords("Random", kwds))
504         return NULL;
505 
506     self = (RandomObject *)type->tp_alloc(type, 0);
507     if (self == NULL)
508         return NULL;
509     tmp = random_seed(self, args);
510     if (tmp == NULL) {
511         Py_DECREF(self);
512         return NULL;
513     }
514     Py_DECREF(tmp);
515     return (PyObject *)self;
516 }
517 
518 static PyMethodDef random_methods[] = {
519     _RANDOM_RANDOM_RANDOM_METHODDEF
520     _RANDOM_RANDOM_SEED_METHODDEF
521     _RANDOM_RANDOM_GETSTATE_METHODDEF
522     _RANDOM_RANDOM_SETSTATE_METHODDEF
523     _RANDOM_RANDOM_GETRANDBITS_METHODDEF
524     {NULL,              NULL}           /* sentinel */
525 };
526 
527 PyDoc_STRVAR(random_doc,
528 "Random() -> create a random number generator with its own internal state.");
529 
530 static PyTypeObject Random_Type = {
531     PyVarObject_HEAD_INIT(NULL, 0)
532     "_random.Random",                   /*tp_name*/
533     sizeof(RandomObject),               /*tp_basicsize*/
534     0,                                  /*tp_itemsize*/
535     /* methods */
536     0,                                  /*tp_dealloc*/
537     0,                                  /*tp_vectorcall_offset*/
538     0,                                  /*tp_getattr*/
539     0,                                  /*tp_setattr*/
540     0,                                  /*tp_as_async*/
541     0,                                  /*tp_repr*/
542     0,                                  /*tp_as_number*/
543     0,                                  /*tp_as_sequence*/
544     0,                                  /*tp_as_mapping*/
545     0,                                  /*tp_hash*/
546     0,                                  /*tp_call*/
547     0,                                  /*tp_str*/
548     PyObject_GenericGetAttr,            /*tp_getattro*/
549     0,                                  /*tp_setattro*/
550     0,                                  /*tp_as_buffer*/
551     Py_TPFLAGS_DEFAULT | Py_TPFLAGS_BASETYPE,           /*tp_flags*/
552     random_doc,                         /*tp_doc*/
553     0,                                  /*tp_traverse*/
554     0,                                  /*tp_clear*/
555     0,                                  /*tp_richcompare*/
556     0,                                  /*tp_weaklistoffset*/
557     0,                                  /*tp_iter*/
558     0,                                  /*tp_iternext*/
559     random_methods,                     /*tp_methods*/
560     0,                                  /*tp_members*/
561     0,                                  /*tp_getset*/
562     0,                                  /*tp_base*/
563     0,                                  /*tp_dict*/
564     0,                                  /*tp_descr_get*/
565     0,                                  /*tp_descr_set*/
566     0,                                  /*tp_dictoffset*/
567     0,                                  /*tp_init*/
568     0,                                  /*tp_alloc*/
569     random_new,                         /*tp_new*/
570     PyObject_Free,                      /*tp_free*/
571     0,                                  /*tp_is_gc*/
572 };
573 
574 PyDoc_STRVAR(module_doc,
575 "Module implements the Mersenne Twister random number generator.");
576 
577 
578 static struct PyModuleDef _randommodule = {
579     PyModuleDef_HEAD_INIT,
580     "_random",
581     module_doc,
582     -1,
583     NULL,
584     NULL,
585     NULL,
586     NULL,
587     NULL
588 };
589 
590 PyMODINIT_FUNC
PyInit__random(void)591 PyInit__random(void)
592 {
593     PyObject *m;
594 
595     if (PyType_Ready(&Random_Type) < 0)
596         return NULL;
597     m = PyModule_Create(&_randommodule);
598     if (m == NULL)
599         return NULL;
600     Py_INCREF(&Random_Type);
601     PyModule_AddObject(m, "Random", (PyObject *)&Random_Type);
602     return m;
603 }
604