Cython¶

Christian Elsasser (che@physik.uzh.ch) & Nicola Chiapolini (nicola.chiapolini@mnf.uzh.ch)

The content of this lecture might be distributed under CC by-sa.

An introductionnary example - Fibonacci numbers¶

Pure Python¶

In [1]:
def fib0(n):
    a,b=1,1
    for i in range(n):
        a,b=a+b,a
    return a
In [2]:
fib0(5)
Out[2]:
13
In [3]:
fib0(6)
Out[3]:
21
In [4]:
%timeit fib0(10000)

1.38 ms ± 13.1 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

Cython: Compile unmodified¶

In [5]:
import cython
In [6]:
%load_ext Cython
In [7]:
%%cython
def helloWorld():
    print("Hello, World!")
In [8]:
helloWorld
Out[8]:
<function _cython_magic_32b5dd934daf90290d8e49fa0caafcbb.helloWorld>
In [9]:
helloWorld()

Hello, World!
In [10]:
%%cython
def fib1(n):
    a,b=1,1
    for i in range(n):
        a,b=a+b,a
    return a
In [11]:
fib1
Out[11]:
<function _cython_magic_64ea3364e12ecfe0d7f2970b838ef6c8.fib1>
In [12]:
fib0
Out[12]:
<function __main__.fib0(n)>
In [13]:
fib1(5)
Out[13]:
13
In [14]:
%timeit fib1(10000)

1.26 ms ± 3.4 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

Cython: Declaring Types¶

In [15]:
%%cython
def fib2(int n):
    cdef int a,b,i
    a,b=1,1
    for i in range(n):
        a,b=a+b,a
    return a
In [16]:
fib2
Out[16]:
<function _cython_magic_87b7995341b8f82b8b733bd8f63a944e.fib2>
In [17]:
fib2(5)
Out[17]:
13
In [18]:
%timeit fib2(10000)

5.6 µs ± 76.1 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)

Impressive, but¶

In [19]:
print("python result: ", fib0(50))
print("cython result: ", fib2(50))

python result:  32951280099
cython result:  -1408458269
In [20]:
print("python result: ", fib0(40))
print("cython result: ", fib2(40))

python result:  267914296
cython result:  267914296
In [21]:
%timeit fib0(40)
%timeit fib2(40)

1.83 µs ± 13.8 ns per loop (mean ± std. dev. of 7 runs, 1000000 loops each)
68.5 ns ± 0.171 ns per loop (mean ± std. dev. of 7 runs, 10000000 loops each)

Conclusion: The speedup is real, but we pay with new potential errors.

Second Example - Numerical integration¶

start with declaring types¶

In [22]:
%%cython -a

from math import sin,exp

def f0(double x):
    return sin(x)*exp(-x)

def integrate0(double a,double b,int N):
    cdef double dx,x,s
    cdef int _
    dx = (b-a)/N
    x = a+dx/2
    s = 0.0
    for _ in range(N):
        s+=f0(x)
        x += dx
    return s*dx
Out[22]:
Cython: _cython_magic_f826ce6229b50973a9ef92e6c0ee0ea0.pyx

Generated by Cython 0.29.30

Yellow lines hint at Python interaction.
Click on a line that starts with a "+" to see the C code that Cython generated for it. 

 01: 
+02: from math import sin,exp
  __pyx_t_1 = PyList_New(2); if (unlikely(!__pyx_t_1)) __PYX_ERR(0, 2, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_t_1);
  __Pyx_INCREF(__pyx_n_s_sin);
  __Pyx_GIVEREF(__pyx_n_s_sin);
  PyList_SET_ITEM(__pyx_t_1, 0, __pyx_n_s_sin);
  __Pyx_INCREF(__pyx_n_s_exp);
  __Pyx_GIVEREF(__pyx_n_s_exp);
  PyList_SET_ITEM(__pyx_t_1, 1, __pyx_n_s_exp);
  __pyx_t_2 = __Pyx_Import(__pyx_n_s_math, __pyx_t_1, 0); if (unlikely(!__pyx_t_2)) __PYX_ERR(0, 2, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_t_2);
  __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
  __pyx_t_1 = __Pyx_ImportFrom(__pyx_t_2, __pyx_n_s_sin); if (unlikely(!__pyx_t_1)) __PYX_ERR(0, 2, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_t_1);
  if (PyDict_SetItem(__pyx_d, __pyx_n_s_sin, __pyx_t_1) < 0) __PYX_ERR(0, 2, __pyx_L1_error)
  __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
  __pyx_t_1 = __Pyx_ImportFrom(__pyx_t_2, __pyx_n_s_exp); if (unlikely(!__pyx_t_1)) __PYX_ERR(0, 2, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_t_1);
  if (PyDict_SetItem(__pyx_d, __pyx_n_s_exp, __pyx_t_1) < 0) __PYX_ERR(0, 2, __pyx_L1_error)
  __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
  __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0;

 03: 
