Isothermal-isobaric Molecular Dynamics Research Articles

Effects of leucine on hydrate formation: A combined experimental and molecular dynamics study

Compared to conventional natural gas storage and transportation methods, hydrates offer a safer, more compact alternative with broad future potential. To enhance hydrate formation, there has been considerable interest in environmentally friendly and cost-effective amino acid promoters, such as leucine. This study utilized a transparent vessel for CH4 hydrate formation experiments, identifying the optimal 0.5–0.6 wt% leucine concentration, which led to enhanced gas consumption, shorter induction time, and maximum gas storage density. In the leucine system, hydrates exhibited distinctive branching growth on container walls, both upward and downward. The formation of numerous microchannels amplified capillary effects, enhancing gas–liquid contact and fostering hydrate formation. Simultaneously, isobaric-isothermal molecular dynamics (MD) simulations integrated system growth configurations, structural parameters, radial distribution functions, energy barriers and mean square displacements, unveiling leucine’s microscopic impact on hydrate formation. Simulation results unequivocally show that leucine significantly accelerates CH4 hydrate growth by momentarily adsorbing at the solid–liquid interface, expediting water molecule ordering into cage-like structures. Leucine also enhances water molecule structural order around the mobile leucine molecule. Moreover, leucine lowers the energy barrier for methane transitioning from gas to liquid, facilitating gas molecule diffusion into the liquid phase and affording plenty interaction time for water and methane molecules, thereby providing favorable conditions for hydrate formation. These findings provide crucial perspectives for future research in the field of eco-friendly hydrate promoters, paving the way for innovative solutions in natural gas storage and transportation.

Molecular dynamics study on relaxation of supercooled liquid water at different cooling rates

The relaxation in supercooled liquids is a phenomena by which the thermodynamic system starts transforming from one state to other after losing its high end thermodynamic properties (i.e. density, potential energy, enthalpy, etc.) with a significant change. Since the relaxation phenomena is quite challenging to observe, specifically near the transition region of liquid water. Hence, the thermodynamic analysis for understanding the relaxation phenomena of liquid water in supercooled region from the high density liquid (HDL) state to low density amorphous (LDA) state are important from both fundamental and as well as applied perspectives. Therefore, we focused our study on performing thermodynamic analysis on the relaxation of super-cooled liquid water by using a monatomic water (mW) potential model with an isothermal-isobaric (NPT) molecular dynamics (MD) simulation technique at zero pressure under different cooling rates. In this work, we examine the relaxation in terms of root mean square (RMS) fluctuation of two and three-body potential energy of the maximum size of four coordinated cluster under different cooling rates. The predicted data over this study shows the independency of the cooling rate, where the relaxation of the liquid state of mW water results in a sharp irreversible drop in potential energy and large fluctuation in two and three body potential energy of the maximum size of four coordinated cluster. Further, we found the relaxation of liquid state of mW water associated with the sharp and continuous increase in the fraction of four coordinated particles as a function of temperature for different cooling rates. The predicted results of this study are consistent with the literature, which helps to understand the relaxation phenomena of mW water and other such thermodynamic systems more significantly.

Active sparse Bayesian committee machine potential for isothermal-isobaric molecular dynamics simulations.

Recent advancements in machine learning potentials (MLPs) have significantly impacted the fields of chemistry, physics, and biology by enabling large-scale first-principles simulations. Among different machine learning approaches, kernel-based MLPs distinguish themselves through their ability to handle small datasets, quantify uncertainties, and minimize over-fitting. Nevertheless, their extensive computational requirements present considerable challenges. To alleviate these, sparsification methods have been developed, aiming to reduce computational scaling without compromising accuracy. In the context of isothermal and isobaric ML molecular dynamics (MD) simulations, achieving precise pressure estimation is crucial for reproducing reliable system behavior under constant pressure. Despite progress, sparse kernel MLPs struggle with precise pressure prediction. Here, we introduce a virial kernel function that significantly enhances the pressure estimation accuracy of MLPs. Additionally, we propose the active sparse Bayesian committee machine (BCM) potential, an on-the-fly MLP architecture that aggregates local sparse Gaussian process regression (SGPR) MLPs. The sparse BCM potential overcomes the steep computational scaling with the kernel size, and a predefined restriction on the size of kernel allows for fast and efficient on-the-fly training. Our advancements facilitate accurate and computationally efficient machine learning-enhanced MD (MLMD) simulations across diverse systems, including ice-liquid coexisting phases, Li10Ge(PS6)2 lithium solid electrolyte, and high-pressure liquid boron nitride.

