test avx512

libmandel/CMakeLists.txt

@@ -39,7 +39,7 @@ SET(MandelSources
 FILE(GLOB MandelHeaders include/*.h)
 
 if (ARCH STREQUAL "X86_64" OR ARCH STREQUAL "X86")
-    list(APPEND MandelSources src/CpuGeneratorsAVX.cpp src/CpuGeneratorsAVXFMA.cpp src/CpuGeneratorsSSE2.cpp)
+    list(APPEND MandelSources src/CpuGeneratorsAVX512.cpp  src/CpuGeneratorsAVX.cpp src/CpuGeneratorsAVXFMA.cpp src/CpuGeneratorsSSE2.cpp)
 elseif(ARCH STREQUAL "ARM")
     list(APPEND MandelSources src/CpuGeneratorsNeon.cpp)
 endif()
@@ -90,6 +90,12 @@ endif(Boost_FOUND)
 
 if (ARCH STREQUAL "X86_64" OR ARCH STREQUAL "X86")
     if (MSVC)
+        set_source_files_properties(src/CpuGeneratorsAVX512.cpp PROPERTIES COMPILE_FLAGS /arch:AVX512F)
+    else()
+        set_source_files_properties(src/CpuGeneratorsAVX512.cpp PROPERTIES COMPILE_FLAGS -mavx512f)
+    endif(MSVC)
+
+    if (MSVC)
         set_source_files_properties(src/CpuGeneratorsAVX.cpp PROPERTIES COMPILE_FLAGS /arch:AVX)
         set_source_files_properties(src/CpuGeneratorsAVXFMA.cpp PROPERTIES COMPILE_FLAGS /arch:AVX)
     else()

libmandel/include/CpuGenerators.h

@@ -11,6 +11,7 @@ namespace mnd
         X86_SSE2,
         X86_AVX,
         X86_AVX_FMA,
+        X86_AVX_512,
         ARM_NEON,
     };
 
@@ -126,11 +127,11 @@ public:
 
 
 template<bool parallel>
-class mnd::CpuGenerator<double, mnd::X86_AVX_FMA, parallel> : public Generator
+class mnd::CpuGenerator<float, mnd::X86_AVX_FMA, parallel> : public Generator
 {
 public:
     inline CpuGenerator(void) :
-        Generator{ mnd::getPrecision<double>() }
+        Generator{ mnd::getPrecision<float>() }
     {
     }
     virtual void generate(const MandelInfo& info, float* data);
@@ -138,11 +139,11 @@ public:
 
 
 template<bool parallel>
-class mnd::CpuGenerator<float, mnd::X86_AVX_FMA, parallel> : public Generator
+class mnd::CpuGenerator<double, mnd::X86_AVX_FMA, parallel> : public Generator
 {
 public:
     inline CpuGenerator(void) :
-        Generator{ mnd::getPrecision<float>() }
+        Generator{ mnd::getPrecision<double>() }
     {
     }
     virtual void generate(const MandelInfo& info, float* data);
@@ -161,4 +162,27 @@ public:
 };
 
 
+template<bool parallel>
+class mnd::CpuGenerator<float, mnd::X86_AVX_512, parallel> : public Generator
+{
+public:
+    inline CpuGenerator(void) :
+        Generator{ mnd::getPrecision<float>() }
+    {
+    }
+    virtual void generate(const MandelInfo& info, float* data);
+};
+
+
+template<bool parallel>
+class mnd::CpuGenerator<double, mnd::X86_AVX_512, parallel> : public Generator
+{
+public:
+    inline CpuGenerator(void) :
+        Generator{ mnd::getPrecision<double>() }
+    {
+    }
+    virtual void generate(const MandelInfo& info, float* data);
+};
+
 #endif // MANDEL_CPUGENERATORS_H

libmandel/src/CpuGeneratorsAVX512.cpp

@@ -5,44 +5,49 @@
 
 #include <memory>
 
-using mnd::CpuGeneratorAvx512Float;
-using mnd::CpuGeneratorAvx512Double;
+using mnd::CpuGenerator;
 
-void CpuGeneratorAvx512Float::generate(const mnd::MandelInfo& info, float* data)
+namespace mnd
+{
+    template class CpuGenerator<float, mnd::X86_AVX_512, false>;
+    template class CpuGenerator<float, mnd::X86_AVX_512, true>;
+
+    template class CpuGenerator<double, mnd::X86_AVX_512, false>;
+    template class CpuGenerator<double, mnd::X86_AVX_512, true>;
+}
+
+template<bool parallel>
+void CpuGenerator<float, mnd::X86_AVX_512, parallel>::generate(const mnd::MandelInfo& info, float* data)
 {
     using T = float;
     const MandelViewport& view = info.view;
-    omp_set_num_threads(2 * omp_get_num_procs());
-#pragma omp parallel for
+
+    const float dppf = float(view.width / info.bWidth);
+    const float viewxf = float(view.x);
+    __m512 viewx = { viewxf, viewxf, viewxf, viewxf, viewxf, viewxf, viewxf, viewxf,
+                     viewxf, viewxf, viewxf, viewxf, viewxf, viewxf, viewxf, viewxf };
+    __m512 dpp = { dppf, dppf, dppf, dppf, dppf, dppf, dppf, dppf,
+                   dppf, dppf, dppf, dppf, dppf, dppf, dppf, dppf };
+
+    if constexpr(parallel)
+        omp_set_num_threads(omp_get_num_procs());
+#pragma omp parallel for schedule(static, 1) if (parallel)
     for (long j = 0; j < info.bHeight; j++) {
-        T y = T(view.y) + T(j) * T(view.height / info.bHeight);
+        T y = T(view.y + double(j) * view.height / info.bHeight);
+        __m512 ys = { y, y, y, y, y, y, y, y, y, y, y, y, y, y, y, y };
         long i = 0;
         for (i; i < info.bWidth; i += 16) {
-            __m512 xs = {
-                float(view.x + double(i) * view.width / info.bWidth),
-                float(view.x + double(i + 1) * view.width / info.bWidth),
-                float(view.x + double(i + 2) * view.width / info.bWidth),
-                float(view.x + double(i + 3) * view.width / info.bWidth),
-                float(view.x + double(i + 4) * view.width / info.bWidth),
-                float(view.x + double(i + 5) * view.width / info.bWidth),
-                float(view.x + double(i + 6) * view.width / info.bWidth),
-                float(view.x + double(i + 7) * view.width / info.bWidth),
-                float(view.x + double(i + 8) * view.width / info.bWidth),
-                float(view.x + double(i + 9) * view.width / info.bWidth),
-                float(view.x + double(i + 10) * view.width / info.bWidth),
-                float(view.x + double(i + 11) * view.width / info.bWidth),
-                float(view.x + double(i + 12) * view.width / info.bWidth),
-                float(view.x + double(i + 13) * view.width / info.bWidth),
-                float(view.x + double(i + 14) * view.width / info.bWidth),
-                float(view.x + double(i + 15) * view.width / info.bWidth)
-            };
+            __m512 pixc = { float(i), float(i + 1), float(i + 2), float(i + 3), float(i + 4), float(i + 5), float(i + 6), float(i + 7),
+                            float(i + 8), float(i + 9), float(i + 10), float(i + 11), float(i + 12), float(i + 13), float(i + 14), float(i + 15) };
+            __m512 xs = _mm512_fmadd_ps(dpp, pixc, viewx);
 
-            __m512 counter = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
-            __m512 adder = {1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1};
+            __m512 counter = { 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 };
+            __m512 adder = { 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1 };
 
-            __m512 threshold = {16.0f, 16.0f, 16.0f, 16.0f, 16.0f, 16.0f, 16.0f, 16.0f, 16.0f, 16.0f, 16.0f, 16.0f, 16.0f, 16.0f, 16.0f, 16.0f};
+            __m512 threshold = { 16.0f, 16.0f, 16.0f, 16.0f, 16.0f, 16.0f, 16.0f, 16.0f,
+                                 16.0f, 16.0f, 16.0f, 16.0f, 16.0f, 16.0f, 16.0f, 16.0f };
 
-            __m512 ys = {y, y, y, y, y, y, y, y, y, y, y, y, y, y, y};
+            __m512 ys = { y, y, y, y, y, y, y, y, y, y, y, y, y, y, y };
             __m512 a = xs;
             __m512 b = ys;
 
@@ -77,6 +82,67 @@ void CpuGeneratorAvx512Float::generate(const mnd::MandelInfo& info, float* data)
 }
 
 
-void CpuGeneratorAvx512Double::generate(const mnd::MandelInfo& info, float* data)
+template<bool parallel>
+void CpuGenerator<double, mnd::X86_AVX_512, parallel>::generate(const mnd::MandelInfo& info, float* data)
 {
+    using T = float;
+    const MandelViewport& view = info.view;
+
+    const float dppf = float(view.width / info.bWidth);
+    const float viewxf = float(view.x);
+    __m512 viewx = { viewxf, viewxf, viewxf, viewxf, viewxf, viewxf, viewxf, viewxf,
+                     viewxf, viewxf, viewxf, viewxf, viewxf, viewxf, viewxf, viewxf };
+    __m512 dpp = { dppf, dppf, dppf, dppf, dppf, dppf, dppf, dppf,
+                   dppf, dppf, dppf, dppf, dppf, dppf, dppf, dppf };
+
+    if constexpr(parallel)
+        omp_set_num_threads(omp_get_num_procs());
+#pragma omp parallel for schedule(static, 1) if (parallel)
+    for (long j = 0; j < info.bHeight; j++) {
+        T y = T(view.y + double(j) * view.height / info.bHeight);
+        __m512 ys = { y, y, y, y, y, y, y, y, y, y, y, y, y, y, y, y };
+        long i = 0;
+        for (i; i < info.bWidth; i += 16) {
+            __m512 pixc = { float(i), float(i + 1), float(i + 2), float(i + 3), float(i + 4), float(i + 5), float(i + 6), float(i + 7),
+                            float(i + 8), float(i + 9), float(i + 10), float(i + 11), float(i + 12), float(i + 13), float(i + 14), float(i + 15) };
+            __m512 xs = _mm512_fmadd_ps(dpp, pixc, viewx);
+
+            __m512 counter = { 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 };
+            __m512 adder = { 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1 };
+
+            __m512 threshold = { 16.0f, 16.0f, 16.0f, 16.0f, 16.0f, 16.0f, 16.0f, 16.0f,
+                                 16.0f, 16.0f, 16.0f, 16.0f, 16.0f, 16.0f, 16.0f, 16.0f };
+
+            __m512 ys = { y, y, y, y, y, y, y, y, y, y, y, y, y, y, y };
+            __m512 a = xs;
+            __m512 b = ys;
+
+            for (int k = 0; k < info.maxIter; k++) {
+                __m512 aa = _mm512_mul_ps(a, a);
+                __m512 bb = _mm512_mul_ps(b, b);
+                __m512 abab = _mm512_mul_ps(a, b); abab = _mm512_add_ps(abab, abab);
+                a = _mm512_add_ps(_mm512_sub_ps(aa, bb), xs);
+                b = _mm512_add_ps(abab, ys);
+                __mmask16 cmp = _mm512_cmp_ps_mask(_mm512_add_ps(aa, bb), threshold, _CMP_LE_OQ);
+                counter = _mm512_mask_add_ps(counter, cmp, counter, adder);
+                if (cmp == 0) {
+                    break;
+                }
+            }
+
+            auto alignVec = [](float* data) -> float* {
+                void* aligned = data;
+                ::size_t length = 64;
+                std::align(32, 16 * sizeof(float), aligned, length);
+                return static_cast<float*>(aligned);
+            };
+
+            float resData[32];
+            float* ftRes = alignVec(resData);
+
+            _mm512_store_ps(ftRes, counter);
+            for (int k = 0; k < 8 && i + k < info.bWidth; k++)
+                data[i + k + j * info.bWidth] = ftRes[k] > 0 ? ftRes[k] : info.maxIter;
+        }
+    }
 }