On the performance of molecular polarization methods. I. Water and carbon tetrachloride close to a point charge

Exploring optimization strategies for improving explicit water models: Rigid n-point model and polarizable model based on Drude oscillator

Rigid n-point water models are widely used in atomistic simulations, but have known accuracy drawbacks. Increasing the number of point charges, as well as adding electronic polarizability, are two common strategies for accuracy improvements. Both strategies come at considerable computational cost, which weighs heavily against modest possible accuracy improvements in practical simulations. In an effort to provide guidance for model development, here we have explored the limiting accuracy of "electrostatically globally optimal" n-point water models in terms of their ability to reproduce properties of water dimer-a mimic of the condensed state of water. For a given n, each model is built upon a set of reference multipole moments (e.g. ab initio) and then optimized to reproduce water dimer total dipole moment. The models are then evaluated with respect to the accuracy of reproducing the geometry of the water dimer. We find that global optimization of the charge distribution alone can deliver high accuracy of the water model: for n = 4 or n = 5, the geometry of the resulting water dimer can be almost within 50 of the ab initio reference, which is half that of the experimental error margin. Thus, global optimization of the charge distribution of classical n-point water models can lead to high accuracy models. We also find that while the accuracy improvement in going from n = 3 to n = 4 is substantial, the additional accuracy increase in going from n = 4 to n = 5 is marginal. Next, we have explored accuracy limitations of the standard practice of adding electronic polarizability (via a Drude particle) to a "rigid base"-pre-optimization rigid n-point water model. The resulting model (n = 3) shows a relatively small improvement in accuracy, suggesting that the strategy of merely adding the polarizability to an inferior accuracy water model used as the base cannot fix the defects of the latter. An alternative strategy in which the parameters of the rigid base model are globally optimized along with the polarizability parameter is much more promising: the resulting 3-point polarizable model out-performs even the 5-point optimal rigid model by a large margin. We suggest that future development efforts consider 3- and 4-point polarizable models where global optimization of the "rigid base" is coupled to optimization of the polarizability to deliver globally optimal solutions.

Capturing Many-Body Interactions with Classical Dipole Induction Models

The nonadditive many-body interactions are significant for structural and thermodynamic properties of condensed phase systems. In this work we examined the many-body interaction energy of a large number of common organic/biochemical molecular clusters, which consist of 18 chemical species and cover nine common organic elements, using the Møller-Plesset perturbation theory to the second order (MP2) [ Møller et al. Phys. Rev. 1934 , 46 , 618 . ]. We evaluated the capability of Thole-based dipole induction models to capture the many-body interaction energy. Three models were compared: the original model and parameters used by the AMOEBA force field, a variation of this original model where the damping parameters have been reoptimized to MP2 data, and a third model where the damping function form applied to the permanent electric field is modified. Overall, we find the simple classical atomic dipole models are able to capture the 3- and 4-body interaction energy across a wide variety of organic molecules in various intermolecular configurations. With modified Thole models, it is possible to further improve the agreement with MP2 results. These models were also tested on systems containing metal/halogen ions to examine the accuracy and transferability. This work suggests that the form of damping function applied to the permanent electrostatic field strongly affects the distance dependence of polarization energy at short intermolecular separations.

A transferable ab initio based force field for aqueous ions

We present a new polarizable force field for aqueous ions (Li(+), Na(+), K(+), Rb(+), Cs(+), Mg(2 +), Ca(2 +), Sr(2 +), and Cl(-)) derived from condensed phase ab initio calculations. We use maximally localized Wannier functions together with a generalized force and dipole-matching procedure to determine the whole set of parameters. Experimental data are then used only for validation purposes and a good agreement is obtained for structural, dynamic, and thermodynamic properties. The same procedure applied to crystalline phases allows to parametrize the interaction between cations and the chloride anion. Finally, we illustrate the good transferability of the force field to other thermodynamic conditions by investigating concentrated solutions.

Physical properties of small water clusters in low and moderate electric fields

Likely candidates for the lowest minima of water clusters (H(2)O)(N) for N ≤ 20 interacting with a uniform electric field strength in the range E ≤ 0.6 V/Å have been identified using basin-hopping global optimization. Two water-water model potentials were considered, namely TIP4P and the polarizable Dang-Chang potential. The two models produce some consistent results but also exhibit significant differences. The cluster internal energy and dipole moment indicate two varieties of topological transition in the structure of the global minimum as the field strength is increased. The first takes place at low field strengths (0.1 V/Å<E < 0.2 V/Å) and reorganizes the hydrogen-bonds to orient the water permanent dipoles along the field. The second type of transition occurs at larger field strengths (0.3 V/Å<E < 0.5 V/Å) and corresponds to an extensive structural reorganization, where several hydrogen-bonds break as the cluster stretches along the field direction, the larger clusters (N > 10) usually forming helical structures.

Polarizable Atomic Multipole-based Molecular Mechanics for Organic Molecules

An empirical potential based on permanent atomic multipoles and atomic induced dipoles is reported for alkanes, alcohols, amines, sulfides, aldehydes, carboxylic acids, amides, aromatics and other small organic molecules. Permanent atomic multipole moments through quadrupole moments have been derived from gas phase ab initio molecular orbital calculations. The van der Waals parameters are obtained by fitting to gas phase homodimer QM energies and structures, as well as experimental densities and heats of vaporization of neat liquids. As a validation, the hydrogen bonding energies and structures of gas phase heterodimers with water are evaluated using the resulting potential. For 32 homo- and heterodimers, the association energy agrees with ab initio results to within 0.4 kcal/mol. The RMS deviation of hydrogen bond distance from QM optimized geometry is less than 0.06 Å. In addition, liquid self-diffusion and static dielectric constants computed from molecular dynamics simulation are consistent with experimental values. The force field is also used to compute the solvation free energy of 27 compounds not included in the parameterization process, with a RMS error of 0.69 kcal/mol. The results obtained in this study suggest the AMOEBA force field performs well across different environments and phases. The key algorithms involved in the electrostatic model and a protocol for developing parameters are detailed to facilitate extension to additional molecular systems.

The polarizable point dipoles method with electrostatic damping: implementation on a model system

Recently, the use of polarizable force fields in Molecular Dynamics simulations has been gaining importance, since they allow a better description of heterogeneous systems compared to simple point charges force fields. Among the various techniques developed in the last years the one based on polarizable point dipoles represents one of the most used. In this paper, we review the basic technical issues of the method, illustrating the way to implement intramolecular and intermolecular damping of the electrostatic interactions, either with and without the Ewald summation method. We also show how to reduce the computational overhead for evaluating the dipoles, introducing to the state-of-the-art methods: the extended Lagrangian method and the always stable predictor corrector method. Finally we discuss the importance of screening the electrostatic interactions at short range, defending this technique against simpler approximations usually made. We compare results of density functional theory and classical force field-based Molecular Dynamics simulations of chloride in water.

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 modelling 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 Langrangian 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(4), 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.

Folding of gas-phase polyalanines in a static electric field: alignment, deformations, and polarization effects

Monte Carlo simulations of the temperature-induced unfolding of small gas-phase polyalanines in a static, homogeneous electric field are reported, based on the AMBER ff96 force field. The peptides exhibit a structural transition from the native alpha-helix state to entropically favored beta-sheet conformations, before eventually turning to extended coil at higher temperatures. Upon switching the electric field, the molecules undergo preferential alignment of their dipole moment vector toward the field axis and a shift of the alpha-beta transition to higher temperatures. At higher field strengths (>10(8) V/m) the molecules stretch and the alpha-beta and beta-coil transitions merge. A simple three-state model is shown to account for the observed behavior. Under even higher fields, density functional theory calculations and a polarizable force field both show that electronic rearrangements tend to further increase the dipole moment, polarization effects being approximately half in magnitude with respect to stretching effect. Finally a tentative (temperature, field-strength) phase diagram is sketched.

Simulation of Ca2+ and Mg2+ solvation using polarizable atomic multipole potential

The alkaline earth metals calcium and magnesium are critically involved in many biomolecular processes. To understand the hydration thermodynamics of these ions, we have performed molecular dynamics simulations using a polarizable potential. Particle-mesh Ewald for point multipoles has been applied to the calculation of electrostatic interactions. The parameters in this model have been determined from an ab initio quantum mechanical calculation of dimer interactions between ions and water. Two methods for ion solvation free energy calculation, free energy perturbation, and the Bennett acceptance ratio have been compared. Both predict results consistent with other theoretical estimations while the Bennett approach leads to a much smaller statistical error. Based on the Born theory and the ion-oxygen radial distribution functions, we estimate the effective size of the ions in solution, concluding that K(+) > Na(+) congruent with Ca(2+) > Mg(2+). There appears to be much stronger perturbation in water structure, dynamics, and dipole moment around the divalent cations than the monovalent K(+) and Na(+). The average water coordination numbers for Ca(2+) and Mg(2+) are 7.3 and 6, respectively. The lifetime of water molecules in the first solvation shell of Mg(2+) is on the order of hundreds of picoseconds, in contrast to only few picoseconds for Ca(2+), K(+), or Na(+).

Key Role of the Polarization Anisotropy of Water in Modeling Classical Polarizable Force Fields

We have evaluated the extent to which classical polarizable force fields, based either on the chemical potential equalization principle or on distributed polarizabilities in the framework of the Sum of Interactions Between Fragments Ab initio computed (SIBFA), can reproduce the ab initio polarization energy and the dipole moment of three distinct water oligomers: bifurcated chains, transverse hydrogen-bonded chains, and longitudinal hydrogen-bonded chains of helical shape. To analyze the many-body polarization effect, chains of different size, i.e., from 2 to 12 water monomers, have been considered. Although the dipole moment is a well-defined quantity in both classical polarizable models and quantum mechanical methods, polarization energy can be defined unequivocally only in the former type of approaches. In this study we have used the Kitaura-Morokuma (KM) procedure. Although the KM approach is on the one hand known to overestimate the polarization energy for strongly interacting molecules, on the other hand it can account for the many-body polarization effectively, whereas some other procedures do not. Our data show that, if off-centered lone pair polarizabilities are explicitly represented, classical polarizable force fields can afford a close agreement with the ab initio results, both in terms of polarization energy and in terms of dipole moment.

Behavior of polarizable models in presence of strong electric fields. I. Origin of nonlinear effects in water point-charge systems

In the current opinion, the inclusion of polarization response in classical computer simulations is considered as one of the most important and urgent improvements to be implemented in modern empirical potential models. In this work we focus on the capability of polarizable models, based on the pairwise Coulomb interactions, to model systems where strong electric fields enter into play. As shown by Masia, Probst, and Rey (MPR) [in J. Chem. Phys. 121, 7362 (2004)], when a molecule interacts with point charges, polarizable models show underpolarization with respect to ab initio methods. We prove that this underpolarization, clearly related to nonlinear polarization effects, cannot be simply ascribed to the lack of hyperpolarization in the polarizable models, as suggested by MPR. Analysis of the electron-density rearrangement induced on a water molecule by a point charge reveals a twofold level of polarization response. One level involves intramolecular charge transfer on the whole molecular volume, with the related polarization exhibiting a seemingly linear behavior with the external electric field. The other nonlinear polarization level occurs only at strong electric fields and is found to be strictly correlated to the quantum-mechanical nature of the water molecule. The latter type of polarization has a local character, being limited to the space region of the water lone pairs.

On the performance of molecular polarization methods. II. Water and carbon tetrachloride close to a cation

Our initial study on the performance of molecular polarization methods close to a positive point charge [M. Masia, M. Probst, and R. Rey, J. Chem. Phys. 121, 7362 (2004)] is extended to the case in which a molecule interacts with a real cation. Two different methods (point dipoles and shell model) are applied to both the ion and the molecule. The results are tested against high-level ab initio calculations for a molecule (water or carbon tetrachloride) close to Li+, Na+, Mg2+, and Ca2+. The monitored observable is in all cases the dimer electric dipole as a function of the ion-molecule distance for selected molecular orientations. The moderate disagreement previously obtained for point charges at intermediate distances, and attributed to the linearity of current polarization methods (as opposed to the nonlinear effects evident in ab initio calculations), is confirmed for real cations as well. More importantly, it is found that at short separations the phenomenological polarization methods studied here substantially overestimate the dipole moment induced if the ion is described quantum chemically as well, in contrast to the dipole moment induced by a point-charge ion, for which they show a better degree of accord with ab initio results. Such behavior can be understood in terms of a decrease of atomic polarizabilities due to the repulsion between electronic charge distributions at contact separations. It is shown that a reparametrization of the Thole method for damping of the electric field, used in conjunction with any polarization scheme, allows to satisfactorily reproduce the dimer dipole at short distances. In contrast with the original approach (developed for intramolecular interactions), the present reparametrization is ion and method dependent, and corresponding parameters are given for each case.

Calculation of the free energy of polarization: quantifying the effect of explicitly treating electronic polarization on the transferability of force-field parameters

The lack of an explicit description of electronic polarization in nonpolarizable force fields usually results in an incomplete transferability of force-field parameter sets when applied in simulations of the system of interest in either a polar or an apolar environment. For example, the use of nonpolarizable parameter sets optimized to reproduce experimental data on properties of pure liquids of polar compounds commonly yields too low solubilities in water for the corresponding compounds. The reason is that the fixed charge distributions calibrated for the pure liquid might correspond to too low molecular dipole moments in case of hydration. In the current study, we quantitatively show that explicit inclusion of electronic polarization can improve the transferability of biomolecular force-field parameter sets. With this aim, free energies of polarization, DeltaGpola, have been calculated, with DeltaGpola corresponding to the free energy difference between identical systems described by a polarizable and a nonpolarizable model. Using a nonpolarizable model and a polarizable one (based on the charge-on-spring approach) for dimethyl ether (DME), which were both parametrized to reproduce experimental values for pure liquid properties, small values were found for DeltaGpola for the pure liquid or when a DME solute was solvated in the apolar solvent cyclohexane. For the solute hydrated in water, however, DeltaGpola was found to be of the same order of magnitude as the discrepancy between the free energy of hydration from simulation using a nonpolarizable solute model and the experimental value. Thus, introducing polarizabilities clearly improves the transferability of the parameter set. Additionally, in calculations of an anion solvated in DME, DeltaGpola for the solvent adopted relatively large values. From an estimation of the errors in the calculated free energy differences, it was furthermore shown that the calculation of DeltaGpola offers an effective and accurate method to obtain differences in solvation (or excess) free energies between systems described by polarizable and nonpolarizable models when compared to a direct calculation of solvation (or excess) free energies.

Development of Many−Body Polarizable Force Fields for Li-Battery Components: 1. Ether, Alkane, and Carbonate-Based Solvents

Classical many-body polarizable force fields were developed for n-alkanes, perflouroalkanes, polyethers, ketones, and linear and cyclic carbonates on the basis of quantum chemistry dimer energies of model compounds and empirical thermodynamic liquid-state properties. The dependence of the electron correlation contribution to the dimer binding energy on basis-set size and level of theory was investigated as a function of molecular separation for a number of alkane, ether, and ketone dimers. Molecular dynamics (MD) simulations of the force fields accurately predicted structural, dynamic, and transport properties of liquids and unentangled polymer melts. On average, gas-phase dimer binding energies predicted with the force field were between those from MP2/aug-cc-pvDz and MP2/aug-cc-pvTz quantum chemistry calculations.

On the coupling between molecular diffusion and solvation shell exchange

The connection between diffusion and solvent exchanges between first and second solvation shells is studied by means of molecular dynamics simulations and analytic calculations, with detailed illustrations for water exchange for the Li(+) and Na(+) ions, and for liquid argon. First, two methods are proposed which allow, by means of simulation, to extract the quantitative speed-up in diffusion induced by the exchange events. Second, it is shown by simple kinematic considerations that the instantaneous velocity of the solute conditions to a considerable extent the character of the exchanges. Analytic formulas are derived which quantitatively estimate this effect, and which are of general applicability to molecular diffusion in any thermal fluid. Despite the simplicity of the kinematic considerations, they are shown to well describe many aspects of solvent exchange/diffusion coupling features for nontrivial systems.