A Simple Additive Potential Model for Simulating ... - ACS Publications

A Simple Additive Potential Model for Simulating Hydrogen Peroxide in Chemical and Biological Systems. Esam A. Orabi* and Ann M. English*. Center for ...
0 downloads 0 Views 2MB Size
Subscriber access provided by Eastern Michigan University | Bruce T. Halle Library

Structure Prediction

A Simple Additive Potential Model for Simulating Hydrogen Peroxide in Chemical and Biological Systems Esam A. Orabi, and Ann M. English J. Chem. Theory Comput., Just Accepted Manuscript • DOI: 10.1021/acs.jctc.8b00246 • Publication Date (Web): 09 Apr 2018 Downloaded from http://pubs.acs.org on April 11, 2018

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

A Simple Additive Potential Model for Simulating Hydrogen Peroxide in Chemical and Biological Systems

Esam A. Orabi* [a,†] and Ann M. English*[a,b]

[a] Department of Chemistry and Biochemistry, Concordia University, 7141 Sherbrooke Street West, Montréal, Québec H4B 1R6, Canada [b] Center for Research in Molecular Modeling (CERMM) and Quebec Network for Research on Protein Function, Engineering, and Applications (PROTEO)

ACS Paragon Plus Environment

1

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 2 of 47

Abstract Hydrogen peroxide (H2O2) has numerous industrial, environmental, medical, cosmetic and biological applications. Given its importance, we provide a simple model as an alternative to experiment for studying the properties of pure liquid H2O2 and its concentrated aqueous solutions, which are hazardous, and for understanding the biological roles of H2O2 at the molecular level. A four-site additive model is calibrated for H2O2 based on the ab initio and experimental properties of the gaseous monomer and the density and heat of vaporization of liquid H2O2 at 0 °C. Our model together with the TIP3P water model reproduce the ab initio binding energies of (H2O2)m, H2O2.(H2O)n and H2O.(H2O2)n clusters (m = 2, 3 and n = 1, 2) calculated at the MP2 level using the 6-311++G(d,p) or the 6-311++G(3df,3pd) basis set. It yields structure, the self-diffusion coefficient, heat capacity and densities at temperatures up to 200 ºC of the pure liquid in good agreement with experiment. The models correctly predicts the hydration free energy of H2O2 and reproduces the experimental density of aqueous H2O2 solutions at 0‒96 °C. Investigation of the solvation of H2O2 and H2O in aqueous H2O2 solutions reveals that, as in the gas phase, H2O2 is a better H-bond donor but poorer acceptor than H2O and the bonding stability follows the order: Op‒Hp···Ow > Ow‒Hw···Ow ≥ Op‒Hp···Op > Ow‒Hw···Op. Stronger H-bonding in H2O/H2O2 mixtures than in the pure liquids is consistent with exothermic heats of mixing and explains why the observed density and vapor pressure of the aqueous solutions are higher and lower, respectively, than expected from ideal mixing. Results also show that H2O2 adopts a skewed equilibrium geometry in gas and liquid phases but more polar cis and nonpolar trans conformations also are accessible and will stabilize H2O2 in environments of different polarity. In sum, our simple model presents a reliable tool for simulating H2O2 in chemistry and biology.

ACS Paragon Plus Environment

2

Page 3 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

1. Introduction Hydrogen peroxide was first prepared in 1818 by reacting barium peroxide with nitric acid.1 Photochemical H2O2 production occurs in surface waters and ozone and aldehyde photolysis is a source of atmospheric H2O2, while anthrahydroquinone autoxidation is the most widely used commercial production method.2 Aerobic organisms produce H2O2 mostly3 via the spontaneous or catalytic breakdown of superoxide anions4 and from NADPH oxidase-independent mechanisms.5 Industrial uses of H2O2 include the bleaching of wood pulp, cotton, wool, silk and solid surfaces.2,6,7 Environmentally, it is used in the treatment of domestic and industrial effluents and in the removal of toxic or malodorous compounds from industrial gas streams.2 Also, it is a major oxidant of sulfur dioxide to sulfuric acid in clouds. In the chemical industry H2O2 serves as a reagent in the synthesis of organic and inorganic peroxides, the epoxidation of oils and unsaturated esters, the hydroxylation of phenols and the preparation of detergents.2,7 It also is employed as a propellant in the polishing and treatment of metal surfaces. H2O2 concentrations in excess of 10% are hazardous2,8 but dilute solutions (< 6%) have numerous medical and cosmetic applications because of their disinfecting, antiseptic and mild bleaching properties. H2O2 was considered a toxic byproduct of aerobic metabolism in biological systems but we now know that it vitally regulates eukaryotic signal transduction3,4 and is essential in the innate immune system.9 It becomes cytotoxic at elevated concentrations and biological antioxidant defenses include enzymes, catalases and peroxidases, that efficiently scavenge H2O2.3,4 The properties of H2O2 have been extensively studied because of its importance. Available experimental data include crystal structures,7,10–14 vibrational spectra,15–22 vapor pressure,23 density and viscosity,24,25 heat of vaporization26,27 and heat capacity.27 Aqueous H2O2

ACS Paragon Plus Environment

3

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 4 of 47

solutions also have been investigated experimentally and crystal structure,28 vapor pressure,29 density24,30,31 and freezing point32 measurements have been reported but these have been limited to temperatures below 100 °C due to the dangers associated with explosions above 150 ºC.8 Various theoretical studies have been carried out in the gas phase on the H2O2 monomer,33–42 H2O2 clusters36,38,40,41 and H2O2-H2O clusters.43–49 However, only a few conformers have been reported for these clusters and complexes between H2O with more than one H2O2 have not been reported to date. Optimization of a force field (FF) for H2O2 is required to study pure liquid H2O2 and its concentrated solutions at the high temperatures found in industry and the environment. Reliable FFs also are needed for elucidating differences in H2O2 vs. H2O association and for investigating and understanding the biological roles of H2O2 at the molecular level. Fortunately, the literature provides a wealth of experimental and theoretical data for accurate FF calibration. Although the properties of liquid37,50 and aqueous H2O246,50 and of H2O2 in biology51–54 have been studied using a number of potential models, these are not calibrated for H2O2 interactions with itself or with H2O. An exception is the ABEEM/MM fluctuating charge model developed by Yu and Yang.40 This model was optimized based on the properties of (H2O2)n clusters (n = 1‒6) and was used to investigate the properties of liquid H2O2 at 0 ºC. It consists of 11 charge sites (4 atomic, 4 lone pairs and three along the O‒O and O‒H bonds), which make it computationally expensive relative to the model presented here. Here we demonstrate that a simple additive four-atomic-sites model for H2O2 reliably predicts the structural and thermodynamic properties of pure and aqueous H2O2 liquids. Our model is calibrated based on the experimental and ab initio properties of the H2O2 monomer and the density and heat of vaporization of liquid H2O2 at 0 °C. We use our model to investigate the

ACS Paragon Plus Environment

4

Page 5 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

properties of (H2O2)m, H2O2.(H2O)n and H2O.(H2O2)n clusters (m = 2, 3 and n = 1, 2) and the structural and thermodynamic properties of pure and aqueous H2O2 solutions under a wide range of temperatures (−20 to 200 ºC) and pressures (0.069 to 3132.41 mm Hg). Our model will be valuable for simulating industrial, environmental and atmospheric applications and for exploring the biological functions of H2O2.

2. Computational Methods 2.1 Ab initio Calculations These were performed using Gaussian 09 software.55 The H2O2 monomer and small H2O2 complexes are amenable to high-level ab initio quantum mechanical calculations. Correlation effects are important for describing the geometry of H2O233 and its experimental geometry21 can be

closely

reproduced

using

MP2/6-311G(3d,2p),34

B3LYP/6-311++G(3df,3pd)

and

CCSD(T)/aug-cc-pVTZ model chemistry.39 Geometry optimization of (H2O2)m, H2O2.(H2O)n and H2O.(H2O2)n clusters (m = 1−3 and n = 1, 2) were performed here at the MP2(full) level using the 6-311++G(d,p) and the 6-311++G(3df,3pd) basis sets. The structures of the complexes were optimized starting from various initial binding conformations generated by placing the interacting molecules in linear, bent or cage arrangements, while considering different relative positions and numbers of H-bonds in each arrangement. The structures were first optimized using the smaller basis set and the resulting structures were further optimized using the larger basis set. Unless otherwise stated, no geometry constraints were imposed during optimization and all reported structures are energy-minima (no imaginary frequencies). The geometry of the monomer is also optimized at the CCSD(T)/6-311++G(3df,3pd) level. The counterpoise (CP) procedure of Boys and Bernardi56 was used to correct binding energies for basis set superposition error (BSSE) and both corrected (ECP) and uncorrected (E) binding energies are reported.

ACS Paragon Plus Environment

5

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 47

The barriers associated with the rotation of the O‒H bonds relative to the O‒O bond are investigated from a relaxed potential energy scan of the HOOH dihedral angle from 0° to 180° in 5° increments at the MP2(full) level using both basis sets. To understand differences in the structural and binding energies of H2O2- vs H2O-complexes, electrostatic potential surfaces are calculated for H2O2 and H2O at the MP2(full)/6-311++G(3df,3pd) level. 2.2 Molecular Mechanics Calculations High level ab initio calculations remain computationally prohibitive for studying bulk liquids so molecular mechanical (MM) calculations are typically used. MM calculations were performed with the program CHARMM57 using our model for H2O2 and the CHARMM-compatible TIP3P model for H2O.58,59 The MM-calculated binding energy (EMM) is the difference in energy between a complex and its isolated ligands with all species constrained at their ab initiooptimized geometry. 2.2.1 Development of a FF for H2O2

The potential energy function U(r) used in the all-atom additive CHARMM FF is given by:57    =    −  +   −  



+ + +



 



(,--



  −   +

* + − +  +



7

4,(,01 ./01 23 6 01



()

∅ 1 + cos&∅ − '

9

4,(,01 ;0 ;1 −23 6 : + < 1 01 01

U(r) depends on bond distance (b), valence angle (θ), Urey‒Bradley 1,3-distance (S), dihedral angle (Φ) and improper dihedral angle (ω), where the subscript zero indicates the equilibrium

ACS Paragon Plus Environment

6

Page 7 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

value. Kb, Kθ, KS, KΦ and Kω are the corresponding force constants. Parameters n and δ in the dihedral term represent the multiplicity and phase angle. The last term in eq 1 describes nonbonded interactions between atoms i and j via van der Waals or LJ (6‒12) and electrostatic or Coulomb terms. The partial charge on atom i is qi and rij is the distance between nonbonded atoms i and j. The minimum interaction radius Rmin,ij and the well-depth ϵij in the LJ (6‒12) potential are calculated from the individual LJ parameters of atoms i and j using the Lorentz‒ Berthelot combination rules: /01 = =/( × /? and 4,(,(? =

4,(,( + 4,(,? 2 2

To simplify the H2O2 model, no Urey‒Bradley or improper potential-energy terms are used. The b0 for OH and OO bonds and θ0 for the HOO angle were obtained from geometry optimization of the H2O2 monomer at the MP2(full)/6-311++G(3df,3pd) level. Kb and Kθ were adjusted to reproduce the experimental bond-stretch and angle-bending vibrational frequencies. The atomic charges (;C and ;D ) and the dihedral parameters (n, δ, and KΦ) were adjusted to yield dipole moment and energy as functions of ΦHOOH in agreement with ab initio calculations. The atomic LJ parameters of H are set to those of the TIP3P water model in CHARMM (ϵH = 0.046 kcal/mol and Rmin,H = 0.449 Å)59 and the LJ parameters for O are optimized to reproduce the density and enthalpy of vaporization of liquid H2O2 at 0 °C. We use the LJ parameters of O in the TIP3P model as a first guess to simulate a system of 500 H2O2 molecules for 2 ns and the parameters were varied until the experimental density and heat of vaporization were reproduced. 2.2.2 Molecular Dynamics

Unless otherwise stated, MD simulations were performed with cubic periodic boundary conditions (PBC) in an isothermal-isobaric ensemble (NpT). Simulations of liquid H2O2 were performed in a box of 500 H2O2 molecules at various temperatures using vapor pressures

ACS Paragon Plus Environment

7

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 47

estimated from the reported relationship between p and T.23 In addition, simulations of 250 and 1000 H2O2 molecules were performed at 0 ºC and 0.37 mm Hg to check if the self-diffusion coefficient depends on system size.60 The solvation structure of H2O2 in liquid H2O and of H2O in liquid H2O2 was simulated by adding one solute molecule to 500 solvent molecules at 25 ºC and 760.00 mm Hg. The properties of aqueous H2O2 solutions were investigated from MD simulations of 500 molecules (H2O2+H2O) where the H2O2 mole fraction (EFG HG ) was increased from 0.0 to 1.0 in increments of 0.1. The aqueous H2O2 solutions were simulated at 0, 25, 50 and 96 ºC using the corresponding vapor pressure of pure water (4.59, 23.79, 92.69 and 658.58 mm Hg). All simulations were run for 17 ns and the first 2 ns are considered the equilibration period and discarded. Electrostatic interactions are computed using the particle-mesh Ewald method61 with I = 0.33 for charge screening and a 1.0 Å grid spacing with fourth-order splines for mesh interpolation. Real-space interactions (LJ and electrostatic) are cut off at 15 Å. A Nosé−Hoover thermostat62 and an Andersen−Hoover barostat63 maintained the temperature and pressure at the preset values with relaxation times of 0.1 and 0.2 ps, respectively. The equations of motion are integrated in 1fs steps. All bonds involving H atoms are kept at their reference lengths using the RATTLE/Roll algorithm.64 2.2.3. Free Energy Calculations

The transferability of our H2O2 model and the TIP3 water model58,59 for aqueous H2O2 is examined by calculating the hydration free energy (∆KLMN ) of H2O2 in liquid H2O at 25 ºC and

760.00 mm Hg. ∆KLMN is decomposed into electrostatic (∆KO ) plus attractive (∆K(- ) and

repulsive components P∆K- Q of the LJ interaction:65

∆KLMN = ∆KO + ∆K(- + ∆K- 3

ACS Paragon Plus Environment

8

Page 9 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

∆KO and ∆K(- are computed using thermodynamic integration (λ = 0, 0.1, 0.2, …1) and ∆K- is computed (λ = 0, 0.05, 0.1, 0.15, 0.2, 0.25, 0.3, 0.35, 0.4, 0.5, …1) using a soft-core

scheme65 and unbiased using the weighted histogram analysis method (WHAM).66

3. Results and Discussion 3.1 Ab initio Geometries and Binding Energies Experimental values of 119.1° and 111.5° have been estimated for the equilibrium dihedral angle of H2O2 (ΦHOOH) based on microwave19,20 and far-infrared spectra,16 respectively. Fitting the spectral transitions to a large-amplitude Hamiltonian to account for vibration-torsion-rotation interaction, yielded ΦHOOH = 111.5°‒111.8º.21 Thus, H2O2 adopts a skewed structure in its global-minimum geometry, making it one of the smallest molecules with chirality and capable of internal rotation (i.e., rotation of the OH bonds around the peroxide bond). Because its chirality depends on internal rotation, H2O2 is of interest in chiral separations and in studies of chirality exchange during molecular collisions.67 Geometry optimization of the H2O2 monomer with different model chemistry yields ΦHOOH = 111.3°‒121.5° (Table 1) and the larger basis correctly predicts ΦHOOH. Given the computational cost of CCSD(T)/6-311++G(3df,3pd) calculations, the geometries of (H2O2)2, (H2O2)3, H2O2.H2O, H2O2.2H2O and 2H2O2.H2O clusters were optimized at the MP2 level and results show that H-bonds are shorter and binding energies are larger (more negative E and ECP) for the larger basis set. Geometries optimized at the MP2(full)/6-311++G(3df,3pd) level are reported and their atomic coordinates are given in the Supporting Information (SI). Average ΦHOOH and binding energies (E, ECP) calculated at the MP2(full)/6-311++G(d,p) level are reported for comparison.

ACS Paragon Plus Environment

9

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 47

3.1.1 Molecular Properties of Gaseous H2O2 and H2O Monomers

Geometry optimization of H2O2 yields one energy minima (1a) and transition-state cis (1b) and trans (1c) structures (Figure 1A, see SI for coordinates). Except for overestimating ΦHOOH at the MP2(full)/6-311++G(d,p) level, the ab initio calculations yield structural parameters for conformer 1a in good agreement with experiment (Table 1) but overestimate the dipole moment and vibration frequencies of 1a. Structures 1b and 1c possess the same rOH as 1a but rOO is 0.010 Å longer and θHOO is either 1.5° smaller (1b) or 4.5° larger (1c) than that of 1a (data not shown).

Table 1. Molecular properties of gaseous H2O2 Property a

MP2(full)/ 6-311++G(d,p)

MP2(full)/ 6-311++G(3df,3pd)

CCSD(T)/ 6-311++G(3df,3pd)

Additive model

(Ref)

rOH (Å)

0.964

0.963

0.964

0.963

0.96 (17)

rOO (Å)

1.449

1.442

1.452

1.454

1.47 (17)

Experiment

1.467 (10) rHH (Å)

2.426

2.368

2.377

2.452

θHOO (deg)

99.6

99.9

100.0

102.6

97

ΦHOOH (deg)

121.5

111.7

111.3

114.5

111.5 (16)

Dipole (D)

1.81

1.91

1.92

2.00

1.76 (22,42)

ν1 (symmetric O‒H stretch)

3851

3852

3823

3597

3599 (15)

ν2 (symmetric angle bending)

1459

1436

1431

1368

1394 (17)

ν3 (O‒O stretch)

922

939

906

877

880 (15)

ν4 (torsional oscillation)

394

410

402

371

371 (16)

ν5 (antisymmetric O‒H stretch)

3852

3852

3824

3609

3608 (15)

ν6 (antisymmetric angle bending)

1302

1334

1329

1299

1266 (15)

a

(17)

Vibrational frequencies in cm-1.

ACS Paragon Plus Environment

10

Page 11 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

A plot of energy vs. ΦHOOH reveals that 1b and 1c are 1.2 and 7.8 kcal/mol less stable than 1a at the MP2(full)/6-311++G(3df,3pd) level or 1.0 and 9.0 kcal/mol at the MP2(full)/6311++G(d,p) level (Figure 1B). These values bracket the experimental barrier for rotation around the OO bond to give the trans conformer (0.9‒1.1 kcal/mol)15,16,21 but the larger basis set better reproduces the experimental barrier for the cis conformer (7.0‒7.1 kcal/mol).16,21 This is also the case at the CCSD(T)/6-311++G(3df,3pd) level, which predicts that 1b and 1c are 1.2 and 7.5 kcal/mol less stable than 1a. Rotation around the OO bond also is accompanied by a dramatic change in the dipole moment (Figure 1C) and the larger basis set predicts lower dipole moments as ΦHOOH approaches 0°. Optimizing the geometry of H2O with MP2(full)/6-311++G(3df,3pd) reveals that its OH bond length (rOH = 0.958 Å) is slightly less than in H2O2. However, the HOH angle in H2O (104.1°) is greater than the HOO angle in H2O2, which is likely compressed by electrostatic attraction of H with the second O atom. To understand the strength of noncovalent interactions in the two molecules, we show in Figure 1D the electrostatic potential maps for H2O and H2O2. Overlap of the electron density from the O lone pairs of H2O results in a single Vs,min (−30.5 kcal/mol) centered on the O atom in the HOH plane while overlap between the Vs,min on each O atom of H2O2 gives a single Vs,min (−26.4 kcal/mol) located at the center of the O‒O bond. Its smaller Vs,min makes H2O2 a weaker H-bond acceptor than H2O. However, since the H atoms of H2O2 are more positive than those of H2O (Vs,max = 40.0 vs 35.7 kcal/mol), H2O2 is a stronger Hbond donor as reflected in the pKa values of 11.6 (H2O2) and 14.0 (H2O).68

ACS Paragon Plus Environment

11

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 47

Figure 1. Molecular properties of gaseous H2O2 and electrostatic maps of H2O and H2O2. (A) Optimized geometry at the MP2(full)/6-311++G(3df,3pd) level of the energy minimum (1a) and two transition-state structures (1b and 1c) of H2O2 (see SI for coordinates) and their ΦHOOH and point group. Variation in (B) energy and (C) dipole moment of H2O2 as function of ΦHOOH. (D) Electrostatic potential maps of H2O (top) and H2O2 (bottom) oriented to show the charge on the O (left) and H (right) atoms (see geometries below each surface). The most negative (Vs,min, red) and positive (Vs,max, blue) values of Vs(r) are around the O and H atoms, respectively.

ACS Paragon Plus Environment

12

Page 13 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

3.1.2 Structures and Binding Energies of (H2O2)n (n = 2, 3) Clusters

The energy-minima structures of the H2O2 dimers and trimers obtained from geometry optimization at the MP2(full)/6-311++G(3df,3pd) level are shown in Figure 2 (see SI for coordinates). Optimization at the MP2(full)/6-311++G(d,p) level yield similar geometries but Hbonds are up to 0.2 Å longer (structures not shown). The shorter H-bonds at the larger basis set indicate strengthening of the interactions as reflected in the 15‒25% larger (more negative) binding energies (ECP, Table 2). Table S1 provides BSSE-uncorrected binding energies (E) and the average ΦHOOH calculated at both basis sets. BSSEs calculated with the 6-311++G(3df,3pd) basis set are 20‒30% smaller than those calculated with the smaller basis set for H2O2 clusters. As for the H2O2 monomer, the smaller basis set systematically predicts larger average ΦHOOH values. The four energy-minima conformers found for (H2O2)2 (2pa‒2pd, Figure 2A) are stabilized by two H-bonds with an additional weak H-bond in 2pc. The six-membered cyclic structures of 2pa and 2pb are ~2 kcal/mol more stable than the five-membered cyclic structures of 2pc and 2pd because H-bonding is more linear in larger rings (θO‒H···O ~160° vs ~135°). The ECP values of −7.98 and −7.83 kcal/mol for 2pa and 2pb, which differ only in the trans vs cis orientation of their non-H-bonded H atoms, are in good agreement with previous ab initio calculations at the MP2 level with various basis sets (−7.37 to −8.34 kcal/mol for 2pa and −7.15 to −8.17 kcal/mol for 2pb).36,38,40 To the best of our knowledge, 2pc and 2pd have not been reported previously and the additional weak H-bond in 2pc stabilizes it by 0.5 kcal/mol relative to 2pd. Structure 2pe in Figure 2A was obtained from a constrained optimization to force formation of a single linear H-bond between the two H2O2 molecules for comparison with (H2O)2 (2w, Figure 2B). This H-bond is 0.05 Å shorter than that in 2w. However, the fact that

ACS Paragon Plus Environment

13

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 47

2pe is stabilized by two additional weak and bent (θO‒H···O = 72‒76°) H-bonds (Figure 2) suggests that 2pe and 2w exhibit similar stability (Table 2) probably because the more negative electrostatic potential on O of H2O is counterbalanced by the larger electropositive potential on H of H2O2 (Figure 1D).

Figure 2. Optimized structures of all stable conformers of (A) (H2O2)2 (2p), (C) (H2O2)3 (3p) and (B) (H2O)2 (2w) in the gas phase at the MP2(full)/6-311++G(3df,3pd) level. The conformations of each complex are presented in decreasing order of stability. H-bonding is indicated by dotted lines and numbers indicate H-bond distances in Å. The atomic coordinates are provided in the SI.

ACS Paragon Plus Environment

14

Page 15 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Table 2. Binding energies (kcal/mol) for the (H2O2)2 (2p), (H2O2)3 (3p) and (H2O)2 (2w) conformers Conformer a ECP b

EMM c

ECP d

EMM e

2pa

−6.83

−6.83

−7.98

−6.87

2pb

−6.61

−6.48

−7.83

−6.45

2pc

−5.25

−6.20

−6.44

−6.29

2pd

−5.27

−4.83

−5.91

−5.27

2pe f

−4.31

−5.40

−4.76

−5.81

2w

−4.49

−6.29

−4.62

−6.09

3pa

−15.87 −17.78 −18.84 −18.29

3pb

−14.72 −15.70 −17.38 −16.26

3pc

−14.35 −14.10 −17.26 −14.22

3pd

−14.38 −14.24 −17.12 −14.25

3pe

−13.80 −14.41 −16.91 −14.29

3pf

−14.03 −13.84 −16.70 −13.62

3pg

−14.26 −15.45 −16.59 −16.62

3ph

−13.70 −14.53 −16.40 −15.39

3pi

−13.85 −15.31 −16.01 −16.08

3pj

−12.53 −13.93 −15.66 −13.71

3pk

−12.76 −14.85 −15.37 −14.98

3pl

−12.15 −14.03 −15.18 −14.13

3pm

−12.61 −13.23 −14.98 −13.72

3pn

−10.91 −12.02 −13.31 −12.59

a

Structures are presented in Figure 2.

b

Calculated with MP2/6-311++G(d,p) and corrected for BSSE

c

Calculated using the additive model with MP2/6-311++G(d,p) optimized geometries

d

Calculated with MP2/6-311++G(3df,3pd) and corrected for BSSE

e

Calculated using the additive model with MP2/6-311++G(3df,3pd) optimized geometries

f

This structure was obtained by constraining the geometry so that a single H-bond forms

ACS Paragon Plus Environment

15

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 47

Compared to the maximum of five structures found in previous ab initio investigations,36,38,40 we report 14 stable conformers for (H2O2)3 that are stabilized by three to five H-bonds (Figure 2C). The engagement of O in accepting an H-bond weakens its ability to act as H-bond donor or as a second H-bond acceptor. This together with the geometry of the Hbond (length and angle) explains the observed trend in stability of the various structures. For example, only one O acts as both H-bond donor and acceptor in the global minimum conformer (3pa), while two O atoms act as bi-acceptor in the least stable conformer (3pn). Notably, many of the trimers (e.g., 3pa, 3pn) possess cage-like structures. 3.1.3 Structures and Binding Energies of H2O2.H2O, H2O2.2H2O and 2H2O2.H2O Clusters

The MP2(full)/6-311++G(3df,3pd)-optimized geometries of these binary clusters are reported in Figure 3 and are characterized by up to 0.4 Å shorter H-bonds compared to those optimized with the 6-311++G(d,p) basis set (structures not shown), which is reflected in larger ECP (Table 3) and E (Table S2) binding energies at the larger basis set. BSSEs calculated at the larger basis set are 30‒40% lower than those at the smaller basis set. Table S2 also shows that the smaller basis set systematically predicts larger average ΦHOOH. Three energy-minima structures are obtained for the H2O2.H2O pair in addition to pwc (Figure 3A), where H2O2 is constrained to donate a single H-bond to H2O. The shorter H-bond and the ~2 kcal/mol higher binding energy in pwc compared to pwd (Table 3) reveal that H2O2 is a better H-bond donor than H2O. This is consistent with the larger electropositive potential on H of H2O2 and the larger electronegative potential on O of H2O (Figure 1D). Both pwa and pwb are stabilized by two H-bonds and the free H atoms are oriented trans in pwa but cis in pwb. Previous ab initio studies report structures similar to pwa and pwd43–48 with ECP of −6.68 and −6.40 kcal/mol for pwa and −3.69 and −3.86 kcal/mol for pwd with MP2/6-

ACS Paragon Plus Environment

16

Page 17 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

311++G(2d,2p)(6d,10f)43 and MP2/aug-cc-pVTZ//MP2/aug-cc-pVDZ,46 respectively, values that bracket ours (Table 3). Of the six structures optimized for H2O2.2H2O (Figure 3B), the p2wa/p2wb, p2wc/p2wd and p2we/p2wf pairs differ only in the relative orientation of the free H of the H2O molecules. Conformers p2wa/p2wb are characterized by three H-bonds and a planar cyclic arrangement of the four O atoms. One peroxide O is exocyclic in p2wc/p2wd, and p2we/p2wf possess four Hbonds and a bent arrangement of the three molecules making these pairs are ~1.5 and 3.0 kcal/mol less stable than p2wa/p2wb. The H-bond distance (Hp···Ow < Hw···Ow < Hw···Op; p = peroxide and w = water) varies inversely with H-bond strength, which reflects the electrostatic potentials on the atoms of H2O2 and H2O (Figure 1D). Previous investigations on H2O2.2H2O reported structures similar to p2wa, p2wc and p2wd with ECP = −15.81 and −16.71, −14.61 and −15.35, and −14.60 and −15.22 kcal/mol with MP2/6-311++G(2d,2p)(6d,10f)43 and MP2/augcc-pVTZ//MP2/aug-cc-pVDZ level.46 A structure similar to p2wb was reported as the global minimum conformer with B3LYP/cc-pVTZ44 and structures similar to p2wa and p2we were reported with PBE/6-31+G(d)45 but not their binding energies.

ACS Paragon Plus Environment

17

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 47

Figure 3. Optimized geometries of the (A) H2O2.H2O (pw), (B) H2O2.2H2O (p2w) and (C) 2H2O2.H2O (2pw) gaseous complexes with MP2/6-311++G(3df,3pd). H-bonding is presented as dotted lines and the numbers (Å) indicate the H-bond distances. To distinguish H2O from H2O2, the O atoms of H2O molecules are labeled with W. The atomic coordinates are provided in the SI.

ACS Paragon Plus Environment

18

Page 19 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Table 3. Binding energies (kcal/mol) of the gas-phase H2O2.H2O (pw), H2O2.2H2O (p2w) and 2H2O2.H2O (2pw) conformers Conformer a ECP b

EMM c

ECP d

EMM e

pwa

−6.10

−6.52

−6.68

−6.37

pwb

−5.71

−5.91

−6.25

−5.32

pwc f

−5.51

−7.20

−5.92

−7.09

pwd

−3.55

−5.06

−3.79

−5.03

p2wa

−15.06 −17.68 −16.65 −17.54

p2wb

−15.06 −17.52 −16.61 −17.36

p2wc

−13.81 −16.91 −15.20 −17.13

p2wd

−13.74 −16.17 −15.09 −16.06

p2we

−12.57 −13.21 −14.00 −12.95

p2wf

−12.35 −12.77 −13.81 −12.32

2pwa

−15.40 −18.05 −17.77 −18.61

2pwb

−15.09 −17.29 −17.44 −17.65

2pwc

−14.89 −17.38 −16.81 −17.66

2pwd

−14.78 −16.92 −16.64 −17.22

2pwe

−14.47 −17.16 −16.55 −17.11

2pwf

−14.14 −16.07 −16.37 −16.18

2pwg

−13.70 −14.31 −15.82 −14.10

2pwh

−13.24 −15.11 −15.80 −14.85

2pwi

−13.54 −14.08 −15.64 −13.36

2pwj

−13.15 −13.51 −15.12 −12.83

2pwk

−11.86 −12.25 −13.45 −12.46

2pwl

−11.55 −11.70 −13.15 −11.73

2pwm

−11.67 −12.01 −12.91 −11.77

2pwn

−11.33 −11.76 −12.48 −11.54

2pwo

−11.07 −11.68 −12.19 −10.90

a

Structures are presented in Figure 3.

b‒f

See footnotes of Table 2

ACS Paragon Plus Environment

19

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 47

Optimization of 2H2O2.H2O gives 15 stable structures with 3 or 4 H-bonds (Figure 3C). In all structures, H2O acts as both an H-bond donor and acceptor and the most stable structures are those that form strong H-bonds with H2O as seen in 2pwa and 2pwb, which adopts a cage conformation and differ mainly in the orientation of the free H atom of H2O. Three H-bonds stabilize the cyclic arrangements of the five O atoms in 2pwc. A single O atom is exocyclic in 2pwd and 2pwe. H-bonded atoms in 2pwf and 2pwh occupy perpendicular planes. Conformers 2pwi‒2pwo have bent arrangements of the molecules with structures 2pwi/2pwj, 2pwk/2pwl and 2pwn/2pwo being different in the position of the free H atoms. Compared to the equilibrium ΦHOOH of 111.7° for the H2O2 monomer, the average value is 102.6°‒114.1° for (H2O2)2, 98.4°‒145.5° for (H2O2)3, 111.0°‒114.8° for H2O2.H2O, 98.2°‒ 114.4° for H2O2.2H2O and 83.0°‒117.5° for 2H2O2.H2O (Tables S1, S2). Notably, the range of equilibrium ΦHOOH values suggests that the molecule flexibility increases with cluster size. 3.2 Optimized Force Field for H2O2 The optimized FF parameters for H2O2 are listed in Table 4. The additive potential model yields structural parameters and dipole moment of the equilibrium gas-phase H2O2 monomer in good agreement with experiment (Table 1). It predicts that the trans and cis H2O2 conformers are 1.6 and 9.7 kcal/mol less stable than the global minimum skewed structure (Figure 1A,B), which are in fair agreement with the experimental differences of 0.9‒1.1 kcal/mol15,16,21 and 7.0‒7.1 kcal/mol,16,21 respectively. It yields dipole moments as function of ΦHOOH in agreement with MP2/6-311++G(d,p) (Figure 1C), and also predicts a smaller θHOO angle for the trans conformer (101.9°) but larger angle for the cis conformer (106.1°) compared to 102.6° of the equilibrium structure of H2O2.

ACS Paragon Plus Environment

20

Page 21 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Table 4. Parameters of the optimized H2O2 model Parameter

STU (Å)

Value 0.963 2

521.0

Kb(OH) (kcal/mol/Å )

STT (Å)

1.442 2

285.5

Kb(OO) (kcal/mol/Å )

VUTT (°)

99.92 2

Kθ(HOO) (kcal/mol/rad )

60.4

δ (°)

0.0

n

2

KΦ(HOOH) (kcal/mol)

2.02

WT (e)

–0.41

XZ[\ (kcal/mol) Y

0.04600

XZ[\ (kcal/mol) ^

0.20384

WU (e)

]Z[\ (Å) Y ]Z[\ (Å) ^

0.41

0.44900

3.34846

Despite not being optimized to reproduce the ab initio properties of the clusters, the model together with the TIP3P water model yield binding energies for the various conformers of the (H2O2)2, (H2O2)3, H2O2.H2O, H2O2.2H2O and 2H2O2.H2O clusters in good agreement with ab initio ECP values (Tables 2, 3). The average absolute percent error in binding energies calculated on the MP2(full)/6-311++G(d,p) geometries of the pure and binary clusters is 8.6% and 12.1% and is 9.8% and 9.2% for binding energies calculated using geometries optimized with the larger basis set.

ACS Paragon Plus Environment

21

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 22 of 47

3.3 Thermodynamic Properties of Liquid H2O2 The H2O2 model is optimized to reproduce the density and vaporization enthalpy of liquid H2O2 at 0 ºC and 0.37 mm Hg. The molar enthalpy of vaporization (∆Hvap, i.e., the quantity of heat required to convert one mole of a substance from a liquid to a gas at constant temperature) is given by:69 Δ`a- = 4b − Δc = 4b − P〈c〉 − 〈c〉 Q 4 R is the gas constant and T is the absolute temperature. The average potential energy per mole in the liquid phase 〈c〉 is obtained from the NpT MD simulations of liquid H2O2, and that in the gas

phase 〈c〉 is obtained from vacuum MD simulations of 500 H2O2 molecules in the NVT ensemble using infinite nonbonded cutoffs at a volume calculated using the ideal gas law. The model yields an average potential energy for gaseous and liquid H2O2 of 23.70 and 11.73 kcal/mol, which gives (eq. 4) ∆Hvap = 12.51 kcal/mol (Table 5). Yu and Yang40 neglected the potential energy of gaseous H2O2 and calculated a value of ∆Hvap = 12.82 kcal/mol using the ABEEM/MM model. Simulation of 500 H2O2 molecules at 0 ºC and 0.37 mm Hg yields an average molecular volume of 38.52 Å3 and hence a density of 1.466 ± 0.001 g/cm3. In comparison, values of 28.85 Å3 and 1.037 g/cm3 are obtained for simulating 500 H2O molecules under the same conditions. Replacing H in H2O by OH in H2O2 increases the mass by 90% but the molecular volume by only 34%. The relatively small increase in volume accounts for the larger density of H2O2 and suggests stronger association between H2O2 vs H2O molecules. The self-diffusion coefficient of an H2O2 molecule in the liquid state is obtained from the long-time limit of the mean-square displacement of its oxygen atoms:69 gh

i

x

1 1

= lim 〈 stu,0 q − tu,0 0w 〉 5 m→o 6q r 0y7

ACS Paragon Plus Environment

22

Page 23 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

The resulting diffusion coefficient is corrected for system-size dependence using the formula:60 g = gh

i

+

2.837297  b 6 6€‚

where  is Boltzmann constant,  the shear viscosity of the solvent and L the average length of the cubic simulation box. Simulation of 250, 500 and 1000 H2O2 molecules at 0 °C results in values of the self-diffusion coefficient DPBC of 0.71, 0.72, and 0.75 × 10−9 m2/s, respectively. Using the experimental viscosity of liquid H2O2 at 0 °C (1.819 cp)25 results in D = 0.84 × 10−9 m2/s, which is closer to the experimental value (1.0 × 10−9 m2/s)70 than D from the ABEEM/MM model (1.62 × 10−9 m2/s).40 The dielectric constant of the liquid (ƒ) is calculated from:69 ƒ = 1 +

4€ 〈† 〉 − 〈†〉  7 3〈„〉… b

Here M is the total dipole moment of the box. The convergence of ƒ is monitored by plotting its value as a function of time, and is averaged over the last 15 ns of 17-ns simulations. Our model overestimates the dielectric constant of liquid H2O2 (Table 5) probably because it overestimates the dipole moment of the gaseous monomer (Table 1). The specific heat capacity at constant pressure (Cp) is calculated from five simulations at 278.15, 283.15, 288.15, 293.15 and 298.15 K with a constant pressure of 2.25 mm Hg (the vapor pressure of liquid H2O2 at 298.15 K). The total energy 〈m‡m 〉 and volume 〈„〉 are averaged over time and a linear fit of 〈m‡m 〉 + ˆ〈„〉 vs. T is used to calculate Cp.69 ‰Š =

1 ‹〈m‡m 〉 + ˆ〈„〉 3 6 8 ‹b r Œy Ž.7 , Šy . ‘’ ““ F”

The heat capacity derived from our model (19.4 cal mol–1 K–1 at 25 °C) agrees well with the experimental value (21.35 cal mol–1 K–1).27

ACS Paragon Plus Environment

23

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 24 of 47

Table 5. Calculated and experimental properties of liquid H2O2 a,b Method

ρ (g/cm3)

∆Hvap (kcal/mol)

D × 10–9 m2/s

Model

1.466 ± 0.001

12.51 ± 0.01

0.84 ± 0.1

110 ± 5

19.4 ± 0.1

Expt.

1.463 (24) 1.471 (31)

12.59 (26) 12.62 (27)

1.0 (70)

74.6 (68)

21.35 (27; at 25 °C)

ε

Cp (cal mol–1 K–1)

a

All properties were calculated at 0 °C and 0.37 mm Hg except for Cp, which is calculated at 25 °C and 2.25 mm Hg. b References for the experimental data are in parentheses.

3.4 Structure of Liquid H2O2 and the Equilibrium Geometry of H2O2 in the Liquid State The microscopic structure of liquid H2O2 is analyzed from the •H– H– , •H– F–  and •F– F–  radial distribution functions (RDFs) derived from the MD simulations (Figure 4A‒C). The •H–H–  function calculated from our additive potential model shows an asymmetric peak with a maximum at 2.84 Å, a shoulder at ~3.4 Å and a minimum at 4.05 Å, together with weak peaks at 4.3 and 5.5 Å. Integration to 4.05 Å yields 14.7 Op (1 intramolecular + 14 intermolecular Op) per given Op, indicating that each H2O2 molecule is surrounded by a minimum of 7 H2O2 molecules in its first solvation shell. The ABEEM/MM model40 predicts a more symmetric and sharp peak with a maximum at 2.74 Å in addition to three other peaks centered at 3.8, 5.2 and 7.6 Å (Figure 4A). The predicted maximum in the closest intermolecular O···O distance from both models (2.85 Å and 2.74 Å) brackets the experimental values (2.758 − 2.799 Å) from X-ray and neutron diffraction measurements at −20 to −163 °C.11–13 The •H– F–  function displays a sharp maximum centered at 1.95 Å that corresponds to both intramolecular and intermolecular O···H contacts and is in good agreement with the values (1.786 − 1.98 Å) from the low-temperature X-ray and Neutron diffraction studies.12,13 Integrating the •H– F–  function up to 2.45 Å yields 3.7 H atoms, corresponding to 2 intramolecular H atoms plus an average of 1.7 H-bonds per Op. By comparison, the ABEEM model predicts

ACS Paragon Plus Environment

24

Page 25 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

slightly a more intense peak centered at 2.01, and both models predict additional weaker and shoulders peaks between 3−6 Å that correspond to second and third shell H atoms. Calculation of the •H— F—  function (Figure 4B) allows us to compare H-bonding in

liquid H2O and liquid H2O2. The intermolecular O···H peak in •H— F—  at 1.85 Å is 0.10 Å shorter than that in liquid H2O2, indicating stronger H-bonding in liquid water (similar to gasphase). Consistent with our finding, a vibrational study suggests that intermolecular H-bonds in crystalline (or frozen) H2O2 are 10‒15% weaker than in frozen (or crystalline) H2O.18 It should be noted that the •H– F–  function is more intense due to intramolecular O···H interactions, which are found at 2.37 Å in the optimized isolated monomer (Table 1). Integration of the •H— F—  function up to 2.45 Å yields 4.0 H atoms, indicating 2 H-bonds/Ow. While H-bonding is slightly stronger in liquid H2O, fewer H-bonds are formed per molecule (4.0) than in liquid H2O2 (6.8). Thus, association between H2O2 molecules is larger than between H2O molecules. The first peak centered at 2.36 Å in the •F– F–  function corresponds to intramolecular H···H distances while less intense peaks at 2.69, 4.9 and 7.3 Å are due to intermolecular H···H pairs. In comparison, the ABEEM/MM model displays more intense peaks centered at 2.41, 2.86, 4.5 and 6.9 Å.40 To investigate the equilibrium geometry of the H2O2 molecule in the liquid state, we examined the distribution of rOO, θHOO and ϕHOOH during the simulations with the OH bond at its equilibrium value of 0.963 Å (Figure 4D‒F). Notably, rOO (1.45 Å) and θHOO (103º) maintain their equilibrium values found in the gaseous monomer (Table 1), but adopt values in the range 1.3−1.6 Å and 80−120º. In contrast, the equilibrium value of ϕHOOH (93º) is 22º smaller in the liquid vs gaseous monomer. This is likely because the molecule is more polar at small ϕHOOH (Figure 1C), which would allow for stronger H-bonds between the molecules. It should be noted

ACS Paragon Plus Environment

25

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 26 of 47

that values between 0º and 180º are found for ϕHOOH (not shown in Figure 4F) but given the much lower stability of cis (ϕHOOH = 0º) vs. trans H2O2 (ϕHOOH = 180º), only 270 structures possess ϕHOOH ≤ 10º compared to 2200 structures with ϕHOOH ≥ 170º of the 1.5×107 structures examined during the simulation. Finally, our additive potential model predicts equilibrium bond length and angles that match the X-ray and neutron diffraction values for H2O2 at −192 to −20 ºC (avg rOO = 1.453–1.461 Å; avg θHOO = 97−102.7º; and avg ϕHOOH = 90.2−99.0º).11–15

Figure 4. Structure of liquid H2O2 and equilibrium geometry of H2O2 in the liquid state. (A‒C) •H– H– ,

•H– F–  and •F– F–  radial distribution functions calculated from 15-ns MD simulation of 500 H2O2

molecules at 0 °C and 0.37 mm Hg using our optimized additive model for H2O2 (solid lines) vs the

ABEEM/MM-derived functions (dotted lines).40 The •H— F—  function calculated for liquid H2O from

15-ns MD simulation of 500 TIP3P molecules at 0 °C and 0.37 mm Hg is also shown in panel B (red line). Note that peaks at rOO = 1.45 Å and rOH 0.96 Å due to intramolecular covalent bonds are omitted

ACS Paragon Plus Environment

26

Page 27 of 47

from panels A and B. (D‒F) Probability distribution during the simulations of rOO, θHOO and ϕHOOH, respectively, for H2O2 in the liquid state with the OH bond constrained at its equilibrium value of 0.963 Å.

3.5 Liquid H2O2 at Various p and T The transferability of the H2O2 model to different thermodynamic conditions is validated by calculating the density of the liquid in the temperature range −20.0 to 200.0 °C at the corresponding vapor pressure (0.069‒3132.41 mm Hg).23 Since H2O2 freezes at −0.46 °C,32 it is at the supercooled state during the simulation at −20.0 °C. The list of temperatures and pressures used for the simulations and the corresponding calculated densities are shown in Table S3. Figure 5 plots the calculated as well as the experimental density24 of liquid H2O2 vs. temperature and shows that calculated densities are in excellent agreement with available experimental data.

1.5

Additive model Experiment

Density (g/cm3)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

1.4

1.3

1.2 200

250

300

350

400

450

500

T (K)

Figure 5. Comparison of calculated (black squares) and experimental24 (red circles) density of liquid H2O2 vs temperature. Data from Table S3 of the SI.

3.6 Solvation of H2O2 in Liquid Water Our model for H2O2 and the TIP3P model for H2O yield H2O2 hydration free energy components: ∆KO = −9.3 kcal/mol, ∆K(- = −5.2 kcal/mol and ∆K- = 6.4 kcal/mol. Inputting

ACS Paragon Plus Environment

27

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 28 of 47

these values into eq 3 gives ∆K) = −8.1 kcal/mol, in very good agreement with experiment (−8.5 to −9.0 kcal/mol).50 We next analyze the hydration structure of H2O2 using the gXY(r) RDFs where X = Op, Hp; and Y = Ow, Hw (Figure 6A). The •D˜D™  function displays overlapping peaks centered at 2.76 and 3.35 Å. The former is in excellent agreement with X-ray measurements on crystals of H2O2.2H2O (2.69 and 2.76 Å).28 Integration up to the minima at 3.15 and 4.15 Å yields 3.0 and 9.4 water molecules in the micro49 and entire first hydration shells of H2O2, respectively. In comparison, a Monte Carlo simulation on a single H2O2 molecule in 300 water molecules found twice as many water molecules in the micro (6) and entire first solvation (17) shells of H2O249 and combined quantum mechanics and molecular mechanics (QM/MM) MD simulations reveals 6 water molecules in the entire first solvation shell of H2O2.47 The •D˜C™  function is characterized by a low intensity peak centered at 1.95 Å, a minimum at 2.45 Å, a second peak at 3.30 Å with a shoulder at 4.03 Å and a minimum at 4.9 Å (Figure 6A). Integration up to 2.45 Å yields 1.5 Hw per Op, consistent with the calculated 3 water molecules in close vicinity with H2O2. The peak at 1.95 Å is within the range (1.85 - 2.1 Å) reported from previous MD,46 Monte Carlo49 and QM/MM MD simulations.47,71 Also, previous investigations found that each Op acts as H-bond acceptor to an average of 1.047 or 1.449 H2O molecules. The •C˜ D™  function possesses a sharp peak centered at 1.75 Å, a minimum at 2.40 Å followed by a second peak at 3.65 Å and a minimum at 5.0 Å (Figure 6A). Integration up to 2.40 Å yields a single Ow per Hp. In comparison, from previous simulations the first peak for •C˜D™  is centered at 1.75,46 1.81,49 1.6671 and 1.8 Å.47

ACS Paragon Plus Environment

28

Page 29 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

The •C˜ C™  function displays two peaks centered at 2.55 and 4.4 Å with minima at 3.15 and 5.6 Å, respectively. Integration of the function to 3.15 Å yields 6.6 Hw, consistent with 3 nearest water neighbors. Thus, from the four RDFs examined, we conclude that the first hydration shell of H2O2 is composed of around nine H2O molecules, with two acting as H-bond acceptors. The remaining seven water molecules act mainly as transient H-bond donors, three at a time. 3.7 Solvation of H2O in Liquid Hydrogen Peroxide We also investigated the solvation structure of H2O in liquid hydrogen peroxide (Figure 6B). Integration up to 4.0 Å of •D™ D˜ , which displays a maximum at 2.84 Å with a shoulder at

3.24 Å, reveals 14 H2O2 around each H2O. Integration up to 2.40 Å of •D™ C˜ , which has the

same characteristics as •C˜D™  (Figure 6A), yields 1.9 Hp atoms for each Ow. Also, similar to •D˜C™  (Figure 6A), •C™ D˜  displays a peak at 1.95 Å and integration up to 2.45 Å gives

1.4 Op atoms per Hw. The •C™ C˜  function displays two peaks centered at 2.55 and 4.4 Å with minima at 3.35 and 5.6 Å. In summary, H2O is solvated by 14 peroxide molecules, with up to four acting as H-bond donors and the remaining acting as transient H-bond acceptors, three at a time. The relative stability of the various types of H-bonds in aqueous solution can be determined by comparing H-bond distances and the number of H-bond calculated for the first intermolecular peak in the •Dš C›  functions (x, y = w, p; Figure 4B and 6). The •D™ C˜  and

•D™ C™  functions yield the largest number of H-bonds (~2Hp/w per Ow) but the 0.1 Å shorter H-bond distance in •D™ C˜  indicates that Op‒Hp···Ow H-bonding is the strongest in solution.

While •D˜C˜  and •D˜C™  possess similar H-bond distances (1.95 Å), the larger number of

ACS Paragon Plus Environment

29

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 30 of 47

Op‒Hp···Op H-bonds indicates that Op‒Hp···Ow H-bonds are the weakest. Thus, as in the gasphase, the strength of H-bonding in aqueous H2O2 solution is as follows: Op‒Hp···Ow > Ow‒ Hw···Ow ≥ Op‒Hp···Op > Ow‒Hw···Op.

Figure 6. Solvation of H2O2 in liquid water and of H2O in liquid hydrogen peroxide. Radial distribution functions gXY(r) between: (A) a single H2O2 molecule (X = Op, Hp; and Y = Ow, Hw) solvated in 500 water molecules and (B) a single H2O molecule (X = Ow, Hw; and Y = Op, Hp) solvated in 500 hydrogen peroxide molecules at T = 25.0 ºC and p = 760.00 mm Hg. See caption to Figure 4 for details.

3.8 Density and Structure of Aqueous H2O2 H2O2 is completely miscible with H2O at all concentrations32 allowing solution density to be measured over a wide range of compositions.24,30,31 To further test the performance of the H2O2 model, we calculated solution density at an H2O2 mole fraction (EFG HG ) of 0.0 to 1.0. Table S4 lists the simulated T and p values and the calculated densities. Plots of density vs. EFG HG expose

ACS Paragon Plus Environment

30

Page 31 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

the excellent agreement between theory and experiment30 (Figure 7), with the maximum error (3.7%) at EFG HG = 0.0 at 0 ºC, due to an overestimation of the density of liquid water. Drawing a

straight line on the plots between EFG HG = 0.0 and 1.0 reveals that both the calculated and

experimental densities are higher than expected from ideal mixing. This we attribute to the replacement of the Op‒Hp···Op and Ow‒Hw···Ow H-bonds in the pure liquids by the shorter and stronger Op‒Hp···Ow H-bonds in the binary mixtures. A reduction in volume ensues although mixing also allows weaker Ow‒Hw···Op H-bonds to form. However, the influence of the stronger Op‒Hp···Ow H-bonds dominates since adding H2O2 to H2O is exothermic27 and the partial vapor pressures of H2O2-H2O mixtures are lower than expected for ideal solutions.29 Additionally, shock-frozen H2O2-H2O mixtures exhibit O‒H stretching vibrations that are downshifted in the dihydrate (H2O2.2H2O) compared to the pure components, indicating stronger H-bonding between H2O and H2O2 than between like molecules.7

Figure 7. Comparison of calculated and experimental Densities of H2O2/H2O mixtures. Densities were calculated from the MD simulations of 500 (H2O2 + H2O) molecules (solid lines). Experimental densities30 (dotted lines) were measured at 0 ºC (black), 25 ºC (red), 50 ºC (blue) and 96 ºC (green). The calculated densities are listed in Table S4.

ACS Paragon Plus Environment

31

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 32 of 47

Changes in the structure of aqueous H2O2 with temperature were examined for EFG HG = 0.0 to 1.0. The results show that the number of H-bonds of each type decreases with increasing temperature and H-bond-donor concentration (Figure 8). Also, at any given EFG HG , the number of Op‒Hp···Ow or Ow‒Hw···Ow H-bonds per Ow is larger than the number of Op‒Hp···Op or Ow‒ Hw···Op per Op. The data also show that more Op‒Hp···Ow or Op‒Hp···Op H-bonds are present at EFG HG than Ow‒Hw···Ow or Ow‒Hw···Op H-bonds at (1−EFG HG ). These results reveal that the relative strength of the different types of H-bonds remains the same, irrespective of the mixture composition, which is consistent with the miscibility of the two compounds over the entire EFG HG range. The conformational flexibility of H2O2 in water (50 H2O2 + 450 H2O) is investigated by monitoring the dihedral angle during the 15-ns MD simulations (Figure S1). Similar to that liquid H2O2 at 0 ᵒC, the average equilibrium value of ΦHOOH is 93ᵒ from 0 − 96 ᵒC, which is in good agreement with QM/MM MD simulations of a single H2O2 molecule in water (ΦHOOH = 90 ‒ 95ᵒ).47,48 The decrease in ΦHOOH compared to its equilibrium gas-phase value (114.5ᵒ, Table 1) allows for stronger intermolecular H-bonding. The similar behavior of H2O2 in its pure and aqueous solutions is likely due to the comparable polarity of both compounds. Figure S1 also shows that the probability distribution of ΦHOOH expands as the temperature increases, revealing that H2O2 becomes conformationally more flexible. This also is evidenced from the increase between 0 to 96 ᵒC in the number (24 to 92) of cis (ΦHOOH < 10ᵒ) and the number (408 to 2421) of trans (ΦHOOH < 170ᵒ) structures out of the 1.5×106 structures analyzed at each temperature.

ACS Paragon Plus Environment

32

Page 33 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Figure 8. The number of H-bonds per O atom of H2O or H2O2 vs. EFG HG . The maximum and minimum of

the first intermolecular peak in the RDFs, •H— F– , •H— F— , •H– F–  and •H– F— , are similar at the four temperatures and the number of H-bonds are calculated by integrating the first function up to 2.40 Å and other functions up to 2.45 Å. Table S5 lists the calculated number of H-bonds.

ACS Paragon Plus Environment

33

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 34 of 47

4. Conclusions Gaseous H2O2 is characterized by a global minimum conformer with skewed geometry (ΦHOOH = 111.5º) plus transition-state cis (ΦHOOH = 0.0º) and trans (ΦHOOH = 180.0º) conformers. Geometry optimization of the monomer with MP2(full)/6-311++G(d,p) overestimates ΦHOOH of the equilibrium conformer and the energy barrier for its conversion to the cis conformer but use of the larger basis set, 6-311++G(3df,3pd), corrects for these discrepancies. Molecular electrostatic potentials predict that, compared to H2O, H2O2 is a better H-bond donor but a poorer H-bond acceptor. Also, internal rotation allows H2O2 to behave differently from water. For example, the cis conformer should stabilize strong interactions with anions while the zero dipole moment of the trans conformer should stabilize H2O2 in nonpolar environments that exclude polar H2O. Our additive potential model for H2O2 is optimized based on the properties of the gaseous monomer determined with MP2(full)/6-311++G(3df,3pd), as well as the experimental vibrational frequencies of isolated H2O2 and the density and enthalpy of vaporization of liquid H2O2 at 0 ºC. The model yields energy barriers for the cis and trans conformers in fair agreement with experiment. It yields binding energies for (H2O2)m clusters (m = 2, 3) in very good agreement with those calculated with MP2(full)/6-311++G(3df,3pd), and predicts the experimental structure, self-diffusion coefficient and heat capacity for liquid H2O2 at 0 ºC. Furthermore, it was used to simulate liquid H2O2 up to 200 ºC and analysis of the structure at 0 ºC shows that each H2O2 molecule forms ~7 H-bonds, by donating and accepting 3.4 H-bonds. Our H2O2 model together with the TIP3P water model reproduce the ab initio-calculated binding energies of H2O2.(H2O)n and H2O.(H2O2)n clusters (n = 1, 2) and yield hydration free energy of H2O2 and densities for aqueous H2O2 solutions in agreement with experiment. The ab initio gas-phase data as well as the results from the MD simulations of the aqueous solutions

ACS Paragon Plus Environment

34

Page 35 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

indicate that H-bonds between H2O2 and H2O are stronger or weaker than those between H2O or H2O2 molecules. The net effect, however, is stronger interaction between H2O2 and H2O in aqueous solutions than between like molecules in the pure liquids. This explains why adding H2O2 to H2O is exothermic and why the binary mixtures have lower vapor pressure than the ideal mixtures. H2O2 can be viewed as the product of substituting a single H atom in water by an OH group. Thus, O‒H···O σ-type H-bonding is the main packing force in both liquid H2O2 and liquid H2O but the two liquids display quite different physical properties. Our results indicate that H-bonding in liquid H2O is slightly stronger than in liquid H2O2, but the larger number of Hbonds per H2O2 (7) vs H2O (4) indicates that molecules in liquid H2O2 are more associated. This together with the larger molar mass of H2O2 (34.0147 g/mol) compared to H2O (18.0153 g/mol) results in the observed higher density, viscosity, boiling point and heat of vaporization of liquid H2O2.72 Since our simple model for H2O2 accurately reproduces the experimental properties of pure H2O2 liquid under various thermodynamic conditions, it can be reliably used as an alternative to experiment for investigating pure or concentrated H2O2 solutions, which are hazardous. Our model also will serve to simulate H2O2 in different solvents to establish the response of H2O2 structure and flexibility to solvent polarity. To validate our model for biological applications, we have tested it, together with the CHARMM36 FF, for binding energies of H2O2 with model compounds for the neutral and charged amino acid side chains. Results show that, except for the sulfur-containing models, the FF reproduces the ab initio-calculated binding energies for the complexes (manuscript in preparation). Following calibration for sulfur-H2O2 interactions, our model will be implemented

ACS Paragon Plus Environment

35

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 36 of 47

to expose differences in H2O2 vs. H2O binding and diffusion through proteins. This will provide a better understanding at the atomic level of the biological roles of H2O2, a critical signaling molecule in cells.

ASSOCIATED CONTENT Supporting Information. Tables showing the BSSE-uncorrected binding energies (E) and average ΦHOOH calculated at the MP2(full)/6-311++G(d,p) and at the MP2(full)/6-311++G(3df,3pd) levels for the (H2O2)m, H2O2.(H2O)n and H2O.(H2O2)n clusters (m = 2, 3 and n = 1, 2), Tables listing the calculated densities of pure and aqueous H2O2 solutions and the number of the four different H-bonds in aqueous H2O2, a figure showing the probability distribution of ΦHOOH during the simulation of χCG DG = 0.1 aqueous H2O2 solution at different temperatures and the atomic coordinates of the structures presented in Figures 1‒3. This information is available free of charge via the Internet at http://pubs.acs.org.

AUTHOR INFORMATION Corresponding Author * E-mail: [email protected] and [email protected] Tel.: +1-514-848-2424, extension 5835; Fax: +1-514-848-2868 * E-mail: [email protected]. Tel.: +1-514-848-2424, extension 3338; Fax: +1-514-848-2868

Notes The authors declare no competing financial interest

ACS Paragon Plus Environment

36

Page 37 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation



On leave from Department of Chemistry, Faculty of Science, Assiut University, Assiut 71516,

Egypt.

ACKNOWLEDGMENT The authors thank Professor Guillaume Lamoureux (Concordia) for his valuable comments. This work was supported by a research grant from the Natural Sciences and Engineering Research Council (NSERC) of Canada awarded to AME. Computations were enabled by the support provided by the Centre for Research in Molecular Modeling (CERMM) at Concordia and by Compute Canada (www.computecanada.ca).

References (1)

Thénard, L. J. Observations Sur Des Combinaisons Nouvelles Entre L’oxigène et Divers Acides. Ann. Chim. Phys. 1818, 8, 306–312.

(2)

Eul, W.; Moeller, A.; Steiner, N. Hydrogen Peroxide. Kirk-Othmer Encycl. Chem. Technol. 2001, 13, 1–58.

(3)

Marinho, H. S.; Real, C.; Cyrne, L.; Soares, H.; Antunes, F. Hydrogen Peroxide Sensing, Signaling and Regulation of Transcription Factors. Redox Biol. 2014, 2, 535–562.

(4)

Veal, E. A.; Day, A. M.; Morgan, B. A. Hydrogen Peroxide Sensing and Signaling. Mol. Cell 2007, 26, 1–14.

(5)

Fay, A. J.; Qian, X.; Jan, Y. N.; Jan, L. Y. SK Channels Mediate NADPH OxidaseIndependent Reactive Oxygen Species Production and Apoptosis in Granulocytes. Proc. Natl. Acad. Sci. U. S. A. 2006, 103, 17548–17553.

ACS Paragon Plus Environment

37

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(6)

Page 38 of 47

Dence, C. W.; Reeve, D. W. Pulp Bleaching: Principles and Practice; Tappi Press: Atlanta, 1996.

(7)

Albers, P. W.; Glenneberg, J.; Refson, K.; Parker, S. F. Inelastic Incoherent Neutron Scattering Study of the Molecular Properties of Pure Hydrogen Peroxide and Its Water Mixtures of Different Concentration. J. Chem. Phys. 2014, 140, 164504.

(8)

Matheson, G. L.; Maass, O. The Properties of Pure Hydrogen Peroxide VI. J. Am. Chem. Soc. 1929, 51, 674–687.

(9)

Niethammer, P.; Grabher, C.; Look, A. T.; Mitchison, T. J. A Tissue-Scale Gradient of Hydrogen Peroxide Mediates Rapid Wound Detection in Zebrafish. Nature 2009, 459, 996–999.

(10)

Giguère, P. A.; Schomaker, V. An Electron Diffraction Study of Hydrogen Peroxide and Hydrazine. J. Am. Chem. Soc. 1943, 65, 2025–2029.

(11)

Abrahams, S. C.; Collin, R. L.; Lipscomb, W. N. The Crystal Structure of Hydrogen Peroxide. Acta Crystallogr. 1951, 4, 15–20.

(12)

Busing, W. R.; Levy, H. A. Crystal and Molecular Structure of Hydrogen Peroxide: A Neutron-Diffraction Study. J. Chem. Phys. 1965, 42, 3054–3059.

(13)

Savariault, J. M.; Lehmann, M. S. Experimental Determination of the Deformation Electron Density in Hydrogen Peroxide by Combination of X-Ray and Neutron Diffraction Measurements. J. Am. Chem. Soc. 1980, 102, 1298–1303.

(14)

Fritchie Jr. C. J., Mcmullan, R. K. Neutron Diffraction Study of the 1 : 1 Urea : Hydrogen Peroxide Complex at 81 K. Acta Crystallogr. 1981, B37, 1086–1091.

ACS Paragon Plus Environment

38

Page 39 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

(15)

Redington, R. L.; Olson, W. B.; Cross, P. C. Studies of Hydrogen Peroxide: The Infrared Spectrum and the Internal Rotation Problem. J. Chem. Phys. 1962, 36, 1311–1326.

(16)

Hunt, R. H.; Leacock, R. A.; Peters, C. W.; Hecht, K. T. Internal-Rotation in Hydrogen Peroxide: The Far-Infrared Spectrum and the Determination of the Hindering Potential. J. Chem. Phys. 1965, 42, 1931–1946.

(17)

Giguére, P. A.; Srinivasan, T. K. K. A Raman Study of H2O2 and D2O2 Vapor. J. Raman Spectrosc. 1974, 2, 125–132.

(18)

Arnau, J. L.; Giguère, P. A.; Abe, M.; Taylor, R. C. Vibrational Spectra and Normal Coordinate Analysis of Crystalline H2O2, D2O2, and HDO2. Spectrochim. Acta A 1974, 30, 777–796.

(19)

Khachkuruzov, G. A.; Przhevalskii, I. N. Molecular Constants of Hydrogen Peroxide. IV. Structural Parameters. Opt. Spectrosc. 1974, 36, 172–174.

(20)

Khachkuruzov, G. A.; Przhevalskii, I. N. Molecular Constants of Hydrogen Peroxide. 7: Parameters of the Potential Energy of Internal Rotation. Opt. Spectrosc. 1978, 44, 112– 114.

(21)

Koput, J. On the r0* Structure and the Torsional Potential Function of Hydrogen Peroxide. J. Mol. Spectrosc. 1986, 115, 438–441.

(22)

Perrin, A.; Flaud, J.-M.; Camy-Peyret, C.; Schermaul, R.; Winnewisser, M.; Mandin, J.Y.; Dana, V.; Badaoui, M.; Koput, J. Line Intensities in the Far-Infrared Spectrum of H2O2. J. Mol. Spectrosc. 1996, 176, 287–296.

(23)

Maass, O.; Hiebert, P. G. The Properties of Pure Hydrogen Peroxide V: Vapor Pressure. J.

ACS Paragon Plus Environment

39

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 40 of 47

Am. Chem. Soc. 1924, 46, 2693–2700. (24)

Maass, O.; Hatcher, W. H. The Properties of Pure Hydrogen Peroxide I. J. Am. Chem. Soc. 1920, 42, 2548–2569.

(25)

Phibbs, M. K.; Giguére, P. A. Hydrogen Peroxide and Its Analogues. I. Density, Refractive Index, Viscosity, and Surface Tension of Deuterium Peroxide-Deuterium Oxide Solutions. Can. J. Chem. 1951, 29, 173–181.

(26)

Foley, W. T.; Giguére, P. A. Hydrogen Peroxide and Its Analogues: IV. Some Thermal Properties of Hydrogen Peroxide. Can. J. Chem. 1951, 29, 895–903.

(27)

Giguére, P. A.; Morissette, B. G.; Olmos, A. W.; Knop, O. Hydrogen Peroxide and Its Analogues: VII. Calorimetric Properties of the Systems H2O-H2O2 and D2O-D2O2. Can. J. Chem. 1955, 33, 804–820.

(28)

Olovsson, I.; Templeton, D. H.; Rundqvist, S.; Varde, E.; Westin, G. The Crystal Structure of Hydrogen Peroxide Dihydrate. Acta Chem. Scand. 1960, 14, 1325–1332.

(29)

Scatchard, G.; Kavanagh, G. M.; Ticknor, L. B. Vapor-Liquid Equilibrium. VIII. Hydrogen Peroxide-Water Mixtures. J. Am. Chem. Soc. 1952, 74, 3715–3720.

(30)

Easton, M. F.; Mitchell, A. G.; Wynne-Jones, W. F. K. The Behaviour of Mixtures of Hydrogen Peroxide and Water. Part 1.—Determination of the Densities of Mixtures of Hydrogen Peroxide and Water. Trans. Faraday Soc. 1952, 48, 796–801.

(31)

Huckaba, C. E.; Keyes, F. G. The Density of Aqueous Hydrogen Peroxide Solutions. J. Am. Chem. Soc. 1948, 70, 2578–2581.

ACS Paragon Plus Environment

40

Page 41 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

(32)

Foley, W. T.; Giguère, P. A. Hydrogen Peroxide and Its Analogues: II. Phase Equilibrium in the System Hydrogen Peroxide – Water. Can. J. Chem. 1951, 29, 123–132.

(33)

Cremer, D. Theoretical Determination of Molecular Structure and Conformation. I. The Role of Basis Set and Correlation Effects in Calculations on Hydrogen Peroxide. J. Chem. Phys. 1978, 69, 4440–4455.

(34)

Carpenter, J. E.; Weinhold, F. Torsion-Vibration Interactions in Overtone Excited States of Hydrogen Peroxide. J. Phys. Chem. 1986, 90, 6405–6408.

(35)

Koput, J.; Carter, S.; Handy, N. C. Potential Energy Surface and Vibrational−Rotational Energy Levels of Hydrogen Peroxide. J. Phys. Chem. A 1998, 102, 6325–6330.

(36)

Kulkarni, S. A.; Bartolotti, L. J.; Pathak, R. K. Ab Initio Investigations on Neutral Hydrogen Peroxide Clusters: (H2O2)n (n=2–4). Chem. Phys. Lett. 2003, 372, 620–626.

(37)

Praprotnik, M.; Janežič, D. Molecular Dynamics Integration and Molecular Vibrational Theory. II. Simulation of Nonlinear Molecules. J. Chem. Phys. 2005, 122, 174102.

(38)

Elango, M.; Parthasarathi, R.; Subramanian, V.; Ramachandran, C. N.; Sathyamurthy, N. Hydrogen Peroxide Clusters: The Role of Open Book Motif in Cage and Helical Structures. J. Phys. Chem. A 2006, 110, 6294–6300.

(39)

Maciel, G. S.; Bitencourt, A. C. P.; Ragni, M.; Aquilanti, V. Studies of the Dynamics around the O-O Bond: Orthogonal Local Modes of Hydrogen Peroxide. Chem. Phys. Lett. 2006, 432, 383–390.

(40)

Yu, C.-Y.; Yang, Z.-Z. A Systemic Investigation of Hydrogen Peroxide Clusters (H2O2)n (n = 1-6) and Liquid-State Hydrogen Peroxide: Based on Atom-Bond Electronegativity

ACS Paragon Plus Environment

41

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 42 of 47

Equalization Method Fused into Molecular Mechanics and Molecular Dynamics. J. Phys. Chem. A 2011, 115, 2615–2626. (41)

Sánchez-Sanz, G.; Alkorta, I.; Elguero, J. Theoretical Study of the HXYH Dimers (X, Y = O, S, Se). Hydrogen Bonding and Chalcogen–chalcogen Interactions. Mol. Phys. 2011, 109, 2543–2552.

(42)

McGrath, M. P. Torsional Electric Dipole Moment Functions Calculated for HOOH and ClOOCl. J. Chem. Phys. 2013, 138, 94305.

(43)

Kulkarni, A. D.; Pathak, R. K.; Bartolotti, L. J. Structures, Energetics, and Vibrational Spectra of H2O2···(H2O)n, n = 1−6 Clusters: Ab Initio Quantum Chemical Investigations. J. Phys. Chem. A 2005, 109, 4583–4590.

(44)

Ferreira, C.; Martiniano, H. F. M. C.; Cabral, B. J. C.; Aquilanti, V. Electronic Excitation and Ionization of Hydrogen Peroxide-Water Clusters: Comparison with Water Clusters. Int. J. Quantum Chem. 2011, 111, 1824–1835.

(45)

Thürmer, S.; Seidel, R.; Winter, B.; Ončák, M.; Slavíček, P. Flexible H2O2 in Water: Electronic Structure from Photoelectron Spectroscopy and Ab Initio Calculations. J. Phys. Chem. A 2011, 115, 6239–6249.

(46)

Yu, C.; Gong, L.; Yang, Z. Theoretical Study on the Hydration of Hydrogen Peroxide in Terms of Ab Initio Method and Atom-Bond Electronegativity Equalization Method Fused into Molecular Mechanics. Front. Chem. China 2011, 6, 287–299.

(47)

Moin, S. T.; Hofer, T. S.; Randolf, B. R.; Rode, B. M. An Ab Initio Quantum Mechanical Charge Field Molecular Dynamics Simulation of Hydrogen Peroxide in Water. Comput.

ACS Paragon Plus Environment

42

Page 43 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Theor. Chem. 2012, 980, 15–22. (48)

Fedorov, D. G.; Sugita, Y.; Choi, C. H. Efficient Parallel Implementations of QM/MMREMD (Quantum Mechanical/Molecular Mechanics-Replica-Exchange MD) and Umbrella Sampling: Isomerization of H2O2 in Aqueous Solution. J. Phys. Chem. B 2013, 117, 7996–8002.

(49)

Caputo, M.; Provasi, P. F.; Benitez, L.; Georg, H. C.; Canuto, S.; Coutinho, K. Monte Carlo−Quantum Mechanics Study of Magnetic Properties of Hydrogen Peroxide in Liquid Water. J. Phys. Chem. A 2014, 6239–6247.

(50)

Vácha, R.; Slavíček, P.; Mucha, M.; Finlayson-Pitts, B. J.; Jungwirth, P. Adsorption of Atmospherically Relevant Gases at the Air/Water Interface: Free Energy Profiles of Aqueous Solvation of N2, O2, O3, OH, H2O, HO2, and H2O2. J. Phys. Chem. A 2004, 108, 11573–11579.

(51)

Domínguez, L.; Sosa-Peinado, A.; Hansberg, W. Catalase Evolved to Concentrate H2O2 at Its Active Site. Arch. Biochem. Biophys. 2010, 500, 82–91.

(52)

Chung, Y. H.; Xia, J.; Margulis, C. J. Diffusion and Residence Time of Hydrogen Peroxide and Water in Crowded Protein Environments. J. Phys. Chem. B 2007, 111, 13336–13344.

(53)

Domínguez, L.; Sosa-Peinado, A.; Hansberg, W. How Catalase Recognizes H2O2 in a Sea of Water. Proteins Struct. Funct. Bioinforma. 2014, 82, 45–56.

(54)

Cordeiro, R. M. Reactive Oxygen Species at Phospholipid Bilayers: Distribution, Mobility and Permeation. Biochim. Biophys. Acta 2014, 1838, 438–444.

ACS Paragon Plus Environment

43

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(55)

Page 44 of 47

Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Petersson, G. A.; Nakatsuji, H.; Li, X.; Caricato, M.; Marenich, A. V.; Bloino, J.; Janesko, B. G.; Gomperts, R.; Mennucci, B.; Hratchian, H. P.; Ortiz, J. V.; Izmaylov, A. F.; Sonnenberg, J. L.; Williams-Young, D.; Ding, F.; Lipparini, F.; Egidi, F.; Goings, J.; Peng, B.; Petrone, A.; Henderson, T.; Ranasinghe, D.; Zakrzewski, V. G.; Gao, J.; Rega, N.; Zheng, G.; Liang, W.; Hada, M.; Ehara, M.; Toyota, K.; Fukuda, R.; Hasegawa, J.; Ishida, M.; Nakajima, T.; Honda, Y.; Kitao, O.; Nakai, H.; Vreven, T.; Throssell, K.; Montgomery, J. A., Jr.; Peralta, J. E.; Ogliaro, F.; Bearpark, M. J.; Heyd, J. J.; Brothers, E. N.; Kudin, K. N.; Staroverov, V. N.; Keith, T. A.; Kobayashi, R.; Normand, J.; Raghavachari, K.; Rendell, A. P.; Burant, J. C.; Iyengar, S. S.; Tomasi, J.; Cossi, M.; Millam, J. M.; Klene, M.; Adamo, C.; Cammi, R.; Ochterski, J. W.; Martin, R. L.; Morokuma, K.; Farkas, O.; Foresman, J. B.; Fox, D. J. Gaussian 16, revision A.03; Gaussian Inc.: Wallingford, CT, 2016.

(56)

Boys, S. F.; Bernardi, F. The Calculation of Small Molecular Interactions by the Differences of Separate Total Energies. Some Procedures with Reduced Errors. Mol. Phys. 1970, 19, 553–566.

(57)

Brooks, B. R.; Brooks, C. L.; Mackerell, A. D.; Nilsson, L.; Petrella, R. J.; Roux, B.; Won, Y.; Archontis, G.; Bartels, C.; Boresch, S.; Caflisch, A.; Caves, L.; Cui, Q.; Dinner, A. R.; Feig, M.; Fischer, S.; Gao, J.; Hodoscek, M.; Im, W.; Kuczera, K.; Lazaridis, T.; Ma, J.; Ovchinnikov, V.; Paci, E.; Pastor, R. W.; Post, C. B.; Pu, J. Z.; Schaefer, M.; Tidor, B.; Venable, R. M.; Woodcock, H. L.; Wu, X.; Yang, W.; York, D. M.; Karplus, M. CHARMM: The Biomolecular Simulation Program. J. Comput. Chem. 2009, 30, 1545–1614.

ACS Paragon Plus Environment

44

Page 45 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

(58)

Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L. Comparison of Simple Potential Functions for Simulating Liquid Water. J. Chem. Phys. 1983, 79, 926–935.

(59)

MacKerell, A. D.; Bashford, D.; Bellott, M.; Dunbrack, R. L.; Evanseck, J. D.; Field, M. J.; Fischer, S.; Gao, J.; Guo, H.; Ha, S.; Joseph-McCarthy, D.; Kuchnir, L.; Kuczera, K.; Lau, F. T. K.; Mattos, C.; Michnick, S.; Ngo, T.; Nguyen, D. T.; Prodhom, B.; Reiher, W. E.; Roux, B.; Schlenkrich, M.; Smith, J. C.; Stote, R.; Straub, J.; Watanabe, M.; Wiórkiewicz-Kuczera, J.; Yin, D.; Karplus, M. All-Atom Empirical Potential for Molecular Modeling and Dynamics Studies of Proteins. J. Phys. Chem. B 1998, 102, 3586–3616.

(60)

Yeh, I.-C.; Hummer, G. System-Size Dependence of Diffusion Coefficients and Viscosities from Molecular Dynamics Simulations with Periodic Boundary Conditions. J. Phys. Chem. B 2004, 108, 15873–15879.

(61)

Essmann, U.; Perera, L.; Berkowitz, M. L.; Darden, T.; Lee, H.; Pedersen, L. G. A Smooth Particle Mesh Ewald Method. J. Chem. Phys. 1995, 103, 8577–8593.

(62)

Hoover, W. G. Canonical Dynamics: Equilibrium Phase-Space Distributions. Phys. Rev. A, Gen. Phys. 1985, 31, 1695–1697.

(63)

Martyna, G. J.; Tobias, D. J.; Klein, M. L. Constant Pressure Molecular Dynamics Algorithms. J. Chem. Phys. 1994, 101, 4177–4189.

(64)

Martyna, G. J.; Tuckerman, M. E.; Tobias, D. J.; Klein, M. L. Explicit Reversible Integrators for Extended Systems Dynamics. Mol. Phys. 1996, 87, 1117–1157.

ACS Paragon Plus Environment

45

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(65)

Page 46 of 47

Deng, Y.; Roux, B. Hydration of Amino Acid Side Chains: Nonpolar and Electrostatic Contributions Calculated from Staged Molecular Dynamics Free Energy Simulations with Explicit Water Molecules. J. Phys. Chem. B. 2004, 108, 16567–16576.

(66)

Kumar, S.; Rosenberg, J. M.; Bouzida, D.; Swendsen, R. H.; Kollman, P. A. THE Weighted Histogram Analysis Method for Free-Energy Calculations on Biomolecules. I. The Method. J. Comput. Chem. 1992, 13, 1011–1021.

(67)

Lombardi, A.; Palazzetti, F.; Maciel, G. S.; Aquilanti, V.; Sevryuk, M. B. Simulation of Oriented Collision Dynamics of Simple Chiral Molecules. Int. J. Quantum Chem. 2011, 111, 1651–1658.

(68)

Lide, D. R. CRC Handbook of Chemistry and Physics, 90th ed.; CRC Press/Taylor and Francis: Boca Raton, FL, 2009.

(69)

Orabi, E. A.; Lamoureux, G. Simulation of Liquid and Supercritical Hydrogen Sulfide and of Alkali Ions in the Pure and Aqueous Liquid. J. Chem. Theory Comput. 2014, 10, 3221– 3235.

(70)

Henri, V. Gesetze Der Enzymwirkung Und Heterogene Katalyse. Zeitschrift für Elektrotechnik und Elektrochemie 1905, 11, 790–794.

(71)

Martins-Costa, M. T. C.; Ruiz-López, M. F. Molecular Dynamics of Hydrogen Peroxide in Liquid Water Using a Combined Quantum/classical Force Field. Chem. Phys. 2007, 332, 341–347.

(72)

Giguère, P. A. Molecular Association and Structure of Hydrogen Peroxide. J. Chem. Educ. 1983, 60, 399–401.

ACS Paragon Plus Environment

46

Page 47 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

TOC Graphic

ACS Paragon Plus Environment

47