diff --git a/.github/workflows/build-sphinx.yml b/.github/workflows/build-sphinx.yml index f6664f66c1f0..f4e3b74c6237 100644 --- a/.github/workflows/build-sphinx.yml +++ b/.github/workflows/build-sphinx.yml @@ -17,6 +17,7 @@ jobs: env: python-ver: '3.9' + CHANNELS: '-c dppy/label/dev -c intel -c conda-forge --override-channels' steps: - name: Cancel Previous Runs @@ -74,10 +75,10 @@ jobs: - name: Install dpnp dependencies run: | conda install dpctl mkl-devel-dpcpp onedpl-devel tbb-devel dpcpp_linux-64 \ - cmake cython pytest ninja scikit-build -c dppy/label/dev -c intel -c conda-forge + cmake cython pytest ninja scikit-build sysroot_linux-64">=2.28" ${{ env.CHANNELS }} - name: Install cuPy dependencies - run: conda install -c conda-forge cupy cudatoolkit=10.0 + run: conda install cupy cudatoolkit=10.0 - name: Conda info run: conda info diff --git a/.github/workflows/generate_coverage.yaml b/.github/workflows/generate_coverage.yaml index 4d054274502f..fd38dde9a5d3 100644 --- a/.github/workflows/generate_coverage.yaml +++ b/.github/workflows/generate_coverage.yaml @@ -15,6 +15,7 @@ jobs: env: python-ver: '3.10' + CHANNELS: '-c dppy/label/dev -c intel -c conda-forge --override-channels' steps: - name: Cancel Previous Runs @@ -34,7 +35,6 @@ jobs: python-version: ${{ env.python-ver }} miniconda-version: 'latest' activate-environment: 'coverage' - channels: intel, conda-forge - name: Install Lcov run: | @@ -42,7 +42,7 @@ jobs: - name: Install dpnp dependencies run: | conda install cython llvm cmake scikit-build ninja pytest pytest-cov coverage[toml] \ - dppy/label/dev::dpctl dpcpp_linux-64 mkl-devel-dpcpp tbb-devel onedpl-devel + dpctl dpcpp_linux-64 sysroot_linux-64">=2.28" mkl-devel-dpcpp tbb-devel onedpl-devel ${{ env.CHANNELS }} - name: Conda info run: | conda info diff --git a/CMakeLists.txt b/CMakeLists.txt index cdecc3cefd72..efa35ac50869 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -57,6 +57,9 @@ set(CYTHON_FLAGS "-t -w \"${CMAKE_SOURCE_DIR}\"") find_package(Cython REQUIRED) find_package(Dpctl REQUIRED) +message(STATUS "Dpctl_INCLUDE_DIRS=" ${Dpctl_INCLUDE_DIRS}) +message(STATUS "Dpctl_TENSOR_INCLUDE_DIR=" ${Dpctl_TENSOR_INCLUDE_DIR}) + if(WIN32) string(CONCAT WARNING_FLAGS "-Wall " diff --git a/dpnp/backend/extensions/lapack/CMakeLists.txt b/dpnp/backend/extensions/lapack/CMakeLists.txt index 8040fb433cb9..e54de4068c01 100644 --- a/dpnp/backend/extensions/lapack/CMakeLists.txt +++ b/dpnp/backend/extensions/lapack/CMakeLists.txt @@ -45,6 +45,7 @@ target_include_directories(${python_module_name} PRIVATE ${CMAKE_CURRENT_SOURCE_ target_include_directories(${python_module_name} PRIVATE ${CMAKE_CURRENT_SOURCE_DIR}/../../src) target_include_directories(${python_module_name} PUBLIC ${Dpctl_INCLUDE_DIRS}) +target_include_directories(${python_module_name} PUBLIC ${Dpctl_TENSOR_INCLUDE_DIR}) if (WIN32) target_compile_options(${python_module_name} PRIVATE diff --git a/dpnp/backend/extensions/lapack/heevd.cpp b/dpnp/backend/extensions/lapack/heevd.cpp index 8c943646ff0a..f99fb94c18ec 100644 --- a/dpnp/backend/extensions/lapack/heevd.cpp +++ b/dpnp/backend/extensions/lapack/heevd.cpp @@ -26,7 +26,12 @@ #include +// dpctl tensor headers +#include "utils/memory_overlap.hpp" +#include "utils/type_utils.hpp" + #include "heevd.hpp" +#include "types_matrix.hpp" #include "dpnp_utils.hpp" @@ -42,19 +47,34 @@ namespace lapack namespace mkl_lapack = oneapi::mkl::lapack; namespace py = pybind11; +namespace type_utils = dpctl::tensor::type_utils; + +typedef sycl::event (*heevd_impl_fn_ptr_t)(sycl::queue, + const oneapi::mkl::job, + const oneapi::mkl::uplo, + const std::int64_t, + char*, + char*, + std::vector&, + const std::vector&); + +static heevd_impl_fn_ptr_t heevd_dispatch_table[dpctl_td_ns::num_types][dpctl_td_ns::num_types]; template -static sycl::event call_heevd(sycl::queue exec_q, +static sycl::event heevd_impl(sycl::queue exec_q, const oneapi::mkl::job jobz, const oneapi::mkl::uplo upper_lower, const std::int64_t n, - T* a, - RealT* w, + char* in_a, + char* out_w, std::vector& host_task_events, const std::vector& depends) { - validate_type_for_device(exec_q); - validate_type_for_device(exec_q); + type_utils::validate_type_for_device(exec_q); + type_utils::validate_type_for_device(exec_q); + + T* a = reinterpret_cast(in_a); + RealT* w = reinterpret_cast(out_w); const std::int64_t lda = std::max(1UL, n); const std::int64_t scratchpad_size = mkl_lapack::heevd_scratchpad_size(exec_q, jobz, upper_lower, n, lda); @@ -163,13 +183,11 @@ std::pair heevd(sycl::queue exec_q, throw py::value_error("Execution queue is not compatible with allocation queues"); } - // check that arrays do not overlap, and concurrent access is safe. - // TODO: need to be exposed by DPCTL headers - // auto const &overlap = dpctl::tensor::overlap::MemoryOverlap(); - // if (overlap(eig_vecs, eig_vals)) - // { - // throw py::value_error("Arrays index overlapping segments of memory"); - // } + auto const& overlap = dpctl::tensor::overlap::MemoryOverlap(); + if (overlap(eig_vecs, eig_vals)) + { + throw py::value_error("Arrays with eigenvectors and eigenvalues are overlapping segments of memory"); + } bool is_eig_vecs_f_contig = eig_vecs.is_f_contiguous(); bool is_eig_vals_c_contig = eig_vals.is_c_contiguous(); @@ -182,38 +200,51 @@ std::pair heevd(sycl::queue exec_q, throw py::value_error("An array with output eigenvalues must be C-contiguous"); } - int eig_vecs_typenum = eig_vecs.get_typenum(); - int eig_vals_typenum = eig_vals.get_typenum(); - auto const& dpctl_capi = dpctl::detail::dpctl_capi::get(); + auto array_types = dpctl_td_ns::usm_ndarray_types(); + int eig_vecs_type_id = array_types.typenum_to_lookup_id(eig_vecs.get_typenum()); + int eig_vals_type_id = array_types.typenum_to_lookup_id(eig_vals.get_typenum()); - sycl::event heevd_ev; - std::vector host_task_events; + heevd_impl_fn_ptr_t heevd_fn = heevd_dispatch_table[eig_vecs_type_id][eig_vals_type_id]; + if (heevd_fn == nullptr) + { + throw py::value_error("No heevd implementation defined for a pair of type for eigenvectors and eigenvalues"); + } + + char* eig_vecs_data = eig_vecs.get_data(); + char* eig_vals_data = eig_vals.get_data(); const std::int64_t n = eig_vecs_shape[0]; const oneapi::mkl::job jobz_val = static_cast(jobz); const oneapi::mkl::uplo uplo_val = static_cast(upper_lower); - if ((eig_vecs_typenum == dpctl_capi.UAR_CDOUBLE_) && (eig_vals_typenum == dpctl_capi.UAR_DOUBLE_)) - { - std::complex* a = reinterpret_cast*>(eig_vecs.get_data()); - double* w = reinterpret_cast(eig_vals.get_data()); + std::vector host_task_events; + sycl::event heevd_ev = + heevd_fn(exec_q, jobz_val, uplo_val, n, eig_vecs_data, eig_vals_data, host_task_events, depends); - heevd_ev = call_heevd(exec_q, jobz_val, uplo_val, n, a, w, host_task_events, depends); - } - else if ((eig_vecs_typenum == dpctl_capi.UAR_CFLOAT_) && (eig_vals_typenum == dpctl_capi.UAR_FLOAT_)) - { - std::complex* a = reinterpret_cast*>(eig_vecs.get_data()); - float* w = reinterpret_cast(eig_vals.get_data()); + sycl::event args_ev = dpctl::utils::keep_args_alive(exec_q, {eig_vecs, eig_vals}, host_task_events); + return std::make_pair(args_ev, heevd_ev); +} - heevd_ev = call_heevd(exec_q, jobz_val, uplo_val, n, a, w, host_task_events, depends); - } - else +template +struct HeevdContigFactory +{ + fnT get() { - throw py::value_error("Unexpected types of either eigenvectors or eigenvalues"); + if constexpr (types::HeevdTypePairSupportFactory::is_defined) + { + return heevd_impl; + } + else + { + return nullptr; + } } +}; - sycl::event args_ev = dpctl::utils::keep_args_alive(exec_q, {eig_vecs, eig_vals}, host_task_events); - return std::make_pair(args_ev, heevd_ev); +void init_heevd_dispatch_table(void) +{ + dpctl_td_ns::DispatchTableBuilder contig; + contig.populate_dispatch_table(heevd_dispatch_table); } } } diff --git a/dpnp/backend/extensions/lapack/heevd.hpp b/dpnp/backend/extensions/lapack/heevd.hpp index 93ce6fe560e1..85696d147f66 100644 --- a/dpnp/backend/extensions/lapack/heevd.hpp +++ b/dpnp/backend/extensions/lapack/heevd.hpp @@ -45,6 +45,8 @@ namespace lapack dpctl::tensor::usm_ndarray eig_vecs, dpctl::tensor::usm_ndarray eig_vals, const std::vector& depends); + + extern void init_heevd_dispatch_table(void); } } } diff --git a/dpnp/backend/extensions/lapack/lapack_py.cpp b/dpnp/backend/extensions/lapack/lapack_py.cpp index ea7506308032..eaa3e6873b6a 100644 --- a/dpnp/backend/extensions/lapack/lapack_py.cpp +++ b/dpnp/backend/extensions/lapack/lapack_py.cpp @@ -33,25 +33,45 @@ #include "heevd.hpp" #include "syevd.hpp" +namespace lapack_ext = dpnp::backend::ext::lapack; namespace py = pybind11; +// populate dispatch vectors +void init_dispatch_vectors(void) +{ + lapack_ext::init_syevd_dispatch_vector(); +} + +// populate dispatch tables +void init_dispatch_tables(void) +{ + lapack_ext::init_heevd_dispatch_table(); +} + PYBIND11_MODULE(_lapack_impl, m) { + init_dispatch_vectors(); + init_dispatch_tables(); + m.def("_heevd", - &dpnp::backend::ext::lapack::heevd, + &lapack_ext::heevd, "Call `heevd` from OneMKL LAPACK library to return " "the eigenvalues and eigenvectors of a complex Hermitian matrix", py::arg("sycl_queue"), - py::arg("jobz"), py::arg("upper_lower"), - py::arg("eig_vecs"), py::arg("eig_vals"), + py::arg("jobz"), + py::arg("upper_lower"), + py::arg("eig_vecs"), + py::arg("eig_vals"), py::arg("depends") = py::list()); m.def("_syevd", - &dpnp::backend::ext::lapack::syevd, + &lapack_ext::syevd, "Call `syevd` from OneMKL LAPACK library to return " "the eigenvalues and eigenvectors of a real symmetric matrix", py::arg("sycl_queue"), - py::arg("jobz"), py::arg("upper_lower"), - py::arg("eig_vecs"), py::arg("eig_vals"), + py::arg("jobz"), + py::arg("upper_lower"), + py::arg("eig_vecs"), + py::arg("eig_vals"), py::arg("depends") = py::list()); } diff --git a/dpnp/backend/extensions/lapack/syevd.cpp b/dpnp/backend/extensions/lapack/syevd.cpp index a4dded7543ab..d03c2dff372c 100644 --- a/dpnp/backend/extensions/lapack/syevd.cpp +++ b/dpnp/backend/extensions/lapack/syevd.cpp @@ -26,7 +26,12 @@ #include +// dpctl tensor headers +#include "utils/memory_overlap.hpp" +#include "utils/type_utils.hpp" + #include "syevd.hpp" +#include "types_matrix.hpp" #include "dpnp_utils.hpp" @@ -42,18 +47,33 @@ namespace lapack namespace mkl_lapack = oneapi::mkl::lapack; namespace py = pybind11; +namespace type_utils = dpctl::tensor::type_utils; + +typedef sycl::event (*syevd_impl_fn_ptr_t)(sycl::queue, + const oneapi::mkl::job, + const oneapi::mkl::uplo, + const std::int64_t, + char*, + char*, + std::vector&, + const std::vector&); + +static syevd_impl_fn_ptr_t syevd_dispatch_vector[dpctl_td_ns::num_types]; template -static sycl::event call_syevd(sycl::queue exec_q, +static sycl::event syevd_impl(sycl::queue exec_q, const oneapi::mkl::job jobz, const oneapi::mkl::uplo upper_lower, const std::int64_t n, - T* a, - T* w, + char* in_a, + char* out_w, std::vector& host_task_events, const std::vector& depends) { - validate_type_for_device(exec_q); + type_utils::validate_type_for_device(exec_q); + + T* a = reinterpret_cast(in_a); + T* w = reinterpret_cast(out_w); const std::int64_t lda = std::max(1UL, n); const std::int64_t scratchpad_size = mkl_lapack::syevd_scratchpad_size(exec_q, jobz, upper_lower, n, lda); @@ -162,13 +182,11 @@ std::pair syevd(sycl::queue exec_q, throw py::value_error("Execution queue is not compatible with allocation queues"); } - // check that arrays do not overlap, and concurrent access is safe. - // TODO: need to be exposed by DPCTL headers - // auto const& overlap = dpctl::tensor::overlap::MemoryOverlap(); - // if (overlap(eig_vecs, eig_vals)) - // { - // throw py::value_error("Arrays index overlapping segments of memory"); - // } + auto const& overlap = dpctl::tensor::overlap::MemoryOverlap(); + if (overlap(eig_vecs, eig_vals)) + { + throw py::value_error("Arrays with eigenvectors and eigenvalues are overlapping segments of memory"); + } bool is_eig_vecs_f_contig = eig_vecs.is_f_contiguous(); bool is_eig_vals_c_contig = eig_vals.is_c_contiguous(); @@ -181,43 +199,56 @@ std::pair syevd(sycl::queue exec_q, throw py::value_error("An array with output eigenvalues must be C-contiguous"); } - int eig_vecs_typenum = eig_vecs.get_typenum(); - int eig_vals_typenum = eig_vals.get_typenum(); - auto const& dpctl_capi = dpctl::detail::dpctl_capi::get(); + auto array_types = dpctl_td_ns::usm_ndarray_types(); + int eig_vecs_type_id = array_types.typenum_to_lookup_id(eig_vecs.get_typenum()); + int eig_vals_type_id = array_types.typenum_to_lookup_id(eig_vals.get_typenum()); - sycl::event syevd_ev; - std::vector host_task_events; + if (eig_vecs_type_id != eig_vals_type_id) + { + throw py::value_error("Types of eigenvectors and eigenvalues are missmatched"); + } + + syevd_impl_fn_ptr_t syevd_fn = syevd_dispatch_vector[eig_vecs_type_id]; + if (syevd_fn == nullptr) + { + throw py::value_error("No syevd implementation defined for a type of eigenvectors and eigenvalues"); + } + + char* eig_vecs_data = eig_vecs.get_data(); + char* eig_vals_data = eig_vals.get_data(); const std::int64_t n = eig_vecs_shape[0]; const oneapi::mkl::job jobz_val = static_cast(jobz); const oneapi::mkl::uplo uplo_val = static_cast(upper_lower); - if (eig_vecs_typenum != eig_vals_typenum) - { - throw py::value_error("Types of eigenvectors and eigenvalues aare missmatched"); - } - else if (eig_vecs_typenum == dpctl_capi.UAR_DOUBLE_) - { - double* a = reinterpret_cast(eig_vecs.get_data()); - double* w = reinterpret_cast(eig_vals.get_data()); + std::vector host_task_events; + sycl::event syevd_ev = + syevd_fn(exec_q, jobz_val, uplo_val, n, eig_vecs_data, eig_vals_data, host_task_events, depends); - syevd_ev = call_syevd(exec_q, jobz_val, uplo_val, n, a, w, host_task_events, depends); - } - else if (eig_vecs_typenum == dpctl_capi.UAR_FLOAT_) - { - float* a = reinterpret_cast(eig_vecs.get_data()); - float* w = reinterpret_cast(eig_vals.get_data()); + sycl::event args_ev = dpctl::utils::keep_args_alive(exec_q, {eig_vecs, eig_vals}, host_task_events); + return std::make_pair(args_ev, syevd_ev); +} - syevd_ev = call_syevd(exec_q, jobz_val, uplo_val, n, a, w, host_task_events, depends); - } - else +template +struct SyevdContigFactory +{ + fnT get() { - throw py::value_error("Unexpected types with num=" + std::to_string(eig_vecs_typenum) + - " for eigenvectors and eigenvalues"); + if constexpr (types::SyevdTypePairSupportFactory::is_defined) + { + return syevd_impl; + } + else + { + return nullptr; + } } +}; - sycl::event args_ev = dpctl::utils::keep_args_alive(exec_q, {eig_vecs, eig_vals}, host_task_events); - return std::make_pair(args_ev, syevd_ev); +void init_syevd_dispatch_vector(void) +{ + dpctl_td_ns::DispatchVectorBuilder contig; + contig.populate_dispatch_vector(syevd_dispatch_vector); } } } diff --git a/dpnp/backend/extensions/lapack/syevd.hpp b/dpnp/backend/extensions/lapack/syevd.hpp index 14d167ec02a7..c5f0bc1b1531 100644 --- a/dpnp/backend/extensions/lapack/syevd.hpp +++ b/dpnp/backend/extensions/lapack/syevd.hpp @@ -45,6 +45,8 @@ namespace lapack dpctl::tensor::usm_ndarray eig_vecs, dpctl::tensor::usm_ndarray eig_vals, const std::vector& depends = {}); + + extern void init_syevd_dispatch_vector(void); } } } diff --git a/dpnp/backend/extensions/lapack/types_matrix.hpp b/dpnp/backend/extensions/lapack/types_matrix.hpp new file mode 100644 index 000000000000..4175873b541f --- /dev/null +++ b/dpnp/backend/extensions/lapack/types_matrix.hpp @@ -0,0 +1,80 @@ +//***************************************************************************** +// Copyright (c) 2023, Intel Corporation +// All rights reserved. +// +// Redistribution and use in source and binary forms, with or without +// modification, are permitted provided that the following conditions are met: +// - Redistributions of source code must retain the above copyright notice, +// this list of conditions and the following disclaimer. +// - Redistributions in binary form must reproduce the above copyright notice, +// this list of conditions and the following disclaimer in the documentation +// and/or other materials provided with the distribution. +// +// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" +// AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE +// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE +// ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE +// LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR +// CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF +// SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS +// INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN +// CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) +// ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF +// THE POSSIBILITY OF SUCH DAMAGE. +//***************************************************************************** + +#pragma once + +#include + +// dpctl tensor headers +#include "utils/type_dispatch.hpp" + +// dpctl namespace for operations with types +namespace dpctl_td_ns = dpctl::tensor::type_dispatch; + +namespace dpnp +{ +namespace backend +{ +namespace ext +{ +namespace lapack +{ +namespace types +{ +/** + * @brief A factory to define pairs of supported types for which + * MKL LAPACK library provides support in oneapi::mkl::lapack::heevd function. + * + * @tparam T Type of array containing input matrix A and an output array with eigenvectors. + * @tparam RealT Type of output array containing eigenvalues of A. + */ +template +struct HeevdTypePairSupportFactory +{ + static constexpr bool is_defined = std::disjunction, RealT, double>, + dpctl_td_ns::TypePairDefinedEntry, RealT, float>, + // fall-through + dpctl_td_ns::NotDefinedEntry>::is_defined; +}; + +/** + * @brief A factory to define pairs of supported types for which + * MKL LAPACK library provides support in oneapi::mkl::lapack::syevd function. + * + * @tparam T Type of array containing input matrix A and an output arrays with eigenvectors and eigenvectors. + */ +template +struct SyevdTypePairSupportFactory +{ + static constexpr bool is_defined = std::disjunction, + dpctl_td_ns::TypePairDefinedEntry, + // fall-through + dpctl_td_ns::NotDefinedEntry>::is_defined; +}; +} +} +} +} +}