Loading...
Searching...
No Matches
chill Namespace Reference
Chill

CHILL and CHILL+ structure classification. More...

Classes

struct  BondClassifier
 One rule set for classifying bond correlations \(c_{ij}\). More...
struct  LinearClassifier
 One-versus-rest ridge classifier on a dense feature matrix. More...
struct  QlmAtom
 This is the local orientational bond order parameter \(q_{lm}\), of length \(2l+1\). More...
struct  SteinhardtQl
 Per-particle Steinhardt order parameters of a single degree \(l\). More...
struct  TemplateHit
struct  VoronoiWeights
 The Voronoi facet neighbours of one particle and their area weights. More...
struct  YlmAtom
 This contains a complex vector of length \(2l+1\). More...

Enumerations

enum class  CrystalKind {
  other = 0 , sc , fcc , hcp ,
  bcc
}

Functions

BondClassifier chillRule ()
BondClassifier chillPlusRule ()
BondClassifier bondClassifier (const std::string &name)
void registerBondClassifier (const std::string &name, const BondClassifier &rule)
std::vector< std::string > bondClassifierNames ()
 Names of every registered rule set.
void classifyBonds (molSys::PointCloud< molSys::Point< double >, double > &yCloud, const std::vector< std::vector< int > > &nList, const BondClassifier &rule, bool isSlice=false)
void getCorrel (molSys::PointCloud< molSys::Point< double >, double > &yCloud, const std::vector< std::vector< int > > &nList, bool isSlice=false, int coordinationNumber=4)
 Function for getting the bond order correlations \(c_{ij}\) (or \(a_{ij}\) in some treatments) according to the CHILL algorithm.
void getIceTypeNoPrint (molSys::PointCloud< molSys::Point< double >, double > &yCloud, const std::vector< std::vector< int > > &nList, bool isSlice=false)
 Function that classifies every particle's molSys::atom_state_type ice type, according to the CHILL algorithm.
void getIceType (molSys::PointCloud< molSys::Point< double >, double > &yCloud, const std::vector< std::vector< int > > &nList, std::string path, int firstFrame, bool isSlice=false, std::string outputFileName="chill.txt")
 Function that classifies every particle's molSys::atom_state_type ice type, according to the CHILL algorithm.
void getCorrelPlus (molSys::PointCloud< molSys::Point< double >, double > &yCloud, const std::vector< std::vector< int > > &nList, bool isSlice=false, int coordinationNumber=4)
void getIceTypePlus (molSys::PointCloud< molSys::Point< double >, double > &yCloud, const std::vector< std::vector< int > > &nList, std::string path, int firstFrame, bool isSlice=false, std::string outputFileName="chillPlus.txt")
 Classifies each atom according to the CHILL+ algorithm.
void getIceTypePlusNoPrint (molSys::PointCloud< molSys::Point< double >, double > &yCloud, const std::vector< std::vector< int > > &nList, bool isSlice=false)
 CHILL+ ice types on the cloud. Does not write a file.
std::vector< double > getq6 (molSys::PointCloud< molSys::Point< double >, double > &yCloud, const std::vector< std::vector< int > > &nList, bool isSlice=false)
void reclassifyWater (molSys::PointCloud< molSys::Point< double >, double > &yCloud, std::vector< double > &q6)
int printIceType (molSys::PointCloud< molSys::Point< double >, double > &yCloud, std::string path, int firstFrame, bool isSlice=false, std::string outputFileName="superChill.txt")
 Prints out the iceType for a particular frame onto the terminal.
bool isInterfacial (molSys::PointCloud< molSys::Point< double >, double > &yCloud, const std::vector< std::vector< int > > &nList, int iatom, int num_staggrd, int num_eclipsd, bool chillPlus=false)
int numStaggered (molSys::PointCloud< molSys::Point< double >, double > &yCloud, const std::vector< std::vector< int > > &nList, int jatom)
 Finds the number of staggered bonds for a given atom of index jatom.
