QMCPACK
CrystalLattice< T, D > Struct Template Reference

a class that defines a supercell in D-dimensional Euclean space. More...

+ Inheritance diagram for CrystalLattice< T, D >:
+ Collaboration diagram for CrystalLattice< T, D >:

Public Types

enum  { DIM = D }
 enumeration for the dimension of the lattice More...
 
using Base = LRBreakupParameters< T, D >
 alias to the base class More...
 
using Scalar_t = T
 the type of scalar More...
 
using SingleParticlePos = TinyVector< T, D >
 the type of a D-dimensional position vector More...
 
using SingleParticleIndex = TinyVector< int, D >
 the type of a D-dimensional index vector More...
 
using Tensor_t = Tensor< T, D >
 the type of a D-dimensional Tensor More...
 

Public Member Functions

 CrystalLattice ()
 default constructor, assign a huge supercell More...
 
SingleParticlePos a (int i) const
 Provide interfaces familiar to fotran users. More...
 
SingleParticlePos b (int i) const
 Provide interfaces familiar to fotran users. More...
 
template<class T1 >
SingleParticlePos toUnit (const TinyVector< T1, D > &r) const
 Convert a cartesian vector to a unit vector. More...
 
template<class T1 >
SingleParticlePos toUnit_floor (const TinyVector< T1, D > &r) const
 
template<class T1 >
SingleParticlePos toCart (const TinyVector< T1, D > &c) const
 Convert a unit vector to a cartesian vector. More...
 
bool isValid (const TinyVector< T, D > &u) const
 return true if all the open direction of reduced coordinates u are in the range [0,1) More...
 
bool outOfBound (const TinyVector< T, D > &u) const
 return true if any direction of reduced coordinates u goes larger than 0.5 More...
 
void applyMinimumImage (TinyVector< T, D > &c) const
 
Dot (const SingleParticlePos &ra, const SingleParticlePos &rb) const
 evaluate the cartesian distance More...
 
SingleParticlePos k_cart (const SingleParticlePos &kin) const
 conversion of a reciprocal-vector More...
 
SingleParticlePos k_unit (const SingleParticlePos &kin) const
 conversion of a caresian reciprocal-vector to unit k-vector More...
 
ksq (const SingleParticlePos &kin) const
 evaluate $k^2$ More...
 
template<typename T1 >
CrystalLattice< T, D > & operator= (const CrystalLattice< T1, D > &rhs)
 assignment operator More...
 
CrystalLattice< T, D > & operator*= (T sc)
 scale the lattice vectors by sc. More...
 
template<class TT >
void set (const Tensor< TT, D > &lat)
 set the lattice vector from the command-line options More...
 
void reset ()
 Evaluate the reciprocal vectors, volume and metric tensor. More...
 
void print (std::ostream &, int level=2) const
 Print out CrystalLattice Data. More...
 

Public Attributes

bool DiagonalOnly
 true, if off-diagonal elements are zero so that other classes can take advantage of this More...
 
int SuperCellEnum
 supercell enumeration More...
 
TinyVector< int, D > BoxBConds
 The boundary condition in each direction. More...
 
VacuumScale
 The scale factor for adding vacuum. More...
 
SingleParticlePos ABC
 
bool explicitly_defined
 true, the lattice is defined by the input instead of an artificial default More...
 
Scalar_t Volume
 Physical properties of a supercell. More...
 
Scalar_t WignerSeitzRadius
 Wigner-Seitz cell radius. More...
 
Scalar_t SimulationCellRadius
 simulation cell radii More...
 
Scalar_t CellRadiusSq
 SimulationCellRadius*SimulationCellRadius. More...
 
Scalar_t WignerSeitzRadius_G
 Wigner-Seitz cell radius in reciprocal space. More...
 
Tensor_t R
 Real-space unit vectors. R(i,j) i=vector and j=x,y,z. More...
 
Tensor_t G
 Reciprocal unit vectors. G(j,i) i=vector and j=x,y,z. More...
 
Tensor_t Gt
 Transpose of reciprocal unit vectors: More...
 
Tensor_t M
 Metric tensor. More...
 
Tensor_t Mg
 Metric tensor for G vectors. More...
 
