Eigen benchmark for 3D vector multiplication

benchmark.cpp

#include <eigen3/Eigen/StdVector>
#include <iostream>
#include <random>
#include <vector>

#include <benchmark/benchmark.h>

template<typename Vector_>
using AlignedVector = std::vector<Vector_, Eigen::aligned_allocator<Vector_>>;

template<typename V_, int N_>
static void BM_VectorMulti(benchmark::State& state) {
  size_t s = state.range(0);

  using Vec_ = Eigen::Matrix<V_, N_, 1>;

  auto vec_a = AlignedVector<Vec_>(s);
  auto vec_b = AlignedVector<Vec_>(s);
  auto vec_c = std::vector<V_>(s);

  auto rng = std::mt19937 {std::random_device{}()};
  auto dist = std::uniform_real_distribution<V_>{ -1.0, 1.0 };

  auto vecgen = [&dist, &rng]() -> Vec_ {
    if constexpr(N_ == 4) {
      return { dist(rng), dist(rng), dist(rng), 0 };
    } else if constexpr(N_ == 3) {
      return { dist(rng), dist(rng), dist(rng) };
    } else if constexpr(N_ == Eigen::Dynamic) {
      Vec_ v; 
      v.resize(3, 1);
      v << dist(rng), dist(rng), dist(rng);
      return v;
    }
  };

  std::generate(vec_a.begin(), vec_a.end(), vecgen);
  std::generate(vec_b.begin(), vec_b.end(), vecgen);

  for (auto _ : state) {
    // This code gets timed
    //
    // The benchmark is to do a vector multiplication of element i in 
    // vector A with element in vector B, then rotate vector B one element
    // and repeat. So, we will have N^2 vector multiplications
    for (size_t i = 0; i < s; i++) {
      for (size_t j = 0; j < s; j++) {
        asm("# here is Eigen math");
        vec_c[j] = vec_a[(j+1)%s].transpose() * vec_b[j];
      }

      benchmark::DoNotOptimize(vec_c);
    }

    benchmark::ClobberMemory();
  }
}

// Register the function as a benchmark
BENCHMARK_TEMPLATE(BM_VectorMulti, float, 4)->RangeMultiplier(4)->Range(16, 256);
BENCHMARK_TEMPLATE(BM_VectorMulti, double, 4)->RangeMultiplier(4)->Range(16, 256);
BENCHMARK_TEMPLATE(BM_VectorMulti, float, 3)->RangeMultiplier(4)->Range(16, 256);
BENCHMARK_TEMPLATE(BM_VectorMulti, double, 3)->RangeMultiplier(4)->Range(16, 256);
BENCHMARK_TEMPLATE(BM_VectorMulti, float, Eigen::Dynamic)->RangeMultiplier(4)->Range(16, 256);
BENCHMARK_TEMPLATE(BM_VectorMulti, double, Eigen::Dynamic)->RangeMultiplier(4)->Range(16,256);

// Run the benchmark
BENCHMARK_MAIN();