Additive empirical force field for hexopyranose monosaccharides

We present an all-atom additive empirical force field for the hexopyranose monosaccharide form of glucose and its diastereomers allose, altrose, galactose, gulose, idose, mannose, and talose. The model is developed to be consistent with the CHARMM all-atom biomolecular force fields, and the same parameters are used for all diastereomers, including both the alpha- and beta-anomers of each monosaccharide. The force field is developed in a hierarchical manner and reproduces the gas-phase and condensed-phase properties of small-molecule model compounds corresponding to fragments of pyranose monosaccharides. The resultant parameters are transferred to the full pyranose monosaccharides, and additional parameter development is done to achieve a complete hexopyranose monosaccharide force field. Parametrization target data include vibrational frequencies, crystal geometries, solute-water interaction energies, molecular volumes, heats of vaporization, and conformational energies, including those for over 1800 monosaccharide conformations at the MP2/cc-pVTZ//MP2/6-31G(d) level of theory. Although not targeted during parametrization, free energies of aqueous solvation for the model compounds compare favorably with experimental values. Also well-reproduced are monosaccharide crystal unit cell dimensions and ring pucker, densities of concentrated aqueous glucose systems, and the thermodynamic and dynamic properties of the exocyclic torsion in dilute aqueous systems. The new parameter set expands the CHARMM additive force field to allow for simulation of heterogeneous systems that include hexopyranose monosaccharides in addition to proteins, nucleic acids, and lipids.     

On geometries of stacked and H-bonded nucleic acid base pairs determined at various DFT, MP2, and CCSD(T) levels up to the CCSD(T)/complete basis set limit level

The geometries and interaction energies of stacked and hydrogen-bonded uracil dimers and a stacked adeninecdots, three dots, centeredthymine pair were studied by means of high-level quantum chemical calculations. Specifically, standard as well as counterpoise-corrected optimizations were performed at second-order Moller-Plesset (MP2) and coupled cluster level of theory with single, double, and perturbative triple excitations [CCSD(T)] levels with various basis sets up to the complete basis set limit. The results can be summarized as follows: (i) standard geometry optimization with small basis set (e.g., 6-31G(*)) provides fairly reasonable intermolecular separation; (ii) geometry optimization with extended basis sets at the MP2 level underestimates the intermolecular distances compared to the reference CCSD(T) results, whereas the MP2/cc-pVTZ counterpoise-corrected optimization agrees well with the reference geometries and, therefore, is recommended as a next step for improving MP2/cc-pVTZ geometries; (iii) the stabilization energy of stacked nucleic acids base pairs depends considerably on the method used for geometry optimization, so the use of reliable geometries, such as counterpoise-corrected MP2/cc-pVTZ ones, is recommended; (iv) the density functional theory methods fail completely in locating the energy minima for stacked structures and when the geometries from MP2 calculations are used, the resulting stabilization energies are strongly underestimated; (v) the self-consistent charges-density functional tight binding method, with inclusion of the empirical dispersion energy, accurately reproduces interaction energies and geometries of dispersion-bonded (stacked) complexes; this method can thus be recommended for prescanning the potential energy surfaces of van der Waals complexes.     

Force field development for poly(dimethylsilylenemethylene) with the aid of ab initio calculations

The molecular geometries, conformational energies, and zero-point energies of di(trimethylsilylene)methylene have been determined from high-level quantum chemistry calculations. The results are further used in the parametrization of a classical potential energy function suitable for performing simulations of the corresponding polymer, namely, poly(dimethylsilylenemethylene). Di(trimethylsilylene)methylene geometrical parameter optimizations for a proper location of the global minimum and other local minima, constrained at certain dihedral and bond angles, were performed at both the B3LYP/6-311G and MP2(full)/6-311G levels of theory. The global minimum configuration is slightly displaced from a perfectly staggered geometry, approximately by 16.0 degrees, at both levels of theory. Molecular mechanics and Monte Carlo calculations for isolated polymer chains together with molecular dynamics runs for the modeled dimer provide very good results in terms of conformational and thermodynamic properties.     

Microwave spectrum and conformational composition of 3-fluoropropionitrile (FCH2CH2CN)

The microwave spectrum of 3-fluoropropionitrile, FCH(2)CH(2)C≡N, has been investigated in the whole 17-75 GHz spectral region. Selected portions of the spectrum in the 75-95 GHz have also been recorded. The microwave spectra of the ground state as well as of three vibrationally excited states of each of two conformers have been assigned. The spectra of the vibrationally excited states belong to the lowest torsional and bending vibrations. The F-C-C-C chain of atoms is exactly antiperiplanar in one of these rotamers and synclinal in the second conformer. The F-C-C-C dihedral angle is 65(2)° in the synclinal form. The energy difference between the two forms has been obtained from relative intensity measurements performed on microwave transitions. It was found that the antiperiplanar conformer is more stable than the synclinal form by 1.4(5) kJ/mol. It is argued that the gauche effect is a significant force in this compound. Quantum chemical calculations at the high CCSD(full)/cc-pVTZ, MP2(full)/cc-pVTZ, and B3LYP/cc-pVTZ levels of theory have been performed. Most, but not all, of the theoretical predictions are in good agreement with experiment.     

Conformational Analysis of Free and Bound Retinoic Acid

The conformational profiles of unbound all-trans and 9-cis retinoic acid (RA) have been determined using classical and quantum mechanical calculations. Sixty-six all-trans-RA (ATRA) and forty-eight 9-cis-RA energy minimum conformers were identified via HF/6-31G* geometry optimizations in vacuo. Their relative conformational energies were estimated utilizing the M06, M06-2x and MP2 methods combined with the 6-311+G(d,p), aug-cc-pVDZ and aug-cc-pVTZ basis sets, as well as complete basis set MP2 extrapolations using the latter two basis sets. Single-point energy calculations performed with the M06-2x density functional were found to yield similar results to MP2/CBS for the low-energy retinoic acid conformations. Not unexpectedly, the conformational propensities of retinoic acid were governed by the orientation and arrangement of the torsion angles associated with the polyene tail. We also used previously reported QM/MM X-ray refinement results on four ATRA-protein crystal structures plus one newly refined 9-cis-RA complex (PDB ID 1XDK) in order to investigate the conformational preferences of bound retinoic acid. In the re-refined RA conformers the conjugated double bonds are nearly coplanar, which is consistent with the global minimum identified by the Omega/QM method rather than the corresponding crystallographically determined conformations given in the PDB. Consequently, a 91.3% average reduction of the local strain energy in the gas phase, as well as 92.1% in PCM solvent, was observed using the QM/MM refined structures versus the PDB deposited RA conformations. These results thus demonstrate that our QM/MM X-ray refinement approach can significantly enhance the quality of X-ray crystal structures refined by conventional refinement protocols, thereby providing reliable drug-target structural information for use in structure-based drug discovery applications.     

Potential energy surfaces of SimOn cluster formation and isomerization

The reaction paths for formation and isomerization of a set of silica SimOn (m = 2,3, n = 1-5) nanoclusters have been investigated using second-order perturbation theory (MP2) with the 6-31G(d) basis set. The MP2/6-31G(d) calculations have predicted singlet ground states for all clusters excluding Si3O2. The total energies of the most important points on the potential energy surfaces (PES) have been determined using the completely renormalized (CR) singles and doubles coupled cluster method including perturbative triples, CR-CCSD(T) with the cc-pVTZ basis set. Although transition states have been located for many isomerization reactions, only for Si3O3 and Si3O4 have some transition states been found for the formation of a cluster from the separated reactants. In all other cases, the process of formation of SimOn clusters appears to proceed without potential energy barriers.     

Ab initio/GIAO-CCSD(T) study of structures, energies, and 13C NMR chemical shifts of C4H7(+) and C5H9(+) ions: relative stability and dynamic aspects of the cyclopropylcarbinyl vs bicyclobutonium ions

The structures and energies of the carbocations C 4H 7 (+) and C 5H 9 (+) were calculated using the ab initio method. The (13)C NMR chemical shifts of the carbocations were calculated using the GIAO-CCSD(T) method. The pisigma-delocalized bisected cyclopropylcarbinyl cation, 1 and nonclassical bicyclobutonium ion, 2 were found to be the minima for C 4H 7 (+) at the MP2/cc-pVTZ level. At the MP4(SDTQ)/cc-pVTZ//MP2/cc-pVTZ + ZPE level the structure 2 is 0.4 kcal/mol more stable than the structure 1. The (13)C NMR chemical shifts of 1 and 2 were calculated by the GIAO-CCSD(T) method. Based on relative energies and (13)C NMR chemical shift calculations, an equilibrium involving the 1 and 2 in superacid solutions is most likely responsible for the experimentally observed (13)C NMR chemical shifts, with the latter as the predominant equilibrating species. The alpha-methylcyclopropylcarbinyl cation, 4, and nonclassical bicyclobutonium ion, 5, were found to be the minima for C 5H 9 (+) at the MP2/cc-pVTZ level. At the MP4(SDTQ)/cc-pVTZ//MP2/cc-pVTZ + ZPE level ion 5 is 5.9 kcal/mol more stable than the structure 4. The calculated (13)C NMR chemical shifts of 5 agree rather well with the experimental values of C 5H 9 (+).     

Quantum mechanical calculations of conformationally relevant 1H and 13C NMR chemical shifts of calixarene systems

[graphs: see text] QM GIAO calculations of 13C and 1H chemical shift values of the ArCH2Ar group have been performed, using the hybrid DFT functional MPW1PW91 and the 6-31G(d,p) basis set, on some representative calixarenes and on a series of simplified calixarene models allowing derivation of chemical shift surfaces versus phi and chi dihedral angles. A good reproduction of experimental data was obtained. The applicability of chemical shift surfaces in the study of calixarene conformational features is illustrated.     

Quantum mechanical calculations of conformationally relevant 1H and 13C NMR chemical shifts of N-, O-, and S-substituted calixarene systems

QM GIAO calculations of (13)C and (1)H chemical shift values of the ArCH(2)Ar group in N-, O-, and S-substituted calixarene systems were performed with a hybrid DFT functional MPW1PW91 and 6-31G(d,p) basis set. A good reproduction of experimental data was obtained for some representative calixarenes and for a series of simplified calixarene models. This allowed the derivation of chemical shift surfaces versus phi and chi dihedral angles. The applicability of chemical shift surfaces in the study of calixarene conformational features is illustrated.     

An ab initio study of the Cl(P)+C2H6→C2H5+HCl abstraction reaction

Electronic energies, geometries, and harmonic vibration frequencies for the reactants, products, and transition state for the Cl(³P)+C₂H₆→C₂H₅+HCl abstraction reaction were evaluated at the HF and MP2 levels using several correlation consistent polarized-valence basis sets. Single-point calculations at PMP2, MP4, QCISD(T), and CCSD(T) levels were also carried out. The values of the forward activation energies obtained at the MP4/cc-pVTZ, QCISD(T)/cc-pVTZ, and CCSD(T)/cc-pVTZ levels using the MP2/cc-pVTZ structures are equal to −0.1, −0.4, and −0.3 kcal/mol, respectively. The experimental value is equal to 0.3±0.2 kcal/mol. We found that the MP2/aug-cc-pVTZ adiabatic vibration energy for the reaction (−2.4 kcal/mol) agrees well with the experimental value −(2.2–2.6) kcal/mol. Rate constants calculated with the zeroth-order interpolated variational transition state (IVTST-0) method are in good agreement with experiment. In general, the theoretical rate constants differ from experiment by, at most, a factor of 2.6.     

Extending the treatment of backbone energetics in protein force fields: limitations of gas-phase quantum mechanics in reproducing protein conformational distributions in molecular dynamics simulations