Crystal Structure Predictions for 4-Amino-2,3,6-trinitrophenol Using a Tailor-Made First-Principles-Based Force Field

Predictions of crystal structures from first-principles electronic structure calculations and molecular simulations have been performed for an energetic molecule, 4-amino-2,3,6-trinitrophenol. This physics-based approach consists of a series of steps. First, a tailor-made two-body potential energy surface (PES) was constructed with recently developed software, autoPES, using symmetry-adapted perturbation theory based on a density-functional theory description of monomers [SAPT(DFT)]. The fitting procedure ensures asymptotic correctness of the PES by employing a rigorous asymptotic multipole expansion, which seamlessly integrates with SAPT(DFT) interaction energies. Next, crystal structure prediction (CSP) was performed by generating possible crystal structures with rigid molecules, minimizing these structures using the SAPT(DFT) force field, and running isothermal–isobaric molecular dynamics (MD) simulations with flexible molecules based on the tailor-made SAPT(DFT) intermolecular force field and a generic/SAPT(DFT) intramolecular one. This workflow led to the experimentally observed structure being identified as one of the forms with the lowest lattice energy, demonstrating the success of a first-principles, bottom-up approach to CSP. Importantly, we argue that the accuracy of the intermolecular potential, here the SAPT(DFT)-based potential, is determinative of the crystal structure, while generic/SAPT(DFT) force fields can be used to represent the intramolecular potential. This force field approach simplifies the CSP workflow, without significantly compromising the accuracy of the prediction.

Revisiting the composition dependence of the properties of water-dimethyl sulfoxide liquid mixtures. Molecular dynamics computer simulations

We have revisited the composition dependence of principal properties of liquid water-DMSO mixtures by using the isobaric-isothermal molecular dynamics computer simulations. A set of non-polarizable semi-flexible models for the DMSO molecule combined with the TIP4P-2005 and TIP4P/ε water models is considered. We restrict calculations to atmospheric pressure, 0.1013 MPa, and room temperature, 298.15 K. Composition trends of density, excess mixing volume and excess mixing enthalpy, partial molar volumes and partial molar enthalpies of species, apparent molar volumes are reported. Besides, we explore composition trends of the self-diffusion of species, the static dielectric constant and the surface tension. Evolution of the microscopic structure of the mixture with composition is analyzed in terms of radial distributions functions, coordination numbers and the fractions of hydrogen-bonded molecules. We intend to capture the peculiarities of mixing the species in the mixture upon the DMSO molar fraction and the anomalous behaviors, if manifested in each of the properties under study. The quality of several combinations of the models for species is evaluated in detail to establish the possibility of necessary improvements.

Group-based evaluation of temperature and pressure for molecular dynamics simulation with a large time step

