|
Highly Efficient FFT for Exascale: HeFFTe v2.3
|
Classes | |
| struct | heffte_plan_options |
| Equivalent to heffte::plan_options but defined for the C API. More... | |
| struct | heffte_fft_plan |
| Generic wrapper around some fft3d/fft3d_r2c instance, use heffte_plan instead of this. More... | |
Macros | |
| #define | Heffte_BACKEND_STOCK 0 |
| Set the use of the stock backend. | |
| #define | Heffte_BACKEND_FFTW 1 |
| Set the use of the FFTW backend. | |
| #define | Heffte_BACKEND_MKL 2 |
| Set the use of the MKL backend. | |
| #define | Heffte_BACKEND_CUFFT 10 |
| Set the use of the cuFFT backend. | |
| #define | Heffte_BACKEND_ROCFFT 11 |
| Set the use of the rocfft backend. | |
| #define | Heffte_SUCCESS 0 |
| The success return of a HeFFTe method. | |
| #define | Heffte_RESHAPE_ALGORITHM_ALLTOALLV 0 |
| Corresponds to heffte::reshape_algorithm::alltoallv. | |
| #define | Heffte_RESHAPE_ALGORITHM_P2P_PLINED 1 |
| Corresponds to heffte::reshape_algorithm::p2p_plined. | |
| #define | Heffte_RESHAPE_ALGORITHM_P2P 2 |
| Corresponds to heffte::reshape_algorithm::p2p. | |
| #define | Heffte_SCALE_NONE 0 |
| Indicate no scaling, see heffte::scale::none. | |
| #define | Heffte_SCALE_FULL 1 |
| Indicate full scaling, see heffte::scale::full. | |
| #define | Heffte_SCALE_SYMMETRIC 2 |
| Indicate symmetric scaling, see heffte::scale::symmetric. | |
Typedefs | |
| typedef heffte_fft_plan * | heffte_plan |
| C-style wrapper around an instance of heffte::fft3d or heffte::fft3d_r2c using some backend. | |
Functions | |
| int | heffte_set_default_options (int backend, heffte_plan_options *options) |
| Populates the fields of the options with the ones default for the backend. More... | |
| int | heffte_plan_create (int backend, int const inbox_low[3], int const inbox_high[3], int const *inbox_order, int const outbox_low[3], int const outbox_high[3], int const *outbox_order, MPI_Comm const comm, heffte_plan_options const *options, heffte_plan *plan) |
| Creates a new heffte_plan with the given backend, boxes, communicator, and options. More... | |
| int | heffte_plan_create_r2c (int backend, int const inbox_low[3], int const inbox_high[3], int const *inbox_order, int const outbox_low[3], int const outbox_high[3], int const *outbox_order, int r2c_direction, MPI_Comm const comm, heffte_plan_options const *options, heffte_plan *plan) |
| Creates a new heffte_plan for an r2c operation with the given backend, boxes, communicator, and options. More... | |
| int | heffte_plan_destroy (heffte_plan plan) |
| Destory a heffte_plan, call the destructor of the internal object. More... | |
| int | heffte_size_inbox (heffte_plan const plan) |
| Wrapper around heffte::fft3d::size_inbox() and heffte::fft3d_r2c::size_inbox() | |
| int | heffte_size_outbox (heffte_plan const plan) |
| Wrapper around heffte::fft3d::size_outbox() and heffte::fft3d_r2c::size_outbox() | |
| int | heffte_size_workspace (heffte_plan const plan) |
| Wrapper around heffte::fft3d::size_workspace() and heffte::fft3d_r2c::size_workspace() | |
| int | heffte_get_backend (heffte_plan const plan) |
| Returns the backend used in the create command. | |
| int | heffte_is_r2c (heffte_plan const plan) |
| Returns 1 if heffte_plan_create_r2c() was called and 0 otherwise. | |
| void | heffte_forward_s2c (heffte_plan const plan, float const *input, void *output, int scale) |
| Wrapper around heffte::fft3d::forward() and heffte::fft3d_r2c::forward() More... | |
| void | heffte_forward_c2c (heffte_plan const plan, void const *input, void *output, int scale) |
| Forward transform, single precision, complex to complex, see heffte_forward_s2c() | |
| void | heffte_forward_d2z (heffte_plan const plan, double const *input, void *output, int scale) |
| Forward transform, double precision, real to complex, see heffte_forward_s2c() | |
| void | heffte_forward_z2z (heffte_plan const plan, void const *input, void *output, int scale) |
| Forward transform, double precision, complex to complex, see heffte_forward_s2c() | |
| void | heffte_forward_s2c_buffered (heffte_plan const plan, float const *input, void *output, void *workspace, int scale) |
| Wrapper around heffte::fft3d::forward() and heffte::fft3d_r2c::forward() with user allocated workspace. More... | |
| void | heffte_forward_c2c_buffered (heffte_plan const plan, void const *input, void *output, void *workspace, int scale) |
| Forward transform, buffered, double precision, complex to complex, see heffte_forward_s2c_buffered() | |
| void | heffte_forward_d2z_buffered (heffte_plan const plan, double const *input, void *output, void *workspace, int scale) |
| Forward transform, buffered, double precision, complex to complex, see heffte_forward_s2c_buffered() | |
| void | heffte_forward_z2z_buffered (heffte_plan const plan, void const *input, void *output, void *workspace, int scale) |
| Forward transform, buffered, double precision, complex to complex, see heffte_forward_s2c_buffered() | |
| void | heffte_backward_c2s (heffte_plan const plan, void const *input, float *output, int scale) |
| Wrapper around heffte::fft3d::backward() and heffte::fft3d_r2c::backward() More... | |
| void | heffte_backward_c2c (heffte_plan const plan, void const *input, void *output, int scale) |
| Backward transform, single precision, complex to complex, see heffte_backward_c2s() | |
| void | heffte_backward_z2d (heffte_plan const plan, void const *input, double *output, int scale) |
| Backward transform, double precision, complex to real, see heffte_backward_c2s() | |
| void | heffte_backward_z2z (heffte_plan const plan, void const *input, void *output, int scale) |
| Backward transform, double precision, complex to complex, see heffte_backward_c2s() | |
| void | heffte_backward_c2s_buffered (heffte_plan const plan, void const *input, float *output, void *workspace, int scale) |
| Wrapper around heffte::fft3d::backward() and heffte::fft3d_r2c::backward() with user allocated workspace. More... | |
| void | heffte_backward_c2c_buffered (heffte_plan const plan, void const *input, void *output, void *workspace, int scale) |
| Backward transform, buffered, single precision, complex to complex, see heffte_backward_c2s_buffered() | |
| void | heffte_backward_z2d_buffered (heffte_plan const plan, void const *input, double *output, void *workspace, int scale) |
| Backward transform, buffered, double precision, complex to real, see heffte_backward_c2s_buffered() | |
| void | heffte_backward_z2z_buffered (heffte_plan const plan, void const *input, void *output, void *workspace, int scale) |
| Backward transform, buffered, double precision, complex to complex, see heffte_backward_c2s_buffered() | |
The templated classes heffte::fft3d and heffte::fft3d_r2c are replaced by opaque plan handles of type heffte_plan. The plan has to be created and later destroyed with
Queries of the problem size, backend and r2c reduction of indexes can be performed with:
Forward and backward transforms are performed using variants of:
The buffered versions of the methods accept an extra parameter that is a user allocated workspace.
Complex arrays are accepted as void-pointers which removes type safety but works around the issues with multiple complex types. The methods use different suffixes to indicate the input/output types using standard BLAS naming conventions:
| int heffte_set_default_options | ( | int | backend, |
| heffte_plan_options * | options | ||
| ) |
Populates the fields of the options with the ones default for the backend.
| int heffte_plan_create | ( | int | backend, |
| int const | inbox_low[3], | ||
| int const | inbox_high[3], | ||
| int const * | inbox_order, | ||
| int const | outbox_low[3], | ||
| int const | outbox_high[3], | ||
| int const * | outbox_order, | ||
| MPI_Comm const | comm, | ||
| heffte_plan_options const * | options, | ||
| heffte_plan * | plan | ||
| ) |
Creates a new heffte_plan with the given backend, boxes, communicator, and options.
Creates a new heffte_plan with the given set of parameters.
| backend | must be Heffte_BACKEND_FFTW, Heffte_BACKEND_MKL or Heffte_BACKEND_CUFFT and the corresponding backend must be enabled in CMake |
| inbox_low | is the low indexes of the inbox in heffte::fft3d() |
| inbox_high | is the high indexes of the inbox in heffte::fft3d() |
| inbox_order | is the order indexes of the inbox in heffte::fft3d(), if set to NULL the default order (0, 1, 2) will be used |
| outbox_low | is the low indexes of the outbox in heffte::fft3d() |
| outbox_high | is the high indexes of the outbox in heffte::fft3d() |
| outbox_order | is the order indexes of the outbox in heffte::fft3d(), if set to NULL the default order (0, 1, 2) will be used |
| comm | is the MPI communicator holding all the boxes |
| options | is either NULL, which will force the use of the default backend options, of an instance of heffte_plan_options that will translate to the options used by the heffte::fft3d() |
| plan | will be overwritten with a new heffte_plan must be destroyed by heffte_plan_destroy() |
| int heffte_plan_create_r2c | ( | int | backend, |
| int const | inbox_low[3], | ||
| int const | inbox_high[3], | ||
| int const * | inbox_order, | ||
| int const | outbox_low[3], | ||
| int const | outbox_high[3], | ||
| int const * | outbox_order, | ||
| int | r2c_direction, | ||
| MPI_Comm const | comm, | ||
| heffte_plan_options const * | options, | ||
| heffte_plan * | plan | ||
| ) |
Creates a new heffte_plan for an r2c operation with the given backend, boxes, communicator, and options.
Identical to heffte_plan_create() with the addition of the r2c_direction and the internal class will be heffte::fft3d_r2c
| int heffte_plan_destroy | ( | heffte_plan | plan | ) |
Destory a heffte_plan, call the destructor of the internal object.
| plan | to be destoryed |
| void heffte_forward_s2c | ( | heffte_plan const | plan, |
| float const * | input, | ||
| void * | output, | ||
| int | scale | ||
| ) |
Wrapper around heffte::fft3d::forward() and heffte::fft3d_r2c::forward()
Performs forward Fourier transform from real input data to complex output data, both in single precision.
| plan | is an already created heffte_plan |
| input | has size at least heffte_size_inbox() and contains the input data |
| output | has size at least heffte_size_outbox() of complex-single precision numbers |
| scale | is Heffte_SCALE_NONE, Heffte_SCALE_FULL or Heffte_SCALE_SYMMETRIC, see heffte::scale for details |
The _s2c suffix indicates the real/complex types and precision using standard BLAS conventions:
| Input | Output | Method |
| float | std::complex<float> | heffte_forward_s2c() |
| std::complex<float> | std::complex<float> | heffte_forward_c2c() |
| double | std::complex<double> | heffte_forward_d2z() |
| std::complex<double> | std::complex<double> | heffte_forward_z2z() |
| void heffte_forward_s2c_buffered | ( | heffte_plan const | plan, |
| float const * | input, | ||
| void * | output, | ||
| void * | workspace, | ||
| int | scale | ||
| ) |
Wrapper around heffte::fft3d::forward() and heffte::fft3d_r2c::forward() with user allocated workspace.
Accepts a user allocated workspace buffer, but otherwise identical to heffte_forward_s2c()
| plan | is an already created heffte_plan |
| input | has size at least heffte_size_inbox() and contains the input data |
| output | has size at least heffte_size_outbox() of complex-single precision numbers |
| workspace | has size at least heffte_size_workspace() of complex-single precision numbers |
| scale | is Heffte_SCALE_NONE, Heffte_SCALE_FULL or Heffte_SCALE_SYMMETRIC, see heffte::scale for details |
The method has four variants corresponding to different precision and real/complex types:
| Input | Output | Workspace | Method |
| float | std::complex<float> | std::complex<float> | heffte_forward_s2c_buffered() |
| std::complex<float> | std::complex<float> | std::complex<float> | heffte_forward_c2c_buffered() |
| double | std::complex<double> | std::complex<double> | heffte_forward_d2z_buffered() |
| std::complex<double> | std::complex<double> | std::complex<double> | heffte_forward_z2z_buffered() |
| void heffte_backward_c2s | ( | heffte_plan const | plan, |
| void const * | input, | ||
| float * | output, | ||
| int | scale | ||
| ) |
Wrapper around heffte::fft3d::backward() and heffte::fft3d_r2c::backward()
Performs backward Fourier transform from complex input data to real output data, both in single precision.
| plan | is an already created heffte_plan |
| input | has size at least heffte_size_inbox() of complex-single precision numbers |
| output | has size at least heffte_size_outbox() will be overwritten with the output |
| scale | is Heffte_SCALE_NONE, Heffte_SCALE_FULL or Heffte_SCALE_SYMMETRIC, see heffte::scale for details |
The _s2c suffix indicates the real/complex types and precision using standard BLAS conventions:
| Input | Output | Method |
| std::complex<float> | float | heffte_backward_c2s() |
| std::complex<float> | std::complex<float> | heffte_backward_c2c() |
| std::complex<double> | double | heffte_backward_z2d() |
| std::complex<double> | std::complex<double> | heffte_backward_z2z() |
| void heffte_backward_c2s_buffered | ( | heffte_plan const | plan, |
| void const * | input, | ||
| float * | output, | ||
| void * | workspace, | ||
| int | scale | ||
| ) |
Wrapper around heffte::fft3d::backward() and heffte::fft3d_r2c::backward() with user allocated workspace.
Accepts a user allocated workspace buffer, but otherwise identical to heffte_backward_s2c()
| plan | is an already created heffte_plan |
| input | has size at least heffte_size_inbox() of complex-single precision numbers |
| output | has size at least heffte_size_outbox() will be overwritten with the output |
| workspace | has size at least heffte_size_workspace() of complex-single precision numbers |
| scale | is Heffte_SCALE_NONE, Heffte_SCALE_FULL or Heffte_SCALE_SYMMETRIC, see heffte::scale for details |
The method has four variants corresponding to different precision and real/complex types:
| Input | Output | Workspace | Method |
| std::complex<float> | float | std::complex<float> | heffte_backward_c2s_buffered() |
| std::complex<float> | std::complex<float> | std::complex<float> | heffte_backward_c2c_buffered() |
| std::complex<double> | double | std::complex<double> | heffte_backward_z2d_buffered() |
| std::complex<double> | std::complex<double> | std::complex<double> | heffte_backward_z2z_buffered() |