Quantum Mechanical Molecular Dynamics Research Articles

Fast method for quantum mechanical molecular dynamics

As the processing power available for scientific computing grows, first-principles Born-Oppenheimer molecular dynamics simulations are becoming increasingly popular for the study of a wide range of problems in materials science, chemistry, and biology. Nevertheless, the computational cost of Born-Oppenheimer molecular dynamics still remains prohibitively large for many potential applications. Here we show how to avoid a major computational bottleneck: the self-consistent-field optimization prior to force calculations. The optimization-free quantum mechanical molecular dynamics method gives trajectories that are almost indistinguishable from an ``exact'' microcanonical Born-Oppenheimer molecular dynamics simulation even when low-prefactor linear scaling sparse matrix algebra is used. Our findings show that the computational gap between classical and quantum mechanical molecular dynamics simulations can be significantly reduced.

Effects of Molecular Dynamics Thermostats on Descriptions of Chemical Nonequilibrium.

The performance of popular molecular dynamics (MD) thermostat algorithms in constant temperature simulations of equilibrium systems is well-known. This is not the case, however, in the context of nonequilibrium chemical systems, such as chemical reactions or nanoscale self-assembly processes. In this work, we investigate the effect of popular thermostat algorithms on the "natural" (i.e., Hamiltonian) dynamics of a nonequilibrium, chemically reacting system. By comparing constant-temperature quantum mechanical MD (QM/MD) simulations of carbon vapor condensation using velocity scaling, Berendsen, Andersen, Langevin, and Nosé-Hoover chain thermostat algorithms with natural NVE simulations, we show that efficient temperature control and reliable reaction dynamics are mutually exclusive in such a system. This problem may be circumvented, however, by placing the reactive system in an inert He atmosphere, which is itself described using NVT MD. We demonstrate that both realistic temperature control and dynamics consistent with natural NVE dynamics can then be obtained simultaneously. In essence, the thermal energy created by the natural dynamics of the NVE subsystem is drained by the thermostat acting on the NVT atmosphere, without adversely affecting the dynamics of the reactive system itself.

The Stability of Bisulfite and Sulfonate Ions in Aqueous Solution Characterized by Hydration Structure and Dynamics

The aqueous solutions of bisulfite (SO(3)H(-)) and sulfonate (HSO(3)(-)) were simulated by the ab initio quantum mechanical charge field molecular dynamics (QMCF MD) formalism. All superimposed trajectories for the atomic coordinates of solutes with three-dimensional alignment here illustrated the reactivities of the ions. Power spectra were evaluated on the basis of the velocity autocorrelation functions (VACFs) with the normal-mode analysis, presenting a higher frequency of the symmetric SO(3) deformation (δ(s)(SO(3))) than the asymmetric SO(3) deformation (δ(as)(SO(3))) modes for the sulfonate ion. The different influence of solvent on the frequency of the O-H and S-H stretching suggests a higher stability of hydrated sulfonate ion. The bisulfite shows a slightly stronger molecular hydration shell than the sulfonate ion with the average number of ion-solvent hydrogen bonds (H-bonds) of 5.3 and 5.0, respectively. Extra water molecules within the molecular hydration shell are found for bisulfite (1.2) and for sulfonate (1.6). The mean residence times for the water ligands classify each ion as a structure maker, while the S-H bond within the sulfonate ion displays a hydrophobic behavior. No tautomerization was observed within the simulation period.

A theoretical study of N–H ··· π H-bond interaction of pyrrole: from clusters to the liquid

The N–H ··· π H-bond interactions of clusters and liquid formed by pyrrole molecules were investigated by Quantum Mechanics (QM) and classical Molecular Dynamics (MD), respectively. Based on the optimized geometry at the B97-D/aug-cc-pVTZ level of theory including dispersion correction, the nature and the origin of N–H ··· π H-bond interactions were unveiled by atoms in molecules (AIM), natural bond orbital (NBO), and energy decomposition analysis (EDA). Among them, the AIM analysis gives evidence to the presence of N–H ··· π H-bond interactions, the NBO examination reveals that π → σ* donor-acceptor orbital interaction is of great importance. EDA study indicates that N–H ··· π interactions are governed by the electrostatic and dispersion term. Meanwhile, MD simulation with OPLS-AA (optimized potentials for liquid simulations all-atom) was applied to study the pure liquid pyrrole at different temperature. The results confirm the existence of the N–H ··· π H-bond in the pure liquid pyrrole, and further characterized the structures of this H-bond which is somewhat different to the clusters.

