Bond Orientational Order Parameters at 3D

The class static.boo.boo_3d() calculates the bond orientational order (BOO) parameters in three-dimensions. All calculations start from measuring the \(l\)-fold symmetry BOO as a \(2l+1\) complex vector \(q_{lm}\):

\[ q_{{lm}}(i) = \frac{1}{N_i} \sum_i^{N_i} Y_{{lm}} (\theta(\vec r_{ij}), \phi(\vec r_{ij})) \tag{1} \]

where \(Y_{lm}\) are the spherical harmonics and \(N_i\) the number of bonds of particles \(i\). This is the original one with all (\(N_i\)) neighboring particles treated the same. However, in some cases, we need to treat these neighboring particles differently based on their properties with respect to the center particle \(i\). Then we can obtain the weighted \(2l+1\) complex vector \(q_{lm}\) as

\[ q_{{lm}}(i) = \sum_i^{N_i} \frac{A_j}{\sum A_j} Y_{{lm}} (\theta(\vec r_{ij}), \phi(\vec r_{ij})) \tag{2} \]

where \(A_j\) is the property of the \(j-\)th neighbor. Therefore, the normalized weights \(\frac{A_j}{\sum A_j}\) are used in the calculation. We can treat \(A_j=1\) for the case shown in equation (1). The weights are provided in the same way as for the neighbors.

Furthermore, dealing with \(q_{lm}\) by coarse-graining is sometimes helpful by

\[ Q_{lm}(i)=\frac {1}{N_i+1} \left(q_{lm}(i) + \sum_{j=1} ^{N_i} q_{lm}(j) \right) \tag{3} \]

Therefore, we can define the local complex vectors in a flexible way for different purposes. In principal, we name below the local quantities with small letters (associate with \(q_{lm}\)) and coarse-grained ones with large letters (associate with \(Q_{lm}\)). In the calculations one can set coarse_graining=True to calculate for the latter case.

By either using the complex vetcor \(q_{lm}\) or \(Q_{lm}\), various BOO parameters can be calculated. For example,

  • the rotational invariants q_l (local, \(q_l\)) and Q_l (coarse-grained, \(Q_l\)) by boo_3d.ql_Ql(),

\[ q_l=\sqrt{\frac{4\pi}{2l+1} \sum_{m=-l}^{l} \left| q_{lm} \right|^2} \tag{4} \]

  • local coherence of the complex vectors sij for \(q_{lm}\) and Sij for \(Q_{lm}\) by boo_3d.sij_ql_Ql(),

\[ s(i,j)=\frac {\sum_{m=-l}^l q_{lm}(i) q_{lm}(j)}{\sqrt{\sum_{m=-l}^l \left| q_{lm}(i) \right|^2} \sqrt{\sum_{m=-l}^l \left| q_{lm}(j) \right|^2}} \tag{5} \]

  • w_l/w_cap_l (local, \(w_l\), \(\hat w_l\)) and W_l/W_cap_l (coarse-grained, \(W_l\), \(\hat W_l\)) by boo_3d.w_W_cap()

\[\begin{split} w_{l} = \sum_{m_{1} + m_{2} + m_{3} = 0}^{} \left( \begin{matrix} l & l & l \\ m_{1} & m_{2} & m_{3} \ \end{matrix} \right) q_{lm_{1}}q_{lm_{2}}q_{lm_{3}} \end{split}\]
\[ \hat{w_l}=w_l \left(\sum_{m=-l}^l \left| q_{lm} \right|^2 \right) ^{-\frac {3}{2}} \tag{7} \]

  • spatial correlation of the complex vectors by boo_3d.spatial_corr()

\[ G_l(r)=\frac{4\pi}{2l+1} \frac{\sum_{ij} \sum_{m=-l}^l Q_{lm}(i) Q_{lm}^*(j) \delta(\vec r_{ij} - \vec r)} {\sum_{ij} \delta(\vec r_{ij} - \vec r)} \tag{8} \]

  • time correlation of the complex vectors by boo_3d.time_corr()

