Article pubs.acs.org/JPCA
Ab Initio Calculations of the Interaction between CO2 and the Acetate Ion Janice A. Steckel* National Energy Technology Laboratory, U.S. Department of Energy, Pittsburgh, Pennsylvania 15236, United States S Supporting Information *
ABSTRACT: A series of ab initio calculations designed to investigate the interaction of CO2 with acetate are presented. The lowest energy structure, AC−CO2-η2, is predicted by CCSD(T)/aVTZ to be bound by −10.6 kcal/ mol. Six of the bound complexes have binding energies on the order of −8 kcal/mol, but analysis shows that the η1-CT complex is fundamentally different from the others. The η1-CT complex is characterized by geometric distortion, large polarization and induction effects and charge transfer whereas the other five complexes have little geometric distortion and negligible charge transfer. The amount of charge that is transferred from the anion to the CO2 in the η1-CT complex is estimated to be about half an electron by NPA, DMA, CHELPG, and Mulliken analyses, whereas the EDA-ALMO-CTA (B3LYP) approach predicts a charge transfer of 75 me−. However, the transfer of this small amount of charge leads to an energy lowering of −56 kcal/mol, without which the complex would not be bound. The RI-MP2 geometries closely approximate those resulting from the CCSD optimizations, and the optimized second-order opposite spin (O2) method performs well for all the complexes except for the η1CT complex. DFT methods do not reproduce all the ab initio geometries, binding energies and/or energy ordering of these complexes although the range-separated hybrid meta-GGA (M11) and nonlocal (VV10 and vdwDF10) functionals are shown to yield results significantly better than other functionals considered for this system. The fact that there is such variation among DFT methods has implications for DFT-based ab initio molecular dynamics simulations and for the parametrization of classical force fields based on DFT calculations.
■
INTRODUCTION Ionic liquids (ILs) are a class of liquid compounds with appealing features such as low melting point, low vapor pressure, excellent thermal stability, good solubility, and low flammability.1−3 ILs have been defined as salts with melting temperatures below 100 °C.4 Modification of the anion or the cation has a significant impact on the gas solubility and other properties of the IL. As a consequence of the number of possible substitutions, the number of potential IL combinations is on the order of 1019 and the properties of ILs can be made to vary widely.5 Understanding the origins of the properties of ILs and synthesizing new ILs with desired properties would provide a wide array of opportunities to enhance their potential applications. Among other potential applications for which ILs have been considered (such as solvents, batteries, dye-sensitized solar cells, lubricants, catalysis), there is intense interest in the use of ILs for the separation of CO2 from flue gas exhaust.6−13 Various types of supported ionic liquid membranes and IL-impregnated polymeric membranes have been proposed for this purpose.14 Much interest has been focused on CO2 solubility in ILs,15 and it has been noted that the effective capture of CO2 from flue gas exhaust will require strong (6−10 kcal/mol) absorption.16 The prediction of CO2 solubility in novel ILs via molecular simulation methods would be a valuable contribution, but typical predictions of the Henry’s law constant do not compare well with experimental values.5 This article not subject to U.S. Copyright. Published 2012 by the American Chemical Society
In this paper, ab initio calculations of the binding energy of CO2 with acetate are presented and features of the acetate− CO2 interaction are analyzed. We are interested in modeling the interaction of CO2 with acetate as a prototypical example of an anion that binds CO2 strongly. There have been previous studies of the acetate−CO2 complex; Carvalho et. al carried out MP2 calculations on the acetate and trifluoroacetate−CO2 complexes and Kim et al. carried out MP2 and B3LYP/631+G(d) calculations on acetate and CO2 and characterized the acetate anion as an electron donor to the CO2.17−19 In the current work, additional configurations and accurate coupledcluster results with a large basis set are presented. The details of this interaction have implications for the use of DFT-based molecular dynamics as well as force field development.
■
AB INITIO CALCULATIONS The purpose of the ab initio calculations was to obtain an accurate description of the interaction of CO2 with the acetate anion. Coupled-cluster singles and doubles (CCSD) geometry optimizations with the aug-cc-VTZ (aVTZ) basis set were carried out on eight acetate−CO2 configurations to characterize the interaction of CO2 with the acetate anion. For each optimized structure, a subsequent single point energy was Received: June 29, 2012 Revised: October 19, 2012 Published: October 26, 2012 11643
dx.doi.org/10.1021/jp306446d | J. Phys. Chem. A 2012, 116, 11643−11650
The Journal of Physical Chemistry A
Article
Figure 1. CCSD/aVTZ optimized structures of acetate−CO2.
calculated adding the noniterative triples correction (CCSD(T)/aVTZ). These calculations were carried out using MOLPRO.20 The counterpoise corrected binding energy, (EbindCPC), and the interaction energy, (EintCPC), of each acetate−CO2 complex were calculated according to the following definition:
Resolution of the identity second order Møller−Plesset perturbation theory (RI-MP2)23 calculations were also carried out on the acetate−CO2 configurations. The optimized fitting basis sets of Weigend and co-workers have been employed for the RI-MP2 calculations.24,25 RI-MP2 results were compared with MP2 results for a number of the complexes. The interactions differ by ≤0.01 kcal/mol and may be obtained at up to 90% reduction in the computational time. RI-MP2 geometry optimizations were carried out with aug-cc-pVXZ (X = 2, 3, 4).26−29 Use of the correlation-consistent bases enabled extrapolation to the complete basis set (CBS) limit using the Klopper-type30 expression:
CPC E bind = [EAB − EA − E B] + [EA − 0 − EA − X ]
+ [E B − 0 − E B − X ] CPC A B = E int + Edef + Edef
(1)
in which EAB is the energy of the optimized complex, EA is the energy of relaxed monomer A calculated with its own basis set, EA−X is the energy of monomer A at the complex geometry using the basis set of the complex with ghost atoms replacing monomer B, EA−0 is the energy of monomer A at the complex geometry using only its own basis sets, and EB, EB−X, and EB−0 have definitions analogous to those for monomer A. EAdef and EBdef are the energies required to deform the monomers from their equilibrium geometries to the geometries they assume in the complex. The BSSE is estimated according to E BSSE = [EA − 0 − EA − X ] + [E B − 0 − E B − X ]
X CBS E bind = E bind +
a b + 4 3 X X
(3)
where EXbind was calculated with basis set X, a and b are fitting parameters, and EbindCBS is the extrapolated CBS energy. The binding energies used to extrapolate to the CBS limit were not corrected for BSSE. Frequency calculations were carried out at the RI-MP2/aVTZ level of theory to characterize each configuration as a minimum, transition state or saddle point. There have been several methods based on scaled opposite spin MP2 (SOS-MP2) developed recently to improve the scaling and accuracy of ab initio descriptions of nonbonded interactions.31−34 The optimized second-order opposite-spin (O2/aVTZ) method introduces optimized orbitals and takes the opposite spin component of the correlation energy (times an empirical scaling factor) as the correlation energy. We have applied the O2 method (with the resolution of the identity (RI)
(2)
CPC
Eint determined via eq 1 may be compared to the interaction energies determined via symmetry adapted perturbation theory (SAPT). Furthermore, it is Eint (rather than Ebind) that is employed in the determination of gas solubility via the Widom insertion method in the context of Monte Carlo-based determinations of gas solubility.21,22 11644
dx.doi.org/10.1021/jp306446d | J. Phys. Chem. A 2012, 116, 11643−11650
The Journal of Physical Chemistry A
Article
Table 1. Calculated Charge Transfer and Binding Energies for Eight Acetate−CO2 Complexesa complex
type
ΔQEDA HF
ΔQEDA B3LYP
ΔQNPA
ΔQDMA
ΔQCHELPG
ΔQMUL
EBSSE
ECPC bind
ECBS bind
η η1-CT twist TS-par η1 TS-perp SP M
min min min TS min TS saddle min
5 56 5 5 6 4 2 0
17 75 18 19 20 16 7 2
18 473 17 18 19 15 7 1
93 621 104 113 113 95 54 90
117 491 109 113 109 103 101 47
26 633 22 29 11 5 13 29
0.7 2.6 0.8 0.8 0.7 0.7 0.6 0.4
−10.6 −8.3 −8.5 −8.4 −8.4 −8.0 −7.9 −2.6
−10.1 −7.9 −7.8 −7.5 −7.8 −7.4 −7.2 −2.2
2
a
Characterization as minimum (min), transition state (TS), or saddle point (saddle) is based on RI-MP2/aVTZ optimization and frequency calculations. Charge transfer calculated via several methods (ΔQ) is reported in me−. Counterpoise-corrected CCSD(T)/aVTZ binding energies CBC (ECPC bind ) (CCSD/aVTZ-optimized geometries) are compared to the complete basis set limit result (Ebind ) based on RI-MP2/aVXZ, X = 2, 3, 4 calculations, as described in the text. Energies are expressed in kcal/mol.
Table 2. Calculated Interaction Energy Data for Eight Acetate−CO2 Complexesa complex
ESAPT exch
ESAPT pol
ESAPT disp
ESAPT ind
δHF
ESAPT int
ECPC int
EAdef
EBdef
η η1-CT twist TS-par η1 TS-perp SP M
17.3 250.9b 17.1 17.6 18.0 15.8 11.1 3.3
−17.2 −150.8b −15.5 −15.7 −16.2 −14.5 −11.3 −3.0
−6.3 −29.8b −5.8 −5.9 −5.9 −5.3 −4.2 −2.4
−5.8 −141.8b −5.7 −5.9 −5.8 −5.3 −4.4 −0.7
−1.3 −34.7b −1.5 −1.5 −1.5 −1.3 −0.9 −0.1
−11.9 −71.5b −9.9 −9.8 −9.9 −9.3 −8.8 −2.9
−12.2 −48.6 −10.2 −10.1 −10.1 −9.5 −8.9 −2.7
0.1 4.6 0.1 0.1 0.1 0.1 0.1 0.0
1.5 35.7 1.5 1.6 1.6 1.4 1.0 0.1
2
The DF-DFT-SAPT/aVQZ exchange, polarization, dispersion, induction, δHF, and interaction energies were calculated on the CCSD/aVTZ A geometries. ECPC int is the counterpoise-corrected CCSD(T) interaction energy, as defined in the text. Edef is the deformation energy of the anion, and EBdef is the deformation energy of CO2. Energies are expressed in kcal/mol. bDF-DFT-SAPT results must be considered possibly inaccurate for this complex. See text. a
were carried out on the CCSD/aVTZ-optimized geometries using the Q-Chem and MOLPRO programs.20,53 Symmetry adapted perturbation theory (SAPT) is a powerful and accurate method ideally suited to directly calculate weak or nonbonding interactions between molecules or ions.54 When the interacting monomers are expressed in terms of the Kohn− Sham orbitals, the method is known as DFT-SAPT.55−62 For additional computational speedup without significant loss of accuracy, the method has also been combined with density fitting (DF).63,64 The DF-DFT-SAPT/aVQZ calculations were carried out (on the CCSD/aVTZ optimized geometries) to analyze the interaction energy in terms of physically meaningful components such as electrostatic, induction, dispersion, and exchange-repulsion energies. The δ Hartree−Fock term (δHF) has been calculated and included with the SAPT induction energy. These calculations were carried out using MOLPRO.20
approximation) to the acetate−CO2 complexes with a scaling factor of 1.2, to determine if this method provides an accurate description of the acetate CO2 interaction. Starting with the CCSD-optimized structures, density functional theory (DFT) geometry reoptimizations and evaluation of the binding energies were carried out out using the aVTZ basis and a number of popular functionals. The functionals used included PBE, BLYP, and B3LYP.35−39 The Grimme dispersion correction (D3) in conjunction with the PBE functional was examined as well.40 To better account for the effects of dispersion, calculations were undertaken with the nonlocal correlation functional of Vydrov and Van Voorhis (VV10)41 and that of Lundqvist and Langreth (vdW-DF-10).42 The VV10 nonlocal correlation functional was used with the parameters b = 5.9 and C = 0.0093 together with Perdew and Wang local correlation (PW92).35 With both the nonlocal correlation functionals, use was made of the exchange functional of Perdew and Wang (rPW86) refitted for use with the vdW-DF-10 and VV10 dispersion functionals.43 Finally, calculations were carried out using the range-separated hybrid meta-GGA (M11) of Truhlar’s group44 and the rangeseparated functional LRC-ωPBEh with an attenuation parameter of ω = 0.2 bohr−1.45 The amount of charge transferred from the anion to the CO2 molecule was estimated using several methods including the charge transfer analysis based on energy decomposition analysis with absolutely localized molecular orbitals (EDA-ALMOCTA),46,47 natural population analysis (NPA),48 distributed multipole analysis (DMA),49 charges from electrostatic potentials using a grid based method (CHELPG),50,51 and Mulliken population analysis (MULL).52 These calculations
■
RESULTS AND DISCUSSION Eight stationary points on the CCSD/aVTZ potential energy surface for the acetate−CO2 system are illustrated in Figure 1, and CCSD(T) counterpoise-corrected binding and interaction energies for these complexes are presented in Tables 1 and 2. (CCSD(T) single point energies were calculated for the CCSD-optimized structures.) The values for the T1 diagnostic65 for the CCSD(T) calculations were all
