The QM‐MM interface for CHARMM‐deMon

We present a new QM/MM interface for fast and efficient simulations of organic and biological molecules. The CHARMM/deMon interface has been developed and tested to perform minimization and atomistic simulations for multi‐particle systems. The current features of this QM/MM interface include readability for molecular dynamics, tested compatibility with Free Energy Perturbation simulations (FEP) using the dual topology/single coordinate method. The current coupling scheme uses link atoms, but further extensions of the code to incorporate other available schemes are planned. We report the performance of different levels of theory for the treatment of the QM region, while the MM region was represented by a classical force‐field (CHARMM27) or a polarizable force‐field based on a simple Drude model. The current QM/MM implementation can be coupled to the dual‐thermostat method and the VV2 integrator to run molecular dynamics simulations. © 2009 Wiley Periodicals, Inc. J Comput Chem, 2010     

Interaction of organic molecules with the TiO2 (110) surface: ab inito calculations and classical force fields

We have studied the adsorption of a number of organic molecules consisting of methyl, benzyl, and carboxylic groups on the rutile TiO2 (110) surface using both ab initio and atomistic simulation techniques. We have tested the applicability of a simple embedded cluster model to studying the adsorption of small organic molecules on the perfect rutile TiO2 (110) surface, and used this model to develop a classical force field for the interactions of a wide class of organic molecules consisting of these groups with the rutile TiO2 (110) surface. The force field accounts for physisorption and ionic bonding of organic molecules at the surface. It allows the reproduction of adsorption energies and of geometries of organic molecules on the rutile surface. It should be useful for studying diffusion of these molecules and their manipulation with use of AFM and STM tips.     

Strategies for multiscale modeling and simulation of organic materials: polymers and biopolymers

Advances in theory and methods are making it practical to consider fully first principles (de novo) predictions of structures, properties and processes for organic materials. However, despite the progress there remains an enormous challenge in bridging the vast range of distances and time scales between de novo atomistic simulations and the quantitative continuum models for the macroscopic systems essential in industrial design and operations. Recent advances relevant to such developments include: quantum chemistry including continuum solvation and force field embedding, de novo force fields to describe phase transitions, molecular dynamics (MD) including continuum solvent, non equilibrium MD for rheology and thermal conductivity and mesoscale simulations. To provide some flavor for the opportunities we will illustrate some of the progress and challenges by summarizing some recent developments in methods and their applications to polymers and biopolymers. Four different topics will be covered: (1) hierarchical modeling approach applied to modeling olfactory receptors, (2) stabilization of leucine zipper coils by introduction of trifluoroleucine, (3) modeling response of polymers sensors for electronic nose, and (4) diffusion of gases in amorphous polymers.     

ForceFit: A code to fit classical force fields to quantum mechanical potential energy surfaces

The ForceFit program package has been developed for fitting classical force field parameters based upon a force matching algorithm to quantum mechanical gradients of configurations that span the potential energy surface of the system. The program, which runs under UNIX and is written in C++, is an easy‐to‐use, nonproprietary platform that enables gradient fitting of a wide variety of functional force field forms to quantum mechanical information obtained from an array of common electronic structure codes. All aspects of the fitting process are run from a graphical user interface, from the parsing of quantum mechanical data, assembling of a potential energy surface database, setting the force field, and variables to be optimized, choosing a molecular mechanics code for comparison to the reference data, and finally, the initiation of a least squares minimization algorithm. Furthermore, the code is based on a modular templated code design that enables the facile addition of new functionality to the program. © 2010 Wiley Periodicals, Inc. J Comput Chem, 2010     

Computer simulation of polypeptide adsorption on model biomaterials

When biomaterials are inserted in a biological environment, for instance in a body implant, proteins do quickly adsorb on the exposed surface. Such process is of fundamental importance, since it directs the subsequent cell adhesion. Here we review recent advances in this field obtained with molecular simulations. While coarse-grained models can provide important general results, as it has long been recognized in polymer science, the hierarchical structure of a very complex copolymer such as a protein, together with the nature of the biomaterial surface suggest that atomistic models are better suited to investigate these phenomena. Thus, after briefly mentioning some common features of coarse-grained and atomistic force fields, we first discuss early theoretical and coarse-grained simulation results about protein adsorption, and then we highlight the main results recently obtained by us with atomistic models. In particular, we discuss some conformational and energetic aspects of the adsorption of protein fragments with different secondary structure on surfaces of different wettability, including hydrophobic graphite and hydrophilic poly(vinylalcohol). We also consider other features, such as the simulation of the materials wettability, the hydration of the adsorbed fragments, their kinetics of spreading, and the sequential adsorption of two protein fragments on top of each other, highlighting the results of general interest.     

Molecular dynamics simulations of aqueous ions at the liquid–vapor interface accelerated using graphics processors

Molecular dynamics (MD) simulations are a vital tool in chemical research, as they are able to provide an atomistic view of chemical systems and processes that is not obtainable through experiment. However, large‐scale MD simulations require access to multicore clusters or supercomputers that are not always available to all researchers. Recently, scientists have returned to exploring the power of graphics processing units (GPUs) for various applications, such as MD, enabled by the recent advances in hardware and integrated programming interfaces such as NVIDIA's CUDA platform. One area of particular interest within the context of chemical applications is that of aqueous interfaces, the salt solutions of which have found application as model systems for studying atmospheric process as well as physical behaviors such as the Hoffmeister effect. Here, we present results of GPU‐accelerated simulations of the liquid–vapor interface of aqueous sodium iodide solutions. Analysis of various properties, such as density and surface tension, demonstrates that our model is consistent with previous studies of similar systems. In particular, we find that the current combination of water and ion force fields coupled with the ability to simulate surfaces of differing area enabled by GPU hardware is able to reproduce the experimental trend of increasing salt solution surface tension relative to pure water. In terms of performance, our GPU implementation performs equivalent to CHARMM running on 21 CPUs. Finally, we address possible issues with the accuracy of MD simulaions caused by nonstandard single‐precision arithmetic implemented on current GPUs. © 2010 Wiley Periodicals, Inc. J Comput Chem, 2011     

A discrete interaction model/quantum mechanical method for simulating surface-enhanced Raman spectroscopy

We have derived and implemented analytical gradients for the discrete interaction model/quantum mechanics (DIM/QM) method. DIM/QM combines an atomistic electrodynamics model with time-dependent density functional theory and thus enables modeling of the optical properties for a molecule while taking into account the local environment of a nanoparticle's surface. The DIM/QM analytical gradients allow for geometry optimizations, vibrational frequencies, and Raman spectra to be simulated for molecules interacting with metal nanoparticles. We have simulated the surface-enhanced Raman scattering (SERS) spectra for pyridine adsorbed on different sites of icosahedral nanoparticles with diameters between 1 and 8 nm. To describe the adsorption of the pyridine molecule onto the metal surface, we have implemented a coordination-dependent force field to differentiate the various local surface environments. We find that the DIM/QM method predicts geometries and frequencies that are in good agreement with full QM simulations and experiments. For the simulated SERS spectra of pyridine, we find a significant dependence on the adsorption site and the size of the metal nanoparticle. This illustrates the importance of accounting for the local environment around the molecule. The Raman enhancement factors are shown to roughly mirror the magnitude of the nanoparticle's local field about the molecule. Because the simulated nanoparticles are small, the plasmon peaks are quite broad which results in weak local electric fields and thus modest Raman enhancement factors.     

Revealing the Transient Concentration of CO2 in a Mixed‐Matrix Membrane by IR Microimaging and Molecular Modeling

Through IR microimaging the spatially and temporally resolved development of the CO₂ concentration in a ZIF‐8@6FDA‐DAM mixed matrix membrane (MMM) was visualized during transient adsorption. By recording the evolution of the CO₂ concentration, it is observed that the CO₂ molecules propagate from the ZIF‐8 filler, which acts as a transport “highway”, towards the surrounding polymer. A high‐CO₂‐concentration layer is formed at the MOF/polymer interface, which becomes more pronounced at higher CO₂ gas pressures. A microscopic explanation of the origins of this phenomenon is suggested by means of molecular modeling. By applying a computational methodology combining quantum and force‐field based calculations, the formation of microvoids at the MOF/polymer interface is predicted. Grand canonical Monte Carlo simulations further demonstrate that CO₂ tends to preferentially reside in these microvoids, which is expected to facilitate CO₂ accumulation at the interface.     

Atomistic insight into chondroitin‐6‐sulfate glycosaminoglycan chain through quantum mechanics calculations and molecular dynamics simulation

Chondroitin‐6‐sulfate (C6S) is a glycosaminoglycan (GAG) constituent in the extracellular matrix, which participates actively in crucial biological processes, as well as in various pathological conditions, such as atherosclerosis and cancer. Molecular interactions involving the C6S chain are therefore of considerable interest. A computational model for atomistic simulation was built. This work describes the design and validation of a force field for a C6S dodecasaccharide chain. The results of an extensive molecular dynamics simulation performed with the new force field provide a novel insight into the structure and dynamics of the C6S chain. The intramolecular H‐bonds in the disaccharide linkage region are suggested to play a major role in determining the chain structural dynamics. Moreover, the unravelling of an additional H‐bond involving the sulfate groups in C6S is interesting as changes in sulfation have been claimed to be an important factor in several diseases. The force field will prove useful for future studies of crucial interactions between C6S and various nanoassemblies. It can also be used as a basis for modeling of other GAGs. © 2010 Wiley Periodicals, Inc. J Comput Chem, 2010     