SingleParticlePos Length
 Length[idim] length of the idim-th lattice vector. More...
 
SingleParticlePos OneOverLength
 OneOverLength[idim] 1/length of the idim-th lattice vector. More...
 
SingleParticlePos Center
 Center of the cell sum(Rv[i])/2. More...
 
TinyVector< SingleParticlePos, D > Rv
 Real-space unit vectors. More...
 
TinyVector< SingleParticlePos, D > Gv
 Reciprocal unit vectors. More...
 

Detailed Description

template<class T, unsigned D>
struct qmcplusplus::CrystalLattice< T, D >

a class that defines a supercell in D-dimensional Euclean space.

CrystalLattice handles the physical properties of a supercell, such as lattice vectors, reciprocal vectors and metric tensors and provides interfaces to access the lattice properties and convert units of position vectors or a single-particle position from Cartesian to Lattice Unit vice versa.

Definition at line 55 of file CrystalLattice.h.

Member Typedef Documentation

◆ Base

using Base = LRBreakupParameters<T, D>

alias to the base class

Definition at line 58 of file CrystalLattice.h.

◆ Scalar_t

using Scalar_t = T

the type of scalar

Definition at line 67 of file CrystalLattice.h.

◆ SingleParticleIndex

the type of a D-dimensional index vector

Definition at line 71 of file CrystalLattice.h.

◆ SingleParticlePos

the type of a D-dimensional position vector

Definition at line 69 of file CrystalLattice.h.

◆ Tensor_t

using Tensor_t = Tensor<T, D>

the type of a D-dimensional Tensor

Definition at line 73 of file CrystalLattice.h.

Member Enumeration Documentation

◆ anonymous enum

anonymous enum

enumeration for the dimension of the lattice

Enumerator
DIM 

Definition at line 61 of file CrystalLattice.h.

62  {
63  DIM = D
64  };

Constructor & Destructor Documentation

◆ CrystalLattice()

default constructor, assign a huge supercell

Default constructor. Initialized to a 1x1x1 cubic supercell.

Definition at line 29 of file CrystalLattice.cpp.

30 {
31  explicitly_defined = false;
32  BoxBConds = 0;
33  VacuumScale = 1.0;
34  R.diagonal(1e10);
35  G = R;
36  M = R;
37  Volume = 1;
38  reset();
39 }
void reset()
Evaluate the reciprocal vectors, volume and metric tensor.
T VacuumScale
The scale factor for adding vacuum.
Tensor_t G
Reciprocal unit vectors. G(j,i) i=vector and j=x,y,z.
Tensor_t M
Metric tensor.
Scalar_t Volume
Physical properties of a supercell.
TinyVector< int, D > BoxBConds
The boundary condition in each direction.
void diagonal(const T &rhs)
Definition: Tensor.h:205
bool explicitly_defined
true, the lattice is defined by the input instead of an artificial default
Tensor_t R
Real-space unit vectors. R(i,j) i=vector and j=x,y,z.

Member Function Documentation

◆ a()

SingleParticlePos a ( int  i) const
inline

Provide interfaces familiar to fotran users.

Parameters
ithe index of the directional vector, $i\in [0,D)$
Returns
The lattice vector of the ith direction

Definition at line 137 of file CrystalLattice.h.

Referenced by DTD_BConds< T, 3, PPPG >::DTD_BConds(), DTD_BConds< T, 3, PPNG >::DTD_BConds(), DTD_BConds< T, 3, PPPG+SOA_OFFSET >::DTD_BConds(), and DTD_BConds< T, 3, PPNG+SOA_OFFSET >::DTD_BConds().

137 { return Rv[i]; }
TinyVector< SingleParticlePos, D > Rv
Real-space unit vectors.

◆ applyMinimumImage()

void applyMinimumImage ( TinyVector< T, D > &  c) const
inline

Definition at line 194 of file CrystalLattice.h.

195  {
196  if (SuperCellEnum)
197  {
198  TinyVector<T, D> u = dot(c, G);
199  for (int i = 0; i < D; ++i)
200  u[i] = u[i] - round(u[i]);
201  c = dot(u, R);
202  }
203  }
int SuperCellEnum
supercell enumeration
Tensor_t G
Reciprocal unit vectors. G(j,i) i=vector and j=x,y,z.
Tensor< typename BinaryReturn< T1, T2, OpMultiply >::Type_t, D > dot(const AntiSymTensor< T1, D > &lhs, const AntiSymTensor< T2, D > &rhs)
Tensor_t R
Real-space unit vectors. R(i,j) i=vector and j=x,y,z.