Car-Parinello molecular dynamics simulations of thionitroxide and S-nitrosothiol in the gas phase, methanol, and water—A theoretical study of S-nitrosylation

A dilemma about whether thionitroxide radical (RSNHO) or S-nitrosothiol (RSNO) is observed in protein S-nitrosylation has arisen recently. To illustrate the effect of chemical environment on these structures, this paper presents quantum mechanical molecular dynamics of thionitroxide, and cis- and trans-S-nitrosothiols in the gas phase, methanol, and water. By using Car-Parrinello molecular dynamics (CPMD), we have observed that there is free rotation about the S-N bond at 300 K in thionitroxide, but no such rotation is observed for S-nitrosothiol. The C-S-N-O torsion angle distribution in thionitroxide is significantly dependent upon the surrounding environment, leading to either gauche-, cis-, or trans-conformation. In the case of S-nitrosothiol the C-S-N-O plane is twisted slightly by 5°–15° in the cis-isomer, while the periplanar structure is well-retained in the trans-isomer. The calculated results are in agreement with the X-ray crystallographic data of small molecular RSNO species. Interestingly, for both compounds, the CPMD simulations show that solvation can cause a decrease in the S-N bond length. Moreover, the oxygen atom of thionitroxide is found to be a good hydrogen-bond acceptor, forming an oxyanionhole-like hydrogen bonding network.

The Reversible Opening of Water Channels in Cytochrome c Modulates the Heme Iron Reduction Potential

Dynamic protein-solvent interactions are fundamental for life processes, but their investigation is still experimentally very demanding. Molecular dynamics simulations up to hundreds of nanoseconds can bring to light unexpected events even for extensively studied biomolecules. This paper reports a combined computational/experimental approach that reveals the reversible opening of two distinct fluctuating cavities in Saccharomyces cerevisiae iso-1-cytochrome c. Both channels allow water access to the heme center. By means of a mixed quantum mechanics/molecular dynamics (QM/MD) theoretical approach, the perturbed matrix method (PMM), that allows to reach long simulation times, changes in the reduction potential of the heme Fe(3+)/Fe(2+) couple induced by the opening of each cavity are calculated. Shifts of the reduction potential upon changes in the hydration of the heme propionates are observed. These variations are relatively small but significant and could therefore represent a tool developed by cytochrome c for the solvent driven, fine-tuning of its redox functionality.

Ab Initio Quantum Mechanical Charge Field Molecular Dynamics Simulation (QMCF-MD) of Bi3+ in Water

The hydration of the Bi(III) ion was determined via an ab initio quantum mechanical charge field molecular dynamics (QMCF-MD) simulation. Ten picosecond sampling was carried out to determine structural and dynamical properties of the Bi(III) ion in aqueous solution. In the first hydration shell, the ion is 9-fold coordinated with a maximum probability of the Bi-O distance at 2.51 Å. In total, 11 exchanges were observed in the first-shell showing associative, dissociative, and interexchange character. As with the dominant existence of 9-fold coordination, the geometry of the Bi(III) ion is in between the tricapped trigonal prism and the capped square antiprism.

Quantum mechanical/effective fragment potential molecular dynamics (QM/EFP-MD) study on intra-molecular proton transfer of glycine in water

We show that the hybrid quantum mechanical/effective fragment potential (QM/EFP) can be a very effective and practical quantum mechanical molecular dynamics method, when it is properly combined with well-developed traditional molecular dynamics (MD) techniques. QM/EFP-MD simulations on intra-molecular proton transfer of glycine with 290 EFP waters yielded accurate free energy change and reaction barrier of the zwitterion→neutral form conversion. Water rearrangements turned out to be the main driving force of the proton transfer.

