old.pygsl.c
fdc645e5b8f2master

old.pygsl.c
View Options

			#include <Python.h>
			#include <math.h>
			#include <numpy/libnumarray.h>


			#include <stdio.h>
			#include <gsl/gsl_qrng.h>

			/*********************************/
			/* tests */
			/*********************************/

			static PyObject *
			pygsl_sobol_sequence(self, args)
			PyObject *self;
			PyObject *args;
			{

			PyArrayObject *x;
			int d,n;
			int dim[2];

			int i,j;


			/* parse arguments */
			if (!PyArg_ParseTuple(args, "ii", &n ,&d))
			return NULL;

			if (d > 40){
			PyErr_SetString(PyExc_ValueError,"second argument must not be larger than 40.");
			return NULL;
			}

			/* allocate */
			gsl_qrng * q = gsl_qrng_alloc (gsl_qrng_sobol, d);

			/* create output */
			dim[0]=n;
			dim[1]=d;
			x = (PyArrayObject *) PyArray_FromDims(2,dim,PyArray_DOUBLE);


			/* compute sequence */

			for (i = 0; i < n; i++)
			{
			double v[d];
			gsl_qrng_get (q, v);

			for (j = 0; j < d; j++)
			{
			(double )(x->data + i(x->strides[0]) + j(x->strides[1])) = v[j];
			}

			}

			gsl_qrng_free (q);












			return PyArray_Return(x);
			}



			/* definition of the method table */

			static PyMethodDef pygslMethods[] = {

			{"sobol_sequence", pygsl_sobol_sequence, METH_VARARGS,
			"Return a sobol sequence of len n and of dimension d."},

			{NULL, NULL, 0, NULL} /* Sentinel */
			};





			void initpygsl(void)
			{
			(void) Py_InitModule("pygsl", pygslMethods);

			import_array();
			}