Clustering Algorithms to Analyze Molecular Dynamics Simulation Trajectories for Complex Chemical and Biological Systems

Molecular dynamics （MD） simulation has become a powerful tool to investigate the structure- function relationship of proteins and other biological macromolecules at atomic resolution and biologically relevant timescales. MD simulations often produce massive datasets con- taining millions of snapshots describing proteins in motion. Therefore, clustering algorithms have been in high demand to be developed and applied to classify these MD snapshots and gain biological insights. There mainly exist two categories of clustering algorithms that aim to group protein conformations into clusters based on the similarity of their shape （geometric clustering） and kinetics （kinetic clustering）. In this paper, we review a series of frequently used clustering algorithms applied in MD simulations, including divisive algorithms, ag- glomerative algorithms （single-linkage, complete-linkage, average-linkage, centroid-linkage and ward-linkage）, center-based algorithms （K-Means, K-Medoids, K-Centers, and APM）, density-based algorithms （neighbor-based, DBSCAN, density-peaks, and Robust-DB）, and spectral-based algorithms （PCCA and PCCA＋）. In particular, differences between geomet- ric and kinetic clustering metrics will be discussed along with the performances of diflhrent clustering algorithms. We note that there does not exist a one-size-fits-all algorithm in the classification of MD datasets. For a specific application, the right choice of clustering algo- rithm should be based on the purpose of clustering, and the intrinsic properties of the MD conformational ensembles. Therefore, a main focus of our review is to describe the merits and limitations of each clustering algorithm. We expect that this review would be helpful to guide researchers to choose appropriate clustering algorithms for their own MD datasets.

Molecular dynamics simulation of an argon cluster filled inside carbon nanotubes

The effects of the diameters of single-walled carbon nanotubes （SWCNTs） （7.83A to 27.40A） and temperature （20 K-45 K） on the equilibrium structure of an argon cluster are systematically studied by molecular dynamics simulation with consideration of the SWCNTs to be fixed. Since the diameters of SWCNTs with different chiralities increase when temperature is fixed at 20 K, the equilibrium structures of the argon cluster transform from monoatomic chains to helical and then to multishell coaxial cylinders. Chirality has almost no noticeable influence on these cylindrosymmetric structures. The effects of temperature and a non-equilibrium sudden heating process on the structures of argon clusters in SWCNTs are also studied by molecular dynamics simulation.

Coalescence between Cu_(57) and Cu_(58) clusters at a room temperature:molecular dynamics simulations

Three coalescence processes of Cu57-Cu57, Cu57-Cu58, and Cu5s Cu58 clusters at 300 K are investigated by employing molecular dynamics simulations. According to the evolutions of mean square displacement and local atom packing, the coalescence process can be separated into three stages including an approaching stage, a coalescing stage, and a coalesced stage. The simulations show that the coalescence processes and the formed products are sensitive to the respective initial structures of, and the relative configuration between, the two coalescing icosahedron-based clusters.

MOLECULAR DYNAMICS SIMULATIONS OF FILLED AND EMPTY CAGE-LIKE WATER CLUSTERS IN LIQUID WATER AND THEIR SIGNIFICANCE TO GAS HYDRATE FORMATION MECHANISMS

Molecular dynamics simulations are performed to observe the evolutions of 512 and 51262 cage-like water clusters filled with or without a methane molecule immersed in bulk liquid water at 250 K and 230 K. The lifetimes of these clusters are calculated according to their Lindemann index δ (t) using the criteria of δ≥0.07. For both the filled and empty clusters, we find the dynamics of bulk water determines the lifetimes of cage-like water clusters, and that the lifetime of 512 62 cage-like cluster is the same as that of 512 cage-like cluster. Although the methane molecule indeed makes the filled cage-like cluster more stable than the empty one, the empty cage-like cluster still has chance to be long-lived compared with the filled clusters. These observations support the labile cluster hypothesis on the formation mechanisms of gas hydrates.

Solid−solid phase transition of tungsten induced by high pressure:A molecular dynamics simulation

The phase transition of tungsten(W)under high pressures was investigated with molecular dynamics simulation.The structure was characterized in terms of the pair distribution function and the largest standard cluster analysis(LSCA).It is found that under 40−100 GPa at a cooling rate of 0.1 K/ps a pure W melt first crystallizes into the body-centred cubic(BCC)crystal,and then transfers into the hexagonal close-packed(HCP)crystal through a series of BCC−HCP coexisting states.The dynamic factors may induce intermediate stages during the liquid−solid transition and the criss-cross grain boundaries cause lots of indistinguishable intermediate states,making the first-order BCC−HCP transition appear to be continuous.Furthermore,LSCA is shown to be a parameter-free method that can effectively analyze both ordered and disordered structures.Therefore,LSCA can detect more details about the evolution of the structure in such structure transition processes with rich intermediate structures.

Molecular dynamics simulations of cascade damage near the Y_2Ti_2O_7 nanocluster/ferrite interface in nanostructured ferritic alloys

A comparative study of cascades in nanostructured ferritic alloys and pure Fe is performed to reveal the influence of Y_2Ti_2O_7 nanocluster on cascades by molecular dynamics simulations. The cascades with energies of primary knock-on atom（PKA） ranging from 0.5 keV to 4.0 keV and PKA＇s distances to the interface from 0 to 50 are simulated. It turns out that the Y_2Ti_2O_7 nanocluster can absorb the kinetic energy of cascade atoms, prevent the cascade from extending and reduce the defect production significantly when the cascades overlap with the nanocluster. When the cascade affects seriously the nanocluster, the weak sub-cascade collisions are rebounded by the nanocluster and thus leave more interstitials in the matrix. On the contrary, when the cascade contacts weakly the nanocluster, the interface can pin the arrived interstitials and this leaves more vacancies in the matrix. Moreover, the results indicate that the Y_2Ti_2O_7 nanocluster keeps stable upon the displacement cascade damage.

Molecular Dynamics Studies of Solid-Liquid Phase Transitions Ⅱ——Homogeneous Nucleation of (KI)_(108) Clusters

在一系列分子动力学计算中，（ＫＩ）１０８离子簇冷却时自发地由液相变为立方面心相．在３５０Ｋ和４００Ｋ的成核速率大于１０３６核／ｍ３·ｓ．由成核速率得到的固－液界面自由能约为５０ｍＪ／ｍ２．

Theoretical insights into oxygen reduction reaction on Au-based single-atom alloy cluster catalysts

Developing highly active alloy catalysts that surpass the performance of platinum group metals in the oxygen reduction reaction(ORR)is critical in electrocatalysis.Gold-based single-atom alloy(AuSAA)clusters are gaining recognition as promising alternatives due to their potential for high activity.However,enhancing its activity of AuSAA clusters remains challenging due to limited insights into its actual active site in alkaline environments.Herein,we studied a variety of Au_(54)M_(1) SAA cluster catalysts and revealed the operando formed MO_(x)(OH)_(y) complex acts as the crucial active site for catalyzing the ORR under the basic solution condition.The observed volcano plot indicates that Au_(54)Co_(1),Au_(54)M_(1),and Au_(54)Ru_(1) clusters can be the optimal Au_(54)M_(1) SAA cluster catalysts for the ORR.Our findings offer new insights into the actual active sites of AuSAA cluster catalysts,which will inform rational catalyst design in experimental settings.

Clusters Deposition on Surface an Atomic Scale Study by Computer Simulation Method

The investigation is generalized to clusters with sizes up to 3000 atoms, covering this way the range of sizes experimentally available for low energy cluster beam deposition. The atomic scale modeling is carried on by both Molecular Dynamics and Metropolis Monte Carlo. This represents a huge series of simulations (175 cases) to which further calculations are added by spot when finer tuning of the parameters is necessary. Analyzing the results is a major task which is still in progress. This way, not only a realistic range of sizes is covered, but also the whole range of compositions and the temperature range relevant to the solid and the liquid states.

First-Principles Calculations for Structures and Melting Temperature of Si6 Clusters

We investigate the structures and the melting temperature of the Si6 cluster by using the first-principles pseudopotential method in real space and Langevin molecular dynamics. It is shown that the ground structure of the Si6 cluster is a square bipyramid, and the corresponding melting temperature is about 1923 K. In the heating procedure, the structures of the Si6 cluster change from high symmetry structures containing 5-8 bonds, via prolate structures containing 3-4 bonds, to oblate structures containing 1-2 bonds.

Simulation of molecular dynamics of silver subcritical nuclei and crystal clusters during solidification

Molecular dynamics simulation is carried out to investigate the effects of cooling rate on the final configurations of silver after rapid solidification. The cooling rate for the formation of a silver amorphous phase is determined by analyzing its pair distribution function, H-A bond index, and the largest crystal cluster. Further, the equilibrium structures of the subcritical nuclei and crystal clusters are studied. The results show that the solidified microstructure is composed of a mixture of crystal clusters and amorphous phases at a certain cooling rate range. The size of the largest crystal cluster decreases with the increasing cooling rate, and it completely disappears when the cooling rate exceeds a critical value. The structures of the subcritical nuclei and the largest crystal cluster are composed of lamellar structures of fcc and hcp atoms, indicating that the lamellar structure of fcc and hcp atoms in the silver crystal originates from nucleation, and not from the growth of crystals.

Local structure changes of Cu55 cluster during heating

The structural relaxation of a cluster containing 55 atoms at elevated temperatures is simulated by molecular dynamics. The interatomic interactions are given by using the embedded atom method （EAM） potential. By decomposing the peaks of the radial distribution functions （RDFs） according to the pair analysis technique, the local structural patterns are identified for this cluster. During increasing temperature, structural changes of different shells determined by atom density profiles result in an abrupt increase in internal energy. The simulations show how local structural changes can strongly cause internal energy to change accordingly.

Modification to solution-diffusion model for performance prediction of nanofiltration of long-alkyl-chain ionic liquids aqueous solutions based on ion cluster

Mathematical modeling for nanofiltration of ionic liquids(ILs) solutions could assist to understand transfer mechanism and predict experimental values. In this work, modeling by solution-diffusion model for nanofiltration of long-alkyl-chain ILs aqueous solutions was proposed. Molecular simulations were performed to validate the existence of ion cluster in long-alkyl-chain ILs aqueous solution. Based on the results of simulations, parameters used in the solution-diffusion model were modified, such as concentration of ILs and diameter of ion cluster.The modeling process was developed for three long-alkyl-chain ILs aqueous solutions with different concentrations(1-alkyl-3-methylimidazolium chloride: [C6 mim]Cl, [C8 mim]Cl, [C10 mim]Cl). The calculated values obtained from modified solution-diffusion model could well match the experimental values.

Size effect of He clusters on the interactions with self-interstitial tungsten atoms at different temperatures

The behaviors of helium clusters and self-interstitial tungsten atoms at different temperatures are investigated with the molecular dynamics method. The self-interstitial tungsten atoms prefer to form crowdions which can tightly bind the helium cluster at low temperature. The crowdion can change its position around the helium cluster by rotating and slipping at medium temperatures, which leads to formation of combined crowdions or dislocation loop locating at one side of a helium cluster. The combined crowdions or dislocation loop even separates from the helium cluster at high temperature. It is found that a big helium cluster is more stable and its interaction with crowdions or dislocation loop is stronger.

Formation and evolution properties of clusters in liquid metal copper during rapid cooling processes

Based on the quantum Sutton-Chen many-body potential,a molecular dynamics simulation was performed to investigate the formation and evolution properties of clusters in liquid Cu with 50 000 atoms.The cluster-type index method(CTIM)was used to describe the complex microstructure transitions.It is demonstrated that the amorphous structures are mainly formed with the three bond-types of 1551,1541 and 1431 in the system,and the icosahedral cluster(12 0 12 0)and other basic polyhedron clusters of(12 2 8 2),(13 1 10 2),(13 3 6 4),(14 1 10 3),(14 2 8 4)and(14 3 6 5)play a critical and leading role in the transition from liquid to glass.The nano-clusters formed in the system consist of some basic clusters and middle cluster configurations by connecting to each other,and distinguish from those obtained by gaseous deposition and ionic spray.From the results of structural parameter pair distribution function g(r),bond-types and basic cluster-types,it is found that the glass transition temperature Tg for liquid metal Cu is about 673 K at the cooling rate of 1.0×1014 K/s.

Molecular dynamics simulations of displacement cascades in Fe-10%Cr systems

Molecular dynamics simulations of the displacement cascades in Fe 10%Cr systems are used to sinmlate the primary knocked-on atom events of the irradiation damage at temperatures 300, 600, and 750 K with primary knockedon atom energies between 1 and 15 keV. The results indicate that the vacancies produced by the cascade are all in the central region of the displacement cascade. During the cascade, all recoil Fe and Cr atoms combine with each other to form Fe Cr or Fe Fe interstitial dumbbells as well as interstitial clusters. The number and the size of interstitial clusters increase with the energy of the primary knocked-on atom and the temperature. A few large clusters consist of a large number of lee interstitials with a few Cr atoms, the rest are lee Cr clusters with small and medium sizes. The interstitial dumbbells of Fe lee and Fe-Cr are in the （111）and （110） series directions, respectively.

Atomic mobility in energetic cluster deposition

This paper tries to outline the influence of atomic mobility on the initial fabrication of thin films formed by LECBD. Based on our recent studies on low-energy cluster beam deposition (LECBD) by molecular dynamics simulation, two examples, the deposition of small carbon clusters on Si and diamond surfaces and Al clusters on Ni substrate, were mainly discussed. The impact energy of the cluster ranges from 0.1 eV to 100 eV. In the former case, the mobility and the lateral migration of surface atoms, especially the recoil atoms, are enhanced with increasing the impact energy, which promote the film to be smoother and denser. For the latter case, the transverse kinetic energy of cluster atoms, caused mainly by the collision between moving cluster atoms, dominates the lateral spread of cluster atoms on the surface, which is contributive to layer-by-layer growth of thin films. Our result is consistent with the experimental observations that the film structure is strongly dependent on the impact energy. In addition, it elucidates that the atomic mobility takes a leading role in the structure characteristic of films formed by LECBD.

Molecular Dynamics Studies of Collapse Stages in the chain folding process of Polyethylene

Molecular dynamics simulations of model polyethylene chains with various lengths up to 4000 CH2 units were performed. Our results verified that the modification of force fields to demonstrate the collapse stages in the folding process is not necessary. The transition between the stages of chain folding became sharp, and the clusters became stable as the chain became longer.

MD and OKMC simulations of the displacement cascades in nickel

The molecular dynamics(MD) method was used to investigate the displacement cascades with primary knock-on atom(PKA) energies of 2-40 keV at 100 and600 K.The migration energy of defects and their clusters was calculated by nudged elastic band(NEB) method.Object kinetic Monte Carlo(OKMC) was used to simulate the evolution of defects in Ni under annealing.In each annealing stage,the recombination mechanism was discussed and evolution of the defects under different cascade conditions was compared.It was found that the defects generated in high-temperature cascades are more stable than those in the low-temperature cascades.In addition,almost all the defects are annihilated during annealing process at low PKA energy.At PKA energy of 20-40 keV,however,a large number of defects would remain after annealing.

Molecular Dynamics Simulations of Displacement Cascade in Ni-Based Concentrated Solid Solution Alloys

Single-phase concentrated solid solution alloys(SP-CSAs),including high-entropy alloys,have received extensive attention due to their excellent irradiation resistance.In this work,displacement cascade simulations are conducted using the molecular dynamics method to study the evolution of defects in Ni-based SP-CSAs.Compared with pure Ni,the NiCr,NiCo,and NiCu alloys exhibit a larger number of Frankel pairs(FPs)in the thermal peak stage,but a smaller number of surviving FPs.However,the NiFe alloy displays the opposite phenomenon.To explain these different observations for NiFe and other alloys,the formation energy and migration energy of interstitials/vacancies are calculated.In the NiFe alloy,both the formation energy and migration energy barrier are higher.On the other hand,in NiCr and other alloys,the formation energy of interstitials/vacancies is lower,as is the migration energy barrier of interstitials.The energy analysis agrees well with previous observations.The present work provides new insights into the mechanism behind the irradiation resistance of binary Ni-based SP-CSAs.