◆ b()

SingleParticlePos b ( int  i) const
inline

Provide interfaces familiar to fotran users.

Parameters
ithe index of the directional vector, $i\in [0,D)$
Returns
The reciprocal vector of the ith direction

Definition at line 143 of file CrystalLattice.h.

143 { return Gv[i]; }
TinyVector< SingleParticlePos, D > Gv
Reciprocal unit vectors.

◆ Dot()

T Dot ( const SingleParticlePos ra,
const SingleParticlePos rb 
) const
inline

evaluate the cartesian distance

Parameters
raa vector in the supercell unit
rba vector in the supercell unit
Returns
Cartesian distance with two vectors in SC unit
Note
The distance between two cartesian vectors are handled by dot function defined in OhmmsPETE/TinyVector.h

Definition at line 213 of file CrystalLattice.h.

213 { return dot(ra, dot(M, rb)); }
Tensor_t M
Metric tensor.
Tensor< typename BinaryReturn< T1, T2, OpMultiply >::Type_t, D > dot(const AntiSymTensor< T1, D > &lhs, const AntiSymTensor< T2, D > &rhs)

◆ isValid()

bool isValid ( const TinyVector< T, D > &  u) const
inline

return true if all the open direction of reduced coordinates u are in the range [0,1)

Definition at line 177 of file CrystalLattice.h.

Referenced by ParticleSet::makeMoveAllParticles(), ParticleSet::makeMoveAllParticlesWithDrift(), ParticleSet::makeMoveAndCheck(), and qmcplusplus::TEST_CASE().

178  {
179  bool inside = true;
180  for (int dim = 0; dim < D; dim++)
181  inside &= (BoxBConds[dim] || (u[dim] >= T(0) && u[dim] < T(1)));
182  return inside;
183  }
TinyVector< int, D > BoxBConds
The boundary condition in each direction.

◆ k_cart()

SingleParticlePos k_cart ( const SingleParticlePos kin) const
inline

conversion of a reciprocal-vector

Parameters
kinan input reciprocal vector in the Reciprocal-vector unit
Returns
k(reciprocal vector) in cartesian unit

Definition at line 219 of file CrystalLattice.h.

Referenced by EinsplineSetBuilder::AnalyzeTwists2(), BsplineReader::check_twists(), Gvectors< ST, LT >::Gvectors(), BsplineReader::initialize_spo2band(), MomentumDistribution::MomentumDistribution(), MomentumEstimator::putSpecial(), PWBasis::readbasis(), and PWBasis::trimforecut().

219 { return TWOPI * dot(G, kin); }
Tensor_t G
Reciprocal unit vectors. G(j,i) i=vector and j=x,y,z.
Tensor< typename BinaryReturn< T1, T2, OpMultiply >::Type_t, D > dot(const AntiSymTensor< T1, D > &lhs, const AntiSymTensor< T2, D > &rhs)

◆ k_unit()

SingleParticlePos k_unit ( const SingleParticlePos kin) const
inline

conversion of a caresian reciprocal-vector to unit k-vector

Parameters
kinan input reciprocal vector in cartesian form
Returns
k(reciprocal vector) as unit vector

Definition at line 225 of file CrystalLattice.h.

225 { return dot(R, kin) / TWOPI; }
Tensor< typename BinaryReturn< T1, T2, OpMultiply >::Type_t, D > dot(const AntiSymTensor< T1, D > &lhs, const AntiSymTensor< T2, D > &rhs)
Tensor_t R
Real-space unit vectors. R(i,j) i=vector and j=x,y,z.

◆ ksq()

T ksq ( const SingleParticlePos kin) const
inline

evaluate $k^2$

Parameters
kinan input reciprocal vector in reciprocal-vector unit
Returns
$k_{in}^2$

Definition at line 232 of file CrystalLattice.h.

