Showing papers in "Journal of Computational Chemistry in 1995"

IMOMM: A new integrated ab initio + molecular mechanics geometry optimization scheme of equilibrium structures and transition states

Abstract: A new computational scheme integrating ab initio and molecular mechanics descriptions in different parts of the same molecule is presented. In contrast with previous approaches, this method is especially designed to allow the introduction of molecular mechanics corrections in full geometry optimizations concerning problems usually studied through ab initio calculations on model systems. The scheme proposed in this article intends to solve some of the systematic error associated with modeling through the use of molecular mechanics corrections. This method, which does not require any new parameter, evaluates explicitly the energy derivatives with respect to geometrical parameters and therefore has a straightforward application to geometry optimization. Examples of its performance on two simple cases are provided: the equilibrium geometry of cyclopropene and the energy barriers on SN2 reactions of alkyl chloride systems. Results are in satisfactory agreement with those of full ab initio calculations in both cases. © 1995 by John Wiley & Sons, Inc.

Remarks on the use of the apparent surface charges (ASC) methods in solvation problems: Iterative versus matrix‐inversion procedures and the renormalization of the apparent charges

Abstract: We present a formal and numerical comparison between the iterative and matrix-inversion approaches of the polarizable continuum model. The formal analysis shows completely the equivalence of the two approaches. Numerical equivalence is also recovered, introducing in both methods the proper boundary conditions on the apparent charge distribution. © 1995 John Wiley & Sons, Inc.

Application of the multimolecule and multiconformational RESP methodology to biopolymers: Charge derivation for DNA, RNA, and proteins

Abstract: We present the derivation of charges of ribo‐ and deoxynucleosides, nucleotides, and peptide fragments using electrostatic potentials obtained from ab initio calculations with the 6‐31G* basis set. For the nucleic acid fragments, we used electrostatic potentials of the four deoxyribonucleosides (A, G, C, T) and four ribonucleosides (A, G, C, U) and dimethylphosphate. The charges for the deoxyribose nucleosides and nucleotides are derived using multiple‐molecule fitting and restrained electrostatic potential (RESP) fits,1,2 with Lagrangian multipliers ensuring a net charge of 0 or ± 1. We suggest that the preferred approach for deriving charges for nucleosides and nucleotides involves allowing only C1′ and H1′ of the sugar to vary as the nucleic acid base, with the remainder of sugar and backbone atoms forced to be equivalent. For peptide fragments, we have combined multiple conformation fitting, previously employed by Williams3 and Reynolds et al.,4 with the RESP approach1,2 to derive charges for blocked dipeptides appropriate for each of the 20 naturally occuring amino acids. Based on our results for propyl amine,1,2 we suggest that two conformations for each peptide suffice to give charges that represent well the conformationally dependent electrostatic properties of molecules, provided that these two conformations contain different values of the dihedral angles that terminate in heteroatoms or hydrogens attached to heteroatoms. In these blocked dipeptide models, it is useful to require equivalent N—H and CO charges for all amino acids with a given net charge (except proline), and this is accomplished in a straightforward fashion with multiple‐molecule fitting. Finally, the application of multiple Lagrangian constraints allows for the derivation of monomeric residues with the appropriate net charge from a chemically blocked version of the residue. The multiple Lagrange constraints also enable charges from two or more molecules to be spliced together in a well‐defined fashion. Thus, the combined use of multiple molecules, multiple conformations, multiple Lagrangian constraints, and RESP fitting is shown to be a powerful approach to deriving electrostatic charges for biopolymers. © 1995 John Wiley & Sons, Inc.

The double cubic lattice method: Efficient approaches to numerical integration of surface area and volume and to dot surface contouring of molecular assemblies

