Quantum-classical Molecular Dynamics Simulations Research Articles

Infrared-laser-pulse-enhanced ultrafast fragmentation of N22+ following Auger decay: Mixed quantum-classical simulations

We employ mixed quantum-classical molecular dynamics simulations to investigate the fragmentation of ${\mathrm{N}}_{2}$ molecules after core-level photoionization by an x-ray laser, subsequent Auger decay, and followed by a femtosecond IR pulse that interacts with ${{\mathrm{N}}_{2}}^{2+}$. The delayed IR pulse couples the dissociative electronic states of ${{\mathrm{N}}_{2}}^{2+}$ with electronic states that can support long-lived vibrational resonances. We compare our simulations with previous quantum dynamics calculations in a quasidiabatic representation, which employed a small number of electronic states. Good agreement for both the Auger spectrum as well as the influence of the delayed IR pulse is found. By employing the mixed quantum-classical treatment, we can greatly reduce the computational cost to simulate the fragmentation dynamics compared to the quantum dynamics simulations. Furthermore, we reinvestigate the title process by employing an extended set of adiabatic potential energy surfaces and also investigate the role of nonadiabatic coupling in the process. The use of the full set of adiabatic potentials increases the dissociation probability and changes the details of the interaction with the IR pulse, but no effect due to the nonadiabatic coupling is found.

Low-energy hydrogen uptake by small-cage Cn and Cn-1B fullerenes

We present a theoretical study of the hydrogen uptake capability of carbon fullerene cages Cn and their boron-doped heterofullerene equivalents Cn-1B, with n = 20, 40, and 60, irradiated by hydrogen atoms in an impact energy range of 0.1–100 eV. In order to predict exohedral and endohedral hydrogen captures as well as the scattering probability of hydrogen for various cage types and sizes, we perform quantum-classical molecular dynamics (QCMD) calculations using the self-consistent-charge density-functional tight-binding (SCC-DFTB) method. Maximum endohedral hydrogen capture probabilities of 20% for n = 60 and 14% for n = 40 are found at impact energies close to 15 eV for both Cn and Cn-1B systems. For n = 20, however, endohedral capture is observed at a maximum of 2%, while the exohedral capture reaches a maximum of 5% both at 15 eV. Similar results for the hydrogen capture are obtained by classical molecular dynamics based on the ReaxFF potential. Finally, the stopping cross section per carbon atom from the QCMD simulations for all cage sizes displays a linear dependence on the projectile velocity with a threshold at 0.8 eV, and extrapolates well to the available theoretical data.

A path for synthesis of boron-nitride nanostructures in volume of arc plasma

We find a possible channel for direct nanosynthesis of boron-nitride (BN) nanostructures, including growth of BN nanotubes from a mixture of BN diatomic molecules by quantum–classical molecular dynamics simulations. No catalyst or boron nanoparticle is needed for this synthesis, however the conditions for the synthesis of each of the nanostructures, such as temperature and flux of the BN feedstock are identified and are compatible with the conditions in an electric arc at high pressure. We also find that BN nanostructures can be synthetized by feeding a boron nanoparticle by BN diatomic molecules, however if hydrogen rich molecules like NH3 or HBNH are used as a feedstock, two-dimensional nanoflake stable structures are formed.

Excess electrons in methanol clusters: Beyond the one-electron picture.

We performed a series of comparative quantum chemical calculations on various size negatively charged methanol clusters, CH3OHn-. The clusters are examined in their optimized geometries (n = 2-4), and in geometries taken from mixed quantum-classical molecular dynamics simulations at finite temperature (n = 2-128). These latter structures model potential electron binding sites in methanol clusters and in bulk methanol. In particular, we compute the vertical detachment energy (VDE) of an excess electron from increasing size methanol cluster anions using quantum chemical computations at various levels of theory including a one-electron pseudopotential model, several density functional theory (DFT) based methods, MP2 and coupled-cluster CCSD(T) calculations. The results suggest that at least four methanol molecules are needed to bind an excess electron on a hydrogen bonded methanol chain in a dipole bound state. Larger methanol clusters are able to form stronger interactions with an excess electron. The two simulated excess electron binding motifs in methanol clusters, interior and surface states, correlate well with distinct, experimentally found VDE tendencies with size. Interior states in a solvent cavity are stabilized significantly stronger than electron states on cluster surfaces. Although we find that all the examined quantum chemistry methods more or less overestimate the strength of the experimental excess electron stabilization, MP2, LC-BLYP, and BHandHLYP methods with diffuse basis sets provide a significantly better estimate of the VDE than traditional DFT methods (BLYP, B3LYP, X3LYP, PBE0). A comparison to the better performing many electron methods indicates that the examined one-electron pseudopotential can be reasonably used in simulations for systems of larger size.

