Intermolecular Terms Research Articles

A Variational Linear-Scaling Framework to Build Practical, Efficient Next-Generation Orbital-Based Quantum Force Fields

We introduce a new hybrid molecular orbital/density-functional modified divide-and-conquer (mDC) approach that allows the linear-scaling calculation of very large quantum systems. The method provides a powerful framework from which linear-scaling force fields for molecular simulations can be developed. The method is variational in the energy, and has simple, analytic gradients and essentially no break-even point with respect to the corresponding full electronic structure calculation. Furthermore, the new approach allows intermolecular forces to be properly balanced such that non-bonded interactions can be treated, in some cases, to much higher accuracy than the full calculation. The approach is illustrated using the second-order self-consistent charge density-functional tight-binding model (DFTB2). Using this model as a base Hamiltonian, the new mDC approach is applied to a series of water systems, where results show that geometries and interaction energies between water molecules are greatly improved relative to full DFTB2. In order to achieve substantial improvement in the accuracy of intermolecular binding energies and hydrogen bonded cluster geometries, it was necessary to extend the DFTB2 model to higher-order atom-centered multipoles for the second-order self-consistent intermolecular electrostatic term. Using generalized, linear-scaling electrostatic methods, timings demonstrate that the method is able to calculate a water system of 3000 atoms in less than half of a second, and systems of up to one million atoms in only a few minutes using a conventional desktop workstation.

Effects of the dispersion interaction in liquid water

The properties of liquid water obtained from molecular simulations using a density functional theory (B3LYP) based force field [O. Akin-Ojo, F. Wang, J. Comput. Chem. 32 (2011) 453.] with and without an intermolecular dispersion interaction term are compared. Although it is small, the dispersion interaction is found to significantly affect some properties such as the density and the heat of vaporization but not others, for example, the atom–atom radial distribution functions and the infrared spectra. The inclusion of the dispersion interaction in the force field produces an increase in the computed self diffusivity comparable to that resulting from a ∼5 K rise in temperature.

Thermodynamics of liquids: standard molar entropies and heat capacities of common solvents from 2PT molecular dynamics

We validate here the Two-Phase Thermodynamics (2PT) method for calculating the standard molar entropies and heat capacities of common liquids. In 2PT, the thermodynamics of the system is related to the total density of states (DoS), obtained from the Fourier Transform of the velocity autocorrelation function. For liquids this DoS is partitioned into a diffusional component modeled as diffusion of a hard sphere gas plus a solid component for which the DoS(υ) → 0 as υ→ 0 as for a Debye solid. Thermodynamic observables are obtained by integrating the DoS with the appropriate weighting functions. In the 2PT method, two parameters are extracted from the DoS self-consistently to describe diffusional contributions: the fraction of diffusional modes, f, and DoS(0). This allows 2PT to be applied consistently and without re-parameterization to simulations of arbitrary liquids. We find that the absolute entropy of the liquid can be determined accurately from a single short MD trajectory (20 ps) after the system is equilibrated, making it orders of magnitude more efficient than commonly used perturbation and umbrella sampling methods. Here, we present the predicted standard molar entropies for fifteen common solvents evaluated from molecular dynamics simulations using the AMBER, GAFF, OPLS AA/L and Dreiding II forcefields. Overall, we find that all forcefields lead to good agreement with experimental and previous theoretical values for the entropy and very good agreement in the heat capacities. These results validate 2PT as a robust and efficient method for evaluating the thermodynamics of liquid phase systems. Indeed 2PT might provide a practical scheme to improve the intermolecular terms in forcefields by comparing directly to thermodynamic properties.

Proton Transport in Nafion117 Membranes by Molecular Dynamics

A molecular dynamics model was performed to study the proton transport of Nafion series membrane which is often used in low temperature fuel cells. The simulations intents to investigate the microstructure and the phenomenon of the proton transport processing. The model includes all-atom of main and side chains. The force field includes intermolecular Coulomb and Lennard-Jones terms and intra-molecular terms for harmonic bond stretch potentials, harmonic angle bend potentials and cosine terms for the torsions. The simulations were carried out in two cube systems at different temperature where water content differed from 5 and 10 water molecules per acid group in the polymer, respectively. The results showed that proton transport affected with water content by analyses of snap- shots of the MD simulation, the radial distribution functions between the sulfur atoms of SO3− groups, between the oxygen atoms of H3O+ ions, between oxygen atoms of water molecules and Nafion atoms at various stages.

Large-scale compensation of errors in pairwise-additive empirical force fields: comparison of AMBER intermolecular terms with rigorous DFT-SAPT calculations

The intermolecular interaction energy components for several molecular complexes were calculated using force fields available in the AMBER suite of programs and compared with Density Functional Theory-Symmetry Adapted Perturbation Theory (DFT-SAPT) values. The extent to which such comparison is meaningful is discussed. The comparability is shown to depend strongly on the intermolecular distance, which means that comparisons made at one distance only are of limited value. At large distances the coulombic and van der Waals 1/r(6) empirical terms correspond fairly well with the DFT-SAPT electrostatics and dispersion terms, respectively. At the onset of electronic overlap the empirical values deviate from the reference values considerably. However, the errors in the force fields tend to cancel out in a systematic manner at equilibrium distances. Thus, the overall performance of the force fields displays errors an order of magnitude smaller than those of the individual interaction energy components. The repulsive 1/r(12) component of the van der Waals expression seems to be responsible for a significant part of the deviation of the force field results from the reference values. We suggest that further improvement of the force fields for intermolecular interactions would require replacement of the nonphysical 1/r(12) term by an exponential function. Dispersion anisotropy and its effects are discussed. Our analysis is intended to show that although comparing the empirical and non-empirical interaction energy components is in general problematic, it might bring insights useful for the construction of new force fields. Our results are relevant to often performed force-field-based interaction energy decompositions.

Polarizable Intermolecular Potentials for Water and Benzene Interacting with Halide and Metal Ions

A complete derivation of polarizable intermolecular potentials based on high-level, gas-phase quantum-mechanical calculations is proposed. The importance of appreciable accuracy together with inherent simplicity represents a significant endeavor when enhancement of existing force fields for biological systems is sought. Toward this end, symmetry-adapted perturbation theory (SAPT) can provide an expansion of the total interaction energy into physically meaningful e.g. electrostatic, induction and van der Waals terms. Each contribution can be readily compared with its counterpart in classical force fields. Since the complexity of the different intermolecular terms cannot be fully embraced using a minimalist description, it is necessary to resort to polyvalent expressions capable of encapsulating overlooked contributions from the quantum-mechanical expansion. This choice results in consistent force field components that reflect the underlying physical principles of the phenomena. This simplified potential energy function is detailed and definitive guidelines are drawn. As a proof of concept, the methodology is illustrated through a series of test cases that include the interaction of water and benzene with halide and metal ions. In each case considered, the total energy is reproduced accurately over a range of biologically relevant distances.

Crystal structure studies using ab-initio potential functions from partitioned ab-initio MODPOT/VRDDO SCF energy calculations. I. N2 and CO2 test cases. II. Nitromethane, CH3NO2

Several years ago we embarked on a project to calculate optimal crystal-packing and crystal-structure parameters based on potential functions from energy-partitioned ab-initio intermolecular SCF calculations plus calculations of dispersion energy contributions. Our approach, aimed at evaluating the intermolecular interactions, such as in polymer propagation and in molecular crystals, and at analyzing the optimal crystal packing, is based on nonempirical ab-initio calculations for smaller molecular aggregates (monomers, dimers, trimers, etc.), partitioning the total SCF interaction energies into the different components, and then fitting these components individually to functional forms or when necessary recalculating or estimating explicitly for certain interaction components for each different unit cell dimension change. We derived and implemented a complete program, including higher order terms and special features necessary only for calculations on crystals. We also derived and implemented a complete CRYSTAL program which, given the crystal symmetry, allows us to vary and optimize the crystal-structure parameters. We carried out test crystal optimization calculations on N2 and on CO2 both to check out the entire CRYSTAL program package and to explore the extensive programs written for various ab-initio based methods for calculating reliable intermolecular interactions. We investigated such aspects as the comparison of the various intermolecular interaction terms calculated from the potential functions compared with those calculated by ab-initio wave functions. We also investigated for the CRYSTAL program how many surrounding unit cells had to be included for energy convergence. We investigated the crystal CH3NO2. We calculated the intermolecular interactions at more than 50 geometries, partitioned the interaction energies, and fit these to functional forms. The agreement of our calculated unit cell dimensions with experiment was excellent, within 1 to 2.8%. Furthermore, our method fits the partitioned ΔESCF results to atom-class-atom-class potential functions, which enables a general library of such potential functions to be built up for the molecular classes of interest.