Classical Reactive Molecular Dynamics Implementations: State of the Art

Reactive molecular dynamics (RMD) implementations equipped with force field approaches to simulate both the time evolution as well as chemical reactions of a broad class of materials are reviewed herein. We subdivide the RMD approaches developed during the last decade as well as older ones already reviewed in 1995 by Srivastava and Garrison and in 2000 by Brenner into two classes. The methods in the first RMD class rely on the use of a reaction cutoff distance and employ a sudden transition from the educts to the products. Due to their simplicity these methods are well suited to generate equilibrated atomistic or material‐specific coarse‐grained polymer structures. In connection with generic models they offer useful qualitative insight into polymerization reactions. The methods in the second RMD class are based on empirical reactive force fields and implement a smooth and continuous transition from the educts to the products. In this RMD class, the reactive potentials are based on many‐body or bond‐order force fields as well as on empirical standard force fields, such as CHARMM, AMBER or MM3 that are modified to become reactive. The aim with the more sophisticated implementations of the second RMD class is the investigation of the reaction kinetics and mechanisms as well as the evaluation of transition state geometries. Pure or hybrid ab initio, density functional, semi‐empirical, molecular mechanics, and Monte Carlo methods for which no time evolution of the chemical systems is achieved are excluded from the present review. So are molecular dynamics techniques coupled with quantum chemical methods for the treatment of the reactive regions, such as Car–Parinello molecular dynamics.     

Paramfit: Automated optimization of force field parameters for molecular dynamics simulations

The generation of bond, angle, and torsion parameters for classical molecular dynamics force fields typically requires fitting parameters such that classical properties such as energies and gradients match precalculated quantum data for structures that scan the value of interest. We present a program, Paramfit, distributed as part of the AmberTools software package that automates and extends this fitting process, allowing for simplified parameter generation for applications ranging from single molecules to entire force fields. Paramfit implements a novel combination of a genetic and simplex algorithm to find the optimal set of parameters that replicate either quantum energy or force data. The program allows for the derivation of multiple parameters simultaneously using significantly fewer quantum calculations than previous methods, and can also fit parameters across multiple molecules with applications to force field development. Paramfit has been applied successfully to systems with a sparse number of structures, and has already proven crucial in the development of the Assisted Model Building with Energy Refinement Lipid14 force field. © 2014 Wiley Periodicals, Inc.     

Efficient parametrization of complex molecule–surface force fields

We present an efficient scheme for parametrizing complex molecule–surface force fields from ab initio data. The cost of producing a sufficient fitting library is mitigated using a 2D periodic embedded slab model made possible by the quantum mechanics/molecular mechanics scheme in CP2K. These results were then used in conjunction with genetic algorithm (GA) methods to optimize the large parameter sets needed to describe such systems. The derived potentials are able to well reproduce adsorption geometries and adsorption energies calculated using density functional theory. Finally, we discuss the challenges in creating a sufficient fitting library, determining whether or not the GA optimization has completed, and the transferability of such force fields to similar molecules. © 2015 Wiley Periodicals, Inc.     

A Complete Multiscale Modelling Approach for Polymer–Clay Nanocomposites

We present an innovative, multiscale computational approach to probe the behaviour of polymer–clay nanocomposites (PCNs). Our modeling recipe is based on 1) quantum/force‐field‐based atomistic simulation to derive interaction energies among all system components; 2) mapping of these values onto mesoscopic bead–field (MBF) hybrid‐method parameters; 3) mesoscopic simulations to determine system density distributions and morphologies (i.e., intercalated versus exfoliated); and 4) simulations at finite‐element levels to calculate the relative macroscopic properties. The entire computational procedure has been applied to two well‐known PCN systems, namely Nylon 6/Cloisite 20A and Nylon 6/Cloisite 30B, as test materials, and their mechanical properties were predicted in excellent agreement with the available experimental data. Importantly, our methodology is a truly bottom‐up approach, and no “learning from experiment” was needed in any step of the entire procedure.     

Structure and energetics of diphenylalanine self-assembling on Cu(110)

