Add vectorized version of filter

Author: ProExpertProg

cpp/CMakeLists.txt

@@ -31,7 +31,14 @@ FetchContent_Declare(
         GIT_TAG e593f6695c6065e6b345fe2862f04a519ed484e0
 )
 
-FetchContent_MakeAvailable(argparse mdspan cnpy spdlog)
+FetchContent_Declare(
+        eve
+        GIT_REPOSITORY https://github.com/jfalcou/eve.git
+        GIT_TAG b2d8b637e71d132654c52480549e9b79944d1f74
+)
+
+option(EVE_BUILD_TEST "Build EVE tests" OFF)
+FetchContent_MakeAvailable(argparse mdspan cnpy spdlog eve)
 target_include_directories(cnpy PUBLIC ${CNPY_SOURCE_DIR})
 
 # ==========================
@@ -44,6 +51,7 @@ target_link_libraries(fftw-cpp INTERFACE mdspan)
 option(ENABLE_AVX512 "Enable AVX512 instructions (e.g. might want to turn off for the sake of Valgrind)" OFF)
 if (ENABLE_AVX512)
     add_compile_options("-mavx512f" "-march=native")
+    add_compile_definitions(AVX512_ENABLED)
 endif ()
 
 option(ENABLE_CILK "Enable Cilk parallelism" ON)
@@ -66,7 +74,8 @@ add_library(src-lib OBJECT lib/Naive.cpp lib/Naive.hpp lib/HermiteRunner.cpp lib
         lib/Brackets.hpp lib/Brackets.cpp
 )
 target_include_directories(src-lib PUBLIC lib/)
-target_link_libraries(src-lib mdspan fftw-cpp cnpy spdlog::spdlog)
+# TODO give eve only to vectorized targets
+target_link_libraries(src-lib mdspan fftw-cpp cnpy spdlog::spdlog eve::eve)
 
 add_subdirectory(bench)
 

cpp/bench/hl-filter.cpp

@@ -35,3 +35,7 @@ BENCHMARK(BM_HouLiFilter<HouLiFilterCached>)
 BENCHMARK(BM_HouLiFilter<HouLiFilterCached1D>)
     ->ArgsProduct({{2048, 4096, 8192}, {2048, 4096, 8192}})
     ->Unit(benchmark::kMillisecond);
+
+BENCHMARK(BM_HouLiFilter<HouLiFilterCached1DVector>)
+    ->ArgsProduct({{2048, 4096, 8192}, {2048, 4096, 8192}})
+    ->Unit(benchmark::kMillisecond);

cpp/lib/Filter.cpp

@@ -3,6 +3,9 @@
 //
 
 #include "Filter.hpp"
+
+#include <eve/module/core.hpp>
+
 namespace ahr {
 void HouLiFilter::operator()(Grid::View::C_XY view) const {
   grid.for_each_kxky([&](Dim kx, Dim ky) {
@@ -39,4 +42,33 @@ void HouLiFilterCached1D::operator()(Grid::View::C_XY view) const {
     view(kx, ky) *= factors_x[kx] * factors_y[ky];
   });
 }
+
+HouLiFilterCached1DVector::HouLiFilterCached1DVector(Grid const &grid)
+    : HouLiFilterCached1D(grid), factors_x_duped(2 * grid.KX) {
+  assert(grid.KY % C_WIDTH == 0);
+  for (Dim kx = 0; kx < grid.KX; ++kx) {
+    factors_x_duped[kx * 2] = factors_x[kx];
+    factors_x_duped[kx * 2 + 1] = factors_x[kx];
+  }
+}
+
+void HouLiFilterCached1DVector::operator()(Grid::View::C_XY view) const {
+  for (int ky = 0; ky < grid.KY; ++ky) {
+    int kx = 0;
+    for (; kx <= grid.KX - C_WIDTH; kx += C_WIDTH) {
+      Real *view_addr = (Real *)&view(kx, ky);
+      VReal input{view_addr};
+
+      VReal vfx{&factors_x_duped[kx * 2]};
+      VReal vfy{factors_y[ky]};
+
+      eve::store(input * vfx * vfy, view_addr);
+    }
+
+    // tail
+    for (; kx < grid.KX; ++kx) {
+      view(kx, ky) *= factors_x[kx] * factors_y[ky];
+    }
+  }
+}
 } // namespace ahr

cpp/lib/Filter.hpp

@@ -1,6 +1,7 @@
 #pragma once
 #include "constants.hpp"
 #include "grid.hpp"
+#include <eve/wide.hpp>
 
 namespace ahr {
 
@@ -33,4 +34,19 @@ class HouLiFilterCached1D : protected HouLiFilter {
 protected:
   std::vector<Real> factors_x, factors_y;
 };
+
+class HouLiFilterCached1DVector : HouLiFilterCached1D {
+public:
+  explicit HouLiFilterCached1DVector(Grid const &grid);
+  void operator()(Grid::View::C_XY view) const;
+
+private:
+  using VReal = eve::wide<Real>;
+  static auto constexpr R_WIDTH = VReal::size();
+  static auto constexpr C_WIDTH = VReal::size() / 2;
+
+  /// A pre-expanded vector of 2d factors
+  std::vector<Real> factors_x_duped;
+};
+
 } // namespace ahr

cpp/test/filter.cpp

@@ -14,7 +14,7 @@ template <typename TestedFilter> class TestFilter : public ::testing::Test {
   TestedFilter filter_t{grid};
 };
 
-using Types = ::testing::Types<HouLiFilterCached, HouLiFilterCached1D>;
+using Types = ::testing::Types<HouLiFilterCached, HouLiFilterCached1D, HouLiFilterCached1DVector>;
 TYPED_TEST_SUITE(TestFilter, Types);
 
 TYPED_TEST(TestFilter, Filter) {