Understanding the Structure of the Hydrogen Bonds Network and its

6 days ago - Analysis of the hydrogen bond network in aprotic ionic liquid 1-ethyl-3-methylimidazolium bis(trifluoromethylsulfonyl)imide has been perf...
0 downloads 0 Views 27MB Size
Subscriber access provided by University of South Dakota

B: Liquids, Chemical and Dynamical Processes in Solution, Spectroscopy in Solution

Understanding the Structure of the Hydrogen Bonds Network and its Influence on Vibrational Spectra in a Prototypical Aprotic Ionic Liquid Mateusz Zbigniew Brela, Piotr Kubisiak, and Andrzej Eilmes J. Phys. Chem. B, Just Accepted Manuscript • DOI: 10.1021/acs.jpcb.8b05839 • Publication Date (Web): 21 Sep 2018 Downloaded from http://pubs.acs.org on September 22, 2018

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

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

Page 1 of 44 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

Understanding the Structure of the Hydrogen Bonds Network and its Influence on Vibrational Spectra in a Prototypical Aprotic Ionic Liquid

Mateusz Z. Brela, Piotr Kubisiak, Andrzej Eilmes*

Jagiellonian University, Faculty of Chemistry, Gronostajowa 2, 30-387, Kraków, Poland

*

[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

Analysis of the hydrogen bond network in aprotic ionic liquid 1-ethyl-3methylimidazolium bis(trifluoromethylsulfonyl)imide has been performed based on structures obtained from ab initio or classical molecular dynamics simulations. Statistics of different donor and acceptor atoms and the amount of chelating or bifurcated bonds has been presented. Most of hydrogen bonds in EMIM–TFSI are formed with oxygen atom as hydrogen acceptors; and the most probable bifurcated bonds are those with a mixed pair of oxygen and nitrogen acceptors. Spectral graph analysis has shown that the cations may form hydrogen bonds to up to five different anions and the connectivity of the whole hydrogen bond network is supported mainly by H-O bonds. In the structures of the liquid simulated via force fieldbased dynamics the number of hydrogen bond is smaller and fluorine atoms are the most favored hydrogen acceptors. One-dimensional potential energy profiles for hydrogen atom displacements and corresponding vibrational frequencies have been calculated for selected CH bonds. Individual C-H stretching frequencies vary by 200-300 cm-1 indicating differences in local environment of hydrogen atoms forming C-H···O hydrogen bonds.

2

ACS Paragon Plus Environment

Page 2 of 44

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

1. Introduction Room-temperature ionic liquids (ILs) are salts, typically with bulky organic cations, molten at relatively low temperatures (below 100 °C). Large scientific attention gained by these systems stems from their physicochemical properties and prospective use as alternative solvents or solvents tailored for specific applications. Ionic liquids are composed entirely of ions; therefore strong electrostatic interactions between charged species occur in the solvent and dominate properties of ILs. Nevertheless, there is a possibility of hydrogen bonds (HBs) formation in ILs composed of cations with hydrogen atoms and anions with appropriate H-atom acceptors. Although HBs effects are weak compared to electrostatic interactions between ions, there is no doubt that hydrogen bonding is present in ILs and its importance for the solvation properties and reactivity of ionic liquids has been recognized.1,2 Experimental studies on hydrogen bonding in ionic liquids have been supported by theoretical methods. Some recent examples of quantum chemical calculations include investigations of geometries and binding energies of ion pairs,3,4,5,6,7,8 interaction energy decomposition via Symmetry Adapted Perturbation Theory (SAPT),4,9 electrostatic potential analysis,10 gas phase or solid state calculations of vibrational frequencies and intensities for isolated ion pairs or crystals,4,6,8 Natural Bond Orbital (NBO)4 or Quantum Theory of Atoms In Molecules (QTAIM) analysis for ion pairs7,8,10,11 and crystals6,8 and Non-Covalent Interaction index (NCI) calculations.10,11

Information on the structure of the liquid,

distribution of HB donor and acceptor atoms and the number of possible hydrogen bonds is obtainable from molecular dynamics (MD) simulations either employing classical force fields12,13,14,15 or based on ab initio methodology.14,16,17 Infrared spectra of bulk ILs can be calculated from ab initio molecular dynamics (AIMD) simulations by means of Fourier transform of the autocorrelation function of the dipole moment.14,16 Combined experimental and theoretical studies have been extensively used for better understanding effects of hydrogen bonding on physicochemical properties of various aprotic ionic liquids.18 Several recent works study the H-bonding correlations with the spectroscopic features of selected systems. Gu et al.19 discussed in details the HBs in the alpha-hydroxy carboxylic acid complexes reviewing current knowledge. Analysis of molecular dynamics simulations in relation to infrared (IR) spectroscopic results provides a great opportunity to gain deeper insight into physical properties of considered systems; e.g. hydrogen bond network effects on aggregation of water molecules in mixed solvents.20

3

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

Page 4 of 44

A class of commonly investigated prototypical ionic liquids is based on imidazolium derivatives. A summary of vibrational spectroscopy studies on such systems can be found in a recent review.21 Interactions in the liquid, its structure and the possibility of hydrogen bonding is usually analyzed based on quantum-chemical calculations for ion pairs7,22 and MD simulations.23,24 In

our

recent

work25

we

performed

AIMD

simulations

for

1-ethyl-3-

methylimidazolium bis(trifluoromethylsulfonyl)imide (EMIM-TFSI) liquid to study solvent effects on UV-vis and vibrational spectra. Infrared spectrum of bulk EMIM-TFSI was calculated in the range of 600 – 3500 cm-1. In this paper we used the structures or the IL recorded in the previous study to get a better understanding of the network of hydrogen bonds in EMIM-TFSI. We produced statistics of individual bonds with different donors and acceptors and determined the connectivity of the total network of HBs. For the latter purpose the methodology of spectral graph analysis was used. We also checked the structure of the IL for the presence of special cases of HBs, such as bifurcated or chelating bonds. Ab initio molecular dynamics simulations are computationally challenging and despite increasing power of computational resources, classical, inexpensive, MD simulations will continue to be a method of choice for many applications. Therefore, we performed a similar analysis for EMIM-TFSI structures obtained from force field-based simulations to compare the predictions of ab initio and classical molecular dynamics. Finally, using the geometries from AIMD simulations to construct the potential energy profiles, we calculated vibrational frequencies of one-dimensional quantum-mechanical oscillators modeling the effects of hydrogen bonds on the IR spectra.

2. Computational details 1-ethyl-3-methylimidazolium bis(trifluoromethylsulfonyl)imide is a prototypical aprotic ionic liquid. Hydrogen bonds in EMIM-TFSI are formed between cations and ions and therefore they may be classified as doubly ionic HBs.1 To test the force field-based molecular dynamics for this system we started with classical MD simulations of EMIM-TFSI using NAMD v. 2.12 simulation package.26 A nonpolarizable force field was based on the OPLS parameterization27 with values of bonded parameters adapted from the Lopes/Pádua field28,29 (including corrections for the cation from Ref. 30) and non-bonded terms from the Köddermann’s paper.31 In addition to the non4

ACS Paragon Plus Environment

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

polarizable field we used another parametrization in which polarization effects were taken into account via Drude oscillators.32 In our simulations we set the force constant of Drude oscillators to kD = 1000 kcal·mol-1·Å-2. Drude particles were attached to all non-hydrogen atoms. For the polarizable parameterization polarizabilities of atoms and partial charges were adapted from the APPLE&P polarizable force field for liquids and electrolytes.33 As a test of correctness of the FF parameters we performed additional simulations in the NPT ensemble. The calculated density of the liquid was 1.515 and 1.521 g/cm3 in non-polarizable and polarizable FF, respectively; therefore it deviated less than 0.01 g/cm3 from the experimental value of 1.518 g/cm3.34 We performed classical MD simulations for systems of 500 EMIM-TFSI ion pairs in the NVT ensemble at T=298 K with Langevin dynamics and time step of 1.0 fs used to integrate the equations of motion. To assess the possible size effects, another series of simulations with 15 ion pairs in the simulation box was performed; the smaller size corresponds to the size of the systems used in ab initio computations. Periodic boundary conditions were applied to the system, with the periodic box size set to reproduce the experimental density of the neat IL. Electrostatic interactions were taken into account via particle mesh Ewald algorithm;35 cutoff for van der Waals interactions was set to 9 Å (this value was limited by the side length of the periodic simulation box for smaller systems). The length of classical MD simulations comprised 120-150 ns of equilibration stage and the production part of the trajectory with frames saved every 1 ps (small systems) or 5 ps (large systems) for another 40 ns. As found during analysis of the results, differences between systems with 15 and 500 ion pairs are practically negligible, except for better statistics (smaller standard deviations, smoother radial distribution function plots) obtained for large systems. Therefore in the main paper we will show only the results of classical MD simulations for 500 EMIM-TFSI pairs; selected parameters are compared for both system sizes in the Supporting Information. The AIMD structures of the liquid obtained using CP2K v. 2.6.1 software36,37 were reused from the Born-Oppenheimer molecular dynamics simulations performed in our recent paper on solvent effects in EMIM-TFSI.25 Density functional theory methodology with Pade functional38 and molecularly optimized DZVP basis set was used in simulations; other details can be found in Ref. 25. In our previous work25 this setup of simulations resulted in good reproduction of the experimental IR spectrum of EMIM-TFSI,22,39 suggesting that the structure of the liquid and ion-ion interaction have been properly described.

5

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

For the analysis of HBs reported in this work we used 40 ps of the AIMD trajectory with frames saved in 1 fs intervals. The results were averaged over the whole production part of the trajectory, therefore the total number of frames used in averaging was of the same order for classical and ab initio MD, although the classical simulations spanned 1000 times longer total time. Selected frames of the AIMD trajectory were used in a series of single point energy calculations at the same level of DFT methodology as applied in MD simulations. Based on these computations potential energy profiles for displacements of hydrogen atoms in EMIM cations were constructed. This treatment allowed us to analyze the C—H stretching motions (for C2—H, C4—H and C5—H bonds) and their spectroscopic features in details. These studies were performed based on the one-dimensional (1D) potential curves for protons where displacement vectors for hydrogen atoms were along the axis of the considered C—H bonds. Snapshot structures of the analyzed system were extracted every 500 MD steps (equivalent to 0.5 ps). The corresponding internal coordinate (x) was defined as the distance between the hydrogen and the carbon atom. In each scan, position of the hydrogen atom was changed in steps of 0.05 Å in the range from 0.80 Å to 2.30 Å with respect to the hydrogen-donating carbon atom; positions of all other atoms in the structure were kept fixed. Each potential for proton motion was created based upon 31 single point calculations. In each structure of the MD simulation box, five selected hydrogen bonds were analyzed. Our procedure involved totally 37 200 single point calculations. This number looks discouragingly large; however, original calculations of 40 ps ab initio MD trajectory with 1 fs time step required 40000 calculations, therefore computational effort for potential energy calculations, although significant, was not larger than that for AIMD simulations. It should be stressed that the big advantage of such a post-MD quantization methodology is the possibility of full parallelization because each single point energy can be computed independently. This methodology was originally postulated by Stare et al.40,41 and successfully applied to several crystalline systems.42,43,44

3. Results and discussion 3.1. Structure of the liquid Before investigating hydrogen bond network in simulated structures we performed some preliminary analysis of relevant distributions of atoms. In Fig. 1 we present radial distribution functions (RDFs) for pairs of imidazolium hydrogen atom and oxygen or fluorine atom from the anion. The HIm-O and Him-F RDFs for classical simulations generally agree 6

ACS Paragon Plus Environment

Page 6 of 44

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

with the MD results reported for EMIM-TFSI liquid.45 In the HIm-O RDFs there is a well pronounced peak just above 2 Å, followed by other, lower maxima above 4 Å. Position of the first peak depends slightly on the simulation protocol: in the classical, non-polarizable simulations the maximum is located approximately at 2.25 Å and shifts to about 2.2 Å in polarizable FF. In the structures from AIMD the first maximum appears at about 2.15 Å. Interactions of imidazolium hydrogen atoms with TFSI oxygens in the ab initio simulations are therefore slightly stronger than in the non-polarizable force field-based dynamics; inclusion of polarization effects in classical simulations improves agreement with ab initio results. The main maximum in all HIm-F RDFs is located at about 5 Å. In classical MD simulations there is another, weak maximum at 2 Å, indicating that the probability of finding fluorine atom close to imidazolium hydrogen in the structures obtained from force field-based simulations is larger than in AIMD trajectories.

Fig. 1. Radial distribution functions for HIm-O and HIm-F atom pairs in EMIM-TFSI structures simulated in ab initio and classical molecular dynamics. The possibility of formation of the X-H···Y hydrogen bonds depends not only on the H-Y distance but also on the angle formed by the three atoms involved in the HB. In Fig. 2 we visualized the probability of finding the Y atom at the rXY distance from the donor atom X depending on the angle αXHY between XH and XY vectors. The differences between nonpolarizable and polarizable classical MD are barely noticeable, therefore in Fig. 2 we compared only structures from AIMD and classical non-polarizable simulations. In AIMD structures there is an increased probability of finding an O atom at the distance about 3 Å from the C atom of the imidazolium ring with the αCHO angle about 30 – 40 degrees. In classical NP simulations this probability maximum moves toward larger angles and longer distances. Finding an oxygen atom at short distance from the carbon atom of the 7

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

CH3 group is less probable than in the vicinity of imidazolium C atom, even in ab initio simulations, nevertheless, in AIMD structures there are some oxygen atoms found between 3 and 3.5 Å from the CH3 carbon atom and with the αCHO angle less than 60 degrees. Only few O atoms can be found within these limits in classical simulations.

Fig.2. Spatial distributions of H donor (X) – H acceptor (Y) pairs in structures obtained in ab initio and classical non-polarizable MD. Finding an F atom close to imidazolium carbon atom with small αCHF angle is much less probable than finding of an oxygen atom in similar location; this probability increases in classical NP simulations compared to AIMD (as already indicated by Him-F RDFs in Fig. 1).

8

ACS Paragon Plus Environment

Page 8 of 44

Page 9 of 44 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

On the other hand, there is little difference between distributions of fluorine atoms around CH3 groups obtained from ab initio and classical simulations. From this part of the analysis we can conclude that in AIMD structures most HBs involving imidazolium or CH3 hydrogen atoms will be formed to oxygen acceptor atom. In classical NP simulations number of such bonds will be greatly reduced, in particular in the case of CH3 donor atoms. Distribution of hydrogen bonds from CH3 hydrogens to fluorine atoms is expected to be much less affected by the simulation protocol and the probability of forming a HB between imidazolium H and fluorine atom, although small, will increase in classical MD structures.

3.2. Statistics of hydrogen bonds In order to check our simulated structures for the presence of hydrogen bonds we have to set the criteria of existence of a X—H···Y bond between donor X and acceptor atom Y. In this work we set the following thresholds for HB detection: 1) the X—Y distance less than or equal 3.5 Å and 2) the angle αXHY between XH and XY vectors less than or equal 40 degrees. These values were set based on previous analysis of hydrogen bonds parameters46,47 and are in agreement with geometries of EMIM-TFSI ion pairs calculated in Ref. 22.

Fig. 3. Distributions of H-Y distances in H-Y pairs selected according to the geometrical thresholds. Vertical broken lines mark the sum of the van der Waals radii of H and Y atoms. As stated in Ref. 47, setting the angle limit as geometrical criterion for the HB detection is problematic; therefore we focused on the condition that the distance between the proton and the acceptor atom must be shorter than the sum of their van der Waals radii47 and checked whether it is fulfilled when the above thresholds are used. In Fig. 3 we present histograms of H-Y distances between atom pairs selected according to our criterions; the H-O 9

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

distribution was decomposed into contributions arising from different types of hydrogen atoms (from the imidazolium ring, CH2 and CH3 groups, respectively). Maxima for H-F and H-N distributions appear between 2.3-2.4 Å and the maximum for H-O distribution is located at 2.2 Å. From its decomposition one may note that the distances between oxygen atoms and CH2/CH3 hydrogen atoms (maximum at about 2.4 Å) are significantly longer than distances between O and imidazolium H atoms (maximum of the distribution around 2.1 Å) indication larger strength of the latter bonds. The sum of H and Y Pauling radii48 equals 2.6, 2.55 and 2.7 Å for Y=O, F and N, respectively. It is readily noticeable that our thresholds for the X—Y distances and the αXHY angles lead to H-Y distances shorter than the sum of atomic radii, and therefore may serve as geometrical indicators of hydrogen bonds. We begin our analysis of the hydrogen bond network in EMIM-TFSI with the statistics of HBs obtained from AIMD and classical MD simulations. In Fig. 4 we show total number of HBs per pair of IL ions and its breakdown into contributions from HBs to different acceptor atoms. In AIMD structures there are about 3.6 hydrogen bonds per cation-anion pair, most of which are H-O bonds (2.6 bonds per pair). Bonds to F or N acceptors are less probable, with about 0.7 and 0.3 H-F and H-N bonds per ion pair, respectively.

