General Force Field Research Articles

Nanoindentation Response of Self-Assembled Monolayers on Gold with Conical Carbon Indenters

Nanoindentation characteristics of 1-hexadecanethiol self-assembled monolayers (SAMs) on an Au surface are studied using molecular dynamics (MD) simulation. The effects of different conical angles of indenter and indentation depths on the mechanical properties of the SAMs are simulated. The interaction of SAM atoms is described by a general universal force field (UFF), the secondmoment approximation for tight-binding (TB-SMA) is used for Au substrate, and Lennard–Jones potential function is employed to describe interaction among the indenter, the SAMs, and the Au substrate atoms. The simulation results show that the adhesion force and relaxation force increase with increasing the indentation depth and conical angle of the indenter. In addition, the plastic energy of the SAMs is an order larger than the elastic energy.

Accurate Prediction of Methane Adsorption in a Metal–Organic Framework with Unsaturated Metal Sites by Direct Implementation of an ab Initio Derived Potential Energy Surface in GCMC Simulation

Whereas grand-canonical Monte Carlo (GCMC) simulations based on generic force fields provide good predictions of adsorption isotherms in metal–organic frameworks (MOFs), especially at higher temperature, they fail to correctly describe the adsorption mechanism in MOFs with coordinatively unsaturated sites (cus's) at low temperatures, even for nonpolar fluids such as methane. To address this problem, we directly implemented the potential energy surface calculated by a hybrid DFT/ab inito method in the GCMC simulations using the adsorption of methane on CuBTC as an example. A comparison with previously published in situ experiments shows that our approach not only quantitatively predicts adsorption isotherms for a wide range of temperatures and pressures but also provides the correct description of the adsorption mechanism, including adsorption on the cus's. We also show that care must be taken when selecting the ab initio method to be coupled with GCMC simulations to obtain accurate predictions.

MATCH: An atom‐typing toolset for molecular mechanics force fields

We introduce a toolset of program libraries collectively titled multipurpose atom-typer for CHARMM (MATCH) for the automated assignment of atom types and force field parameters for molecular mechanics simulation of organic molecules. The toolset includes utilities for the conversion of multiple chemical structure file formats into a molecular graph. A general chemical pattern-matching engine using this graph has been implemented whereby assignment of molecular mechanics atom types, charges, and force field parameters are achieved by comparison against a customizable list of chemical fragments. While initially designed to complement the CHARMM simulation package and force fields by generating the necessary input topology and atom-type data files, MATCH can be expanded to any force field and program, and has core functionality that makes it extendable to other applications such as fragment-based property prediction. In this work, we demonstrate the accurate construction of atomic parameters of molecules within each force field included in CHARMM36 through exhaustive cross validation studies illustrating that bond charge increment rules derived from one force field can be transferred to another. In addition, using leave-one-out substitution it is shown that it is also possible to substitute missing intra and intermolecular parameters with ones included in a force field to complete the parameterization of novel molecules. Finally, to demonstrate the robustness of MATCH and the coverage of chemical space offered by the recent CHARMM general force field (Vanommeslaeghe, et al., J Comput Chem 2010, 31, 671), one million molecules from the PubChem database of small molecules are typed, parameterized, and minimized.

On the brachistochrone of a variable mass particle in general force fields

The problem of the brachistochronic motion of a variable mass particle is considered. The particle moves through a resistant medium in the field of arbitrary active forces. Beginning from these general assumptions, and applying Pontryagin’s minimum principle along with singular optimal control theory, a corresponding two-point boundary value problem is obtained and solved. The solution proposed involves an appropriate numerical procedure based upon the shooting method. In this numerical procedure, the evaluation of ranges for unknown values of costate variables is avoided by the choice of a corresponding Cartesian coordinate of the particle as an independent variable. A numerical example assuming the resistance force proportional to the square of the particle speed is presented. A review of existing results for related problems is provided, and it can be shown that these problems may be regarded as special cases of the brachistochrone problem formulated and solved in this paper under very general assumptions by means of optimal control theory.

Multiscale Molecular Dynamics of Protein Aggregation

The 60's gave birth to the practical implementation of classical mechanics to unravel the dynamics and energetics of biomolecules. In the 70's the use of generalized force fields and more advanced integrative solutions to the microscopic understanding of nature (like hybrid QM/MM) were introduced. During the 80's, algorithms to obtain free energy values were further developed and in the 90's practical integration schemes of molecular mechanics force fields with other levels of detail (QM on one extreme and advances in implicit solvation on the other) were implemented in widely spread software. In the first decade of the XXIst century a considerable effort has been put in two seemingly discordant models for the simulation of biomolecules. On the one hand, extraordinary advances in computing technologies (both in terms of processor power and of new efficient parallel and distributed computing schemas) have allowed researchers to deal with bigger systems and longer simulations, reaching molecular processes including millions of particles or lying in the microsecond scale. On the other hand, the realization that the relevant answers to many biomolecular problems are not homogeneously distributed through the molecular structure, something already envisioned by the QM/MM pioneers more than three decades ago, has led researchers to find smart ways of putting different emphases on different ranges of the spatial or system time scale. In this context, e.g., molecular aggregation represents a paradigm for multiscalability, as molecular recognition can be understood with simple (semi-)macroscopic terms when the two fragments are far apart, while the atomic interactions need to be considered in full detail upon close distances. In this manuscript the current status of the techniques that use multiple scale representations of biomolecules are reviewed, and the findings are synthesized in a modular schema that can be extensively used when studying aggregation processes. It is shown that a smart alternative to brute force and massive computation of uninteresting regions in the all atom potential energy surface is the consideration of a simplified reference potential, explored thoroughly in the relevant regions, combined with a free energy perturbation approach that transforms this simple representation to a full atom representation.

In Silico Model Escherichia Coli Membranes: Simulating a Lipid with a Cyclopropane Ring

Bacterial membranes are composed mostly of lipids with phosphoethanolamine (PE) and phosphoglycerol (PG) head groups and a variety of fatty acid chains. A defining characteristic of bacterial membranes is the presence of lipids with a cyclic-containing chain. To our knowledge, a phospholipid with such a chain has never been studied in a molecular dynamics simulation. Thus, novel force field parameters to describe a cyclopropane moiety were developed using parameters from the CHARMM36 (C36) general force field as a basis and high-level quantum mechanical energies. Two simple membranes and one complex membrane for use as realistic model Escherichia coli (E. coli) cytoplasmic membranes were designed for future work with secondary active transporters expressed in E. coli. Compositions of the membranes were based on several different experimental methods using an E. coli K12 strain grown on Luria broth. One simple membrane consisted of 85 % PE 18:1/16:0 (POPE) and 15 % PG 18:1/16:0(POPG), and the other of only PE cy17:0/16:0 (PMPE). The complex membrane consisted of six different phospholipids, the most prevalent being the cyclic-containing lipid, PMPE. NPT simulations were carried out at 310 K for 50 ns using the C36 lipid force field and the developed cyclopropane moiety force field for each membrane. NMR deuterium order parameters (SCD), density profiles, and diffusion constants are compared with experiment. The model membranes will provide a basis to study membrane bound proteins or other bound molecules expressed in E. coli with realistic molecular dynamic simulations.

Automation of the Charmm General Force Field for Drug-Like Molecules

The CHARMM Force Field[1] is widely used for the study of biological macromolecules. Recently, we released the CHARMM General Force Field (CGenFF) which allows simulating drug-like molecules in an environment represented by the CHARMM biomolecular force field.[2] To facilitate the use of CGenFF, a computer program was developed that automatically assigns atom types, bonded parameters and charges to any organic molecule. The atom typing is rule-based and programmable, making it easy to update the atom typing scheme as the force field grows. Assignment of bonded parameters is based on substituting atom types in the definition of the parameter, and returns a “penalty score” as a measure for the accuracy of the approximation. Charges are assigned using an extended bond-charge increment scheme that should be able to capture short- and medium-range inductive and mesomeric effects.[3] This functionality is accessible on the web, which makes it possible to apply CGenFF on arbitrary molecules in a matter of seconds.[4] Although the program should ideally accept any molecule, the quality of the resulting parameters (as indicated by their “penalty scores”) varies depending on whether a similar compound has been parametrized before. Therefore, an effort has been started to build a computational engine with the capability of automatically validating and/or optimizing parameters for drug-like model compounds of the user's interest in the framework of CGenFF, thus providing high-quality force field parameters for computer-aided drug design.[4]1. A. D. MacKerell et al., J. Phys. Chem. B 1998, 102, 3586-36162. K. Vanommeslaeghe et al., J. Comput. Chem. 2010, 31, 671-6903. K. Vanommeslaeghe et al., in preparation4. NSF Award #0823198; https://www.paramchem.org/

The Harmonic Vibration Problem of the Cubic XY8 Molecule

The secular equations for the various irreducible representations of the cubic XY8 molecule (point group Oh) are deduced by means of Wilson's F-G matrix method. The potential function used corresponds to the most general harmonic force field. The equations for the case of the SVFF are also given. The forms of the normal vibrations are shown in schematic diagrams.

Structural and energetic aspects of the differences between real and predicted polymorphs

AbstractIn crystal structure prediction simulations based on lattice energy minimization, usually hundreds of structures within a reasonable range of lattice energy and density are found, whereas in practice, it is very rare to find more than a few polymorphs of the same compound. In the work presented here, this discrepancy is investigated from a structural and energetic point of view. 56 crystal structures of 26 polymorphic mono‐ and disubstituted aromatic compounds, extracted from the Cambridge Structural Database, have been analysed with respect to inter‐polymorphic structural similarity. For comparison, potential crystal packing arrangements of the substances have been predicted with molecular mechanics simulations using a generic force field. The predicted structures are analysed with respect to structural features and similarity, and with respect to the number of structures and their lattice energy. It is found that the real polymorphs studied in this work tend to be structurally quite dissimilar with regard to hydrogen bonding and spatial packing of structural motifs, while many of the predicted structures of a given compound are very similar to each other. The results suggest that structure and lattice energy alone cannot explain why so few polymorphs are found in practice compared to the very large numbers predicted in simulations. (© 2010 WILEY‐VCH Verlag GmbH & Co. KGaA, Weinheim)

Diffusion of tagged particle in an exclusion process

We study the diffusion of tagged hard-core interacting Brownian point particles under the influence of an external force field in one dimension. Using the Jepsen line we map this many-particle problem onto a single particle one. We obtain general equations for the distribution and the mean-square displacement <(xT)2> of the tagged center particle valid for rather general external force fields and initial conditions. The case of symmetric distribution of initial conditions around the initial position of the tagged particle on x=0 and symmetric potential fields V(x)=V(-x) yields zero drift <xT>=0 and is investigated in detail. We find <(xT)2>=R(1-R)/2Nr2 where 2N is the (large) number of particles in the system. R is a single particle reflection coefficient, i.e., the probability that a particle free of collisions starts on x0>0 and remains in x>0 while r is the probability density of noninteracting particles on the origin. We show that this equation is related to the mathematical theory of order statistics and it can be used to find <(xT)2> even when the motion between collision events is not Brownian (e.g., it might be ballistic or anomalous diffusion). As an example we derive the Percus relation for non-Gaussian diffusion. A wide range of physical behaviors emerge which are very different than the classical single file subdiffusion <(xT)2> approximately t1/2 found for uniformly distributed particles in an infinite space and in the absence of force fields.

A multivariable hierarchical finite element for static and vibration analysis of beams

Formulations of a multivariable hierarchical beam element for static and vibration analysis are presented based on the generalized variational principle with two kinds of variables. Two forms of shifted Legendre hierarchical polynomials are used as interpolating basis functions of displacement and generalized force field functions for the beam element respectively, which will simplify the computations of the relevant matrices. The multivariable hierarchical beam element formulations, in which the displacement and generalized force field functions are independently constructed, are derived by applying the generalized variational principle with two kinds of variables. Since differential operations to obtain stress fields in conventional displacement based finite element methods are not required, the present method has very high accuracy for the two kinds of independent variables simultaneously, especially for the generalized forces. Static and vibration numerical examples demonstrate the applicability of the proposed method. The proposed method can be easily extended to deal with structural analysis of shells or plates.

Molecular Modeling of Crystalline Alkylthiophene Oligomers and Polymers

We present the results of a thorough molecular modeling study of several alkylthiophene-based oligomers and polymers. In particular, we consider two polymers whose limit-ordered crystal structures have been recently reported by our group, on the basis of powder X-ray data analysis: poly(3-(S)-2-methylbutylthiophene) (P3MBT) and form I' of poly(3-butylthiophene) (P3BT). We first describe the development of a series general purpose force fields for the simulation of these and related systems. The force fields incorporate the results of ab initio calculations of the bond torsion energies of selected oligomers and differ in the set of atomic charges used to represent the electrostatic interactions. We then present the results of an extensive validation of these force fields, by means of molecular mechanics (MM) energy minimizations and molecular dynamics (MD) simulations of the crystal structures of these oligomers and polymers. While our "best" force field does not outperform the others on each of the investigated systems, it provides a balanced description of their overall structure and energetics. Finally, our MM minimizations and MD simulations confirm that the reported crystal structures of P3MBT and P3BT are stable and correspond to well-defined energetic minima. The room-temperature MD simulations reveal a certain degree of side-chain disorder, even in our virtually defect-free polymer crystal models.

Validity of current force fields for simulations on boron nitride nanotubes

Past molecular dynamics (MD) studies of boron nitride nanotubes (BNNTs) have used van der Waals parameters from generic force fields, combined with various values for the partial charges on the boron and nitrogen atoms. This study explores the validity of these parameters by first using quantum chemical packages Car-Parrinello molecular dynamics (CPMD) and Gaussian to compute partial charges for isolated and periodic BNNTs, both with and without water. Then in order to test the accuracy of the molecular mechanics force field using our computed charges, the authors calculate the interaction energy between each water molecule in a hydrated nanotube with the nanotube itself using two methods: first using a quantum chemical calculation, and secondly using the molecular mechanics force field. The authors show that in order to obtain satisfactory agreement in the interaction energies the boron and nitrogen Lennard-Jones parameters must be adjusted from their usual values. Modified Lennard-Jones parameters and partial charges, obtained by fitting, are presented as candidates for future MD simulations of hydrated BNNTs.

CHARMM general force field: A force field for drug-like molecules compatible with the CHARMM all-atom additive biological force fields.

The widely used CHARMM additive all-atom force field includes parameters for proteins, nucleic acids, lipids, and carbohydrates. In the present article, an extension of the CHARMM force field to drug-like molecules is presented. The resulting CHARMM General Force Field (CGenFF) covers a wide range of chemical groups present in biomolecules and drug-like molecules, including a large number of heterocyclic scaffolds. The parametrization philosophy behind the force field focuses on quality at the expense of transferability, with the implementation concentrating on an extensible force field. Statistics related to the quality of the parametrization with a focus on experimental validation are presented. Additionally, the parametrization procedure, described fully in the present article in the context of the model systems, pyrrolidine, and 3-phenoxymethylpyrrolidine will allow users to readily extend the force field to chemical groups that are not explicitly covered in the force field as well as add functional groups to and link together molecules already available in the force field. CGenFF thus makes it possible to perform "all-CHARMM" simulations on drug-target interactions thereby extending the utility of CHARMM force fields to medicinally relevant systems.

Molecular dynamics of contact behavior of self-assembled monolayers on gold using nanoindentation

Molecular dynamics simulation is used to study nanoindentation of the self-assembled monolayers (SAMs) on an Au surface. The interaction of SAM atoms is described by a general universal force field (UFF), the tight-binding second-moment approximation (TB-SMA) is used for Au substrate, and the Lennard-Jones potential function is employed to describe interaction among the indenter, the SAMs, and the Au substrate atoms. The model consists of a planar Au substrate with n-hexadecanethiol SAM chemisorbed to the substrate. The simulation results show that the contact pressure increases as the SAMs temperature increases. In addition, the contact pressure also increases as the depth and velocity of indentation increase.

Simulation of Water Adsorption on Kaolinite under Atmospheric Conditions

Grand canonical Monte Carlo calculations are employed to investigate water adsorption on kaolinite at 298 and 235 K. Both basal planes (the Al and Si surfaces) as well as two edge-like surfaces are considered. The general force field CLAYFF is used together with the SPC/E and TIP5P-E models for water. Problems that occur in single slab simulations due to arbitrary truncation of the point charge lattice are identified, and a working remedy is discussed. The edges and the Al surface adsorb water at subsaturation in the atmospherically relevant pressure range. The Si surface remains dry up to saturation. Both edges have a very strong affinity for water and adsorb continuously up to monolayer coverage. The Al surface has a weaker affinity for water but forms a subsaturation monolayer. On the Al surface, the monolayer is formed in an essentially sharp transition, and strong hysteresis is observed upon desorption. This indicates collective behavior among the water molecules which is not present for the edges. Binding energies of singly adsorbed water molecules at 10 K were determined to understand the differences in water uptake by the four kaolinite surfaces. Binding energies (SPC/E) of -21.6, -46.4, -73.5, and -94.1 kJ/mol, were determined for the Si surface, Al surface, unprotonated edge, and protonated edge, respectively. The water monolayer on the Al surface, particularly at 235 K, exhibits hexagonal patterns. However, the associated lattice parameters are not compatible with ice Ih. Water density and hydrogen bonding in the monolayers at both 298 and 235 K were also determined to better understand the structure of the adsorbed water.

Adsorption and Separation of Noble Gases by IRMOF-1: Grand Canonical Monte Carlo Simulations

The gas storage capacity of metal−organic frameworks (MOFs) is well-known and has been investigated using both experimental and computational methods. Previous Monte Carlo computer simulations of gas adsorption by MOFs have made several questionable approximations regarding framework−framework and framework−adsorbate interactions: potential parameters from general force fields have been used, and framework atoms were fixed at their crystallographic coordinates (rigid framework). We assess the validity of these approximations with grand canonical Monte Carlo simulations for a well-known Zn-based MOF (IRMOF-1), using potential parameters specifically derived for IRMOF-1. Our approach is validated by comparison with experimental results for hydrogen and xenon adsorption at room temperature. The effects of framework flexibility on the adsorption of noble gases and hydrogen are described, as well as the selectivity of IRMOF-1 for xenon versus other noble gases. At both low temperature (78 K) and room temperatu...

Electrodynamics as a particular case of the most general relativistic force field

Electrodynamics as a particular case of the most general relativistic force field

On the accuracy of predicting shear viscosity of molecular liquids using the periodic perturbation method

Non-equilibrium molecular dynamics simulations are performed to calculate shear viscosities of 16 representative molecular liquids using the periodic perturbation method (PPM). A perturbation index is defined to measure the strength of the perturbation. It is identified that the predictions are systematically underestimated using PPM. The origin of the underestimate is the acoustic wave in the liquid density, which is persistent in the simulation box unless the perturbation is completely removed. However, there is a linear correlation between the perturbation indexes and the apparent viscosities, which can be utilized to accurately predict the shear viscosities. Finally, it is demonstrated that general force fields derived based on equilibrium properties can be used to predict the shear viscosities of small molecular liquids with relative errors less than 10%.

REACH Coarse-Grained Biomolecular Simulation: Transferability between Different Protein Structural Classes

Coarse graining of protein interactions provides a means of simulating large biological systems. The REACH (Realistic Extension Algorithm via Covariance Hessian) coarse-graining method, in which the force constants of a residue-scale elastic network model are calculated from the variance-covariance matrix obtained from atomistic molecular dynamics (MD) simulation, involves direct mapping between scales without the need for iterative optimization. Here, the transferability of the REACH force field is examined between protein molecules of different structural classes. As test cases, myoglobin (all α), plastocyanin (all β), and dihydrofolate reductase ( α/ β) are taken. The force constants derived are found to be closely similar in all three proteins. An MD version of REACH is presented, and low-temperature coarse-grained (CG) REACH MD simulations of the three proteins are compared with atomistic MD results. The mean-square fluctuations of the atomistic MD are well reproduced by the CGMD. Model functions for the CG interactions, derived by averaging over the three proteins, are also shown to produce fluctuations in good agreement with the atomistic MD. The results indicate that, similarly to the use of atomistic force fields, it is now possible to use a single, generic REACH force field for all protein studies, without having first to derive parameters from atomistic MD simulation for each individual system studied. The REACH method is thus likely to be a reliable way of determining spatiotemporal motion of a variety of proteins without the need for expensive computation of long atomistic MD simulations.