Combining Dyad Protonation and Active Site Plasticity in BACE-1 Structure-Based Drug Design

The ability of the BACE-1 catalytic dyad to adopt multiple protonation states and the conformational flexibility of the active site have hampered the reliability of computational screening campaigns carried out on this drug target for Alzheimer's disease. Here, we propose a protocol that, for the first time, combining quantum mechanical calculations, molecular dynamics, and conformational ensemble virtual ligand screening addresses these issues simultaneously. The encouraging results prefigure this approach as a valuable tool for future drug discovery campaigns.

Hydration of Mg2+ and its influence on the water hydrogen bonding network via ab initio QMCF MD

An ab initio quantum mechanical charge field molecular dynamics (QMCF MD) simulation has been carried out in order to obtain detailed insight into the hydration properties of Mg2+. The ion, being a ‘kosmotrope’, considerably influences the water hydrogen bond network in its vicinity. Certain types of hydrogen bonds are strengthened at the cost of weakening others. The cationic influence on the water hydrogen bonding network have been investigated primarily employing time correlation functions of H-bond dynamics. The frequencies of the MgO6 skeletal normal modes calculated using the normal coordinate analysis shows a close resemblance to the experimental data, thus validating the accuracy of ab initio quantum chemical simulations.

A Theoretical Reappraisal of Polylysine in the Investigation of Secondary Structure Sensitivity of Infrared Spectra

Infrared spectroscopy has long provided a means to estimate the secondary structure of proteins and peptides. In particular, the vibrational spectra of the amide I' band have been widely used for this purpose as the frequency positions of the amide I' bands are related to the presence of specific secondary structures. Here, we calculate the amide I' IR spectra of polylysine in aqueous solution in its three secondary structure states, i.e., α-helix, β-sheet, and random coil, by means of a mixed quantum mechanics/molecular dynamics (QM/MD) theoretical-computational methodology based on the perturbed matrix method (PMM). The computed spectra show a good agreement with the experimental ones. Although our calculations confirm the importance of the excitonic coupling in reproducing important spectral features (e.g., the width of the absorption band), the frequency shift due to secondary-structure changes is also well reproduced without the inclusion of the excitonic coupling, pointing to a role played by the local environment. Concerning the β-conformation spectrum, which is characterized by a double-peak amide I' band due to excitonic coupling, our results indicate that it does not correspond to a generic antiparallel β-sheet (e.g., of the typical size present in native proteins) but is rather representative of extended β-structures, which are common in β-aggregates. Moreover, we also show that the solvent has a crucial role in the shape determination of the β-conformation amide I' band and in particular in the disappearance of the high-frequency secondary peak in the case of small sheets (e.g., 6-stranded).

On the structure and dynamics of the hydrated sulfite ion in aqueous solution – an ab initio QMCF MD simulation and large angle X-ray scattering study

Theoretical ab initio quantum mechanical charge field molecular dynamics (QMCF MD) formalism has been applied in conjunction to experimental large angle X-ray scattering to study the structure and dynamics of the hydrated sulfite ion in aqueous solution. The results show that there is a considerable effect of the lone electron-pair on sulfur concerning structure and dynamics in comparison with the sulfate ion with higher oxidation number and symmetry of the hydration shell. The S-O bond distance in the hydrated sulfite ion has been determined to 1.53(1) Å by both methods. The hydrogen bonds between the three water molecules bound to each sulfite oxygen are only slightly stronger than those in bulk water. The sulfite ion can therefore be regarded as a weak structure maker. The water exchange rate is somewhat slower for the sulfite ion than for the sulfate ion, τ(0.5) = 3.2 and 2.6 ps, respectively. An even more striking observation in the angular radial distribution (ARD) functions is that the for sulfite ion the water exchange takes place in close vicinity of the lone electron-pair directed at its sides, while in principle no water exchange did take place of the water molecules hydrogen bound to sulfite oxygens during the simulation time. This is also confirmed when detailed pathway analysis is conducted. The simulation showed that the water molecules hydrogen bound to the sulfite oxygens can move inside the hydration shell to the area outside the lone electron-pair and there be exchanged. On the other hand, for the hydrated sulfate ion in aqueous solution one can clearly see from the ARD that the distribution of exchange events is symmetrical around the entire hydration sphere.

