CASM
AClustersApproachtoStatisticalMechanics
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules
CASM::Lattice Class Reference

#include <Lattice.hh>

+ Inheritance diagram for CASM::Lattice:

Detailed Description

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...
 
Latticemake_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 Latticederived () 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...
 

Member Typedef Documentation

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.

Constructor & Destructor Documentation

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.

Member Function Documentation

bool CASM::Lattice::_eq ( const Lattice RHS) const
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.

void CASM::Lattice::_generate_voronoi_table ( ) const
private

populate voronoi information.

Definition at line 563 of file Lattice.cc.

bool CASM::Comparisons< Lattice >::_ne ( const Lattice B) const
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.

Lattice CASM::Lattice::bcc ( )
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.

Array< int > CASM::Lattice::calc_kpoints ( Array< int >  prim_kpoints,
Lattice  prim_lat 
)

Calculates the kpoint mesh for a supercell lattice given the kpoint mesh for the primitive lattice

Definition at line 161 of file Lattice.cc.

Lattice CASM::Lattice::canonical_form ( double  tol = TOL) const

Returns the canonical equivalent Lattice, using the point group of the Lattice.

  • To specify different symmetry, see canonical_equivalent_lattice

Definition at line 289 of file Lattice.cc.

Lattice CASM::Lattice::canonical_form ( const SymGroup pg,
double  tol = TOL 
) const

Returns the canonical equivalent Lattice, using the provided point group.

Definition at line 298 of file Lattice.cc.

Lattice CASM::Lattice::cubic ( )
static

Construct simple cubic primitive cell of unit volume.

Definition at line 69 of file Lattice.cc.

const Lattice & CASM::Comparisons< Lattice >::derived ( ) const
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.

Lattice CASM::Lattice::fcc ( )
static

Construct FCC primitive cell of unit volume.

Definition at line 45 of file Lattice.cc.

template<typename SymOpIterator , typename SymOpOutputIterator >
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

Parameters
sub_groupwith subset of
super_groupthat leaves this lattice invariant
sub_groupshould be empty

Definition at line 222 of file Lattice.cc.

SymOp CASM::Lattice::from_canonical ( double  tol = TOL) const

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.

SymOp CASM::Lattice::from_canonical ( const SymGroup pg,
double  tol = TOL 
) const

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.

void CASM::Lattice::generate_point_group ( SymGroup point_group,
double  pg_tol = TOL 
) const

Populate.

Parameters
point_groupwith the point group of this lattice
point_groupshould be empty
pg_tolcan 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.

Parameters
supercellwith symmetrically distinct supercells of this lattice Superlattices are enumerated with volumes
min_prim_vol<= volume <=
max_prim_vol
effective_pgis 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.

template<typename CoordType , typename CoordType2 >
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.

Parameters
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.

Lattice CASM::Lattice::hexagonal ( )
static

Construct cubic primitive cell of unit volume.

Definition at line 75 of file Lattice.cc.

double CASM::Lattice::inner_voronoi_radius ( ) const
inline

Radius of the largest sphere that fits wholly within the voronoi cell.

Definition at line 98 of file Lattice.hh.

const Eigen::Matrix3d& CASM::Lattice::inv_lat_column_mat ( ) const
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.

bool CASM::Lattice::is_canonical ( const SymGroup pg,
double  tol = TOL 
) const

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.

bool CASM::Lattice::is_supercell_of ( const Lattice tile,
double  _tol = TOL 
) const

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.

const Eigen::Matrix3d& CASM::Lattice::lat_column_mat ( ) const
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.

bool CASM::Comparisons< Lattice >::operator!= ( const Lattice B) const
inlineinherited

Definition at line 41 of file Comparisons.hh.

bool CASM::Lattice::operator< ( const Lattice RHS) const

Compare two Lattice.

  • First compares is_niggli(*this, TOL) with is_niggli(RHS, TOL)
  • Then compares via standard_orientation_compare(this->lat_column_mat(), RHS.lat_column_mat(), TOL)

