A polarizable empirical force field for molecular dynamics simulation of liquid hydrocarbons

Electronic polarizability is usually treated implicitly in molecular simulations, which may lead to imprecise or even erroneous molecular behavior in spatially electronically inhomogeneous regions of systems such as proteins, membranes, interfaces between compounds, or mixtures of solvents. The majority of available molecular force fields and molecular dynamics simulation software packages does not account explicitly for electronic polarization. Even the simplest charge‐on‐spring (COS) models have only been developed for few types of molecules. In this work, we report a polarizable COS model for cyclohexane, as this molecule is a widely used solvent, and for linear alkanes, which are also used as solvents, and are the precursors of lipids, amino acid side chains, carbohydrates, or nucleic acid backbones. The model is an extension of a nonpolarizable united‐atom model for alkanes that had been calibrated against experimental values of the density, the heat of vaporization and the Gibbs free energy of hydration for each alkane. The latter quantity was used to calibrate the parameters governing the interaction of the polarizable alkanes with water. Subsequently, the model was tested for other structural, thermodynamic, dielectric, and dynamic properties such as trans/gauche ratios, excess free energy, static dielectric permittivity, and self‐diffusion. A good agreement with the experimental data for a large set of properties for each considered system was obtained, resulting in a transferable set of polarizable force‐field parameters for CH₂, CH₃, and CH₄ moieties. © 2014 Wiley Periodicals, Inc.     

Molecular modeling and dynamics studies with explicit inclusion of electronic polarizability: theory and applications

A current emphasis in empirical force fields is on the development of potential functions that explicitly treat electronic polarizability. In the present article, the commonly used methodologies for modeling electronic polarization are presented along with an overview of selected application studies. Models presented include induced point-dipoles, classical Drude oscillators, and fluctuating charge methods. The theoretical background of each method is followed by an introduction to extended Lagrangian integrators required for computationally tractable molecular dynamics simulations using polarizable force fields. The remainder of the review focuses on application studies using these methods. Emphasis is placed on water models, for which numerous examples exist, with a more thorough discussion presented on the recently published models associated with the Drude-based CHARMM and the AMOEBA force fields. The utility of polarizable models for the study of ion solvation is then presented followed by an overview of studies of small molecules (e.g., CCl₄, alkanes, etc.) and macromolecule (proteins, nucleic acids and lipid bilayers) application studies. The review is written with the goal of providing a general overview of the current status of the field and to facilitate future application and developments.     

Properties of alkali-halide salt solutions about polarizable nanoparticle solutes for different ion models

We investigate the distributions of various salts about large hydrophobic polarizable solutes in aqueous electrolyte solutions. The solutes are modeled as nanometer-sized cylindrical objects, a scale relevant to biomolecules and nanomaterials, and particularly high aspect ratio nanoparticles. Interactions, including image charge forces arising from the finite polarizability of the solute, between explicit solvent/ions and the solute are computed explicitly using a molecular dynamics simulation methodology we have recently introduced. Comparisons are made between several salt species and different models of the force fields for each ionic component of the salt. We find evidence that both small cations, Li(+), and large anions, I(-), adsorb at hydrophobic interfaces. Our results indicate that the ion structure about the solute is strongly dependent on the force field investigated, suggesting that ion selectivity is quite sensitive to the respective parameters defining the ion's size and binding energy as well as to the polarizability of the solute.     

Polarizable empirical force field for sulfur‐containing compounds based on the classical Drude oscillator model

Condensed‐phase computational studies of molecules using molecular mechanics approaches require the use of force fields to describe the energetics of the systems as a function of structure. The advantage of polarizable force fields over nonpolarizable (or additive) models lies in their ability to vary their electronic distribution as a function of the environment. Toward development of a polarizable force field for biological molecules, parameters for a series of sulfur‐containing molecules are presented. Parameter optimization was performed to reproduce quantum mechanical and experimental data for gas phase properties including geometries, conformational energies, vibrational spectra, and dipole moments as well as for condensed phase properties such as heats of vaporization, molecular volumes, and free energies of hydration. Compounds in the training set include methanethiol, ethanethiol, propanethiol, ethyl methyl sulfide, and dimethyl disulfide. The molecular volumes and heats of vaporization are in good accordance with experimental values, with the polarizable model performing better than the CHARMM22 nonpolarizable force field. Improvements with the polarizable model were also obtained for molecular dipole moments and in the treatment of intermolecular interactions as a function of orientation, in part due to the presence of lone pairs and anisotropic atomic polarizability on the sulfur atoms. Significant advantage of the polarizable model was reflected in calculation of the dielectric constants, a property that CHARMM22 systematically underestimates. The ability of this polarizable model to accurately describe a range of gas and condensed phase properties paves the way for more accurate simulation studies of sulfur‐containing molecules including cysteine and methionine residues in proteins. © 2010 Wiley Periodicals, Inc. J Comput Chem, 2010     

Implementation of extended Lagrangian dynamics in GROMACS for polarizable simulations using the classical Drude oscillator model

Explicit treatment of electronic polarization in empirical force fields used for molecular dynamics simulations represents an important advancement in simulation methodology. A straightforward means of treating electronic polarization in these simulations is the inclusion of Drude oscillators, which are auxiliary, charge‐carrying particles bonded to the cores of atoms in the system. The additional degrees of freedom make these simulations more computationally expensive relative to simulations using traditional fixed‐charge (additive) force fields. Thus, efficient tools are needed for conducting these simulations. Here, we present the implementation of highly scalable algorithms in the GROMACS simulation package that allow for the simulation of polarizable systems using extended Lagrangian dynamics with a dual Nosé–Hoover thermostat as well as simulations using a full self‐consistent field treatment of polarization. The performance of systems of varying size is evaluated, showing that the present code parallelizes efficiently and is the fastest implementation of the extended Lagrangian methods currently available for simulations using the Drude polarizable force field.     

Revisiting the hexane-water interface via molecular dynamics simulations using nonadditive alkane-water potentials

We present results addressing properties of a polarizable force field for hexane based on the fluctuating charge (FQ) formalism and developed in conjunction with the Chemistry at Harvard Molecular Mechanics (CHARMM) potential function. Properties of bulk neat hexane, its liquid-vapor interface, and its interface with a polarizable water model (TIP4P-FQ) are discussed. The FQ model is compared to a recently modified alkane model, C27r, also based on the CHARMM potential energy function. With respect to bulk properties, both models predict bulk density within 1%; the FQ model predicts the liquid vaporization enthalpy within 2%, while the C27r force field underestimates the property by roughly 20% (and in this sense reflects the quality of the C27r force field across the spectrum of linear and branched alkanes). The FQ hexane model realistically captures the dielectric properties of the bulk in terms of a dielectric constant of 1.94, in excellent agreement with experimental values in the range of 1.9-2.02. This behavior is also in conformity with a recent polarizable alkane model based on Drude oscillators. Furthermore, the bulk dielectric is essentially captured in the infinite frequency, or optical, dielectric contribution. The FQ model is in this respect a more realistic force field for modeling lipid bilayer interiors for which most current state-of-the-art force fields do not accurately capture the dielectric environment. The molecular polarizability of the FQ model is 11.79 A3, in good agreement with the range of experimental and ab initio values. In contrast to FQ models of polar solvents such as alcohols and water, there was no need to scale gas-phase polarizabilities in order to avoid polarization catastrophes in the pure bulk. In terms of the liquid-vapor and liquid-liquid interfaces, the FQ model displays a rich orientational structure of alkane and water in the respective interfacial systems, in general conforming with earlier simulation studies of such interfaces. The FQ force field shows a marked deviation in the interfacial dipole potentials computed from the charge densities averaged over simulation trajectories. At the liquid-vapor interface, the FQ model predicts a potential drop of -178.71 mV in contrast to the C27r estimate of -433.80 mV. For the hexane-water interface, the FQ force field predicts a dipole potential drop of -379.40 mV in contrast to the C27r value of -105.42 mV. Although the surface dipole potential predicted by the FQ model is roughly 3.5 times that predicted by the C27r potential, it is consistent with reported experimental potentials across solvated lipid bilayers in the range of 400-600 mV.     

Computer simulation of electron transfer processes across the electrode|electrolyte interface: a treatment of solvent and electrode polarizability

A polarizable molecular dynamics model for adiabatic electron transfer across the electrode|electrolyte interface is presented. The electronic polarizability of the water and of the metal electrode is accounted for by a dynamical fluctuating charge algorithm, image charges, and the Ewald summation adapted for a conducting interface. The effects of the solvent electronic polarizability are studied by computing the diabatic and adiabatic free energy curves for both polarizable and non-polarizable water models. This represents the first effort to compute the adiabatic free energy curves from simulation for a fully polarizable electrochemical system.     

Polarizable empirical force field for aromatic compounds based on the classical drude oscillator

The polarizable empirical CHARMM force field based on the classical Drude oscillator has been extended to the aromatic compounds benzene and toluene. Parameters were optimized for benzene and then transferred directly to toluene, with parameters for the methyl moiety of toluene taken from the previously published work on the alkanes. Optimization of all parameters was performed against an extensive set of quantum mechanical and experimental data. Ab initio data was used for determination of the electrostatic parameters, for the vibrational analysis, and in the optimization of the relative magnitudes of the Lennard-Jones parameters. The absolute values of the Lennard-Jones parameters were determined by comparing computed and experimental heats of vaporization, molecular volumes, free energies of hydration, and dielectric constants. The newly developed parameter set was extensively tested against additional experimental data such as diffusion constants, heat capacities at constant pressure, and isothermal compressibilities including data as a function of temperature. Moreover, the structures of liquid benzene, liquid toluene, and solutions of each in water were studied. In the case of benzene, the computed and experimental total distribution function were compared, with the developed model shown to be in excellent agreement with experiment.     

Understanding the dielectric properties of liquid amides from a polarizable force field

The role played by electronic polarization in the dielectric properties of liquid N-methyl acetamide (NMA) is examined using molecular dynamics simulations with a polarizable force field based on classical Drude oscillators. The model presented is the first force field shown to reproduce the anomalously large dielectric constant of liquid NMA. Details of the molecular polarizability are found to be important. For instance, all elements of the polarizability tensor, rather then just the trace, impact on the condensed phase properties. Two factors related to electronic polarizability are found to contribute to this large dielectric constant. First is the significant enhancement of the mean amide molecular dipole magnitude, which is 50% larger in the liquid than in the gas phase. Second is the consequent strong hydrogen bonding between molecular neighbors that enhances the orientational alignment of the molecular dipoles. Polarizable models of amide compounds that have two (acetamide) and zero (N,N-dimethyl acetamide) polar hydrogen-bond donor atoms are also investigated. Experimentally, the neat liquid dielectric constants at 373 K are 100 for NMA, 66 for acetamide and 26 for N,N-dimethyl acetamide. The polarizable models replicate this trend, predicting a dielectric constant of 92+/-5 for NMA, 66+/-3 for acetamide and 23+/-1 for N,N-dimethyl acetamide.     

A polarizable force field for water using an artificial neural network

A force field for liquid water including polarization effects has been constructed using an artificial neural network (ANN). It is essential to include a many-body polarization effect explicitly into a potential energy function in order to treat liquid water which is dense and highly polar. The new potential energy function is a combination of empirical and nonempirical potentials. The TIP4P model was used for the empirical part of the potential. For the nonempirical part, an ANN with a back-propagation of error algorithm (BPNN) was introduced to reproduce the complicated many-body interaction energy surface from ab initio quantum mechanical calculations. BPNN, described in terms of a matrix, provides enough flexibility to describe the complex potential energy surface (PES). The structural and thermodynamic properties, calculated by isobaric-isothermal (constant-NPT) Monte Carlo simulations with the new polarizable force field for water, are compatible with experimental results. Thus, the simulation establishes the validity of using our estimated PES with a polarization effect for accurate predictions of liquid state properties. Applications of this approach are simple and systematic so that it can easily be applied to the development of other force fields besides the water-water system.     

Analytic energy gradient in combined time-dependent density functional theory and polarizable force field calculation

Formulas for evaluating analytic energy gradient are derived for combined time-dependent density functional theory (TDDFT) and polarizable force field methods that incorporate dipole polarizability tensors and linearly induced point dipoles. The Z-vector method for determining relaxed one-particle difference density matrix in regular TDDFT methods is extended to include induced dipoles. The analytic gradient of the mutual polarization energy of the force field and the TDDFT excited state can be formulated by using the TDDFT difference density-induced dipoles and the transition state density-induced dipoles. All the forces and torques involving induced dipoles can be efficiently evaluated using standard electrostatic formulas as if the induced dipoles were permanent dipoles. The formulas are given in the most general form and are applicable to various flavors of polarizable force fields. Implementation and tests with a polarizable five-point water model show that the formulas are rigorous. The carbonyl vibration modes and infrared spectrum intensities of a cluster formed by acetone and two water molecules are studied.     

The anharmonic potential energy surface of methyl fluoride

A large experimental spectroscopic data set sensitive to the cubic anharmonic potential energy surface (PES) of methyl fluoride has been compiled from the literature for six symmetric and asymmetric top isotopomers of methyl fluoride: 12CH3F, 13CH3FF, 12CD3F, 13CD3F, 12CHD2F and 12CH2DF. This empirical data set has been used to critically assess the best available literature ab initio cubic force field and various 'improved' theoretical force fields. A perturbation-resonance approach to the calculation of the observables from the force constants has been utilized and existing PESs were found to reproduce the data poorly. The careful treatment required for the correct theoretical reproduction of empirical spectroscopic constants arising from highly correlated least-squares fits to the original data is discussed. A new fit to the data has been performed (optimising 19 of the 38 cubic force constants) resulting in a much improved PES. The latter has been used to predict currently unknown spectroscopic constants for the least well-characterised isotopomer 12CH2DF. The prospects for a future empirical determination of the complete cubic force field of methyl fluoride are discussed and new data most likely to yield new information on the PES identified.     

Accounting for electronic polarization in non-polarizable force fields

The issues of electronic polarizability in molecular dynamics simulations are discussed. We argue that the charges of ionized groups in proteins, and charges of ions in conventional non-polarizable force fields such as CHARMM, AMBER, GROMOS, etc should be scaled by a factor about 0.7. Our model explains why a neglect of electronic solvation energy, which typically amounts to about a half of total solvation energy, in non-polarizable simulations with un-scaled charges can produce a correct result; however, the correct solvation energy of ions does not guarantee the correctness of ion-ion pair interactions in many non-polarizable simulations. The inclusion of electronic screening for charged moieties is shown to result in significant changes in protein dynamics and can give rise to new qualitative results compared with the traditional non-polarizable force field simulations. The model also explains the striking difference between the value of water dipole μ～ 3D reported in recent ab initio and experimental studies with the value μ(eff)～ 2.3D typically used in the empirical potentials, such as TIP3P or SPC/E. It is shown that the effective dipole of water can be understood as a scaled value μ(eff) = μ/√ε(el), where ε(el) = 1.78 is the electronic (high-frequency) dielectric constant of water. This simple theoretical framework provides important insights into the nature of the effective parameters, which is crucial when the computational models of liquid water are used for simulations in different environments, such as proteins, or for interaction with solutes.     

Molecule‐specific determination of atomic polarizabilities with the polarizable atomic multipole model

Recently, many polarizable force fields have been devised to describe induction effects between molecules. In popular polarizable models based on induced dipole moments, atomic polarizabilities are the essential parameters and should be derived carefully. Here, we present a parameterization scheme for atomic polarizabilities using a minimization target function containing both molecular and atomic information. The main idea is to adopt reference data only from quantum chemical calculations, to perform atomic polarizability parameterizations even when relevant experimental data are scarce as in the case of electronically excited molecules. Specifically, our scheme assigns the atomic polarizabilities of any given molecule in such a way that its molecular polarizability tensor is well reproduced. We show that our scheme successfully works for various molecules in mimicking dipole responses not only in ground states but also in valence excited states. The electrostatic potential around a molecule with an externally perturbing nearby charge also exhibits a near‐quantitative agreement with the reference data from quantum chemical calculations. The limitation of the model with isotropic atoms is also discussed to examine the scope of its applicability. © 2012 Wiley Periodicals, Inc.     

Parametrizing PCM to obtain solvation free energies from group contributions

A parametrization methodology for evaluating the solvation free energy, using the polarizable continuum model implemented in Gamess software, is presented in a formulation which makes use of a group contribution conception to construct the cavities. The systems studied include alkanes, alcohols, aldehydes and ketones embeded in a continuous medium simulating the water as the solvent. For each family, the CH2, OH, and C=O moieties of atoms are put together in single spheres forming a group. The cavities are constructed in two different ways, one for the electrostatic component and the other for nonelectrostatic contributions, i.e., the cavitation, dispersion, and repulsion components of free energy of solvation. A multivariate analysis is performed to obtain an assembly of variables, for each homologous series, able to give the results which are close to experiment. The analysis is addressed in order to (i) compare the theoretical free energy of solvation with the experimental trends of the solutes in aqueous media, when the chain is increased, (ii) compare the behavior of each component of free energy with the increasing CH2 number, (iii) investigate the influence of the oxygen atom on the components, and (iv) quantify the relative contribution of each component to the final free energy of solvation for some homologous series.     

Effect of polarization on the opsin shift in rhodopsins. 2. Empirical polarization models for proteins

The explicit treatment of polarization as a many-body interaction in condensed-phase systems represents a current problem in empirical force-field development. Although a variety of efficient models for molecular polarization have been suggested, polarizable force fields are still far from common use nowadays. In this work, we consider interactive polarization models employing Thole's short-range damping scheme and assess them for application on polypeptides. Despite the simplicity of the model, we find mean polarizabilities and anisotropies of amino acid side chains in excellent agreement with MP2/cc-pVQZ benchmark calculations. Combined with restrained electrostatic potential (RESP) derived atomic charges, the models are applied in a quantum-mechanical/molecular-mechanical (QM/MM) approach. An iterative scheme is used to establish a self-consistent mutual polarization between the QM and MM moieties. This ansatz is employed to study the influence of the protein polarizability on calculated optical properties of the protonated Schiff base of retinal in rhodopsin (Rh), bacterio-rhodopsin (bR), and pharaonis sensory rhodopsin II (psRII). The shifts of the excitation energy due to the instantaneous polarization response of the protein to the charge transfer on the retinal chromophore are quantified using the high level ab initio multireference spectroscopy-oriented configuration interaction (SORCI) method. The results are compared with those of previously published QM1/QM2/MM models for bR and psRII.