+04: def f0(double x):
/* Python wrapper */
static PyObject *__pyx_pw_46_cython_magic_f826ce6229b50973a9ef92e6c0ee0ea0_1f0(PyObject *__pyx_self, PyObject *__pyx_arg_x); /*proto*/
static PyMethodDef __pyx_mdef_46_cython_magic_f826ce6229b50973a9ef92e6c0ee0ea0_1f0 = {"f0", (PyCFunction)__pyx_pw_46_cython_magic_f826ce6229b50973a9ef92e6c0ee0ea0_1f0, METH_O, 0};
static PyObject *__pyx_pw_46_cython_magic_f826ce6229b50973a9ef92e6c0ee0ea0_1f0(PyObject *__pyx_self, PyObject *__pyx_arg_x) {
  double __pyx_v_x;
  PyObject *__pyx_r = 0;
  __Pyx_RefNannyDeclarations
  __Pyx_RefNannySetupContext("f0 (wrapper)", 0);
  assert(__pyx_arg_x); {
    __pyx_v_x = __pyx_PyFloat_AsDouble(__pyx_arg_x); if (unlikely((__pyx_v_x == (double)-1) && PyErr_Occurred())) __PYX_ERR(0, 4, __pyx_L3_error)
  }
  goto __pyx_L4_argument_unpacking_done;
  __pyx_L3_error:;
  __Pyx_AddTraceback("_cython_magic_f826ce6229b50973a9ef92e6c0ee0ea0.f0", __pyx_clineno, __pyx_lineno, __pyx_filename);
  __Pyx_RefNannyFinishContext();
  return NULL;
  __pyx_L4_argument_unpacking_done:;
  __pyx_r = __pyx_pf_46_cython_magic_f826ce6229b50973a9ef92e6c0ee0ea0_f0(__pyx_self, ((double)__pyx_v_x));
  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_f826ce6229b50973a9ef92e6c0ee0ea0_f0(CYTHON_UNUSED PyObject *__pyx_self, double __pyx_v_x) {
  PyObject *__pyx_r = NULL;
  __Pyx_RefNannyDeclarations
  __Pyx_RefNannySetupContext("f0", 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_XDECREF(__pyx_t_5);
  __Pyx_AddTraceback("_cython_magic_f826ce6229b50973a9ef92e6c0ee0ea0.f0", __pyx_clineno, __pyx_lineno, __pyx_filename);
  __pyx_r = NULL;
  __pyx_L0:;
  __Pyx_XGIVEREF(__pyx_r);
  __Pyx_RefNannyFinishContext();
  return __pyx_r;
}
/* … */
  __pyx_tuple_ = PyTuple_Pack(2, __pyx_n_s_x, __pyx_n_s_x); if (unlikely(!__pyx_tuple_)) __PYX_ERR(0, 4, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_tuple_);
  __Pyx_GIVEREF(__pyx_tuple_);
/* … */
  __pyx_t_2 = PyCFunction_NewEx(&__pyx_mdef_46_cython_magic_f826ce6229b50973a9ef92e6c0ee0ea0_1f0, NULL, __pyx_n_s_cython_magic_f826ce6229b50973a9); if (unlikely(!__pyx_t_2)) __PYX_ERR(0, 4, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_t_2);
  if (PyDict_SetItem(__pyx_d, __pyx_n_s_f0, __pyx_t_2) < 0) __PYX_ERR(0, 4, __pyx_L1_error)
  __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0;
  __pyx_codeobj__2 = (PyObject*)__Pyx_PyCode_New(1, 0, 2, 0, CO_OPTIMIZED|CO_NEWLOCALS, __pyx_empty_bytes, __pyx_empty_tuple, __pyx_empty_tuple, __pyx_tuple_, __pyx_empty_tuple, __pyx_empty_tuple, __pyx_kp_s_home_nchiapol_cache_ipython_cyt, __pyx_n_s_f0, 4, __pyx_empty_bytes); if (unlikely(!__pyx_codeobj__2)) __PYX_ERR(0, 4, __pyx_L1_error)

+05:     return sin(x)*exp(-x)
  __Pyx_XDECREF(__pyx_r);
  __Pyx_GetModuleGlobalName(__pyx_t_2, __pyx_n_s_sin); if (unlikely(!__pyx_t_2)) __PYX_ERR(0, 5, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_t_2);
  __pyx_t_3 = PyFloat_FromDouble(__pyx_v_x); if (unlikely(!__pyx_t_3)) __PYX_ERR(0, 5, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_t_3);
  __pyx_t_4 = NULL;
  if (CYTHON_UNPACK_METHODS && unlikely(PyMethod_Check(__pyx_t_2))) {
    __pyx_t_4 = PyMethod_GET_SELF(__pyx_t_2);
    if (likely(__pyx_t_4)) {
      PyObject* function = PyMethod_GET_FUNCTION(__pyx_t_2);
      __Pyx_INCREF(__pyx_t_4);
      __Pyx_INCREF(function);
      __Pyx_DECREF_SET(__pyx_t_2, function);
    }
  }
  __pyx_t_1 = (__pyx_t_4) ? __Pyx_PyObject_Call2Args(__pyx_t_2, __pyx_t_4, __pyx_t_3) : __Pyx_PyObject_CallOneArg(__pyx_t_2, __pyx_t_3);
  __Pyx_XDECREF(__pyx_t_4); __pyx_t_4 = 0;
  __Pyx_DECREF(__pyx_t_3); __pyx_t_3 = 0;
  if (unlikely(!__pyx_t_1)) __PYX_ERR(0, 5, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_t_1);
  __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0;
  __Pyx_GetModuleGlobalName(__pyx_t_3, __pyx_n_s_exp); if (unlikely(!__pyx_t_3)) __PYX_ERR(0, 5, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_t_3);
  __pyx_t_4 = PyFloat_FromDouble((-__pyx_v_x)); if (unlikely(!__pyx_t_4)) __PYX_ERR(0, 5, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_t_4);
  __pyx_t_5 = NULL;
  if (CYTHON_UNPACK_METHODS && unlikely(PyMethod_Check(__pyx_t_3))) {
    __pyx_t_5 = PyMethod_GET_SELF(__pyx_t_3);
    if (likely(__pyx_t_5)) {
      PyObject* function = PyMethod_GET_FUNCTION(__pyx_t_3);
      __Pyx_INCREF(__pyx_t_5);
      __Pyx_INCREF(function);
      __Pyx_DECREF_SET(__pyx_t_3, function);
    }
  }
  __pyx_t_2 = (__pyx_t_5) ? __Pyx_PyObject_Call2Args(__pyx_t_3, __pyx_t_5, __pyx_t_4) : __Pyx_PyObject_CallOneArg(__pyx_t_3, __pyx_t_4);
  __Pyx_XDECREF(__pyx_t_5); __pyx_t_5 = 0;
  __Pyx_DECREF(__pyx_t_4); __pyx_t_4 = 0;
  if (unlikely(!__pyx_t_2)) __PYX_ERR(0, 5, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_t_2);
  __Pyx_DECREF(__pyx_t_3); __pyx_t_3 = 0;
  __pyx_t_3 = PyNumber_Multiply(__pyx_t_1, __pyx_t_2); if (unlikely(!__pyx_t_3)) __PYX_ERR(0, 5, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_t_3);
  __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
  __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0;
  __pyx_r = __pyx_t_3;
  __pyx_t_3 = 0;
  goto __pyx_L0;

 06: 
+07: def integrate0(double a,double b,int N):
/* Python wrapper */
static PyObject *__pyx_pw_46_cython_magic_f826ce6229b50973a9ef92e6c0ee0ea0_3integrate0(PyObject *__pyx_self, PyObject *__pyx_args, PyObject *__pyx_kwds); /*proto*/
static PyMethodDef __pyx_mdef_46_cython_magic_f826ce6229b50973a9ef92e6c0ee0ea0_3integrate0 = {"integrate0", (PyCFunction)(void*)(PyCFunctionWithKeywords)__pyx_pw_46_cython_magic_f826ce6229b50973a9ef92e6c0ee0ea0_3integrate0, METH_VARARGS|METH_KEYWORDS, 0};
static PyObject *__pyx_pw_46_cython_magic_f826ce6229b50973a9ef92e6c0ee0ea0_3integrate0(PyObject *__pyx_self, PyObject *__pyx_args, PyObject *__pyx_kwds) {
  double __pyx_v_a;
  double __pyx_v_b;
  int __pyx_v_N;
  PyObject *__pyx_r = 0;
  __Pyx_RefNannyDeclarations
  __Pyx_RefNannySetupContext("integrate0 (wrapper)", 0);
  {
    static PyObject **__pyx_pyargnames[] = {&__pyx_n_s_a,&__pyx_n_s_b,&__pyx_n_s_N,0};
    PyObject* values[3] = {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  3: values[2] = PyTuple_GET_ITEM(__pyx_args, 2);
        CYTHON_FALLTHROUGH;
        case  2: values[1] = PyTuple_GET_ITEM(__pyx_args, 1);
        CYTHON_FALLTHROUGH;
        case  1: values[0] = PyTuple_GET_ITEM(__pyx_args, 0);
        CYTHON_FALLTHROUGH;
        case  0: break;
        default: goto __pyx_L5_argtuple_error;
      }
      kw_args = PyDict_Size(__pyx_kwds);
      switch (pos_args) {
        case  0:
        if (likely((values[0] = __Pyx_PyDict_GetItemStr(__pyx_kwds, __pyx_n_s_a)) != 0)) kw_args--;
        else goto __pyx_L5_argtuple_error;
        CYTHON_FALLTHROUGH;
        case  1:
        if (likely((values[1] = __Pyx_PyDict_GetItemStr(__pyx_kwds, __pyx_n_s_b)) != 0)) kw_args--;
        else {
          __Pyx_RaiseArgtupleInvalid("integrate0", 1, 3, 3, 1); __PYX_ERR(0, 7, __pyx_L3_error)
        }
        CYTHON_FALLTHROUGH;
        case  2:
        if (likely((values[2] = __Pyx_PyDict_GetItemStr(__pyx_kwds, __pyx_n_s_N)) != 0)) kw_args--;
        else {
          __Pyx_RaiseArgtupleInvalid("integrate0", 1, 3, 3, 2); __PYX_ERR(0, 7, __pyx_L3_error)
        }
      }
      if (unlikely(kw_args > 0)) {
        if (unlikely(__Pyx_ParseOptionalKeywords(__pyx_kwds, __pyx_pyargnames, 0, values, pos_args, "integrate0") < 0)) __PYX_ERR(0, 7, __pyx_L3_error)
      }
    } else if (PyTuple_GET_SIZE(__pyx_args) != 3) {
      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);
    }
    __pyx_v_a = __pyx_PyFloat_AsDouble(values[0]); if (unlikely((__pyx_v_a == (double)-1) && PyErr_Occurred())) __PYX_ERR(0, 7, __pyx_L3_error)
    __pyx_v_b = __pyx_PyFloat_AsDouble(values[1]); if (unlikely((__pyx_v_b == (double)-1) && PyErr_Occurred())) __PYX_ERR(0, 7, __pyx_L3_error)
    __pyx_v_N = __Pyx_PyInt_As_int(values[2]); if (unlikely((__pyx_v_N == (int)-1) && PyErr_Occurred())) __PYX_ERR(0, 7, __pyx_L3_error)
  }
  goto __pyx_L4_argument_unpacking_done;
  __pyx_L5_argtuple_error:;
  __Pyx_RaiseArgtupleInvalid("integrate0", 1, 3, 3, PyTuple_GET_SIZE(__pyx_args)); __PYX_ERR(0, 7, __pyx_L3_error)
  __pyx_L3_error:;
  __Pyx_AddTraceback("_cython_magic_f826ce6229b50973a9ef92e6c0ee0ea0.integrate0", __pyx_clineno, __pyx_lineno, __pyx_filename);
  __Pyx_RefNannyFinishContext();
  return NULL;
  __pyx_L4_argument_unpacking_done:;
  __pyx_r = __pyx_pf_46_cython_magic_f826ce6229b50973a9ef92e6c0ee0ea0_2integrate0(__pyx_self, __pyx_v_a, __pyx_v_b, __pyx_v_N);
  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_f826ce6229b50973a9ef92e6c0ee0ea0_2integrate0(CYTHON_UNUSED PyObject *__pyx_self, double __pyx_v_a, double __pyx_v_b, int __pyx_v_N) {
  double __pyx_v_dx;
  double __pyx_v_x;
  double __pyx_v_s;
  CYTHON_UNUSED int __pyx_v__;
  PyObject *__pyx_r = NULL;
  __Pyx_RefNannyDeclarations
  __Pyx_RefNannySetupContext("integrate0", 0);
/* … */
  /* function exit code */
  __pyx_L1_error:;
  __Pyx_XDECREF(__pyx_t_5);
  __Pyx_XDECREF(__pyx_t_6);
  __Pyx_XDECREF(__pyx_t_7);
  __Pyx_XDECREF(__pyx_t_8);
  __Pyx_XDECREF(__pyx_t_9);
  __Pyx_AddTraceback("_cython_magic_f826ce6229b50973a9ef92e6c0ee0ea0.integrate0", __pyx_clineno, __pyx_lineno, __pyx_filename);
  __pyx_r = NULL;
  __pyx_L0:;
  __Pyx_XGIVEREF(__pyx_r);
  __Pyx_RefNannyFinishContext();
  return __pyx_r;
}
/* … */
  __pyx_tuple__4 = PyTuple_Pack(7, __pyx_n_s_a, __pyx_n_s_b, __pyx_n_s_N, __pyx_n_s_dx, __pyx_n_s_x, __pyx_n_s_s, __pyx_n_s__3); if (unlikely(!__pyx_tuple__4)) __PYX_ERR(0, 7, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_tuple__4);
  __Pyx_GIVEREF(__pyx_tuple__4);
/* … */
  __pyx_t_2 = PyCFunction_NewEx(&__pyx_mdef_46_cython_magic_f826ce6229b50973a9ef92e6c0ee0ea0_3integrate0, NULL, __pyx_n_s_cython_magic_f826ce6229b50973a9); if (unlikely(!__pyx_t_2)) __PYX_ERR(0, 7, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_t_2);
  if (PyDict_SetItem(__pyx_d, __pyx_n_s_integrate0, __pyx_t_2) < 0) __PYX_ERR(0, 7, __pyx_L1_error)
  __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0;

 08:     cdef double dx,x,s
 09:     cdef int _
+10:     dx = (b-a)/N
  __pyx_t_1 = (__pyx_v_b - __pyx_v_a);
  if (unlikely(__pyx_v_N == 0)) {
    PyErr_SetString(PyExc_ZeroDivisionError, "float division");
    __PYX_ERR(0, 10, __pyx_L1_error)
  }
  __pyx_v_dx = (__pyx_t_1 / ((double)__pyx_v_N));

+11:     x = a+dx/2
  __pyx_v_x = (__pyx_v_a + (__pyx_v_dx / 2.0));

+12:     s = 0.0
  __pyx_v_s = 0.0;

+13:     for _ in range(N):
  __pyx_t_2 = __pyx_v_N;
  __pyx_t_3 = __pyx_t_2;
  for (__pyx_t_4 = 0; __pyx_t_4 < __pyx_t_3; __pyx_t_4+=1) {
    __pyx_v__ = __pyx_t_4;

+14:         s+=f0(x)
    __pyx_t_5 = PyFloat_FromDouble(__pyx_v_s); if (unlikely(!__pyx_t_5)) __PYX_ERR(0, 14, __pyx_L1_error)
    __Pyx_GOTREF(__pyx_t_5);
    __Pyx_GetModuleGlobalName(__pyx_t_7, __pyx_n_s_f0); if (unlikely(!__pyx_t_7)) __PYX_ERR(0, 14, __pyx_L1_error)
    __Pyx_GOTREF(__pyx_t_7);
    __pyx_t_8 = PyFloat_FromDouble(__pyx_v_x); if (unlikely(!__pyx_t_8)) __PYX_ERR(0, 14, __pyx_L1_error)
    __Pyx_GOTREF(__pyx_t_8);
    __pyx_t_9 = NULL;
    if (CYTHON_UNPACK_METHODS && unlikely(PyMethod_Check(__pyx_t_7))) {
      __pyx_t_9 = PyMethod_GET_SELF(__pyx_t_7);
      if (likely(__pyx_t_9)) {
        PyObject* function = PyMethod_GET_FUNCTION(__pyx_t_7);
        __Pyx_INCREF(__pyx_t_9);
        __Pyx_INCREF(function);
        __Pyx_DECREF_SET(__pyx_t_7, function);
      }
    }
    __pyx_t_6 = (__pyx_t_9) ? __Pyx_PyObject_Call2Args(__pyx_t_7, __pyx_t_9, __pyx_t_8) : __Pyx_PyObject_CallOneArg(__pyx_t_7, __pyx_t_8);
    __Pyx_XDECREF(__pyx_t_9); __pyx_t_9 = 0;
    __Pyx_DECREF(__pyx_t_8); __pyx_t_8 = 0;
    if (unlikely(!__pyx_t_6)) __PYX_ERR(0, 14, __pyx_L1_error)
    __Pyx_GOTREF(__pyx_t_6);
    __Pyx_DECREF(__pyx_t_7); __pyx_t_7 = 0;
    __pyx_t_7 = PyNumber_InPlaceAdd(__pyx_t_5, __pyx_t_6); if (unlikely(!__pyx_t_7)) __PYX_ERR(0, 14, __pyx_L1_error)
    __Pyx_GOTREF(__pyx_t_7);
    __Pyx_DECREF(__pyx_t_5); __pyx_t_5 = 0;
    __Pyx_DECREF(__pyx_t_6); __pyx_t_6 = 0;
    __pyx_t_1 = __pyx_PyFloat_AsDouble(__pyx_t_7); if (unlikely((__pyx_t_1 == (double)-1) && PyErr_Occurred())) __PYX_ERR(0, 14, __pyx_L1_error)
    __Pyx_DECREF(__pyx_t_7); __pyx_t_7 = 0;
    __pyx_v_s = __pyx_t_1;

+15:         x += dx
    __pyx_v_x = (__pyx_v_x + __pyx_v_dx);
  }

+16:     return s*dx
  __Pyx_XDECREF(__pyx_r);
  __pyx_t_7 = PyFloat_FromDouble((__pyx_v_s * __pyx_v_dx)); if (unlikely(!__pyx_t_7)) __PYX_ERR(0, 16, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_t_7);
  __pyx_r = __pyx_t_7;
  __pyx_t_7 = 0;
  goto __pyx_L0;
In [23]:
from math import pi
In [24]:
integrate0(0.0,pi,1000000)
Out[24]:
0.5216069591343814
In [25]:
%timeit integrate0(0.0,pi,1000000)

106 ms ± 476 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)

