Nonbonded Forces Research Articles

Accelerating molecular dynamics simulations by linear prediction of time series

We present a molecular dynamics simulation scheme which allows to speed up molecular dynamics simulations by linear prediction of force time series. The explicit calculation of nonbonding forces is periodically replaced by linear prediction from past values. Applying our method to liquid oxygen consisting of flexible molecules we obtained real speedups between 5.4 and 6.5, compared to conventional molecular dynamics simulations. Here only the bond-stretching forces were calculated at each time step. We demonstrate that essential dynamical quantities, such as the mean-square displacement and the velocity autocorrelation function, are preserved.

Influence of microscopic interactions on the spectra of polyacetylene

The influence of interatomic forces on the spectra of polyacetylene has been analyzed in detail by means of molecular dynamics simulations. Several simulations have been undertaken in which different terms have been added to the force field. Stretching, bending, torsional, and nonbonded forces have been subsequently incorporated in the simulations, in order to point out the specific spectral effects of every interaction. It has been obtained that the influence of all interatomic forces is qualitatively the same in both the trans and cis systems except for the torsions involving hydrogen atoms, which are relevant for the spectra parallel to the chain axis in the cis system, but not in the trans. Stretching and bending forces have a strong influence on the spectra parallel to the chain axis. The interchain forces do not significantly affect any of the spectra. Moreover, the results obtained in our simulations have been used to assign the experimental frequencies of polyacetylene to particular microscopic interactions.

An electrochemically controllable nanomechanical molecular system utilizing edge-to-face and face-to-face aromatic interactions.

[formula: see text] A new molecular system, 2,11-dithio[4,4]metametaquinocyclophane containing a quinone moiety, was designed and synthesized. As the quinone moiety can readily be converted into an aromatic pi-system (hydroquinone) upon reduction, the nanomechanical molecular cyclophane system exhibits a large flapping motion like a molecular flipper from the electrochemical redox process. The conformational changes upon reduction and oxidation are caused by changes of nonbonding interaction forces (devoid of bond formation/breaking) from the edge-to-face to face-to-face aromatic interactions and vice versa, respectively.

Quantum Mechanical Study of the Nonbonded Forces in Water−Methanol Complexes

The water-methanol dimer can adopt two possible configurations (WdM or MdW) depending on whether the water or the methanol acts as the hydrogen bond donor. The relative stability between the two configurations is less than 1 kcal/mol, and as a result, this dimer has been a challenging system to investigate using either theoretical or experimental techniques. In this paper, we present a systematic study of the dependence of the geometries, interaction energies, and harmonic frequencies on basis sets and treatment of electron correlation for the two configurations. At the highest theory level, MP2/aug-cc-pVQZ//MP2/aug-cc-pVTZ, interaction energies of -5.72 and -4.95 kcal/mol were determined for the WdM and MdW configurations, respectively, after correcting for basis set superposition error using the Boys-Bernardi counterpoise scheme. Extrapolating to the complete basis set limit resulted in interaction energies of -5.87 for WdM and -5.16 kcal/mol for MdW. The energy difference between the two configurations is larger than the majority of previously reported values, confirming that the WdM complex is preferred, in agreement with experimental observations. The effects that electron correlation have on the geometry were investigated by performing optimization at the MP2(full), MP4, and CCSD levels of theory. The approach trajectories for the formation of each dimer configuration are presented and the importance of these trajectories in the development of parameters for use in classical force fields is discussed.

Restrained electrostatic potential atomic partial charges for condensed-phase simulations of carbohydrates

