gistfile1.py

import numpy as np
cimport numpy as np

from cpython cimport PyCObject_AsVoidPtr
from scipy.linalg.blas import cblas

# The following doesn' work (compiles, but returns incorrect values).
# ctypedef float (*sdot_ptr) (const int *N, const float *X, const int *incX, const float *Y, const int *incY) nogil
# cdef sdot_ptr sdot=<sdot_ptr>PyCObject_AsVoidPtr(cblas.sdot._cpointer)

# Works. Why does the return value have to be "double"?
ctypedef double (*sdot_ptr) (const int *N, const float *X, const int *incX, const float *Y, const int *incY) nogil
cdef sdot_ptr sdot=<sdot_ptr>PyCObject_AsVoidPtr(cblas.sdot._cpointer)


def test_sdot():
    cdef int one = 1
    cdef int size = 100

    x = np.arange(1, size + 1).astype(np.float32)
    y = 1.0 / x
    return sdot(&size, <float *>(np.PyArray_DATA(x)), &one, <float *>(np.PyArray_DATA(y)), &one)


print test_sdot()  # 100.0 when "double", but 0.0 when "float"