shuuji3.xyz
[READ-ONLY] Mirror of https://github.com/shuuji3/python-extension-twister. Experimental Python C extension to generate random numbers by using Mersenne twister
1つの乱数を生成する twister.random() メソッドを実装。乱数生成に mt19937ar を使用するように修正。
author
TAKAHASHI Shuuji
date
(Aug 5, 2018, 5:34 PM +0900)
commit
078f9f25
078f9f25c4f4edb2bcd9542d313b42034846c14b
parent
ea0198d8
ea0198d8552c8a2d9a156c0ecef5fd33ccc9a084
+489
-114
-8
CMakeLists.txt
-8
CMakeLists.txt
···
1
-
cmake_minimum_required(VERSION 3.12)
2
-
project(python_extension C)
3
-
4
-
set(CMAKE_C_STANDARD 99)
5
-
include_directories(SYSTEM /Users/shuuji/.pyenv/versions/3.6.0/include/python3.6m/)
6
-
include_directories(./SFMT-src-1.5.1)
7
-
8
-
add_executable(python_extension twistermodule.c twister.h twister.c)
+1
Pipfile
+1
Pipfile
+136
Pipfile.lock
+136
Pipfile.lock
···
1
+
{
2
+
"_meta": {
3
+
"hash": {
4
+
"sha256": "844c79ae5c316e32fb91de6ac11b768b695eb20562a5de3e3c54d1694467b917"
5
+
},
6
+
"pipfile-spec": 6,
7
+
"requires": {
8
+
"python_version": "3.6"
9
+
},
10
+
"sources": [
11
+
{
12
+
"name": "pypi",
13
+
"url": "https://pypi.org/simple",
14
+
"verify_ssl": true
15
+
}
16
+
]
17
+
},
18
+
"default": {},
19
+
"develop": {
20
+
"appnope": {
21
+
"hashes": [
22
+
"sha256:5b26757dc6f79a3b7dc9fab95359328d5747fcb2409d331ea66d0272b90ab2a0",
23
+
"sha256:8b995ffe925347a2138d7ac0fe77155e4311a0ea6d6da4f5128fe4b3cbe5ed71"
24
+
],
25
+
"markers": "sys_platform == 'darwin'",
26
+
"version": "==0.1.0"
27
+
},
28
+
"backcall": {
29
+
"hashes": [
30
+
"sha256:38ecd85be2c1e78f77fd91700c76e14667dc21e2713b63876c0eb901196e01e4",
31
+
"sha256:bbbf4b1e5cd2bdb08f915895b51081c041bac22394fdfcfdfbe9f14b77c08bf2"
32
+
],
33
+
"version": "==0.1.0"
34
+
},
35
+
"decorator": {
36
+
"hashes": [
37
+
"sha256:2c51dff8ef3c447388fe5e4453d24a2bf128d3a4c32af3fabef1f01c6851ab82",
38
+
"sha256:c39efa13fbdeb4506c476c9b3babf6a718da943dab7811c206005a4a956c080c"
39
+
],
40
+
"version": "==4.3.0"
41
+
},
42
+
"ipython": {
43
+
"hashes": [
44
+
"sha256:007dcd929c14631f83daff35df0147ea51d1af420da303fd078343878bd5fb62",
45
+
"sha256:b0f2ef9eada4a68ef63ee10b6dde4f35c840035c50fd24265f8052c98947d5a4"
46
+
],
47
+
"index": "pypi",
48
+
"version": "==6.5.0"
49
+
},
50
+
"ipython-genutils": {
51
+
"hashes": [
52
+
"sha256:72dd37233799e619666c9f639a9da83c34013a73e8bbc79a7a6348d93c61fab8",
53
+
"sha256:eb2e116e75ecef9d4d228fdc66af54269afa26ab4463042e33785b887c628ba8"
54
+
],
55
+
"version": "==0.2.0"
56
+
},
57
+
"jedi": {
58
+
"hashes": [
59
+
"sha256:b409ed0f6913a701ed474a614a3bb46e6953639033e31f769ca7581da5bd1ec1",
60
+
"sha256:c254b135fb39ad76e78d4d8f92765ebc9bf92cbc76f49e97ade1d5f5121e1f6f"
61
+
],
62
+
"version": "==0.12.1"
63
+
},
64
+
"parso": {
65
+
"hashes": [
66
+
"sha256:35704a43a3c113cce4de228ddb39aab374b8004f4f2407d070b6a2ca784ce8a2",
67
+
"sha256:895c63e93b94ac1e1690f5fdd40b65f07c8171e3e53cbd7793b5b96c0e0a7f24"
68
+
],
69
+
"version": "==0.3.1"
70
+
},
71
+
"pexpect": {
72
+
"hashes": [
73
+
"sha256:2a8e88259839571d1251d278476f3eec5db26deb73a70be5ed5dc5435e418aba",
74
+
"sha256:3fbd41d4caf27fa4a377bfd16fef87271099463e6fa73e92a52f92dfee5d425b"
75
+
],
76
+
"markers": "sys_platform != 'win32'",
77
+
"version": "==4.6.0"
78
+
},
79
+
"pickleshare": {
80
+
"hashes": [
81
+
"sha256:84a9257227dfdd6fe1b4be1319096c20eb85ff1e82c7932f36efccfe1b09737b",
82
+
"sha256:c9a2541f25aeabc070f12f452e1f2a8eae2abd51e1cd19e8430402bdf4c1d8b5"
83
+
],
84
+
"version": "==0.7.4"
85
+
},
86
+
"prompt-toolkit": {
87
+
"hashes": [
88
+
"sha256:1df952620eccb399c53ebb359cc7d9a8d3a9538cb34c5a1344bdbeb29fbcc381",
89
+
"sha256:3f473ae040ddaa52b52f97f6b4a493cfa9f5920c255a12dc56a7d34397a398a4",
90
+
"sha256:858588f1983ca497f1cf4ffde01d978a3ea02b01c8a26a8bbc5cd2e66d816917"
91
+
],
92
+
"version": "==1.0.15"
93
+
},
94
+
"ptyprocess": {
95
+
"hashes": [
96
+
"sha256:923f299cc5ad920c68f2bc0bc98b75b9f838b93b599941a6b63ddbc2476394c0",
97
+
"sha256:d7cc528d76e76342423ca640335bd3633420dc1366f258cb31d05e865ef5ca1f"
98
+
],
99
+
"version": "==0.6.0"
100
+
},
101
+
"pygments": {
102
+
"hashes": [
103
+
"sha256:78f3f434bcc5d6ee09020f92ba487f95ba50f1e3ef83ae96b9d5ffa1bab25c5d",
104
+
"sha256:dbae1046def0efb574852fab9e90209b23f556367b5a320c0bcb871c77c3e8cc"
105
+
],
106
+
"version": "==2.2.0"
107
+
},
108
+
"simplegeneric": {
109
+
"hashes": [
110
+
"sha256:dc972e06094b9af5b855b3df4a646395e43d1c9d0d39ed345b7393560d0b9173"
111
+
],
112
+
"version": "==0.8.1"
113
+
},
114
+
"six": {
115
+
"hashes": [
116
+
"sha256:70e8a77beed4562e7f14fe23a786b54f6296e34344c23bc42f07b15018ff98e9",
117
+
"sha256:832dc0e10feb1aa2c68dcc57dbb658f1c7e65b9b61af69048abc87a2db00a0eb"
118
+
],
119
+
"version": "==1.11.0"
120
+
},
121
+
"traitlets": {
122
+
"hashes": [
123
+
"sha256:9c4bd2d267b7153df9152698efb1050a5d84982d3384a37b2c1f7723ba3e7835",
124
+
"sha256:c6cb5e6f57c5a9bdaa40fa71ce7b4af30298fbab9ece9815b5d995ab6217c7d9"
125
+
],
126
+
"version": "==4.3.2"
127
+
},
128
+
"wcwidth": {
129
+
"hashes": [
130
+
"sha256:3df37372226d6e63e1b1e1eda15c594bca98a22d33a23832a90998faa96bc65e",
131
+
"sha256:f4ebe71925af7b40a864553f761ed559b43544f8f71746c2d756c7fe788ade7c"
132
+
],
133
+
"version": "==0.1.7"
134
+
}
135
+
}
136
+
}
+3
-2
setup.py
+3
-2
setup.py
···
1
1
from distutils.core import setup, Extension
2
2
3
-
twistermodule = Extension('twister', sources=['twistermodule.c'])
3
+
twistermodule = Extension('twister',
4
+
sources=['src/twistermodule.c', 'src/mt19937ar/mt19937ar.c'])
4
5
5
6
setup(name='twister',
6
-
version='1.0.0',
7
+
version='1.1.0',
7
8
author='TAKAHASHI Shuuji',
8
9
author_email='shuuji3@gmail.com',
9
10
url='https://github.com/shuuji3/python-extension-twister/',
+72
src/mt19937ar/mt19937ar.h
+72
src/mt19937ar/mt19937ar.h
···
1
+
/*
2
+
A C-program for MT19937, with initialization improved 2002/1/26.
3
+
Coded by Takuji Nishimura and Makoto Matsumoto.
4
+
5
+
Before using, initialize the state by using init_genrand(seed)
6
+
or init_by_array(init_key, key_length).
7
+
8
+
Copyright (C) 1997 - 2002, Makoto Matsumoto and Takuji Nishimura,
9
+
All rights reserved.
10
+
Copyright (C) 2005, Mutsuo Saito
11
+
All rights reserved.
12
+
13
+
Redistribution and use in source and binary forms, with or without
14
+
modification, are permitted provided that the following conditions
15
+
are met:
16
+
17
+
1. Redistributions of source code must retain the above copyright
18
+
notice, this list of conditions and the following disclaimer.
19
+
20
+
2. Redistributions in binary form must reproduce the above copyright
21
+
notice, this list of conditions and the following disclaimer in the
22
+
documentation and/or other materials provided with the distribution.
23
+
24
+
3. The names of its contributors may not be used to endorse or promote
25
+
products derived from this software without specific prior written
26
+
permission.
27
+
28
+
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
29
+
"AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
30
+
LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
31
+
A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR
32
+
CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
33
+
EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
34
+
PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
35
+
PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
36
+
LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
37
+
NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
38
+
SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
39
+
40
+
41
+
Any feedback is very welcome.
42
+
http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/emt.html
43
+
email: m-mat @ math.sci.hiroshima-u.ac.jp (remove space)
44
+
*/
45
+
46
+
/* initializes mt[N] with a seed */
47
+
void init_genrand(unsigned long s);
48
+
49
+
/* initialize by an array with array-length */
50
+
/* init_key is the array for initializing keys */
51
+
/* key_length is its length */
52
+
/* slight change for C++, 2004/2/26 */
53
+
void init_by_array(unsigned long init_key[], int key_length);
54
+
55
+
/* generates a random number on [0,0xffffffff]-interval */
56
+
unsigned long genrand_int32(void);
57
+
58
+
/* generates a random number on [0,0x7fffffff]-interval */
59
+
long genrand_int31(void);
60
+
61
+
/* These real versions are due to Isaku Wada, 2002/01/09 added */
62
+
/* generates a random number on [0,1]-real-interval */
63
+
double genrand_real1(void);
64
+
65
+
/* generates a random number on [0,1)-real-interval */
66
+
double genrand_real2(void);
67
+
68
+
/* generates a random number on (0,1)-real-interval */
69
+
double genrand_real3(void);
70
+
71
+
/* generates a random number on [0,1) with 53-bit resolution*/
72
+
double genrand_res53(void);
+65
src/mt19937ar/mtTest.c
+65
src/mt19937ar/mtTest.c
···
1
+
/*
2
+
A C-program for MT19937, with initialization improved 2002/1/26.
3
+
Coded by Takuji Nishimura and Makoto Matsumoto.
4
+
5
+
Before using, initialize the state by using init_genrand(seed)
6
+
or init_by_array(init_key, key_length).
7
+
8
+
Copyright (C) 1997 - 2002, Makoto Matsumoto and Takuji Nishimura,
9
+
All rights reserved.
10
+
Copyright (C) 2005, Mutsuo Saito,
11
+
All rights reserved.
12
+
13
+
Redistribution and use in source and binary forms, with or without
14
+
modification, are permitted provided that the following conditions
15
+
are met:
16
+
17
+
1. Redistributions of source code must retain the above copyright
18
+
notice, this list of conditions and the following disclaimer.
19
+
20
+
2. Redistributions in binary form must reproduce the above copyright
21
+
notice, this list of conditions and the following disclaimer in the
22
+
documentation and/or other materials provided with the distribution.
23
+
24
+
3. The names of its contributors may not be used to endorse or promote
25
+
products derived from this software without specific prior written
26
+
permission.
27
+
28
+
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
29
+
"AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
30
+
LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
31
+
A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR
32
+
CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
33
+
EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
34
+
PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
35
+
PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
36
+
LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
37
+
NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
38
+
SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
39
+
40
+
41
+
Any feedback is very welcome.
42
+
http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/emt.html
43
+
email: m-mat @ math.sci.hiroshima-u.ac.jp (remove space)
44
+
*/
45
+
46
+
#include <stdio.h>
47
+
#include "mt19937ar.h"
48
+
49
+
int main(void)
50
+
{
51
+
int i;
52
+
unsigned long init[4]={0x123, 0x234, 0x345, 0x456}, length=4;
53
+
init_by_array(init, length);
54
+
printf("1000 outputs of genrand_int32()\n");
55
+
for (i=0; i<1000; i++) {
56
+
printf("%10lu ", genrand_int32());
57
+
if (i%5==4) printf("\n");
58
+
}
59
+
printf("\n1000 outputs of genrand_real2()\n");
60
+
for (i=0; i<1000; i++) {
61
+
printf("%10.8f ", genrand_real2());
62
+
if (i%5==4) printf("\n");
63
+
}
64
+
return 0;
65
+
}
+79
src/mt19937ar/readme-mt.txt
+79
src/mt19937ar/readme-mt.txt
···
1
+
This is a Mersenne Twister pseudorandom number generator
2
+
with period 2^19937-1 with improved initialization scheme,
3
+
modified on 2002/1/26 by Takuji Nishimura and Makoto Matsumoto.
4
+
modified on 2005/4/26 by Mutsuo Saito
5
+
6
+
Contents of this tar ball:
7
+
readme-mt.txt this file
8
+
mt19937ar.c the C source (ar: initialize by ARray)
9
+
mt19937ar.h the C header file for mt19937ar
10
+
mtTest.c the C test main program of mt19937ar.c
11
+
mt19937ar.out Test outputs of six types generators. 1000 for each
12
+
13
+
1. Initialization
14
+
The initialization scheme for the previous versions of MT
15
+
(e.g. 1999/10/28 version or earlier) has a tiny problem, that
16
+
the most significant bits of the seed is not well reflected
17
+
to the state vector of MT.
18
+
19
+
This version (2002/1/26) has two initialization schemes:
20
+
init_genrand(seed) and init_by_array(init_key, key_length).
21
+
22
+
init_genrand(seed) initializes the state vector by using
23
+
one unsigned 32-bit integer "seed", which may be zero.
24
+
25
+
init_by_array(init_key, key_length) initializes the state vector
26
+
by using an array init_key[] of unsigned 32-bit integers
27
+
of length key_kength. If key_length is smaller than 624,
28
+
then each array of 32-bit integers gives distinct initial
29
+
state vector. This is useful if you want a larger seed space
30
+
than 32-bit word.
31
+
32
+
2. Generation
33
+
After initialization, the following type of pseudorandom numbers
34
+
are available.
35
+
36
+
genrand_int32() generates unsigned 32-bit integers.
37
+
genrand_int31() generates unsigned 31-bit integers.
38
+
genrand_real1() generates uniform real in [0,1] (32-bit resolution).
39
+
genrand_real2() generates uniform real in [0,1) (32-bit resolution).
40
+
genrand_real3() generates uniform real in (0,1) (32-bit resolution).
41
+
genrand_res53() generates uniform real in [0,1) with 53-bit resolution.
42
+
43
+
Note: the last five functions call the first one.
44
+
if you need more speed for these five functions, you may
45
+
suppress the function call by copying genrand_int32() and
46
+
replacing the last return(), following to these five functions.
47
+
48
+
3. main()
49
+
main() is an example to initialize with an array of length 4,
50
+
then 1000 outputs of unsigned 32-bit integers,
51
+
then 1000 outputs of real [0,1) numbers.
52
+
53
+
4. The outputs
54
+
The output of the mt19937ar.c is in the file mt19937ar.out.
55
+
If you revise or translate the code, check the output
56
+
by using this file.
57
+
58
+
5. Cryptography
59
+
This generator is not cryptoraphically secure.
60
+
You need to use a one-way (or hash) function to obtain
61
+
a secure random sequence.
62
+
63
+
6. Correspondence
64
+
See:
65
+
URL http://www.math.keio.ac.jp/matumoto/emt.html
66
+
email matumoto@math.keio.ac.jp, nisimura@sci.kj.yamagata-u.ac.jp
67
+
68
+
7. Reference
69
+
M. Matsumoto and T. Nishimura,
70
+
"Mersenne Twister: A 623-Dimensionally Equidistributed Uniform
71
+
Pseudo-Random Number Generator",
72
+
ACM Transactions on Modeling and Computer Simulation,
73
+
Vol. 8, No. 1, January 1998, pp 3--30.
74
+
75
+
-------
76
+
Copyright (C) 1997 - 2002, Makoto Matsumoto and Takuji Nishimura,
77
+
All rights reserved.
78
+
Copyright (C) 2005, Mutsuo Saito
79
+
All rights reserved.
+108
src/twistermodule.c
+108
src/twistermodule.c
···
1
+
#include <Python.h>
2
+
3
+
#include <stdio.h>
4
+
#include <time.h>
5
+
#include "mt19937ar/mt19937ar.h"
6
+
7
+
/**
8
+
* 整数の乱数を num 個 print する関数
9
+
* @param num print する乱数の数を指定する
10
+
*/
11
+
void print_rand_int(int num) {
12
+
// Print header
13
+
printf("# print %d random int\n", num);
14
+
fflush(stdout);
15
+
16
+
// Print rand
17
+
for (int i = 0; i < num; i++) {
18
+
printf("%10lu ", genrand_int32());
19
+
if (i % 5 == 4) printf("\n");
20
+
}
21
+
printf("\n");
22
+
fflush(stdout);
23
+
}
24
+
25
+
/**
26
+
* [0, 1) の乱数を num 個 print する関数
27
+
* @param num print する乱数の数を指定する
28
+
*/
29
+
void print_rand(int num) {
30
+
// Print header
31
+
printf("# print %d random numbers\n", num);
32
+
fflush(stdout);
33
+
34
+
// Print rands
35
+
for (int i = 0; i < num; i++) {
36
+
printf("%10.8f ", genrand_real2());
37
+
if (i % 5 == 4) printf("\n");
38
+
}
39
+
printf("\n");
40
+
fflush(stdout);
41
+
}
42
+
43
+
/**
44
+
* print_rand_int() を実行する Python メソッド
45
+
*/
46
+
static PyObject *
47
+
twister_print_rand_int(PyObject *self, PyObject *args) {
48
+
int num = 10;
49
+
if (!PyArg_ParseTuple(args, "|l", &num)) {
50
+
return NULL;
51
+
}
52
+
print_rand_int(num);
53
+
Py_RETURN_NONE;
54
+
}
55
+
56
+
/**
57
+
* print_rand() を実行する Python メソッド
58
+
*/
59
+
static PyObject *
60
+
twister_print_rand(PyObject *self, PyObject *args) {
61
+
int num = 10;
62
+
if (!PyArg_ParseTuple(args, "|l", &num)) {
63
+
return NULL;
64
+
}
65
+
print_rand(num);
66
+
Py_RETURN_NONE;
67
+
}
68
+
69
+
/**
70
+
* 乱数を1つ生成する Python メソッド
71
+
*/
72
+
static PyObject *
73
+
twister_random(PyObject *self, PyObject *args) {
74
+
return PyFloat_FromDouble(genrand_real2());
75
+
}
76
+
77
+
/**
78
+
* 拡張モジュールのメソッドの定義
79
+
*/
80
+
static PyMethodDef TwisterMethods[] = {
81
+
{"print_rand_int", twister_print_rand_int, METH_VARARGS, "Print random int num times."},
82
+
{"print_rand", twister_print_rand, METH_VARARGS, "Print random x in [0, 1) num times."},
83
+
{"random", twister_random, METH_VARARGS, "Generate a random number x in [0, 1)."},
84
+
{NULL, NULL},
85
+
};
86
+
87
+
/**
88
+
* 拡張モジュールの定義
89
+
*/
90
+
static struct PyModuleDef twistermodule = {
91
+
PyModuleDef_HEAD_INIT,
92
+
"twister",
93
+
"Generate random numbers by using Mersenne twister.",
94
+
-1,
95
+
TwisterMethods,
96
+
};
97
+
98
+
/**
99
+
* モジュール読み込み時に行われる初期化処理
100
+
* @return
101
+
*/
102
+
PyMODINIT_FUNC
103
+
PyInit_twister(void) {
104
+
// 現在時刻をシードに与えて乱数を初期化する
105
+
init_genrand((unsigned long) time(NULL));
106
+
107
+
return PyModule_Create(&twistermodule);
108
+
}
+25
-22
twister.c
src/mt19937ar/mt19937ar.c
+25
-22
twister.c
src/mt19937ar/mt19937ar.c
···
1
-
/*
1
+
/*
2
2
A C-program for MT19937, with initialization improved 2002/1/26.
3
3
Coded by Takuji Nishimura and Makoto Matsumoto.
4
4
5
-
Before using, initialize the state by using init_genrand(seed)
5
+
Before using, initialize the state by using init_genrand(seed)
6
6
or init_by_array(init_key, key_length).
7
7
8
8
Copyright (C) 1997 - 2002, Makoto Matsumoto and Takuji Nishimura,
9
-
All rights reserved.
9
+
All rights reserved.
10
+
Copyright (C) 2005, Mutsuo Saito,
11
+
All rights reserved.
10
12
11
13
Redistribution and use in source and binary forms, with or without
12
14
modification, are permitted provided that the following conditions
···
19
21
notice, this list of conditions and the following disclaimer in the
20
22
documentation and/or other materials provided with the distribution.
21
23
22
-
3. The names of its contributors may not be used to endorse or promote
23
-
products derived from this software without specific prior written
24
+
3. The names of its contributors may not be used to endorse or promote
25
+
products derived from this software without specific prior written
24
26
permission.
25
27
26
28
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
···
42
44
*/
43
45
44
46
#include <stdio.h>
47
+
#include "mt19937ar.h"
45
48
46
-
/* Period parameters */
49
+
/* Period parameters */
47
50
#define N 624
48
51
#define M 397
49
52
#define MATRIX_A 0x9908b0dfUL /* constant vector a */
···
58
61
{
59
62
mt[0]= s & 0xffffffffUL;
60
63
for (mti=1; mti<N; mti++) {
61
-
mt[mti] =
62
-
(1812433253UL * (mt[mti-1] ^ (mt[mti-1] >> 30)) + mti);
64
+
mt[mti] =
65
+
(1812433253UL * (mt[mti-1] ^ (mt[mti-1] >> 30)) + mti);
63
66
/* See Knuth TAOCP Vol2. 3rd Ed. P.106 for multiplier. */
64
67
/* In the previous versions, MSBs of the seed affect */
65
68
/* only MSBs of the array mt[]. */
···
73
76
/* init_key is the array for initializing keys */
74
77
/* key_length is its length */
75
78
/* slight change for C++, 2004/2/26 */
76
-
void init_by_array(const unsigned long init_key[], int key_length)
79
+
void init_by_array(unsigned long init_key[], int key_length)
77
80
{
78
81
int i, j, k;
79
82
init_genrand(19650218UL);
···
81
84
k = (N>key_length ? N : key_length);
82
85
for (; k; k--) {
83
86
mt[i] = (mt[i] ^ ((mt[i-1] ^ (mt[i-1] >> 30)) * 1664525UL))
84
-
+ init_key[j] + j; /* non linear */
87
+
+ init_key[j] + j; /* non linear */
85
88
mt[i] &= 0xffffffffUL; /* for WORDSIZE > 32 machines */
86
89
i++; j++;
87
90
if (i>=N) { mt[0] = mt[N-1]; i=1; }
···
89
92
}
90
93
for (k=N-1; k; k--) {
91
94
mt[i] = (mt[i] ^ ((mt[i-1] ^ (mt[i-1] >> 30)) * 1566083941UL))
92
-
- i; /* non linear */
95
+
- i; /* non linear */
93
96
mt[i] &= 0xffffffffUL; /* for WORDSIZE > 32 machines */
94
97
i++;
95
98
if (i>=N) { mt[0] = mt[N-1]; i=1; }
96
99
}
97
100
98
-
mt[0] = 0x80000000UL; /* MSB is 1; assuring non-zero initial array */
101
+
mt[0] = 0x80000000UL; /* MSB is 1; assuring non-zero initial array */
99
102
}
100
103
101
104
/* generates a random number on [0,0xffffffff]-interval */
···
124
127
125
128
mti = 0;
126
129
}
127
-
130
+
128
131
y = mt[mti++];
129
132
130
133
/* Tempering */
···
145
148
/* generates a random number on [0,1]-real-interval */
146
149
double genrand_real1(void)
147
150
{
148
-
return genrand_int32()*(1.0/4294967295.0);
149
-
/* divided by 2^32-1 */
151
+
return genrand_int32()*(1.0/4294967295.0);
152
+
/* divided by 2^32-1 */
150
153
}
151
154
152
155
/* generates a random number on [0,1)-real-interval */
153
156
double genrand_real2(void)
154
157
{
155
-
return genrand_int32()*(1.0/4294967296.0);
158
+
return genrand_int32()*(1.0/4294967296.0);
156
159
/* divided by 2^32 */
157
160
}
158
161
159
162
/* generates a random number on (0,1)-real-interval */
160
163
double genrand_real3(void)
161
164
{
162
-
return (((double)genrand_int32()) + 0.5)*(1.0/4294967296.0);
165
+
return (((double)genrand_int32()) + 0.5)*(1.0/4294967296.0);
163
166
/* divided by 2^32 */
164
167
}
165
168
166
169
/* generates a random number on [0,1) with 53-bit resolution*/
167
-
double genrand_res53(void)
168
-
{
169
-
unsigned long a=genrand_int32()>>5, b=genrand_int32()>>6;
170
-
return(a*67108864.0+b)*(1.0/9007199254740992.0);
171
-
}
170
+
double genrand_res53(void)
171
+
{
172
+
unsigned long a=genrand_int32()>>5, b=genrand_int32()>>6;
173
+
return(a*67108864.0+b)*(1.0/9007199254740992.0);
174
+
}
172
175
/* These real versions are due to Isaku Wada, 2002/01/09 added */
-13
twister.h
-13
twister.h
···
1
-
#ifndef PYTHON_EXTENSION_TWISTER_H
2
-
#define PYTHON_EXTENSION_TWISTER_H
3
-
4
-
void init_genrand(unsigned long s);
5
-
void init_by_array(const unsigned long init_key[], int key_length);
6
-
unsigned long genrand_int32(void);
7
-
long genrand_int31(void);
8
-
double genrand_real1(void);
9
-
double genrand_real2(void);
10
-
double genrand_real3(void);
11
-
double genrand_res53(void);
12
-
13
-
#endif //PYTHON_EXTENSION_TWISTER_H
-69
twistermodule.c
-69
twistermodule.c
···
1
-
#include <Python.h>
2
-
3
-
#include <time.h>
4
-
#include "twister.h"
5
-
6
-
// 整数の乱数を20個 print する関数
7
-
void print_rand_int(int num) {
8
-
printf("# print %d random int\n", num);
9
-
for (int i = 0; i < num; i++) {
10
-
printf("%10lu ", genrand_int32());
11
-
if (i % 5 == 4) printf("\n");
12
-
}
13
-
}
14
-
15
-
// 0 <= r < 1 の乱数を20個 print する関数
16
-
void print_rand(int num) {
17
-
printf("# print %d random number\n", num);
18
-
for (int i = 0; i < num; i++) {
19
-
printf("%10.8f ", genrand_real2());
20
-
if (i % 5 == 4) printf("\n");
21
-
}
22
-
}
23
-
24
-
static PyObject *
25
-
twister_print_rand_int(PyObject *self, PyObject *args) {
26
-
int num = 1;
27
-
if (!PyArg_ParseTuple(args, "|l", &num)) {
28
-
return NULL;
29
-
}
30
-
print_rand_int(num);
31
-
Py_RETURN_NONE;
32
-
}
33
-
34
-
static PyObject *
35
-
twister_print_rand(PyObject *self, PyObject *args) {
36
-
int num = 1;
37
-
if (!PyArg_ParseTuple(args, "|l", &num)) {
38
-
return NULL;
39
-
}
40
-
print_rand(num);
41
-
Py_RETURN_NONE;
42
-
}
43
-
44
-
static PyMethodDef TwisterMethods[] = {
45
-
{"print_rand_int", twister_print_rand_int, METH_VARARGS, "Print random int num times."},
46
-
{"print_rand", twister_print_rand, METH_VARARGS, "Print random r (0 <= r < 1) num times."},
47
-
{NULL, NULL},
48
-
};
49
-
50
-
static struct PyModuleDef twistermodule = {
51
-
PyModuleDef_HEAD_INIT,
52
-
"twister",
53
-
"Generate random numbers by using Mersenne twister.",
54
-
-1,
55
-
TwisterMethods,
56
-
};
57
-
58
-
PyMODINIT_FUNC
59
-
PyInit_twister(void) {
60
-
// モジュール読み込み時に現在時刻で乱数を初期化する
61
-
init_genrand((unsigned long) time(NULL));
62
-
63
-
return PyModule_Create(&twistermodule);
64
-
}
65
-
66
-
int main(void) {
67
-
print_rand(10);
68
-
print_rand_int(20);
69
-
}