Fig. 4. Average total number of HBs per cation and its breakdown according to different acceptor atoms. The number of hydrogen bonds is significantly smaller in classical MD simulations: the total number of bonds is 1.8 and 2.0 bonds per pair in non-polarizable and polarizable simulations, respectively. The difference to ab initio structures is mainly due to much smaller number of H-O bonds (0.8 – 1.0 bonds/pair). Classical simulations predict slightly more H-F bonds and less H-N bonds than ab initio MD. Inclusion of polarizability in classical force field via Drude oscillators has a rather limited effect shifting the results toward AIMD values by 0.1 – 0.2 bonds/ion pair. 10

ACS Paragon Plus Environment

Page 10 of 44

Page 11 of 44 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

Fig. 5. Average numbers of HBs per cation depending on the donor/acceptor atom and the type of MD used in simulations. In the next step we try to analyze the relative abundance of different hydrogen bonds presenting for each X—H···Y bond a breakdown into contributions arising from different hydrogen donor atoms X: carbon atoms of the imidazolium ring, carbon atoms from CH3 groups and carbon atom from CH2 unit (Fig. 5). In AIMD structures regardless of the donor atom most HBs are formed to oxygen atoms: on average there are 1.2, 1.0 and 0.4 bonds/pair from imidazolium, CH3 and CH2 carbon donor atom, respectively. Taking into account different number of imidazolium, CH3 and CH2 hydrogen atoms in the cation (3, 6 and 2,

11

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

Page 12 of 44

respectively) we can estimate that the probability of a H atom to form a HB is the largest for imidazolium hydrogens, and is about 2 times smaller for CH3 or CH2 hydrogen atoms. The average number of 1.22 HBs per imidazolium atom may be further decomposed into contributions of 0.47, 0.43 and 0.32 HBs per C2, C4 and C5 atom, respectively; in accord with expected larger acidity of the C2 atom (numbering of imidazolium atoms is shown in Fig. S1 in Supporting Information). Probability to form a HB to fluorine atom is similar for all types of hydrogen atoms in the EMIM cation, therefore most H-F bonds (0.4 bond per pair of ions) are formed by CH3 hydrogens. In the structures obtained from classical non-polarizable simulations (NP MD) practically there are no HBs between CH3 or CH2 hydrogen atoms and oxygen acceptor. The number of H-O bonds from imidazolium ring hydrogens is reduced to 0.7 bonds/ion pair. Compared to AIMD results there is small increase of the number of H-F bonds, mainly those involving imidazolium H atoms. In classical polarizable simulations the abundance of different types of hydrogen bonds is similar to non-polar MD; the most noticeable difference is the increase of the number of imidazolium H-O bonds by 0.2 to 0.9 bonds/pair. Data shown in Fig. 5 confirm the conclusions of the analysis of radial and spatial distribution of acceptor atoms (Figs. 1 and 2) that the main difference between ab initio and classical MD is decrease in the number of H-O bonds and the lack of H-O bonds from CH3 and CH2 groups in classical simulations.