232 { return dot(kin, dot(Mg, kin)); }
Tensor< typename BinaryReturn< T1, T2, OpMultiply >::Type_t, D > dot(const AntiSymTensor< T1, D > &lhs, const AntiSymTensor< T2, D > &rhs)
Tensor_t Mg
Metric tensor for G vectors.

◆ operator*=()

CrystalLattice< T, D > & operator*= ( sc)

scale the lattice vectors by sc.

Rescale this supercell by a scalar.

All the internal data are reset.

Parameters
scthe scaling value
Returns
a new CrystalLattice
Parameters
scA scaling factor.

Definition at line 95 of file CrystalLattice.cpp.

96 {
97  R *= sc;
98  reset();
99  return *this;
100 }
void reset()
Evaluate the reciprocal vectors, volume and metric tensor.
Tensor_t R
Real-space unit vectors. R(i,j) i=vector and j=x,y,z.

◆ operator=()

CrystalLattice<T, D>& operator= ( const CrystalLattice< T1, D > &  rhs)
inline

assignment operator

Definition at line 236 of file CrystalLattice.h.

237  {
238  Base::LR_dim_cutoff = rhs.LR_dim_cutoff;
239  Base::LR_kc = rhs.LR_kc;
240  Base::LR_rc = rhs.LR_rc;
241  Base::num_ewald_grid_points = rhs.num_ewald_grid_points;
242  Base::ndim = rhs.ndim;
243 
244  explicitly_defined = rhs.explicitly_defined;
245  BoxBConds = rhs.BoxBConds;
246  VacuumScale = rhs.VacuumScale;
247  R = rhs.R;
248  reset();
249  return *this;
250  }
void reset()
Evaluate the reciprocal vectors, volume and metric tensor.
T VacuumScale
The scale factor for adding vacuum.
TinyVector< int, D > BoxBConds
The boundary condition in each direction.
bool explicitly_defined
true, the lattice is defined by the input instead of an artificial default
Tensor_t R
Real-space unit vectors. R(i,j) i=vector and j=x,y,z.

◆ outOfBound()

bool outOfBound ( const TinyVector< T, D > &  u) const
inline

return true if any direction of reduced coordinates u goes larger than 0.5

Definition at line 186 of file CrystalLattice.h.

Referenced by ParticleSet::makeMoveAllParticles(), ParticleSet::makeMoveAllParticlesWithDrift(), ParticleSet::makeMoveAndCheck(), and qmcplusplus::TEST_CASE().

187  {
188  for (int i = 0; i < D; ++i)
189  if (std::abs(u[i]) > 0.5)
190  return true;
191  return false;
192  }
MakeReturn< UnaryNode< FnFabs, typename CreateLeaf< Vector< T1, C1 > >::Leaf_t > >::Expression_t abs(const Vector< T1, C1 > &l)

◆ print()

void print ( std::ostream &  os,
int  level = 2 
) const

Print out CrystalLattice Data.

Definition at line 103 of file CrystalLattice.cpp.

Referenced by SimulationCell::resetLRBox().

104 {
105  /*\note level == 0: print only the lattice vectors
106  * level == 1: lattice vectors, boundary conditions, grid
107  * level == 2: + all the internal values
108  */
109  std::string unit_name = "bohr";
110 
111  std::string lattice_name = " Lattice (" + unit_name + "):";
112  std::string pad(lattice_name.length(), ' ');
113  os << lattice_name;
114  for (int i = 0; i < D; ++i)
115  {
116  if (i > 0)
117  {
118  os << pad;
119  }
120  os << Rv[i] << std::endl;
121  }
122  if (level > 0)
123  {
124  os << std::endl;
125  os << " Boundary Conditions: ";
126  for (int i = 0; i < D; ++i)
127  {
128  if (BoxBConds[i])
129  os << " p ";
130  else
131  os << " n ";
132  }
133  os << std::endl;
134  if (VacuumScale != 1.0)
135  os << " Vacuum scale: " << VacuumScale << std::endl;
136  }
137  if (level > 1)
138  {
139  os << std::endl;
140  os << " Volume (bohr^3) = " << Volume << std::endl;
141  os << std::endl;
142  os << " Reciprocal vectors without 2*pi.\n";
143  for (int i = 0; i < D; ++i)
144  os << " g_" << i + 1 << " = " << Gv[i] << std::endl;
145  os << std::endl;
146  os << " Metric tensor in real-space.\n";
147  for (int i = 0; i < D; ++i)
148  {
149  os << " h_" << i + 1 << " = ";
150  for (int j = 0; j < D; ++j)
151  {
152  os << M(i, j) << " ";
153  }
154  os << std::endl;
155  }
156  os << std::endl;
157  os << " Metric tensor in g-space.\n";
158  for (int i = 0; i < D; ++i)
159  {
160  os << " h_" << i + 1 << " = ";
161  for (int j = 0; j < D; ++j)
162  {
163  os << Mg(i, j) << " ";
164  }
165  os << std::endl;
166  }
167  }
168 }
T VacuumScale
The scale factor for adding vacuum.
Tensor_t M
Metric tensor.
TinyVector< SingleParticlePos, D > Rv
Real-space unit vectors.
Scalar_t Volume
Physical properties of a supercell.
TinyVector< int, D > BoxBConds
The boundary condition in each direction.
TinyVector< SingleParticlePos, D > Gv
Reciprocal unit vectors.
Tensor_t Mg
Metric tensor for G vectors.

