JDFTx  1.7.0
GridInfo Class Reference

Simulation grid descriptor. More...

#include <GridInfo.h>

Public Types

enum  LatticeType {
  Manual , Triclinic , Monoclinic , Orthorhombic ,
  Tetragonal , Rhombohedral , Hexagonal , Cubic
}
 
enum  LatticeModification { Simple , BodyCentered , BaseCentered , FaceCentered }
 
enum  PlanType {
  PlanForward , PlanInverse , PlanForwardInPlace , PlanInverseInPlace ,
  PlanRtoC , PlanCtoR
}
 

Public Member Functions

void update ()
 Update the grid information on changing the lattice vectors.
 
void printLattice ()
 Print the lattice vectors.
 
void printReciprocalLattice ()
 Print the reciprocal lattice vectors.
 
void setLatticeVectors ()
 Set lattice vectors based on lattice type, modification and parameters above.
 
void initialize (bool skipHeader=false, const std::vector< SpaceGroupOp > sym=std::vector< SpaceGroupOp >(1, SpaceGroupOp()))
 
fftw_plan getPlan (PlanType planType, int nThreads) const
 
vector3< int > wrapGcoords (const vector3< int > iG) const
 < wrap negative G-indices to the positive side
 
int fullRindex (const vector3< int > iR) const
 < Index into the full real-space box:
 
int fullGindex (const vector3< int > iG) const
 < Index into the full reciprocal-space box:
 
int halfGindex (const vector3< int > iG) const
 < Index into the half-reduced reciprocal-space box:
 

Public Attributes

enum GridInfo::LatticeType latticeType
 
enum GridInfo::LatticeModification latticeModification
 
double a
 
double b
 
double c
 
double alpha
 
double beta
 
double gamma
 lattice specified by type, modification and standard side lengths, angles in degrees
 
vector3 lattScale
 Remember lattice scale specified at input (Note R always includes scale, once latt-scale has processed)
 
matrix3 R
 directly specified lattice vectors
 
double Gmax
 radius of wavefunction G-sphere, whole density sphere (double the radius) must be inscribable within the FFT box
 
double GmaxRho
 if non-zero, override the FFT box inscribable sphere radius
 
vector3< int > S
 sample points in each dimension (if 0, will be determined automatically based on Gmax)
 
double detR
 cell volume
 
matrix3 RT
 
matrix3 RTR
 
matrix3 invR
 
matrix3 invRT
 
matrix3 invRTR
 various combinations of lattice-vectors
 
matrix3 G
 
matrix3 GT
 
matrix3 GGT
 
matrix3 invGGT
 various combinations of reciporcal lattice-vectors
 
double dV
 volume per grid point
 
vector3 h [3]
 real space sample vectors
 
int nr
 position space grid count = S[0]*S[1]*S[2]
 
int nG
 reciprocal lattice count = S[0]*S[1]*(S[2]/2+1) [on account of using r2c and c2r ffts]
 
double dGradial
 recommended spacing of radial G functions
 
double GmaxSphere
 recommended maximum G-vector for radial functions for the wavefunction sphere
 
double GmaxGrid
 recommended maximum G-vector for radial functions for the density grid
 
int irStart
 
int irStop
 
int iGstart
 
int iGstop
 
cufftHandle planZ2Z
 CUFFT plan for all the complex transforms.
 
cufftHandle planD2Z
 CUFFT plan for R -> G.
 
cufftHandle planZ2D
 CUFFT plan for G -> R.
 

Static Public Attributes

static const double maxAllowedStrain
 maximum strain allowed when updating lattice vectors and using update()
 

Detailed Description

Simulation grid descriptor.

To setup a simulation grid, create a blank GridInfo object, set the public data members S and R, and call initialize().

This sets up all the auxilliary grid information, and a bunch of shared utilities such as a random number generator, plans for Fourier transforms etc.

Member Enumeration Documentation

◆ LatticeType

Enumerator
Manual 

Directly specify R.

◆ PlanType

Enumerator
PlanForward 

Forward complex transform.

PlanInverse 

Inverse complex transform.

PlanForwardInPlace 

Forward in-place complex transform.

PlanInverseInPlace 

Inverse in-place complex transform.

PlanRtoC 

Real to complex transform.

PlanCtoR 

Complex to real transform.

Member Function Documentation

◆ initialize()

void GridInfo::initialize ( bool  skipHeader = false,
const std::vector< SpaceGroupOp sym = std::vector< SpaceGroupOp >(1, SpaceGroupOp()) 
)

Initialize the dependent quantities below. If S is specified and is too small for the given Gmax, the call will abort. If skipHeader=true, the "initializing the grid" banner will be supressed (useful for auxiliary grids in calculation) If sym is specified, then auto-computed S will be ensured commensurate with those symmetries


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