turn f into a c-function¶

In [26]:
%%cython -a

from math import sin,exp

cdef double f1(double x):
    return sin(x)*exp(-x)

def integrate1(double a,double b,int N):
    cdef double dx,s
    cdef int i
    dx = (b-a)/N
    s = 0.0
    for i in range(N):
        s+=f1(a+(i+0.5)*dx)
    return s*dx
Out[26]:
Cython: _cython_magic_91f39410cca3e4dc6d7a18d3a7581e6f.pyx

Generated by Cython 0.29.30

Yellow lines hint at Python interaction.
Click on a line that starts with a "+" to see the C code that Cython generated for it. 

 01: 
+02: from math import sin,exp
  __pyx_t_1 = PyList_New(2); if (unlikely(!__pyx_t_1)) __PYX_ERR(0, 2, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_t_1);
  __Pyx_INCREF(__pyx_n_s_sin);
  __Pyx_GIVEREF(__pyx_n_s_sin);
  PyList_SET_ITEM(__pyx_t_1, 0, __pyx_n_s_sin);
  __Pyx_INCREF(__pyx_n_s_exp);
  __Pyx_GIVEREF(__pyx_n_s_exp);
  PyList_SET_ITEM(__pyx_t_1, 1, __pyx_n_s_exp);
  __pyx_t_2 = __Pyx_Import(__pyx_n_s_math, __pyx_t_1, 0); if (unlikely(!__pyx_t_2)) __PYX_ERR(0, 2, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_t_2);
  __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
  __pyx_t_1 = __Pyx_ImportFrom(__pyx_t_2, __pyx_n_s_sin); if (unlikely(!__pyx_t_1)) __PYX_ERR(0, 2, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_t_1);
  if (PyDict_SetItem(__pyx_d, __pyx_n_s_sin, __pyx_t_1) < 0) __PYX_ERR(0, 2, __pyx_L1_error)
  __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
  __pyx_t_1 = __Pyx_ImportFrom(__pyx_t_2, __pyx_n_s_exp); if (unlikely(!__pyx_t_1)) __PYX_ERR(0, 2, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_t_1);
  if (PyDict_SetItem(__pyx_d, __pyx_n_s_exp, __pyx_t_1) < 0) __PYX_ERR(0, 2, __pyx_L1_error)
  __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
  __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0;

 03: 