◆ reset()

void reset ( )

Evaluate the reciprocal vectors, volume and metric tensor.

Definition at line 51 of file CrystalLattice.cpp.

Referenced by OneBodyDensityMatrices::OneBodyDensityMatrices(), CrystalLattice< ST, 3 >::operator=(), SimulationCell::resetLRBox(), and qmcplusplus::TEST_CASE().

52 {
53  G = inverse(R); //G = transpose(Inverse(R));
54  Gt = transpose(G);
55  Volume = std::abs(det(R));
56  //M = dot(transpose(R),R);
57  M = dot(R, transpose(R));
58  T t = TWOPI * TWOPI;
59  Mg = t * dot(transpose(G), G);
60  for (int i = 0; i < D; ++i)
61  for (int j = 0; j < D; ++j)
62  Rv[i][j] = R(i, j);
63  for (int i = 0; i < D; ++i)
64  for (int j = 0; j < D; ++j)
65  Gv[i][j] = G(j, i);
66  for (int i = 0; i < D; ++i)
67  {
68  Length[i] = std::sqrt(dot(Rv[i], Rv[i]));
69  OneOverLength[i] = 1.0 / Length[i];
70  }
71  Center = 0.0;
72  for (int i = 0; i < D; ++i)
73  Center += Rv[i];
74  Center *= .5;
75  //analysis of a lattice using LatticeAnalyzer
76  LatticeAnalyzer<T, D> ldesc;
77  SuperCellEnum = ldesc(BoxBConds);
78  DiagonalOnly = ldesc.isDiagonalOnly(R);
79  ABC = ldesc.calcSolidAngles(Rv, OneOverLength);
80  WignerSeitzRadius = ldesc.calcWignerSeitzRadius(Rv);
81  WignerSeitzRadius_G = ldesc.calcWignerSeitzRadius(Gv);
82  SimulationCellRadius = ldesc.calcSimulationCellRadius(Rv);
83  // set equal WignerSeitzRadius and SimulationCellRadius when they are very close.
85  WignerSeitzRadius - SimulationCellRadius <= WignerSeitzRadius * std::numeric_limits<float>::epsilon() * 2)
88 }
int SuperCellEnum
supercell enumeration
SingleParticlePos Center
Center of the cell sum(Rv[i])/2.
Tensor< T, D >::Type_t det(const Tensor< T, D > &a)
Definition: TensorOps.h:838
Scalar_t WignerSeitzRadius_G
Wigner-Seitz cell radius in reciprocal space.
MakeReturn< UnaryNode< FnFabs, typename CreateLeaf< Vector< T1, C1 > >::Leaf_t > >::Expression_t abs(const Vector< T1, C1 > &l)
Tensor_t G
Reciprocal unit vectors. G(j,i) i=vector and j=x,y,z.
Tensor_t M
Metric tensor.
TinyVector< SingleParticlePos, D > Rv
Real-space unit vectors.
SingleParticlePos Length
Length[idim] length of the idim-th lattice vector.
Tensor_t Gt
Transpose of reciprocal unit vectors:
Scalar_t Volume
Physical properties of a supercell.
TinyVector< int, D > BoxBConds
The boundary condition in each direction.
SingleParticlePos OneOverLength
OneOverLength[idim] 1/length of the idim-th lattice vector.
AntiSymTensor< T, D > transpose(const AntiSymTensor< T, D > &rhs)
MakeReturn< UnaryNode< FnSqrt, typename CreateLeaf< Vector< T1, C1 > >::Leaf_t > >::Expression_t sqrt(const Vector< T1, C1 > &l)
Tensor< T, D > inverse(const Tensor< T, D > &a)
Definition: TensorOps.h:879
TinyVector< SingleParticlePos, D > Gv
Reciprocal unit vectors.
Tensor< typename BinaryReturn< T1, T2, OpMultiply >::Type_t, D > dot(const AntiSymTensor< T1, D > &lhs, const AntiSymTensor< T2, D > &rhs)
Scalar_t SimulationCellRadius
simulation cell radii
bool DiagonalOnly
true, if off-diagonal elements are zero so that other classes can take advantage of this ...
Scalar_t WignerSeitzRadius
Wigner-Seitz cell radius.
Scalar_t CellRadiusSq
SimulationCellRadius*SimulationCellRadius.
Tensor_t Mg
Metric tensor for G vectors.
Tensor_t R
Real-space unit vectors. R(i,j) i=vector and j=x,y,z.