Study of water clusters in the n = 2–34 size regime, based on the ABEEM/MM model

Various properties of water clusters in the n = 2–34 size regime with the change of cluster size have been systemically explored based on the newly developed flexible-body and charge-fluctuating ABEEM/MM water potential model. The ABEEM/MM water model is to take ABEEM charges of all atoms, bonds, and lone-pairs of water molecules into the intermolecular electrostatic interaction term in molecular mechanics. The computed correlating properties characterizing water clusters (H2O)n (n = 2–34) include optimal structures, structural parameters, ABEEM charge distributions, binding energies, hydrogen bonds, dipole moments, and so on. The study of optimal structures shows that the ABEEM/MM model can correctly predict the following important structural features, such as the transition from two-dimensional (from dimer to pentamer) to three-dimensional (for clusters larger than the hexamer) structures at hexamer region, the transition from cubes to cages at dodecamer (H2O)12, the transition from all-surface (all water molecules on the surface of the cluster) to one water-centered (one water molecule at the center of the cluster, fully solvated) structures at (H2O)17, the transition from one to two internal molecules in the cage at (H2O)33, and so on. The first three structural transitions are in good agreement with those obtained from previous work, while the fourth transition is different from that identified by Hartke. Subsequently, a systematic investigation of structural parameters, ABEEM charges, energetic properties, and dipole moments of water clusters with increasing cluster size can provide important reference for describing the objective trait of hydrogen bonds in water cluster system, and also provide a strong impetus toward understanding how the water clusters approach the bulk limit.

A graph theoretical approach for assessing bio-macromolecular complex structural stability

Fast and proper assessment of bio macro-molecular complex structural rigidity as a measure of structural stability can be useful in systematic studies to predict molecular function, and can also enable the design of rapid scoring functions to rank automatically generated bio-molecular complexes. Based on the graph theoretical approach of Jacobs et al. [Jacobs DJ, Rader AJ, Kuhn LA, Thorpe MF (2001) Protein flexibility predictions using graph theory. Proteins: Struct Funct Genet 44:150-165] for expressing molecular flexibility, we propose a new scheme to analyze the structural stability of bio-molecular complexes. This analysis is performed in terms of the identification in interacting subunits of clusters of flappy amino acids (those constituting regions of potential internal motion) that undergo an increase in rigidity at complex formation. Gains in structural rigidity of the interacting subunits upon bio-molecular complex formation can be evaluated by expansion of the network of intra-molecular inter-atomic interactions to include inter-molecular inter-atomic interaction terms. We propose two indices for quantifying this change: one local, which can express localized (at the amino acid level) structural rigidity, the other global to express overall structural stability for the complex. The new system is validated with a series of protein complex structures reported in the protein data bank. Finally, the indices are used as scoring coefficients to rank automatically generated protein complex decoys.

Origins of Resistance Conferred by the R292K Neuraminidase Mutation via Molecular Dynamics and Free Energy Calculations

