Fix Emscripten build with SIMD

.github/workflows/c-cpp.yml

@@ -15,6 +15,19 @@ env:
   BUILD_TYPE: Release
 
 jobs:
+  build_emscripten:
+    runs-on: ubuntu-latest
+
+    steps:
+    - uses: actions/checkout@v4
+    - uses: emscripten-core/setup-emsdk@v16
+    - name: cmake configure
+      run: emcmake cmake -S . -B build_wasm -DCMAKE_BUILD_TYPE=$BUILD_TYPE -DCMAKE_CROSSCOMPILING_EMULATOR=node -DPFFFT_USE_BENCH_GREEN=OFF -DPFFFT_USE_BENCH_KISS=OFF -DPFFFT_USE_BENCH_POCKET=OFF -DPFFFT_BUILD_TESTS=ON -DPFFFT_BUILD_BENCHMARKS=OFF -DPFFFT_BUILD_EXAMPLES=OFF
+    - name: build
+      run: cmake --build build_wasm
+    - name: test
+      run: cd build_wasm && ctest --output-on-failure
+
   build_w_mipp_ubuntu-amd64:
     runs-on: ubuntu-latest
 

README.md

@@ -136,7 +136,7 @@ sudo apt-get install cmake-curses-gui
 Some of the options:
 * `PFFFT_USE_TYPE_FLOAT` to activate single precision 'float' (default: ON)
 * `PFFFT_USE_TYPE_DOUBLE` to activate 'double' precision float (default: ON)
-* `PFFFT_USE_SIMD` to use SIMD (SSE/AVX/NEON/ALTIVEC) CPU features? (default: ON)
+* `PFFFT_USE_SIMD` to use SIMD (SSE/AVX/NEON/ALTIVEC/WASM SIMD) CPU features? (default: ON)
 * `DISABLE_SIMD_AVX` to disable AVX CPU features (default: OFF)
 * `PFFFT_USE_SIMD_NEON` to force using NEON on ARM (requires PFFFT_USE_SIMD) (default: OFF)
 * `PFFFT_USE_SCALAR_VECT` to use 4-element vector scalar operations (if no other SIMD) (default: ON)
