Long-range Electrostatic Interactions Research Articles

Charge compensation at the interface between the polar NaCl(111) surface and a NaCl aqueous solution.

Periodic supercell models of electric double layers formed at the interface between a charged surface and an electrolyte are subject to serious finite size errors and require certain adjustments in the treatment of the long-range electrostatic interactions. In a previous publication Zhang and Sprik [Phys. Rev. B 94, 245309 (2016)], we have shown how this can be achieved using finite field methods. The test system was the familiar simple point charge model of a NaCl aqueous solution confined between two oppositely charged walls. Here this method is extended to the interface between the (111) polar surface of a NaCl crystal and a high concentration NaCl aqueous solution. The crystal is kept completely rigid and the compensating charge screening the polarization can only be provided by the electrolyte. We verify that the excess electrolyte ionic charge at the interface conforms to the Tasker 1/2 rule for compensating charge in the theory of polar rock salt (111) surfaces. The interface can be viewed as an electric double layer with a net charge. We define a generalized Helmholtz capacitance CH which can be computed by varying the applied electric field. We find CH=8.23 μF cm-2, which should be compared to the 4.23 μF cm-2 for the (100) non-polar surface of the same NaCl crystal. This is rationalized by the observation that compensating ions shed their first solvation shell adsorbing as contact ions pairs on the polar surface.

Role of Electrostatics in Protein-RNA Binding: The Global vs the Local Energy Landscape.

U1A protein-stem loop 2 RNA association is a basic step in the assembly of the spliceosomal U1 small nuclear ribonucleoprotein. Long-range electrostatic interactions due to the positive charge of U1A are thought to provide high binding affinity for the negatively charged RNA. Short range interactions, such as hydrogen bonds and contacts between RNA bases and protein side chains, favor a specific binding site. Here, we propose that electrostatic interactions are as important as local contacts in biasing the protein-RNA energy landscape toward a specific binding site. We show by using molecular dynamics simulations that deletion of two long-range electrostatic interactions (K22Q and K50Q) leads to mutant-specific alternative RNA bound states. One of these states preserves short-range interactions with aromatic residues in the original binding site, while the other one does not. We test the computational prediction with experimental temperature-jump kinetics using a tryptophan probe in the U1A-RNA binding site. The two mutants show the distinct predicted kinetic behaviors. Thus, the stem loop 2 RNA has multiple binding sites on a rough RNA-protein binding landscape. We speculate that the rough protein-RNA binding landscape, when biased to different local minima by electrostatics, could be one way that protein-RNA interactions evolve toward new binding sites and novel function.

Deciphering the promiscuous interactions between intrinsically disordered transactivation domains and the KIX domain.

The kinase-inducible domain interacting (KIX) domain of the transcriptional coactivator CBP protein carries 2 isolated binding sites (designated as the c-Myb site and the MLL site) and is capable of binding numerous intrinsically disordered transactivation domains (TADs), including c-Myb and pKID via the c-Myb site, and MLL, E2A and c-Jun via the MLL site. In this study we compared the kinetics for binding of various disordered TADs to the KIX domain via computational biophysical analyses. We found that the binding rates are heavily affected by long-range electrostatic interactions. The basal rate constants for forming the encounter complexes are similar for different KIX binding peptides, favorable electrostatic interactions between the MLL site and the peptides result in greater association rates when peptides bind to the MLL site. FOXO3a and p53 TAD each contains 2 copies of KIX binding motif and each motif interacts with both the MLL site and the c-Myb site. Our kinetics studies suggest that binding of FOXO3a or p53 TAD to the KIX domain is via a sequential mechanism, where one KIX binding motif binds to the MLL site first and then the other KIX binding motif binds to the c-Myb site. Considering the promiscuous interactions between FOXO3a and KIX, and p53 TAD and KIX, electrostatic steering simplifies the binding mechanism. This study highlights the importance of long-range electrostatic interactions in molecular recognition process involving multi-motif intrinsically disordered proteins and promiscuous interactions.

Detection of long-range electrostatic interactions between charged molecules by means of fluorescence correlation spectroscopy.

In the present paper, an experimental feasibility study on the detection of long-range intermolecular interactions through three-dimensional molecular diffusion in solution is performed. This follows recent theoretical and numerical analyses reporting that long-range electrodynamic forces between biomolecules could be identified through deviations from Brownian diffusion. The suggested experimental technique was fluorescence correlation spectroscopy (FCS). By considering two oppositely charged molecular species in aqueous solution, namely, lysozymes and fluorescent dye molecules (Alexa488), the diffusion coefficient of the dyes has been measured for different values of the concentration of lysozyme, that is, for different average distances between the oppositely charged molecules. For our model, long-range interactions are of electrostatic origin, suggesting that their action radius can be varied by changing the ionic strength of the solution. The experimental outcomes clearly prove the detectability of long-range intermolecular interactions by means of the FCS technique. Molecular dynamics simulations provide a clear and unambiguous interpretation of the experimental results.

Toward chemical accuracy in the description of ion-water interactions through many-body representations. Alkali-water dimer potential energy surfaces.

This study presents the extension of the MB-nrg (Many-Body energy) theoretical/computational framework of transferable potential energy functions (PEFs) for molecular simulations of alkali metal ion-water systems. The MB-nrg PEFs are built upon the many-body expansion of the total energy and include the explicit treatment of one-body, two-body, and three-body interactions, with all higher-order contributions described by classical induction. This study focuses on the MB-nrg two-body terms describing the full-dimensional potential energy surfaces of the M+(H2O) dimers, where M+ = Li+, Na+, K+, Rb+, and Cs+. The MB-nrg PEFs are derived entirely from "first principles" calculations carried out at the explicitly correlated coupled-cluster level including single, double, and perturbative triple excitations [CCSD(T)-F12b] for Li+ and Na+ and at the CCSD(T) level for K+, Rb+, and Cs+. The accuracy of the MB-nrg PEFs is systematically assessed through an extensive analysis of interaction energies, structures, and harmonic frequencies for all five M+(H2O) dimers. In all cases, the MB-nrg PEFs are shown to be superior to both polarizable force fields and ab initio models based on density functional theory. As previously demonstrated for halide-water dimers, the MB-nrg PEFs achieve higher accuracy by correctly describing short-range quantum-mechanical effects associated with electron density overlap as well as long-range electrostatic many-body interactions.

Nanoscale Investigation of Generation 1 PAMAM Dendrimers Interaction with a Protein Nanopore

Herein, we describe at uni-molecular level the interactions between poly(amidoamine) (PAMAM) dendrimers of generation 1 and the α-hemolysin protein nanopore, at acidic and neutral pH, and ionic strengths of 0.5 M and 1 M KCl, via single-molecule electrical recordings. The results indicate that kinetics of dendrimer-α-hemolysin reversible interactions is faster at neutral as compared to acidic pH, and we propose as a putative explanation the fine interplay among conformational and rigidity changes on the dendrimer structure, and the ionization state of the dendrimer and the α-hemolysin. From the analysis of the dendrimer’s residence time inside the nanopore, we posit that the pH- and salt-dependent, long-range electrostatic interactions experienced by the dendrimer inside the ion-selective α-hemolysin, induce a non-Stokesian diffusive behavior of the analyte inside the nanopore. We also show that the ability of dendrimer molecules to adapt their structure to nanoscopic spaces, and control the flow of matter through the α-hemolysin nanopore, depends non-trivially on the pH- and salt-induced conformational changes of the dendrimer.

A New Mixed All-Atom/Coarse-Grained Model: Application to Melittin Aggregation in Aqueous Solution.

We introduce a new mixed resolution, all-atom/coarse-grained approach (AACG), for modeling peptides in aqueous solution and apply it to characterizing the aggregation of melittin. All of the atoms in peptidic components are represented, while a single site is used for each water molecule. With the full flexibility of the peptide retained, our AACG method achieves speedups by a factor of 3–4 for CPU time reduction and another factor of roughly 7 for diffusion. An Ewald treatment permits the inclusion of long-range electrostatic interactions. These characteristics fit well with the requirements for studying peptide association and aggregation, where the system sizes and time scales require considerable computational resources with all-atom models. In particular, AACG is well suited for biologics since changes in peptide shape and long-range electrostatics may play an important role. The application of AACG to melittin, a 26-residue peptide with a well-known propensity to aggregate in solution, serves as an initial demonstration of this technology for studying peptide aggregation. We observed the formation of melittin aggregates during our simulations and characterized the time-evolution of aggregate size distribution, buried surface areas, and residue contacts. Key interactions including π-cation and π-stacking involving TRP19 were also examined. Our AACG simulations demonstrated a clear salt effect and a moderate temperature effect on aggregation and support the molten globule model of melittin aggregates. As a showcase, this work illustrates the useful role for AACG in investigations of peptide aggregation and its potential to guide formulation and design of biologics.

Combining frozen-density embedding with the conductor-like screening model using Lagrangian techniques for response properties.

We present the explicit derivation of an approach to the multiscale description of molecules in complex environments that combines frozen-density embedding (FDE) with continuum solvation models, in particular the conductor-like screening model (COSMO). FDE provides an explicit atomistic description of molecule-environment interactions at reduced computational cost, while the outer continuum layer accounts for the effect of long-range isotropic electrostatic interactions. Our treatment is based on a variational Lagrangian framework, enabling rigorous derivations of ground- and excited-state response properties. As an example of the flexibility of the theoretical framework, we derive and discuss FDE + COSMO analytical molecular gradients for excited states within the Tamm-Dancoff approximation (TDA) and for ground states within second-order Møller-Plesset perturbation theory (MP2) and a second-order approximate coupled cluster with singles and doubles (CC2). It is shown how this method can be used to describe vertical electronic excitation (VEE) energies and Stokes shifts for uracil in water and carbostyril in dimethyl sulfoxide (DMSO), respectively. In addition, VEEs for some simplified protein models are computed, illustrating the performance of this method when applied to larger systems. The interaction terms between the FDE subsystem densities and the continuum can influence excitation energies up to 0.3 eV and, thus, cannot be neglected for general applications. We find that the net influence of the continuum in presence of the first FDE shell on the excitation energy amounts to about 0.05 eV for the cases investigated. The present work is an important step toward rigorously derived ab initio multilayer and multiscale modeling approaches. © 2017 Wiley Periodicals, Inc.

Ordering and Dynamics for the Formation of Two-Dimensional Molecular Crystals on Black Phosphorene

Phosphorene, a monolayer of elemental phosphorus, has emerged as the most important 2D atomic crystal post graphene. As a consequence of its suitable direct band gap, it has high carrier mobility, thereby making it suitable for potential applications in nanoelectronic devices. In this article, we have performed molecular dynamics simulations for the first time to unveil the potential of black phosphorene (BP) as a platform to fabricate two-dimensional crystals through van der Waals epitaxy of organic molecules. Eight molecules (benzene, 1,3,5-trimethylbenzene, 1,3,5-trifluorobenzene, 1,3,5-trihydroxybenzene, 7,7′,8,8′-tetracyanoquinodimethane (TCNQ), pentacene, coronene, and fullerene) of various size, shape, and polarity are shown to form highly periodic 2D monolayers on BP. We show that self-assembly proceeds through a hierarchy of steps, beginning with adsorption on the surface followed by interplay between the kinetic energy and intermolecular interactions. Crystalline 2D layers are obtained through the minimization of van der Waals interactions for nonpolar and electrostatic interactions for polar molecules. Long-range attractive electrostatic interactions contribute significantly toward the stabilities of the superlattices, thereby making the melting points of polar crystals much higher compared to their nonpolar analogs (Tm (polar) > Tm (nonpolar)). Interestingly, although nonpolar fullerene is shown to form highly stable superlattices due to extensive intermolecular van der Waals interactions, the underlying BP monolayer helps to orient the molecules and, hence, indirectly assists molecular ordering. It also contributes to the overall stability of the molecular crystals through interfacial adsorption. A comparison of self-assembly patterns over BP vis-a-vis graphene and hexagonal boron nitride shows that the nonaromatic nature and undulation of black phosphorene do not affect the intrinsic stability of the assemblies. These findings offer new avenues for generating new molecular doped crystalline superlattices for semiconductor industry.

Electrostatic Origin of the Red Solvatochromic Shift of DFHBDI in RNA Spinach.

Interactions with the environment tune the spectral properties of biological chromophores, e.g., fluorescent proteins. Understanding the relative contribution of the various types of noncovalent interactions in the spectral shifts can provide rational design principles toward developing new fluorescent probes. In this work, we investigate the origin of the red shift in the absorption spectra of the difluoro hydroxybenzylidene dimethyl imidazolinone (DFHBDI) chromophore in RNA spinach as compared to the aqueous solution. We systematically decompose the effects of various components of interactions, namely, stacking, hydrogen bonding, and long-range electrostatics, in order to elucidate the relative role of these interactions in the observed spectral behavior. We find that the absorption peak of DFHBDI is red-shifted by ∼0.35 eV in RNA relative to the aqueous solution. Earlier proposals from Huang and co-workers have implicated the stacking interactions between DFHBDI and nucleic acid bases to be the driving force behind the observed red shift. In contrast, our findings reveal that the long-range electrostatic interactions between DFHBDI and negatively charged RNA make the most significant contribution. Moreover, we notice that the opposing electrostatic fields due to the RNA backbone and the polarized water molecules around the RNA give rise to the resultant red shift. Our results emphasize the effect of strong heterogeneity in the various environmental factors that might be competing with each other.

Accurate size determination of polystyrene latex nanoparticles in aqueous media using a particle tracking analysis method

The size determination of nanoparticles in aqueous solution is important in bio- and nanotechnology. The sizes of particles in aqueous solution are determined by particle tracking analysis (PTA) with the Stokes–Einstein assumption from diffusion coefficient measurements. However, the interactions between particles affect the particle size determined by PTA. In this work, we observed particle/electrolyte concentration-dependent change in the apparent size of polystyrene (PS)-latex particles. Based on the Derjaguin, Landau, Verwey, and Overbeek theory, the change in the observed size of PS-latex particles estimated by PTA in various aqueous electrolyte solutions is caused by the long-range electrostatic interactions between the particles. Accurate size determination of PS-latex particles in aqueous solution was performed after reduction of the electrostatic interaction between particles by varying both the concentrations of the PS-latex particles and the electrolyte in aqueous solution. Zeta potential measurement and DLVO assessment gave a good indicator to reduce the electrostatic interactions between the particles for determination of accurate size of particle by PTA method. The proposed approach makes it possible to compare the determined sizes of various nanoparticles in the liquid-phase by PTA to those obtained by differential mobility analysers in the gas-phase. Although the determined sizes of nanoparticles using different sizing methods are difficult to compare, our investigation would allow researchers to achieve some degree of concordance on the size of nanoparticles using different sizing methods.

Au Nanoparticles in Lipid Bilayers: A Comparison between Atomistic and Coarse-Grained Models

The computational study of the interaction between charged, ligand-protected metal nanoparticles and model lipid membranes has been recently addressed at both atomistic and coarse-grained levels. Here we compare the performance of three versions of the coarse-grained Martini force field at describing the nanoparticle–membrane interaction. The three coarse-grained models differ in terms of treatment of long-range electrostatic interactions and water polarizability. The NP–membrane interaction consists of the transition from a metastable NP–membrane complex, in which the NP is only partially embedded in the membrane, to a configuration in which the NP is anchored to both membrane leaflets. All three of the coarse-grained models provide a description of the metastable NP–membrane complex that is consistent with that obtained using an atomistic force field. As for the anchoring transition, the polarizable-water Martini correctly describes the molecular mechanisms and the energetics of the transition. The standard version of the Martini model, instead, underestimates the free-energy barriers for anchoring and does not completely capture the membrane deformations involved in the transition process.

A Unified Theory for the Blue- and Red-Shifting Phenomena in Hydrogen and Halogen Bonds.

Typical hydrogen and halogen bonds exhibit red-shifts of their vibrational frequencies upon the formation of hydrogen and halogen bonding complexes (denoted as D···Y-A, Y = H and X). The finding of blue-shifts in certain complexes is of significant interest, which has led to numerous studies of the origins of the phenomenon. Because charge transfer mixing (i.e., hyperconjugation in bonding systems) has been regarded as one of the key forces, it would be illuminating to compare the structures and vibrational frequencies in bonding complexes with the charge transfer effect "turned on" and "turned off". Turning off the charge transfer mixing can be achieved by employing the block-localized wave function (BLW) method, which is an ab initio valence bond (VB) method. Further, with the BLW method, the overall stability gained in the formation of a complex can be analyzed in terms of a few physically meaningful terms. Thus, the BLW method provides a unified and physically lucid way to explore the nature of red- and blue-shifting phenomena in both hydrogen and halogen bonding complexes. In this study, a direct correlation between the total stability and the variation of the Y-A bond length is established based on our BLW computations, and the consistent roles of all energy components are clarified. The n(D) → σ*(Y-A) electron transfer stretches the Y-A bond, while the polarization due to the approach of interacting moieties reduces the HOMO-LUMO gap and results in a stronger orbital mixing within the YA monomer. As a consequence, both the charge transfer and polarization stabilize bonding systems with the Y-A bond stretched and red-shift the vibrational frequency of the Y-A bond. Notably, the energy of the frozen wave function is the only energy component which prefers the shrinking of the Y-A bond and thus is responsible for the associated blue-shifting. The total variations of the Y-A bond length and the corresponding stretching vibrational frequency are thus determined by the competition between the frozen-energy term and the sum of polarization and charge transfer energy terms. Because the frozen energy is composed of electrostatic and Pauli exchange interactions and frequency shifting is a long-range phenomenon, we conclude that long-range electrostatic interaction is the driving force behind the frozen energy term.

Stability of Polar Vortex Lattice in Ferroelectric Superlattices.

A novel mesoscale state comprising of an ordered polar vortex lattice has been demonstrated in ferroelectric superlattices of PbTiO3/SrTiO3. Here, we employ phase-field simulations, analytical theory, and experimental observations to evaluate thermodynamic conditions and geometric length scales that are critical for the formation of such exotic vortex states. We show that the stability of these vortex lattices involves an intimate competition between long-range electrostatic, long-range elastic, and short-range polarization gradient-related interactions leading to both an upper and a lower bound to the length scale at which these states can be observed. We found that the critical length is related to the intrinsic domain wall width, which could serve as a simple intuitive design rule for the discovery of novel ultrafine topological structures in ferroic systems.

A refined polarizable water model for the coarse-grained MARTINI force field with long-range electrostatic interactions.

We present a refined version of the polarizable Martini water model - coined refPOL - designed specifically for the use with long-range electrostatics. The refPOL model improves the agreement with the experimentally measured dielectric constant and the mass density of water at room temperature compared to the original polarizable Martini water force field when particle mesh Ewald electrostatics are employed. Our study reveals that the model remains applicable with various commonly used settings for the non-bonded interactions, including reaction field electrostatics. The oil/water partitioning behavior of uncharged Martini bead types is thoroughly investigated: Lennard-Jones interactions between the refPOL model and the remaining Martini beads are adjusted to reproduce the hydration free energies obtained with the original polarizable water model, while free energies of solvation in apolar media remain unchanged. The cross-interactions with charged bead types are parameterized to agree with the experimentally observed area per lipid of a fully solvated dipalmitoylphosphatidylcholine bilayer. We additionally verify the model by analyzing the potentials of mean force between different sample pairs in refPOL water and comparing the results to reference data obtained using the original polarizable Martini water model as well as fully atomistic simulations. Based on the results, we suggest to replace the original polarizable Martini water model with the new refPOL model for future applications.

The Competition between Electrostatic-Steering and Conformational Dynamics in the Diffusion-Limited Association of Calcineurin and Calmodulin

Calcineurin (CaN) is a serine/threonine phosphatase that regulates a variety of physiological and pathophysiological processes in most mammalian tissue. In its inactive state, the CaN catalytic domain is inhibited by the auto-inhibitory domain (AID); the active state is obtained upon CaM binding to the CaN regulatory domain (RD). It has been established that the RD is highly disordered when inhibiting CaN, yet it undergoes a disorder-to-order transition upon binding calmodulin (CaM) to activate the phosphatase. Given that the RD is richly populated with polar and charged amino acids, arguably electrostatic interactions may influence the rate of CaM association. However, it is likely that properties of the RD conformational ensemble, such as its ‘effective volume’ and accessibility of its CaM binding motif, influence the CaM/CaN association rate. In present study, we investigated via computational modeling the extent to which electrostatics and structural disorder co-facilitate or hinder CaM/CaN binding kinetics. We examined several peptides containing the CaM binding motif, for which lengths and amino acid charge distributions were varied, to isolate the contributions of electrostatics versus conformational diversity to predicted, diffusion-limited association rates. These rates were predicted using Replica Exchange Molecular Dynamics (REMD) and Brownian Dynamics (BD) simulations. Our results indicate that association rates vary as a function of increasing CaN RD length (beyond the required CaM recognition sequence), thus indicating that RD conformational ensemble properties influence CaM binding. Second, we found that increasing the solvent ionic strength generally depressed CaM/CaN association rates, owing to the attenuation of long-range electrostatic interactions that would normally accelerate protein-protein association. Finally, CaN peptides with positively-charged amino acids substituted at native negatively-charged sites had complex effects on the predicted association rate, owing to ‘off-target’ interactions that in some cases competed with the intending binding site. Our findings detail the interplay between conformational diversity and electrostatically-driven protein-protein association involving CaN, which are likely to extend to wide-ranging processes regulated by intrinsically-disordered proteins.

Discovery and Mechanistic Characterization of Selective Inhibitors of H2S-producing Enzyme: 3-Mercaptopyruvate Sulfurtransferase (3MST) Targeting Active-site Cysteine Persulfide

Very recent studies indicate that sulfur atoms with oxidation state 0 or −1, called sulfane sulfurs, are the actual mediators of some physiological processes previously considered to be regulated by hydrogen sulfide (H2S). 3-Mercaptopyruvate sulfurtransferase (3MST), one of three H2S-producing enzymes, was also recently shown to produce sulfane sulfur (H2Sn). Here, we report the discovery of several potent 3MST inhibitors by means of high-throughput screening (HTS) of a large chemical library (174,118 compounds) with our H2S-selective fluorescent probe, HSip-1. Most of the identified inhibitors had similar aromatic ring-carbonyl-S-pyrimidone structures. Among them, compound 3 showed very high selectivity for 3MST over other H2S/sulfane sulfur-producing enzymes and rhodanese. The X-ray crystal structures of 3MST complexes with two of the inhibitors revealed that their target is a persulfurated cysteine residue located in the active site of 3MST. Precise theoretical calculations indicated the presence of a strong long-range electrostatic interaction between the persulfur anion of the persulfurated cysteine residue and the positively charged carbonyl carbon of the pyrimidone moiety of the inhibitor. Our results also provide the experimental support for the idea that the 3MST-catalyzed reaction with 3-mercaptopyruvate proceeds via a ping-pong mechanism.

DNA packaging in viral capsids with peptide arms.

Strong chain rigidity and electrostatic self-repulsion of packed double-stranded DNA in viruses require a molecular motor to pull the DNA into the capsid. However, what is the role of electrostatic interactions between different charged components in the packaging process? Though various theories and computer simulation models were developed for the understanding of viral assembly and packaging dynamics of the genome, long-range electrostatic interactions and capsid structure have typically been neglected or oversimplified. By means of molecular dynamics simulations, we explore the effects of electrostatic interactions on the packaging dynamics of DNA based on a coarse-grained DNA and capsid model by explicitly including peptide arms (PAs), linked to the inner surface of the capsid, and counterions. Our results indicate that the electrostatic interactions between PAs, DNA, and counterions have a significant influence on the packaging dynamics. We also find that the packed DNA conformations are largely affected by the structure of the PA layer, but the packaging rate is insensitive to the layer structure.

The effect of electrostatic boundaries in molecular simulations: symmetry matters.

Artifacts arise when the long-ranged electrostatic interaction is inappropriately treated in molecular simulations of electrolytes. When the usual Ewald3D sum method with the tinfoil boundary condition (e3dtf) is used for simulations of an interfacial liquid under an external electric field, a straightforward analysis of the liquid structure often suggests unphysical dielectric properties as a consequence of the inaccurate treatment of the electrostatics. In order to understand the underlying mechanism that leads to this apparent violation of thermodynamics, we now derive a new equation in the weak-field limit that, in a mean field view, accounts for the average effect arising from the difference between e3dtf and the sophisticated Ewald2D sum method (e2d). Numerical simulations of a water system in slab geometry confirm the validity of the weak-field limit equation for a series of parameter setup associated with e3dtf. Moreover, a similar procedure applied to a spherically confined water system suggests that corrections to the seemingly inappropriate treatment of the electrostatics in fact vanish. This cancellation of the boundary effect due to symmetry immediately sheds light on the long-lasting problem of the validity of the ad hoc application of e3dtf for bulk systems. In total, we argue that artifacts arising from e3dtf are often predictable and analytical corrections to the straightforward analysis might be applied to reveal consistent thermodynamic properties in liquid simulations.

Perturbation Approach for Computing Infrared Spectra of the Local Mode of Probe Molecules.

Linear and two-dimensional infrared (IR) spectroscopy of site-specific probe molecules provides an opportunity to gain a molecular-level understanding of the local hydrogen-bonding network, conformational dynamics, and long-range electrostatic interactions in condensed-phase and biological systems. A challenge in computation is to determine the time-dependent vibrational frequencies that incorporate explicitly both nuclear quantum effects of vibrational motions and an electronic structural representation of the potential energy surface. In this paper, a nuclear quantum vibrational perturbation (QVP) method is described for efficiently determining the instantaneous vibrational frequency of a chromophore in molecular dynamics simulations. Computational efficiency is achieved through the use of (a) discrete variable representation of the vibrational wave functions, (b) a perturbation theory to evaluate the vibrational energy shifts due to solvent dynamic fluctuations, and (c) a combined QM/MM potential for the systems. It was found that first-order perturbation is sufficiently accurate, enabling time-dependent vibrational frequencies to be obtained on the fly in molecular dynamics. The QVP method is illustrated in the mode-specific linear and 2D-IR spectra of the H-Cl stretching frequency in the HCl-water clusters and the carbonyl stretching vibration of acetone in aqueous solution. To further reduce computational cost, a hybrid strategy was proposed, and it was found that the computed vibrational spectral peak position and line shape are in agreement with experimental results. In addition, it was found that anharmonicity is significant in the H-Cl stretching mode, and hydrogen-bonding interactions further enhance anharmonic effects. The present QVP method complements other computational approaches, including path integral-based molecular dynamics, and represents a major improvement over the electrostatics-based spectroscopic mapping procedures.