The large quadrupole of water molecules

Many quantum mechanical calculations indicate water molecules in the gas and liquid phase have much larger quadrupole moments than any of the common site models of water for computer simulations. Here, comparisons of multipoles from quantum mechanical∕molecular mechanical (QM∕MM) calculations at the MP2∕aug-cc-pVQZ level on a B3LYP∕aug-cc-pVQZ level geometry of a waterlike cluster and from various site models show that the increased square planar quadrupole can be attributed to the p-orbital character perpendicular to the molecular plane of the highest occupied molecular orbital as well as a slight shift of negative charge toward the hydrogens. The common site models do not account for the p-orbital type electron density and fitting partial charges of TIP4P- or TIP5P-type models to the QM∕MM dipole and quadrupole give unreasonable higher moments. Furthermore, six partial charge sites are necessary to account reasonably for the large quadrupole, and polarizable site models will not remedy the problem unless they account for the p-orbital in the gas phase since the QM calculations show it is present there too. On the other hand, multipole models by definition can use the correct multipoles and the electrostatic potential from the QM∕MM multipoles is much closer than that from the site models to the potential from the QM∕MM electron density. Finally, Monte Carlo simulations show that increasing the quadrupole in the soft-sticky dipole-quadrupole-octupole multipole model gives radial distribution functions that are in good agreement with experiment.     

A combined molecular dynamics+quantum mechanics method for investigation of dynamic effects on local surface structures

A combined molecular dynamics (MD)+quantum mechanics (QM) method for studying processes on ionic surfaces is presented. Through the combination of classical MD and ab initio embedded-cluster calculations, this method allows the modeling of surface processes involving both the structural and dynamic features of the substrate, even for large-scale systems. The embedding approach used to link the information from the MD simulation to the cluster calculation is presented, and rigorous tests have been carried out to ensure the feasibility of the method. The electrostatic potential and electron density resulting from our embedded-cluster model have been compared with periodic slab results, and confirm the satisfying quality of our embedding scheme as well as the importance of applying embedding in our combined MD+QM approach. We show that a highly accurate representation of the Madelung potential becomes a prerequisite when the embedded-cluster approach is applied to temperature-distorted surface snapshots from the MD simulation.     

Multipole moments for embedding potentials: Exploring different atomic allocation algorithms

Polarizable quantum mechanical (QM)/molecular mechanics (MM)‐embedding methods are currently among the most promising methods for computationally feasible, yet reliable, production calculations of localized excitations and molecular response properties of large molecular complexes, such as proteins and RNA/DNA, and of molecules in solution. Our aim is to develop a computational methodology for distributed multipole moments and their associated multipole polarizabilities which is accurate, computationally efficient, and with smooth convergence with respect to multipole order. As the first step toward this goal, we herein investigate different ways of obtaining distributed atom‐centered multipole moments that are used in the construction of the electrostatic part of the embedding potential. Our objective is methods that not only are accurate and computationally efficient, but which can be consistently extended with site polarizabilities including internal charge transfer terms. We present a new way of dealing with well‐known problems in relation to the use of basis sets with diffuse functions in conventional atomic allocation algorithms, avoiding numerical integration schemes. Using this approach, we show that the classical embedding potential can be systematically improved, also when using basis sets with diffuse functions, and that very accurate embedding potentials suitable for QM/MM embedding calculations can be acquired. © 2016 Wiley Periodicals, Inc.     

An extensible interface for QM/MM molecular dynamics simulations with AMBER

We present an extensible interface between the AMBER molecular dynamics (MD) software package and electronic structure software packages for quantum mechanical (QM) and mixed QM and classical molecular mechanical (MM) MD simulations within both mechanical and electronic embedding schemes. With this interface, ab initio wave function theory and density functional theory methods, as available in the supported electronic structure software packages, become available for QM/MM MD simulations with AMBER. The interface has been written in a modular fashion that allows straight forward extensions to support additional QM software packages and can easily be ported to other MD software. Data exchange between the MD and QM software is implemented by means of files and system calls or the message passing interface standard. Based on extensive tests, default settings for the supported QM packages are provided such that energy is conserved for typical QM/MM MD simulations in the microcanonical ensemble. Results for the free energy of binding of calcium ions to aspartate in aqueous solution comparing semiempirical and density functional Hamiltonians are shown to demonstrate features of this interface. © 2013 Wiley Periodicals, Inc.     

Quantum mechanics/molecular mechanics minimum free-energy path for accurate reaction energetics in solution and enzymes: sequential sampling and optimization on the potential of mean force surface

To accurately determine the reaction path and its energetics for enzymatic and solution-phase reactions, we present a sequential sampling and optimization approach that greatly enhances the efficiency of the ab initio quantum mechanics/molecular mechanics minimum free-energy path (QM/MM-MFEP) method. In the QM/MM-MFEP method, the thermodynamics of a complex reaction system is described by the potential of mean force (PMF) surface of the quantum mechanical (QM) subsystem with a small number of degrees of freedom, somewhat like describing a reaction process in the gas phase. The main computational cost of the QM/MM-MFEP method comes from the statistical sampling of conformations of the molecular mechanical (MM) subsystem required for the calculation of the QM PMF and its gradient. In our new sequential sampling and optimization approach, we aim to reduce the amount of MM sampling while still retaining the accuracy of the results by first carrying out MM phase-space sampling and then optimizing the QM subsystem in the fixed-size ensemble of MM conformations. The resulting QM optimized structures are then used to obtain more accurate sampling of the MM subsystem. This process of sequential MM sampling and QM optimization is iterated until convergence. The use of a fixed-size, finite MM conformational ensemble enables the precise evaluation of the QM potential of mean force and its gradient within the ensemble, thus circumventing the challenges associated with statistical averaging and significantly speeding up the convergence of the optimization process. To further improve the accuracy of the QM/MM-MFEP method, the reaction path potential method developed by Lu and Yang [Z. Lu and W. Yang, J. Chem. Phys. 121, 89 (2004)] is employed to describe the QM/MM electrostatic interactions in an approximate yet accurate way with a computational cost that is comparable to classical MM simulations. The new method was successfully applied to two example reaction processes, the classical SN2 reaction of Cl-+CH3Cl in solution and the second proton transfer step of the reaction catalyzed by the enzyme 4-oxalocrotonate tautomerase. The activation free energies calculated with this new sequential sampling and optimization approach to the QM/MM-MFEP method agree well with results from other simulation approaches such as the umbrella sampling technique with direct QM/MM dynamics sampling, demonstrating the accuracy of the iterative QM/MM-MFEP method.     

QM:QM electronic embedding using Mulliken atomic charges: energies and analytic gradients in an ONIOM framework

An accurate first-principles treatment of chemical reactions for large systems remains a significant challenge facing electronic structure theory. Hybrid models, such as quantum mechanics:molecular mechanics (QM:MM) and quantum mechanics:quantum mechanics (QM:QM) schemes, provide a promising avenue for such studies. For many chemistries, including important reactions in materials science, molecular mechanics or semiempirical methods may not be appropriate, or parameters may not be available (e.g., surface chemistry of compound semiconductors such as indium phosphide or catalytic chemistry of transition metal oxides). In such cases, QM:QM schemes are of particular interest. In this work, a QM:QM electronic embedding model within the ONIOM (our own N-layer integrated molecular orbital molecular mechanics) extrapolation framework is presented. To define the embedding potential, we choose the real-system low-level Mulliken atomic charges. This results in a set of well-defined and unique embedding charges. However, the parametric dependence of the charges on molecular geometry complicates the energy gradient that is necessary for the efficient exploration of potential energy surfaces. We derive an efficient form for the forces where a single set of self-consistent field response equations is solved. Initial tests of the method and key algorithmic issues are discussed.     

Combined linear response quantum mechanics and classical electrodynamics (QM/ED) method for the calculation of surface-enhanced Raman spectra

A multiscale method is presented that allows for evaluation of plasmon-enhanced optical properties of nanoparticle/molecule complexes with no additional cost compared to standard electrodynamics (ED) and linear response quantum mechanics (QM) calculations for the particle and molecule, respectively, but with polarization and orientation effects automatically described. The approach first calculates the total field of the nanoparticle by ED using the finite difference time domain (FDTD) method. The field intensity in the frequency domain as a function of distance from the nanoparticle is calculated via a Fourier transform. The molecular optical properties are then calculated with QM in the frequency domain in the presence of the total field of the nanoparticle. Back-coupling due to dipolar reradiation effects is included in the single-molecule plane wave approximation. The effects of polarization and partial orientation averaging are considered. The QM/ED method is evaluated for the well-characterized test case of surface-enhanced Raman scattering (SERS) of pyridine bound to silver, as well as for the resonant Raman chromophore rhodamine 6G. The electromagnetic contribution to the enhancement factor is 10(4) for pyridine and 10(2) for rhodamine 6G.     

Contributions of pauli repulsions to the energetics and physical properties computed in QM/MM methods

Conventional combined quantum mechanical/molecular mechanical (QM/MM) methods lack explicit treatment of Pauli repulsions between the quantum‐mechanical and molecular‐mechanical subsystems. Instead, classical Lennard‐Jones (LJ) potentials between QM and MM nuclei are used to model electronic Pauli repulsion and long‐range London dispersion, despite the fact that the latter two are inherently of quantum nature. Use of the simple LJ potential in QM/MM methods can reproduce minimal geometries and energies of many molecular clusters reasonably well, as compared to full QM calculations. However, we show here that the LJ potential cannot correctly describe subtle details of the electron density of the QM subsystem because of the neglect of Pauli repulsions between the QM and MM subsystems. The inaccurate electron density subsequently affects the calculation of electronic and magnetic properties of the QM subsystem. To explicitly consider Pauli interactions with QM/MM methods, we propose a method to use empirical effective potentials on the MM atoms. The test case of the binding energy and magnetic properties of a water dimer shows promising results for the general application of effective potentials to mimic Pauli repulsions in QM/MM calculations. © 2013 Wiley Periodicals, Inc.     

Solvent effect on the absorption spectra of coumarin 120 in water: A combined quantum mechanical and molecular mechanical study

The solvent effect on the absorption spectra of coumarin 120 (C120) in water was studied utilizing the combined quantum mechanical∕molecular mechanical (QM∕MM) method. In molecular dynamics (MD) simulation, a new sampling scheme was introduced to provide enough samples for both solute and solvent molecules to obtain the average physical properties of the molecules in solution. We sampled the structure of the solute and solvent molecules separately. First, we executed a QM∕MM MD simulation, where we sampled the solute molecule in solution. Next, we chose random solute structures from this simulation and performed classical MD simulation for each chosen solute structure with its geometry fixed. This new scheme allowed us to sample the solute molecule quantum mechanically and sample many solvent structures classically. Excitation energy calculations using the selected samples were carried out by the generalized multiconfigurational perturbation theory. We succeeded in constructing the absorption spectra and realizing the red shift of the absorption spectra found in polar solvents. To understand the motion of C120 in water, we carried out principal component analysis and found that the motion of the methyl group made the largest contribution and the motion of the amino group the second largest. The solvent effect on the absorption spectrum was studied by decomposing it in two components: the effect from the distortion of the solute molecule and the field effect from the solvent molecules. The solvent effect from the solvent molecules shows large contribution to the solvent shift of the peak of the absorption spectrum, while the solvent effect from the solute molecule shows no contribution. The solvent effect from the solute molecule mainly contributes to the broadening of the absorption spectrum. In the solvent effect, the variation in C-C bond length has the largest contribution on the absorption spectrum from the solute molecule. For the solvent effect on the absorption spectrum from the solvent molecules, the solvent structure around the amino group of C120 plays the key role.     

Reaction path potential for complex systems derived from combined ab initio quantum mechanical and molecular mechanical calculations