Recently, we proposed novel temperature and pressure evaluations in molecular dynamics (MD) simulations to preserve the accuracy up to the third order of a time step, δt [J. Jung, C. Kobayashi, and Y. Sugita, J. Chem. Theory Comput. 15, 84-94 (2019); J. Jung, C. Kobayashi, and Y. Sugita, J. Chem. Phys. 148, 164109 (2018)]. These approaches allow us to extend δt of MD simulations under an isothermal-isobaric condition up to 5 fs with a velocity Verlet integrator. Here, we further improve the isothermal-isobaric MD integration by introducing the group-based evaluations of system temperature and pressure to our previous approach. The group-based scheme increases the accuracy even using inaccurate temperature and pressure evaluations by neglecting the high-frequency vibrational motions of hydrogen atoms. It also improves the overall performance by avoiding iterations in thermostat and barostat updates and by allowing a multiple time step integration such as r-RESPA (reversible reference system propagation algorithm) with our proposed high-precision evaluations of temperature and pressure. Now, the improved integration scheme conserves physical properties of lipid bilayer systems up to δt = 5 fs with velocity Verlet as well as δt = 3.5 fs for fast motions in r-RESPA, respectively.

On the apparent molar volume of methanol in water-methanol mixtures. Composition and temperature effects from molecular dynamics study

We explore the effects of composition and temperature on the apparent molar volumes of species of water-methanol mixtures. Isothermal-isobaric molecular dynamics simulations are used with this purpose. Several combinations of models for water and for methanol are explored. Validity of predictions concerned with a puzzling minimum of apparent molar volume of methanol in water-rich solutions is tested against experimental results.

Second critical point in two realistic models of water.

The hypothesis that water has a second critical point at deeply supercooled conditions was formulated to provide a thermodynamically consistent interpretation of numerous experimental observations. A large body of work has been devoted to verifying or falsifying this hypothesis, but no unambiguous experimental proof has yet been found. Here, we use histogram reweighting and large-system scattering calculations to investigate computationally two molecular models of water, TIP4P/2005 and TIP4P/Ice, widely regarded to be among the most accurate classical force fields for this substance. We show that both models have a metastable liquid-liquid critical point at deeply supercooled conditions and that this critical point is consistent with the three-dimensional Ising universality class.

Liquid-liquid critical point of water

Chemical Physics Despite a broad range of experimental observations that indirectly point out the possible existence of a liquid-liquid phase transition in deeply supercooled water, no unambiguous experiment has shown this yet. The challenges of performing experiments under such conditions are associated with the inevitable rapid crystallization of the metastable liquid state, and therefore computational simulations become a crucial alternative. Using two accurate classical models of water and long (more than 40 microseconds) isobaricisothermal molecular dynamics simulations, Debenedetti et al. provide strong computational evidence of the presence of a second metastable liquid-liquid critical point in the deeply supercooled region that has critical behavior consistent with the three-dimensional Ising universality class. Science this issue p. [289][1] [1]: /lookup/doi/10.1126/science.abb9796

The cocrystal mechanism of HMX and LLM-105 by theoretical simulations

Isothermal-isobaric molecular dynamics (NPT-MD) simulations and density functional theory (DFT) calculations were performed to explore the cocrystal formation mechanism of energetic materials. Based on the morphologically dominant crystal habits of HMX (octahydro-1,3,5,7-tetranitro-1,3,5,7-tetrazocine) and LLM-105 (2,6-diamino-3,5-dinitropyrazine-l-oxide), six interface structures of HMX with LLM-105 have been established to obtain the interaction energy, radial distribution function (RDF), electrostatic potential and mechanical properties. It revealed that the interaction energies in different models are different, owing to the various molecular orientation and stackingmanner of the contacted crystal face of HMX and LLM-105. The (1 1 –1)/(1 1 0) interface is mainly governed by the O···H hydrogen bond and strong van der Waals forces and this model possesses the largest interaction energy. Electrostatic potentials indicate (1 1 –1)/(1 1 0) interface blending the strongest, which reveals that the cocrystal process is stronglyrelevant to the interface behavior of (1 1 –1)/(1 1 0) among the 6 interfaces. Moreover, the mechanical parameters show that the (1 1 –1)/(1 1 0) interface structure (VI) exhibits excellent ductility and fracture strength. Additionally, this work presents a novel insight, which from the crystal face of HMX and LLM-105 as well as the interaction behavior to research the cocrystal mechanism.