Abstract: The double cubic lattice method (DCLM) is an accurate and rapid approach for computing numerically molecular surface areas (such as the solvent accessible or van der Waals surface) and the volume and compactness of molecular assemblies and for generating dot surfaces. The algorithm has no special memory requirements and can be easily implemented. The computation speed is extremely high, making interactive calculation of surfaces, volumes, and dot surfaces for systems of 1000 and more atoms possible on single-processor workstations. The algorithm can be easily parallelized. The DCLM is an algorithmic variant of the approach proposed by Shrake and Rupley (J. Mol. Biol., 79, 351–371, 1973). However, the application of two cubic lattices—one for grouping neighboring atomic centers and the other for grouping neighboring surface dots of an atom—results in a drastic reduction of central processing unit (CPU) time consumption by avoiding redundant distance checks. This is most noticeable for compact conformations. For instance, the calculation of the solvent accessible surface area of the crystal conformation of bovine pancreatic trypsin inhibitor (entry 4PTI of the Brookhaven Protein Data Bank, 362-point sphere for all 454 nonhydrogen atoms) takes less than 1 second (on a single R3000 processor of an SGI 4D/480, about 5 MFLOP). The DCLM does not depend on the spherical point distribution applied. The quality of unit sphere tesselations is discussed. We propose new ways of subdivision based on the icosahedron and dodecahedron, which achieve constantly low ratios of longest to shortest arcs over the whole frequency range. The DCLM is the method of choice, especially for large molecular complexes and high point densities. Its speed has been compared to the fastest techniques known to the authors, and it was found to be superior, especially when also taking into account the small memory requirement and the flexibility of the algorithm. The program text may be obtained on request. © 1995 by John Wiley & Sons, Inc.

Multidimensional free-energy calculations using the weighted histogram analysis method

Abstract: The recently formulated weighted histogram analysis method (WHAM)1 is an extension of Ferrenberg and Swendsen's multiple histogram technique for free-energy and potential of mean force calculations. As an illustration of the method, we have calculated the two-dimensional potential of mean force surface of the dihedrals gamma and chi in deoxyadenosine with Monte Carlo simulations using the all-atom and united-atom representation of the AMBER force fields. This also demonstrates one of the major advantages of WHAM over umbrella sampling techniques. The method also provides an analysis of the statistical accuracy of the potential of mean force as well as a guide to the most efficient use of additional simulations to minimize errors. © 1995 John Wiley & Sons, Inc.

Density functional theory and molecular clusters

Abstract: Density functional theory (DFT) methods, including nonlocal density gradient terms in the exchange and correlation energy functionals, were applied to various types of molecular clusters: H-bonded, ionic, electrostatic, and London. Reliable results on the structure and stabilization energy were obtained for the first two types of cluster as long as Becke3LYP and Becke3P86 functionals and basis sets of at least DZ + P quality were used. DFT methods with currently available functionals failed completely, however, for London-type clusters, for which no minimum was found on the potential energy surfaces. DFT interaction energy exhibits the same basis set extension dependence as the Hartree-Fock (HF) interaction energy. Therefore, the Boys-Bernardi function counterpoise procedure should be employed for elimination of the DFT basis set superposition error. © 1995 John Wiley & Sons, Inc.

Harmonic analysis of large systems. I. Methodology

Abstract: Methods have been developed for the determination of vibrational frequencies and normal modes of large systems in the full conformational space (including all degrees of freedom) and in a reduced conformational space (reducing the number of degrees of freedom). The computational method, which includes Hessian generation and storage, full and iterative diagonalization techniques, and the refinement of the results, is presented. A method is given for the quasiharmonic analysis and the reduced basis quasiharmonic analysis. The underlying principle is that from the atomic fluctuations, an effective harmonic force field can be determined relative to the dynamic average structure. Normal mode analysis tools can be used to characterize quasiharmonic modes of vibration. These correspond to conventional normal modes except that anharmonic effects are included. Numerous techniques for the analyses of vibrational frequencies and normal modes are described. Criteria for the analysis of the similarity of low-frequency normal modes is presented. The approach to determining the natural frequencies and normal modes of vibration described here is general and applicable to any large system. © 1995 John Wiley & Sons, Inc.

Poling: Promoting conformational variation

Abstract: This article introduces several methods of assessing the extent to which a collection of conformations represents or covers conformational space. It also describes poling: a novel technique for promoting conformational variation that can be applied to any method of conformational analysis that locally minimizes a penalty or energy function. The function being minimized is modified to force similar conformers away from each other. The method is independent of the origin of the initial conformers and of the particular minimization method used. It is found that, with the modification of the penalty function, clustering of the resulting conformers is generally unnecessary because the conformers are forced to be dissimilar. The functional form of the poling function is presented, and the merits are discussed with reference to (1) efficacy at promoting variation and (2) perturbation of the unmodified function. Results will be presented using conformers obtained from distance geometry with and without poling. It will be shown that the addition of poling eliminates much redundancy in conformer generation and improves the coverage of the conformational space. © 1995 by John Wiley & Sons, Inc.

Numerical solution of the nonlinear Poisson–Boltzmann equation: Developing more robust and efficient methods

