Hipace
|
Wrapper around multiple FFT libraries. More...
Classes | |
struct | DSTplan |
This struct contains the vendor FFT plan and additional metadata. More... | |
Typedefs | |
using | DSTplans = amrex::LayoutData< DSTplan > |
Enumerations | |
enum class | direction { forward , backward } |
Functions | |
DSTplan | CreatePlan (const amrex::IntVect &real_size, amrex::FArrayBox *position_array, amrex::FArrayBox *fourier_array) |
create FFT plan for the backend FFT library. More... | |
void | DestroyPlan (DSTplan &dst_plan) |
Destroy library FFT plan. More... | |
template<AnyDST::direction d> | |
void | Execute (DSTplan &dst_plan) |
Perform FFT with backend library. More... | |
template void | Execute< direction::forward > (DSTplan &dst_plan) |
template void | Execute< direction::backward > (DSTplan &dst_plan) |
void | ExpandR2R (amrex::FArrayBox &dst, amrex::FArrayBox &src) |
Extend src into a symmetrized larger array dst. More... | |
void | ShrinkC2R (amrex::FArrayBox &dst, amrex::BaseFab< amrex::GpuComplex< amrex::Real >> &src) |
Extract symmetrical src array into smaller array dst. More... | |
void | ToComplex (const amrex::Real *const AMREX_RESTRICT in, amrex::GpuComplex< amrex::Real > *const AMREX_RESTRICT out, const int n_data, const int n_batch) |
Make Complex array out of Real array to prepare for fft. out[idx] = -in[2*idx-2] + in[2*idx] + i*in[2*idx-1] for each column with in[-1] = 0; in[-2] = -in[0]; in[n_data] = 0; in[n_data+1] = -in[n_data-1]. More... | |
void | C2Rfft (AnyFFT::VendorFFTPlan &plan, amrex::GpuComplex< amrex::Real > *AMREX_RESTRICT in, amrex::Real *const AMREX_RESTRICT out) |
Complex to Real fft for every column of the input matrix. The output Matrix has its indexes reversed compared to some other libraries. More... | |
void | ToSine (const amrex::Real *const AMREX_RESTRICT in, amrex::Real *const AMREX_RESTRICT out, const int n_data, const int n_batch) |
Make Sine-space Real array out of array from fft. out[idx] = 0.5 *(in[n_data-idx] - in[idx+1] + (in[n_data-idx] + in[idx+1])/ (2*sin((idx+1)*pi/(n_data+1)))) for each column. More... | |
void | Transpose (const amrex::Real *const AMREX_RESTRICT in, amrex::Real *const AMREX_RESTRICT out, const int n_data, const int n_batch) |
Transpose input matrix out[idy][idx] = in[idx][idy]. More... | |
void | C2Rfft (AnyFFT::VendorFFTPlan &plan, amrex::GpuComplex< amrex::Real > *AMREX_RESTRICT in, amrex::Real *const AMREX_RESTRICT out, rocfft_execution_info execinfo) |
Variables | |
cufftType | VendorR2C = CUFFT_D2Z |
cufftType | VendorC2R = CUFFT_Z2D |
const auto | VendorCreatePlanR2R2D = fftw_plan_r2r_2d |
Wrapper around multiple FFT libraries.
The header file defines the API and the base types (Complex and VendorFFTPlan), and the implementation for different FFT libraries is done in different cpp files. This wrapper only depends on the underlying FFT library AND on AMReX (There is no dependence on WarpX).
using AnyDST::DSTplans = typedef amrex::LayoutData<DSTplan> |
Collection of FFT plans, one FFTplan per box
|
strong |
Direction in which the FFT is performed.
Enumerator | |
---|---|
forward | |
backward |
void AnyDST::C2Rfft | ( | AnyFFT::VendorFFTPlan & | plan, |
amrex::GpuComplex< amrex::Real > *AMREX_RESTRICT | in, | ||
amrex::Real *const AMREX_RESTRICT | out | ||
) |
Complex to Real fft for every column of the input matrix. The output Matrix has its indexes reversed compared to some other libraries.
[in] | plan | cuda fft plan for transformation |
[in] | in | input complex array |
[out] | out | output real array |
void AnyDST::C2Rfft | ( | AnyFFT::VendorFFTPlan & | plan, |
amrex::GpuComplex< amrex::Real > *AMREX_RESTRICT | in, | ||
amrex::Real *const AMREX_RESTRICT | out, | ||
rocfft_execution_info | execinfo | ||
) |
DSTplan AnyDST::CreatePlan | ( | const amrex::IntVect & | real_size, |
amrex::FArrayBox * | position_array, | ||
amrex::FArrayBox * | fourier_array | ||
) |
create FFT plan for the backend FFT library.
[in] | real_size | Size of the real array, along each dimension. |
[out] | position_array | Real array from/to where R2R DST is performed |
[out] | fourier_array | Real array to/from where R2R DST is performed |
void AnyDST::DestroyPlan | ( | DSTplan & | dst_plan | ) |
Destroy library FFT plan.
[out] | dst_plan | plan to destroy |
void AnyDST::Execute | ( | DSTplan & | dst_plan | ) |
Perform FFT with backend library.
[out] | dst_plan | plan for which the FFT is performed |
template void AnyDST::Execute< direction::backward > | ( | DSTplan & | dst_plan | ) |
template void AnyDST::Execute< direction::forward > | ( | DSTplan & | dst_plan | ) |
void AnyDST::ExpandR2R | ( | amrex::FArrayBox & | dst, |
amrex::FArrayBox & | src | ||
) |
Extend src into a symmetrized larger array dst.
[in,out] | dst | destination array, odd symmetry around 0 and the middle points in x and y |
[in] | src | source array |
void AnyDST::ShrinkC2R | ( | amrex::FArrayBox & | dst, |
amrex::BaseFab< amrex::GpuComplex< amrex::Real >> & | src | ||
) |
Extract symmetrical src array into smaller array dst.
[in,out] | dst | destination array |
[in] | src | destination array, symmetric in x and y |
void AnyDST::ToComplex | ( | const amrex::Real *const AMREX_RESTRICT | in, |
amrex::GpuComplex< amrex::Real > *const AMREX_RESTRICT | out, | ||
const int | n_data, | ||
const int | n_batch | ||
) |
Make Complex array out of Real array to prepare for fft. out[idx] = -in[2*idx-2] + in[2*idx] + i*in[2*idx-1] for each column with in[-1] = 0; in[-2] = -in[0]; in[n_data] = 0; in[n_data+1] = -in[n_data-1].
[in] | in | input real array |
[out] | out | output complex array |
[in] | n_data | number of (contiguous) rows in position matrix |
[in] | n_batch | number of (strided) columns in position matrix |
void AnyDST::ToSine | ( | const amrex::Real *const AMREX_RESTRICT | in, |
amrex::Real *const AMREX_RESTRICT | out, | ||
const int | n_data, | ||
const int | n_batch | ||
) |
Make Sine-space Real array out of array from fft. out[idx] = 0.5 *(in[n_data-idx] - in[idx+1] + (in[n_data-idx] + in[idx+1])/ (2*sin((idx+1)*pi/(n_data+1)))) for each column.
[in] | in | input real array |
[out] | out | output real array |
[in] | n_data | number of (contiguous) rows in position matrix |
[in] | n_batch | number of (strided) columns in position matrix |
void AnyDST::Transpose | ( | const amrex::Real *const AMREX_RESTRICT | in, |
amrex::Real *const AMREX_RESTRICT | out, | ||
const int | n_data, | ||
const int | n_batch | ||
) |
Transpose input matrix out[idy][idx] = in[idx][idy].
[in] | in | input real array |
[out] | out | output real array |
[in] | n_data | number of (contiguous) rows in input matrix |
[in] | n_batch | number of (strided) columns in input matrix |
cufftType AnyDST::VendorC2R = CUFFT_Z2D |
const auto AnyDST::VendorCreatePlanR2R2D = fftw_plan_r2r_2d |
cufftType AnyDST::VendorR2C = CUFFT_D2Z |