WarpX
AnyFFT.H
Go to the documentation of this file.
1 /* Copyright 2019-2023
2  *
3  * This file is part of ABLASTR.
4  *
5  * License: BSD-3-Clause-LBNL
6  */
7 
8 #ifndef ABLASTR_ANYFFT_H_
9 #define ABLASTR_ANYFFT_H_
10 
11 #ifdef ABLASTR_USE_FFT
12 # include <AMReX_Config.H>
13 # include <AMReX_LayoutData.H>
14 
15 # if defined(AMREX_USE_CUDA)
16 # include <cufft.h>
17 # elif defined(AMREX_USE_HIP)
18 # if __has_include(<rocfft/rocfft.h>) // ROCm 5.3+
19 # include <rocfft/rocfft.h>
20 # else
21 # include <rocfft.h>
22 # endif
23 # else
24 # include <fftw3.h>
25 # endif
26 #endif
27 
28 
36 {
37 
41  void setup();
42 
46  void cleanup();
47 
48 #ifdef ABLASTR_USE_FFT
49 
50  // First, define library-dependent types (complex, FFT plan)
51 
53 # if defined(AMREX_USE_CUDA)
54 # ifdef AMREX_USE_FLOAT
55  using Complex = cuComplex;
56 # else
57  using Complex = cuDoubleComplex;
58 # endif
59 # elif defined(AMREX_USE_HIP)
60 # ifdef AMREX_USE_FLOAT
61  using Complex = float2;
62 # else
63  using Complex = double2;
64 # endif
65 # else
66 # ifdef AMREX_USE_FLOAT
67  using Complex = fftwf_complex;
68 # else
69  using Complex = fftw_complex;
70 # endif
71 # endif
72 
76 # if defined(AMREX_USE_CUDA)
77  using VendorFFTPlan = cufftHandle;
78 # elif defined(AMREX_USE_HIP)
79  using VendorFFTPlan = rocfft_plan;
80 # else
81 # ifdef AMREX_USE_FLOAT
82  using VendorFFTPlan = fftwf_plan;
83 # else
84  using VendorFFTPlan = fftw_plan;
85 # endif
86 # endif
87 
88  // Second, define library-independent API
89 
91  enum struct direction {R2C, C2R};
92 
95  struct FFTplan
96  {
97  amrex::Real* m_real_array;
98  Complex* m_complex_array;
99  VendorFFTPlan m_plan;
100  direction m_dir;
101  int m_dim;
102  };
103 
105  using FFTplans = amrex::LayoutData<FFTplan>;
106 
115  FFTplan CreatePlan(const amrex::IntVect& real_size, amrex::Real* real_array,
116  Complex* complex_array, direction dir, int dim);
117 
121  void DestroyPlan(FFTplan& fft_plan);
122 
126  void Execute(FFTplan& fft_plan);
127 
128 #endif
129 
130 }
131 
132 #endif // ABLASTR_ANYFFT_H_
int m_dir
Definition: AnyFFT.H:36
FFTplan CreatePlan(const amrex::IntVect &real_size, amrex::Real *const real_array, Complex *const complex_array, const direction dir, const int dim)
Definition: WrapCuFFT.cpp:30
void DestroyPlan(FFTplan &fft_plan)
Definition: WrapCuFFT.cpp:72
void cleanup()
Definition: WrapCuFFT.cpp:18
void Execute(FFTplan &fft_plan)
Definition: WrapCuFFT.cpp:78
void setup()
Definition: WrapCuFFT.cpp:16