SolOrder

A module to calculate various structural order parameter related to water molecules. This project initiated with the aim to distinguish between different polymorphs of ice.
Explore the docs »

View Demo · Report Bug · Request Feature

Table of Contents

1. About The Project
2. Getting Started
   • Prerequisites
   • Installation
3. Usage
4. Contributing
5. License
6. Contact
7. Acknowledgments

About The Project

SolOrder is purposed to calculate simple structural property based order parameter of water molecules in MD simulation output. This module is validated for gromacs output files.
Currently it can calculate Four body structural order parameter ( $F_4$ ), Orientational and translational tetrahedral order parameter(OTO & TTO), Local structure index(LSI), Local entropy fingerprints etc by follwing the protocols specified in the literature so far. There are a lot more order parameters related to water structure and we wish to have them under one roof. Our current interest is to include bond orientational order parameters and SOAP(Smooth Overlap of Atomic Positions) related order paramter to the SolOrder module. There's a lot we can actually do to incrase number of ways to accurately detect phases, of which immediates are invoking ideas of molecular graph and ML.

(back to top)

Getting Started

Prerequisites

This module uses few optimized modules such as MDAnalysis, pandas, numba, numpy etc. to fast-track the calculation of OPs. They are listed in requirements.txt file. You can install them using either pip or conda. Like below...

• MDAnalysis

  pip install MDAnalysis

Installation

You may not need to install the module. Just clone the repository and import them locally. To install the SolOrder module globally you can do either of the following ways .

1. Clone the repo

   git clone https://github.com/dmighty007/SolOrder.git
   cd SolOrder
   python setup.py install

2. Or you can install directly from pypi store..

   pip install SolOrder

(back to top)

Usage

For more examples, please refer to the Documentation

######## Import two classes from the module...
from SolOrder import OrderParameter as OP
from SolOrder import eFingerprint as E

########## Basic modules....
import numpy as np
import MDAnalysis as mda
import seaborn as sns
import matplotlib.pyplot as plt

########## The MDAnalysis universe object...
u = mda.Universe("sH.gro")
u

<Universe with 7344 atoms>

Dihedral angle distribution and other order parameters

######## Create the object..
obj = OP(u)

########## Dihedral angle distribution...
obj.get_dihedrals()
sns.kdeplot(obj.dihedrals, bw_method = 0.05)

F4 order parameter is a mathematical function of the total number of water molecules and the torsion angle between the O-H bond of two water molecules, such that oxygen atoms of both water molecules have a distance that is less than 0.3 nm. This particular parameter can discriminate between bulk liquid water, ice, and sI crystal. It has a value of -0.4, -0.04, and 0.7 for bulk liquid water, regular ice, and sI hydrate crystals respectively.It is given as..

$$F_4 = \frac{1}{N}\sum(cos3\phi)$$ The torsion angle $\phi$ is defined as H-O...O-H, for two adjacent water molecules, where the hydrogens have been defined as the outer-most hydrogens in the water dimer.

obj.F4()
sns.kdeplot(obj.f4, bw_method = 0.05)

The below is the most common type of tetrahedral order paramater that is being used. This order parameter uses the four closest water Oxygen neighbors for the calculation. The value of q can rang from 0 to 1, where 0 is for an ideal gas and 1 is for a regular tetrahedron. $$q = 1− \frac{3}{8}\sum_{j=1}^3\sum_{k=j+1}^4(cos\phi_{jk}+\frac{1}{3})^2$$

$q =$ Orientational tetrahedral order parameter.
$\phi_{jk}=$ Angle formed by the Oxygen atom under consideration & the two nearest neighbor atoms j & k.

obj.OTO()
sns.kdeplot(obj.tetra_orient,bw_method = 0.05)

Similar to orientational tetrahedral order. But here, the variance of the radial distance between the central water oxygen atom and the four nearest neighbors are calculated. THe value of $S_k$ is close to 1 and equals to 1 for the perfect tetrahedron. As the local tetrahedral nature increases, $S_k$ becomes more closer to 1.

$$ S_k = 1−\frac{1}{3}\sum_{k=1}^4\frac{(r_k−\bar{r})^2}{4\bar{r}^2}$$

$S_k$ = Translational tetrahedral order parameter

$r_k$ = Radial distance from the cental Oxygen atom to the k th peripheral closest neighbor

$\bar{r}$ = Arithmatic mean of the four radial distances

obj.TTO()
sns.kdeplot(obj.tetra_trans,bw_method = 0.05)

The local-structure-index (LSI) is defined for each molecule i by ordering the nearest neighbors j according to increasing distance to reference molecule i as $r_1 < r_2 < r_3 < ... < r_{n(i)} < 0.37 nm < r_{n(i)+1}$ where $n(i)$ is the number of molecules that are within $0.37 nm $ from molecule i (oxygen atom positions are used). The LSI distinguishes molecules with well separated first and second coordination shells from molecules with disordered environment, containing molecules in interstitial positions, through the parameter $I(i)$ defined by

$$ I(i) = \frac{1}{n(i)} \sum_{j=1}^{n(i)}[\Delta(j;i) - \bar{\Delta}(i)]^2 $$

Here $\Delta(j;i) = r_{j+1}- r_j$ and $\bar{\Delta}(i)$ is the average of $\Delta(j;i)$ over all neighbors j of molecule i within the cutoff. A low LSI corresponds to a disordered local environment (HDL) while a high LSI indicates a highly structured, tetrahedral coordination (LDL).

obj.LSI()
sns.kdeplot(obj.lsi,bw_method = 0.05)

obj.MinimumAngle()
sns.kdeplot(np.cos(obj.minAngles),bw_method = 0.05)

Entropy Fingerprints

Entropy is a global property, hence, to define alocal order parameter, it is projected onto each particle i as :

$$ e(i) = -2 \pi \rho k_b \int_0^{r_m} [{g_m^i(r)}\ln(g_m^i(r)) - g_m^i(r) + 1 ] r^2 dr $$

where $g_m^i(r)$ is the radial distribution function centred on particle i, which can be smoothed with Gaussians of finite width (here we choose 0.15 $\unicode{x212B}$) for obtaining a continuous and differentiable order parameter.In eqn $r_m$ represents an upper integration limit which, in a rigorous definition of entropy, should approach infinity; here $r_m$ is set to 5.0 $\unicode{x212B}$.\ The modifier function stores the computed pair entropy values, which are always negative, in the output particle property Entropy. Lower entropy values correspond to more ordered structural environments. The modifier function does not take into account the types of the input particles; it assumes the model is a single-component system. To improve the resolution between different environments,a locally averaged entropy, le, can be defined as:

$$ le(i) = \frac{1}{(N+1)} \sum_{i=1}^N e(i) $$

entropy = E(u)

entropy.Entropy()

  0%|          | 0/1836 [00:00<?, ?it/s]

entropy.localEntropy()

  0%|          | 0/1836 [00:00<?, ?it/s]

sns.kdeplot(entropy.le, bw_method = 0.05)

(back to top)

Contributing

Contributions are what make the open source community such an amazing place to learn, inspire, and create. Any contributions you make are greatly appreciated.

If you have a suggestion that would make this better, please fork the repo and create a pull request. You can also simply open an issue with the tag "enhancement". Don't forget to give the project a star! Thanks again!

1. Fork the Project
2. Create your Feature Branch (git checkout -b feature/AmazingFeature)
3. Commit your Changes (git commit -m 'Add some AmazingFeature')
4. Push to the Branch (git push origin feature/AmazingFeature)
5. Open a Pull Request

(back to top)

License

Distributed under the MIT License. See LICENSE for more information.

(back to top)

Contact

Dibyendu Maity - @MAITYDIBYENDU99 - dibyendumaity1999@bose.res.in

Project Link: https://github.com/dmighty007/SolOrder

(back to top)

Acknowledgments

Below are the references, basically from where the ideas and algorithms are taken. We might want to benchmark with their results!!

• Investigating the quasi-liquid layer on ice surfaces: a comparison of order parameters
• Mining of Effective Local Order Parameters to Classify Ice Polymorphs
• Characterization of the Local Structure in Liquid Water by Various Order Parameters
• Simulations of the methane hydrate/methane gas interface near hydrate forming conditions conditions

(back to top)