+04: cdef double f1(double x):
static double __pyx_f_46_cython_magic_91f39410cca3e4dc6d7a18d3a7581e6f_f1(double __pyx_v_x) {
  double __pyx_r;
  __Pyx_RefNannyDeclarations
  __Pyx_RefNannySetupContext("f1", 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_XDECREF(__pyx_t_5);
  __Pyx_WriteUnraisable("_cython_magic_91f39410cca3e4dc6d7a18d3a7581e6f.f1", __pyx_clineno, __pyx_lineno, __pyx_filename, 1, 0);
  __pyx_r = 0;
  __pyx_L0:;
  __Pyx_RefNannyFinishContext();
  return __pyx_r;
}

+05:     return sin(x)*exp(-x)
  __Pyx_GetModuleGlobalName(__pyx_t_2, __pyx_n_s_sin); if (unlikely(!__pyx_t_2)) __PYX_ERR(0, 5, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_t_2);
  __pyx_t_3 = PyFloat_FromDouble(__pyx_v_x); if (unlikely(!__pyx_t_3)) __PYX_ERR(0, 5, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_t_3);
  __pyx_t_4 = NULL;
  if (CYTHON_UNPACK_METHODS && unlikely(PyMethod_Check(__pyx_t_2))) {
    __pyx_t_4 = PyMethod_GET_SELF(__pyx_t_2);
    if (likely(__pyx_t_4)) {
      PyObject* function = PyMethod_GET_FUNCTION(__pyx_t_2);
      __Pyx_INCREF(__pyx_t_4);
      __Pyx_INCREF(function);
      __Pyx_DECREF_SET(__pyx_t_2, function);
    }
  }
  __pyx_t_1 = (__pyx_t_4) ? __Pyx_PyObject_Call2Args(__pyx_t_2, __pyx_t_4, __pyx_t_3) : __Pyx_PyObject_CallOneArg(__pyx_t_2, __pyx_t_3);
  __Pyx_XDECREF(__pyx_t_4); __pyx_t_4 = 0;
  __Pyx_DECREF(__pyx_t_3); __pyx_t_3 = 0;
  if (unlikely(!__pyx_t_1)) __PYX_ERR(0, 5, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_t_1);
  __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0;
  __Pyx_GetModuleGlobalName(__pyx_t_3, __pyx_n_s_exp); if (unlikely(!__pyx_t_3)) __PYX_ERR(0, 5, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_t_3);
  __pyx_t_4 = PyFloat_FromDouble((-__pyx_v_x)); if (unlikely(!__pyx_t_4)) __PYX_ERR(0, 5, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_t_4);
  __pyx_t_5 = NULL;
  if (CYTHON_UNPACK_METHODS && unlikely(PyMethod_Check(__pyx_t_3))) {
    __pyx_t_5 = PyMethod_GET_SELF(__pyx_t_3);
    if (likely(__pyx_t_5)) {
      PyObject* function = PyMethod_GET_FUNCTION(__pyx_t_3);
      __Pyx_INCREF(__pyx_t_5);
      __Pyx_INCREF(function);
      __Pyx_DECREF_SET(__pyx_t_3, function);
    }
  }
  __pyx_t_2 = (__pyx_t_5) ? __Pyx_PyObject_Call2Args(__pyx_t_3, __pyx_t_5, __pyx_t_4) : __Pyx_PyObject_CallOneArg(__pyx_t_3, __pyx_t_4);
  __Pyx_XDECREF(__pyx_t_5); __pyx_t_5 = 0;
  __Pyx_DECREF(__pyx_t_4); __pyx_t_4 = 0;
  if (unlikely(!__pyx_t_2)) __PYX_ERR(0, 5, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_t_2);
  __Pyx_DECREF(__pyx_t_3); __pyx_t_3 = 0;
  __pyx_t_3 = PyNumber_Multiply(__pyx_t_1, __pyx_t_2); if (unlikely(!__pyx_t_3)) __PYX_ERR(0, 5, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_t_3);
  __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
  __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0;
  __pyx_t_6 = __pyx_PyFloat_AsDouble(__pyx_t_3); if (unlikely((__pyx_t_6 == (double)-1) && PyErr_Occurred())) __PYX_ERR(0, 5, __pyx_L1_error)
  __Pyx_DECREF(__pyx_t_3); __pyx_t_3 = 0;
  __pyx_r = __pyx_t_6;
  goto __pyx_L0;

 06: 
+07: def integrate1(double a,double b,int N):
/* Python wrapper */
static PyObject *__pyx_pw_46_cython_magic_91f39410cca3e4dc6d7a18d3a7581e6f_1integrate1(PyObject *__pyx_self, PyObject *__pyx_args, PyObject *__pyx_kwds); /*proto*/
static PyMethodDef __pyx_mdef_46_cython_magic_91f39410cca3e4dc6d7a18d3a7581e6f_1integrate1 = {"integrate1", (PyCFunction)(void*)(PyCFunctionWithKeywords)__pyx_pw_46_cython_magic_91f39410cca3e4dc6d7a18d3a7581e6f_1integrate1, METH_VARARGS|METH_KEYWORDS, 0};
static PyObject *__pyx_pw_46_cython_magic_91f39410cca3e4dc6d7a18d3a7581e6f_1integrate1(PyObject *__pyx_self, PyObject *__pyx_args, PyObject *__pyx_kwds) {
  double __pyx_v_a;
  double __pyx_v_b;
  int __pyx_v_N;
  PyObject *__pyx_r = 0;
  __Pyx_RefNannyDeclarations
  __Pyx_RefNannySetupContext("integrate1 (wrapper)", 0);
  {
    static PyObject **__pyx_pyargnames[] = {&__pyx_n_s_a,&__pyx_n_s_b,&__pyx_n_s_N,0};
    PyObject* values[3] = {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  3: values[2] = PyTuple_GET_ITEM(__pyx_args, 2);
        CYTHON_FALLTHROUGH;
        case  2: values[1] = PyTuple_GET_ITEM(__pyx_args, 1);
        CYTHON_FALLTHROUGH;
        case  1: values[0] = PyTuple_GET_ITEM(__pyx_args, 0);
        CYTHON_FALLTHROUGH;
        case  0: break;
        default: goto __pyx_L5_argtuple_error;
      }
      kw_args = PyDict_Size(__pyx_kwds);
      switch (pos_args) {
        case  0:
        if (likely((values[0] = __Pyx_PyDict_GetItemStr(__pyx_kwds, __pyx_n_s_a)) != 0)) kw_args--;
        else goto __pyx_L5_argtuple_error;
        CYTHON_FALLTHROUGH;
        case  1:
        if (likely((values[1] = __Pyx_PyDict_GetItemStr(__pyx_kwds, __pyx_n_s_b)) != 0)) kw_args--;
        else {
          __Pyx_RaiseArgtupleInvalid("integrate1", 1, 3, 3, 1); __PYX_ERR(0, 7, __pyx_L3_error)
        }
        CYTHON_FALLTHROUGH;
        case  2:
        if (likely((values[2] = __Pyx_PyDict_GetItemStr(__pyx_kwds, __pyx_n_s_N)) != 0)) kw_args--;
        else {
          __Pyx_RaiseArgtupleInvalid("integrate1", 1, 3, 3, 2); __PYX_ERR(0, 7, __pyx_L3_error)
        }
      }
      if (unlikely(kw_args > 0)) {
        if (unlikely(__Pyx_ParseOptionalKeywords(__pyx_kwds, __pyx_pyargnames, 0, values, pos_args, "integrate1") < 0)) __PYX_ERR(0, 7, __pyx_L3_error)
      }
    } else if (PyTuple_GET_SIZE(__pyx_args) != 3) {
      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);
    }
    __pyx_v_a = __pyx_PyFloat_AsDouble(values[0]); if (unlikely((__pyx_v_a == (double)-1) && PyErr_Occurred())) __PYX_ERR(0, 7, __pyx_L3_error)
    __pyx_v_b = __pyx_PyFloat_AsDouble(values[1]); if (unlikely((__pyx_v_b == (double)-1) && PyErr_Occurred())) __PYX_ERR(0, 7, __pyx_L3_error)
    __pyx_v_N = __Pyx_PyInt_As_int(values[2]); if (unlikely((__pyx_v_N == (int)-1) && PyErr_Occurred())) __PYX_ERR(0, 7, __pyx_L3_error)
  }
  goto __pyx_L4_argument_unpacking_done;
  __pyx_L5_argtuple_error:;
  __Pyx_RaiseArgtupleInvalid("integrate1", 1, 3, 3, PyTuple_GET_SIZE(__pyx_args)); __PYX_ERR(0, 7, __pyx_L3_error)
  __pyx_L3_error:;
  __Pyx_AddTraceback("_cython_magic_91f39410cca3e4dc6d7a18d3a7581e6f.integrate1", __pyx_clineno, __pyx_lineno, __pyx_filename);
  __Pyx_RefNannyFinishContext();
  return NULL;
  __pyx_L4_argument_unpacking_done:;
  __pyx_r = __pyx_pf_46_cython_magic_91f39410cca3e4dc6d7a18d3a7581e6f_integrate1(__pyx_self, __pyx_v_a, __pyx_v_b, __pyx_v_N);
  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_91f39410cca3e4dc6d7a18d3a7581e6f_integrate1(CYTHON_UNUSED PyObject *__pyx_self, double __pyx_v_a, double __pyx_v_b, int __pyx_v_N) {
  double __pyx_v_dx;
  double __pyx_v_s;
  int __pyx_v_i;
  PyObject *__pyx_r = NULL;
  __Pyx_RefNannyDeclarations
  __Pyx_RefNannySetupContext("integrate1", 0);
/* … */
  /* function exit code */
  __pyx_L1_error:;
  __Pyx_XDECREF(__pyx_t_5);
  __Pyx_AddTraceback("_cython_magic_91f39410cca3e4dc6d7a18d3a7581e6f.integrate1", __pyx_clineno, __pyx_lineno, __pyx_filename);
  __pyx_r = NULL;
  __pyx_L0:;
  __Pyx_XGIVEREF(__pyx_r);
  __Pyx_RefNannyFinishContext();
  return __pyx_r;
}
/* … */
  __pyx_tuple_ = PyTuple_Pack(6, __pyx_n_s_a, __pyx_n_s_b, __pyx_n_s_N, __pyx_n_s_dx, __pyx_n_s_s, __pyx_n_s_i); if (unlikely(!__pyx_tuple_)) __PYX_ERR(0, 7, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_tuple_);
  __Pyx_GIVEREF(__pyx_tuple_);