Abstract: We present a robust and efficient numerical method for solution of the nonlinear Poisson-Boltzmann equation arising in molecular biophysics. The equation is discretized with the box method, and solution of the discrete equations is accomplished with a global inexact-Newton method, combined with linear multilevel techniques we have described in an article appearing previously in this journal. A detailed analysis of the resulting method is presented, with comparisons to other methods that have been proposed in the literature, including the classical nonlinear multigrid method, the nonlinear conjugate gradient method, and nonlinear relaxation methods such as successive overrelaxation. Both theoretical and numerical evidence suggests that this method will converge in the case of molecules for which many of the existing methods will not. In addition, for problems which the other methods are able to solve, numerical experiments show that the new method is substantially more efficient, and the superiority of this method grows with the problem size. The method is easy to implement once a linear multilevel solver is available and can also easily be used in conjunction with linear methods other than multigrid. © 1995 by John Wiley & Sons, Inc.

Accurate Modeling of the Intramolecular Electrostatic Energy of Proteins

Abstract: The (ϕ, ψ) energy surface of blocked alanine (N-acetyl–N′-methyl alanineamide) was calculated at the Hartree-Fock (HF)/6-31G* level using ab initio molecular orbital theory. A collection of six electrostatic models was constructed, and the term electrostatic model was used to refer to (1) a set of atomic charge densities, each unable to deform with conformation; and (2) a rule for estimating the electrostatic interaction energy between a pair of atomic charge densities. In addition to two partial charge and three multipole electrostatic models, this collection includes one extremely detailed model, which we refer to as nonspherical CPK. For each of these six electrostatic models, parameters—in the form of partial charges, atomic multipoles, or generalized atomic densities—were calculated from the HF/6-31G* wave functions whose energies define the ab initio energy surface. This calculation of parameters was complicated by a problem that was found to originate from the locking in of a set of atomic charge densities, each of which contains a small polarization-induced deformation from its idealized unpolarized state. It was observed that the collective contribution of these small polarization-induced deformations to electrostatic energy differences between conformations can become large relative to ab initio energy differences between conformations. For each of the six electrostatic models, this contribution was reduced by an averaging of atomic charge densities (or electrostatic energy surfaces) over a large collection of conformations. The ab initio energy surface was used as a target with respect to which relative accuracies were determined for the six electrostatic models. A collection of 42 more complete molecular mechanics models was created by combining each of our six electrostatic models with a collection of seven models of repulsion + dispersion + intrinsic torsional energy, chosen to provide a representative sample of functional forms and parameter sets. A measure of distance was defined between model and ab initio energy surfaces; and distances were calculated for each of our 42 molecular mechanics models. For most of our 12 standard molecular mechanics models, the average error between model and ab initio energy surfaces is greater than 1.5 kcal/mol. This error is decreased by (1) careful treatment of the nonspherical nature of atomic charge densities, and (2) accurate representation of electrostatic interaction energies of types 1—2 and 1—3. This result suggests an electrostatic origin for at least part of the error between standard model and ab initio energy surfaces. Given the range of functional forms that is used by the current generation of protein potential functions, these errors cannot be corrected by compensating for errors in other energy components. © 1995 by John Wiley & Sons, Inc.

A molecular mechanics/grid method for evaluation of ligand–receptor interactions

Abstract: We present a computational method for prediction of the conformation of a ligand when bound to a macromolecular receptor. The method is intended for use in systems in which the approximate location of the binding site is known and no large-scale rearrangements of the receptor are expected upon formation of the complex. The ligand is initially placed in the vicinity of the binding site and the atomic motions of the ligand and binding site are explicitly simulated, with solvent represented by an implicit solvation model and using a grid representation for the bulk of the receptor protein. These two approximations make the method computationally efficient and yet maintain accuracy close to that of an all-atom calculation. For the benzamidine/trypsin system, we ran 100 independent simulations, in many of which the ligand settled into the low-energy conformation observed in the crystal structure of the complex. The energy of these conformations was lower than and well-separated from that of others sampled. Extensions of this method are also discussed. © 1995 by John Wiley & Sons, Inc.

The relaxation of molecular crystal structures using a distributed multipole electrostatic model

