Hydrogen Bonding in Protic Ionic Liquids: Structural Correlations, Vibrational Spectroscopy, and Rotational Dynamics of Liquid Ethylammonium Nitrate
Tobias Zentel Viviane Overbeck Dirk Michalik Oliver Kühn Ralf Ludwig
physics.chem-ph


Hydrogen Bonding in Protic Ionic Liquids]Hydrogen Bonding in Protic Ionic Liquids: Structural Correlations, Vibrational Spectroscopy, and  Rotational Dynamics of Liquid Ethylammonium Nitrate



^1Institut für Physik, Universität Rostock, Albert-Einstein-Str. 23-24, D-18059 Rostock, Germany
^2Abteilung für Physikalische Chemie, Institut für Chemie, Universität Rostock, Dr.-Lorenz-Weg 1, D-18059 Rostock, Germany

oliver.kuehn@uni-rostock.de 



  July 26, 2017




The properties of the hydrogen bonds in ethylammonium nitrate are analyzed by using molecular dynamics simulations and infrared as well as nuclear magnetic resonance experiments.  Ethylammonium nitrate features a flexible three-dimensional network of hydrogen bonds with moderate strengths, which makes it distinct from related triethylammonium-based ionic liquids. First, the network's flexibility is manifested in a not very pronounced correlation of the hydrogen bond geometries, which is caused by rapid interchanges of bonding partners. The large flexibility of the network leads to a substantial broadening of the mid-IR absorption band, with the contributions due to N-H stretching motions ranging from 2800 to 3250 . Finally, the different dynamics are also seen in the rotational correlation of the N-H bond vector, where a correlation time as short as 16.1 ps is observed.



























§ INTRODUCTION


Ethylammonium nitrate (EAN) has been synthesized more than hundred years ago by Paul Walden <cit.>. EAN is formed by the neutralization of ethylamine with concentrated nitric acid. It consists solely of ions, namely ethylammonium and nitrate, and melts significantly below 100^∘C, i.e. at 13-14^∘C. Thus, EAN can be regarded as the first room temperature ionic liquid (IL) <cit.>. EAN belongs to the subset of protic ionic liquids (PILs) <cit.>, which are of practical interest due to their ongoing development for  high-performance Li-ion cells in electric vehicles, instead of the flammable carbonate electrolytes <cit.>, as well as electrolytes for advanced fuel cells <cit.>. Recent studies mainly addressed the structure and thermodynamics of EAN. This archetypical PIL forms a three-dimensional hydrogen bonded (H-bonded) network, similar to that of water. In principle, these structural features were supported by molecular dynamics (MD) simulations <cit.>. In addition, symmetric and asymmetric stretching modes of the H-bonds (HBs) were observed in far infrared (FIR) spectra <cit.>.  EAN evaporates as contact ion pairs leading to low vapor pressures and an enthalpy of vaporization of about 105.3 kJmol^-1 <cit.>. This stability makes EAN attractive as an electrolyte for new generation fuel cells up to 419 K. In contrast, EAN evaporates as molecules near the boiling point at 513 K. Although structure and thermodynamics of EAN are widely studied, a detailed analysis of the HB network and the related dynamical information is still lacking. This is somewhat surprising because the “flickering HB network”, similar to that of water, suggests strongly enhanced dynamics of molecular vectors that are involved in H-bonding compared to those in the one-dimensional H-bonded PIL analogues <cit.>.       

For EAN, wherein the proton exchange is fast compared to the time scale of nuclear magnetic resonance (NMR) spectroscopy, the reorientational correlation times for the N-H molecular vector, τ_ NH, can  only be obtained from ^15N-enhanced proton relaxation time experiments. To avoid the demanding and costly synthesis of isotopic substituted compounds, one can instead measure the deuteron relaxation times of the N-D deuterons. This relaxation mechanism is strong and purely intramolecular in nature. However, proper correlation times, τ_ NH, can be obtained from the relaxation times, T_1, only, when reliable deuteron quadrupole coupling constants, χ_ D, are known. Unfortunately, this is not the case for ILs, neither for the solid, the liquid nor the gas phase. For this reason we recently developed a method for deriving χ_ D from a relation between density functional theory (DFT) calculated χ_ D (for N-D) and proton chemical shifts δ^1H (for N-H) within clusters of ILs <cit.>. Thus, with simple measurements of the proton chemical shifts,  one can then determine accurate χ_ D values for the liquid phase.  

MD simulations  are an established tool to study ILs <cit.>. To get a reliable description of the complex atomic motions intrinsic to HBs, methods beyond simple force field-based MD  are needed. In principle, DFT-based ab initio MD would be the method of choice, but it is hampered by the high computational costs <cit.>. An alternative is provided by the self-consistent charge DFT-based tight-binding approach (DFTB) <cit.>. Although applications of DFTB to the dynamics of ILs are scarce, available simulations provided evidence, e.g., for reliable prediction of  X-ray diffraction data <cit.> and vibrational spectra as well as geometric correlations related to the HB motion of PILs <cit.>. 

In the present work, we continue our investigation of the molecular dynamics of HBs in alkylammonium nitrate PILs <cit.> by focussing on liquid EAN. Using a combined experimental and theoretical approach, the N-H⋯O HBs will be analyzed in terms of their geometric correlations,   FIR and mid infrared (MIR) signatures as well as the rotational correlation of the N-H bond  vectors. It is the purpose of this work to provide detailed information about the structure and dynamics of the very first IL EAN, with particular focus on the HB network, which resembles that of water. Section <ref> starts with an outline of the DFTB MD simulation protocol. Next, the method for obtaining the reorientational correlation time, τ_ NH, for N-H/D bonds in EAN is sketched, which is based on the δ^1H  versus χ_ D relation calculated for EAN clusters. Further, experimental details for NMR and FIR/MIR measurements are provided. In Section <ref> results for   structure, infrared (IR) spectra, and  rotational correlation time are discussed. A summary is provided in Section <ref>.



