Highly Efficient FFT for Exascale: HeFFTe v2.4
Loading...
Searching...
No Matches
heffte_backend_stock.h
1/*
2 -- heFFTe --
3 Univ. of Tennessee, Knoxville
4 @date
5*/
6
7#ifndef HEFFTE_BACKEND_STOCK_FFT_H
8#define HEFFTE_BACKEND_STOCK_FFT_H
9
10#include "heffte_r2r_executor.h"
11
12#include "stock_fft/heffte_stock_tree.h"
13
21namespace heffte{
22
23namespace backend{
28 template<> struct is_enabled<stock> : std::true_type{};
33 template<> struct is_enabled<stock_cos> : std::true_type{};
38 template<> struct is_enabled<stock_sin> : std::true_type{};
43 template<> struct is_enabled<stock_cos1> : std::true_type{};
44
45// Specialization is not necessary since the default behavior assumes CPU parameters.
46// template<>
47// struct buffer_traits<stock>{
48// using location = tag::cpu;
49// template<typename T> using container = std::vector<T>;
50// };
51}
52
54template<> struct is_ccomplex<stock::Complex<float, 1>> : std::true_type{};
55template<> struct is_zcomplex<stock::Complex<double, 1>> : std::true_type{};
56#ifdef Heffte_ENABLE_AVX
61template<> struct is_ccomplex<stock::Complex<float, 4>> : std::true_type{};
62template<> struct is_ccomplex<stock::Complex<float, 8>> : std::true_type{};
63#ifdef Heffte_ENABLE_AVX512
64template<> struct is_ccomplex<stock::Complex<float, 16>> : std::true_type{};
65#endif // Heffte_ENABLE_AVX512
66
71template<> struct is_zcomplex<stock::Complex<double, 2>> : std::true_type{};
72template<> struct is_zcomplex<stock::Complex<double, 4>> : std::true_type{};
73#ifdef Heffte_ENABLE_AVX512
74template<> struct is_zcomplex<stock::Complex<double, 8>> : std::true_type{};
75#endif // Heffte_ENABLE_AVX512
76#endif // Heffte_ENABLE_AVX
77
78#ifdef Heffte_ENABLE_AVX
79
81template<typename F, int L>
82stock::Complex<F, L> copy_pad(std::complex<F> const *c, int c_len, int i_stride) {
83 if(c_len > L/2) {
84 throw std::runtime_error("Invalid padding requested!\n");
85 }
86 std::complex<F> ret[L];
87 for(int i = 0; i < c_len; i++) ret[i] = c[i*i_stride];
88 return stock::Complex<F,L> (ret);
89}
91template<typename F, int L>
92stock::Complex<F,L> copy_pad(F const *c, int c_len, int i_stride) {
93 if(c_len > L/2) {
94 throw std::runtime_error("Invalid padding requested!\n");
95 }
96 std::complex<F> ret[L];
97 for(int i = 0; i < c_len; i++) ret[i] = std::complex<F>(c[i*i_stride], 0.0);
98 return stock::Complex<F,L> (ret);
99}
100
101template<typename F> struct pack_size { };
102#ifdef Heffte_ENABLE_AVX512
103template<> struct pack_size<float> {constexpr static int size = 16;};
104template<> struct pack_size<double>{constexpr static int size = 8;};
105#else
106template<> struct pack_size<float> {constexpr static int size = 8;};
107template<> struct pack_size<double>{constexpr static int size = 4;};
108#endif // Heffte_ENABLE_AVX512
109
119template<typename F, direction dir>
120struct plan_stock_fft{
130 plan_stock_fft(int size, int howmanyffts, int stride, int rdist, int cdist):
131 N(size), num_ffts(howmanyffts), stride_sz(stride), real_d(rdist), comp_d(cdist) {
132 numNodes = stock::getNumNodes(N);
133 root = std::unique_ptr<stock::biFuncNode<F,L>[]>(new stock::biFuncNode<F,L>[numNodes]);
134 init_fft_tree(root.get(), N);
135 }
137 int N, num_ffts, stride_sz, real_d, comp_d, numNodes;
138 constexpr static int L = pack_size<F>::size;
139 std::unique_ptr<stock::biFuncNode<F,L>[]> root;
141 void execute(std::complex<F> const idata[], F odata[]) {
142 // Allocate input and output temporary arrays
143 stock::complex_vector<F,L> inp (N);
144 stock::complex_vector<F,L> out (N);
145
146 // Perform batch transform on everything save for the remainder
147 for(int p = 0; p+((L/2)-1) < num_ffts; p += (L/2)) {
148 // Convert types
149 inp[0] = stock::Complex<F,L> (&idata[p*comp_d], comp_d);
150 for(int i = 1; i < (N+2)/2; i++) {
151 int idx = p*comp_d + i*stride_sz;
152 inp[i] = stock::Complex<F,L> (&idata[idx], comp_d);
153 inp[N-i] = inp[i].conjugate();
154 }
155 // Perform fft
156 root[0].fptr(inp.data(), out.data(), 1, 1, root.get(), dir);
157 // Convert type back
158 for(int i = 0; i < N; i++) {
159 int idx = p*real_d + i*stride_sz;
160 stock::Complex<F,L> arr = out[i];
161 for(int j = 0; j < L/2; j++) odata[idx+j*real_d] = arr[j].real();
162 }
163 }
164
165 // Handle remainder
166 if(num_ffts % (L/2) > 0) {
167 int rem = num_ffts % (L/2);
168 // Init p for ease of use
169 int p = num_ffts - rem;
170 inp[0] = copy_pad<F,L>(&idata[p*comp_d], rem, comp_d);
171 for(int i = 1; i < (N+2)/2; i++) {
172 int idx = p*comp_d + i*stride_sz;
173 inp[i] = copy_pad<F,L>(&idata[idx], rem, comp_d);
174 inp[N-i] = inp[i].conjugate();
175 }
176 root[0].fptr(inp.data(), out.data(), 1, 1, root.get(), dir);
177 for(int i = 0; i < N; i++) {
178 int idx = p*real_d + i*stride_sz;
179 stock::Complex<F,L> arr = out[i];
180 // Only need first k columns
181 for(int k = 0; k < rem; k++) {
182 odata[idx + k*real_d] = arr[k].real();
183 }
184 }
185 }
186 }
188 void execute(F const idata[], std::complex<F> odata[]) {
189 // Allocate input and output temporary arrays
190 stock::complex_vector<F,L> inp (N);
191 stock::complex_vector<F,L> out (N);
192
193 // Perform batch transform on everything save for the remainder
194 for(int p = 0; p+((L/2)-1) < num_ffts; p += L/2) {
195 // Convert types
196 for(int i = 0; i < N; i++) {
197 int idx = p*real_d + i*stride_sz;
198 inp[i] = copy_pad<F,L>(&idata[idx], L/2, real_d);
199 }
200 // Perform fft
201 root[0].fptr(inp.data(), out.data(), 1, 1, root.get(), dir);
202 // Convert type back
203 for(int i = 0; i < (N+2)/2; i++) {
204 int idx = p*comp_d + i*stride_sz;
205 stock::Complex<F,L> arr = out[i];
206 for(int j = 0; j < L/2; j++) odata[idx+j*comp_d] = arr[j];
207 }
208 }
209
210 // Handle remainder
211 if(num_ffts % (L/2) > 0) {
212 int rem = num_ffts % (L/2);
213 // Init p for ease of use
214 int p = num_ffts - rem;
215 for(int i = 0; i < N; i++) {
216 int idx = p*real_d + i*stride_sz;
217 // remainder columns are all zeros
218 inp[i] = copy_pad<F,L>(&idata[idx], rem, real_d);
219 }
220 root[0].fptr(inp.data(), out.data(), 1, 1, root.get(), dir);
221 for(int i = 0; i < (N+2)/2; i++) {
222 int idx = p*comp_d + i*stride_sz;
223 stock::Complex<F,L> arr = out[i];
224 // Only need first k columns
225 for(int k = 0; k < rem; k++) {
226 odata[idx + k*comp_d] = arr[k];
227 }
228 }
229 }
230 }
231};
232
239template<typename F, direction dir>
240struct plan_stock_fft<std::complex<F>, dir>{
249 plan_stock_fft(int size, int howmanyffts, int stride, int dist):
250 N(size), num_ffts(howmanyffts), stride_sz(stride), idist(dist), odist(dist) {
251 numNodes = stock::getNumNodes(N);
252 root = std::unique_ptr<stock::biFuncNode<F,L>[]>(new stock::biFuncNode<F,L>[numNodes]);
253 init_fft_tree(root.get(), N);
254 }
255
256 int N, num_ffts, stride_sz, idist, odist, numNodes;
257 constexpr static int L = pack_size<F>::size;
258 std::unique_ptr<stock::biFuncNode<F, L>[]> root;
260 void execute(std::complex<F> data[]) {
261 // Allocate input and output temporary arrays
262 stock::complex_vector<F,L> inp (N);
263 stock::complex_vector<F,L> out (N);
264
265 // Perform batch transform on everything save for the remainder
266 for(int p = 0; p + (L/2 - 1) < num_ffts; p += L/2) {
267 // Convert types
268 for(int i = 0; i < N; i++) {
269 int idx = p*idist + i*stride_sz;
270 inp[i] = stock::Complex<F,L> (&data[idx], idist);
271 }
272 // Perform fft
273 root[0].fptr(inp.data(), out.data(), 1, 1, root.get(), dir);
274 // Convert type back
275 for(int i = 0; i < N; i++) {
276 int idx = p*odist + i*stride_sz;
277 stock::Complex<F,L> arr = out[i];
278 for(int j = 0; j < L/2; j++) data[idx+j*odist] = arr[j];
279 }
280 }
281
282 // Handle remainder
283 if(num_ffts % (L/2) > 0) {
284 int rem = num_ffts % (L/2);
285 // Init p for ease of use
286 int p = num_ffts - rem;
287 for(int i = 0; i < N; i++) {
288 int idx = p*idist + i*stride_sz;
289 // remainder columns are all zeros
290 inp[i] = copy_pad<F,L>(&data[idx], rem, idist);
291 }
292 root[0].fptr(inp.data(), out.data(), 1, 1, root.get(), dir);
293 for(int i = 0; i < N; i++) {
294 int idx = p*odist + i*stride_sz;
295 stock::Complex<F,L> arr = out[i];
296 // Only need first k columns
297 for(int k = 0; k < rem; k++) {
298 data[idx + k*odist] = arr[k];
299 }
300 }
301 }
302 }
303};
304
305#else
306
307// In the case the user does not have AVX installed
308
313template<typename F, direction dir>
324 plan_stock_fft(int size, int howmanyffts, int stride, int rdist, int cdist):
325 N(size), num_ffts(howmanyffts), stride_sz(stride), real_d(rdist), comp_d(cdist) {
326 numNodes = stock::getNumNodes(N);
327 root = std::unique_ptr<stock::biFuncNode<F,1>[]>(new stock::biFuncNode<F,1>[numNodes]);
328 init_fft_tree(root.get(), N);
329 }
331 int N, num_ffts, stride_sz, real_d, comp_d, numNodes;
332 std::unique_ptr<stock::biFuncNode<F, 1>[]> root;
334 void execute(std::complex<F> const idata[], F odata[]) {
335 // Allocate input and output temporary arrays
336 stock::complex_vector<F,1> inp (N);
337 stock::complex_vector<F,1> out (N);
338
339 // Perform batch transform on everything save for the remainder
340 for(int p = 0; p < num_ffts; p++) {
341 // Convert types
342 inp[0] = stock::Complex<F,1> {idata[p*comp_d]};
343 for(int i = 1; i < (N+2)/2; i++) {
344 int idx = p*comp_d + i*stride_sz;
346 inp[N-i] = inp[i].conjugate();
347 }
348 // Perform fft
349 root[0].fptr(inp.data(), out.data(), 1, 1, root.get(), dir);
350 // Convert type back
351 for(int i = 0; i < N; i++) {
352 int idx = p*real_d + i*stride_sz;
354 odata[idx+0*real_d] = arr[0].real();
355 }
356 }
357 }
359 void execute(F const idata[], std::complex<F> odata[]) {
360 // Allocate input and output temporary arrays
361 stock::complex_vector<F,1> inp (N);
362 stock::complex_vector<F,1> out (N);
363
364 // Perform batch transform on everything save for the remainder
365 for(int p = 0; p < num_ffts; p++) {
366 // Convert types
367 for(int i = 0; i < N; i++) {
368 int idx = p*real_d + i*stride_sz;
369 inp[i] = stock::Complex<F,1> {idata[idx+0*real_d], 0};
370 }
371 // Perform fft
372 root[0].fptr(inp.data(), out.data(), 1, 1, root.get(), dir);
373 // Convert type back
374 for(int i = 0; i < (N+2)/2; i++) {
375 int idx = p*comp_d + i*stride_sz;
377 odata[idx+0*comp_d] = arr[0];
378 }
379 }
380 }
381};
382
389template<typename F, direction dir>
390struct plan_stock_fft<std::complex<F>, dir>{
399 plan_stock_fft(int size, int howmanyffts, int stride, int dist):
400 N(size), num_ffts(howmanyffts), stride_sz(stride), idist(dist), odist(dist) {
401 numNodes = stock::getNumNodes(N);
402 root = std::unique_ptr<stock::biFuncNode<F,1>[]>(new stock::biFuncNode<F,1>[numNodes]);
403 init_fft_tree(root.get(), N);
404 }
405
406 int N, num_ffts, stride_sz, idist, odist, numNodes;
407 std::unique_ptr<stock::biFuncNode<F,1>[]> root;
409 void execute(std::complex<F> data[]) {
410 // Allocate input and output temporary arrays
411 stock::complex_vector<F,1> inp (N);
412 stock::complex_vector<F,1> out (N);
413
414 // Perform batch transform on everything save for the remainder
415 for(int p = 0; p < num_ffts; p++) {
416 // Convert types
417 for(int i = 0; i < N; i++) {
418 int idx = p*idist + i*stride_sz;
419 inp[i] = stock::Complex<F, 1> {data[idx+0*idist]};
420 }
421 // Perform fft
422 root[0].fptr(inp.data(), out.data(), 1, 1, root.get(), dir);
423 // Convert type back
424 for(int i = 0; i < N; i++) {
425 int idx = p*odist + i*stride_sz;
427 data[idx+0*odist] = arr[0];
428 }
429 }
430 }
431};
432#endif // Heffte_ENABLE_AVX512
433
447public:
455 template<typename index>
457 size(box.size[dimension]),
458 num_ffts(fft1d_get_howmany(box, dimension)),
460 dist((dimension == box.order[0]) ? size : 1),
461 blocks((dimension == box.order[1]) ? box.osize(2) : 1),
462 block_stride(box.osize(0) * box.osize(1)),
463 total_size(box.count())
464 {}
466 template<typename index> stock_fft_executor(void*, box3d<index> const, int, int){
467 throw std::runtime_error("2D transforms for the stock backend are not available yet!");
468 }
470 template<typename index> stock_fft_executor(void*, box3d<index> const){
471 throw std::runtime_error("3D transforms for the stock backend are not available yet!");
472 }
473
475 void forward(std::complex<float> data[], std::complex<float>*) const override{
476 make_plan(cforward);
477 for(int i=0; i<blocks; i++) {
478 cforward->execute(data + i * block_stride);
479 }
480 }
482 void backward(std::complex<float> data[], std::complex<float>*) const override{
483 make_plan(cbackward);
484 for(int i=0; i<blocks; i++) {
485 cbackward->execute(data + i * block_stride);
486 }
487 }
489 void forward(std::complex<double> data[], std::complex<double>*) const override{
490 make_plan(zforward);
491 for(int i=0; i<blocks; i++) {
492 zforward->execute(data + i * block_stride);
493 }
494 }
496 void backward(std::complex<double> data[], std::complex<double>*) const override{
497 make_plan(zbackward);
498 for(int i=0; i<blocks; i++) {
499 zbackward->execute(data + i * block_stride);
500 }
501 }
502
504 void forward(float const indata[], std::complex<float> outdata[], std::complex<float> *workspace) const override{
506 forward(outdata, workspace);
507 }
509 void backward(std::complex<float> indata[], float outdata[], std::complex<float> *workspace) const override{
510 backward(indata, workspace);
511 for(int i=0; i<total_size; i++) outdata[i] = std::real(indata[i]);
512 }
514 void forward(double const indata[], std::complex<double> outdata[], std::complex<double> *workspace) const override{
516 forward(outdata, workspace);
517 }
519 void backward(std::complex<double> indata[], double outdata[], std::complex<double> *workspace) const override{
520 backward(indata, workspace);
521 for(int i=0; i<total_size; i++) outdata[i] = std::real(indata[i]);
522 }
523
525 int box_size() const override{ return total_size; }
526
528 size_t workspace_size() const override{ return 0; }
529
530private:
532 template<typename scalar_type, direction dir>
533 void make_plan(std::unique_ptr<plan_stock_fft<scalar_type, dir>> &plan) const{
534 if (!plan) plan = std::unique_ptr<plan_stock_fft<scalar_type, dir>>(new plan_stock_fft<scalar_type, dir>(size, num_ffts, stride, dist));
535 }
536
537 int size, num_ffts, stride, dist, blocks, block_stride, total_size;
538 mutable std::unique_ptr<plan_stock_fft<std::complex<float>, direction::forward>> cforward;
539 mutable std::unique_ptr<plan_stock_fft<std::complex<float>, direction::backward>> cbackward;
540 mutable std::unique_ptr<plan_stock_fft<std::complex<double>, direction::forward>> zforward;
541 mutable std::unique_ptr<plan_stock_fft<std::complex<double>, direction::backward>> zbackward;
542};
543
553public:
563 template<typename index>
565 size(box.size[dimension]),
566 num_ffts(fft1d_get_howmany(box, dimension)),
568 blocks((dimension == box.order[1]) ? box.osize(2) : 1),
569 rdist((dimension == box.order[0]) ? size : 1),
570 cdist((dimension == box.order[0]) ? size/2 + 1 : 1),
571 rblock_stride(box.osize(0) * box.osize(1)),
572 cblock_stride(box.osize(0) * (box.osize(1)/2 + 1)),
573 rsize(box.count()),
574 csize(box.r2c(dimension).count())
575 {}
576
578 void forward(float const indata[], std::complex<float> outdata[], std::complex<float>*) const override{
579 make_plan(sforward);
580 for(int i=0; i<blocks; i++) {
581 sforward->execute(indata + i*rblock_stride, outdata + i*cblock_stride);
582 }
583 }
585 void backward(std::complex<float> indata[], float outdata[], std::complex<float>*) const override{
586 make_plan(sbackward);
587 for(int i=0; i<blocks; i++) {
588 sbackward->execute(indata + i*cblock_stride, outdata + i*rblock_stride);
589 }
590 }
592 void forward(double const indata[], std::complex<double> outdata[], std::complex<double>*) const override{
593 make_plan(dforward);
594 for(int i=0; i<blocks; i++) {
595 dforward->execute(indata + i*rblock_stride, outdata + i*cblock_stride);
596 }
597 }
599 void backward(std::complex<double> indata[], double outdata[], std::complex<double>*) const override{
600 make_plan(dbackward);
601 for(int i=0; i<blocks; i++) {
602 dbackward->execute(indata + i*cblock_stride, outdata + i*rblock_stride);
603 }
604 }
605
607 int box_size() const override{ return rsize; }
609 int complex_size() const override{ return csize; }
611 size_t workspace_size() const override{ return 0; }
612
613private:
615 template<typename scalar_type, direction dir>
616 void make_plan(std::unique_ptr<plan_stock_fft<scalar_type, dir>> &plan) const{
617 if (!plan) plan = std::unique_ptr<plan_stock_fft<scalar_type, dir>>(new plan_stock_fft<scalar_type, dir>(size, num_ffts, stride, rdist, cdist));
618 }
619
620 int size, num_ffts, stride, blocks;
621 int rdist, cdist, rblock_stride, cblock_stride, rsize, csize;
622 mutable std::unique_ptr<plan_stock_fft<float, direction::forward>> sforward;
623 mutable std::unique_ptr<plan_stock_fft<double, direction::forward>> dforward;
624 mutable std::unique_ptr<plan_stock_fft<float, direction::backward>> sbackward;
625 mutable std::unique_ptr<plan_stock_fft<double, direction::backward>> dbackward;
626};
633template<> struct one_dim_backend<backend::stock>{
638};
675
680template<> struct default_plan_options<backend::stock>{
682 static const bool use_reorder = true;
683};
688template<> struct default_plan_options<backend::stock_cos>{
690 static const bool use_reorder = true;
691};
696template<> struct default_plan_options<backend::stock_sin>{
698 static const bool use_reorder = true;
699};
704template<> struct default_plan_options<backend::stock_cos1>{
706 static const bool use_reorder = true;
707};
708
709}
710
711#endif /* HEFFTE_BACKEND_STOCK_FFT_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
Custom complex type taking advantage of vectorization A Complex Type intrinsic to HeFFTe that takes a...
Definition heffte_stock_complex.h:29
Wrapper to stock API for real-to-complex transform with shortening of the data.
Definition heffte_backend_stock.h:552
int box_size() const override
Returns the size of the box with real data.
Definition heffte_backend_stock.h:607
int complex_size() const override
Returns the size of the box with complex coefficients.
Definition heffte_backend_stock.h:609
size_t workspace_size() const override
Return the size of the needed workspace.
Definition heffte_backend_stock.h:611
void backward(std::complex< float > indata[], float outdata[], std::complex< float > *) const override
Backward transform, single precision.
Definition heffte_backend_stock.h:585
void forward(float const indata[], std::complex< float > outdata[], std::complex< float > *) const override
Forward transform, single precision.
Definition heffte_backend_stock.h:578
stock_fft_executor_r2c(void *, box3d< index > const box, int dimension)
Constructor defines the box and the dimension of reduction.
Definition heffte_backend_stock.h:564
void forward(double const indata[], std::complex< double > outdata[], std::complex< double > *) const override
Forward transform, double precision.
Definition heffte_backend_stock.h:592
void backward(std::complex< double > indata[], double outdata[], std::complex< double > *) const override
Backward transform, double precision.
Definition heffte_backend_stock.h:599
Wrapper around the Stock FFT API.
Definition heffte_backend_stock.h:446
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_stock.h:504
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_stock.h:509
size_t workspace_size() const override
Return the size of the needed workspace.
Definition heffte_backend_stock.h:528
void forward(std::complex< double > data[], std::complex< double > *) const override
Forward fft, double-complex case.
Definition heffte_backend_stock.h:489
void forward(std::complex< float > data[], std::complex< float > *) const override
Forward fft, float-complex case.
Definition heffte_backend_stock.h:475
void backward(std::complex< double > data[], std::complex< double > *) const override
Backward fft, double-complex case.
Definition heffte_backend_stock.h:496
stock_fft_executor(void *, box3d< index > const, int, int)
Placeholder, unimplemented.
Definition heffte_backend_stock.h:466
void backward(std::complex< float > data[], std::complex< float > *) const override
Backward fft, float-complex case.
Definition heffte_backend_stock.h:482
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_stock.h:514
stock_fft_executor(void *, box3d< index > const)
Placeholder, unimplemented.
Definition heffte_backend_stock.h:470
int box_size() const override
Returns the size of the box.
Definition heffte_backend_stock.h:525
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_stock.h:519
stock_fft_executor(void *, box3d< index > const box, int dimension)
Constructor, specifies the box and dimension.
Definition heffte_backend_stock.h:456
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
@ backward
Inverse DFT transform.
@ forward
Forward DFT transform.
Namespace containing all HeFFTe methods and classes.
Definition heffte_backend_cuda.h:38
Allows to define whether a specific backend interface has been enabled.
Definition heffte_common.h:226
Type-tag for the Cosine Transform type 1 using the stock FFT backend.
Definition heffte_common.h:140
Type-tag for the Cosine Transform using the stock FFT backend.
Definition heffte_common.h:130
Type-tag for the Sine Transform using the stock FFT backend.
Definition heffte_common.h:135
Type-tag for the stock FFT backend.
Definition heffte_common.h:125
Defines a set of default plan options for a given backend.
Definition heffte_common.h:761
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
plan_stock_fft(int size, int howmanyffts, int stride, int dist)
Constructor to plan out an FFT using the stock backend.
Definition heffte_backend_stock.h:399
void execute(std::complex< F > data[])
Execute an FFT inplace on std::complex<F> data.
Definition heffte_backend_stock.h:409
Specialization for r2c single precision.
Definition heffte_backend_stock.h:314
void execute(std::complex< F > const idata[], F odata[])
Execute C2R FFT.
Definition heffte_backend_stock.h:334
void execute(F const idata[], std::complex< F > odata[])
Execute R2C FFT.
Definition heffte_backend_stock.h:359
plan_stock_fft(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_stock.h:324
int N
Identical to the F-complex specialization.
Definition heffte_backend_stock.h:331
Class to represent the call-graph nodes.
Definition heffte_stock_algos.h:78