Charges derived from fitting a classical Coulomb model to quantum mechanical molecular electrostatic potentials (so called ESP-charges) are frequently used in simulations of macromolecules. Simulational methods that use ESP-charges generally reproduce the geometries of hydrogen bonded complexes, despite the fact that these charges are known to overestimate the strengths of these interactions. Through the use of a restraint function during the fitting of the partial charges to the electrostatic potentials the magnitudes of the charges may be attenuated (so called RESP-charges). For the AMBER force field RESP-charges have been proposed for proteins and nucleic acids. Here we examine a novel approach for determining the RESP-charges for carbohydrates based on molecular dynamics (MD) simulations of crystal structures. During a simulation, the crystallographic unit cell geometry is sensitive to both inter-molecular non-bonded forces and internal torsional rotations. However, for polar molecules, and specifically carbohydrates, the crystal geometries are particularly sensitive to the set of partial atomic charges employed in the simulation. Thus, given a force field in which the van der Waals and torsion terms are well parameterized, it is possible to assess the suitability of a set of partial charges by monitoring the properties of the crystal during an MD simulation. We have examined several charge sets for use with the GLYCAM parameters for carbohydrate and glycoprotein simulations and found that a restraint weight of 0.01 gives the best agreement with the neutron diffraction structure of α- d-glucopyranose. Unrestrained ESP-charges performed poorly as did the charges obtained from Mulliken and distributed multipole analyses of the quantum mechanical HF/6-31G ∗ wavefunctions.

Interactions of dissolved organic matter with xenobiotic compounds: Molecular modeling in water

AbstractAn hypothesis for the structure of dissolved organic matter (DOM) in water is proposed. It is based on previously published humic acid and soil organic matter (SOM) models. Personal computer (PC)‐based molecular modeling and geometry optimization of DOM and humic/xenobiotic complexes in vacuo and water were performed using modern PC software in order to determine low energy conformations and to simulate site‐specific processes such as trapping and binding of biological and anthropogenic substances. Nanochemistry (10−9 m level) allows the evaluation of atomic and molecular space requirements, voids, inter‐ and intramolecular hydrogen bonds, and interactions with water, metal cations, and xenobiotics. The described modeling approach in general allows hydrophilic and hydrophobic reactions to be examined. Structural, molecular, and environmental properties of DOM and its xenobiotic complexes were determined by quantitative structure‐activity relationship software. Focal points were molecular properties, such as solvent accessibility as well as van der Waals surface areas and volumes, partial charges, hydration energy (peptides), hydrophobicity (log P), refractivity, and polarizabilities of humic/xenobiotic complexes were determined. Molecular mechanics calculations show that nonbonded forces (e.g., van der Waals) and hydrogen bonds were the main reasons for temporary immobility of xenobiotic substances retained in DOM. Preliminary experiments to simulate the acidity of water molecules by protonation‐enhanced reactions with polar xenobiotics (e.g., hydroxyatrazine) but left nonpolar substances (e.g., DDT) unchanged.

Derivation of Class II Force Fields. 7. Nonbonded Force Field Parameters for Organic Compounds

A set of nonbonded force field parameters consisting of atomic partial charges and van der Waals parameters has been derived by fitting experimental data for a broad set of organic compounds. The compounds in the fit spanned 11 functional groups: alcohols, aldehydes, amides, amines, carboxylic acids, esters, ethers, N-heterocycles, hydrocarbons, ketones, and sulfur compounds. The data consist of 136 crystal structures, 34 sublimation energies, and 63 gas-phase dipole moments. The nonbonded potential function is of the “9-6” form and is used in our class II CFF force field. This paper describes the fitting procedure and the quality of the fit as measured by comparing computed and experimental gas-phase dipole moments and crystal lattice energies. As a further test, each crystal structure was optimized. We report the accuracy of the computed crystal structures, including results obtained with force fields derived using two different combination rules. A summary is presented of the accuracy and consistency of predicted crystal lattice vectors and the lattice energies at the computed crystal structures. The results are discussed separately for each of the functional groups.

Molecular modelling of ( η6-arene)–Cr(CO) 3 complex systems : The barrier to rotation about the arene–metal bond axis

