Block-Structured AMR Software Framework
amrex::FFT::LocalR2C< T, D, M > Class Template Reference

Local Discrete Fourier Transform. More...

#include <AMReX_FFT_LocalR2C.H>

Public Member Functions

 LocalR2C (IntVectND< M > const &fft_size, T *p_fwd=nullptr, GpuComplex< T > *p_bwd=nullptr, bool cache_plan=true)
 Constructor. More...
 
 ~LocalR2C ()
 
 LocalR2C ()=default
 
 LocalR2C (LocalR2C &&) noexcept
 
LocalR2Coperator= (LocalR2C &&) noexcept
 
 LocalR2C (LocalR2C const &)=delete
 
LocalR2Coperator= (LocalR2C const &)=delete
 
template<Direction DIR = D, std::enable_if_t< DIR==Direction::forward||DIR==Direction::both, int > = 0>
void forward (T const *indata, GpuComplex< T > *outdata)
 Forward transform. More...
 
void clear ()
 
template<Direction DIR = D, std::enable_if_t< DIR==Direction::backward||DIR==Direction::both, int > = 0>
void backward (GpuComplex< T > const *indata, T *outdata)
 Backward transform. More...
 
scalingFactor () const
 
IntVectND< M > const & spectralSize () const
 Spectral domain size. More...
 

Private Attributes

Plan< T > m_fft_fwd
 
Plan< T > m_fft_bwd
 
T * m_p_fwd = nullptr
 
GpuComplex< T > * m_p_bwd = nullptr
 
IntVectND< M > m_real_size
 
IntVectND< M > m_spectral_size
 
bool m_cache_plan = false
 

Detailed Description

template<typename T, FFT::Direction D = FFT::Direction::both, int M = AMREX_SPACEDIM>
class amrex::FFT::LocalR2C< T, D, M >

Local 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.

For more details, we refer the users to https://amrex-codes.github.io/amrex/docs_html/FFT_Chapter.html.

Constructor & Destructor Documentation

◆ LocalR2C() [1/4]

template<typename T , FFT::Direction D, int M>
amrex::FFT::LocalR2C< T, D, M >::LocalR2C ( IntVectND< M > const &  fft_size,
T *  p_fwd = nullptr,
GpuComplex< T > *  p_bwd = nullptr,
bool  cache_plan = true 
)
explicit

Constructor.

Given the diverse interfaces of FFT libraries we use, this constructo has a number of optional arguments.

The user can provide the data pointers to the constructor. They are only needed by FFTW because its plan creation requires the input and output arrays. If they are null, we will delay the plan creation for FFTW until the forward or backward function is called.

The cache_plan option is only used when we use cufft, rocfft and onemkl, but not FFTW.

Parameters
fft_sizeThe forward domain size (i.e., the domain of the real data)
p_fwdForward domain data pointer (optional)
p_bwdBackward domain data pointer (optional)
cache_planTry to cache the plan or not (optionl)

◆ ~LocalR2C()

template<typename T , FFT::Direction D, int M>
amrex::FFT::LocalR2C< T, D, M >::~LocalR2C

◆ LocalR2C() [2/4]

template<typename T , FFT::Direction D = FFT::Direction::both, int M = AMREX_SPACEDIM>
amrex::FFT::LocalR2C< T, D, M >::LocalR2C ( )
default

◆ LocalR2C() [3/4]

template<typename T , FFT::Direction D, int M>
amrex::FFT::LocalR2C< T, D, M >::LocalR2C ( LocalR2C< T, D, M > &&  rhs)
noexcept

◆ LocalR2C() [4/4]

template<typename T , FFT::Direction D = FFT::Direction::both, int M = AMREX_SPACEDIM>
amrex::FFT::LocalR2C< T, D, M >::LocalR2C ( LocalR2C< T, D, M > const &  )
delete

Member Function Documentation

◆ backward()

template<typename T , FFT::Direction D, int M>
template<Direction DIR, std::enable_if_t< DIR==Direction::backward||DIR==Direction::both, int > >
void amrex::FFT::LocalR2C< T, D, M >::backward ( GpuComplex< T > const *  indata,
T *  outdata 
)

Backward transform.

This function is not available when this class template is instantiated for forward-only transform. For GPUs, this function is synchronous on the host.

Parameters
indatainput data
outdataoutput data

◆ clear()

template<typename T , FFT::Direction D, int M>
void amrex::FFT::LocalR2C< T, D, M >::clear

◆ forward()

template<typename T , FFT::Direction D, int M>
template<Direction DIR, std::enable_if_t< DIR==Direction::forward||DIR==Direction::both, int > >
void amrex::FFT::LocalR2C< T, D, M >::forward ( T const *  indata,
GpuComplex< T > *  outdata 
)

Forward transform.

This function is not available when this class template is instantiated for backward-only transform. For GPUs, this function is synchronous on the host.

Parameters
indatainput data
outdataoutput data

◆ operator=() [1/2]

template<typename T , FFT::Direction D, int M>
LocalR2C< T, D, M > & amrex::FFT::LocalR2C< T, D, M >::operator= ( LocalR2C< T, D, M > &&  rhs)
noexcept

◆ operator=() [2/2]

template<typename T , FFT::Direction D = FFT::Direction::both, int M = AMREX_SPACEDIM>
LocalR2C& amrex::FFT::LocalR2C< T, D, M >::operator= ( LocalR2C< T, D, M > const &  )
delete

◆ scalingFactor()

template<typename T , FFT::Direction D, int M>
T amrex::FFT::LocalR2C< T, D, M >::scalingFactor

Scaling factor. If the data goes through forward and then backward, the result multiplied by the scaling factor is equal to the original data.

◆ spectralSize()

template<typename T , FFT::Direction D = FFT::Direction::both, int M = AMREX_SPACEDIM>
IntVectND<M> const& amrex::FFT::LocalR2C< T, D, M >::spectralSize ( ) const
inline

Spectral domain size.

Member Data Documentation

◆ m_cache_plan

template<typename T , FFT::Direction D = FFT::Direction::both, int M = AMREX_SPACEDIM>
bool amrex::FFT::LocalR2C< T, D, M >::m_cache_plan = false
private

◆ m_fft_bwd

template<typename T , FFT::Direction D = FFT::Direction::both, int M = AMREX_SPACEDIM>
Plan<T> amrex::FFT::LocalR2C< T, D, M >::m_fft_bwd
private

◆ m_fft_fwd

template<typename T , FFT::Direction D = FFT::Direction::both, int M = AMREX_SPACEDIM>
Plan<T> amrex::FFT::LocalR2C< T, D, M >::m_fft_fwd
private

◆ m_p_bwd

template<typename T , FFT::Direction D = FFT::Direction::both, int M = AMREX_SPACEDIM>
GpuComplex<T>* amrex::FFT::LocalR2C< T, D, M >::m_p_bwd = nullptr
private

◆ m_p_fwd

template<typename T , FFT::Direction D = FFT::Direction::both, int M = AMREX_SPACEDIM>
T* amrex::FFT::LocalR2C< T, D, M >::m_p_fwd = nullptr
private

◆ m_real_size

template<typename T , FFT::Direction D = FFT::Direction::both, int M = AMREX_SPACEDIM>
IntVectND<M> amrex::FFT::LocalR2C< T, D, M >::m_real_size
private

◆ m_spectral_size

template<typename T , FFT::Direction D = FFT::Direction::both, int M = AMREX_SPACEDIM>
IntVectND<M> amrex::FFT::LocalR2C< T, D, M >::m_spectral_size
private

The documentation for this class was generated from the following file: