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 private:
77 // static constexpr double ROUNDOFF_MARGIN_MULTIPLIER = 1.000001;
78 
80  unsigned int getMajorCellIndexInternal(Vector3d coordinates) const;
81 
84  Vector3d getFirstCellCoordinates(Vector3d coordinate) const;
85 };
86 
87 inline unsigned int WignerSeitzCell::getMajorCellIndexInternal(
88  Vector3d coordinates
89 ) const{
90  const std::vector<Vector3d> &basisVectors = getBasisVectors();
91  double coordinatesNorm = coordinates.norm();
92 
93  unsigned int cellIndex = 0;
94  switch(getNumDimensions()){
95  case 1:
96  {
97  TBTKExit(
98  "WignerSeitzCell::getMajorCellIndexInternal()",
99  "Not yet implemented for 1D.",
100  ""
101  );
102  }
103  case 2:
104  {
105  const Vector3d &b0 = basisVectors.at(0);
106  const Vector3d &b1 = basisVectors.at(1);
107  Vector3d b2 = (b0*b1).unit();
108 
109  double minDistanceLine0 = Vector3d::dotProduct(b0, (b1*b2).unit());
110  double minDistanceLine1 = Vector3d::dotProduct(b1, (b2*b0).unit());
111  double X = abs(2*coordinatesNorm/minDistanceLine0);
112  double Y = abs(2*coordinatesNorm/minDistanceLine1);
113 
114  for(int x = -X-1; x < X+1; x++){
115  for(int y = -Y-1; y < Y+1; y++){
116  if(x == 0 && y == 0)
117  continue;
118 
119  Vector3d latticePoint
120  = x*basisVectors.at(0)
121  + y*basisVectors.at(1);
122 
123  if(
125  latticePoint,
126  coordinates
128  latticePoint,
129  latticePoint
130  ) > 1/2.
131  ){
132  cellIndex++;
133  }
134  }
135  }
136 
137  break;
138  }
139  case 3:
140  {
141  const Vector3d &b0 = basisVectors.at(0);
142  const Vector3d &b1 = basisVectors.at(1);
143  const Vector3d &b2 = basisVectors.at(2);
144 
145  double minDistancePlane0 = Vector3d::dotProduct(b0, (b1*b2).unit());
146  double minDistancePlane1 = Vector3d::dotProduct(b1, (b2*b0).unit());
147  double minDistancePlane2 = Vector3d::dotProduct(b2, (b0*b1).unit());
148  double X = abs(2*coordinatesNorm/minDistancePlane0);
149  double Y = abs(2*coordinatesNorm/minDistancePlane1);
150  double Z = abs(2*coordinatesNorm/minDistancePlane2);
151 
152  for(int x = -X-1; x < X+1; x++){
153  for(int y = -Y-1; y < Y+1; y++){
154  for(int z = -Z-1; z < Z+1; z++){
155  if(x == 0 && y == 0 && z == 0)
156  continue;
157 
158  Vector3d latticePoint
159  = x*basisVectors.at(0)
160  + y*basisVectors.at(1)
161  + z*basisVectors.at(2);
162 
163  if(
165  latticePoint,
166  coordinates
168  latticePoint,
169  latticePoint
170  ) > 1/2.
171  ){
172  cellIndex++;
173  }
174  }
175  }
176  }
177 
178  break;
179  }
180  default:
181  TBTKExit(
182  "WignerSeitzCell::getCellIndex()",
183  "Only coordinates with 1-3 componenents supported.",
184  ""
185  );
186  }
187 
188  return cellIndex;
189 }
190 
191 inline Vector3d WignerSeitzCell::getFirstCellCoordinates(Vector3d coordinates) const{
192  const std::vector<Vector3d> &basisVectors = getBasisVectors();
193  Vector3d b[3];
194  for(unsigned int n = 0; n < 3; n++)
195  b[n] = basisVectors.at(n);
196 
197  for(unsigned int n = 0; n < getNumDimensions(); n++){
198  double coordinatesProjection = Vector3d::dotProduct(
199  coordinates,
200  (b[(n+1)%3]*b[(n+2)%3]).unit()
201  );
202  double bProjection = Vector3d::dotProduct(
203  b[n],
204  (b[(n+1)%3]*b[(n+2)%3]).unit()
205  );
206 
207  coordinates = coordinates - ((int)((coordinatesProjection - bProjection)/bProjection + 3/2.))*b[n];
208  }
209 
210  bool done = false;
211  while(!done){
212  done = true;
213  for(int x = -1; x <= 1; x++){
214  for(int y = -1; y <= 1; y++){
215  for(int z = -1; z <= 1; z++){
216  const Vector3d v = x*b[0] + y*b[1] + z*b[2];
217 
218  if((coordinates + v).norm() < coordinates.norm()){
219  coordinates = coordinates + v;
220  done = false;
221  }
222  }
223  }
224  }
225  }
226 
227  return coordinates;
228 }
229 
230 }; //End namespace TBTK
231 
232 #endif
const std::vector< Vector3d > & getBasisVectors() const
Definition: SpacePartition.h:205
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:153
double norm() const
Definition: Vector3d.h:184
virtual std::vector< std::vector< double > > getMinorMesh(const std::vector< unsigned int > &numMeshPoints) const
Definition: ParallelepipedCell.h:36
MeshType
Definition: SpacePartition.h:46
Flexible physical index.
Definition: Index.h:70
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.