§ METHODS




 §.§ Molecular Dynamics


A simulation box consisting of 60 EAN ion pairs and 900 atoms with periodic boundary conditions was built from 10 randomly positioned hexamer clusters. 
To obtain a bulk structure, the system was equilibrated  using the force field parameters from <cit.>. To converge the box size a  1.5 ns simulation in the NPT ensemble at 300 K was performed using a  Parrinello-Rahman barostat. The resulting density of 1.244 g/cm^3 is very close to the literature value of 1.216 g/cm^3 <cit.>. To control the temperature in the subsequent NVT simulation (at 300 K), a Nosé-Hover thermostat with chain length 10 is used. The length of the NVT trajectory has been  5 ps. The force field-based MD simulations were performed with the Gromacs 4.5.5 software <cit.>.

The final equilibrated force field structure  was used as input for subsequent DFTB simulations with the DFTB+ software package <cit.>. Here, 3rd order as well as dispersion  corrections have been used together with the  3-ob  Slater-Koster parameter set <cit.>.  A 25 ps NVT trajectory has been simulated (time step 0.5 fs), which was used to sample starting points  for NVE trajectories every 5 ps. The NVE trajectories had a length of 30 ps with a time step of 0.5 ps.

The IR spectrum, I(ω), is obtained by Fourier transformation of  the dipole moment trajectory, μ⃗(t), using a standard Kaiser window function, κ(t), with Kaiser parameter 10,  i.e. <cit.>
  
    I(ω) = ω^2 ⟨ | ∫_0^Texp(-i ω t)   μ⃗(t) κ(t)   dt  |^2 ⟩  .

Here, the canonical ensemble average denoted by < ⋯> is obtained by sampling initial configurations for  NVE trajectories as described above. The window function ensures that the dipole is damped to zero at the edges. The Kaiser parameter controls the time scale of the decay and is chosen large enough to avoid artificial broadening of  the spectra.
The Mulliken charges from DFTB simulations are used to calculate the dipoles of the individual ions and  of the NH_3 and CH_3 groups. The total box dipole is obtained as the sum of the dipoles of the individual ions.
Note that the  dipole of an ionic system is not uniquely defined, as discussed in Ref. <cit.>. Here, we employed the center of mass of the considered subsystem as the reference for the dipole calculations.
The spectra presented below are smoothed by convolution with a Gaussian function with width parameter σ= 25  cm^-1 in the MIR and σ= 12.5  cm^-1 in the FIR region.




 §.§ Reorientational Correlation Time


For the determination of the reorientational correlation time of the N-H molecular vectors, τ_ NH,  we measured the deuteron relaxation time (T_1)_ D of the deuterated species [EtND_3][NO_3] at 303 K. We have chosen quadrupolar relaxation because it is the strongest relaxation mechanism and presents a purely intramolecular process. Deuteron nuclear magnetic relaxation is driven by the interaction of the electrostatic quadrupole moment, eQ, of the deuteron nucleus with the main component of the electric field gradient (EFG) tensor at the nucleus, q_zz, generated by the electron distribution surrounding the nucleus along the N-H bonds. The relaxation rate (1/T_1)_ D is given by <cit.>

    (1/T_1)_ D=3/10π^2 ( 1+ η_ D^2/3)χ_ D^2( τ_ NH/1+ω_0^2τ_ NH^2 +4τ_ NH/1+4ω_0^2τ_ NH^2)   ,

where χ_ D is the deuteron nuclear quadrupole coupling constant and η_ D=(q_xx-q_yy)/q_zz is the related asymmetry parameter. If the relaxation process does not depend on frequency and the extreme narrowing condition, ω_0τ_ NH≪ 1, holds, this equation simplifies to


    (1/T_1)_ D=3/2π^2 ( 1+ η_ D^2/3)χ_ D^2 τ_ NH  .

Note that we tacitly assumed that the reorientational correlation time of the principal axis of the deuterium EFG, actually entering these equations, is identical to τ_ NH <cit.>. Additionally, experimental data indicate that there is a negligible isotope effect for the rotational diffusion <cit.>.

Before the reorientational correlation times  can be determined from the quadrupolar relaxation rates, the deuteron quadrupole coupling constant  must be known. Despite what the term `coupling constant' implies, χ_ D has been shown to be temperature and solvent dependent 
<cit.>.
Thus χ_ D is a sensitive probe for H-bonding. Unfortunately, for EAN the deuteron quadrupole coupling constant is not known at all. Thus, we used an approach for the determination of the proper deuteron quadrupole coupling constant in the liquid phase, which has been successfully applied for N-D coupling constants in other PILs   and is briefly described here <cit.>.

The relation between the principal component of the EFG tensor, eq_zz, and the deuteron nuclear quadrupole coupling constant, χ_ D, is given by (h is Planck's constant)


    χ_ D=(eQeq_zz/h)   .


In principle, the deuteron quadrupole coupling constant  can be now obtained by multiplying calculated EFG  by a calibrated nuclear quadrupole moment, eQ. The latter  is obtained by plotting the measured gas phase quadrupole coupling constants from microwave spectroscopy versus the calculated EFGs for small molecules such as CD_4, CD_3OH, ND_3 etc. as described by Huber et al.
<cit.>.
 The slope gives a reasonable value of eQ =295.5 fm^2. This approach holds regardless of whether gas phase molecules, H-bonded clusters or IL complexes are investigated. 

Recently, we could show for a large set of DFT-calculated clusters of trialkylammonium-based PILs that there exists a linear correlation between the calculated proton chemical shifts  and the calculated  deuteron quadrupole coupling constants <cit.>. The advantage now is that δ^1H can be easily measured in the liquid phase, and thus, by virtue of the computationally established linear correlation, provides access to  the  χ_ D values. Thus, we applied this approach to the present EAN. First, we obtained the proton chemical shifts, δ^1H, and the quadrupole coupling constants, χ_ D, from DFT calculated clusters including n=2, 4, 6, 8 ion-pairs. 
Throughout, the geometries of the clusters as well as the proton chemical shifts and the EFGs  of the deuterons were calculated using Gaussian 09 <cit.> at the B3LYP-D3/6-31+G* level of theory, including D3 dispersion correction as introduced by Grimme et al. <cit.>.
 The proton chemical shifts were referenced against  TMS (tetramethylsilane), as it was done in the experiment (see below). The χ_ D values were derived by multiplying the calculated main components of the EFG tensor q_zz with the calibrated nuclear quadrupole moment eQ =295.5 fm^2. 
The asymmetry parameter of the EFG  for N-D, η_ D, was found to be negligible and did not have to be considered in (<ref>). This has been shown for molecular systems such as ammonia, formamide, N-methyl formamide and N-methyl acetamide 
<cit.>.
 
 The reorientational correlation times, τ_ NH, thus obtained from (<ref>) will be compared with   simulation results. To this end, the rotational time correlation function, C(t), is computed according to

    C(t) = < P_2( u⃗(0) u⃗(t)) >    .


Here, P_2 is the second Legendre polynomial and u⃗(t) the unit vector along the N-H bond <cit.>. 

From C(t) the rotational correlation time  follows as

    τ_ NH = ∫_0^∞< P_2( u⃗(0) u⃗(t) > dt   .



 


 §.§ Experimental Methods


The  EAN samples were of commercial origin (IOLITEC) with purity degrees of  >97%. Prior to experiments, the sample was subjected to vacuum evaporation at 333 K for more than 72 hours to remove possible traces of solvents and moisture. Once purified, the sample was transferred under inert argon atmosphere and stored in a hermetically sealed bottle. The purity of the EAN sample was confirmed with ^1H and ^13C NMR (see Supplement). The residual water concentration in the sample for combustion calorimetry was determined by Karl Fischer titration and was lower than 100 ppm. The ^1H NMR spectra of the protonated EAN were recorded on a Bruker Avance 500 MHz spectrometer using 5 mm probes. The δ^1H chemical shifts were measured versus TMS.  The EAN sample was deuterated by H/D exchange in D_2O and dried again. Longitudinal magnetic relaxation times T_1 were measured with the same spectrometer at a resonance frequency of ν_0=ω_0/2π= 76.7 MHz employing the inversion recovery (180^∘ - τ - 90^∘) pulse sequence <cit.> (see Supplement). T_1 is estimated to be accurate to within ±2%. Temperature calibrations were carried out using methanol and  ethylene glycol NMR thermometers.	

The  FIR measurements were performed with a Bruker Vertex 70 FTIR spectrometer. The instrument was equipped with an extension for measurements in this spectral region. This equipment consists of a multilayer mylar beam splitter, a room temperature DLATGS detector with preamplifier and polyethylene   windows for the internal optical path. The accessible spectral region for this configuration lies between 30 and 680 . The spectrometer was purged continuously with dry air during the experiments in order to minimize contributions from atmospheric water vapor. Further reduction of these signals was achieved by utilization of telescope tubes with polyethylene windows reducing the optical path in the sample chamber to a minimum. In this way, fluctuations in the atmospheric water content can be prevented. The L.O.T.-Oriel cell was equipped with polyethylene windows. The 0.012 mm path length was realized by tin spacers. For each spectrum 100 scans at a spectral resolution of 1  were recorded. The experimental spectrum shown in Figure <ref> has been published elsewhere <cit.>. 
The MIR measurements were performed with a Bruker Vector 22 FTIR spectrometer. An L.O.T.-Oriel variable-temperature cell equipped with CaF_2 windows was used. For each spectrum 128 scans were recorded at a spectral resolution of 1 .



§ RESULTS AND DISCUSSION




 §.§ H-Bond Geometries

















EAN facilitates multiple HB opportunities and thus is capable of forming a complex network structure, which has been discussed to be comparable to that of water <cit.>. There exist three acceptor oxygen atoms per anion and three donor positions per cation, so that in principle a fully saturated HB network is possible. Note that the trajectory simulations discussed below did not exhibit any proton transfer events leading to a stable product of two neutral molecules. This might be related to the relatively short propagation time.



Nitrate is a strongly interacting anion and we assume all HB donors to be involved in H-bonding. In the following a HB is counted always between an ammonium hydrogen and its closest oxygen.  The average HB distances are r_NO=2.91 Å, r_NH=1.04 Å, and    r_OH=2.02 Å. These values are in the range  expected  for moderate strong HBs <cit.>. Further, EAN exhibits a bent HB geometry with an average angle, α, between the NH and OH vector of 34^∘ (see also figure <ref>). Only 16% of the HBs are linear (α <15^∘). Hayes et al. investigated the structure of various PILs including EAN with Empirical Potential Structure Refinement (EPSR) and Monte Carlo methods <cit.>. They found 12% of the HBs being linear (α <15^∘) and 1.6  oxygen atoms per ammonium hydrogen. In total this classifies the  HBs in EAN as being substantially bifurcated. 
 
The continuous lifetime of a H-bonded ion pair is   (average lifetime of non-interrupted HBs) and the intermittent lifetime is  (total average lifetime along trajectories). It is not unusual that the two lifetimes do not match <cit.>. However, it should be noted that in the current analysis only one HB per ammonium hydrogen is possible. Thus, thermal motion of a bifurcated HB might, according to the distance criterion, be counted as bond breaking and reforming. In this case, the HB lifetimes would be  underestimated.




The HBs of one ammonium group are almost exclusively formed with three different anions (less than 1% of the HBs of one ammonium group are to the same counter anion). In 27% of the cases an oxygen atom is acceptor of two HBs, leaving some of the oxygen atoms unpaired. Furthermore, there exist a number of ring-like structures, where two cations are  H-bonded to the same two anions, an exemplary structure is shown in figure <ref>. This can be either via different oxygen atoms of the same anion or via a single oxygen involved in two HBs (as shown). In the present bulk MD simulations containing 60 ion pairs, on average 18.4 rings exist, consisting of 73.6 ions. Since there are three HB possibilities per ion  it is possible for ions to be involved in more than one ring. Especially interesting are 'double rings' that form a HB network cluster consisting of two cations that share all three anions. On average 2.6 cluster structures are found. 
We can conclude that the  ring-like geometries appear to be  an abundant structural motif in EAN and can be taken as evidence for the generally assumed three-dimensional network character of the HBs in EAN.




In order to quantify geometric correlations of the HBs in EAN, the valence bond model of Pauling is used in the following <cit.>. Here, starting from the HB distances, r_NH and r_OH, bond orders are defined as p_i = exp(-(r_i - r_i^eq )/b_i) with i = { NH, OH} and r_i^eq being the gas phase equilibrium bond distances. Under the constraint that the sum of the two bond orders must equal one, the two coordinates depend on each other and the HB geometry change along a proton transfer reaction can be described by a single coordinate <cit.>. This path is commonly plotted into a graph showing the correlation between the coordinates r_1=0.5(r_ NH-r_ OH) and r_2=r_ NH+r_ OH.  For the case of EAN the gas phase equilibrium distances r_i^eq are obtained from single molecule gas phase DFTB optimizations. The NH bond length is set to the average of the NH bond lengths of the optimized structure. 



The bond order decay parameters were changed simultaneously starting from the values reported for HBs in crystal structures <cit.> until a good visual fit with the DFTB data sampled along the trajectory was observed. This results in b_OH=0.3 Å and b_NH=0.33 Å.  All parameters of the Pauling model are collected in table <ref>.  In passing we note that the same values were found for EAN gas phase cluster simulations <cit.>. 

The resulting  reaction path is drawn in figure <ref> (blue line) together with the distribution of geometries according to the trajectory. If the distances r_ NH and r_  OH were fully uncorrelated, the geometries would lie on a straight line with  slope 2, plotted for reference as a gray line. From the apparent deviation from linearity one can conclude on a weak correlation between the bond distances across the HB.

This analysis was previously performed for triethylammonium nitrate (TEAN) using a similar setup <cit.>. The resulting parameters of the valence bond order model are also given in table <ref>. In particular the  decay parameters were found to be b_NH =0.35 Å and b_OH=0.32 Å. Thus, judging from the geometric correlations the HBs in EAN are slightly weaker than in TEAN. This can be attributed to the HB angle. Compared with TEAN, in EAN not only the average angle is larger (26^∘ versus 34^∘), but also the distribution is broader (cf. figure <ref>b).





 §.§ IR Spectra



In figure <ref> the FIR spectrum calculated via (<ref>) is presented together with an experimental spectrum taken from  <cit.>. The calculated spectrum (blue) shows three peaks in the FIR region, i.e. at 92 cm^-1, 201 cm^-1 (with a shoulder at 260 cm^-1), and 439 cm^-1. To understand the origin of these peaks, the spectrum has been additionally computed from the dipoles of the anions, μ⃗_anion (red), and cations, μ⃗_cation (green), only. In general the intensities of the cation and anion spectra  do not add up to the full spectrum, because cross-terms μ⃗_cation·μ⃗_anion originating from the square of the total dipole,  μ⃗_cation + μ⃗_anion, in (<ref>) are not included in the separate spectra. 
The peak at the lowest frequency  stems almost exclusively from cationic motion and is attributed to be essentially due to dispersion interactions, which usually are dominated by alkyl tails, with smaller contributions due to interionic HB motion.  Note that this peak was not present in previous hexamer gas phase cluster simulations, since in this case the alkyl chains are pointing outwards what reduces their interaction <cit.>. 

The peak at 201 cm^-1 has a low intensity in cationic as well as anionic spectra. Hence, it can be attributed to cross-terms, which suggests that the band is due to a correlated cation-anion motion likely originating from interionic HBs.
The peak around 440 cm^-1 was previously assigned by normal mode calculations to a N-C-C bending motion of the cation <cit.>. 

These assignments are in line with previous FIR experimental investigations of EAN <cit.>.
Overall the simulation results reproduce the experimental spectra very well. The largest difference is the absence of a clear  low-frequency dispersion band around 92 cm^-1 in the experiment. This could be due to a larger broadening of this band in the experimental data. In fact related bands around 100  are discernible in less flexible TEA-like systems <cit.>.






The MIR spectral region is shown with spectra from experiment (black) and simulation (blue) in figure <ref>. The experimental signal covers a wide range from 2700 to 3300 cm^-1 and shows a left-skewed distribution. The maximum intensity is located at 3092 cm^-1 with a side peak at 3159 cm^-1 and a shoulder at 3250 cm^-1. Further, there are several local maxima on the low-frequency side. In general, the MIR spectrum results from the extended three-dimensional HB network, including possible  proton transfers between anions and cations.
 Due to H-bonding the symmetric and anti-symmetric NH stretching vibrations are strongly red-shifted and overlap with the corresponding C-H vibrational bands of the ethyl group. To assign the  observed structure, the separate contributions to the calculated spectrum of the NH_3 (green) and CH_3 (red) dipoles are also shown in figure <ref>. For better visual comparison the maximum of the CH_3-only simulation is used as reference to scale all theoretical spectra by aligning it to the position found for the respective modes in experimental Raman spectra <cit.>. This results in a scaling of the frequencies by a factor of 0.978.

Overall, the agreement between simulation and experiment is very good, i.e. position and general shape of the spectra match very well. As expected the contribution of  the polar NH_3 groups has a much higher intensity than those of the non-polar CH_3  groups.  The IR signal above 3050 cm^-1 stems exclusively from NH_3 motion and the maximum coincides with that of the full spectrum. Note the perfect match of the global maximum at 3092 (unscaled 3162) cm^-1, mainly originating from the NH_3 dipole, which is achieved even though only the CH_3 peak was used to obtain the scaling factor. Additionally, the contribution from the NH_3 dipoles shows the same left-skewed behavior as the full spectrum  and it is very broad, i.e. spanning the range  from 2800 to 3250 cm^-1. The latter indicates a substantial distribution of HB strengths. The alkyl CH stretching vibrations have a much weaker intensity and span a narrower range in between two maxima at 2923 (3060) cm^-1 and 2863 (2927) cm^-1.





We note that the spectrum according to hexamer gas phase simulations using the DFTB method <cit.> deviates from that of the present bulk case. Especially, the cluster shows a larger HB strength, as indicated by a larger blue shift in FIR and a larger  red shift in MIR of related spectral features. The results of the present simulations can be further compared with bulk phase calculation for TEAN reported in <cit.>. One finds that for EAN the red shift of the NH-stretching vibrations in the MIR region is smaller, but the blue shift in the FIR region is more pronounced. At first sight, this is a surprising result since the magnitudes of both shifts are associated with the HB strength. According to the above analysis of HB geometries, EAN should feature a weaker HB than TEAN. Thus, the comparison of shifts in the FIR region is misleading insofar as the interionic motion in EAN is influenced by three possible HBs per ion, whereas there is only one for TEAN. On average this may cause a tighter potential and hence a larger blue-shift.













 §.§ Rotational Correlation


In figure <ref> we show the  calculated deuteron quadrupole coupling constants versus calculated proton chemical shifts for HBs taken from different sized EAN clusters. Apparently, there is  an almost linear dependence. From linear regression we obtain the relation χ_ D=274.86 kHz - 12.342 δ^1H kHz/ppm. 
We observe that the proton chemical shifts vary between 2 and 11 ppm depending on the strength of the NH⋯O cation-anion interaction at the different  positions within the clusters. This range for the chemical shifts corresponds to deuteron quadrupole coupling constants varying from 250 down to 140 kHz.   The weaker HBs are characterized by smaller δ^1H  shifts, and larger χ_ D values, and vice versa. It does not matter, whether these pairs of properties are calculated for N-H/N-D in varying sized clusters or for different configurations within these clusters, they all show a linear dependence. This is not at all surprising, because both properties, δ^1H  as well as χ_ D, are sensitive to local and directional interactions, such as H-bonding. 

The ultimate test for the reliability of the linear relationship between both properties in PILs is carried out by interpolating the chemical shifts to 0 ppm, indicating the absence of intermolecular interactions.  For this case we expect the calculated χ_ D values to be similar to those measured for the gas phase of similar molecules. And indeed, the estimated value of 275 kHz for χ_ D in EAN is only slightly lower than the measured gas phase and calculated monomer χ_ D values for ammonia (290.6 kHz), formamide (292 kHz) and N-methylformamide (294 kHz), respectively
 <cit.>.

We now use this relation for deriving the deuteron quadrupole coupling constant of EAN in the liquid phase. From the measured proton chemical shift δ^1H=7.375 ppm at 303 K we obtain a deuteron quadrupole constant of χ_ D=184 kHz (see table <ref>). In figure <ref> we compare the χ_ D value for EAN with those of triethylammonium-based PILs, including bis(trifluoromethylsulfonyl)-imide [NTf_2], trifluoromethyl sulfonate [CF_3SO_3], and methylsulfonate [CH_3SO_3] as anions (PILs I-III). At a first glance, it seems to be surprising that EAN, which includes the strongest interacting nitrate anion, does not give the lowest χ_ D value. 
However, in the TEA-based PILs we have single, one-dimensional HBs between cation and anion, whereas in EAN a three-dimensional HB network is formed. Thus the charge transfer from the anions to N-H bonds of the cations is distributed on three bonds rather than one as it is the case for the other PILs.    Clearly, the quadrupole coupling constant is a sensitive probe for H-bonding in ILs, as it has been previously shown for molecular liquids <cit.>.







The obtained χ_ D value of 184 kHz can now be used for the determination of the reorientational correlation time  using (<ref>), which gives τ_ NH= 16.1 ps at 303 K. Compared to the other PILs, which are described by single N-H bonds and weaker interacting anions, the correlation time of EAN are surprisingly low. The τ_ NH values for PILs I-III are 154.2 ps, 190.2 ps and 316.9 ps, respectively and thus one order of magnitude larger than those of EAN (see table <ref>). In PILs I-III the molecular motion of the cations was slowed down in the order [NTf_2] < [CF_3SO_3] < [CH_3SO_3]. For the strongly interacting [NO_3] anion even larger correlation times should be expected. The strongly enhanced dynamics in EAN can ony be explained by the high flexibility of the HB network and fast proton exchange. That becomes obvious if we compare the correlation times for all four PILs with the corresponding viscosities.


We recently discussed the relation between reorientational correlation times   and viscosities, η, for PILs I-III. The observed linear behavior is expected from the Stokes-Einstein-Debye (SED) relation <cit.>, τ_ NH= η V_ eff/k_ BT if the volume/size of the solute is similar for all PILs. Here, V_ eff is the effective volume obtained by multiplying the volume V with the so-called Gierer-Wirtz factor <cit.>, which for neat liquids has a value of about 0.16. 
 From the slope of the plot in figure <ref>, we can estimate the effective volume, V_ eff, to be about 0.1217 nm^3. If we assume a spherical shape of the solute particles, we obtain an effective radius, R_ eff, of about 3.07 Å. This radius is in reasonable agreement with the one of the TEA cation, which has been calculated to be 3.22 Å <cit.>. Since we use reorientational correlation times of the cation only, which are identical in all three PILs, a linear behavior between τ_ NH  and  η follows. That the straight line goes through zero further supports the idea that the microscopic correlation times and macroscopic viscosities describe the dynamical behavior of similar molecular species.  If we add τ_ NH and η for EAN to figure <ref>, we observe that this relation does not hold for the smaller ethylammonium cation.      Even if we take into account, that the volume of ethylammonium is half as large as that of the TEA cation, the measured correlation time is too small to be described by the SED relation.  Obviously, the correlation time reflects the shorter life times of the HBs in the three-dimensional network of EAN.  In contrast, in PILs I-III τ_ NH reflects the rotational dynamics of the strong HB within the ion pair and is thus more suitable to describe the simple hydrodynamic basis of the SED relation.  





 

Further insight into this conclusion comes from the MD simulations.  The rotational correlation function C(t), calculated according to (<ref>), is presented in figure <ref>. C(t) shows an initially fast decay, followed by a much slower decay, which doesn't reach zero in the presented time interval. For the following analysis C(t) has been modeled by a tri-exponential function, i.e. C(t) = ∑_i=1^3 a_i exp(- t/τ_i ). The resulting three time scales (amplitudes) are τ_1=0.167 ps (0.24), τ_2=1.69 ps (0.34) and τ_3=13.4 ps (0.33). The rotational correlation time according to (<ref>) for this tri-exponential fit follows as 4.4 ps, i.e. it is about four times smaller than the experimental one of 16.1 ps. However, the simulations clearly confirm the order of magnitude for EAN as compared to the TEA-based systems. We note in passing, that force field MD simulations  of PILs I-III  by Strauch et al. yielded correlation times, which for II and III substantially exceeded the experimental values <cit.>.  

The three different decay times contributing to C(t) can be discussed in the context of the findings by Hunger et al., who reported on dielectric relaxation and femtosecond IR experiments <cit.>.  Focussing on the transient dynamics of the initially excited ND-stretching vibration, they observed three time scales at 298 K, i.e. 0.2 ps, 2.5 ps, and 15.7 ps. The former two have been related to vibrational energy relaxation via intermediate vibrational states.  The slowest one was attributed to reorientational dynamics of the HBs.  By comparison with the time constant obtained from dielectric relaxation, they proposed a cooperative reorientation through large angle jumps around the CN bond axis <cit.>. Comparison with the present correlation function analysis yields a good agreement of the slowest component of the tri-exponential fit (13.4 ps) with the above time scale (15.7 ps). However, we must note that the large angle jumps, which were considered to be at the origin for this time scale in <cit.>, have not been clearly observed in the present simulations.


Relating the two shorter time scales of <cit.> to the present findings is not that straightforward, although the agreement is striking. Of course, one can argue that any change in HB geometry leading to the loss of orientational correlation is connected to anharmonic couplings, which could be responsible for the observed vibrational relaxation. This is supported by the fact that the τ_1=0.167 ps decay time is equal to the continuous lifetime τ_cont of a HB in EAN. Hence, this time scale reflects the simultaneous effect of HB widening/shortening and NH-bending, which is driven by anharmonic couplings. In passing we note that the disruption of the HBs can also be viewed as being due to librational motion. 




§ SUMMARY



We provide valuable information about the structure and dynamics of the HB network in EAN, which can be regarded as archetype of protic ionic liquids.  The properties of the HBs have been analyzed using simulations and experiments within a multi-dimensional approach, which focussed on geometric correlations, IR, and NMR spectroscopy. DFTB-based molecular dynamics provided the frame for unraveling molecular details behind the macroscopic observables. EAN features a flexible three-dimensional network of HBs with moderate strength. This makes EAN  distinct from related TEA-based PILs, which have a one-dimensional network only. Geometric correlations between HB  and NH bond lengths are not very prononced, which can be seen as a consequence of the rapid interchange of HB partners, expressed by the rather short continuous lifetime. Signatures of H-bonding are also observed in the MIR and FIR spectra as red and blue-shifted bands, respectively. DFTB has the advantage of providing a means to calculate IR spectra due to molecular fragments in a very efficient way, which allowed to unravel to origin of broad spectral bands. This way we could show that the width of NH vibrational contributions spans the range from 2800 to 3250 . Rotational correlation of the NH bond vectors within the HB provides an additional means for scrutinizing the HB network. Whereas the deuteron quadrupole coupling constants for EAN and TEA-based PILs suggest comparable cation-anion interaction strengths, the dynamics is strikingly different as is manifest in the order of magnitude difference of the respective rotational correlation times. 



§ ACKNOWLEDGMENTS

The authors thank the Deutsche Forschungsgemeinschaft (DFG) for financial support through the Sfb 652. 



§ REFERENCES



10
url<#>1#1urlprefixURL 



walden14_405
Walden P 1914 Bull. Acad. Imper. Sci. (St. Petersburg) 8 405–422

plechkova08_123
Plechkova N and Seddon K 2008 Chem. Soc. Rev. 37 123–150

greaves08_206
Greaves T L and Drummond C J 2008 Chem. Rev. 108 206–237

angell07_1228
Angell C A, Byrne N and Belieres J P 2007 Acc. Chem. Res. 40
  1228–1236

fumino09_3184
Fumino K, Wulf A and Ludwig R 2009 Angew. Chem. Int. Ed. 48
  3184–3186

kruger10_101101
Krüger M, Bründermann E, Funkner S, Weingärtner H and Havenith M
  2010 J. Chem. Phys. 132 101101

hayes12_7468
Hayes R, Imberti S, Warr G G and Atkin R 2012 Angew. Chem. Int. Ed.
  51 7468–7471

menne13_39
Menne S, Pires J, Anouti M and Balducci A 2013 Electrochem. Commun.
  31 39–41

markusson07_8717
Markusson H, Belières J P, Johansson P, Angell C A and Jacobsson P 2007
  J. Phys. Chem. A 111 8717–8723

emelyanenko14_11640
Emel'yanenko V N, Boeck G, Verevkin S P and Ludwig R 2014 Chem. Eur. J.
  20 11640–11645

strauch16_17788
Strauch M, Bonsa A M, Golub B, Overbeck V, Michalik D, Paschek D and Ludwig R
  2016 PhysChemChemPhys 18 17788–17794

ludwig98_9312
Ludwig R, Reis O, Winter R, Weinhold F and Farrar T C 1998 J. Phys. Chem.
  B 102 9312–9318

Maginn09_373101
Maginn E J 2009 J. Phys. Condens. Matt. 21 373101–373118

kirchner09_213
Kirchner B 2009 Top. Curr. Chem. 290 213–262

kirchner15_202
Kirchner B, Hollóczki O, Lopes J N C and Pádua A A H 2015 Wiley
  Interdiscip. Rev. Comput. Mol. Sci. 5 202–214

zahn10_124506
Zahn S, Thar J and Kirchner B 2010 J. Chem. Phys. 132 124506

Elstner98_7260
Elstner M, Porezag D, Jungnickel G, Elsner J, Haugk M, Frauenheim T, Suhai S
  and Seifert G 1998 Phys. Rev. B 58 7260–7268

Addicoat14_4633
Addicoat M A, Stefanovic R, Webber G B, Atkin R and Page A J 2014 J. Chem.
  Theory Comput. 10 4633–4643

zentel16_234504
Zentel T and Kühn O 2016 J. Chem. Phys. 145 234504

zentel17_a
Zentel T and Kühn O 2017 Theor. Chem. Acc. 136 87

zentel17_56
Zentel T and Kühn O 2017 J. Mol. Liq. 226 56–62

canongialopes04_16893
Canongia Lopes J N and Pádua A A 2004 J. Phys. Chem. B 108
  16893–16898

canongialopes04_2038
Canongia Lopes J N, Deschamps J and Pádua A A H 2004 J. Phys. Chem.
  B 108 2038–2047

hess08_435
Hess B, Kutzner C, van der Spoel D and Lindahl E 2008 J. Chem. Theory
  Comput. 4 435–447

aradi07_5678
Aradi B, Hourahine B and Frauenheim T 2007 J. Phys. Chem. A 111
  5678–5684

gaus11_931
Gaus M, Cui Q and Elstner M 2011 J. Chem. Theory Comput. 7
  931–948

yangy07_10861
Yang, Y, Yu H, York D, Cui Q and Elstner M 2007 J. Phys. Chem. A 111 10861–10873

gaus13_338
Gaus M, Goez A and Elstner M 2013 J. Chem. Theory Comput. 9
  338–354

elstner01_5149
Elstner M, Hobza P, Frauenheim T, Suhai S and Kaxiras E 2001 J. Chem.
  Phys. 114 5149–5155

kaiser80_105
Kaiser J F and Schafer R W 1980 IEEE Trans. 
Acoustics, Speech and Signal Processing 28 105–107

Brehm12_5030
Brehm M, Weber H, Pensado A S, Stark A and Kirchner B 2012 PhysChemChemPhys 14 5030

abragam61_
Abragam A 1961 The Principle of Nuclear Magnetism (Oxford: Oxford
  Clarendon Press)

farrar87_
Farrar T C 1987 An Introduction to Pulse NMR Spectroscopy (Chicago:
  Farragut Press)

huber85_4591
Huber H 1985 J. Chem. Phys. 83 4591–4598

wulf07_323
Wulf A, Ludwig R, Sasisanker P and Weingärtner H 2007 Chem. Phys.
  Lett. 439 323–326

wulf07_2265
Wulf A, Fumino K, Michalik D and Ludwig R 2007 ChemPhysChem 8
  2265–2269

ludwig94_6684
Ludwig R, Farrar T C and Zeidler M D 1994 J. Phys. Chem. 98
  6684–6687

ludwig95_6941
Ludwig R, Weinhold F and Farrar T C 1995 J. Chem. Phys. 103
  6941–6950

ludwig95_329
Ludwig R 1995 Chem. Phys. 195 329–337

ludwig91_89
Ludwig R, Gill D S and Zeidler M D 1991 Z. Naturforsch. A 46a
  89–94

ludwig94_313
Ludwig R and Zeidler M D 1994 Mol. Phys. 82 313–323

ludwig95_19
Ludwig R, Zeidler M D and Farrar T C 1995 Zeitschr. Phys. Chem. 189 19–27

gordalla86_817
Gordalla B C and Zeidler M D 1986 Mol. Phys. 59 817–828

gordalla91_975
Gordalla B C and Zeidler M D 1991 Mol. Phys. 74 975–984

eggenberger92_5898
Eggenberger R, Gerber S, Huber H, Searles D and Welker M 1992 J. Chem.
  Phys. 97 5898–5904

searles_
Searles D and Huber H Encyclopedia of Nuclear Magnetic Resonance,
  vol 9 ed Grant D and Harris R (New York: Wiley)

frisch09
Frisch M J, Trucks G W, Schlegel H B, Scuseria G E, Robb M A, Cheeseman J R,
  Scalmani G, Barone V, Mennucci B, Petersson G A, Nakatsuji H, Caricato M, Li
  X, Hratchian H P, Izmaylov A F, Bloino J, Zheng G, Sonnenberg J L, Hada M,
  Ehara M, Toyota K, Fukuda R, Hasegawa J, Ishida M, Nakajima T, Honda Y, Kitao
  O, Nakai H, Vreven T, Montgomery J J A, Peralta J E, Ogliaro F, Bearpark M,
  Heyd J J, Brothers E, Kudin K N, Staroverov V N, Kobayashi R, Normand J,
  Raghavachari K, Rendell A, Burant J C, Iyengar S S, Tomasi J, Cossi M, Rega
  N, Millam J M, Klene M, Knox J E, Cross J B, Bakken V, Adamo C, Jaramillo J,
  Gomperts R, Stratmann R E, Yazyev O, Austin A J, Cammi R, Pomelli C,
  Ochterski J W, Martin R L, Morokuma K, Zakrzewski V G, Voth G A, Salvador P,
  Dannenberg J J, Dapprich S, Daniels A D, Farkas O, Foresman J B, Ortiz J V,
  Cioslowski J and Fox D J 2009 Gaussian 09 (Wallingfort: Gaussian
  Inc.)

grimme10_154104
Grimme S, Antony J, Ehrlich S and Krieg H 2010 J. Chem. Phys. 132
  154104

ehrlich11_3414
Ehrlich S, Moellmann J, Reckien W, Bredow T and Grimme S 2011 ChemPhysChem 12 3414–3420

grimme16_5105
Grimme S, Hansen A, Brandenburg J G and Bannwarth C 2016 Chem. Rev.
  116 5105–5154

ludwig95_5118
Ludwig R, Weinhold F and Farrar T C 1995 J. Chem. Phys. 102
  5118–5125

ludwig95_3636
Ludwig R, Weinhold F and Farrar T C 1995 J. Chem. Phys. 103
  3636–3642

ludwig97_8861
Ludwig R, Weinhold F and Farrar T C 1997 J. Phys. Chem. A 101
  8861–8870

ludwig97_499
Ludwig R, Weinhold F and Farrar T C 1997 J. Chem. Phys. 107
  499–507

ludwig98_197
Ludwig R, Weinhold F and Farrar T C 1998 Ber. Bunsen Ges. Phys. Chem.
  102 197–204

ludwig98_205
Ludwig R, Weinhold F and Farrar T C 1998 Ber. Bunsen Ges. Phys. Chem.
  102 205–212

lawrence03_264
Lawrence C P and Skinner J L 2003 J. Chem. Phys. 118 264

evans81_481
Evans D F, Chen S H, Schriver G W and Arnett E M 1981 J. Am. Chem. Soc.
  103 481–482

Hunt15_1257
Hunt P A, Ashworth C R and Matthews R P 2015 Chem. Soc. Rev. 44
  1257–1288

Hayes13_4623
Hayes R, Imberti S, Warr G G and Atkin R 2013 Angew. Chem. Int. Ed.
  Engl. 52 4623–4627

pauling47_542
Pauling L 1947 J. Am. Chem. Soc. 69 542

Limbach04_5195
Limbach H H, Pietrzak M, Sharif S, Tolstoy P M, Shenderovich I G, Smirnov S N,
  Golubev N S and Denisov G S 2004 Chem. Eur. J. 10 5195–5204

yan10_15695
Yan Y and Kühn O 2010 PhysChemChemPhys 12 15695–15703
yan11_1
Yan Y and Kühn O 2011 in Hydrogen Bonding and Transfer in the Excited State,
  vol 1 ed Han K-L and Zhao G J R (New York: John Wiley & Sons), p 1
steiner98_7041
Steiner T 1998 J. Phys. Chem. A 102 7041–7052

fumino12_6236
Fumino K, Reichert E, Wittler K, Hempelmann R and Ludwig R 2012 Angew.
  Chem. Int. Ed. 51 6236–6240

Bodo13_144309
Bodo E, Sferrazza A, Caminiti R, Mangialardo S and Postorino P 2013 J.
  Chem. Phys. 139 144309

seipelt97_1501
Seipelt C and Zeidler M 1997 Ber. Bunsen Ges. Phys. Chem. 101
  1501–1508

wendt98_1077
Wendt M and Farrar T C 1998 Mol. Phys. 95 1077–1081

wendt99_753
Wendt M, Zeidler M D and Farrar T C 1999 Mol. Phys. 97 753–756

wulf06_266
Wulf A and Ludwig R 2006 ChemPhysChem 7 266–272

einstein56_
Einstein A 1956 Investigations on the Theory of Brownian Motion (New
  York: Dover)

debye29_
Debye P 1929 Polar Molecules (New York: Dover)

gierer53_532
Gierer A and Wirtz K 1953 Z. Naturforsch. A 8 532–538

zarrougui15_686
Zarrougui R, Dhahbi M and Lemordant D 2015 J. Solution Chem. 44
  686–702

hunger12_3034
Hunger J, Sonnleitner T, Liu L, Buchner R, Bonn M and Bakker H J 2012 J.
  Phys. Chem. Lett. 3 3034–3038