Point mutations in the influenza virus enzyme neuraminidase (NA) have been reported that lead to dramatic loss of activity for known NA inhibitors including the FDA approved sialic acid mimics zanamivir and oseltamivir. A more complete understanding of the molecular basis for such resistance is a critical component toward development of improved next-generation drugs. In this study, we have used explicit solvent all-atom molecular dynamics simulations, free energy calculations (MM-GBSA), and residue-based decomposition to model binding of four ligands with NA from influenza virus subtype N9. The goal is to elucidate which structural and energetic properties change as a result of a mutation at position R292K. Computed binding free energies show strong correlation with experiment (r(2) = 0.76), and an examination of individual energy components reveal that changes in intermolecular Coulombic terms (ΔEcoul) best describe the variation in affinity with structure (r(2) = 0.93). H-bond populations also parallel the experimental ordering (r = -0.96, r(2) = 0.86) reinforcing the view that electrostatics modulate binding in this system. Notably, in every case, the simulation results correctly predict that loss of binding occurs as a result of the R292K mutation. Per-residue binding footprints reveal that changes in ΔΔEcoul for R292K-wildtype at position 292 parallel the change in experimental fold resistance energies (ΔΔGR292K-WT) with S03 < S00 < S02 < S01. The footprints also reveal that the most potent ligands have (1) less reliance on R292 for intrinsic affinity, (2) enhanced binding via residues E119, E227, and E277, and (3) flatter ΔEcoul and ΔH-bond profiles. Improved resistance for S03 appears to be a function of the ligand's larger guanidinium group which leads to an increased affinity for wildtype NA while at the same time a reduction in favorable interactions localized to R292. Overall, the computational results significantly enhance experimental observations through quantification of specific interactions which govern molecular recognition along the N9-ligand binding interface.

X-ray Diffraction Study on the Structure of Concentrated Aqueous Solutions Involving Alanine Molecules with Different Optical Activities

Abstract X-ray diffraction measurements on aqueous 2.5 mol % Dl-, l-, and d-alanine solutions in D2O were carried out at 26±2 °C in order to obtain information concerning the difference in the hydrogen-bonded structure between aqueous solutions involving amino acid molecules with different optical activities. The difference function, Δiinter(Q), between intermolecular interference term observed for Dl- and l-alanine and between Dl- and d-alanine solutions both exhibited a first peak at Q=1.6 Å−1, followed by oscillatory features extending to higher-Q region, implying that there is a difference in the intermolecular structure is present between these solutions. The difference distribution function, Δginter(r), obtained from the Fourier transform of the Δiinter(Q) between Dl- and l-, and between Dl- and d-alanine solutions showed an obvious negative peak at r=2.8 Å, which was attributed to the nearest neighbor hydrogen-bonded O···O interaction. The least squares fitting analysis of the observed Δiinter(Q) showed that the intermolecular O···O distance and the difference in the coordination number between Dl- and l-, and between Dl- and d-alanine solutions are 2.76(2) Å and −0.18(1), and 2.81(3) Å and −0.18(1), respectively. It was concluded that the intermolecular hydrogen-bonded network in aqueous l- and d-alanine solutions is stronger than that in the Dl-alanine solution.

Structure of amorphous NaCl–glucose system studied by X-ray diffraction and IR methods

X-ray diffraction and diffuse reflection IR spectroscopic measurements have been carried out on amorphous NaCl–glucose mixtures, a-(NaCl) x (glucose) 1− x , with x = 0, 0.05, 0.1 and 0.15, in order to obtain structural information on the intermolecular hydrogen-bonded interaction between glucose molecules affected by the presence of NaCl. The difference distribution function Δ g inter( r) was derived from the Fourier transform of the difference intermolecular interference term Δ i inter( Q) between X-ray intermolecular interference terms observed from amorphous NaCl-glucose and pure glucose samples. A negative peak appears at r ≈ 2.5 Å in the Δ g inter( r) observed for the 10 mol% NaCl sample, while, the Δ g inter( r) for the 15 mol% NaCl sample does not show pronounced negative peak. On the other hand, the IR spectra for the O–H stretching region (2300 ⩽ ν ⩽ 3800 cm −1) indicate that hydrogen bonds between glucose molecules are significantly collapsed in the samples containing 5–15 mol% NaCl. These results imply that the contribution from the Na +⋯Cl − contact ion pair is dominated in the 15 mol% NaCl sample.

Theoretical Studies of Energy Transfer Rates of Secondary Explosives

