J. Chem. Sci. (2018) 130:3
https://doi.org/10.1007/s12039-017-1404-1© Indian Academy of Sciences

REGULAR ARTICLE

Nonpolar solvation dynamics for a nonpolar solute in room
temperature ionic liquid: a nonequilibrium molecular dynamics
simulation study
SANDIPA INDRAa and SNEHASIS DASCHAKRABORTYb,∗
a 19, Nischintapur Lane, Santipur, Nadia, West Bengal 741 404, Indiab Department of Chemistry, Indian Institute of Technology Patna, Patna, Bihar 801 103, India

E-mail: snehasis@iitp.ac.in
MS received 11 September 2017; revised 4 November 2017; accepted 5 November 2017

Abstract. Nonpolar solvation dynamics of a nonpolar diatomic solute in a room temperature ionic liquid (RTIL)
has been followed via nonequilibrium molecular dynamics (MD) simulation. Frank-Condon type excitation of
the solute, previously in equilibrium in RTIL solvent, has been modelled by abruptly changing the Lennard-Jones(LJ) diameter of the solute atoms and thereby disrupting the equilibrium situation. The rearrangement of the
RTIL solvent molecules, which has been seen to be mostly contributed by the solute’s first solvation shell, aroundthe excited solute results overall spectral narrowing and biphasic decay of the solvation energy; a dominant andvery rapid process having sub-100 fs relaxation time, followed by a slower one relaxing at a timescale of ∼5 ps.
A mode-coupling theory based calculation is also used to obtain the nonpolar solvation relaxation function fora model nonpolar solute dissolved in model RTIL solvent. The theoretical relaxation decay is not in very good
agreement with the simulated nonequilibrium solvation response function; the theory predicts the short time
relaxation component slower and the longtime component faster than those of the simulated nonequilibrium
relaxation. We have also checked the validity of the linear response theory (LRT) for nonpolar solvation in
RTIL by looking at the equilibrium solvation energy correlation in the RTIL solvent in presence of the groundstate (GS) and the excited state (ES) solute. Apparent breakdown of the LRT in the present case elucidates
the probable disagreement between the theoretical and simulated nonequilibrium nonpolar solvation response
functions.
Keywords. Solvation dynamics; nonpolar solvation; ionic liquid; molecular dynamics; linear response theory.

1. Introduction
Room temperature ionic liquids (RTILs) are interesting
new materials with multiple fascinating properties. An
increasing number of studies are putting forward newer
applications of these RTILs; as the reaction medium, 1,2
in energy applications, 3 as CO2 capturing material, 4,5etc. These applications are directly connected to thekey macroscopic properties of RTILs, which too originate from the structural and dynamical properties at themolecular level. This is the rationale of growing interestof understanding molecular structure and dynamics of
RTILs using experiment, theory and computer simulation.