SteinhardtQl steinhardtQl (const molSys::PointCloud< molSys::Point< double >, double > &yCloud, const std::vector< std::vector< int > > &nList, int orderL)
std::vector< TemplateHit > classifyTemplates (const molSys::PointCloud< molSys::Point< double >, double > &yCloud, const std::vector< std::vector< int > > &nList, int kNeigh=12)
std::vector< double > soapSpectrum (const molSys::PointCloud< molSys::Point< double >, double > &yCloud, int iatom, const std::vector< std::vector< int > > &nList, int nMax, int lMax, double rcut)
std::vector< std::vector< double > > soapSpectrumAll (const molSys::PointCloud< molSys::Point< double >, double > &yCloud, const std::vector< std::vector< int > > &nList, int nMax, int lMax, double rcut)
std::vector< double > voronoiFeature (const molSys::PointCloud< molSys::Point< double >, double > &yCloud, int iatom, double candidateCutoff)
 [q4, q6, q8] from the Voronoi-weighted Steinhardt path (Mickel).
std::vector< std::vector< double > > voronoiFeatures (const molSys::PointCloud< molSys::Point< double >, double > &yCloud, double candidateCutoff)
 [q4, q6, q8] for every particle. One tessellation per order (l = 4, 6, 8).
std::vector< VoronoiWeights > voronoiFacetWeights (const molSys::PointCloud< molSys::Point< double >, double > &yCloud, double candidateCutoff)
 Facet neighbours and area weights for every particle, with an a posteriori exactness certificate.
SteinhardtQl steinhardtQlVoronoi (const molSys::PointCloud< molSys::Point< double >, double > &yCloud, double candidateCutoff, int orderL)

Detailed Description

CHILL and CHILL+ structure classification.

This namespace contains functions that are used in the CHILL/CHILL+ classification scheme, as well as a yodaCloud struct to hold all the information.

Both the CHILL and CHILL+ methods are based on the local bond order parameter method, developed by ten Wolde et al. for the identification of crystal nuclei in Lennard-Jones systems. The local order parameter method is based on order parameters introduced by Steinhardt et al.

The local environment of the water molecules is classified using an algorithm based on spherical harmonics, which is independent of the specific crystal structure and does not require the definition of a reference frame. The local order around each water molecule \(i\) is described by a local orientational bond order parameter \(q_{lm}(i)\), with \( 2l+1 \) complex components:

\[ q_{lm}(i) = @frac{1}{N_b(i)} \Sigma_{j=1}^{N_b(i)} Y_{lm}(r_{ij}) \]

Here, \(N_b(i)=4\) is the number of nearest neighbours for water molecule \(i\), \(l\) is a free integer parameter, \(m\) is an integer such that \(m \in [-l,l]\). The functions \(Y_{lm}(r_{ij})\) are the spherical harmonics and \(r_{ij}\) is the distance vector between molecule \(i\) and each one of its four nearest neighbours \(j\).

The alignment of the orientation of the local structures is measured by the normalized dot product of the local orientational bond order parameter, \(q_l(i)\), between each pair of neighbor molecules \(i\) and \(j\), given by:

\[ a(i,j) = \frac{\sum_{m=-l}^{l} q_{lm}(i)\, q_{lm}^*(j)} {|q_l(i)|\,|q_l(j)|} \]

where, \(q_{lm}^*\) is the complex conjugate of \(q_{lm}\).

Depending on the values of \(a(i,j)\), the bond between each pair of molecules \(i\) and \(j\) is classified as being either staggered or eclipsed.

Algorithm Eclipsed Bonds Staggered Bonds
CHILL \(-0.05 >= a(i,j) >= -0.2\) \(a(i,j) <= -0.8\)
CHILL+ \(0.25 >= a(i,j) >= -0.35\) \(a(i,j) <= -0.8\)

Since each molecule \(i\) in deeply supercooled water has four nearest neighbours, the type of the resultant four bonds is used to identify phases according to the CHILL and CHILL+ algorithms.

The CHILL algorithm can classify water molecules as belonging to the cubic, hexagonal, interfacial or liquid (amorphous) phase:

| Phase | E | S | Neighbours | Description | |----------—|--—|--—|---------—|------------------------------------------------------—| | Cubic | 0 | 4 | 4 | All 4 neighbours are staggered. | | Hexagonal | 1 | 3 | 4 | 3 staggered and 1 eclipsed bond. | | Interfacial | any | 2 | 4 | Belongs to the ice phase, but does not satisfy strict criteria of cubic and hexagonal | | | 0 | 3 | 4 | | | Liquid | N/A | N/A | any | Classification by exclusion. |

Here, E and S refer to eclipsed and staggered bonds. The neighbours column contains the number of nearest neighbours.

The CHILL+ algorithm, which is a modified version of CHILL, additionally identifies clathrate and interfacial clathrate phases. The criteria are enumerated in the table below:

Phase E S Neighbours Description
Cubic 0 4 4 no change was made to

staggered bond criterion from CHILL | | Hexagonal | 1 | 3 | 4 | wider eclipsed range compared to CHILL; identifies 99% hexagonal ice up to 270 K | | Interfacial | any | 2 | 4 | must have at least one first neighbor water with more than two staggered bonds | | | 0 | 3 | 4 | must have at least one first neighbor water with more than one staggered bond | | Clathrate | 4 | 0 | 4 | bulk clathrate and part of the interface of clathrates can be identified with four eclipsed bonds | | Interfacial clathrate | 3 | any | 4 | partial clathrate cages and threads of clathrate-like order | | Liquid | N/A | N/A | any | classifies as liquid if none of the above criteria are fulfilled |

Although both the CHILL and CHILL+ classification schemes take into account the local order of the environment of each particle, both schemes output a per-particle classification.

Changelog

Enumeration Type Documentation

◆ CrystalKind

enum class chill::CrystalKind
strong
Enumerator
other 
sc 
fcc 
hcp 
bcc 

Definition at line 21 of file structure_desc.hpp.

21{ other = 0, sc, fcc, hcp, bcc };

Function Documentation

◆ classifyTemplates()

std::vector< TemplateHit > chill::classifyTemplates ( const molSys::PointCloud< molSys::Point< double >, double > & yCloud,
const std::vector< std::vector< int > > & nList,
int kNeigh = 12 )
nodiscard

Overlay the k nearest neighbours of each particle onto FCC, HCP, BCC and SC shells. Uses IRA when linked; otherwise Horn on a distance-sorted correspondence (high-symmetry shells only).

◆ soapSpectrum()

std::vector< double > chill::soapSpectrum ( const molSys::PointCloud< molSys::Point< double >, double > & yCloud,
int iatom,
const std::vector< std::vector< int > > & nList,
int nMax,
int lMax,
double rcut )
nodiscard

Bartok SOAP power spectrum for one particle (nMax radial Gaussians, spherical harmonics through lMax). Length is nMax*nMax*(lMax+1).

◆ soapSpectrumAll()

std::vector< std::vector< double > > chill::soapSpectrumAll ( const molSys::PointCloud< molSys::Point< double >, double > & yCloud,
const std::vector< std::vector< int > > & nList,
int nMax,
int lMax,
double rcut )
nodiscard

Bartok SOAP power spectrum for every particle. Length is nop; each row is nMax*nMax*(lMax+1).

◆ steinhardtQlVoronoi()

SteinhardtQl chill::steinhardtQlVoronoi ( const molSys::PointCloud< molSys::Point< double >, double > & yCloud,
double candidateCutoff,
int orderL )
nodiscard

Steinhardt parameters of degree orderL (3, 4, 6 or 8) with the bond sum weighted by Voronoi facet areas. qlBar averages the weighted q_lm over the particle and its facet neighbours, after Lechner and Dellago.

◆ voronoiFacetWeights()

std::vector< VoronoiWeights > chill::voronoiFacetWeights ( const molSys::PointCloud< molSys::Point< double >, double > & yCloud,
double candidateCutoff )
nodiscard

Facet neighbours and area weights for every particle, with an a posteriori exactness certificate.

candidateCutoff seeds the bisector search; the cutoff is enlarged automatically until the certificate holds or the growth cap is reached.

Certificate. The cell of particle i clipped against every candidate within cutoff c is the true Voronoi cell whenever every vertex of the clipped cell lies within c/2 of i: a particle k beyond the cutoff has |d_k| > c, its bisector plane sits at distance |d_k|/2 > c/2 from i, and a half-space whose boundary is farther than c/2 cannot cut a region contained in the ball of radius c/2. An open or under-clipped cell keeps vertices on the seeding square at distance ~c and fails the certificate, which is exactly the failure the enlargement retries.

◆ voronoiFeature()

std::vector< double > chill::voronoiFeature ( const molSys::PointCloud< molSys::Point< double >, double > & yCloud,
int iatom,
double candidateCutoff )
nodiscard

[q4, q6, q8] from the Voronoi-weighted Steinhardt path (Mickel).

◆ voronoiFeatures()

std::vector< std::vector< double > > chill::voronoiFeatures ( const molSys::PointCloud< molSys::Point< double >, double > & yCloud,
double candidateCutoff )
nodiscard

[q4, q6, q8] for every particle. One tessellation per order (l = 4, 6, 8).