Understanding the mechanism of shock-induced chemical reaction in secondary explosives is necessary to pursue the development and the safe use of new explosives having high performance and low sensitivity. In an effort to understand the mechanism, the energy transfer rates of such secondary explosives as PETN(I), PETN(II), delta-HMX, alpha-HMX, beta-HMX, RDX, ANTA, DMN, and NM have been evaluated based on the formula derived by Fried and Ruggiero [Fried, L. E.; Ruggiero, A. J. J. Phys. Chem. 1994, 98, 9786]. The energy transfer rates were determined in terms of the density of vibrational states and the unharmonic vibron-phonon coupling term, which were calculated by using a flexible potential containing both intra- and intermolecular terms. For the secondary explosives, a good correlation was found between the energy transfer rates and the impact sensitivity. The energy transfer rates are several times faster for the explosives with higher sensitivity such as PETN, HMX, and RDX than those with lower sensitivity such as ANTA, DMN, and NM. The calculations presented suggest the energy transfer rate in secondary explosive crystals is a significant factor in their sensitivity and introduction of double bond, or hydrogen bonds, or caged structure into secondary explosives is expected to decrease the sensitivity.

Inter- and intra-molecular dynamics in liquid Se 2Br 2

Liquid Se 2Br 2 is well known to have a strong orientational correlation among the neighbouring molecules. Molecular dynamics of this liquid has been investigated by means of quasi-elastic neutron scattering. The temperature dependence of elastic component of the dynamic structure factor reveals that the dynamic structure consists of two scattering contributions in the wave number region where the first peak in the static structure factor lies. One is the non-correlating inter-molecular term and the other is the orientational correlation between strongly coupled neighbouring molecules.

Neutron Diffraction Studies on the Hydrogen-Bonded Structure of Water Molecules in Ion Exchange Resins

Abstract Neutron diffraction measurements were carried out for anion and cation exchange resins (AER and CER, respectively) absorbing H/D isotopically substituted water molecules in order to obtain structural information on the hydrogen-bonded network among water molecules within the resins. Partial structure factors, aHH(Q), aXH(Q), and aXX(Q) (X: atoms involved in the sample except for water-hydrogen atom), for AER and CER samples were respectively deduced from combined analyses of intermolecular interference terms observed for resins containing D2O (99.8% D), 0H2O (35.9% D), and 0–2H2O (67.9% D). The nearest neighbor distances, r(O···H) = 1.95(3) Å and r(H···H) = 2.40(3) Å, determined for the AER sample agree well with those reported for pure liquid water. On the other hand, considerably shorter intermolecular distances, r(O···H) = 1.78(3) Å and r(H···H) = 2.28(3) Å, were found for the CER sample, implying that the hydrogen-bonded network among water molecules in CER would be strongly affected by a highly acidic environment.

Intertube interactions in carbon nanotube bundles

Energetics of pairs and seven membered bundles formed by achiral as well as chiral carbon nanotubes is explored. The applied model Hamiltonian consists of a tight binding intratube part, an intermolecular hopping term, and a Lennard-Jones pair potential. Heterochiral aligned tube pairs are found the most stable energetically. Rotation of tubes around the tube axis in an aligned pair of identical tubes requires the smallest energy input if a cog-wheel rotation (disrotation) is exercised for heterochiral and achiral pairs or anti-cog-wheel rotation (conrotation) for a homochiral pair. Energetic preference for disrotation in heterochiral identical pairs manifests even for short (several nanometer long) tubes. Conrotation of homochiral pairs is hindered by end effects and becomes favorable only if modeling long (infinite) identical tubes. Among achiral tubes, those that possess a ${C}_{3n}$ axis with odd $n$ are shown to form close-packed bundles in an energetically optimal way. Similarly, among achiral tubes those having a ${C}_{3n}$ axis, $n$ either even or odd can form an energetically favorable close-packed bundle.

Free-energy force-field three-dimensional quantitative structure–activity relationship analysis of a set of p38-mitogen activated protein kinase inhibitors

