A unified theoretical framework for fluctuating-charge models in atom-space and in bond-space

Completely Multipolar Model as a General Framework for Many-Body Interactions as Illustrated for Water

We introduce a general framework for many-body force fields, the Completely Multipolar Model (CMM), that utilizes multipolar electrical moments modulated by exponential decay of electron density as a common functional form for all terms of an energy decomposition analysis of intermolecular interactions. With this common functional form, the CMM model establishes well-formulated damped tensors that reach the correct asymptotes at both long- and short-range while formally ensuring no short-range catastrophes. CMM describes the separable EDA terms of dispersion, exchange polarization, and Pauli repulsion with short-ranged anisotropy, polarization as intramolecular charge fluctuations and induced dipoles, while charge transfer describes explicit movement of charge between molecules, and naturally describes many-body charge transfer by coupling into the polarization equations. We also utilize a new one-body potential that accounts for intramolecular polarization by including an electric field-dependent correction to the Morse potential to ensure that CMM reproduces all physically relevant monomer properties including the dipole moment, molecular polarizability, and dipole and polarizability derivatives. The quality of CMM is illustrated through agreement of individual terms of the EDA and excellent extrapolation to energies and geometries of an extensive validation set of water cluster data.

Effect of Electric Fields on the Decomposition of Phosphate Esters

Phosphate esters decompose on metal surfaces and form protective polyphosphate films. For many applications, such as in lubricants for electric vehicles and wind turbines, an understanding of the effect of electric fields on molecular decomposition is urgently required. Experimental investigations have yielded contradictory results, with some suggesting that electric fields improve tribological performance, while others have reported the opposite effect. Here, we use nonequilibrium molecular dynamics (NEMD) simulations to study the decomposition of tri-n-butyl phosphate (TNBP) molecules nanoconfined between ferrous surfaces (iron and iron oxide) under electrostatic fields. The reactive force field (ReaxFF) method is used to model the effects of chemical bonding and molecular dissociation. We show that the charge transfer with the polarization current equalization (QTPIE) method gives more realistic behavior compared to the standard charge equilibration (QEq) method under applied electrostatic fields. The rate of TNBP decomposition via carbon-oxygen bond dissociation is faster in the nanoconfined systems than that in the bulk due to the catalytic action of the surfaces. In all cases, the application of an electric field accelerates TNBP decomposition. When electric fields are applied to the confined systems, the phosphate anions are pulled toward the surface with high electric potential, while the alkyl cations are pulled to the surface with lower potential, leading to asymmetric film growth. Analysis of the temperature- and electric field strength-dependent dissociation rate constants using the Arrhenius equation suggests that, on reactive iron surfaces, the increased reactivity under an applied electric field is driven mostly by an increase in the pre-exponential factor, which is linked to the number of molecule-surface collisions. Conversely, the accelerated decomposition of TNBP on iron oxide surfaces can be attributed to a reduction in the activation energy with increasing electric field strength. Single-molecule nudged-elastic band (NEB) calculations also show a linear reduction in the energy barrier for carbon-oxygen bond breaking with electric field strength, due to stabilization of the charged transition state. The simulation results are consistent with experimental observations of enhanced and asymmetric tribofilm growth under electrostatic fields.

Unifying Charge-Flow Polarization Models

We show that several models where electric polarization in molecular systems is modeled by charge-flow between atoms can all be considered as different manifestations of a general underlying mathematical structure. The models can be classified according to whether they employ atomic or bond parameters and whether they employ atom/bond hardness or softness. We show that an ab initio calculated charge response kernel can be considered as the inverse screened Coulombic matrix projected onto the zero-charge subspace, and this may provide a method for deriving charge screening functions to be used in force fields. The analysis suggests that some models contain redundancies, and we argue that a parameterization of charge-flow models in terms of bond softness is preferable as it depends on local quantities and decay to zero upon bond dissociation, while bond hardness depends on global quantities and increases toward infinity upon bond dissociation.

Assessing the quality of QM/MM approaches to describe vacuo-to-water solvatochromic shifts

The performance of different quantum mechanics/molecular mechanics embedding models to compute vacuo-to-water solvatochromic shifts is investigated. In particular, both nonpolarizable and polarizable approaches are analyzed and computed results are compared to reference experimental data. We show that none of the approaches outperform the others and that errors strongly depend on the nature of the molecular transition to be described. Thus, we prove that the best choice of embedding model highly depends on the molecular system and that the use of a specific approach as a black box can lead to significant errors and, sometimes, totally wrong predictions.

Using atomic charges to model molecular polarization

We review different models for introducing electric polarization in force fields, with special focus on methods where polarization is modelled at the atomic charge level. While electric polarization has been included in several force fields, the common approach has been to focus on atomic dipole polarizability. Several approaches allow modelling electric polarization by using charge-flow between charge sites instead, but this has been less exploited, despite that atomic charges and charge-flow is expected to be more important than atomic dipoles and dipole polarizability. A number of challenges are required to be solved for charge-flow models to be incorporated into polarizable force fields, for example how to parameterize the models and how to make them computational efficient.

Water in an External Electric Field: Comparing Charge Distribution Methods Using ReaxFF Simulations

The growing interest in the effects of external electric fields on reactive processes requires predictive methods that can reach longer length and time scales than quantum mechanical simulations. Recently, many studies have included electric fields in ReaxFF, a widely used reactive molecular dynamics method. In the case of modeling an external electric field, the charge distribution method used in ReaxFF is critical. The most common charge distribution method used in previous studies of electric fields is the charge equilibration (QEq) method, which assumes that the system is a contiguous conductor and that charge transfer can occur across any distance. In contrast, many systems of interest are insulators or semiconductors, and long-distance charge transfer should not occur in response to a small difference in potential. This study focuses on the limitations of the QEq method in the context of water in an external electric field. We demonstrate that QEq can predict unphysical charge distributions and exhibits properties that do not converge as a function of system size. Furthermore, we show that electric fields within the recently developed atom-condensed Kohn-Sham density functional theory (DFT) approximated to the second-order (ACKS2) approach address the major limitations of electric fields in QEq. With ACKS2, we observe more physical charge distributions and properties that converge as a function of system size. We do not suggest that ACKS2 is perfect in all circumstances but rather show specific cases where it addresses the major shortcomings of QEq in the context of an external electric field.

Generalized Method for Charge-Transfer Equilibration in Reactive Molecular Dynamics

Variable charge models (e.g., electronegativity equalization method (EEM), charge equilibration (QEq), electrostatic plus (ES+)) used in reactive molecular dynamics simulations often inherently impose a global charge transfer between atoms (approximating each system as an ideal metal). Consequently, most surface processes (e.g., adsorption, desorption, deposition, sputtering) are affected, potentially causing dubious dynamics. This issue has been addressed by certain split charge variants (i.e., split charge equilibration (SQE), redoxSQE) through a distance-dependent bond hardness, by the atomic charge ACKS2 and QTPIE models, which are based on the Kohn-Sham density functional theory, as well as by an electronegativity screening extension to the QEq model (approximating each system as an ideal insulator). In a brief review of the QEq and the QTPIE model, their applicability for studying surface interactions is assessed in this work. Following this evaluation, a revised generalization of the QEq and QTPIE models is proposed and formulated, called the charge-transfer equilibration model or in short the QTE model. This method is based on the equilibration of charge-transfer variables, which locally constrain the split charge transfer per unit time (i.e., due to overlapping orbitals) without any kind of bond hardness specification. Furthermore, a formalism relying solely on atomic charges is obtained by a respective transformation, employing an extended Lagrangian method. We moreover propose a mirror boundary condition and its implementation to accelerate surface investigations. The models proposed in this work facilitate reactive molecular dynamics simulations, which describe various materials and surface phenomena appropriately.

Quantum Mechanical and Molecular Mechanics Modeling of Membrane-Embedded Rhodopsins

Computational chemistry provides versatile methods for studying the properties and functioning of biological systems at different levels of precision and at different time scales. The aim of this article is to review the computational methodologies that are applicable to rhodopsins as archetypes for photoactive membrane proteins that are of great importance both in nature and in modern technologies. For each class of computational techniques, from methods that use quantum mechanics for simulating rhodopsin photophysics to less-accurate coarse-grained methodologies used for long-scale protein dynamics, we consider possible applications and the main directions for improvement.

Molecular Dynamics Simulations of Ionic Liquids and Electrolytes Using Polarizable Force Fields

Many applications in chemistry, biology, and energy storage/conversion research rely on molecular simulations to provide fundamental insight into structural and transport properties of materials with high ionic concentrations. Whether the system is comprised entirely of ions, like ionic liquids, or is a mixture of a polar solvent with a salt, e.g., liquid electrolytes for battery applications, the presence of ions in these materials results in strong local electric fields polarizing solvent molecules and large ions. To predict properties of such systems from molecular simulations often requires either explicit or mean-field inclusion of the influence of polarization on electrostatic interactions. In this manuscript, we review the pros and cons of different treatments of polarization ranging from the mean-field approaches to the most popular explicit polarization models in molecular dynamics simulations of ionic materials. For each method, we discuss their advantages and disadvantages and emphasize key assumptions as well as their adjustable parameters. Strategies for the development of polarizable models are presented with a specific focus on extracting atomic polarizabilities. Finally, we compare simulations using polarizable and nonpolarizable models for several classes of ionic systems, discussing the underlying physics that each approach includes or ignores, implications for implementation and computational efficiency, and the accuracy of properties predicted by these methods compared to experiments.

Describing Molecular Polarizability by a Bond Capacity Model

We propose a bond capacity model for describing molecular polarization in force field energy functions at the charge-only level. Atomic charges are calculated by allowing charge to flow between atom pairs according to a bond capacity and a difference in electrostatic potential. The bond capacity is closely related to the bond order and decays to zero as the bond distance is increased. The electrostatic potential is composed of an intrinsic potential, identified as the electronegativity, and a screened Coulomb potential from all other charges. The bond capacity model leads to integer fragment charges upon bond dissociation and displays linear scaling of the polarizability with system size. Bond capacity parameters can be derived from reference molecular polarizabilities, while electronegativity parameters can be derived from reference atomic charges or a reference molecular electrostatic potential. Out-of-plane polarization for planar systems is modeled by off-nuclei charge sites. The model is shown to be able to reproduce anisotropic reference molecular polarizabilities with an accuracy of ∼10% using a limited set of bond capacity parameters and can describe both inter- and intramolecular polarization.

Mechanisms of Iodide⁻Triiodide Exchange Reactions in Ionic Liquids: A Reactive Molecular-Dynamics Exploration

Efficient charge transport has been observed in iodine-doped, iodide-based room-temperature ionic liquids, yielding high ionic conductivity. To elucidate preferred mechanistic pathways for the iodide ( I - )-to-triiodide ( I 3 - ) exchange reactions, we have performed 10 ns reactive molecular-dynamics calculations in the liquid state for 1-butyl-3-methylimidazolium iodide ([BMIM][I]) at 450 to 750 K. Energy-barrier distributions for the iodine-swapping process were determined as a function of temperature, employing a charge-reassignment scheme drawn in part from electronic-structure calculations. Bond-exchange events were observed with rate-determining energy barriers ranging from ~0.19 to 0.23 ± 0.06 eV at 750 and 450 K, respectively, with an approximately Arrhenius temperature dependence for iodine self-diffusivity and reaction kinetics, although diffusion dominates/limits the bond-exchange events. This charge transfer is not dissimilar in energetics to those in solid-state superionic conductors.

An atomistic model for the charge distribution in layered MoS2

We present an atomistic model for predicting the distribution of doping electric charges in layered molybdenum disulfide (MoS₂). This model mimics the charge around each ion as a net Gaussian-spatially distributed charge plus an induced dipole, and is able to predict the distribution of doping charges in layered MoS₂ in a self-consistent scheme. The profiles of doping charges in monolayer MoS₂ flakes computed by this charge-dipole model are in good agreement with those obtained by density-functional-theory calculations. Using this model, we quantitatively predict the charge enhancement in MoS₂ monolayer nanoribbons, with which strong ionic charge-localization effects are shown.

Study of interactions between metal ions and protein model compounds by energy decomposition analyses and the AMOEBA force field

The interactions between metal ions and proteins are ubiquitous in biology. The selective binding of metal ions has a variety of regulatory functions. Therefore, there is a need to understand the mechanism of protein-ion binding. The interactions involving metal ions are complicated in nature, where short-range charge-penetration, charge transfer, polarization, and many-body effects all contribute significantly, and a quantitative description of all these interactions is lacking. In addition, it is unclear how well current polarizable force fields can capture these energy terms and whether these polarization models are good enough to describe the many-body effects. In this work, two energy decomposition methods, absolutely localized molecular orbitals and symmetry-adapted perturbation theory, were utilized to study the interactions between Mg²⁺/Ca²⁺ and model compounds for amino acids. Comparison of individual interaction components revealed that while there are significant charge-penetration and charge-transfer effects in Ca complexes, these effects can be captured by the van der Waals (vdW) term in the AMOEBA force field. The electrostatic interaction in Mg complexes is well described by AMOEBA since the charge penetration is small, but the distance-dependent polarization energy is problematic. Many-body effects were shown to be important for protein-ion binding. In the absence of many-body effects, highly charged binding pockets will be over-stabilized, and the pockets will always favor Mg and thus lose selectivity. Therefore, many-body effects must be incorporated in the force field in order to predict the structure and energetics of metalloproteins. Also, the many-body effects of charge transfer in Ca complexes were found to be non-negligible. The absorption of charge-transfer energy into the additive vdW term was a main source of error for the AMOEBA many-body interaction energies.

A charge equilibration formalism for treating charge transfer effects in MD simulations: Application to water clusters

Conventional classical force fields by construction do not explicitly partition intermolecular interactions to include polarization and charge transfer effects, whereas fully quantum mechanical treatments allow a means to effect this dissection (although not uniquely due to the lack of a charge transfer operator). Considering the importance of polarization in a variety of systems, a particular class of classical models, charge equilibration models, have been extensively developed to study those systems; since these types of interaction models are inherently based on movement of charge throughout a system, they are natural platform for including polarization and charge transfer effects within the context of molecular simulations. Here, we present two bond-space charge equilibration models we term as QE2 and mixed QE2 treat charge transfer in classical molecular mechanical calculations those provide practical solutions to two major drawbacks of charge equilibration models: (a) a nonvanishing amount of charge transfer between two heteroatoms at large separations, and (b) superlinear polarizability scaling during bond dissociation due to charge transfer over unphysical, large distances. To control charge transfer during dissociation of a bond in a molecular system, we introduce a distance-dependent scaling function (QE2 model) which, controls and recovers physical behavior of the homonuclear and heteronuclear charge transfer between two atoms at small and large values of internuclear separation; and the mixed QE2 model in which we combine the QE2 model under allow and disallow charge transfer situations that describe both charge transfer and polarizability in a distance-dependent manner. We demonstrate the utility of both models in the case of a water dimer, and compare the results with other existing models, and further, we perform short molecular dynamics simulations for few water clusters with the QE2 model to show the charge transfer and internuclear separation are correlated in dynamics. © 2017 Wiley Periodicals, Inc.

A New Split Charge Equilibration Model and REPEAT Electrostatic Potential Fitted Charges for Periodic Frameworks with a Net Charge

Periodic frameworks that possess a net charge, such as zeolites, are an important class of materials in wide use. For guest-host interactions to be simulated in these materials, partial atomic charges are often used. In this work, we investigate two methods for the generation of partial atomic charges in periodic systems having a net framework charge. We first examine the validity of generating REPEAT electrostatic potential fitted charges derived from periodic electronic structure calculations, where a constant background charge is added to neutralize the net charge on the framework. The constant background charge obviates the need to add neutralizing counterions, which may induce artifacts such as polarization in the infinite periodic system. The second method we explore is the split charge equilibration (SQE) method for the rapid generation of partial atomic charges. The original formulation of the SQE method cannot be applied to systems with a net charge. In this work, we reformulate the SQE method by transforming the split charges into an atomic charge basis that allows for non-neutral systems to be treated. The new SQE model, which we call SQE_AB (for atomic basis), was validated with a series of tests using both charged and neutral metal organic frameworks and zeolites. It was shown that SQE_AB gives equivalent results to those of the original SQE model for neutral systems. We then demonstrated that the SQE_AB method is able to "capture" the chemical structure of a charged framework better than that of the charge equilibration model by comparing to REPEAT electrostatic potential fitted charges.

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.

An Estimation of Hybrid Quantum Mechanical Molecular Mechanical Polarization Energies for Small Molecules Using Polarizable Force-Field Approaches

In this work, we report two polarizable molecular mechanics (polMM) force field models for estimating the polarization energy in hybrid quantum mechanical molecular mechanical (QM/MM) calculations. These two models, named the potential of atomic charges (PAC) and potential of atomic dipoles (PAD), are formulated from the ab initio quantum mechanical (QM) response kernels for the prediction of the QM density response to an external molecular mechanical (MM) environment (as described by external point charges). The PAC model is similar to fluctuating charge (FQ) models because the energy depends on external electrostatic potential values at QM atomic sites; the PAD energy depends on external electrostatic field values at QM atomic sites, resembling induced dipole (ID) models. To demonstrate their uses, we apply the PAC and PAD models to 12 small molecules, which are solvated by TIP3P water. The PAC model reproduces the QM/MM polarization energy with a R² value of 0.71 for aniline (in 10,000 TIP3P water configurations) and 0.87 or higher for other 11 solute molecules, while the PAD model has a much better performance with R² values of 0.98 or higher. The PAC model reproduces reference QM/MM hydration free energies for 12 solute molecules with a RMSD of 0.59 kcal/mol. The PAD model is even more accurate, with a much smaller RMSD of 0.12 kcal/mol, with respect to the reference. This suggests that polarization effects, including both local charge distortion and intramolecular charge transfer, can be well captured by induced dipole type models with proper parametrization.

Local electric field factors by a combined charge-transfer and point–dipole interaction model

A model for the local electric field as a linear response to a frequency-dependent external electric field is presented based on a combined charge-transfer and point–dipole interaction force-field model.

Symmetry laws improve electronegativity equalization by orders of magnitude and call for a paradigm shift in conceptual density functional theory

The strict Wigner-Witmer symmetry constraints on chemical bonding are shown to determine the accuracy of electronegativity equalization (ENE) to a high degree. Bonding models employing the electronic chemical potential, μ, as the negative of the ground-state electronegativity, χ(GS), frequently collide with the Wigner-Witmer laws in molecule formation. The violations are presented as the root of the substantially disturbing lack of chemical potential equalization (CPE) in diatomic molecules. For the operational chemical potential, μ(op), the relative deviations from CPE fall between -31% ≤ δμ(op) ≤ +70%. Conceptual density functional theory (cDFT) cannot claim to have operationally (not to mention, rigorously) proven and unified the CPE and ENE principles. The solution to this limitation of cDFT and the symmetry violations is found in substituting μ(op) (i) by Mulliken's valence-state electronegativity, χ(M), for atoms and (ii) its new generalization, the valence-pair-affinity, α(VP), for diatomic molecules. Mulliken's χ(M) is equalized into the α(VP) of the bond, and the accuracy of ENE is orders of magnitude better than that of CPE using μ(op). A paradigm shift replacing the dominance of ground states by emphasizing valence states seems to be in order for conceptual DFT.

Systematic Parametrization of Polarizable Force Fields from Quantum Chemistry Data

We introduce ForceBalance, a method and free software package for systematic force field optimization with the ability to parametrize a wide variety of functional forms using flexible combinations of reference data. We outline several important challenges in force field development and how they are addressed in ForceBalance, and present an example calculation where these methods are applied to develop a highly accurate polarizable water model. ForceBalance is available for free download at https://simtk.org/home/forcebalance.