We investigate the dynamical features of the adsorption of diphenylalanine molecules on the Cu(110) surface and of their assembling into supramolecular structures by a combination of quantum and classical atomistic modeling with dynamic scanning tunneling microscopy and spectroscopic experiments. Our results reveal a self-assembling mechanism in which isolated adsorbed molecules change their conformation and adsorption mode as a consequence of their mutual interactions. In particular, the formation of zwitterions after proton transfer between initially neutral molecules is found to be the key event of the assembling process, which stabilizes the supramolecular structures. Because of the constraints on the intermolecular bonds exerted by the surface-molecule interactions, the assembly process is strictly stereoselective, and may suggest a general model for patterning and functionalization of bare metal surfaces with short chiral peptides.     

A polarizable ellipsoidal force field for halogen bonds

The anisotropic effects and short‐range quantum effects are essential characters in the formation of halogen bonds. Since there are an array of applications of halogen bonds and much difficulty in modeling them in classical force fields, the current research reports solely the polarizable ellipsoidal force field (PEff) for halogen bonds. The anisotropic charge distribution was represented with the combination of a negative charged sphere and a positively charged ellipsoid. The polarization energy was incorporated by the induced dipole model. The resulting force field is “physically motivated,” which includes separate, explicit terms to account for the electrostatic, repulsion/dispersion, and polarization interaction. Furthermore, it is largely compatible with existing, standard simulation packages. The fitted parameters are transferable and compatible with the general AMBER force field. This PEff model could correctly reproduces the potential energy surface of halogen bonds at MP2 level. Finally, the prediction of the halogen bond properties of human Cathepsin L (hcatL) has been found to be in excellent qualitative agreement with the cocrystal structures. © 2013 Wiley Periodicals, Inc.     

QuanPol: A full spectrum and seamless QM/MM program

The quantum chemistry polarizable force field program (QuanPol) is implemented to perform combined quantum mechanical and molecular mechanical (QM/MM) calculations with induced dipole polarizable force fields and induced surface charge continuum solvation models. The QM methods include Hartree–Fock method, density functional theory method (DFT), generalized valence bond theory method, multiconfiguration self‐consistent field method, Møller–Plesset perturbation theory method, and time‐dependent DFT method. The induced dipoles of the MM atoms and the induced surface charges of the continuum solvation model are self‐consistently and variationally determined together with the QM wavefunction. The MM force field methods can be user specified, or a standard force field such as MMFF94, Chemistry at Harvard Molecular Mechanics (CHARMM), Assisted Model Building with Energy Refinement (AMBER), and Optimized Potentials for Liquid Simulations‐All Atom (OPLS‐AA). Analytic gradients for all of these methods are implemented so geometry optimization and molecular dynamics (MD) simulation can be performed. MD free energy perturbation and umbrella sampling methods are also implemented. © 2013 Wiley Periodicals, Inc.     

A hybrid all-atom/coarse grain model for multiscale simulations of DNA

Hybrid simulations of molecular systems, which combine all-atom (AA) with simplified (or coarse grain, CG) representations, propose an advantageous alternative to gain atomistic details on relevant regions while getting profit from the speedup of treating a bigger part of the system at the CG level. Here we present a reduced set of parameters derived to treat a hybrid interface in DNA simulations. Our method allows us to forthrightly link a state-of-the-art force field for AA simulations of DNA with a CG representation developed by our group. We show that no modification is needed for any of the existing residues (neither AA nor CG). Only the bonding parameters at the hybrid interface are enough to produce a smooth transition of electrostatic, mechanic and dynamic features in different AA/CG systems, which are studied by molecular dynamics simulations using an implicit solvent. The simplicity of the approach potentially permits us to study the effect of mutations/modifications as well as DNA binding molecules at the atomistic level within a significantly larger DNA scaffold considered at the CG level. Since all the interactions are computed within the same classical Hamiltonian, the extension to a quantum/classical/coarse-grain multilayer approach using QM/MM modules implemented in widely used simulation packages is straightforward.     

The IR spectrum of supercritical water: Combined molecular dynamics/quantum mechanics strategy and force field for cluster sampling

Supercritical water was analyzed recently as a gas of small clusters of waters linked to each other by intermolecular hydrogen‐bonds, but unexpected “linear” conformations of clusters are required to reproduce the infra‐red (IR) spectra of the supercritical state. Aiming at a better understanding of clusters in supercritical water, this work presents a strategy combining classical molecular dynamics to explore the potential energy landscape of water clusters with quantum mechanical calculation of their IR spectra. For this purpose, we have developed an accurate and flexible force field of water based on the TIP5P 5‐site model. Water dimers and trimers obtained with this improved force field compare well with the quantum mechanically optimized clusters. Exploration by simulated annealing of the potential energy surface of the classical force field reveals a new trimer conformation whose IR response determined from quantum calculations could play a role in the IR spectra of supercritical water. © 2011 Wiley Periodicals, Inc. Int J Quantum Chem, 2011