Definition at line 1013 of file Lattice.cc.

bool CASM::Comparisons< Lattice >::operator<= ( const Lattice B) const
inlineinherited

Definition at line 29 of file Comparisons.hh.

bool CASM::Comparisons< Lattice >::operator== ( const Lattice B) const
inlineinherited

Definition at line 37 of file Comparisons.hh.

bool CASM::Comparisons< Lattice >::operator> ( const Lattice B) const
inlineinherited

Definition at line 25 of file Comparisons.hh.

bool CASM::Comparisons< Lattice >::operator>= ( const Lattice B) const
inlineinherited

Definition at line 33 of file Comparisons.hh.

LatVec CASM::Lattice::operator[] ( Index  i)
inline

Get i'th lattice vector as column expression.

Definition at line 56 of file Lattice.hh.

ConstLatVec CASM::Lattice::operator[] ( Index  i) const
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.

Parameters
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.

Parameters
_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.

SymOp CASM::Lattice::to_canonical ( double  tol = TOL) const

Returns the operation that applied to *this returns the canonical form.

Definition at line 249 of file Lattice.cc.

SymOp CASM::Lattice::to_canonical ( const SymGroup pg,
double  tol = TOL 
) const

Returns the operation that applied to *this returns the canonical form.

Definition at line 261 of file Lattice.cc.

std::tuple<LatVec, LatVec, LatVec> CASM::Lattice::vectors ( )
inline

Definition at line 65 of file Lattice.hh.

std::tuple<ConstLatVec, ConstLatVec, ConstLatVec> CASM::Lattice::vectors ( ) const
inline

Definition at line 69 of file Lattice.hh.

double CASM::Lattice::vol ( ) const
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.

Eigen::MatrixXd const& CASM::Lattice::voronoi_table ( ) const
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.

Friends And Related Function Documentation

Eigen::Matrix3i canonical_hnf ( const Eigen::Matrix3i &  T,
const SymGroup effective_pg,
const Lattice ref_lattice 
)
related

Return canonical hermite normal form of the supercell matrix.

Returns
Eigen::Matrix3i of H in canonical form
Parameters
TA supercell matrix (Eigen::Matrix3i), such that S = U*T, where S is the superlattice and U the unit lattice, as column vector matrices
ref_latticeThe lattice the transformation matrix T is meant to be acting on
effective_pgGroup 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]
form the highest lexicographic order when considering equivalent superlattices by point group operations.
  • Equivalent superlattices can be obtained using point group operations: S' = op*S = U*H*V, where V is integer and has determinant +/- 1
  • Substituting S = U*T, we have op*U*T = U*H*V.
  • Or H*V = U.inverse*op*U*T, which is the hermite normal form of U.inverse*op*U*T
  • So T is canonical if it is in hermite normal form and for all point group operations it has a higher lexicographic order than the resulting H

Definition at line 653 of file SupercellEnumerator.cc.

Lattice replace_vector ( const Lattice lat,
const Eigen::Vector3d &  new_vector,
double  tol 
)
related

Returns a minimum volume Lattice obtainable by replacing one Lattice vector with tau.

  • No guarantee on the result being canonical in any way

Definition at line 1214 of file Lattice.cc.

Member Data Documentation

friend CASM::Lattice::Comparisons< Lattice >
private

Definition at line 253 of file Lattice.hh.

double CASM::Lattice::m_inner_voronoi_radius
mutableprivate

Definition at line 261 of file Lattice.hh.

Eigen::Matrix3d CASM::Lattice::m_inv_lat_mat
private

Definition at line 269 of file Lattice.hh.

Eigen::Matrix3d CASM::Lattice::m_lat_mat
private

Definition at line 269 of file Lattice.hh.

Eigen::MatrixXd CASM::Lattice::m_voronoi_table
mutableprivate

Definition at line 262 of file Lattice.hh.


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