bench-gemv.cc 4.8 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137
  1. // -*- mode: c++; coding: utf-8 -*-
  2. /// @file bench-gemv.hh
  3. /// @brief Benchmark for BLAS-2 type ops
  4. // (c) Daniel Llorens - 2017
  5. // This library is free software; you can redistribute it and/or modify it under
  6. // the terms of the GNU Lesser General Public License as published by the Free
  7. // Software Foundation; either version 3 of the License, or (at your option) any
  8. // later version.
  9. // These operations aren't really part of the ET framework, just standalone
  10. // functions.
  11. // Cf bench-gemm.cc for BLAS-3 type ops.
  12. #include <iostream>
  13. #include <iomanip>
  14. #include "ra/test.hh"
  15. #include "ra/complex.hh"
  16. #include "ra/ra.hh"
  17. #include "ra/bench.hh"
  18. using std::cout, std::endl, std::setw, std::setprecision, ra::TestRecorder;
  19. using ra::Small, ra::View, ra::Unique;
  20. using real = double;
  21. // -------------------
  22. // variants of the defaults, should be slower if the default is well picked.
  23. // TODO compare with external GEMV/GEVM
  24. // -------------------
  25. enum trans_t { NOTRANS, TRANS };
  26. int main()
  27. {
  28. TestRecorder tr(std::cout);
  29. auto gemv_i = [&](auto const & a, auto const & b)
  30. {
  31. int const M = a.size(0);
  32. ra::Big<decltype(a(0, 0)*b(0)), 1> c({M}, ra::none);
  33. for (int i=0; i<M; ++i) {
  34. c(i) = dot(a(i), b);
  35. }
  36. return c;
  37. };
  38. auto gemv_j = [&](auto const & a, auto const & b)
  39. {
  40. int const M = a.size(0);
  41. int const N = a.size(1);
  42. ra::Big<decltype(a(0, 0)*b(0)), 1> c({M}, 0.);
  43. for (int j=0; j<N; ++j) {
  44. c += a(ra::all, j)*b(j);
  45. }
  46. return c;
  47. };
  48. auto gevm_j = [&](auto const & b, auto const & a)
  49. {
  50. int const N = a.size(1);
  51. ra::Big<decltype(b(0)*a(0, 0)), 1> c({N}, ra::none);
  52. for (int j=0; j<N; ++j) {
  53. c(j) = dot(b, a(ra::all, j));
  54. }
  55. return c;
  56. };
  57. auto gevm_i = [&](auto const & b, auto const & a)
  58. {
  59. int const M = a.size(0);
  60. int const N = a.size(1);
  61. ra::Big<decltype(b(0)*a(0, 0)), 1> c({N}, 0.);
  62. for (int i=0; i<M; ++i) {
  63. c += b(i)*a(i);
  64. }
  65. return c;
  66. };
  67. auto bench_all = [&](int k, int m, int n, int reps)
  68. {
  69. auto bench_mv = [&tr, &m, &n, &reps](auto && f, char const * tag, trans_t t)
  70. {
  71. ra::Big<real, 2> aa({m, n}, ra::_0-ra::_1);
  72. auto a = t==TRANS ? transpose<1, 0>(aa) : aa();
  73. ra::Big<real, 1> b({a.size(1)}, 1-2*ra::_0);
  74. ra::Big<real, 1> ref = gemv(a, b);
  75. ra::Big<real, 1> c;
  76. auto bv = Benchmark().repeats(reps).runs(3).run([&]() { c = f(a, b); });
  77. tr.info(std::setw(5), std::fixed, Benchmark::avg(bv)/(m*n)/1e-9, " ns [",
  78. Benchmark::stddev(bv)/(m*n)/1e-9 ,"] ", tag, t==TRANS ? " [T]" : " [N]").test_eq(ref, c);
  79. };
  80. auto bench_vm = [&tr, &m, &n, &reps](auto && f, char const * tag, trans_t t)
  81. {
  82. ra::Big<real, 2> aa({m, n}, ra::_0-ra::_1);
  83. auto a = t==TRANS ? transpose<1, 0>(aa) : aa();
  84. ra::Big<real, 1> b({a.size(0)}, 1-2*ra::_0);
  85. ra::Big<real, 1> ref = gevm(b, a);
  86. ra::Big<real, 1> c;
  87. auto bv = Benchmark().repeats(reps).runs(4).run([&]() { c = f(b, a); });
  88. tr.info(std::setw(5), std::fixed, Benchmark::avg(bv)/(m*n)/1e-9, " ns [",
  89. Benchmark::stddev(bv)/(m*n)/1e-9 ,"] ", tag, t==TRANS ? " [T]" : " [N]").test_eq(ref, c);
  90. };
  91. tr.section(m, " x ", n, " times ", reps);
  92. // some variants are way too slow to check with larger arrays.
  93. if (k>0) {
  94. bench_mv(gemv_i, "mv i", NOTRANS);
  95. bench_mv(gemv_i, "mv i", TRANS);
  96. bench_mv(gemv_j, "mv j", NOTRANS);
  97. bench_mv(gemv_j, "mv j", TRANS);
  98. bench_mv([&](auto const & a, auto const & b) { return gemv(a, b); }, "mv default", NOTRANS);
  99. bench_mv([&](auto const & a, auto const & b) { return gemv(a, b); }, "mv default", TRANS);
  100. bench_vm(gevm_i, "vm i", NOTRANS);
  101. bench_vm(gevm_i, "vm i", TRANS);
  102. bench_vm(gevm_j, "vm j", NOTRANS);
  103. bench_vm(gevm_j, "vm j", TRANS);
  104. bench_vm([&](auto const & a, auto const & b) { return gevm(a, b); }, "vm default", NOTRANS);
  105. bench_vm([&](auto const & a, auto const & b) { return gevm(a, b); }, "vm default", TRANS);
  106. }
  107. };
  108. bench_all(3, 10, 10, 10000);
  109. bench_all(3, 100, 100, 100);
  110. bench_all(3, 500, 500, 1);
  111. bench_all(3, 10000, 1000, 1);
  112. bench_all(3, 1000, 10000, 1);
  113. bench_all(3, 100000, 100, 1);
  114. bench_all(3, 100, 100000, 1);
  115. return tr.summary();
  116. }