Combined ab initio quantum mechanical and molecular mechanical calculations have been widely used for modeling chemical reactions in complex systems such as enzymes, with most applications being based on the determination of a minimum energy path connecting the reactant through the transition state to the product in the enzyme environment. However, statistical mechanics sampling and reaction dynamics calculations with a combined ab initio quantum mechanical (QM) and molecular mechanical (MM) potential are still not feasible because of the computational costs associated mainly with the ab initio quantum mechanical calculations for the QM subsystem. To address this issue, a reaction path potential energy surface is developed here for statistical mechanics and dynamics simulation of chemical reactions in enzymes and other complex systems. The reaction path potential follows the ideas from the reaction path Hamiltonian of Miller, Handy and Adams for gas phase chemical reactions but is designed specifically for large systems that are described with combined ab initio quantum mechanical and molecular mechanical methods. The reaction path potential is an analytical energy expression of the combined quantum mechanical and molecular mechanical potential energy along the minimum energy path. An expansion around the minimum energy path is made in both the nuclear and the electronic degrees of freedom for the QM subsystem internal energy, while the energy of the subsystem described with MM remains unchanged from that in the combined quantum mechanical and molecular mechanical expression and the electrostatic interaction between the QM and MM subsystems is described as the interaction of the MM charges with the QM charges. The QM charges are polarizable in response to the changes in both the MM and the QM degrees of freedom through a new response kernel developed in the present work. The input data for constructing the reaction path potential are energies, vibrational frequencies, and electron density response properties of the QM subsystem along the minimum energy path, all of which can be obtained from the combined quantum mechanical and molecular mechanical calculations. Once constructed, it costs much less for its evaluation. Thus, the reaction path potential provides a potential energy surface for rigorous statistical mechanics and reaction dynamics calculations of complex systems. As an example, the method is applied to the statistical mechanical calculations for the potential of mean force of the chemical reaction in triosephosphate isomerase.     

Simulation of Electronic Circular Dichroism of Nucleic Acids: From the Structure to the Spectrum

We present a quantum mechanical (QM) simulation of the electronic circular dichroism (ECD) of nucleic acids (NAs). The simulation combines classical molecular dynamics, to obtain the structure and its temperature‐dependent fluctuations, with a QM excitonic model to determine the ECD. The excitonic model takes into account environmental effects through a polarizable embedding and uses a refined approach to calculate the electronic couplings in terms of full transition densities. Three NAs with either similar conformations but different base sequences or similar base sequences but different conformations have been investigated and the results were compared with experimental observations; a good agreement was seen in all cases. A detailed analysis of the nature of the ECD bands in terms of their excitonic composition was also carried out. Finally, a comparison between the QM and the DeVoe models clearly revealed the importance of including fluctuations of the excitonic parameters and of accurately determining the electronic couplings. This study demonstrates the feasibility of the ab initio simulation of the ECD spectra of NAs, that is, without the need of experimental structural or electronic data.     

The merits of the frozen-density embedding scheme to model solvatochromic shifts

We investigate the usefulness of a frozen-density embedding scheme within density-functional theory [J. Phys. Chem. 97, 8050 (1993)] for the calculation of solvatochromic shifts. The frozen-density calculations, particularly of excitation energies have two clear advantages over the standard supermolecule calculations: (i) calculations for much larger systems are feasible, since the time-consuming time-dependent density functional theory (TDDFT) part is carried out in a limited molecular orbital space, while the effect of the surroundings is still included at a quantum mechanical level. This allows a large number of solvent molecules to be included and thus affords both specific and nonspecific solvent effects to be modeled. (ii) Only excitations of the system of interest, i.e., the selected embedded system, are calculated. This allows an easy analysis and interpretation of the results. In TDDFT calculations, it avoids unphysical results introduced by spurious mixings with the artificially too low charge-transfer excitations which are an artifact of the adiabatic local-density approximation or generalized gradient approximation exchange-correlation kernels currently used. The performance of the frozen-density embedding method is tested for the well-studied solvatochromic properties of the n-->pi(*) excitation of acetone. Further enhancement of the efficiency is studied by constructing approximate solvent densities, e.g., from a superposition of densities of individual solvent molecules. This is demonstrated for systems with up to 802 atoms. To obtain a realistic modeling of the absorption spectra of solvated molecules, including the effect of the solvent motions, we combine the embedding scheme with classical molecular dynamics (MD) and Car-Parrinello MD simulations to obtain snapshots of the solute and its solvent environment, for which then excitation energies are calculated. The frozen-density embedding yields estimated solvent shifts in the range of 0.20-0.26 eV, in good agreement with experimental values of between 0.19 and 0.21 eV.     

