CASM
AClustersApproachtoStatisticalMechanics
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules
Niggli.hh
Go to the documentation of this file.
1 #ifndef NIGGLI_HH
2 #define NIGGLI_HH
3 
4 #include "casm/external/Eigen/Core"
5 
6 namespace CASM {
7  class Lattice;
8  class SymGroup;
9 
46  class NiggliRep {
47  public:
48 
49  NiggliRep(const Lattice &init_lat);
50 
51  NiggliRep(const Eigen::Matrix3d &init_lat_col_mat);
52 
54  double A() const;
55 
57  double B() const;
58 
60  double C() const;
61 
63  double ksi() const;
64 
66  double eta() const;
67 
69  double zeta() const;
70 
71  const Eigen::Matrix3d &metrical_matrix() const;
72 
74  bool meets_criteria_1(double compare_tol) const;
75 
77  bool meets_criteria_2(double compare_tol) const;
78 
80  bool meets_criteria_3(double compare_tol) const;
81 
83  bool meets_criteria_4(double compare_tol) const;
84 
86  bool meets_criteria_5(double compare_tol) const;
87 
89  bool meets_criteria_6(double compare_tol) const;
90 
92  bool meets_criteria_7(double compare_tol) const;
93 
96  bool meets_criteria_8(double compare_tol) const;
97 
99  bool is_niggli(double compare_tol) const;
100 
102  bool is_niggli_type1(double compare_tol) const;
103 
105  bool is_niggli_type2(double compare_tol) const;
106 
107  void debug_criteria(double compare_tol) const;
108 
109  private:
110 
113 
114  };
115 
116 
118  Lattice canonical_equivalent_lattice(const Lattice &in_lat, const SymGroup &point_grp, double compare_tol);
119 
121  Lattice niggli(const Lattice &in_lat, double compare_tol, bool keep_handedness = false);
122 
124  bool is_niggli(const Eigen::Matrix3d &test_lat_mat, double compare_tol);
125 
127  bool is_niggli(const Lattice &test_lat, double compare_tol);
128 
130  Eigen::VectorXd spatial_unroll(const Eigen::Matrix3d &lat_mat, double compare_tol);
131 
133  bool standard_orientation_spatial_compare(const Eigen::Matrix3d &low_score_lat_mat, Eigen::Matrix3d &high_score_lat_mat, double compare_tol);
134 
136  bool standard_orientation_compare(const Eigen::Matrix3d &low_score_lat_mat, const Eigen::Matrix3d &high_score_lat_mat, double compare_tol);
137 
139 }
140 
141 #endif
bool meets_criteria_2(double compare_tol) const
B<=C OR (B==C, |eta| <= |zeta|)
Definition: Niggli.cc:59
bool standard_orientation_spatial_compare(const Eigen::Matrix3d &low_score_lat_mat, Eigen::Matrix3d &high_score_lat_mat, double compare_tol)
Compare the spatial orientation (ignoring matrix symmetry) and determine which one is oriented more s...
double A() const
Square of lattice length a.
Definition: Niggli.cc:15
bool meets_criteria_7(double compare_tol) const
|zeta|<=A OR (zeta==A, eta<=2*ksi) OR (zeta==-A, eta==0)
Definition: Niggli.cc:139
bool is_niggli_type2(double compare_tol) const
True if all conditions except 3 are true.
Definition: Niggli.cc:182
bool standard_orientation_compare(const Eigen::Matrix3d &low_score_lat_mat, const Eigen::Matrix3d &high_score_lat_mat, double compare_tol)
Determine whether high_score has a more standard format than low_score.
Definition: Niggli.cc:391
Main CASM namespace.
Definition: complete.cpp:8
double C() const
Square of lattice length c.
Definition: Niggli.cc:23
SymGroup is a collection of symmetry operations that satisfy the group property The symmetry operatio...
Definition: SymGroup.hh:33
bool is_niggli_type1(double compare_tol) const
True if all conditions except 4 are true.
Definition: Niggli.cc:178
Eigen::VectorXd spatial_unroll(const Eigen::Matrix3d &lat_mat, double compare_tol)
Generate a vector whose lexicographical value determines how well it's oriented in space...
Definition: Niggli.cc:319
bool meets_criteria_4(double compare_tol) const
For type II: ksi<=0 && eta<=0 && zeta<=0 (all angles >= 90)
Definition: Niggli.cc:86
double zeta() const
2ab*cos(gamma)
Definition: Niggli.cc:35
bool is_niggli(const Eigen::Matrix3d &test_lat_mat, double compare_tol)
Check whether the given lattice (represented as a matrix) is in niggli TYPE ?? reduced form (does not...
Definition: Niggli.cc:218
double B() const
Square of lattice length b.
Definition: Niggli.cc:19
Lattice canonical_equivalent_lattice(const Lattice &in_lat, const SymGroup &point_grp, double compare_tol)
Find the niggli, most standard oriented version of the given orbit (defined by the given SymGroup) of...
Definition: Niggli.cc:276
NiggliRep(const Lattice &init_lat)
Definition: Niggli.cc:11
const Eigen::Matrix3d & metrical_matrix() const
Definition: Niggli.cc:39
Eigen::VectorXd VectorXd
double ksi() const
2bc*cos(alpha)
Definition: Niggli.cc:27
bool meets_criteria_1(double compare_tol) const
A<=B OR (A==B, |ksi| <= |eta|)
Definition: Niggli.cc:43
bool meets_criteria_8(double compare_tol) const
Definition: Niggli.cc:160
const Eigen::Matrix3d m_metrical_matrix
Transpose of initialization lattice dotted with itself.
Definition: Niggli.hh:112
Eigen::Matrix3d Matrix3d
bool is_niggli(double compare_tol) const
True if all conditions are true, and either 4 OR 3 is false.
Definition: Niggli.cc:186
void debug_criteria(double compare_tol) const
Definition: Niggli.cc:202
bool meets_criteria_3(double compare_tol) const
For type I: ksi>0 && eta>0 && zeta>0 (all angles < 90)
Definition: Niggli.cc:75
double eta() const
2ac*cos(beta)
Definition: Niggli.cc:31
bool meets_criteria_6(double compare_tol) const
|eta|<=A OR (eta==A, zeta<=2*ksi) OR (eta==-A, zeta==0)
Definition: Niggli.cc:118
Lattice niggli(const Lattice &in_lat, double compare_tol, bool keep_handedness=false)
Convert the given lattice into it's niggli reduced form, with the most standard orientation possilbe...
Definition: Niggli.cc:240
bool meets_criteria_5(double compare_tol) const
|ksi|<=B OR (ksi==B, zeta<=2*eta) OR (ksi==-B, zeta==0)
Definition: Niggli.cc:97