◆ set()

void set ( const Tensor< TT, D > &  lat)

set the lattice vector from the command-line options

Parameters
lata tensor representing a supercell

Definition at line 43 of file CrystalLattice.cpp.

Referenced by InitMolecularSystem::initWithVolume(), LatticeParser::put(), SpinDensity::put(), OrbitalImages::put(), EinsplineSetBuilder::ReadOrbitalInfo_ESHDF(), SpinDensityInput::readXML(), and EinsplineSetBuilder::set_metadata().

44 {
45  explicitly_defined = true;
46  R = lat;
47  reset();
48 }
void reset()
Evaluate the reciprocal vectors, volume and metric tensor.
bool explicitly_defined
true, the lattice is defined by the input instead of an artificial default
Tensor_t R
Real-space unit vectors. R(i,j) i=vector and j=x,y,z.

◆ toCart()

SingleParticlePos toCart ( const TinyVector< T1, D > &  c) const
inline

Convert a unit vector to a cartesian vector.

Boundary conditions are not applied.

Definition at line 171 of file CrystalLattice.h.

Referenced by MomentumDistribution::accumulate(), MomentumEstimator::evaluate(), OrbitalImages::evaluate(), OneBodyDensityMatrices::generateUniformGrid(), OneBodyDensityMatrices::generateUniformSamples(), and OneBodyDensityMatrices::normalizeBasis().

172  {
173  return dot(c, R);
174  }
Tensor< typename BinaryReturn< T1, T2, OpMultiply >::Type_t, D > dot(const AntiSymTensor< T1, D > &lhs, const AntiSymTensor< T2, D > &rhs)
Tensor_t R
Real-space unit vectors. R(i,j) i=vector and j=x,y,z.

◆ toUnit()

SingleParticlePos toUnit ( const TinyVector< T1, D > &  r) const
inline

Convert a cartesian vector to a unit vector.

Boundary conditions are not applied.

Definition at line 149 of file CrystalLattice.h.

Referenced by SpinDensityNew::accumulate(), MagnetizationDensity::computeBin(), SplineR2R< ST >::convertPos(), SpinDensity::evaluate(), ParticleSet::makeMoveAllParticles(), ParticleSet::makeMoveAllParticlesWithDrift(), ParticleSet::makeMoveAndCheck(), SpinDensity::test_evaluate(), and CrystalLattice< ST, 3 >::toUnit_floor().

150  {
151  return dot(r, G);
152  }
Tensor_t G
Reciprocal unit vectors. G(j,i) i=vector and j=x,y,z.
Tensor< typename BinaryReturn< T1, T2, OpMultiply >::Type_t, D > dot(const AntiSymTensor< T1, D > &lhs, const AntiSymTensor< T2, D > &rhs)

◆ toUnit_floor()