/* … */
  __pyx_t_2 = PyCFunction_NewEx(&__pyx_mdef_46_cython_magic_91f39410cca3e4dc6d7a18d3a7581e6f_1integrate1, NULL, __pyx_n_s_cython_magic_91f39410cca3e4dc6d); if (unlikely(!__pyx_t_2)) __PYX_ERR(0, 7, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_t_2);
  if (PyDict_SetItem(__pyx_d, __pyx_n_s_integrate1, __pyx_t_2) < 0) __PYX_ERR(0, 7, __pyx_L1_error)
  __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0;

 08:     cdef double dx,s
 09:     cdef int i
+10:     dx = (b-a)/N
  __pyx_t_1 = (__pyx_v_b - __pyx_v_a);
  if (unlikely(__pyx_v_N == 0)) {
    PyErr_SetString(PyExc_ZeroDivisionError, "float division");
    __PYX_ERR(0, 10, __pyx_L1_error)
  }
  __pyx_v_dx = (__pyx_t_1 / ((double)__pyx_v_N));

+11:     s = 0.0
  __pyx_v_s = 0.0;

+12:     for i in range(N):
  __pyx_t_2 = __pyx_v_N;
  __pyx_t_3 = __pyx_t_2;
  for (__pyx_t_4 = 0; __pyx_t_4 < __pyx_t_3; __pyx_t_4+=1) {
    __pyx_v_i = __pyx_t_4;

+13:         s+=f1(a+(i+0.5)*dx)
    __pyx_v_s = (__pyx_v_s + __pyx_f_46_cython_magic_91f39410cca3e4dc6d7a18d3a7581e6f_f1((__pyx_v_a + ((__pyx_v_i + 0.5) * __pyx_v_dx))));
  }