Time-dependent quantum dynamical simulations of C2condensation under extreme conditions

We report theoretical studies of the initial phase of bulk C(2) condensation into carbon nano-structures by means of Born-Oppenheimer and time-dependent quantum mechanical Liouville-von Neumann molecular dynamics based on the density-functional tight-binding (DFTB) framework for electrons. We observe that the time-dependent quantum mechanical approach leads to faster formation of carbon nanostructures than analogous Born-Oppenheimer simulations. Our results suggest that the condensation of bulk carbon is nonadiabatic in nature, with the critical role of electronic stopping as in ion-irradiation of materials. Contrary to time-dependent quantum mechanical simulations, Born-Oppenheimer dynamics incorrectly predict that the short carbon chains obtained from initial reactive collisions between C(2) quickly evaporate, leading to much lower probability of secondary collisions and condensation. We also discuss some deficiencies in Born-Oppenheimer dynamics that lead to unphysical charge polarization and electron transfer.

Simulation of electronic excitation in the liquid state by quantum mechanical charge field molecular dynamics

The ab initio quantum mechanical charge field molecular dynamics simulation technique was applied for the first time to an electronically excited system, namely hydrated Li(I). This approach does not only provide data for the influence of electronic excitation on structure and composition of the hydrate, but also leads to a spectrum for the excitation band including rotational, translational and vibrational energy transfers. The first excitation band for Li(I) in water peaks at 13.1eV. The excitation leads to reduced coordination numbers of first and second hydration shell and stabilises the first and labilises the second shell.

An ab initio quantum mechanical charge field molecular dynamics simulation of hydrogen peroxide in water

An ab initio quantum mechanical charge field molecular dynamics (QMCF MD) simulation was performed on a single molecule of hydrogen peroxide immersed in water to investigate its stability in aqueous solution, since pure hydrogen peroxide is very unstable. The structural parameters such as radial distribution functions (RDFs), coordination number distributions (CNDs) and angular distribution functions (ADFs) indicate the existence of ∼4 hydrogen bonds between hydrogen peroxide and water molecules, with both molecules acting as hydrogen bond donors and hydrogen bond acceptors. The overall hydration shell consists of ∼6 water molecules surrounding the hydrogen peroxide molecule. The analysis of the hydrogen bond dynamics verified the presence of strong hydrogen bonds compared to pure water, thus confirming the stabilization of hydrogen peroxide in aqueous solution.

Symmetry Breaking and Hydration Structure of Carbonate and Nitrate in Aqueous Solutions: A Study by Ab Initio Quantum Mechanical Charge Field Molecular Dynamics

The ab initio quantum mechanical charge field molecular dynamics (QMCF MD) formalism was applied to simulate carbonate and nitrate anions in aqueous solution. The out-of-plane (ν(2)) spectra obtained from the velocity autocorrelation functions (VACFs) and the torsion angle-time functions indicate that the symmetry of carbonate is reduced from D(3h) to a lower degree by breaking up the molecular plane, whereas the planarity of nitrate anion is retained. The calculated frequencies are in good agreement with the Raman and IR data. Carbonate shows a stronger molecular hydration shell than the nitrate anion with the average molecular coordination numbers of 8.9 and 7.9, respectively. A comparison with the average number of ion-solvent hydrogen bonds (H-bonds) indicates the extra water molecules within the hydration shell of carbonate (∼2) and nitrate (∼3), readily migrating from one coordinating site to another. The mean residence times for water ligands in general classify carbonate and nitrate as moderate and weak structure-making anions, while the specific values for individual sites of nitrate reveal local weak structure-breaking properties.

Adaptive-Partitioning Redistributed Charge and Dipole Schemes for QM/MM Dynamics Simulations: On-the-fly Relocation of Boundaries that Pass through Covalent Bonds.

