In [35]:
%load_ext cythonmagic
import numpy as np
import numpy.random as rn
import matplotlib.pyplot as plt
The cythonmagic extension is already loaded. To reload it, use:
%reload_ext cythonmagic
In [42]:
%%cython -a
from libc.math cimport log
cdef extern from "stdlib.h":
int c_rand "rand" ()
int RAND_MAX
#@cython.cdivision(True)
cpdef double crand():
return c_rand() / (RAND_MAX + 0.0)
#@cython.cdivision(True)
cpdef double cexp(float scale=0.1):
return -log(1-crand())*scale
Out[42]:
Generated by Cython 0.20.1 on Mon Aug 11 12:28:14 2014
1: from libc.math cimport log
2: cdef extern from "stdlib.h":
3: int c_rand "rand" ()
4: int RAND_MAX
5:
6: #@cython.cdivision(True)
7: cpdef double crand():
static PyObject *__pyx_pw_46_cython_magic_df8f9ba600a6940eaea087072e1b1942_1crand(PyObject *__pyx_self, CYTHON_UNUSED PyObject *unused); /*proto*/
static double __pyx_f_46_cython_magic_df8f9ba600a6940eaea087072e1b1942_crand(CYTHON_UNUSED int __pyx_skip_dispatch) {
double __pyx_r;
__Pyx_RefNannyDeclarations
__Pyx_RefNannySetupContext("crand", 0);
/* … */
/* function exit code */
__pyx_L1_error:;
__Pyx_WriteUnraisable("_cython_magic_df8f9ba600a6940eaea087072e1b1942.crand", __pyx_clineno, __pyx_lineno, __pyx_filename, 0);
__pyx_r = 0;
__pyx_L0:;
__Pyx_RefNannyFinishContext();
return __pyx_r;
}
/* Python wrapper */
static PyObject *__pyx_pw_46_cython_magic_df8f9ba600a6940eaea087072e1b1942_1crand(PyObject *__pyx_self, CYTHON_UNUSED PyObject *unused); /*proto*/
static PyObject *__pyx_pw_46_cython_magic_df8f9ba600a6940eaea087072e1b1942_1crand(PyObject *__pyx_self, CYTHON_UNUSED PyObject *unused) {
PyObject *__pyx_r = 0;
__Pyx_RefNannyDeclarations
__Pyx_RefNannySetupContext("crand (wrapper)", 0);
__pyx_r = __pyx_pf_46_cython_magic_df8f9ba600a6940eaea087072e1b1942_crand(__pyx_self);
/* function exit code */
__Pyx_RefNannyFinishContext();
return __pyx_r;
}
static PyObject *__pyx_pf_46_cython_magic_df8f9ba600a6940eaea087072e1b1942_crand(CYTHON_UNUSED PyObject *__pyx_self) {
PyObject *__pyx_r = NULL;
__Pyx_RefNannyDeclarations
__Pyx_RefNannySetupContext("crand", 0);
__Pyx_XDECREF(__pyx_r);
__pyx_t_1 = PyFloat_FromDouble(__pyx_f_46_cython_magic_df8f9ba600a6940eaea087072e1b1942_crand(0)); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 7; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_1);
__pyx_r = __pyx_t_1;
__pyx_t_1 = 0;
goto __pyx_L0;
/* function exit code */
__pyx_L1_error:;
__Pyx_XDECREF(__pyx_t_1);
__Pyx_AddTraceback("_cython_magic_df8f9ba600a6940eaea087072e1b1942.crand", __pyx_clineno, __pyx_lineno, __pyx_filename);
__pyx_r = NULL;
__pyx_L0:;
__Pyx_XGIVEREF(__pyx_r);
__Pyx_RefNannyFinishContext();
return __pyx_r;
}
8: return c_rand() / (RAND_MAX + 0.0)
__pyx_t_1 = rand();
__pyx_t_2 = (RAND_MAX + 0.0);
if (unlikely(__pyx_t_2 == 0)) {
#ifdef WITH_THREAD
PyGILState_STATE __pyx_gilstate_save = PyGILState_Ensure();
#endif
PyErr_SetString(PyExc_ZeroDivisionError, "float division");
#ifdef WITH_THREAD
PyGILState_Release(__pyx_gilstate_save);
#endif
{__pyx_filename = __pyx_f[0]; __pyx_lineno = 8; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
}
__pyx_r = (__pyx_t_1 / __pyx_t_2);
goto __pyx_L0;
9:
10: #@cython.cdivision(True)
11: cpdef double cexp(float scale=0.1):
static PyObject *__pyx_pw_46_cython_magic_df8f9ba600a6940eaea087072e1b1942_3cexp(PyObject *__pyx_self, PyObject *__pyx_args, PyObject *__pyx_kwds); /*proto*/
static double __pyx_f_46_cython_magic_df8f9ba600a6940eaea087072e1b1942_cexp(CYTHON_UNUSED int __pyx_skip_dispatch, struct __pyx_opt_args_46_cython_magic_df8f9ba600a6940eaea087072e1b1942_cexp *__pyx_optional_args) {
float __pyx_v_scale = ((float)0.1);
double __pyx_r;
__Pyx_RefNannyDeclarations
__Pyx_RefNannySetupContext("cexp", 0);
if (__pyx_optional_args) {
if (__pyx_optional_args->__pyx_n > 0) {
__pyx_v_scale = __pyx_optional_args->scale;
}
}
/* … */
/* function exit code */
__pyx_L0:;
__Pyx_RefNannyFinishContext();
return __pyx_r;
}
/* Python wrapper */
static PyObject *__pyx_pw_46_cython_magic_df8f9ba600a6940eaea087072e1b1942_3cexp(PyObject *__pyx_self, PyObject *__pyx_args, PyObject *__pyx_kwds); /*proto*/
static PyObject *__pyx_pw_46_cython_magic_df8f9ba600a6940eaea087072e1b1942_3cexp(PyObject *__pyx_self, PyObject *__pyx_args, PyObject *__pyx_kwds) {
float __pyx_v_scale;
PyObject *__pyx_r = 0;
__Pyx_RefNannyDeclarations
__Pyx_RefNannySetupContext("cexp (wrapper)", 0);
{
static PyObject **__pyx_pyargnames[] = {&__pyx_n_s_scale,0};
PyObject* values[1] = {0};
if (unlikely(__pyx_kwds)) {
Py_ssize_t kw_args;
const Py_ssize_t pos_args = PyTuple_GET_SIZE(__pyx_args);
switch (pos_args) {
case 1: values[0] = PyTuple_GET_ITEM(__pyx_args, 0);
case 0: break;
default: goto __pyx_L5_argtuple_error;
}
kw_args = PyDict_Size(__pyx_kwds);
switch (pos_args) {
case 0:
if (kw_args > 0) {
PyObject* value = PyDict_GetItem(__pyx_kwds, __pyx_n_s_scale);
if (value) { values[0] = value; kw_args--; }
}
}
if (unlikely(kw_args > 0)) {
if (unlikely(__Pyx_ParseOptionalKeywords(__pyx_kwds, __pyx_pyargnames, 0, values, pos_args, "cexp") < 0)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 11; __pyx_clineno = __LINE__; goto __pyx_L3_error;}
}
} else {
switch (PyTuple_GET_SIZE(__pyx_args)) {
case 1: values[0] = PyTuple_GET_ITEM(__pyx_args, 0);
case 0: break;
default: goto __pyx_L5_argtuple_error;
}
}
if (values[0]) {
__pyx_v_scale = __pyx_PyFloat_AsFloat(values[0]); if (unlikely((__pyx_v_scale == (float)-1) && PyErr_Occurred())) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 11; __pyx_clineno = __LINE__; goto __pyx_L3_error;}
} else {
__pyx_v_scale = ((float)0.1);
}
}
goto __pyx_L4_argument_unpacking_done;
__pyx_L5_argtuple_error:;
__Pyx_RaiseArgtupleInvalid("cexp", 0, 0, 1, PyTuple_GET_SIZE(__pyx_args)); {__pyx_filename = __pyx_f[0]; __pyx_lineno = 11; __pyx_clineno = __LINE__; goto __pyx_L3_error;}
__pyx_L3_error:;
__Pyx_AddTraceback("_cython_magic_df8f9ba600a6940eaea087072e1b1942.cexp", __pyx_clineno, __pyx_lineno, __pyx_filename);
__Pyx_RefNannyFinishContext();
return NULL;
__pyx_L4_argument_unpacking_done:;
__pyx_r = __pyx_pf_46_cython_magic_df8f9ba600a6940eaea087072e1b1942_2cexp(__pyx_self, __pyx_v_scale);
int __pyx_lineno = 0;
const char *__pyx_filename = NULL;
int __pyx_clineno = 0;
/* function exit code */
__Pyx_RefNannyFinishContext();
return __pyx_r;
}
static PyObject *__pyx_pf_46_cython_magic_df8f9ba600a6940eaea087072e1b1942_2cexp(CYTHON_UNUSED PyObject *__pyx_self, float __pyx_v_scale) {
PyObject *__pyx_r = NULL;
__Pyx_RefNannyDeclarations
__Pyx_RefNannySetupContext("cexp", 0);
__Pyx_XDECREF(__pyx_r);
__pyx_t_2.__pyx_n = 1;
__pyx_t_2.scale = __pyx_v_scale;
__pyx_t_1 = __pyx_f_46_cython_magic_df8f9ba600a6940eaea087072e1b1942_cexp(0, &__pyx_t_2);
__pyx_t_3 = PyFloat_FromDouble(__pyx_t_1); if (unlikely(!__pyx_t_3)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 11; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_3);
__pyx_r = __pyx_t_3;
__pyx_t_3 = 0;
goto __pyx_L0;
/* function exit code */
__pyx_L1_error:;
__Pyx_XDECREF(__pyx_t_3);
__Pyx_AddTraceback("_cython_magic_df8f9ba600a6940eaea087072e1b1942.cexp", __pyx_clineno, __pyx_lineno, __pyx_filename);
__pyx_r = NULL;
__pyx_L0:;
__Pyx_XGIVEREF(__pyx_r);
__Pyx_RefNannyFinishContext();
return __pyx_r;
}
12: return -log(1-crand())*scale
__pyx_r = ((-log((1.0 - __pyx_f_46_cython_magic_df8f9ba600a6940eaea087072e1b1942_crand(0)))) * __pyx_v_scale);
goto __pyx_L0;
In [50]:
%matplotlib inline
N=10000
value = np.zeros((1,N))
mu = 1./2.
for ii in range(0,N):
value[0,ii] = cexp(1./mu)
h,bins = np.histogram(value,bins=np.arange(0,20,0.1),density=True)
plt.figure()
plt.bar(bins[0:-1],h,width=0.1)
plt.show()
In [36]:
print "check uniform dist."
%timeit rn.random_sample()
%timeit crand()
print
print "check exp distr."
%timeit rn.exponential(scale=1.0)
%timeit cexp(1.0)
check uniform dist.
1000000 loops, best of 3: 869 ns per loop
1000000 loops, best of 3: 440 ns per loop
check exp distr.
100000 loops, best of 3: 1.97 µs per loop
1000000 loops, best of 3: 739 ns per loop
In [41]:
update_time = 0.1
bond = 1
system_time = 1
p1 = 0.03
p2 = 0.02
p = p2/(p1+p2)
mu = p1 + p2
#%timeit update_single_bond(p, mu, bond, update_time)
#%timeit update_bond(p, mu, bond, update_time)
bond, time = update_bond(p, mu, bond, update_time)
print bond
print time
<type 'bool'>
0.181702651229
In [45]:
def update_single_bond(p,mu,bond,update_time):
if bond:
bond = rn.random_sample() > p
else:
bond = rn.random_sample() < 1-p
update_time += rn.exponential(scale=1/mu)
return bond, update_time
In [33]:
%%cython -a
import numpy as np
cimport numpy as np
cimport cython
from libc.math cimport log
cdef extern from "stdlib.h":
int c_rand "rand" ()
int RAND_MAX
#@cython.cdivision(True)
cpdef inline double crand():
return c_rand() / (RAND_MAX + 0.0)
#@cython.cdivision(True)
cpdef inline double cexp(double scale):
return log(1-crand())/(-scale)
@cython.boundscheck(False)
@cython.cdivision(True)
def update_bond(double p, double mu, unsigned int bond, double update_time):
if bond == 1:
bond = <int>(crand() > p)
else:
bond = <int>(crand() < 1-p)
update_time += cexp(1.0/mu)
return bond, update_time
Out[33]:
Generated by Cython 0.20.1 on Mon Aug 11 12:28:14 2014
1: import numpy as np
__pyx_t_1 = __Pyx_Import(__pyx_n_s_numpy, 0, -1); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 1; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_1);
if (PyDict_SetItem(__pyx_d, __pyx_n_s_np, __pyx_t_1) < 0) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 1; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
/* … */
__pyx_t_1 = PyDict_New(); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 1; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_1);
if (PyDict_SetItem(__pyx_d, __pyx_n_s_test, __pyx_t_1) < 0) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 1; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
2: cimport numpy as np
3: cimport cython
4: from libc.math cimport log
5: cdef extern from "stdlib.h":
6: int c_rand "rand" ()
7: int RAND_MAX
8:
9: #@cython.cdivision(True)
10: cpdef inline double crand():
static PyObject *__pyx_pw_46_cython_magic_c0f97b3e9793bd8868e8aa700276015a_1crand(PyObject *__pyx_self, CYTHON_UNUSED PyObject *unused); /*proto*/
static CYTHON_INLINE double __pyx_f_46_cython_magic_c0f97b3e9793bd8868e8aa700276015a_crand(CYTHON_UNUSED int __pyx_skip_dispatch) {
double __pyx_r;
__Pyx_RefNannyDeclarations
__Pyx_RefNannySetupContext("crand", 0);
/* … */
/* function exit code */
__pyx_L1_error:;
__Pyx_WriteUnraisable("_cython_magic_c0f97b3e9793bd8868e8aa700276015a.crand", __pyx_clineno, __pyx_lineno, __pyx_filename, 0);
__pyx_r = 0;
__pyx_L0:;
__Pyx_RefNannyFinishContext();
return __pyx_r;
}
/* Python wrapper */
static PyObject *__pyx_pw_46_cython_magic_c0f97b3e9793bd8868e8aa700276015a_1crand(PyObject *__pyx_self, CYTHON_UNUSED PyObject *unused); /*proto*/
static PyObject *__pyx_pw_46_cython_magic_c0f97b3e9793bd8868e8aa700276015a_1crand(PyObject *__pyx_self, CYTHON_UNUSED PyObject *unused) {
PyObject *__pyx_r = 0;
__Pyx_RefNannyDeclarations
__Pyx_RefNannySetupContext("crand (wrapper)", 0);
__pyx_r = __pyx_pf_46_cython_magic_c0f97b3e9793bd8868e8aa700276015a_crand(__pyx_self);
/* function exit code */
__Pyx_RefNannyFinishContext();
return __pyx_r;
}
static PyObject *__pyx_pf_46_cython_magic_c0f97b3e9793bd8868e8aa700276015a_crand(CYTHON_UNUSED PyObject *__pyx_self) {
PyObject *__pyx_r = NULL;
__Pyx_RefNannyDeclarations
__Pyx_RefNannySetupContext("crand", 0);
__Pyx_XDECREF(__pyx_r);
__pyx_t_1 = PyFloat_FromDouble(__pyx_f_46_cython_magic_c0f97b3e9793bd8868e8aa700276015a_crand(0)); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 10; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_1);
__pyx_r = __pyx_t_1;
__pyx_t_1 = 0;
goto __pyx_L0;
/* function exit code */
__pyx_L1_error:;
__Pyx_XDECREF(__pyx_t_1);
__Pyx_AddTraceback("_cython_magic_c0f97b3e9793bd8868e8aa700276015a.crand", __pyx_clineno, __pyx_lineno, __pyx_filename);
__pyx_r = NULL;
__pyx_L0:;
__Pyx_XGIVEREF(__pyx_r);
__Pyx_RefNannyFinishContext();
return __pyx_r;
}
11: return c_rand() / (RAND_MAX + 0.0)
__pyx_t_1 = rand();
__pyx_t_2 = (RAND_MAX + 0.0);
if (unlikely(__pyx_t_2 == 0)) {
#ifdef WITH_THREAD
PyGILState_STATE __pyx_gilstate_save = PyGILState_Ensure();
#endif
PyErr_SetString(PyExc_ZeroDivisionError, "float division");
#ifdef WITH_THREAD
PyGILState_Release(__pyx_gilstate_save);
#endif
{__pyx_filename = __pyx_f[0]; __pyx_lineno = 11; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
}
__pyx_r = (__pyx_t_1 / __pyx_t_2);
goto __pyx_L0;
12:
13: #@cython.cdivision(True)
14: cpdef inline double cexp(double scale):
static PyObject *__pyx_pw_46_cython_magic_c0f97b3e9793bd8868e8aa700276015a_3cexp(PyObject *__pyx_self, PyObject *__pyx_arg_scale); /*proto*/
static CYTHON_INLINE double __pyx_f_46_cython_magic_c0f97b3e9793bd8868e8aa700276015a_cexp(double __pyx_v_scale, CYTHON_UNUSED int __pyx_skip_dispatch) {
double __pyx_r;
__Pyx_RefNannyDeclarations
__Pyx_RefNannySetupContext("cexp", 0);
/* … */
/* function exit code */
__pyx_L1_error:;
__Pyx_WriteUnraisable("_cython_magic_c0f97b3e9793bd8868e8aa700276015a.cexp", __pyx_clineno, __pyx_lineno, __pyx_filename, 0);
__pyx_r = 0;
__pyx_L0:;
__Pyx_RefNannyFinishContext();
return __pyx_r;
}
/* Python wrapper */
static PyObject *__pyx_pw_46_cython_magic_c0f97b3e9793bd8868e8aa700276015a_3cexp(PyObject *__pyx_self, PyObject *__pyx_arg_scale); /*proto*/
static PyObject *__pyx_pw_46_cython_magic_c0f97b3e9793bd8868e8aa700276015a_3cexp(PyObject *__pyx_self, PyObject *__pyx_arg_scale) {
double __pyx_v_scale;
PyObject *__pyx_r = 0;
__Pyx_RefNannyDeclarations
__Pyx_RefNannySetupContext("cexp (wrapper)", 0);
assert(__pyx_arg_scale); {
__pyx_v_scale = __pyx_PyFloat_AsDouble(__pyx_arg_scale); if (unlikely((__pyx_v_scale == (double)-1) && PyErr_Occurred())) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 14; __pyx_clineno = __LINE__; goto __pyx_L3_error;}
}
goto __pyx_L4_argument_unpacking_done;
__pyx_L3_error:;
__Pyx_AddTraceback("_cython_magic_c0f97b3e9793bd8868e8aa700276015a.cexp", __pyx_clineno, __pyx_lineno, __pyx_filename);
__Pyx_RefNannyFinishContext();
return NULL;
__pyx_L4_argument_unpacking_done:;
__pyx_r = __pyx_pf_46_cython_magic_c0f97b3e9793bd8868e8aa700276015a_2cexp(__pyx_self, ((double)__pyx_v_scale));
int __pyx_lineno = 0;
const char *__pyx_filename = NULL;
int __pyx_clineno = 0;
/* function exit code */
__Pyx_RefNannyFinishContext();
return __pyx_r;
}
static PyObject *__pyx_pf_46_cython_magic_c0f97b3e9793bd8868e8aa700276015a_2cexp(CYTHON_UNUSED PyObject *__pyx_self, double __pyx_v_scale) {
PyObject *__pyx_r = NULL;
__Pyx_RefNannyDeclarations
__Pyx_RefNannySetupContext("cexp", 0);
__Pyx_XDECREF(__pyx_r);
__pyx_t_1 = PyFloat_FromDouble(__pyx_f_46_cython_magic_c0f97b3e9793bd8868e8aa700276015a_cexp(__pyx_v_scale, 0)); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 14; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_1);
__pyx_r = __pyx_t_1;
__pyx_t_1 = 0;
goto __pyx_L0;
/* function exit code */
__pyx_L1_error:;
__Pyx_XDECREF(__pyx_t_1);
__Pyx_AddTraceback("_cython_magic_c0f97b3e9793bd8868e8aa700276015a.cexp", __pyx_clineno, __pyx_lineno, __pyx_filename);
__pyx_r = NULL;
__pyx_L0:;
__Pyx_XGIVEREF(__pyx_r);
__Pyx_RefNannyFinishContext();
return __pyx_r;
}
15: return log(1-crand())/(-scale)
__pyx_t_1 = log((1.0 - __pyx_f_46_cython_magic_c0f97b3e9793bd8868e8aa700276015a_crand(0)));
__pyx_t_2 = (-__pyx_v_scale);
if (unlikely(__pyx_t_2 == 0)) {
#ifdef WITH_THREAD
PyGILState_STATE __pyx_gilstate_save = PyGILState_Ensure();
#endif
PyErr_SetString(PyExc_ZeroDivisionError, "float division");
#ifdef WITH_THREAD
PyGILState_Release(__pyx_gilstate_save);
#endif
{__pyx_filename = __pyx_f[0]; __pyx_lineno = 15; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
}
__pyx_r = (__pyx_t_1 / __pyx_t_2);
goto __pyx_L0;
16:
17: @cython.boundscheck(False)
18: @cython.cdivision(True)
19: def update_bond(double p, double mu, unsigned int bond, double update_time):
/* Python wrapper */
static PyObject *__pyx_pw_46_cython_magic_c0f97b3e9793bd8868e8aa700276015a_5update_bond(PyObject *__pyx_self, PyObject *__pyx_args, PyObject *__pyx_kwds); /*proto*/
static PyMethodDef __pyx_mdef_46_cython_magic_c0f97b3e9793bd8868e8aa700276015a_5update_bond = {__Pyx_NAMESTR("update_bond"), (PyCFunction)__pyx_pw_46_cython_magic_c0f97b3e9793bd8868e8aa700276015a_5update_bond, METH_VARARGS|METH_KEYWORDS, __Pyx_DOCSTR(0)};
static PyObject *__pyx_pw_46_cython_magic_c0f97b3e9793bd8868e8aa700276015a_5update_bond(PyObject *__pyx_self, PyObject *__pyx_args, PyObject *__pyx_kwds) {
double __pyx_v_p;
double __pyx_v_mu;
unsigned int __pyx_v_bond;
double __pyx_v_update_time;
PyObject *__pyx_r = 0;
__Pyx_RefNannyDeclarations
__Pyx_RefNannySetupContext("update_bond (wrapper)", 0);
{
static PyObject **__pyx_pyargnames[] = {&__pyx_n_s_p,&__pyx_n_s_mu,&__pyx_n_s_bond,&__pyx_n_s_update_time,0};
PyObject* values[4] = {0,0,0,0};
if (unlikely(__pyx_kwds)) {
Py_ssize_t kw_args;
const Py_ssize_t pos_args = PyTuple_GET_SIZE(__pyx_args);
switch (pos_args) {
case 4: values[3] = PyTuple_GET_ITEM(__pyx_args, 3);
case 3: values[2] = PyTuple_GET_ITEM(__pyx_args, 2);
case 2: values[1] = PyTuple_GET_ITEM(__pyx_args, 1);
case 1: values[0] = PyTuple_GET_ITEM(__pyx_args, 0);
case 0: break;
default: goto __pyx_L5_argtuple_error;
}
kw_args = PyDict_Size(__pyx_kwds);
switch (pos_args) {
case 0:
if (likely((values[0] = PyDict_GetItem(__pyx_kwds, __pyx_n_s_p)) != 0)) kw_args--;
else goto __pyx_L5_argtuple_error;
case 1:
if (likely((values[1] = PyDict_GetItem(__pyx_kwds, __pyx_n_s_mu)) != 0)) kw_args--;
else {
__Pyx_RaiseArgtupleInvalid("update_bond", 1, 4, 4, 1); {__pyx_filename = __pyx_f[0]; __pyx_lineno = 19; __pyx_clineno = __LINE__; goto __pyx_L3_error;}
}
case 2:
if (likely((values[2] = PyDict_GetItem(__pyx_kwds, __pyx_n_s_bond)) != 0)) kw_args--;
else {
__Pyx_RaiseArgtupleInvalid("update_bond", 1, 4, 4, 2); {__pyx_filename = __pyx_f[0]; __pyx_lineno = 19; __pyx_clineno = __LINE__; goto __pyx_L3_error;}
}
case 3:
if (likely((values[3] = PyDict_GetItem(__pyx_kwds, __pyx_n_s_update_time)) != 0)) kw_args--;
else {
__Pyx_RaiseArgtupleInvalid("update_bond", 1, 4, 4, 3); {__pyx_filename = __pyx_f[0]; __pyx_lineno = 19; __pyx_clineno = __LINE__; goto __pyx_L3_error;}
}
}
if (unlikely(kw_args > 0)) {
if (unlikely(__Pyx_ParseOptionalKeywords(__pyx_kwds, __pyx_pyargnames, 0, values, pos_args, "update_bond") < 0)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 19; __pyx_clineno = __LINE__; goto __pyx_L3_error;}
}
} else if (PyTuple_GET_SIZE(__pyx_args) != 4) {
goto __pyx_L5_argtuple_error;
} else {
values[0] = PyTuple_GET_ITEM(__pyx_args, 0);
values[1] = PyTuple_GET_ITEM(__pyx_args, 1);
values[2] = PyTuple_GET_ITEM(__pyx_args, 2);
values[3] = PyTuple_GET_ITEM(__pyx_args, 3);
}
__pyx_v_p = __pyx_PyFloat_AsDouble(values[0]); if (unlikely((__pyx_v_p == (double)-1) && PyErr_Occurred())) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 19; __pyx_clineno = __LINE__; goto __pyx_L3_error;}
__pyx_v_mu = __pyx_PyFloat_AsDouble(values[1]); if (unlikely((__pyx_v_mu == (double)-1) && PyErr_Occurred())) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 19; __pyx_clineno = __LINE__; goto __pyx_L3_error;}
__pyx_v_bond = __Pyx_PyInt_As_unsigned_int(values[2]); if (unlikely((__pyx_v_bond == (unsigned int)-1) && PyErr_Occurred())) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 19; __pyx_clineno = __LINE__; goto __pyx_L3_error;}
__pyx_v_update_time = __pyx_PyFloat_AsDouble(values[3]); if (unlikely((__pyx_v_update_time == (double)-1) && PyErr_Occurred())) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 19; __pyx_clineno = __LINE__; goto __pyx_L3_error;}
}
goto __pyx_L4_argument_unpacking_done;
__pyx_L5_argtuple_error:;
__Pyx_RaiseArgtupleInvalid("update_bond", 1, 4, 4, PyTuple_GET_SIZE(__pyx_args)); {__pyx_filename = __pyx_f[0]; __pyx_lineno = 19; __pyx_clineno = __LINE__; goto __pyx_L3_error;}
__pyx_L3_error:;
__Pyx_AddTraceback("_cython_magic_c0f97b3e9793bd8868e8aa700276015a.update_bond", __pyx_clineno, __pyx_lineno, __pyx_filename);
__Pyx_RefNannyFinishContext();
return NULL;
__pyx_L4_argument_unpacking_done:;
__pyx_r = __pyx_pf_46_cython_magic_c0f97b3e9793bd8868e8aa700276015a_4update_bond(__pyx_self, __pyx_v_p, __pyx_v_mu, __pyx_v_bond, __pyx_v_update_time);
int __pyx_lineno = 0;
const char *__pyx_filename = NULL;
int __pyx_clineno = 0;
/* function exit code */
__Pyx_RefNannyFinishContext();
return __pyx_r;
}
static PyObject *__pyx_pf_46_cython_magic_c0f97b3e9793bd8868e8aa700276015a_4update_bond(CYTHON_UNUSED PyObject *__pyx_self, double __pyx_v_p, double __pyx_v_mu, unsigned int __pyx_v_bond, double __pyx_v_update_time) {
PyObject *__pyx_r = NULL;
__Pyx_RefNannyDeclarations
__Pyx_RefNannySetupContext("update_bond", 0);
/* … */
/* function exit code */
__pyx_L1_error:;
__Pyx_XDECREF(__pyx_t_2);
__Pyx_XDECREF(__pyx_t_3);
__Pyx_XDECREF(__pyx_t_4);
__Pyx_AddTraceback("_cython_magic_c0f97b3e9793bd8868e8aa700276015a.update_bond", __pyx_clineno, __pyx_lineno, __pyx_filename);
__pyx_r = NULL;
__pyx_L0:;
__Pyx_XGIVEREF(__pyx_r);
__Pyx_RefNannyFinishContext();
return __pyx_r;
}
/* … */
__pyx_tuple__7 = PyTuple_Pack(4, __pyx_n_s_p, __pyx_n_s_mu, __pyx_n_s_bond, __pyx_n_s_update_time); if (unlikely(!__pyx_tuple__7)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 19; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_tuple__7);
__Pyx_GIVEREF(__pyx_tuple__7);
/* … */
__pyx_t_1 = PyCFunction_NewEx(&__pyx_mdef_46_cython_magic_c0f97b3e9793bd8868e8aa700276015a_5update_bond, NULL, __pyx_n_s_cython_magic_c0f97b3e9793bd8868); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 19; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_1);
if (PyDict_SetItem(__pyx_d, __pyx_n_s_update_bond, __pyx_t_1) < 0) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 19; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
20: if bond == 1:
__pyx_t_1 = ((__pyx_v_bond == 1) != 0);
if (__pyx_t_1) {
21: bond = crand() > p
__pyx_v_bond = (__pyx_f_46_cython_magic_c0f97b3e9793bd8868e8aa700276015a_crand(0) > __pyx_v_p);
goto __pyx_L3;
}
/*else*/ {
22: else:
23: bond = crand() < 1-p
__pyx_v_bond = (__pyx_f_46_cython_magic_c0f97b3e9793bd8868e8aa700276015a_crand(0) < (1.0 - __pyx_v_p));
}
__pyx_L3:;
24: update_time += cexp(1.0/mu)
__pyx_v_update_time = (__pyx_v_update_time + __pyx_f_46_cython_magic_c0f97b3e9793bd8868e8aa700276015a_cexp((1.0 / __pyx_v_mu), 0));
25: return bond, update_time
__Pyx_XDECREF(__pyx_r);
__pyx_t_2 = __Pyx_PyInt_From_unsigned_int(__pyx_v_bond); if (unlikely(!__pyx_t_2)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 25; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_2);
__pyx_t_3 = PyFloat_FromDouble(__pyx_v_update_time); if (unlikely(!__pyx_t_3)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 25; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_3);
__pyx_t_4 = PyTuple_New(2); if (unlikely(!__pyx_t_4)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 25; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_4);
PyTuple_SET_ITEM(__pyx_t_4, 0, __pyx_t_2);
__Pyx_GIVEREF(__pyx_t_2);
PyTuple_SET_ITEM(__pyx_t_4, 1, __pyx_t_3);
__Pyx_GIVEREF(__pyx_t_3);
__pyx_t_2 = 0;
__pyx_t_3 = 0;
__pyx_r = __pyx_t_4;
__pyx_t_4 = 0;
goto __pyx_L0;
In [30]:
dim = 1000
p1 = 0.0003
p2 = 0.0005
p = p2/(p1+p2)
mu = p1 + p2
N_runs = 100
bonds = (rn.rand(2,dim,dim) < p).astype(int)
times = rn.exponential(scale=1.0/mu, size=(2,dim,dim))
system_time = 0
In [31]:
%%timeit
global times
global bonds
global system_time
for run in xrange(N_runs):
update_lattice(bonds,times,system_time,p,mu,dim)
system_time +=1
10 loops, best of 3: 172 ms per loop
In [27]:
%%cython
cimport cython
import numpy as np
cimport numpy as np
from libc.math cimport log
cdef extern from "stdlib.h":
int c_rand "rand" ()
int RAND_MAX
@cython.cdivision(True)
cdef inline double crand():
return c_rand() / (RAND_MAX + 0.0)
@cython.cdivision(True)
cdef inline int crandint(int MAX):
return <int>(c_rand() / (RAND_MAX + 0.0) * MAX)
@cython.cdivision(True)
cdef inline double cexp(double scale):
return -log(1-crand())*scale
@cython.boundscheck(False)
@cython.wraparound(False)
@cython.cdivision(True)
cdef void update_bond(long[:,:,::1] bond, double[:,:,::1] update_time,
double system_time, double p, double mu,
np.intp_t ii, np.intp_t jj, np.intp_t kk):
while update_time[ii,jj,kk] < system_time:
if bond[ii,jj,kk] == 1:
bond[ii,jj,kk] = <int>(crand() > p)
else:
bond[ii,jj,kk] = <int>(crand() < 1-p)
update_time[ii,jj,kk] += cexp(1.0/mu)
return
@cython.boundscheck(False)
@cython.wraparound(False)
cpdef update_lattice(long[:,:,::1] bonds, double[:,:,::1] times,
double system_time,
double p, double mu, int dim):
cdef np.intp_t ii,jj,kk
for ii in range(2):
for jj in range(dim):
for kk in range(dim):
if times[ii,jj,kk] < system_time:
update_bond(bonds,times,system_time, p, mu,ii,jj,kk)
return
In [20]:
%%cython -a
# cython: profile=False
cimport cython
import numpy as np
import numpy.random as rn
cimport numpy as np
import pyximport
pyximport.install()
from lattice import update_lattice, crandint
@cython.boundscheck(False)
@cython.wraparound(False)
def runsimulation(long N, long dim, double p, double mu):
cdef double [::1] dist = np.empty((N,),dtype=float)
cdef long [:,:,::1] bonds = (rn.rand(2,dim,dim) < p).astype(int)
cdef double [:,:,::1] times = rn.exponential(scale=1.0/mu, size=(2,dim,dim))
cdef unsigned int x, y, x0, y0, choice, system_time
x = int(dim/2)
y = int(dim/2)
x0 = int(dim/2)
y0 = int(dim/2)
system_time = 0
cdef np.intp_t run
#cdef long [:,::1] neighbours = np.empty((4,2),dtype=long)
for run in xrange(N):
choice = crandint(4)
if choice == 0:
x = <unsigned int>(x + bonds[0,x,y])
elif choice == 1:
x = <unsigned int>(x-bonds[0,x-1,y])
elif choice == 2:
y = <unsigned int>(y-bonds[1,x,y-1])
elif choice == 3:
y = <unsigned int>(y+bonds[1,x,y])
update_lattice(bonds,times,system_time,p,mu,dim)
system_time += 1
dist[run] = (x-x0)*(x-x0) + (y-y0)*(y-y0)
return np.asarray(dist)
Out[20]:
Generated by Cython 0.20.1 on Tue Aug 12 10:30:01 2014
1: # cython: profile=False
__pyx_t_1 = PyDict_New(); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 1; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_1);
if (PyDict_SetItem(__pyx_d, __pyx_n_s_test, __pyx_t_1) < 0) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 1; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
2:
3: cimport cython
4: import numpy as np
__pyx_t_1 = __Pyx_Import(__pyx_n_s_numpy, 0, -1); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 4; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_1);
if (PyDict_SetItem(__pyx_d, __pyx_n_s_np, __pyx_t_1) < 0) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 4; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
5: import numpy.random as rn
__pyx_t_1 = PyList_New(1); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 5; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_1);
__Pyx_INCREF(__pyx_n_s__21);
PyList_SET_ITEM(__pyx_t_1, 0, __pyx_n_s__21);
__Pyx_GIVEREF(__pyx_n_s__21);
__pyx_t_2 = __Pyx_Import(__pyx_n_s_numpy_random, __pyx_t_1, -1); if (unlikely(!__pyx_t_2)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 5; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_2);
__Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
if (PyDict_SetItem(__pyx_d, __pyx_n_s_rn, __pyx_t_2) < 0) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 5; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0;
6: cimport numpy as np
7: import pyximport
__pyx_t_2 = __Pyx_Import(__pyx_n_s_pyximport, 0, -1); if (unlikely(!__pyx_t_2)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 7; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_2);
if (PyDict_SetItem(__pyx_d, __pyx_n_s_pyximport, __pyx_t_2) < 0) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 7; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0;
8: pyximport.install()
__pyx_t_2 = __Pyx_GetModuleGlobalName(__pyx_n_s_pyximport); if (unlikely(!__pyx_t_2)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 8; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_2);
__pyx_t_1 = __Pyx_PyObject_GetAttrStr(__pyx_t_2, __pyx_n_s_install); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 8; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_1);
__Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0;
__pyx_t_2 = __Pyx_PyObject_Call(__pyx_t_1, __pyx_empty_tuple, NULL); if (unlikely(!__pyx_t_2)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 8; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_2);
__Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
__Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0;
9: from lattice import update_lattice, crandint
__pyx_t_2 = PyList_New(2); if (unlikely(!__pyx_t_2)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 9; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_2);
__Pyx_INCREF(__pyx_n_s_update_lattice);
PyList_SET_ITEM(__pyx_t_2, 0, __pyx_n_s_update_lattice);
__Pyx_GIVEREF(__pyx_n_s_update_lattice);
__Pyx_INCREF(__pyx_n_s_crandint);
PyList_SET_ITEM(__pyx_t_2, 1, __pyx_n_s_crandint);
__Pyx_GIVEREF(__pyx_n_s_crandint);
__pyx_t_1 = __Pyx_Import(__pyx_n_s_lattice, __pyx_t_2, -1); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 9; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_1);
__Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0;
__pyx_t_2 = __Pyx_ImportFrom(__pyx_t_1, __pyx_n_s_update_lattice); if (unlikely(!__pyx_t_2)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 9; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_2);
if (PyDict_SetItem(__pyx_d, __pyx_n_s_update_lattice, __pyx_t_2) < 0) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 9; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0;
__pyx_t_2 = __Pyx_ImportFrom(__pyx_t_1, __pyx_n_s_crandint); if (unlikely(!__pyx_t_2)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 9; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_2);
if (PyDict_SetItem(__pyx_d, __pyx_n_s_crandint, __pyx_t_2) < 0) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 9; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0;
__Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
10:
11: @cython.boundscheck(False)
12: @cython.wraparound(False)
13: def runsimulation(long N, long dim, double p, double mu):
/* Python wrapper */
static PyObject *__pyx_pw_46_cython_magic_2cacec797f4cd0125184bd8a61a91229_1runsimulation(PyObject *__pyx_self, PyObject *__pyx_args, PyObject *__pyx_kwds); /*proto*/
static PyMethodDef __pyx_mdef_46_cython_magic_2cacec797f4cd0125184bd8a61a91229_1runsimulation = {__Pyx_NAMESTR("runsimulation"), (PyCFunction)__pyx_pw_46_cython_magic_2cacec797f4cd0125184bd8a61a91229_1runsimulation, METH_VARARGS|METH_KEYWORDS, __Pyx_DOCSTR(0)};
static PyObject *__pyx_pw_46_cython_magic_2cacec797f4cd0125184bd8a61a91229_1runsimulation(PyObject *__pyx_self, PyObject *__pyx_args, PyObject *__pyx_kwds) {
long __pyx_v_N;
long __pyx_v_dim;
double __pyx_v_p;
double __pyx_v_mu;
PyObject *__pyx_r = 0;
__Pyx_RefNannyDeclarations
__Pyx_RefNannySetupContext("runsimulation (wrapper)", 0);
{
static PyObject **__pyx_pyargnames[] = {&__pyx_n_s_N,&__pyx_n_s_dim,&__pyx_n_s_p,&__pyx_n_s_mu,0};
PyObject* values[4] = {0,0,0,0};
if (unlikely(__pyx_kwds)) {
Py_ssize_t kw_args;
const Py_ssize_t pos_args = PyTuple_GET_SIZE(__pyx_args);
switch (pos_args) {
case 4: values[3] = PyTuple_GET_ITEM(__pyx_args, 3);
case 3: values[2] = PyTuple_GET_ITEM(__pyx_args, 2);
case 2: values[1] = PyTuple_GET_ITEM(__pyx_args, 1);
case 1: values[0] = PyTuple_GET_ITEM(__pyx_args, 0);
case 0: break;
default: goto __pyx_L5_argtuple_error;
}
kw_args = PyDict_Size(__pyx_kwds);
switch (pos_args) {
case 0:
if (likely((values[0] = PyDict_GetItem(__pyx_kwds, __pyx_n_s_N)) != 0)) kw_args--;
else goto __pyx_L5_argtuple_error;
case 1:
if (likely((values[1] = PyDict_GetItem(__pyx_kwds, __pyx_n_s_dim)) != 0)) kw_args--;
else {
__Pyx_RaiseArgtupleInvalid("runsimulation", 1, 4, 4, 1); {__pyx_filename = __pyx_f[0]; __pyx_lineno = 13; __pyx_clineno = __LINE__; goto __pyx_L3_error;}
}
case 2:
if (likely((values[2] = PyDict_GetItem(__pyx_kwds, __pyx_n_s_p)) != 0)) kw_args--;
else {
__Pyx_RaiseArgtupleInvalid("runsimulation", 1, 4, 4, 2); {__pyx_filename = __pyx_f[0]; __pyx_lineno = 13; __pyx_clineno = __LINE__; goto __pyx_L3_error;}
}
case 3:
if (likely((values[3] = PyDict_GetItem(__pyx_kwds, __pyx_n_s_mu)) != 0)) kw_args--;
else {
__Pyx_RaiseArgtupleInvalid("runsimulation", 1, 4, 4, 3); {__pyx_filename = __pyx_f[0]; __pyx_lineno = 13; __pyx_clineno = __LINE__; goto __pyx_L3_error;}
}
}
if (unlikely(kw_args > 0)) {
if (unlikely(__Pyx_ParseOptionalKeywords(__pyx_kwds, __pyx_pyargnames, 0, values, pos_args, "runsimulation") < 0)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 13; __pyx_clineno = __LINE__; goto __pyx_L3_error;}
}
} else if (PyTuple_GET_SIZE(__pyx_args) != 4) {
goto __pyx_L5_argtuple_error;
} else {
values[0] = PyTuple_GET_ITEM(__pyx_args, 0);
values[1] = PyTuple_GET_ITEM(__pyx_args, 1);
values[2] = PyTuple_GET_ITEM(__pyx_args, 2);
values[3] = PyTuple_GET_ITEM(__pyx_args, 3);
}
__pyx_v_N = __Pyx_PyInt_As_long(values[0]); if (unlikely((__pyx_v_N == (long)-1) && PyErr_Occurred())) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 13; __pyx_clineno = __LINE__; goto __pyx_L3_error;}
__pyx_v_dim = __Pyx_PyInt_As_long(values[1]); if (unlikely((__pyx_v_dim == (long)-1) && PyErr_Occurred())) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 13; __pyx_clineno = __LINE__; goto __pyx_L3_error;}
__pyx_v_p = __pyx_PyFloat_AsDouble(values[2]); if (unlikely((__pyx_v_p == (double)-1) && PyErr_Occurred())) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 13; __pyx_clineno = __LINE__; goto __pyx_L3_error;}
__pyx_v_mu = __pyx_PyFloat_AsDouble(values[3]); if (unlikely((__pyx_v_mu == (double)-1) && PyErr_Occurred())) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 13; __pyx_clineno = __LINE__; goto __pyx_L3_error;}
}
goto __pyx_L4_argument_unpacking_done;
__pyx_L5_argtuple_error:;
__Pyx_RaiseArgtupleInvalid("runsimulation", 1, 4, 4, PyTuple_GET_SIZE(__pyx_args)); {__pyx_filename = __pyx_f[0]; __pyx_lineno = 13; __pyx_clineno = __LINE__; goto __pyx_L3_error;}
__pyx_L3_error:;
__Pyx_AddTraceback("_cython_magic_2cacec797f4cd0125184bd8a61a91229.runsimulation", __pyx_clineno, __pyx_lineno, __pyx_filename);
__Pyx_RefNannyFinishContext();
return NULL;
__pyx_L4_argument_unpacking_done:;
__pyx_r = __pyx_pf_46_cython_magic_2cacec797f4cd0125184bd8a61a91229_runsimulation(__pyx_self, __pyx_v_N, __pyx_v_dim, __pyx_v_p, __pyx_v_mu);
int __pyx_lineno = 0;
const char *__pyx_filename = NULL;
int __pyx_clineno = 0;
/* function exit code */
__Pyx_RefNannyFinishContext();
return __pyx_r;
}
static PyObject *__pyx_pf_46_cython_magic_2cacec797f4cd0125184bd8a61a91229_runsimulation(CYTHON_UNUSED PyObject *__pyx_self, long __pyx_v_N, long __pyx_v_dim, double __pyx_v_p, double __pyx_v_mu) {
__Pyx_memviewslice __pyx_v_dist = { 0, 0, { 0 }, { 0 }, { 0 } };
__Pyx_memviewslice __pyx_v_bonds = { 0, 0, { 0 }, { 0 }, { 0 } };
__Pyx_memviewslice __pyx_v_times = { 0, 0, { 0 }, { 0 }, { 0 } };
unsigned int __pyx_v_x;
unsigned int __pyx_v_y;
unsigned int __pyx_v_x0;
unsigned int __pyx_v_y0;
unsigned int __pyx_v_choice;
unsigned int __pyx_v_system_time;
__pyx_t_5numpy_intp_t __pyx_v_run;
PyObject *__pyx_r = NULL;
__Pyx_RefNannyDeclarations
__Pyx_RefNannySetupContext("runsimulation", 0);
/* … */
/* function exit code */
__pyx_L1_error:;
__Pyx_XDECREF(__pyx_t_1);
__Pyx_XDECREF(__pyx_t_2);
__Pyx_XDECREF(__pyx_t_3);
__Pyx_XDECREF(__pyx_t_4);
__PYX_XDEC_MEMVIEW(&__pyx_t_5, 1);
__PYX_XDEC_MEMVIEW(&__pyx_t_6, 1);
__Pyx_XDECREF(__pyx_t_7);
__PYX_XDEC_MEMVIEW(&__pyx_t_8, 1);
__Pyx_XDECREF(__pyx_t_23);
__Pyx_XDECREF(__pyx_t_24);
__Pyx_XDECREF(__pyx_t_25);
__Pyx_AddTraceback("_cython_magic_2cacec797f4cd0125184bd8a61a91229.runsimulation", __pyx_clineno, __pyx_lineno, __pyx_filename);
__pyx_r = NULL;
__pyx_L0:;
__PYX_XDEC_MEMVIEW(&__pyx_v_dist, 1);
__PYX_XDEC_MEMVIEW(&__pyx_v_bonds, 1);
__PYX_XDEC_MEMVIEW(&__pyx_v_times, 1);
__Pyx_XGIVEREF(__pyx_r);
__Pyx_RefNannyFinishContext();
return __pyx_r;
}
/* … */
__pyx_tuple__22 = PyTuple_Pack(14, __pyx_n_s_N, __pyx_n_s_dim, __pyx_n_s_p, __pyx_n_s_mu, __pyx_n_s_dist, __pyx_n_s_bonds, __pyx_n_s_times, __pyx_n_s_x, __pyx_n_s_y, __pyx_n_s_x0, __pyx_n_s_y0, __pyx_n_s_choice, __pyx_n_s_system_time, __pyx_n_s_run); if (unlikely(!__pyx_tuple__22)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 13; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_tuple__22);
__Pyx_GIVEREF(__pyx_tuple__22);
/* … */
__pyx_t_1 = PyCFunction_NewEx(&__pyx_mdef_46_cython_magic_2cacec797f4cd0125184bd8a61a91229_1runsimulation, NULL, __pyx_n_s_cython_magic_2cacec797f4cd01251); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 13; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_1);
if (PyDict_SetItem(__pyx_d, __pyx_n_s_runsimulation, __pyx_t_1) < 0) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 13; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
__pyx_codeobj__23 = (PyObject*)__Pyx_PyCode_New(4, 0, 14, 0, 0, __pyx_empty_bytes, __pyx_empty_tuple, __pyx_empty_tuple, __pyx_tuple__22, __pyx_empty_tuple, __pyx_empty_tuple, __pyx_kp_s_users_stud_koher_cache_ipython, __pyx_n_s_runsimulation, 13, __pyx_empty_bytes); if (unlikely(!__pyx_codeobj__23)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 13; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
14:
15: cdef double [::1] dist = np.empty((N,),dtype=float)
__pyx_t_1 = __Pyx_GetModuleGlobalName(__pyx_n_s_np); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 15; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_1);
__pyx_t_2 = __Pyx_PyObject_GetAttrStr(__pyx_t_1, __pyx_n_s_empty); if (unlikely(!__pyx_t_2)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 15; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_2);
__Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
__pyx_t_1 = __Pyx_PyInt_From_long(__pyx_v_N); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 15; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_1);
__pyx_t_3 = PyTuple_New(1); if (unlikely(!__pyx_t_3)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 15; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_3);
PyTuple_SET_ITEM(__pyx_t_3, 0, __pyx_t_1);
__Pyx_GIVEREF(__pyx_t_1);
__pyx_t_1 = 0;
__pyx_t_1 = PyTuple_New(1); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 15; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_1);
PyTuple_SET_ITEM(__pyx_t_1, 0, __pyx_t_3);
__Pyx_GIVEREF(__pyx_t_3);
__pyx_t_3 = 0;
__pyx_t_3 = PyDict_New(); if (unlikely(!__pyx_t_3)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 15; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_3);
if (PyDict_SetItem(__pyx_t_3, __pyx_n_s_dtype, ((PyObject *)((PyObject*)(&PyFloat_Type)))) < 0) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 15; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__pyx_t_4 = __Pyx_PyObject_Call(__pyx_t_2, __pyx_t_1, __pyx_t_3); if (unlikely(!__pyx_t_4)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 15; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_4);
__Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0;
__Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
__Pyx_DECREF(__pyx_t_3); __pyx_t_3 = 0;
__pyx_t_5 = __Pyx_PyObject_to_MemoryviewSlice_dc_double(__pyx_t_4);
if (unlikely(!__pyx_t_5.memview)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 15; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_DECREF(__pyx_t_4); __pyx_t_4 = 0;
__pyx_v_dist = __pyx_t_5;
__pyx_t_5.memview = NULL;
__pyx_t_5.data = NULL;
16: cdef long [:,:,::1] bonds = (rn.rand(2,dim,dim) < p).astype(int)
__pyx_t_4 = __Pyx_GetModuleGlobalName(__pyx_n_s_rn); if (unlikely(!__pyx_t_4)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 16; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_4);
__pyx_t_3 = __Pyx_PyObject_GetAttrStr(__pyx_t_4, __pyx_n_s_rand); if (unlikely(!__pyx_t_3)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 16; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_3);
__Pyx_DECREF(__pyx_t_4); __pyx_t_4 = 0;
__pyx_t_4 = __Pyx_PyInt_From_long(__pyx_v_dim); if (unlikely(!__pyx_t_4)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 16; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_4);
__pyx_t_1 = __Pyx_PyInt_From_long(__pyx_v_dim); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 16; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_1);
__pyx_t_2 = PyTuple_New(3); if (unlikely(!__pyx_t_2)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 16; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_2);
__Pyx_INCREF(__pyx_int_2);
PyTuple_SET_ITEM(__pyx_t_2, 0, __pyx_int_2);
__Pyx_GIVEREF(__pyx_int_2);
PyTuple_SET_ITEM(__pyx_t_2, 1, __pyx_t_4);
__Pyx_GIVEREF(__pyx_t_4);
PyTuple_SET_ITEM(__pyx_t_2, 2, __pyx_t_1);
__Pyx_GIVEREF(__pyx_t_1);
__pyx_t_4 = 0;
__pyx_t_1 = 0;
__pyx_t_1 = __Pyx_PyObject_Call(__pyx_t_3, __pyx_t_2, NULL); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 16; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_1);
__Pyx_DECREF(__pyx_t_3); __pyx_t_3 = 0;
__Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0;
__pyx_t_2 = PyFloat_FromDouble(__pyx_v_p); if (unlikely(!__pyx_t_2)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 16; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_2);
__pyx_t_3 = PyObject_RichCompare(__pyx_t_1, __pyx_t_2, Py_LT); __Pyx_XGOTREF(__pyx_t_3); if (unlikely(!__pyx_t_3)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 16; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
__Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0;
__pyx_t_2 = __Pyx_PyObject_GetAttrStr(__pyx_t_3, __pyx_n_s_astype); if (unlikely(!__pyx_t_2)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 16; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_2);
__Pyx_DECREF(__pyx_t_3); __pyx_t_3 = 0;
__pyx_t_3 = PyTuple_New(1); if (unlikely(!__pyx_t_3)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 16; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_3);
__Pyx_INCREF(((PyObject *)((PyObject*)(&PyInt_Type))));
PyTuple_SET_ITEM(__pyx_t_3, 0, ((PyObject *)((PyObject*)(&PyInt_Type))));
__Pyx_GIVEREF(((PyObject *)((PyObject*)(&PyInt_Type))));
__pyx_t_1 = __Pyx_PyObject_Call(__pyx_t_2, __pyx_t_3, NULL); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 16; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_1);
__Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0;
__Pyx_DECREF(__pyx_t_3); __pyx_t_3 = 0;
__pyx_t_6 = __Pyx_PyObject_to_MemoryviewSlice_d_d_dc_long(__pyx_t_1);
if (unlikely(!__pyx_t_6.memview)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 16; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
__pyx_v_bonds = __pyx_t_6;
__pyx_t_6.memview = NULL;
__pyx_t_6.data = NULL;
17: cdef double [:,:,::1] times = rn.exponential(scale=1.0/mu, size=(2,dim,dim))
__pyx_t_1 = __Pyx_GetModuleGlobalName(__pyx_n_s_rn); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 17; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_1);
__pyx_t_3 = __Pyx_PyObject_GetAttrStr(__pyx_t_1, __pyx_n_s_exponential); if (unlikely(!__pyx_t_3)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 17; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_3);
__Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
__pyx_t_1 = PyDict_New(); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 17; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_1);
if (unlikely(__pyx_v_mu == 0)) {
#ifdef WITH_THREAD
PyGILState_STATE __pyx_gilstate_save = PyGILState_Ensure();
#endif
PyErr_SetString(PyExc_ZeroDivisionError, "float division");
#ifdef WITH_THREAD
PyGILState_Release(__pyx_gilstate_save);
#endif
{__pyx_filename = __pyx_f[0]; __pyx_lineno = 17; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
}
__pyx_t_2 = PyFloat_FromDouble((1.0 / __pyx_v_mu)); if (unlikely(!__pyx_t_2)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 17; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_2);
if (PyDict_SetItem(__pyx_t_1, __pyx_n_s_scale, __pyx_t_2) < 0) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 17; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0;
__pyx_t_2 = __Pyx_PyInt_From_long(__pyx_v_dim); if (unlikely(!__pyx_t_2)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 17; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_2);
__pyx_t_4 = __Pyx_PyInt_From_long(__pyx_v_dim); if (unlikely(!__pyx_t_4)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 17; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_4);
__pyx_t_7 = PyTuple_New(3); if (unlikely(!__pyx_t_7)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 17; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_7);
__Pyx_INCREF(__pyx_int_2);
PyTuple_SET_ITEM(__pyx_t_7, 0, __pyx_int_2);
__Pyx_GIVEREF(__pyx_int_2);
PyTuple_SET_ITEM(__pyx_t_7, 1, __pyx_t_2);
__Pyx_GIVEREF(__pyx_t_2);
PyTuple_SET_ITEM(__pyx_t_7, 2, __pyx_t_4);
__Pyx_GIVEREF(__pyx_t_4);
__pyx_t_2 = 0;
__pyx_t_4 = 0;
if (PyDict_SetItem(__pyx_t_1, __pyx_n_s_size, __pyx_t_7) < 0) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 17; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_DECREF(__pyx_t_7); __pyx_t_7 = 0;
__pyx_t_7 = __Pyx_PyObject_Call(__pyx_t_3, __pyx_empty_tuple, __pyx_t_1); if (unlikely(!__pyx_t_7)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 17; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_7);
__Pyx_DECREF(__pyx_t_3); __pyx_t_3 = 0;
__Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
__pyx_t_8 = __Pyx_PyObject_to_MemoryviewSlice_d_d_dc_double(__pyx_t_7);
if (unlikely(!__pyx_t_8.memview)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 17; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_DECREF(__pyx_t_7); __pyx_t_7 = 0;
__pyx_v_times = __pyx_t_8;
__pyx_t_8.memview = NULL;
__pyx_t_8.data = NULL;
18: cdef unsigned int x, y, x0, y0, choice, system_time
19: x = int(dim/2)
__pyx_v_x = ((unsigned int)__Pyx_div_long(__pyx_v_dim, 2));
20: y = int(dim/2)
__pyx_v_y = ((unsigned int)__Pyx_div_long(__pyx_v_dim, 2));
21: x0 = int(dim/2)
__pyx_v_x0 = ((unsigned int)__Pyx_div_long(__pyx_v_dim, 2));
22: y0 = int(dim/2)
__pyx_v_y0 = ((unsigned int)__Pyx_div_long(__pyx_v_dim, 2));
23: system_time = 0
__pyx_v_system_time = 0;
24: cdef np.intp_t run
25: #cdef long [:,::1] neighbours = np.empty((4,2),dtype=long)
26:
27: for run in xrange(N):
__pyx_t_9 = __pyx_v_N;
for (__pyx_t_10 = 0; __pyx_t_10 < __pyx_t_9; __pyx_t_10+=1) {
__pyx_v_run = __pyx_t_10;
28:
29: choice = crandint(4)
__pyx_t_7 = __Pyx_GetModuleGlobalName(__pyx_n_s_crandint); if (unlikely(!__pyx_t_7)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 29; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_7);
__pyx_t_1 = __Pyx_PyObject_Call(__pyx_t_7, __pyx_tuple_, NULL); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 29; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_1);
__Pyx_DECREF(__pyx_t_7); __pyx_t_7 = 0;
__pyx_t_11 = __Pyx_PyInt_As_unsigned_int(__pyx_t_1); if (unlikely((__pyx_t_11 == (unsigned int)-1) && PyErr_Occurred())) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 29; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
__pyx_v_choice = __pyx_t_11;
/* … */
__pyx_tuple_ = PyTuple_Pack(1, __pyx_int_4); if (unlikely(!__pyx_tuple_)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 29; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_tuple_);
__Pyx_GIVEREF(__pyx_tuple_);
30: if choice == 0:
case 0:
31: x = <unsigned int>(x + bonds[0,x,y])
__pyx_t_12 = 0;
__pyx_t_11 = __pyx_v_x;
__pyx_t_13 = __pyx_v_y;
__pyx_v_x = ((unsigned int)(__pyx_v_x + (*((long *) ( /* dim=2 */ ((char *) (((long *) ( /* dim=1 */ (( /* dim=0 */ (__pyx_v_bonds.data + __pyx_t_12 * __pyx_v_bonds.strides[0]) ) + __pyx_t_11 * __pyx_v_bonds.strides[1]) )) + __pyx_t_13)) )))));
break;
32: elif choice == 1:
case 1:
33: x = <unsigned int>(x-bonds[0,x-1,y])
__pyx_t_14 = 0;
__pyx_t_15 = (__pyx_v_x - 1);
__pyx_t_16 = __pyx_v_y;
__pyx_v_x = ((unsigned int)(__pyx_v_x - (*((long *) ( /* dim=2 */ ((char *) (((long *) ( /* dim=1 */ (( /* dim=0 */ (__pyx_v_bonds.data + __pyx_t_14 * __pyx_v_bonds.strides[0]) ) + __pyx_t_15 * __pyx_v_bonds.strides[1]) )) + __pyx_t_16)) )))));
break;
34: elif choice == 2:
case 2:
35: y = <unsigned int>(y-bonds[1,x,y-1])
__pyx_t_17 = 1;
__pyx_t_18 = __pyx_v_x;
__pyx_t_19 = (__pyx_v_y - 1);
__pyx_v_y = ((unsigned int)(__pyx_v_y - (*((long *) ( /* dim=2 */ ((char *) (((long *) ( /* dim=1 */ (( /* dim=0 */ (__pyx_v_bonds.data + __pyx_t_17 * __pyx_v_bonds.strides[0]) ) + __pyx_t_18 * __pyx_v_bonds.strides[1]) )) + __pyx_t_19)) )))));
break;
36: elif choice == 3:
switch (__pyx_v_choice) {
/* … */
case 3:
37: y = <unsigned int>(y+bonds[1,x,y])
__pyx_t_20 = 1;
__pyx_t_21 = __pyx_v_x;
__pyx_t_22 = __pyx_v_y;
__pyx_v_y = ((unsigned int)(__pyx_v_y + (*((long *) ( /* dim=2 */ ((char *) (((long *) ( /* dim=1 */ (( /* dim=0 */ (__pyx_v_bonds.data + __pyx_t_20 * __pyx_v_bonds.strides[0]) ) + __pyx_t_21 * __pyx_v_bonds.strides[1]) )) + __pyx_t_22)) )))));
break;
default: break;
}
38:
39: update_lattice(bonds,times,system_time,p,mu,dim)
__pyx_t_1 = __Pyx_GetModuleGlobalName(__pyx_n_s_update_lattice); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 39; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_1);
__pyx_t_7 = __pyx_memoryview_fromslice(__pyx_v_bonds, 3, (PyObject *(*)(char *)) __pyx_memview_get_long, (int (*)(char *, PyObject *)) __pyx_memview_set_long, 0);; if (unlikely(!__pyx_t_7)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 39; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_7);
__pyx_t_3 = __pyx_memoryview_fromslice(__pyx_v_times, 3, (PyObject *(*)(char *)) __pyx_memview_get_double, (int (*)(char *, PyObject *)) __pyx_memview_set_double, 0);; if (unlikely(!__pyx_t_3)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 39; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_3);
__pyx_t_4 = __Pyx_PyInt_From_unsigned_int(__pyx_v_system_time); if (unlikely(!__pyx_t_4)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 39; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_4);
__pyx_t_2 = PyFloat_FromDouble(__pyx_v_p); if (unlikely(!__pyx_t_2)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 39; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_2);
__pyx_t_23 = PyFloat_FromDouble(__pyx_v_mu); if (unlikely(!__pyx_t_23)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 39; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_23);
__pyx_t_24 = __Pyx_PyInt_From_long(__pyx_v_dim); if (unlikely(!__pyx_t_24)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 39; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_24);
__pyx_t_25 = PyTuple_New(6); if (unlikely(!__pyx_t_25)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 39; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_25);
PyTuple_SET_ITEM(__pyx_t_25, 0, __pyx_t_7);
__Pyx_GIVEREF(__pyx_t_7);
PyTuple_SET_ITEM(__pyx_t_25, 1, __pyx_t_3);
__Pyx_GIVEREF(__pyx_t_3);
PyTuple_SET_ITEM(__pyx_t_25, 2, __pyx_t_4);
__Pyx_GIVEREF(__pyx_t_4);
PyTuple_SET_ITEM(__pyx_t_25, 3, __pyx_t_2);
__Pyx_GIVEREF(__pyx_t_2);
PyTuple_SET_ITEM(__pyx_t_25, 4, __pyx_t_23);
__Pyx_GIVEREF(__pyx_t_23);
PyTuple_SET_ITEM(__pyx_t_25, 5, __pyx_t_24);
__Pyx_GIVEREF(__pyx_t_24);
__pyx_t_7 = 0;
__pyx_t_3 = 0;
__pyx_t_4 = 0;
__pyx_t_2 = 0;
__pyx_t_23 = 0;
__pyx_t_24 = 0;
__pyx_t_24 = __Pyx_PyObject_Call(__pyx_t_1, __pyx_t_25, NULL); if (unlikely(!__pyx_t_24)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 39; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_24);
__Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
__Pyx_DECREF(__pyx_t_25); __pyx_t_25 = 0;
__Pyx_DECREF(__pyx_t_24); __pyx_t_24 = 0;
40: system_time += 1
__pyx_v_system_time = (__pyx_v_system_time + 1);
41: dist[run] = (x-x0)*(x-x0) + (y-y0)*(y-y0)
__pyx_t_26 = __pyx_v_run;
*((double *) ( /* dim=0 */ ((char *) (((double *) __pyx_v_dist.data) + __pyx_t_26)) )) = (((__pyx_v_x - __pyx_v_x0) * (__pyx_v_x - __pyx_v_x0)) + ((__pyx_v_y - __pyx_v_y0) * (__pyx_v_y - __pyx_v_y0)));
}
42:
43: return np.asarray(dist)
__Pyx_XDECREF(__pyx_r);
__pyx_t_24 = __Pyx_GetModuleGlobalName(__pyx_n_s_np); if (unlikely(!__pyx_t_24)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 43; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_24);
__pyx_t_25 = __Pyx_PyObject_GetAttrStr(__pyx_t_24, __pyx_n_s_asarray); if (unlikely(!__pyx_t_25)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 43; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_25);
__Pyx_DECREF(__pyx_t_24); __pyx_t_24 = 0;
__pyx_t_24 = __pyx_memoryview_fromslice(__pyx_v_dist, 1, (PyObject *(*)(char *)) __pyx_memview_get_double, (int (*)(char *, PyObject *)) __pyx_memview_set_double, 0);; if (unlikely(!__pyx_t_24)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 43; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_24);
__pyx_t_1 = PyTuple_New(1); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 43; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_1);
PyTuple_SET_ITEM(__pyx_t_1, 0, __pyx_t_24);
__Pyx_GIVEREF(__pyx_t_24);
__pyx_t_24 = 0;
__pyx_t_24 = __Pyx_PyObject_Call(__pyx_t_25, __pyx_t_1, NULL); if (unlikely(!__pyx_t_24)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 43; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_24);
__Pyx_DECREF(__pyx_t_25); __pyx_t_25 = 0;
__Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
__pyx_r = __pyx_t_24;
__pyx_t_24 = 0;
goto __pyx_L0;
In [19]:
import pyximport
pyximport.install()
import simulation
simulation.runsimulation(N = 100, dim = 400, p = 0.3, mu = 1./100.)
Out[19]:
array([ 0., 0., 1., 1., 1., 1., 1., 0., 1., 1., 1.,
1., 2., 2., 2., 2., 2., 2., 1., 1., 2., 1.,
1., 0., 1., 1., 1., 1., 1., 1., 1., 1., 1.,
1., 1., 1., 1., 1., 1., 0., 1., 1., 1., 1.,
1., 0., 1., 1., 1., 1., 1., 1., 0., 1., 1.,
0., 1., 0., 0., 1., 2., 1., 2., 1., 1., 4.,
4., 9., 10., 17., 10., 13., 13., 20., 13., 13., 20.,
13., 20., 17., 16., 17., 20., 13., 18., 13., 13., 20.,
13., 13., 13., 13., 13., 20., 17., 16., 16., 16., 17.,
17.])
In [18]:
%timeit runsimulation(N = 100, dim = 400, p = 0.3, mu = 1./100.)
10 loops, best of 3: 78.6 ms per loop
In [ ]:
Content source: andreasko/adlershof-at-work
Similar notebooks: