Theoretical Investigation of CâH Vibrational Spectroscopy. 1

Subscriber access provided by UNIV OF NEWCASTLE

Article

Theoretical Investigation of C-H Vibrational Spectroscopy I. Modeling of Methyl and Methylene Groups of Ethanol with Different Conformers Lin Wang, Tatsuya Ishiyama, and Akihiro Morita J. Phys. Chem. A, Just Accepted Manuscript • DOI: 10.1021/acs.jpca.7b05320 • Publication Date (Web): 17 Jul 2017 Downloaded from http://pubs.acs.org on July 18, 2017

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 free 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 accessible to all readers and 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.

The Journal of Physical Chemistry A 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 51

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

The Journal of Physical Chemistry

Theoretical Investigation of C-H Vibrational Spectroscopy I. Modeling of Methyl and Methylene Groups of Ethanol with Different Conformers Lin Wang,†,‡ Tatsuya Ishiyama,¶ and Akihiro Morita∗,†,‡ Department of Chemistry, Graduate School of Science, Tohoku University, Aoba-ku, Sendai 980-8578, Japan, Elements Strategy Initiative for Catalysts and Batteries (ESICB), Kyoto University, Kyoto 615-8520, Japan, and Department of Applied Chemistry, Graduate School of Science and Engineering, University of Toyama, Toyama 930-8555, Japan E-mail: [email protected]

Abstract A flexible and polarizable molecular model of ethanol is developed to extend our investigation of thermodynamic, structural and vibrational properties of liquid and interface. Molecular dynamics (MD) simulation with the present model confirmed that this model well reproduces a number of properties of liquid ethanol, including density, heat of vaporization, surface tension, molecular dipole moment, and trans/gauche ratio. In particular, the present model can describe vibrational IR, Raman and sum frequency generation (SFG) spectra of ethanol and partially deuterated analogues with ∗

To whom correspondence should be addressed Tohoku University ‡ Kyoto University ¶ University of Toyama †

1 ACS Paragon Plus Environment


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

reliable accuracy. The improved accuracy is largely attributed to proper modeling of the conformational dependence and the intramolecular couplings including Fermi resonance in C-H vibrations. Precise dependence of torsional motions is found to be critical in representing vibrational spectra of the C-H bending. This model allows for further vibrational analysis of complicated alkyl groups widely observed in various organic molecules with MD simulation.

1

Introduction

In modern chemistry and related disciplines, vibrational spectroscopy is a fundamental tool for studying molecular structure. Vibrational spectra of molecules and materials are quite useful to extract detailed information on chemical moieties and structure, since various chemical moieties have their own fingerprint bands in the frequency domain. Infrared (IR) absorption and Raman scattering have been commonly used as optical techniques of the vibrational spectroscopy. 1 Recently sum frequency generation (SFG) spectroscopy is widely used to selectively probe the vibrational spectra at interfaces. 2–6 To investigate the details of molecular structure from these vibrational spectroscopies, correct assignment of the spectra is of indispensable importance. Among various bands studied by the vibrational spectroscopy, the C-H stretching band of alkyl moieties is one of the most widely observed bands to date. The C-H band in 2800 ∼ 3000 cm−1 region is fairly easy to measure by conventional IR or Raman spectroscopy. The C-H band is regarded as a fingerprint to investigate polymers, 7,8 membranes, 9,10 proteins 11–13 , ionic liquids 14–16 etc., as alkyl groups are ubiquitously included in those organic and biological molecules. However, the C-H band is often troublesome to be interpreted, because the band structure tends to be complicated by various overlapping components. Analogous alkyl moieties, such as methyl (CH3 ) and methylene (CH2 ) groups, have their vibrational modes consisting of symmetric stretching, asymmetric stretching and their Fermi resonance components. While the IR or Raman spectra in the C-H stretching region have 2 ACS Paragon Plus Environment

Page 2 of 51

Page 3 of 51

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


been interpreted rather empirically for a long time, the assignment of the C-H spectra for large organic molecules is often still confusing. The vibrational SFG spectra of alkyl groups are particularly challenging to interpret in case where no additional information on the corresponding feature is available in conventional IR or Raman spectra of bulk materials. Therefore, reliable theoretical support is quite desirable to help interpreting these complicated band structure. The present study aims at establishing a theoretical method for interpreting the complicated C-H bands in a unified manner. Since the C-H stretching band is associated with the Fermi resonance to the C-H bending vibrations, modeling of C-H bending vibrations has to be also required for this purpose. In our previous works, 17,18 we carried out detailed analysis of methyl C-H band of methanol, and thereby presented the unified interpretation of the IR, Raman, and SFG spectra. While the C-H band shapes of the three spectra are apparently similar, which is rather a source of confusion, their assignment should be distinct in terms of the roles of asymmetric stretching and Fermi resonance. We clearly showed that the differences of the these band structures are understood in a unified viewpoint with the help of molecular dynamics (MD) simulation. 17,18 Encouraged by the success, we extend the development of molecular modeling toward more complicated organic molecules. Here we treat ethanol as a suitable example for the purpose. Ethanol is a good model system including both the methyl and methylene groups, and also has conformational (trans and gauche) isomers. These features are common to other organic molecules including alkyl chains. Ethanol in itself has been widely used as solvent in modern industry. Various spectroscopic studies of IR, 19–21 Raman, 19–26 and SFG 27–29 as well as theoretical calculations 30,31 have been carried out for gas and liquid ethanol. Partially deuterated ethanol, such as CD3 CH2 OH and CH3 CD2 OH, have been also experimentally investigated, 21,22,27 which provide useful information to disentangle the overlapping band structure. Nevertheless, the C-H stretching band of ethanol has not yet been fully understood. For example, the vibrational band of CD3 CH2 OH at ∼2981 cm−1 has been assigned to be either asymmetric stretching of



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

CH2 21 or Fermi resonance, 22,27 while a recent Raman study 24 pointed out a possible role of different conformers. The MD simulation readily allows for investigating a specific conformer selectively, in contrast to usual experimental observations, and thus could offer the decisive assignment if its computational accuracy is confirmed. In reliable MD calculations of vibrational spectra, molecular modeling is of key importance. The model should describe instantaneous dipole moment and polarizability to calculate IR and Raman spectra, respectively, using the time correlation function formulas. 17,32 The calculation of vibrational SFG spectra requires both the dipole and polarizability of the system. 33 These vibrational spectra are commonly represented with the molecular modeling based on the Charge Response Kernel (CRK) theory. 34,35 The CRK theory is derived from ab initio or density functional theory calculations in fully non-empirical manner, and the derived CRK model is guaranteed to represent instantaneous dipole moment and polarizability tensor in the equivalent accuracy with the ab initio or density functional calculations. We have applied the CRK model to successfully reproduce the SFG spectra in the C-H stretching region of methanol 17 and benzene. 36 However, the present modeling of ethanol is somewhat different from those of methanol and benzene in a sense that ethanol has multiple conformational isomers. In this work, we extend our modeling to allow for different conformers. The different conformers should have distinct molecular parameters to be accurately described, while they are interchangeable by some large amplitude motion(s). Our modeling of ethanol is developed to satisfy these requirements, and successfully represents the various vibrational spectra of ethanol and partially deuterated ones in a comprehensive manner. The model is utilized to establish reliable assignment of the complicated C-H band of ethanol. The remainder of this paper is constructed as follows. Section 2 presents the detailed modeling method of ethanol, and Sec. 3 describes conditions of MD simulation. Section 4 presents the results of the thorough validation of our model by MD simulation, including calculations of structural, thermodynamic, and spectroscopic properties. We will argue that the present treatment of multiple conformers as well as intramolecular vibrational couplings


Page 4 of 51

Page 5 of 51

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


lead to improved performance in a number of aspects. Brief summary and conclusion follow in Sec. 5.

2

Molecular Modeling

This section summarizes the development of flexible and polarizable model of ethanol using the CRK theory. 34,35 Figure 1 shows the trans and gauche conformers of ethanol molecule and their atomic labels. The present modeling has following two features beyond the previous CRK modeling, i.e. (i) unified treatment of methyl and methylene groups and (ii) the large amplitude coordinates which connect different conformers. Besides the two fundamental features of modeling, the present modeling of methyl and methylene groups also considers the Fermi resonance mechanisms, and abides by the permutation symmetry of equivalent hydrogens. The total potential energy of liquid ethanol is divided into two parts,

U = Uintra + Uinter ,

(1)

where Uintra and Uinter are intramolecular and intermolecular potentials, respectively. The following subsections describe these potential functions that satisfy the above requirements of modeling.

2.1

Intramolecular Potential

Uintra in Eq. (1) is represented by the sum of molecular contributions, Uintra =

∑ i

uintra (i),

where uintra (i) is the intramolecular potential of i-th ethanol molecule. Here we treat the intramolecular potential of one molecule uintra , and omit the suffix i in the following. uintra is represented as a function of natural internal coordinates (NICs). 37 The ethanol molecule has 21 (= 3 × 9 − 6) NICs, S1 ∼ S21 , whose definitions are shown in Table 1. Note that



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 51

S3 (H-O-C-C torsion) and S21 (methyl rotation) are associated to large amplitude motions, with multiple minima along each of the coordinates.

Potential function: The functional form of uintra is given by 1,2,4−20 3 ∑ √ 1 ∑ ki,j ∆Si ∆Sj + ) + c3(n) cos(nS3 ) + c21 cos( 3S21 ) } 2 | {z } n=1 i,j (iv) {z } | {z } |

−α∆S1 2

uintra = D(1 − e | {z (i)

(iii)

(ii)

+ ct |

2 ∑ i=1

CH3 CH2 ∆ri2 cos(τi ) + UFR + UFR | {z } | {z } (vi) (vii) {z }

(2)

(v)

where ∆Si denotes the displacement of Si from its equilibrium value, ∆Si = Si − Sieq . The right-hand side of Eq. (2) consists of seven terms, namely (i)-(vii). We briefly explain these terms in the following. The first term (i) is a Morse potential along S1 (O-H stretching) including anharmonicity. The second term (ii) is a general harmonic potential with a force constant matrix ki,j . Note that k1,1 is set to be 0, since the S1 motion (O-H stretching) has been treated in the first term. The general harmonic potential (ii) excludes the two coordinates, S3 and S21 , associated to the large amplitude motions. The potentials along S3 and S21 are treated in the third (iii) and fourth (iv) terms, respectively. The coordinates of large amplitude motions are not properly treated with their displacement from equilibrium, ∆Si , but with their cosine terms to account for the periodic boundary. The torsional angles τi in (v) are treated in the similar manner. The last three terms (v)-(vii) describe anharmonic couplings associated to the C-H stretching in the methyl and methylene groups. The fifth term (v) describes the coupling in the methylene group illustrated in Figure 1 (b), where ∆ri = ri − rieq (i = 1, 2) is the displacement of the Cα -Hi distance from its equilibrium value and τi is the Hi -Cα -O-H dihedral angle. The sixth (vi) and seventh (vii) terms represent the anharmonic couplings related to the Fermi resonance for methyl and methylene groups, respectively, as detailed in Eqs. (3)


Page 7 of 51

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


and (4) below. The parameters in the intramolecular potential are obtained using quantum mechanical (QM) calculations in the following procedure. The equilibrium configuration is obtained by geometry optimization of the trans conformer, and the Hessian is calculated at the configuration. The initial values of ki,j are given from the calculated Hessian. Then the diagonal elements ki,i and ct are optimized to reproduce the experimental IR, Raman and SFG spectra of liquid ethanol in the C-H stretching region. D and α are determined by the second and third derivatives of the QM energy with respect to S1 . c3(n) and c21 are obtained by the least square fitting of scanned potential surface along S3 and S21 , respectively, by QM calculation. All the QM calculations are performed at the B3LYP 38,39 level and with aug-cc-pVTZ basis set 40 using Gaussian09 package. 41

Anharmonic couplings: To describe the C-H stretching vibrations of alkyl moieties, it is crucial to consider key anharmonic couplings (v)-(vii) in the methyl and methylene groups. Although these terms are minor in terms of energetics, they are significant to describe the C-H stretching vibrations. The coupling mechanism (v) considers the fact that the conformation of OH group has some influence on the adjacent Cα group. This mechanism elucidates that equivalence of two methylene hydrogens H1 and H2 is broken in the gauche conformer, and thus distinguishes the vibrations of methylene group in the trans and gauche conformers. Since the influence of OH conformation is not negligible in the adjacent Cα group, we have incorporated this mechanism in the previous modeling of methanol about the Cα methyl group. 17 The other two terms (vi) and (vii) take account of the Fermi resonances. The modeling of the Fermi resonance has been described in details in our previous work, 17 and here we briefly show the forms of the potential functions. The methyl group involves the resonance couplings between three C-H stretching modes (S13 , S14 , S15 ) and overtones or combination of three bending modes (S16 , S17 , S18 ), while the methylene group involves the couplings be-



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 51

tween two C-H stretching modes (S7 , S8 ) and two C-H bending modes (S9 , S11 ). Considering the permutation symmetry among equivalent hydrogen sites, these coupling terms take the following forms,

CH3 2 2 UFR = κ1 (2∆S15 ∆S17 ∆S18 + ∆S14 ∆S18 − ∆S14 ∆S17 ) 2 2 2 +κ2 (∆S13 ∆S17 + ∆S13 ∆S18 ) + κ3 ∆S13 ∆S16

+κ4 (2∆S14 ∆S16 ∆S17 + 2∆S15 ∆S16 ∆S18 ), CH2 2 UFR + κ7 ∆S7 ∆S9 ∆S11 , = κ5 ∆S7 ∆S92 + κ6 ∆S7 ∆S11

(3) (4)

where the coupling parameters κ1 ∼ κ7 are determined by B3LYP/aug-cc-pVTZ calculations with quantum corrections. 17 All the intramolecular parameters are shown in Supporting Information (SI).

2.2

Intermolecular Potential

Uinter in Eq. (1) consists of Lennard-Jones (LJ) interaction ULJ and electrostatic interaction UC ; Uinter = ULJ + UC . 2.2.1

LJ interaction

The LJ interaction ULJ is calculated by

ULJ =

molecules sites ∑ ∑ i>j

[( 4εab

a,b

σab rai,bj

)12

( −

σab rai,bj

)6 ] (5)

where rai,bj is the distance between the site a of molecule i and the site b of molecule j. The LJ parameters are taken from the OPLS force field 42 listed in Table 2. The LJ parameters for √ distinct atom species are given by the Lorentz-Berthelot rule: σab = (σa +σb )/2, εab = εa εb .


Page 9 of 51

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


2.2.2

CRK model

The electrostatic interaction UC is given by

UC =

1 ∑∑ 1 ∑∑ Qai Vai − Kabi Vai Vbi , 2 i a 2 i a,b

(6)

where Qai and Vai are the partial charge and electrostatic potential of site a in molecule i. Kabi = ∂Qai /∂Vbi is the CRK of molecule i. The first term of the right-hand side of Eq. (6) indicates the intermolecular site-site Coulomb interactions, while the second term indicates the electronic reorganization energy. Vai in Eq. (6) is given with the partial charges of surrounding molecules j(̸= i) by

Vai =

∑ ∑ Qbj fai,bj , rai,bj b

(7)

j(̸=i)

where fai,bj is a damping function of short-range Coulomb interaction to avoid polarization catastrophe. 43,44 In this work, we employed the damping function determined in our previous work 45 at the hydroxyl O and H sites, where the width parameter ξ of the Gaussian charge distribution was optimized to be 0.503 ˚ A.

2.2.3

Conformational dependence of Q and K

The parameters Q and K in Eq. (6) depend on the internal coordinates of the ethanol molecule. The dependence on the ordinary coordinates is described in (a) in the following, while the internal coordinates associated with the large amplitude motions are modeled with special care in (b) and (c). (a) Ordinary coordinates St (t = 1, 2, 4 ∼ 20): The ordinary coordinates have a single minimum along each of the coordinates. Thus the conformational dependence of Qa and Kab



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 51

with respect to St (t = 1, 2, 4 ∼ 20) is represented by the first-order Taylor expansion: 1,2,4∼20 (

Qa Kab

)0 ∂Qa = + ∆St , ∂St t )0 1,2,4∼20 ( ∑ ∂Kab 0 = Kab + ∆St , ∂St t Q0a

∑

(8a)

(8b)

0 , (∂Qa /∂St )0 and (∂Kab /∂St )0 are determined at the equilibrium with respect where Q0a , Kab

to these coordinates, {∆St = 0} (t = 1, 2, 4 ∼ 20). We denote these quantities by f 0 0 , (∂Qa /∂St )0 or (∂Kab /∂St )0 ) in the following. These quantities f 0 depend on (= Q0a , Kab

the coordinates of large amplitude motions which are not considered in Eq. (8), and the dependence is described in (b) and (c) as follows.

(b) Trans-gauche torsion S3 : The S3 coordinate distinguishes the trans and gauche conformations of ethanol, where S3 = 180◦ corresponds to trans, S3 = 60◦ to gauche1, and S3 = −60◦ to gauche2. The S3 dependence of f 0 is represented as f 0 (S3 ) = f (0) + f (c) cos S3 + f (s) sin S3 .

(9)

The coefficients f (0) , f (c) , and f (s) are determined from the f 0 values at trans and gauche forms. Suppose that f 0,trans , f 0,gauche1 , and f 0,gauche2 are the parameters obtained at the optimized geometries of trans, gauche1, and gauche2 by the QM calculations, respectively. Then the coefficients in Eq. (9) are determined from these values by

f (0) = (f 0,trans + f 0,gauche1 + f 0,gauche2 )/3,

(10a)

f (c) = (−2f 0,trans + f 0,gauche1 + f 0,gauche2 )/3 √ f (s) = (f 0,gauche1 − f 0,gauche2 )/ 3

(10b)


(10c)

Page 11 of 51

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


The parameters in Eq. (10) naturally interpolate arbitrary S3 values between the trans and gauche conformers with Eq. (9).

(c) Methyl rotation τ3 , τ4 , τ5 : The remaining degree of freedom is the methyl rotation, which is described with S21 . In our previous model of methanol, 17 the dependence of the methyl rotation was neglected for simplicity, and thus the three hydrogens of methyl group possessed same parameters due to the permutation symmetry. The present model explicitly accounts for the methyl rotation, which allows for nonequivalent parameters of methyl hydrogens according to the different torsional position. This improvement is important for accurate description of C-H vibrations, as will be discussed in Sec. 4.3. The methyl rotation is treated with three torsional angles, τ3 , τ4 , τ5 , instead of S21 . As depicted in Figure 1 (c), τ3 , τ4 , τ5 are the dihedral angles of H3 -Cβ -Cα -O, H4 -Cβ -Cα -O, H5 Cβ -Cα -O, respectively. We use the angles τ3 , τ4 , τ5 in the methyl modeling for convenience to distinguish the three methyl hydrogens H3 , H4 , and H5 . Accordingly, the parameters included in Eq. (10), i.e. (f (0) , f (c) , f (s) ) or (f 0,trans , f 0,gauche1 , f 0,gauche2 ), are treated as dependent on the methyl torsional angles τ3 , τ4 , τ5 . In what follows, we discuss the dependence with an example of partial charges, f 0,X = Q0,X in Eq. (10), for each conformer X (= trans, gauche1, gauche2 ). The other parameters of f 0 = K 0 , (∂Q/∂S)0 , (∂K/∂S)0 are treated in an analogous manner, and their treatments are summarized in SI. The torsional dependence of the partial charges of methyl hydrogens at the conformation X, Q0,X Ha , is represented by X X X Q0,X Ha (τa ) = t0 + tc cos τa + ts sin τa

(a = 3, 4, 5)

(11)

where QH3 , QH4 , QH5 are the partial charges of methyl hydrogen sites H3 , H4 , H5 , respectively.



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 51

X X The coefficients tX 0 , tc , and ts are determined from the calculated partial charges by

  eq,X eq,X  = (Qeq,X tX  0 H4 + QH5 + QH3 )/3   eq,X eq,X tX = (−2Qeq,X c H4 + QH5 + QH3 )/3    √   tX = (Qeq,X − Qeq,X )/ 3 s H5 H3

(12)

eq,X eq,X where Qeq,X H3 , QH4 and QH5 in Eq. (12) are obtained by the QM calculations at the equilib-

rium X conformation shown in Figure 1. Equation (11) satisfies the permutation symmetry of three methyl hydrogens, and nevertheless reproduces the different partial charges at each eq,X eq,X conformation X (Qeq,X H3 , QH4 , QH5 ) with τ3 = π/3, τ4 = π, and τ5 = −π/3. The in-

stantaneous partial charge of the methyl carbon QCβ is determined to maintain the charge neutrality in response to Eq. (11) as a function of τ3 , τ4 and τ5 . The partial charges of the other sites are independent of the methyl rotation. The conformational dependence in the above cases (b) and (c) is summarized as follows. In the case of f 0 = Q0a , Eq. (9) indicates that f 0 is a function of the trans-gauche torsion S3 . The ingredient parameters in Eq. (9), f (0) , f (c) , and f (s) , are determined from f 0,X = Q0,X a at different conformers X by Eq. (10), and f 0,X at each conformer X is a function of the methyl angles τa (a = 3, 4, 5) by Eq. (11). Detailed formulations in the case of the other 0 parameters f 0 (= Kab , (∂Qa /∂St )0 , (∂Kab /∂St )0 ) are described in SI.

2.2.4

Calculation of parameters

The partial charges and CRK are calculated for both the trans and gauche conformers of ethanol at the B3LYP/aug-cc-pVTZ level of theory. The QM calculations were carried out using Gaussian09 package 41 with our extension to include the CRK calculation. The partial charges are determined with the ChelpG method, 46 where the damping treatment is applied 47 to the electrostatic potential (ESP) fitting with a damping parameter of λ = 0.40 to avoid ill-defined behavior of ESP charges of buried sites. The derivative quantities, ∂Qa /∂St


Page 13 of 51

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


and ∂Kab /∂St , for the trans and gauche conformers are calculated by 4-point numerical differentiation at each of the equilibrium configurations. The calculated partial charges of trans and gauche ethanol are given in Table 2. Other parameters and detailed calculation procedures are summarized in SI.

3

MD Procedures

3.1

MD conditions

MD simulations of liquid ethanol were carried out with three types of isotopes, CH3 CH2 OH, CH3 CD2 OH and CD3 CH2 OH. For each species, we employed two cell geometries for investigating bulk liquid and vapor-liquid interface, respectively. For the bulk simulation, N = 300 molecules were placed in a cubic cell of 30 ˚ A. The simulation was carried out under isothermal-isobaric (N P T ) ensemble with three-dimensional periodic boundary conditions. Initial configuration was prepared by the Packmol package, 48 where the molecules were placed randomly with a distance tolerance of 2.0 ˚ A between two atoms of different molecules. After sufficient equilibration was performed for more than 0.3 ns, a total production run was carried out for sampling time of 1.8 ns for each isotope species. The mean dimensions of the cell during the N P T run were Lx × Ly × Lz = 30.7 ˚ A × 30.7 ˚ A × 30.7 ˚ A. For the interface simulation, N = 370 molecules were placed in a rectangular cell of Lx × Ly × Lz = 30 ˚ A × 30 ˚ A × 150 ˚ A. The molecules were initially placed within the layer of −20 ˚ A < z < 20 ˚ A as illustrated in blue in Figure 2(a). Thereby two vapor-liquid interfaces are formed normal to the z axis in both sides of the liquid slab. The simulation was carried out under canonical (N V T ) ensemble with three-dimensional periodic boundary conditions. 24 independent initial configurations were prepared in parallel. The system was first equilibrated for 100 ps by restricting the molecules in the region of −25 ˚ A < z < 25 ˚ A. Then the restriction is removed and the system was further equilibrated for 100 ps. Production run of 1.2 ns was carried out for each trajectory, which gives a total production 13 ACS Paragon Plus Environment


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 51

time of 28.8 ns for each isotope species. Besides the bulk and interface simulation mentioned above, the vapor phase simulation of ethanol was also carried out for a reference. N = 6 ethanol molecules were randomly placed in a cubic cell of Lx × Ly × Lz = 40 ˚ A × 40 ˚ A × 40 ˚ A. The simulation was carried out under canonical (N V T ) ensemble with three-dimensional periodic boundary conditions. The system was first equilibrated for 300 ps and then sampled for 360 ps. For all the above simulations, the temperature was set to 298.15 K using the Nosè-Hoover thermostat 49,50 with a coupling constant of 0.4 ps. The pressure of the N P T ensemble in the bulk simulation was set to 1 bar using the Hoover barostat 51 with a coupling constant of 2 ps. The equations of motion were integrated using the velocity Verlet algorithm with a time step of 0.61 fs. The Ewald summation method 52 with a separation parameter of 0.242 ˚ A−1 was used for treating the long-range electrostatic interactions. The LJ interactions were evaluated using a 14 ˚ A cutoff, with adopting tail corrections in energy and pressure calculation. 53

3.2

Vibrational Spectra

The vibrational spectra of IR, Raman and SFG were calculated on the basis of time correlation functions. Here we briefly summarize calculated properties of the three vibrational spectroscopies by MD simulation. IR Spectra: The vibrational IR spectrum is represented by the lineshape function: 32,54 1 IIR (ω) = 2π

∫

∞

−∞

dt ⟨M(t) · M(0)⟩ e−iωt ,

(13)

where M is the dipole moment of the system, and ⟨ ⟩ denotes ensemble average. We note that the time correlation functions treated in Eq. (13) and hereafter are defined in the classical mechanics so as to be calculated by classical MD simulations. 33 The absorption coefficient α(ω) is given using the lineshape function IIR (ω) with the harmonic quantum 14 ACS Paragon Plus Environment

Page 15 of 51

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


correction by 55 α(ω)n(ω) =

4π 2 ω 2 IIR (ω), 3V ckB T

(14)

where V is the volume of the cell, c the speed of light in vacuo, kB the Boltzmann constant, T the temperature, and n(ω) is the refractive index of ethanol. 56

Raman Spectra: The isotropic and anisotropic parts of the vibrational Raman spectra are represented by the following lineshape functions: 32,54 ∫ ∞ 1 ¯ A(0)⟩ ¯ = dt ⟨A(t) e−iωt , 2π −∞ ∫ ∞ 1 aniso IRaman (ω) = dt ⟨Tr(B(t)B(0))⟩ e−iωt , 2π −∞ iso IRaman (ω)

(15) (16)

¯ is the where A¯ = (1/3)TrA is the isotropic part of polarizability tensor A, and B = A − AI traceless anisotropic part of the polarizability. Then the polarized Raman scattering, where the polarization of the Raman scattering is parallel to that of the incident radiation, is given by ∥ IRaman (ω)

)( ) ]4 ℏω ( 1 2 aniso iso = (ω0 − ω) IRaman (ω) + IRaman (ω) . c kB T 1 − e−ℏω/kB T 15 [n

(17)

In Eq. (17) ω0 is the frequency of incident radiation, which is adopted from an experimental wavelength of 532 nm. 22

SFG Spectra: The vibrational SFG spectra are represented with the second-order sus(2)

ceptibility tensor χpqr , where the suffixes p, q, r denote the space-fixed coordinates x ∼ z. (2)

(2)

(2),res

χpqr consists of vibrational resonant and nonresonant terms; χpqr = χpqr

(2),nonres

+ χpqr

. The

former is calculated with the time correlation function between the polarizability A and the dipole moment M of the system, 33

χ(2),res pqr (ω2 )

iω2 = kB T

∫

∞

dt⟨Apq (t)Mr (0)⟩ eiω2 t ,

0


(18)


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 51

in electronically off-resonant conditions. In the SFG formulas, ω1 , ω2 and ω3 = ω1 +ω2 denote the frequencies of visible, IR and SFG lights, respectively. The vibrationally nonresonant (2),nonres

term χpqr

is real and constant over the IR frequency range, and its value is assumed to

be consistent to the experimental spectra in the present SFG calculations. The calculated χ(2) tensor is utilized to interpret experimental SFG spectra measured in the optical geometry of Figure 2(b). Measured SFG spectra depend on the polarizations (s or p) of the incident lights and the SFG signal. The combination of polarizations is represented by three letters, such as ssp or sps, which denote the polarizations of SFG, visible and IR lights, respectively. The SFG intensity is proportional to the square of the effective nonlinear (2)

(2)

susceptibility, I ∝ |χeff |2 , and χeff is represented for different polarization combinations by 4 (2)

(19)

(2)

(20)

χeff,ssp = Lyy (ω3 )Lyy (ω1 )Lzz (ω2 ) sin β2 χ(2) yyz , χeff,sps = Lyy (ω3 )Lzz (ω1 )Lyy (ω2 ) sin β1 χ(2) yzy , (2)

χeff,ppp = −Lxx (ω3 )Lxx (ω1 )Lzz (ω2 ) cos β3 cos β1 sin β2 χ(2) xxz −Lxx (ω3 )Lzz (ω1 )Lxx (ω2 ) cos β3 sin β1 cos β2 χ(2) xzx +Lzz (ω3 )Lxx (ω1 )Lxx (ω2 ) sin β3 cos β1 cos β2 χ(2) zxx +Lzz (ω3 )Lzz (ω1 )Lzz (ω2 ) sin β3 sin β1 sin β2 χ(2) zzz .

(21)

(2)

Equations (19)-(21) include the nonlinear susceptibility χpqr as well as the Fresnel factors L(ωi ) and the incident/outgoing angles β1 , β2 and β3 . The Fresnel factors L(ωi ) at frequency ωi (i = 1, 2, 3) are given by, 4 2n1 (ωi ) cos γi , n1 (ωi ) cos γi + n2 (ωi ) cos βi 2n1 (ωi ) cos βi , Lyy (ωi ) = n1 (ωi ) cos βi + n2 (ωi ) cos γi ( )2 n1 (ωi ) 2n2 (ωi ) cos βi , Lzz (ωi ) = n1 (ωi ) cos γi + n2 (ωi ) cos βi n′ (ωi )

Lxx (ωi ) =


(22) (23) (24)

Page 17 of 51

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


where n1 (ωi ), n2 (ωi ) and n′ (ωi ) are the refractive indexes of vapor, liquid and interface, respectively, at the frequency ωi . γi is the refractive angle into liquid ethanol defined by the Snell’s law, n1 (ωi ) sin βi = n2 (ωi ) sin γi . The values of ωi , βi and ni (i = 1, 2, 3) are taken from experimental data and conditions. 27,56 The refractive index of the interface n′ (ωi ) is not an experimentally observable quantity by itself, and its value has been estimated by various ways though it may include some uncertainty. 4,57,58 The previous works 4,57,58 showed that the orientational analysis of ppp spectra is often sensitive to the assumed value of n′ in the three-layer model. The present (2)

MD calculation also found that the relative amplitude of χeff,ppp in Eq. (21) significantly varies with the assumed value of n′ . To deal with this uncertainty, we evaluated the n′ (ωi ) value so that the calcualted relative SFG intensities between ssp and ppp are consistent to the experimental data. The optimum value of n′ is thereby determined to be n′ (ωi ) ≃ 0.2 n1 (ωi ) + 0.8 n2 (ωi ).

(25)

This result suggests that n′ (ωi ) is close to the refractive index of liquid phase than that of vapor. While Eq. (25) is obtained rather empirically to reproduce the relative intensity of ppp polarization, it turned out to be quite consistent to the result of our previous theoretical estimate. 58 In that work, we proposed a microscopic theory of n′ to offer a well-defined, microscopic basis to the phenomenological three-layer model. The theory argued that the proper value of n′ should be defined so that the observed/calculated amplitude of SFG signal becomes equivalent to that of the three-layer model, and thereby derived the effective dielectric constant of the interface layer, ϵ′ (ωi ) = (n′ (ωi ))2 , to be 1.6 in the case of water/vapor interface. 58 That value is close to that of liquid phase (1.72) than that of vapor phase (1.0). It is intriguing to find that both the ethanol/vapor and water/vapor interfaces have the optimum n′ or ϵ′ values close to those of liquid phase. The reason is briefly elucidated based on the detailed interfacial profiles of dielectric properties. 58 In a microscopic resolution, the



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 51

dielectric properties at interface vary from those in gas to in liquid over a monolayer region of interface, and their transient profiles are nearly in parallel with that of the density. The theory dictates that the optimum n′ (ϵ′ ) is given as a weighted average of the dielectric properties with various depth. The weight factor over the interface region is given with the depth profile of the induced nonlinear (SFG/SHG) polarization. In the liquid interfaces, the nonlinear polarization is generated in a few monolayers, where the dielectric properties are already close to those of the bulk liquid. This is the essential mechanism that the calculated dielectric constant of interface becomes close to that of the bulk liquid as a consequence of the weighted average. 58 Detailed mechanism for the optimum refractive indexes (dielectric constants) of various interfaces should be further investigated.

4

Results and Discussion

4.1

Liquid Properties

In this subsection, we examine the calculated properties of ethanol in bulk liquid and interface. The results are compared with available experimental and/or calculated data to validate the performance of the present CRK model. Thermodynamic properties: Table 3 summarizes the calculated properties of ethanol in comparison to available experimental data 59–62 and previous calculation results using polarizable (PIPF, 63 FQ 64 ) and non-polarizable (C22, 64 OPLS 65,66 ) force fields. The heat of vaporization ∆Hv was calculated by ∆Hv = NA (⟨U vapor ⟩ − ⟨U liquid ⟩ + kB T ),

(26)

where ⟨U vapor ⟩ and ⟨U liquid ⟩ represent the potential energies per molecule in vapor and liquid phases, respectively. NA is the Avogadro constant to represent the quantity on molar basis. The present results of the density and heat of vaporization fairly well reproduce the exper18 ACS Paragon Plus Environment

Page 19 of 51

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


imental values. 59 The agreement with the experimental values could be further improved by optimizing the LJ parameters to the present CRK model, though we employed the OPLS parameters here to facilitate transferability to other organic molecules. The surface tension γ was calculated from the slab MD simulation by the anisotropy of pressure tensor: 67 [ ] Lz 1 γ= ⟨Pzz ⟩ − (⟨Pxx ⟩ + ⟨Pyy ⟩) , 2 2

(27)

where Pxx , Pyy , and Pzz is the xx, yy, and zz components of the pressure tensor, respectively. The present model well reproduces the experimental surface tension 60 as well as other models. 64

Dipole moment: The calculated and experimental dipole moments of ethanol in vapor (µv ) and liquid (µl ) phases are also summarized in Table 3. The calculated µv of 1.67 D and µl of 2.64 D are in good agreement with experimental values of 1.68 D and 3.04 D, 61 respectively. The calculated difference ∆µ = µl − µv = 0.97 D indicates an enhanced polarization of ethanol molecule in liquid. This enhanced polarization is mostly attributed to the electronic polarization in the present MD simulation, and the nuclear polarization has a minor contribution to this difference, about 0.05 D. Comparing the calculated values of the dipole moment with other force fields in Table 3, all the force fields somewhat underestimate the experimental value, 64 µl = 3.04 D. It is noteworthy that the dipoles µl = µv of the non-polarizable models (C22, OPLS) are also underestimated, though the non-polarizable models are not capable of reproducing the induced dipole moment in liquid, ∆µ = 0, and thus the dipole moment µv = µl is given as a phenomenological parameter. The present CRK model yields a slightly larger µl value (2.64 D) than the others, though it is still smaller than the experimental one. The underestimated dipole in the present CRK model may be partly due to the short-range damping treatment of Coulomb interaction, as described in Sec. 2.2, which is necessary to avoid divergence of polarization and to make the polarizable MD trajectories stable during liquid simulations. 19 ACS Paragon Plus Environment


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 51

Ratio of conformers: The trans-gauche ratio of ethanol is represented with the population ratio of gauche conformer ω g , ⟨ g

ω =

N (gauche) N (trans) + N (gauche)

⟩ ,

(28)

where N (trans) and N (gauche) are the numbers of trans and gauche molecules, respectively. In the MD simulation of ethanol, the gauche conformer is defined by −120◦ < S3 < 120◦ while the trans is otherwise. We note that the equilibrium trans conformer has a slightly lower energy than gauche by 0.14 kcal/mol with B3LYP/aug-cc-pVTZ calculation. Meanwhile, the potential barrier from trans to gauche is 1.02 kcal/mol, and thus the conformers are readily equilibrated by the MD simulation at room temperature. The calculated and experimental ratio in vapor (ωvg ) and liquid (ωlg ) phases are shown in Table 3. The present result in the vapor phase, ωvg = 0.59, agrees with experimental result of 0.62, 62 quantum chemical calculation of 0.62 68 at B3LYP/aug-cc-pVTZ level, and the experimental result of 0.60 69 in xenon solution. In liquid ethanol, the calculated ratio ωlg = 0.81 indicates more preference to the gauche conformer than in the vapor phase, essentially because the gauche conformer is more polar than trans. It is evidenced from the fact that the dipole moments of an isolated ethanol molecule calculated by B3LYP/aug-cc-pVTZ are 1.58 D and 1.71 D for the trans and gauche ethanol, respectively. The calculated dipole moments in liquid phase by the present MD simulation are 2.44 D and 2.66 D for trans and gauche ethanol, respectively. On the other hand, the OPLS force field showed nearly unchanged ratios of the gauche conformer, ωvg = 0.475 in vapor and ωlg = 0.500 in liquid phase. 65 This is understandable as the non-polarizable force field cannot reproduce the induced dipole moment. At the interface region with a few monolayers of ethanol, the trans-gauche ratio is calculated to be 0.80, which is quite close to that in bulk liquid, ωlg = 0.81. Further analysis of molecular orientations and related SFG spectra of ethanol surface will be discussed in the subsequent paper.


Page 21 of 51

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


Radial distribution function: Figures 3 (a)-(c) show the calculated site-site radial distribution functions (RDFs) between (a) H and H, (b) O and H, and (c) O and O sites of different ethanol molecules in liquid, where H denotes the hydroxyl hydrogen site and O denotes the hydroxyl oxygen site.

The sharp first peaks of three RDFs indicate strong

hydrogen bond network in liquid ethanol. The positions of the first peak are 2.45, 1.75 and 2.65 ˚ A for H-H, O-H and O-O, respectively, which are in agreement with the previous calculated results with PIPF and FQ models. 63,64 The integrated area of the first peak of O-H RDF up to the first minimum (2.55 ˚ A) indicates an average number of 2.0 hydrogen bonds per molecule in the liquid ethanol. In order to further validate the calculated results, we calculated the structure factors and compared them with the experimental neutron scattering results. 70 The structure factors are calculated by the Fourier transform of the site-site RDFs, ∫ Hαβ (Q) = 4πρ

∞

r2 [gαβ (r) − 1]

0

sin Qr dr, Qr

(29)

where α, β are atom sites, Q is the momentum transfer and ρ is atomic number density of total system. Figures 3 (d)-(f) show the calculated and experimental intermolecular structure factors of H-H, X-H and X-X, respectively. Here H denotes the hydroxyl hydrogen and X the heavy atoms (C and O). It is shown that the calculated structure factors are in good agreement with experimental data in lineshape and position of peaks and minima. Slightly larger discrepancies are found in the region of Q → 0, though they should not be used to evaluate the computational accuracy. It is noted that the neutron experiments involve uncertainties in Q → 0, and the H(0) values were estimated from simulations to determine the structure factors in the region of low momentum transfer. 70



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

4.2

Vibrational Spectra

In this subsection, we examine the vibrational IR, Raman and SFG spectra of ethanol as a non-deuterated (CH3 CH2 OH) and partially deuterated (CH3 CD2 OH, CD3 CH2 OH) molecular species to further validate the performance of the present model. Whole IR and Raman spectra: Figure 4 shows the calculated and experimental 71 IR spectra of liquid ethanol. The calculated spectrum well reproduces the lineshape and peak position of the experimental spectrum over the whole range of frequency, including the O-H stretching (∼ 3400 cm−1 ), C-H stretching (∼ 2900 cm−1 ), and C-H bending (∼ 1400 cm−1 ) bands. For example, the calculated peak frequency of the O-H vibrational band at 3330 cm−1 is in good agreement with the experimental peak at 3327 cm−1 . 71 The peak of the O-H stretching frequency is significantly red shifted from that in the gas phase at around 3665 cm−1 , 72 and the amount of the red-shift in the liquid phase is properly reproduced by the present model. The red shift is induced by the proper modeling of anharmonicity in the O-H stretching and electronic polarization effect on the basis of the QM calculations. On the other hand, the intensity of the C-O stretching band (∼ 1100 cm−1 ) is obviously underestimated in the calculated spectrum. This difference is due to the lack of conformational dependence of C-O stretching on H-O-C-C torsion in Eq. (2), which is not our focus in the present modeling of C-H vibrations. Figure 5 shows the calculated and experimental 25 polarized Raman spectra of liquid ethanol. The calculated Raman spectrum captures qualitative features of experimental bands, and it suffices for analyzing the character of the C-H stretching band in comparison to the experiment. In particular, we note that the intensity of C-H bending band (1300 ∼ 1500 cm−1 ) is weaker than that of the C-H stretching (2800 ∼ 3000 cm−1 ) in almost one order of magnitude, and the relative intensity is well reproduced in the calculated Raman spectrum in panel (a). The relative intensity of the two bands provides a critical test for the accuracy of the intramolecular potential model, and this issue will be further discussed in


Page 22 of 51

Page 23 of 51

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


Sec. 4.3.

IR and Raman spectra in C-H stretching region: The present study aims at analyzing the C-H stretching bands of alkyl groups, and we further focus on the C-H stretching region of ethanol and partially deuterated analogues. Figure 6 shows the calculated IR (red lines) and polarized Raman bands (blue lines) of ethanol CH3 CH2 OH and partially deuterated ones, CH3 CD2 OH and CD3 CH2 OH, together with available experimental data of these species. 21,22 Panels (a, d) of Figure 6 show the calculated IR and Raman bands of CH3 CH2 OH, which well describe the experimental spectra with three sub-bands at 2876, 2930 and 2973 cm−1 . Panels (b, e) show the calculated IR and Raman results of CH3 CD2 OH, which indicate the methyl C-H stretching vibrations. Both the IR (red lines) and Raman (blue lines) spectra show two peaks at 2933 and 2973 cm−1 , and the higher-frequency (2973 cm−1 ) band is stronger in the IR and vice versa in Raman. These band shapes are in good agreement with the experimental spectra, as shown in panels (b, e). In previous experimental studies, the higher-frequency (2973 cm−1 ) band was assigned to the asymmetric stretching of CH3 group, 22,23 whereas the assignment of the lower-frequency (2933 cm−1 ) band has been confusing. It was assigned to either symmetric stretching 23 or Fermi resonance. 22 The detailed analysis of these bands will be given in the subsequent paper. 73 We also notice two other peaks in the low-frequency tail at 2871 and 2900 cm−1 in the experimental spectra, while calculated spectra show a broad band at 2850 ∼ 2900 cm−1 instead. We confirmed that there are corresponding band components in the calculated spectra, though their width is fairly broad. The calculated spectra will be used to discuss these components in the subsequent paper. 73 Panels (c, f) shows the C-H stretching spectra of the methylene group in CD3 CH2 OH. The calculated IR and Raman spectra exhibit two peaks at 2876 and 2960 cm−1 , and their lineshapes capture those of the experiment quite well. The former peak at 2876 cm−1 was



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

assigned to symmetric stretching of CH2 group. 22,23 However, the latter peak at 2960 cm−1 has not been established to date. It was assigned to be either Fermi resonance, 22 asymmetric stretching of CH2 , 23 or the C-H stretching of the gauche conformer. 24 The assignment of this band will be also discussed in the subsequent paper. 73 In summary, the present model well reproduces the qualitative features of the experimental IR and Raman spectra of ethanol and partially deuterated analogues using the common molecular model. The subsequent paper will focus on the detailed analysis of these band structures.

SFG intensity: Figure 7 shows the calculated and experimental SFG spectra of ethanol and partially deuterated ones in the C-H stretching region in the ssp, sps and ppp polarization combinations. In general, the calculated spectra well reproduce the peak position and lineshape of experimental data. 27 Some apparent discrepancies are seen in the CH3 CD2 OH spectra in the 2800 - 2900 cm−1 region, which is of the same origin with that in the IR and Raman spectra, as discussed in Figure 6(b) and (e). The band components in this frequency region are present but broadened in the calculated spectra. Besides, the relative intensities between ssp and ppp polarizations are also well reproduced in the calculated SFG spectra, while the calculated intensity of sps spectra shows a qualitative agreement with experimental data. The present model well captures qualitative features of SFG spectra in ssp, ppp and sps polarizations. The detailed assignment will be fully examined by MD simulations in the subsequent paper.

4.3

Modeling of Large Amplitude Motions

The ethanol molecule involves multiple conformations along the H-O-C-C torsion and methyl rotation, and the present model of ethanol takes account of the dependence of their coordinates on molecular properties, as described in Sec. 2.2.3. This is an important extension of previous CRK models 54 where the dependence of small amplitude vibrations alone has been


Page 24 of 51

Page 25 of 51

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


considered. This subsection discusses the improvement of the model accuracy by explicit treatment of the large amplitude motions.

Influence of trans-gauche torsion: The trans and gauche conformers of ethanol are connected by the S3 coordinate, the H-O-C-C torsion. The molecular parameters of trans and gauche conformers are distinguished by Eq. (9) as a function of S3 . The effect of the S3 dependence can be clarified by turning off the S3 dependence on molecular parameters and examining the consequence in calculated MD results. This tentative model by turning off the S3 dependence is called Model A in the following. The Model A is readily prepared by setting 0 , f (c) = 0 and f (s) = 0 in Eq. (9), and consequently the molecular properties f 0 (= Q0a , Kab

(∂Qa /∂Si )0 , (∂Kab /∂Si )0 ) coincide with f (0) in Eq. (9) irrespective of the S3 coordinate. f (0) is the average parameter of those for trans, gauche1, and gauche2 conformers in Eq. (10a). Accordingly, the Model A adopts the average common parameters to all the conformers by neglecting the difference of conformers. This modeling is common to conventional force fields with considering no explicit torsional dependence of molecular parameters. The calculated IR and Raman spectra with Model A are shown in Figure 8 with the red lines. By comparing the red lines (Model A) to the black lines (original model), the O-H stretching band of Model A is more emphasized and red shifted than that of the orignal model. These features indicate that the polarity of the liquid is overestimated by the Model A. This is supported by the fact that the calculated dipole moment for trans and gauche ethanol are 3.15 D and 3.38 D by Model A while 2.44 D and 2.66 D by orginal model. The increased dipole moments indicate an overestimate of polarity with Model A. Another noticeable feature is that the Model A remarkably overestimates the Raman intensity of the C-H bending band (∼ 1400 cm−1 ) in Figure 8 (b). This error implies that the C-H bending band could not be analyzed by conventional molecular models with no explicit torsional dependence. This overestimation of the C-H bending band is also seen with another Model B discussed below, which is shown in the blue lines in Figure 8. This



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 51

error will be further examined in the following.

Influence of methyl rotation: In the present model, the methyl rotation affects molecular parameters through the dependence of the torsional angles τa (a = 3, 4, 5), as indicated in Eq. (11). Here we examine the effect of methyl rotation on liquid properties by comparing the tentative model with no torsional dependence. The effect of methyl rotation is turned X off by setting tX c = 0 and ts = 0 in Eq. (11) in an analogous manner as above in the trans-

gauche torsion. This model is called Model B, which eliminates the effect of methyl rotation on molecular parameters and employs average parameters among three methyl hydrogen sites. Consequently, the Model B necessarily adopts common parameters to three methyl hydrogen sites, irrespective of their conformations with respect to the methyl rotation. The calculated IR and Raman spectra with Model B are shown with the blue lines in Figure 8. Comparing the blue line (Model B ) in Panel (a) with the black one of the original model, we found that the difference between two calculated results is small in the IR spectra. In the Raman spectra in panel (b), however, the Model B remarkably overestimates the peak intensity of the C-H bending region (∼1400 cm−1 ) several times as much as the original model. Such overestimation in the C-H bending band was encountered in our previous study on vibrational analysis of propylene carbonate and dimethyl carbonate, 74 where the dependence of the methyl rotation has not been incorporated. In order to understand the influence of methyl rotation on the vibrational spectra, we examine the IR and Raman intensities of a single molecule by calculating the transition dipole and polarizability. The IR intensity at a certain normal mode qi is determined by the square of the transition dipole moment, which is proportional to the square of derivative of the dipole moment with respect to qi , i.e. ( mol IIR

≡

∂µx ∂qi

(

)2 +

∂µy ∂qi

(

)2 +

∂µz ∂qi

)2 ,

(30)

where the superscript “mol” denotes property of a single molecule. µx , µy , µz are the x, 26 ACS Paragon Plus Environment

Page 27 of 51

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


y, z components of the dipole moment of the isolated molecule, respectively. The polarized mol is represented analogously, Raman intensity IRaman

( mol IRaman

≡

∂α ¯ ∂qi

)2

2 + 15

(

∂β 2 ∂qi

) ,

(31)

where (∂ α ¯ /∂qi ) and (∂β 2 /∂qi ) are given using three principal values of polarizability, α1 , α2 , α3 , by ) ∂α1 ∂α2 ∂α3 + + , ∂qi ∂qi ∂qi [( )2 ( )2 ( )2 ] 1 ∂α1 ∂α2 ∂α2 ∂α3 ∂α3 ∂α1 = − + − + − . 3 ∂qi ∂qi ∂qi ∂qi ∂qi ∂qi

∂α ¯ 1 = ∂qi 3 ∂β 2 ∂qi

(

mol mol Table 4 summarizes the calculated IIR and IRaman of an isolated trans ethanol molecule

by Eqs. (30) and (31), respectively. The table also shows the results of B3LYP/aug-cc-pVTZ calculations, labeled “QM”. The results of the present model, labeled “Present”, show nearly quantitative agreement to those of QM, since the present CRK model is derived from the B3LYP/aug-cc-pVTZ calculations. Table 4 also lists the control results of Models A and B. mol The Models A and B yield comparable results of IIR in the C-H stretching modes, with

the ratio to the QM values in the range of 0.8 – 1.5. However, for C-H bending modes, the mol Models A and B yield remarkable overestimation of IRaman . This table shows the results of

trans ethanol, and the same tendency is also found in the gauche ethanol. These results indicate that the remarkably overestimated Raman intensity in the C-H bending band in Figure 8 is attributed to the inaccurate transition polarizability with the Models A and B. Adequate modeling of the torsional dependence is necessary for analyzing the C-H bending bands. Further analysis of the overestimation mechanism is provided in SI, which supports the importance of the conformational dependence of molecular properties.



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

5

Conclusion

In this work, we have developed a flexible and polarizable model for ethanol on the basis of the charge response kernel (CRK) theory. The present modeling properly takes account of the large amplitude motions along the H-O-C-C torsion and methyl rotation, which greatly improve the accuracy and reliability of the model. As a consequence, the present model of ethanol is able to distinguish different properties of trans and gauche conformers. The performance of the model is examined by calculating various properties of liquid ethanol by MD simulation, including density, heat of vaporization, surface tension, radial distribution functions, molecular dipole moment, and the population ratio of trans and gauche conformers. The calculated results generally well agree with experimental data as well as other simulation results. In addition to these structural and thermodynamic properties, the present model is capable of describing the vibrational spectra of IR, Raman and SFG on the same footing with the MD simulation. The reliability of the results is examined with normal ethanol and partially deuterated ones, and the calculated spectra of the normal and partially deuterated species well agree with available experimental data in their lineshapes and peak positions. The success of the vibrational modeling stems from the proper modeling of anharmonic couplings of methyl and methylene groups and the treatment of methyl rotation. We also found that the accurate treatment of methyl rotation is crucial for the C-H bending vibrations, particularly in the Raman spectra. If we neglect the effect of methyl rotation and simply treat equivalent methyl hydrogens, that approximation would lead to drastic overestimation in the C-H bending intensity. The present work demonstrates that a careful modeling of inter- and intramolecular force field of ethanol molecule can lead to reliable description of various thermodynamic and vibrational properties of condensed phase by MD simulation. While the ab initio MD could lift the load of molecular modeling in principle, we think that a reliable model with sufficient accuracy will be of great use to explore various properties of liquids and interfaces 28 ACS Paragon Plus Environment

Page 28 of 51

Page 29 of 51

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


with modest cost of computation. A remarkable feature of this work is that the present model incorporates different conformers of alkyl molecules and can be easily transferred to other organic molecules which contain different conformers. In the subsequent paper, we will apply this model to the analysis of vibrational spectra of normal ethanol and partially deuterated ones, and thereby provide unified assignment of IR, Raman and SFG spectra.

Supporting Information Available This material is available free of charge via the Internet at http://pubs.acs.org/.

Acknowledgement This work was supported by the Grants-in-Aids (JP25104003, JP26288003) by the Japan Society for the Promotion of Science (JSPS), Japan. The numerical computations were performed with the supercomputers in Research Center for Computational Science, Okazaki, Japan.



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

References (1) Silverstein, R. M.; Webster, F. X.; Kiemle, D. J.; Bryce, D. L. Spectroscopic Identification of Organic Compound, 8th Edition; John Wiley and Sons: Hoboken, NJ, 2015. (2) Shen, Y. R.; Ostroverkhov, V. Sum-Frequency Vibrational Spectroscopy on Water Interfaces: Polar Orientation of Water Molecules at Interfaces. Chem. Rev. 2006, 106, 1140. (3) Richmond, G. L. Molecular Bonding and Interactions at Aqueous Surfaces as Probed by Vibrational Sum Frequency Spectroscopy. Chem. Rev. 2002, 102, 2693–2724. (4) Wang, H. F.; Gan, W.; Lu, R.; Rao, Y.; Wu, B. H. Quantitative Spectral and Ireintational Analysis in Surface Sum Frequency Generation Vibrational Spectroscopy (SFGVS). Int. Rev. Phys. Chem. 2005, 24, 191–256. (5) Ho, J.; Psciuk, B. T.; Chase, H. M.; Rudshteyn, B.; Upshur, M. A.; Fu, L.; Thomson, R. J.; Wang, H.-F.; Geiger, F. M.; Batista, V. S. Sum Frequency Generation Spectroscopy and Molecular Dynamics Simulations Reveal a Rotationally Fluid Adsorption State of α-Pinene on Silica. J. Phys. Chem. C 2016, 120, 12578–12589. (6) McGeachy, A. C.; Olenick, L. L.; Troiano, J. M.; Lankone, R. S.; Melby, E. S.; Kuech, T. R.; Ehimiaghe, E.; Fairbrother, D. H.; Pedersen, J. A.; Geiger, F. M. Resonantly Enhanced Nonlinear Optical Probes of Oxidized Multiwalled Carbon Nanotubes at Supported Lipid Bilayers. J. Phys. Chem. B 2017, 121, 1321–1329. (7) Chen, Z.; Shen, Y. R.; Somorjai, G. A. Studies of Polymer Surfaces by Sum Frequency Generation Vibrational Spectroscopy. Annu. Rev. Phys. Chem. 2002, 53, 437–465. (8) Gerrard, D. L.; Maddam, W. F. Polymer Characterization by Raman Spectroscopy. Appl. Spectrosc. Rev. 1986, 22, 251–334.


Page 30 of 51

Page 31 of 51

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


(9) Gaber, B. P.; Peticolas, W. L. On the Quantitative Interpretation of Biomembrane Structure by Raman Spectroscopy. Biochim. Biophys. Acta 1977, 465, 260–274. (10) Liu, W.; Wang, Z.; Fu, L.; Leblanc, R. M.; Yan, E. C. Y. Lipid Compositions Modulate Fluidity and Stability of Bilayers: Characterization by Surface Pressure and Sum Frequency Generation Spectroscopy. Langmuir 2013, 29, 15022–15031. (11) Rygula, A.; Majzner, K.; Marzec, K. M.; Kaczor, A.; Pilarczyk, M.; Baranska, M. Raman Spectroscopy of Proteins: a review. J. Raman Spectrosc. 2013, 44, 1061–1076. (12) Nevin, A.; Osticioli, I.; Anglos, D.; Burnstock, A.; Cather, S.; Castellucci, E. Raman Spectra of Proteinaceous Materials Used in Paintings: A Multivariate Analytical Approach for Classification and Identification. Anal. Chem. 2007, 79, 6143–6151. (13) Yan, E. C. Y.; Wang, Z.; Fu, L. Proteins at Interfaces Probed by Chiral Vibrational Sum Frequency Generation Spectroscopy. J. Phys. Chem. B 2015, 119, 2769–2785. (14) Paschoal, V. H.; Faria, L. F. O.; Ribeiro, M. C. C. Vibrational Spectroscopy of Ionic Liquids. Chem. Rev. 2017, 117, 7053–7112. (15) Xu, S.; Xing, S.; Pei, S.; Baldelli, S. Sum Frequency Generation Spectroscopy Study of an Ionic Liquid at a Graphene-BaF2 (111) Interface. J. Phys. Chem. B 2014, 118, 5203–5210. (16) Penalber, C. Y.; Baker, G. A.; Baldelli, S. Sum Frequency Generation Spectroscopy of Imidazolium-Based Ionic Liquids with Cyano-Functionalized Anions at the Solid SaltCLiquid Interface. J. Phys. Chem. B 2013, 117, 5939–5949. (17) Ishiyama, T.; Sokolov, V. V.; Morita, A. Molecular Dynamics Simulation of Liquid Methanol. I. Molecular Modeling Including C-H Vibration and Fermi Resonance. J. Chem. Phys. 2011, 134, 024509.



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

(18) Ishiyama, T.; Sokolov, V. V.; Morita, A. Molecular dynamics simulatiojn of liquid methanol. II. Unified Assignments of infrared, Raman, and Sum Frequency Generation Vibrational Spectra in Methyl C-H stretching region. J. Chem. Phys. 2011, 134, 024510. (19) Burikov, S.; Dolenko, T.; Patsaeva, S.; Starokurov, Y.; Yuzhakov, V. Raman and IR spectroscopy research on hydrogen bonding in waterCethanol systems. Mol. Phys. 2010, 108, 2427–2436. (20) Perchard, J. P.; Josien, M. L. Etude des spectres vibrationnels de douze especes isotopioues d’ethanols monomeres. J. Chim. Phys. 1968, 65, 1834–1855. (21) Perchard, J. P.; Josien, M. L. Etude des spectres vibrationnels de douze especes isotopioues d’ethanols autoassocies. J. Chim. Phys. 1968, 65, 1856–1875. (22) Yu, Y.; Lin, K.; Zhou, X.; Wang, H.; Liu, S.; Ma, X. New C-H Stretching Vibrational Spectral Features in the Raman Spectra of Gaseous and Liquid Ethanol. J. Phys. Chem. C 2007, 111, 8971–8978. (23) Kamogawa, K.; Kaminaka, S.; Kitagawa, T. Behavior of Ethanol in Various Binary Solutions: Difference Raman Spectroscopy on the C-H Stretching Vibrations. J. Phys. Chem. 1987, 91, 222–226. (24) Chen, L.; Zhu, W.; Lin, K.; Hu, N.; Yu, Y.; Zhou, X.; Yuan, L.; Hu, S.; Luo, Y. Identification of Alcohol Conformers by Raman Spectra in the CCH Stretching Region. J. Phys. Chem. A 2015, 119, 3209–3217. (25) Dolenko, T. A.; Burikov, S. A.; Dolenko, S. A.; Efitorov, A. O.; Plastinin, I. V.; Yuzhakov, V. I.; Patsaeva, S. V. Raman Spectroscopy of Water-Ethanol Solutions: The Estimation of Hydrogen Bonding Energy and the Appearance of Clathrate-like Structures in Solutions. J. Phys. Chem. A 2015, 119, 10806–10815.


Page 32 of 51

Page 33 of 51

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


(26) Ratajska-Gadomska, B.; Gadomski, W. Influence of Confinement on Solvation of Ethanol in Water Studied by Raman Spectroscopy. J. Chem. Phys. 2010, 133, 234505. (27) Gan, W.; Zhang, Z.; Feng, R.; Wang, H. Identification of Overlapping Features in the Sum Frequency Generation Vibrational Spectra of Air/Ethanol Interface. Chem. Phys. Lett. 2006, 423, 261–265. (28) Gan, W.; Zhang, Z.; Feng, R.; Wang, H. Spectral Interference and Molecular Conformation at Liquid Interface with Sum Frequency Generation Vibrational Spectroscopy (SFG-VS). J. Phys. Chem. C 2007, 111, 8726–8738. (29) Stanners, C.; Du, Q.; Chin, R.; Cremer, P.; Somorjai, G.; Shen, Y. R. Polar Ordering at the Liquid-vapor Interface of n-alcohols (C1-C8). Chem. Phys. Lett. 1995, 232, 407– 413. (30) Pitsevich, G.; Doroshenko, I.; Pogorelov, V.; Shablinskas, V.; Balevichus, V.; Kozlovskaya, E. Nonempiric Anharmonic Computations of IR Spectra of Ethanol Conformers in B3LYP/cc-pVQZ Approximation (Stretch C-H Vibrations). Am. J. Chem. 2012, 2, 218–227. (31) Pitsevich, G. A.; Doroshenko, I. Y.; Pogorelov, V. Y.; Sablinskas, V.; Balevicius, V. Structure and Vibrational Spectra of gauche and trans conformers of Ethanol: Nonempirical Anharmonic Calculations and FTIR Spectra in Argon Matrices. Low Temp. Phys. 2013, 39, 389–400. (32) McQuarrie, D. A. Statistical Mechanics; University Science Books: Sausalito, 2000. (33) Morita, A.; Ishiyama, T. Recent Progress in Theoretical Analysis of Vibrational Sum Frequency Generation Spectroscopy. Phys. Chem. Chem. Phys. 2008, 10, 5801–5816. (34) Morita, A.; Kato, S. Ab Initio Molecular Orbital Theory on Intramolecular Charge



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

Polarization: Effect of Hydrogen Abstraction on the Charge Sensitivity of Aromatic and Nonaromatic Species. J. Am. Chem. Soc. 1997, 119, 4021–4032. (35) Ishida, T.; Morita, A. Extended Treatment of Charge Response Kernel Comprising the Density Functional Theory and Charge Regulation Procedures. J. Chem. Phys 2006, 125, 074112. (36) Kawaguchi, T.; Shiratori, K.; Henmi, Y.; Ishiyama, T.; Morita, A. Mechanisms of Sum Frequency Generation from Liquid Benzene: Symmetry Breaking at Interface and Bulk Contribution. J. Phys. Chem. C 2012, 116, 13169–13182. (37) Pulay, P.; Fogarasi, G.; Pang, F.; Boggs, J. E. Systematic Ab Initio Gradient Calculation of Molecular Geometries, Force Constants, and Dipole Moment Derivatives. J. Am. Chem. Soc. 1979, 101, 2550–2560. (38) Becke, A. D. Density-Functional Thermochemistry. III. The Role of Exact Exchange. J. Chem. Phys. 1993, 98, 5648–5652. (39) Lee, C.; Yang, W.; Parr, R. G. Development of the Colle-Salvetti Correlation-Energy Formula into a Functional of the Electron Density. Phys. Rev. B 1988, 37, 785–789. (40) Dunning Jr., T. H. Gaussian Basis Sets for Use in Correlated Molecular Calculations. I. The Atoms Boron through Neon and Hydrogen. J. Chem. Phys. 1989, 90, 1007–1023. (41) Frisch, M. J. et al. Gaussian 09, revision C.01; Gaussian, Inc.: Wallingford CT, 2010. (42) Jorgensen, W. L.; Maxwell, D. S.; Tirado-Rives, J. Development and Testing of the OPLS All-Atom Force Field on Conformational Energetics and Properties of Organic Liquids. J. Am. Chem. Soc. 1996, 118, 11225–11236. (43) Thole, B. T. Molecular Polarizabilities Calculated with a Modified Dipole Interaction. Chem. Phys. 1981, 59, 341–350.


Page 34 of 51

Page 35 of 51

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


(44) Elking, D.; Darden, T.; Woods, R. J. Gaussian Induced Dipole Polarization Model. J. Comp. Chem. 2007, 28, 1261–1276. (45) Ishiyama, T.; Morita, A. Analysis of Anisotropic Local Field in Sum Frequency Generation Spectroscopy with the Charge Response Kernel Water Model. J. Chem. Phys. 2009, 131, 244714. (46) Heinz, H.; Suter, U. W. Atomic Charges for Classical Simulations of Polar Systems. J. Phys. Chem. B 2004, 108, 18341–18352. (47) Morita, A.; Kato, S. The Charge Response Kernel with Modified Electrostatic Potential Charge Model. J. Phys. Chem. A 2002, 106, 3909–3916. (48) Martinez, L.; Andrade, R.; Birgin, E. G.; Martinez, J. M. PACKMOL: A Package for Building Initial Configurations for Molecular Dynamics Simulations. J. Comp. Chem. 2009, 30, 2157–2164. (49) Nose, S. A Unified Formulation of the Constant Temperature Molecular Dynamics Methods. J. Chem. Phys. 1984, 81, 511–519. (50) Hoover, W. G. Canonical Dynamics: Equilibrium Phase-Space Distributions. Phys. Rev. A 1985, 31, 1695–1697. (51) Martyna, G. J.; Tobias, D. J.; Klein, M. L. Constant Pressure Molecular Dynamics Algorithms. J. Chem. Phys. 1994, 101, 4177–4189. (52) Fincham, D. Optimisation of the Ewald Sum for Large Systems. Mol. Sim. 1994, 13, 1–9. (53) Allen, M. P.; Tildesley, D. J. Computer Simulation of Liquids; Clarendon Press: Oxford, 1987.



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

(54) Ishiyama, T.; Morita, A. Molecular Dynamics Study of Gas-Liquid Aqueous Sodium Halide Interfaces. I. Flexible and Polarizable Molecular Modeling and Interfacial Properties. J. Phys. Chem. C 2007, 111, 721–737. (55) Ishiyama, T.; Morita, A. Molecular Dynamics Study of Gas-Liquid Aqueous Sodium Halide Interfaces. II. Analysis of Vibrational Sum Frequency Generation Spectra. J. Phys. Chem. C 2007, 111, 738–748. (56) Sani, E.; Dell’Oro, A. Spectral Optical Constants of Ethanol and Isopropanol from Ultraviolet to Far Infrared. Opt. Mater. 2016, 60, 137–141. (57) Zhuang, X.; Miranda, P. B.; Kim, D.; Shen, Y. R. Mapping Molecular Orientation and Conformation at Interfaces by Surface Nonlinear Optics. Phys. Rev. B 1999, 59, 12632–12640. (58) Shiratori, K.; Morita, A. Molecular Theory on Dielectric Constant at Interfaces: A Molecular Dynamics Study of the Water/Vapor Interface. J. Chem. Phys. 2011, 134, 234705. (59) Wilhoit, R. C.; Zwolinski, B. J. Physical and Thermodynamic Properties of Aliphatic Alcohols. J. Phys. Chem. Ref. Data. Suppl. 1973, 2, 2. (60) Dass, N.; Singh, O. Temperature Dependence of Surface Tension in Associated Liquids. J. Phys. Soc. Jpn. 1970, 28, 806–806. (61) McClellan, A. L. Tables of Experimental Dipole Moments; Rahara Enterprises: El Cerrito, 1989. (62) Kakar, R. K.; Quade, C. R. Microwave Rotational Spectrum and Internal Rotation in gauche Ethyl Alcohol. J. Chem. Phys. 1980, 72, 4300–4307. (63) Gao, J.; Habibollazadeh, D.; Shao, L. A Polarizable Intermolecular Potential Function for Simulation of Liquid Alcohols. J. Phys. Chem. 1995, 99, 16460. 36 ACS Paragon Plus Environment

Page 36 of 51

Page 37 of 51

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


(64) Patel, S.; Brooks III, C. L. Structure, Thermodynamics, and Liquid-vapor Equilibrium of Ethanol from Molecular-dynamics Simulations Using Nonadditive Interactions. J. Chem. Phys. 2005, 123, 164502. (65) Saiz, L.; Padro, J. A.; Guardia, E. Structure and Dynamics of Liquid Ethanol. J. Phys. Chem. B 1997, 101, 78–86. (66) Jorgensen, W. L. Optimized Intermolecular Potential Functions for Liquid Alcohols. J. Phys. Chem. 1986, 90, 1276–1284. (67) Walton, J.; Tildesley, D.; Rowlinson, J.; Henderson, J. The Pressure Tensor at the Planar Surface of a Liquid. Mol. Phys. 1983, 48, 1357–1368. (68) Garden, A. L.; Paulot, F.; Crounse, J. D.; Maxwell-Cameron, I. J.; Wennberg, P. O.; Kjaergaard, H. G. Calculation of Conformationally Weighted Dipole Moments Useful in IonCmolecule Collision Rate Estimates. Chem. Phys. Lett. 2009, 474, 45–50. (69) Durig, J. R.; Deeb, H.; Darkhalil, I. D.; Klaassen, J. J.; Gounev, T. K.; Ganguly, A. The r0 Structural Parameters, Conformational Stability, Barriers to Internal Rotation, and Vibrational Assignments for trans and gauche Ethanol. J. Mol. Struct. 2011, 985, 202–210. (70) Benmore, C. J.; Loh, Y. L. The Structure of Liquid Ethanol: A Neutron Diffraction and Molecular Dynamics Study. J. Chem. Phys. 2000, 112, 5877–5883. (71) Doroshenko, I.; Pogorelov, V.; Sablinskas, V. Infrared Absorption Spectra of Monohydric Alcohols. Dataset Pap. Chem. 2013, 2013, 329406. (72) NIST Chemistry WebBook, NIST Standard Reference Database, Number 69; National Institute of Standards and Technology: Gaithersburg, MD, 2005.



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

(73) Wang, L.; Ishiyama, T.; Morita, A. Theoretical Investigation of C-H Stretching Vibrational Spectroscopy. II. Unified Assignment of IR, Raman and SFG Spectra of Ethanol. 2017, subsequent paper. (74) Wang, L.; Peng, Q.; Ye, S.; Morita, A. Surface Structure of Organic Carbonate Liquids Investigated by Molecular Dynamics Simulation and Sum Frequency Generation Spectroscopy. J. Phys. Chem. C 2016, 120, 15185–15197.


Page 38 of 51

Page 39 of 51

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


Table 1: Natural internal coordinates S1 ∼ S21 of ethanol. NICs Skeletal:

Methylene: (a)

Methyl: (b)

(a) (b)

S1 = rOH S2 = θHOCα S3 = τHOCα Cβ S4 = rCα O S5 = rCα Cβ S6 = θOCα Cβ S7 = √12 (r1 + r2 ) S8 = √12 (r1 − r2 ) S9 = α S10 = 12 (β1 − β2 + β3 − β4 ) S11 = 12 (β1 + β2 − β3 − β4 ) S12 = 12 (β1 − β2 − β3 + β4 ) S13 = √13 (r3 + r4 + r5 ) S14 = √16 (2r4 − r3 − r5 ) S15 = √12 (r3 − r5 ) S16 = √16 (α3 + α4 + α5 − β3 − β4 − β5 ) S17 = √16 (2α4 − α3 − α5 ) S18 = √12 (α3 − α5 ) S19 = √16 (2β4 − β3 − β5 ) S20 = √12 (β3 − β5 ) S21 = √13 (τ3 + τ4 + τ5 )

Character O-H str. O-H-C def. H-O-C-C torsion C-O str. C-C str. O-C-C def. CH2 sym. str. CH2 asym. str. CH2 bend. CH2 rock. CH2 wag. CH2 twist. CH3 sym. str. CH3 asym. str. 1 CH3 asym. str. 2 CH3 sym. def. CH3 asym. def. 1 CH3 asym. def. 2 CH3 rock. 1 CH3 rock. 2 torsion

ri (i = 1, 2) of the methylene group denotes the Cα -Hi bond length. ri (i = 3, 4, 5) of the methyl group denotes the Cβ -Hi bond length, and τi (i = 3, 4, 5) the Hi -Cβ -Cα -O dihedral angle as illustrated in Figure 1(c).



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 51

Table 2: ESP charges and LJ parameters, σ and ε, of trans and gauche ethanol. The notation of each atom type is shown in Figure 1. Atom Cβ H3 H4 H5 Cα H1 H2 O H

charge [e] trans gauche -0.017 0.014 0.036 0.014 -0.003 -0.012 0.036 0.015 0.031 0.055 0.055 0.093 0.055 0.032 -0.537 -0.571 0.344 0.359

σ [˚ A] 3.50 2.50 2.50 2.50 3.50 2.50 2.50 3.12 0.00


ε [J/mol] 276.14 125.52 125.52 125.52 276.14 125.52 125.52 711.28 0.00

Page 41 of 51

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


Table 3: The calculated density (ρ), heat of vaporization (∆Hv ), surface tension (γ) of liquid ethanol, average dipole moment of vapor (µv ) and liquid (µd ) ethanol, and the ratio of gauche conformer in vapor (ωvg ) and liquid (ωlg ) ethanol. The results in the present work are compared to other calculated and experimental results. Error bars are shown in parentheses. this work PIPFa ρ (g/mL) 0.793 (0.002) 0.759 ∆Hv (kcal/mol) 9.22 (0.11) 10.08 γ (mN/m) 21.9 (2.1) µv (D) 1.67 (0.01) 1.87 µl (D) 2.64 (0.01) 2.44 g ωv 0.59 ωlg 0.81 c b a Ref. 63 Ref. 64 Ref. 65 d Ref. 66 e Ref. 59 f Ref. 60 g Ref. 61 h Ref. 62

FQb C22b 0.787 0.791 10.24 10.30 22.74 23.70 1.63 2.36 2.207 2.36


OPLS Exp. 0.787c 0.785e 10.23c 10.11e 22.8f d 2.2 1.68g 2.2d 3.04g c 0.475 0.62h 0.500c


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 51

Table 4: Calculated IR and Raman intensities of an isolated trans ethanol by Eqs. (30) and (31), respectively. The derivatives are taken with respect to the normal mode coordinates qi with unit reduced mass. ”QM” denotes the calculated results by B3LYP/aug-cc-pVTZ, while “Present”, “Model A”, and “Model B ” the results of the present original model, Model A, and Model B. The values in parentheses indicate the ratio to the QM value. ω /cm−1

QM

C-H stretching CH2 sym. str. CH2 asym. str. CH3 sym. str. CH3 asym. str. 1 CH3 asym. str. 2

2980 3003 3033 3097 3101

36.9 25.2 8.86 15.5 17.5

C-H bending CH2 twist. CH3 sym. def. CH2 wag. CH3 asym. def. 1 CH3 asym. def. 2 CH2 bend.

1300 1406 1448 1483 1500 1526

0.00 0.779 6.92 3.36 1.58 0.855

Normal mode

mol /(10−6 atomic units) IIR Present Model A

37.0 25.1 8.99 15.6 17.5

(1.00) (1.00) (1.01) (1.01) (1.00)

0.00 (—) 0.774 (0.99) 6.88 (1.00) 3.39 (1.01) 1.54 (0.98) 0.845 (0.99)

33.0 21.8 10.9 18.3 20.4

(0.90) (0.87) (1.23) (1.18) (1.16)

0.0619 (—) 2.61 (3.35) 12.0 (1.74) 10.9 (3.24) 1.14 (0.72) 0.200 (0.23)

Model B 38.5 30.6 9.17 12.9 16.2

QM

mol IRaman /(10−4 atomic units) Present Model A

Model B

(1.04) (1.22) (1.03) (0.83) (0.92)

213 89.5 281 50.2 35.8

213 89.6 283 50.5 35.8

(1.00) (1.00) (1.01) (1.01) (1.00)

240 112 345 53.9 33.4

(1.13) (1.25) (1.23) (1.07) (0.93)

239 118 269 76.8 49.6

(1.12) (1.32) (0.96) (1.53) (1.39)

0.147 (—) 0.794 (1.02) 6.92 (1.00) 1.40 (0.42) 0.614 (0.39) 1.80 (2.10)

3.18 0.0837 1.20 4.02 5.96 2.38

3.19 0.0799 1.18 3.95 6.04 2.37

(1.00) (0.95) (0.98) (0.98) (1.01) (0.99)

4.75 5.40 4.05 2.18 6.92 5.95

(1.49) (64.6) (3.38) (0.54) (1.16) (2.49)

11.2 2.83 12.7 17.3 25.1 5.95

(3.51) (33.8) (10.6) (4.31) (4.21) (2.49)


Page 43 of 51

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


Figure 1: (a) Schematic pictures of ethanol molecule in trans (left) and gauche (right) forms, where the atomic labels (Cα , Cβ , O, H, H1 ∼ H5 ) and two natural internal coordinates (S3 , S21 ) are illustrated. (b, c) Newman projections of trans ethanol along (b) Cα O bond and (c) Cβ Cα bond. The torsional angles (τ1 ∼ τ5 ) are also shown.



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

Figure 2: (a) The illustration of interfacial simulation cell in this work. The blue slab denotes the initial area where ethanol molecules were placed. (b) Light geometry of SFG measurement.


Page 44 of 51

Page 45 of 51

8

(a) H−H

6

(d) H−H

4

4

0

2

−4 −8

0 6

(e) X−H

(b) O−H 0

4

H(Q)

g(r)

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


−4

2 −8 0 6

(f) X−X

(c) O−O 0

4

−4

Cal. Exp.

2 −8 0

0

2

4

6

8

10

0

2

4

6

8

10

−1

r(Å)

Q(Å )

Figure 3: The calculated site-site radial distribution functions (RDFs) between (a) H and H, (b) O and H, and (c) O and O sites of different ethanol molecules in liquid, where H denotes the hydroxyl hydrogen site and O denotes the hydroxyl oxygen site. Calculated (solid line) and experimental 70 (dash line) structure factors of (d) H-H, (e) X-H, (f) X-X in liquid ethanol are shown in the right panels, where X denotes heavy atoms (C and O) of ethanol.



(a) Cal.

Absorbance (arb. unit)

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

(b) Exp.

1000 1500 2000 2500 3000 3500

Wavenumber(cm-1)

Figure 4: (a) Calculated and (b) experimental 71 IR spectra of liquid ethanol.


Page 46 of 51

Page 47 of 51

(a) Cal.

Intensity (arb. unit)

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


(b) Exp.

1000 1500 2000 2500 3000 3500

Wavenumber(cm-1)

Figure 5: (a) Calculated and (b) experimental 25 polarized Raman spectra of liquid ethanol.



CH3CH2OH (b)

CD3CH2OH (c)

Absorbance (arb. unit)

(a)

CH3CD2OH

IR Cal.

IR Exp.

(d)

(e)

(f)

Intensity (arb. unit)

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 48 of 51

Raman Cal.

Raman Exp.

2800 2850 2900 2950 2800 2850 2900 2950 2800 2850 2900 2950

Wavenumber(cm-1) Figure 6: Calculated and experimental 21,22 IR (red) and polarized Raman (blue) spectra of CH3 CH2 OH (a, d), CH3 CD2 OH (b, e) and CD3 CH2 OH (c, f) in C-H stretching region.


Page 49 of 51

CH3CH2OH

SFG Intensity (arb. unit)

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


CH3CD2OH

CD3CH2OH

(a) Cal.

(b)

(c)

(d) Exp.

(e)

(f)

ssp ppp sps

2800 2850 2900 2950 2800 2850 2900 2950 2800 2850 2900 2950

Wavenumber(cm-1) Figure 7: Calculated (a-c) and experimental 27 (d-f) SFG spectra of CH3 CH2 OH (a, d), CH3 CD2 OH (b, e) and CD3 CH2 OH (c, f) under ssp (red), ppp (blue) and sps (pink) polarization combinations. The ssp and ppp spectra are offset for clarity.



Absorbance

(a) IR

Present Model Model A Model B

(b) Raman

Intensity

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

1000 1500 2000 2500 3000 3500

Wavenumber(cm-1)

Figure 8: The calculated IR (a) and polarized Raman (b) spectra of liquid ethanol using different models. The scale of y axis is arbitrary but comparable among different models.


Page 50 of 51

Page 51 of 51

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


Figure 9: TOC graphic. 8 cm × 4.5 cm.


Theoretical Investigation of CâH Vibrational Spectroscopy. 1

Recommend Documents

Theoretical Investigation of CâH Vibrational Spectroscopy. 1