Based on the 1H-NMR calculated percentage of population of mono- and disubstituted arene–tricarbonyl–chromium complexes, the preferred conformation with regard to the internal rotation about the arene–chromium axis has been theoretically studied. Analytical force field expressions have also been developed to calculate the arene–(CO) 3 potential barrier and arene–Cr bond stretching. The torsional and stretching force constants are given. The new expressions to calculate force field were derived from reactivity constants and vibrational spectroscopy. The arene–Cr(CO) 3 intermolecular forces which are responsible for the change in equilibrium ratio of stable conformations were divided into two terms, the potential barrier of the arene–(CO) 3 rotating groups, and the arene–chromium stretching energy. The stable forms determined by these forces can be changed by great attractive and/or repulsive electrostatic and steric non-bonded interactions. These non-bonded forces are also depending on the arene–metal distance, which change linearly following the resonance effect of substituents on the aromatic ring. Arene–chromium bonding forces play the most important role in the arene–tricarbonyl–chromium studied here. There are appreciable changes in the arene–chromium distances from one complex to another, and it was found to be dependent of electron density on the aromatic ring. Electron-donating or electron-withdrawing substituent characters affect directly this distance and the energy difference between stable forms.

Longer time steps for molecular dynamics

Simulations of the dynamics of biomolecules have been greatly accelerated by the use of multiple time-stepping methods, such as the Verlet-I/r-RESPA (reversible reference system propagator algorithms) method, which is based on approximating “slow” forces as widely separated impulses. Indeed, numerical experiments have shown that time steps of 4 fs are possible for these slow forces but unfortunately have also shown that a long time step of 5 fs results in a dramatic energy drift. To overcome this instability, a symplectic modification of the impulsive Verlet-I/r-RESPA method has been proposed, called the mollified impulse method. The idea is that one modifies the slow part of the potential energy so that it is evaluated at “time averaged” values of the positions, and one uses the gradient of this modified potential for the slow part of the force. By filtering out excitations to the fastest motions, these averagings allow the use of longer time steps than does the impulse method. We introduce a new mollified method, Equilibrium, that avoids instability in a more effective manner than previous averaging mollified methods. Our experiments show that Equilibrium with a time step of 6 fs is as stable as the impulsive Verlet-I/r-RESPA method with a time step of 4 fs. We show that it may be necessary to include the effect of nonbonded forces in the averaging to make yet longer time steps possible. We also show that the slight modification of the potential has little effect on accuracy. For this purpose we compare self-diffusion coefficients and radial distribution functions against the Leapfrog method with a short time step (0.5 fs).

Acceleration of molecular mechanic simulation by parallelization and fast multipole techniques

Simulations of classical molecular dynamic (MD) systems can be sped up considerably by parallelizing the existing codes for distributed memory machines. In classical MD the CPU time is typically a function of the square of the number of atoms. The size of the molecular system which can be solved is therefore often limited by the CPU available. There are different approaches for reducing computation time. One consists in parallelizing sequential O( N 2) algorithms. The other is replacing the calculation of non-bonding forces by a less complex algorithm which can then be parallelized. We have generated a code (MEGADYN) for the simulation of MD of large simulation ensembles (up to 10 6 atoms) on the basis of classical force field methods. A reduction of complexity of the calculation of forces and energy down to O( N) was achieved by employing Greengards fast multipole method (FMM) to the Coulomb interaction. Within the framework of FMM the periodic boundary conditions are realized in a minimum image convention type manner. Thus MEGADYN can be used to simulate NVT as well as NPT ensembles.

Development of MD Engine: High-speed accelerator with parallel processor design for molecular dynamics simulations

Application of molecular dynamics (MD) simulations to large systems, such as biological macromolecules, is severely limited by the availability of computer resources. As the size of the system increases, the number of nonbonded forces (Coulombic and van der Waals interactions) to be evaluated increases as 𝒪(N2), where N is the number of particles in the system. The force evaluation consumes more than 99% of the CPU time in an MD simulation involving over 10,000 particles. Hence, the major target for reduction of the CPU time should be acceleration of the calculation of nonbonded forces. For this purpose, we developed a custom processor for calculating nonbonded interactions and a scalable plug-in machine (to a workstation), the MD Engine, in which numbers of the custom processors work in parallel. The processor has a pipeline architecture to calculate the total nonbonded force using the coordinates, electric charge, and species of each particle broadcast by the host computer. The force is calculated with sufficient accuracy for practical MD simulations. The processor also calculates virials simultaneously with forces for use in the calculation of pressure, accommodates periodic boundary conditions, and can be used in Ewald summations. An MD Engine system consisting of 76 processors calculates nonbonded interactions about 50 times faster than an UltraSPARC-I processor (Sun Ultra-2, 200 MHz) or an R10000 processor (SGI Origin 200, 180 MHz). On a Sun Ultra-2 workstation with a single UltraSPARC-I processor an MD simulation of a Ras p21 protein molecule immersed in a water sphere (13,258 particles) was accelerated by a factor of 48 using the MD Engine system. ©1999 John Wiley & Sons, Inc. J Comput Chem 20: 185–199, 1999

INTERACTIONS OF DISSOLVED ORGANIC MATTER WITH XENOBIOTIC COMPOUNDS: MOLECULAR MODELING IN WATER

An hypothesis for the structure of dissolved organic matter (DOM) in water is proposed. It is based on previously published humic acid and soil organic matter (SOM) models. Personal computer (PC)-based molecular modeling and geometry optimization of DOM and humic/xenobiotic complexes in vacuo and water were performed using modern PC software in order to determine low energy conformations and to simulate site-specific processes such as trapping and binding of biological and anthropogenic substances. Nanochemistry (10−9 m level) allows the evaluation of atomic and molecular space requirements, voids, inter- and intramolecular hydrogen bonds, and interactions with water, metal cations, and xenobiotics. The described modeling approach in general allows hydrophilic and hydrophobic reactions to be examined. Structural, molecular, and environmental properties of DOM and its xenobiotic complexes were determined by quantitative structure-activity relationship software. Focal points were molecular properties, such as solvent accessibility as well as van der Waals surface areas and volumes, partial charges, hydration energy (peptides), hydrophobicity (log P), refractivity, and polarizabilities of humic/xenobiotic complexes were determined. Molecular mechanics calculations show that nonbonded forces (e.g., van der Waals) and hydrogen bonds were the main reasons for temporary immobility of xenobiotic substances retained in DOM. Preliminary experiments to simulate the acidity of water molecules by protonation-enhanced reactions with polar xenobiotics (e.g., hydroxyatrazine) but left nonpolar substances (e.g., DDT) unchanged.

Changes in the Carbon−Carbon and Carbon−Nitrogen Bond Orders in the Macrocycle of Chlorophyll a upon Singlet and Triplet Excitation As Probed by Resonance Raman Spectroscopy of Natural-Abundance and Singly and Doubly Labeled Species with 15N, 13C, and 2H Isotopes

The Raman spectra of chlorophyll a (1700−1000 cm-1) in the S0, S1, and T1 states were recorded for the species of natural-abundance (NA) isotopic composition and those species that were totally labeled with the 15N, 13C, 13C + 15N, 2H, and 2H + 15N isotopes. An empirical normal-coordinate analysis using the S0 Raman lines of all the isotope species was performed to establish their assignments and to determine a set of force constants including the stretching, bending, nonbonded repulsive force constants and their cross-terms. The Ca−Cm and Ca−N stretching force constants in the S1 and T1 states were also determined, on the basis of the assumption that other force constants were the same as those in the S0 state, to achieve the classification of the observed Raman lines of the NA-, 13C-, and 2H-Chl a into two groups: those which shift largely upon the 15N-substitution and are associated with the Ca−N stretching vibrations and those which do not and are associated with the carbon−oxygen and carbon−carbon s...

Resonance‐Assisted O‐H ⃛ O Hydrogen Bonding: Its Role in the Crystalline Self‐Recognition of β‐Diketone Enols and its Structural and IR Characterization