Computational studies of proteins based on empirical force fields represent a powerful tool to obtain structure-function relationships at an atomic level, and are central in current efforts to solve the protein folding problem. The results from studies applying these tools are, however, dependent on the quality of the force fields used. In particular, accurate treatment of the peptide backbone is crucial to achieve representative conformational distributions in simulation studies. To improve the treatment of the peptide backbone, quantum mechanical (QM) and molecular mechanical (MM) calculations were undertaken on the alanine, glycine, and proline dipeptides, and the results from these calculations were combined with molecular dynamics (MD) simulations of proteins in crystal and aqueous environments. QM potential energy maps of the alanine and glycine dipeptides at the LMP2/cc-pVxZ//MP2/6-31G* levels, where x = D, T, and Q, were determined, and are compared to available QM studies on these molecules. The LMP2/cc-pVQZ//MP2/6-31G* energy surfaces for all three dipeptides were then used to improve the MM treatment of the dipeptides. These improvements included additional parameter optimization via Monte Carlo simulated annealing and extension of the potential energy function to contain peptide backbone phi, psi dihedral crossterms or a phi, psi grid-based energy correction term. Simultaneously, MD simulations of up to seven proteins in their crystalline environments were used to validate the force field enhancements. Comparison with QM and crystallographic data showed that an additional optimization of the phi, psi dihedral parameters along with the grid-based energy correction were required to yield significant improvements over the CHARMM22 force field. However, systematic deviations in the treatment of phi and psi in the helical and sheet regions were evident. Accordingly, empirical adjustments were made to the grid-based energy correction for alanine and glycine to account for these systematic differences. These adjustments lead to greater deviations from QM data for the two dipeptides but also yielded improved agreement with experimental crystallographic data. These improvements enhance the quality of the CHARMM force field in treating proteins. This extension of the potential energy function is anticipated to facilitate improved treatment of biological macromolecules via MM approaches in general.     

烷基亚硝胺非α-位代谢产物反应活性从头算研究

甲基正烷基亚硝胺通过非α-位代谢产生的亲电活性中心在化学致癌过程中具有重要作用. 采用从头计算方法比较了甲基乙基亚硝胺和甲基丙基亚硝胺β-位或γ-位代谢形成的硫酸酯、磷酸酯和葡萄糖醛酸酯在亚硝基的邻基参与作用下的反应活化能以及过渡态和中间体的分子结构. 计算结果表明,硫酸酯的反应活性明显高于磷酸酯和葡萄糖醛酸酯,而β-位代谢物受邻基参与的促进作用更大.     

Mechanisms for the deamination reaction of cytosine with H2O/OH(-) and 2H2O/OH(-): a computational study

Mechanisms for the deamination reaction of cytosine with H 2O/OH (-) and 2H 2O/OH (-) to produce uracil were investigated using ab initio calculations. Optimized geometries of reactants, transition states, intermediates, and products were determined at MP2 and B3LYP using the 6-31G(d) basis set and at B3LYP/6-31+G(d) levels of theory. Single point energies were also determined at MP2/G3MP2Large and G3MP2 levels of theory. Thermodynamic properties (Delta E, Delta H, and Delta G), activation energies, enthalpies, and free energies of activation were calculated for each reaction pathway investigated. Intrinsic reaction coordinate (IRC) analysis was performed to characterize the transition states on the potential energy surface. Seven pathways for the deamination reaction were found. All pathways produce an initial tetrahedral intermediate followed by several conformational changes. The final intermediate for all pathways dissociates to product via a 1-3 proton shift. The activation energy for the rate-determining step, the formation of the tetrahedral intermediate for pathway D, the only pathway that can lead to uracil, is 115.3 kJ mol (-1) at the G3MP2 level of theory, in excellent agreement with the experimental value (117 +/- 4 kJ mol (-1)).     

Computational study of noncovalent complexes between formamide and formic acid