Theoretical study of N2O reduction by CO in Fe-BEA Zeolite.

Quantum mechanical (QM) and QM/molecular mechanics (MM) studies of the full catalytic cycle of N(2)O reduction by CO in Fe-BEA zeolite, that is, oxidation of BEA-Fe by N(2)O and reduction of BEA-Fe-alphaO by CO, is presented. A large QM cluster, representing half of the channel of the BEA zeolite, is used. The contribution of the MM embedding to the calculated activation energies is found to be negligible. The minimum-energy paths for N(2)O decomposition and reduction with CO are calculated using the nudged elastic band (NEB) method. Calculated and experimental activation energies are in good agreement. The two possible orientations for the gaseous molecules adsorbing on the Fe site that are found lead to different activation energies.     

Incorporation of a QM/MM buffer zone in the variational double self-consistent field method

The explicit polarization (X-Pol) potential is an electronic-structure-based polarization force field, designed for molecular dynamics simulations and modeling of biopolymers. In this approach, molecular polarization and charge transfer effects are explicitly treated by a combined quantum mechanical and molecular mechanical (QM/MM) scheme, and the wave function of the entire system is variationally optimized by a double self-consistent field (DSCF) method. In the present article, we introduce a QM buffer zone for a smooth transition from a QM region to an MM region. Instead of using the Mulliken charge approximation for all QM/MM interactions, the Coulombic interactions between the adjacent fragments are determined directly by electronic structure theory. The present method is designed to accelerate the speed of convergence of the total energy and charge density of the system.     

A novel quantum mechanical/molecular mechanical approach to the free energy calculation for isomerization of glycine in aqueous solution

The free energy change associated with the isomerization reaction of glycine in water solution has been studied by a hybrid quantum mechanical/molecular mechanical (QM/MM) approach combined with the theory of energy representation (QM/MM-ER) recently developed. The solvation free energies for both neutral and zwitterionic form of glycine have been determined by means of the QM/MM-ER simulation. The contributions of the electronic polarization and the fluctuation of the QM solute to the solvation free energy have been investigated. It has been found that the contribution of the density fluctuation of the zwitterionic solute is estimated as -4.2 kcal/mol in the total solvation free energy of -46.1 kcal/mol, while that of the neutral form is computed as -3.0 kcal/mol in the solvation free energy of -15.6 kcal/mol. The resultant free energy change associated with the isomerization of glycine in water has been obtained as -7.8 kcal/mol, in excellent agreement with the experimental data of -7.3 or -7.7 kcal/mol, implying the accuracy of the QM/MM-ER approach. The results have also been compared with those computed by other methodologies such as the polarizable continuum model and the classical molecular simulation. The efficiency and advantage of the QM/MM-ER method has been discussed.     

Uptake of phenol on aerosol particles

We present a study of the interaction between a phenol molecule and an aerosol particle. The aerosol particle is represented by a cluster of 128 water molecules. Using a classical approach, we present interaction energy surfaces for different relative distances and for three orientations of phenol relative to the particle. From the energy surfaces we find the reaction pathways with the largest interaction between the molecule and the particle. We use a quantum mechanics/molecular mechanics (QM/MM) method to calculate a potential energy curve for each reaction path. Coupled cluster methods are used for the part of the system described by quantum mechanics, while the part described by molecular mechanics is represented by a polarizable force field. We compare results obtained from the classical approach with the QM/MM results. Furthermore, we use the QM/MM results to calculate mass accommodation coefficients using a quantum-statistical (QM-ST) model and show how the mass accommodation coefficient depends on the relative orientation of phenol with respect to the aerosol particle.