High-performance FFT Library zldsp::fft

By zsliu98 |
Tags:

zldsp::fft is a header-only C++ Fast Fourier Transform (FFT) library built with Google Highway. zldsp::fft is licensed under Apache-2.0 license. You can obtain the source code at GitHub. zldsp::fft supports SSE2/SSE4/AVX2/NEON SIMD targets and both AoS/SoA layouts for complex numbers.

Usage

  1. C++ Standard: C++20 or higher
  2. Google Highway: You must include and link Google Highway in your project. The headers in this library expect the following includes to be resolvable:
#include <hwy/aligned_allocator.h>
#include <hwy/highway.h>

Static Dispatch

target_compile_definitions(my_static_target PRIVATE HWY_COMPILE_ONLY_STATIC)

Use the compiler’s architecture option to select an SSE, AVX2, or NEON static target.

SIMD TargetGCC/ClangMSVC
SSE2-march=x86-64no flag required
SSE4-march=x86-64-v2 -maes -mpclmulnot supported
AVX2-march=x86-64-v3 -maes -mpclmul/arch:AVX2
NEON-march=armv8-a+simd/arch:armv8.0

See static_dispatch_caller for CFFT and RFFT static dispatch examples with AoS and SoA layouts.

Caller-owned Dynamic Dispatch

For example, an x86 application can enable only SSE2 and AVX2:

target_compile_definitions(my_dynamic_target PRIVATE "HWY_DISABLED_TARGETS=~(HWY_SSE2|HWY_AVX2)")

Compile this target for the oldest supported baseline (for example, -march=x86-64 for SSE2).

See wrapper, its interface, and dynamic_dispatch_caller for CFFT and RFFT dynamic dispatch examples with AoS and SoA layouts.

Benchmark

Here are some benchmarks results. The benchmark code is available at GitHub. I have tried my best to configure other libraries correctly. All libraries are built with LLVM/Clang or Apple-Clang. Google Benchmark is used for benchmark.

Included libraries are (they are subjected to their own licenses):

The benchmarks shown below are Real FFT benchmarks from order 5 to order 25.

Apple M4 Pro 14-core (local)

NEON float32
Plot 1
NEON float64
Plot 2

Intel Core i7-8850H (local)

AVX2 float32
Plot 1
AVX2 float64
Plot 2

AMD EPYC 7763 (GitHub Actions runner images)

AVX2 float32
Plot 1
AVX2 float64
Plot 2

AMD AMD EPYC 9V74 (GitHub Actions runner images)

AVX2 float32
Plot 1
AVX2 float64
Plot 2

Ampere Altra Q80-30 (Oracle Cloud VM)

NEON float32
Plot 1
NEON float64
Plot 2

Design

Low Order

For low orders, 0 <= order <= 5, zldsp_fft uses specialized CFFT kernels. Each transform size therefore has a fixed implementation without general stage iteration or stage dispatch to keep values in SIMD registers where practical.

Medium Order

For medium orders, 6 <= order <= switch_order - 1, zldsp_fft uses out-of-place Stockham DIT design. It is scheduled as follows:

  • an even order begins with Radix-4 and following by Radix-4 stages
  • an odd order begins with Radix-8 and following by Radix-4 stages

Intermediate data uses an SIMD-friendly AoSoA layout. The first/final stage is fused with conversion from/to the requested AoS or SoA input/output.

High Order

For high orders, switch_order <= order, zldsp_fft uses a hybrid Cooley–Tukey & Stockham design. It uses Radix-4 Cooley–Tukey DIF macro stages to split the transform into micro CFFTs sized to fit approximately within L1 cache. Then, each micro CFFT is completed with the medium-order Stockham DIT design. The resulting matrix is placed in natural order using a reusable terminal buffer and SIMD register-transpose tiles.

Switch Order

The boundary between medium and high orders is derived from the detected L1 cache size. The implementation estimates the largest transform orders that fit its working sets and selects:

switch_order = maximum_L1_order + 4

Reference