The geometries and binding energies of 1:1, 1:2, and 1:4 formic acid-formamide complexes (FA-FMA) are calculated by quantum chemical procedures. Vibrational spectra and intermolecular distances of the most stable FA-FMA dimers as well as the influence of the basis set superposition error (BSSE) on the geometries and energies of the dimers are also discussed. All FA-FMA dimers are optimized at the B3LYP/cc-pVTZ, the MP2/cc-pVDZ, aug-cc-pVDZ, cc-pVTZ, and aug-cc-pVTZ levels of theory to study the influence of the level of theory on the calculated geometries and energies. CCSD(T)/cc-pVTZ single-point calculations at the MP2/aug-cc-pVTZ-optimized geometries were performed as reference for estimating the quality of lower level calculations. These calculations allow us to qualitatively describe the competition between different types of hydrogen-bonding interactions in FA-FMA complexes. FA-FMA dimers are compared to other formamide complexes and to the FA-FMA crystal structure.     

Conventional strain energies of 1,2-dihydroazete, 2,3-dihydroazete, 1,2-dihydrophosphete, and 2,3-dihydrophosphete

The conventional strain energies of 1,2-dihydroazete, 2,3-dihydroazete, 1,2-dihydrophosphete, and 2,3-dihydrophosphete are determined within the isodesmic, homodesmotic, and hyperhomodesmotic models. Optimum equilibrium geometries, harmonic vibrational frequencies, and corresponding electronic energies and zero-point vibrational energies are computed for all pertinent molecular systems using SCF theory, second-order perturbation theory, and density functional theory and employing the correlation consistent basis sets cc-pVDZ, cc-pVTZ, and cc-pVQZ. Single-point fourth-order perturbation theory, CCSD, and CCSD(T) calculations employing the cc-pVTZ and the cc-pVQZ basis sets are computed using the MP2/cc-pVTZ and MP2/cc-pVQZ optimized geometries, respectfully, to ascertain the contribution of higher order correlation. Three DFT functionals, B3LYP, wB97XD, and M06-2X, are employed to determine whether they can yield results similar to those obtained at the CCSD(T) level.     

Microwave spectrum and conformational composition of 1-vinylimidazole

The microwave spectrum of 1-vinylimidazole has been investigated in the 21-80 GHz spectral region. The spectra of two conformers have been assigned. One of these forms is planar, while the other is nonplanar with the imidazole ring and the vinyl group forming an angle of 15(4)° from coplanarity. The planar form is found to be 5.7(7) kJ/mol more stable than the nonplanar rotamer by relative intensity measurements. The spectra of 10 vibrationally excited states of the planar form and one excited-state spectrum of the nonplanar form were assigned. The vibrational frequencies of several of these states were determined by relative intensity measurements. The microwave work has been augmented by quantum chemical calculations at the CCSD/cc-pVTZ, MP2/cc-pVTZ, and B3LYP/cc-pVTZ levels of theory. The B3LYP calculations predict erroneously that both forms of 1-vinylimidazole are planar, whereas the MP2 and CCSD calculations correctly predict the existence of a planar and a nonplanar conformer of this compound.     

Conformational study of the structure of free 18-crown-6

A conformational search was performed for 18-crown-6 using the CONLEX method at the MM3 level. To have a more accurate energy order of the predicted conformations, the predicted conformations were geometry optimized at the HF/STO-3G level and the 198 lowest energy conformations, according to the HF/STO-3G energy order, were geometry optimized at the HF/6-31+G level. In addition, the 47 nonredundant lowest energy conformations, according to the MP2/6-31+G energy order at the HF/6-31+G optimized geometry, hereafter the MP2/6-31+G//HF/6-31+G energy order, were geometry optimized at the B3LYP/6-31+G level. According to the MP2/6-31+G//B3LYP/6-31+G energy order, three conformations had energies lower than the experimentally known Ci conformation of 18c6. At the MP2/6-31+G//B3LYP/6-31+G level, the S6 lowest energy conformation is more stable by 1.96 kcal/mol than this Ci conformation. This was confirmed by results at the MP2/6-31+G level with an energy difference of 1.84 kcal/mol. Comparison between the structure of the S6 conformation of 18c6 and the S4 lowest energy conformation of 12-crown-4, as well as other important conformations of both molecules, is made. It is concluded that the correlation energy is necessary to have an accurate energy order of the predicted conformations. A rationalization of the conformational energy order in terms of the hydrogen bonding and conformational dihedral angles is given. It is also suggested that to have a better energy order of the predicted conformations at the MM3 level, better empirical force fields corresponding to the hydrogen bond interactions are needed.     