AbstractThe fact that hydrogen bonding is normally stronger than other nonbonding attractive forces can be exploited for the rational design of molecular crystals with known packing features and specific physical properties (crystal engineering). In the present paper the problem of obtaining homodromous molecular chains controlled by strong O–H ⃛ O interactions is investigated, particular attention being paid to β‐chains, that is, infinite hydrogen‐bonded chains of β‐diketone enol fragments ⃛ OCCCOH ⃛, which are linked by stronger‐than‐usual resonance‐assisted hydrogen bonds (RAHBs) and are intrinsically interesting as prototypes of a large family of switching proton bistate molecular devices. Accordingly, the crystal and molecular structures of thirteen new compounds containing the 1,3‐cyclopentanedione and 1,3‐cyclohexanedione fragment (or their heterocyclic analogues) were determined, and most of them were found to give the expected β‐chain packing pattern. Comparison with literature data makes it possible to identify seven fundamental β‐chain patterns, which can be shown to be selected by reason of the relative encumbrances of the substituents. Furthermore, a general analysis of all functional groups able to form strong OH ⃛ O bonds reveals a semiquantitative correspondence between the OH ⃛ O measurable parameters (O ⃛ O, H ⃛ O and OH distances, and ṽ(OH) IR stretching frequencies) and the hydrogen bond energy EHB, and a hierarchy of chemical functionalities that are well characterized by limited EHB ranges and that, in decreasing order of energy, can direct the crystal packing process.

Protein-protein interfaces: architectures and interactions in protein-protein interfaces and in protein cores. Their similarities and differences.

Protein structures generally consist of favorable folding motifs formed by specific arrangements of secondary structure elements. Similar architectures can be adopted by different amino acids sequences, although the details of the structures vary. It has long been known that despite the sequence variability, there is a striking preferential conservation of the hydrophobic character of the amino acids at the buried positions of these folding motifs. Differences in the sizes of the side-chains are accommodated by movements of the secondary structure elements with respect to each other, leading to compact packing. Scanning protein-protein interfaces reveals that similar architectures are also observed at and around their interacting surfaces, with preservation of the hydrophobic character, although not to the same extent. The general forces that determine the origin of the native structures of proteins have been investigated intensively. The major non-bonded forces operating on a protein chain as it folds into a three-dimensional structure are likely to be packing, the hydrophobic effect, and electrostatic interactions. While the substantial hydrophobic forces lead to a compact conformation, they are also nonspecific and cannot serve as a guide to a conformationally unique structure. For the general folding problem, it thus appears that packing is a prime candidate for determining a particular fold. Specific hydrogen-bonding patterns and salt-bridges have also been proposed to play a role. Inspection of protein-protein interfaces reveals that the hallmarks governing single chain protein structures also determine their interactions, suggesting that similar principles underlie protein folding and protein-protein associations. This review focuses on some aspects of protein-protein interfaces, particularly on the architectures and their interactions. These are compared with those present in protein monomers. This task is facilitated by the recently compiled, non-redundant structural dataset of protein-protein interfaces derived from the crystallographic database. In particular, although current view holds that protein-protein interfaces and interactions are similar to those found in the conformations of single-chain proteins, this review brings forth the differences as well. Not only is it logical that such differences would exist, it is these differences that further illuminate protein folding on the one hand and protein-protein recognition on the other. These are also particularly important in considering inhibitor (ligand) design.

Validation of a reparameterized MM3 non-bonded force field for hydrocarbons: crystal lattice studies of C 60 and C 70 and adsorption of hydrocarbons onto graphite

Four hydrocarbon van der Waals potentials have been employed to obtain crystal unit cell dimensions and enthalpy sublimation values for two fullerenes, C 60 and C 70. The same force fields have also been employed in computing adsorption energies between a series of hydrocarbon molecules and a model graphite basal plane. The four hydrocarbon force fields are MM2, MM3, MM3hc, which is a new reparametrization of MM3 by the authors in earlier studies, and the potential of Williams and Starr. The MM3 and MM3hc force fields have been shown to offer the best agreement between the experimental and calculated crystal data for C 60 and C 70. In the graphite-adsorbate molecule computations, both charged and uncharged models for the four force fields have been tested to establish which model provides the most accurate structural and thermodynamic results. These calculations indicate that models with no charges are superior despite their deficiencies. For the adsorption of aliphatic molecules onto graphite, MM3hc yields closer agreement between experimental and calculated energies of adsorption than its counterparts, while its predictions for aromatic compounds are slightly less accurate than the MM2 and MM3 results. All force fields are similar in their predictions of structural arrangements of the adsorbate molecules on graphite. However, the computed structures disagree with experimental findings in several cases, underlining the inadequacy of simple analytical potential functions for graphite-hydrocarbon interactions.

Error evaluation in the design of a special‐purpose processor that calculates nonbonded forces in molecular dynamics simulations

AbstractSpecial‐purpose parallel machines that are plugged into a workstation to accelerate molecular dynamics (MD) simulations are attracting a considerable amount of interest. These machines comprise scalable homogeneous multiprocessors for calculating nonbonded forces (Coulombic and van der Waals forces), which consume more than 99% of the central processing unit (CPU) time in standard MD simulations. Each processor element in the machine has a pipeline architecture to calculate the total nonbonded force exerted on a particle by all of the other particles using information regarding the coordinates, the electric charge, and the species of each particle broadcast by the host computer. The processor then sends the calculated force back to the host computer. This article addresses the precision of the calculated nonbonded forces in the design of a processor LSI with minimal complexity. The precision of the arithmetic inside the processor that is required to calculate forces for MD simulations using Verlet's procedure was critically evaluated. Forward and backward error analysis, coupled with numerical MD experiments on one‐dimensional systems, was performed, and the following results were obtained: (1) Each element of the position vector which the processor receives from the host computer should have a precision of at least 25 bits; and (2) the pairwise forces should be calculated using floating point numbers with at least 29 bits of mantissa in the processor. Calculation of a pairwise force, which involves second‐order polynomial interpolation using a table‐driven algorithm, requires a key which contains a duplicate of at least 11 most significant bits of mantissa of the squared pairwise distance. The final result was that (3) the total force that acts on a particle, which is obtained by summing the forces exerted by all of the other particles, should be calculated using an accumulator that has a mantissa of at least 48 bits. © 1995 by John Wiley & Sons, Inc.

PVM‐AMBER: A parallel implementation of the AMBER molecular mechanics package for workstation clusters

AbstractA parallel version of the popular molecular mechanics package AMBER suitable for execution on workstation clusters has been developed. Computer‐intensive portions of molecular dynamics or free‐energy perturbation computations, such as nonbonded pair list generation or calculation of nonbonded energies and forces, are distributed across a collection of Unix workstations linked by Ethernet or FDDI connections. This parallel implementation utilizes the message‐passing software PVM (Parallel Virtual Machine) from Oak Ridge National Laboratory to coordinate data exchange and processor synchronization. Test simulations performed for solvated peptide, protein, and lipid bilayer systems indicate that reasonable parallel efficiency (70–90%) and computational speedup (2–5 × serial computer runtimes) can be achieved with small workstation clusters (typically six to eight machines) for typical biomolecular simulation problems. PVM‐AMBER is also easily and rapidly portable to different hardware platforms due to the availability of PVM for numerous computers. The current version of PVM‐AMBER has been tested successfully on Silicon Graphics, IBM RS6000, DEC ALPHA, and HP 735 workstation clusters and heterogeneous clusters of these machines, as well as on CRAY T3D and Kendall Square KSR2 parallel supercomputers. Thus, PVM‐AMBER provides a simple and cost‐effective mechanism for parallel molecular dynamics simulations on readily available hardware platforms. Factors that affect the efficiency of this approach are discussed. © 1995 by John Wiley & Sons, Inc.

Critical Evaluation of Benzene Analytical Nonbonded Force Fields. Reparametrization of the MM3 Potential

Critical Evaluation of Benzene Analytical Nonbonded Force Fields. Reparametrization of the MM3 Potential

An Efficient, Box Shape Independent Non-Bonded Force and Virial Algorithm for Molecular Dynamics

A notation is introduced and used to transform a conventional specification of the non-bonded force and virial algorithm in the case of periodic boundary conditions into an alternative specification. The implementation of the transformed specification is simpler and typically a factor of 1.5 faster than a conventional implementation. Moreover, it is generic with respect to the shape of the simulated system, i.e. the same routines can be used to handle triclinic boxes, truncated octahedron boxes etc. An implementation of this method is presented, and the speed achieved on various machines is given. Essence of the new method is that the number of calculations of image particle positions is strongly reduced during non-bonded force calculations.