1
2#include <iostream>
3#include <Eigen/Core>
4#include <bench/BenchTimer.h>
5
6using namespace Eigen;
7using namespace std;
8
9#define END 9
10
11template<int S> struct map_size { enum { ret = S }; };
12template<>  struct map_size<10> { enum { ret = 20 }; };
13template<>  struct map_size<11> { enum { ret = 50 }; };
14template<>  struct map_size<12> { enum { ret = 100 }; };
15template<>  struct map_size<13> { enum { ret = 300 }; };
16
17template<int M, int N,int K> struct alt_prod
18{
19  enum {
20    ret = M==1 && N==1 ? InnerProduct
21        : K==1 ? OuterProduct
22        : M==1 ? GemvProduct
23        : N==1 ? GemvProduct
24        : GemmProduct
25  };
26};
27
28void print_mode(int mode)
29{
30  if(mode==InnerProduct) std::cout << "i";
31  if(mode==OuterProduct) std::cout << "o";
32  if(mode==CoeffBasedProductMode) std::cout << "c";
33  if(mode==LazyCoeffBasedProductMode) std::cout << "l";
34  if(mode==GemvProduct) std::cout << "v";
35  if(mode==GemmProduct) std::cout << "m";
36}
37
38template<int Mode, typename Lhs, typename Rhs, typename Res>
39EIGEN_DONT_INLINE void prod(const Lhs& a, const Rhs& b, Res& c)
40{
41  c.noalias() += typename ProductReturnType<Lhs,Rhs,Mode>::Type(a,b);
42}
43
44template<int M, int N, int K, typename Scalar, int Mode>
45EIGEN_DONT_INLINE void bench_prod()
46{
47  typedef Matrix<Scalar,M,K> Lhs; Lhs a; a.setRandom();
48  typedef Matrix<Scalar,K,N> Rhs; Rhs b; b.setRandom();
49  typedef Matrix<Scalar,M,N> Res; Res c; c.setRandom();
50
51  BenchTimer t;
52  double n = 2.*double(M)*double(N)*double(K);
53  int rep = 100000./n;
54  rep /= 2;
55  if(rep<1) rep = 1;
56  do {
57    rep *= 2;
58    t.reset();
59    BENCH(t,1,rep,prod<CoeffBasedProductMode>(a,b,c));
60  } while(t.best()<0.1);
61
62  t.reset();
63  BENCH(t,5,rep,prod<Mode>(a,b,c));
64
65  print_mode(Mode);
66  std::cout << int(1e-6*n*rep/t.best()) << "\t";
67}
68
69template<int N> struct print_n;
70template<int M, int N, int K> struct loop_on_m;
71template<int M, int N, int K, typename Scalar, int Mode> struct loop_on_n;
72
73template<int M, int N, int K>
74struct loop_on_k
75{
76  static void run()
77  {
78    std::cout << "K=" << K << "\t";
79    print_n<N>::run();
80    std::cout << "\n";
81
82    loop_on_m<M,N,K>::run();
83    std::cout << "\n\n";
84
85    loop_on_k<M,N,K+1>::run();
86  }
87};
88
89template<int M, int N>
90struct loop_on_k<M,N,END> { static void run(){} };
91
92
93template<int M, int N, int K>
94struct loop_on_m
95{
96  static void run()
97  {
98    std::cout << M << "f\t";
99    loop_on_n<M,N,K,float,CoeffBasedProductMode>::run();
100    std::cout << "\n";
101
102    std::cout << M << "f\t";
103    loop_on_n<M,N,K,float,-1>::run();
104    std::cout << "\n";
105
106    loop_on_m<M+1,N,K>::run();
107  }
108};
109
110template<int N, int K>
111struct loop_on_m<END,N,K> { static void run(){} };
112
113template<int M, int N, int K, typename Scalar, int Mode>
114struct loop_on_n
115{
116  static void run()
117  {
118    bench_prod<M,N,K,Scalar,Mode==-1? alt_prod<M,N,K>::ret : Mode>();
119
120    loop_on_n<M,N+1,K,Scalar,Mode>::run();
121  }
122};
123
124template<int M, int K, typename Scalar, int Mode>
125struct loop_on_n<M,END,K,Scalar,Mode> { static void run(){} };
126
127template<int N> struct print_n
128{
129  static void run()
130  {
131    std::cout << map_size<N>::ret << "\t";
132    print_n<N+1>::run();
133  }
134};
135
136template<> struct print_n<END> { static void run(){} };
137
138int main()
139{
140  loop_on_k<1,1,1>::run();
141
142  return 0;
143}
144