Abstract: We describe a method for minimizing the lattice energy of molecular crystal structures, using a realistic anisotropic atom-atom model for the intermolecular forces. Molecules are assumed to be rigid, and the structure is described by the center of mass positions and orientational parameters for each molecule in the unit cell, as well as external strain parameters used to optimize the cell geometry. The resulting program uses a distributed multipole description of the electrostatic forces, which consists of sets of atomic multipoles (charge, dipole, quadrupole, etc.) to represent the lone pair, pi electron density, and other nonspherical features in the atomic charge distribution. Such ab initio based, electrostatic models are essential for describing the orientation dependence of the intermolecular forces, including hydrogen bonding, between polar molecules. Studies on a range of organic crystals containing hydrogen bonds are used to illustrate the use of this new crystal structure relaxation program, DMAREL, and show that it provides a promising new approach to studying the crystal packing of polar molecules. (C) 1995 by John Wiley and Sons, Inc.

Algorithms for constrained molecular dynamics

Abstract: In molecular dynamics simulations, the fastest components of the potential field impose severe restrictions on the stability and hence the speed of computational methods. One possibility for treating this problem is to replace the fastest components with algebraic length constraints. In this article the resulting systems of mixed differential and algebraic equations are studied. Commonly used discretization schemes for constrained Hamiltonian systems are discussed. The form of the nonlinear equations is examined in detail and used to give convergence results for the traditional nonlinear solution technique SHAKE iteration and for a modification based on successive overrelaxation (SOR). A simple adaptive algorithm for finding the optimal relaxation parameter is presented. Alternative direct methods using sparse matrix techniques are discussed. Numerical results are given for the new techniques, which have been implemented in the molecular modeling software package CHARMM and show as much as twofold improvement over SHAKE iteration. © 1995 John Wiley & Sons, Inc.

The fast multipole boundary element method for molecular electrostatics: An optimal approach for large systems

Abstract: We propose a fast implementation of the boundary element method for solving the Poisson equation, which approximately determines the electrostatic field around solvated molecules of arbitrary shape. The method presented uses computational resources of order O(N) only, where N is the number of elements representing the dielectric boundary at the molecular surface. The method is based on the Fast Multipole Algorithm by Rokhlin and Greengard, which is used to calculate the Coulombic interaction between surface elements in linear time. We calculate the solvation energies of a sphere, a small polar molecule, and a moderately sized protein. The values obtained by the boundary element method agree well with results from finite difference calculations and show a higher degree of consistency due to the absence of grid dependencies. The boundary element method can be taken to a much higher accuracy than is possible with finite difference methods and can therefore be used to verify their validity. © 1995 by John Wiley & Sons, Inc.

A simple yet accurate boundary element method for continuum dielectric calculations

Abstract: A simple yet accurate method for calculating electrostatic potentials using the boundary element continuum dielectric method is presented. It is shown that the limiting factor in accuracy is not the evaluation of integrals involving the interaction between boundary elements but rather a proper estimation of the self‐polarization of a patch upon itself. We derive a sum rule that allows us to calculate this important self‐polarization term in a self‐consistent and simple way. Intricate integration schemes used in previous treatments are consequently rendered unnecessary while concurrently achieving at least comparable accuracy over earlier methods. In some model systems for which analytic solutions are available, the computed surface polarization charge and reaction field energy are correct to better than six significant figures. An application of the method to the calculation of hydration free energies is presented. Good agreement with experimental values is obtained.

Parameterization and evaluation of a flexible water model

Abstract: The thermodynamic, dielectric, and dynamic properties of a newly parameterized flexible water model are studied using molecular dynamics simulations. The potential function developed is based on the popular simple point charge (SPC) rigid model with the addition of appropriate harmonic and anharmonic energy terms for stretching and bending. Care was taken to account for the self-polarization and gas-phase monomer energy corrections during the parameterization, which have typically been ignored in past studies. The results indicate that an increased Lennard-Jones repulsive coefficient and slightly scaled partial charges are required when adding flexibility to the rigid model potential to reliably reproduce the experimental density, energy, and O O radial distribution function of water at 298 K and 1 atm. Analysis of the power spectrum derived from the H-velocity autocorrelation function allowed the water potential to be evaluated further and refined by adjusting the valence forces to fit the vibrational frequencies of the gas and liquid. Once a consistent set of parameters was determined, the static dielectric properties of the water model were calculated at two temperatures using the reaction field method to treat long-range forces and correlations. The dielectric constant of 75 ± 7 calculated at 300 K is in good agreement with the experimental value of 78.5. The Kirkwood g factor was also examined for temperature dependence and showed the correct increasing behavior with decreasing T. As a final check of the water potential, the free energies of solvation of a flexible water molecule and neon were predicted using thermodynamic perturbation methods. The calculated solvation energies of −7.0 ± 0.8 for water and 2.7 ± 0.7 for neon are both consistent with the experimental values of −6.3 and 2.7 kcal/mol. Comparisons are made throughout the study with the results of previous rigid and flexible model simulations. © 1995 by John Wiley & Sons, Inc.

