Parallel Discrete Fourier Transform. More...
#include <AMReX_FFT_R2C.H>
Public Types | |
using | cMF = FabArray< BaseFab< GpuComplex< T > > > |
using | MF = std::conditional_t< C, cMF, std::conditional_t< std::is_same_v< T, Real >, MultiFab, FabArray< BaseFab< T > > > > |
Public Member Functions | |
R2C (Box const &domain, Info const &info=Info{}) | |
Constructor. | |
R2C (std::array< int, AMREX_SPACEDIM > const &domain_size, Info const &info=Info{}) | |
Constructor. | |
~R2C () | |
R2C (R2C const &)=delete | |
R2C (R2C &&)=delete | |
R2C & | operator= (R2C const &)=delete |
R2C & | operator= (R2C &&)=delete |
void | setLocalDomain (std::array< int, AMREX_SPACEDIM > const &local_start, std::array< int, AMREX_SPACEDIM > const &local_size) |
Set local domain. | |
std::pair< std::array< int, AMREX_SPACEDIM >, std::array< int, AMREX_SPACEDIM > > | getLocalDomain () const |
Get local domain. | |
void | setLocalSpectralDomain (std::array< int, AMREX_SPACEDIM > const &local_start, std::array< int, AMREX_SPACEDIM > const &local_size) |
Set local spectral domain. | |
std::pair< std::array< int, AMREX_SPACEDIM >, std::array< int, AMREX_SPACEDIM > > | getLocalSpectralDomain () const |
Get local spectral domain. | |
template<typename F , Direction DIR = D, std::enable_if_t< DIR==Direction::both, int > = 0> | |
void | forwardThenBackward (MF const &inmf, MF &outmf, F const &post_forward, int incomp=0, int outcomp=0) |
Forward and then backward transform. | |
template<Direction DIR = D, std::enable_if_t< DIR==Direction::forward||DIR==Direction::both, int > = 0> | |
void | forward (MF const &inmf, int incomp=0) |
Forward transform. | |
template<Direction DIR = D, std::enable_if_t< DIR==Direction::forward||DIR==Direction::both, int > = 0> | |
void | forward (MF const &inmf, cMF &outmf, int incomp=0, int outcomp=0) |
Forward transform. | |
template<typename RT , typename CT , Direction DIR = D, bool CP = C, std::enable_if_t<(DIR==Direction::forward||DIR==Direction::both) &&((sizeof(RT) *2==sizeof(CT) &&!CP)||(sizeof(RT)==sizeof(CT) &&CP)), int > = 0> | |
void | forward (RT const *in, CT *out) |
Forward transform. | |
template<Direction DIR = D, std::enable_if_t< DIR==Direction::both, int > = 0> | |
void | backward (MF &outmf, int outcomp=0) |
Backward transform. | |
template<Direction DIR = D, std::enable_if_t< DIR==Direction::backward||DIR==Direction::both, int > = 0> | |
void | backward (cMF const &inmf, MF &outmf, int incomp=0, int outcomp=0) |
Backward transform. | |
template<typename CT , typename RT , Direction DIR = D, bool CP = C, std::enable_if_t<(DIR==Direction::backward||DIR==Direction::both) &&((sizeof(RT) *2==sizeof(CT) &&!CP)||(sizeof(RT)==sizeof(CT) &&CP)), int > = 0> | |
void | backward (CT const *in, RT *out) |
Backward transform. | |
T | scalingFactor () const |
template<Direction DIR = D, std::enable_if_t< DIR==Direction::forward||DIR==Direction::both, int > = 0> | |
std::pair< cMF *, IntVect > | getSpectralData () const |
Get the internal spectral data. | |
std::pair< BoxArray, DistributionMapping > | getSpectralDataLayout () const |
Get BoxArray and DistributionMapping for spectral data. | |
template<typename F > | |
void | post_forward_doit_0 (F const &post_forward) |
template<typename F > | |
void | post_forward_doit_1 (F const &post_forward) |
template<Direction DIR, std::enable_if_t< DIR==Direction::forward||DIR==Direction::both, int > > | |
std::pair< typename R2C< T, D, C >::cMF *, IntVect > | getSpectralData () const |
Private Member Functions | |
void | prepare_openbc () |
void | backward_doit (MF &outmf, IntVect const &ngout=IntVect(0), Periodicity const &period=Periodicity::NonPeriodic(), int outcomp=0) |
void | backward_doit (cMF const &inmf, MF &outmf, IntVect const &ngout=IntVect(0), Periodicity const &period=Periodicity::NonPeriodic(), int incomp=0, int outcomp=0) |
std::pair< Plan< T >, Plan< T > > | make_c2c_plans (cMF &inout, int ndims) const |
template<typename FA , typename RT > | |
std::pair< std::unique_ptr< char, DataDeleter >, std::size_t > | install_raw_ptr (FA &fa, RT const *p) |
Static Private Member Functions | |
static Box | make_domain_x (Box const &domain) |
static Box | make_domain_y (Box const &domain) |
static Box | make_domain_z (Box const &domain) |
static std::pair< BoxArray, DistributionMapping > | make_layout_from_local_domain (std::array< int, AMREX_SPACEDIM > const &local_start, std::array< int, AMREX_SPACEDIM > const &local_size) |
Friends | |
template<typename U > | |
class | OpenBCSolver |
template<typename U > | |
class | Poisson |
template<typename U > | |
class | PoissonHybrid |
Parallel Discrete Fourier Transform.
This class supports Fourier transforms between real and complex data. The name R2C indicates that the forward transform converts real data to complex data, while the backward transform converts complex data to real data. It should be noted that both directions of transformation are supported, not just from real to complex. The scaling follows the FFTW convention, where applying the forward transform followed by the backward transform scales the original data by the size of the input array.
The arrays are assumed to be in column-major, which is different from FFTW's row-major layout. Because the complex domain data have the Hermitian symmetry, only half of the data in the complex domain are stored. If the real domain size is nx * ny * nz, the complex domain's size will be (nx/2+1) * ny * nz.
For more details, we refer the users to https://amrex-codes.github.io/amrex/docs_html/FFT_Chapter.html.
using amrex::FFT::R2C< T, D, C >::cMF = FabArray<BaseFab<GpuComplex<T> > > |
using amrex::FFT::R2C< T, D, C >::MF = std::conditional_t <C, cMF, std::conditional_t<std::is_same_v<T,Real>, MultiFab, FabArray<BaseFab<T> > > > |
|
explicit |
Constructor.
domain | the forward domain (i.e., the domain of the real data) |
info | optional information |
|
explicit |
Constructor.
If AMREX_SPACEDIM is 3 and you want to do 2D FFT, you just need to set the size of one of the dimensions to 1.
domain_size | size of the forward domain (i.e., the real data domain) |
info | optional information |
amrex::FFT::R2C< T, D, C >::~R2C | ( | ) |
|
delete |
|
delete |
void amrex::FFT::R2C< T, D, C >::backward | ( | cMF const & | inmf, |
MF & | outmf, | ||
int | incomp = 0 , |
||
int | outcomp = 0 |
||
) |
Backward transform.
This function is not available when this class template is instantiated for forward-only transform.
inmf | input data in FabArray<BaseFab<GpuComplex<T>>> |
outmf | output data in MultiFab or FabArray<BaseFab<float>> |
void amrex::FFT::R2C< T, D, C >::backward | ( | CT const * | in, |
RT * | out | ||
) |
Backward transform.
This raw pointer version of backward requires setLocalDomain/getLocalDomain and setLocalSpectralDomain/getLocalSpectralDomain have been called already. Note that one is allowed to call this function multiple times after the set/get domain functions are called only once unless the domain decomposition changes. In fact, that is the preferred way because it has better performance. All processes need to call this function even if their local size is zero. If the local size is zero, one can pass nullptrs.
void amrex::FFT::R2C< T, D, C >::backward | ( | MF & | outmf, |
int | outcomp = 0 |
||
) |
Backward transform.
This function is available only when this class template is instantiated for transforms in both directions.
outmf | output data in MultiFab or FabArray<BaseFab<float>> |
|
private |
|
private |
void amrex::FFT::R2C< T, D, C >::forward | ( | MF const & | inmf, |
cMF & | outmf, | ||
int | incomp = 0 , |
||
int | outcomp = 0 |
||
) |
Forward transform.
This function is not available when this class template is instantiated for backward-only transform.
inmf | input data in MultiFab or FabArray<BaseFab<float>> |
outmf | output data in FabArray<BaseFab<GpuComplex<T>>> |
void amrex::FFT::R2C< T, D, C >::forward | ( | MF const & | inmf, |
int | incomp = 0 |
||
) |
Forward transform.
The output is stored in this object's internal data. This function is not available when this class template is instantiated for backward-only transform.
inmf | input data in MultiFab or FabArray<BaseFab<float>> |
void amrex::FFT::R2C< T, D, C >::forward | ( | RT const * | in, |
CT * | out | ||
) |
Forward transform.
This raw pointer version of forward requires setLocalDomain/getLocalDomain and setLocalSpectralDomain/getLocalSpectralDomain have been called already. Note that one is allowed to call this function multiple times after the set/get domain functions are called only once, unless the domain decomposition changes. In fact, that is the preferred way because it has better performance. All processes need to call this function even if their local size is zero. If the local size is zero, one can pass nullptrs.
|
inline |
Forward and then backward transform.
This function is available only when this class template is instantiated for transforms in both directions. It's more efficient than calling the forward function that stores the spectral data in a caller provided container followed by the backward function, because this can avoid parallel communication between the internal data and the caller's data container.
inmf | input data in MultiFab or FabArray<BaseFab<float>> |
outmf | output data in MultiFab or FabArray<BaseFab<float>> |
post_forward | a callable object for processing the post-forward data before the backward transform. Its interface is (int,int,int,GpuComplex<T>&) , where the integers are indices in the spectral space, and the reference to the complex number allows for the modification of the spectral data at that location. |
std::pair< std::array< int, AMREX_SPACEDIM >, std::array< int, AMREX_SPACEDIM > > amrex::FFT::R2C< T, D, C >::getLocalDomain | ( | ) | const |
Get local domain.
This function returns the domain decomposition chosen by AMReX. The first part of the pair is the local starting indices, and the second part is the local domain size.
This is needed, only if one uses the raw pointer interfaces, not the amrex::MulitFab interfaces. Only one of the functions, setLocalDomain and getLocalDomain, should be called.
std::pair< std::array< int, AMREX_SPACEDIM >, std::array< int, AMREX_SPACEDIM > > amrex::FFT::R2C< T, D, C >::getLocalSpectralDomain | ( | ) | const |
Get local spectral domain.
This function returns the domain decomposition chosen by AMReX for the complex data spectral domain. The returned pair contains the local starting indices and the local domain size.
This is needed, only if one uses the raw pointer interfaces, not the amrex::MulitFab interfaces. Only one of the functions, setLocalSpectralDomain and getLocalSpectralDomain, should be called. Note that one could use this function together with setLocalDomain. That is the user is allowed to choose their own real domain decomposition, while let AMReX choose the spectral data domain decomposition. Also note that the entire spectral domain has the size of (nx+1)/2 * ny * nz, if the real domain is nx * ny * nz.
std::pair< cMF *, IntVect > amrex::FFT::R2C< T, D, C >::getSpectralData | ( | ) | const |
Get the internal spectral data.
This function is not available when this class template is instantiated for backward-only transform. For performance reasons, the returned data array does not have the usual ordering of (x,y,z)
. The order is specified in the second part of the return value.
std::pair< typename R2C< T, D, C >::cMF *, IntVect > amrex::FFT::R2C< T, D, C >::getSpectralData | ( | ) | const |
std::pair< BoxArray, DistributionMapping > amrex::FFT::R2C< T, D, C >::getSpectralDataLayout | ( | ) | const |
Get BoxArray and DistributionMapping for spectral data.
The returned BoxArray and DistributionMapping can be used to build FabArray<BaseFab<GpuComplex<T>>> for spectral data. The returned BoxArray has the usual order of (x,y,z)
.
|
private |
|
private |
|
inlinestaticprivate |
|
inlinestaticprivate |
|
inlinestaticprivate |
|
staticprivate |
|
delete |
|
delete |
void amrex::FFT::R2C< T, D, C >::post_forward_doit_0 | ( | F const & | post_forward | ) |
void amrex::FFT::R2C< T, D, C >::post_forward_doit_1 | ( | F const & | post_forward | ) |
|
private |
T amrex::FFT::R2C< T, D, C >::scalingFactor | ( | ) | const |
Scaling factor. If the data goes through forward and then backward, the result multiplied by the scaling factor is equal to the original data.
void amrex::FFT::R2C< T, D, C >::setLocalDomain | ( | std::array< int, AMREX_SPACEDIM > const & | local_start, |
std::array< int, AMREX_SPACEDIM > const & | local_size | ||
) |
Set local domain.
This is needed, only if one uses the raw pointer interfaces, not the amrex::MulitFab interfaces. This may contain collective MPI calls. So all processes even if their local size is zero should call. This function informs AMReX the domain decomposition chosen by the user for the forward domain (i.e., real data domain). There is no constraint on the domain decomposition strategy. One can do 1D, 2D or 3D domain decomposition. Alternatively, one could also let AMReX choose for you by calling getLocalDomain. The latter could potentially reduce data communication.
Again, this is needed, only if one uses the raw pointer interfaces. Only one of the functions, setLocalDomain and getLocalDomain, should be called.
This should only be called once unless the domain decomposition changes.
local_start starting indices of the local domain local_size size of the local domain
void amrex::FFT::R2C< T, D, C >::setLocalSpectralDomain | ( | std::array< int, AMREX_SPACEDIM > const & | local_start, |
std::array< int, AMREX_SPACEDIM > const & | local_size | ||
) |
Set local spectral domain.
This is needed, only if one uses the raw pointer interfaces, not the amrex::MulitFab interfaces. This may contain collective MPI calls. So all processes even if their local size is zero should call. This function informs AMReX the domain decomposition chosen by the user for the complex data domain. There is no constraint on the domain decomposition strategy. One can do 1D, 2D or 3D domain decomposition. Alternatively, one could also let AMReX choose for you by calling getLocalSpectralDomain. The latter could potentially reduce data communication.
Again, this is needed, only if one uses the raw pointer interfaces. Only one of the functions, setLocalSpectralDomain and getLocalSpectralDomain, should be called. Note that one could use this function together with getLocalDomain. That is the user is allowed to choose their own spectral domain decomposition, while let AMReX choose the real data domain decomposition. Also note that the entire spectral domain has the size of (nx+1)/2 * ny * nz, if the real domain is nx * ny * nz.
This should only be called once unless the domain decomposition changes.
local_start starting indices of the local domain local_size size of the local domain
|
friend |
|
friend |
|
friend |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
mutableprivate |
|
mutableprivate |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |