Reaction analysis and visualization of ReaxFF molecular dynamics simulations

Highlights

• This article presents the algorithms and applications of VARxMD.

• VARxMD allows analysis of chemical reactions for ReaxFF MD for the first time.

• Reactions generated directly from 3D coordinates and bond orders by bonding analysis.

• VARxMD has been applied in pyrolysis of large scale coal models and an HDPE model.

Abstract

ReaxFF MD (Reactive Force Field Molecular Dynamics) is a promising method for investigating complex chemical reactions in relatively larger scale molecular systems. The existing analysis tools for ReaxFF MD lack the capability of capturing chemical reactions directly by analyzing the simulation trajectory, which is critical in exploring reaction mechanisms. This paper presents the algorithms, implementation strategies, features, and applications of VARxMD, a tool for Visualization and Analysis of Reactive Molecular Dynamics. VARxMD is dedicated to detailed chemical reaction analysis and visualization from the trajectories obtained in ReaxFF MD simulations. The interrelationships among the atoms, bonds, fragments, species and reactions are analyzed directly from the three-dimensional (3D) coordinates and bond orders of the atoms in a trajectory, which are accomplished by determination of atomic connectivity for recognizing connected molecular fragments, perception of bond types in the connected fragments for molecules or radicals, indexing of all these molecules or radicals (chemical species) based on their 3D coordinates and recognition of bond breaking or forming in the chemical species for reactions. Consequently, detailed chemical reactions taking place between two sampled frames can be generated automatically. VARxMD is the first tool specialized for reaction analysis and visualization in ReaxFF MD simulations. Applications of VARxMD in ReaxFF MD simulations of coal and HDPE (high-density polyethylene) pyrolysis show that VARxMD provides the capabilities in exploring the reaction mechanism in large systems with complex chemical reactions involved that are difficult to access manually.

Introduction

The combination of classical molecular dynamics (MD) and the Reactive Force Field (ReaxFF) proposed by van Duin et al. [1], known as reactive molecular dynamics simulation, is a useful approach for molecular simulation of a complex system with chemical reactions [2], [3], [4], [5]. ReaxFF MD uses energy terms based on a general bond-order potential that fully addresses the chemistry of dynamic bonds. In addition to non-bonded interactions, ReaxFF MD employs the dynamic charge equilibration using Mortier's electronegativity equalization method (EEM) [6] at each time-step to account for polarization effects. Thus, ReaxFF MD allows the description of the formation, transition, and dissociation of chemical bonds in a molecular system with accuracy close to that of DFT [7] and with reduced computational costs [8]. Particularly, the potential functions in ReaxFF can automatically handle coordination changes associated with reactions; thus, pre-definitions of reactive sites or reaction pathways are unnecessary. Moreover, the recent progress of high performance computing programs for ReaxFF MD has enabled simulations of large molecular systems on supercomputers [9], [10], clusters [11] and desktop workstations [12]. The ReaxFF parallelization scheme in the more recent ADF platform has been rewritten to remove the global limit on the number of atoms in simulation systems [13]. The first Graphic Processing Unit (GPU) enabled ReaxFF molecular dynamics program is available and offers significantly increased computational capability on a single PC with a C2050 GPU attached, which makes performing ReaxFF MD simulations for larger system sizes and longer time scales possible on desktop workstations [12]. The progress that has been made on ReaxFF MD provides a promising approach for studying reaction mechanisms in processes with very complex chemical reactions at an atomic scale, where the reactions sites are difficult to predefine manually.

However, challenges arise in the analysis of ReaxFF MD simulation results for large molecular systems with complex chemical reactions. In addition to three-dimensional (3D) atomic coordinates, the trajectory from a ReaxFF MD simulation is in the form of bond orders without clear bond types between bonded atoms, in contrast to the trajectory obtained from a classical MD simulation. The existing analysis tools for ReaxFF MD simulations can only provide chemical formula-based analysis results that provide a time-step evolution of the system for the number of molecules based on these formulas. As shown in Fig. 1, the reax tool in LAMMPS can only produce a text list file, while the ADF platform can provide similar information in the form of graphs. Because one chemical formula may represent many chemical species in terms of 3D structure of a chemically connected molecular fragment in simulated molecular systems, the existing analysis tools are not capable of distinguishing between the species in the reaction products accurately. More importantly, the information on the chemical reactions in the simulated systems is not available in the analysis results, which is critical in understanding the reaction mechanisms via ReaxFF MD simulation. As a result, the analysis of the trajectory has to be performed manually. Manual analysis is time consuming, yet affordable, for relatively small and simple systems simulated with ReaxFF MD, but this approach is impractical for large systems that involve complex chemical reactions, especially for the investigation of coal pyrolysis. A pyrolysis simulation of a bituminous coal model with 4976 atoms using the ReaxFF MD program revealed that more than 900 reactions might occur at 2000 K within a 250 ps simulation period when the trajectory output interval is 12.5 ps [14]. Even larger numbers of reactions can be expected if the trajectory output interval becomes smaller or for the simulation of larger molecular systems, such as in a pyrolysis simulation of the Liulin coal model with 28,351 atoms [15]. Therefore, VARxMD (Visualization and Analysis of Reactive Molecular Dynamics) was created to provide a new tool specifically designed to facilitate the analysis and visualization of detailed chemical reactions for the trajectory of ReaxFF MD simulation. This paper will present the methods, algorithms and applications of VARxMD.