Benchmark-Quality Semiexperimental Structural Parameters of van

Subscriber access provided by READING UNIV

Article

Benchmark-Quality Semi-Experimental Structural Parameters of Van Der Waals Complexes Peter Kraus, Daniel A Obenchain, and Irmgard Frank J. Phys. Chem. A, Just Accepted Manuscript • DOI: 10.1021/acs.jpca.7b10797 • Publication Date (Web): 04 Jan 2018 Downloaded from http://pubs.acs.org on January 4, 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 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 35 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

Benchmark-Quality Semi-Experimental Structural Parameters of van der Waals Complexes P. Kraus,∗ D. A. Obenchain, and I. Frank Institut f¨ ur Physikalische Chemie und Elektrochemie, Leibniz Universität Hannover, Callinstraße 3A, 30165 Hannover, Germany E-mail: [email protected]

1

ACS Paragon Plus Environment

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

Abstract Accurate datasets including non-covalent interactions have become essential for benchmarking computational methods. However, while there is much focus on obtaining an accurate description of relative energies, reliable prediction of accurate equilibrium geometries is also important. To facilitate the benchmarking of computed geometries, the current work includes an accurate dataset of semi-experimental equilibrium geometries of non-covalent complexes that can be directly compared to ab initio data. The structures are based on high-accuracy spectroscopic data, combined with vibrational corrections at the double-hybrid density functional level. The current work is designed to complement available datasets of semi-experimental geometries of small rigid molecules, and ab initio geometries of complexes. The benchmark-quality data comprises 16 complexes, and includes dispersion interactions, hydrogen bonding, CH/π · ·π interactions and trimers. In addition to the reference data, accurate counterpoise-corrected geometries have been obtained up to the CCSD level, along with interaction energies. A short overview of the performance of computational methods, including dispersion-corrected B3LYP and B2PLYP functionals, is also included.

1

1

2

In recent years, a large part of computational chemistry research focuses on obtaining an ac-

3

curate description of non-covalent interactions. Post Hartree-Fock (HF) wavefunction theory

4

(WFT) methods, including MP2 and the current ”gold standard” coupled-cluster CCSD(T)

5

method, are able to capture the non-covalent interaction energy more accurately than most

6

single-hybrid density functionals, thanks to a more complete description of electron corre-

7

lation. 1 However, a routine application of CCSD(T) on systems above ∼10 non-hydrogen

8

atoms is currently prohibitively expensive, 2 especially when second or third energy deriva-

9

tives are required. An alternative, considerably less costly approach includes (i) the use of

10

double-hybrid density functional theory (DFT) methods (e.g. B2PLYP 3 or PBE0-2 4 ), where

11

a part of the MP2 correlation is mixed into the exchange-correlation functional in addition

12

to the HF exchange, and/or (ii) using an empirical dispersion correction scheme with DFT.

13

In fact, the latter method has become so popular and widely used, that Grimme’s DFT-D2

14

paper 5 is the most cited paper in chemistry from the last decade. 6

15

16

Introduction

However,

most of the recent research,

reviews, 2 and benchmark databases

(e.g. NCCE31 7,8 ) are focused on reproducing accurate interaction energies (∆E int ) of non2


Page 2 of 35

Page 3 of 35 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


17

covalent complexes at pre-determined geometries. 2 From a benchmark perspective, this is

18

a necessary cost-saving measure, as interaction energies can be simply computed from the

19

difference in the energy obtained from a single-point calculation of the monomers and the

20

complex. However, to compare with experimentally accessible dissociation energies (De ),

21

costly geometry optimisations to consider monomer deformation upon complexation, and

22

even more costly anharmonic frequency calculations to correct for changes in zero point

23

energy are required.

24

On the other hand, benchmark databases focused on accurate geometries (e.g. CS20, 8

25

CCse22 9 ) generally only include small organic molecules, or perhaps lattice constants

26

(e.g. PS47 8 ). While these databases are very useful to determine the performance of a given

27

DFT or WFT method in calculating covalent bond lengths, the transferability of those results

28

to non-covalent complexes may be limited. Additionally, theoretical equilibrium structures

29

are often used to predict experimental lines in rotational spectra, as well as a starting point

30

for fitting experimental rotational constants. Therefore, a good initial structure helps obtain

31

experimental data faster, and can significantly improve the quality of the final structure.

32

To our best knowledge, the only database of accurate equilibrium geometries

33

of non-covalent complexes is the A24 database, obtained with counterpoise-corrected

34

ˇ aˇc and Hobza. 10 While the set of geometries is indeed of CCSD(T)/CBS calculations by Rez´

35

state-of-the-art accuracy from theoretical point of view, the presented optimised structures

36

are often not the experimentally observed conformers. A reliable database of benchmark-

37

quality structural parameters for non-covalent interactions based on experimental data does

38

not exist. Therefore, the goal of the current work is to produce a set of benchmark struc-

39

tural parameters covering a wide range of non-covalent interactions. Four sets of van der

40

Waals complexes are considered in this study: (i) rare gas (Rg ) dimers with OCS and ethane

41

to capture the effect of increased polarisability, (ii) dimers with hydrogen halides to in-

42

clude hydrogen bonding, (iii) a selection of trimers to consider higher order structures, and

43

(iv) π · ·π and CH· · π interactions with unsaturated compounds. Complexes of water have 3



Ar

N

Page 4 of 35

N

HCl

3

Ne

C

HCCH

C

C

C

C C

Figure 1: Overview of the four classes of the studied complexes: complexes with predominantly dispersion interactions (yellow), hydrogen-bonded complexes (blue), complexes with CH/π · ·π interactions (red), and trimers (white). 44

been explicitly excluded from this study, as they have been subjected to a large amount of

45

computational and experimental work, summarised in a recent review. 11 Halogen-bonding

46

complexes are also omitted from the current database; this type of bonding has been recently

47

reviewed using computational methods including symmetry-adapted perturbation theory by

48

Koláˇr and Hobza. 12

49

2

50

The framework used here to calculate semi-experimental structures by combining experimen-

51

tal rotational constants ((B0β )E ) and theoretical equilibrium corrections ((∆Beβ )T ) has been

52

described by Piccardo et al. 9 With the exception of the level of theory used, it is applied here

53

analogously to previous work. To summarise, the equilibrium bond distances (re ), accessible

Theoretical methods

4


Page 5 of 35 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


E

(Beβ)SE = (B0β)E - (∆Beβ)T correction between r0 and re eg. vibrational effects

r

v2 v1 v0

state average r0 : ground experimental data geometry re : equilibrium ab initio calculations

Figure 2: A diagram of the equilibrium and ground-state average bond lengths superimposed over a bond potential. 54

by eg. ab initio methods, can differ significantly from the vibrational ground-state averaged

55

bond lengths (r0 ) accessible in experiments. This is predominantly due to zero-point effects

56

and anharmonicity of the bond potential (see Figure 2). However, a semi-experimental equi-

57

librium rotational constant ((Beβ )SE ) along axis β, that is corrected for such effects, can be

58

approximated from Eq. (1),

(Beβ )SE = (B0β )E − (∆Beβ )T

(1)

59

where (B0β )E is obtained directly from experiment, and (∆Beβ )T corresponds to the sum of

60

the quantum-mechanical electronic and vibrational corrections to the rotational constant:

β (∆Beβ )T = ∆Belβ + ∆Bvib

(2)

61

The electronic term is calculated using Eq. (3) from the ratio of electronic and protonic

62

masses (me− and mp+ respectively), the diagonal element of the rotational g tensor (gββ )

63

and the theoretical rotational constant ((Beβ )T ). It is usually negligible, but for consistency

64

with Piccardo et al. 9 it is systematically included. The vibrational term, calculated using 5



Page 6 of 35

65

Eq. (4), is a sum of the anharmonic vibration-rotation coupling constants (αiβ ) over all

66

vibrational modes (i), with di accounting for their degeneracy.

me− ββ β T g (Be ) mp+ 1X β =− α di 2 i i

∆Belβ = β ∆Bvib

(3) (4)

67

The ∆Belβ values are obtained from NMR susceptibility calculations at CP-B3LYP-

68

D3BJ/AVTZ level of theory, combining the B3LYP hybrid functional 13,14 with an empirical

69

dispersion correction term (-D3BJ) 15,16 and a correlation-consistent triple-ζ basis set aug-

70

mented by diffuse functions 17 (AVTZ). The prefix ”CP-” denotes the use of counterpoise

71

correction of the basis set superposition error (BSSE) during geometry optimisation.

72

β The ∆Bvib values are calculated at CP-B2PLYP-D3BJ/AVTZ level of theory. 3,15–17 Cal-

73

culations carried out with other methods were used for validation and comparison purposes.

74

β Piccardo et al. 9 have shown that ∆Bvib obtained with B3LYP/SNSD and B3LYP/AVTZ are

75

in a good agreement with reference CCSD(T) data for small organic molecules; here the level

76

of theory is further adapted for van der Waals complexes by the inclusion of MP2-correlation

77

in the B2PLYP functional, empirical dispersion correction, and counterpoise correction for

78

BSSE.

79

All frequency calculations are performed following a geometry optimisation at the

80

same level of theory, with tightened convergence and integration criteria,† using Gaussian

81

G09.E01. 18 The anharmonic vibration-rotation coupling constants (αiβ ) are obtained from

82

a second-order vibrational perturbation analysis. 19 However, in some cases, where the po-

83

tential energy surface around the rare gas atom is very flat, the second-order vibrational

84

perturbation theory may lead to unreliable results. 19,20 †

All calculations were performed with a pruned (250, 974) grid with SCF convergence of 0.1 nHartree. In geometry optimisations, maximum and RMS thresholds of 2 and 1 µHartree and 6 and 4 µBohr were imposed for forces and displacements respectively.

6


Page 7 of 35 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


85

The (B0β )E values for each studied molecule are obtained from microwave spectroscopy

86

data available in the literature. The experimental rotational constants are first fitted to

87

obtain an r0 structure using our own in-house code, by minimizing the root mean square

88

(RMS) of the residuals of the fitted rotational constants (”RMS residuals”) by a brute-force

89

method. This initial structure is then processed by Strfit ver. 8a.X.2016 21 to confirm con-

90

vergence and obtain uncertainty parameters. As the χ2 values are used as the minimisation

91

target of the fitting procedure in Strfit, they are reported with each r0 and reSE structure,

92

in addition to the error estimates in Section 3.7. The monomeric r0 and re bond lengths for

93

OCS, 22 C2 H4 , 23 H2 CO, 23 HF, 24 HCl, 25 HCN, 9 N2 O, 26 NH3 9 and HCCH 9 are also obtained

94

from literature. As a consequence, the starting r0 structure might differ from the originally

95

reported r0 structures; care is taken so that this difference is minimised. The re structures

96

are obtained using the same procedure and the same fit parameters as for the r0 structures,

97

with the monomeric re bond lengths used instead. Strfit-compatible input files for all

98

structures are included in the Supporting information.

99

3

Results and discussion

100

All calculated vibrational and electronic corrections to the rotational constants are listed,

101

for each isotopologue, in the Supporting information. An overview of all studied complexes

102

is shown in Fig. 1. The structural dataset comprises 16 complexes, 45 bond lengths, and

103

uses experimental rotational constants obtained from literature and calculated electronic

104

and anarmonic vibrational-rotational corrections for 72 isotopologues. As the current work

105

is focused on the accurate geometries of non-covalent interactions, only selected bond lengths

106

are reported in the following text. Lengths of bonds involving hydrogen are omitted from the

107

dataset, as the uncertainty in hydrogen positions is larger than in positions of heavier atoms,

108

as a consequence of the fitting process. Furthermore, for the HCl··H2 CO and CH2 ClF··HCCH

109

dimers, the H/D substitution data is unavailable. However, all r0 and reSE structures are

7



C

Ne Ar

C

Figure 3: A diagram of the reSE geometry of the studied Rg··C2 H4 complexes. 110

available in the Supporting information in xyz format.

111

3.1

112

The ethylene complexes of rare gas (Rg) atoms have been included in this database to

113

represent the weaker (∆E int ≤2 kJ/mol 27 ) non-covalent interactions arising purely due to

114

London forces. Selected semi-experimental equilibrium bond lengths (reSE ), their r0 counter-

115

parts, calculated equilibrium CP-B2PLYP-D3BJ/AVTZ and CP-CCSD/AVTZ bond lengths

116

(reT and reCC respectively) are listed in Table 1, with the reSE geometries shown in Fig. 3. The

117

resulting dataset is composed of the bond lengths between the Rg atom and the C atoms in

118

ethylene.

Rare gas complexes with C2 H4

119

The rotational constants used to fit the two Rg··C2 H4 complexes have been obtained by

120

Liu and Jäger, and include eight Ne··C2 H4 isotopologues and four Ar··C2 H4 isotopologues. 27

121

In the original publication, Liu and Jäger fit the Rg distance from ethylene center of mass

122

using the (B+C)/2 rotational constant. Here, the r0 and reSE structures have been obtained

123

by fitting directly to the respective B and C rotational constants, with the ethylene molecule

124

constrained to the structural parameters from Duncan. 23 The structures were further con-

125

strained to a T-shape planar geometry. 27 Despite these constraints and the somewhat high

126

χ2 value of the fit (see Table 1), the difference between the original and current r0 structures

127

is less than 0.02 ˚ A.

128

In addition to the spectroscopic data, calculations to obtain the potential energy surfaces

129

(PES) for the Rg··C2 H4 system as a function of the dimer separation and two angles have

130

been carried out by Liu and Jäger. The PES were obtained at CP-CCSD(T) level with

8


Page 8 of 35

Page 9 of 35 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: The fitted ground-state (r0 ) and semi-experimental (reSE ) bond lengths of the Rg··C2 H4 complexes in ˚ A and χ2 values. Calculated equilibrium CP-B2PLYPT D3BJ/AVTZ (re ) and CP-CCSD/AVTZ (reCC ) values included for comparison. Bond Ne ·· C1,2 χ2 Ar ·· C1,2 χ2

r0 3.830 292.16 3.953 108.61

reSE 3.513 120.82 3.842 1.07

reT 3.609 – 3.887 –

reCC 3.706 – 3.984 –

131

the AVTZ basis set supplemented by additional bond functions to describe the non-covalent

132

interaction. 27 The re values at the global minima of the three-dimensional PES are 3.612

133

and 3.947 ˚ A for the two dimers respectively. The results of our CP-B2PLYP-D3BJ/AVTZ

134

calculations are in a good agreement with the literature values (within 0.01 and 0.06 ˚ A).

135

Furthermore, our CP-CCSD/AVTZ results are in a good agreement with the literature values

136

for the Ar complex, but for the Ne complex a significantly larger separation is predicted. The

137

discrepancy between our CP-CCSD/AVTZ data and the CCSD(T) results can be attributed

138

to the additional bond functions used by Liu and Jäger 27 and inclusion of higher-order

139

correlation in CCSD(T). For the argon complex, the separation distances at the ground-

140

state average and equilibrium geometries are very similar, while for the Ne··C2 H4 complex

141

the difference is above 0.3 ˚ A, consistent with the work of Liu and Jäger. 27 This large difference

142

is due to the vibrational corrections, which reduce the B and C rotational constants by ∼18 %

143

for the Ne complex and ∼6% for the Ar complex. The calculated vibrational corrections are

144

consistent between CP-B3LYP-D3BJ, CP-MP2 and CP-B2PLYP-D3BJ (all with the AVTZ

145

basis set) for the Ar complex. However, for the Ne complex, significantly lower absolute

146

values of αiβ are obtained with CP-MP2. In both cases, the χ2 values are reduced when the

147

vibrational corrections are applied.

148

The interaction energies at the CP-CCSD(T) level were 0.97 kJ/mol and 2.10 kJ/mol

149

for the Ne and Ar complex respectively, 27 45% higher than our CP-CCSD/AVTZ results

150

(0.67 and 1.55 kJ/mol respectively). The Ar··C2 H4 complex has also been included in the

9



S

Ne Ar Kr C

O

Figure 4: A diagram of the reSE geometry of the studied Rg··OCS complexes. 151

A24 set of geometries obtained with CP-CCSD(T)/CBS. 10 However, the perpendicular T-

152

shaped isomer has been studied instead of the planar isomer. The T-shaped isomer was

153

experimentally found to be 0.9 kJ/mol higher in energy than the planar isomer. 27

154

3.2

155

The interaction between the polar OCS molecule and inert rare gas atoms enables us to

156

investigate the trends in non-covalent interactions arising due to the increasing polarisability.

157

This set of dimers covers non-covalent interactions in the 1–4 kJ/mol range. Selected semi-

158

experimental equilibrium bond lengths (reSE ), their r0 counterparts, calculated equilibrium

159

CP-B2PLYP-D3BJ/AVTZ (reT ) and CP-CCSD/AVTZ (reCC ) values are listed in Table 2,

160

with the reSE geometries shown in Fig. 4. The dataset for this class of dimers consists of the

161

distances of the Rg atom to each atom in the OCS molecule.

Rare gas complexes with OCS

162

The experimental rotational constants used to fit the Rg··OCS structures have been

163

obtained by Xu and Gerry for seven Ne··OCS isotopologues; 28 by Xu, Jäger and Gerry for

164

five Ar··OCS isotopologues; 29 and by Lovas and Suenram for three Kr··OCS isotopologues. 30

165

The r0 structures fitted as part of the current work assume a rigid, linear O=C=S geometry,

166

with the bond lengths fixed to the values of Morino and Matsumura. 22 Two structural

167

parameters – the Rg··C bond length and the O–C–Rg angle – were fitted to match the A

168

and C rotational constants. The maximum deviation between the original and current r0

169

values is 0.01 ˚ A.

170

Two-dimensional counterpoise-corrected PES for each of the three complexes have been

171

calculated by Zhu et al. 31,32 and Feng et al. 33 with CCSD(T) using at least a triple-ζ basis 10


Page 10 of 35

Page 11 of 35

2 A ∆Bvib /B0 [%]

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 C ×0.2

A

B

C

A

B C ×0.5

18 OCS

1 0 −1 −2

Ne··OCS

Ar··OCS

Kr··OCS

Figure 5: Normalised vibrational corrections for all studied isotopologues of the Rg··OCS complexes. Symbols correspond to different levels of theory (×: PBE0-D3BJ, ◦: CP-PBE0D3BJ, 4: CP-B3LYP-D3BJ, : CP-B2PLYP-D3BJ, : CP-MP2) while colours denote the basis set used (red: def2-TZVPPD, green: aug-cc-pVTZ, blue: aug-cc-pVQZ). The Ne··OCS and Kr··OCS data are scaled by a factor of 0.2 and 0.5 respectively. 172

set, augmented by additional bond functions, similarly to the work of Liu and Jäger for the

173

Ne··C2 H4 complex. 27 In each PES the OCS molecule is kept rigid. The interaction energy

174

was shown to increase with the mass of the Rg atom, with He and Ne complexes having a

175

well depth below 1 kJ/mol, while the Ar and Kr complexes have a global well of 2.6 kJ/mol

176

and 3.2 kJ/mol respectively. 33 Our calculated CP-CCSD/AVTZ interaction energies of 0.7,

177

1.7 and 2.0 kJ/mol for the Ne, Ar and Kr complexes are again significantly lower than the

178

CP-CCSD(T) results. 33

179

The normalised vibrational corrections to the rotational constants obtained with various

180

methods are shown in Fig. 5. For the case of Ne··OCS, the results obtained for some iso-

181

topologues with CP-B3LYP-D3BJ/AVTZ (◦) and for all isotopologues with CP-B2PLYP-

182

D3BJ/AVQZ () are significantly different from the remaining results, in extreme cases

183

changing sign. Furthermore, the reSE structure obtained with the CP-B2PLYP-D3BJ/AVQZ

184

data shows a much worse degree of fit with a χ2 of 0.64. On the other hand, the theoretical

185

re structures obtained with CP-B2PLYP-D3BJ and the AVTZ, AVQZ or atomic natural

186

orbital basis sets show only small differences (< 0.01 ˚ A) in the Ne··OCS geometry. As such,

187

the inconsistency in the αiβ values is attributed to the shallow well in the Ne··OCS poten-

188

tial energy surface. With increased mass of the rare gas atom, the agreement between the 11



Page 12 of 35

Table 2: The fitted ground-state (r0 ) and semi-experimental (reSE ) bond lengths of the Rg··OCS complexes in ˚ A and χ2 values. Calculated equilibrium CP-B2PLYPT D3BJ/AVTZ (re ) and CP-CCSD/AVTZ (reCC ) values included for comparison. Bond Ne ·· O Ne ·· C Ne ·· S χ2 Ar ·· O Ar ·· C Ar ·· S χ2 Kr ·· O Kr ·· C Kr ·· S χ2

r0 3.375 3.400 4.004 0.18 3.607 3.590 4.119 0.05 3.727 3.701 4.206

Benchmark-Quality Semiexperimental Structural Parameters of van

Recommend Documents