A periodic charge-dipole electrostatic model: parametrization for silver slabs

We present an extension of the charge-dipole model for the description of periodic systems. This periodic charge-dipole electrostatic model (PCDEM) allows one to describe the linear response of periodic structures in terms of charge- and dipole-type gaussian basis functions. The long-range electrostatic interaction is efficiently described by means of the continuous fast multipole method. As a first application, the PCDEM method is applied to describe the polarizability of silver slabs. We find that for a correct description of the polarizability of the slabs both charges and dipoles are required. However a continuum set of parametrizations, i.e., different values of the width of charge- and dipole-type gaussians, leads to an equivalent and accurate description of the slabs polarizability but a completely unphysical description of induced charge-density inside the slab. We introduced the integral squared density measure which allows one to obtain a unique parametrization which accurately describes both the polarizability and the induced density profile inside the slab. Finally the limits of the electrostatic approximations are also pointed out.

A wave-function based approach for polarizable charge model: Systematic comparison of polarization effects on protic, aprotic, and ionic liquids

We first describe a wave-function based formalism of polarizable charge model by starting from the Hartree product ansatz for the total wave function and making the second-order expansion of individual molecular energies with the use of partial charge operators. The resulting model is shown to be formally equivalent to the charge response kernel model that starts from the linear-response approximation to partial charges, and also closely related to a family of fluctuating charge models that are based on the electronegativity equalization principle. We then apply the above model to a systematic comparison of polarization effects on qualitatively different liquids, namely, protic solvents (water and methanol), an aprotic polar solvent (acetonitrile), and imidazolium-based ionic liquids. Electronic polarization is known to decelerate molecular motions in conventional solvents while it accelerates them in ionic liquids. To obtain more insights into these phenomena, we consider an effective decomposition of total polarization energy into molecular contributions, and show that their statistical distribution is well-correlated with the acceleration/deceleration of molecular motions. In addition, we perform effective nonpolarizable simulations based on mean polarized charges, and compare them with fully polarizable simulations. The result shows that the former can reproduce structural properties of conventional solvents rather accurately, while they fail qualitatively to reproduce acceleration of molecular motions in ionic liquids.

Nonmetallic electronegativity equalization and point-dipole interaction model including exchange interactions for molecular dipole moments and polarizabilities

The electronegativity equalization model (EEM) has been combined with a point-dipole interaction model to obtain a molecular mechanics model consisting of atomic charges, atomic dipole moments, and two-atom relay tensors to describe molecular dipole moments and molecular dipole-dipole polarizabilities. The EEM has been phrased as an atom-atom charge-transfer model allowing for a modification of the charge-transfer terms to avoid that the polarizability approaches infinity for two particles at infinite distance and for long chains. In the present work, these shortcomings have been resolved by adding an energy term for transporting charges through individual atoms. A Gaussian distribution is adopted for the atomic charge distributions, resulting in a damping of the electrostatic interactions at short distances. Assuming that an interatomic exchange term may be described as the overlap between two electronic charge distributions, the EEM has also been extended by a short-range exchange term. The result is a molecular mechanics model where the difference of charge transfer in insulating and metallic systems is modeled regarding the difference in bond length between different types of system. For example, the model is capable of modeling charge transfer in both alkanes and alkenes with alternating double bonds with the same set of carbon parameters only relying on the difference in bond length between carbon sigma- and pi-bonds. Analytical results have been obtained for the polarizability of a long linear chain. These results show that the model is capable of describing the polarizability scaling both linearly and nonlinearly with the size of the system. Similarly, a linear chain with an end atom with a high electronegativity has been analyzed analytically. The dipole moment of this model system can either be independent of the length or increase linearly with the length of the chain. In addition, the model has been parametrized for alkane and alkene chains with data from density functional theory calculations, where the polarizability behaves differently with the chain length. For the molecular dipole moment, the same two systems have been studied with an aldehyde end group. Both the molecular polarizability and the dipole moment are well described as a function of the chain length for both alkane and alkene chains demonstrating the power of the presented model.