Energetics of Zn2+ binding to a series of biologically relevant ligands: A molecular mechanics investigation grounded on ab initio self-consistent field supermolecular computations

Abstract: Detailed investigations are performed of the binding energetics of Zn2+ to a series of neutral and anionic ligands making up the sidechains of amino acid residues of proteins, as well as ligands which can be involved in Zn2+ binding during enzymatic activation: imidazole, formamide, methanethiol, methanethiolate, methoxy, and hydroxy. The computations are performed using the SIBFA molecular mechanics procedure (SMM), which expresses the interaction energy under the form of four separate contributions related to the corresponding ab initio supermolecular ones: electrostatic, short-range repulsion, polarization, and charge transfer. Recent refinements to this procedure are first exposed. To test the reliability of this procedure in large-scale simulations of inhibitor binding to metalloenzyme cavities, we undertake systematic comparisons of the SMM results with those of recent large basis set ab initio self-consistent field (SCF) supermolecule computations, in which a decomposition of the total ΔE into its four corresponding components is done (N. Gresh, W. Stevens, and M. Krauss, J. Comp. Chem., 16, 843, 1995). For each complex, the evolution of each individual SMM energy component as a function of radial and in- and out-of-plane angular variations of the Zn2+ position reproduces with good accuracy the behavior of the corresponding SCF term. Computations performed subsequently on di- and oligoligated complexes of Zn2+ show that the SIBFA molecular mechanics (SMM) functionals, Epol and Ect, closely account for the nonadditive behaviors of the corresponding second-order energy contributions determined from the ab initio SCF calculations on these complexes and their nonlinear dependence on the number of ligands. Thus, the total intermolecular interaction energies computed with this procedure reproduce, with good accuracy, the corresponding SCF ones without the need for additional, extraneous terms in the intermolecular potential of polyligated complexes of divalent cations. © 1995 by John Wiley & Sons, Inc.

Improved methods for semiempirical solvation models

Abstract: We present improved algorithms for the SMx (x = 1, 1a, 2, 3) solvation models presented previously [see the overview in C. J. Cramer and D. G. Truhlar, J. Comp.-Aided Mol. Design, 6, 629 (1992)]. These models estimate the free energy of solvation by augmenting a semiempirical Hartree-Fock calculation on the solute with the generalized Born (GB) model for electric polarization of the solvent and a surface tension term based on solvent-accessible surface area. This article presents three improvements in the algorithms used to carry out such calculations, namely (1) an analytical accessible surface area algorithm, (2) a more efficient radial integration scheme for the dielectric screening computation in the GB model, and (3) a damping algorithm for updating the GB contribution to the Fock update during the iterations to achieve a self-consistent field. Improvements (1) and (2) decrease the computer time, and improvement (3) leads to more stable convergence. Improvement (2) removes a small systematic numerical error that was explicitly absorbed into the parameterization in the SMx models. Therefore, we have adjusted the parameters for one of the previous models to yield essentially identical performance as was obtained originally while simultaneously taking advantage of improvement (2). The resulting model is called SM2.1. The fact that we obtain similar results after removing the systematic quadrature bias attests to the robustness of the original parameterization. © 1995 by John Wiley & Sons, Inc.

Ab initio calculations of the ultraviolet resonance Raman spectra of uracil

Abstract: An equation been derived to calculate, ab initio, the frequencies and intensities of a resonant Raman spectrum from the transform theory of resonance Raman scattering. This equation has been used to calculate the intensities of the ultraviolet resonance Raman spectra from the first π-π* excited state of uracil and 1,3-dideuterouracil. The protocol for this calculation is as follows: (1) The force constant matrix elements in Cartesian coordinate space, the vibrational frequencies, and the minimum energy ground and excited state geometries of the molecule are calculated ab initio using the molecular orbital program Gaussian 92, (2) the force constants in Cartesian coordinates are transformed into force constants in the space of a set of 3N – 6 nonredundant symmetrized internal coordinates, (3) the G matrix is constructed from the energy minimized ground state Cartesian coordinates and the GFL = LΛ eigenvalue equation is solved in internal coordinate space, (4) the elements of the L and L−1 matrices are calculated, (5) the changes in all of the internal coordinates in going from the ground to the excited state are calculated, and (6) these results are used in combination with the transform theory of resonance Raman scattering to calculate the relative intensities of each of the 3N – 6 vibrations as a function of the exciting laser frequency. There are no adjustable parameters in this calculation, which reproduces the experimental frequencies and intensities with remarkable fidelity. This indicates that the Dushinsky rotation of the modes in the excited state of these molecules is not important and that the simplest form of the transform theory is adequate. © 1995 John Wiley & Sons, Inc.

