Highly Efficient FFT for Exascale: HeFFTe v2.4
Loading...
Searching...
No Matches
heffte_backend_fftw.h
1/*
2 -- heFFTe --
3 Univ. of Tennessee, Knoxville
4 @date
5*/
6
7#ifndef HEFFTE_BACKEND_FFTW_H
8#define HEFFTE_BACKEND_FFTW_H
9
10#include "heffte_r2r_executor.h"
11
12#ifdef Heffte_ENABLE_FFTW
13
14#include "fftw3.h"
15
27namespace heffte{
28
29namespace backend{
34 template<> struct is_enabled<fftw> : std::true_type{};
35
40 template<> struct is_enabled<fftw_cos> : std::true_type{};
45 template<> struct is_enabled<fftw_sin> : std::true_type{};
50 template<> struct is_enabled<fftw_cos1> : std::true_type{};
55 template<> struct is_enabled<fftw_sin1> : std::true_type{};
56
57// Specialization is not necessary since the default behavior assumes CPU parameters.
58// template<>
59// struct buffer_traits<fftw>{
60// using location = tag::cpu;
61// template<typename T> using container = std::vector<T>;
62// };
63#ifndef Heffte_ENABLE_MKL
68 template<> struct default_backend<tag::cpu> {
69 using type = fftw;
70 };
71#endif
72}
73
78template<> struct is_ccomplex<fftwf_complex> : std::true_type{};
83template<> struct is_zcomplex<fftw_complex> : std::true_type{};
84
93template<typename, direction> struct plan_fftw{};
94
101template<direction dir>
102struct plan_fftw<std::complex<float>, dir>{
111 plan_fftw(int size, int howmanyffts, int stride, int dist) :
112 plan(fftwf_plan_many_dft(1, &size, howmanyffts, nullptr, nullptr, stride, dist,
113 nullptr, nullptr, stride, dist,
114 (dir == direction::forward) ? FFTW_FORWARD : FFTW_BACKWARD, FFTW_ESTIMATE
115 ))
116 {}
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};
129
130 if (embed[0] == 0 and embed[1] == 0){
131 plan = fftwf_plan_many_dft(2, size.data(), howmanyffts, nullptr, nullptr, stride, dist,
132 nullptr, nullptr, stride, dist,
133 (dir == direction::forward) ? FFTW_FORWARD : FFTW_BACKWARD, FFTW_ESTIMATE
134 );
135 }else{
136 plan = fftwf_plan_many_dft(2, size.data(), howmanyffts, nullptr, embed.data(), stride, dist,
137 nullptr, embed.data(), stride, dist,
138 (dir == direction::forward) ? FFTW_FORWARD : FFTW_BACKWARD, FFTW_ESTIMATE
139 );
140 }
141 }
143 plan_fftw(int size1, int size2, int size3){
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,
146 (dir == direction::forward) ? FFTW_FORWARD : FFTW_BACKWARD, FFTW_ESTIMATE);
147 }
149 ~plan_fftw(){ fftwf_destroy_plan(plan); }
151 operator fftwf_plan() const{ return plan; }
153 fftwf_plan plan;
154};
159template<direction dir>
160struct plan_fftw<std::complex<double>, dir>{
162 plan_fftw(int size, int howmanyffts, int stride, int dist) :
163 plan(fftw_plan_many_dft(1, &size, howmanyffts, nullptr, nullptr, stride, dist,
164 nullptr, nullptr, stride, dist,
165 (dir == direction::forward) ? FFTW_FORWARD : FFTW_BACKWARD, FFTW_ESTIMATE
166 ))
167 {}
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};
171
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,
175 (dir == direction::forward) ? FFTW_FORWARD : FFTW_BACKWARD, FFTW_ESTIMATE
176 );
177 }else{
178 plan = fftw_plan_many_dft(2, size.data(), howmanyffts, nullptr, embed.data(), stride, dist,
179 nullptr, embed.data(), stride, dist,
180 (dir == direction::forward) ? FFTW_FORWARD : FFTW_BACKWARD, FFTW_ESTIMATE
181 );
182 }
183 }
185 plan_fftw(int size1, int size2, int size3){
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,
188 (dir == direction::forward) ? FFTW_FORWARD : FFTW_BACKWARD, FFTW_ESTIMATE);
189 }
191 ~plan_fftw(){ fftw_destroy_plan(plan); }
193 operator fftw_plan() const{ return plan; }
195 fftw_plan plan;
196};
197
211public:
219 template<typename index>
220 fftw_executor(void*, box3d<index> const box, int dimension) :
221 size(box.size[dimension]), size2(0),
222 howmanyffts(fft1d_get_howmany(box, dimension)),
223 stride(fft1d_get_stride(box, dimension)),
224 dist((dimension == box.order[0]) ? size : 1),
225 blocks((dimension == box.order[1]) ? box.osize(2) : 1),
226 block_stride(box.osize(0) * box.osize(1)),
227 total_size(box.count()),
228 embed({0, 0})
229 {}
231 template<typename index>
232 fftw_executor(void*, box3d<index> const box, int dir1, int dir2) :
233 size(box.size[std::min(dir1, dir2)]), size2(box.size[std::max(dir1, dir2)]),
234 blocks(1), block_stride(0), total_size(box.count()), embed({0, 0})
235 {
236 int odir1 = box.find_order(dir1);
237 int odir2 = box.find_order(dir2);
238
239 if (std::min(odir1, odir2) == 0 and std::max(odir1, odir2) == 1){
240 stride = 1;
241 dist = size * size2;
242 howmanyffts = box.size[2];
243 }else if (std::min(odir1, odir2) == 1 and std::max(odir1, odir2) == 2){
244 stride = box.size[0];
245 dist = 1;
246 howmanyffts = box.size[0];
247 }else{ // case of directions (0, 2)
248 stride = 1;
249 dist = size;
250 embed = {static_cast<int>(box.size[2]), static_cast<int>(box.size[1] * box.size[0])};
251 howmanyffts = box.size[1];
252 }
253 }
255 template<typename index>
256 fftw_executor(void*, box3d<index> const box) :
257 size(box.size[0]), size2(box.size[1]), howmanyffts(box.size[2]),
258 stride(0), dist(0),
259 blocks(1), block_stride(0),
260 total_size(box.count()),
261 embed({0, 0})
262 {}
263
265 void forward(std::complex<float> data[], std::complex<float>*) const override{
266 make_plan(cforward);
267 for(int i=0; i<blocks; i++){
268 fftwf_complex* block_data = reinterpret_cast<fftwf_complex*>(data + i * block_stride);
269 fftwf_execute_dft(*cforward, block_data, block_data);
270 }
271 }
273 void backward(std::complex<float> data[], std::complex<float>*) const override{
274 make_plan(cbackward);
275 for(int i=0; i<blocks; i++){
276 fftwf_complex* block_data = reinterpret_cast<fftwf_complex*>(data + i * block_stride);
277 fftwf_execute_dft(*cbackward, block_data, block_data);
278 }
279 }
281 void forward(std::complex<double> data[], std::complex<double>*) const override{
282 make_plan(zforward);
283 for(int i=0; i<blocks; i++){
284 fftw_complex* block_data = reinterpret_cast<fftw_complex*>(data + i * block_stride);
285 fftw_execute_dft(*zforward, block_data, block_data);
286 }
287 }
289 void backward(std::complex<double> data[], std::complex<double>*) const override{
290 make_plan(zbackward);
291 for(int i=0; i<blocks; i++){
292 fftw_complex* block_data = reinterpret_cast<fftw_complex*>(data + i * block_stride);
293 fftw_execute_dft(*zbackward, block_data, block_data);
294 }
295 }
296
298 void forward(float const indata[], std::complex<float> outdata[], std::complex<float> *workspace) const override{
299 for(int i=0; i<total_size; i++) outdata[i] = std::complex<float>(indata[i]);
300 forward(outdata, workspace);
301 }
303 void backward(std::complex<float> indata[], float outdata[], std::complex<float> *workspace) const override{
304 backward(indata, workspace);
305 for(int i=0; i<total_size; i++) outdata[i] = std::real(indata[i]);
306 }
308 void forward(double const indata[], std::complex<double> outdata[], std::complex<double> *workspace) const override{
309 for(int i=0; i<total_size; i++) outdata[i] = std::complex<double>(indata[i]);
310 forward(outdata, workspace);
311 }
313 void backward(std::complex<double> indata[], double outdata[], std::complex<double> *workspace) const override{
314 backward(indata, workspace);
315 for(int i=0; i<total_size; i++) outdata[i] = std::real(indata[i]);
316 }
317
319 int box_size() const override{ return total_size; }
321 size_t workspace_size() const override{ return 0; }
322
323private:
325 template<typename scalar_type, direction dir>
326 void make_plan(std::unique_ptr<plan_fftw<scalar_type, dir>> &plan) const{
327 if (not plan){
328 if (dist == 0)
329 plan = std::unique_ptr<plan_fftw<scalar_type, dir>>(new plan_fftw<scalar_type, dir>(size, size2, howmanyffts));
330 else if (size2 == 0)
331 plan = std::unique_ptr<plan_fftw<scalar_type, dir>>(new plan_fftw<scalar_type, dir>(size, howmanyffts, stride, dist));
332 else
333 plan = std::unique_ptr<plan_fftw<scalar_type, dir>>(new plan_fftw<scalar_type, dir>(size, size2, embed, howmanyffts, stride, dist));
334 }
335 }
336
337 int size, size2, howmanyffts, stride, dist, blocks, block_stride, total_size;
338 std::array<int, 2> embed;
339 mutable std::unique_ptr<plan_fftw<std::complex<float>, direction::forward>> cforward;
340 mutable std::unique_ptr<plan_fftw<std::complex<float>, direction::backward>> cbackward;
341 mutable std::unique_ptr<plan_fftw<std::complex<double>, direction::forward>> zforward;
342 mutable std::unique_ptr<plan_fftw<std::complex<double>, direction::backward>> zbackward;
343};
344
349template<direction dir>
350struct plan_fftw<float, dir>{
360 plan_fftw(int size, int howmanyffts, int stride, int rdist, int cdist) :
361 plan((dir == direction::forward) ?
362 fftwf_plan_many_dft_r2c(1, &size, howmanyffts, nullptr, nullptr, stride, rdist,
363 nullptr, nullptr, stride, cdist,
364 FFTW_ESTIMATE
365 ) :
366 fftwf_plan_many_dft_c2r(1, &size, howmanyffts, nullptr, nullptr, stride, cdist,
367 nullptr, nullptr, stride, rdist,
368 FFTW_ESTIMATE
369 ))
370 {}
372 ~plan_fftw(){ fftwf_destroy_plan(plan); }
374 operator fftwf_plan() const{ return plan; }
376 fftwf_plan plan;
377};
382template<direction dir>
383struct plan_fftw<double, dir>{
385 plan_fftw(int size, int howmanyffts, int stride, int rdist, int cdist) :
386 plan((dir == direction::forward) ?
387 fftw_plan_many_dft_r2c(1, &size, howmanyffts, nullptr, nullptr, stride, rdist,
388 nullptr, nullptr, stride, cdist,
389 FFTW_ESTIMATE
390 ) :
391 fftw_plan_many_dft_c2r(1, &size, howmanyffts, nullptr, nullptr, stride, cdist,
392 nullptr, nullptr, stride, rdist,
393 FFTW_ESTIMATE
394 ))
395 {}
397 ~plan_fftw(){ fftw_destroy_plan(plan); }
399 operator fftw_plan() const{ return plan; }
401 fftw_plan plan;
402};
403
413public:
423 template<typename index>
424 fftw_executor_r2c(void*, box3d<index> const box, int dimension) :
425 size(box.size[dimension]),
426 howmanyffts(fft1d_get_howmany(box, dimension)),
427 stride(fft1d_get_stride(box, dimension)),
428 blocks((dimension == box.order[1]) ? box.osize(2) : 1),
429 rdist((dimension == box.order[0]) ? size : 1),
430 cdist((dimension == box.order[0]) ? size/2 + 1 : 1),
431 rblock_stride(box.osize(0) * box.osize(1)),
432 cblock_stride(box.osize(0) * (box.osize(1)/2 + 1)),
433 rsize(box.count()),
434 csize(box.r2c(dimension).count())
435 {}
436
438 void forward(float const indata[], std::complex<float> outdata[], std::complex<float>*) const override{
439 make_plan(sforward);
440 for(int i=0; i<blocks; i++){
441 float *rdata = const_cast<float*>(indata + i * rblock_stride);
442 fftwf_complex* cdata = reinterpret_cast<fftwf_complex*>(outdata + i * cblock_stride);
443 fftwf_execute_dft_r2c(*sforward, rdata, cdata);
444 }
445 }
447 void backward(std::complex<float> indata[], float outdata[], std::complex<float>*) const override{
448 make_plan(sbackward);
449 for(int i=0; i<blocks; i++){
450 fftwf_complex* cdata = const_cast<fftwf_complex*>(reinterpret_cast<fftwf_complex const*>(indata + i * cblock_stride));
451 fftwf_execute_dft_c2r(*sbackward, cdata, outdata + i * rblock_stride);
452 }
453 }
455 void forward(double const indata[], std::complex<double> outdata[], std::complex<double>*) const override{
456 make_plan(dforward);
457 for(int i=0; i<blocks; i++){
458 double *rdata = const_cast<double*>(indata + i * rblock_stride);
459 fftw_complex* cdata = reinterpret_cast<fftw_complex*>(outdata + i * cblock_stride);
460 fftw_execute_dft_r2c(*dforward, rdata, cdata);
461 }
462 }
464 void backward(std::complex<double> indata[], double outdata[], std::complex<double>*) const override{
465 make_plan(dbackward);
466 for(int i=0; i<blocks; i++){
467 fftw_complex* cdata = const_cast<fftw_complex*>(reinterpret_cast<fftw_complex const*>(indata + i * cblock_stride));
468 fftw_execute_dft_c2r(*dbackward, cdata, outdata + i * rblock_stride);
469 }
470 }
471
473 int box_size() const override{ return rsize; }
475 int complex_size() const override{ return csize; }
477 size_t workspace_size() const override{ return 0; }
478
479private:
481 template<typename scalar_type, direction dir>
482 void make_plan(std::unique_ptr<plan_fftw<scalar_type, dir>> &plan) const{
483 if (!plan) plan = std::unique_ptr<plan_fftw<scalar_type, dir>>(new plan_fftw<scalar_type, dir>(size, howmanyffts, stride, rdist, cdist));
484 }
485
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;
492};
493
501template<typename>
504 using plan_type = fftw_plan;
506 template<typename... vars> static plan_type plan_many(vars... args){ return fftw_plan_many_r2r(args...); }
508 static void plan_destroy(plan_type p){ fftw_destroy_plan(p); }
509};
517template<> struct fftw_r2r_types<float> {
519 using plan_type = fftwf_plan;
521 template<typename... vars> static plan_type plan_many(vars... args){ return fftwf_plan_many_r2r(args...); }
523 static void plan_destroy(plan_type p){ fftwf_destroy_plan(p); }
524};
532template<typename scalar_type, typename preprocessor, direction dir>
535 plan_fftw_r2r(int size, int howmanyffts, int stride, int dist){
536 auto kind = make_kind_array<1>();
537 plan = fftw_r2r_types<scalar_type>::plan_many(1, &size, howmanyffts, nullptr, nullptr, stride, dist,
538 nullptr, nullptr, stride, dist, kind.data(), FFTW_ESTIMATE);
539 }
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};
543 auto kind = make_kind_array<2>();
544
545 if (embed[0] == 0 and embed[1] == 0){
546 plan = fftw_r2r_types<scalar_type>::plan_many(2, size.data(), howmanyffts, nullptr, nullptr, stride, dist,
547 nullptr, nullptr, stride, dist, kind.data(), FFTW_ESTIMATE);
548 }else{
549 plan = fftw_r2r_types<scalar_type>::plan_many(2, size.data(), howmanyffts, nullptr, embed.data(), stride, dist,
550 nullptr, embed.data(), stride, dist, kind.data(), FFTW_ESTIMATE);
551 }
552 }
554 plan_fftw_r2r(int size1, int size2, int size3){
555 std::array<int, 3> size = {size3, size2, size1};
556 auto kind = make_kind_array<3>();
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);
558 }
560 template<size_t dims> static std::array<fftw_r2r_kind, dims> make_kind_array() {
561 std::array<fftw_r2r_kind, dims> kind;
562 if (std::is_same<preprocessor, cpu_cos_pre_pos_processor>::value) {
563 kind[0] = (dir == direction::forward) ? FFTW_REDFT10 : FFTW_REDFT01;
564 } else if (std::is_same<preprocessor, cpu_sin_pre_pos_processor>::value) { // sin transform
565 kind[0] = (dir == direction::forward) ? FFTW_RODFT10 : FFTW_RODFT01;
566 } else if (std::is_same<preprocessor, cpu_cos1_pre_pos_processor>::value) {
567 kind[0] = FFTW_REDFT00;
568 } else if (std::is_same<preprocessor, cpu_sin1_pre_pos_processor>::value) { // sin transform
569 kind[0] = FFTW_RODFT00;
570 }
571 for(size_t i=1; i<kind.size(); i++) kind[i] = kind[0];
572 return kind;
573 }
577 operator typename fftw_r2r_types<scalar_type>::plan_type() const{ return plan; }
580};
581template<typename prepost_processor> struct real2real_executor<backend::fftw, prepost_processor> : public executor_base{
583 template<typename index>
584 real2real_executor(void*, box3d<index> const box, int dimension) :
585 size(box.size[dimension]), size2(0),
586 howmanyffts(fft1d_get_howmany(box, dimension)),
587 stride(fft1d_get_stride(box, dimension)),
588 dist((dimension == box.order[0]) ? size : 1),
589 blocks((dimension == box.order[1]) ? box.osize(2) : 1),
590 block_stride(box.osize(0) * box.osize(1)),
591 total_size(box.count()),
592 embed({0, 0})
593 {}
595 template<typename index>
596 real2real_executor(void*, box3d<index> const box, int dir1, int dir2) :
597 size(box.size[std::min(dir1, dir2)]), size2(box.size[std::max(dir1, dir2)]),
598 blocks(1), block_stride(0), total_size(box.count()), embed({0, 0})
599 {
600 int odir1 = box.find_order(dir1);
601 int odir2 = box.find_order(dir2);
602
603 if (std::min(odir1, odir2) == 0 and std::max(odir1, odir2) == 1){
604 stride = 1;
605 dist = size * size2;
606 howmanyffts = box.size[2];
607 }else if (std::min(odir1, odir2) == 1 and std::max(odir1, odir2) == 2){
608 stride = box.size[0];
609 dist = 1;
610 howmanyffts = box.size[0];
611 }else{ // case of directions (0, 2)
612 stride = 1;
613 dist = size;
614 embed = {static_cast<int>(box.size[2]), static_cast<int>(box.size[1] * box.size[0])};
615 howmanyffts = box.size[1];
616 }
617 }
619 template<typename index>
621 size(box.size[0]), size2(box.size[1]), howmanyffts(box.size[2]),
622 stride(0), dist(0),
623 blocks(1), block_stride(0),
624 total_size(box.count()),
625 embed({0, 0})
626 {}
627
629 int box_size() const override{ return total_size; }
631 size_t workspace_size() const override{
632 return 0;
633 }
635 virtual void forward(float data[], float*) const override{
636 make_plan(sforward);
637 for(int i=0; i<blocks; i++){
638 fftwf_execute_r2r(*sforward, data + i * block_stride, data + i * block_stride);
639 }
640 }
642 virtual void forward(double data[], double*) const override{
643 make_plan(dforward);
644 for(int i=0; i<blocks; i++){
645 fftw_execute_r2r(*dforward, data + i * block_stride, data + i * block_stride);
646 }
647 }
649 virtual void backward(float data[], float*) const override{
650 make_plan(sbackward);
651 for(int i=0; i<blocks; i++){
652 fftwf_execute_r2r(*sbackward, data + i * block_stride, data + i * block_stride);
653 }
654 }
656 virtual void backward(double data[], double*) const override{
657 make_plan(dbackward);
658 for(int i=0; i<blocks; i++){
659 fftw_execute_r2r(*dbackward, data + i * block_stride, data + i * block_stride);
660 }
661 }
662
663private:
665 template<typename scalar_type, direction dir>
666 void make_plan(std::unique_ptr<plan_fftw_r2r<scalar_type, prepost_processor, dir>> &plan) const{
667 if (not plan){
668 if (dist == 0)
669 plan = std::unique_ptr<plan_fftw_r2r<scalar_type, prepost_processor, dir>>(new plan_fftw_r2r<scalar_type, prepost_processor, dir>(size, size2, howmanyffts));
670 else if (size2 == 0)
671 plan = std::unique_ptr<plan_fftw_r2r<scalar_type, prepost_processor, dir>>(new plan_fftw_r2r<scalar_type, prepost_processor, dir>(size, howmanyffts, stride, dist));
672 else
673 plan = std::unique_ptr<plan_fftw_r2r<scalar_type, prepost_processor, dir>>(new plan_fftw_r2r<scalar_type, prepost_processor, dir>(size, size2, embed, howmanyffts, stride, dist));
674 }
675 }
676
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;
683};
684
691template<> struct one_dim_backend<backend::fftw>{
696};
697
728template<> struct one_dim_backend<backend::fftw_cos1>{
732 using executor_r2c = void;
733};
740template<> struct one_dim_backend<backend::fftw_sin1>{
744 using executor_r2c = void;
745};
746
751template<> struct default_plan_options<backend::fftw>{
753 static const bool use_reorder = true;
754};
755
760template<> struct default_plan_options<backend::fftw_cos>{
762 static const bool use_reorder = true;
763};
768template<> struct default_plan_options<backend::fftw_sin>{
770 static const bool use_reorder = true;
771};
776template<> struct default_plan_options<backend::fftw_cos1>{
778 static const bool use_reorder = true;
779};
784template<> struct default_plan_options<backend::fftw_sin1>{
786 static const bool use_reorder = true;
787};
788
789}
790
791#endif
792
793#endif /* HEFFTE_BACKEND_FFTW_H */
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
A generic container that describes a 3d box of indexes.
Definition heffte_geometry.h:67
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
fftwf_plan plan_type
Plan type.
Definition heffte_backend_fftw.h:519
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
void executor_r2c
There is no real-to-complex variant.
Definition heffte_backend_fftw.h:732
void executor_r2c
There is no real-to-complex variant.
Definition heffte_backend_fftw.h:708
void executor_r2c
There is no real-to-complex variant.
Definition heffte_backend_fftw.h:744
void executor_r2c
There is no real-to-complex variant.
Definition heffte_backend_fftw.h:720
Indicates the structure that will be used by the fft backend.
Definition heffte_common.h:663
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