\[ C_l(t)=\frac{4\pi}{2l+1} \frac {\left< \sum_i^N \sum_{m=-l}^l Q_{lm}(i, t) Q_{lm}^*(i, 0) \right>} {\left< \sum_i^N \sum_{m=-l}^l \left| Q_{lm}(i, 0) \right|^2 \right>} \tag{9} \]

In equations (4)-(9), one can replace \(q_{lm}\) by \(Q_{lm}\) to calculate the corresponding coarse-grained values, or vice versa.

One example of using \(s(i,j)\) with \(q_{lm}\) for \(l=6\) in computer simulations is to detect the bcc crystal nuclei. The bond between the center particle \(i\) and a neighbor \(j\) is considered crystalline if \(s(i,j)>c\) (such as c=0.7). A particle is crystalline if the number of crystalline bond is larger than a threshold (usually half of the coordination number).

In these calculations, one important input is the neighbor list of each particle. Currently, the module accepts data generated by the neighbors module. One can also reweight BOO (see equation (2)) by using Voronoi polyhedron face area. Other weighting methods are acceptable as long as the correct format to account for the center-neighbor property is given, similar to the neighboring particles.

Note that this module accounts for both orthogonal and triclinic cells.

1. boo_3d class

Input Arguments

  • snapshots (reader.reader_utils.Snapshots): snapshot object of input trajectory (returned by reader.dump_reader.DumpReader)

  • l (int): degree of spherical harmonics

  • neighborfile (str): file name of particle neighbors (see module neighbors)

  • weightsfile (str): file name of particle-neighbor weights (see module neighbors), default None. one typical example is Voronoi face area of the polyhedron; this file should be consistent with neighborfile. If provided, weightsfile will be read to calculate the weighted BOO, otherwise calculate the original BOO.

  • ppp (npt.NDArray): the periodic boundary conditions, setting 1 for yes and 0 for no, default np.array([1,1,1]),

  • Nmax (int): maximum number for neighbors, default 30

Return:

  • None

Example

from PyMatterSim.reader.dump_reader import DumpReader
from PyMatterSim.reader.reader_utils import DumpFileType
from PyMatterSim.neighbors.voropp_neighbors import cal_voro
from PyMatterSim.static.boo import boo_3d

filename = 'dump.atom'
readdump = DumpReader(filename, ndim=3, filetype=DumpFileType.LAMMPS, moltypes=None)
readdump.read_onefile()

cal_voro(readdump.snapshots, outputfile='test')
boo = boo_3d(readdump.snapshots,
             l=6,
             neighborfile='test.neighbor.dat'
             weightsfile='test.facearea.dat')

2. qlm_Qlm()

qlm_Qlm() method gives \(q_{lm}\) and \(Q_{lm}\) values of different particles in complex vectors as numpy array. The results in different snapshots are stored in an array for each of them, separately. Both of them are returned as (\(q_{lm}\), \(Q_{lm}\)). Note that this method will be called automatically when initiallizing the boo_3d() class.

Input Arguments

  • None

Return

  • smallqlm (npt.NDArray): \(q_{lm}\) in vector complex number, in numpy array with shape [nsnapshots, nparticle, 2l+1]

  • largeQlm (npt.NDArray): \(Q_{lm}\) in vector complex number, in numpy array with shape [nsnapshots, nparticle, 2l+1]

Example

qlm, Qlm = boo.qlm_Qlm()

3. ql_Ql()

ql_Ql() method calculates \(q_l\) (local) or \(Q_l\) (coarse-grained).

Input Arguments

  • coarse_graining (bool): whether use coarse-grained \(Q_{lm}\) or local \(q_{lm}\), default False

  • outputfile (str): file name for \(q_l\) or \(Q_l\) results, default None. To reduce storage size and ensure loading speed, save npy file as default with extension “.npy”. If the file extension is “.dat” or “.txt”, also saved a text file.

