Nitriles are important solvents not just for bulk reactions but also for interfacial processes such as separations, heterogeneous catalysis, and electrochemistry. Although nitriles have a polar end and a lipophilic end, the cyano group is not hydrophilic enough for these substances to be thought of as prototypical amphiphiles. This picture is now changing, as research is revealing that at a silica surface nitriles can organize into structures that, in many ways, resemble lipid bilayers. This unexpected organization may be a key component of unique interfacial behavior of nitriles that make them the solvents of choice for so many applications.The first hints of this lipid-bilayer-like (LBL) organization of nitriles at silica interfaces came from optical Kerr effect (OKE) experiments on liquid acetonitrile confined in the pores of sol–gel glasses. The orientational dynamics revealed by OKE spectroscopy suggested that the confined liquid is composed of a relatively immobile sublayer of molecules that accept hydrogen bonds from the surface silanol groups and an interdigitated, antiparallel layer that is capable of exchanging into the centers of the pores.This picture of acetonitrile has been borne out by molecular dynamics simulations and vibrational sum-frequency generation (VSFG) experiments. Remarkably, these simulations further indicate that the LBL organization is repeated with increasing disorder at least 20 Å into the liquid from a flat silica surface. Simulations and VSFG and OKE experiments indicate that extending the alkyl chain to an ethyl group leads to the formation of even more tightly packed LBL organization featuring entangled alkyl tails. When the alkyl portion of the molecule is a bulky t-butyl group, packing constraints prevent well-ordered LBL organization of the liquid. In each case, the surface-induced organization of the liquid is reflected in its interfacial dynamics.Acetonitrile/water mixtures are favored solvent systems for separations technologies such as hydrophilic interaction chromatography. Simulations had suggested that although a monolayer of water partitions to the silica surface in such mixtures, acetonitrile tends to associate with this monolayer. VSFG experiments reveal that, even at high water mole fractions, patches of well-ordered acetonitrile bilayers remain at the silica surface. Due to its ability to donate and accept hydrogen bonds, methanol also partitions to a silica surface in acetonitrile/methanol mixtures and can serve to take the place of acetonitrile in the sublayer closest to the surface.These studies reveal that liquid nitriles can exhibit an unexpected wealth of new organizational and dynamic behaviors at silica surfaces, and presumably at the surfaces of other chemically important materials as well. This behavior cannot be predicted from the bulk organization of these liquids. Our new understanding of the interfacial behavior of these liquids will have important implications for optimizing a wide range of chemical processes in nitrile solvents.

Mutations in the gatekeeper residue of kinases have emerged as a key way through which cancer cells develop resistance to treatment. As such, the design of gatekeeper mutation resistant kinase inhibitors is a crucial way forward in increasing the efficacy of a broad range of anticancer drugs. In this work we use atomistic simulations to provide detailed thermodynamic and structural insight into how two inhibitors of cSrc kinase, namely, a commercial drug and type I kinase inhibitor Dasatinib and the type II inhibitor RL45, respectively fail and succeed in being effective against the T338M gatekeeper residue mutation in the kinase binding site. Given the well-known limitations of atomistic simulations in sampling biomolecular systems, we use an enhanced sampling technique called free energy perturbation with replica exchange solute tempering (FEP/REST). Our calculations find that the type I inhibitor Dasatinib binds favorably to the wild type but unfavorably to T338M mutated kinase, while RL45 binds favorably to both. The predicted relative binding free energies are well within 1 kcal/mol accuracy compared to experiments. We find that Dasatinib’s impotency against gatekeeper residue mutations arises from a loss of ligand–kinase hydrogen bonding due to T338M mutation and from steric hindrance due to the presence of an inflexible phenyl ring close to the ligand. On the other hand, in the type II binding RL45, the central phenyl ring has very pronounced flexibility. This leads to the inhibitor overcoming effects of steric clashes on mutation and maintaining an electrostatically favorable “edge-to-face” orientation with a neighboring phenylalanine residue. Our work provides useful insight into the mechanisms of mutation resistant kinase inhibitors and demonstrates the usefulness of enhanced sampling techniques in computational drug design.

