QMCPACK
FftContainer Class Reference
+ Collaboration diagram for FftContainer:

Public Member Functions

 FftContainer (int nx, int ny, int nz)
 
double getRsValue (int x, int y, int z, int cplex) const
 
void setRsValue (int x, int y, int z, int cplex, double value)
 
double getKsValue (int x, int y, int z, int cplex) const
 
void setKsValue (int x, int y, int z, int cplex, double value)
 
int getNx () const
 
int getNy () const
 
int getNz () const
 
int getIndex (int x, int y, int z) const
 
int getQboxIndex (int x, int y, int z) const
 
 ~FftContainer ()
 
void executeFFT ()
 
void fixRsNorm (double factor)
 
void fixKsNorm (double factor)
 
double getL2NormRS () const
 
double getL2NormKS () const
 

Public Attributes

int fullSize
 
fftw_complex * rspace
 
fftw_complex * kspace
 

Private Attributes

fftw_plan plan_
 
int nx_
 
int ny_
 
int nz_
 

Detailed Description

Definition at line 16 of file FftContainer.h.

Constructor & Destructor Documentation

◆ FftContainer()

FftContainer ( int  nx,
int  ny,
int  nz 
)

Definition at line 65 of file FftContainer.h.

References fullSize, kspace, nx_, ny_, nz_, plan_, and rspace.

66 {
67  nx_ = nx;
68  ny_ = ny;
69  nz_ = nz;
70  fullSize = nx_ * ny_ * nz_;
71  rspace = (fftw_complex*)fftw_malloc(fullSize * sizeof(fftw_complex));
72  kspace = (fftw_complex*)fftw_malloc(fullSize * sizeof(fftw_complex));
73  plan_ = fftw_plan_dft_3d(nx_, ny_, nz_, rspace, kspace, -1, FFTW_ESTIMATE);
74 }
fftw_plan plan_
Definition: FftContainer.h:19
fftw_complex * kspace
Definition: FftContainer.h:25
fftw_complex * rspace
Definition: FftContainer.h:24

◆ ~FftContainer()

Definition at line 76 of file FftContainer.h.

References kspace, plan_, and rspace.

77 {
78  fftw_destroy_plan(plan_);
79  fftw_free(rspace);
80  fftw_free(kspace);
81 }
fftw_plan plan_
Definition: FftContainer.h:19
fftw_complex * kspace
Definition: FftContainer.h:25
fftw_complex * rspace
Definition: FftContainer.h:24

Member Function Documentation

◆ executeFFT()

void executeFFT ( )

Definition at line 83 of file FftContainer.h.

References plan_.

Referenced by EshdfFile::readInEigFcn().

83 { fftw_execute(plan_); }
fftw_plan plan_
Definition: FftContainer.h:19

◆ fixKsNorm()

void fixKsNorm ( double  factor)
inline

Definition at line 51 of file FftContainer.h.

References fullSize, and kspace.

Referenced by EshdfFile::readInEigFcn().

52  {
53  for (int i = 0; i < fullSize; i++)
54  {
55  kspace[i][0] *= factor;
56  kspace[i][1] *= factor;
57  }
58  }
fftw_complex * kspace
Definition: FftContainer.h:25

◆ fixRsNorm()

void fixRsNorm ( double  factor)
inline

Definition at line 43 of file FftContainer.h.

References fullSize, and rspace.

44  {
45  for (int i = 0; i < fullSize; i++)
46  {
47  rspace[i][0] *= factor;
48  rspace[i][1] *= factor;
49  }
50  }
fftw_complex * rspace
Definition: FftContainer.h:24

◆ getIndex()

int getIndex ( int  x,
int  y,
int  z 
) const
inline

Definition at line 38 of file FftContainer.h.

References ny_, and nz_.

Referenced by getKsValue(), getRsValue(), setKsValue(), and setRsValue().

38 { return z + y * nz_ + x * ny_ * nz_; }

◆ getKsValue()

double getKsValue ( int  x,
int  y,
int  z,
int  cplex 
) const

Definition at line 112 of file FftContainer.h.