+14:     return s*dx
  __Pyx_XDECREF(__pyx_r);
  __pyx_t_5 = PyFloat_FromDouble((__pyx_v_s * __pyx_v_dx)); if (unlikely(!__pyx_t_5)) __PYX_ERR(0, 14, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_t_5);
  __pyx_r = __pyx_t_5;
  __pyx_t_5 = 0;
  goto __pyx_L0;
In [27]:
%timeit integrate1(0.0,pi,1000000)

73.1 ms ± 150 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)

f is now not visible in Python anymore

In [28]:
f1(5)

---------------------------------------------------------------------------
NameError                                 Traceback (most recent call last)
/tmp/ipykernel_82657/1424691929.py in <module>
----> 1 f1(5)

NameError: name 'f1' is not defined

use cpdef to keep f in Python

use math functions from c¶

In [29]:
%%cython -a

from libc.math cimport sin,exp

cpdef double f2(double x):
    return sin(x)*exp(-x)

def integrate2(double a,double b,int N):
    cdef double dx,s
    cdef int i
    dx = (b-a)/N
    s = 0.0
    for i in range(N):
        s+=f2(a+(i+0.5)*dx)
    return s*dx
Out[29]:
Cython: _cython_magic_ba08a20ea0b9f9dcb944b03202da1d8d.pyx

Generated by Cython 0.29.30

Yellow lines hint at Python interaction.
Click on a line that starts with a "+" to see the C code that Cython generated for it. 

 01: 
 02: from libc.math cimport sin,exp
 03: 
+04: cpdef double f2(double x):
static PyObject *__pyx_pw_46_cython_magic_ba08a20ea0b9f9dcb944b03202da1d8d_1f2(PyObject *__pyx_self, PyObject *__pyx_arg_x); /*proto*/
static double __pyx_f_46_cython_magic_ba08a20ea0b9f9dcb944b03202da1d8d_f2(double __pyx_v_x, CYTHON_UNUSED int __pyx_skip_dispatch) {
  double __pyx_r;
  __Pyx_RefNannyDeclarations
  __Pyx_RefNannySetupContext("f2", 0);
/* … */
  /* function exit code */
  __pyx_L0:;
  __Pyx_RefNannyFinishContext();
  return __pyx_r;
}

/* Python wrapper */
static PyObject *__pyx_pw_46_cython_magic_ba08a20ea0b9f9dcb944b03202da1d8d_1f2(PyObject *__pyx_self, PyObject *__pyx_arg_x); /*proto*/
static PyObject *__pyx_pw_46_cython_magic_ba08a20ea0b9f9dcb944b03202da1d8d_1f2(PyObject *__pyx_self, PyObject *__pyx_arg_x) {
  double __pyx_v_x;
  PyObject *__pyx_r = 0;
  __Pyx_RefNannyDeclarations
  __Pyx_RefNannySetupContext("f2 (wrapper)", 0);
  assert(__pyx_arg_x); {
    __pyx_v_x = __pyx_PyFloat_AsDouble(__pyx_arg_x); if (unlikely((__pyx_v_x == (double)-1) && PyErr_Occurred())) __PYX_ERR(0, 4, __pyx_L3_error)
  }
  goto __pyx_L4_argument_unpacking_done;
  __pyx_L3_error:;
  __Pyx_AddTraceback("_cython_magic_ba08a20ea0b9f9dcb944b03202da1d8d.f2", __pyx_clineno, __pyx_lineno, __pyx_filename);
  __Pyx_RefNannyFinishContext();
  return NULL;
  __pyx_L4_argument_unpacking_done:;
  __pyx_r = __pyx_pf_46_cython_magic_ba08a20ea0b9f9dcb944b03202da1d8d_f2(__pyx_self, ((double)__pyx_v_x));
  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_ba08a20ea0b9f9dcb944b03202da1d8d_f2(CYTHON_UNUSED PyObject *__pyx_self, double __pyx_v_x) {
  PyObject *__pyx_r = NULL;
  __Pyx_RefNannyDeclarations
  __Pyx_RefNannySetupContext("f2", 0);
  __Pyx_XDECREF(__pyx_r);
  __pyx_t_1 = PyFloat_FromDouble(__pyx_f_46_cython_magic_ba08a20ea0b9f9dcb944b03202da1d8d_f2(__pyx_v_x, 0)); if (unlikely(!__pyx_t_1)) __PYX_ERR(0, 4, __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_ba08a20ea0b9f9dcb944b03202da1d8d.f2", __pyx_clineno, __pyx_lineno, __pyx_filename);
  __pyx_r = NULL;
  __pyx_L0:;
  __Pyx_XGIVEREF(__pyx_r);
  __Pyx_RefNannyFinishContext();
  return __pyx_r;
}

+05:     return sin(x)*exp(-x)
  __pyx_r = (sin(__pyx_v_x) * exp((-__pyx_v_x)));
  goto __pyx_L0;

 06: 
+07: def integrate2(double a,double b,int N):
/* Python wrapper */
static PyObject *__pyx_pw_46_cython_magic_ba08a20ea0b9f9dcb944b03202da1d8d_3integrate2(PyObject *__pyx_self, PyObject *__pyx_args, PyObject *__pyx_kwds); /*proto*/
static PyMethodDef __pyx_mdef_46_cython_magic_ba08a20ea0b9f9dcb944b03202da1d8d_3integrate2 = {"integrate2", (PyCFunction)(void*)(PyCFunctionWithKeywords)__pyx_pw_46_cython_magic_ba08a20ea0b9f9dcb944b03202da1d8d_3integrate2, METH_VARARGS|METH_KEYWORDS, 0};
static PyObject *__pyx_pw_46_cython_magic_ba08a20ea0b9f9dcb944b03202da1d8d_3integrate2(PyObject *__pyx_self, PyObject *__pyx_args, PyObject *__pyx_kwds) {
  double __pyx_v_a;
  double __pyx_v_b;
  int __pyx_v_N;
  PyObject *__pyx_r = 0;
  __Pyx_RefNannyDeclarations
  __Pyx_RefNannySetupContext("integrate2 (wrapper)", 0);
  {
    static PyObject **__pyx_pyargnames[] = {&__pyx_n_s_a,&__pyx_n_s_b,&__pyx_n_s_N,0};
    PyObject* values[3] = {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  3: values[2] = PyTuple_GET_ITEM(__pyx_args, 2);
        CYTHON_FALLTHROUGH;
        case  2: values[1] = PyTuple_GET_ITEM(__pyx_args, 1);
        CYTHON_FALLTHROUGH;
        case  1: values[0] = PyTuple_GET_ITEM(__pyx_args, 0);
        CYTHON_FALLTHROUGH;
        case  0: break;
        default: goto __pyx_L5_argtuple_error;
      }
      kw_args = PyDict_Size(__pyx_kwds);
      switch (pos_args) {
        case  0:
        if (likely((values[0] = __Pyx_PyDict_GetItemStr(__pyx_kwds, __pyx_n_s_a)) != 0)) kw_args--;
        else goto __pyx_L5_argtuple_error;
        CYTHON_FALLTHROUGH;
        case  1:
        if (likely((values[1] = __Pyx_PyDict_GetItemStr(__pyx_kwds, __pyx_n_s_b)) != 0)) kw_args--;
        else {
          __Pyx_RaiseArgtupleInvalid("integrate2", 1, 3, 3, 1); __PYX_ERR(0, 7, __pyx_L3_error)
        }
        CYTHON_FALLTHROUGH;
        case  2:
        if (likely((values[2] = __Pyx_PyDict_GetItemStr(__pyx_kwds, __pyx_n_s_N)) != 0)) kw_args--;
        else {
          __Pyx_RaiseArgtupleInvalid("integrate2", 1, 3, 3, 2); __PYX_ERR(0, 7, __pyx_L3_error)
        }
      }
      if (unlikely(kw_args > 0)) {
        if (unlikely(__Pyx_ParseOptionalKeywords(__pyx_kwds, __pyx_pyargnames, 0, values, pos_args, "integrate2") < 0)) __PYX_ERR(0, 7, __pyx_L3_error)
      }
    } else if (PyTuple_GET_SIZE(__pyx_args) != 3) {
      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);
    }
    __pyx_v_a = __pyx_PyFloat_AsDouble(values[0]); if (unlikely((__pyx_v_a == (double)-1) && PyErr_Occurred())) __PYX_ERR(0, 7, __pyx_L3_error)
    __pyx_v_b = __pyx_PyFloat_AsDouble(values[1]); if (unlikely((__pyx_v_b == (double)-1) && PyErr_Occurred())) __PYX_ERR(0, 7, __pyx_L3_error)
    __pyx_v_N = __Pyx_PyInt_As_int(values[2]); if (unlikely((__pyx_v_N == (int)-1) && PyErr_Occurred())) __PYX_ERR(0, 7, __pyx_L3_error)
  }
  goto __pyx_L4_argument_unpacking_done;
  __pyx_L5_argtuple_error:;
  __Pyx_RaiseArgtupleInvalid("integrate2", 1, 3, 3, PyTuple_GET_SIZE(__pyx_args)); __PYX_ERR(0, 7, __pyx_L3_error)
  __pyx_L3_error:;
  __Pyx_AddTraceback("_cython_magic_ba08a20ea0b9f9dcb944b03202da1d8d.integrate2", __pyx_clineno, __pyx_lineno, __pyx_filename);
  __Pyx_RefNannyFinishContext();
  return NULL;
  __pyx_L4_argument_unpacking_done:;
  __pyx_r = __pyx_pf_46_cython_magic_ba08a20ea0b9f9dcb944b03202da1d8d_2integrate2(__pyx_self, __pyx_v_a, __pyx_v_b, __pyx_v_N);
  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_ba08a20ea0b9f9dcb944b03202da1d8d_2integrate2(CYTHON_UNUSED PyObject *__pyx_self, double __pyx_v_a, double __pyx_v_b, int __pyx_v_N) {
  double __pyx_v_dx;
  double __pyx_v_s;
  int __pyx_v_i;
  PyObject *__pyx_r = NULL;
  __Pyx_RefNannyDeclarations
  __Pyx_RefNannySetupContext("integrate2", 0);
/* … */
  /* function exit code */
  __pyx_L1_error:;
  __Pyx_XDECREF(__pyx_t_5);
  __Pyx_AddTraceback("_cython_magic_ba08a20ea0b9f9dcb944b03202da1d8d.integrate2", __pyx_clineno, __pyx_lineno, __pyx_filename);
  __pyx_r = NULL;
  __pyx_L0:;
  __Pyx_XGIVEREF(__pyx_r);
  __Pyx_RefNannyFinishContext();
  return __pyx_r;
}
/* … */
  __pyx_tuple_ = PyTuple_Pack(6, __pyx_n_s_a, __pyx_n_s_b, __pyx_n_s_N, __pyx_n_s_dx, __pyx_n_s_s, __pyx_n_s_i); if (unlikely(!__pyx_tuple_)) __PYX_ERR(0, 7, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_tuple_);
  __Pyx_GIVEREF(__pyx_tuple_);
/* … */
  __pyx_t_1 = PyCFunction_NewEx(&__pyx_mdef_46_cython_magic_ba08a20ea0b9f9dcb944b03202da1d8d_3integrate2, NULL, __pyx_n_s_cython_magic_ba08a20ea0b9f9dcb9); if (unlikely(!__pyx_t_1)) __PYX_ERR(0, 7, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_t_1);
  if (PyDict_SetItem(__pyx_d, __pyx_n_s_integrate2, __pyx_t_1) < 0) __PYX_ERR(0, 7, __pyx_L1_error)
  __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;

 08:     cdef double dx,s
 09:     cdef int i
+10:     dx = (b-a)/N
  __pyx_t_1 = (__pyx_v_b - __pyx_v_a);
  if (unlikely(__pyx_v_N == 0)) {
    PyErr_SetString(PyExc_ZeroDivisionError, "float division");
    __PYX_ERR(0, 10, __pyx_L1_error)
  }
  __pyx_v_dx = (__pyx_t_1 / ((double)__pyx_v_N));

+11:     s = 0.0
  __pyx_v_s = 0.0;

+12:     for i in range(N):
  __pyx_t_2 = __pyx_v_N;
  __pyx_t_3 = __pyx_t_2;
  for (__pyx_t_4 = 0; __pyx_t_4 < __pyx_t_3; __pyx_t_4+=1) {
    __pyx_v_i = __pyx_t_4;

+13:         s+=f2(a+(i+0.5)*dx)
    __pyx_v_s = (__pyx_v_s + __pyx_f_46_cython_magic_ba08a20ea0b9f9dcb944b03202da1d8d_f2((__pyx_v_a + ((__pyx_v_i + 0.5) * __pyx_v_dx)), 0));
  }

+14:     return s*dx
  __Pyx_XDECREF(__pyx_r);
  __pyx_t_5 = PyFloat_FromDouble((__pyx_v_s * __pyx_v_dx)); if (unlikely(!__pyx_t_5)) __PYX_ERR(0, 14, __pyx_L1_error)
  __Pyx_GOTREF(__pyx_t_5);
  __pyx_r = __pyx_t_5;
  __pyx_t_5 = 0;
  goto __pyx_L0;
In [30]:
%timeit integrate2(0.0,pi,1000000)

17.1 ms ± 31.5 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)

