TBTK
WignerSeitzCell.h
Go to the documentation of this file.
1 /* Copyright 2016 Kristofer Bj√∂rnson
2  *
3  * Licensed under the Apache License, Version 2.0 (the "License");
4  * you may not use this file except in compliance with the License.
5  * You may obtain a copy of the License at
6  *
7  * http://www.apache.org/licenses/LICENSE-2.0
8  *
9  * Unless required by applicable law or agreed to in writing, software
10  * distributed under the License is distributed on an "AS IS" BASIS,
11  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
12  * See the License for the specific language governing permissions and
13  * limitations under the License.
14  */
15 
23 #ifndef COM_DAFER45_TBTK_WIGNER_SEITZ_CELL
24 #define COM_DAFER45_TBTK_WIGNER_SEITZ_CELL
25 
26 #include "TBTK/Index.h"
28 #include "TBTK/Vector3d.h"
29 
30 #include <initializer_list>
31 #include <vector>
32 
33 namespace TBTK{
34 
37 public:
40  std::initializer_list<std::initializer_list<double>> basisVectors,
41  MeshType meshType
42  );
43 
46  const std::vector<std::vector<double>> &basisVectors,
47  MeshType meshType
48  );
49 
51  virtual ~WignerSeitzCell();
52 
54  virtual Index getMajorCellIndex(
55  const std::vector<double> &coordinates
56  ) const;
57 
59  virtual Index getMinorCellIndex(
60  const std::vector<double> &coordinates,
61  const std::vector<unsigned int> &numMeshPoints
62  ) const;
63 
65  virtual std::vector<std::vector<double>> getMajorMesh(
66  const std::vector<unsigned int> &numMeshPoints
67  ) const;
68 
70  virtual std::vector<std::vector<double>> getMinorMesh(
71  const std::vector<unsigned int> &numMeshPoints
72  ) const;
73 
75  virtual std::vector<double> getMinorMeshPoint(
76  const std::vector<unsigned int> &meshPoint,
77  const std::vector<unsigned int> &numMeshPoints
78  ) const;
79 private:
83 // static constexpr double ROUNDOFF_MARGIN_MULTIPLIER = 1.000001;
84 
86  unsigned int getMajorCellIndexInternal(Vector3d coordinates) const;
87 
90  Vector3d getFirstCellCoordinates(Vector3d coordinate) const;
91 };
92 
93 inline unsigned int WignerSeitzCell::getMajorCellIndexInternal(
94  Vector3d coordinates
95 ) const{
96  const std::vector<Vector3d> &basisVectors = getBasisVectors();
97  double coordinatesNorm = coordinates.norm();
98 
99  unsigned int cellIndex = 0;
100  switch(getNumDimensions()){
101  case 1:
102  {
103  TBTKExit(
104  "WignerSeitzCell::getMajorCellIndexInternal()",
105  "Not yet implemented for 1D.",
106  ""
107  );
108  }
109  case 2:
110  {
111  const Vector3d &b0 = basisVectors.at(0);
112  const Vector3d &b1 = basisVectors.at(1);
113  Vector3d b2 = (b0*b1).unit();
114 
115  double minDistanceLine0 = Vector3d::dotProduct(b0, (b1*b2).unit());
116  double minDistanceLine1 = Vector3d::dotProduct(b1, (b2*b0).unit());
117  double X = abs(2*coordinatesNorm/minDistanceLine0);
118  double Y = abs(2*coordinatesNorm/minDistanceLine1);
119 
120  for(int x = -X-1; x < X+1; x++){
121  for(int y = -Y-1; y < Y+1; y++){
122  if(x == 0 && y == 0)
123  continue;
124 
125  Vector3d latticePoint
126  = x*basisVectors.at(0)
127  + y*basisVectors.at(1);
128 
129  if(
131  latticePoint,
132  coordinates
134  latticePoint,
135  latticePoint
136  ) > 1/2.
137  ){
138  cellIndex++;
139  }
140  }
141  }
142 
143  break;
144  }
145  case 3:
146  {
147  const Vector3d &b0 = basisVectors.at(0);
148  const Vector3d &b1 = basisVectors.at(1);
149  const Vector3d &b2 = basisVectors.at(2);
150 
151  double minDistancePlane0 = Vector3d::dotProduct(b0, (b1*b2).unit());
152  double minDistancePlane1 = Vector3d::dotProduct(b1, (b2*b0).unit());
153  double minDistancePlane2 = Vector3d::dotProduct(b2, (b0*b1).unit());
154  double X = abs(2*coordinatesNorm/minDistancePlane0);
155  double Y = abs(2*coordinatesNorm/minDistancePlane1);
156  double Z = abs(2*coordinatesNorm/minDistancePlane2);
157 
158  for(int x = -X-1; x < X+1; x++){
159  for(int y = -Y-1; y < Y+1; y++){
160  for(int z = -Z-1; z < Z+1; z++){
161  if(x == 0 && y == 0 && z == 0)
162  continue;
163 
164  Vector3d latticePoint
165  = x*basisVectors.at(0)
166  + y*basisVectors.at(1)
167  + z*basisVectors.at(2);
168 
169  if(
171  latticePoint,
172  coordinates
174  latticePoint,
175  latticePoint
176  ) > 1/2.
177  ){
178  cellIndex++;
179  }
180  }
181  }
182  }
183 
184  break;
185  }
186  default:
187  TBTKExit(
188  "WignerSeitzCell::getCellIndex()",
189  "Only coordinates with 1-3 componenents supported.",
190  ""
191  );
192  }
193 
194  return cellIndex;
195 }
196 
197 inline Vector3d WignerSeitzCell::getFirstCellCoordinates(Vector3d coordinates) const{
198  const std::vector<Vector3d> &basisVectors = getBasisVectors();
199  Vector3d b[3];
200  for(unsigned int n = 0; n < 3; n++)
201  b[n] = basisVectors.at(n);
202 
203  for(unsigned int n = 0; n < getNumDimensions(); n++){
204  double coordinatesProjection = Vector3d::dotProduct(
205  coordinates,
206  (b[(n+1)%3]*b[(n+2)%3]).unit()
207  );
208  double bProjection = Vector3d::dotProduct(
209  b[n],
210  (b[(n+1)%3]*b[(n+2)%3]).unit()
211  );
212 
213  coordinates = coordinates - ((int)((coordinatesProjection - bProjection)/bProjection + 3/2.))*b[n];
214  }
215 
216  bool done = false;
217  while(!done){
218  done = true;
219  for(int x = -1; x <= 1; x++){
220  for(int y = -1; y <= 1; y++){
221  for(int z = -1; z <= 1; z++){
222  const Vector3d v = x*b[0] + y*b[1] + z*b[2];
223 
224  if((coordinates + v).norm() < coordinates.norm()){
225  coordinates = coordinates + v;
226  done = false;
227  }
228  }
229  }
230  }
231  }
232 
233  return coordinates;
234 }
235 
236 }; //End namespace TBTK
237 
238 #endif
const std::vector< Vector3d > & getBasisVectors() const
Definition: SpacePartition.h:216
Flexible physical index.
WignerSeitzCell(std::initializer_list< std::initializer_list< double >> basisVectors, MeshType meshType)
Three-dimensional vector with components of double type.
unsigned int getNumDimensions() const
Definition: SpacePartition.h:164
double norm() const
Definition: Vector3d.h:184
virtual std::vector< std::vector< double > > getMinorMesh(const std::vector< unsigned int > &numMeshPoints) const
virtual std::vector< double > getMinorMeshPoint(const std::vector< unsigned int > &meshPoint, const std::vector< unsigned int > &numMeshPoints) const
Definition: ParallelepipedCell.h:36
MeshType
Definition: SpacePartition.h:46
Flexible physical index.
Definition: Index.h:69
Definition: Vector3d.h:33
Definition: ModelFactory.h:35
virtual Index getMajorCellIndex(const std::vector< double > &coordinates) const
static double dotProduct(const Vector3d &lhs, const Vector3d &rhs)
Definition: Vector3d.h:188
virtual Index getMinorCellIndex(const std::vector< double > &coordinates, const std::vector< unsigned int > &numMeshPoints) const
virtual std::vector< std::vector< double > > getMajorMesh(const std::vector< unsigned int > &numMeshPoints) const
Definition: WignerSeitzCell.h:36
virtual ~WignerSeitzCell()
Parallelepiped cell.