Force field parameters for sulfates and sulfamates based on ab initio calculations: Extensions of AMBER and CHARMm fields

Abstract: Ab initio self-consistent field (SCF) Hartree-Fock calculations of sulfates ROSO3(−1) (R = Me, Et, i-Pr) and sulfamates RNHSO3(−1) (R = H, Me, Et, i-Pr) were performed at the 4-31G(*S*N) //3-21G(*S*N) basis set levels, where asterisks indicate d functions on sulfur and nitrogen atoms. These standard levels were determined by comparing calculation results with several basis sets up to MP2/6-31G*//6-31G*. Several conformations per compound were studied to obtain molecular geometries, rotational barriers, and potential derived point charges. In methyl sulfate, the rotational barrier around the CO bond is 1.6 kcal/mol at the MP2 level and 1.4 kcal/mol at the standard level. Its ground state has one of three HCOS torsion angles trans and one of three COSO torsion angles trans. Rotation over 60° around the single OS bond in the sulfate group costs 2.5 kcal/mol at the MP2 and 2.1 kcal/mol at the standard level. For ethyl sulfate, the calculated rotational barrier in going from the ground state, which has its CCOS torsion angle trans, to the syn-periplanar conformation (CCOS torsion angle cis) is 4.8 kcal/mol. However, a much lower barrier of 0.7 kcal/mol leads to a secondary gauchelike conformation about 0.4 kcal/mol above the ground state, with the CCOS torsion angle at 87.6°. Again, one of the COSO torsion angles is trans in the ground state, and the rotational barrier for a 60° rotation of the sulfate group amounts to 1.8 kcal/mol. For methyl sulfamate, the rotational barriers are 2.5 kcal/mol around the CN bond and 3.3 kcal/mol around the NS bond. This is noteworthy because sulfamate itself has a calculated rotational barrier around the NS bond of only 1.7 kcal/mol. These and other data were used to parameterize the well-known empirical force fields AMBER and CHARMm. When the new fields were tested by means of vibrational frequency calculations at the 6-31G*//6-31G* level for methyl sulfate, sulfamate, and methyl sulfamate ground states, the frequencies compared favorably with the AMBER and CHARMm calculated frequencies. The transferability of the force parameters to β-D-glucose-6-sulfate and isopropyl sulfate appears to be better than to isopropyl sulfamate. © 1995 by John Wiley & Sons, Inc.

An efficient approach to calculation of zero‐flux atomic surfaces and generation of atomic integration data

Abstract: A new robust method for variational determination of atomic zero-flux surfaces is presented. The zero-flux surface sheets are expressed in terms of variational trial functions in prolate spheroidal coordinates. The trial functions are optimized with a Newton procedure to satisfy the zero-flux condition on a grid. The data required for radial integrations are generated by an adaptive quadrature procedure that employs model electron densities and utilizes an original third-order algorithm for linear search. Results of test calculations involving variational determination of atomic surfaces are presented for a representative set of 20 molecules. The new approach is both less time consuming and substantially more accurate than the previously published algorithms. © 1995 John Wiley & Sons, Inc.

Calculation of molecular geometries, relative conformational energies, dipole moments, and molecular electrostatic potential fitted charges of small organic molecules of biochemical interest by density functional theory