Line 10 is yellow, due to checking for ZeroDivisionError. We can deactivate this with

import cython

@cython.cdivision(True)
def integrate(...

But:

  • we would again open another door for errors
  • this line is only executed once

Was it worth it?¶

In [31]:
import numpy as np

def f3(x):
    return np.sin(x)*np.exp(-x)

def integrate3(a,b,N):
    dx = (b-a)/N
    x = (np.arange(N)+0.5)*dx+a
    return np.sum(f3(x))*dx
In [32]:
%timeit integrate3(0.0,pi,1000000)

27.5 ms ± 81 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)

Cython and STL containers¶

In [33]:
%%cython -+
#distutils: language = c++
from libcpp.vector cimport vector
def manipulateVector(vector[double] v):
    v.push_back(42.0) # Equivalent to append for list in python
    return v
In [34]:
manipulateVector([1,2,3,4])
Out[34]:
[1.0, 2.0, 3.0, 4.0, 42.0]

Cython and exceptions¶

In [35]:
%%cython
cpdef int raiseError():
    raise RuntimeError("A problem!")
    return 1
In [36]:
raiseError()
print("Still going strong!")

---------------------------------------------------------------------------
RuntimeError                              Traceback (most recent call last)
RuntimeError: A problem!
Exception ignored in: '_cython_magic_b93eae444b061f456eb12c45e2102e26.raiseError'
Traceback (most recent call last):
  File "/tmp/ipykernel_82657/3788160389.py", line 1, in <module>
RuntimeError: A problem!

Still going strong!
In [37]:
%%cython
cpdef int raiseErrorBetter() except *:
    raise RuntimeError("A problem!")
    return 1
In [38]:
raiseErrorBetter()
print("Still going strong!")

---------------------------------------------------------------------------
RuntimeError                              Traceback (most recent call last)
/tmp/ipykernel_82657/3461535642.py in <module>
----> 1 raiseErrorBetter()
      2 print("Still going strong!")

_cython_magic_9274737b2f24bc094c22bd49becb1e29.pyx in _cython_magic_9274737b2f24bc094c22bd49becb1e29.raiseErrorBetter()

_cython_magic_9274737b2f24bc094c22bd49becb1e29.pyx in _cython_magic_9274737b2f24bc094c22bd49becb1e29.raiseErrorBetter()

RuntimeError: A problem!
In [39]:
try:
    raiseErrorBetter()
except RuntimeError:
    print("Error caught")

Error caught
In [40]:
%%cython
cpdef int raiseException(int n) except -1:
    # Assuming C function returning an integer; if -1 --> Error appeared
    return n
In [41]:
raiseException(1)
Out[41]:
1
In [42]:
raiseException(-1)

---------------------------------------------------------------------------
SystemError                               Traceback (most recent call last)
/tmp/ipykernel_82657/4102492971.py in <module>
----> 1 raiseException(-1)

SystemError: <built-in function raiseException> returned NULL without setting an exception

Back to integration - Cython and Classes¶

In [43]:
%%cython

from libc.math cimport sin,exp

cdef class Integrand(object):
    cpdef double evaluate(self,double x):
        raise NotImplementedError()
        
cdef class SinExpFunction(Integrand):
    cpdef double evaluate(self,double x):
        return sin(x)*exp(-x)
    
def integrate(Integrand f, double a, double b, int N):
    cdef int i
    cdef double dx,s
    dx = (b-a)/N
    s = 0.0
    for i in range(N):
        s += f.evaluate(a+(i+0.5)*dx)
    return s*dx
In [44]:
f = SinExpFunction()
In [45]:
import math
In [46]:
%timeit integrate(f,0,math.pi,1000000)

19.8 ms ± 81.9 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)
In [47]:
class Poly(Integrand):
    def evaluate(self,x):
        return 3*x-x**3
In [48]:
g = Poly()
integrate(g,0,math.pi,10000)
Out[48]:
-9.5478660351052
In [49]:
%timeit integrate(g,0,math.pi,1000000)

220 ms ± 1.64 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)
In [ ]: