setuptools example working

Author: oleksandr-pavlyk

common_ext/cy_kde.pyx

@@ -0,0 +1,61 @@
+# distutils: language = c++
+# cython: language_level=3
+
+cimport cython
+from libc.stdint cimport uint16_t
+
+cdef extern from "CL/sycl.hpp" namespace "sycl":
+    cdef cppclass dpcpp_queue "sycl::queue":
+        pass
+
+
+cdef extern from "kde.hpp" namespace "example":
+    void kernel_density_estimate[T](
+        dpcpp_queue,   # execution queue
+        uint16_t,      # datadimensionality
+        const T *,     # evaluation points
+        T *,           # output for kde values
+        size_t,        # number of evaluation points
+        const T *,     # observation data
+        size_t,        # number of observations
+        T              # smoothing parameter
+    )
+
+cimport dpctl as c_dpctl
+import numpy as np
+
+def kde_eval(
+        c_dpctl.SyclQueue py_q,
+        cython.floating[:, :] x,
+        cython.floating[:, :] x_data,
+        cython.floating h
+):
+    cdef cython.floating[:] f
+
+    cdef c_dpctl.DPCTLSyclQueueRef qref = py_q.get_queue_ref()
+    cdef dpcpp_queue* sycl_queue = <dpcpp_queue *>qref
+
+    if x.shape[1] != x_data.shape[1]:
+        raise ValueError("Evaluation data and observation data have different dimensions")
+
+    if (x.shape[0] == 0 or x.shape[1] == 0 or
+        x_data.shape[0] == 0 or x_data.shape[1] == 0):
+        raise ValueError("Evaluation and observation data must be non-empty")
+
+    if cython.floating is float:
+        f = np.empty(x.shape[0], dtype=np.float)
+    else:
+        f = np.empty(x.shape[0], dtype=np.double)
+
+    kernel_density_estimate(
+        sycl_queue[0],
+        <uint16_t>x.shape[1],
+        &x[0, 0],
+        &f[0],
+        f.size,
+        &x_data[0,0],
+        x_data.shape[0],
+        h
+    )
+
+    return np.asarray(f)

common_ext/py_kde.cpp

@@ -0,0 +1,90 @@
+#include <CL/sycl.hpp>
+#include <cstdint>
+#include <limits>
+#include <pybind11/pybind11.h>
+#include <pybind11/numpy.h>
+#include "dpctl4pybind11.hpp"
+#include "kde.hpp"
+
+namespace py = pybind11;
+
+template <typename T>
+py::array_t<T>
+py_kde_eval_t(
+    sycl::queue q,
+    py::array_t<T, py::array::c_style|py::array::forcecast> x,
+    py::array_t<T, py::array::c_style|py::array::forcecast> data,
+    T h)
+{
+    py::buffer_info x_pybuf = x.request();
+    py::buffer_info data_pybuf = data.request();
+
+    if (x_pybuf.ndim != 2 || data_pybuf.ndim != 2) {
+	throw py::value_error("Input arrays must be matrices");
+    }
+
+    auto dim = x_pybuf.shape[1];
+    constexpr auto dim_max = static_cast<decltype(dim)>(std::numeric_limits<std::uint16_t>::max());
+
+    if (dim != data_pybuf.shape[1] || dim > dim_max) {
+	throw py::value_error("Inputs have inconsistent functionality or too large a width");
+    }
+
+    auto n = x_pybuf.shape[0];
+    auto n_data = data_pybuf.shape[0];
+
+    T *x_ptr = reinterpret_cast<T *>(x_pybuf.ptr);
+    T *data_ptr = reinterpret_cast<T *>(data_pybuf.ptr);
+
+    py::array_t<T, py::array::c_style> f({n}, {sizeof(T)});
+
+    example::kernel_density_estimate<T>(
+	q, static_cast<std::uint16_t>(dim),
+	const_cast<const T*>(x_ptr),
+	f.mutable_data(0),
+	n,
+	const_cast<const T*>(data_ptr),
+	n_data,
+	h);
+
+    return f;
+}
+
+py::array
+py_kde_eval(
+    sycl::queue exec_q,
+    py::array x,
+    py::array data,
+    py::object h)
+{
+    if (py::isinstance<py::array_t<double>>(x) && py::isinstance<py::array_t<double>>(data)) {
+	std::cout << "Dealing with doubles" << std::endl;
+	return py_kde_eval_t<double>(
+	    exec_q,
+	    py::cast<py::array_t<double>>(x),
+	    py::cast<py::array_t<double>>(data),
+	    py::cast<double>(h)
+	    );
+    } else if (py::isinstance<py::array_t<float>>(x) && py::isinstance<py::array_t<float>>(data)) {
+	return py_kde_eval_t<float>(
+	    exec_q,
+	    py::cast<py::array_t<float>>(x),
+	    py::cast<py::array_t<float>>(data),
+	    py::cast<float>(h)
+	    );
+    } else {
+	throw py::type_error("Both arrays must be either single or double precision floating point types");
+    }
+}
+
+
+PYBIND11_MODULE(_pybind11_kde, m) {
+    import_dpctl();
+
+    m.def("kde_eval", &py_kde_eval,
+	  "Evaluate kernel density estimation function for every argument for the dataset",
+	  py::arg("exec_q"),
+	  py::arg("x"),
+	  py::arg("data"),
+	  py::arg("h") );
+}

kde_setuptools/.gitignore

@@ -0,0 +1,6 @@
+__pycache__
+build
+*.so
+src/cy_kde.cpp
+*~
+kde_setuptools.egg-info

kde_setuptools/README.md

@@ -0,0 +1,8 @@
+# Building this extension
+
+Assuming `Cython` and `pybind11` are installed and DCP++ has been activated:
+
+```bash
+CC=dpcpp LDSHARED="dpcpp --shared" python setup.py develop
+pytest -m tests
+```

kde_setuptools/kde_setuptools/__init__.py

@@ -0,0 +1,4 @@
+from ._cython_kde import kde_eval as cython_kde_eval
+from ._pybind11_kde import kde_eval as pybind11_kde_eval
+
+__all__ = ["cyton_kde_eval", "pybind11_kde_eval"]

kde_setuptools/setup.py

@@ -0,0 +1,33 @@
+from setuptools import setup, Extension
+from pybind11.setup_helpers import Pybind11Extension
+import dpctl
+
+ext_modules = [
+    Pybind11Extension(
+        "kde_setuptools._pybind11_kde",
+        ["./src/pybind11_kde.cpp",],
+        include_dirs=["./src", dpctl.get_include(),]
+    ),
+    Extension(
+        "kde_setuptools._cython_kde",
+        ["./src/cy_kde.pyx",],
+        include_dirs=["./src", dpctl.get_include(),],
+        language="c++"
+    ),
+]
+
+setup(
+    name="kde_setuptools",
+    author="Intel Corporation",
+    version="0.0.1",
+    description="An example of data-parallel extensions built with oneAPI",
+    long_description="""
+    Example of using oneAPI to build data-parallel extension using setuptools.
+
+    Part of oneAPI for Scientific Python community virtual poster.
+    Also see README.md
+    """,
+    license="Apache 2.0",
+    url="https://intelpython.github.io/oneAPI-for-SciPy",
+    ext_modules=ext_modules
+)

kde_setuptools/src/cy_kde.pyx

@@ -0,0 +1 @@
+../../common_ext/cy_kde.pyx

kde_setuptools/src/kde.hpp

@@ -0,0 +1 @@
+../../common_src/kde.hpp

kde_setuptools/src/pybind11_kde.cpp

@@ -0,0 +1 @@
+../../common_ext/py_kde.cpp

kde_setuptools/tests/test_cy_kde.py

@@ -0,0 +1,76 @@
+import numpy as np
+import dpctl
+
+from kde_setuptools import cython_kde_eval as kde_eval
+
+import pytest
+
+def ref_kde(x, data, h):
+    """
+    Reference NumPy implementation for KDE evaluation
+    """
+    assert x.ndim == 2 and data.ndim == 2
+    assert x.shape[1] == data.shape[1]
+    dim = x.shape[1]
+    n_data = data.shape[0]
+    return np.exp(
+        np.square(x[:, np.newaxis, :]-data).sum(axis=-1)/(-2*h*h)
+    ).sum(axis=1)/(np.sqrt(2*np.pi)*h)**dim / n_data
+
+
+def test_1d():
+    try:
+        q = dpctl.SyclQueue()
+    except dpctl.SyclQueueCreationError:
+        pytest.skip("Execution queue could not be created, skipping...")
+
+    x = np.linspace(0.1, 0.9, num=15).reshape((-1, 1))
+    data = np.random.rand(128,1)
+
+    f_cy = kde_eval(q, x, data, 0.1)
+    f_ref = ref_kde(x, data, 0.1)
+
+    assert np.allclose(f_cy, f_ref)
+
+
+def test_2d():
+    try:
+        q = dpctl.SyclQueue()
+    except dpctl.SyclQueueCreationError:
+        pytest.skip("Execution queue could not be created, skipping...")
+
+    x = np.dstack(
+        np.meshgrid(
+            np.linspace(0.1, 0.9, num=7),
+            np.linspace(0.1, 0.9, num=7)
+        )
+    ).reshape(-1, 2)
+
+    data = np.random.rand(128*128, 2)
+
+    f_cy = kde_eval(q, x, data, 0.05)
+    f_ref = ref_kde(x, data, 0.05)
+
+    assert np.allclose(f_cy, f_ref)
+
+
+def test_3d():
+    try:
+        q = dpctl.SyclQueue()
+    except dpctl.SyclQueueCreationError:
+        pytest.skip("Execution queue could not be created, skipping...")
+
+    x = np.dstack(
+        np.meshgrid(
+            np.linspace(0.1, 0.9, num=7),
+            np.linspace(0.1, 0.9, num=7),
+            np.linspace(0.1, 0.9, num=7)
+        )
+    ).reshape(-1, 3)
+
+    data = np.random.rand(16000,3)
+
+    f_cy = kde_eval(q, x, data, 0.01)
+    f_ref = ref_kde(x, data, 0.01)
+
+    assert np.allclose(f_cy, f_ref)