Fig. 6. Schematic representation of bifurcated and chelating hydrogen bonds.

Table 1. Statistics of bifurcated and chelating HBs in EMIM-TFSI structures simulated via different kinds of MD. n(bifurcated), n(chelating), bonds/pair

bonds/pair

chelating HBs acceptors, in % O

F

N

AIMD

0.052

0.146

93.1

4.9

2.0

classical NP

0.037

0.022

30.9

65.7

3.4

classical polariz.

0.053

0.028

37.8

48.9

13.3

12

ACS Paragon Plus Environment

Page 13 of 44 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

Network of hydrogen bonds in the liquid may be further complicated due to possibility of formation of HBs with more than one hydrogen donor or acceptor, i.e. chelating or bifurcated bonds.1 Presence of such bonds in EMIM-TFSI is facilitated by the fact that several atoms of EMIM cations may serve as donors and, respectively, anions of this IL contain several atoms able to act as proton acceptors. In this work we considered two kinds of such bonds schematically depicted in Fig. 6. Bifurcated hydrogen bond is a pair of HBs from one X—H donor to two different acceptor atoms (from the same or different TFSI anions). The other type is the chelating bond – a pair of HBs to one acceptor atom from two hydrogen atoms within the same EMIM cation; the hydrogen atoms may be covalently bound to two different or to same donor atoms X.

Fig. 7. Percentage of HBs bifurcated to different acceptor pairs. 13

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

Information on the number of bifurcated and chelating bonds in the structures obtained from different types of MD simulations is collected in Table 1. As expected, chelating or bifurcated bonds are much less probable than single HBs. The largest number of chelating bonds was found in AIMD structures (0.15 bonds/ion pair). In more than 90% of these bonds the chelated acceptor Y is an oxygen atom. The number of chelating bonds is reduced to 0.02 bonds/pair in classical NP simulations and the most probable acceptor is the F atom (66 % of chelating bonds). The latter change originates from the fact that most of these bonds are formed between fluorine acceptor and two hydrogen atoms from within the same CH3 group. Polarizable classical simulations yield the ratio of O and F acceptors in chelating HBs similar to non-polarizable MD. Unlike two other parameterizations, percentage of chelating bonds to nitrogen atom is non-negligible only in classical polarizable simulations (13 % of chelating HBs). According to Table 1, bifurcated HBs are even less probable than chelating bonds. Although their number does not depend much on the simulation protocol (about 0.05 bonds/pair for AIMD or classical polarizable simulations), their structure significantly differs. In Fig. 7 we present detailed breakdown of the number of bifurcated HBs depending on the donor atom and the pair of acceptors. With O, F and N HB acceptor atoms present in the anion, there are five possible combinations of acceptor atoms (for steric reasons there are no HBs bifurcated to two nitrogen atoms). In AIMD structures the most probable is the bifurcation of a HB to the mixed pair of oxygen and nitrogen atoms as acceptors. There exist also bifurcated bonds to two oxygen or two fluorine atoms as well as to mixed pair of oxygen and fluorine acceptors. Bonding pattern is similar for imidazolium ring and CH3 hydrogen atoms. Conversely, in classical simulations bifurcated bonds involving O or N acceptor atom are formed only by hydrogen atoms from the imidazolium ring. Regardless of the donor atom X, most bifurcated HBs in classical MD is formed to two fluorine acceptors. Polarization in classical simulations increases the probability of ON or OO acceptor pairs for imidazolium hydrogens, nevertheless the bifurcated bonds still are predominantly to fluorine atoms.

3.3. Hydrogen bond network and spectral graph analysis Until now we discussed only the average numbers of different HBs in the ionic liquid. Now we turn our attention to the whole network of hydrogen bonds and its connectivity. We start with the analysis of the number of HBs formed by an EMIM cation (Fig. 8).

14

ACS Paragon Plus Environment

Page 14 of 44

Page 15 of 44 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

Fig. 8. Distribution of average numbers (total and to given acceptor atoms) of HBs formed by EMIM cations. As seen in the top panel of Fig. 8, in ab initio structures cations can form up to 6 - 7 hydrogen bonds. The distribution is rather flat and the most probable are 2 or 3 bonds. Probability that there is no hydrogen bond formed by the cation at given moment is less than 10 %. Looking at probability distributions of different acceptors one can see that this result is due to wide distribution obtained for H-O bonds. The most probable is the case when each cation forms one or two HBs to oxygen atoms. On the other hand, for F or N acceptors the

15

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

probability distribution has maximum at zero, showing that only few cations form one or sometimes two hydrogen bonds to fluorine or nitrogen atom. In classical simulations all probability distributions (the total and those for different acceptors) have maximum at zero and decay monotonically with increasing number of bonds. At given time about 40 % of cations are not involved in hydrogen bonding at all. Distribution of bonds for O and F acceptors is similar, only up to 20 % of cations form a hydrogen bond to such an acceptor atom.

Fig. 9. Representation of hydrogen bonds in EMIM-TFSI as simple graph or multigraph with corresponding adjacency matrices.

16

ACS Paragon Plus Environment

Page 16 of 44

Page 17 of 44 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

Some more insight into the connectivity of hydrogen bond network in simulated structures of EMIM-TFSI liquid may be obtained by means of graph theory; such a methodology was used recently in a study of ion aggregates in salt solutions.49 For this purpose we represent the liquid as an undirected graph consisting of a set of vertices (EMIM and TFSI ions) connected by edges (hydrogen bonds). We will consider simple graphs as well as more general multigraph representation. In a simple graph vertices are connected with single edges, thus we encode only the fact that there is at least one HB between given ions (vertices), without information on the actual number of HBs. In a multigraph the number of edges between vertices equals the number of hydrogen bonds between corresponding ions. A schematic example of two graph representations of hydrogen bonds between ions is shown in Fig. 9. In this work we analyze only HBs between EMIM+ and TFSI– ions, therefore our graphs are bipartite, with anions and cations forming two sets of vertices (white and black circles in Fig. 9). Connectivity of a graph with N vertices is represented as its N×N adjacency matrix A. Diagonal elements of adjacency matrix Aii are equal 0 and off-diagonal elements Aij equal n if there are n edges between vertices i and j; otherwise Aij = 0. For a simple graph representation the off-diagonal elements of the adjacency matrix are therefore either 0 or 1. Adjacency matrices corresponding to sample graphs are shown in Fig. 9. Normalized eigenvalue spectrum of the adjacency matrix yields the spectrum of the graph. Let us recall some properties of such spectra. For a set of disconnected vertices (an edgeless graph) the spectrum has only one eigenvalue λ = 0. If there is a single isolated pair of vertices connected by a single edge, a pair of eigenvalues λ = ±1 appears in addition to the singularity at 0. When the number of edges in the graph increases and they start to form connected network, the eigenvalue spectrum of the connectivity matrix widens and becomes dense. We should also note that the a pair of eigenvalues ±n can appear in the spectrum only if there is at least one vertex with degree (the total number of edges connecting it to other vertices) at least n , where n is the least integer greater than or equal n. Graph spectra for the hydrogen bond network in simulated IL structures treated either as simple graphs or multigraphs are displayed in Fig. 10. Simple graph representation will give information on the number of the counterions linked via at least one hydrogen bond to given ion, while multigraph representation supplements it by information on the actual number of bonds. We produced spectra for the graphs spanned by all HBs, as well as separate graphs constructed by different types of H—Y bonds.

17

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

Fig. 10. Spectra of graphs formed by HBs in EMIM-TFSI (for different types of HBs and for the whole network of all bonds). As one can see from the top panels of Fig. 10, the whole network of all kinds of HBs in AIMD structures yields a continuous graph spectrum spanning from -4.5 to 4.5. Therefore, there are some ions in the liquid which form hydrogen bonds to five other ions. From the multigraph spectrum we can deduce that (in agreement with Fig. 8) the maximum number of HBs formed by an ion reaches 7. It is also clear that the connectivity and continuity of the HB network is provided by H—O bonds. Some ions are connected to 4 other ions via HBs to oxygen atoms. Connectivity of the network of H—F bonds is smaller. From the discrete set of eigenvalues obtained for the graph of H—N bonds we conclude that HBs to nitrogen acceptors are isolated and do not form a continuous network.

18

ACS Paragon Plus Environment

Page 18 of 44

Page 19 of 44 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

The spectra based on the structures obtained from classical MD are narrower and with some visible gaps. The maximum number of counterions linked to given cation is reduced to 4. It may be also noted that there are much more disruptions in the network of H—O bonds and the total network of HBs is supported more-less equally by H—F and H—O bonds. The contribution of H—F bonds is larger in non-polarizable simulations; polarizability increases a little the connectivity of H—O bond network. The H—N hydrogen bonds are in classical simulations even more isolated than in ab initio structures as seen from few isolated peaks in the eigenvalue spectrum. Graph spectral analysis supports therefore the conclusions drawn in preceding sections.

3.4. Vibrational features of hydrogen bonds Presence of hydrogen bonds in the liquid is supposed to affect its vibrational spectrum. In EMIM-TFSI hydrogen donors are carbon atoms, therefore the relevant part of the IR or Raman spectrum is the region of C-H stretches at about 3000 cm-1. In our recent work25 we used AIMD simulations to reproduce the experimental IR spectrum of EMIM-TFSI liquid.22,39 The spectrum was calculated as the Fourier transform of the autocorrelation function of the total dipole moment of simulated system, therefore it corresponded to the spectrum of bulk liquid and we have not traced contributions from individual ions. In this study we want to look at selected C-H bonds in order to assess to what extent different environments of vibrating C-H groups affect the frequencies of individual transitions and their scatter along the MD trajectory.

Fig.11. Sample one-dimensional potential for proton motion during the stretching vibration of the C2-H bond with first two eigenvalues and eigenfunctions as well as the corresponding vibrational frequency.

19

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

To this end we constructed potential energy profiles for C—HIm stretches (as described in Sec. 2) and then we obtained corresponding frequencies of 0-1 vibrational transitions by numerically solving Schrödinger equation for calculated 1D potentials of proton displacements. Five EMIM cations from AIMD simulations were chosen for this purpose and only C—H stretches involving imidazolium C2, C4 and C5 carbon atoms were analyzed. An example of a proton potential for C2—H stretch is depicted in Fig.11; all potentials calculated for one selected ion pair are presented as movie in the Supporting Information. All proton potentials have single minimum and for small displacement their character is close to the harmonic potential. This kind of proton potentials is characteristic for typical moderate and week hydrogen bonds. No double minima potentials were observed.

Fig. 12. Frequencies of vibrational transitions calculated for C2-H bond stretches in five EMIM cations at different steps of the MD trajectory (bottom panel) and corresponding band contours calculated as superpositions of Gaussian functions with a half width of 50 cm-1 (top panel). The black line is the average of all five bands (colour lines). The calculated individual frequencies of vibrational transitions and overall band contours for each analyzed hydrogen bond involving H atom at C2 carbon atom are shown in Fig. 12; similar plots for C4—H and C5—H bonds are shown in Supporting Information.

20

ACS Paragon Plus Environment

Page 20 of 44

Page 21 of 44 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

Band contours were obtained as superpositions of Gaussian functions with half width of 50 cm-1. As readily seen in Fig. 12, frequency of a given C2—H stretch can vary between different frames of the trajectory by 200 – 300 cm-1. These changes reflect the changes in the local environment of EMIM cation in the timescale of MD simulation, i.e. about 40 ps. We may also note that for all bonds the highest calculated frequencies are very similar (about 3150 cm-1) whereas the low-energy end of the interval varies much more (2700 – 2900 cm-1). This feature is related to the shape of potentials at different snapshots of the trajectory: high frequencies are obtained in symmetric potentials, close to harmonic and this constitutes an upper limit for the vibrational frequency. On the other hand, increasing anharmonicity of the potential decreases the frequency of the vibration. Therefore different lower limits of calculated frequencies obtained for different bonds show that the asymmetry of the potential for the C—H stretch experienced along the MD trajectory by different anions may vary quite significantly.

Fig. 13. Averaged frequencies of C—H stretches calculated from 1D potentials for three imidazolium C—H bonds. Frequencies for all five analyzed bonds span the range from 2700 to 3150 cm–1. Differences between individual C2—H vibrations are also noticeable in simulated vibrational band. In the spectrum averaged over all considered bonds the maximum is located at about 3075 cm-1 and the maxima of contributing bands appear between 3020 and 3100 cm-1. The bands are asymmetric with longer tail on the low-frequency side. For C4—H and C5—H stretches the differences in frequencies calculated along the MD trajectory are smaller and the frequencies are on average shifted to slightly higher values. This can be seen in Fig. 13, where we compare averaged spectra for all three imidazolium hydrogen atoms. Although the average over only five anions may be statistically insufficient (but even for such a small number it was

21

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

computationally demanding), the C5—H and C4—H stretches are apparently shifted to energies 15-25 cm-1 higher than the C2—H vibration. This indicates that the HBs formed by hydrogen atoms at C4 and C5 imidazolium carbons are weaker than in the case of the hydrogen atom at C2, consistently with larger acidity of the C2 carbon. In the experimental spectrum of EMIM-TFSI liquid39 the C—H stretching bands of the imidazole ring are located at 3126 and 3165 cm-1. Our averaged band is at frequency slightly lower than experimental. To some extent it may be related to the DFT functional used in this work, which for some C-HIm stretches gives frequencies lower than other functionals (Cf. Table S1 in Supporting Information). Nevertheless, in the IR spectrum of the IL calculated in Ref. 25 from the same AIMD trajectories, the C—H stretching frequencies were slightly overestimated, thus the frequencies calculated for a single ion pair in gas phase apparently are not a good measure of the effect in a bulk liquid. On the other hand, it has been shown that the 1D model of the potential underestimates the frequencies comparing to the experimental values.44 We conclude therefore that the possible origin of the difference between IR spectrum obtained via Fourier transform for the AIMD structures and the frequencies obtained from the model potential can be attributed to the 1D potential used here, which certainly is a rather simplified description of the real energy surface for C—H bond involved in hydrogen bonding. Nevertheless, the approach proved its usefulness in resolving the overall spectrum into contributions due to individual bonds and assessing the range of frequency shifts resulting form changing environment of C—H oscillator.

4. Conclusions We performed extensive analysis of the hydrogen bond network in a typical aprotic IL, EMIM-TFSI, using quantum (AIMD) and classical MD simulations. MD results were analyzed in terms of spatial distribution of HB donor and acceptor atoms, statistics of different types of bonds and spectral graph analysis of connectivity of hydrogen bonds. In ab initio structures most HBs are formed between hydrogen atom from the imidazolium ring or the methyl group of EMIM cation and the oxygen atom of TFSI anion. The most probable pair of acceptors in bifurcated bonds is the pair of oxygen and nitrogen atom. EMIM cations may form HBs up to five different ions and the connectivity of hydrogen bond network is supported by H-O bonds. We showed that in structures of the liquid obtained from classical MD simulations pattern of HBs is essentially different. HBs with oxygen acceptor are formed only by hydrogen atoms from the imidazolium ring and the probability of fluorine being an acceptor 22

ACS Paragon Plus Environment

Page 22 of 44

Page 23 of 44 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

atom increases both in simple and in bifurcated bonds. Average number of hydrogen bonds is reduced and the network of HBs has more disruptions. In the polarizable force field the number of HBs with oxygen acceptor atom increases slightly, but the results are still far from those obtained in ab initio MD. We conclude that the simulation protocol used in MD may affect significantly the number of hydrogen bonds and the structure of their network. Appropriately parameterized classical MD simulations will continue to be useful tools in the investigations of hydrogen bonds in ionic liquids, owing to their relatively low computational cost making possible calculations for large systems and/or long timescales. Nevertheless, in the view of the results presented in this work, it is advisable to supplement such studies with ab initio MD simulations wherever possible, and to compare predictions of classical and ab initio approaches. Finally, we performed a posteriori quantum chemical calculations of 1D potential energy profiles for hydrogen motions in selected HBs obtaining information on differences between individual bonds as well as on fluctuations of the potential and resulting scatter of vibrational transition energies during the simulation. Results have proven applicability of this methodology to investigate vibrational frequencies of individual bonds and to assess to which extent they may differ. Results presented in this work provide some insights into the structure of hydrogen bond network in EMIM-TFSI liquid at atomistic resolution. We believe that methodological findings regarding the type of molecular dynamics (ab initio vs. classical) and usefulness of frequency calculations based on potential energy profiles extracted from the snapshots of the MD trajectory are more general and may applied to similar systems.

Acknowledgements This research was supported by statutory funding from the Faculty of Chemistry of the Jagiellonian University. M.Z.B. thanks the Foundation for Polish Science (FNP) for support. PL-Grid infrastructure was used in computations.

Supporting Information The Supporting Information is available free of charge on the ACS Publications website. Numbering of EMIM+ atoms, CIm-H stretching frequencies for EMIM-TFSI pair calculated in different DFT functionals, frequencies of vibrational transitions calculated for C4-H and C5-H bond stretches from the 1D potential, animation of 1D potential curves

23

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

calculated for selected C2-H bond, comparison of data obtained for two system sizes used in classical MD simulations.

24

ACS Paragon Plus Environment

Page 24 of 44

Page 25 of 44 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

References (1) Hunt, P. A.; Ashworth, C. R.; Matthews, R. P. Hydrogen Bonding in Ionic Liquids. Chem. Soc. Rev. 2015, 44, 1257-1288. (2) Hayes, R.; Warr, G. G.; Atkin, R. Structure and Nanostructure in Ionic Liquids. Chem. Rev. 2015, 115, 6357-6426. (3) Lehmann, S. B. C.; Roatsch, M.; Schöppke, M.; Kirchner, B. On the Physical Origin of the Cation–Anion Intermediate Bond in Ionic Liquids Part I. Placing a (Weak) Hydrogen Bond between Two Charges. Phys. Chem. Chem. Phys. 2010, 12, 7473-7486. (4) Dong, K.; Song, Y.; Liu, X.; Cheng, W.; Yao, X.; Zhang, S. Understanding Structures and Hydrogen Bonds of Ionic Liquids at the Electronic Level. J. Phys. Chem. B 2012, 116, 1007-1017. (5) Matthews, R. P.; Welton, T.; Hunt, P.A. Competitive Pi Interactions and Hydrogen Bonding within Imidazolium Ionic Liquids. Phys. Chem. Chem. Phys. 2014, 16, 3238-3253. (6) Katsyuba, S. A.; Vener, M. V.; Zvereva, E. E.; Fei, Z.; Scopelliti, R.; Laurenczy, G.; Yan, N.; Paunescu, E.; Dyson, P. J. How Strong Is Hydrogen Bonding in Ionic Liquids? Combined X‐ray Crystallographic, Infrared/Raman Spectroscopic, and Density Functional Theory Study. J. Phys. Chem. B 2013, 117, 9094-9105. (7) Housaindokht, M. R.; Hosseini, H. E.; Googheri, M. S. S.; Monhemi, H.; Najafabadi, R. I.; Ashraf, N.; Gholizadeh, M. Hydrogen Bonding Investigation in 1-Ethyl-3methylimidazolium Based Ionic Liquids from Density Functional Theory and Atoms-inMolecules Methods. J. Mol. Liq. 2013, 177, 94-101. (8) Katsyuba, S. A.; Vener, M. V.; Zvereva, E. E.; Fei, Z.; Scopelliti, R.; Brandenburg, J. G.; Siankevich, S.; Dyson, P. J. Quantification of Conventional and Nonconventional ChargeAssisted Hydrogen Bonds in the Condensed and Gas Phases. J. Phys. Chem. Lett. 2015, 6, 4431-4436. (9) Izgorodina, E. I.; MacFarlane, D. R. Nature of Hydrogen Bonding in Charged HydrogenBonded Complexes and Imidazolium-Based Ionic Liquids. J. Phys. Chem. B 2011, 115, 14659-14667. (10) Matthews, R. P.; Welton, T.; Hunt, P. A. Hydrogen Bonding and Pi–pi Interactions in Imidazolium-chloride Ionic Liquid Clusters. Phys. Chem. Chem. Phys., 2015 17, 1443714453.