Quantitative Prediction of the Structure and Viscosity of Aqueous Micellar Solutions of Ionic Surfactants: A Combined Approach Based on Coarse-Grained MARTINI Simulations Followed by Reverse-Mapped All-Atom Molecular Dynamics Simulations.

We address the problem of the quantitative prediction of micelle formation in dilute aqueous solutions of ionic surfactants using sodium dodecyl sulfate (SDS) as a model system through a computational approach that involves three steps: (a) execution of coarse-grained simulations based on the MARTINI force field (with slightly modified parameters to afford the formation of large micelles); (b) reverse mapping of the final self-assembled coarse-grained configuration into an all-atom configuration; and (c) final relaxation of this all-atom configuration through short-time (on the order of a few tens of nanoseconds), detailed isothermal-isobaric molecular dynamics simulations using the CHARMM36 force field. For a given concentration of the solution in SDS molecules, the modified MARTINI-based coarse-grained simulations lead to the formation of large micelles characterized by mean aggregation numbers above the experimentally observed ones. However, by reintroducing the detailed chemical structure through a strategy that solves a well-defined geometric problem and re-equilibrating, these large micellar aggregates quickly dissolve to smaller ones and equilibrate to sizes that perfectly match the average micelle size measured experimentally at the given surfactant concentration. From the all-atom molecular dynamics simulations, we also deduce the surfactant diffusivity DSDS and the zero-shear rate viscosity, η0, of the solution, which are observed to compare very favorably with the few experimental values that we were able to find in the literature.

Functional uncertainty quantification for isobaric molecular dynamics simulations and defect formation energies

Functional uncertainty quantification (FunUQ) was recently proposed to quantify uncertainties in models and simulations that originate from input functions, as opposed to parameters. This paper extends FunUQ to quantify uncertainties originating from interatomic potentials in isothermal-isobaric molecular dynamics (MD) simulations and to the calculation of defect formation energies. We derive and verify a computationally inexpensive expression to compute functional derivatives in MD based on perturbation theory. We show that this functional derivative of the quantities of interest (average internal energy, volume, and defect energies in our case) with respect to the interatomic potential can be used to predict those quantities for a different interatomic potential, without re-running the simulation. The codes and scripts to perform FunUQ in MD are freely available for download. In addition, to facilitate reproducibility and to enable use of best practices for the approach, we created Jupyter notebooks to perform FunUQ analysis on MD simulations and made them available for online simulation in nanoHUB. The tool uses cloud computing resources and users can view, edit, and run end-to-end workflows from a standard web-browser without the need to download or install any software.

A polymerization scenario of the liquid–liquid transition of GeI4

SnI4 and GeI4 have been shown to exhibit similar polyamorphic nature. We examine the microscopic nature of the liquid–liquid transition in GeI4 by conducting an isothermal–isobaric molecular dynamics simulation for the system composed of rigid tetrahedral molecules. The model allows us to semiquantitatively discuss the structural properties of liquid GeI4 below 1 GPa. We define a physical bond between the nearest intermolecular iodine sites satisfying the conditions of forming the metallic I2 bond. We then focus on the formation of molecular clusters in dynamic networks of the bonds. The clusters are mainly formed by molecules whose nearest pairs are in edge-to-edge, face-to-edge, and vertex-to-edge orientations. The clusters grow as pressure increases, and the onset of percolation is observed below 1 GPa. The finite-size scaling analysis for the percolation probability identifies that the threshold pressure is GPa for the present model, which is near the extension of the boundary between the two liquid phases. We thus speculate that the liquid–liquid transition of GeI4 is attained by polymerization, i.e. percolation of molecular networks. The same percolation scenario with a slight modification such as introducing a bootstrap mechanism is expected to apply to the transition in liquid SnI4.

Molecular dynamics simulations of the properties of water-methanol mixtures. Effects of force fields