The solvation dynamics of RTILs for an electronically excited solute molecule—which is generally
fluorescent—have been studied over the years using
the Stokes shift dynamics method 6–22 in order to detect
the solvation behavior, the knowledge of which is the
key to judge a particular RTIL as a potential reaction
medium. Due to the fact that the dipole moment of
the electronically excited state of the solute generally
increases substantially (e.g., the dipole moment of C153
probe increases by ∼7D upon electronic excitation by
a particular wavelength radiation), 23 the solute primarily probes the solute-RTIL dipole-dipole and ion-dipole
interactions guided solvation of RTILs. The dipolarsolvation response function calculated by following
the time-dependent shift of the fluorescence emission* For correspondenceElectronic supplementary material: The online version of this article (https:// doi.org/ 10.1007/ s12039-017-1404-1) contains
supplementary material, which is available to authorized users.

3

Page 2 of 14

spectrum of the excited probe molecule for imidazolium and pyrrolidinium ILs shows biphasic decay:
the faster component (∼20% of the total response)
relaxes exponentially with sub-picosecond time constant (∼100–700 fs), and the slower component relaxes
stretched-exponentially (stretched exponent  1) withtime constant ranging over a few picoseconds to several nanoseconds. 10–23 While it is somewhat settled thatthe origin of the slower relaxation is usually due to theslower rotation of the dipolar cations and/or the translation of the ions, the origin of the ultrafast time constantis still debated. 24–42 Multiple modes of molecular movements can be responsible for this ultrafast relaxation,
including the hydrogen-bond (H-bond) libration, 25–31
faster rotation of the dipolar cation, 24–29 fast translationof small constitute anion, 10,37–42 etc. Recent classicaldensity functional based theoretical calculations 24–31
supported the experimental findings and thereby put forward possible origins of different solvation relaxation
timescales.
More recent three-pulse photon echo peak shift(3PEPS) measurements using oxazine-4 probe in several
RTILs have reported unexpectedly fast solvation. 43,44
The measured relaxation time scales are way fasterthan those obtained in the time-resolved fluorescence
Stokes shift measurements (TRFSS). 10–13 It’s highly
unlikely that the traditional dynamical modes, discussed
in the preceding paragraph, can completely explain
the origin of this ultrafast solvation. This anomaly
reminds the experimental report of ultrafast solvation
using transient hole-burning technique for a wide varietyof viscous solvents. 45–47 It was observed that the solvation of dimethyl-s-tetrazine—a nonpolar solute—isgoverned by an unprecedented fast relaxation mechanism. The relaxation was sometimes more than twoorders of magnitude faster than the traditional diffusive re-orientational relaxation of polar solvationnear dipolar solutes. The key origin for this very fast
relaxation was later understood to be nonpolar interaction between the excited solute and the solvent andtherefore the associated solvation was termed as nonpolar solvation. 45–53 This opens up the possibility ofidentifying nonpolar solvation as the key relaxation
mechanism found in the recent 3PEPS measurementfor RTILs. Furthermore, the reported Stokes shift valueis ∼200−500 cm−1 , 43,44 which is much smaller thanthat found in TRFSS measurement using C153. 10–13
This is probably due to the low polarity of both theground and excited state oxazine-4 probe molecule.
This motivated one of us to investigate the problem
using mode coupling based theory on non-polar solvation dynamics of different RTILs at a model nonpolar
solute. 30 The theoretical relaxation timescales were in

J. Chem. Sci. (2018) 130:3

semi-quantitative agreement with the measured ones
and thus the theoretical study drew an important conclusion that the 3PEPS and TRFSS experiments areprobably reconnoitring different components of the totalsolvation energy relaxation of a laser-excited dye in agiven RTIL.
Although the above mode-coupling theory based calculation was somewhat successful in understanding the
nonpolar solvation relaxation of the oxazine-4 probe in
RTILs, the theory is not absolutely free from problems,particularly when one uses it for a complex solvent like
RTILs. The key approximations in this theory are as follows: (i) the solvent molecules are perfectly sphericalin shape, (ii) the molecules are uniformly distributedin space, and (iii) the LRT holds perfectly so that the
non-equilibrium solvation response function measured
experimentally can be written as equilibrium solvationenergy correlation function. The first approximation isvalid when the theory is used for smaller molecularliquids and therefore the shape can be well approximated as spherical. But the known complexities of RTIL
molecules somewhat question the spherical approximation. This approximation for RTIL molecules is more ofa concern, particularly for nonpolar solvation. This isdue to the fact that the nonpolar solvation is dominatedby the nonpolar interaction—e.g., LJ interaction, whichis short-ranged and thereby primarily contributed by thenearest neighbour solvent IL molecules. Interestingly,
the similar spherical approximation of RTIL molecule
is proven to be valid for dipolar solvation as the theoretically predicted solvation response semi-quantitatively
agrees with measured one. 24,25,30 This can be understoodfrom the fact that dipolar solvation, unlike the nonpolar
one, takes into account largely the longer ranged soluteIL dipole-dipole interaction and therefore average out
the local shape factor significantly. Similarly, the second
approximation is also questionable due to the fact that
the RTIL is not a microscopically homogeneous at all
and possesses spatial and temporal heterogeneities. 32,54
The importance of local interaction in nonpolar solvation makes the approximation even more questionable.
The approximation, that the LRT holds correctly, has
generally proven to be valid for solvation of a dipolar
solute in common dipolar solvents and RTILs. 39,55 But
the theory severely breaks down for nonpolar solvation
response function in common dipolar solvents when thesolute is subjected to sudden perturbation in size parameters, as done in computer simulation for obtaining
nonequilibrium solvation response. 56 These limitationsof the theory motivate us to study the nonpolar solvation
dynamics for explicit all-atom model of RTIL molecule
near an excited nonpolar solute using nonequilibrium
molecular dynamics simulation technique and address

J. Chem. Sci. (2018) 130:3

the possible contrast between the mode coupling theory
and all-atom molecular dynamics simulation.
In this paper, we determine by nonequilibrium simulation the nonpolar solvation energy relaxation of a
diatomic solute in [Bmim][PF6 ]—an imidazolium IL—
by suddenly changing the LJ parameters of the solute,
the way to mimic the Frank-Condon type excitation.
Then the rearrangement of IL molecules around theexcited solute and thereby the solvation energy relaxation has been monitored over time by looking at the
time-dependent Stokes shift value. We have also compared the mode coupling based theoretical solvationenergy relaxation with the simulated one in order tocheck the rationality of the theory in this RTIL. The
virtue of using all-atom model of IL in the present simulation has given us the opportunity of looking at theindividual contributions of different segments of the IL
to the total Stokes shift and thereby understanding theirrelative significance in the whole process. We have also
calculated the equilibrium solvation energy correlation
function from equilibrium trajectory for both GS and
ES of the solute in order to compare with the nonequilibrium solvation response function and thus check thevalidity of LRT for nonpolar solvation in RTIL.

2. Model description and simulation details
In this section, we first provide the simulation details and theequilibration of the system. The section ends with a discussion
on how the trajectories are generated for analyses.

2.1 Simulation details and system equilibration
The simulated system consists of 128 pairs of [Bmim]+ cation(see Scheme 1a) and [PF6 ]− anion (Scheme 1b) as the solventand one nonpolar solute, modelled as a simple rigid diatomic
molecule (Scheme 1c). Traditionally, this type of solute model
has been used for decades in MD simulation for probing solvation behavior in dipolar liquid. 55–57 In addition, the solutemodel having rotational degrees of freedom is more realistic than a single atomic solute model and takes into accountthe solute’s rotational contribution towards the total solvation
dynamics.
Initial gas-phase structures of the [Bmim]+ cation and[PF6 ]− anion are obtained by optimizing at B3LYP/6−311+
G∗∗ level of theory and basis set with the GAMESS package. 58 128 pairs of the optimized ions and the diatomic solute
molecule, the structure of which is created by the bond distance parameter between the two atoms, are then placed ina cubic box; the configuration is generated using PACKMOL package. 59 The required box length is calculated fromthe previously simulated equilibrium density of this RTIL
(4.705 mol dm−3 ) for the same force field and the temperature (298 K) we use here. 60

Page 3 of 143

For the RTIL solvent, we use the force field developed
by Canongia-Lopes et al., based on the OPLS-AA framework. 60 Detailed description of the potential form including
the meaning of the terms is presented elsewhere. 61 The rationale of using this force field, apart from the fact that the useof this force field results into much less computational costcompared to that of the polarizable model, is its proven predictability of the measured solvation and dielectric relaxation
behavior of a dipolar solute dissolved in RTILs. 33,34,38,39 An
alternate model reducing the charges of the RTIL ions toenhance the mobility also does exist. This model is basedon the understanding that due to the induced polarization
between the ions of RTIL the columbic charges of the ions areless than whole numbers (±1). However, a study on applicability of various models of RTIL in two different dynamical
limits, local and collective, shows that the reduced charge
model is not at all a good choice when the local interaction isimportant, for example in single molecule rotational motionof the dipole. 62 Nonpolar solvation being dominated by theinteraction between the solute and the nearest neighbor RTIL
solvent molecules is predominantly a local phenomenon. Thisjustifies the use of this simpler yet useful model for the
RTIL ions. The solute is modelled as a rigid molecule having
two identical LJ atoms of diameter (σ ) 3.61 Å and LJ welldepth (ε) 1.41 KJ/mol. These values represent the GS of thesolute. Therefore, the solute and the RTIL ions interact only126
L J = 4εvia LJ interaction: V1212 (σ12 /r12 ) − (σ12 /r12 ) ;the cross terms—σ12 and ε12 —are calculated using LorentzBerthelot rule. We use GROMACS package 63 for simulating
the entire system in the box. Periodic boundary condition has
been employed on all sides of the box. We optimize the entire
system using steepest descent algorithm for 1000 steps before
10 ns equilibrium simulation run, of which the first 5 ns simulation is run using NPT ensemble and the next 5 ns simulationis continued with NVT ensemble. This 10 ns simulation is carried out in order to get the proper equilibrated structure. In
both cases, we use Verlet leap-frog algorithm for integrating Newton’s equations of motion with 1 fs time step. While
Nose-Hoover thermostat 64,65 with coupling constant 0.5 ps isused for controlling the temperature at 298 K and Berendsenbarostat 66 , with the same coupling constant of 0.5 ps is used
for controlling the pressure at 1 bar. The long-ranged electrostatic interactions are handled using Particle Mesh Ewald
summation technique. 67
The equilibrated system is then continued to run for another1 ns under the same conditions, discussed above before
we finally start running multiple nonequilibrium trajectories,
which will be discussed now.

2.2 Simulation trajectory details
During the 1 ns equilibrium trajectory with the GS solute,
we save configurations every after 10 fs, the time which is
required to clearly resolve the presumably very fast nonpolar solvation dynamics (see Sect 3.2). Using this trajectory
we calculate inter-atomic radial distribution functions (g(r ))of the RTIL solvent molecules, presented in Section S1 of3Page 4 of 14

J. Chem. Sci. (2018) 130:3

Scheme 1.

The all-atom model of the RTIL, [Bmim]+ [PF6 ]− and the diatomic solute.

Supplementary Information (SI). The simulated g(r )s are inexcellent agreement with those presented in earlier studies,indicating proper equilibration of the system. 68–70 This trajectory is also used to get the equilibrium solvent configuration
around the GS solute, which also assists to make a reference
for the configuration after the nonequilibrium situation is created by suddenly changing the LJ parameter(s) of the solute.
100 nonequilibrium trajectories are initiated at a regularinterval of 30 ps from a long equilibrium trajectory of duration 3.0 ns starting from the end point of the previous 1 ns
equilibrium trajectory. The first nonequilibrium trajectory isthus initiated from the final configuration of the previous
1 ns equilibrium trajectory. All the nonequilibrium trajectories are started by increasing LJ diameters of the soluteatoms by 20% of its GS value, keeping LJ well depth (ε)
unchanged. The new LJ diameter for the solute is thus (σ )4.33 Å. We call the solute with the new LJ parameters as ES
solute. The 100 different nonequilibrium trajectories serve
two purposes. First, they provide a sampling of different
initial solvent environments, a characteristic decisive for aset of such nonequilibrium trajectories. Also, they gradually increase the length of the equilibration process, which
is desired for a slower solvent like [Bmim][PF6 ] RTIL. The
total duration of the equilibrium trajectory in presence of the
GS solute is thus equal to 4 ns, which is used for the equilibrium solvation analysis. We extend the equilibrium trajectory
to 20 ns for obtaining reasonably smooth equilibrium structure of RTIL molecules around the GS solute. We calculate
the time-dependent quantities in Section 3 by averaging over
all these 100 nonequilibrium trajectories at a common reference time at t = 0, the starting point of the trajectories. The
lengths of all the trajectories are 30 ps, during which temperature constraining is lifted so that the natural dynamics are

not affected anyway. We will see in Section 3 that the 30 ps
duration of the nonequilibrium trajectories is sufficient forcompleting the nonpolar solvation around the ES solute. The
configurations have been saved in the interval of 5 fs, which isessential for getting fine details of the time-dependent results,
especially for ultrafast nonpolar solvation. By convergencestudy (Section S2 of SI) we have shown that the set of 100
nonequilibrium trajectories is enough for present sampling
purpose.
One of the nonequilibrium trajectories is then extended foranother 20 ns in order to obtain the equilibrium trajectory
including new RTIL solvent configurations around the solutein its ES. Configurations from both the equilibrium trajectories in presence of GS and ES solute have been saved every
after 10 fs. The equilibrium solvation picture around the GS
and the ES solute are extracted from the equilibrium trajectories of 20 ns each.

3. Results and Discussion
In this section, we first compare the equilibrium solvation structure of the solute between its GS andthe ES. The equilibrium average peak shift is shown
next. The time-dependent analysis is started by looking at the time-dependent spectral width function andpeak shift. We then focus on the relative contributions
from the various RTIL sites and the solute’s solvation shells towards the total time-dependent peak shift.
Simulated nonequilibrium nonpolar solvation response
function is then presented and compared with the solvation energy-energy correlation function, obtained from

J. Chem. Sci. (2018) 130:3

Page 5 of 14

mode-coupling theory based calculation. The section
ends with a discussion on examining the validity of LRT
for nonpolar solvation of the solute in the RTIL solvent.
3.1 Structural change
Figure 1 presents g(r )s between the solute and the centerof mass of the RTIL cation, [Bmim]+ and the anion,[PF6 ]− both at the GS and the ES of the solute. Clearly,
the g(r ) peaks shift towards longer r values for ES,
which is simply due to the increase of the σ value of thesolute atoms. Also, the RTIL ions are more structurallyorganized around the ES solute than that around the GS
one. This stems from the fact that the increase of σ valueof the solute atoms increases the attractive part of thepotential and thereby increases the structural order of the
RTIL solvent molecules around the solute. The changeof the solvation structure around the solute molecule
is associated with the solvation dynamics, the primaryobjective of the study. A closer comparison between
Figure 1a and 1b reveals that the change of g(r ) peakheight is more pronounced for solute-[PF6 ]− than thatfor solute-[Bmim]+ .
1.5

(a)

GS
ES

S-[Bmim]+

g(r)

1.0

0.5

0.0
2
1.5

4

6

(b)

8

10

12

14

S-[PF6]-

GS
ES

g(r)

1.0

0.5

0.0
2

4

6

8

10

12

14

r (Å)
Figure 1. Comparison between the solute-solvent g(r ) at the GS and the ES of the solute:g(r ) between (a) COM of the solute and COM of[Bmim]+ cation and (b) COM of the solute andCOM of [PF6 ]− anion of the RTIL (See Scheme 1
for the atom labelling).

3

In order to get further insight about the equilibrium
solvation structure of the solute in [Bmim][PF6 ], wehave shown in Figure 2 the g(r ) between the soluteatoms and some specific atoms of the RTIL. g(r )sbetween the solute and the remaining alkyl carbon atomsof the alkyl chain of [Bmim]+ cation are presented in
Figure S3 of SI. Visibly the solvation structure aroundthe solute changes appreciably upon the excitation of
the solute. Generally, the first peaks shift towards longer
interatomic distance and the overall solvation structures
look more ordered ones. This reiterates the direct proportionality between the increase of σ value of the soluteatoms with the increase of both repulsive and attractive
part of the potential. Comparative analysis of the solvent’s site-specific g(r ) shows that at GS the solute stays
more in the hydrophobic domain (near CT ), as evidenced
by the most intense g(r ) peak for the solute-RTIL. This
is in accordance with the fact that the solute is nonpolar
and therefore prefers the nonpolar domain of the RTIL,
which is the hydrophobic side chain. Comparison of theintensities of the g(r ) peaks between the solute and thealkyl side chains of the [Bmim]+ cation (Figure 2c andFigure S3 of SI) shows that the GS solute stays closer
to CT than C2 or CS; the order of the peak intensities
is CT > CS > C2. Interestingly, similar behavior isobserved for solvation structure of longer chain RTILs([C10 mim][Cl] and [C10 mim][PF6 ]) around 129 Xe soluteprobe. 71
The solvent environment around the excited solutechanges appreciably. The sharp decrease of the firstg(r ) peak height between the ES solute and RTIL
CT indicates that the ES solute, unlike the GS solute,
is not overwhelmingly situated near the hydrophobic
domain. Also, the first peak heights between the ES
solute and the cationic head group atoms (N A , C R )and [PF6 ]− anion are more intense than those with GS
solute. Therefore, the solvent environment is somewhat
diverse around the ES solute as indicated by the similar first peak positions and heights of different g(r )sin Figure 2. The diverse solvent environment aroundthe ES solute corroborates with the features, one generally expects from a solvent like RTIL, which possesses
great extent of spatiotemporal microscopic heterogeneity. 72–77
We now focus on the average dynamic red shift ΔE(t)(see Section S4 of SI for the definition), the value of
E(t) =using
 the equation,
which is calculated
LJJ
VE S,I L (t) − VG S,I L (t) , where VG S,I L and VELS,I
L are
the LJ interactions between the GS and ES solute respectively with the RTIL solvent. The normalized probability
distributions of ΔE values in the presence of GS and
ES solute are presented in Figure 3. These two distributions are calculated from the two equilibrium trajectories

3

Page 6 of 14

J. Chem. Sci. (2018) 130:3
4

2.0

(a)

GS
ES

g(r)

1.5

(c)

S-NA
3

1.0

2

0.5

1

S-CT

GS
ES

0

0.0
2

4

6

8

r (Å)

10

12

2

14

4

6

8

10

r (Å)

14

2.5

2.0

(b)

(d)

S-CR

GS
ES

1.5

S-P

2.0

GS
ES

1.5

g(r)

12

1.0

1.0

0.5

0.5
0.0

0.0
2

4

6

8

r (Å)

10

12

2

14

4

6

8

r (Å)

10

12

14

(4 ns), one involving the GS solute and the other involving the ES solute. Clearly, neither of the distributions
are Gaussian in nature indicating the presence of local
structural heterogeneity around the solute. It’s worth
mentioning that the same distributions for the dipolarsolvation of a diatomic solute in a model polar solvent
were seen to be Gaussian. 55 The GS distribution is set tobe centered at ΔE = 0 and therefore the peak positionof the ES distribution measures the equilibrium averaged shift. The ES distribution peak is shifted towards
the red end by ∼3000 cm−1 , indicating the average
red shift. Note we have obtained much larger Stokes
shift compared to the reported shift in 3PEPS measurement. 43,44 This is clearly due to the high increase of theσ values by 20% of its initial value. Despite the factthat smaller increase of σ value might result into thesmaller shift, we increased the value of σ of the solute
so high. The rationale is following. Figure 3 clearly
demonstrates that the distribution of shift in the presenceof GS solute is far wide (the width at half maximumis ∼4000 cm−1 ). This distribution in presence of the
ES solute should be amply different from that of the
GS solute in order to get a clear time-dependent shift(see Section 3.2) with minimum error especially when

Normalized Probablity

Figure 2. g(r ) between solute atoms and the various atoms of the RTIL ions. The g(r ) for GS solute islabelled as blue and that for the ES solute is labelled as red. See Scheme 1 for the atom-labelling.

0.04

ES

Red Shift

GS

0.02

0.00
-6000

-3000

E (cm-1)

0

3000

Figure 3. Normalized equilibrium probability distributionsof ΔE values in the presence of the GS solute (blue line) andthe ES solute (red line). The sampling is done after every 10 fs
over 4 ns equilibrium trajectories. The vertical dashed lines
indicate the maxima of the peaks.

we have limited number of nonequilibrium trajectories.
Choice of smaller increment of σ value of the solutewould result into significant overlapping between thetwo distributions in Figure 3 for GS and ES solute andJ. Chem. Sci. (2018) 130:3

Page 7 of 14

thereby the time-dependent shift would have low signal/noise ratio. However, the solvation energy relaxation
time scales, the primary objective of this paper, should
not have much effect on this extent of change of the σvalue.
Further inspection of Figure 3 also shows that theinitial E probability distribution appropriate for the
GS solute is much wider than the final E probability distribution appropriate for the ES solute. Thissuggests TRFSS spectrum should show an overall timedependent narrowing, which we will examine in the next
subsection. This spectral narrowing must be directly
related to the larger solvent force constant near the minima of the solvent-free energy curve for the ES solute.
The larger solvent force constant for the ES solute isagain due to the increased interaction as indicated by
the more ordered solvation structure of the ES solutethan that of the GS solute.
3.2 Nonpolar solvation dynamics
We now examine the time-dependent width function
followed by simulated nonequilibrium nonpolar solvation dynamics. We then discuss another well-establishedapproach namely, the mode-coupling theory based calculation to obtain the nonpolar solvation response function and finally compare the theoretical relaxation timescales with the simulated ones in order to understandthe origin of the simulated relaxation timescales.
The time-dependent width function, W (t), is calculated by looking at the standard deviation of the ΔE
sample for 100 different
non-equilibrium trajectories;
W (t) = (E)2 − E2 . Figure 4 presents timedependent W (t), which shows very fast decay as it
almost reaches equilibrium in 150 fs. This ultrafast
decay must be seen also in time-dependent Stokes shift,
as we will see now.
3.0

2.0

6

-2

W(t) /(10 cm )

2.5

1.5
1.0
0.5
0.0
0.0

0.2

0.4

0.6

0.8

1.0

time (ps)Figure 4.

The width function W (t) as a function of time.

3

The time-dependent shift E (t) (see Eq. 4 in Section S4 of SI) has been plotted in Figure 5a as a functionof time. The overall dynamic Stokes shift is ∼3580 cm−1
(obtained by time averaging at long time decay), whichis slightly higher than the equilibrium average red shift(3000 cm−1 ). The probable reason for this difference islesser sampling in the time-resolved study than that inaverage equilibrium ones calculated in Figure 3. Notethat while sampling of only 100 different configurations
is considered for nonequilibrium calculation, normalized equilibrium probability distributions of ΔE values
in the presence of the GS and ES solute are calculatedusing ∼5 × 105 configurations (once in every 10 fs during 5 ns equilibrium trajectory). Furthermore, the verybroad distribution of E in presence of the GS solute(width at the half maximum ∼4000 cm−1 ) introduces
uncertainty in E(t = 0) no matter what. All thesefactors cause the observed difference between the equilibrium and time-dependent shift values. Note that thelarge dynamic Stokes shift is obtained as in the case of
equilibrium Stokes shift, detailed in Section 3.1, where
we have also discussed the possible reason for obtaining
this large Stokes shift.
Now, we focus on the relative contributions of thetotal dynamic Stokes shift from different solvation shells
and different segments of the RTIL molecule in order togain further knowledge about the origin of this dynamic
Stokes shift value. In Figure 5b, we have plotted thetwo time-dependent shifts: (i) that due to the interaction between the solute and the solvent RTIL molecules
which are only in the 1st solvation shell of the soluteand (ii) the one due to the interaction between the soluteand the solvent RTIL molecules which are outside the
1st solvation shell. The separation between the solvation shells is entirely based on the distance criteria (thesolute – the center of mass of RTIL molecule) which
was reflected in Figure 1. For example, the first solvation shell of the solute is considered up to the distance
where the first minimum appears after the first peak ofg(r ).
Figure 5b very clearly shows that the total dynamic
Stokes shift ∼ 3580 cm−1 is remarkably contributed by
the first solvation shell of the solvent RTIL; the value is∼ 3534 cm−1 (obtained by time averaging at long time
decay), which is almost ∼ 99% of the total shift. The
rest of the shift (∼ 46 cm−1 ) is contributed by the solvation shells beyond the first one. This is consistent withthe fact that the LJ interaction—the only solute-solvent
interaction in the system—is short ranged and decays
very rapidly as a function of internuclear separation(∝ 1/r −6 ).
We now investigate the relative contributions of four
key segments of the RTIL—the methyl group (Me), thePage 8 of 140

J. Chem. Sci. (2018) 130:3

1.0

(a)

1.0

0.8
SNP(t)< E> cm-1

-1000

-2000

0.6

Long time decay

SNP(t)

3

0.5

0.4
0.0
0

0.2

-3000

2

4

6

8

10

12

14

time (ps)

0.0
-4000
0

1

2

3

4

5

0.0

0

(b)

st

Beyond the 1 Hydration Shell< E> cm-1< E> cm-1

-2000

1.0
time (ps)

1.5

2.0

Figure 6. Nonequilibrium nonpolar solvation response
N eqm
decay, S N P (t), of the nonpolar solute in [Bmim][PF6 ]N eqm
RTIL. The longtime decay of S N P (t) is shown in the inset.
The blue horizontal dashed line indicates the value “0” forN eqm
S N P (t).

0

-1000

0.5

-50
-100
-150
-200
0

1

2

3

4

5

time (ps)

-3000
st

1 Hydration Shell-4000
0< E> cm-1

0

1

2

3

4

5

(c)

-500

Me
Bu

-1000

Interestingly, the [PF6 ]− anion has the highest contribution of all. The order of increasing contribution from thesegments towards the total shift is as follows: [PF6 ]− >Im > Bu ≈ Me
Figure 5a clearly shows that the shift is very much
rapid as it is almost over by ∼ 150 fs. However, avery subtle shift occurs at longer times up to ∼ 15 ps,
which we examine now. The time-dependent nonequilibrium solvation response function SNN Peqm (t) has been
displayed in Figure 6. This has been calculated using
the following equation:
SNN Peqm (t) =ImE (t) − E (∞),E (0) − E (∞)

(1)

-

[PF6]
-1500
0

1

2

3

4

5

time (ps)Figure 5. (a) Trajectory-averaged time-dependent total
Stokes shift E (t) as a function of time. (b) Separate contributions of the shift from 1st solvation shell (blue line) andfrom the outer shells (red line). The inset shows the latter
contribution more clearly. (c) Division of the total shift incontributions from separate segments of [Bmim][PF6 ] RTIL
molecule. The blue, red, black and pink lines represent thecontributions from the methyl group, the butyl group, the imidazolium ring, and the [PF6 ]− anion.

imidazolium ring (Im), the butyl group (Bu), and the[PF6 ]− anion—towards the total shift. This particularsegmentation is done on the basis of an established foursite model of RTIL, proposed by Roy and Maroncelli. 78

where “<>” indicates the quantity averaged over various initial conditions. Note that the E(t) values during
the last 5 ps of the nonequilibrium trajectories areused for obtaining the E (∞) values. The very fast
dynamics involved in this relaxation process is again
quite apparent. In addition, there is a slower component as well, relaxing up to ∼ 15 ps. The slower part of
the relaxation is clearly visible in the inset of Figure 6.
In order to elucidate further, we have fitted the relaxation function with a biexponential equation of the form,
a1 exp (−t/τ1 ) + a2 exp (−t/τ2 ), where a1 + a2 = 1.
While the faster component having the contributionof ∼ 94% of the total shift relaxes at a timescale of∼ 50 fs, the slower component with only ∼ 6% contribution has the relaxation time constant of ∼ 5.07 ps.
Note that 3PEPS measurement using oxazine-4 soluteprobe molecule, which has a very small change of dipole
moment upon photo-excitation, reported the biphasic

J. Chem. Sci. (2018) 130:3

Page 9 of 14

solvation relaxation with time constants of ∼ 150 fs
and ∼ 2.20 ps, which are in the same order of magnitudes with the simulated time scales. 43,44 This indicatesthat the measured solvation can instigate largely from
nonpolar solvation due to a sudden change of nonpolar interaction between the solute and the RTIL solvent,which can be caused by the Frank-Condon type excitation of the oxazine-4 solute molecule. Note that 3PEPS
measurements are associated with the line-broadening
function, determined by the total solvation energy correlation function, 43 both polar and nonpolar components
of solute-solvent interactions contribute to the total peakshift.
Equilibrium nonpolar solvation energy correlation
function can also be obtained using mode-couplingbased semi-molecular level theory. 30,48,49,79 Recently,
this theory has been used for this system in order to
get an insight on the origin of the relaxation timescales,
found in the 3PEPS experiment. 43,44 We now discuss thetheory briefly and compare with the present simulated
nonequilibrium solvation response function.
The nonequilibrium solvation response function
within the realm of LRT can be equated to the normalized nonpolar solvation energy correlation function as,y
SNT heor(t) =P

C N P (t),C N P (t = 0)(2)where C N P (t)—the un-normalized nonpolar solvation energy-energy time correlation function—can be
expressed as the sum of a short time binary part and along time collective part assuming complete decoupling
between them. 30,48,49 Therefore, one can write:
C N P (t) =C NB P (t = 0) exp − (t/τ B )2 + C NρρP (t) . (3)The first and the second term of the right-hand side of
Eq. 3 are due to the interactions at binary (solute-IL) andcollective limit. The faster binary component relaxes at
a timescale of τ EB , which can be written as 30,48,49
−2C NB P (t = 0)
,
(4)
τ B=
C̈ NB P (t = 0)
d 2 C B (t=0)

NP
where C̈ NB P (t = 0) =. The binary part is basidt 2cally determined by the solute-IL solvent LJ interaction
potential (v12 (r )), where the solvent distribution aroundthe solute can be approximated by the solute-IL solvent
pair correlation function, g12 (r ), weighted by the solvent
density (ρ). Therefore, C NB P (t = 0) and C̈ NB P (t = 0)
can be calculated using the following equations respectively: 30,48,49
 ∞
B
(5)
C N P (t = 0) = 4πρ
drr 2 [v12 (r )]2 g12 (r ),

0

3

and
C̈ NB P (t = 0)4πρ ∞drr 2 [v12 (r )]2 ∇r g12 (r )∇r v12 (r ) .
= ∗
m
0

(6)

Where m ∗ denotes the reduced mass of the solutesolvent composite system.
The slower part which has the contribution from thecollective density fluctuation can be expressed as, 30,48,49(k B T )2 ρ ∞ρρ
2(7)
C N P (t) =dkk 2 c12(k)S (k, t) ,6π 20
where c12 (k) is wave-number, k dependent two particles
direct correlation function and S(k, t) is the dynamic
solvent structure factor of RTIL solvent. Note that thedetails of the calculation of g12 (r ), c12 (k), and S(k, t) in
Eq. 5–7 for this system are discussed elsewhere. 30
We now use the above theoretical formulation to caly
culate the NP solvation response function SNtheor(t) for aPnondipolar solute dissolved in [Bmim][PF6 ] IL assuming the size of the solute molecule same as the RTIL
solvent. The parameters necessary for the calculation
like the sizes of the ions, the density of the RTIL, thedielectric constant of the solvent RTIL, etc., are availy
able in the previous paper. 30 SNT heor(t) is plotted as aPfunction of time in Figure 7. However, for comparison,
we also include the simulated nonequilibrium solvation
response function SNNPeqm (t), already shown in Figure 6.
yVisibly, the initial part of SNT heor(t) decays slower thanP
N eqm
that of SN P (t) but at the later time, during the slower
relaxation, both decay reasonably at a similar rate. Fitting with a function, sum of Gaussian and exponentialfunctions having the following form:(8)
S (t) = a1 exp − (t/τ1 )2 + a2 exp (−t/τ2 )we get the timescale—associated with the initial Gaussian relaxation—(τ1 ) ∼ 225 fs with a contribution∼ 89% and the final slower exponential part having only
∼ 11% contribution relaxing at a timescale of ∼ 2.76 ps.
Clearly, the faster component is ∼ 5 times slower andthe slower component is ∼ 2 times faster than thesimulated relaxation. However, the theoretical and simulated amplitudes of the two modes of the relaxation are
in close agreement with each other. The disagreement
between the timescales probably stems from multiple
reasons. Firstly, the calculation—being based on modecoupling theory—has considered the RTIL molecule asa spherical particle, a very simplified model for the RTIL
molecule, which has a complex structure altogether.
The local complexity in the solute-RTIL solvent interaction is therefore completely ignored. This may leadto serious problem in predicting particularly the initial

3

Page 10 of 14

J. Chem. Sci. (2018) 130:3

and equilibrium solvation energy correlation function. 55
These are equilibrium solvation energy correlation forneutral pair and ion pair. Along the same line we can also
connect the nonequilibrium solvation response function
with the equilibrium solvation energy correlation funcStion either for the GS (SNEqm,G(t)) or for the ES soluteP
Eqm,E S
(SN P(t)).
The equilibrium solvation energy correlation functions for the GS and the ES solute are calculated using
the following equations: 55

1.0
a1

1(fs)

a2

2(ps)

< S>(fs)

t

0.94

50

0.06

5.07

351

S theory
t
NP

0.89

225

0.11

2.76

504

0.8
Neqm

SNP(t)

S NP

0.6
0.4

Neqm
S NP
t

S theory
t
NP

0.2
0.0
0

1

2
3
time (ps)4

5

Figure 7. Time-dependent nonequilibrium nonpolar solvaN eqmtion response function, S N P (t), of the nonpolar solute in[Bmim][PF6 ] RTIL (black solid line) and the equilibrium solvation response function from the theoretical calculation (blueN eqm
dashed line) using Eq. 7. The fitting parameters of S N P (t)T heor yand of S N P(t) are listed in the table (inset). Note that while
T heor y
theS N P(t) is fitted with a function, sum of Gaussian andN eqm
exponential functions (see Eq. 8), S N P (t) is fitted with asimple biexponential function.

response, which is due to the binary interaction between
the solute and the RTIL solvent. Secondly, the theory isbased on the approximation that the RTIL is microscopically homogeneous. This approximation may lead toserious error in estimating the dynamic structure factor(S(k, t)) of the RTIL, the key factor in regulating thecollective density fluctuation of the solvent. Therefore,ythe slower component of SNT heor(t) can deviate fromP
the simulated one simply due to an inaccuracy of thecalculated S(k, t). Note that the present simulation—
being based on the well-established all-atom model of[Bmim][PF6 ]—takes into account the detailed structural and dynamical behavior of the RTIL. Finally, by
assuming the validity of LRT for nonpolar solvation in
RTIL the theory considers the solvation energy-energy
correlation to be same as the nonequilibrium solvation
response. Therefore, the invalidity of the LRT may leadto inaccuracy in calculated solvation response function
and the associated relaxation timescales. In the next subsection, we discuss this in more details.

3.3 Validity of LRTIn the study of solvation dynamics of an ion pair in apolar solvent by Carter and Hynes it was shown that two
forms of LRT can be employed to make the connection
between the nonequilibrium solvation response function

S
SNEqm,G(t) =P

δE (0) δE (t)G S
,(δE (0))2 G S(9)δE (0) δE (t) E S,(δE (0))2 E S

(10)

and
S
SNEqm,E(t) =P

where δE (t) = E (t) − E (t), the fluctuation of
the time-dependent shift from time-averaged value, and“<>” denotes the equilibrium average obtained fromequilibrium simulation trajectories in presence of GS
and ES solute.
Therefore the validity of LRT for NP solvation in[Bmim][PF6 ] RTIL implies thatδE (0) δE (t)G S
(δE (0))2 G S
δE (0) δE (t) E SS= SNEqm,E. (11)(t) =P(δE (0))2 E SSSNN Peqm (t) = SNEqm,G(t) =P

Figure 8 displays the three response functions, SNN Peqm (t),S
Eqm,E S
SNEqm,G(t), and SN P(t) as a function of time.
P
S
Although, two equilibrium correlations—SNEqm,G(t)P
Eqm,E Sand SN P(t)—decay almost similarly with similar
timescales
andamplitudes,
SNN Peqm (t)—most importantly—does not match either withS
Eqm,E S
SNEqm,G(t) or SN P(t). Apparently, the equilibriumP
correlations decay slower than the nonequilibrium correlation. In order to view this with more clarity, we fit thetwo equilibrium relaxation decays with biexponentialfunctions of the form, a1 exp (−t/τ1 ) + a2 exp (−t/τ2 ),where a1 + a2 = 1. The fitting parameters, including the
average solvation time (τ S ) are described in Table 1.
Although the relaxation time constants of the equilibrium correlations are not very different from those of
SNN Peqm (t), the amplitude of the slower part of the relaxation of the formers is almost double than that of thelatter. This makes the overall equilibrium relaxations
much slower than that of SNN Peqm (t). Clearly, this indicates a breakdown of LRT for nonpolar solvation of a
nonpolar solute in this RTIL.

J. Chem. Sci. (2018) 130:3

Page 11 of 14

Lastly, the theoretical and simulated timescales are
in semi-quantitative agreement with the experimentally measured timescales. This corroborates with theunderstanding that the origin of the ultrafast timescales,
observed in 3PEPS measurement in RTILs, might be thenonpolar interaction between excited solute and RTIL
molecules.

1.0
0.8
SNP(t)

Neqm
S NP
t

Eqm ,ES
S NP
t

0.6

Eqm ,GSS NP
t
TheoryS NP t

0.4
0.2
0.0
0.0

4. Conclusions
0.2

0.4
0.6
time (ps)

0.8

1.0

Figure 8. Time-dependent NP solvation response functions
from various routes. The representation is as follows: the simN eqm
ulated nonequilibrium solvation response S N P (t) (blacksolid line), the equilibrium solvation energy correlation funcEqm,G S
tions for the GS (S N P(t) (blue solid line) and for the
Eqm,E S
ES solute (S N P(t)) (red solid line) and the theoretical
T heor ysolvation response function, S N P(t)(green dashed line).
Table 1. The fitting parameters—including the
average solvation time (τ S )—for the four solvation
response functions, shown in Figure 8. The solvation
response timescales, measured by 3PEPS experiment,
are also listed.
S N P (t)N eqm

S N P (t)Eqm,G SSN P
(t)
Eqm,E SSN P(t)
T heor y
SN P(t)
S E x pt (t)3τ1 (fs)a2

τ2 (ps) τsolv  (fs)0.94 50
0.90 93
0.90 89
0.89 225
150

0.06
0.10
0.10
0.11

5.07
6.33
6.18
2.76
2.26

a1

351
716
698
504

In order to elucidate further we have also includedin Figure 8 the theoretical solvation response funcy
tion SNT heor(t) and in Table 1 the corresponding fittingPyparameters. Note that SNT heor(t) is fitted with a function
P—sum of Gaussian and exponential functions—having
the form of Eq. 8. Apparent from Figure 8 and Table 1
is that the theory predicts the initial part of the solvation slowest and the later part of the solvation fastest ofall. Interestingly, the amplitudes of the two relaxationycomponents of SNT heor(t) match very well with the equiP
librium relaxations, indicating the better predictabilityof the theory for the simulated equilibrium solvationyrelaxation. The overall relaxation of, SNT heor(t), decaysP
at a rate in between the equilibrium and nonequilibrium
solvent relaxation.

We have presented here results for nonequilibrium nonpolar solvation dynamics of an ES solute dissolvedin [Bmim][PF6 ] RTIL. The Frank-Condon type excitation of the solute due to its sudden change of the
LJ diameter creates a perturbation in the equilibrium
solvation of the GS solute in the RTIL. The solvent
molecules —located primarily in the first solvation
shell—therefore needs to rearrange until a new equilibrium is attained for the solvation of the ES solute.
This causes an overall time-dependent spectral narrowing and shifts towards smaller frequency. The process
is biphasic with the initial dominant component relaxing at a timescale of ∼ 50 fs followed by a very weakslower decay having a timescale of ∼ 5 ps. These
timescales are not very different from the ones observedin a 3PEPS measurement using Oxazine-4 solute as
the probe molecule. 43,44 This indicates a possibility thatthe 3PEPS measurement probes predominantly nonpolar part of the solvation. We have also used a simple
mode-coupling based approach to calculate solvationenergy-energy correlation function and then comparedwith the simulated nonequilibrium relaxation in orderto verify the theoretical prediction. 30,48,49 Breakdown of
the LRT has been viewed for nonpolar solvation in RTIL
while comparing the simulated nonequilibrium solvent
response function with equilibrium response functions,calculated using equilibrium trajectories in presence ofGS and ES solute.
Note that in this study we have excited the solute bychanging only σ of the solute atoms without affecting
their ε parameter. This is due to the fact that the LJ
potential—being a function of the 12th and 6th orderof σ respectively for repulsive and attractive part of thepotential—is much more sensitive towards the change of
σ value compared to that with the change of ε parameter,
which linearly changes the LJ interaction.
It’s worth mentioning that temperature dependentstudies for the solvation dynamics can be very important in the sense that the viscosity—a function oftemperature—can in principally affect the nonpolar solvation dynamics in general and the slower part of the
decay to be more specific. Note that while the faster part3

Page 12 of 14

of the relaxation is governed by binary solute-solvent
interaction, the longtime decay is due to the collective
fluctuation of the solvent cage, which must be coupledwith the viscosity of the RTIL solvent.
Supplementary Information (SI)
Information regarding the radial distribution functions (g(r ))between [Bmim][PF6 ] RTIL atoms in order to verify that theproper solvent equilibration, convergence of nonequilibriumand equilibrium solvation response functions, g(r ) between
solute atoms and different carbon atoms of the side chain of
the [Bmim]+ cation, and the outline of the origin of nonpolar solvation dynamics of an excited solute in a solvent areavailable at www.ias.ac.in/chemsci.

Acknowledgements
SD acknowledges financial and computational support from
Indian Institute of Technology, Patna, India. We acknowledge
the anonymous reviewers for their sincere effort in improving
the quality of the manuscript.

References
1. Welton T 1999 Room-Temperature Ionic Liquids. Solvents for Synthesis and Catalysis Chem. Rev. 99 2071
2. Hallett J P and Welton T 2011 Room-Temperature Ionic
Liquids. Solvents for Synthesis and Catalysis 2 Chem.
Rev. 111 3508
3. MacFarlane D, Tachikawa N, Forsvth M, Pringle J M,
Howlett P C, Elliott G D, Davis Jr. J H, Watanable M,
Simon P and Angell C A 2014 Energy applications ofionic liquids Energy Environ. Sci. 7 232
4. Bates E D, Mayton R D, Ntai I and Davis Jr. J H
2002 CO2 Capture by a Task-Specific Ionic Liquid J.
Am. Chem. Soc. 124 926
5. Ramdin M, de Loos T W and Vlugt T J H 2012 Stateof-the-Art of CO2 Capture with Ionic Liquids Ind. Eng.
Chem. Res. 51 8149
6. Mandal P K, Sarkar M and Samanta A 2004 ExcitationWavelength-Dependent Fluorescence Behavior of Some
Dipolar Molecules in Room-Temperature Ionic Liquids J. Phys. Chem. A 108 9048
7. Saha S, Mandal P K and Samanta A 2004 Solvation
dynamics of Nile Red in a room temperature ionic liquidusing streak camera Phys. Chem. Chem. Phys. 6 3106
8. Samanta A 2010 Solvation dynamics in ionic liquids:
What we have learned from the dynamic fluorescence
Stokes shift studies J. Phys. Chem. Lett. 1 1557
9. Samanta A 2006 Dynamic Stokes shift and excitation
wavelength dependent fluorescence of dipolar molecules
in room temperature ionic liquids J. Phys. Chem. B 110
13704
10. Arzhantsev S, Jin H, Baker G A and Maroncelli M 2007
Measurements of the Complete Solvation Response in
Ionic Liquids J. Phys. Chem. B 111 4978
11. Maroncelli M, Zhang X-X, Liang M, Roy D and Ernsting N P 2012 Measurements of the complete solvation

J. Chem. Sci. (2018) 130:3
response of coumarin 153 in ionic liquids and the accuracy of simple dielectric continuum predictions Faraday
Discuss. 154 409
12. Zhang X-X, Liang M, Ernsting N P and Maroncelli M
2013 Complete Solvation Response of Coumarin 153 inIonic Liquids J. Phys. Chem. B 117 4291
13. Zhang X-X, Liang M, Ernsting N P and Maroncelli M
2013 Conductivity and Solvation Dynamics in Ionic Liquids J. Phys. Chem. Lett. 4 1205
14. Zhang X-X, Liang M, Hunger J, Buchner R andMaroncelli M 2013 Dielectric Relaxation and Solvation
Dynamics in a Prototypical Ionic Liquid + Dipolar Protic
Liquid Mixture: 1-Butyl-3-Methylimidazolium Tetrafluoroborate + Water J. Phys. Chem. B 117 15356
15. Halder M, Headley L S, Mukherjee P, Song X and PetrichJ W 2006 Experimental and Theoretical Investigationsof Solvation Dynamics of Ionic Fluids: Appropriatenessof Dielectric Theory and the Role of DC Conductivity J.
Phys. Chem. A 110 8623
16. Chakrabarty D, Chakraborty A, Seth D, Hazra P
and Sarkar N 2004 Dynamics of solvation androtational relaxation of Coumarin 153 in 1-butyl-3methylimidazolium hexafluorophosphate [bmim][PF6 ]–water mixtures Chem. Phys. Lett. 397 469
17. Chakrabarty D, Chakraborty A, Seth D and Sarkar N
2005 Effect of Water, Methanol, and Acetonitrile on Solvent Relaxation and Rotational Relaxation of Coumarin
153 in Neat 1-Hexyl-3-methylimidazolium Hexafluorophosphate J. Phys. Chem. A 109 1764
18. Chakrabarty D, Seth D, Chakrabarty A and Sarkar N
2005 Dynamics of solvation and rotational relaxation of
coumarin 153 in ionic liquid confined nanometer-sized
microemulsions J. Phys. Chem. B 109 5753
19. Adhikari A, Sahu K, Dey S, Ghosh S, Mandal U and
Bhattacharyya K 2007 Femtosecond Solvation Dynamics in a Neat Ionic Liquid and Ionic Liquid Microemulsion: Excitation Wavelength Dependence J. Phys. Chem.
B 111 12809
20. Das S K, Sahu P K and Sarkar M 2013 Diffusion–
Viscosity Decoupling in Solute Rotation and Solvent
Relaxation of Coumarin153 in Ionic Liquids Containing
Fluoroalkylphosphate (FAP) Anion: A Thermophysicaland Photophysical Study J. Phys. Chem. B 117 636
21. Das S K and Sarkar M 2012 Studies on the Solvation Dynamics of Coumarin 153 in 1-Ethyl-3Methylimidazolium Alkylsulfate Ionic Liquids: Dependence on Alkyl Chain Length ChemPhysChem 13 2761
22. Das S K, Sahu P K and Sarkar M 2013 Probing theMicroscopic Aspects of 1-Butyl-3-Methylimidazolium
Trifluoroacetate Ionic Liquid and Its Mixture with Water
and Methanol: A Photophysical and Theoretical (DFT)
Study J. Fluores. 23 1217
23. Horng M L, Gardecki J A, Papazyan A and Maroncelli M
1995 Subpicosecond Measurements of Polar Solvation
Dynamics: Coumarin 153 Revisited J. Phys. Chem. 99
17311
24. Kashyap H K and Biswas R 2008 Dipolar Solvation
Dynamics in Room Temperature Ionic Liquids: An
Effective Medium Calculation Using Dielectric Relaxation Data J. Phys. Chem. B 112 12431
25. Kashyap H K and Biswas R 2010 Solvation Dynamics
of Dipolar Probes in Dipolar Room Temperature Ionic

J. Chem. Sci. (2018) 130:3

26.
27.
28.

29.

30.

31.

32.
33.

34.
35.
36.
37.
38.
39.

40.
41.

42.
43.

Liquids: Separation of Ion− Dipole and Dipole− Dipole
Interaction Contributions J. Phys. Chem. B 114 254
Kashyap H K and Biswas R 2010 Solvation Dynamics
in Imidazolium and Phosphonium Ionic Liquids: Effectsof Solute Motion Ind. J. Chem. 49A 685
Kashyap H K and Biswas R 2010 Stokes’ Shift Dynamics in Imidazolium Ionic Liquids: Temperature Dependence J. Phys. Chem. B 114 16811
Daschakraborty S and Biswas R 2011 Stokes’ ShiftDynamics in Alkylimidazolium Aluminate Ionic Liquids: Domination of Solute-IL Dipole-Dipole Interaction Chem. Phys. Lett. 510 202
Daschakraborty S and Biswas R 2012 Stokes shift
dynamics of [Na][TOTO] - a new class of ionic liquids:
A comparative study with more common imidazolium
analogs Chem. Phys. Lett. 545 54
Daschakraborty S and Biswas R 2012 Ultrafast Solvation
Response in Room Temperature Ionic Liquids: Possible
Origin, and Importance of the Collective and the Nearest
Neighbour Solvent Modes J. Chem. Phys. 137 114501
Daschakraborty S, Pal T and Biswas R 2013 Stokes shift
dynamics of ionic liquids: Solute probe dependence, andeffects of self-motion, dielectric relaxation frequency
window, and collective intermolecular solvent modes J.
Chem. Phys. 139 164503
Pal T and Biswas R 2014 Slow solvation in ionic liquids:
Connections to non-Gaussian moves and multi-point correlations J. Chem. Phys. 141 104501
Shim Y and Kim H J 2008 Dielectric Relaxation, Ion
Conductivity, Solvent Rotation, and Solvation Dynamics in a Room-Temperature Ionic Liquid J. Phys. Chem.
B 112 11028
Shim Y and Kim H J 2013 Dielectric Relaxation and Solvation Dynamics in a Room-Temperature Ionic Liquid:
Temperature Dependence. J. Phys. Chem. B 117 11743
Terranova Z L and Corcelli S A 2013 On the Mechanismof Solvation Dynamics in Imidazolium-Based Ionic Liquids. J. Phys. Chem. B 117 15659
Song X 2009 Solvation dynamics in ionic fluids: An
extended Debye–Hückel dielectric continuum model J.
Chem. Phys. 131 044503
Kobrak M N 2007 A comparative study of solvation
dynamics in room-temperature ionic liquids J. Chem.
Phys. 127 184507
Shim Y, Duan J, Choi M Y and Kim H J 2003 Solvationin Molecular Ionic Liquids J. Chem. Phys. 119 6411
Shim Y, Choi M Y and Kim H J 2005 A Molecular Dynamics Computer Simulation Study of RoomTemperature Ionic Liquids. II. Equilibrium and Nonequilibrium Solvation Dynamics J. Chem. Phys. 122 044511
Kobrak M N 2006 Characterization of the Solvation
Dynamics of an Ionic Liquid via Molecular Dynamics
Simulation J. Chem. Phys. 125 064502
Schmollngruber M, Schröder C and Steinhauser O
2013 Polarization Effects on the Solvation Dynamics of
Coumarin C153 in Ionic Liquids: Components and Their
Cross-Correlations J. Chem. Phys. 138 204504
Wu E C and Kim H J 2016 MD Study of Stokes Shifts inIonic Liquids: Temperature Dependence J. Phys. Chem.
B 120 4644
Muramatsu M, Nagasawa Y and Miyasaka Y 2011 Ultrafast Solvation Dynamics in Room Temperature Ionic

Page 13 of 14

44.
45.

46.
47.
48.
49.

50.

51.

52.
53.

54.

55.

56.
57.
58.

59.

60.

3

Liquids Observed by Three-Pulse Photon Echo Peak
Shift Measurements J. Phys. Chem. A 115 3886
Nagasawa Y and Miyasaka Y 2014 Ultrafast solvation
dynamics and charge transfer reactions in room temperature ionic liquids Phys. Chem. Chem. Phys. 16 13008
Ma J, Bout D V and Berg M 1995 Transient hole burningof s-tetrazine in propylene carbonate: A comparison of
mechanical and dielectric theories of solvation J. Chem.
Phys. 103 9146
Fourkas J T, Benigno A and Berg M 1993 Time-resolved
nonpolar solvation dynamics in supercooled and low viscosity n-butylbenzene J. Chem. Phys. 99 8552
Fourkas J T and Berg M 1993 Temperature-dependent
ultrafast solvation dynamics in a completely nonpolarsystem J. Chem. Phys. 98 7773
Bagchi B 1994 Molecular theory of nonpolar solvation
dynamics J. Chem. Phys. 100 6658
Biswas R, Nandi N and Bagchi B 1997 Solvation dynamics in monohydroxy alcohols: Agreement between theory and different experiments J. Phys. Chem. B 101
2968
Reynolds L, Gardecki J A, Frankland S J V, Horng M L
and Maroncelli M 1996 Dipole solvation in nondipolar
solvents: experimental studies of reorganization energies
and solvation dynamics J. Phys. Chem. 100 10337
Saven J S and Skinner J L 1993 A molecular theoryof the line shape: Inhomogeneous and homogeneouselectronic spectra of dilute chromophores in nonpolar
fluids J. Chem. Phys. 99 4391
Ladanyi B M and Stratt R M 1996 Short-Time Dynamics
of Solvation: Relationship between Polar and Nonpolar
Solvation J. Phys. Chem. 100 1266
Ladanyi B M and Maroncelli M 1998 Mechanismsof solvation dynamics of polyatomic solutes in polarand nondipolar solvents: A simulation study J. Chem.
Phys. 109 3204
Pal T and Biswas R 2015 Composition dependenceof dynamic heterogeneity time-and length scales in[Omim][BF4]/water binary mixtures: Molecular dynamics simulation study J. Phys. Chem. B 119 15683
Carter E A and Hynes J T 1991 Solvation dynamics for an
ion pair in a polar solvent: Time-dependent fluorescenceand photochemical charge transfer J. Chem. Phys. 94
5961
Bernardi E, Martins M M and Stassen H 2005 The breakdown of linear response theory in non-polar solvation
dynamics Chem. Phys. Lett. 407 171
Faeder J and Ladanyi B M 2001 Solvation Dynamics
in Aqueous Reverse Micelles: A Computer Simulation
Study J. Phys. Chem. B 105 11148
Schmidt M W, Baldridge K K, Boatz J A, Elbert S T, Gordon M S, Jensen J H, Koseki S, Matsunaga N, Nguyen
K A, Su S, Windus T L, Dupuis M and Montgomery
Jr. J A 1993 General Atomic and Molecular Electronic
Structure System J. Comput. Chem. 14 1347
Martínez L, Andrade R, Birgin E G and Martínez J M
2009 Packmol: A package for building initial configurations for molecular dynamics simulations J. Comput.
Chem. 30 2157
Lopes J N C, Deschamps J and Pádua A A H 2004 Modeling Ionic Liquids Using a Systematic All-Atom Force
Field J. Phys. Chem. B 108 2038

3

Page 14 of 14

61. Jorgensen W L, Maxwell D S and Tirado-Rives J
1996 Development and Testing of the OPLS All-AtomForce Field on Conformational Energetics and Properties of Organic Liquids J. Am. Chem. Soc. 118
11225
62. Schroder C 1996 Comparing reduced partial charge models with polarizable simulations of ionic liquids Phys.
Chem. Chem. Phys. 118 11225
63. van der Spoel D, Lindahl E, Hess B, Groenhof G, MarkA E and Berendsen H J C GROMACS: Fast, Flexible,
and Free J. Comput. Chem. 26 1701
64. Nosé S 1984 A unified formulation of the constant temperature molecular-dynamics methods J. Chem. Phys. 81
511
65. Hoover W G 1985 Canonical dynamics: Equilibriumphase-space distributions Phys. Rev. 31 1695
66. Berendsen H J C, Postma J P M, van Gunsteren W
F, DiNola A and Haak J R 1984 Molecular-Dynamics
with Coupling to an External Bath J. Chem. Phys. 81
3684
67. Darden T, York D and Pedersen L 1993 Particle mesh
Ewald: An N ·log(N ) method for Ewald sums in large
systems J. Chem. Phys. 98 10089
68. Maginn E J 2009 Molecular simulation of ionic liquids:
current status and future opportunities J. Phys.: Condens.
Matter 21 373101
69. Pal T and Biswas R 2013 Rank-dependent orientational
relaxation in an ionic liquid: an all-atom simulation
study Theo. Chem. Acc. 132 1348

J. Chem. Sci. (2018) 130:3
70. Lopes J N C, Deschamps J and Pádua A A H 2006 Nanostructural organization in ionic liquids J. Phys. Chem.
B 110 3330
71. Morgado P, Shimizu K, José M S S, Esperança J M S S,
Reis P M, Rebelo L P N, Lopes J N C and Filipe E J M
2013 Using 129 Xe NMR to Probe the Structure of Ionic
Liquids J. Phys. Chem. Lett. 4 2758
72. Hu Z and Margulis C J 2006 Heterogeneity in a roomtemperature ionic liquid: Persistent local environments
and the red-edge effect Proc. Natl. Acad. Sci. U.S.A. 103
831
73. Richert R 2002 Heterogeneous liquids: fluctuations inspace and time J. Phys. Condens. Matter. 14 R703
74. Jin H, Li X and Maroncelli M 2007 Heterogeneous
Solute Dynamics in Room Temperature Ionic Liquids J.
Phys. Chem. B 111 13473
75. Wang Y and Voth G A 2005 Unique Spatial Heterogeneity in Ionic Liquids J. Am. Chem. Soc. 127 12192
76. Triolo A, Russina O, Bleif H-J and Cola E D 2007
Nanoscale Segregation in Room Temperature Ionic Liquids J. Phys. Chem. B 111 4641
77. Zheng Z-P, Gun W-H, Roy S, Mazur K, Nazet A, Buchner R, Bonn M and Hunger J 2015 Ionic Liquids:
Not only Structurally but also Dynamically Heterogeneous Angew. Chem. Int. Ed. 54 687
78. Roy D and Maroncelli M 2010 An improved four-siteionic liquid model. J. Phys. Chem. B 114 12629
79. Egorov S A 2003 Nonpolar solvation dynamics in supercritical fluids J. Chem. Phys. 118 10643