25

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

(11) Marekha, B. A.; Kalugin, O. N.; Idrissi, A. Non-covalent Interactions in Ionic Liquid Ion Pairs and Ion Pair Dimers: a Quantum Chemical Calculation Analysis. Phys. Chem. Chem. Phys. 2015, 17, 16846-16857. (12) Kohagen, M.; Brehm, M.; Lingscheid, Y.; Giernoth, R.; Sangoro, J.; Kremer, F.; Naumov, S.; Iacob, C.; Kärger, J.; Valiullin, R.; Kirchner, B. How Hydrogen Bonds Influence the Mobility of Imidazolium-Based Ionic Liquids. A Combined Theoretical and Experimental Study of 1-n-Butyl-3-methylimidazolium Bromide. J. Phys. Chem. B 2011, 115, 15280-15288. (13) Skarmoutsos, I.; Dellis, D.; Matthews, R. P.; Welton, T.; Hunt, P. A. Hydrogen Bonding in 1-Butyl- and 1-Ethyl-3-methylimidazolium Chloride Ionic Liquids. J. Phys. Chem. B 2012, 116, 4921-4983. (14) Mondal, A.; Balasubramanian, S. Vibrational Signatures of Cation−Anion Hydrogen Bonding in Ionic Liquids: A Periodic Density Functional Theory and Molecular Dynamics Study. J. Phys. Chem. B 2015, 119, 1994-2002. (15) Paschek, D.; Golub, B.; Ludwig, R. Hydrogen Bonding in a Mixture of Protic Ionic Liquids: a Molecular Dynamics Simulation Study. Phys. Chem. Chem. Phys. 2015, 17, 8431-8440. (16) Zentel, T.; Kühn, O. Hydrogen Bonding in the Protic Ionic Liquid Triethylammonium Nitrate Explored by Density Functional Tight Binding Simulations. J. Chem. Phys. 2016, 145, 234504. (17) Campetella, M.; Montagna, M.; Gontrani, L.; Scarpellini, E.; Bodo, E. Unexpected Proton Mobility in the Bulk Phase of Cholinium-based Ionic Liquids: New Insights from Theoretical Calculations. Phys. Chem. Chem. Phys. 2017, 19, 11869-11880. (18) Fumino, K.; Peppel T.; Geppert-Rybczyńska, M.; Zaitsau, D. H.; Lehmann, J. K.; Verevkin, S. P.; Köckerling, M.; Ludwig, R. The Influence of Hydrogen Bonding on the Physical Properties of Ionic Liquids. Phys. Chem. Chem. Phys. 2011, 13, 14064-14075. (19) Gu, Q.; Su, P.; Xia, Y.; Yang, Z.; Trindle, C. O.; Knee, J. L. Quantitative Probing of Subtle Interactions Among H-bonds in Alpha Hydroxy Carboxylic Acid Complexes. Phys. Chem. Chem. Phys. 2017, 19, 24399-24411. (20) Mondal, S.; Biswas, B.; Nandya, T.; Singh, P. C. Hydrophobic Fluorine Mediated Switching of the Hydrogen Bonding Site as well as Orientation of Water Molecules in the Aqueous Mixture of Monofluoroethanol: IR, Molecular Dynamics and Quantum Chemical Studies. Phys. Chem. Chem. Phys. 2017, 19, 24667-24677.

26

ACS Paragon Plus Environment

Page 26 of 44

Page 27 of 44 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

(21) Paschoal, V. H.; Faria, L. F. O.; Ribeiro, M. C. C. Vibrational Spectroscopy of Ionic Liquids. Chem. Rev. 2017, 117, 7053-7112. (22) Dhumal, N. R.; Noack, K.; Kiefer, J.; Kim, H. J. Molecular Structure and Interactions in the Ionic Liquid 1-Ethyl-3-methylimidazolium Bis(Trifluoromethylsulfonyl)imide. J. Phys. Chem. A 2014, 118, 2547-2557. (23) Fujii, K.; Soejima, Y.; Kyoshoin, Y.; Fukuda, S.; Kanzaki, R.; Umebayashi, Y.; Yamaguchi, T.; Ishiguro, S.; Takamuku, T. Liquid Structure of Room-Temperature Ionic Liquid, 1-Ethyl-3-methylimidazolium Bis-(trifluoromethanesulfonyl) Imide. J. Phys. Chem. B 2008, 112, 4329-4336. (24) Kempter, V.; Kirchner, B. The Role of Hydrogen Atoms in Interactions Involving Imidazolium-Based Ionic Liquids. J. Mol. Struct. 2010, 972, 22-34. (25) Eilmes, A.; Kubisiak, P.; Brela, M. Explicit Solvent Modeling of IR and UV−Vis Spectra of 1-Ethyl-3-methylimidazolium Bis(trifluoromethylsulfonyl)imide Ionic Liquid. J. Phys. Chem. B 2016, 120, 11026-11034. (26) Phillips, J. C.; Braun, R.; Wang, W.; Gumbart, J.; Tajkhorshid, E.; Villa, E.; Chipot, Ch.; Skeel, R. D.; Kalé, L.; Schulten, K. Scalable Molecular Dynamics with NAMD. J. Comput. Chem. 2005, 26, 1781-1802. (27) 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. (28) Canongia Lopes, J. N.; Deschamps, J.; Pádua, A. A. H. Modeling Ionic Liquids Using a Systematic All-Atom Force Field. J. Phys. Chem. B 2004, 108, 2038-2047 (29) Canongia Lopes, J. N.; Pádua, A. A. H. Molecular Force Field for Ionic Liquids Composed of Triflate or Bistriflylimide Anions. J. Phys. Chem. B 2004, 108, 16893-16898. (30) Canongia-Lopes, J. N.; Deschamps, J.; Pádua, A. A. H. Correction: Modeling Ionic Liquids Using a Systematic All-Atom Force Field. J. Phys. Chem. B 2004, 108, 1125011250 (31) Köddermann, T.; Paschek, D.; Ludwig, R. Molecular Dynamic Simulations of Ionic Liquids: A Reliable Description of Structure, Thermodynamics and Dynamics. ChemPhysChem 2007, 8, 2464-2470. (32) Lamoreux, G.; Roux, B. Modeling Induced Polarization with Classical Drude Oscillators: Theory and Molecular Dynamics Simulation Algorithm. J. Chem. Phys. 2003, 119, 30253039.

27

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

(33) Borodin, O. Polarizable Force Field Development and Molecular Dynamics Simulations of Ionic Liquids. J. Phys. Chem. B 2009, 113, 11463-11478. (34) Noda, A.; Hayamizu, K.; Watanabe, M. Pulsed-Gradient Spin−Echo 1H and 19F NMR Ionic Diffusion Coefficient, Viscosity, and Ionic Conductivity of Non-Chloroaluminate Room Temperature Ionic Liquids. J. Phys. Chem. B 2001, 105, 4603−4610. (35) Darden, T.; York, D.; Pedersen, L. Particle Mesh Ewald: An Nlog(N) Method for Ewald Sums in Large Systems. J. Chem. Phys. 1993, 98, 10089−10092. (36) VandeVondele, J.; Krack, M.; Mohamed, F.; Parrinello, M.; Chassaing, T.; Hutter, J. QUICKSTEP: Fast and Accurate Density Functional Calculations Using a Mixed Gaussian and Plane Waves Approach. Comp. Phys. Comm. 2005, 167, 103-128. (37) Hutter, J.; Iannuzzi, M.; Schiffmann, F.; VandeVodele, J. CP2K: Atomistic Simulations of Condensed Matter Systems. WIREs Comput. Mol. Sci. 2014, 4, 15–25. (38) Goedecker, S.; Teter, M.; Hutter, J. Separable Dual-Space Gaussian Pseudopotentials. Phys. Rev. B 1996, 54, 1703-1710. (39) Vyas, S.; Dreyer, C.; Slingsby, J.; Bicknase, D.; Porter, J. M.; Maupin, C. M. Electronic Structure and Spectroscopic Analysis of 1-Ethyl-3-methylimidazolium Bis(trifluoromethylsulfonyl)imide Ion Pair. J. Phys. Chem. A 2014, 118, 6873-6882. (40) Pirc, G.; Mavri, J.; Stare, J. Program Package for Numerical Solving Time-Independent Schrödinger Equation for Vibrational Problems: Inclusion of Coordinate Dependent Reduced Masses. Vib. Spectrosc. 2012, 58, 153-162. (41) Stare, J.; Mavri, J.; Grdadolnik, J.; Zidar, J.; Maksić, Z. B.; Vianello, R. Hydrogen Bond Dynamics of Histamine Monocation in Aqueous Solution: Car-Parrinello Molecular Dynamics and Vibrational Spectroscopy Study. J. Phys. Chem. B 2011, 115, 5999-6010. (42) Brela, M. Z.; Wójcik, M. J.; Boczar, M.; Hashim, R. Car-Parrinello Simulation of the Vibrational Spectra of Strong Hydrogen Bonds with Isotopic Substitution Effects: Application to Oxalic Acid Dehydrate. Chem. Phys. Lett. 2013, 558, 88–-92. (43) Brela, M.; Stare, J.; Pirc, G.; Sollner Dolenc, M.; Boczar, M.; Wójcik, M. J.; Mavri, J. Car-Parrinello Simulation of the Vibrational Spectrum of a Medium Strong Hydrogen Bond by Two-Dimensional Quantization of the Nuclear Motion: Application to 2-hydroxy-5nitrobenzamide. J. Phys. Chem. B 2012, 116, 4510– 4518. (44) Brela, M. Z.; Wójcik, M. J.; Boczar, M.; Witek, Ł.;Yasuda, M.; Ozaki, Y. CarParrinello Molecular Dynamics Simulations of Infrared Spectra of Crystalline Vitamin C

28

ACS Paragon Plus Environment

Page 28 of 44

Page 29 of 44 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

with Analysis of Double Minimum Proton Potentials for Medium-Strong Hydrogen Bonds. J. Phys. Chem. B 2015, 119, 7922-7930. (45) Porter, A. R.; Liem, S. Y.; Popelier, P. L. A. Room Temperature Ionic Liquids Containing Low Water Concentrations - A Molecular Dynamics Study. Phys. Chem. Chem. Phys. 2008, 10, 4240-4248. (46) Desiraju, G. R.; Steiner, T. The Weak Hydrogen Bond in Structural Chemistry and Biology; Oxford University Press: New York, 1999. (47) Grabowski, S. J. Hydrogen Bonding Strength - Measures Based on Geometric and Topological Parameters. J. Phys. Org. Chem. 2004, 17, 18–31 (48) Batsanov, S. S. Van der Waals Radii of Elements. Inorg. Materials 2001, 37, 871-885. (49) Choi, J.-H.; Cho, M. Ion Aggregation in High Salt Solutions. VI. Spectral Graph Analysis of Chaotropic Ion Aggregates. J. Chem. Phys. 2016, 145, 174501.

29

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

TOC graphic

30

ACS Paragon Plus Environment

Page 30 of 44

Page 31 of 44 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

55x37mm (300 x 300 DPI)

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

156x316mm (300 x 300 DPI)

ACS Paragon Plus Environment

Page 32 of 44

Page 33 of 44 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

59x46mm (300 x 300 DPI)

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

55x39mm (300 x 300 DPI)

ACS Paragon Plus Environment

Page 34 of 44

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

151x318mm (300 x 300 DPI)

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

30x12mm (300 x 300 DPI)

ACS Paragon Plus Environment

Page 36 of 44

Page 37 of 44 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

151x318mm (300 x 300 DPI)

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

169x358mm (300 x 300 DPI)

ACS Paragon Plus Environment

Page 38 of 44

Page 39 of 44 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

173x389mm (300 x 300 DPI)

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

143x262mm (300 x 300 DPI)

ACS Paragon Plus Environment

Page 40 of 44

Page 41 of 44 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

65x51mm (300 x 300 DPI)

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

121x178mm (300 x 300 DPI)

ACS Paragon Plus Environment

Page 42 of 44

Page 43 of 44 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

57x44mm (300 x 300 DPI)

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

TOC graphic 39x19mm (300 x 300 DPI)

ACS Paragon Plus Environment

Page 44 of 44