我正在尝试从Python调用SFMT Mersenne Twister实现(可在 http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/SFMT/中找到)。我这样做是因为我希望能够以 4 种概率从离散 pdf 中快速采样。我正在编写一些模拟，这是迄今为止我代码中最慢的部分。
我已经设法编写了一些有效的 C 代码，并使用使用 SFMT 算法创建的 [0,1] 中的随机数对输入 PDF 进行采样。
但是，我不知道从 Python 调用时如何正确初始化 SFMT 随机数生成器。当然，我只想初始化它一次，然后我需要将用于初始化它的结构的地址(sfmt)传递到我对sfmt_genrand_real1的调用中。
所以一些示例代码将是：
// Create a struct which stores the Mersenne Twister's state
sfmt_t sfmt;
// Initialise the MT with a random seed from the system time
// NOTE: Only want to do this once
int seed = time(NULL);
sfmt_init_gen_rand(&sfmt, seed);
// Get a random number
double random_number = sfmt_genrand_real1(&sfmt);
问题是我不知道在从 Python 调用此代码时如何只播种一次 SFMT 随机数生成器。如果我只是用 C 编写所有内容，我会在 main() 函数中进行初始化，然后将&sfmt参数传递给所有后续sfmt_genrand_real1()调用。
但是因为我正在编译这个 C 代码，然后从 Python 调用它，所以我不能初始化该变量一次。现在，我已经将初始化卡在sfmt_genrand_real1()调用之前，因为这是我获得代码甚至编译并能够访问调用随机数生成器中的sfmt变量的唯一方法。我所有试图以某种方式使sfmt变量"全局"的尝试都适得其反。
所以我的问题是：有没有办法只初始化一次 C SFMT 随机数生成器，然后在我所有后续从 Python 调用c_random_sample中访问用于该初始化的sfmt结构？
提前感谢任何可以帮助解决这个问题的人。
这是我的完整 C 代码。(要编译，您需要将所有 SMFT.c和.h文件放在同一个文件夹中，然后使用python setup.py build_ext --inplace进行编译)
#include "Python.h"
#include <stdio.h>
#include <time.h>
#include "SFMT.h"
static double*
get_double_array(PyObject* data)
{
int i, size;
double* out;
PyObject* seq;
seq = PySequence_Fast(data, "expected a sequence");
if (!seq)
return NULL;
size = PySequence_Size(seq);
if (size < 0)
return NULL;
out = (double*) PyMem_Malloc(size * sizeof(double));
if (!out) {
Py_DECREF(seq);
PyErr_NoMemory();
return NULL;
}
for (i=0; i<size; i++) {
out[i] = PyFloat_AsDouble(PySequence_Fast_GET_ITEM(seq, i));
}
Py_DECREF(seq);
if (PyErr_Occurred()) {
PyMem_Free(out);
out = NULL;
}
return out;
}
static PyObject*
c_random_sample(PyObject* self, PyObject* args)
{
int i;
double* pdf;
PyObject* pdf_in;
if (!PyArg_ParseTuple(args, "O:c_random_sample", &pdf_in))
return NULL;
pdf = get_double_array(pdf_in);
if (!pdf)
return NULL;
// Create a struct which stores the Mersenne Twister's state
sfmt_t sfmt;
int seed = time(NULL);
// Initialise the MT with a random seed from the system time
sfmt_init_gen_rand(&sfmt, seed);
// NOTE: This simply re-initialises the random number generator
// on every call. We need to only initialise it once...
double r = sfmt_genrand_real1(&sfmt);
for (i=0; i<4; i++) {
r -= pdf[i];
if (r < 0) {
break;
}
}
PyMem_Free(pdf);
return PyInt_FromLong(i);
}
static PyMethodDef functions[] = {
{"c_random_sample", c_random_sample, METH_VARARGS},
{NULL, NULL}
};
PyMODINIT_FUNC initc_random_sample(void)
{
Py_InitModule4(
"c_random_sample", functions, "Trying to implement random number sampling in C", NULL, PYTHON_API_VERSION
);
}
sfmt在C语言中处于全局范围：
static sfmt_t sfmt; /* outside any functions */
初始化进入你的模块init函数 - 当模块首次导入 Python 时，它被调用一次：
PyMODINIT_FUNC initc_random_sample(void)
{
/* Py_InitModule4 to initialize the module as before */
int seed = time(NULL);
sfmt_init_gen_rand(&sfmt, seed);
}
这不是你要找的，但如果你不介意依赖关系(主要是 Numpy 和 Cython)，你可能会发现 ng-numpy-randomstate 库很有用，它包含许多快速随机数生成器(特别是 dSFMT)。
如果你能使用它，它将节省你编写C包装器的工作。
此外，这个模块(或至少某些功能？)将来可能会合并到 Numpy 中(参见 Numpy 开放问题 #6967)。