Merge 10df422 into 92faa20

Author: vlad-perevezentsev

dpnp/backend/extensions/lapack/CMakeLists.txt

@@ -32,6 +32,7 @@ set(_module_src
     ${CMAKE_CURRENT_SOURCE_DIR}/gesv.cpp
     ${CMAKE_CURRENT_SOURCE_DIR}/gesv_batch.cpp
     ${CMAKE_CURRENT_SOURCE_DIR}/gesvd.cpp
+    ${CMAKE_CURRENT_SOURCE_DIR}/gesvd_batch.cpp
     ${CMAKE_CURRENT_SOURCE_DIR}/getrf.cpp
     ${CMAKE_CURRENT_SOURCE_DIR}/getrf_batch.cpp
     ${CMAKE_CURRENT_SOURCE_DIR}/getri_batch.cpp

dpnp/backend/extensions/lapack/gesvd.cpp

@@ -26,21 +26,20 @@
 #include <pybind11/pybind11.h>
 
 // dpctl tensor headers
-#include "utils/memory_overlap.hpp"
 #include "utils/type_utils.hpp"
 
+#include "common_helpers.hpp"
 #include "gesvd.hpp"
+#include "gesvd_common_utils.hpp"
 #include "types_matrix.hpp"
 
-#include "dpnp_utils.hpp"
-
 namespace dpnp::extensions::lapack
 {
 namespace mkl_lapack = oneapi::mkl::lapack;
 namespace py = pybind11;
 namespace type_utils = dpctl::tensor::type_utils;
 
-typedef sycl::event (*gesvd_impl_fn_ptr_t)(sycl::queue,
+typedef sycl::event (*gesvd_impl_fn_ptr_t)(sycl::queue &,
                                            const oneapi::mkl::jobsvd,
                                            const oneapi::mkl::jobsvd,
                                            const std::int64_t,
@@ -58,26 +57,8 @@ typedef sycl::event (*gesvd_impl_fn_ptr_t)(sycl::queue,
 static gesvd_impl_fn_ptr_t gesvd_dispatch_table[dpctl_td_ns::num_types]
                                                [dpctl_td_ns::num_types];
 
-// Converts a given character code (ord) to the corresponding
-// oneapi::mkl::jobsvd enumeration value
-static oneapi::mkl::jobsvd process_job(std::int8_t job_val)
-{
-    switch (job_val) {
-    case 'A':
-        return oneapi::mkl::jobsvd::vectors;
-    case 'S':
-        return oneapi::mkl::jobsvd::somevec;
-    case 'O':
-        return oneapi::mkl::jobsvd::vectorsina;
-    case 'N':
-        return oneapi::mkl::jobsvd::novec;
-    default:
-        throw std::invalid_argument("Unknown value for job");
-    }
-}
-
 template <typename T, typename RealT>
-static sycl::event gesvd_impl(sycl::queue exec_q,
+static sycl::event gesvd_impl(sycl::queue &exec_q,
                               const oneapi::mkl::jobsvd jobu,
                               const oneapi::mkl::jobsvd jobvt,
                               const std::int64_t m,
@@ -102,16 +83,14 @@ static sycl::event gesvd_impl(sycl::queue exec_q,
 
     const std::int64_t scratchpad_size = mkl_lapack::gesvd_scratchpad_size<T>(
         exec_q, jobu, jobvt, m, n, lda, ldu, ldvt);
-    T *scratchpad = nullptr;
+
+    T *scratchpad = helper::alloc_scratchpad<T>(scratchpad_size, exec_q);
 
     std::stringstream error_msg;
-    std::int64_t info = 0;
     bool is_exception_caught = false;
 
     sycl::event gesvd_event;
     try {
-        scratchpad = sycl::malloc_device<T>(scratchpad_size, exec_q);
-
         gesvd_event = mkl_lapack::gesvd(
             exec_q,
             jobu,  // Character specifying how to compute the matrix U:
@@ -138,26 +117,7 @@ static sycl::event gesvd_impl(sycl::queue exec_q,
             scratchpad_size, depends);
     } catch (mkl_lapack::exception const &e) {
         is_exception_caught = true;
-        info = e.info();
-        if (info < 0) {
-            error_msg << "Parameter number " << -info
-                      << " had an illegal value.";
-        }
-        else if (info == scratchpad_size && e.detail() != 0) {
-            error_msg
-                << "Insufficient scratchpad size. Required size is at least "
-                << e.detail();
-        }
-        else if (info > 0) {
-            error_msg << "The algorithm computing SVD failed to converge; "
-                      << info << " off-diagonal elements of an intermediate "
-                      << "bidiagonal form did not converge to zero.\n";
-        }
-        else {
-            error_msg << "Unexpected MKL exception caught during gesvd() "
-                         "call:\nreason: "
-                      << e.what() << "\ninfo: " << e.info();
-        }
+        gesvd_utils::handle_lapack_exc(scratchpad_size, e, error_msg);
     } catch (sycl::exception const &e) {
         is_exception_caught = true;
         error_msg << "Unexpected SYCL exception caught during gesvd() call:\n"
@@ -182,7 +142,7 @@ static sycl::event gesvd_impl(sycl::queue exec_q,
 }
 
 std::pair<sycl::event, sycl::event>
-    gesvd(sycl::queue exec_q,
+    gesvd(sycl::queue &exec_q,
           const std::int8_t jobu_val,
           const std::int8_t jobvt_val,
           dpctl::tensor::usm_ndarray a_array,
@@ -191,103 +151,26 @@ std::pair<sycl::event, sycl::event>
           dpctl::tensor::usm_ndarray out_vt,
           const std::vector<sycl::event> &depends)
 {
-    const int a_array_nd = a_array.get_ndim();
-    const int out_u_array_nd = out_u.get_ndim();
-    const int out_s_array_nd = out_s.get_ndim();
-    const int out_vt_array_nd = out_vt.get_ndim();
-
-    if (a_array_nd != 2) {
-        throw py::value_error(
-            "The input array has ndim=" + std::to_string(a_array_nd) +
-            ", but a 2-dimensional array is expected.");
-    }
-
-    if (out_s_array_nd != 1) {
-        throw py::value_error("The output array of singular values has ndim=" +
-                              std::to_string(out_s_array_nd) +
-                              ", but a 1-dimensional array is expected.");
-    }
-
-    if (jobu_val == 'N' && jobvt_val == 'N') {
-        if (out_u_array_nd != 0) {
-            throw py::value_error(
-                "The output array of the left singular vectors has ndim=" +
-                std::to_string(out_u_array_nd) +
-                ", but it is not used and should have ndim=0.");
-        }
-        if (out_vt_array_nd != 0) {
-            throw py::value_error(
-                "The output array of the right singular vectors has ndim=" +
-                std::to_string(out_vt_array_nd) +
-                ", but it is not used and should have ndim=0.");
-        }
-    }
-    else {
-        if (out_u_array_nd != 2) {
-            throw py::value_error(
-                "The output array of the left singular vectors has ndim=" +
-                std::to_string(out_u_array_nd) +
-                ", but a 2-dimensional array is expected.");
-        }
-        if (out_vt_array_nd != 2) {
-            throw py::value_error(
-                "The output array of the right singular vectors has ndim=" +
-                std::to_string(out_vt_array_nd) +
-                ", but a 2-dimensional array is expected.");
-        }
-    }
-
-    // check compatibility of execution queue and allocation queue
-    if (!dpctl::utils::queues_are_compatible(
-            exec_q, {a_array.get_queue(), out_s.get_queue(), out_u.get_queue(),
-                     out_vt.get_queue()}))
-    {
-        throw std::runtime_error(
-            "USM allocations are not compatible with the execution queue.");
-    }
-
-    auto const &overlap = dpctl::tensor::overlap::MemoryOverlap();
-    if (overlap(a_array, out_s) || overlap(a_array, out_u) ||
-        overlap(a_array, out_vt) || overlap(out_s, out_u) ||
-        overlap(out_s, out_vt) || overlap(out_u, out_vt))
-    {
-        throw py::value_error("Arrays have overlapping segments of memory");
-    }
-
-    bool is_a_array_f_contig = a_array.is_f_contiguous();
-    if (!is_a_array_f_contig) {
-        throw py::value_error("The input array must be F-contiguous");
-    }
-
-    bool is_out_u_array_f_contig = out_u.is_f_contiguous();
-    bool is_out_vt_array_f_contig = out_vt.is_f_contiguous();
-
-    if (!is_out_u_array_f_contig || !is_out_vt_array_f_contig) {
-        throw py::value_error("The output arrays of the left and right "
-                              "singular vectors must be F-contiguous");
-    }
-
-    bool is_out_s_array_c_contig = out_s.is_c_contiguous();
-    bool is_out_s_array_f_contig = out_s.is_f_contiguous();
-
-    if (!is_out_s_array_c_contig || !is_out_s_array_f_contig) {
-        throw py::value_error("The output array of singular values "
-                              "must be contiguous");
-    }
+    constexpr int expected_a_u_vt_ndim = 2;
+    constexpr int expected_s_ndim = 1;
+
+    gesvd_utils::common_gesvd_checks(exec_q, a_array, out_s, out_u, out_vt,
+                                     jobu_val, jobvt_val, expected_a_u_vt_ndim,
+                                     expected_s_ndim);
+
+    // // Ensure `m` and 'n' are non-zero, otherwise return empty
+    // // events
+    // if (gesvd_utils::check_zeros_shape_gesvd(a_array, out_s, out_u, out_vt,
+    //                                          jobu_val, jobvt_val))
+    // {
+    //     // nothing to do
+    //     return std::make_pair(sycl::event(), sycl::event());
+    // }
 
     auto array_types = dpctl_td_ns::usm_ndarray_types();
     int a_array_type_id =
         array_types.typenum_to_lookup_id(a_array.get_typenum());
-    int out_u_type_id = array_types.typenum_to_lookup_id(out_u.get_typenum());
     int out_s_type_id = array_types.typenum_to_lookup_id(out_s.get_typenum());
-    int out_vt_type_id = array_types.typenum_to_lookup_id(out_vt.get_typenum());
-
-    if (a_array_type_id != out_u_type_id || a_array_type_id != out_vt_type_id) {
-        throw py::type_error(
-            "Input array, output left singular vectors array, "
-            "and outpuy right singular vectors array must have "
-            "the same data type");
-    }
 
     gesvd_impl_fn_ptr_t gesvd_fn =
         gesvd_dispatch_table[a_array_type_id][out_s_type_id];
@@ -311,8 +194,8 @@ std::pair<sycl::event, sycl::event>
     const std::int64_t ldvt =
         std::max<std::size_t>(1UL, jobvt_val == 'S' ? (m > n ? n : m) : n);
 
-    const oneapi::mkl::jobsvd jobu = process_job(jobu_val);
-    const oneapi::mkl::jobsvd jobvt = process_job(jobvt_val);
+    const oneapi::mkl::jobsvd jobu = gesvd_utils::process_job(jobu_val);
+    const oneapi::mkl::jobsvd jobvt = gesvd_utils::process_job(jobvt_val);
 
     std::vector<sycl::event> host_task_events;
     sycl::event gesvd_ev =

dpnp/backend/extensions/lapack/gesvd.hpp

@@ -33,7 +33,7 @@
 namespace dpnp::extensions::lapack
 {
 extern std::pair<sycl::event, sycl::event>
-    gesvd(sycl::queue exec_q,
+    gesvd(sycl::queue &exec_q,
           const std::int8_t jobu_val,
           const std::int8_t jobvt_val,
           dpctl::tensor::usm_ndarray a_array,
@@ -42,5 +42,16 @@ extern std::pair<sycl::event, sycl::event>
           dpctl::tensor::usm_ndarray out_vt,
           const std::vector<sycl::event> &depends);
 
+extern std::pair<sycl::event, sycl::event>
+    gesvd_batch(sycl::queue &exec_q,
+                const std::int8_t jobu_val,
+                const std::int8_t jobvt_val,
+                dpctl::tensor::usm_ndarray a_array,
+                dpctl::tensor::usm_ndarray out_s,
+                dpctl::tensor::usm_ndarray out_u,
+                dpctl::tensor::usm_ndarray out_vt,
+                const std::vector<sycl::event> &depends);
+
 extern void init_gesvd_dispatch_table(void);
+extern void init_gesvd_batch_dispatch_table(void);
 } // namespace dpnp::extensions::lapack