Vibrational Quantum Decoherence in Liquid Water.

Traditional descriptions of vibrational energy transfer consider a quantum oscillator interacting with a classical environment. However, a major limitation of this simplified description is the neglect of quantum decoherence induced by the different interactions between two distinct quantum states and their environment, which can strongly affect the predicted energy-transfer rate and vibrational spectra. Here, we use quantum-classical molecular dynamics simulations to determine the vibrational quantum decoherence time for an OH stretch vibration in liquid heavy water. We show that coherence is lost on a sub-100 fs time scale due to the different responses of the first shell neighbors to the ground and excited OH vibrational states. This ultrafast decoherence induces a strong homogeneous contribution to the linear infrared spectrum and suggests that resonant vibrational energy transfer in H2O may be more incoherent than previously thought.

Ultrafast phosphate hydration dynamics in bulk H2O.

Phosphate vibrations serve as local probes of hydrogen bonding and structural fluctuations of hydration shells around ions. Interactions of H2PO4(-) ions and their aqueous environment are studied combining femtosecond 2D infrared spectroscopy, ab-initio calculations, and hybrid quantum-classical molecular dynamics (MD) simulations. Two-dimensional infrared spectra of the symmetric (νS(PO2(-))) and asymmetric (νAS(PO2(-))) PO2(-) stretching vibrations display nearly homogeneous lineshapes and pronounced anharmonic couplings between the two modes and with the δ(P-(OH)2) bending modes. The frequency-time correlation function derived from the 2D spectra consists of a predominant 50 fs decay and a weak constant component accounting for a residual inhomogeneous broadening. MD simulations show that the fluctuating electric field of the aqueous environment induces strong fluctuations of the νS(PO2(-)) and νAS(PO2(-)) transition frequencies with larger frequency excursions for νAS(PO2(-)). The calculated frequency-time correlation function is in good agreement with the experiment. The ν(PO2(-)) frequencies are mainly determined by polarization contributions induced by electrostatic phosphate-water interactions. H2PO4(-)/H2O cluster calculations reveal substantial frequency shifts and mode mixing with increasing hydration. Predicted phosphate-water hydrogen bond (HB) lifetimes have values on the order of 10 ps, substantially longer than water-water HB lifetimes. The ultrafast phosphate-water interactions observed here are in marked contrast to hydration dynamics of phospholipids where a quasi-static inhomogeneous broadening of phosphate vibrations suggests minor structural fluctuations of interfacial water.

Utilizing fast multipole expansions for efficient and accurate quantum-classical molecular dynamics simulations.

Recently, a novel approach to hybrid quantum mechanics/molecular mechanics (QM/MM) molecular dynamics (MD) simulations has been suggested [Schwörer et al., J. Chem. Phys. 138, 244103 (2013)]. Here, the forces acting on the atoms are calculated by grid-based density functional theory (DFT) for a solute molecule and by a polarizable molecular mechanics (PMM) force field for a large solvent environment composed of several 10(3)-10(5) molecules as negative gradients of a DFT/PMM hybrid Hamiltonian. The electrostatic interactions are efficiently described by a hierarchical fast multipole method (FMM). Adopting recent progress of this FMM technique [Lorenzen et al., J. Chem. Theory Comput. 10, 3244 (2014)], which particularly entails a strictly linear scaling of the computational effort with the system size, and adapting this revised FMM approach to the computation of the interactions between the DFT and PMM fragments of a simulation system, here, we show how one can further enhance the efficiency and accuracy of such DFT/PMM-MD simulations. The resulting gain of total performance, as measured for alanine dipeptide (DFT) embedded in water (PMM) by the product of the gains in efficiency and accuracy, amounts to about one order of magnitude. We also demonstrate that the jointly parallelized implementation of the DFT and PMM-MD parts of the computation enables the efficient use of high-performance computing systems. The associated software is available online.

