dot, dotc: add customization points and tests kokkos#96

Author: fnrizzi

examples/kokkos-based/CMakeLists.txt

@@ -1,5 +1,6 @@
 
-linalg_add_example(dot_kokkos)
 linalg_add_example(add_kokkos)
+linalg_add_example(dot_kokkos)
+linalg_add_example(dotc_kokkos)
 linalg_add_example(simple_scale_kokkos)
 linalg_add_example(matrix_vector_product_kokkos)

examples/kokkos-based/add_kokkos.cpp

@@ -15,34 +15,37 @@ void print_elements(const T1 & v, const std::vector<ScalarType> & gold)
 
 int main(int argc, char* argv[])
 {
-  std::cout << "running add example calling custom kokkos" << std::endl;
-  int N = 50;
+  std::cout << "add example: calling kokkos-kernels" << std::endl;
+
+  std::size_t N = 50;
   Kokkos::initialize(argc,argv);
   {
-    Kokkos::View<double*> x_view("x",N);
-    Kokkos::View<double*> y_view("y",N);
-    Kokkos::View<double*> z_view("z",N);
+    using value_type = double;
+
+    Kokkos::View<value_type*> x_view("x",N);
+    Kokkos::View<value_type*> y_view("y",N);
+    Kokkos::View<value_type*> z_view("z",N);
 
-    double* x_ptr = x_view.data();
-    double* y_ptr = y_view.data();
-    double* z_ptr = z_view.data();
+    value_type* x_ptr = x_view.data();
+    value_type* y_ptr = y_view.data();
+    value_type* z_ptr = z_view.data();
 
     using dyn_1d_ext_type = std::experimental::extents<std::experimental::dynamic_extent>;
-    using mdspan_type  = std::experimental::mdspan<double, dyn_1d_ext_type>;
+    using mdspan_type  = std::experimental::mdspan<value_type, dyn_1d_ext_type>;
     mdspan_type x(x_ptr,N);
     mdspan_type y(y_ptr,N);
     mdspan_type z(z_ptr,N);
 
-    std::vector<double> gold(N);
+    std::vector<value_type> gold(N);
     for(int i=0; i<x.extent(0); i++){
       x(i) = i;
-      y(i) = i + (double)10;
+      y(i) = i + (value_type)10;
       z(i) = 0;
       gold[i] = x(i) + y(i);
     }
 
     namespace stdla = std::experimental::linalg;
-    const double init_value = 2.0;
+    const value_type init_value = 2.0;
 
     {
       // This goes to the base implementation
@@ -51,7 +54,7 @@ int main(int argc, char* argv[])
 
     {
       // reset z since it is modified above
-      for(int i=0; i<z.extent(0); i++){ z(i) = 0; }
+      for(std::size_t i=0; i<z.extent(0); i++){ z(i) = 0; }
 
       // This forwards to KokkosKernels
       stdla::add(KokkosKernelsSTD::kokkos_exec<>(), x,y,z);

examples/kokkos-based/dot_kokkos.cpp

@@ -4,35 +4,37 @@
 
 int main(int argc, char* argv[])
 {
-  std::cout << "running dot example calling custom kokkos" << std::endl;
-  int N = 50;
+  std::cout << "dot example: calling kokkos-kernels" << std::endl;
+
+  std::size_t N = 50;
   Kokkos::initialize(argc,argv);
   {
-    Kokkos::View<double*> a_view("A",N);
-    Kokkos::View<double*> b_view("B",N);
-    double* a_ptr = a_view.data();
-    double* b_ptr = b_view.data();
+    using value_type = double;
+
+    Kokkos::View<value_type*> a_view("A",N);
+    Kokkos::View<value_type*> b_view("B",N);
+    value_type* a_ptr = a_view.data();
+    value_type* b_ptr = b_view.data();
 
     using dyn_1d_ext_type = std::experimental::extents<std::experimental::dynamic_extent>;
-    using mdspan_type  = std::experimental::mdspan<double, dyn_1d_ext_type>;
+    using mdspan_type  = std::experimental::mdspan<value_type, dyn_1d_ext_type>;
     mdspan_type a(a_ptr,N);
     mdspan_type b(b_ptr,N);
-    for(int i=0; i<a.extent(0); i++){
+    for(std::size_t i=0; i<a.extent(0); i++){
       a(i) = i;
       b(i) = i;
     }
 
     namespace stdla = std::experimental::linalg;
-    const double init_value = 2.0;
+    const value_type init_value(2.0);
 
     // This goes to the base implementation
     const auto res_seq = stdla::dot(std::execution::seq, a, b, init_value);
+    printf("Seq result    = %lf\n", res_seq);
 
     // This forwards to KokkosKernels
     auto res_kk = stdla::dot(KokkosKernelsSTD::kokkos_exec<>(), a, b, init_value);
-
     printf("Kokkos result = %lf\n", res_kk);
-    printf("Seq result    = %lf\n", res_seq);
   }
   Kokkos::finalize();
 }

examples/kokkos-based/dotc_kokkos.cpp

@@ -0,0 +1,42 @@
+
+#include <experimental/linalg>
+#include <iostream>
+
+int main(int argc, char* argv[])
+{
+  std::cout << "dotc example: calling kokkos-kernels" << std::endl;
+
+  std::size_t N = 10;
+  Kokkos::initialize(argc,argv);
+  {
+    using value_type = std::complex<double>;
+    using view_t = Kokkos::View<value_type*>;
+    view_t a_view("A",N);
+    view_t b_view("B",N);
+    value_type* a_ptr = a_view.data();
+    value_type* b_ptr = b_view.data();
+
+    using dyn_1d_ext_type = std::experimental::extents<std::experimental::dynamic_extent>;
+    using mdspan_type  = std::experimental::mdspan<value_type, dyn_1d_ext_type>;
+    mdspan_type a(a_ptr,N);
+    mdspan_type b(b_ptr,N);
+    for(std::size_t i=0; i<a.extent(0); i++){
+      const value_type a_i(double(i) + 1.0, double(i) + 1.0);
+      const value_type b_i(double(i) - 2.0, double(i) - 2.0);
+      a(i) = a_i;
+      b(i) = b_i;
+    }
+
+    namespace stdla = std::experimental::linalg;
+    const value_type init_value(2., 3.);
+
+    // This goes to the base implementation
+    const auto res_seq = stdla::dotc(std::execution::seq, a, b, init_value);
+    std::cout << "Seq result    = " << res_seq << "\n";
+
+    // This forwards to KokkosKernels
+    auto res_kk = stdla::dotc(KokkosKernelsSTD::kokkos_exec<>(), a, b, init_value);
+    std::cout << "Kokkos result = " << res_kk << "\n";
+  }
+  Kokkos::finalize();
+}

include/experimental/__p1673_bits/blas1_dot.hpp

@@ -78,12 +78,6 @@ struct is_custom_dot_avail<
 
 } // end anonymous namespace
 
-
-// ------------
-// PUBLIC API:
-// ------------
-
-// dot, with init value
 template<class ExecutionPolicy,
          class ElementType1,
          extents<>::size_type ext1,
@@ -100,6 +94,9 @@ Scalar dot(
   std::experimental::mdspan<ElementType2, std::experimental::extents<ext2>, Layout2, Accessor2> v2,
   Scalar init)
 {
+  static_assert(v1.static_extent(0) == dynamic_extent ||
+                v2.static_extent(0) == dynamic_extent ||
+                v1.static_extent(0) == v2.static_extent(0));
 
   constexpr bool use_custom = is_custom_dot_avail<
     decltype(execpolicy_mapper(exec)), decltype(v1), decltype(v2), Scalar
@@ -132,7 +129,6 @@ Scalar dot(std::experimental::mdspan<ElementType1, std::experimental::extents<ex
   return dot(std::experimental::linalg::impl::default_exec_t(), v1, v2, init);
 }
 
-// Conjugated dot, with init value
 template<class ElementType1,
          extents<>::size_type ext1,
          class Layout1,
@@ -150,7 +146,6 @@ Scalar dotc(
   return dot(conjugated(v1), v2, init);
 }
 
-// conjugated dot: with policy, with init value
 template<class ExecutionPolicy,
          class ElementType1,
          extents<>::size_type ext1,
@@ -162,12 +157,12 @@ template<class ExecutionPolicy,
          class Accessor2,
          class Scalar>
 Scalar dotc(
-  ExecutionPolicy&& /* exec */,
+  ExecutionPolicy&& exec,
   std::experimental::mdspan<ElementType1, std::experimental::extents<ext1>, Layout1, Accessor1> v1,
   std::experimental::mdspan<ElementType2, std::experimental::extents<ext2>, Layout2, Accessor2> v2,
   Scalar init)
 {
-  return dotc(v1, v2, init);
+  return dot(exec, conjugated(v1), v2, init);
 }
 
 namespace dot_detail {
@@ -190,7 +185,6 @@ namespace dot_detail {
   -> decltype(x(0) * y(0));
 } // namespace dot_detail
 
-// dot, without init value
 template<class ElementType1,
          extents<>::size_type ext1,
          class Layout1,
@@ -240,7 +234,7 @@ auto dotc(
   std::experimental::mdspan<ElementType2, std::experimental::extents<ext2>, Layout2, Accessor2> v2)
   -> decltype(dot_detail::dot_return_type_deducer(conjugated(v1), v2))
 {
-  using return_t = decltype(dot_detail::dot_return_type_deducer(v1, v2));
+  using return_t = decltype(dot_detail::dot_return_type_deducer(conjugated(v1), v2));
   return dotc(v1, v2, return_t{});
 }