In modern-day simulations of many-body systems, much of the computational complexity is shifted to the identification of slowly changing molecular order parameters called collective variables (CVs) or reaction coordinates. A vast array of enhanced-sampling methods are based on the identification and biasing of these low-dimensional order parameters, whose fluctuations are important in driving rare events of interest. Here, we describe a new algorithm for finding optimal low-dimensional CVs for use in enhanced-sampling biasing methods like umbrella sampling, metadynamics, and related methods, when limited prior static and dynamic information is known about the system, and a much larger set of candidate CVs is specified. The algorithm involves estimating the best combination of these candidate CVs, as quantified by a maximum path entropy estimate of the spectral gap for dynamics viewed as a function of that CV. The algorithm is called spectral gap optimization of order parameters (SGOOP). Through multiple practical examples, we show how this postprocessing procedure can lead to optimization of CV and several orders of magnitude improvement in the convergence of the free energy calculated through metadynamics, essentially giving the ability to extract useful information even from unsuccessful metadynamics runs.

A key factor influencing a drug’s efficacy is its residence time in the binding pocket of the host protein. Using atomistic computer simulation to predict this residence time and the associated dissociation process is a desirable but extremely difficult task due to the long timescales involved. This gets further complicated by the presence of biophysical factors such as steric and solvation effects. In this work, we perform molecular dynamics (MD) simulations of the unbinding of a popular prototypical hydrophobic cavity–ligand system using a metadynamics-based approach that allows direct assessment of kinetic pathways and parameters. When constrained to move in an axial manner, the unbinding time is found to be on the order of 4,000 s. In accordance with previous studies, we find that the cavity must pass through a region of sharp wetting transition manifested by sudden and high fluctuations in solvent density. When we remove the steric constraints on ligand, the unbinding happens predominantly by an alternate pathway, where the unbinding becomes 20 times faster, and the sharp wetting transition instead becomes continuous. We validate the unbinding timescales from metadynamics through a Poisson analysis, and by comparison through detailed balance to binding timescale estimates from unbiased MD. This work demonstrates that enhanced sampling can be used to perform explicit solvent MD studies at timescales previously unattainable, to our knowledge, obtaining direct and reliable pictures of the underlying physiochemical factors including free energies and rate constants.

It is currently the consensus belief that protective osmolytes such as trimethylamine N-oxide (TMAO) favor protein folding by being excluded from the vicinity of a protein, whereas denaturing osmolytes such as urea lead to protein unfolding by strongly binding to the surface. Despite there being consensus on how TMAO and urea affect proteins as a whole, very little is known as to their effects on the individual mechanisms responsible for protein structure formation, especially hydrophobic association. In the present study, we use single-molecule atomic force microscopy and molecular dynamics simulations to investigate the effects of TMAO and urea on the unfolding of the hydrophobic homopolymer polystyrene. Incorporated with interfacial energy measurements, our results show that TMAO and urea act on polystyrene as a protectant and a denaturant, respectively, while complying with Tanford–Wyman preferential binding theory. We provide a molecular explanation suggesting that TMAO molecules have a greater thermodynamic binding affinity with the collapsed conformation of polystyrene than with the extended conformation, while the reverse is true for urea molecules. Results presented here from both experiment and simulation are in line with earlier predictions on a model Lennard–Jones polymer while also demonstrating the distinction in the mechanism of osmolyte action between protein and hydrophobic polymer. This marks, to our knowledge, the first experimental observation of TMAO-induced hydrophobic collapse in a ternary aqueous system.

Computer simulations are used to determine the free energy landscape for the binding of the anticancer drug Dasatinib to its src kinase receptor and show that before settling into a free energy basin the ligand must surmount a free energy barrier. An analysis based on using both the ligand-pocket separation and the pocket-water occupancy as reaction coordinates shows that the free energy barrier is a result of the free energy cost for almost complete desolvation of the binding pocket. The simulations further show that the barrier is not a result of the reorganization free energy of the binding pocket. Although a continuum solvent model gives the location of free energy minima, it is not able to reproduce the intermediate free energy barrier. Finally, it is shown that a kinetic model for the on rate constant in which the ligand diffuses up to a doorway state and then surmounts the desolvation free energy barrier is consistent with published microsecond time-scale simulations of the ligand binding kinetics for this system [Shaw, D. E. et al. J. Am. Chem. Soc. 2011, 133, 9181−9183].

Single-molecule force spectroscopies are remarkable tools for studying protein folding and unfolding, but force unfolding explores protein configurations that are potentially very different from the ones traditionally explored in chemical or thermal denaturation. Understanding these differences is crucial because such configurations serve as starting points of folding studies, and thus can affect both the folding mechanism and the kinetics. Here we provide a detailed comparison of both chemically induced and force-induced unfolded state ensembles of ubiquitin based on extensive, all-atom simulations of the protein either extended by force or denatured by urea. As expected, the respective unfolded states are very different on a macromolecular scale, being fully extended under force with no contacts and partially extended in urea with many nonnative contacts. The amount of residual secondary structure also differs: A significant population of α-helices is found in chemically denatured configurations but such helices are absent under force, except at the lowest applied force of 30 pN where short helices form transiently. We see that typical-size helices are unstable above this force, and β-sheets cannot form. More surprisingly, we observe striking differences in the backbone dihedral angle distributions for the protein unfolded under force and the one unfolded by denaturant. A simple model based on the dialanine peptide is shown to not only provide an explanation for these striking differences but also illustrates how the force dependence of the protein dihedral angle distributions give rise to the worm-like chain behavior of the chain upon force.

Heterogeneous ice growth exhibits a maximum in freezing rate arising from the competition between kinetics and the thermodynamic driving force between the solid and liquid states. Here, we use molecular dynamics simulations to elucidate the atomistic details of this competition, focusing on water properties in the interfacial region along the secondary prismatic direction. The crystal growth velocity is maximized when the efficiency of converting interfacial water molecules to ice, collectively known as the attachment kinetics, is greatest. We find water molecules that contact the intermediate ice layer in concave regions along the atomistically roughened surface are more likely to freeze directly. An increased roughening of the solid surface at large undercoolings consequently plays an important limiting role in the rate of ice growth, as water molecules are unable to integrate into increasingly deeper surface pockets. These results provide insight into the molecular mechanisms for self-assembly of solid phases that are important in many biological and atmospheric processes.

Longstanding mechanistic questions about the role of protecting osmolyte trimethylamine N-oxide (TMAO) that favors protein folding and the denaturing osmolyte urea are addressed by studying their effects on the folding of uncharged polymer chains. Using atomistic molecular dynamics simulations, we show that 1 M TMAO and 7 M urea solutions act dramatically differently on these model polymer chains. Their behaviors are sensitive to the strength of the attractive dispersion interactions of the chain with its environment: when these dispersion interactions are sufficiently strong, TMAO suppresses the formation of extended conformations of the hydrophobic polymer as compared to water while urea promotes the formation of extended conformations. Similar trends are observed experimentally for real protein systems. Quite surprisingly, we find that both protecting and denaturing osmolytes strongly interact with the polymer, seemingly in contrast with existing explanations of the osmolyte effect on proteins. We show that what really matters for a protective osmolyte is its effective depletion as the polymer conformation changes, which leads to a negative change in the preferential binding coefficient. For TMAO, there is a much more favorable free energy of insertion of a single osmolyte near collapsed conformations of the polymer than near extended conformations. By contrast, urea is preferentially stabilized next to the extended conformation and thus has a denaturing effect.

Force spectroscopies have emerged as a powerful and unprecedented tool to study and manipulate biomolecules directly at a molecular level. Usually, protein and DNA behavior under force is described within the framework of the worm-like chain (WLC) model for polymer elasticity. Although it has been surprisingly successful for the interpretation of experimental data, especially at high forces, the WLC model lacks structural and dynamical molecular details associated with protein relaxation under force that are key to the understanding of how force affects protein flexibility and reactivity. We use molecular dynamics simulations of ubiquitin to provide a deeper understanding of protein relaxation under force. We find that the WLC model successfully describes the simulations of ubiquitin, especially at higher forces, and we show how protein flexibility and persistence length, probed in the force regime of the experiments, are related to how specific classes of backbone dihedral angles respond to applied force. Although the WLC model is an average, backbone model, we show how the protein side chains affect the persistence length. Finally, we find that the diffusion coefficient of the protein’s end-to-end distance is on the order of 108 nm2/s, is position and side-chain dependent, but is independent of the length and independent of the applied force, in contrast with other descriptions.

