From f65b97b96ec4111807222fd95057813a454a2bea Mon Sep 17 00:00:00 2001
From: Jaron Kent-Dobias <jaron@kent-dobias.com>
Date: Sat, 5 Apr 2025 20:50:12 -0300
Subject: Enable multithreaded FFTW with openmp

---
 Makefile    | 6 +++---
 fourier.cpp | 2 ++
 fourier.hpp | 4 ++++
 3 files changed, 9 insertions(+), 3 deletions(-)

diff --git a/Makefile b/Makefile
index d4c15bf..547e336 100644
--- a/Makefile
+++ b/Makefile
@@ -1,6 +1,6 @@
 all: walk correlations integrator fourier_integrator get_energy
 
-CC := clang++ -std=c++17 -flto -Wno-mathematical-notation-identifier-extension -O3 -march=native -mtune=native -fopenmp
+CC := clang++ -std=c++17 -flto -Wno-mathematical-notation-identifier-extension -O3 -march=native -mtune=native -fopenmp -DFFTW_THREADS=1
 
 walk: walk.cpp
 	$(CC) walk.cpp -o walk
@@ -15,7 +15,7 @@ fourier.o: fourier.cpp fourier.hpp
 	$(CC) fourier.cpp -c -o fourier.o
 
 fourier_integrator: fourier_integrator.cpp fourier.hpp fourier.o
-	$(CC) fourier_integrator.cpp fourier.o -lfftw3 -o fourier_integrator
+	$(CC) fourier_integrator.cpp fourier.o -lfftw3 -lfftw3_omp -o fourier_integrator
 
 get_energy: get_energy.cpp fourier.hpp fourier.o
-	$(CC) get_energy.cpp fourier.o -lfftw3 -o get_energy
+	$(CC) get_energy.cpp fourier.o -lfftw3 -lfftw3_omp -o get_energy
diff --git a/fourier.cpp b/fourier.cpp
index 0c22637..c2ed600 100644
--- a/fourier.cpp
+++ b/fourier.cpp
@@ -25,6 +25,8 @@ Real ddf(Real λ, unsigned p, unsigned s, Real q) {
 }
 
 FourierTransform::FourierTransform(unsigned n, Real Δω, Real Δτ, unsigned flags) : a(2 * n), â(n + 1), Δω(Δω), Δτ(Δτ) {
+  fftw_init_threads();
+  fftw_plan_with_nthreads(FFTW_THREADS);
   fftw_import_wisdom_from_filename("fftw.wisdom");
   plan_r2c = fftw_plan_dft_r2c_1d(2 * n, a.data(), reinterpret_cast<fftw_complex*>(â.data()), flags);
   plan_c2r = fftw_plan_dft_c2r_1d(2 * n, reinterpret_cast<fftw_complex*>(â.data()), a.data(), flags);
diff --git a/fourier.hpp b/fourier.hpp
index 406f2f2..1ebb7bc 100644
--- a/fourier.hpp
+++ b/fourier.hpp
@@ -5,6 +5,10 @@
 #include <vector>
 #include <tuple>
 
+#ifndef FFTW_THREADS
+#define FFTW_THREADS 1
+#endif
+
 using Real = double;
 using Complex = std::complex<Real>;
 using namespace std::complex_literals;
-- 
cgit v1.2.3-70-g09d2