Abstract: Density functional theory is tested on a large ensemble of model compounds containing a wide variety of functional groups to understand better its ability to reproduce experimental molecular geometries, relative conformational energies, and dipole moments. We find that gradient-corrected density functional methods with triple-ζ plus polarization basis sets reproduce geometries well. Most bonds tend to be approximately 0.015 A longer than the experimental results. Bond angles are very well reproduced and most often fall within a degree of experiment. Torsions are, on average, within 4 degrees of the experimental values. For relative conformational energies, comparisons with Hartree-Fock calculations and correlated conventional ab initio methods indicate that gradient-corrected density functionals easily surpass the Hartree-Fock approximation and give results which are nearly as accurate as MP2 calculations. For the 35 comparisons of conformational energies for which experimental data was available, the root mean square (rms) deviation for gradient-corrected functionals was approximately 0.5 kcal mol−1. Without gradient corrections, the rms deviation is 0.8 kcal mol−1, which is even less accurate than the Hartree-Fock calculations. Calculations with extended basis sets and with gradient corrections incorporated into the self-consistent procedure generate dipole moments with an rms deviation of 5%. Dipole moments from local density functional calculations, with more modest basis sets, can be scaled down to achieve roughly the same accuracy. In this study, all density functional geometries were generated by local density functional self-consistent calculations with gradient corrections added in a perturbative fashion. Such an approach generates results that are almost identical to the self-consistent gradient-corrected calculations, which require significantly more computer time. Timings on scalar and vector architectures indicate that, for moderately sized systems, our density functional implementation requires only slightly less computer resources than established Hartree-Fock programs. However, our density functional calculations scale much better and are significantly faster than their MP2 counterparts, whose results they approach. © 1995 John Wiley & Sons, Inc.

Flexible ligand docking without parameter adjustment across four ligand–receptor complexes

Abstract: Understanding molecular recognition is one of the fundamental problems in molecular biology. Computationally, molecular recognition is formulated as a docking problem. Ideally, a molecular docking algorithm should be computationally efficient, provide reasonably thorough search of conformational space, obtain solutions with reasonable consistency, and not require parameter adjustments. With these goals in mind, we developed DIVALI (Docking wIth eVolutionary AlgorIthms), a program which efficiently and reliably searches for the possible binding modes of a ligand within a fixed receptor. We use an AMBER-type potential function and search for good ligand conformations using a genetic algorithm (GA). We apply our system to study the docking of both rigid and flexible ligands in four different complexes. Our results indicate that it is possible to find diverse binding modes, including structures like the crystal structure, all with comparable potential function values. To achieve this, certain modifications to the standard GA recipe are essential. © 1995 John Wiley & Sons, Inc.

Molecular dynamics simulation with a continuum electrostatic model of the solvent

Abstract: The accuracy and simplicity of the Poisson-Boltzmann electrostatics model has led to the suggestion that it might offer an efficient solvent model for use in molecular mechanics calculations on biomolecules. We report a successful merger of the Poisson-Boltzmann and molecular dynamics approaches, with illustrative calculations on the small solutes dichloroethane and alanine dipeptide. The algorithm is implemented within the program UHBD. Computational efficiency is achieved by the use of rather coarse finite difference grids to solve the Poisson-Boltzmann equation. Nonetheless, the conformational distributions generated by the new method agree well with reference distributions obtained as Boltzmann distributions from energies computed with fine finite difference grids. The conformational distributions also agree well with the results of experimental measurements and conformational analyses using more detailed solvent models. We project that when multigrid methods are used to solve the finite difference problem and the algorithm is implemented on a vector supercomputer, the computation of solvent electrostatic forces for a protein of modest size will add only about 0.1 s computer time per simulation step relative to a vacuum calculation. © 1995 by John Wiley & Sons, Inc.

Ab initio calculation of molar volumes: Comparison with experiment and use in solvation models

Abstract: A simple and efficient procedure of calculating molecular volume (V) based on the Monte Carlo method is presented. The volume of a molecule is defined by the volume occupied by the 0.001‐au electron density envelope. We have employed this method to compute the molecular volumes (V) of a large selection of organic molecules and compare them with the corresponding molar volumes (V) measured in the liquid state. A strong correlation is found to exist between the V and V values (V/V ≈ 0.75). Using this linear relationship, the calculated molecular volume may provide an estimate of the cavity‐volume parameter in solvent‐effect calculations. As a chemical application of molecular volume, we have investigated the conformational equilibrium of 1,2‐difluoroethane in the liquid state using the self‐consistent reaction field theory. © 1995 by John Wiley & Sons, Inc. Copyright

On the accuracy of density functionals and their basis set dependence: an extensive study on the main group homonuclear diatomic molecules li2 to br2

Abstract: The equilibrium bond distances, harmonic frequencies, and bond dissociation energies of the 21 homonuclear diatomics Li2—F2, Na2—Cl2, and K2—Br2 have been determined using approximate density functional theory (DFT) employing various widely used functionals and basis sets ranging from single zeta to triple zeta plus polarization quality. The results are in general much less sensitive to the size of the basis set as in conventional ab initio molecular orbital (MO) theory, while the choice of the functional is of much more significance. For one basis set (6-311G*), the performance of the DFT-based calculations has been compared and found to be superior to Hartree-Fock (HF) Moller Plesset second order perturbation theory (MP2), or configuration interaction with single and double excitations (CISD) calculations. Particularly, no pathological cases, such as the group 2 dimers (Be2, Mg2, Ca2), are observed. © 1995 by John Wiley & Sons, Inc.

