Python code coverage for Modules/_randommodule.c

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