7#ifndef HEFFTE_BACKEND_MKL_H
8#define HEFFTE_BACKEND_MKL_H
10#include "heffte_r2r_executor.h"
12#ifdef Heffte_ENABLE_MKL
83template<
typename scalar_type>
94 static_assert(std::is_same<scalar_type, std::complex<float>>::value
95 or std::is_same<scalar_type, std::complex<double>>::value,
96 "plan_mkl requires std::complex scalar_type with float or double, see plan_mkl_r2c for real types");
120 static_assert(std::is_same<scalar_type, std::complex<float>>::value
121 or std::is_same<scalar_type, std::complex<double>>::value,
122 "plan_mkl requires std::complex scalar_type with float or double, see plan_mkl_r2c for real types");
184 template<
typename index>
191 block_stride(
box.osize(0) *
box.osize(1)),
192 total_size(
box.count()),
196 template<
typename index>
199 blocks(1), block_stride(0), total_size(
box.count())
208 howmanyffts =
box.size[2];
210 stride =
box.size[0];
213 howmanyffts =
box.size[0];
218 howmanyffts =
box.size[1];
222 template<
typename index>
224 size(
box.size[0]), size2(
box.size[1]), howmanyffts(
box.size[2]),
226 blocks(1), block_stride(0),
227 total_size(
box.count()),
232 void forward(std::complex<float>
data[], std::complex<float>*)
const override{
234 for(
int i=0;
i<blocks;
i++){
240 void backward(std::complex<float>
data[], std::complex<float>*)
const override{
242 for(
int i=0;
i<blocks;
i++){
248 void forward(std::complex<double>
data[], std::complex<double>*)
const override{
250 for(
int i=0;
i<blocks;
i++){
256 void backward(std::complex<double>
data[], std::complex<double>*)
const override{
258 for(
int i=0;
i<blocks;
i++){
265 void forward(
float const indata[], std::complex<float>
outdata[], std::complex<float> *workspace)
const override{
275 void forward(
double const indata[], std::complex<double>
outdata[], std::complex<double> *workspace)
const override{
286 int box_size()
const override{
return total_size; }
292 template<
typename scalar_type>
298 plan = std::unique_ptr<plan_mkl<scalar_type>>(
new plan_mkl<scalar_type>(size, howmanyffts, stride, dist));
300 plan = std::unique_ptr<plan_mkl<scalar_type>>(
new plan_mkl<scalar_type>(size, size2, embed, howmanyffts, dist));
304 int size, size2, howmanyffts, stride, dist, blocks, block_stride, total_size;
305 std::array<MKL_LONG, 2> embed;
306 mutable std::unique_ptr<plan_mkl<std::complex<float>>> cplan;
307 mutable std::unique_ptr<plan_mkl<std::complex<double>>> zplan;
314template<
typename scalar_type, direction dir>
327 static_assert(std::is_same<scalar_type, float>::value
or std::is_same<scalar_type, double>::value,
328 "plan_mkl_r2c requires scalar_type with float or double, see plan_mkl for complex types");
377 template<
typename index>
385 rblock_stride(
box.osize(0) *
box.osize(1)),
386 cblock_stride(
box.osize(0) * (
box.osize(1)/2 + 1)),
394 for(
int i=0;
i<blocks;
i++){
395 float *
rdata =
const_cast<float*
>(
indata +
i * rblock_stride);
402 make_plan(sbackward);
403 for(
int i=0;
i<blocks;
i++){
411 for(
int i=0;
i<blocks;
i++){
412 double *
rdata =
const_cast<double*
>(
indata +
i * rblock_stride);
419 make_plan(dbackward);
420 for(
int i=0;
i<blocks;
i++){
435 template<
typename scalar_type, direction dir>
440 int size, howmanyffts, stride, blocks;
441 int rdist, cdist, rblock_stride, cblock_stride, rsize, csize;
442 mutable std::unique_ptr<plan_mkl_r2c<float, direction::forward>> sforward;
443 mutable std::unique_ptr<plan_mkl_r2c<double, direction::forward>> dforward;
444 mutable std::unique_ptr<plan_mkl_r2c<float, direction::backward>> sbackward;
445 mutable std::unique_ptr<plan_mkl_r2c<double, direction::backward>> dbackward;
491 static const bool use_reorder =
false;
499 static const bool use_reorder =
true;
507 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 mkl API for real-to-complex transform with shortening of the data.
Definition heffte_backend_mkl.h:366
int box_size() const override
Returns the size of the box with real data.
Definition heffte_backend_mkl.h:427
int complex_size() const override
Returns the size of the box with complex coefficients.
Definition heffte_backend_mkl.h:429
void forward(double const indata[], std::complex< double > outdata[], std::complex< double > *) const override
Forward transform, double precision.
Definition heffte_backend_mkl.h:409
void forward(float const indata[], std::complex< float > outdata[], std::complex< float > *) const override
Forward transform, single precision.
Definition heffte_backend_mkl.h:392
size_t workspace_size() const override
Return the size of the needed workspace.
Definition heffte_backend_mkl.h:431
mkl_executor_r2c(void *, box3d< index > const box, int dimension)
Constructor defines the box and the dimension of reduction.
Definition heffte_backend_mkl.h:378
void backward(std::complex< double > indata[], double outdata[], std::complex< double > *) const override
Backward transform, double precision.
Definition heffte_backend_mkl.h:418
void backward(std::complex< float > indata[], float outdata[], std::complex< float > *) const override
Backward transform, single precision.
Definition heffte_backend_mkl.h:401
Wrapper around the MKL API.
Definition heffte_backend_mkl.h:175
void forward(std::complex< double > data[], std::complex< double > *) const override
Forward fft, double-complex case.
Definition heffte_backend_mkl.h:248
void forward(std::complex< float > data[], std::complex< float > *) const override
Forward fft, float-complex case.
Definition heffte_backend_mkl.h:232
size_t workspace_size() const override
Return the size of the needed workspace.
Definition heffte_backend_mkl.h:288
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_mkl.h:275
mkl_executor(void *, box3d< index > const box)
Merges three FFTs into one.
Definition heffte_backend_mkl.h:223
int box_size() const override
Returns the size of the box.
Definition heffte_backend_mkl.h:286
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_mkl.h:280
mkl_executor(void *, box3d< index > const box, int dir1, int dir2)
Merges two FFTs into one.
Definition heffte_backend_mkl.h:197
mkl_executor(void *, box3d< index > const box, int dimension)
Constructor, specifies the box and dimension.
Definition heffte_backend_mkl.h:185
void backward(std::complex< float > data[], std::complex< float > *) const override
Backward fft, float-complex case.
Definition heffte_backend_mkl.h:240
void backward(std::complex< double > data[], std::complex< double > *) const override
Backward fft, double-complex case.
Definition heffte_backend_mkl.h:256
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_mkl.h:265
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_mkl.h:270
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.
void check_error(MKL_LONG status, std::string const &function_name)
Checks the status of a call to the MKL backend.
Definition heffte_backend_mkl.h:69
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
Allows to define whether a specific backend interface has been enabled.
Definition heffte_common.h:226
Type-tag for the Cosine Transform using the MKL FFT backend.
Definition heffte_common.h:151
Type-tag for the Sine Transform using the MKL FFT backend.
Definition heffte_common.h:156
Type-tag for the MKL backend.
Definition heffte_common.h:146
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
Unlike the C2C plan R2C is non-symmetric and it requires that the direction is built into the plan.
Definition heffte_backend_mkl.h:315
DFTI_DESCRIPTOR_HANDLE plan
Identical to the float-complex specialization.
Definition heffte_backend_mkl.h:355
plan_mkl_r2c(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_mkl.h:325
~plan_mkl_r2c()
Identical to the float-complex specialization.
Definition heffte_backend_mkl.h:351
Base plan for backend::mkl, works only for float and double complex.
Definition heffte_backend_mkl.h:84
plan_mkl(int size, int howmanyffts, int stride, int dist)
Constructor, takes inputs identical to MKL FFT descriptors.
Definition heffte_backend_mkl.h:93
~plan_mkl()
Destructor, deletes the plan.
Definition heffte_backend_mkl.h:156
plan_mkl(size_t size1, size_t size2, std::array< MKL_LONG, 2 > const &embed, size_t howmanyffts, size_t dist)
Constructor, takes inputs identical to MKL FFT descriptors.
Definition heffte_backend_mkl.h:119
DFTI_DESCRIPTOR_HANDLE plan
The MKL opaque structure (pointer to struct).
Definition heffte_backend_mkl.h:160
plan_mkl(int size1, int size2, int size3)
Constructor, takes inputs identical to MKL FFT descriptors.
Definition heffte_backend_mkl.h:145