Bio


My research focuses on method development and multiscale simulation of chemical reactions in large-scale systems, especially in proteins of biomedical importance. The chemical reactions I studied include both electronically adiabatic ones, where only one electronic state is involved, as well as electronically non-adiabatic ones, where the coupled electron-nuclei motion leads to transitions among a manifold of electronic states.

Proton transport is a ubiquitous chemical reaction in biomolecules. Characterization of the free energy profile for proton transport events is key for understanding the structure-function relationship of a variety of proteins. However, this is a significant challenge for traditional computational methods because they are unable to accurately and efficiently treat the reactive nature of proton transport. To overcome this challenge, in my Ph.D. research under Prof. Gregory Voth at University of Chicago, I developed the multiscale reactive molecular dynamics (MS-RMD) method for accurate and efficient calculation of free energy profiles for proton transport in biomolecules. In addition, I applied this method to reveal the acid activation mechanism of influenza virus A M2 channel and proton pumping mechanism in cytochrome c oxidase (CcO). My work not only expanded the simulation tool for studying proton transport in biomolecules, but also provided new insights into the design principles of more efficient flu drugs and artificial proton pumps.

Electronically non-adiabatic chemical reactions are ubiquitous in the photochemistry of light-sensitive proteins. Understanding the dynamical aspect of these reactions are essential for designing new tools for photocontrol of biological systems. The main challenge is to accurately and efficiently simulate the coupled motion of nuclear and electronic degrees of freedoms in such reactions. To overcome this challenge, in my postdoctoral research under Prof. William Miller at UC Berkeley, I developed the “symmetrical quasi-classical windowing model applied to Meyer-Miller Hamiltonian” (SQC/MM) method for accurate and efficient simulation of electronically non-adiabatic dynamic processes, and applied it to studying the excitation energy transfer in organic semiconducting polymers. In addition, in my current research under Prof. Todd Martinez at Stanford University, I applied ab initio multiple spawning (AIMS) method to studying the photochemical reaction in channelrhodopsin 2 (ChR2), an important tool for photocontrol of neural signals in optogenetics. My work not only expanded the simulation tool for studying non-adiabatic dynamics of photochemical reactions in large-scale biomolecular systems, but also provided new insights into the functional mechanism of widely-used optogenetic tools.

Honors & Awards


  • Outstanding Young Researcher Award, From Computational Biophysics to Systems symposium (05/2015)
  • The Windt Memorial Fund Graduate Fellowship, The University of Chicago (05/2015)
  • Chemical Computing Group Excellence Award, American Chemical Society (08/2015)

Professional Education


  • Doctor of Philosophy, University of Chicago (2016)
  • Master of Science, University of Chicago, Chemistry (2012)
  • Bachelor of Science, Tsinghua University, Chemical Biology (2011)

Stanford Advisors


Current Research and Scholarly Interests


Channelrhodopsin 2 (ChR2) is a light-gated ion channel and an important tool in optogenetics. Photoisomerization of retinal protonated Schiff base (RPSB) in ChR2 triggers channel activation. Despite the importance of ChR2 in optogenetics, the detailed mechanism for photoisomerization and channel activation is still not fully understood. Here, we report on computer simulations to investigate the photoisomerization mechanism and its effect on the activation of ChR2. Nonadiabatic dynamics simulation of ChR2 was carried out using the ab initio multiple spawning (AIMS) method and QM/MM with a restricted ensemble Kohn-Sham (REKS) treatment of the QM region. Our results agree well with spectroscopic measurements, and reveal that the RPSB isomerization is highly specific around the C13=C14 bond and follows the “aborted bicycle-pedal” mechanism. In addition, RPSB photoisomerization facilitates its deprotonation and partially increases the hydration level in the channel, which could trigger subsequent channel opening and ion conduction. This work presents the first simulation of the photodynamics of RPSB isomerization in ChR2, and provides possible design principles for improving the optogenetic tools.

Projects


  • Multiscale simulation of the photodynamics of channelrhodopsin 2, Stanford University

    Channelrhodopsin 2 (ChR2) is a light-gated ion channel and an important tool in optogenetics. Photoisomerization of retinal protonated Schiff base (RPSB) in ChR2 triggers channel activation. Despite the importance of ChR2 in optogenetics, the detailed mechanism for photoisomerization and channel activation is still not fully understood. Here, we report on computer simulations to investigate the photoisomerization mechanism and its effect on the activation of ChR2. Nonadiabatic dynamics simulation of ChR2 was carried out using the ab initio multiple spawning (AIMS) method and QM/MM with a restricted ensemble Kohn-Sham (REKS) treatment of the QM region. Our results agree well with spectroscopic measurements, and reveal that the RPSB isomerization is highly specific around the C13=C14 bond and follows the “aborted bicycle-pedal” mechanism. In addition, RPSB photoisomerization facilitates its deprotonation and partially increases the hydration level in the channel, which could trigger subsequent channel opening and ion conduction. This work presents the first simulation of the photodynamics of RPSB isomerization in ChR2, and provides possible design principles for improving the optogenetic tools.

    Location

    Keck Building 185, 380 Roth Way, Stanford, CA 94305

Lab Affiliations


All Publications


  • First-Principles Characterization of the Elusive I Fluorescent State and the Structural Evolution of Retinal Protonated Schiff Base in Bacteriorhodopsin. Journal of the American Chemical Society Yu, J. K., Liang, R., Liu, F., Martinez, T. J. 2019

    Abstract

    The conversion of light energy into work is essential to life on earth. Bacteriorhodopsin (bR), a light-activated proton pump in Archae, has served for many years as a model system for the study of this process in photoactive proteins. Upon absorption of a photon, its chromophore, the retinal protonated Schiff base (RPSB), isomerizes from its native all-trans form to a 13-cis form and pumps a proton out of the cell in a process that is coupled to eventual ATP synthesis. Despite numerous time-resolved spectroscopic studies over the years, the details of the photodynamics of bR on the excited state, particularly the characterization of the I fluorescent state, the time-resolved reaction mechanism, and the role of the counterion cluster of RPSB, remain uncertain. Here, we use ab initio multiple spawning (AIMS) with spin-restricted ensemble Kohn-Sham (REKS) theory to simulate the nonadiabatic dynamics of the ultrafast photoreaction in bR. The excited state dynamics can be partitioned into three distinct phases: (1) relaxation away from the Franck-Condon region dominated by changes in retinal bond length alternation, (2) dwell time on the excited state in the I fluorescent state featuring an untwisted, bond length inverted RPSB, and (3) rapid torsional evolution to the conical intersection after overcoming a small excited state barrier. We fully characterize the I fluorescent state and the excited state barrier that hinders direct evolution to the conical intersection following photoexcitation. We also find that photoisomerization is accompanied by weakening of the interaction between RPSB and its counterion cluster. However, in contradiction with a recent time-resolved X-ray experiment, hydrogen bond cleavage is not necessary to reproduce the observed photoisomerization dynamics.

    View details for DOI 10.1021/jacs.9b08941

    View details for PubMedID 31621314

  • Nonadiabatic Photodynamics of Retinal Protonated Schiff Base in Channelrhodopsin 2 JOURNAL OF PHYSICAL CHEMISTRY LETTERS Liang, R., Liu, F., Martinez, T. J. 2019; 10 (11): 2862–68
  • Proton-Induced Conformational and Hydration Dynamics in the Influenza A M2 Channel. Journal of the American Chemical Society Watkins, L. C., Liang, R., Swanson, J. M., DeGrado, W. F., Voth, G. A. 2019; 141 (29): 11667–76

    Abstract

    The influenza A M2 protein is an acid-activated proton channel responsible for acidification of the inside of the virus, a critical step in the viral life cycle. This channel has four central histidine residues that form an acid-activated gate, binding protons from the outside until an activated state allows proton transport to the inside. While previous work has focused on proton transport through the channel, the structural and dynamic changes that accompany proton flux and enable activation have yet to be resolved. In this study, extensive Multiscale Reactive Molecular Dynamics simulations with explicit Grotthuss-shuttling hydrated excess protons are used to explore detailed molecular-level interactions that accompany proton transport in the +0, + 1, and +2 histidine charge states. The results demonstrate how the hydrated excess proton strongly influences both the protein and water hydrogen-bonding network throughout the channel, providing further insight into the channel's acid-activation mechanism and rectification behavior. We find that the excess proton dynamically, as a function of location, shifts the protein structure away from its equilibrium distributions uniquely for different pH conditions consistent with acid-activation. The proton distribution in the xy-plane is also shown to be asymmetric about the channel's main axis, which has potentially important implications for the mechanism of proton conduction and future drug design efforts.

    View details for DOI 10.1021/jacs.9b05136

    View details for PubMedID 31264413

  • Structural Coupling Throughout the Active Site Hydrogen Bond Networks of Ketosteroid Isomerase and Photoactive Yellow Protein JOURNAL OF THE AMERICAN CHEMICAL SOCIETY Pinney, M. M., Natarajan, A., Yabukarski, F., Sanchez, D. M., Liu, F., Liang, R., Doukov, T., Schwans, J. P., Martinez, T. J., Herschlag, D. 2018; 140 (31): 9827–43

    Abstract

    Hydrogen bonds are fundamental to biological systems and are regularly found in networks implicated in folding, molecular recognition, catalysis, and allostery. Given their ubiquity, we asked the fundamental questions of whether, and to what extent, hydrogen bonds within networks are structurally coupled. To address these questions, we turned to three protein systems, two variants of ketosteroid isomerase and one of photoactive yellow protein. We perturbed their hydrogen bond networks via a combination of site-directed mutagenesis and unnatural amino acid substitution, and we used 1H NMR and high-resolution X-ray crystallography to determine the effects of these perturbations on the lengths of the two oxyanion hole hydrogen bonds that are donated to negatively charged transition state analogs. Perturbations that lengthened or shortened one of the oxyanion hole hydrogen bonds had the opposite effect on the other. The oxyanion hole hydrogen bonds were also affected by distal hydrogen bonds in the network, with smaller perturbations for more remote hydrogen bonds. Across 19 measurements in three systems, the length change in one oxyanion hole hydrogen bond was propagated to the other, by a factor of -0.30 ± 0.03. This common effect suggests that hydrogen bond coupling is minimally influenced by the remaining protein scaffold. The observed coupling is reproduced by molecular mechanics and quantum mechanics/molecular mechanics (QM/MM) calculations for changes to a proximal oxyanion hole hydrogen bond. However, effects from distal hydrogen bonds are reproduced only by QM/MM, suggesting the importance of polarization in hydrogen bond coupling. These results deepen our understanding of hydrogen bonds and their networks, providing strong evidence for long-range coupling and for the extent of this coupling. We provide a broadly predictive quantitative relationship that can be applied to and can be further tested in new systems.

    View details for DOI 10.1021/jacs.8b01596

    View details for Web of Science ID 000441475800011

    View details for PubMedID 29990421

  • The symmetrical quasi-classical approach to electronically nonadiabatic dynamics applied to ultrafast exciton migration processes in semiconducting polymers JOURNAL OF CHEMICAL PHYSICS Liang, R., Cotton, S. J., Binder, R., Hegger, R., Burghardt, I., Miller, W. H. 2018; 149 (4): 044101

    Abstract

    In the last several years, a symmetrical quasi-classical (SQC) windowing model applied to the classical Meyer-Miller (MM) vibronic Hamiltonian has been shown to be a simple, efficient, general, and quite-accurate method for treating electronically nonadiabatic processes at the totally classical level. Here, the SQC/MM methodology is applied to ultrafast exciton dynamics in a Frenkel/site-exciton model of oligothiophene (OT) as a model of organic semiconductor polymers. In order to keep the electronic representation as compact and efficient as possible, the adiabatic version of the MM Hamiltonian was employed, with dynamical calculations carried out in the recently developed "kinematic momentum" representation, from which site/monomer-specific (diabatic) excitation probabilities were extracted using a new procedure developed in this work. The SQC/MM simulation results are seen to describe coherent exciton transport driven by planarization of a central torsion defect in the OT oligomer as well as to capture exciton self-trapping effects in good agreement with benchmark quantum calculations using the multi-layer multiconfiguration time-dependent Hartree approach. The SQC/MM calculations are also seen to significantly outperform the standard Ehrenfest approach, which shows serious discrepancies. These results are encouraging, not only because they illustrate a significant further application of the SQC/MM approach and its utility, but because they strongly suggest that classical mechanical simulations (with the potential for linear scaling efficiency) can be used to capture, quantitatively, important dynamical features of electronic excitation energy transfer in semiconducting polymers.

    View details for DOI 10.1063/1.5037815

    View details for Web of Science ID 000440586200005

    View details for PubMedID 30068189

  • On the adiabatic representation of Meyer-Miller electronic-nuclear dynamics JOURNAL OF CHEMICAL PHYSICS Cotton, S. J., Liang, R., Miller, W. H. 2017; 147 (6): 064112

    Abstract

    The Meyer-Miller (MM) classical vibronic (electronic + nuclear) Hamiltonian for electronically non-adiabatic dynamics-as used, for example, with the recently developed symmetrical quasiclassical (SQC) windowing model-can be written in either a diabatic or an adiabatic representation of the electronic degrees of freedom, the two being a canonical transformation of each other, thus giving the same dynamics. Although most recent applications of this SQC/MM approach have been carried out in the diabatic representation-because most of the benchmark model problems that have exact quantum results available for comparison are typically defined in a diabatic representation-it will typically be much more convenient to work in the adiabatic representation, e.g., when using Born-Oppenheimer potential energy surfaces (PESs) and derivative couplings that come from electronic structure calculations. The canonical equations of motion (EOMs) (i.e., Hamilton's equations) that come from the adiabatic MM Hamiltonian, however, in addition to the common first-derivative couplings, also involve second-derivative non-adiabatic coupling terms (as does the quantum Schrödinger equation), and the latter are considerably more difficult to calculate. This paper thus revisits the adiabatic version of the MM Hamiltonian and describes a modification of the classical adiabatic EOMs that are entirely equivalent to Hamilton's equations but that do not involve the second-derivative couplings. The second-derivative coupling terms have not been neglected; they simply do not appear in these modified adiabatic EOMs. This means that SQC/MM calculations can be carried out in the adiabatic representation, without approximation, needing only the PESs and the first-derivative coupling elements. The results of example SQC/MM calculations are presented, which illustrate this point, and also the fact that simply neglecting the second-derivative couplings in Hamilton's equations (and presumably also in the Schrödinger equation) can cause very significant errors.

    View details for DOI 10.1063/1.4995301

    View details for Web of Science ID 000407746400015

    View details for PubMedID 28810754

  • Understanding the essential proton-pumping kinetic gates and decoupling mutations in cytochrome c oxidase PROCEEDINGS OF THE NATIONAL ACADEMY OF SCIENCES OF THE UNITED STATES OF AMERICA Liang, R., Swanson, J. J., Wikstrom, M., Voth, G. A. 2017; 114 (23): 5924–29

    Abstract

    Cytochrome c oxidase (CcO) catalyzes the reduction of oxygen to water and uses the released free energy to pump protons against the transmembrane proton gradient. To better understand the proton-pumping mechanism of the wild-type (WT) CcO, much attention has been given to the mutation of amino acid residues along the proton translocating D-channel that impair, and sometimes decouple, proton pumping from the chemical catalysis. Although their influence has been clearly demonstrated experimentally, the underlying molecular mechanisms of these mutants remain unknown. In this work, we report multiscale reactive molecular dynamics simulations that characterize the free-energy profiles of explicit proton transport through several important D-channel mutants. Our results elucidate the mechanisms by which proton pumping is impaired, thus revealing key kinetic gating features in CcO. In the N139T and N139C mutants, proton back leakage through the D-channel is kinetically favored over proton pumping due to the loss of a kinetic gate in the N139 region. In the N139L mutant, the bulky L139 side chain inhibits timely reprotonation of E286 through the D-channel, which impairs both proton pumping and the chemical reaction. In the S200V/S201V double mutant, the proton affinity of E286 is increased, which slows down both proton pumping and the chemical catalysis. This work thus not only provides insight into the decoupling mechanisms of CcO mutants, but also explains how kinetic gating in the D-channel is imperative to achieving high proton-pumping efficiency in the WT CcO.

    View details for DOI 10.1073/pnas.1703654114

    View details for Web of Science ID 000402703800048

    View details for PubMedID 28536198

    View details for PubMedCentralID PMC5468613

  • Acid activation mechanism of the influenza A M2 proton channel PROCEEDINGS OF THE NATIONAL ACADEMY OF SCIENCES OF THE UNITED STATES OF AMERICA Liang, R., Swanson, J. J., Madsen, J. J., Hong, M., DeGrado, W. F., Voth, G. A. 2016; 113 (45): E6955–E6964

    Abstract

    The homotetrameric influenza A M2 channel (AM2) is an acid-activated proton channel responsible for the acidification of the influenza virus interior, an important step in the viral lifecycle. Four histidine residues (His37) in the center of the channel act as a pH sensor and proton selectivity filter. Despite intense study, the pH-dependent activation mechanism of the AM2 channel has to date not been completely understood at a molecular level. Herein we have used multiscale computer simulations to characterize (with explicit proton transport free energy profiles and their associated calculated conductances) the activation mechanism of AM2. All proton transfer steps involved in proton diffusion through the channel, including the protonation/deprotonation of His37, are explicitly considered using classical, quantum, and reactive molecular dynamics methods. The asymmetry of the proton transport free energy profile under high-pH conditions qualitatively explains the rectification behavior of AM2 (i.e., why the inward proton flux is allowed when the pH is low in viral exterior and high in viral interior, but outward proton flux is prohibited when the pH gradient is reversed). Also, in agreement with electrophysiological results, our simulations indicate that the C-terminal amphipathic helix does not significantly change the proton conduction mechanism in the AM2 transmembrane domain; the four transmembrane helices flanking the channel lumen alone seem to determine the proton conduction mechanism.

    View details for DOI 10.1073/pnas.1615471113

    View details for Web of Science ID 000388073300010

    View details for PubMedID 27791184

    View details for PubMedCentralID PMC5111692

  • Multiscale simulations reveal key features of the proton-pumping mechanism in cytochrome c oxidase PROCEEDINGS OF THE NATIONAL ACADEMY OF SCIENCES OF THE UNITED STATES OF AMERICA Liang, R., Swanson, J. J., Peng, Y., Wikstrom, M., Voth, G. A. 2016; 113 (27): 7420–25

    Abstract

    Cytochrome c oxidase (CcO) reduces oxygen to water and uses the released free energy to pump protons across the membrane. We have used multiscale reactive molecular dynamics simulations to explicitly characterize (with free-energy profiles and calculated rates) the internal proton transport events that enable proton pumping during first steps of oxidation of the fully reduced enzyme. Our results show that proton transport from amino acid residue E286 to both the pump loading site (PLS) and to the binuclear center (BNC) are thermodynamically driven by electron transfer from heme a to the BNC, but that the former (i.e., pumping) is kinetically favored whereas the latter (i.e., transfer of the chemical proton) is rate-limiting. The calculated rates agree with experimental measurements. The backflow of the pumped proton from the PLS to E286 and from E286 to the inside of the membrane is prevented by large free-energy barriers for the backflow reactions. Proton transport from E286 to the PLS through the hydrophobic cavity and from D132 to E286 through the D-channel are found to be strongly coupled to dynamical hydration changes in the corresponding pathways and, importantly, vice versa.

    View details for DOI 10.1073/pnas.1601982113

    View details for Web of Science ID 000379021700049

    View details for PubMedID 27339133

    View details for PubMedCentralID PMC4941487

  • Computationally Efficient Multiscale Reactive Molecular Dynamics to Describe Amino Acid Deprotonation in Proteins JOURNAL OF CHEMICAL THEORY AND COMPUTATION Lee, S., Liang, R., Voth, G. A., Swanson, J. J. 2016; 12 (2): 879–91

    Abstract

    An important challenge in the simulation of biomolecular systems is a quantitative description of the protonation and deprotonation process of amino acid residues. Despite the seeming simplicity of adding or removing a positively charged hydrogen nucleus, simulating the actual protonation/deprotonation process is inherently difficult. It requires both the explicit treatment of the excess proton, including its charge defect delocalization and Grotthuss shuttling through inhomogeneous moieties (water and amino residues), and extensive sampling of coupled condensed phase motions. In a recent paper (J. Chem. Theory Comput. 2014, 10, 2729-2737), a multiscale approach was developed to map high-level quantum mechanics/molecular mechanics (QM/MM) data into a multiscale reactive molecular dynamics (MS-RMD) model in order to describe amino acid deprotonation in bulk water. In this article, we extend the fitting approach (called FitRMD) to create MS-RMD models for ionizable amino acids within proteins. The resulting models are shown to faithfully reproduce the free energy profiles of the reference QM/MM Hamiltonian for PT inside an example protein, the ClC-ec1 H(+)/Cl(-) antiporter. Moreover, we show that the resulting MS-RMD models are computationally efficient enough to then characterize more complex 2-dimensional free energy surfaces due to slow degrees of freedom such as water hydration of internal protein cavities that can be inherently coupled to the excess proton charge translocation. The FitRMD method is thus shown to be an effective way to map ab initio level accuracy into a much more computationally efficient reactive MD method in order to explicitly simulate and quantitatively describe amino acid protonation/deprotonation in proteins.

    View details for DOI 10.1021/acs.jctc.5b01109

    View details for Web of Science ID 000370112900038

    View details for PubMedID 26734942

    View details for PubMedCentralID PMC4750100

  • Multiscale simulation reveals a multifaceted mechanism of proton permeation through the influenza A M2 proton channel PROCEEDINGS OF THE NATIONAL ACADEMY OF SCIENCES OF THE UNITED STATES OF AMERICA Liang, R., Li, H., Swanson, J. J., Voth, G. A. 2014; 111 (26): 9396–9401

    Abstract

    The influenza A virus M2 channel (AM2) is crucial in the viral life cycle. Despite many previous experimental and computational studies, the mechanism of the activating process in which proton permeation acidifies the virion to release the viral RNA and core proteins is not well understood. Herein the AM2 proton permeation process has been systematically characterized using multiscale computer simulations, including quantum, classical, and reactive molecular dynamics methods. We report, to our knowledge, the first complete free-energy profiles for proton transport through the entire AM2 transmembrane domain at various pH values, including explicit treatment of excess proton charge delocalization and shuttling through the His37 tetrad. The free-energy profiles reveal that the excess proton must overcome a large free-energy barrier to diffuse to the His37 tetrad, where it is stabilized in a deep minimum reflecting the delocalization of the excess charge among the histidines and the cost of shuttling the proton past them. At lower pH values the His37 tetrad has a larger total charge that increases the channel width, hydration, and solvent dynamics, in agreement with recent 2D-IR spectroscopic studies. The proton transport barrier becomes smaller, despite the increased charge repulsion, due to backbone expansion and the more dynamic pore water molecules. The calculated conductances are in quantitative agreement with recent experimental measurements. In addition, the free-energy profiles and conductances for proton transport in several mutants provide insights for explaining our findings and those of previous experimental mutagenesis studies.

    View details for DOI 10.1073/pnas.1401997111

    View details for Web of Science ID 000338118900029

    View details for PubMedID 24979779

    View details for PubMedCentralID PMC4084430

  • Benchmark Study of the SCC-DFTB Approach for a Biomolecular Proton Channel JOURNAL OF CHEMICAL THEORY AND COMPUTATION Liang, R., Swanson, J. J., Voth, G. A. 2014; 10 (1): 451–62

    Abstract

    The self-consistent charge density functional tight binding (SCC-DFTB) method has been increasingly applied to study proton transport (PT) in biological environments. However, recent studies revealing some significant limitations of SCC-DFTB for proton and hydroxide solvation and transport in bulk aqueous systems call into question its accuracy for simulating PT in biological systems. The current work benchmarks the SCC-DFTB/MM method against more accurate DFT/MM by simulating PT in a synthetic leucine-serine channel (LS2), which emulates the structure and function of biomolecular proton channels. It is observed that SCC-DFTB/MM produces over-coordinated and less structured pore water, an over-coordinated excess proton, weak hydrogen bonds around the excess proton charge defect and qualitatively different PT dynamics. Similar issues are demonstrated for PT in a carbon nanotube, indicating that the inaccuracies found for SCC-DFTB are not due to the point charge based QM/MM electrostatic coupling scheme, but rather to the approximations of the semiempirical method itself. The results presented in this work highlight the limitations of the present form of the SCC-DFTB/MM approach for simulating PT processes in biological protein or channel-like environments, while providing benchmark results that may lead to an improvement of the underlying method.

    View details for DOI 10.1021/ct400832r

    View details for Web of Science ID 000330142400043

    View details for PubMedID 25104919

    View details for PubMedCentralID PMC4120842

  • Application of the SCC-DFTB Method to Hydroxide Water Clusters and Aqueous Hydroxide Solutions JOURNAL OF PHYSICAL CHEMISTRY B Choi, T., Liang, R., Maupin, C., Voth, G. A. 2013; 117 (17): 5165–79

    Abstract

    The self-consistent charge density functional tight binding (SCC-DFTB) method has been applied to hydroxide water clusters and a hydroxide ion in bulk water. To determine the impact of various implementations of SCC-DFTB on the energetics and dynamics of a hydroxide ion in gas phase and condensed phase, the DFTB2, DFTB2-γ(h), DFTB2-γ(h)+gaus, DFTB3-diag, DFTB3-diag+gaus, DFTB3-Full+gaus, and DFTB3-3OB implementations have been tested. Energetic stabilities for small hydroxide clusters, OH(-)(H2O)n, where n = 4-7, are inconsistent with the results calculated with the B3LYP and second order Møller-Plesset (MP2) levels of ab initio theory. The condensed phase simulations, OH(-)(H2O)127, using the DFTB2, DFTB2-γ(h), DFTB2-γ(h)+gaus, DFTB3-diag, DFTB3-diag+gaus, DFTB3-Full+gaus and DFTB3-3OB methods are compared to Car-Parrinello molecular dynamics (CPMD) simulations using the BLYP functional. The SCC-DFTB method including a modified O-H repulsive potential and the third order correction (DFTB3-diag/Full+gaus) is shown to poorly reproduce the CPMD computational results, while the DFTB2 and DFTB2-γ(h) method somewhat more closely describe the structural and dynamical nature of the hydroxide ion in condensed phase. The DFTB3-3OB outperforms the MIO parameter set but is no more accurate than DFTB2. It is also shown that the overcoordinated water molecules lead to an incorrect bulk water density and result in unphysical water void formation. The results presented in this paper point to serious drawbacks for various DFTB extensions and corrections for a hydroxide ion in aqueous environments.

    View details for DOI 10.1021/jp400953a

    View details for Web of Science ID 000318536700047

    View details for PubMedID 23566052