Return

  • calculated \(q_l\) or \(Q_l\) in npt.NDArray with shape [nsnapshots, nparticle]

Example

Ql = boo.ql_Ql(coarse_graining=True, outputfile='./results/Ql.npy')

4. sij_ql_Ql()

sij_ql_Ql() method calculates the orientation correlation of \(q_{lm}\) or \(Q_{lm}\), named as \(s_{ij}\) or \(S_{ij}\).

Input Arguments

  • coarse_graining (bool): whether use coarse-grained \(Q_{lm}\) or local \(q_{lm}\), default False

  • c (float): cutoff defining bond property, such as solid or not, default 0.7

  • outputqlQl (str): npy file name of \(q_l\) or \(Q_l\) from particle-level sij, default None; The header is “id sum_sij num_neighbors”, available from logger

  • outputsij (str): npy file name for original \(s_{ij}\) of \(q_l\) or \(Q_l\), default None; The header is “id CN sij”, available from logger

Return

  • aggregated \(s_{ij}\) or \(S_{ij}\) to particle-level in npt.NDArray with shape [nsnapshots*nparticle, 3]. This result is also stored in outputqlQl file for each particle.

  • calculated \(s_{ij}\) or \(S_{ij}\) in npt.NDArray with shape [nsnapshots*nparticle, 2+max_neighbors], where max_neighbors is the maximum number of neighbors and the first two columns are particle index and coordination number. The returned result is also stored in outputsij file for each pair

In outputqlQl file, the first column is particle index, the second is the sum of \(s_{ij}\), and third one is number of neighbors.

outputsij file stores the \(s_{ij}\) value of each particle with each of its neighbors. Results in different snapshots are written in sequence, with one header “id cn sij”.

Example

sij_Ql = boo.sij_ql_Ql(coarse_graining=True, 
                       outputqlQl='sum_sij.npy',
                       outputsij='sij_Ql.npy')

5. w_W_cap()

w_W_cap() calculates the Wigner 3-j symbol boo based on qlm or Qlm.

Input Arguments

  • coarse_graining (bool): whether use coarse-grained \(Q_{lm}\) or local \(q_{lm}\), default False

  • outputw (str): txt file name for w (original) based on \(q_{lm}\) or \(Q_{lm}\), default None

  • outputwcap (str): file name for wcap (normalized) based on qlm or Qlm, default None. To reduce storage size and ensure loading speed, save npy file as default with extension “.npy”. If the file extension is “.dat” or “.txt”, also saved a text file.

Return

  • w_W (np.adarray), calculated \(w\) or \(W\) with shape [nsnapshots, nparticle]

  • w_W_cap (np.adarray), calculated \(w_{cap}\) or \(W_{cap}\) with shape [nsnapshots, nparticle]

Example

W, W_cap = boo.w_W_cap(coarse_graining=True, 
                       outputw='W.dat',
                       outputwcap='W_cap.npy')

6. spatial_corr()

spatial_corr() method calculates spatial correlation function of \(q_{lm}\) or \(Q_{lm}\)

Input Arguments

  • coarse_graining (bool): whether use coarse-grained \(Q_{lm}\) or local \(q_{lm}\), default False

  • rdelta (float): bin size in calculating g(r) and Gl(r), default 0.01

  • outputfile (str): csv file name for \(g_l\), default None

Example

Gl_Q = boo.spatial_corr(coarse_graining=True, outputfile='Gl_Q.csv')

7. time_corr()

time_corr() method calculates time correlation of \(q_{lm}\) or \(Q_{lm}\)

Input Arguments

  • coarse_graining (bool): whether use coarse-grained \(Q_{lm}\) or local \(q_{lm}\), default False

  • dt (float): timestep used in user simulations, default 0.002

  • outputfile (str): csv file name for time correlation results, default None

Example

Gl_time = boo.time_corr(coarse_graining=True, outputfile='Gl_time.csv')