CASM
AClustersApproachtoStatisticalMechanics
|
#include <Lattice.hh>
Definition at line 31 of file Lattice.hh.
Public Types | |
typedef Eigen::Matrix3d::ColXpr | LatVec |
typedef Eigen::Matrix3d::ConstColXpr | ConstLatVec |
Public Member Functions | |
Lattice (const Eigen::Vector3d &vec1, const Eigen::Vector3d &vec2, const Eigen::Vector3d &vec3) | |
Lattice (const Eigen::Ref< const Eigen::Matrix3d > &lat_mat=Eigen::Matrix3d::Identity()) | |
LatVec | operator[] (Index i) |
Get i'th lattice vector as column expression. More... | |
ConstLatVec | operator[] (Index i) const |
Get i'th lattice vector as column expression. More... | |
std::tuple< LatVec, LatVec, LatVec > | vectors () |
std::tuple< ConstLatVec, ConstLatVec, ConstLatVec > | vectors () const |
Lattice | scaled_lattice (double scale) const |
Return scaled copy of this lattice (Note: Volume will be scaled by scale^3) More... | |
double | angle (Index i) const |
Return angle between lattice vectors (*this)[(i+1)%3] and (*this)[(i+2)%3], in degrees. More... | |
double | length (Index i) const |
Return length of i'th lattice vector. More... | |
double | vol () const |
Return signed volume of this lattice. More... | |
Eigen::MatrixXd const & | voronoi_table () const |
Return voronoi table, which specifies outward-pointing normals of Lattice Voronoi cell. outward-pointing normals are given as rows of the matrix, and defined such that if (voronoi_table()*coord.cart()).maxCoeff()>1, then 'coord' is outside of the voronoi cell. More... | |
double | inner_voronoi_radius () const |
Radius of the largest sphere that fits wholly within the voronoi cell. More... | |
const Eigen::Matrix3d & | lat_column_mat () const |
3x3 matrix with lattice vectors as its columne More... | |
const Eigen::Matrix3d & | inv_lat_column_mat () const |
Inverse of Lattice::lat_column_mat() It is the transformation matrix 'C2F', such that f = C2F * c where 'c' is Cartesian coordinate/vector and 'f' and 'f' is fractional coordinate/vector. More... | |
Array< int > | calc_kpoints (Array< int > prim_kpoints, Lattice prim_lat) |
Lattice | get_reciprocal () const |
Return reciprocal lattice. More... | |
void | generate_point_group (SymGroup &point_group, double pg_tol=TOL) const |
Populate. More... | |
template<typename SymOpIterator , typename SymOpOutputIterator > | |
SymOpOutputIterator | find_invariant_subgroup (SymOpIterator begin, SymOpIterator end, SymOpOutputIterator result, double pg_tol=TOL) const |
Output the SymOp that leave this lattice invariant. More... | |
void | find_invariant_subgroup (const SymGroup &super_group, SymGroup &sub_group, double pg_tol=TOL) const |
populate More... | |
bool | is_canonical (double tol=TOL) const |
Check if Lattice is in the canonical form. More... | |
bool | is_canonical (const SymGroup &pg, double tol=TOL) const |
Check if Lattice is in the canonical form. More... | |
SymOp | to_canonical (double tol=TOL) const |
Returns the operation that applied to *this returns the canonical form. More... | |
SymOp | to_canonical (const SymGroup &pg, double tol=TOL) const |
Returns the operation that applied to *this returns the canonical form. More... | |
SymOp | from_canonical (double tol=TOL) const |
Returns the operation that applied to the canonical form returns *this. More... | |
SymOp | from_canonical (const SymGroup &pg, double tol=TOL) const |
Returns the operation that applied to the canonical form returns *this. More... | |
Lattice | canonical_form (double tol=TOL) const |
Returns the canonical equivalent Lattice, using the point group of the Lattice. More... | |
Lattice | canonical_form (const SymGroup &pg, double tol=TOL) const |
Returns the canonical equivalent Lattice, using the provided point group. More... | |
void | generate_supercells (Array< Lattice > &supercell, const SymGroup &effective_pg, const ScelEnumProps &enum_props) const |
Populate. More... | |
template<typename T > | |
Lattice | make_supercell (const Eigen::Matrix< T, 3, 3 > &trans_mat) const |
make a supercell of this lattice. Equivalent to Lattice(lat_column_mat()*trans_mat) More... | |
Lattice | get_reduced_cell () const |
Find the lattice vectors which give most compact unit cell Compactness is measured by how close lat_column_mat().transpose()*lat_column_mat() is to a diagonal matrix. More... | |
void | print_voronoi_table (std::ostream &stream) const |
double | min_voronoi_radius () const |
Radius of largest sphere that totally fits inside the voronoi cell. More... | |
double | max_voronoi_measure (const Eigen::Vector3d &pos, Eigen::Vector3d &lattice_trans) const |
int | voronoi_number (const Eigen::Vector3d &pos) const |
Eigen::Vector3i | enclose_sphere (double radius) const |
template<typename CoordType , typename CoordType2 > | |
Array< CoordType > | gridstruc_build (double max_radius, double min_radius, Array< CoordType > basis, CoordType2 lat_point) |
Make a grid of lattice sites such that min_radius <= distance <= max_radius from. More... | |
void | read (std::istream &stream) |
void | print (std::ostream &stream, int _prec=8) const |
bool | is_equivalent (const Lattice &RHS, double tol) const |
Are two lattices the same, even if they have different lattice vectors. More... | |
bool | operator< (const Lattice &RHS) const |
Compare two Lattice. More... | |
bool | is_supercell_of (const Lattice &tile, Eigen::Matrix3d &multimat, double _tol=TOL) const |
Matrix that relates two lattices (e.g., strain or slat) More... | |
bool | is_supercell_of (const Lattice &tile, const Array< SymOp > &symlist, Eigen::Matrix3d &multimat, double _tol=TOL) const |
bool | is_supercell_of (const Lattice &tile, double _tol=TOL) const |
Checks if lattice is a supercell of tile, applying operations from symlist. More... | |
bool | is_supercell_of (const Lattice &tile, const Array< SymOp > &symlist, double _tol=TOL) const |
Lattice | box (const Lattice &prim, const Lattice &scel, bool verbose=false) const |
Return a lattice with diagonal matrix that fits around starting lattice. More... | |
Lattice & | make_right_handed () |
Flip c vector if it's on the wrong side of a-b plane – return (*this) More... | |
bool | is_right_handed () const |
Check if the lattice is right handed. More... | |
Eigen::Vector3i | get_millers (Eigen::Vector3d plane_normal, double tolerance=TOL) const |
Given a normal vector, a Vector3 containing the miller indeces for the lattice is generated. More... | |
Lattice | get_lattice_in_plane (Eigen::Vector3i millers, int max_vol=20) const |
Generates a lattice with vectors a and b parallel to the plane described by the miller indeces. More... | |
Array< double > | pg_converge (double large_tol) |
void | pg_converge (double small_tol, double large_tol, double increment) |
void | symmetrize (const SymGroup &relaxed_pg) |
Force this lattice to have symmetry of group. More... | |
void | symmetrize (double _tol) |
Force this lattice to have symmetry of point group calculated based on tolerance. More... | |
bool | operator> (const Lattice &B) const |
bool | operator<= (const Lattice &B) const |
bool | operator>= (const Lattice &B) const |
bool | operator== (const Lattice &B) const |
bool | operator!= (const Lattice &B) const |
Static Public Member Functions | |
static Lattice | fcc () |
Construct FCC primitive cell of unit volume. More... | |
static Lattice | bcc () |
Construct BCC primitive cell of unit volume. More... | |
static Lattice | cubic () |
Construct simple cubic primitive cell of unit volume. More... | |
static Lattice | hexagonal () |
Construct cubic primitive cell of unit volume. More... | |
Protected Member Functions | |
const Lattice & | derived () const |
bool | _ne (const Lattice &B) const |
Private Member Functions | |
bool | _eq (const Lattice &RHS) const |
Are lattice vectors identical for two lattices, within TOL. More... | |
void | _generate_voronoi_table () const |
populate voronoi information. More... | |
Private Attributes | |
friend | Comparisons< Lattice > |
double | m_inner_voronoi_radius |
Eigen::MatrixXd | m_voronoi_table |
Eigen::Matrix3d | m_lat_mat |
Eigen::Matrix3d | m_inv_lat_mat |
Related Functions | |
(Note that these are not member functions.) | |
Lattice | replace_vector (const Lattice &lat, const Eigen::Vector3d &new_vector, double tol) |
Returns a minimum volume Lattice obtainable by replacing one Lattice vector with tau. More... | |
Eigen::Matrix3i | canonical_hnf (const Eigen::Matrix3i &T, const SymGroup &effective_pg, const Lattice &ref_lattice) |
Return canonical hermite normal form of the supercell matrix. More... | |
typedef Eigen::Matrix3d::ConstColXpr CASM::Lattice::ConstLatVec |
Definition at line 34 of file Lattice.hh.
typedef Eigen::Matrix3d::ColXpr CASM::Lattice::LatVec |
Definition at line 33 of file Lattice.hh.
CASM::Lattice::Lattice | ( | const Eigen::Vector3d & | vec1, |
const Eigen::Vector3d & | vec2, | ||
const Eigen::Vector3d & | vec3 | ||
) |
Definition at line 27 of file Lattice.cc.
CASM::Lattice::Lattice | ( | const Eigen::Ref< const Eigen::Matrix3d > & | lat_mat = Eigen::Matrix3d::Identity() | ) |
Construct Lattice from a matrix of lattice vectors, where lattice vectors are columns (e.g., lat_mat is equivalent to coord_trans[FRAC])
Construct Lattice from a matrix of lattice vectors, where lattice vectors are columns (e.g., lat_mat is equivalent to lat_column_mat())
Definition at line 38 of file Lattice.cc.
|
private |
Are lattice vectors identical for two lattices, within TOL.
Are lattice vectors identical for two lattices.
Definition at line 1024 of file Lattice.cc.
|
private |
populate voronoi information.
Definition at line 563 of file Lattice.cc.
|
inlineprotectedinherited |
Definition at line 56 of file Comparisons.hh.
double CASM::Lattice::angle | ( | Index | i | ) | const |
Return angle between lattice vectors (*this)[(i+1)%3] and (*this)[(i+2)%3], in degrees.
Definition at line 99 of file Lattice.cc.
|
static |
Construct BCC primitive cell of unit volume.
Definition at line 57 of file Lattice.cc.
Lattice CASM::Lattice::box | ( | const Lattice & | prim, |
const Lattice & | scel, | ||
bool | verbose = false |
||
) | const |
Return a lattice with diagonal matrix that fits around starting lattice.
Calculates the kpoint mesh for a supercell lattice given the kpoint mesh for the primitive lattice
Definition at line 161 of file Lattice.cc.
Returns the canonical equivalent Lattice, using the point group of the Lattice.
Definition at line 289 of file Lattice.cc.
Returns the canonical equivalent Lattice, using the provided point group.
Definition at line 298 of file Lattice.cc.
|
static |
Construct simple cubic primitive cell of unit volume.
Definition at line 69 of file Lattice.cc.
|
inlineprotectedinherited |
Definition at line 48 of file Comparisons.hh.
Eigen::Vector3i CASM::Lattice::enclose_sphere | ( | double | radius | ) | const |
Gives a scaling of the lattice vectors such that after scaling, The eight parallelipipeds touching the origin enclose a sphere of given radius
This function, given a linearly independent set of lattice vectors, finds the dimensions along the unit cell vectors such that a sphere of given radius fits within a uniform grid of 2dim[1]x2dim[2]x2dim[3] lattice points centered at the origin.
The algorithm works by getting the normal (e.g. n1) to each pair of lattice vectors (e.g. a2, a3), scaling this normal to have length radius and then projecting this normal parallel to the a2,a3 plane onto the remaining lattice vector a1. This will tell us the number of a1 vectors needed to make a grid to encompass the sphere.
Definition at line 626 of file Lattice.cc.
|
static |
Construct FCC primitive cell of unit volume.
Definition at line 45 of file Lattice.cc.
SymOpOutputIterator CASM::Lattice::find_invariant_subgroup | ( | SymOpIterator | begin, |
SymOpIterator | end, | ||
SymOpOutputIterator | result, | ||
double | pg_tol = TOL |
||
) | const |
Output the SymOp that leave this lattice invariant.
Definition at line 77 of file Lattice_impl.hh.
void CASM::Lattice::find_invariant_subgroup | ( | const SymGroup & | super_group, |
SymGroup & | sub_group, | ||
double | pg_tol = TOL |
||
) | const |
populate
sub_group | with subset of |
super_group | that leaves this lattice invariant |
sub_group | should be empty |
Definition at line 222 of file Lattice.cc.
Returns the operation that applied to the canonical form returns *this.
Returns the operation that applied to *this returns the canonical form.
Definition at line 273 of file Lattice.cc.
Returns the operation that applied to the canonical form returns *this.
Returns the operation that applied to *this returns the canonical form.
Definition at line 280 of file Lattice.cc.
Populate.
point_group | with the point group of this lattice |
point_group | should be empty |
pg_tol | can be increased to find point group of lattice vectors that are slightly distorted due to numerical noise |
Definition at line 304 of file Lattice.cc.
void CASM::Lattice::generate_supercells | ( | Array< Lattice > & | supercell, |
const SymGroup & | effective_pg, | ||
const ScelEnumProps & | enum_props | ||
) | const |
Populate.
Generate super Lattice.
supercell | with symmetrically distinct supercells of this lattice Superlattices are enumerated with volumes |
min_prim_vol | <= volume <= |
max_prim_vol | |
effective_pg | is a group that should either be equivalent to the full point group of this lattice or be a subgroup of that full point group |
Use SupercellEnumerator to enumerate possible HNF transformation matrices. Unique supercells are identified by applying point group operations and keeping the supercell if the HNF is 'canonical', meaning that the HNF indices in order H00, H11, H22, H12, H02, H01 are the lexicographically greatest.
The supercell that is inserted in the 'supercell' container is the niggli cell, rotated to a standard orientation (see standard_orientation function).
See PrimcClex::generate_supercells for information on dims and G.
Definition at line 410 of file Lattice.cc.
Lattice CASM::Lattice::get_lattice_in_plane | ( | Eigen::Vector3i | millers, |
int | max_vol = 20 |
||
) | const |
Generates a lattice with vectors a and b parallel to the plane described by the miller indeces.
Using miller indeces, the intercept of the plane is calculated. All intercepts are multiplied by a factor such that their values are integers. This corresponds to the plane being shifted so that it lands on three lattice sites. Using these lattice sites new vectors A and B that lie on the desired plane can be constructed. The choice of vector C is somewhat arbitrary. get_lattice_in_plane will first construct C to be the shortest most orthogonal vector possible. The user is then given the option to make C more orthogonal to A and B in exchange for a greater length. The returned lattice will be the one containing the most orthogonal C that is still smaller than max_vol
Definition at line 717 of file Lattice.cc.
Eigen::Vector3i CASM::Lattice::get_millers | ( | Eigen::Vector3d | plane_normal, |
double | tolerance = TOL |
||
) | const |
Given a normal vector, a Vector3 containing the miller indeces for the lattice is generated.
Definition at line 694 of file Lattice.cc.
Lattice CASM::Lattice::get_reciprocal | ( | ) | const |
Return reciprocal lattice.
Definition at line 149 of file Lattice.cc.
Lattice CASM::Lattice::get_reduced_cell | ( | ) | const |
Find the lattice vectors which give most compact unit cell Compactness is measured by how close lat_column_mat().transpose()*lat_column_mat() is to a diagonal matrix.
This function finds the reduced cell from the given primitive cell.
First, the translation vectors of the primitive cell are determined, then reduced. The Bravais lattice is determined using Niggli's transformations. The method used to determine the basis vectors are to find body and face diagonals across the cell. The vectors satisfying these main conditions can be found by ensuring: i) the shortest face diagonal is not shorter than the longest edge of the same face ii) the shortest body diagonal is not shorter than the longest edge of the cell
This is accomplished by enumerating the skew matrices (identity matrix with one of the off-diagonal terms =1 or =-1). This will replace one of the lattice vectors with the face diagonal corresponding to that linear combination. If the length of any of the new skewed lattice vectors is shorter than the length of the original cell, it is replaced. All skew matrices have determinants equal to 1, which preserves the volume of the original cell.
Definition at line 442 of file Lattice.cc.
Array< CoordType > CASM::Lattice::gridstruc_build | ( | double | max_radius, |
double | min_radius, | ||
Array< CoordType > | basis, | ||
CoordType2 | lat_point | ||
) |
Make a grid of lattice sites such that min_radius <= distance <= max_radius from.
lat_point |
This function generates a grid of points between max_radius and min_radius. Additionally, it also fills up the points with a basis
Definition at line 15 of file Lattice_impl.hh.
|
static |
Construct cubic primitive cell of unit volume.
Definition at line 75 of file Lattice.cc.
|
inline |
Radius of the largest sphere that fits wholly within the voronoi cell.
Definition at line 98 of file Lattice.hh.
|
inline |
Inverse of Lattice::lat_column_mat() It is the transformation matrix 'C2F', such that f = C2F * c where 'c' is Cartesian coordinate/vector and 'f' and 'f' is fractional coordinate/vector.
Definition at line 113 of file Lattice.hh.
bool CASM::Lattice::is_canonical | ( | double | tol = TOL | ) | const |
Check if Lattice is in the canonical form.
Definition at line 235 of file Lattice.cc.
Check if Lattice is in the canonical form.
Definition at line 242 of file Lattice.cc.
bool CASM::Lattice::is_equivalent | ( | const Lattice & | RHS, |
double | tol | ||
) | const |
Are two lattices the same, even if they have different lattice vectors.
Are lattice vectors identical for two lattices.
Definition at line 1002 of file Lattice.cc.
bool CASM::Lattice::is_right_handed | ( | ) | const |
Check if the lattice is right handed.
Definition at line 1096 of file Lattice.cc.
bool CASM::Lattice::is_supercell_of | ( | const Lattice & | tile, |
Eigen::Matrix3d & | multimat, | ||
double | _tol = TOL |
||
) | const |
Matrix that relates two lattices (e.g., strain or slat)
Checks if lattice is a supercell of tile, acting on multiplication matrix. Check is performed applying operations from symlist
Definition at line 655 of file Lattice.cc.
bool CASM::Lattice::is_supercell_of | ( | const Lattice & | tile, |
const Array< SymOp > & | symlist, | ||
Eigen::Matrix3d & | multimat, | ||
double | _tol = TOL |
||
) | const |
Definition at line 645 of file Lattice.cc.
Checks if lattice is a supercell of tile, applying operations from symlist.
Definition at line 664 of file Lattice.cc.
bool CASM::Lattice::is_supercell_of | ( | const Lattice & | tile, |
const Array< SymOp > & | symlist, | ||
double | _tol = TOL |
||
) | const |
Definition at line 670 of file Lattice.cc.
|
inline |
3x3 matrix with lattice vectors as its columne
Definition at line 104 of file Lattice.hh.
double CASM::Lattice::length | ( | Index | i | ) | const |
Return length of i'th lattice vector.
Definition at line 93 of file Lattice.cc.
Lattice & CASM::Lattice::make_right_handed | ( | ) |
Flip c vector if it's on the wrong side of a-b plane – return (*this)
A lattice is considered right handed when the determinant of the lattice vector matrix is positive. This routine will flip the sign of all the lattice vectors if it finds that the determinant is negative
Definition at line 681 of file Lattice.cc.
double CASM::Lattice::max_voronoi_measure | ( | const Eigen::Vector3d & | pos, |
Eigen::Vector3d & | lattice_trans | ||
) | const |
Given cartesian vector 'pos', returns double v_dist and populates 'lattice_trans'. v_dist is the fractional number of half lattice-planes along 'lattice_trans' between 'pos' and the lattice plane that passes through the origin. lattice_trans is a lattice translation that is normal to a face of the voronoi cell and translating pos to 'pos-lattice_trans' will yield a translationally-equivalent vector between the -1 and +1 half-plane. v_dist is maximized over all faces of the voronoi cell.
Definition at line 530 of file Lattice.cc.
double CASM::Lattice::min_voronoi_radius | ( | ) | const |
Radius of largest sphere that totally fits inside the voronoi cell.
|
inlineinherited |
Definition at line 41 of file Comparisons.hh.
bool CASM::Lattice::operator< | ( | const Lattice & | RHS | ) | const |
Compare two Lattice.
Definition at line 1013 of file Lattice.cc.
|
inlineinherited |
Definition at line 29 of file Comparisons.hh.
|
inlineinherited |
Definition at line 37 of file Comparisons.hh.
|
inlineinherited |
Definition at line 25 of file Comparisons.hh.
|
inlineinherited |
Definition at line 33 of file Comparisons.hh.
Get i'th lattice vector as column expression.
Definition at line 56 of file Lattice.hh.
|
inline |
Get i'th lattice vector as column expression.
Definition at line 61 of file Lattice.hh.
Array< double > CASM::Lattice::pg_converge | ( | double | large_tol | ) |
Definition at line 351 of file Lattice.cc.
void CASM::Lattice::pg_converge | ( | double | small_tol, |
double | large_tol, | ||
double | increment | ||
) |
Definition at line 369 of file Lattice.cc.
void CASM::Lattice::print | ( | std::ostream & | stream, |
int | _prec = 8 |
||
) | const |
Definition at line 129 of file Lattice.cc.
void CASM::Lattice::print_voronoi_table | ( | std::ostream & | stream | ) | const |
void CASM::Lattice::read | ( | std::istream & | stream | ) |
Definition at line 116 of file Lattice.cc.
Lattice CASM::Lattice::scaled_lattice | ( | double | scale | ) | const |
Return scaled copy of this lattice (Note: Volume will be scaled by scale^3)
Definition at line 87 of file Lattice.cc.
void CASM::Lattice::symmetrize | ( | const SymGroup & | relaxed_pg | ) |
Force this lattice to have symmetry of group.
relaxed_pg |
Applies all operations of given SymGroup to the lattice and averages out the lattice vectors, changing your lattice to perfectly match with the SymGroup.
Definition at line 1038 of file Lattice.cc.
void CASM::Lattice::symmetrize | ( | double | tol | ) |
Force this lattice to have symmetry of point group calculated based on tolerance.
_tol |
Same as the other symmetrize routine, except this one assumes that the symmetry group you mean to use is the point group of your lattice within a certain tolerance.
Definition at line 1087 of file Lattice.cc.
Returns the operation that applied to *this returns the canonical form.
Definition at line 249 of file Lattice.cc.
Returns the operation that applied to *this returns the canonical form.
Definition at line 261 of file Lattice.cc.
Definition at line 65 of file Lattice.hh.
|
inline |
Definition at line 69 of file Lattice.hh.
|
inline |
Return signed volume of this lattice.
Definition at line 83 of file Lattice.hh.
int CASM::Lattice::voronoi_number | ( | const Eigen::Vector3d & | pos | ) | const |
return number of voronoi cell faces that 'pos' is on, within +/- TOL 0 indicates that 'pos' is within the voronoi cell, and the origin is the nearest lattice site -1 indicates that 'pos' is outside the voronoi cell, and there is a lattice site closer than the origin Values of 1<=n<=7 indicate that there n lattice sites equally as close as the origin
Definition at line 541 of file Lattice.cc.
|
inline |
Return voronoi table, which specifies outward-pointing normals of Lattice Voronoi cell. outward-pointing normals are given as rows of the matrix, and defined such that if (voronoi_table()*coord.cart()).maxCoeff()>1, then 'coord' is outside of the voronoi cell.
Definition at line 90 of file Lattice.hh.
|
related |
Return canonical hermite normal form of the supercell matrix.
T | A supercell matrix (Eigen::Matrix3i), such that S = U*T, where S is the superlattice and U the unit lattice, as column vector matrices |
ref_lattice | The lattice the transformation matrix T is meant to be acting on |
effective_pg | Group of symmetry operations to use for determining whether a canonical HNF matrix has been found. This is probably the point group of a structure or configuration Canonical form is such that T is in Hermite normal form (as from CASM::hermite_normal_form), and the unrolled coefficients [a f e]
[0 b d] -> abcdef
[0 0 c]
|
Definition at line 653 of file SupercellEnumerator.cc.
|
related |
Returns a minimum volume Lattice obtainable by replacing one Lattice vector with tau.
Definition at line 1214 of file Lattice.cc.
|
private |
Definition at line 253 of file Lattice.hh.
|
mutableprivate |
Definition at line 261 of file Lattice.hh.
|
private |
Definition at line 269 of file Lattice.hh.
|
private |
Definition at line 269 of file Lattice.hh.
|
mutableprivate |
Definition at line 262 of file Lattice.hh.