The p38-mitogen-activated protein kinases (p38-MAPKs) belong to a family of serine-threonine kinases activated by pro-inflammatory or stressful stimuli that are known to be involved in several diseases. Their biological importance, related to the release of inflammatory pro-cytokines such as tumor necrosis factor-alpha (TNF-alpha) and interleukin-1 (IL-1), has generated many studies aiming at the development of selective inhibitors for the treatment of inflammatory diseases. In this work, we developed receptor-based three dimensional (3D) quantitative structure-activity relationship (QSAR) models for a series of 33 pyridinyl imidazole compounds [Liverton et al. (1999) 42:2180], using a methodology named free-energy force-field (FEFF) [Tokarski and Hopfinger (1997) 37:792], in which scaled intra- and intermolecular energy terms of the Assisted Model Building Energy Refinement (AMBER) force field combined with a hydration-shell solvation model are the independent variables used in the QSAR studies. Multiple temperature molecular-dynamics simulations (MDS) of ligand-protein complexes and genetic-function approximation (GFA) were employed using partial least squares (PLS) as the fitting functions to develop FEFF-3D-QSAR models for the binding process. The best model obtained in the FEFF-3D-QSAR receptor-dependent (RD) method shows the importance of the van der Waals energy change upon binding and the electrostatic energy in the interaction of ligands with the receptor. The QSAR equations described here show good predictability and may be regarded as representatives of the binding process of ligands to p38-MAPK. Additionally, we have compared the top FEFF-3D-QSAR model with receptor independent (RI) 4D-QSAR models developed in a recent study [Romeiro et al. (2005) 19:385].

Study of Lithium Cation in Water Clusters: Based on Atom-Bond Electronegativity Equalization Method Fused into Molecular Mechanics

We present a potential model for Li(+)-water clusters based on a combination of the atom-bond electronegativity equalization and molecular mechanics (ABEEM/MM) that is to take ABEEM charges of the cation and all atoms, bonds, and lone pairs of water molecules into the intermolecular electrostatic interaction term in molecular mechanics. The model allows point charges on cationic site and seven sites of an ABEEM-7P water molecule to fluctuate responding to the cluster geometry. The water molecules in the first sphere of Li(+) are strongly structured and there is obvious charge transfer between the cation and the water molecules; therefore, the charge constraint on the ionic cluster includes the charged constraint on the Li(+) and the first-shell water molecules and the charge neutrality constraint on each water molecule in the external hydration shells. The newly constructed potential model based on ABEEM/MM is first applied to ionic clusters and reproduces gas-phase state properties of Li(+)(H(2)O)(n) (n = 1-6 and 8) including optimized geometries, ABEEM charges, binding energies, frequencies, and so on, which are in fair agreement with those measured by available experiments and calculated by ab initio methods. Prospects and benefits introduced by this potential model are pointed out.

Intermolecular interaction from biorthogonal orbitals: new developments

New developments to compute the intermolecular interaction between two subsystems are proposed. They are based on biorthogonal orbitals and second quantization techniques. Test calculations have been performed on some geometrical configurations of a water dimer and on the strongly H-bonded OH−(H2O) system. All second-order intermolecular and mixed (inter–intra) terms are included. The results are compared with previously published SAPT calculations and with MP2 supermolecule treatment. They are all generally in good agreement. However, the description of the repulsive wall is sensitive to the quality of the basis set. The choice of a suitable basis set for the strongly H-bonded system is especially delicate.

A stable discretization of the lattice Boltzmann equation for simulation of incompressible two-phase flows at high density ratio

A stable discretization of the lattice Boltzmann equation (LBE) for non-ideal gases is presented for simulation of incompressible two-phase flows having high density and viscosity ratios. The stiffness of the discretized forcing terms in LBE for non-ideal gases is known to trigger severe numerical instability and restrict practical application of the LBE method. Use of a proper pressure updating scheme is also crucial to the stability of the LBE method because of non-negligible pressure variation across the phase interface. To deal with these issues, we propose a stable discretization scheme, which assumes the low Mach number approximation, and utilizes the stress and potential forms of the surface tension force, the incompressible transformation, and the consistent discretization of the intermolecular forcing terms. The proposed stable discretization scheme is applied to simulate 1-D advection equation with a source term, a stationary droplet, droplet oscillation and droplet splashing and deposition on a thin film at a density ratio of 1000 with varying Reynolds numbers. The numerical solutions of stationary and oscillatory droplets agree well with analytic solutions including the Laplace’s law. The time history of the spread factor of the liquid sheet emitted after the droplet impact also follows the known spreading power law.