@@ -169,6 +169,17 @@ ctest -C Release
 
 see [https://cmake.org/cmake/help/v3.15/manual/cmake-generators.7.html#visual-studio-generators](https://cmake.org/cmake/help/v3.15/manual/cmake-generators.7.html#visual-studio-generators)
 
+With [Emscripten](https://emscripten.org/) for WebAssembly (requires the [Emscripten SDK](https://emscripten.org/docs/getting_started/downloads.html)):
+
+```sh
+mkdir build
+cd build
+emcmake cmake ..
+cmake --build .
+ctest
+```
+
+WASM SIMD is enabled automatically. Emscripten provides NEON-to-WASM SIMD translation via [SIMDe](https://github.com/simd-everywhere/simde) (SIMD Everywhere) compatibility headers, so pffft's NEON code paths are reused for WebAssembly.
 
 ## Using pffft in your CMake project
 

benchmarks/bench_conv.cpp

@@ -17,7 +17,7 @@
 
 #define MIPP_VECTOR  mipp::vector
 #else
-#define MIPP_VECTOR  std::vector
+#define MIPP_VECTOR  ::std::vector
 #endif
 
 #include "pf_conv_dispatcher.h"
@@ -30,7 +30,7 @@
 MIPP_VECTOR<float> generate_rng_vec(int M, int N = -1, int seed_value = 1)
 {
     MIPP_VECTOR<float> v(N < 0 ? M : N);
-    std::mt19937 g;
+    ::std::mt19937 g;
     g.seed(seed_value);
     constexpr float scale = 1.0F / (1.0F + float(INT_FAST32_MAX));
     for (int k = 0; k < M; ++k)
@@ -107,7 +107,7 @@ int bench_oop(
     {
         move_rest(buffer, &state);
         //memcpy(buffer+state.size, &s[off], B * sizeof(s[0]));
-        std::copy(&signal[off], &signal[off+blockLen], buffer+state.size);
+        ::std::copy(&signal[off], &signal[off+blockLen], buffer+state.size);
         state.size += blockLen;
         int n_out = conv_oop(buffer, &state, filter, sz_filter, &y[n_out_sum]);
         n_out_sum += n_out;
@@ -139,7 +139,7 @@ int bench_cx_real_oop(
     {
         move_rest(buffer, &state);
         //memcpy(buffer+state.size, &s[off], B * sizeof(s[0]));
-        std::copy(&signal[off], &signal[off+blockLen], &buffer[state.size]);
+        ::std::copy(&signal[off], &signal[off+blockLen], &buffer[state.size]);
         state.size += blockLen;
         int n_out = conv_oop(buffer, &state, filter, sz_filter, &y[n_out_sum]);
         n_out_sum += n_out;
@@ -234,7 +234,7 @@ int main(int argc, char *argv[])
         //else
         //    float_simd_size[a] = 0;
     }
-    //const int max_simd_size = *std::max_element( &float_simd_size[0], &float_simd_size[num_arch] );
+    //const int max_simd_size = *::std::max_element( &float_simd_size[0], &float_simd_size[num_arch] );
     if (verbose)
         fprintf(stderr, "max float simd size: %d\n", max_simd_size);
 
@@ -315,12 +315,12 @@ int main(int argc, char *argv[])
         fprintf(stderr, "y[%2d] = %g  %+g * i\n", k, y[2*k], y[2*k+1]);
     fprintf(stderr, "\n");
 
-    const std::complex<float> * sc = reinterpret_cast< std::complex<float>* >( s.data() );
+    const ::std::complex<float> * sc = reinterpret_cast< ::std::complex<float>* >( s.data() );
     const int Nc = N /2;
     fprintf(stderr, "reference with std::complex<float>:\n");
     for (int off = 0; off +filterLen <= Nc; ++off )
     {
-        std::complex<float> sum(0.0F, 0.0F);
+        ::std::complex<float> sum(0.0F, 0.0F);
         for (int k=0; k < filterLen; ++k)
             sum += sc[off+k] * filter[k];
         fprintf(stderr, "yv[%2d] = %g  %+g * i\n", off, sum.real(), sum.imag() );

benchmarks/pf_conv.cpp

@@ -56,7 +56,7 @@ void ARCHFUNCNAME(conv_float_move_rest)(float * RESTRICT s, conv_buffer_state *
     if (R > 0)
     {
         // memmove(s, &s[state->offset], R * sizeof(s[0]));   // move them to the begin
-        std::copy(&s[state->offset], &s[state->size], s);
+        ::std::copy(&s[state->offset], &s[state->size], s);
     }
     else
         R = 0;
@@ -71,7 +71,7 @@ void ARCHFUNCNAME(conv_cplx_move_rest)(complexf * RESTRICT s, conv_buffer_state
     if (R > 0)
     {
         // memmove(s, &s[state->offset], R * sizeof(s[0]));   // move them to the begin
-        std::copy(&s[state->offset], &s[state->size], s);
+        ::std::copy(&s[state->offset], &s[state->size], s);
     }
     else
         R = 0;

cmake/target_optimizations.cmake

@@ -62,6 +62,11 @@ elseif (CMAKE_SYSTEM_PROCESSOR MATCHES "armv7l")
     set(GCC_MARCH_DESC "native/ARMwNEON:armv7-a")
     set(GCC_MARCH_VALUES "none;native;armv7-a" CACHE INTERNAL "List of possible architectures")
     set(GCC_EXTRA_VALUES "none;neon_vfpv4;neon_rpi3_a53;neon_rpi4_a72" CACHE INTERNAL "List of possible additional options")
+elseif (EMSCRIPTEN)
+    # Emscripten WASM SIMD is handled automatically in target_set_c/cxx_arch_flags
+    set(GCC_MARCH_DESC "wasm-simd")
+    set(GCC_MARCH_VALUES "none" CACHE INTERNAL "List of possible architectures")
+    set(GCC_EXTRA_VALUES "" CACHE INTERNAL "List of possible additional options")
 else()
     message(WARNING "unsupported CMAKE_SYSTEM_PROCESSOR '${CMAKE_SYSTEM_PROCESSOR}'")
     # other PROCESSORs could be "ppc", "ppc64",  "arm" - or something else?!
@@ -117,6 +122,13 @@ endif()
 ######################################################
 
 function(target_set_c_arch_flags target)
+    # Emscripten WASM SIMD via NEON emulation
+    if (EMSCRIPTEN)
+        message(STATUS "Emscripten detected: enabling WASM SIMD with NEON emulation for C target ${target}")
+        target_compile_options(${target} PRIVATE "-msimd128")
+        target_compile_definitions(${target} PRIVATE PFFFT_ENABLE_NEON=1)
+        return()
+    endif()
     if ( ("${TARGET_C_ARCH}" STREQUAL "") OR ("${TARGET_C_ARCH}" STREQUAL "none") )
         message(STATUS "C ARCH for target ${target} is not set!")
     else()
@@ -150,6 +162,13 @@ function(target_set_c_arch_flags target)
 endfunction()
 
 function(target_set_cxx_arch_flags target)
+    # Emscripten WASM SIMD via NEON emulation
+    if (EMSCRIPTEN)
+        message(STATUS "Emscripten detected: enabling WASM SIMD with NEON emulation for C++ target ${target}")
+        target_compile_options(${target} PRIVATE "-msimd128")
+        target_compile_definitions(${target} PRIVATE PFFFT_ENABLE_NEON=1)
+        return()
+    endif()
     if ( ("${TARGET_CXX_ARCH}" STREQUAL "") OR ("${TARGET_CXX_ARCH}" STREQUAL "none") )
         message(STATUS "C++ ARCH for target ${target} is not set!")
     else()

examples/example_cpp11_cplx_dbl_fwd.cpp

@@ -7,25 +7,25 @@
 
 void cxx11_forward_complex_double(const int transformLen)
 {
-  std::cout << "running " << __FUNCTION__ << "()" << std::endl;
+  ::std::cout << "running " << __FUNCTION__ << "()" << ::std::endl;
 
   // first check - might be skipped
-  using FFT_T = pffft::Fft< std::complex<double> >;
+  using FFT_T = pffft::Fft< ::std::complex<double> >;
   if (transformLen < FFT_T::minFFtsize())
   {
-    std::cerr << "Error: minimum FFT transformation length is " << FFT_T::minFFtsize() << std::endl;
+    ::std::cerr << "Error: minimum FFT transformation length is " << FFT_T::minFFtsize() << ::std::endl;
     return;
   }
 
   // instantiate FFT and prepare transformation for length N
-  pffft::Fft< std::complex<double> > fft(transformLen);
+  pffft::Fft< ::std::complex<double> > fft(transformLen);
 
   // one more check
   if (!fft.isValid())
   {
-    std::cerr << "Error: transformation length " << transformLen << " is not decomposable into small prime factors. "
+    ::std::cerr << "Error: transformation length " << transformLen << " is not decomposable into small prime factors. "
               << "Next valid transform size is: " << FFT_T::nearestTransformSize(transformLen)
-              << "; next power of 2 is: " << FFT_T::nextPowerOfTwo(transformLen) << std::endl;
+              << "; next power of 2 is: " << FFT_T::nextPowerOfTwo(transformLen) << ::std::endl;
     return;
   }
 
@@ -34,26 +34,26 @@ void cxx11_forward_complex_double(const int transformLen)
   auto Y = fft.spectrumVector();
 
   // alternative access: get raw pointers to aligned vectors
-  std::complex<double> *Xs = X.data();
-  std::complex<double> *Ys = Y.data();
+  ::std::complex<double> *Xs = X.data();
+  ::std::complex<double> *Ys = Y.data();
 
   // prepare some input data
   for (int k = 0; k < transformLen; k += 2)
   {
-    X[k] = std::complex<double>(k, k&1);        // access through AlignedVector<double>
-    Xs[k+1] = std::complex<double>(-1-k, k&1);  // access through raw pointer
+    X[k] = ::std::complex<double>(k, k&1);        // access through AlignedVector<double>
+    Xs[k+1] = ::std::complex<double>(-1-k, k&1);  // access through raw pointer
   }
 
   // do the forward transform; write complex spectrum result into Y
   fft.forward(X, Y);
 
   // print spectral output
-  std::cout << "output should be complex spectrum with " << fft.getSpectrumSize() << " bins" << std::endl;
-  std::cout << "output vector has size " << Y.size() << " (complex bins):" << std::endl;
+  ::std::cout << "output should be complex spectrum with " << fft.getSpectrumSize() << " bins" << ::std::endl;
+  ::std::cout << "output vector has size " << Y.size() << " (complex bins):" << ::std::endl;
   for (unsigned k = 0; k < Y.size(); k += 2)
   {
-    std::cout << "Y[" << k << "] = " << Y[k] << std::endl;
-    std::cout << "Y[" << k+1 << "] = " << Ys[k+1] << std::endl;
+    ::std::cout << "Y[" << k << "] = " << Y[k] << ::std::endl;
+    ::std::cout << "Y[" << k+1 << "] = " << Ys[k+1] << ::std::endl;
   }
 }
 

examples/example_cpp11_real_dbl_fwd.cpp

@@ -7,13 +7,13 @@
 
 void cxx11_forward_real_double(const int transformLen)
 {
-  std::cout << "running " << __FUNCTION__ << "()" << std::endl;
+  ::std::cout << "running " << __FUNCTION__ << "()" << ::std::endl;
 
   // first check - might be skipped
   using FFT_T = pffft::Fft<double>;
   if (transformLen < FFT_T::minFFtsize())
   {
-    std::cerr << "Error: minimum FFT transformation length is " << FFT_T::minFFtsize() << std::endl;
+    ::std::cerr << "Error: minimum FFT transformation length is " << FFT_T::minFFtsize() << ::std::endl;
     return;
   }
 
@@ -23,9 +23,9 @@ void cxx11_forward_real_double(const int transformLen)
   // one more check
   if (!fft.isValid())
   {
-    std::cerr << "Error: transformation length " << transformLen << " is not decomposable into small prime factors. "
+    ::std::cerr << "Error: transformation length " << transformLen << " is not decomposable into small prime factors. "
               << "Next valid transform size is: " << FFT_T::nearestTransformSize(transformLen)
-              << "; next power of 2 is: " << FFT_T::nextPowerOfTwo(transformLen) << std::endl;
+              << "; next power of 2 is: " << FFT_T::nextPowerOfTwo(transformLen) << ::std::endl;
     return;
   }
 
@@ -35,7 +35,7 @@ void cxx11_forward_real_double(const int transformLen)
 
   // alternative access: get raw pointers to aligned vectors
   double *Xs = X.data();
-  std::complex<double> *Ys = Y.data();
+  ::std::complex<double> *Ys = Y.data();
 
   // prepare some input data
   for (int k = 0; k < transformLen; k += 2)
@@ -48,12 +48,12 @@ void cxx11_forward_real_double(const int transformLen)
   fft.forward(X, Y);
 
   // print spectral output
-  std::cout << "output should be complex spectrum with " << fft.getSpectrumSize() << " bins" << std::endl;
-  std::cout << "output vector has size " << Y.size() << " (complex bins):" << std::endl;
+  ::std::cout << "output should be complex spectrum with " << fft.getSpectrumSize() << " bins" << ::std::endl;
+  ::std::cout << "output vector has size " << Y.size() << " (complex bins):" << ::std::endl;
   for (unsigned k = 0; k < Y.size(); k += 2)
   {
-    std::cout << "Y[" << k << "] = " << Y[k] << std::endl;
-    std::cout << "Y[" << k+1 << "] = " << Ys[k+1] << std::endl;
+    ::std::cout << "Y[" << k << "] = " << Y[k] << ::std::endl;
+    ::std::cout << "Y[" << k+1 << "] = " << Ys[k+1] << ::std::endl;
   }
 }
 

examples/example_cpp98_cplx_flt_fwd.cpp

@@ -7,53 +7,53 @@
 
 void cxx98_forward_complex_float(const int transformLen)
 {
-  std::cout << "running " << __FUNCTION__ << "()" << std::endl;
+  ::std::cout << "running " << __FUNCTION__ << "()" << ::std::endl;
 
   // first check - might be skipped
-  typedef pffft::Fft< std::complex<float> > FFT_T;
+  typedef pffft::Fft< ::std::complex<float> > FFT_T;
   if (transformLen < FFT_T::minFFtsize())
   {
-    std::cerr << "Error: minimum FFT transformation length is " << FFT_T::minFFtsize() << std::endl;
+    ::std::cerr << "Error: minimum FFT transformation length is " << FFT_T::minFFtsize() << ::std::endl;
     return;
   }
 
   // instantiate FFT and prepare transformation for length N
-  pffft::Fft< std::complex<float> > fft(transformLen);
+  pffft::Fft< ::std::complex<float> > fft(transformLen);
 
   // one more check
   if (!fft.isValid())
   {
-    std::cerr << "Error: transformation length " << transformLen << " is not decomposable into small prime factors. "
+    ::std::cerr << "Error: transformation length " << transformLen << " is not decomposable into small prime factors. "
               << "Next valid transform size is: " << FFT_T::nearestTransformSize(transformLen)
-              << "; next power of 2 is: " << FFT_T::nextPowerOfTwo(transformLen) << std::endl;
+              << "; next power of 2 is: " << FFT_T::nextPowerOfTwo(transformLen) << ::std::endl;
     return;
   }
 
   // allocate aligned vectors for input X and output Y
-  pffft::AlignedVector< std::complex<float> > X = fft.valueVector();
-  pffft::AlignedVector< std::complex<float> > Y = fft.spectrumVector();
+  pffft::AlignedVector< ::std::complex<float> > X = fft.valueVector();
+  pffft::AlignedVector< ::std::complex<float> > Y = fft.spectrumVector();
 
   // alternative access: get raw pointers to aligned vectors
-  std::complex<float> *Xs = X.data();
-  std::complex<float> *Ys = Y.data();
+  ::std::complex<float> *Xs = X.data();
+  ::std::complex<float> *Ys = Y.data();
 
   // prepare some input data
   for (int k = 0; k < transformLen; k += 2)
   {
-    X[k] = std::complex<float>(k, k&1);        // access through AlignedVector<float>
-    Xs[k+1] = std::complex<float>(-1-k, k&1);  // access through raw pointer
+    X[k] = ::std::complex<float>(k, k&1);        // access through AlignedVector<float>
+    Xs[k+1] = ::std::complex<float>(-1-k, k&1);  // access through raw pointer
   }
 
   // do the forward transform; write complex spectrum result into Y
   fft.forward(X, Y);
 
   // print spectral output
-  std::cout << "output should be complex spectrum with " << fft.getSpectrumSize() << " bins" << std::endl;
-  std::cout << "output vector has size " << Y.size() << " (complex bins):" << std::endl;
+  ::std::cout << "output should be complex spectrum with " << fft.getSpectrumSize() << " bins" << ::std::endl;
+  ::std::cout << "output vector has size " << Y.size() << " (complex bins):" << ::std::endl;
   for (unsigned k = 0; k < Y.size(); k += 2)
   {
-    std::cout << "Y[" << k << "] = " << Y[k] << std::endl;
-    std::cout << "Y[" << k+1 << "] = " << Ys[k+1] << std::endl;
+    ::std::cout << "Y[" << k << "] = " << Y[k] << ::std::endl;
+    ::std::cout << "Y[" << k+1 << "] = " << Ys[k+1] << ::std::endl;
   }
 }