Recently, Heyden, Lin, and Truhlar (J. Phys. Chem. B2007, 111, 2231-2241) formularized the adaptive-partitioning schemes for quantum mechanical and molecular mechanical (QM/MM) molecular dynamics simulations. The adaptive-partitioning schemes permit on-the-fly reclassification of atoms/groups as part of the QM or MM subsystems during dynamics simulations. Test simulations of argon atoms in a periodic box with dual-level MM potentials in the microcanonical ensemble demonstrated that the adaptive-partitioning schemes conserved energy and momentum, which is critical to ensure correct sampling of configuration spaces of desired ensembles. In this work, we reported the extension of the adaptive-partitioning schemes to deal with groups that are molecular fragments. The newly developed adaptive-partitioning redistributed charge scheme and adaptive-partitioning redistributed charge and dipole schemes allow on-the-fly relocation of the QM/MM boundaries that cut through covalent bonds during dynamics simulations. Test QM/MM simulations with a variety of QM levels of theory in the microcanonical ensembles demonstrated that the new schemes conserve energy and momentum.

Absorption and Fluorescence of PRODAN in Phospholipid Bilayers: A Combined Quantum Mechanics and Classical Molecular Dynamics Study

Absorption and fluorescence spectra of PRODAN (6-propionyl-2-dimethylaminonaphthalene) were studied by means of the time-dependent density functional theory and the algebraic diagrammatic construction method. The influence of environment, a phosphatidylcholine lipid bilayer and water, was taken into account employing a combination of quantum chemical calculations with empirical force-field molecular dynamics simulations. Additionally, experimental absorption and emission spectra of PRODAN were measured in cyclohexane, water, and lipid vesicles. Both planar and twisted configurations of the first excited state of PRODAN were taken into account. The twisted structure is stabilized in both water and a lipid bilayer, and should be considered as an emitting state in polar environments. Orientation of the excited dye in the lipid bilayer significantly depends on configuration. In the bilayer, the fluorescence spectrum can be regarded as a combination of emission from both planar and twisted structures.

Hydration of highly charged ions

Based on a series of ab initio quantum mechanical charge field molecular dynamics (QMCF MD) simulations, the broad spectrum of structural and dynamical properties of hydrates of trivalent and tetravalent ions is presented, ranging from extreme inertness to immediate hydrolysis. Main group and transition metal ions representative for different parts of the periodic system are treated, as are 2 threefold negatively charged anions. The results show that simple predictions of the properties of the hydrates appear impossible and that an accurate quantum mechanical simulation in cooperation with sophisticated experimental investigations seems the only way to obtain conclusive results.

Acetylene⋅⋅⋅Furan Trimer Formation at 0.37 K as a Model for Ultracold Aggregation of Non‐ and Weakly Polar Molecules

We have studied the aggregation process of (C(2)H(2))⋅⋅⋅furan trimers at ultracold temperatures (0.37 K) in helium nanodroplets. Computational sampling of the potential energy surface using the multiple-minima-hypersurface (MMH) approach yielded seven possible minimum structures, optimized at the MP2 level of theory with the cc-pVTZ and 6-311++G(d,p) basis sets. Experimentally, we could assign five transitions in the IR spectrum of acetylene-furan aggregates in the acetylene C-H(asym) stretch region between 3240 and 3300 cm(-1) to vibrational bands of the 2:1 acetylene-furan trimer. The transitions were assigned to three ring structures that all contain the T-shaped acetylene dimer as structural sub-motif. Two of the structures form a nonplanar ring involving a C-H(Ac) ⋅⋅⋅π(Fu) bond, the third is a nearly planar ring containing a C-H(Ac) ⋅⋅⋅O(Fu) bond. This assignment was corroborated by quantum mechanical/molecular dynamics (QM/MD) simulations mimicking in detail the aggregation process of precooled monomers. The simulations provided evidence for a transition from a higher level local minimum to the global minimum state over a small barrier during the aggregation process. The experimentally observed structures can be explained by a step-by-step aggregation of moieties pre-cooled to 0.37 K that are steered by intermediate and short-range electrostatic interactions. Thus, we are able to unravel a special aggregation mechanism which differs from aggregation of molecules with large dipole moments where this aggregation process is dominated by long range 1/r(3) dipole-dipole interaction ("electrostatic steering"). This mechanism is expected to be a general mechanism in ultracold chemistry.