Computational study of the deamination reaction of cytosine with H2O and OH-

The mechanism for the deamination reaction of cytosine with H(2)O and OH(-) to produce uracil was investigated using ab initio calculations. Optimized geometries of reactants, transition states, intermediates, and products were determined at RHF/6-31G(d), MP2/6-31G(d), and B3LYP/6-31G(d) levels and for anions at the B3LYP/6-31+G(d) level. Single-point energies were also determined at B3LYP/6-31+G(d), MP2/GTMP2Large, and G3MP2 levels of theory. Thermodynamic properties (DeltaE, DeltaH, and DeltaG), activation energies, enthalpies, and free energies of activation were calculated for each reaction pathway that was investigated. Intrinsic reaction coordinate analysis was performed to characterize the transition states on the potential energy surface. Two pathways for deamination with H(2)O were found, a five-step mechanism (pathway A) and a two-step mechanism (pathway B). The activation energy for the rate-determining steps, the formation of the tetrahedral intermediate for pathway A and the formation of the uracil tautomer for pathway B, are 221.3 and 260.3 kJ/mol, respectively, at the G3MP2 level of theory. The deamination reaction by either pathway is therefore unlikely because of the high barriers that are involved. Two pathways for deamination with OH(-) were also found, and both of them are five-step mechanisms. Pathways C and D produce an initial tetrahedral intermediate by adding H(2)O to deprotonated cytosine which then undergoes three conformational changes. The final intermediate dissociates to product via a 1-3 proton shift. Deamination with OH(-), through pathway C, resulted in the lowest activation energy, 148.0 kJ/mol, at the G3MP2 level of theory.     

A theoretical study on three conformational structures of 2,6‐distyrylpyridine

Ab initio methods at the levels HF/cc‐pVDZ, HF/6‐31G(d,p), MP2/cc‐pVDZ, and MP2/6‐31G(d,p), as well as methods based on density functional theory (DFT) employing the hybrid functional B3LYP with the basis sets cc‐pVDZ and 6‐31G(d,p), have been applied to study the conformers of 2,6‐distyrylpyridine. Bond distances, bond angles, and dihedral angles have been calculated at the B3LYP level. The calculated values were in good agreement with those measured by X‐ray diffraction analysis of 2,6‐distyrylpyridine. The values calculated using the Hartree‐Fock method and second‐order perturbation theory (MP2) were inconsistent. The optimized lowest‐energy geometries were calculated from the reported X‐ray structural data by the B3LYP/cc‐pVDZ method. Three conformations, A, B, and C, were proposed for 2,6‐distyrylpyridine. Calculations at the three levels of theory indicated that conformation A was the most stable structure, with conformations C and B being higher in energy by 1.10 and 2.57 kcal/mol, respectively, using the same method and basis function. The same trend in the relative energies of the three possible conformations was observed at the two levels of theory and with the different basis sets employed. The reported X‐ray data were utilized to optimize total molecular energy of conformation A at the different calculation levels. The bond lengths, bond angles, and dihedral angles were then obtained from the optimized geometries by ab initio methods and by applying DFT using the two basis functions cc‐pVDZ and 6‐31G(d,p). The values were analyzed and compared. The calculated total energies, the relative energies of the molecular orbitals, the gap between them, and the dipole moment for each conformational structure proposed for 2,6‐distyrylpyridine are also reported. © 2009 Wiley Periodicals, Inc. Int J Quantum Chem, 2010     

Potential energy surfaces for low-lying electronic states of SO2

The potential energy surfaces for the nine low-lying electronic states of SO₂ have been constructed by using the multi-reference second order perturbation theory (MRPT2) with the basis set cc-pVTZ. The optimized geometries and the adiabatic excitation energies of these states are in good agreement with experiments and previous calculations. The crossings and avoided crossings displayed in the potential energy surfaces are expounded.