Isothermal-isobaric molecular dynamics simulations are used to examine the microscopic structure and some properties of water-methanol liquid mixture. The TIP4P/2005 and SPC/E water models are combined with the united atom TraPPE and the all-atom force field model for methanol. Our principal focus is to evaluate the quality of predictions of different combinations of model force fields concerning the composition dependence of basic properties of this system. Specifically, we explored the composition effects on density, excess molar volume and excess entropy, as well as on the surface tension and static dielectric constant. In addition, the structural properties are described in terms of the coordination numbers and the average number of hydrogen bonds between molecules of constituent species. Finally, the composition dependence of self-diffusion coefficients of the species is evaluated. All theoretical predictions are tested with respect to experimental data.

Successful test of the classical nucleation theory by molecular dynamic simulations of BaS

We used a recently developed two-body interatomic potential for barium sulfide (BaS), and performed molecular dynamics (MD) simulations with 36,000 particles to determine the kinetics of spontaneous, homogeneous nucleation and growth of BaS crystals in the supercooled liquid. Isothermal-isobaric MD simulations were accomplished at three temperatures. The calculated pair correlation function, along with several snapshots, allowed us to quantify the nucleation times, their crystal growth rates, and the time evolution of overall crystallization. Nucleation was spontaneously achieved in the supercooled liquid state; therefore, we computed the average birth times of the critical nuclei (also known as average onset time of the first nucleus) from 15 samples at each temperature, from which we computed the steady-state nucleation rates, MD Jss(T). Then, we independently obtained by MD the diffusion coefficients, D(T), the melting point, Tm, and the enthalpy of melting, ΔHm. Thus, the MD Jss could be compared with the predictions of the Classical Nucleation Theory (CNT) for homogeneous nucleation using only one fitting parameter, the nucleus/liquid interfacial free energy, σ. The calculated critical nucleus is made of only 2 to 3-unit cells, which is consistent with the critical size observed in the simulations, approximately 10–15 atoms. Using a constant (fitted value of) σ, and the D(T) and thermodynamic parameters from the simulations, we found that the MD pre-exponential factor has the same order of magnitude as the theoretical value predicted by the CNT. To the best of our knowledge, this agreement of the predictions of CNT and MD simulations has seldom been reported. Therefore, our results corroborate the validity of the CNT for simple supercooled liquids.

First principles isothermal-isobaric centroid molecular dynamics simulation of high pressure ices

To assess the performance of our implemented first principles isothermal-isobaric centroid molecular dynamics (CMD) we performed Born-Oppenheimer (BO) CMD simulations in the NPT ensemble for high pressure ices (HPIs) VIII, VII, and X at 270 K. By elevating pressure the crystal system is found to transform from tetragonal to cubic at ∼30 GPa, accompanying the proton order (ice VIII)-disorder (ice VII) transition observed to occur in our BO-CMD simulations. Moreover, our simulations are proved to be able to reproduce well both infrared and Raman spectra for HPIs in the frequency range above ∼500 cm−1.

Quantumness and state boundaries hidden in supercritical helium-4: A path integral centroid molecular dynamics study.

Isothermal-isobaric path integral centroid molecular dynamics simulations were conducted for fluid 4He at more than 600 state points in the proximity of the critical point to reveal the detailed states and relevant quantumness underlying the supercritical state. Through intensive analyses of the thermodynamic, dynamic, and quantum properties, we revealed the hidden state boundaries that separate the liquid-like and gas-like states in the supercritical region of this fluid. The Widom line, defined as the locus of the maxima of isobaric heat capacity C P , is also the quantum boundary at which there are changes in the isobaric temperature-dependence of the quantum wavelength, λ quantum, i.e., maximum amplitude of the Feynman imaginary-time paths (necklaces) of individual atoms. The Frenkel line, the famous dynamic state boundary, was observed to start from nearly the same point, 0.73-0.76 T c, on the P-T plane as observed for classical fluids. Several state boundaries based on the new criteria were found to emanate from the critical point or its vicinity on the P-T plane and are discussed in comparison with these boundaries. The quantumness of this fluid was expressed as (a) non-classical significant depression of C P observed in the liquid-like state; (b) the depression of the slopes dP/dT of the Widom line and the liquid-gas coexistence line near the critical point; and (c) the depression of the heat of pseudo-boiling across the Widom line. This is explained in terms of the decreasing kinetic energy with temperature observed in the liquid-like state below the Widom temperature T Widom, or alternatively in terms of the lattice model heat capacity, including the λ quantum.