A model of protein–ligand binding kinetics, in which slow solvent dynamics results from hydrophobic drying transitions, is investigated. Molecular dynamics simulations show that solvent in the receptor pocket can fluctuate between wet and dry states with lifetimes in each state that are long enough for the extraction of a separable potential of mean force and wet-to-dry transitions. We present a diffusive surface hopping model that is represented by a 2D Markovian master equation. One dimension is the standard reaction coordinate, the ligand–pocket separation, and the other is the solvent state in the region between ligand and binding pocket which specifies whether it is wet or dry. In our model, the ligand diffuses on a dynamic free-energy surface which undergoes kinetic transitions between the wet and dry states. The model yields good agreement with results from explicit solvent molecular dynamics simulation and an improved description of the kinetics of hydrophobic assembly. Furthermore, it is consistent with a “non-Markovian Brownian theory” for the ligand–pocket separation coordinate alone.

When two phases of water are at equilibrium, the ratio of hydrogen isotopes in each is slightly altered because of their different phase affinities. This isotopic fractionation process can be utilized to analyze water’s movement in the world’s climate. Here we show that equilibrium fractionation ratios, an entirely quantum mechanical property, also provide a sensitive probe to assess the magnitude of nuclear quantum fluctuations in water. By comparing the predictions of a series of water models, we show that those describing the OH chemical bond as rigid or harmonic greatly overpredict the magnitude of isotope fractionation. Models that account for anharmonicity in this coordinate are shown to provide much more accurate results because of their ability to give partial cancellation between inter- and intramolecular quantum effects. These results give evidence of the existence of competing quantum effects in water and allow us to identify how this cancellation varies across a wide-range of temperatures. In addition, this work demonstrates that simulation can provide accurate predictions and insights into hydrogen fractionation.

The elastic restoring force of tissues must be able to operate over the very wide range of loading rates experienced by living organisms. It is surprising that even the fastest events involving animal muscle tissues do not surpass a few hundred hertz. We propose that this limit is set in part by the elastic dynamics of tethered proteins extending and relaxing under a changing load. Here we study the elastic dynamics of tethered proteins using a fast force spectrometer with sub-millisecond time resolution, combined with Brownian and Molecular Dynamics simulations. We show that the act of tethering a polypeptide to an object, an inseparable part of protein elasticity in vivo and in experimental setups, greatly reduces the attempt frequency with which the protein samples its free energy. Indeed, our data shows that a tethered polypeptide can traverse its free-energy landscape with a surprisingly low effective diffusion coefficient Deff ∼ 1,200 nm2/s. By contrast, our Molecular Dynamics simulations show that diffusion of an isolated protein under force occurs at Deff ∼ 108 nm2/s. This discrepancy is attributed to the drag force caused by the tethering object. From the physiological time scales of tissue elasticity, we calculate that tethered elastic proteins equilibrate in vivo with Deff ∼ 104–106 nm2/s which is two to four orders magnitude smaller than the values measured for untethered proteins in bulk.

We apply a free energy perturbation simulation method, free energy perturbation/replica exchange with solute tempering, to two modifications of protein–ligand complexes that lead to significant conformational changes, the first in the protein and the second in the ligand. The approach is shown to facilitate sampling in these challenging cases where high free energy barriers separate the initial and final conformations and leads to superior convergence of the free energy as demonstrated both by consistency of the results (independence from the starting conformation) and agreement with experimental binding affinity data. The second case, consisting of two neutral thrombin ligands that are taken from a recent medicinal chemistry program for this interesting pharmaceutical target, is of particular significance in that it demonstrates that good results can be obtained for large, complex ligands, as opposed to relatively simple model systems. To achieve quantitative agreement with experiment in the thrombin case, a next generation force field, Optimized Potentials for Liquid Simulations 2.0, is required, which provides superior charges and torsional parameters as compared to earlier alternatives.

Past experience provides a rich source of predictive information about the world that could be used to guide and optimize ongoing perception. However, the neural mechanisms that integrate information coded in long-term memory (LTM) with ongoing perceptual processing remain unknown. Here, we explore how the contents of LTM optimize perception by modulating anticipatory brain states. By using a paradigm that integrates LTM and attentional orienting, we first demonstrate that the contents of LTM sharpen perceptual sensitivity for targets presented at memory-predicted spatial locations. Next, we examine oscillations in EEG to show that memory-guided attention is associated with spatially specific desynchronization of alpha-band activity over visual cortex. Additionally, we use functional MRI to confirm that target-predictive spatial information stored in LTM triggers spatiotopic modulation of preparatory activity in extrastriate visual cortex. Finally, functional MRI results also implicate an integrated cortical network, including the hippocampus and a dorsal frontoparietal circuit, as a likely candidate for organizing preparatory states in visual cortex according to the contents of LTM.

Acetonitrile confined in silica nanopores with surfaces of varying functionality is studied by means of molecular dynamics simulation. The hydrogen-bonding interaction between the surface and the liquid is parametrized by means of first-principles molecular dynamics simulations. It is found that acetonitrile orders into bilayer like structures near the surface, in agreement with prior simulations and experiments. A newly developed method is applied to calculate relevant time correlation functions for molecules in different layers of the pore. This method takes into account the short lifetimes of the molecules in the layers. We compare this method with prior techniques that do not take this lifetime into account and discuss their pitfalls. We show that in agreement with experiment, the dynamics of the system may be described by a two population model that accounts for bulk-like relaxation in the center and frustrated dynamics near the surface of the pore. Specific hydrogen-bonding interactions are found to play a large role in engendering this behavior.

Solvent plays an important role in the relative motion of nanoscopic bodies, and the study of such phenomena can help elucidate the mechanism of hydrophobic assembly, as well as the influence of solvent-mediated effects on in vivo motion in crowded cellular environments. Here we study important aspects of this problem within the framework of Brownian dynamics. We compute the free energy surface that the Brownian particles experience and their hydrodynamic interactions from molecular dynamics simulations in explicit solvent. We find that molecular scale effects dominate at short distances, thus giving rise to deviations from the predictions of continuum hydrodynamic theory. Drying phenomena, solvent layering, and fluctuations engender distinct signatures of the molecular scale. The rate of assembly in the diffusion-controlled limit is found to decrease from molecular scale hydrodynamic interactions, in opposition to the free energy driving force for hydrophobic assembly, and act to reinforce the influence of the free energy surface on the association of more hydrophilic bodies.

We study the kinetics of assembly of two plates of varying hydrophobicity, including cases where drying occurs and water strongly solvates the plate surfaces. The potential of mean force and molecular-scale hydrodynamics are computed from molecular dynamics simulations in explicit solvent as a function of particle separation. In agreement with our recent work on nanospheres [J. Phys. Chem. B2012,116, 378–389], regions of high friction are found to be engendered by large and slow solvent fluctuations. These slow fluctuations can be due to either drying or confinement. The mean first passage times for assembly are computed by means of molecular dynamics simulations in explicit solvent and by Brownian dynamics simulations along the reaction path. Brownian dynamics makes use of the potential of mean force and hydrodynamic profile that we determined. Surprisingly, we find reasonable agreement between full-scale molecular dynamics and Brownian dynamics, despite the role of slow solvent relaxation in the assembly process. We found that molecular-scale hydrodynamic interactions are essential in describing the kinetics of assembly.

Biological processes often depend on protein–ligand binding events, yet accurate calculation of the associated energetics remains as a significant challenge of central importance to structure-based drug design. Recently, we have proposed that the displacement of unfavorable waters by the ligand, replacing them with groups complementary to the protein surface, is the principal driving force for protein–ligand binding, and we have introduced the WaterMap method to account this effect. However, in spite of the adage “nature abhors vacuum,” one can occasionally observe situations in which a portion of the receptor active site is so unfavorable for water molecules that a void is formed there. In this paper, we demonstrate that the presence of dry regions in the receptor has a nontrivial effect on ligand binding affinity, and suggest that such regions may represent a general motif for molecular recognition between the dry region in the receptor and the hydrophobic groups in the ligands. With the introduction of a term attributable to the occupation of the dry regions by ligand atoms, combined with the WaterMap calculation, we obtain excellent agreement with experiment for the prediction of relative binding affinities for a number of congeneric ligand series binding to the major urinary protein receptor. In addition, WaterMap when combined with the cavity contribution is more predictive than at least one specific implementation [Abel R, Young T, Farid R, Berne BJ, Friesner RA (2008) J Am Chem Soc 130:2817–2831] of the popular MM-GBSA approach to binding affinity calculation.

Calculation of protein−ligand binding affinities continues to be a hotbed of research. Although many techniques for computing protein−ligand binding affinities have been introduced—ranging from computationally very expensive approaches, such as free energy perturbation (FEP) theory, to more approximate techniques, such as empirically derived scoring functions, which, although computationally efficient, lack a clear theoretical basis—there remains a pressing need for more robust approaches. A recently introduced technique, the displaced-solvent functional (DSF) method, was developed to bridge the gap between the high accuracy of FEP and the computational efficiency of empirically derived scoring functions. In order to develop a set of reference data to test the DSF theory for calculating absolute protein−ligand binding affinities, we have pursued FEP theory calculations of the binding free energies of a methane ligand with 13 different model hydrophobic enclosures of varying hydrophobicity. The binding free energies of the methane ligand with the various hydrophobic enclosures were then recomputed by DSF theory and compared with the FEP reference data. We find that the DSF theory, which relies on no empirically tuned parameters, shows excellent quantitative agreement with the FEP. We also explored the ability of buried solvent accessible surface area and buried molecular surface area models to describe the relevant physics and find the buried molecular surface area model to offer superior performance over this data set.

Multiple time scale methodologies have gained widespread use in molecular dynamics simulations and are implemented in a variety of ways across numerous packages. However, performance of the algorithms depends upon the details of the implementation. This is particularly important in the way in which the nonbonded interactions are partitioned. In this work, we show why some previous implementations give rise to energy drifts, and how this can be corrected. We also provide a recipe for using multiple time step methods to generate stable trajectories in large scale biomolecular simulations, where long trajectories are needed.

In atomic force spectroscopic studies of the elastomeric protein ubiquitin, the β-strands 1-5 serve as the force clamp. Simulations show how the rupture force in the force-induced unfolding depends on the kinetics of water molecule insertion into positions where they can eventually form hydrogen bonding bridges with the backbone hydrogen bonds in the force-clamp region. The intrusion of water into this region is slowed down by the hydrophobic shielding effect of carbonaceous groups on the surface residues of β-strands 1-5, which thereby regulates water insertion prior to hydrogen bond breakage. The experiments show that the unfolding of the mechanically stressed protein is nonexponential due to static disorder. Our simulations show that different numbers and/or locations of bridging water molecules give rise to a long-lived distribution of transition states and static disorder. We find that slowing down the translational (not rotational) motions of the water molecules by increasing the mass of their oxygen atoms, which leaves the force field and thereby the equilibrium structure of the solvent unchanged, increases the average rupture force; however, the early stages of the force versus time behavior are very similar for our “normal” and fictitious “heavy” water models. Finally, we construct six mutant systems to regulate the hydrophobic shielding effect of the surface residues in the force-clamp region. The mutations in the two termini of β-sheets 1-5 are found to determine a preference for different unfolding pathways and change mutant’s average rupture force.

Because of its fundamental importance to molecular biology, great interest has continued to persist in developing novel techniques to efficiently characterize the thermodynamic and structural features of liquid water. A particularly fruitful approach, first applied to liquid water by Lazaridis and Karplus, is to use molecular dynamics or Monte Carlo simulations to collect the required statistics to integrate the inhomogeneous solvation theory equations for the solvation enthalpy and entropy. We here suggest several technical improvements to this approach, which may facilitate faster convergence and greater accuracy. In particular, we devise a nonparametric kth nearest-neighbors (NN)-based approach to estimate the water−water correlation entropy, and we suggest an alternative factorization of the water−water correlation function that appears to more robustly describe the correlation entropy of the neat fluid. It appears that the NN method offers several advantages over the more common histogram-based approaches, including much faster convergence for a given amount of simulation data; an intuitive error bound that may be readily formulated without resorting to block averaging or bootstrapping; and the absence of empirically tuned parameters, which may bias the results in an uncontrolled fashion.

The phase ordering of 1-bromoeicosane (C20H41Br) at the liquid-graphite and vacuum-graphite interfaces is examined through a joint experimental (part I) and theoretical effort (part II). The stable morphologies under solvent and ultrahigh vacuum conditions are revealed by STM experiments to be the head-to-head structures with 90° and 60° lamella−backbone angles, respectively. At 90° and 60° close packing is attained, independent of the corrugation of the graphite lattice. The potential energy of the minimized 60° structure is slightly lower than that of the 90° structure under vacuum conditions. In addition, the basin of the potential energy surface about the 90° form is very narrow. All-atom molecular dynamics simulations depict a 90°-to-60° phase transition in vacuum. Both morphologies are stable when an explicit solvent model is included. We speculate that the choice of the 90° form under solvent is driven by symmetry considerations and the self-assembly pathway. For example, the 90° structure may serve as a superior template for solvent capping. An implicit solvent model fails to stabilize the 90° form; however, it does lower the potential energy of this structure relative to the 60° geometry.

Scanning tunneling microscopy (STM) images of self-assembled monolayers of close-packed alkane chains on highly oriented pyrolitic graphite often display an alternating bright and dark spot pattern. Classical simulations suggest that a tilt of the alkane backbone is unstable and, therefore, unlikely to account for the contrast variation. First principles calculations based on density functional theory show that an electronic effect can explain the observed alternation. Furthermore, the asymmetric spot pattern associated with the minimum energy alignment is modulated depending on the registry of the alkane adsorbate relative to the graphite surface, explaining the characteristic moiré pattern that is often observed in STM images with close packed alkyl assemblies.

The role of many-body effects in modeling silica was investigated using self-consistent force matching. Both pairwise and polarizable classical force fields were developed systematically from ab initio density functional theory force calculations, allowing for a direct comparison of the role of polarization in silica. It was observed that the pairwise potential performed remarkably well at reproducing the basic silica tetrahedral structure. However, the Si−O−Si angle that links the silica tetrahedra showed small but distinct differences with the polarizable potential, a result of the inability of the pairwise potential to properly account for variations in the polarization of the oxygens. Furthermore, the transferability of the polarizable potential was investigated and suggests that additional forces may be necessary to more completely describe silica annealing.

The structural properties of self-assembled monolayers of short 1-bromoalkanes and n-alkanes on graphite were investigated by a combination of ultrahigh vacuum scanning tunneling microscopy (UHV-STM) at 80 K and theoretical methods. STM images of 1-bromohexane reveal a lamellar packing structure in which the molecules form a 57° ± 3° lamella-molecular backbone angle and a head-to-head assembly of the bromine atoms (Müller, et al. Proc. Natl. Acad. Sci. U.S.A. 2005, 102, 5315). STM images of 1-bromoheptane also show a head-to-head 60° ± 3° lamella-molecular backbone pattern; however, the molecules pack in a herringbone structure. The odd/even chain-length alternation in the monolayer morphologies of 1-bromoalkanes is similar to that observed for the self-assembly of short n-alkanes on graphite, suggesting that the bromine atom acts effectively as an extension of the carbon backbone. The analogy, however, is incomplete. Odd and even short n-alkanes (hexane, heptane, octane) display 60° herringbone and rectangular (not 60°) lamella-molecular backbone configurations, respectively. The balance of intermolecular forces and packing considerations responsible for this odd/even alternation in monolayer morphology for short 1-bromoalkanes on graphite is examined here according to classical molecular dynamics simulations and in light of the structural properties of analogous n-alkane assemblies.

The mechanism of denaturation of proteins by urea is explored by using all-atom microseconds molecular dynamics simulations of hen lysozyme generated on BlueGene/L. Accumulation of urea around lysozyme shows that water molecules are expelled from the first hydration shell of the protein. We observe a 2-stage penetration of the protein, with urea penetrating the hydrophobic core before water, forming a “dry globule.” The direct dispersion interaction between urea and the protein backbone and side chains is stronger than for water, which gives rise to the intrusion of urea into the protein interior and to urea's preferential binding to all regions of the protein. This is augmented by preferential hydrogen bond formation between the urea carbonyl and the backbone amides that contributes to the breaking of intrabackbone hydrogen bonds. Our study supports the “direct interaction mechanism” whereby urea has a stronger dispersion interaction with protein than water.

We propose a mechanism, based on a ≥10-μs molecular dynamics simulation, for the surprising misfolding of hen egg-white lysozyme caused by a single mutation (W62G). Our simulations of the wild-type and mutant lysozymes in 8 M urea solution at biological temperature (with both pH 2 and 7) reveal that the mutant structure is much less stable than that of the wild type, with the mutant showing larger fluctuations and less native-like contacts. Analysis of local contacts reveals that the Trp-62 residue is the key to a cooperative long-range interaction within the wild type, where it acts like a bridge between two neighboring basic residues. Thus, a native-like cluster or nucleation site can form near these residues in the wild type but not in the mutant. The time evolution of the secondary structure also exhibits a quicker loss of the β-sheets in the mutant than in the wild type, whereas some of the α-helices persist during the entire simulation in both the wild type and the mutant in 8 M urea (even though the tertiary structures are basically all gone). These findings, while supporting the general conclusions of a recent experimental study by Dobson and coworkers [Klein-Seetharam J, Oikama M, Grimshaw SB, Wirmer J, Duchardt E, Ueda T, Imoto T, Smith LJ, Dobson CM, Schwalbe H (2002) Science 295:1719–1722], provide a detailed but different molecular picture of the misfolding mechanism.

We examine the role of electronic polarizability in water on short (tens of femtoseconds), intermediate (hundreds of femtoseconds), and long (≈1 ps) time scales by comparing molecular dynamics results to experimental data for vibrational spectroscopy of HOD in liquid D2O. Because the OH absorption frequency is sensitive to the details of the atomic forces experienced in the liquid, our results provide important quantitative comparisons for several popular empirical water potentials. When compared with their fixed-charge counterparts, the polarizable models give similar slower long time constants for the decay of vibrational correlations and reorientational motion that is in better agreement with experiments. Polarizable potentials yield qualitatively dissimilar predictions for frequency fluctuations and transition dipole moment fluctuations at equilibrium. Models that confine the polarizability to the plane of the molecule (i.e., TIP4P–FQ) overestimate the width of the distribution describing frequency fluctuations by more than a factor of two. These models also underestimate the amplitude of the hydrogen-bond stretch at 170 cm–1. A potential that has both an out-of-plane polarization and fluctuating charges, POL5–TZ, compares best with experiments. We interpret our findings in terms of microscopic dynamics and make suggestions that may improve the quality of emerging polarizable force fields for water.

We studied by molecular dynamics simulations the temperature dependence of hydrophobic association and drying transition of large-scale solutes. Similar to the behavior of small solutes, we found the association process to be characterized by a large negative heat capacity change. The origin of this large change in heat capacity is the high fragility of hydrogen bonds between water molecules at the interface with hydrophobic solutes; an increase in temperature breaks more hydrogen bonds at the interface than in the bulk. With increasing temperature, both entropy and enthalpy changes for association strongly decrease, while the change in free energy weakly varies, exhibiting a small minimum at high temperatures. At around T = Ts = 360 K, the change in entropy is zero, a behavior similar to the solvation of small nonpolar solutes. Unexpectedly, we find that at Ts, there is still a substantial orientational ordering of the interfacial water molecules relative to the bulk. Nevertheless, at this point, the change in entropy vanishes due to a compensating contribution of translational entropy. Thus, at Ts, there is rotational order and translational disorder of the interfacial water relative to bulk water. In addition, we studied the temperature dependence of the drying−wetting transition. By calculating the contact angle of water on the hydrophobic surface at different temperatures, we compared the critical distance observed in the simulations with the critical distance predicted by macroscopic theory. Although the deviations of the predicted from the observed values are very small (8−23%), there seems to be an increase in the deviations with an increase in temperature. We suggest that these deviations emerge due to increased fluctuations, characterizing finite systems, as the temperature increases.