An examination of a density functional/molecular mechanical coupled potential

Abstract: In this article we describe the coupling of a density functional (DF) Hamiltonian with the molecular mechanics (MM) potential function AMBER. We examine a series of test cases in which we compare the binding energies and geometries of the complexes predicted by this coupled potential with those predicted by other theoretical methods and experiment to establish the relative accuracy of the DF/MM coupled potential. We find that the DF/MM coupled potential performs well in most cases studied and, in general, outperforms the semiempirical/MM approach. The interaction energies and structures obtained using this method appear to be insensitive to the use of nonlocal (NL) corrections to the DF method. The is fortuitous because the NL treatment is significantly more computationally expensive than the local treatment. However, NL corrections may be required to predict accurately the shape of potential energy surfaces that involve bond breaking and formation. The DF/MM method has also been applied to the determination of the solvation free energy for a series of ions using free-energy perturbation methods. The results obtained are good and can be improved by a simple scaling of the Lennard-Jones parameters for the ion in question. © 1995 by John Wiley & Sons, Inc.

Harmonic analysis of large systems. III. Comparison with molecular dynamics

Abstract: Atomic motions in bovine pancreatic trypsin inhibitor (BPTI), derived from molecular dynamics, harmonic analysis, and quasiharmonic analysis, are compared when a single protein model, energy parameters, and environment are employed. Molecular dynamics (MD) was carried out for 2 nanoseconds. An average structure was determined from the last nanosecond of the MD simulation, when no major structural changes were observed. This structure was used for several harmonic analysis calculations as well as for a reference structure for the quasiharmonic analysis, for both full basis and reduced basis sets. In contrast to the harmonic analysis results, the quasiharmonic reduced basis calculation using a spherical harmonics reduced basis provided good agreement with the full basis calculation, suggesting that when anharmonic effects are considered, BPTI can behave as a homogeneous object. An extensive analysis of the normal modes from a diverse set of 201 minimized MD simulation frames was performed. On only the sub-picosecond time scale were energy minima revisited after a transition to another state. This analysis shows that the dynamics average structure is not representative of the simulation frames in terms of energy and vibrational frequencies. For this model of BPTI, 42% of the motion (mean-squared fluctuation) can be attributed to harmonic limit behavior. A spectral analysis of the correlation function of deformation for a particular normal mode or quasiharmonic mode can be used to determine the time scales of motions which correspond to harmonic vibration, large-scale drift, or sharp transitions between local substrates. © 1995 John Wiley & Sons, Inc.

Development of optimized mst/scrf methods for semiempirical calculations: the mndo and pm3 hamiltonians

Abstract: The self-consistent reaction field (SCRF) method proposed by Miertus, Scrocco, and Tomasi (MST) has been optimized for MNDO and PM3 semiempirical Hamiltonians. Different algorithms used to compute the molecular electrostatic potential (MEP) and different solute cavities have been investigated. The ability of the optimized models to reproduce experimental free energies of solvation and to mimic the solvent effect in several chemical processes has been compared with the ab initio and AM1 versions of the MST method as well as with experimental data. © 1995 by John Wiley & Sons, Inc.

The geometry of pyrazole: A test for ab initio calculations

Abstract: Ab initio calculations on the structure of pyrazole have been carried out at different levels of accuracy. At the Hartree‐Fock (HF) level, the performance of several basis sets, namely 3‐21G, 6‐31G, 6‐31G**, and 6–311G** was investigated. The influence of electron correlation effects also was studied by carrying out geometry optimizations at the MP2, MP4, and QCISD levels. The performance of a density functional method also was evaluated. We have also investigated the possible influence of the frozen core approximation on the final optimized geometry. Three different statistical analyses were considered in determining which geometry is closest to the experimental microwave geometry—namely Paul Curtin's diagrams, cluster analysis, and multidimensional scaling. From these analyses, we conclude that there is no asymptotic approach to the experimental geometry by increasing the quality of the theoretical model, although, as expected, the more reliable structures are those obtained at the MP2, MP4, and QCISD levels, as well as those obtained by the B3LYP density functional method. We have also found that the values of the rotational constants are a tight criterion to define the quality of a molecular geometry. © 1995 by John Wiley & Sons, Inc.