Chris Knight

The water–amorphous silica interface: Analysis of the Stern layer and surface conduction

To explain why dynamical properties of an aqueous electrolyte near a charged surface seem to be governed by a surface charge less than the actual one, the canonical Stern model supposes an interfacial layer of ions and immobile fluid. However, large ion mobilities within the Stern layer are needed to reconcile the Stern model with surface conduction measurements. Modeling the aqueous electrolyte-amorphous silica interface at typical charge densities, a prototypical double layer system, the flow velocity does not vanish until right at the surface. The Stern model is a good effective model away from the surface, but cannot be taken literally near the surface. Indeed, simulations show no ion mobility where water is immobile, nor is such mobility necessary since the surface conductivity in the simulations is comparable to experimental values.

Multiscale reactive molecular dynamics

Many processes important to chemistry, materials science, and biology cannot be described without considering electronic and nuclear-level dynamics and their coupling to slower, cooperative motions of the system. These inherently multiscale problems require computationally efficient and accurate methods to converge statistical properties. In this paper, a method is presented that uses data directly from condensed phase ab initio simulations to develop reactive molecular dynamics models that do not require predefined empirical functions. Instead, the interactions used in the reactive model are expressed as linear combinations of interpolating functions that are optimized by using a linear least-squares algorithm. One notable benefit of the procedure outlined here is the capability to minimize the number of parameters requiring nonlinear optimization. The method presented can be generally applied to multiscale problems and is demonstrated by generating reactive models for the hydrated excess proton and hydroxide ion based directly on condensed phase ab initio molecular dynamics simulations. The resulting models faithfully reproduce the water-ion structural properties and diffusion constants from the ab initio simulations. Additionally, the free energy profiles for proton transfer, which is sensitive to the structural diffusion of both ions in water, are reproduced. The high fidelity of these models to ab initio simulations will permit accurate modeling of general chemical reactions in condensed phase systems with computational efficiency orders of magnitudes greater than currently possible with ab initio simulation methods, thus facilitating a proper statistical sampling of the coupling to slow, large-scale motions of the system.

The Dissociated Amorphous Silica Surface: Model Development and Evaluation

At pH 7, amorphous silica has a characteristic negative charge due to the deprotonation of silanol groups on the surface. Electrokinetic phenomena and transport of biomolecules in devices depend sensitively on the surface morphology, distribution of ions and solvent, and adsorption properties of solutes close to the surface in the electrical double layer region. Hence, simulation of these phenomena requires detailed atomistic models of the double layer region. In this Article, we extend our undissociated silica surface model [J. Phys. Chem. B 2007, 111, 11181-11193] to include dissociated Si-O(-) groups, which interact with both water and salt (Na(+) and Cl(-)). We have also conducted ab initio molecular dynamics (AIMD) simulations of a smaller system consisting of a hydrated silica slab. The radial distribution functions predicted by the empirical model are in qualitative agreement with those from the AIMD simulations. The hydrophobic and hydrophilic nature of silanol-poor and silanol-rich regions of the amorphous silica surface observed in our empirical model is reproduced in the AIMD simulations of the smaller slab. In the initial stages of our AIMD simulations, we observe various chemical processes that represent different hydroxylation mechanisms of the surface.

A reexamination of the ice III/IX hydrogen bond ordering phase transition.

Ice III is a hydrogen bond disordered crystal which when cooled 1 K / min or faster transforms to an antiferroelectric hydrogen bond ordered structure, ice IX. Throughout its region of stability, experiments indicate that the H bonds in ice III are, in fact, partially ordered, i.e., some proton arrangements are preferred. In addition, there has been evidence that the structure of ice IX retains some residual disorder after the transition. Diffraction experiments and calorimetry apparently conflict with regard to the degree of ordering at the ice III/IX transition. Mean field statistical mechanical theories have been used to link partial occupations from diffraction data with thermodynamics. In this work, we investigate the ice III/IX proton ordering phase transition using electronic density functional theory calculations for small unit cells, extended to simulate the phase transition in a large unit cell using graph invariants. In agreement with experiment, we observe partial ordering over a wide range of temperatures as ice III transforms to partially disordered ice IX, near 126 K, which becomes fully ordered at lower temperatures. We compare our results from full statistical mechanical simulations with mean field models, finding small errors for the low-temperature ice IX phase and much larger errors for the high-temperature ice III phase. The failure of mean field theories may explain the apparent conflict between diffraction experiments and calorimetry.

Hydrogen bond ordering in ice V and the transition to ice XIII

The proton ordered version of ice V, ice XIII, was recently identified using Raman spectroscopy and neutron diffraction techniques. The transformation, between 108 and 117 K, only occurred in the presence of a small amount of dopant, similar to the proton ordering transition of ice Ih/XI. In this work, we investigate the hydrogen bond fluctuations in ice V and XIII with statistical mechanical techniques that use results from periodic electronic density functional theory calculations as input. We find a number of low-lying hydrogen bond configurations, approximately 20 within 10 K/water above the ground state state configuration, the structure of which agrees with fully ordered ice XIII. Using an analytic theory, graph invariants, we developed effective spin-lattice Hamiltonians governing hydrogen bond fluctuations to perform statistical mechanical calculations for a large simulation cell containing 6048 water molecules. Two models were constructed, one more elaborate than the first, to gauge the robustness of our methods when the unit cell is very complex and a large number of configurations lie close in energy to the ground state. The predicted proton ordering transitions, 62 and 72 K for the two models, are in qualitative agreement with experiment. Occupation probabilities, obtained from our simulations, compare well with values from recent neutron diffraction experiments and help verify our effective Hamiltonians. In both models, we find that a second order phase transition intervenes about 10 K above the transition to ice XIII, but its effect is negligible on the behavior of thermodynamic functions near the transition to ice XIII.

Exploring the behaviour of the hydrated excess proton at hydrophobic interfaces

The affinity of the excess proton for the aqueous solution-hydrophobic interface was examined for two specific examples, the air-water and hydrophobic wall-water cases, using a multiconfigurational molecular dynamics algorithm. The use of a reactive simulation method is important as it allows for a realistic description of the excess proton, namely, its propensity to hop between water molecules via the Grotthuss mechanism. The free energy profile reveals a minimum at these interfaces due to a favourable enthalpic term that outweighs the entropic penalty. The key factors that contribute to this enthalpic minimum were examined using a generalization of a scheme that decomposes the interaction energy into separate terms arising from various local environments [Otten et al., Proc. Natl. Acad. Sci. USA, 109, 701 (2012)] (coordination shell, bulk, and interface) and the delocalization energy (which allows the proton to hop). For both systems, it was observed that the energetic penalty for loss of coordinating water molecules as the excess proton moves toward the hydrophobic interface is more than compensated by the displacement of unfavourable interfacial water molecules. In addition, the ion becomes more delocalized, more Zundel-like, and therefore possesses a larger effective radius as it moves to the interface. The fluctuations of the instantaneous interface were reduced near the vicinity of the ion, thereby giving rise to an entropic penalty. This paper will discuss the application of energy decomposition schemes to multiconfigurational simulations and the resulting consequences realized for the excess proton at hydrophobic interfaces.

Many-Body Interactions in Ice

Many-body effects in ice are investigated through a systematic analysis of the lattice energies of several proton ordered and disordered phases, which are calculated with different flexible water models, ranging from pairwise additive (q-TIP4P/F) to polarizable (TTM3-F and AMOEBA) and explicit many-body (MB-pol) potential energy functions. Comparisons with available experimental and diffusion Monte Carlo data emphasize the importance of an accurate description of the individual terms of the many-body expansion of the interaction energy between water molecules for the correct prediction of the energy ordering of the ice phases. Further analysis of the MB-pol results, in terms of fundamental energy contributions, demonstrates that the differences in lattice energies between different ice phases are sensitively dependent on the subtle balance between short-range two-body and three-body interactions, many-body induction, and dispersion energy. By correctly reproducing many-body effects at both short range and long range, it is found that MB-pol accurately predicts the energetics of different ice phases, which provides further support for the accuracy of MB-pol in representing the properties of water from the gas to the condensed phase.

Site disorder in ice VII arising from hydrogen bond fluctuations.

Recently, multisite models used for the refinement of neutron diffraction data have suggested that the structure of ice VII is quite unlike that of its ordered counterpart, ice VIII. We investigate the oxygen site disorder by modeling the site displacement, obtained from periodic DFT calculations, as a function of the local hydrogen bond network. Then, using graph invariants to describe hydrogen bond fluctuations in the thermodynamic limit, we perform statistical mechanical calculations using the oxygen site displacement model developed here. We find that the probability distribution of the oxygen atom about its perfect lattice site more closely resembles the 100 model rather than the recently suggested 111 model, although both models represent a simplification of the actual site distribution. We also find a unimodel distribution for the hydrogen bonded oxygen-oxygen distance and a trimodel distribution for the nearest nonbonded oxygen-oxygen distance with a peak separation of approximately 0.1 A.

Bringing Reactivity to the Aggregation-Volume-Bias Monte Carlo Based Simulation Framework: Water Nucleation Induced by a Reactive Proton.

The development of the aggregation-volume-bias Monte Carlo based simulation technique has led to recent success in studying rare nucleation events, but thus far, this simulation method has been limited to nonreactive systems. This work presents the first application of this technique to study a reactive system of relevance to atmospheric chemistry, i.e., formation of water droplets in the presence of a reactive proton, by combining this approach with a multistate empirical valence bond (MSEVB) description of the excess proton (or the hydronium). It was shown that the ability for the hydronium to share its charge with adjacent water molecules changes dramatically with the cluster size, especially when clusters are small and the distribution of the charge is affected by the presence of an interface, emphasizing the need to use this more sophisticated MSEVB model for such a reactive system. In addition, the simulation results obtained from this system are compared to those with nonreactive hard-sphere ions of different sizes. Overall, the presence of a hydronium or ions appeared to dramatically change the free energy landscape of nucleation compared to the pure water system, leading to the formation of a stable precritical cluster. Although the free energy change due to the addition of the first few water molecules was shown to be very sensitive to the ionic details, the later portion of the free energy profile was found to be nearly independent of the nature of the ion.