Hydrated Electrons in Water Clusters: Inside or Outside, Cavity or Noncavity?

In this work, we compare the applicability of three electron–water molecule pseudopotentials in modeling the physical properties of hydrated electrons. Quantum model calculations illustrate that the recently suggested Larsen–Glover–Schwartz (LGS) model and its modified m-LGS version have a too-attractive potential in the vicinity of the oxygen. As a result, LGS models predict a noncavity hydrated electron structure in clusters at room temperature, as seen from mixed one-electron quantum–classical molecular dynamics simulations of water cluster anions, with the electron localizing exclusively in the interior of the clusters. Comparative calculations using the cavity-preferring Turi–Borgis (TB) model predict interior-state and surface-state cluster isomers. The computed associated physical properties are also analyzed and compared to available experimental data. We find that the LGS and m-LGS potentials provide results that appear to be inconsistent with the size dependence of the experimental data. The simulated TB tendencies are qualitatively correct. Furthermore, ab initio calculations on static LGS noncavity structures indicate weak stabilization of the excess electron in regions where the LGS potential preferably and strongly binds the electron. TB calculations give stabilization energies that are in line with the ab initio results. In conclusion, we observe that the cavity-preferring pseudopotential model predicts cluster physical properties in better agreement with experimental data and ab initio calculations than the models predicting noncavity structures for the hydrated electron.

Atomistic Modeling of Two-Dimensional Electronic Spectra and Excited-State Dynamics for a Light Harvesting 2 Complex

The Light Harvesting 2 (LH2) complex is a vital part of the photosystem of purple bacteria. It is responsible for the absorption of light and transport of the resulting excitations to the reaction center in a highly efficient manner. A general description of the chromophores and the interaction with their local environment is crucial to understand this highly efficient energy transport. Here we include this interaction in an atomistic way using mixed quantum-classical (molecular dynamics) simulations of spectra. In particular, we present the first atomistic simulation of nonlinear optical spectra for LH2 and use it to study the energy transport within the complex. We show that the frequency distributions of the pigments strongly depend on their positions with respect to the protein scaffold and dynamics of their local environment. Furthermore, we show that although the pigments are closely packed the transition frequencies of neighboring pigments are essentially uncorrelated. We present the simulated linear absorption spectra for the LH2 complex and provide a detailed explanation of the states responsible for the observed two-band structure. Finally, we discuss the energy transfer within the complex by analyzing population transfer calculations and 2D spectra for different waiting times. We conclude that the energy transfer from the B800 ring to the B850 ring is mediated by intermediate states that are delocalized over both rings, allowing for a stepwise downhill energy transport.

Water Gated Transitions in Proton Pumping of Respiratory Complex I

Respiratory Complex I or NADH:ubiquinone oxidoreductase is a redox-driven proton-pump. Powered by quinone reduction, Complex I drives the energetically uphill translocation of protons across the mitochondrial inner membrane and bacterial cytoplasmic membrane. The established electrochemical proton gradient provides the driving force for active transport and synthesis of ATP and is thus crucial for biological energy conversion. Complex I comprises a membrane domain that includes three antiporter-like subunits, involved in the proton-pumping process, and a soluble domain, responsible for reduction of quinones by electron transfer from NADH. Remarkably, site-directed mutagenesis experiments show that mutations of titratable residues in the antiporter-like subunits, ∼200 A away from site of quinone reduction, inhibit both proton pumping and quinone reduction. To explain this long-range proton-coupled electron transfer mechanism, both indirect and direct coupling models have been suggested. However, despite the recent elucidation of the complete intact structure of Complex I, the molecular principles of the coupling principles remain elusive. We present here results of large-scale classical and hybrid quantum-classical (QM/MM) molecular dynamics (MD) simulations of Complex I, embedded in biologically realistic environments. Our simulations indicate that water molecules provide important elements in the proton-pumping process. Our findings may form a basis for understanding long-range energy transduction in Complex I, and mechanistic similarities to other redox-driven proton-pumps such as cytochrome c oxidase and bacteriorhodopsin.