SingleParticlePos toUnit_floor ( const TinyVector< T1, D > &  r) const
inline

Definition at line 155 of file CrystalLattice.h.

Referenced by EinsplineSetBuilder::ReadOrbitalInfo_ESHDF().

156  {
157  SingleParticlePos val_dot;
158  val_dot = toUnit(r);
159  for (int i = 0; i < D; i++)
160  if (-std::numeric_limits<T1>::epsilon() < val_dot[i] && val_dot[i] < 0)
161  val_dot[i] = T1(0.0);
162  else
163  val_dot[i] -= std::floor(val_dot[i]);
164  return val_dot;
165  }
TinyVector< T, D > SingleParticlePos
the type of a D-dimensional position vector
SingleParticlePos toUnit(const TinyVector< T1, D > &r) const
Convert a cartesian vector to a unit vector.
MakeReturn< UnaryNode< FnFloor, typename CreateLeaf< Vector< T1, C1 > >::Leaf_t > >::Expression_t floor(const Vector< T1, C1 > &l)

Member Data Documentation

◆ ABC

Definition at line 126 of file CrystalLattice.h.

◆ BoxBConds

TinyVector<int, D> BoxBConds

◆ CellRadiusSq

Scalar_t CellRadiusSq

SimulationCellRadius*SimulationCellRadius.

Definition at line 93 of file CrystalLattice.h.

◆ Center

◆ DiagonalOnly

bool DiagonalOnly

true, if off-diagonal elements are zero so that other classes can take advantage of this

Definition at line 77 of file CrystalLattice.h.

◆ explicitly_defined

◆ G

◆ Gt

Transpose of reciprocal unit vectors:

Definition at line 101 of file CrystalLattice.h.

◆ Gv

Reciprocal unit vectors.

Introduced to efficiently return one vector at a time. Gv[i] is D-dim vector of the ith direction.

Definition at line 123 of file CrystalLattice.h.

Referenced by CrystalLattice< ST, 3 >::b(), and MomentumDistribution::MomentumDistribution().

◆ Length

Length[idim] length of the idim-th lattice vector.

Definition at line 107 of file CrystalLattice.h.

◆ M

Metric tensor.

Definition at line 103 of file CrystalLattice.h.

Referenced by CrystalLattice< ST, 3 >::Dot().

◆ Mg

Metric tensor for G vectors.

Definition at line 105 of file CrystalLattice.h.

Referenced by CrystalLattice< ST, 3 >::ksq().

◆ OneOverLength

SingleParticlePos OneOverLength

OneOverLength[idim] 1/length of the idim-th lattice vector.

Definition at line 109 of file CrystalLattice.h.

◆ R

◆ Rv

Real-space unit vectors.

Introduced to efficiently return one vector at a time. Rv[i] is D-dim vector of the ith direction.

Definition at line 117 of file CrystalLattice.h.

Referenced by CrystalLattice< ST, 3 >::a(), SpinDensity::put(), SpinDensity::report(), SpinDensityNew::report(), OrbitalImages::report(), SimulationCell::resetLRBox(), qmcplusplus::TEST_CASE(), and OrbitalImages::write_orbital_xsf().

◆ SimulationCellRadius

Scalar_t SimulationCellRadius

simulation cell radii

Definition at line 91 of file CrystalLattice.h.

Referenced by LatticeParser::put().

◆ SuperCellEnum

◆ VacuumScale

T VacuumScale

The scale factor for adding vacuum.

Definition at line 83 of file CrystalLattice.h.

Referenced by CrystalLattice< ST, 3 >::operator=(), LatticeParser::put(), SimulationCell::resetLRBox(), and qmcplusplus::TEST_CASE().

◆ Volume

◆ WignerSeitzRadius

Scalar_t WignerSeitzRadius

Wigner-Seitz cell radius.

Definition at line 89 of file CrystalLattice.h.

Referenced by LatticeParser::put().

◆ WignerSeitzRadius_G

Scalar_t WignerSeitzRadius_G

Wigner-Seitz cell radius in reciprocal space.

Definition at line 95 of file CrystalLattice.h.

Referenced by MomentumDistribution::MomentumDistribution().


The documentation for this struct was generated from the following files: