Fluctuating Charge Model Research Articles

Parameterization of a Fluctuating Charge Model for Complexes Containing 3d Transition Metals.

Metalloproteins widely exist in biology, playing pivotal roles in diverse life processes. Meanwhile, molecular dynamics (MD) simulations based on classical force fields has emerged as an important tool in scientific research. Partial charges are critical parameters within classical force fields and usually derived from quantum mechanical (QM) calculations. However, QM calculations are often time-consuming and prone to basis set dependence. Alternatively, fluctuating charge (FQ) models offer another avenue for partial charge derivation, which has significant speed advantages and can be used for large-scale screening. Building upon our previous work, which introduced an FQ model for zinc-containing complexes, herein we extend this model to include additional 3d transition metals which are important to the life sciences, namely chromium, manganese, iron, cobalt, and nickel. Employing CM5 charges as target for parametrization, our FQ model accurately reproduces partial charges for 3d metal complexes featuring biologically relevant ligands. Furthermore, by using atomic charges derived by our FQ model, MD simulations have been performed. These charges showed excellent performance in simulating proteomic metal sites housing multiple metal ions, specifically, a metalloprotein containing an iron-sulfur cluster and another containing a dimanganese metal site, showcasing comparable performance to those of RESP charges. We anticipate that our study can accelerate the parametrization of atomic charges for metalloproteins featuring 3d transition metals, thereby facilitating simulations of relevant systems.

Development of a Fluctuating Charge Model for Zinc-Containing Metalloproteins.

Metalloproteins widely exist in biology and play important roles in various processes. To accurately simulate metalloprotein systems, modeling polarization and charge transfer effects is vital. The fluctuating charge (FQ) model can efficiently generate atomic charges and simulate the charge transfer effect; it has been developed for a wide range of applications, but few models have been specifically tailored for metalloproteins. In this study, we present a fluctuating charge model specifically for zinc-containing metalloproteins based on the extended charge equilibration (EQeq) scheme. Our model was parametrized to reproduce CM5 charges instead of RESP/CHELPG charges because the former is less dependent on the conformation or basis set, does not suffer from unphysical charges for buried atoms, and is still being able to well reproduce the molecular dipoles. During our study, we found that adding the Pauling-bond-order-like term (referred to as the "+C term" in a previous study) between the zinc ion and ligating atoms significantly improves the model's performance. Although our model was trained for four-coordinated zinc sites only, our results indicated it can well describe the atomic charges in diverse zinc sites. Morever, our model was used to generate partial charges for the metal sites in three different zinc-containing metalloproteins (with four-, five-, and six-coordinated metal sites, respectively). These charges exhibited performance comparable to that of the RESP charges in molecular dynamics (MD) simulations. Additional tests indicated our model could also well reproduce the CM5 charges when geometric changes were involved. Those results indicate that our model can efficiently calculate the atomic charges for metal sites and well simulate the charge transfer effect, which marks an important step toward developing versatile polarizable force fields for metalloproteins.

Application of a Fluctuating Charge Polarization Model to Large Polyaromatic Hydrocarbons and Graphene Nanoflakes.

We present a polarization model incorporating coupled fluctuating charges and point inducible dipoles that is able to accurately describe the dipole polarizabilities of small hydrocarbons and, for sufficiently large graphene nanoflakes, reproduce the classical image potential of an infinite conducting sheet. When our fluctuating charge model is applied to the hexagonal carbon nanoflake C60000 we attain excellent agreement with the image potential and induced charge distribution of a conducting sheet. With the inclusion of inducible dipole terms, the model predicts an image plane of zim = 1.3334 a0, which falls in line with prior estimates for graphene. We consider the case of two charges placed on opposite sides of C60000 and find that the fluctuating charge model reproduces classical electrostatics once again. By testing opposing and similar signs of the external charges, we conclude that an atomically thin molecule or extended system does not fully screen their interaction.

Development, Validation, and Pilot Application of a Generalized Fluctuating Charge Model for Computational Spectroscopy in Solution.

A general approach enforcing nonperiodic boundary conditions for the computation of spectroscopic properties in solution has been improved including an effective description of charge-transfer contributions and coordination number adjustment for explicit solvent molecules. Both contributions are obtained from a continuous description of intermolecular hydrogen bonds, which has been employed also for an effective clustering of molecular dynamics trajectories. Fine tuning of the model has been performed for several water clusters, and then its efficiency and reliability have been demonstrated by computing the absorption spectra of different creatinine tautomers in aqueous solution.

Introducing the effective polarizable bond (EPB) model in DNA simulations

Electrostatic polarization plays an important role in characterizing non-bonded interactions of molecule dynamics (MD) simulation of biomolecules. In this work, we extended a simple fluctuating-charge model, the effective polarizable bond (EPB) model to DNA simulations. Following the previously proposed EPB method, we re-parametrized the EPB model for DNA systems with the consideration of a realistic local electric-field environment and derived a set of DNA-specific EPB parameters. Two typical B-DNA systems, a pure DNA 1BNA and a ligand-DNA 8BNA systems were simulated respectively in the ff14SB and ff14SBEPB force fields. Results demonstrate that use of the EPB model in the DNA simulations can cause shrinking of DNA double strands along central axis and strengthen base-pair/ligand-DNA hydrogen-bond (HB) interactions. Furthermore, it stabilizes residue-level fluctuations for the 1BNA. However, it is largely different for these two DNA systems to mimic conformation changes of the minor/major grooves by the EPB model. Lastly, it is tested that approximately %9∼%9.9 CPU costs are additionally taken to execute the EPB calculations compared with the conventional ff14SB force field.

Bridging the 12-6-4 Model and the Fluctuating Charge Model.

Metal ions play important roles in various biological systems. Molecular dynamics (MD) using classical force field has become a popular research tool to study biological systems at the atomic level. However, meaningful MD simulations require reliable models and parameters. Previously we showed that the 12-6 Lennard-Jones nonbonded model for ions could not reproduce the experimental hydration free energy (HFE) and ion-oxygen distance (IOD) values simultaneously when ion has a charge of +2 or higher. We discussed that this deficiency arises from the overlook of the ion-induced dipole interaction in the 12-6 model, and this term is proportional to 1/r 4 based on theory. Hence, we developed the 12-6-4 model and showed it could solve this deficiency in a physically meaningful way. However, our previous research also found that the 12-6-4 model overestimated the coordination numbers (CNs) for some highly charged metal ions. And we attributed this artifact to that the current 12-6-4 scheme lacks a correction for the interactions among the first solvation shell water molecules. In the present study, we considered the ion-included dipole interaction by using the 12-6 model with adjusting the atomic charges of the first solvation shell water molecules. This strategy not only considers the ion-induced dipole interaction between ion and the first solvation shell water molecules but also well accounts for the increased repulsion among these water molecules compared to the bulk water molecules. We showed this strategy could well reproduce the experimental HFE and IOD values for Mg2+, Zn2+, Al3+, Fe3+, and In3+ and solve the CN overestimation issue of the 12-6-4 model for Fe3+ and In3+. Moreover, our simulation results showed good agreement with previous ab initio MD simulations. In addition, we derived the physical relationship between the C 4 parameter and induced dipole moment, which agreed well with our simulation results. Finally, we discussed the implications of the present work for simulating metalloproteins. Due to the fluctuating charge model uses a similar concept to the 12-6 model with adjusting atomic charges, we believe the present study builds a bridge between the 12-6-4 model and the fluctuating charge model.

The phase diagram of carbon dioxide from correlation functions and a many-body potential.

The phase stability and equilibria of carbon dioxide are investigated from 125-325K and 1-10 000 atm using extensive molecular dynamics (MD) simulations and the Two-Phase Thermodynamics (2PT) method. We devise a direct approach for calculating phase diagrams, in general, by considering the separate chemical potentials of the isolated phase at specific points on the P-T diagram. The unique ability of 2PT to accurately and efficiently approximate the entropy and Gibbs energy of liquids allows for assignment of phase boundaries from relatively short (∼100ps) MD simulations. We validate our approach by calculating the critical properties of the flexible elementary physical model 2, showing good agreement with previous results. We show, however, that the incorrect description of the short-range Pauli force and the lack of molecular charge polarization lead to deviations from experiments at high pressures. We, thus, develop a many-body, fluctuating charge model for CO2, termed CO2-Fq, from high level quantum mechanics (QM) calculations that accurately capture the condensed phase vibrational properties of the solid (including the Fermi resonance at 1378cm-1) as well as the diffusional properties of the liquid, leading to overall excellent agreement with experiments over the entire phase diagram. This work provides an efficient computational approach for determining phase diagrams of arbitrary systems and underscores the critical role of QM charge reorganization physics in molecular phase stability.

Stochastic PDEs, random fields and exact mean-values

Introducing projector-operator technique and algebra of Terwiel’s cumulants we study stochastic linear partial differential equations with global and local disorder. We present the evolution equation for the mean-value of the field as a series in terms of Terwiel’s cumulant operators. Then, we prove that if we use binary disorder with time exponential-correlated structure, as source of the stochastic perturbation, this series cuts leading to a treatable evolution equation. We apply this approach to find the exact mean-value solution of electromagnetic waves with stochastic absorption of energy in conducting media. This model shows the occurrence of novel time-scale separation phenomena. Local disorder in telegrapher’s equation is also presented. Thus we show that strong disorder leads to anomalous behavior at short and long time regimes. In addition, other physical systems with global disorder are worked out to find exact mean-value solutions: finite-velocity diffusion in the presence of a deterministic force (Smoluchoswki-like process generalizing, in this way, Feynman–Kac’s formula for its numerical solution); Lorentz’ force on a fluctuating charge model (we calculate the diffusion coefficient transverse to the applied magnetic field); and a generalized non-Maxwellian velocity distribution (Ornstein–Uhlenbeck like process showing a noise-induced transition in the stationary distribution).

Mass-zero constrained molecular dynamics for electrode charges in simulations of electrochemical systems.

Classical molecular dynamics simulations have recently become a standard tool for the study of electrochemical systems. State-of-the-art approaches represent the electrodes as perfect conductors, modeling their responses to the charge distribution of electrolytes via the so-called fluctuating charge model. These fluctuating charges are additional degrees of freedom that, in a Born-Oppenheimer spirit, adapt instantaneously to changes in the environment to keep each electrode at a constant potential. Here, we show that this model can be treated in the framework of constrained molecular dynamics, leading to a symplectic and time-reversible algorithm for the evolution of all the degrees of freedom of the system. The computational cost and the accuracy of the new method are similar to current alternative implementations of the model. The advantage lies in the accuracy and long term stability guaranteed by the formal properties of the algorithm and in the possibility to systematically introduce additional kinematic conditions of arbitrary number and form. We illustrate the performance of the constrained dynamics approach by enforcing the electroneutrality of the electrodes in a simple capacitor consisting of two graphite electrodes separated by a slab of liquid water.

A force consistent method for electrostatic energy calculation in fluctuating charge model.

A practical approach to include the polarization effect in a molecular force field is the fluctuating charge method in which atomic charges vary as the configuration of the molecular system changes. However, the use of the Coulomb formula to evaluate energy in a fluctuating charge method is theoretically inconsistent with the forces given by the fluctuating method. In this work, we propose a force-consistent method to correctly calculate electrostatic energies of molecular systems using a fluctuating charge model (Effective Polarizable Bond or EPB). In this protocol, the electrostatic energy is obtained by numerical interaction of the atomic forces along the MD trajectory, rather than using the default Coulomb formula in the EPB model. Test study on the benchmark Barnase-Barstar protein-protein interaction system demonstrates that although the total electrostatic energy of the system shows little deviation due to the averaging effect, specific residue-residue electrostatic interaction energy is affected and the level of the effect depends on the charges of the interacting residues with charged residues showing pronounced differences in calculated energies between using the current protocol and the standard Coulomb formula. It is recommended that the proposed numerical interaction method should be preferred in the calculation of electrostatic energy in fluctuating charge models used in molecular dynamics simulations.

Fully Polarizable QM/Fluctuating Charge Approach to Two-Photon Absorption of Aqueous Solutions.

We present the extension of the quantum/classical polarizable fluctuating charge model to the calculation of single residues of quadratic response functions, as required for the computational modeling of two-photon absorption cross sections. By virtue of a variational formulation of the quantum/classical polarizable coupling, we are able to exploit an atomic orbital-based quasienergy formalism to derive the additional coupling terms in the response equations. Our formalism can be extended to the calculation of arbitrary order response functions and their residues. The approach has been applied to the challenging problem of one- and two-photon spectra of rhodamine 6G (R6G) in aqueous solution. Solvent effects on one- and two-photon spectra of R6G in aqueous solution have been analyzed by considering three different approaches, from a continuum (QM/PCM) to two QM/MM models (nonpolarizable QM/TIP3P and polarizable QM/FQ). Both QM/TIP3P and QM/FQ simulated OPA and TPA spectra show that the inclusion of discrete water solvent molecules is essential to increase the agreement between theory and experiment. QM/FQ has been shown to give the best agreement with experiments.

Higher-Order Extended Lagrangian Born-Oppenheimer Molecular Dynamics for Classical Polarizable Models.

Generalized extended Lagrangian Born-Oppenheimer molecular dynamics (XLBOMD) methods provide a framework for fast iteration-free simulations of models that normally require expensive electronic ground state optimizations prior to the force evaluations at every time step. XLBOMD uses dynamically driven auxiliary degrees of freedom that fluctuate about a variationally optimized ground state of an approximate "shadow" potential which approximates the true reference potential. While the requirements for such shadow potentials are well understood, constructing such potentials in practice has previously been ad hoc, and in this work, we present a systematic development of XLBOMD shadow potentials that match the reference potential to any order. We also introduce a framework for combining friction-like dissipation for the auxiliary degrees of freedom with general-order integration, a combination that was not previously possible. These developments are demonstrated with a simple fluctuating charge model and point induced dipole polarization models.

Structure and Liquid Fragility in Sodium Carbonate.

The relationship between local structure and dynamics is explored for molten sodium carbonate. A flexible fluctuating-charge model, which allows for changes in the shape and charge distribution of the carbonate molecular anion, is developed. The system shows the evolution of highly temperature-dependent complex low-dimensional structures which control the dynamics (and hence the liquid fragility). By varying the molecular anion charge distribution, the key interactions responsible for the formation of these structures can be identified and rationalized. An increase in the mean charge separation within the carbonate ions increases the connectivity of the emerging structures and leads to an increase in the system fragility.

Studies on hydrolysis mechanism of anticancer ruthenium drug ImH[trans-Ru(Im)2Cl4] via ABEEMσπ polarizable force field combined with QM and MD-FEP

We used ABEEMσπ(atom-bond electronegativity equalization method) polarizable force field(ABEEMσπ PFF) method combined with QM and molecular dynamics-free energy perturbation(MD-FEP) methods to investigate the function of water molecules in hydrolysis process of ImH[trans-Ru(Im)2Cl4](ICR). The activation free energies obtained via MD-FEP calculation are in fair agreement with the experimental data. In addition, QM/MM(ABEEM) rationally describes the charge distributions and the electrostatic interaction between molecules. ABEEMσπ fluctuating charge model has the following good characteristics: (1) not only atomic charge regions but also σ, π bond and lone pair charge regions are explicitly represented for a molecule; (2) the region charges are geometry dependent and calculated from time to time in the dynamic simulation without any iterative procedure so that its performance is time-saving compared with the Drude model and induced dipole model.

Electronic absorption spectra of pyridine and nicotine in aqueous solution with a combined molecular dynamics and polarizable QM/MM approach.

The electronic absorption spectra of pyridine and nicotine in aqueous solution have been computed using a multistep approach. The computational protocol consists in studying the solute solvation with accurate molecular dynamics simulations, characterizing the hydrogen bond interactions, and calculating electronic transitions for a series of configurations extracted from the molecular dynamics trajectories with a polarizable QM/MM scheme based on the fluctuating charge model. Molecular dynamics simulations and electronic transition calculations have been performed on both pyridine and nicotine. Furthermore, the contributions of solute vibrational effect on electronic absorption spectra have been taken into account in the so called vertical gradient approximation. © 2016 The Authors. Journal of Computational Chemistry Published by Wiley Periodicals, Inc.

Evaluation of Representations and Response Models for Polarizable Force Fields.

For classical simulations of condensed-phase systems, such as organic liquids and biomolecules, to achieve high accuracy, they will probably need to incorporate an accurate, efficient model of conformation-dependent electronic polarization. Thus, it is of interest to understand what determines the accuracy of a polarizable electrostatics model. This study approaches this problem by breaking polarization models down into two main components: the representation of electronic polarization and the response model used for mapping from an inducing field to the polarization within the chosen representation. Among the most common polarization representations are redistribution of atom-centered charges, such as those used in the fluctuating charge model, and atom-centered point dipoles, such as those used in a number of different polarization models. Each of these representations has been combined with one or more response models. The response model of fluctuating charge, for example, is based on the idea of electronegativity equalization in the context of changing electrostatic potentials (ESPs), whereas point-dipole representations typically use a response model based on point polarizabilities whose induced dipoles are computed based on interaction with other charges and dipoles. Here, we decouple polarization representations from their typical response models to analyze the strengths and weaknesses of various polarization approximations. First, we compare the maximal possible accuracies achievable by the charge redistribution and point-dipole model representations, by testing their ability to replicate quantum mechanical (QM) ESPs around small molecules polarized by external inducing charges. Perhaps not surprisingly, the atom-centered dipole model can yield higher accuracy. Next, we test two of the most commonly used response functions used for the point-dipole representations, self-consistent and direct (or first-order) inducible point polarizabilities, where the polarizabilities are optimized to best fit the full set of polarized QM potentials for each molecule studied. Strikingly, the induced-dipole response model markedly degrades accuracy relative to that obtainable with optimal point dipoles. In fact, the maximal accuracy achievable with this response model is even lower than that afforded by an optimal charge-redistribution representation. This means that, if coupled with a sufficiently accurate response function, the point-charge representation could outperform the standard induced-dipole model. Furthermore, although a key advantage of the point-dipole representation, relative to charge redistribution, is its ability to capture out-of-plane polarization, the inducible dipole response model causes it to be less accurate than optimal charge redistribution for out-of-plane induction of the planar nitrobenzene molecule. Thus, the widely used inducible dipole response function falls short of the full potential accuracy achievable with the point-dipole representation it employs. Additional results reported here bear on the relative accuracy of self-consistent inducible dipoles versus that of the first-order, or direct, approximation and on methods for assigning partial atomic charges for use in conjunction with inducible dipole models. In sum, these results point to the improvement of polarization response models as an important direction for future research aimed at improving the accuracy of molecular simulations.

Polarizable force fields based on physical models and quantum chemical calculations

This perspective gives a brief overview of recent developments of the polarizable force fields (FFs)—a kind of specific FF method that includes polarization effect into conventional molecular mechanics. The commonly adopted polarizable models, that is, the fluctuating charge model, Drude model, and the inducible dipole model are expatiated. Taking advantage of the recent development of computational techniques and fragment‐based low‐scaling quantum mechanics (QM) methods, QM‐based polarizable FFs appeared and particularly aroused great interest in biological systems. Current applications and limitations of several models are discussed. Opportunities and challenges for future development are also addressed. © 2014 Wiley Periodicals, Inc.

An efficient fluctuating charge model for transition metal complexes

A fluctuating charge model for transition metal complexes, based on the Hirshfeld partitioning scheme, spectroscopic energy data from the NIST Atomic Spectroscopy Database and the electronegativity equalization approach, has been developed and parameterized for organic ligands and their high- and low-spin Fe(II) and Fe(III), low-spin Co(III) and Cu(II) complexes, using atom types defined in the Momec force field. Based on large training sets comprising a variety of transition metal complexes, a general parameter set has been developed and independently validated which allows the efficient computation of geometry-dependent charge distributions in the field of transition metal coordination compounds.

Microhydration of alanine in gas phase studied by quantum chemical method and ABEEMσπ/MM fluctuating charge model

A fluctuating charge interaction potential function for alanine-water was constructed in the spirit of newly developed ABEEMσπ/MM(atom-bond electronegativity equalization method at the σπ level fused into molecular mechanics). The properties of gaseous neutral alanine-(H2O)n(n=1–7) clusters were systematically investigated by quantum mechanics(QM) and the constructed ABEEMσπ/MM potential, such as conformations, hydrogen bonds (H-bonds), interaction energies, charge distributions, and so on. The results of ABEEMσπ/MM model are in fair agreement with those of QM and available experimental data. For isolated alanine, compared with those of experimental structure, the average absolute deviations(AAD) of bond length and bond angle are 0.002 nm and 1.4°, respectively. For alanine-water clusters, the AAD of interaction energies and H-bond lengths are only 3.77 kJ/mol and 0.012 nm, respectively, compared to the results of MP2/aug-cc-pVDZ//MP2/6-311+G** method. The ABEEMσπ charges fluctuate with the changing conformation of the system, and can accurately and reasonably reflect the interpolarization between water and alanine. The presented alanine-water potential function may provide a basis for further simulations on related aqueous solutions of biomolecules.

Mg2+/Ca2+ binding to DNA bases: a quantum chemical method and ABEEMσπ/MM fluctuating charge model study

The interactions of Mg2+ and Ca2+ binding to adenine, cytosine, guanine, and thymine at various binding sites were studied by a high-level quantum chemical method and ABEEMσπ/MM fluctuating charge model. The geometries and binding energies of M2+-bases complexes were determined at CCSD(T)/6-311 ++G(2d,2p)//MP2/6-311 ++G(2d,2p) level of theory, with the basis set superposition error corrections for the binding energy calculations. In comparison with the ab initio results, an accurate classical metal cation–base interaction potential was constructed and parameterized in terms of ABEEMσπ/MM model. It is revealed that Mg2+/Ca2+ prefers to bind with bases at the bidentate position (between two nitrogen atoms or oxygen and nitrogen atoms in purine and pyrimidine), where the binding energy is the largest. Moreover, the distance between M2+ and the base increases from Mg2+ to Ca2+, while the binding energy of Mg2+–base is greater than that of Ca2+-base. The ABEEMσπ/MM potential gives reasonable geometries and binding energies compared with the present quantum chemical calculations, and the overall percentage RMSDs are 1.4 and 1.6% for geometries and binding energies, respectively. Furthermore, the transferability of the parameters of the new potential is validated by investigation of Mg2+/Ca2+ binding to tautomer of bases, and results from our potential also show quite good consistency with those of MP2/6-311 ++G(2d,2p)//B3LYP/6-311 ++G(d,p) method, with the overall percentage RMSDs of 2.2 and 4.7% for geometries and binding energies, respectively. This work will serve as a basis for further investigations of the mechanisms of cation effects on the structure and property of nucleic acids.