A Mixed Quantum-Classical Molecular Dynamics Study of anti-Tetrol and syn-Tetrol Dissolved in Liquid Chloroform: Hydrogen-Bond Structure and Its Signature on the Infrared Absorption Spectrum

The intramolecular hydrogen-bond structure of stereoselectively synthesized syn-tetrol and anti-tetrol dissolved in deuterated chloroform is investigated via a mixed quantum-classical molecular dynamics simulation. An extensive conformational analysis is performed in order to determine the dominant conformations, the distributions among them, and their sensitivity to the method for assigning partial charges (RESP vs AM1-BCC). The signature of the conformational distribution and method of assigning partial charges on the infrared absorption spectra is analyzed in detail. The relationship between the spectra and the underlying hydrogen-bond structure is elucidated.

A Mixed Quantum-Classical Molecular Dynamics Study of anti-Tetrol and syn-Tetrol Dissolved in Liquid Chloroform II: Infrared Emission Spectra, Vibrational Excited-State Lifetimes, and Nonequilibrium Hydrogen-Bond Dynamics

The effect of vibrational excitation and relaxation of the hydroxyl stretch on the hydrogen-bond structure and dynamics of stereoselectively synthesized syn-tetrol and anti-tetrol dissolved in deuterated chloroform are investigated via a mixed quantum-classical molecular dynamics simulation. Emphasis is placed on the changes in hydrogen-bond structure upon photoexcitation and the nonequilibrium hydrogen-bond dynamics that follows the subsequent relaxation from the excited to the ground vibrational state. The propensity to form hydrogen bonds is shown to increase upon photoexcitation of the hydroxyl stretch, thereby leading to a sizable red-shift of the infrared emission spectra relative to the corresponding absorption spectra. The vibrational excited state lifetimes are calculated within the framework of Fermi's golden rule and the harmonic-Schofield quantum correction factor, and found to be sensitive reporters of the underlying hydrogen-bond structure. The energy released during the relaxation from the excited to the ground state is shown to break hydrogen bonds involving the relaxing hydroxyl. The spectral signature of this nonequilibrium relaxation process is analyzed in detail.

Computational Study of the One‐ and Two‐Dimensional Infrared Spectra of a Proton‐Transfer Mode in a Hydrogen‐Bonded Complex Dissolved in a Polar Nanocluster

The signatures of nanosolvation on the one- and two-dimensional (1D and 2D) IR spectra of a proton-transfer mode in a hydrogen-bonded complex dissolved in polar solvent molecule nanoclusters of varying size are elucidated by using mixed quantum-classical molecular dynamics simulations. For this particular system, increasing the number of solvent molecules successively from N=7 to N=9 initiates the transition of the system from a cluster state to a bulk-like state. Both the 1D and 2D IR spectra reflect this transition through pronounced changes in their peak intensities and numbers, but the time-resolved 2D IR spectra also manifest spectral features that uniquely identify the onset of the cluster-to-bulk transition. In particular, it is observed that in the 1D IR spectra, the relative intensities of the peaks change such that the number of peaks decreases from three to two as the size of the cluster increases from N=7 to N=9. In the 2D IR spectra, off-diagonal peaks are observed in the N=7 and N=8 cases at zero waiting time, but not in the N=9 case. It is known that there are no off-diagonal peaks in the 2D IR spectrum of the bulk version of this system at zero waiting time, so the disappearance of these peaks is a unique signature of the onset of bulk-like behavior. Through an examination of the trajectories of various properties of the complex and solvent, it is possible to relate the emergence of these off-diagonal peaks to an interplay between the vibrations of the complex and the solvent polarization dynamics.

Metal-Dependent Activity of Fe and Ni Acireductone Dioxygenases: How Two Electrons Reroute the Catalytic Pathway

Two virtually identical acireductone dioxygenases, ARD and ARD′, catalyze completely different oxidation reactions of the same substrate, 1,2-dihydroxy-3-keto-5-(methylthio)pentene, depending exclusively on the nature of the bound metal. Fe2+-dependent ARD′ produces the α-keto acid precursor of methionine and formate and allows for the recycling of methionine in cells. Ni2+-dependent ARD instead produces methylthiopropionate, CO, and formate, and exits the methionine salvage cycle. This mechanistic difference has not been understood to date but has been speculated to be due to the difference in coordination of the substrate to Fe2+versus Ni2+: forming a five-membered ring versus a six-membered ring, respectively, thus exposing different carbon atoms for the attack by O2. Here, using mixed quantum-classical molecular dynamics simulations followed by the density functional theory mechanistic investigation, we show that, contrary to the old hypothesis, both metals preferentially bind the substrate as a six-membered ring, exposing the exact same sites to the attack by O2. It is the electronic properties of the metals that are then responsible for the system following different reaction paths, to yield the respective products. We fully explain the puzzling metal-induced difference in functionality between ARD and ARD′ and, in particular, propose a new mechanism for ARD′. All results are in agreement with available isotopic substitution and other experimental data.

On the active site of mononuclear B1 metallo β-lactamases: a computational study

Metallo-β-lactamases (MβLs) are Zn(II)-based bacterial enzymes that hydrolyze β-lactam antibiotics, hampering their beneficial effects. In the most relevant subclass (B1), X-ray crystallography studies on the enzyme from Bacillus Cereus point to either two zinc ions in two metal sites (the so-called '3H' and 'DCH' sites) or a single Zn(II) ion in the 3H site, where the ion is coordinated by Asp120, Cys221 and His263 residues. However, spectroscopic studies on the B1 enzyme from B. Cereus in the mono-zinc form suggested the presence of the Zn(II) ion also in the DCH site, where it is bound to an aspartate, a cysteine, a histidine and a water molecule. A structural model of this enzyme in its DCH mononuclear form, so far lacking, is therefore required for inhibitor design and mechanistic studies. By using force field based and mixed quantum-classical (QM/MM) molecular dynamics (MD) simulations of the protein in aqueous solution we constructed such structural model. The geometry and the H-bond network at the catalytic site of this model, in the free form and in complex with two common β-lactam drugs, is compared with experimental and theoretical findings of CphA and the recently solved crystal structure of new B2MβL from Serratia fonticola (Sfh-I). These are MβLs from the B2 subclass, which features an experimentally well established mono-zinc form, in which the Zn(II) is located in the DCH site. From our simulations the εεδ and δεδ protomers emerge as possible DCH mono-zinc reactive species, giving a novel contribution to the discussion on the MβL reactivity and to the drug design process.

Molecular Basis of the Mechanism of Thiol Oxidation by Hydrogen Peroxide in Aqueous Solution: Challenging the SN2 Paradigm

The oxidation of cellular thiol-containing compounds, such as glutathione and protein Cys residues, is considered to play an important role in many biological processes. Among possible oxidants, hydrogen peroxide (H(2)O(2)) is known to be produced in many cell types as a response to a variety of extracellular stimuli and could work as an intracellular messenger. This reaction has been reported to proceed through a S(N)2 mechanism, but despite its importance, the reaction is not completely understood at the atomic level. In this work, we elucidate the reaction mechanism of thiol oxidation by H(2)O(2) for a model methanethiolate system using state of the art hybrid quantum-classical (QM-MM) molecular dynamics simulations. Our results show that the solvent plays a key role in positioning the reactants, that there is a significant charge redistribution in the first stages of the reaction, and that there is a hydrogen transfer process between H(2)O(2) oxygen atoms that occurs after reaching the transition state. These observations challenge the S(N)2 mechanism hypothesis for this reaction. Specifically, our results indicate that the reaction is driven by a tendency of the slightly charged peroxidatic oxygen to become even more negative in the product via an electrophilic attack on the negative sulfur atom. This is inconsistent with the S(N)2 mechanism, which predicts a protonated sulfenic acid and hydroxyl anion as stable intermediates. These intermediates are not found. Instead, the reaction proceeds directly to unprotonated sulfenic acid and water.

Sampling the Proton Transfer Reaction Coordinate in Mixed Quantum-Classical Molecular Dynamics Simulations

An umbrella sampling approach based on the vibrational energy gap is presented and examined for exploring the reaction coordinate for a proton transfer (PT) reaction. The technique exploits the fact that for a PT reaction the energy gap between the vibrational ground and excited states of the transferring proton reaches a minimum at the transition state. Umbrella sampling is used within mixed quantum-classical simulations to identify the transition state configurations and explore the reaction free energy curve and vibrationally nonadiabatic coupling. The method is illustrated by application to a model phenol-amine proton transfer reaction complex in a nanoconfined solvent. The results from this new umbrella sampling approach are consistent with those obtained from previous umbrella sampling calculations based on a collective solvent coordinate. This sampling approach further provides insight into the vibrationally nonadiabatic coupling for the proton transfer reaction and has potential for simulating vibrational spectra of PT reaction complexes in solution.

Watching the Solvation of Atoms in Liquids One Solvent Molecule at a Time

We use mixed quantum-classical molecular dynamics simulations and ultrafast transient hole-burning spectroscopy to build a molecular-level picture of the motions of solvent molecules around Na atoms in liquid tetrahydrofuran. We find that even at room temperature, the solvation of Na atoms occurs in discrete steps, with the number of solvent molecules nearest the atom changing one at a time. This explains why the rate of solvent relaxation differs for different initial nonequilibrium states, and reveals how the solvent helps determine the identity of atomic species in liquids.

A new electron-methanol molecule pseudopotential and its application for the solvated electron in methanol

A new electron-methanol molecule pseudopotential is developed and tested in the present paper. The formal development of the potential is based on quantum mechanical calculations on the electron-methanol molecule model in the static exchange approximation. The computational model includes a steep confining potential that keeps the otherwise unbound excess electron in the vicinity of the methanol molecule. Using the Phillips-Kleinman theorem we introduce a smooth pseudowave function of the excess electron with the exact eigenenergy and correct asymptotic behavior. The nonlocal potential energy operator of the model Hamiltonian is then replaced to a local potential that reproduces the ground-state properties of the excess electron satisfactorily. The pseudopotential is then optimized in an analytically simple functional form to fit this approximate local potential in conjunction with the point charges and the geometry of a classical, all-site methanol-methanol interaction potential. Of the adjustable parameters, the parameters for the carbon and the methyl hydrogen atoms are optimized, while those for the oxygen and the hydroxyl hydrogen are taken from a previous electron-water molecule pseudopotential. A polarization term is added to the potential a posteriori. The polarization parameters are chosen to reproduce the experimental position of the optical absorption spectrum of an excess electron in mixed quantum-classical molecular dynamics simulations. The energetic, structural and spectroscopic properties of the solvated electron in a methanol bath are simulated at 300 K and compared with previous solvated electron simulations and available experimental data.

Computational study of the signature of hydrogen-bond strength on the infrared spectra of a hydrogen-bonded complex dissolved in a polar liquid

The signature of hydrogen-bond strength on the one- and two-dimensional infrared spectra of the hydrogen-stretch in a hydrogen-bonded complex dissolved in a polar liquid was investigated via mixed quantum–classical molecular dynamics simulations. Non-Condon effects were found to intensify with increasing hydrogen-bond strength and to shift oscillator strength from the stable configurations that correspond to the ionic and covalent tautomers into unstable configurations that correspond to the transition-state between them. The transition-state peak is observed to blue shift and increase in intensity with increasing hydrogen-bond strength, and to dominate the spectra in the case of a strong hydrogen-bond. It is argued that the application of multidimensional infrared spectroscopy in the region of the transition-state peak can provide a uniquely direct probe of the molecular events underlying breaking and forming of hydrogen-bonds in the condensed phase.