7#ifndef HEFFTE_BACKEND_FFTW_H
8#define HEFFTE_BACKEND_FFTW_H
10#include "heffte_r2r_executor.h"
12#ifdef Heffte_ENABLE_FFTW
63#ifndef Heffte_ENABLE_MKL
101template<direction dir>
111 plan_fftw(
int size,
int howmanyffts,
int stride,
int dist) :
127 plan_fftw(
int size1,
int size2, std::array<int, 2>
const &embed,
int howmanyffts,
int stride,
int dist){
128 std::array<int, 2> size = {size2,
size1};
130 if (embed[0] == 0
and embed[1] == 0){
132 nullptr,
nullptr, stride, dist,
136 plan =
fftwf_plan_many_dft(2, size.data(), howmanyffts,
nullptr, embed.data(), stride, dist,
137 nullptr, embed.data(), stride, dist,
144 std::array<int, 3> size = {
size3, size2,
size1};
145 plan =
fftwf_plan_many_dft(3, size.data(), 1,
nullptr,
nullptr, 1, 1,
nullptr,
nullptr, 1, 1,
159template<direction dir>
162 plan_fftw(
int size,
int howmanyffts,
int stride,
int dist) :
169 plan_fftw(
int size1,
int size2, std::array<int, 2>
const &embed,
int howmanyffts,
int stride,
int dist){
170 std::array<int, 2> size = {size2,
size1};
172 if (embed[0] == 0
and embed[1] == 0){
173 plan =
fftw_plan_many_dft(2, size.data(), howmanyffts,
nullptr,
nullptr, stride, dist,
174 nullptr,
nullptr, stride, dist,
178 plan =
fftw_plan_many_dft(2, size.data(), howmanyffts,
nullptr, embed.data(), stride, dist,
179 nullptr, embed.data(), stride, dist,
186 std::array<int, 3> size = {
size3, size2,
size1};
187 plan =
fftw_plan_many_dft(3, size.data(), 1,
nullptr,
nullptr, 1, 1,
nullptr,
nullptr, 1, 1,
219 template<
typename index>
226 block_stride(
box.osize(0) *
box.osize(1)),
227 total_size(
box.count()),
231 template<
typename index>
234 blocks(1), block_stride(0), total_size(
box.count()), embed({0, 0})
242 howmanyffts =
box.size[2];
243 }
else if (std::min(odir1, odir2) == 1 and std::max(odir1, odir2) == 2){
244 stride = box.size[0];
246 howmanyffts = box.size[0];
250 embed = {
static_cast<int>(box.size[2]),
static_cast<int>(box.size[1] * box.size[0])};
251 howmanyffts = box.size[1];
255 template<
typename index>
257 size(
box.size[0]), size2(
box.size[1]), howmanyffts(
box.size[2]),
259 blocks(1), block_stride(0),
260 total_size(
box.count()),
265 void forward(std::complex<float>
data[], std::complex<float>*)
const override{
267 for(
int i=0;
i<blocks;
i++){
273 void backward(std::complex<float>
data[], std::complex<float>*)
const override{
274 make_plan(cbackward);
275 for(
int i=0;
i<blocks;
i++){
281 void forward(std::complex<double>
data[], std::complex<double>*)
const override{
283 for(
int i=0;
i<blocks;
i++){
289 void backward(std::complex<double>
data[], std::complex<double>*)
const override{
290 make_plan(zbackward);
291 for(
int i=0;
i<blocks;
i++){
298 void forward(
float const indata[], std::complex<float>
outdata[], std::complex<float> *workspace)
const override{
308 void forward(
double const indata[], std::complex<double>
outdata[], std::complex<double> *workspace)
const override{
319 int box_size()
const override{
return total_size; }
325 template<
typename scalar_type, direction dir>
333 plan = std::unique_ptr<plan_fftw<scalar_type, dir>>(
new plan_fftw<scalar_type, dir>(size, size2, embed, howmanyffts, stride, dist));
337 int size, size2, howmanyffts, stride, dist, blocks, block_stride, total_size;
338 std::array<int, 2> embed;
341 mutable std::unique_ptr<plan_fftw<std::complex<double>,
direction::forward>> zforward;
349template<direction dir>
360 plan_fftw(
int size,
int howmanyffts,
int stride,
int rdist,
int cdist) :
382template<direction dir>
385 plan_fftw(
int size,
int howmanyffts,
int stride,
int rdist,
int cdist) :
423 template<
typename index>
431 rblock_stride(
box.osize(0) *
box.osize(1)),
432 cblock_stride(
box.osize(0) * (
box.osize(1)/2 + 1)),
440 for(
int i=0;
i<blocks;
i++){
441 float *
rdata =
const_cast<float*
>(
indata +
i * rblock_stride);
448 make_plan(sbackward);
449 for(
int i=0;
i<blocks;
i++){
457 for(
int i=0;
i<blocks;
i++){
458 double *
rdata =
const_cast<double*
>(
indata +
i * rblock_stride);
465 make_plan(dbackward);
466 for(
int i=0;
i<blocks;
i++){
481 template<
typename scalar_type, direction dir>
483 if (!plan) plan = std::unique_ptr<plan_fftw<scalar_type, dir>>(
new plan_fftw<scalar_type, dir>(size, howmanyffts, stride, rdist, cdist));
486 int size, howmanyffts, stride, blocks;
487 int rdist, cdist, rblock_stride, cblock_stride, rsize, csize;
488 mutable std::unique_ptr<plan_fftw<float, direction::forward>> sforward;
489 mutable std::unique_ptr<plan_fftw<double, direction::forward>> dforward;
490 mutable std::unique_ptr<plan_fftw<float, direction::backward>> sbackward;
491 mutable std::unique_ptr<plan_fftw<double, direction::backward>> dbackward;
532template<
typename scalar_type,
typename preprocessor, direction dir>
541 plan_fftw_r2r(
int size1,
int size2, std::array<int, 2>
const &embed,
int howmanyffts,
int stride,
int dist){
542 std::array<int, 2> size = {size2,
size1};
545 if (embed[0] == 0
and embed[1] == 0){
555 std::array<int, 3> size = {
size3, size2,
size1};
557 plan =
fftw_r2r_types<scalar_type>::plan_many(3, size.data(), 1,
nullptr,
nullptr, 1, 1,
nullptr,
nullptr, 1, 1,
kind.data(),
FFTW_ESTIMATE);
561 std::array<fftw_r2r_kind, dims>
kind;
562 if (std::is_same<preprocessor, cpu_cos_pre_pos_processor>::value) {
564 }
else if (std::is_same<preprocessor, cpu_sin_pre_pos_processor>::value) {
566 }
else if (std::is_same<preprocessor, cpu_cos1_pre_pos_processor>::value) {
568 }
else if (std::is_same<preprocessor, cpu_sin1_pre_pos_processor>::value) {
583 template<
typename index>
590 block_stride(
box.osize(0) *
box.osize(1)),
591 total_size(
box.count()),
595 template<
typename index>
598 blocks(1), block_stride(0), total_size(
box.count()), embed({0, 0})
606 howmanyffts =
box.size[2];
607 }
else if (std::min(odir1, odir2) == 1 and std::max(odir1, odir2) == 2){
608 stride = box.size[0];
610 howmanyffts = box.size[0];
614 embed = {
static_cast<int>(box.size[2]),
static_cast<int>(box.size[1] * box.size[0])};
615 howmanyffts = box.size[1];
619 template<
typename index>
621 size(
box.size[0]), size2(
box.size[1]), howmanyffts(
box.size[2]),
623 blocks(1), block_stride(0),
624 total_size(
box.count()),
629 int box_size()
const override{
return total_size; }
637 for(
int i=0;
i<blocks;
i++){
644 for(
int i=0;
i<blocks;
i++){
650 make_plan(sbackward);
651 for(
int i=0;
i<blocks;
i++){
657 make_plan(dbackward);
658 for(
int i=0;
i<blocks;
i++){
665 template<
typename scalar_type, direction dir>
677 int size, size2, howmanyffts, stride, dist, blocks, block_stride, total_size;
678 std::array<int, 2> embed;
679 mutable std::unique_ptr<plan_fftw_r2r<float, prepost_processor, direction::forward>> sforward;
680 mutable std::unique_ptr<plan_fftw_r2r<double, prepost_processor, direction::forward>> dforward;
681 mutable std::unique_ptr<plan_fftw_r2r<float, prepost_processor, direction::backward>> sbackward;
682 mutable std::unique_ptr<plan_fftw_r2r<double, prepost_processor, direction::backward>> dbackward;
753 static const bool use_reorder =
true;
762 static const bool use_reorder =
true;
770 static const bool use_reorder =
true;
778 static const bool use_reorder =
true;
786 static const bool use_reorder =
true;
Base class for all backend executors.
Definition heffte_common.h:561
virtual int complex_size() const
Return the size of the complex-box (r2c executors).
Definition heffte_common.h:594
virtual void backward(float[], float *) const
Backward r2r, single precision.
Definition heffte_common.h:570
virtual void forward(float[], float *) const
Forward r2r, single precision.
Definition heffte_common.h:566
Wrapper to fftw3 API for real-to-complex transform with shortening of the data.
Definition heffte_backend_fftw.h:412
void forward(float const indata[], std::complex< float > outdata[], std::complex< float > *) const override
Forward transform, single precision.
Definition heffte_backend_fftw.h:438
int box_size() const override
Returns the size of the box with real data.
Definition heffte_backend_fftw.h:473
int complex_size() const override
Returns the size of the box with complex coefficients.
Definition heffte_backend_fftw.h:475
void forward(double const indata[], std::complex< double > outdata[], std::complex< double > *) const override
Forward transform, double precision.
Definition heffte_backend_fftw.h:455
fftw_executor_r2c(void *, box3d< index > const box, int dimension)
Constructor defines the box and the dimension of reduction.
Definition heffte_backend_fftw.h:424
void backward(std::complex< double > indata[], double outdata[], std::complex< double > *) const override
Backward transform, double precision.
Definition heffte_backend_fftw.h:464
void backward(std::complex< float > indata[], float outdata[], std::complex< float > *) const override
Backward transform, single precision.
Definition heffte_backend_fftw.h:447
size_t workspace_size() const override
Return the size of the needed workspace.
Definition heffte_backend_fftw.h:477
Wrapper around the FFTW3 API.
Definition heffte_backend_fftw.h:210
void forward(double const indata[], std::complex< double > outdata[], std::complex< double > *workspace) const override
Converts the real data to complex and performs double-complex forward transform.
Definition heffte_backend_fftw.h:308
void forward(std::complex< double > data[], std::complex< double > *) const override
Forward fft, double-complex case.
Definition heffte_backend_fftw.h:281
void backward(std::complex< float > indata[], float outdata[], std::complex< float > *workspace) const override
Performs backward float-complex transform and truncates the complex part of the result.
Definition heffte_backend_fftw.h:303
void forward(float const indata[], std::complex< float > outdata[], std::complex< float > *workspace) const override
Converts the real data to complex and performs float-complex forward transform.
Definition heffte_backend_fftw.h:298
size_t workspace_size() const override
Return the size of the needed workspace.
Definition heffte_backend_fftw.h:321
fftw_executor(void *, box3d< index > const box, int dimension)
Constructor, specifies the box and dimension.
Definition heffte_backend_fftw.h:220
fftw_executor(void *, box3d< index > const box, int dir1, int dir2)
Merges two FFTs into one.
Definition heffte_backend_fftw.h:232
void backward(std::complex< float > data[], std::complex< float > *) const override
Backward fft, float-complex case.
Definition heffte_backend_fftw.h:273
void forward(std::complex< float > data[], std::complex< float > *) const override
Forward fft, float-complex case.
Definition heffte_backend_fftw.h:265
fftw_executor(void *, box3d< index > const box)
Merges three FFTs into one.
Definition heffte_backend_fftw.h:256
void backward(std::complex< double > data[], std::complex< double > *) const override
Backward fft, double-complex case.
Definition heffte_backend_fftw.h:289
void backward(std::complex< double > indata[], double outdata[], std::complex< double > *workspace) const override
Performs backward double-complex transform and truncates the complex part of the result.
Definition heffte_backend_fftw.h:313
int box_size() const override
Returns the size of the box.
Definition heffte_backend_fftw.h:319
int fft1d_get_howmany(box3d< index > const box, int const dimension)
Return the number of 1-D ffts contained in the box in the given dimension.
Definition heffte_geometry.h:159
int fft1d_get_stride(box3d< index > const box, int const dimension)
Return the stride of the 1-D ffts contained in the box in the given dimension.
Definition heffte_geometry.h:169
direction
Indicates the direction of the FFT (internal use only).
Definition heffte_common.h:652
@ backward
Inverse DFT transform.
@ forward
Forward DFT transform.
Namespace containing all HeFFTe methods and classes.
Definition heffte_backend_cuda.h:38
Defines inverse mapping from the location tag to a default backend tag.
Definition heffte_common.h:430
stock type
Defaults to the stock backend.
Definition heffte_common.h:432
Type-tag for the Cosine Transform type 1 using the FFTW backend.
Definition heffte_common.h:114
Type-tag for the Cosine Transform using the FFTW backend.
Definition heffte_common.h:104
Type-tag for the Sine Transform type 1 using the FFTW backend.
Definition heffte_common.h:119
Type-tag for the Sine Transform using the FFTW backend.
Definition heffte_common.h:109
Type-tag for the FFTW backend.
Definition heffte_common.h:99
Allows to define whether a specific backend interface has been enabled.
Definition heffte_common.h:226
Defines a set of default plan options for a given backend.
Definition heffte_common.h:761
static void plan_destroy(plan_type p)
Destructor.
Definition heffte_backend_fftw.h:523
static plan_type plan_many(vars... args)
Make the plan.
Definition heffte_backend_fftw.h:521
Definition heffte_backend_fftw.h:502
fftw_plan plan_type
Plan type.
Definition heffte_backend_fftw.h:504
static void plan_destroy(plan_type p)
Destructor.
Definition heffte_backend_fftw.h:508
static plan_type plan_many(vars... args)
Make the plan.
Definition heffte_backend_fftw.h:506
Struct to specialize to allow HeFFTe to recognize custom single precision complex types.
Definition heffte_utils.h:252
Struct to specialize to allow HeFFTe to recognize custom double precision complex types.
Definition heffte_utils.h:270
Indicates the structure that will be used by the fft backend.
Definition heffte_common.h:663
Wrapper around cufftHandle plans, set for float or double complex.
Definition heffte_backend_cuda.h:346
fftw_plan plan
Identical to the float-complex specialization.
Definition heffte_backend_fftw.h:401
~plan_fftw()
Identical to the float-complex specialization.
Definition heffte_backend_fftw.h:397
plan_fftw(int size, int howmanyffts, int stride, int rdist, int cdist)
Identical to the float-complex specialization.
Definition heffte_backend_fftw.h:385
~plan_fftw()
Identical to the float-complex specialization.
Definition heffte_backend_fftw.h:372
fftwf_plan plan
Identical to the float-complex specialization.
Definition heffte_backend_fftw.h:376
plan_fftw(int size, int howmanyffts, int stride, int rdist, int cdist)
Constructor taking into account the different sizes for the real and complex parts.
Definition heffte_backend_fftw.h:360
fftw_plan plan
Identical to the float-complex specialization.
Definition heffte_backend_fftw.h:195
plan_fftw(int size1, int size2, std::array< int, 2 > const &embed, int howmanyffts, int stride, int dist)
Identical to the float-complex specialization.
Definition heffte_backend_fftw.h:169
~plan_fftw()
Identical to the float-complex specialization.
Definition heffte_backend_fftw.h:191
plan_fftw(int size, int howmanyffts, int stride, int dist)
Identical to the float-complex specialization.
Definition heffte_backend_fftw.h:162
plan_fftw(int size1, int size2, int size3)
Identical to the float-complex specialization.
Definition heffte_backend_fftw.h:185
plan_fftw(int size1, int size2, int size3)
Identical to the float-complex specialization.
Definition heffte_backend_fftw.h:143
plan_fftw(int size, int howmanyffts, int stride, int dist)
Constructor, takes inputs identical to fftwf_plan_many_dft().
Definition heffte_backend_fftw.h:111
~plan_fftw()
Destructor, deletes the plan.
Definition heffte_backend_fftw.h:149
plan_fftw(int size1, int size2, std::array< int, 2 > const &embed, int howmanyffts, int stride, int dist)
Constructor, takes inputs identical to fftwf_plan_many_dft().
Definition heffte_backend_fftw.h:127
fftwf_plan plan
The FFTW3 opaque structure (pointer to struct).
Definition heffte_backend_fftw.h:153
Wrapper for the r2r plan to allow for RAII style of management.
Definition heffte_backend_fftw.h:533
~plan_fftw_r2r()
Identical to the other specialization.
Definition heffte_backend_fftw.h:575
plan_fftw_r2r(int size1, int size2, std::array< int, 2 > const &embed, int howmanyffts, int stride, int dist)
Identical to the float-complex specialization.
Definition heffte_backend_fftw.h:541
plan_fftw_r2r(int size, int howmanyffts, int stride, int dist)
Identical to the float-complex specialization.
Definition heffte_backend_fftw.h:535
plan_fftw_r2r(int size1, int size2, int size3)
Identical to the float-complex specialization.
Definition heffte_backend_fftw.h:554
static std::array< fftw_r2r_kind, dims > make_kind_array()
Make the array with the kind of the transform.
Definition heffte_backend_fftw.h:560
fftw_r2r_types< scalar_type >::plan_type plan
The actual fftw plan.
Definition heffte_backend_fftw.h:579
Base plan for fftw, using only the specialization for float and double complex.
Definition heffte_backend_fftw.h:93
virtual void backward(double data[], double *) const override
Backward r2r, double precision.
Definition heffte_backend_fftw.h:656
real2real_executor(void *, box3d< index > const box)
Construct a plan for a single 3D transform, not implemented currently.
Definition heffte_backend_fftw.h:620
virtual void forward(double data[], double *) const override
Forward r2r, double precision.
Definition heffte_backend_fftw.h:642
real2real_executor(void *, box3d< index > const box, int dimension)
Construct a plan for batch 1D transforms.
Definition heffte_backend_fftw.h:584
virtual void backward(float data[], float *) const override
Backward r2r, single precision.
Definition heffte_backend_fftw.h:649
int box_size() const override
Returns the size of the box.
Definition heffte_backend_fftw.h:629
virtual void forward(float data[], float *) const override
Forward r2r, single precision.
Definition heffte_backend_fftw.h:635
size_t workspace_size() const override
Returns the size of the box.
Definition heffte_backend_fftw.h:631
real2real_executor(void *, box3d< index > const box, int dir1, int dir2)
Construct a plan for batch 2D transforms, not implemented currently.
Definition heffte_backend_fftw.h:596
Template algorithm for the Sine and Cosine transforms.
Definition heffte_r2r_executor.h:192