References getIndex(), and kspace.

112 { return kspace[getIndex(x, y, z)][cplex]; }
fftw_complex * kspace
Definition: FftContainer.h:25
int getIndex(int x, int y, int z) const
Definition: FftContainer.h:38

◆ getL2NormKS()

double getL2NormKS ( ) const

Definition at line 95 of file FftContainer.h.

References fullSize, and kspace.

96 {
97  double l2norm = 0.0;
98  for (int i = 0; i < fullSize; i++)
99  {
100  l2norm += (kspace[i][0] * kspace[i][0] + kspace[i][1] * kspace[i][1]);
101  }
102  return l2norm;
103 }
fftw_complex * kspace
Definition: FftContainer.h:25

◆ getL2NormRS()

double getL2NormRS ( ) const

Definition at line 85 of file FftContainer.h.

References fullSize, and rspace.

86 {
87  double l2norm = 0.0;
88  for (int i = 0; i < fullSize; i++)
89  {
90  l2norm += (rspace[i][0] * rspace[i][0] + rspace[i][1] * rspace[i][1]);
91  }
92  return l2norm;
93 }
fftw_complex * rspace
Definition: FftContainer.h:24

◆ getNx()

int getNx ( ) const
inline

Definition at line 34 of file FftContainer.h.

References nx_.

Referenced by EshdfFile::readInEigFcn().

34 { return nx_; }

◆ getNy()

int getNy ( ) const
inline

Definition at line 35 of file FftContainer.h.

References ny_.

Referenced by EshdfFile::readInEigFcn().

35 { return ny_; }

◆ getNz()

int getNz ( ) const
inline

Definition at line 36 of file FftContainer.h.

References nz_.

Referenced by EshdfFile::readInEigFcn().

36 { return nz_; }

◆ getQboxIndex()

int getQboxIndex ( int  x,
int  y,
int  z 
) const
inline

Definition at line 39 of file FftContainer.h.

References nx_, and ny_.

Referenced by EshdfFile::readInEigFcn().

39 { return x + y * nx_ + z * ny_ * nx_; }

◆ getRsValue()

double getRsValue ( int  x,
int  y,
int  z,
int  cplex 
) const

Definition at line 105 of file FftContainer.h.

References getIndex(), and rspace.

105 { return rspace[getIndex(x, y, z)][cplex]; }
int getIndex(int x, int y, int z) const
Definition: FftContainer.h:38
fftw_complex * rspace
Definition: FftContainer.h:24

◆ setKsValue()

void setKsValue ( int  x,
int  y,
int  z,
int  cplex,
double  value 
)

Definition at line 114 of file FftContainer.h.

References getIndex(), and kspace.

115 {
116  kspace[getIndex(x, y, z)][cplex] = value;
117 }
fftw_complex * kspace
Definition: FftContainer.h:25
int getIndex(int x, int y, int z) const
Definition: FftContainer.h:38

◆ setRsValue()

void setRsValue ( int  x,
int  y,
int  z,
int  cplex,
double  value 
)

Definition at line 107 of file FftContainer.h.

References getIndex(), and rspace.

108 {
109  rspace[getIndex(x, y, z)][cplex] = value;
110 }
int getIndex(int x, int y, int z) const
Definition: FftContainer.h:38
fftw_complex * rspace
Definition: FftContainer.h:24

Member Data Documentation

◆ fullSize

◆ kspace

fftw_complex* kspace

◆ nx_

int nx_
private

Definition at line 20 of file FftContainer.h.

Referenced by FftContainer(), getNx(), and getQboxIndex().

◆ ny_

int ny_
private

Definition at line 20 of file FftContainer.h.

Referenced by FftContainer(), getIndex(), getNy(), and getQboxIndex().

◆ nz_

int nz_
private

Definition at line 20 of file FftContainer.h.

Referenced by FftContainer(), getIndex(), and getNz().

◆ plan_

fftw_plan plan_
private

Definition at line 19 of file FftContainer.h.

Referenced by executeFFT(), FftContainer(), and ~FftContainer().

◆ rspace

fftw_complex* rspace

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