Effects of ion concentration and solvent composition on the properties of water-methanol solutions of NaCl. NPT molecular dynamics computer simulation results

Isothermal-isobaric molecular dynamics simulations are used to examine the microscopic structure and other properties of a model solution consisting of NaCl salt dissolved in water-methanol mixture. The SPC/E water model and the united atom model for methanol are combined with the force field for ions by Dang [J. Amer. Chem. Soc., 1995, 117, 6954] to describe the entire system. Our principal focus is to study the effects of two variables, namely, the solvent composition and ion concentrations on the solution's density, on the structural properties, self-diffusion coefficients of the species and the dielectric constant. Moreover, we performed a detailed analysis of the first coordination numbers of the species. Trends of the behaviour of the average number of hydrogen bonds between solvent molecules are evaluated.

Thermodynamics of 1-ethyl-3-methylimidazolium methylsulfate co-solvent mixtures using molecular dynamics simulation

In the present work, isothermal-isobaric (NPT) and canonical (NVT) ensemble molecular dynamics (MD) simulations have been used to calculate the binary mixture properties (excess molar heats and volumes) of 1-ethyl-3-methylimidazolium methyl sulfate [EMIM][MeSO4] with six polar covalent molecules (acetic acid, acetone, chloroform, dimethyl sulfoxide, isopropyl alcohol and methanol) over the entire composition range. Excess partial molar properties as a function of mole fraction of co-solvent were calculated by using a Redlich-Kister expansion equation. The excess molar heats and volumes for [EMIM][MeSO4] with all the six co-solvents were all negative over the entire composition range. The size and structure of these molecules, as well as the nature of interaction between the mixture components seems to play a major role in the observed effects. These results have been interpreted in terms of hydrogen bonding and intermolecular interactions between ionic liquid ions and co-solvent molecules via radial distribution functions. The results showed that hydrogen bond interactions between ionic liquid and co-solvents play a major role in the volumetric and energetic properties of the mixtures. Further, radial distribution functions suggest that ionic liquid ions aggregate at higher concentrations of co-solvents.

Molecular mechanism of melting of a helical polymer crystal: Role of conformational order, packing and mobility of polymers

The molecular mechanism of melting of a superheated helical polymer crystal has been investigated using isothermal-isobaric molecular dynamics simulation that allows anisotropic deformation of the crystal lattice. A detailed microscopic analysis of the onset and progression of melting and accompanying changes in the polymer conformational order, translational, and orientation order of the solid along the melting pathway is presented. Upon gradual heating from room temperature to beyond the melting point at ambient pressure, the crystal exhibits signatures of premelting well below the solid-to-liquid melting transition at the melting point. The melting transition is manifested by abrupt changes in the crystal volume, lattice energy, polymer conformation, and dynamical properties. In the premelting stage, the crystal lattice structure and backbone orientation of the polymer chains are retained but with the onset of weakening of long-range helical order and interchain packing of polymers perpendicular to the fibre axis of the crystal. The premelting also marks the onset of conformational defects and anisotropic solid-state diffusion of polymers along the fibre axis. The present study underscores the importance of the interplay between intermolecular packing, interactions, and conformational dynamics at the atomic level in determining the macroscopic melting behavior of polymer crystals.