Multiscale Studies on Ionic Liquids - Chemical Reviews (ACS

May 10, 2017 - (1) a new force field model for ionic liquids and the multiscale simulation ... She was awarded the Nomination Award of the fourth Top ...
3 downloads 0 Views 25MB Size
Review pubs.acs.org/CR

Multiscale Studies on Ionic Liquids Kun Dong, Xiaomin Liu, Haifeng Dong, Xiangping Zhang, and Suojiang Zhang*

Downloaded via UNIV OF TOLEDO on June 30, 2018 at 11:26:40 (UTC). See https://pubs.acs.org/sharingguidelines for options on how to legitimately share published articles.

State Key Laboratory of Multiphase Complex Systems, Beijing Key Laboratory of Ionic Liquids Clean Process, Institute of Process Engineering, Chinese Academy of Sciences, Beijing 100190, China ABSTRACT: Ionic liquids (ILs) offer a wide range of promising applications because of their much enhanced properties. However, further development of such materials depends on the fundamental understanding of their hierarchical structures and behaviors, which requires multiscale strategies to provide coupling among various length scales. In this review, we first introduce the structures and properties of these typical ILs. Then, we introduce the multiscale modeling methods that have been applied to the ILs, covering from molecular scale (QM/MM), to mesoscale (CG, DPD), to macroscale (CFD for unit scale and thermodynamics COSMO-RS model and environmental assessment GD method for process scale). In the following section, we discuss in some detail their applications to the four scales of ILs, including molecular scale structures, mesoscale aggregates and dynamics, and unit scale reactor design and process design and optimization of typical IL applications. Finally, we address the concluding remarks of multiscale strategies in the understanding and predictive capabilities of ILs. The present review aims to summarize the recent advances in the fundamental and application understanding of ILs.

CONTENTS 1. Introduction 2. Ionic Liquids and Multiscale Strategy 2.1. Structures of Typical ILs 2.2. Properties of Typical ILs Transport Properties Surface Property 2.3. Multiscale Strategy for ILs 3. Multiscale Modeling Methods for IL Systems 3.1. Molecular Scale Methods 3.2. Mesoscale Methods 3.2.1. CG Method 3.2.2. DPD Method 3.2.3. MM&Meso Bridging 3.3. CFD Model for Unit Scale 3.4. Process Scale Methods 3.4.1. COSMO-RS Model 3.4.2. Green Degree Method 4. Molecular Scale Structures of ILs 4.1. Single Ions and Ion Pairs 4.1.1. Structures of Single Ions 4.1.2. Structures, Interaction Energies, and Charge Distributions of Ion Pairs 4.1.3. Ion Pair Lifetime 4.2. H-Bond and Network 4.2.1. H-Bond Structures 4.2.2. H-Bonding Network 4.2.3. Effect on Properties of ILs 4.3. QM/MM Illustrating Reaction Mechanism 5. Mesoscale Structures of ILs 5.1. Bulk Structures of ILs 5.1.1. Self-assembly in Neat ILs 5.1.2. Aggregates in Aqueous Solutions © 2017 American Chemical Society

5.1.3. Aggregates in Nonaqueous Solutions 5.2. Interface Structures 5.2.1. Structures of Liquid−Gas Interfaces 5.2.2. Structures of Liquid−Solid Interfaces 5.3. Comparison of Structures to Other Solvents 6. Unit Scale Design Based on ILs Media 6.1. Gas−Liquid Two-Phase Flow Hydrodynamics 6.2. Mass Transfer 6.3. Heat Transfer 6.4. Reactor Design 7. Process Scale Design for ILs Applications 7.1. Thermodynamic Properties Prediction 7.2. Environmental Assessment 7.2.1. GD Analysis 7.2.2. LCA Analysis 7.3. Objective Optimization Analysis 7.4. Process Simulations for IL Key Applications 7.4.1. CO2 Capture and Conversion 7.4.2. Biomass Conversion 7.4.3. Separation Process 8. Concluding Remarks Author Information Corresponding Author ORCID Notes Biographies Acknowledgments Glossary References

6637 6638 6638 6639 6640 6641 6642 6643 6643 6644 6644 6645 6645 6646 6646 6646 6647 6647 6647 6647 6647 6648 6648 6648 6649 6650 6651 6652 6652 6652 6654

6656 6658 6658 6659 6660 6662 6662 6664 6665 6665 6667 6668 6668 6668 6668 6669 6669 6669 6672 6674 6675 6676 6676 6676 6676 6676 6676 6676 6678

Special Issue: Ionic Liquids Received: November 25, 2016 Published: May 10, 2017 6636

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

1. INTRODUCTION Serious drawbacks related to environmental, health, or safety issues have forced our society to create evidence-based solutions from chemistry research, with strong attention to solvents. Ionic liquids (ILs) are a new, innovative class of solvents with huge potential to have an impact across many areas of scientific and engineering research. ILs are fluids that exemplify some excellent properties, and thus, they harvest intense interest and continuing success. More importantly, the sheer number of ILs (estimated at as high as 1018 ILs accessible,1,2 including their binary and ternary mixtures3−6) provides enormous scope for scientific innovation. Before ILs, some inorganic molten salts and eutectics were significantly developed to replace volatile organic compounds (VOCs) in chemical industries to resolve the deteriorating environmental problem.7 In 1914, Walden8 synthesized an ionic salt, ethylammonium nitrate, that displayed a low melting point of 12 °C and was the first defined IL. In the 1970s−80s, the alkylimidazolium halogenoaluminate ILs were reported by mixing aluminum halides with theimidazolium halides, and these were also termed first generation ILs.9,10 Besides being applied in electrochemistry, it was also attempted to apply these liquid salts in catalytic reactions,11 such as Diels−Alder reactions and alkylation of sodium β-naphthoxide.12,13 Until the middle of the 1990s, air-stable second generation IL salts containing tetrafluoroborate ([BF4]−), hexafluorophosphate ([PF6]−), etc.14,15 were synthesized. As these ILs are easily stored and handled outside an air environment, the interest in them followed an exponential growth. According to the ISI Web of Knowledge, the growth rate of published papers related to ILs in 2000−2016 is 130%, and the papers titled by ILs in 2014, 2015, and 2016 are 5156, 4617, and 4752 respectively. In recent years, a great number of novel ILs were synthesized.16,17 Many property databases, such as the frequently used ILTHERMO v2.0,18 INFO thermophysical properties database,19 IUPAC database, UFT Merck database (seeming offline), DDB search, and our IPE IL database,20 have been built. The properties can involve thermodynamic, thermochemical, electrochemistry, and transport aspects. The book titled Ionic Liquids: Physicochemical Properties by Zhang et al.21 summarizes the unique physical properties of the ILs, such as density (ρ), viscosity (η), phase transition temperatures (melting point, Tm, glass transition, Tg, boiling point, Tb), thermal decomposition temperature, polarity, and electrochemical windows. With the unique properties, ILs have been used in diverse ranges of technologies and applications currently, involving making things (catalysis,22−25 medicine,26 plastics27), processes (mineral extraction,28−30 CO2 capture,31−36 cellulose treatment 37−39 ), electrochemistry (battery, 40−42 supercapacitor43−45), materials (metal, nanotube,46−49 graphete48,50−52), and transporting things (lubricants,53 fuels54). Even with the increasing miniaturization of science and devices,55,56 the ILs remain a convenient and widely used medium in which to do chemistry. In the early comprehensive review by Welton,11 the catalysis of chloroaluminate ILs formed by mixing aluminum halides with the corresponding imidazolium or pyridinium cations was introduced. However, the ILs were highly sensitive to atmospheric moisture, which prevented their industrial applications. In the updated review in 2011, Welton et al.23 found that the number of ILs has been proliferated to the end of 2009 with the primary literature surveyed up; in particular, the second generation ILs, including [BF4]−, [PF6]−, and [NTf2]− anions,

have been used in more catalytic reactions. Those ILs combining the quaternary ammonium cations and the coordinating anions have been demonstrated to be ideal immobilizing agents for various transition-metal catalyst precursors in the reactions ranging from the Ziegler−Natta type to the hydroformylation of olefins.22 In some organic reactions, especially in the heterocyclic synthesis, the functional, acid, and zwitterionic ILs were used as organocatalyst.16,57 The important and promising applications are that the ILs were used in homogeneous catalytic reactions, such as in the hydrogenation,58,59 oxidation,60−65 carbonylation,66,67 Heck reaction,68−70 and hydroformylation,71 and the ILs “mill” the catalyst to form nanoparticles, enabling more efficient reactions to take place compared with molecular solvents. For some time, ILs have been beyond the range of solvents and bore the important responsibility for the alternative cleaner and green technologies. In recent years, the potential applications in electrochemistry, gas absorption, and biology (biomass) engineering have promoted their rapid progress. They are very good candidates as electrolyte with higher conductivities in the broad range of 0.1−20 mS·cm−1 and wider electrochemical windows of 4.5−5 V.72 The fundamental understanding of the liquid−solid interface between IL electrolytes and electrodes was reviewed by Fedorov et al.73 In gas separation, especially on CO2 capture, ILs showed potential applications. Lei et al.74 reviewed the experimental and modeling progress based on the ILs gas absorption. In the biomass pretreatment, ILs also showed their advantages with respect to the environment. Pinkert et al.37 specially discussed the interaction between ILs and cellulose in their Chemical Reviews paper to discover the molecular mechanism for the high dissolution and benign recyclability. One key area of research ripe for the investigation linking to properties and applications is structures of the ILs. More recently, Hayes et al.75 made a thorough review on the structures of ILs. The self-assembled nanoaggregates in bulk phases and at solid interfaces were integrally discussed. The similarities and differences between the ILs and the traditional solvents (water, benzene, alcohol, supercritical fluids and melt salts, electrolytes, etc.) were compared. The driving interests for the structural aspects of ILs can be seen even before they were widely appreciated historically. In the past decade or so, advances in theoretical, computational, and experimental techniques have provided unprecedented resolution pictures for structures of ILs. For this review, an overall picture about structures, properties, and applications of ILs is presented. However, when one really plans to design a technology or a process using ILs, the open question is “what is the right IL for my requirement?” The most common way to answer this question likely refers to the structures of the ILs. However, for today’s chemical engineers, understanding structure only plays as a partial role, and the reactor innovation and process design (optimization) assume other necessary constituents. Such deployment embodies a “multiscale philosophy”. A general definition for multiscale study is the bridging of length (spatial) and time (temporal) scales to describe the macroscopic properties and behaviors based on microscopic structures. In fact, the multiscale covers the entire design chain from molecule, to mesoscale, to reactor, even a system, which is believed to be the most promising methodology for handling the complex chemical systems.76,77 Recently, Xu et al.78 reviewed the multiscale application of ILs on fixation of CO2 into cyclic carbonates. Three scales were identified for the whole process: molecule scale, unit scale, and system scale. Such identification is an attempt at multiscale 6637

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

Figure 1. Chemical structures of the typical IL cations and anions. Cations. Neutral: 1,3-dialkylimidazolium ([CnCnim]+), 1,2,3-trialkylimidazolium ([CnCnCnim]+), N-alkylpyridinium ([Cnpy]+), 1,1-dialkylpyrrolidinium ([CnCnPyr]+), 1,2,3-triazolium, 1,2,4-triazolium, 1,2,3-trialkylsulfonium, tetraalkylammonium ([Nn,n,n,n]+), t(i)tetraalkylphosphonium ([Pn,n,n,n]+), 1,1,3,3-tetramethylguanidinium ([TMG]+), and di-imidazolium ([Cn(mim)2]+). Acidic: 3-methyl-1-(4-sulfo-butyl)-benzimidazolium, 3-propylsulfonic-triethylammonium, 3-propylsulfonic-triethylphosphonium, 1(3-propylcarboxyl)-3-alkylimidazolium, 1-(2-benzenesulfonic)-3-alkylimidazolium, 3-propylsulfonic-tris(2-hydroxylethyl)ammonium. Anions. Neutral: halides (X−), nitrate ([NO3]−), tetrafluoroborate ([BF4]−), hexafluorophosphate ([PF6]−), tricyanomethide, methylsulfate, trifluoromethanesulfonates ([OTf]−), thiocyanate, dicyanamide, trifluoroacetyltrifluoromethanesulfonamide, bis(trifluoromethylsulfonyl)imide ([NTf2]−), benzenesulfate. Acidic: hydrogensulfate, hydrogen phosphate, aluminum chloride, tetrachloroferrate, p-toluenesulfonic. Basic: dicyanamide, formate, acetate, phenolic.

In light of the recent successful cases, we imagine that the multiscale study is the tendency for more physicochemical and engineering data, improving theoretical prediction tools, and developing dedicated process units and devices. We hope that this review will attract more engineers to this highly exciting field of research.

utilization on IL-based application. However, much more information is necessary to evaluate the efficiency of ILs in an engineering process generated only in a more advanced state of development in which more attention can be given to the aspects such as recycling efficiency, recovery, or the degradation degree of ILs. Chemical engineering has evolved toward multiscale design,76 and the related knowledge has been extended from focusing on micromechanisms to understanding macroscale material distribution and individual phenomena.79,80 Although there is absolutely no doubt that the unique properties of ILs offer great potential to improve the existing engineering applications, up to date, the multiscale strategy has not been truly applied to the studies of ILs. The studies of ILs are still compartmentalized, and molecule is molecular and the mesoscale structures are little considered, which frustrates the industrial development of ILs. Despite the advances made in the modeling and simulation for the structures, thermodynamics, mechanism, and transport of ILs, a comprehensive review about the multiscale studies on ILs, especially based on industrial interests, has not been reported. In this review, we discussed IL fundamental and application investigations based on the multiscale modeling and the related experimental work. First, the structures and properties of the typical ILs are introduced. Then, we introduce the multiscale modeling methods that have been applied to the ILs, covering from molecular scale (QM/MM), to mesoscale (CG, DPD), to macroscale (CFD for unit scale and thermodynamics COSMORS model and environmental assessment GD method for process scale). In the following section, we discuss in some detail their applications to the four scales of ILs, including molecular scale structures, mesoscale aggregates and dynamics, and unit scale reactor design and process design and optimization of typical IL applications. Finally, we address the concluding remarks of multiscale strategies in the understanding and predictive capabilities of ILs.

2. IONIC LIQUIDS AND MULTISCALE STRATEGY 2.1. Structures of Typical ILs

One advantage of using ILs in an industrially process is that ILs have no detectable vapor pressure at ordinary conditions, and therefore contribute no VOCs to the atmosphere. Another is the great number of ILs possible: that is, at least thousands of neat ILs, a million binary, and 1018 ternary mixtures are potentially possible,1 which enables the ILs to be designed and tuned81 for optimized selectivity, substrate solubility, product separation, and enantioselectivity. Figure 1 shows the structures of the typical IL cations and anions. Though they have occurred in many review papers, here we emphasize the necessity, since these ILs are frequently used in most IL-related applications.82 It is found that the cations are organic with larger molecular volumes and the anions are inorganic with better geometric symmetry in Figure 1. We used their abbreviations in this review. [Bmim]+ and [BMIM]+ both refer to the 1-butyl-3-methylimidazolium cation. Sometimes, ambiguity exits; for example, [Pmim]+ may refer to the 1-propyl3-methylimidazolium cation, but also to 1-pentyl-3-methylimidazolium. To circumvent these ambiguities, an alphabetic abbreviation was used to describe the alkyl chains with the charged center; thus, the 1-butyl-3-methylimidazolium cation is expressed by [C4C1im]+. For functionalized ILs, the type and position of the functional group should be noted; for example, [(HO)4C4C1im]+ denotes that there is an alcohol group on the terminal carbon of the butyl chain. Alkyl chains on the cations are assumed to be saturated. For cyclic cations, such as the 1-butyl2,3-dimethylimidazolium, the alkylation is assumed to be on the 6638

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

Table 1. Properties of the Typical ILsa Typical ILs [C2C1im]X

[C4C1im]X

[C2C1im]Cl/AlCl3(1:1)110 [C2C1im]Cl/AlCl3(1:2)110 [C4C1im][BF4]111 [C4C1im][PF6]111 [C4C1im][CF3CO2]111 [C4C1im][OTf]111 [C4C1im][NTf2]111 [C6C1im][NTf2]115 [C8C1im][NTf2]115 [C4pyr][BF4] [C4pyr][NTf2] [N2226][NTf2] [P66614][NTf2] [Mptrim][OTf]b132 [Mptrim][NTf2]c132 [EtNH3][NO3]83,133 [EtNH3][HSO4]83,133 [EtNH3]lactate83,133 [EA][NO3]d [TMG][OTf]e83,136 [C4C1im][HSO4]137 [(HSO3)3C3C1im][HSO4]16 [(HSO3)4C4C1im][HSO4]16 [(HSO3)4N4][HSO4]138 [(HSO3)3C3C4bim][HSO4]139 [(HSO3)3C3C4bim][CF3SO3]139

Tm (Tg, Td, Tb)

C

ρ

X = Cl 361.75101 X = Br 352.0102 X = I 351.0103 X = Cl 314.1106 X = Br 346.1 X = I 188.1107 280.15 177.15 190.15 283.15 (196.15) 195.15 290.15 270.15 (168.15) 231.15 (192.15) 193.15 288.45 (202.0)118 299.0 (197.0)120 287.1124 206.1128 201.55 201.55 282.15 (181.65)134 313.15 (189.15) 216.15 324.15 (191.15) 314.15 349.00n 579.15p 584.15p 543.15p 682.00p 747.10p

235.0104 264.8f

1.1241105

48.3105

323.0108 317.0108

1.0800106 1.3000107 1.4600107 1.2941 1.3887 1.2035 1.3709 1.2152 1.3015 1.4397 1.6400 1.5902 1.2100j119 1.4360120 1.2886125 1.0661130

30.0105 1500.0107 400.0107 17.80 19.11 100.0 258.0 73.0 83.0 49.0 68.0 91.0 163.0119 61.0121 187.0126 336.7130 129.0 87.0 32.0 128.0 803.0 113.0 167.0 2147.31o 1296.0 1427.0 9450.0o 35.72o 42,895o

363.1 407.6 417.0114 536.0114 631.6116 692.7117 383.0114 596l117 1360m129

205.37135

1.5200 1.2160 1.4380 1.1100 1.2650 1.2900 1.2815o 1.5017 1.4647 1.3532o 1.0753o 1.4545o

η

σ

γ

0.249101 10.6g102 10.0h102 1.0i107

54.70109

22.6 14.7 3.50 1.46 2.80 2.90 3.90 2.20 1.30 2.5k120 3.07122 0.887127

43.84112 37.28113

37.0109

33.40123 34.10126 30.90131

4.10 4.90 26.90 4.40 0.26 9.35 2.68

47.30 56.30 39.30 50.60 48.36o

7.57o

46.79o 35.60o 45.35o

a

Temperature is 298 K and pressure is 1 atm; Tm melting point, Tg glass transition, Tg decomposition temperature and Tb boiling temperature, K; C heat capacity, J·K−1 mol−1; ρ density, g·m−3; η viscosity, cP; σ conductivity, mS·m−1; γ surface tension, mN·m−1. b[Mptrim][OTf] 4-methyl-1propyl-1,2,4-triazolium [OTf]. c[Mptrim][NTf2] 4-methyl-1-propyl-1,2,4-triazolium [NTf2]. d[EA][NO3] ethanolammonium [NO3]. e[TMG][OTf]. f350 K. g335.6 K. h348.1 K. i313.16 K. j312.8 K. k303 K. l333.15 K. m333.15 K. nIt is decomposition temperature (Td) by thermal gravimetric analysis (TGA). o303.5 K. pThe calculated boiling temperature (Tb), the calculated equation adapted from Shereshefsky.140

achiral ILs and chiral ILs.85−87 In addition, magnetic ILs with a paramagnetic atom/group,88 polymeric ILs with a polymerizable ion,89 divalent ILs, or solvate ILs90 are identified as the new class of ILs. We also found the reports about amino acid ILs91,92 and aryl alkyl ILs named from specific functional groups introduced on the ions. The synthesis of the ILs enables the crossfertilization of chemistries and guiding future research into these ILs types.

heteroatom(s) and can be expressed by [C4C1C1im]+. Nonimidazolium cations, such as tetraalkylammonium in Figure 1, are simply expressed as [Nn,n,n,n]+ (n is the number of carbons in each alkyl chain). For the cations with special structures, such as guanidinium ([TMG]+) and dialkylimidazolium ([Cn(mim)]2+), the abbreviation is expressed alphabetically. The anion abbreviation is more structural; for example, the nitrate, tetrafluoroborate, and hexafluorophosphate can be abbreviated as [NO3]−, [BF4]−, and [PF6]−, respectively. ILs are usually classified depending on the functional group considered most important. At present, protic ILs (PILs)17,83 can be identified from the reported ILs based on the well-established division between proton-donating ILs and other ILs that are nonproton-donating, such as aprotic ILs (AILs). Figure 1 shows the most-used PILs, including primary, secondary, or tertiary ammonium, mono- or di-imidazolium ions, triazolium, and guanidinium, and the coupled anions including [NO3]−, [NTf2]−, and methylsulfate. However, the definition is not as rigid as it would seem; Mirjafari et al.84 characterized dicationic ILs with one protic and one aprotic charge center, enabling both functionalities to be expressed. Based on the presence or absence of stereogenic groups, ILs can also be divided into two categories:

2.2. Properties of Typical ILs

ILs present the dramatic varied properties along with different structures of them. The structure−property relationship is the fundamental understanding for the behaviors of ILs.16,17,75,93 A survey of those reviews has shown that the IL fundamentals in relation to chemical structure are reasonably advanced, and some works were investigated by the multiscale modeling,93−97 and approximate98,99 and classical74,100 theories, although ion arrangements are once again ignored in these theories because the liquid is modeled as a continuum of ions of fixed molecular volume and density. Table 1 lists the properties of a series of typical ILs. For these properties, we will then discuss the relationship between these properties and IL structures. 6639

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

[PF6] > [BF4] ≈ [NTf2] > [OTf] > [CF3CO2]. Therefore, the thermal stability of the ILs depends on the anionic species. As a basic thermodynamic property of liquids, the heat capacity (Cp) of some ILs has been measured by calorimetry. The types of calorimeters used include differential scanning calorimeters (DSCs) and modulated DSCs (MDSCs) as the most popular instruments, but their uncertainty is at least 5%. More accurate data can be obtained with adiabatic and Tian− Calvet DSC calorimeters. For the imidazolium- and pyridiniumbased series of ILs, the heat capacity increases with the number of atoms in the alkyl chain, as observed in Table 1. A linear trend is also observed for the numbers of carbon atoms of n = 3−14.114 A 32.1 ± 0.3 J K−1 mol−1 value of ΔCp (H → CH3) at 298.15 K was calculated with the increments per CH2 group for various series of ILs. Similar correlations were found at other temperatures. The constant value of the increment in different series over a wide range of n confirms the validity of additive procedures toward estimating the heat capacity of ILs. The observation means that there is not only the increment with length of the alkyl chain, but also the increment with interaction between the different cations and anions. Transport Properties. Viscosity (η) is the most often mentioned transport property of ILs, involved in electrochemical and separating applications. Most ILs show high viscosities, but generally the viscosity values of ILs vary between 20 and 40000 cP; even for the low viscous ILs, their viscosities are around 20− 30 cP, that is 1 order of magnitude higher than that of conventional solvents, such as water, with a viscosity of ∼1 cP.146 From Table 1, it is seen that first-generation chloroaluminatebased ILs show low viscosity η, such as 17.8 cP for [C2C1im][AlCl4] ([C2C1im]Cl/AlCl3(1:1) mixing) and 19.11 cP for [C2C1im][Al2Cl7] ([C2C1im]Cl/AlCl3(1:2) mixing);110 thus, the ILs have favorable electrochemical properties for the process, such as electrodeposition.147 In comparison, the hydrophobic ILs with [PF6]− and [NTf2]− anions show high viscosity η; even the viscosity of [C4C1im][PF6] reaches 258.0 cP. A general trend is that the viscosity decreases with increasing size of ions, but is more anion-dependent. For cations (e.g., imidazolium and pyrodinium), the viscosity increases with increasing length of the alkyl chains through stronger van der Waals (VDW) interactions. As listed in Table 1, the viscosity η of the [CnC1im][NTf2] increase from 49.0, to 68.0, to 91.0 cP when n increases from 4, to 6, to 8. Okoturo et al.148 found that viscosities are not linear for Arrhenius law but better fitted to the Vogel−Tammann−Fulcher (VTF) equation for the temperature dependency. The ILs with asymmetric cations obey Arrhenius law, and the ILs with small and symmetric cations are not well fitted to the VTF equation; however, the functional ILs with, e.g., −OH and C−O groups do not comply with both. Except for the size and shape of ions, the viscosity is closely related with intermolecular interaction. Zhang et al.149 found that the high interaction energy (electrostatic force and VDW interaction) between anion and cation led to a higher viscosity. In addition, the H-bond is responsible for ion stacking, and has an effect on the viscosity,150,151 and the trend is the opposite of smaller Hbonding interaction leading to a higher viscosity.152 It is also true for the PILs that the viscosity is dependent on the ion−ion interactions, and greater interactions lead to higher viscosities. Alkyl chains increase the viscosity through stronger VDW interactions, while delocalization of the charge on the anion decreases the viscosity by weakening the H-bond. For the PILs with alkylammonium cations, the viscosity increased with increasing alkyl chain length and significantly increased with

Thermal Properties are particularly important to classify the ILs from the conventional salts. Most ILs at normal atmospheric pressure stay in the liquid state up to temperatures as high as 200−300 °C. According to Aparicio et al.141 the thermal properties, for example Tm, Tg, and Td, depend on a balance between cation and anion symmetry, flexibility of chains, and charge accessibility. The halogen-based ILs with C2 and C4 alkyl chains show higher Tm than those non-halogen analogues, but the iodide one shows an exception of low Tm of 188 K. When mixing with aluminum chloride, haloaluminate ILs show decreased Tm. Compared with the Tm and Tg of [CnC1im][NTf2] (n = 4, 6, 8), increasing the length of alkyl chains in the cations generally decreases the Tm. The PILs show similar Tm values to those for the AILs. Nevertheless, in Greaves et al.’s review,83 a great number of PILs, including ammonium, mono- or di-imidazolium ions, caprolactam, and guanidinium cations, exhibit a great gap in Tm. The Tm values of alkylammonium salts are >100 °C, but the PILs that contained alkylammonium cations with small carboxylic acids and most of the salts with an imidazolium cation have low Tm. Among those, particularly the PILs containing [NTf2]− generally have low melting points, which is attributed to the relatively large radius of these fluorinated anions that weakens the electrostatic interaction with the cation.142 Acidic ILs can be classified into the Lewis acidic ILs due to electronic deficiency and Brönsted acidic ILs due to the protons according to the type of acid function on ions.16 Santos et al.143 studied the Tg of a series of systems of the type [C1Him][NTf2], [C2Him][NTf2], [C1C1Him][NTf2], and [C4C1C1im][NTf2] and established a linear correlation between the Tg and the alkyl chain size. The methylation at the C2 position as well as the N−H acidic group at the imidazolium ring contribute to a significant variation in the cation−anion interactions and their dynamics, which is reflected in their charge distribution and polarizability, leading to a significant differentiation of the refractive indices, surface tension, and heat capacities. For these sulfonic acid functionalized ILs, the decomposition temperatures Td were measured and deceased in the order methylimidazolium > triethanolammonium > pyridinium.144 The methylimidazolium ILs such as 1-(4-propylsulfonic)-3-methylimidazolium hydrogensulfate ([(HSO3)3C3C1im][HSO4]) and 1-(4-butylsulfonic)3-methylimidazolium hydrogensulfate ([(HSO3)4C4C1im][HSO4]) in Table 1 show decomposition onset temperatures (air) in the 213−353 °C ranges. Additionally, it is found that the Td of these ILs are highly dependent on the nature of the anion. The thermophysical properties of the IL 1-(4-butylsulfonic)-3methylimidazolium chloride [(HSO3)4C4C1im][Cl] have been compared with the [C4C1im]Cl, and the acidic ILs decomposed in three steps with DTG peaks centered at 269, 373, and 387 °C, whereas the neutral [C4C1im]Cl decomposed in one step. Actually, the glass transition is indicative of the cohesive energy within the salt, which is decreased by repulsive Pauli forces from the overlapping of closed electron shells and increased through the attractive Coulombic and van der Waals interactions.145 According to the DSC measurements,111 a heat capacity change occurs for most ILs, corresponding to a Tg, and the Tg values of the ILs in Table 1 are fairly low; for example, the Tg of [C4C1im][NTf2] is 168.15 K, indicating that the ILs possess very wide liquid temperature. But there are also the exceptions; for example, [C4C1im][OTf] does not exhibit a Tg at the heating scan. Furthermore, the thermogravimetric results show that most of the ILs have excellent thermal stability up to 400 °C. The decomposition temperature (Td) for the ILs follows the order 6640

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

plot shows an excellent linear correlation between ILs’ conductivities and viscosities with an average slope of R = 0.922 ± 0.022 and intercept log(C/S·cm2·mol−1) ranging from −0.02 to −0.14. An important usage of the Walden plot allows accurate calculation of ILs’ viscosity from their conductivity, which is usually far easier and quick with high accuracy. However, Schreiner et al.157 pointed out that the so-called ideal KCl line is at least a bit questionable with respect to their deviation from three considerations. First, there is no theoretical foundation justifying the comparison of pure ILs’ conductivity and viscosity to that of any diluted electrolyte solution; second, the choice of 1 M aqueous KCl as the “ideal” solution to compare to ILs’ data seems to be arbitrary; third, the assumption of ideal behavior for that has not been supported and arguable. Literature data for infinitely diluted aqueous KCl solutions at several temperatures show clearly that in this case the slope α in the Walden plot does not equal 1 but is actually about 0.87 (cf. Figure 2). For PILs, more delocalized charges and ion−ion interactions lead to their high conductivity. The alkylammonium-based PILs show high ionic conductivity; for example, the conductivity of [EtNH3][NO3] is 26.9 m·S cm−1. The conductivity for PILs with alkylammonium and alkylimidazolium cations decreased as the alkyl chain length increased.83 From the Walden plot of PILs, it is found that the PILs generally show poor ionicity.100 MacFarlane et al.158 measured the Walden plot of a series of ammonium acetate ILs and found that the primary amines exhibit a “‘closeto-ideal”’ behavior and produce highly ionized ILs, whereas tertiary amines sit far below the ideal line (vertical deviation, ΔW > 1) and tend to form mixtures with a low degree of proton transfer. The VTF equation is also chosen to analyze the ion conduction behavior of the acidic ILs. It is suggested that the temperature dependence of the ionic conductivity for [(HSO3)4N4][HSO4] is well-fitted by the VTF equation.138 In the Walden plots of the acidic IL, the molar conductivity, Λ (S· cm2·mol−1), was calculated by using the equation Λ = σMmρ−1 (where Mm, σ, and ρ are the molar mass, specific conductivity, and density) of the IL. Compared with the KCl ideal line, the Walden plots of the [(HSO3)4N4][HSO4] fall above the ideal line, suggesting that these SO3H-functionalized Brönsted acidic ILs belong to “superionic liquids”, which means the system obeys some decoupled transport mechanism, for example, a Grotthuss proton transport or small ion penetration mechanism.159 Surface Property. The surface is critical to understanding the IL processes, such as the gas−liquid surface for vaporization and gas-solubility,160 the liquid−liquid surface for separation and extraction,161 and the liquid−solid surface for the electrified interface in electrochemistry.73 Moreover, all chemical reactions are catalyzed at the interfaces between the participating species, and are intensified by transport between the interfaces of different phases. Surface tension (γ) as the basic property is a powerful means to explore the different types of arrangement/orientation of ions at a surface. Such measurements are also vital to anchor theoretical models that are able to describe realistically, at a molecular level, the fluid properties and structures of ILs. Law et al.162 have performed a systematic investigation of the surface tension of imidazolium-based ILs with a tensiometer. They found a linear temperature dependence of the surface tension. The surface energy decreases with increasing length of the alkyl chain. In the case of [CnC1im][NTf2] ILs in Table 1, a slow fall of surface tension γ is shown with n. At 298 K, the values of γ are in the range 35 to 42 mN m−1 for [C2C1im][NTf2], and decrease to

hydroxyl or methyl substitution on the alkyl chain. Similar trends were seen for the PILs with imidazolium cations, with the viscosity affected more by changes to the anion than to the cation.83 In addition, the VFT equation is used to fit the experimental data of viscosity to temperature for the acidic ILs as the estimated viscosity value of [(HSO3)4N4][HSO4] in Table 1.138 Conductivity (σ) is considered as another basic transport property. The high viscosity and low mobility lead to the low conductivity for the ILs in Table 1, which is typically in range of 0.1−30 m·S cm−1.153 The conductivity of ILs is also temperaturedependency and can be well fit by the VFT equation. For example the conductivity of [C2C1im][BF4] is ∼16 m·S cm−1 at 25 °C, ∼ 70 m·S cm-1 at 100 °C and ∼120 m·S cm-1 at 150 °C, which is increased with temperature. At low temperatures, the structures of cations and anions have an effect on the conductivity of ILs. The conductivity usually decreases with the increasing length of alkyl chains of cations, and the anions have more significant effect on the conductivity. For the ILs with [C4C1im]+ cation in Table 1, the conductivities follow the order of [NTf2]− > [BF4]− > [OTf]− > [CF3CO2]− > [PF6]−, which seems not relevant to the molecular volume and symmetry, but maybe more depend on the charge distribution of the anion and the number of charge carrier.154 It is important to take into account the conductivity when comparing the viscosity; such is called the Walden rule. Based on the rule, the Walden plot has been used increasingly in the last few years to illustrate the conductivity−viscosity relationship of ILs.100,134,155,156 Figure 2 shows the Walden plots of several typical imidazolium-based ILs. It can be seen that the Walden

Figure 2. Walden plot correlating conductivities (Y-axis) with viscosities (X-axis) for the imidazolium-based ILs. ●, [C2C1im][BF4]; ▽, [C2C1im][NTf2]; ■, [C2C1im][DCA]; ○, [C4C1im][BF4]. (+) represents conductivities of infinitely diluted aqueous KCl solutions from 0 to 55 °C and their linear extrapolation (dashed line). (×) represents Walden plot of aqueous KCl solutions at 25 °C in the concentration range 0 to 4 mol·L−1. The inset shows numerical values next to data points in the concentration (mol·L−1) for aqueous KCl solutions and temperature (°C) for infinitely diluted aqueous KCl solutions, respectively.157 Reproduced with permission from ref 157. Copyright 2010 American Chemical Society. 6641

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

values in the range 29 to 32 mN m−1 for [C10C1im][NTf2], which suggests that the trend is not linear when n is equal to or greater than six and the γ becomes almost constant.163,164 However, the larger γ occurs in PILs and the values of γ for ILs [EtNH3][NO3] and [EtNH3][HSO4] are 47.3 and 56.3 mN/m, respectively, which may be related with the particular structures of PILs. The structure at the surface is very different from the bulky one. Sloutskin et al.165 discovered a 6−7 Å thick surface layer by X-ray reflectivity (XRR), and the electron density at the surface layer is 10−12% higher than that of the unstructured liquid surface. Information about the ion orientation at the gas−liquid interface has been provided by the direct recoil spectrometry and MD calculation.166,167 We will discuss some details of structures at gas−liquid and solid−liquid surfaces in a later section. The γ of the acidic IL [(HSO3)4N4][HSO4] at T = 303.15 K is 46.79 mN·m−1 in Table 1, which is much higher than those of typical organic solvents. Generally speaking, γ of many liquid exhibits a nearly linear decrease with the increase of temperature and can be expressed by the Eötvos equation.140 Unlike water (where, as temperature increases, there are profound structural modifications due to the weakening of the H-bond network), the electrostatic force and VDW interactions are responsible for most of the internal cohesive energy of ILs, which indicated that the organization level of the ILs structure is high.

Figure 3. Harmonic wavenumbers (in cm−1) for two equilibrium structures of [C2C1im]Cl, denoted by geometry 1 and geometry 2, calculated by pHF-MP2 and DFT.95 The main components of the vibrational motions have been assigned and indicated in the figure. Reproduced with permission from ref 95. Copyright 2010 Elsevier.

one can aim to extend the DFT calculations to larger systems. In the classical MD simulation, each molecule carries its explicit chemical identity, and the electrostatic properties can be expressed by the multipole moments. Thus, classical force field parameters should represent the electric monopole and dipole properties found in DFT calculations which describe the electronic structure explicitly.178 Wang et al.179 performed a similar multiscale simulation to investigate the thermodynamics, including density, volume expansivity and isothermal compressibility, and transport properties, such as zero-shear viscosity of the trihexyltetradecylphosphonium chloride ([P6,6,6,14]Cl) IL, and the results are generally consistent with the experimental data. Although the current multiscale studies have focused on the formation of nanoscale structures and/or molecule-clusters, there are many opportunities to mix, match, and incorporate different IL atoms, functional groups, or ion types to optimize performance. Solvent properties and thus function are strongly influenced by the net organization of ions in the bulk, which requires deeper understanding of the role of IL self-assembly. Nevertheless, satisfying the requirement of a new industrial design, transport phenomena and hydrodynamics at unit scale have to be considered when the new IL solvents were used in reaction columns. The engineering continuation is a parallel development from the fundamental structural investigations.76,78−80,180 Therefore, the “multiscale” proposed in this review is a cross-domain research strategy covering different temporal and spatial scales to satisfy the requirement of a new paradigm of industrial design covering the entire design chain from molecule to system as shown in Figure 4. In this sense, the modeling and simulation are more advanced than the experimental measurement in exploring IL multiscale structures with different spatial scales. For IL systems, four spatial scales are mainly involved in Figure 4. Molecular scale ions, ion pairs, and their H-bonds can be calculated by the sophisticated QM and MM, and the mesoscopic self-assembled aggregates and interfaces can be explored by force field-based MD, CG, and DPD simulation. Many empirical models based on the traditional flow models were also developed to simulate the unit scale flow and mass and heat transfers in IL media. Process optimizing methods that assess the energy and environment of a process are prior to IL industrialization. Some details of the models and methods will be discussed in the next section.

2.3. Multiscale Strategy for ILs

Along with their complex structures and diverse properties, the ILs become rather complex systems. The interactions between anions and cations range from weak VDW, solvophobic,168,169 dispersion forces170 to strong coulomb and anisotropic Hbonding,171 halogen bond,172 dipole−dipole,173 magnetic dipole, electron pair donor/acceptor interactions, which leads to their obviously different properties.75 Additionally, structural changes may cause the dramatic changes in properties; for example, the increasing length of alkyl chains of the cations leads to remarkable decreases in the melting points.141 Sometimes, atomic structures are difficult to correlate with the macroproperties when the alkyl chains become longer or anions become more complicated.174−176 In the bulk phase, the ILs with long alkyl chains usually show self-assembled nanostructures,75 which causes a gap for directly correlating the macroscopic properties with the local structures. ILs structurally show multiscale organization resulting from single ions, ion pairs, and aggregates. In order to discover the structures of the ILs, a state-of-the-art multiscale model has been proposed by Dommert et al.95 The approach involved three levels of calculating techniques acting on different scales, namely quantum chemical (QC) post-Hartree− Fock calculations (pHF), density functional theory (DFT), and molecular dynamics (MD) simulations to go from individual molecules to bulk properties in a setup that is able to bridge the nanosecond time and nanometer length scale gaps. The pHF calculations were performed using the second order Møller− Plesset perturbation theory (MP2) with extended basis sets, coupled cluster theory, and basis set extrapolation techniques.177 For two [C2C1im]Cl equilibrium structures (geometry 1 and 2), pHF second order Møller−Plesset perturbation theory (MP2) and DFT agree reasonably well regarding structure, vibrational frequencies, and dipole moments (Figure 3). The total dipole moments from the DFT and the MP2 theory agreed within 5% for the single cation and one ion pair. But pHF approaches are computationally too expensive to calculate clusters with several ions; the principal agreement was taken as a good indication that 6642

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

Figure 4. Multiscale strategy of ILs satisfying the requirement of a new paradigm of industrial design covering the whole design chain from molecule to process. Molecular scale represents the single ions or ion pairs with 1−10 nm scale, which are modeled by QM and MM methods. Mesoscale represents the self-assembled nano- and mircoaggregates in and their solutions with larger 10 nm scale, which are modeled by CG and DPD methods. Unit scale can be involved in the flow and transportation phenomena of ILs in reaction columns with millimeter and meter scales. Process scale is a process integration from thermodynamic property prediction, to environmental and life assessment analysis, to objective optimization and process simulation.

3. MULTISCALE MODELING METHODS FOR IL SYSTEMS In chemical reactions, the structures are usually complex and vary from the atomistic level of the individual backbone bond of a single chain to the scale of the radius of gyration, which can reach tens of nanometers, even in the reactors of melts, blends, and solutions. Catalyst particles reach microns, millimeters, and even larger. The time span is an even wider range, from femtoseconds to milliseconds, seconds, or even hours for glassy materials, or large scale ordering processes. An approximate solution is the multiscale modelings by bridging of length and time scales and linking.181−183 In the past decades, advances in the multiscale modeling for predicting the structural, thermal, and mechanical properties at the microscopic level have been made.184,185 The models in general refer to the subatomic quantum mechanics (QM) at less than 1 nm scale and the molecular mechanics (MM) at 1−10 nm scales. The multiscale design stemmed from Karplus, Warshel, and Levitt’s works to deal with the biological macromolecular systems,186−188 in which the large systems are divided into small regions where the QM description is essential (bonds broking) to the remainder of the system being represented on a simpler level by MM, that is the QM/MM model. However, the complexity of a chemical process is far beyond the description from a biomolecular system; for example, transport properties in a reactor unit have to be considered in detail. Fermeglia et al.184 proposed the message passing model to bridge the atomistic scale and macroscopic scale passing through the mesoscopic scale. The strategy is based on an overlapping array of successively coarser modeling techniques. At each plateau, the coarse description parameters are based on representative results of the immediately smaller (finer) description. One obtains computational information at a smaller (finer) scale and passes it to a larger model by relaying some parameters.189,190 The ultimate goal is to

predict the macroscopic behavior of an engineering process from the microscopic structures. Although the methods have been developed extensively, we do not aim to provide a detailed description for each method, and only their basic principles are discussed that are used to simulate ILs in this review. With respect to some details of the methods, interested readers can refer to relevant books, reviews, and research articles. We would like to classify these methods according to their length scales from molecular, to mesoscale, to macroscopic prediction. 3.1. Molecular Scale Methods

The methods in general refer to the subatomic QM and the atomistic MM. Even now, it is still computationally difficult to use high-level QM to obtain convergent sampling of the many configurations needed to reliably describe the free energies of even medium sized systems. Even the DFT191−193 can only deal with molecular size less than 1−10 nm. The solution to this problem has emerged when realizing that a description of complex systems does not require the representation of all parts of the system at the same level. The QM/MM model is coupling between the electrostatic effects of the quantum and classical subsystems, and it has eventually become a key to advances in describing the complex system.188 The general quantum self-consistent Hamiltonian is expressed by ii Fμμ = Uμμ +

1 Pμμγμμ − 2

∑ Pννγνμ − ∑ Q iγii′ − ∑ Q jγij ν≠μ

i ′≠ i

j

(1)

where U is the core Hamiltonian, P is the quantum mechanical bond order, Q is the net atomic charge, γ is the electronic repulsive integral, and μ and ν are atomic orbitals on atom i. Now, assigning atom i to the part of the system that should be treated quantum mechanically indicated that the other atoms (denoted 6643

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

by j) can be treated classically, assuming that their charge is constant. This equation can be generalized to cases where the charge distribution of the classical atoms is not fixed and can be polarized by the field of the quantum atoms. Thus, the solute− solvent coupling Hamiltonian is obtained by adding the potential from the solvent atoms to the solute Hamiltonian. The total potential energy is then given by Vtotal = ES(FS) + ESs′ + Ess

Fi ⃗(t ) = m

∑ bond

ki (li − li0)2 + 2

(2)



+

torsion N

∑ i=1

∑ angle

dt 2

(4)

There are many algorithms to integrate the above equation, including the algorithms of verlet, velocity verlet, leapfrog, and Beeman. The verlet algorithm is widely used, in which the position r(⃗ t) and acceleration a⃗(t) at time t and the position r(⃗ t − δt) from the previous step t − δt are used to calculate the new position r(⃗ t + δt) at t + δt in eq 5 r (⃗ t + δt ) = 2 r (⃗ t ) − r (⃗ t − δt ) + a ⃗(t )δt 2

ES (FS) is the energy that is quantum mechanically obtained that includes the given electrostatic potential from the solvent. E′Ss is the nonelectrostatic solute−solvent interaction term, and Ess is the solvent−solvent classical force field. At this level of approximation, the nonelectrostatic term is evaluated by the standard van der Waals potential function. For some supermolecules, the “connection” between the quantum and classical regions is treated by a classical force field (which is included in E′Ss), where the quantum atoms at the boundaries are connected to dummy hydrogen-like atoms in order to balance the electrons in the quantum system. The introduction of a realistic electrostatic model for the large molecule, such as an enzyme and its surrounding water, together with the incorporation of this effect in a quantum Hamiltonian, finally for the first time yielded the energy of heterolytic bond breaking processes. This is the QM/MM approach to treat large systems quantum mechanically combined with the MM method.188 MM generally depends on the force fields. Equation 3 represents schematically a frequently used force field form. The former three terms are bonded interactions including bond, angle, and torsion contributions, and the latter two terms are nonbonded interactions including electrostatic and LennardJones interactions. U (r N ) =

d 2 ri ⃗

(5)

MD simulation can be performed in many ensembles, such as microcanonical (NVE), canonical (NVT), and isothermal− isobarric (NPT). The temperature and pressure are controlled by adding an appropriate thermostat (e.g., Berendsen, Nose, Nose− Hoover, and Nose−Poincare) and barostat (e. g., Andersen, Hoover, and Berendsen), respectively. 3.2. Mesoscale Methods

The length range 10−100 nm or larger is referred to as the mesoscale, and some are the morphology of polymer blend and composite, which often dominates actual aggregate properties or transfer properties. The methods, including coarse-grained (CG), dissipative particle dynamics (DPD),199 lattice Boltzmann (LB), 2 0 0 time-dependent Ginzburg−Landau method (TDGL),201 and dynamic density functional theory,202,203 have been developed.204 3.2.1. CG Method. In order to bridge the length scales by simple models, Warshel et al.205 proposed the CG model to solve protein folding. An improved treatment of electrostatic energies and long time-scale behavior of the simplified model is solved by applying MD simulations of the full model and the same forces in the reduced model by Langevin dynamics.206 The CG approach is determined by the effective CG energy function UCG between CG sites by finding the best fit of all total forces on the underlying atomic groups, under the central force pairwise assumption for the CG forces. The primary challenge is significantly easier to simulate, but it reproduces the same physical behavior as a reference all-atom one. For this, three choices need to be made: (1) CG models typically entail pseudoatom sites that are designed to represent combined groups of multiple atoms. (2) UCG for the CG model defines the interactions between the pseudoatoms and is chosen so as to reproduce the same thermodynamic properties of the reference detailed system. Some optimal UCG values have been proposed,207−210 such as the following equation optimized by relative entropy

ki (θi − θi0) 2

Vn (1 + cos(nω − γ ))+ 2

⎛ ⎡⎛ ⎞12 ⎛ ⎞6 ⎤ qiqj ⎞⎟ σij ⎥ ⎜ ⎢⎜ σij ⎟ ⎜ ⎟ ∑ ⎜4εij⎢⎜ ⎟ − ⎜ ⎟ ⎥ + ⎟ ⎜ r ⎝ rij ⎠ ⎦ 4πε0rij ⎟⎠ j=i+1 ⎝ ⎣⎝ ij ⎠ N

(3)

Different force field parameters can be derived from dissimilar types of experimental data, such as enthalpy of vaporization (OPLS), enthalpy of sublimation, dipole moments, or various spectroscopic parameters. Some functions do not account for electronic polarization of the environment, an effect that can significantly reduce electrostatic interactions of partial atomic charges. This problem was addressed by developing “polarizable force fields”.194,195 CFF/ind and ENZYMIX are the first polarizable force fields that have subsequently been used in many applications for biological systems.195 Some new force fields were developed, such as DRF90, developed by Duijnen et al., PIPF developed by Gao et al.,196 the CHARMM polarizable force field developed by Patel and Brooks III,197 and the new generation AMBER force field developed by Caldwell et al.198 Based on the force fields, MD simulation is developed by solving a set of classic Newtonian equations (eq 4) of motion for all atoms.

UCG(R ) = −kBT ln

∫ e−βU δ[M(r) − R]dr + Const A

(6)

(3) Affective dynamical equations for describing the timeevolution of the CG system. When we remove some degrees of freedom during CG, we remove the time scales associated with them; to properly evolve the CG system, we need to build those time scales back into the dynamical equations to get the dynamics of the CG system (approximately) correct. Voth et al. have applied the multiscale CG model to liquid,211 biomolecular,212,213 and IL97,214−216 systems. The CG of ILs involves a reduction in the number of interaction sites inside of 6644

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

partial charge of the CG sites of kind α; Nβ and K are, respectively, the number of the CG sites of kind β and the total number of kinds of CG sites in the system. Involving the bonded CG interactions, the FM method permits separation of bonded and nonbonded CG forces, which is important for the sites having an overlap in intra- and intermolecular regions, causing the FM force field to be a mixture of bonded and nonbonded components. All bonds might be described by simple harmonic potentials. In principle, the harmonic bonded potentials can then be obtained through least-squares fit of the FM bonded forces to a harmonic force. However, the CG bonded forces for ILs often deviate from being harmonic, and some further generalizations are required. A systematic and accurate way of retrieving effective bonded force field parameters from the saved MD configurations is described in Voth et al.’s work.215 3.2.2. DPD Method. DPD was originally developed by Hoogerbrugge and Koelman.199 It can be used to simulate both Newtonian and non-Newtonian fluids, including polymer melts and blends, on microscopic length and time scales. DPD particles are defined by their mass Mi, position ri, and momentum pi. The force between particles i and j can be described by a sum of conservative FCij , dissipative FDij , and random forces FRij :

each ion. The center of a CG site can be taken to be the center-ofmass (COM) of the underlying atoms. A group of atoms can be most readily CG as one site if there are relatively few degrees of freedom inside that group. As a principle, an aromatic ring can be considered as one CG site because its backbone structure is almost a rigid body and a methyl group can also be one CG site when it rotates along its principal axis symmetrically. The anion of ILs is also generally small and amenable as a single CG site. By contrast, the bulky cation generally has to be grouped into several CG sites. The bonded and nonbonded interactions were mixed together, with the bond force field parameters being decided according to the matched force, and the angles were treated as harmonic bonds. The nonbonded CG sites are handled by a force-matching (FM) procedure.217−219 In the FM methodology, the atoms are first grouped into CG sites with the CG force on a CG site α being the net force on all underlying atoms, FαCG = ∑i f αi , where fiα is the force on the ith underlying atom. The pairwise effective force f p(rij) between CG sites i and j is represented as the sum of a short-range part and a long-range Coulomb part, i. e., ⎛ Q iQ j ⎞ ⎟nij f p (rij) = −⎜⎜f (rij) + rij2 ⎟⎠ ⎝

(7)

where rij is the modulus of the vector connecting the two CG sites and Q is the partial charge of the CG site. The short-range term f(rij) is represented by third-order polynomials (cubic splines) connecting a set of points {rk}, thus preserving continuity of its functions and their first two derivatives across the junction such that

= A(r , {rk})fi + B(r , {rk})fi + 1 + C(r , {rk})f i″

∑ ∑ ∑ γ = nb , b β = 1, k j = 1, Nβ

+

Q αβ

(11)

FijD = −γωD(rij)(eiĵ ·pij )eiĵ

(12)

FijR = σξijωR (rij)eiĵ

(13)

∑ FijC + ∑ FijD + ∑ FijR j≠i

j≠i

j≠i

(14)

Because the forces are pairwise and momentum is conserved, however, energy is not conserved because of the presence of the dissipative and random force terms similar to those of BD, but it incorporates the effects of Brownian motion on larger length scales. DPD has some advantages over MD; mainly, the hydrodynamic behavior is observed with far fewer particles than required in a MD simulation because of its larger particle size. 3.2.3. MM&Meso Bridging. Some mesoscale methods, such as DPD, utilize the CG model to simulate the morphology of inhomogeneous materials or solutions using Flory−Huggins χ parameters between components, and the χ parameters are estimated through the atomistic simulation.202,220 Fermeglia et al. estimated χ parameters from the bulk atomistic simulations of blends and single component systems via the relevant cohesive energies as shown in eq 12.184 The cohesive energies are calculated for pure components 1 and 2, and their mixture, using a MM technique and then combined to obtain the interaction for the mixture.

⎛ ⎜f {rαil , βjl}, {rαβ , γ , k}, {f ″ ,γ ,k ) }, (f αβ ⎜ αβ , γ , k ⎝

⎞ ⎟

δγ , nb⎟nαil , βjl rα2il , βjl ⎠

= FαCG il

FijC = Π 0ωC (rij)eiĵ

Fi(t ) =

(8)

where r ∈ [ri,ri+1], A, B, C, and D are known functions of r, and {rk}, {f k}, and {f ″k } are tabulations of f(r) and its second derivatives on a radial mesh grid {rk}. A key property for the success of this method is that a spline representation depends linearly on its parameters, which are tabulations of f(r) and its second derivative, {f k, f ″k }, on a radial mesh {rk}. The parameters {f k and f ″k } are obtained from the fit. Equalization of the atomistic MD forces FilCG to forces predicted by using the representation in eqs 7 and 8, which act on the ith CG site in the ith configuration sampled along the atomistic trajectories, results in the following set of linear equations −

(10)

Here rij = |ri − rj|, êij = r̂ij/rij, Π0 is a constant related to the fluid compressibility, and ωC, ωD, and ωR are the weight functions for each force. The total force acting on particle i at time t is

f (r , {rk}, {fk }, {f k″ })

+ D(r , {rk})f i″+ 1

Fij = FijC + FijD + FijR

(9)

⎛ ΔE ⎞ χ = ⎜ mix ⎟Vmon ⎝ RT ⎠

with respect to the force parameters subject to the fit {fαβ,γ,k, f″αβ,γ,k, Qαβ}, where α = 1, ..., Nα. In eq 9, {αil} labels the ith atom of kind α in the lth configuration; rαil,βil is the distance between atoms {αi} and {βj} in the lth configuration; Qαβ = Qα − Qβ with Qα the

(15)

where Vnom is the volume of the polymer or cluster. The energy of mixing is defined in terms of the cohesive energy densities for 6645

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

pure components and their blends calculated from the cohesive energy per volume; ΔEmix is an energy per unit volume, and the cohesion energy is related to the solubility parameter δ (δ = Ecoh /V ):

ρ=

⎧+|d| x ∈ the primary phase ⎪ φ (x , t ) = ⎨ 0 x ∈ the interface ⎪ −| | d x ∈ the second phase ⎩

(17)

∂ (φ) + ∇·(φv ⃗) = 0 ∂t

For all kinds of fluid flows, the conservation equations for mass and momentum are solved with the computational fluid dynamics (CFD) method, which can be written as

(25)

The surface tension source term can be expressed as FSf⃗ = σκδ(φ)∇φ

(18)

∂ (ρv ⃗) + ∇(ρvv⃗ ⃗) ∂t

(26)

Compared with the VOF method, the level-set function is smooth and continuous; thus, its spatial gradients and the topologically complex interfaces can be accurately calculated. However, the LS method is found to have a deficiency in preserving volume conservation. Therefore, in order to overcome the deficiencies, the coupled LS-VOF (CLSVOF) method224,225 is generally applied in recent works.

(19)

where v⃗ is the velocity vector and P is the pressure. F⃗S stands for the source term. The finite volume method (FVM) is a common approach to solve these equations, due to its advantage in memory usage and solution speed, especially for large systems, high Reynolds number turbulent flows and source term dominated flows. In FVM, the above partial differential equations are recast in a conservative form, and solved over discrete control volumes, which guarantees the conservation of fluxes through a particular control volume. Multiphase flow plays an important role in chemical engineering; in particular, the gas−liquid flow is common in ILs processes. For the interfacial multiphase flow, the volume of fluid (VOF) method and level set (LS) method are the most popular CFD methods. The VOF method221 is a flexible and efficient approach to the fluid−fluid interface tracking. In the VOF method, the volume fraction equation is solved, which is represented as

3.4. Process Scale Methods

Process simulation is a term in chemical engineering which emphasizes the unity of a chemical process and considers the interactions between different unit operations from the outset, rather than optimizing them separately. For the IL process, the models mainly integrate property prediction, environmental assessment, and objective optimization. 3.4.1. COSMO-RS Model. The conductor-like screening model for realistic solvents (COSMO-RS),226 developed by Klamt et al.227 in 1993, is proven as an effective method for predicting the thermodynamic properties of ILs.228 In this model, each molecule is described by the surface charge densities and the interaction between molecules is characterized by the surface charge densities (σ), as shown in eq 24 and eq 25. Emisfit (σ , σ ′) = αeff emisfit (σ , σ ′) = αeff

(20)

The volume fraction of each phase obeys the following expression:

α′ (σ + σ ′) 2

(27)

Ehb(σ , σ ′) = αeff ehb(σ , σ ′) = αeff chb min{0, σσ ′ + σhb 2} (28)

∑ αq = 1 q

(24)

where d is the distance from the interface. The evolution of the level-set function can be given in a similar fashion as the VOF model:

3.3. CFD Model for Unit Scale

∂ (αq) + ∇·(αqv ⃗) ∂t

(23)

where σ is the surface tension and κ is the curvature, which is ∇α defined as κ = ∇· |∇ α| . In the LS method,223 the level-set function φ is defined as a signed distance to the interface. Accordingly, the interface is the zero level-set, and φ(x,t) can be expressed in a two-phase flow system:

where Z is the coordination number and ΔE12 is the interaction energy between two systems.

= −∇P + ∇(μ(∇·v ⃗ + ∇·v ⃗T )) + ρg ⃗ + FS⃗

(22)

q

FS⃗ = σκ ∇α

(16)

∇·(ρv ⃗) = 0

∑ αqμq

The interaction of surface tension is introduced into the momentum equation by the source term F⃗S. For instance, the continuum surface force (CSF) model222 is usually employed to define the surface tension source term, which can be written as

in which Φ is the volume fraction and Ecoh is the cohesive energy. Accordingly, the values of Ecoh at different temperatures can be obtained, for both mixtures and pure components, from simulation by calculating the difference between the nonbonded energy of the periodic structure and the corresponding value of an isolated parent chain in a vacuum. An alternative method, based on a less precise procedure, is to obtain an estimation of χ using advanced quantitative structure− property relationships (QSPRs). The value of χ can also be determined by MC simulation using the following equation:

Z ΔE12 RT

μ=

q

⎛E ⎞ ⎛E ⎞ ⎛E ⎞ ΔEmix = Φ1⎜ coh ⎟ + Φ2⎜ coh ⎟ − ⎜ coh ⎟ ⎝ V ⎠ pure1 ⎝ V ⎠ pure2 ⎝ V ⎠mix

χ=

∑ αqρq

where αef f is the contacted surface area; α′ is the energy factor; σ and σ′ are the net screening charge density; and Emisf it is the misfit energy of two molecules. Ehb is the hydrogen interaction energy. emisf it and emisf it are the misfit energy and hydrogen interaction energy of the unit area.

(21)

where αq is the volume fraction of phase q. The fluid density ρ and viscosity μ in the above equations are defined as follows 6646

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

As a result, the chemical potential of each fragment μδ(σ) can be calculated as given in eq 26.

C2−H moiety is more acidic than that on C4/5−H. However, ab initio MD simulation232 indicated that the hydrogen atoms on C2−H and C4/5−H carry very similar charges. The contradictory results can be explained when the charges of the hydrogen atom are summed into the heavy atoms; the combined charges on the C−H moiety such as the charges q(C2−H) = +0.536, q(C4−H) = +0.239, and q(C5−H) = +0.242 in [C4C1im]+ the cation show the same order as the experimental acidities.233 For Nalkylpyridinium,234,235 tetraalkylammonium,127,236 and tetraalkylphosphonium237,238 cations, the nitrogen and phosphorus atoms are positively charged, and α-carbon atoms are negatively charged. The α-hydrogen atoms show some acidity and may be important for some properties and reactions of these ILs.237−239 Torsional motions of the alkyl groups can give rise to conformational equilibria in some cations. For example, trans− trans and trans−gauche-conformations of the n-butyl chain in [C4C1im]+ result in different crystalline polymorphs from orthorhombic to monoclinic.231,240 It is also likely that the coexistence of these conformations significantly affects the liquid structure and the bulk properties of ILs, even resulting in nanostructured fluids.241 Anions such as [BF4]− or [PF6]− have been a focus, but in the presence of moisture, these anions hydrolyze, forming HF.242 More complex anions, such as [NTf2]− or [OTf]−, or halogenfree ions, such as dicyanimide, tosylate, or N-alkylsulfates, are now preferred.243 In particular, the [NTf2]− forms liquid salts of low viscosity with high thermal and electrochemical stability.148,153,242 The torsional motions of the anions also have an effect on the properties of ILs. For example, in the liquid form of [C1C1im][NTf2] the anion trans form prevails, whereas in the crystal structure the anion cis conformer is found.244−246 The different steric and electronic environments of the two conformations are the major reason for the bulk properties of these salts. 4.1.2. Structures, Interaction Energies, and Charge Distributions of Ion Pairs. Ion pairs are the simplest repeat unit in ILs and thus are important references in predicting the properties of ILs. QM calculations have depicted the static structures of the ion pairs.174,247−250 In the [C4C1im]Cl ion pair, it was found that the Cl− anion has a favorable position above or below the imidazolium ring; however, other favorable positions are within the plane of the ring, particularly in front of the C2−H group, and form a C2−H···Cl H-bond. The favorable positions illustrate the competition between different interactions. The above or below position is determined by Coulombic force, and the plane position is the H-bonded arrangement. The interaction energy (ΔE) between anion and cation is a measurement of the stability of the ion pair. The ΔE values of the [C2C1im][CF3CO2], [C2C1im][BF4], [C2C1im][CF3SO3], [C2C1im][NTf2], and [C2C1im][PF6] ion pairs are 375.72, 356.48, 345.60, 329.70, and 328.03 kJ/mol, calculated at the MP2/6-311* level, respectively.251 It seems that the energies are not closely relevant to the sizes and symmetries of these anions. Contrarily, the more negative charge on the −CO2 group than that on the −SO3 group may result from the larger ΔE of [C2C1im][CF3CO2] than that of [C2C1im][CF3SO3], indicating the charge distribution is more responsible for the ΔE. It is also found that the conformation with the C2−H···anion H-bond in the ion pairs has the larger ΔE values, such as the C2−H···O Hbonds in [C2C1im][CF3CO2] and [C2C1im][CF3SO3] ion pairs. According to the charge scaling effect, the total charge of a single ion may be reduced in the ion pairs.252,253 Krossing et al.254 found that the ionic charges (q±) in 1,3-dimethylimidazolium

⎡ ⎧ μ (σ ′) − αeff e(σ , σ ′) ⎫ δ ⎬ μδ (σ ) = −RT × ln⎢ pδ (σ ′) exp⎨ RT ⎢⎣ ⎭ ⎩











⎤ dσ ′⎥ ⎥⎦

(29)

where pδ(σ′) is the surface screening charge density distribution, also called the σ profile. e(σ,σ′) is the total molecular interaction energy, which is the sum of Emisf it and Ehb. R is the gas constant, and T is absolute temperature. The activity coefficient γi for the residual term of the solute can be described as in eq 27 ⎛ μ − μ0⎞ i ⎟ γi = exp⎜⎜ i ⎟ RT ⎝ ⎠

(30)

where μi is the chemical potential and μ0i is the chemical potential under the standard state. 3.4.2. Green Degree Method. ILs are honored as green solvents, and environmental assessment is critical to their industrial applications. Recently, Zhang et al.229 proposed the green degree (GD) method to quantitatively evaluate the comprehensive environmental impact of a substance, stream, or a process based on nine environmental impact categories. Based on this method, every process can be considered as a system, which would influence the environment through the exchange of material and energy with the environment. ΔGD =

∑ GDks ,out + ∑ GDke ,out + ∑ GDks ,in k

+

k



GDke , in

k

k

(31)

e,in where GDs.in k and GDk are the input GD values of a material e,out stream and an energy stream into a process. GDs,out k and GDk are the corresponding output values.

4. MOLECULAR SCALE STRUCTURES OF ILS ILs were intuitively thought to be same analogues as classic molten salts. However, the low Tm of most ILs defies this consideration. Molecular scale structures of ILs have been explored by QM, MM, and QM/MM methods.95 A crucial application in QM/MM is treating the electrostatic potentials through the Blöchl approach employed on the DFT level to obtain values for partial charges, which are considerably different from the ones used in the standard force fields. 4.1. Single Ions and Ion Pairs

4.1.1. Structures of Single Ions. Single ions are the threshold to understand the structures of ILs. For 1-alkyl-3methylimidazolium cation, a main feature is π-electron delocalization on the imidazolium ring to lead to a 3-center-4electron configuration across the N1−C2-N3 moiety and a double bond between the C4 and C5 atoms.174 The negative charges are located on the nitrogen atoms. Owing to the electron deficit in the CN bond, the C2 atom is positively charged and the C4/5 carbon atoms remain essentially neutral. The acidity of the hydrogen atoms on the imidazolium ring is related with some chemical reactions and particularly interesting.151,174,175 Some experimental results230,231 suggest that the hydrogen atom on the 6647

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

NMR,257,265,279 X-ray,279−281 neutron scattering,174,282 and calculations248,273,274,283−285) has confirmed the assertion. Although the H-bonds between a neutral molecule and a charged ion, such as NH4+···NH3 and NH4+···H2O, have been discussed,286 the H-bonds between two charged ions, such as in ILs, are not the same case. More recently, a critical review by Hunt et al.175 has highlighted the H-bonds in ILs. QM parameters, including interaction energies, partial charges, and density descriptors encompassed by the AIM methodology, qualitative molecular orbital theory, and NBO analysis, were externalized. To differentiate from the classic H-bond and the ionic H-bond,286 the H-bond in ILs was called a doubly ionic Hbond by the authors. Dong et al.287 in their recent paper gave a new “Z-bond” definition to describe the coupling interaction between the H-bond and the electrostatic force in ILs, which shows some features of the geometric, energetic, electronic, and self-assembly aspects. For these imidazolium-based ILs, multiple H-bonds can form within an ion pair and exhibit bifurcated and chelated structures,175 as shown in Figure 5. It can be seen that there

methylsulfate ([C1C1im][C1SO4]) are reduced to 0.90e by high resolution XRD measurement and to 0.87e by atom in molecule (AIM) calculation, respectively, which suggests the charge transferring between cation and anion. In the [C1C1im]Cl ion pair, the charge of −0.82 upon extrapolation to infinite cluster size was reported for the chloride anion. Schmidt et al.177 calculated an even smaller net charge of q± = 0.63e for the ions of [C1C1im]Cl in the bulk phase. MD simulations gave the charged scaling of 0.9 at the IL−water interfaces, while the scaling by ε−0.5 was used in high- and low-dielectric media, ion−solvent systems, and ionic groups in proteins, where ε is the electronic component or high-frequency limit of the relative permittivity.255−258 4.1.3. Ion Pair Lifetime. The stabilities of an ion pair in bulk ILs are questionable.259 Dielectric spectroscopy260 revealed that the [EtNH3][NO3] is a mixture of ion pairs and free ions, in which only 8% are contact ion pairs and 92% are free ions. The lifetime of the ion pair was found to be on the order of hundreds of picoseconds at 298 K. For the ILs [C2C1im][X] (X = Cl−, Br−, and I−), the contact ion pairs were observed by multinuclear NMR, and in a quasi-molecular state stabilized by strong Hbonds.261 However, ion pair formation has not been observed in the pyrrolidinium-, pyridinium-, tetraalkylammonium-, and triethylsulfonium-based ILs by recent dielectric spectroscopies sensitive to picosecond-to-nanosecond dynamics.262−264 Corresponding NMR does not detect the ion pairs in a micro-tomilliseconds time scale.265 Thus, any potential ion pair maybe be the shorter lifetimes less than a few picoseconds. MD simulations showed proofs for the short lifetime ion pairs. Lee et al.266 found that free ions outnumber ion pairs by 2:1, and pairs are stable only for times of 60−140 ps in the [C4C1im][PF6]. The cation−anion interaction was best described using the concept of an ion association rather than an ion pair, as each ion is not coupled to a single counterion in the ionic atmosphere. Lynden-Bell267 found that the forming anion + cation ion pairs were not important for describing the bulk IL structure. Structural analysis suggests that the destabilization of an ion pair has to do with the overscreening of electrostatic charge in the first solvation shell. Many reports have suggested that ion pairs can be stable when ILs are dissolved in molecular solvents at low concentration. Longer-lived ion pairs were found in [C2C1im][X] ILs + propionitrile and [C4C1im][BF4] + DMSO mixtures by NMR measurements.261,268 MD simulations of [C4C1im][PF6] with naphthalene and other aromatics have also obtained anion + cation ion pairs and cation + aromatic units across a range of solute concentrations.269

Figure 5. H-bonding distances of the typical AIL ion pairs (the distances in angstroms are noted by dashed lines): (A) two C−H···Cl H-bonds in the [C2C1im]Cl ion pair; (B) four C−H···F H-bonds in the [C2C1im][BF4] ion pair; (C) four C−H···O H-bonds in the [C2C1im][CH3CO2] ion pair; (D) C−H···N and C−H···O H-bonds in the [C2C1im][NTf2] ion pair.

4.2. H-Bond and Network

are two C−H···Cl H-bonds in [C2C1im]Cl (Figure 5A) and four C−H···F H-bonds in [C2C1im][BF4] (Figure 2B). For the anions N[CN]2−, [OAc]−, [CF3SO3]−, and [HSO4]−, etc. with the stronger H-bonding acceptors, the anion can “occupy” the multiple positions around the same cation to form multiple but distant (not adjacent) H-bonds. For example, in [C2C1im][OAc] in Figure 5C, the −COO group orients the front position of the cation to form four C−H···O H-bonds. For the complex anions, including linker groups/atoms and internal torsional motions, e.g. [NTf2]− with more than one type of H-bonding acceptor (e.g., O and N of [NTf2]−), there are different types of H-bonds in Figure 5D, which have both an electron-rich central N atom and pendant groups, containing O and N atoms, respectively. Normally, the O atoms are assumed to dominate the H-bonding interactions, but all the electronegative atoms are available in forming the H-bonds.150

Noncovalent interactions play a key role in determining the properties of solvents.75 In the ILs, Coulomb force, H-bonding, and VDW interactions are of great importance for the properties and behaviors of ILs. However, of them, H-bonding interactions are controversial in these Coulomb fluids. The questions are opened: Do H-bonds exist in these Coulomb fluids? To what extent do H-bonds contribute to the overall interaction? And are H-bonds important for the properties of ILs and reactions in ILs? Some researchers, typically Dong, 150,151 Ludwig, 270−273 Hunt,174,175,274 and Izgorodina249 et al., have performed many investigations of H-bonds in ILs. 4.2.1. H-Bond Structures. Nowadays, it is widely accepted that the H-bond is the important interaction existing in ILs, especially in the PILs. A large amount of evidence from experiments (FTIR, 272,275−278 Raman spectroscopy, 278 6648

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

Figure 6. IR spectrum and calculated vibrational modes of [C2C1im][BF4]. The calculated vibrational modes are of the corresponding ion clusters ([C2C1im][BF4])n with n = 2, 3, 4, and 5 at the B3LYP/6-31+G** level and are corrected by the factor 0.964−0.967.285 Reproduced with permission from ref 285. Copyright 2012 American Chemical Society.

Ludwig et al.276,280,288 estimated the H-bonding contribution to be 30 and 50 kJ mol−1, and thus representing between 10% and 16% of the overall interaction energy. The H-bonding strength of ILs can be correlated to the second stabilization energies (ΔE2) between anion and cation. In order to compare the ΔE2, a strategy is to replace the methyl group at the C2 position of [C4C1C1im] by a hydrogen atom obtaining [C4C1im] while the anion [BF4]− is equal. Now the H-bond via C2−H···A can occur in [C4C1im][BF4], and replacing the weakly interacting anion [BF4]− by the strongly interacting anion [NO3]− in the imidazolium-based AIL leads to a stronger H-bonding strength, and finally substituting the imidazolium cation [C4mim]+ by the propyl ammonium cation [PrAm] while retaining the anion [NO3]− forms the PIL [PrAm][NO3] with the strongest Hbonding strength. The calculated ΔE2 of the AILs [C4C1C1im][BF4], [C4C1im][BF4], [C4C1im][NO3], and PIL [PrAm][NO3] is increased gradually,289 and the H-bonding interaction of the PIL is stronger than that of the AIL generally. The liquid structures of the AIL [C1C1im][NTf2] and the PILs [C1Him][NTf2] and [H1H2im][NTf2] were revealed by the high-energy X-ray total scattering (HEXTS) experiments and MD simulations.290 It was found that the N−H···O H-bond causes the closest cation−anion orientation perturbation. MD simulations clearly exhibit that the anion O atom prefers the positively charged NH hydrogen or the imidazolium ring, hydrogen, while the F atom of the anion dominantly distributes right above and right below the imidazolium ring plane. The N− H···O hydrogen bond is short and linear, while the C2−H···O interaction is long and bent. Due to the [NTf2]− anion being quite a weak H-bonding acceptor and/or the electron pair donor, the short and linear N−H···O H-bonds imply that spatial availability, that is, the configuration entropy, could contribute to the formation of the H-bonds in PILs. Sahu et al.291 reported the rotational relaxation behavior of charged and uncharged solutes in the H-bonded [C2C1im][NTf2] IL and non-H-bonded C2methylated [C2C1C1im][NTf2] IL. The results indicated the obvious role the ionic component of the H-bonding interaction has on rotational relaxation of charged solutes. Dispersion forces in ILs are currently studied intensively and thought to become crucial for the structures of the ILs.170,292,293 The experimental results from far-IR spectroscopy and DFT calculations both support the explicit dispersion forces in both

PILs and AILs.294,295 In particular, it is shown that dispersion forces can outbalance the H-bond in ILs with increasing temperature. Ludwig et al. reported that the H-bonded species (+N−H···O3S−) are dominant at low temperatures, and are replaced by the dispersion-interacting species (+N(C6H13)3··· O3S−) with increasing temperature.270 The transition between the two ion-pair species requires about 34 kJ/mol, as obtained from a van’t Hoff analysis. 4.2.2. H-Bonding Network. Recent spectroscopic studies296 found that H-bonds can extend spatially to form the 3D Hbonding network when the imidazolium cations bearing a chiral center or small halide anion, which may stabilize individual ion conformation, even are important for IL properties.152 Dymek et al.297 characterized the hygroscopic crystal structure of [C 2 C 1 im]Cl by XRD and found that the cations are perpendicular to the c axis in four distinct layers and the anions are a layered one in the orthorhmic crystal. They found that each anion interacts with three cations and each cation interacts with three anions and the Cl− anions are situated in an H-bonded position rather than at random, characteristic of a C−H···Cl Hbond network by IR spectroscopy. Dong et al.248 gave a similar H-bonding network of [C 2 C 1 im]Cl calculated by QM calculation. Moreover, Dong et al.285 combined the FTIR spectrum with QM calculation to characterize the H-bonding network of [C2C1im][BF4] (Figure 6). The major vibrations show the good agreement between the experimental and the calculated results. The vibrational peaks between 3000 and 3500 cm−1 are assigned to the stretching vibration of C−H bonds of the imidazolium rings. The strongest peak around 1136 cm−1 is the stretching vibration of the B−F bonds of the anion. Interesting is the far-IR range244,262 in the inset in Figure 6, and the frequencies around 122 cm−1 can be attributed to the asymmetric and symmetric stretching vibrations of the CH···F H-bond, and those around 86 cm−1 can be assigned to the corresponding bending modes of the H-bond. MD simulation found that the anion prefers the positions close to the hydrogen atoms of the imidazolium ring and the H-bonds can be extended into the 3-dimensional network in the bulk [C2C1im][BF4] IL phase.298,299 The relative strength of the Hbonds can be estimated by ESI-MS experiments by dissociating the mixed supramolecular clusters. The relative order of Hbonding strengths to the [C4C1im]+ cation is [CF3CO2]− > 6649

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

Figure 7. (A) Tm, η, and ΔvapH values for the imidazolium-based ILs with a series of anions. [C2C1C1im]+ is C2 methylated cations.306 (A) (a) Tm of [C2C1im]+ (squares) and [C2C1C1im]+ (circles). (b) η of [C2C1im]+ (squares) and [C2C1C1im]+ (circles) with [NTf2]− and [OTf]− anions, and η for methylated ILs at the C5 position (diamonds). (c) ΔvapH of [C4C1im]+ (squares) and [C3C1C1im]+ (circles) with anions [NTf2]− and [OTf]−. (Reproduced with permission from ref 306. Copyright 2009 Wiley-VCH Verlag GmbH & Co. KGaA, Weinheim.) (B) Measured far-IR spectra for 11 imidazolium-based ILs with a series of anions. (Reproduced with permission from ref 311. Copyright 2010 Wiley-VCH Verlag GmbH & Co. KGaA, Weinheim).

[BF4]− > [PF6]− > [InCl4]− > [BPh4]− for the five anions. Neutron diffraction studies also found that charge ordering can endure in the liquid state.244,282 For the PIL [EtNH3][NO3], it was found that the phase transfer of rare gases and hydrocarbons from cyclohexane to the IL was accompanied by negative enthalpy and entropy values, which led to postulating that proton donor and acceptor sites on the ion may form a H-bond network resembling water.300 Ludwig et al.273 measured far-IR spectra of [EtNH3][NO3], propylammonium nitrate ([PrNH3][NO3]), and dimethylammonium nitrate ([DmNH3][NO3]) in the range from 30 to 600 cm−1 that were assigned to the H-bonding stretching and bending vibrations. The DFT calculations deconvoluted the frequencies and assigned specific H-bond interaction. Frequencies between 199 and 224 cm−1 as well as between 134 and 159 cm−1 can be attributed to the asymmetric and symmetric stretching modes of the N−H···O H-bond, and those around 60−78 cm−1 can be assigned to bending modes of these Hbonds. The difference between the asymmetric and symmetric stretches is around 65 cm−1, suggesting comparable interaction strengths for the H-bonds. However, the symmetric stretching and bending modes of the H-bonds can be compared to those found for liquid water and ice, which leads to a conclusion that PILs easily form the H-bonding network. The structure of the Hbonding network was measured by neutron diffraction and empirical potential structure refinement.301 The average H-bond lengths for ethylammonium thiocyanate ([EtNH3][SCN], 1.71 Å) and ethylammonium hydrogen sulfate ([EtNH3][HS], 1.62 Å) possess much shorter H-bonds than other PILs. The H-bond angles show that short H-bond PILs have maxima at 180°, indicating a linear N−H···X arrangement is favored. Conversely, the long H-bonding PILs have maxima at 110°, indicating bent H-bonds are preferred. 4.2.3. Effect on Properties of ILs. Recent works152,272,277,302−304 have highlighted the roles of H-bond on IL properties. A finding is the local and directional H-bonds represent defects in the Coulomb system leading to more fluid rather than more rigid ILs. In this respect, the H-bonds have the

opposite effect on Tm, η, and enthalpy of vaporization (ΔvapH) compared with molecular liquids. X-ray photoelectron spectroscopy (XPS), NMR, and DFT calculations257,305 found that H-bond and charge transfer are independent each other in the imidazolium ILs, which suggests that other interactions, especially electrostatic force, are not changed substantially. Figure 7A shows the changes of Tm, η, and Δ vap H along the [C 2 C 1 im]X ILs and C 2 -methylated [C2C1C1im]X ILs (X is a serial of anions).306 It can be seen that the Tm and η values are increased from [C2C1im]X to [C2C1C1im]X.307 Bonhôte et al.308 found that methylation of [C2C1im][NTf2] at the C2−H position increases the η values substantially from 34 mPa·s to 88 mPa·s, whereas it is only less affected (η = 37 mPa·s) for methylation at the C5−H position. The same behavior is also observed for the Tm, which may be explained by considering the fact that the C2 proton of the [C2C1im]+ cation can engage in stronger H-bonding interaction than the [C2C1C1im]+ cation. However, MacFarlane et al.249,309 believed that the increase in Tm and η was related with the potential energy surface (PES) profiles of the anion moving around the imidazolium cation. In the PES of the [C2C1C1im]X ion pair, the anion was restricted to conformations above/below the imidazolium ring, whereas in the PES of the [C2C1im]X ion pair, the anion can freely move around, and the potential barrier that separates the two conformations exceeds 40 kJ/mol. Consequently, the methylated imidazolium ILs exhibit much higher Tm and η. Maginn et al.310 found that the entropy is likely an important factor in determining the Tm and η values of the ILs by MD simulation for [C2C1im][PF6] and [C2C1C1][PF6]. The ΔvapH of [C3C1C1im]X (X = anion) is increased. Ludwig et al.311 measured FIR spectra of 11 imidazolium-based ILs (Figure 7B) with the [C2C1im]+ cation and a series of anions, and they obtained a linear relationship (eq 29) between ΔvapH values and the vibrational frequencies. ΔvapH = 75.1 kJ/mol + 0.75 kJ/mol·cm−1 × ν

(32)

where ν is the wavenumber. A relationship312 between η and ΔvapH was reported as shown in eq 30, which suggests that for the same class of liquids ΔvapH is 6650

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

Figure 8. Illustrated TS structures corresponding to the QM/MM calculations for different reactions. (A) Diels−Alder reaction by acidic [C2C1im][Al2Cl7] IL and basic [C2C1im][AlCl4] IL.315,316 The distances (angstrom) are given from the 10 million overall Monte Carlo simulation. (B) Elimination of benzisoxazole using piperidine by [C4C1im][PF6] IL.318 (C) β-Elimination with piperidinem by [C4C1im][BF4] IL, and the distance (angstrom) are average values over the final 10 million configurations of the QM/MM simulation.321 (D) SNAr TS structure between morpholine and 2bromo-3-nitrothiophene by [C4C1im][PF6] IL.321 (E) Illustration of first-solvation shell [C4C1im][BF4] (shown as sticks) encapsulation of the TS.321 Reproduced with permission from ref 315, Copyright 2007, ref 316, Copyright 2007, ref 318, Copyright 2008, and ref 321, Copyright 2013, American Chemical Society.

structure requires a description of the bulk phase of supermolecular complexes, even couples with continuum method.314 The reactions including substitutions, eliminations, and creation of pericyclics is investigated by QM/MM to give the molecular origin of the observed rate accelerations and mechanistic changes in IL solution compared to conventional solution. For cyclopentadiene and methyl acrylate Diels−Alder reactions,315 the experimental reaction rate is 10 times faster in 51%-AlCl3 acidic [C2C1im][Al2Cl7] IL than that in water, but the rate is 2.4 times slower in the 48%-AlCl3 basic [C2C1im][AlCl4] IL. Accordingly, the computed free energy (ΔΔG‡, kcal/mol) values are 0.0, + 0.6, and −2.9 in water, [C2C1im][AlCl4], and [C2C1im][Al2Cl7], respectively. From the structures in Figure 8A, it is found that ion-pairing is less dominant in the acidic IL. The transition state (TS) via the H-bond with the [C2C1im]+ cation in the acidic salt is also greater than that afforded by the weaker Lewis-acid effect provided by the H-bond with water molecules.316 The ability of the IL to act as a H-bonding donor (cation effect), moderated by its H-bonding accepting ability (anion effect), may also explain the observed endo/exo ratios.317 The Kemp Elimination reaction of benzisoxazole via piperidine in IL [C4C1im][PF6] solvent was investigated by the QM/MM method. A 2D potential of mean force (PMF) was obtained along two reaction coordinates (proton abstraction by piperidine, RNH−RCH, and ring opening of the benzisoxazole ring, RNO, Figure 8B). The TS structure is followed, and the RNO distance was 2.06 Å, whereas the RNH and RCH distances were 1.10 and 1.75 Å. The ΔΔG‡ value is 25.2 ± 1 kcal/mol compared to the experimental value 22.6 ± 0.5 kcal/mol for the reaction under similar conditions.318 β-Elimination between 1,1,1tribromo-2,2-bis(3,4-dimethoxyphenyl)ethane and the cyclic

strongly correlated with the activation energy for viscosity, and the barrier height that molecules should overcome for diffusion/ rotation. ηVm ∼ RT exp(αHvapRT )

(33)

This effect allows us to tune important properties toward a higher temperature range and expand the range of potential applications of ILs. The effect from H-bonds appears to be markedly different from molecular liquids, possibly by introducing defects in the Coulombic lattice289 to promote directional interactions at the expense of electrostatic force151 favoring the liquid over the solid. This promotes higher order arrangements, for example, ion clusters or amphiphilic structures. Dispersion interactions are of importance for IL properties. Izgorodina et al. showed that calculated energies are closely related to measured melting points if dispersion forces are taken into account.170 For enthalpies of vaporization, it was observed that they increase nearly linearly with the increasing alkyl chain length of the imidazolium cation. A comparison with n-alkanes and n-alcohols showed that the linear increase in intermolecular interaction strength with each methylene group results from dispersion forces only.313 4.3. QM/MM Illustrating Reaction Mechanism

Despite the numerous reactions reported involving IL medium,16,23,24,71,72,83 the reaction thermodynamics and kinetics are rarely systematically reported. One reason is that the molecular factors are largely unknown. The multiscale QM/MM method has been applied to investigate the reaction mechanism of the IL-catalysis systems.94 The rationale that QM/MM is the most appropriate method for describing chemical reactivity in ILs by capturing important phenomena related to global solvent 6651

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

Dupont et al.,328 based on measurements of the diffusion coefficients of three electroactive solutes dissolved in ILs as a function of water content. A large difference was observed between the diffusion of neutral and charged species when comparing “wet” IL samples. The authors concluded the binary IL−water mixtures “may not be regarded as homogenous solvents, but have to be considered as nanostructured with polar and nonpolar regions”. The direct evidence of the ions clustering in ILs came from the ESI-MS experiments.329−332 AILs form cationic [Cq+1Xq]+ and anionic [CqXq+1]− supermolecular aggregates. Typically, the supramolecular aggregates of the type [C2A]+ (C = cation, A = anion) were observed in positiveion mode.329,332 The larger aggregates have been observed in the PILs. In the PILs [EtNH3][NO3] and [PrNH3][NO3], the [C8A7]+ aggregate mass/charge (m/z) peak was dominant.333,334 It is found that two and four ion-paired aggregates show cubiclike arrangement but six and eight ion-paired aggregates already show a high degree of structural disorder.145 The ESI-MS divides fragments violently to contribute to the clusters detection, and the polydispersed aggregates may not be evident for the underlying structure in ILs. The mesoscale aggregates in ILs have attracted much attention when the properties, reactions, and separation were correlated with such structures. Many physical phenomena in ILs were interpreted by the IL aggregates, including the surfactant-like micelles formed in IL−water mixtures and heterogeneous selfdiffusion.335−337 However, the issues so far are as follows: What are the aggregates in the neat ILs and ILs solution systems like? How do the aggregates form? Despite numerous pieces of evidence for the ion aggregates having been exhibited by IR, Raman, 10,265,338,339 surface tension, 340−357 conductivity,340,341,345,347−351,355,357−367 NMR,337,342,343,347,349,355,361−364,368 fluorescence spectroscopy,340,347,362,364 dynamic light scattering,337,354,357,365 and small angle neutron scattering,340,343 we would like to emphasize the results of simulations because they provided the direct images of the mesoscale structures.

amines piperidine and pyrrolidine, as a fundamentally important reaction, has been studied in IL solution.319 It was found that the reaction pathway changed from an irreversible E1cb mechanism to an E2 pathway when [C4C1im][BF4] and [C4C1im][PF6] ILs replaced methanol as solvent.320 In order to distinguish the mechanism between E1 cd and E2, Acevedo et al.321 performed QM/MM calculation. It was found that the beneficial electrostatic interactions and π+−π interactions between the [C4C1im]+ cation and β-phenyl substituents are principal contributors to the mechanism change in the structural configuration of the ions (Figure 8C), although E2 possesses a significant amount of E1cb character and no carbanion intermediate was located. Therefore, in the E2 pathway, the amine abstracts the β-hydrogen with concurrent cleavage of the α-C−Br bond; in contrast, the E1cb process requires two separate steps that figure out the abstraction of the β-hydrogen to be rate limiting and often reversible. Nucleophilic substitutions have been reported to derive sizable rate enhancements from the ILs [C4C1im][BF4] and [C4C1im][PF6] relative to conventional solvents.322 The nucleophilic aromatic substitution (SNAr) reaction between cyclic secondary amines (i.e., piperidine, pyrrolidine, and morpholine) and the 2L-5-nitrothiophene (para-like) and 2-L-3-nitrothiophene (ortholike) isomers has been computationally investigated to determine the mechanism for enhancing reactivity in the salts. The calculated ΔG‡ values for the SNAr reaction by perturbing the RCN distance between the reacting ipso carbon of the solute and the amine nitrogen in increments of 0.02 Å (Figure 8D) using the QM/MM method321 are 25.7, 26.5, and 25.1 kcal/mol in methanol, respectively; however, the ΔG‡ values are 25.8, 27.6, and 23.8 kcal/mol in [C4C1im][BF4] IL when the L is bromo, methoxy, and phenoxy, respectively. Multiple factors lead to beneficial IL effects on SNAr reactivity: (1) π+−π interactions between the [C4C1im]+ cations and the aromatic rings can potentially enforce a coplanarity structure; (2) significant electrostatic enhancements between the solvent and the developing charge separation formed at the TS compared to the case for the reactants; (3) the highly ordered IL structure produces site-specific H-bonding stabilization that, despite an entropy penalty, provided a 5 kcal/mol reduction in ΔH‡ compared to that for methanol. The encapsulation of the addition step TS and the Meisenheimer intermediate finds 8−9 ions that are within 3 Å. Figure 8E illustrates the [C4C1im][BF4] ions forming a cage-like structure to favorably interact with the substrates.321

5.1. Bulk Structures of ILs

5.1.1. Self-assembly in Neat ILs. Ionic aggregation and selfassembly in neat ILs were reported by Voth et al.96,214,215 through CG simulations. Figure 9A shows the CG partition for the [C2C1im][NO3] ion pair, in which the 23 atoms are coarsegrained to be 5 CG sites, and site A represents an imidazolium ring and site D represents the [NO3]− anion. Due to the simplified interactions, the CG simulation is highly efficient, with a larger MD simulation time step and faster ion diffusion. Figure 9B shows one snapshot of ionic aggregates, in which the groups of cation’s tails form nonpolar domains, while the charged head groups of the cations and anions form polar domains with a continuous charged network by retaining their local structures. The simulation found that the tail domains remain relatively stable at a low temperature, but the domains begin to diffuse with increasing temperature; however, when the temperature increases beyond a point, the average density distributed almost uniformly, although tail groups’ instantaneous domains still existed.216 Voth et al.97 found that the finite size effect on the spatial heterogeneity was small. Partial structures displayed a peak at a wave vector smaller than the main peak, indicating the occurrence of an intermediate range order in 1-butyl- and 1octyl- derivatives due to the presence of long alkyl chains.370 Furthermore, the nano-organization was observed in 1,3dialkylimidazolium-based ILs.371,372 An important conclusion is

5. MESOSCALE STRUCTURES OF ILS The liquids, in general, are thought of as the homogeneous bulk. However, the structures of water at nanoscale have been described as “flickering clusters”, linear molecular chains, or larger units.323,324 Even in the mixing solution of water−alcohol, most of the water molecules exist as small H-bonding strings and clusters in a “fluid” of close-packed methyl groups.325 This structures suggest that the anomalous thermodynamics of water−alcohol systems arises from incomplete mixing and from retention of remnants of the three-dimensional H-bonding network structure. ILs bear large molecular volume, strong electrostatic force, and H-bonding direction to prompt them to form heterogeneous structures in bulk phases.326 This spatial heterogeneity is thought to be one of the key features in interpreting many physical phenomena of ILs, including the heterogeneous self-diffusion, surface layering, and surfactant-like micelles formed in IL−water mixtures327 The first mesoscopic IL structure was suggested by 6652

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

and nonpolar domains formed by the color-coding with red and green. An obvious conclusion is that the [C2C1im][Gly] shows less heterogeneity, and [C6C1im][Gly] and [C10C1im][Gly] show obvious heterogeneity. Electrostatic forces are always dominant between ions in the short chain ILs, but VDW interaction increases with longer alkyl chain length in the ions aggregation. In addition, the aggregates were also observed in nonimidazolium guanidinium-based ILs by analyzing the structural and energetic features.381 Optical heterodyne-detected Raman-induced Kerr effect spectroscopy (OHD-RIKES) suggested that heterogeneity exists in these densities of ILs.382 Wide-angle X-ray scattering383 observed the aggregation deduced through the intermolecular partial atom−atom correlation function between terminal methyl carbons of the butyl chain in [C4pyr][NTf2] IL. Diffusion ordered NMR spectroscopy384 was used to determine the aggregate number of the diastereomeric (1R,2S)-ephedrinium (RS)-methoxytrifluorophenylacetate IL. Atmospheric pressure measurements385 found the microheterogenous structures with polar and nonpolar domains in the 1-butyl-3-methylimidazolium octylsulfate IL. ESI-MS experiment333 showed that the formation and size of aggregates are dependent on the nature of the ions and can strongly influence the entropic component of the free energy of amphiphile self-assembly in some PILs. Importantly, the nonpolar region of the IL is commonly manipulated by the length of alkyl side chains, which can be incorporated in either the cation or the anion. Thus, many properties of the ILs could be tailored by combining different kinds of cation and anion, or changing the proportion of polar and nonpolar parts according to the proportion of electrically charged territories versus apolarlike regions.96,386−389 The structural dynamics observed occur on a much faster time scale than the overall orientational relaxation of the ILs measured through OHD-RIKES.390 The results indicate that overall structural relaxation is not required for most local anions from observable to entire relaxation. Such a picture is that of complex layers of structure and dynamics that arise from the spatial heterogeneity and need to be considered when interpreting properties of bulk ILs and optimizing their properties for specific applications. Influence of microscopic structures on bulk properties, such as the effect of H-bond on viscosity152,391,392 or the effect of adding water,393 are different from experiment, and even fail in predicting the conductivity.394 The reason is the properties involve local structure and dynamics or, like conductivity, are intimately tied to local effects; thus, exploring the links between composition, ordering, and dynamics is a necessary step for the identification of the structural handles and

Figure 9. (A) CG scheme for the [C2C1im][NO3] IL.215 The above shows the molecular structure of the[C2C1im][NO3] ion pair. The blue spheres represent nitrogen atoms, cyan is carbon, white is hydrogen, and red is oxygen. The below shows the CG structure. (B) One snapshot illustrating the spatial heterogeneous aggregation.96,369 The white spheres represent the cationic terminal groups, the gold spheres represent the cationic head groups, and the red spheres represent the anions. The ellipses in blue illustrate the approximate positions of the isolated tail domains. Reproduced with permission from ref 215, Copyright 2006, ref 96, Copyright 2007, and ref 369, Copyright 2005, American Chemical Society.

that aggregation of nonpolar domains composed or the alkyl chains was observed for ILs with alkyl side chains longer than or equal to C4. Increasing the length of the alkyl chain caused the nonpolar domains to become larger in a manner analogous to [CnC1im][PF6] (n = 2, 4, 6) systems exhibiting microphase separation (Figure 10).373,374 Spheral aggregates with diameters of 2−8 nm can be observed in the IL [C4C1im][PF6].375 The ions distribute with butyl tail groups of cations protruding outward from the aggregate surface, and the ring centers lying beneath and the number density distribution of imidazolium-ring centers was not uniform near the surface compared to anions. CG simulations for [CnC1im][PF6] were performed376−379 to compare the structures and dynamics with all-atomic models. The cation is generally grouped by three sites and anions by a single site. The structural and energetic properties by the CG model are in reasonable agreement with all-atomic results, the ILs show nano-organized structures whose character depends on the length of the alkyl chains as shown in Figure 10; however, the dynamics are unrealistically slow. A temperature shift of ∼100 K is required to produce agreement between the viscosity and diffusion coefficients of the model and experimental values. A series of amino acid ILs ([CnC1im][Gly], n = 2, 6, 10) were investigated by MD simulations.380 Figure 11 shows the polar

Figure 10. Snapshots of nonpolar domains of [CnC1im][PF6] IL calculated by MD simulation.373 (A) [C2C1im][PF6]; (B) [C4C1im][PF6]; (C) [C6C1im][PF6]. Reproduced with permission from ref 373. Copyright 2006 American Chemical Society. 6653

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

Figure 11. Snapshots of ion aggregates in the imidazolium amino acid ILs by MD simulation.380 (A) [C2C1im][Gly]; (B) [C6C1im][Gly]; (C) [C10C1im][Gly]. Reproduced with permission from ref 380. Copyright 2010 American Chemical Society.

Figure 12. Snapshots of [C8C1im][NO3] at different water concentrations by CG simulation.214 The water mole fractions are 0%, 20%, 50%, 75%, 80%, and 95.2%, respectively. The polar domains, nonpolar domains, and water are colored red, yellow, and dark blue, respectively. Reproduced with permission from ref 214. Copyright 2007 American Chemical Society.

isolated ions or small ion clusters, respectively (XH2O > 0.95).399 Similar structures occur in [C6C1im][NTf2] aqueous solution.400 At low concentration, water molecules were surrounded by several cation−anion pairs. Phase separation could be observed at medium water concentration, while at high water concentration, cationic tails aggregate to form the micelles. It was suggested that the water-induced increased spatial correlation among the polar components lead to the initial increase in density, and the density could rise more obviously in the system of ILs with longer cation tails due to hydrophobic clustering at moderate water concentration.401 Besides the imidazolium ILs, the structures were also found in phosphoniumbased ILs. Wang et al.402 reported that the water molecules are dispersed and embedded in cavities between ions and tend to aggregate to form small clusters, intermediate chain-like structures, large aggregates, and eventually a water network with increasing the water proportion in the mixture of trihexyltetradecylphosphonium bis(oxalato)borate IL ([P6,6,6,14][BOB]) and water at low concentration. In the dilute solution, self-organized micelle-like aggregates can form composed of a hydrophobic core and a hydrophilic shell.403 Triolo et al.404 observed the micellar-like morphology for [CnC1im][PF6] ILs by XRD to explain the linear correlation between the alkyl chain length and the position of the observed bulk length.

tuning parameters that may mitigate impediments to applications. 5.1.2. Aggregates in Aqueous Solutions. The amphiphilic nature of the cations causes ionic aggregation in the solutions.168,395−397 Imidazolium cations with long appended fluorous alkyl chain tails can act as surfactants to facilitate emulsification when added into solvents.398 The aggregates are especially affected by both the solution concentrations and alkyl chain length. First is the effect of concentration. Figure 12 shows the structures of [C8C1im][NO3] at different water concentrations by CG simulation.214 The polar network was broken up by adding water, while both the water network and the micelle exhibited turnover. The mechanism analysis found that the competition between the hydrophobic interactions of the nonpolar groups and the hydrophilic interactions of polar groups results in the turnover for ion micelle aggregation; meanwhile, the water network may result from the competition between the water molecules and water−anion interactions. In aqueous solutions of 1-ethyl-3-methylimidazolium ethylsulfate ([C2C1im][EtSO4]), four concentration ranges can be identified with distinct structural regimes: isolated water molecules (XH2O < 0.5), chain-like water aggregates (0.5 < XH2O < 0.8), a bicontinuous system (0.8 < XH2O < 0.95), and 6654

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

[C8C1im][BF4] aqueous solutions were also simulated and the [C8C1im]+ cations displayed stronger aggregation and slower diffusion.237 Dicationic imidazolium ILs with a tetradecyl chain exhibit a higher thermal stability and potential application.411−414 As a gemini surfactant in aqueous solutions, the dicationic [C5(C10im)2]Br2 IL forms cross-linked micellar aggregates with near-hexagonal arrangement of the hydrophobic cores connected by hydrophilic head groups.415 The aggregates interconnection was mediated by head groups, which could not be observed in monocationic IL solutions. Starting from a random distribution, the hexagonal structure evolved in the 40% (w/w) aqueous mixture as the snapshots shown at various time intervals in Figure 14.416 Xu et al.417 investigated the surface activity and thermodynamic properties of the gemini imidazolium surfactants ([Cn-4-Cnim]Br2, n = 10, 12, 14) through surface tension and electrical conductivity methods. They found418 that the micellization of [Cn-4-Cnim]Br2 is entropy-driven at low temperature, but at high temperature the aggregation is enthalpydriven. The micropolarity of micelles increases, but the aggregation numbers decrease with increasing the spacer length of the IL for the other gemini imidazolium ILs ([C12-sC12im]Br2, s = 2, 4, 6). For the micelle aggregates in the [CnC1im]Br (n = 9, 10, 12, 14, 16) aqueous solution, a linear relationship between the logarithm of the critical micelle concentration (CMC) and the number of carbon atoms of hydrocarbon chain was found.359 Sirieix-Plénet et al.360 found that the [C10C1im]Br behaves like a classical cationic amphiphile at low water concentrations. Particularly, Zheng et al.344,347 found the lower CMC and higher micellar number may be caused by the enhanced hydrophobicity and the π−π interaction in 1-[n-(N-carbazole)alkyl]-3-methylimidazolium bromide ([carbazoleCnC1im]Br (n = 6, 10, and 12)) solution when they incorporated the fluorescent carbazole moiety at the terminal of the hydrocarbon chain by electrical conductivity and NMR spectroscopy. A looser packing of the ILs in the micelles was suggested caused by the hydrophilic head group size and electrostatic repulsion.348 Wang et al.362−364 found that no aggregation occurs for [C4C1im]Br and [C6C1im]Br, forming aggregates above the CMC in aqueous solutions, and the aggregation of [C8C1im]Br, [C10C1im]Br, and [C12C1im]Br is well-defined to form micelles above the CMC. There are also many works on the ILs with Cl− or I− anions.340,341,358 CMC values of [CnC1im]Cl (n = 4, 6, 8, 10, 12, 14, 16, and 18) were obtained. An approximate linear correlation between the retention time of cation variation with Cl− and CMC values was repoted.350,366 Butts et al.343 observed the small near-spherical aggregates in [CnC1im]Cl (n = 2, 4, 6, 8, 10) solution with water at concentrations just above the CMC.

Alkyl chain length is another factor that has an effect on the ionic aggregation. Some works focused on the ILs with halide anions (Cl−, Br−, and I−), especially for Br−.340,341,343−345,347,348,350,356,358−360,362−364,366,405 In the [CnC1im]Br (n = 12, 14, 16) aqueous solutions, the cations spontaneously self-assembly from quasi-spherical polydisperse aggregates, and the aggregate numbers steadily increase with increasing the chain length, and the most ordered structure was found in [C16C1im]Br.406 In the [CnC1im]Br (n = 2, 4, 6, 8) aqueous solutions with shorter chains, the VDW interaction between the chain tails is responsible for the aggregation of the tails of cations. The [C2C1im]Br remains isotropic, [C4C1im]Br weakly associates to form small clusters, [C6C1im]Br forms the definite aggregates, and [C8C1im]Br forms the well-defined aggregates. Figure 13 shows the distribution of cations into an

Figure 13. Distribution of cations into aggregates of size n in [CnC1im] Br aqueous solution.407 Reproduced with permission from ref 407. Copyright 2009 Royal Society of Chemistry.

aggregate with size n in the aqueous solutions.407 The number and the size of the aggregates increased with increasing alkyl chain length, and more spherical aggregates formed in solutions of [C10C1im]Br.408 Additionally, in [CnC1im][PF6] (n = 4, 7, 10) aqueous solutions with longer side chains (n = 7 and 10) a bicontinuous morphology is exhibited.377 Especially, the micellar aggregates were observed in the solution of [C10C1im]Br.409,410 Inside the micelle, the alkyl chains bury themselves and the polar head groups are orientated to water. Besides the growth mechanism, fusion of the micelles to form larger clusters was proposed. To elucidate the effect of alkyl side chain length on the structure and dynamics of the solutions, the [C4C1im][BF4] and

Figure 14. Snapshots of aggregation of a dicationic [C5(C10im)2]Br2 IL in a dilute bulk system at different simulation times.416 (a) 1 ns, (b) 16 ns, (c) 110 ns, and (d) 1 μs. Polar and nonpolar groups are shown in yellow and magenta, respectively. Water and anions are not shown. Reproduced with permission from ref 416. Copyright 2011 American Chemical Society. 6655

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

result in a cooperative character of ion diffusion and even in an appearance of diffusion anisotropy for the cation diffusion.368 The surface activities of 1-butyl-3-methylimidazolium bis(2ethylhexyl)sulfosuccinate ([C4C1im][AOT]) and 1-butyl-3methylimidazolium bis(2-ethylhexyl)phosphate ([C4C1im][DEHP]) are superior to those of their sodium analogs, NaAOT and Na-DEHP.357 When adding NaCl, NaBr, and NaI into the solutions of longchain imidazolium ILs, the CMC values decreased significantly,346 which suggested the hydrophobicity of the anions plays important roles on the aggregation of [C10C1im]Br in aqueous solutions.365 Sadeghi et al.352 studied the surface and micellar properties of [C12C1im]Br in aqueous solutions with different kinds of electrolytes. They found that the surface tensions as well as the CMC values decrease with adding electrolytes. The aggregation promoting ability decreased in the order I− > Br− > Cl−, and (C4H9)4N+ > (C3H7)4N+ > (C2H5)4N+ > (CH3)4N+, respectively. The aggregation behavior in pyridinium-, piperidinium-, and pyrrolidinium-based ILs was also observed.342,351,367 CMC values provided early evidence for aggregation in [CnC1py]Cl, (n = 4, 10, 12, 14, 16, or 18), [C12C1pip]Br, [C12C1py]Br, and [C12C1pyrr]Br by fluorescence technology.342 Aggregation in aqueous solutions of ILs, [C4pyr][BF4], [C4pyr][CF3SO3], [C4C1pyr][N(CN)2], and [C8C1pyr][BF4], was characterized through speed of sound and conductivity measurements.367 Wu et al.424 considered the water effect on the interaction energies and local packing of [C2C1im][Gly]. They indentified that the polar network was continuously broken up and the nonpolar groups aggregation enhanced as the amount of H2O increased. The micellar shuttle pathways of poly(1,2-butadiene)− poly(ethylene-oxide),425 poly(N-isopropylacrylamide-block-ethylene-oxide),426 and poly(ethylene oxide)427 copolymers in a mixture of [C4C1im][PF6] and water was investigated by DPD simulation at different temperatures. The micellar transfer presents four stages in general (Figure 17): (1) formation and stabilization of PB−PEO micelles in the watery phase, (2) migration of PB−PEO micelles toward the water−IL interface, (3) diffusion of PB−PEO micelles into the IL phase, and (4) complete diffusion of PB−PEO micelles into the IL phase. The micelle shrinks a little as the temperature is raised, and it keeps its cuasi-spherical shape as it goes though the water−ionic liquid interface. 5.1.3. Aggregates in Nonaqueous Solutions. The aggregation of ILs in some other solutions, including the supercritical CO2 (scCO2),428 other ILs,429 and organic solvents,405,430−434 was also observed. Takamuku et al.435 found the aggregation of [C2C1im]Cl and [C2C1im][NTf2] in methanol, acetonitrile, and benzene by small-angle neutron scattering and NMR. Polarized optical microscopy and small-

No evident aggregates were found in [C2C1im]Cl and [C4C1im] Cl systems, while oblate aggregates appeared with radius ∼9 Å in the [C6C1im]Cl system. For the ILs with fluorinated anions [PF6]−, [BF4]−, and [(PFBu)SO3]−,340,341,349,354,361,367,419 IR and Raman spectroscopies found that the water can associate with anions at low water contents, while the water exhibits an aggregating network at high concentration. It is proposed that water clusters organized in the polar network, and the hydrophobic alkyl tails aggregated in a micelle-like structure.419 The short-chain [C4C1im][BF4] can be best modeled as a dispersion of polydisperse spherical aggregates that form above a critical aggregation concentration.340 Figure 15 is the proposed schematic models of the aggregates for the different ILs.361

Figure 15. Schematic representation of polydisperse spherical aggregates of different ILs.361 (A) [C4C1im][BF4]; (B) [C4C1im]Cl; (C) [C4C1py]Cl; and (D) [C8C1im]Cl. Reproduced with permission from ref 361. Copyright 2007 American Chemical Society.

Besides spherical, near-spherical, or near-hexagonal aggregates, the vesicle aggregates were also reported in some ILs systems.420 The aggregates of [CnC1im]Br (n = 10, 12, 14) ILs change from spherical to rodlike micelles and then to vesicles with the increasing concentration of ILs. Liu et al.421 revealed a forming process of the unilamellar vesicle in [C12C1im]Br aqueous solution by MD simulations (Figure 16). Interestingly, the inner layer of the polar part is packed more densely than the outer layer, and the H-bond numbers at different concentrations for binding between counterions enhanced with increasing IL concentration.421,422 Additionally, the vesicles were observed in two kinds of trisiloxane-tailed ILs in water through using the MARTINI force field.423 The imidazolium-based ILs composed of pH responsive anions can form the vesicle aggregates without any additives.337 The ILs [CnH2n+1C1im][CmH2m+1SO3] (n = 8, 10, or 12; m = 1 and n = 4 or 8; m = 4 or 8) behave as conventional single chain cationic surfactants with relatively small methylsulfonate anions (n = 8, 10, and 12, m = 1), and the CMC values are lower than those of the corresponding cationic analogues.353 The aggregation ability of the anionic surface active ILs [C4C1im][CnH2n−1O2] (n = 8, 10, 12) is mainly caused by the lower electronegativity of their anions.355 Additionally, aggregates in ILs ([C2C1im][EtSO4] and [C2C1im][TfO]) aqueous solutions

Figure 16. Formation of the unilamellar vesicles in the [C12C1im]Br aqueous solution by MD simulation.421 Green depicts [C12mim]+, and magenta depicts water. Snapshots are shown at (A) 0.6 ns within 6.0 nm, (B) 1.0 ns within 6.0 nm, (C) 20 ns within 4.5 nm, (D) 120 ns within 4.5 nm, and (E) 220 ns within 4.5 nm. Reproduced with permission from ref 421. Copyright 2016 American Chemical Society. 6656

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

guanidinium trifluoroacetate.438 At some pressures, CO2-in-IL microemulsions were first reported when 1,1,3,3-tetramethylguanidinium acetate was mixed with the surfactants [EtNH3][SO2C8F17] and scCO2. Figure 19 illustrates the formation of a CO2-in-IL microemulsion.439 The CO2-swollen micelles can be “tunable” when the micellar size is adjusted easily by changing the pressure of CO2.

Figure 17. Micellar shuttle pathway by DPD simulation.425 (A) Quasispherical PB−PEO copolymer micelle dispersed in the water−IL biphasic environment, T = 298 K. (B) Migration of PB−PEO micelles toward the water−IL interface, T = 303 K. (C) Diffusion into the IL, T = 345 K. (D) Micelle movement into the IL, T = 360 K. (D) Inverse process of micellar shuttle pathway from IL phase to aqueous medium, T = 298 K. (F) At T = 303 K, the IL or water molecules have been deleted in some cases for better visualization. Reproduced with permission from ref 425. Copyright 2012 Elsevier.

Figure 19. Schematic illustration for the formation of a CO2-in-IL microemulsion:439 (a) “dry” micelle dispersed in IL; (b) CO2-bound micelle; (c) CO2-swollen micelle. Reproduced with permission from ref 439. Copyright 2011 Wiley-VCH Verlag GmbH & Co. KGaA.

The structure of guanidium acetate IL in scCO2 was investigated by MD simulations. The result indicated the existence of stable IL droplets within the continuous scCO2 phase. Besides, it was revealed that the IL anion−-head group interactions may be the primary reason for the stability of the microemulsions.440 It is found that hydrophobic IL [C4C1im][PF6] can be dispersed in the hydrophilic IL propylammonium formate by adding the surfactant sodium bis(2-ethylhexyl)sulfosuccinate. The droplets in the microemulsions were spherical, and the sizes of the droplets increase the [C4C1im][PF6] amount.441 Besides, emulsions formed in the mixture of the hydrophilic propylammoniumformate and hydrophobic [C8mim][PF6] with the aid of the surfactant sodium bis(2ethylhexyl)sulfosuccinate.442 Except the supercritical CO2 and the other ILs, some binary and ternary mixtures of [C4C1im][PF6] or [C4C1im][BF4] with organic solvent systems were also investigated by combination of UV−vis spectroscopy and conductivity measurements. The relationship between the degree of aggregation and the dielectric constants of the solvents was clarified.443 Moreover, mixture of ILs with other organic solvents were also investigated. Microemulsions of [C4mim][BF4]/Triton X-100 (TX-100)/cyclohexane were characterized by conductivity measurement, dynamic light scattering, freezefracturing electron microscopy, and UV−vis techniques. In the microemulsions, the heterogeneous structures with IL as the polar core were obtained.444 Mixtures of two kinds of ILs were also studied by multiscale DPD simulation. Zhao et al.436 obtained similar aggregation of [C16C1im]Cl in the [EtNH3][NO3] by DPD. Different topography aggregates, including micelles, normal hexagonal, lamellar, and reverse bicontinuous cubic liquid crystalline phases, were observed. A series of thermodynamic parameters (ΔGm° , ΔHm° , and ΔSm° ) of the micellization of N-alkyl-N-methylpiperidinium bromide ([CnPD]Br, n = 12, 14, 16) in the [EtNH3][NO3]445 were investigated by surface tension and DPD simulation, and the results showed that the micellization was entropy-driven.

angle X-ray scattering also observed the like micelles, normal hexagonal, lamellar, and reverse bicontinuous cubic liquid crystalline phases aggregation of [C16C1im]Cl in [EtNH3][NO3] IL.436 Sphere aggregates with a much larger size than traditional micelles when mixing [CnC1im]Br (n = 10, 12, 14, 16) with [C4C1im][BF4] were probed through measuring surface tension.429 Figure 18 presents the aggregates of S-3-hexadecyl-

Figure 18. Schematic representation of [C16hpim]Br IL aggregates in water (A) and [C4C1im][BF4] IL (B).429 Hollow ellipses and hollow circles represent [C4C1im]+ and Br−, respectively. Reproduced with permission from ref 429. Copyright 2008 Royal Society of Chemistry.

1-(1-hydroxy-propan-2-yl)-imidazolium bromide ([C16hpim]Br) in water and in ILs, respectively, showing the superior capacity for micelle formation compared to traditional ionic surfactants. Han et al.437 found that the properties of the micelles in scCO2 can be affected by the IL [C4C1im][BF4], indicating that the IL enhances the micropolarity and enlarges the size of the micelles considerably, as its amount was very small. They also found that N-ethyl perfluorooctylsulfonamide ([EtNH3][SO2C8F17]) can form reverse micelles in scCO 2 when adding 1,1,3,3tetramethylguanidinium ([(CH 3 ) 2 N] 2 CNH 2 + ) acetate, 1,1,3,3-tetramethylguanidinium lactate, and 1,1,3,3-tetramethyl6657

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

5.2. Interface Structures

segregation of the polar and apolar molecular moieties in surfacefrozen alcohols. The single-box model449,450 fitting the measured XRR (R(qz), qz = (4π/λ) sin α, where α is the grazing angle of incidence of the X-ray beam onto the surface) by a distorted-crystal model was developed to describe a sine-wave-like oscillatory density profile, the amplitude of which decays with increasing depth z below the surface.165 The measured R(qz) by this model invariably yields a very short decay length, such that only the first oscillation is observed, and the second one is already too small to observe above the constant bulk density. Consequently, the density profile obtained using this model comes out to be identical with that obtained from the single-box model. Therefore, in contrast with liquid metals, the ILs do not exhibit an oscillatory surface electron density profile but only a slight enhancement of the electron density in a 6−7 Å thick surface layer. However, MD simulation for the dimethylimidazolium chloride IL yielded a layered surface with an oscillatory density profile, and the reason may lie in more spherical cations than the butyl-tailed cations.167 Further MD simulation found that the number density, mass density, and electron density occur oscillating at the [C4C1im][PF6]−vapor interface.451 However, the oscillations in the electron density are much diminished when compared to the amplitude of oscillations in the number densities. The electron density at the interface is about 12% larger than that in the bulk. This result is in excellent agreement with results of XRR measurement.165 The butyl chain of the cations at the interface exhibits a preference to orient along the surface normal, and the surface has a certain amount of hydrophobic character. The observation of large-amplitude oscillations in the number density at the interface is in disagreement with the measurement from XRR. The underlying reason is the difference in the molecular sizes of the two species, being appropriate to nearly cancel out

5.2.1. Structures of Liquid−Gas Interfaces. IL liquid−air interfaces have been investigated by XRR, sum frequency generation spectroscopy (SFG), and neutral impact collision ion scattering spectroscopy (NICISS).446−448 XRR studies on [C4C1im][PF6] and [C4C1im][BF4] suggested the existence of a surface layer 6−7 Å thick, an electron density 10−12% higher than that of an equivalent simple, and an unstructured liquid surface.165 Two kinds of ions arrangements were proposed at the surface layer, one with the butyl chains normal to the surface (Figure 20A and B) and the other with the chains parallel to the

Figure 20. Possible ion arrangements of [C4C1im][PF6] near the air− liquid surface: (A and B) with the butyl chains normal to the surface, and (C) with the butyl chains parallel to the surface. The six-arm stars symbolize the [PF6]− anions, and the wavy line is the surface. Note that for (B) the imidazole groups are excluded from the surface layer. The obtained density ρbmim above requires dense packing of the alkyl tails if no imidazole groups are present in the surface layer. The incorporation of the spherical anions into such a hydrophobic layer is unlikely, so that a significant incorporation of the anions into the hydrophobic layer of butyl tails would be needed to account for the high electron density obtained from the XRR fits. In (C), a stack of two to three layers of lyingdown ions is presented.

surface (Figure 20C). Both arrangements differ from the surface ordering observed in other systems to date, namely the oscillatory surface-normal density profiles (ρ) observed in the condensed alkyl chain layer in surface-frozen alkanes, and the

Figure 21. Topographies and structures of the AIL [C2C1im][NTf2] and PIL [EtNH3][NO3] near the mica interface. (A) AFM image of the [C2C1im][NTf2] innermost layer adsorbed to mica; the arrows indicate areas of different heights.465 (B and C) Top and side views for the MD simulation structure of the [C2C1im][NTf2] absorbed in the interfacial layer.464 (D) Phase profile along a side axis. (E) Amplitude (dotted) and phase (black) oscillation from the AFM tip approaching the [EtNH3][NO3]−mica interface. Arrows indicate the direction of the axis. (F and G) Topographic AFM images of the innermost ion layer and the first near-surface layer at the [EtNH3][NO3]−mica interfaces.465 Section analyses of the near-surface structure, as indicated by the blue line, are shown as insets to the images. Reproduced with permission from refs 464, Copyright 2012 John Wiley & Sons, and ref 465, Copyright 2013 Royal Society of Chemistry. 6658

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

is consistent with the layer thickness suggested by the height profile in Figure 21C. For PIL [EtNH3][NO3], three different regimes were identified at the mica interface: (1) a sponge-like nanostructure in the bulk (far from the interface), (2) a more lamellar-like arrangement (close to the interface), and (3) rows of worm-like structures in contact with the surface, as identified by phase oscillation in Figure 21E. Subsequent AFM measurements showed that both the adsorbed and the near surface structure determine nanoscale friction of the IL-coated interfaces in Figure 21F and G. Given the structural simplicity of the [EtNH3][NO3]−mica surface, these results are far-reaching and suggest that many ILs will exhibit similar patterns of self-assembly at different solid−liquid interfaces, such as at the surfaces of montmorillonite and kaolinite. However, MD simulations329,330 show the lower charge density, as well as rougher and more porous solid surface, which leads to weaker ion ordering. Carbon−IL interfaces are particularly interesting because they have been fueled by the desire to find new composite materials. The interfaces have been classified by neutral interfaces and an electrified one. The former is important to understand solvation of nanotubes,49 graphene,466 and C60467 in ILs as well as the confined structure of ILs in these materials,47,468,469 and the latter is particularly important for electrochemical processes. Close to neutral surfaces, such as the graphite, AFM,470 and XRR,471,472 studies suggest that the ions are arranged by layer and the number density is increased compared to the bulk of ILs. On graphite, the IL layer varies with applied surface bias and shifts from cation-enriched (negative biases) to anion-enriched (positive biases), similar to that described by Canongia Lopes et al.473 at a silica substrate. Interestingly, the dynamics of ions is related to their relative position near the interface; for example, ions in layers closer to the graphene surface display longer relaxation times and increasingly slaved diffusion.474 Many MD simulations at the different interfaces, such as graphite,475−477 graphene,48,52,471,478 amorphous carbon,479 nanotube,48,469,475 and nanopore interfaces,480 have been performed to obtain the explanation of the experimental phenomena. In more recent years, attention has been paid to electrochemically active interfaces. Several reviews73,481,482 have emerged to understand the structure and electric potential distribution of ILs at the charged interfaces. The conclusion is that ILs exhibit the electrical double layer (EDL) distribution near the charged interface. Numerous experimental studies have been performed to investigate the IL EDL structures. Early reports using SFG459,460,462 revealed a monolayer structure with interfacial cations exhibiting orientational ordering, whereas a structure with an adsorbed ion layer, similar to aqueous electrolytes, was suggested by other experiments.483 Nanometer resolution of the IL EDL structure normal to the interface was probed via AFM for the AILs.455,484 It is found that more than one layer thick ion-pair layers form at the electrode interface, which were quite different from the case in the aqueous electrolytes, as shown by in situ electrochemical AFM force measurements, and increase in number and become more strongly arranged at higher voltages. The AFM data suggested the capacitor-like double layer structure present in ILs oscillates between each plane of charge with a potential decay. In addition, elegant XRR experiments on Au485 and graphene472 probed the structural reorganization of the layers in real time as a function of potential. The results showed that ILs were subject to slow potential-dependent kinetics and structural hysteresis, both of which related to the unusually well-defined layering and slow ion

their contributions to the total electron density. MD simulations at the liquid−vapor interface of [C2C1im][NO3] at 400 K452 using a polarizable model exhibit similar features of an out-ofphase oscillation in the number densities of cations and anions. For the air−liquid surface of PILs, the XRR results indicated a rough or diffuse interface with significant gas interpenetration.446 The interfacial structure was postulated to be dotted with small, dynamic clusters of anions and cations in which the short alkyl moieties surround the charged core in a roughly spherical geometry. The SFG and NICISS data supported this conclusion, as the hydrophobic alkyl groups were shown to shield the charged groups and protrude outward into the air, which suggests that the enormous air−liquid surface area created upon droplet ionization contributes to aggregate detection.447 SFG gave a structural feature at the gas−liquid interface of the water miscible IL.453 When the IL is dry, the cations are oriented with the imidazolium ring parallel to the surface plane for both hydrophilic and hydrophobic ILs. However, when the water is added, the cation reorients to be perpendicular to the surface for the hydrophobic liquid, while the orientation of cations in the hydrophilic liquid is unaffected. In comparison, the water molecules in hydrophilic ILs are stable by intermolecular interaction, such as H-bond and dipole−dipole force; thus, it is less favorable to partition water to the surface. However, the water is less stable in hydrophobic ILs and is forced to the surface. This is evidenced by the appearance of the C−H stretching mode on the imidazolium ring in [BMIM][imide]. 5.2.2. Structures of Liquid−Solid Interfaces. The structures near solid−liquid interfaces are usually very different from that of the bulk phase. The interfacial structures and dynamics have been investigated because the spatial distribution of ions at solid−liquid interfaces is a central scientific problem, as it controls processes as diverse as amphiphile adsorption, colloid stability, lubrication, heterogeneous catalysis, charge transfer, and protein folding.24,53,73,335,454 In the early period, as reviewed by Hayes et al.,455 the solid surfaces often involved mica, silica, quartz, and metals, which are smooth and have a net negative potential.456 Experimental approaches, including SFG (silica,457 quartz,458 Pt,459,460 and TiO461,462), XRR165 and neutron reflectivity,463 have been employed to probe the orientation of the ions adjacent to the solid surfaces. The results suggested that the number density and electron density profiles show the oscillation near the interfaces. The period of the oscillations is consistent with predicted (anion + cation) dimensions, and the interfacial structure propagating into the bulk liquid depends on the structures of the ions. Recent atomic force microscope (AFM) experiments for both AILs and PILs next to a mica surface found the ILs are laterally inhomogeneous and reminiscent of self-assembled surfactant aggregates (AIL [C2C1im][NTf2] and PIL [EtNH3][NO3] in Figure 21).464,465 For the [C2C1im][NTf2], the peaks of the AFM tip approaching the mica surface suggest that the [C2C1im]+ layer adsorbed to the substrate is about 0.4 nm thick, and vibrational sum frequency spectroscopy (VSFS) studies suggest the cation is adsorbed with the imidazolium ring slightly tilted toward the silica surface (between 16° and 32° from surface normal). MD simulations464 provided an agreeable conclusion about the structural distribution at the interface, as shown in Figure 21B. The rings of imidazolium cations tend to be perpendicular to the mica surface, but the ethyl tails of cations tend to parallel with the mica surface, whereas with increasing alkyl chain length, the tendency is changed and the alkyl tails of cations start to be oriented perpendicular to the surface,457 which 6659

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

dynamics. The lateral structure of IL EDL at a charged Au substrate was investigated by neutron reflectivity.486 Little difference is shown in the interfacial structure at potentials more positive than, more negative than, and close to the point of zero charge. Theoretical models of the IL EDL have developed rapidly to support experimental insights. However, the complexity of the IL EDL makes many of the basic assumptions to mean-field models for the aqueous electrolytes unable to be applied. Typically, the IL ions are large with significant charge delocalization, and thus cannot be modeled as point charges. In addition, the ionic composition means that the concentration of ions near the electrode interface simply cannot differ greatly from the bulk, so Debye lengths are smaller than 1 nm. Nevertheless, a simplest mean field theory of EDL for ILs at a flat interface was suggested by Kornyshev in his feature articles.487,488 The capacitance expression through the overall potential drop across the double layer, written for compactness in dimensionless units, u = eU/ kBT, reads

Figure 22. Comparison of AMF observation and MD simulation for a biased [C2C1im][NTf2] IL near an ordered pyrolytic graphite surface. In part D the anion is denoted Tf2N− and the cation is denoted Emim+.496 (A and B) Experimental force separation curves under bias as 2D histogram and MD simulation ion positions for positive and negative potentials. (C) MD simulation ion−electrode interaction at PZC and 1 V bias applied to the electrode. (D) Orientation of the ion layer for −1 and 1 V, respectively. The white line denotes the position of the electrode surface. Reproduced with permission from ref 496. Copyright 2013 American Chemical Society.

( 2u ) ( 2u ) ln⎡⎣1 + 2γ sinh2( 2u )⎤⎦

⎛u⎞ 1 C = CD cosh⎜ ⎟ ⎝ 2 ⎠ 1 + 2γ sinh2

2γ sinh

(34)

The case when γ = 1 is

(u)

2 sinh 2 ⎛u⎞ C = CD cosh⎜ ⎟ ⎝ 2 ⎠ ln[cosh u] cosh u

Here CD =

ε* 4π L D

λD =

1 , 8πLBc0

and LB =

(35) e2 ε * kBT

structure of the EDL formed at a biased electrode can be exhibited. In many cases, simulations are sufficiently advanced such that they often uncover important structural features of the EDL prior to experimental measurement.497,498 They all predict many experimentally testable parameters, which are used to further refine the structural models and so deepen our understanding of the EDL.

are the Debye

capacitance, length, and Bjerrum length. The first two factors in eq 31 determine the Gouy−Chapman capacitance, but the two other factors convert the result into something very different. The second factor in eq 32 makes the voltage dependence of capacitance opposite to the Gouy−Chapman one: the curve takes a bell-shape instead of the U-shape. The behavior of the more general case depends on the value γ. For γ > 1/3, it has a bell-like shape, whereas, for γ < 1/3, it has a double hump shape, called a camel-shape.487 However, as a model for the IL EDL, the mean field theories will not likely continue to have considerable potential in electrochemistry (e.g., capacitors,471,489,490 electrodeposition,147 electrowetting491,492) due to a number of approximations behind them. MD simulations have been developed to trace the evolution of ion structure in the first and subsequent layers with applied bias.493−495 Black et al.496 observed reconfiguration under applied bias and the orientational transitions in the highly ordered pyrolytic graphite surface measured by AFM combined with MD simulation. MD simulation reproduced the results of AFM measurements at either unbiased or biased electrodes, shown as the comparison between the force separation curves of AFM and the ion density of MD simulation in Figure 22A and B. Once a bias is applied to the electrode, Coulombic interaction increases to dominate the interaction between ions and electrode (Figure 22C), which changes the structure of the EDL. Ions are almost oriented parallel to the electrode surface from the unbiased electrode (Figure 22 D). MD simulation provides the ion density and distance between ion layers and the charged surface that can be verified experimentally. By combining the information from AFM and the MD simulation, a comprehensive picture of the

5.3. Comparison of Structures to Other Solvents

As discussed above, ILs show some unusual self-assembled structures in the bulk phases and near air−liquid/solid−liquid interfaces. In these parts, ILs present a preferred separation, orientation of cations and anions, and repeat correlation length (d) between polar and apolar domains (cf. Table 2). Essentially speaking, the formation of polar and apolar domains is dictated by the complicated ion−ion interactions, such as Coulomb force and H-bond, solvophobic, and dispersion forces. Specifically speaking, Coulomb force, as a strong interaction, enforces an electroneutral distribution of positive and negative charges, further leading to the formation of polar domains: the charged groups of cations selectively attract the charged anions, and uncharged groups are expelled from this domain. The localization of charges on the ions leads to strong charge−charge correlations to assist in the formation of a polar domain. Solvophobic interaction becomes extruded in strength with cation alkyl chain length. The effect is to promote formation of larger nonpolar domains and, consequently, more well-defined self-assembled nanostructures. In addition, polar and nonpolar domains may be controlled by the relative dimensions of the polar and nonpolar moieties defined as geometry packing, in which the effect of anions on nanostructure depends on the cation alkyl chain length with amphiphilic cations. It is likely that the polar−nonpolar interface is subject to an area-minimization 6660

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

Table 2. Comparison of IL Structures and Interactions with Other Solventsa

a

Repeat correlation length between domains in the bulk (d). Abbreviations used for interaction types are as follows: C = Coulomb force, D = dipole interaction, Re = repulsive (Born) overlap of electron clouds, HB = H-bonding interaction, VDW = van der Waals interactions, π−π = pi−pi interaction, disp. = dispersion interaction, P = packing geometry, solv. = solvophobic interaction, σ = molecular diameter. [CnC1im] = 1-alkyl-3methylimidazolium cation for alkyl chains length n, [EtNH3] = ethylammonium cation. Note a diversity of structures exists in most liquid types; here, the typical example of solvent structure/key interactions/d are compared.

constraint, although the “tension” of this interface must be extremely low. The H-bonding interaction in the majority of ILs has been discussed above, but it is quite different from that in water. In some respects, H-bonds do not control self-assembly, but form between adjacent donor and acceptor sites in the nanoaggregates. Sometimes H-bonds form though accommodating between ions, but they are not the principal driving force for the aggregate structures. Dispersion forces are usually non-negligible and were evaluated to be 8−15% of the total interaction energy.170,288 Two factors lead to the main contributions from dispersion forces: (1) the presence of the electron density-rich imidazolium ring and (2) bulkier electron-rich anions with multiple interaction sites (e.g., for up to 59 kJ/mol in [NTf2]based ion pairs). The majority of molecular liquids show a favored separation of nearest neighbor molecules.499 The simple atomic liquids, such

as Hg, Ar, and Xe, are well-described by the VDW model of hard spheres with order created over a few diameters (cf. σ in Table 2) oscillating between regions of enhanced and depleted density.500 The diatomic (H2,501 N2,502 O2,502 halogen503) and triatomic (sc-CO2504) liquids, though being VDW model, are nonspherical, and electron densities are inhomogeneous. Classic molecular liquids, such as water,505−508 methanol, methylamine,509 and ammonia,510 possess more well-defined structure via a well-defined H-bonding network, such as the water tetrahedral H-bonding network,511 and other solvents via a nontetrahedral H-bonding network.512 For some molecular liquids (benzene513 in Table 2), a preferred near-neighbor separation and orientation emerges in the bulk. A favored molecular arrangement is exhibited in the first coordination shell (3.5 to ∼7 Å) in the solvents. The molten salts (NaCl in Table 2) display isotropic charge ordering. The bulk structure shows an oscillation around each ion in solution, beginning with the 6661

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

Figure 23. Drag coefficient and aspect ratio (A) and relative errors for the correlation in ILs (B).533 Reproduced with permission from ref 533. Copyright 2010 Elsevier.

species of unlike charge in the first coordination shell, and then regions of successively alternating sign,514 which is qualitatively different from the van der Waals model of the hard-sphere model. The solid−IL interface structures show many similarities with that of the bulk phase. But near the solid surface, the segregation of polar and apolar groups is retained normal to and laterally along the interface, and termed “layered”, “lamellar-like”, or “bilayer” in place of “solvation layers” to describe (1) normal structure and (2) worm-like, micelle-like, etc. for lateral structure (cf. Table 2). The solid−IL interface also exhibits usual features of simple, molecular, and self-assembled liquid systems, and the ions are highly organized in both normal and lateral directions.523−526 Interestingly, the associated push-through forces are typically larger in ILs than for solvation layers of most other solvent types. This indicates that the IL interfacial structure is much more pronounced than those of other solvents, similar to the bulk phase (cf., Table 2) due to strong cohesive interactions. From the comparison, it can be seen that IL structures show quite a difference from the molecular solvents and traditional molten salts, and the interactions result in at least an order of magnitude larger self-assembled structures, and IL thermodynamic properties (enthalpy, entropy, internal energy) are closely related to ion arrangements, which means that many of the models used to simulate traditional solvent organization and behavior527,528 are not fully established.

properties of gas−liquid systems (density, viscosity, surface tension, etc.). Single bubble behavior was investigated in three ILs ([C4C1im][BF4], [C8C1im][BF4], and [C8C1im][PF6]).533 The bubble velocity and bubble shape in the ILs were systematically measured by using a high-speed camera. Two new empirical correlations (eq 31 and eq 32) were proposed and gave satisfactory results in predicting the drag coefficient and aspect ratio of bubbles rising in the ILs (Figure 23 A), showing a small relative error (Figure 23 B). Zhang et al.534 investigated water effects on single CO2 bubbles in IL solutions. They found that the viscosity plays a more important role in determining the bubble behavior than the surface tension. A small amount of water shows a significant influence on the bubble behavior. The emprirical correlations (eq 33 and eq 34) were established to predict the bubble diameter and terminal velocity. Kaji et al.530 investigated different bubble shapes and the coalescences in a small-diameter bubble column by using IL [C2C1im][EtSO4] as a continuous phase. These bubbles are divided into three types depending on their diameter and shape: Taylor bubbles or spherical-cap bubbles, intermediate spherical bubbles of 0.5−2.5 mm, and spherical bubbles with diameters of ≤250 μm. This result showed that coalescence between Taylor and spherical-cap bubbles easily occurred (Figure 24). However, the coalescence between spherical-cap

6. UNIT SCALE DESIGN BASED ON ILS MEDIA Although the self-assembly structures are similar with traditional surfactant, the length-scales are at least an order of magnitude smaller than in the case with the traditional surfactant, which means that many concepts used in simulating traditional surfactant organization cannot be employed for predicting IL thermodynamic properties (enthalpy, entropy, internal energy). Moreover, the dynamics of ILs are not the same as those of molecular solvents, which suggests that the bulk phase nanostructure of ILs exists over much longer lifetimes than observed in other liquids. From the chemical engineering viewpoint, the reactors based on IL media should be redesigned. Related studies regarding the hydrodynamics of bubble columns529−531 and the mass transfer and heat transfer532 were performed in ILs or IL solution.

Figure 24. Sequence of the coalescence of two spherical-cap bubbles.530 Reproduced with permission from ref 530. Copyright 2009 American Chemical Society.

and spherical bubbles is rarely seen. The IL behaves as a noncoalescing liquid with a high viscosity. −0.849 ⎧ Mo 0.020 0.5 ≤ Re ≤ 5 ⎪ 22.73Re CD = ⎨ ⎪ −0.636 Mo 0.046 5 ≤ Re ≤ 50 ⎩ 20.08Re

6.1. Gas−Liquid Two-Phase Flow Hydrodynamics

The hydrodynamics of gas−liquid two-phase flow (bubble size, bubble deformation, bubble velocity, gas holdup, etc.) significantly influence the performance of the gas−liquid reactors unit. The bubble behavior in ILs is mainly controlled by the

E= 6662

1 1 + 0.0187(Re·Eo)0.67

(36)

(37) DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

⎡ d 2(ρ − ρ )2 ⎤1/3 eq l g ⎥ = 0.355Mo−0.172Eo 0.548 VT ⎢ ⎢⎣ ⎥⎦ σμl

(38)

⎡ D σ ⎤−1/3 = 2.675FroMo 0.016 deq⎢ o ⎥ ⎢⎣ ρl g ⎥⎦

(39)

In the correlations above, the dimensionless numbers are defined as follows: Eo =

Fro =

Mo =

Re =

(ρl − ρg )gdb2 μl

Figure 25. Pictures of the bubble column operated at superficial gas velocity, uG = 0.02 m/s (N2, p = 0.1 MPa, T = 473 K) for three liquids. Dibenzyltoluene (DBT): homogeneous regime; polydimethylsiloxane (X-BF) and IL [C4C1C1im][NTf2]: heterogeneous regime.535 Reproduced with permission from ref 535. Copyright 2016 Elsevier.

(40)

Ug 2 dog

(41)

(ρl − ρg )gμl 4 σ 3ρl 2

decreases with an increasing difference between particle and liquid density when adding some solids into the liquid phase. Zhang et al.180,536 proposed a new drag force model (eq 42) considering the unique properties of IL solvents to understand the flow fields of bubbles.

(42)

(ρl − ρg )dbVT μl

∂ (ρv ⃗) + ∇(ρvv⃗ ⃗) ∂t

(43)

The bubble size distribution, Sauter diameter, interfacial area, and total gas holdup are systematically investigated in an IL capturing CO2 system.531 The CO2 bubble diameter is smaller than that of N2 because the physical properties of the gas are different. The bubble diameter and gas holdup are mainly influenced by temperature and gas superficial velocity. A new Sauter diameter correlation531 (eq 41) is proposed for the prediction of in gas−ionic liquid systems with a relative error less than 6%.

= −∇P + ∇(μ(∇·v ⃗ + ∇·v ⃗T )) + ρg ⃗ +

2fρg v |⃗ v |⃗ db

(45)

where f stands for the gas−liquid friction factor

{

0.687 Re ≤ 1000 . The velocity fields and and f = 1 + 0.15Re 0.0183Re Re > 1000 pressure fields were simulated by CFD combined with an improved VOF method. The bubble behavior (velocity, deformation, equivalent diameter, etc.) was better predicted. For the stirred tank, Ouyang et al.537 adopted a new drag coefficient model and population balance model (PBM) to describe the total gas holdup, bubble size distribution, and local Sauter diameter. The predicted results for the total gas holdup and local Sauter diameter agree well with the experimental data. They also investigated the effects of the gas flow rate, agitation speed, and temperature on the hydrodynamics of gas-ILs (local gas holdup, bubble size distribution, and interfacial area). The empirical correlations have better predictions for ILs than the traditional models for molecular solvents. But most of the properties in the models come from the experimental measurement; some are not precise. It is necessary to correlate the ionic fluid from the molecular scale. Recently, Gabl et al.538 found that many properties of IL shows a linear dependence on the inverse box length x = 1/L by MD simulations, provided the system is sufficiently large, and at least 500 ion pairs were needed for the correct handling of dynamics. In the case of translational diffusion coefficients, this linear relation can be rationalized in hydrodynamic terms. The respective formula is essentially determined by viscosity and the inverse box length. The correlation from MD simulation to hydrodynamics can be mediated by the interparticle diffusion tensor shown in eq 43.

0.278 ⎛ ρ ⎞0.103⎛ σ 3ρ ⎞−0.034 d32 ⎛ H ⎞0.175⎛ Ug ⎞ l ⎜ ⎟ ⎜⎜ ⎟⎟ = 33.67⎜⎜ l ⎟⎟ ⎜⎜ 4l ⎟⎟ ⎝ D ⎠ ⎝ gD ⎠ do ⎝ ρg ⎠ ⎝ gμl ⎠

(44)

ILs usually have high viscosity, which leads to high flow resistance. But the viscosity decreases dramatically with increasing temperatures; meanwhile, the surface tension decreases slowly. Goetz et al.535 investigated the hydrodynamic behaviors of gas−liquid and gas−liquid−solid at higher temperatures in a slurry bubble column reactor. By elevating the temperatures, the liquid viscosity drops and, consequently, the regime transition shifts to higher gas velocities. Compared with the bubble behaviors at low temperature, the behavior of ILs at elevated temperature is ruled not only by their viscosity but also by surface tension. This result extends the findings530,534 which show that viscosity is the predominant factor of ILs at low temperatures. Figure 25 shows the bubble population in dibenzyltoluene (DBT), polydimethylsiloxane (X-BF), and IL [C4C1C1im][NTf2]. Small bubbles can be observed for the oil XBF and DBT, but large bubbles and small gas holdup were found in homogeneous IL (IL, Figure 25) when compared to the heat transfer oil DBT. Surprisingly, the oil DBT with nearly the same surface tension as the IL leads to smaller gas bubbles than the IL, and in DBT the regime transition is shifted to higher gas velocities. Small concentrations of fine particles can reduce bubble coalescence and subsequently stabilize the homogeneous regime. The primary bubble size increases and the gas holdup

Dijλδ =

∫0



⟨viγ (0)vjδ(t )⟩dt = Dγδ ( ri ⃗ − rj⃗)

(46)

Here i and j refer to a pair of molecular centers ri⃗ ,rj⃗ , while the superscripts γδ denote the Cartesian tensor components. 6663

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

In hydrodynamics, the velocity flow field u(⃗ r,t) replaces the molecular velocities vi⃗ = p⃗i/mi; the link between both worlds is given by uγ ( r ⃗ , t ) =

V N

piγ



mi

i

δ( r ⃗ − ri ⃗)

carbon monoxide. The influences of different process variables on kL in a batch reactor with a flat gas−liquid interface were assessed. A new correlation Sh = 10−2.64·Re1.07·Sc0.75 is proposed and used for estimating the mass transfer coefficient within the studied range of temperatures, CuCl concentrations, and stirring speeds. Hardacre et al.543 reported the heterogeneous catalytic reduction of phenyl acetylene to styrene over Pd/CaCO3 in IL and heptane. The mass transfer contributions in the system were disentangled for the first time. A modified additivity equation was proposed by introducing the diffusion of the gas, and an improved fit with the experimental data was acquired in the ILs. The lower reaction rates are attributed to the mass transfer of hydrogen within the liquid film. Sharma et al.544 reported volumetric gas−liquid mass transfer coefficient kLa for hydrogen and CO in IL [C4C1im][PF6]. Correlations were proposed using various dimensionless parameters to account for the critical impeller speed. They found that models of the form Sh = 10aReb(Fr − Frcr)cSc0.5 or Sh = 10a(Re − Re)bFrcSc0.5 are the best fit for describing the effects of different operating parameters on ́ et al.539,545 studied the oxygen mass kLa. Torres-Martinez transfer characteristics of a three-phase (gas−water−IL) stirred tank partitioning bioreactor based on the interaction of three operational factors (agitation, aeration rates, and IL volume fraction). The results showed that the IL had a negative impact on kLa with increasing the concentrations of the IL. The empirical correlation kLa = δ(Pg/v)αvsβ(1 − φ)γ was proposed for the three-phase system and showed good agreement between measured and predicted kLa values. For further understanding the mass transfer rules in a gas− liquid system, some models were developed. CFD simulation was carried out to model gas absorption. The coupled CFD-PBM model was first used to predict the mass transfer properties for CO2 absorption with ILs considering the unique properties of ILs.546 The drag force, which adopted a new empirical drag coefficient correlation of bubbles rising in ILs, is taken into account and implemented into FLUENT software by user definition. The mass transfer coefficients were estimated with Higbie’s penetration model. The predicted values of volumetric mass transfer coefficients are 3.4 × 10−5 s−1 at 300 K and 2.2 × 10−5 s−1 at 320 K, while the corresponding experimental values are 2.99 × 10−5 s−1 and 1.841 × 10−5 s−1. Zhang et al.536 developed a CFD method by adopting a drag force equation suitable for IL systems and considering the influence of CO2 concentration in the liquid phase on the viscosity of ILs. The CO2 concentration field in ILs is obtained by simulations for further understanding of the interfacial mass transfer. The results illustrate that the dissolved CO2 is distributed behind the bubble and forms a straight line in pure ILs as shown in Figure 26. The simulation results also showed that the influence of temperature on CO2 mass flux in pure ILs is different from that in aqueous ILs. The CO2 mass transfer is kinetic-controlled in pure ILs and thermodynamic-controlled in aqueous ILs, respectively. A model based on penetration theory was proposed to explore the mechanism of gas−liquid mass transfer of an ILs system in a rotating packed bed (RPB) for CO2 absorption.547 The results showed that high kLa was mainly influenced by the gas−liquid interfacial area and its frequent renewal on the packing surface. The predicted values of the model are in good agreement with the experimental data, with a maximum relative error of less than 15%.

(47)

For a periodic system, one can express the velocity flow field by a Fourier series uγ ( r ⃗ , t ) =



∑ ukγ(t )eik r ⃗ ,

ukγ (t ) =

k⃗

1 N

∑ i

piγ mi



e−ik ri ⃗ (48)

Assembling the relations, the diffusion tensor turns out for a periodic system. In the real space analog, it appears as a multiple of the Oseen tensor (eq 46), referring to a truly infinite system. The merit of the derivation above is the formulation of the viscosity in molecular terms, i.e., via the time correlation function of the Virial, thus preserving the link to the molecular world. D⃡ ( r ⃗) =

kBT Vη

1 − kk̂ ̂ ik ⃗r ⃗ e , k2

∑ k⃗

D⃡ ( r ⃗) =

kBT (1 + rr̂ )̂ 8πηr (49)

The self-diffusion coefficient can be easily obtained from the diffusion tensor by taking the trace, D(r ) =

kBT 1 6πη V

∑ k



4π ik ⃗r ⃗ e k2

(50)

The Ewald potential can be recognized by eq 48 when 2 2 introducing the convergence factor e−k /4α ϕEW ( r ⃗) =

1 V

∑ k⃗

4π −k2 /4α 2 ik ⃗r ⃗ e e k2

(51)

Among the diversity of all methods to calculate Coulomb interactions in a finite system, the Ewald method has been shown to have the minimum size dependence. Actually, one can correct the size dependence caused by the Ewald summation, as long as the box size suffices to reach the linear regime. 6.2. Mass Transfer

Mass transfer reduces with the high viscosity of ILs.539 However, in order to design and operate a new reactor, more attention is paid to the mass transfer rates of gas in some ILs. From this aspects, Liu et al.540 investigated dispersed IL in aqueous solution to enhance the CO2 absorption. They found that CO2 mass transfer rates are increased by dispersing the IL in water. The collaborative action of the hydrodynamic effect and the shuttle effect is proved to be the main mechanism for increasing the mass transfer driving force. For a carbon capture system mixing ILs with an amine, Zhang et al.532 measured the liquid-side mass transfer coefficients kL in a stirred cell reactor by the decreasing pressure method. They found that the kL in IL systems is influenced not only by the viscosity but also by the molecular structures of ILs. Then they proposed a new criterion “absorption degree’’, which combines both thermodynamics and mass transfer properties, to evaluate the performance of solvents for separating CO2 and CH4.541 The results showed that IL−NHD mixtures have high CO2/CH4 selectivity and fast kL for biogas separation. Based on the absorption degree, [C4C1im][NO3] + NHD mixtures were selected as good solvents for separating biogas. Zarca et al.542 evaluated the CO mass transfer rates in the IL imidazolium chlorocuprate (I) for recovery of 6664

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

Figure 26. CO2 mass fraction field in pure [C4C1im][BF4] IL (temperature 313 K, pressure 2 MPa).536 Reproduced with permission from ref 536. Copyright 2015 Elsevier.

Figure 27. Variation of Nusselt number with Reynolds number.556 Reproduced with permission from ref 556. Copyright 2015 Elsevier.

6.3. Heat Transfer

Thermal properties correlating with heat transfer applications are Tm, Tb, liquids range, Cp, heat of fusion, p, and thermal conductivity.548 ILs have good thermal stability and are selected for applications as high-temperature lubricants, thermal fluids, and media for heat transfer and short heat term storage.549 Thermal stabilities of many ILs were measured.550 Salgado et al.551 measured the thermal stability of five imidazolium-based ILs and found 1-butyl-2,3-dimethylimidazolium trifluoromethanesulfonate ([C4C1C1im][OTf]) is the most stable as a longterm heat transfer fluid. Compared with some commercial heat transfer fluids, substituted imidazolium and pyrrolidinium[NTf2] are found to be best suited as heat transfer fluids among their studied ILs.552 Due to the high viscosity and low thermal conductivity, the thermal entrance length of the IL is very large. Ding et al.553 investigated the convective heat transfer coefficient of [C4C1im][NTf2] and found it is much lower than that of distilled water under the laminar flow condition. The convective heat transfer data of [C4C1im][NTf2] were fit well to the convectional Shah’s equation under their experimental conditions. Paul et al.554,555 studied the natural convection heat transfer of [C4C1C1im][NTf2] and [C4C1im][NTf2] in a rectangular cavity with different aspect ratios. The Nusselt number of IL was observed to be higher than that of DI water. The thermal conductivity of IL is approximately 20% of DI water, and the natural convection heat transfer increases with aspect ratio increasing. Forced convective heat transfer of [C6C1im][BF4] was also investigated in smooth and microfin tubes with IL flowing in the inner tube and cooling water flowing in the annulus. The Nusselt number of the microfin tube is 5.4−11.3% higher than those of the smooth tube, and a correlation was developed to calculate the convective heat transfer coefficient for the microfin tube under laminar flow conditions (Figure 27).556 While considering pure ILs as heat transfer fluids, their operating temperatures do not usually exceed 230−250 °C. However, the operating temperature range can be substantially extended by introducing appropriate additives. Zhang et al.51 investigated the thermodynamic properties of nanoparticleenhanced ILs (NEILs) by dispersing graphene. The thermal conductivity, which increases by 15.2%−22.9% containing graphene of 0.06 wt %, is superior to that of the commercial heat transfer fluid Therminol VP-1. The initial decomposition temperature of the graphene-dispersed nanofluids is very close to 440.6 °C, indicating that it has good thermal stability and great

potential for use in medium- and high-temperature systems. The thermal conductivity of nanoparticle-enhanced IL decreased with increasing temperature for adding multiwalled carbon nanotubes (MWNTs) and single-walled carbon nanotubes (SWNTs), while for adding Al2O3 whiskers, an increase in thermal conductivity with increasing temperature is observed (Figure 28).557 NEILs heat-transfer fluids, which have the potential to deliver higher thermal conductivity than pure ILs without significantly increasing viscosity, are an alternative to conventional organic-based heat-transfer fluids. Paul et al.558−560 investigated the heat transfer performance of NEILs by dispersing spherical and whisker Al2O3 nanoparticles into the ILs at different concentrations (0.5, 1.0, and 2.5 wt %). For [C4C1im][NTf2] with 0.5 wt % loading of Al2O3 nanoparticles, the heat transfer coefficient was found to be enhanced by 13% compared with that of the base IL.558 By dispersing 0.5% Al2O3 nanoparticles (spherical and whiskers) in [N4111]][NTf2] IL, the thermal conductivity of NEILs enhanced by ∼5% and heat capacity enhanced by ∼20% compared to the base IL, which also gives 15% enhancement in heat transfer performance.559 The heat capacity of NEILs enhanced by 62% and 45% for 2.5 wt % Al2O3 nanoparticles loading of spherical and whiskers nanoparticles in N-butyl-N-methylpyrrolidinium IL, respectively. Heat transfer was observed to deteriorate for the NEILs compared to the base IL. The observed degradation of the heat transfer performance of the NEILs is not only caused by the change of effective thermophysical properties, but also caused by interaction of particle−fluid and clustering of nanoparticles.560 6.4. Reactor Design

For a chemical process, transport and reaction intensifications are critical for reactor design. Despite the considerable developments made on IL catalytic systems from the molecular to unit levels, there was still a delay in exploring appropriate reactors in IL-based media. Many reactions were carried out in the batch reactors that were designed by experimental reaction conditions (e.g., temperature, pressure, and catalyst lifetime). Recently, continuous flow reactors were also investigated by using IL as solvent and catalyst. A bubbling bed reactor was designed for synthesis of ethylene carbonate (EC) from ethylene oxide (EO) and CO2 by Zhang et al.,561 in which EO was initially pressurized with nitrogen to stay at the bed bottom as a liquid form, which was fed with ILs (5.6 wt % of EO) (Figure 29). CO2 was continuously injected and bubbled through the liquid 6665

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

Figure 28. Thermal conductivity changes with temperature in Al2O3-containing (A) and carbon-containing (B) NEILs.557 Reproduced with permission from ref 557. Copyright 2013 American Chemical Society.

For the gas-phase dimerization/isomerization of ethene in the presence of the cationic Ni-SIL catalysts, the complexes were dissolved in the IL [C2C1im][FAP] (Figure 31) and dispersed as a thin film over a silica gel 100 support. These SILs showed high activity in the first phase of the reaction. Under these conditions, the combination of high catalytic activity and strong exothermicity of the dimerization reaction still resulted in hot spot formation inside the fixed catalyst bed (Figure 31). Reaction intensification can be enhanced by the availability of a reliable and highly efficient microreactor with a large surface-tovolume ratio and short transport path channels to enhance heat and mass transfer dramatically.567,568 Many reactions have been carried out in a microreactor, such as hydrogenation,569,570 oxidation,571 Fischer−Tropsch synthesis,572 and carbonylation.573 IL catalytic microreactors were designed574 to overcome the shortcomings of the fixed bed with low transport and reaction efficiency. A microreactor setup was designed for the cycloaddition of propylene oxide (PO) and CO2 by IL catalyst in Figure 32.575 The microreactor includes a micromixer and a spiral capillary reactor. The micromixer is fabricated in the stainless steel plate, and the gas distributor microchannels are located. The channels are rectangular with a 600 μm (width) × 300 μm (depth) cross section. The internal channel volume is about 0.31 mL. The results indicate that the process can be intensified by increasing the reaction temperature, operating pressure, and catalyst concentration in the microreactor. The mass transfer performance of CO2 and PO is improved in the microreactor compared to conventional reactors. Besides designing new reactors, the IL-self-as solvents can enhance transport and reaction intensifications. The typical reactions are the liquid−liquid biphasic catalytic reactions using

mixture of EO and ILs. Such reactor promotes the homogeneous reaction and increases efficient EO conversion.

Figure 29. Schematic representation of a bubbling bed reactor for EC synthesis from EO and CO2 in ILs.

In the case of heterogeneous reactions, the fixed bed reactor was adopted.561,562 For the above reaction, CO2 and EO were fed to the bottom of the reactor with a flow ratio of 2 to 1. On the top of this reactor, only the EC product (nearly 99.5%) flowed out and was transferred to the next column for alcoholysis (Figure 30A). In the fixed bed reactor, the supported IL (SIL) catalysts were usually used and filled.563−566 Figure 30B shows the preparation of a coconut shell activated carbon tethered SIL catalyst for continuous cycloaddition of CO2 to epichlorohydrin,562 and in these SIL phase systems, a thin film of IL containing the homogeneous catalyst is immobilized on the surface of a high-area, porous support material.

Figure 30. Schematic illustration of SIL phase catalysts in a fixed-bed reactor (A)566 and preparation of coconut shell activated carbon tethered SIL (B).562 Reproduced with permission from refs 562, Copyright 2015 Elsevier, and ref 566, Copyright 2006 Wiley-VCH Verlag GmbH & Co.KGaA. 6666

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

Figure 31. Schematic representation of the applied Ni-based SIL catalyst systems and the reaction/deactivation front moving through the fixed bed containing the SIL catalyst.563 The temperature at the top of the bed (Ttop) and at the bottom of the fixed bed (Tbottom) were measured over time using two thermocouples. At the beginning of the reaction, a rise in temperature to Ttop = 59 °C was observed, while the temperature at the bottom increased to only 42 °C. After that, the temperature at the top rapidly decreased while the temperature at the bottom end slowly increased. With concomitant approach of the maximum temperature at the bottom of Tbottom = 56 °C, a drop in conversion occurred. The observed course of temperature indicates deactivation of the catalyst by a migrating reaction−deactivation front. This front obviously first develops at the top end of the fixed bed reactor where ethene comes in contact with the catalyst bed. Reproduced with permission from ref 563. Copyright 2014 Royal Society of Chemistry.

experiments, and this ligand system provided an n/iso ratio of up to 65. The catalytic efficiency depends on the type of metal complex dissolved in the IL. The metal catalyst in the IL is catalytically active by itself or after self-activation. ILs with weakly coordinating, inert anions (e.g., [NTf2]−, [BF4]−, [SbF6]−, or [PF3(C2F5)3]−) and cations are preferred. The role of the IL is to provide a moderately polar medium that does not compete with the substrate for free coordination sites. Many ILs offer a special combination of polarity and nucleophilicity that provides a superior reaction environment for many electrophilic and cationic transition metal complexes. In addition, chemical reaction between the IL and the catalyst precursor forms the active sites in situ. It is also essential to understand and control exactly the reaction that forms the active catalyst in the IL to make sure that 100% of the precious catalyst precursor ends up in the most active, most selective, and most stable catalyst complex.

Figure 32. Schematic of the microreactor setup for the cycloaddition of PO and CO2.575 A: CO2 cylinder; B: liquid tank (including PO and HETBAB); C: filter; D: pressure regulating valve; E: series II digital pump; F: gas mass flow controller; G: liquid damper; H: one-way valve; I: micromixer; J: spiral capillary in oil bath; K: gas−liquid separator; L: back pressure regulating valve; M: needle valve. Reproduced with permission from ref 575. Copyright 2013 Royal Society of Chemistry.

ILs to dissolve the transition metal complex to intensify the reactions.576 In the process, the products are easily separated from the reactive ionic catalyst phase and the metal catalyst and ILs are recycled without applying any thermal stress to the ionic catalyst solution.577 The typical biphasic reactions are hydrogenation reaction,578−580 dimerization reaction,563,581 and hydroformylation reactions.71,582 In the hydrogenation of cyclohexene by RhCl(TPP)3 catalyst, it was found that the rates of hydrogenation in [C4C1im][SbF6] are almost five times faster than those in acetone, and the ILs can be reused repeatedly and the loss of rhodium leaching into the organic phase was below 0.02%.583 In the hydroformylation of 1-octene reaction, the imidazolium-based dicationic phenoxaphosphino-modified xantphos-type ligand was used in [C4C1im][PF6] IL.584 No loss in catalyst activity or selectivity was observed in seven recycling

7. PROCESS SCALE DESIGN FOR ILS APPLICATIONS Process design and optimization are indispensable in developing new industrial applications. The IL-containing process is usually unpredictable when lacking reliable thermodynamic properties, and the IL performance and screening ILs also need priori prediction. Figure 33 shows a multiscale process design scheme proposed by Zhang et al.,585 as the example of a 1,3-butadiene separation process with IL additive. Some key parameters are given by QM calculation to improve the absorbency and selectivity of the ILs. Then, thermodynamic data, the activity coefficients (γ), equilibrium vapor pressures (P), and Henry constants (x) were produced by a continuum solvation model (CSM) with COSMO-RS. To predict the system over a wide

Figure 33. Schematic diagram using multiscale methodology in process design. 6667

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

7.2. Environmental Assessment

range of concentration and temperature, the group contribution (GC) parameters in the UNIFAC model74,586 are fitted to data from COSMO-RS. Then, process simulation was performed by the commercial Aspen software. In following section, we exhibit the process optimization and simulation of the IL-based industrial level applications. First, thermodynamic properties prediction and environmental and life cycle assessments were introduced. Afterward, we demonstrated several process simulations for IL-based applications.

7.2.1. GD Analysis. The environmental impact should be considered in the chemical process design. The GD method,229 as an accurate and quantitative method, has been designed to assess the environmental impact of IL systems. Based on nine environmental impact categories, this method can evaluate the comprehensive environmental impact of a substance, stream, or process. Xu et al.609 tested three biogas upgradings for removing CO2 technologies, i.e., pressurized water scrubbing (PWS), monoethanolamine aqueous scrubbing (MAS), and ionic liquid scrubbing (ILS). The results showed that high purity CO2 product can be obtained by MAS and ILS methods. Based on GD analysis, the environmental aspect of the three technologies shows that ILS has the highest GD production value, while MAS and PWS produce serious environmental impacts. Yan et al.610 developed an IL green process for the separation of MAL with [C4C1im][BF4]. To make a good comparison, an azetropic distillation process was simulated. The GD analysis for these two processes found that the IL green process exhibits GD superiority. Zhang et al.611 chose GD as the environmental objective to obtain the best parameters to separate 1,3-butadiene using IL-containing acetonitrile. The separation operation is involved only in the process; that is, no material exchange exists between the system and its environment. Therefore, eq 49 could be expressed as

7.1. Thermodynamic Properties Prediction

Although experimental measurement is indispensably important, it is not always practical for screening the right ILs for the new systems. Starting from the molecular structures is fundamentally requisite. Many methods, including the adjustable parameters QM-based QSPR587 and empirical COSMO-RS,227,588 and based on the concepts of corresponding states principle,589 GC,590 equation of state (EOS),591 and artificial neural network (ANN),592 were developed. QSPR593 can predict the melting point,594 heat capacity,595 viscosity,590 density,596 surface tension,597 glass transition temperature,598 and toxicity.599 New QSPR models have been proposed using the linear and nonlinear algorithms to predict the toxicity of ILs. The results provide some theoretical guidance and reference to screen and identify ILs.599 Although the QSPR method reveals the quantitative correlation between structure and property of ILs, this method often needs QM calculation and is considerably time-consuming. The COSMO-RS227 can predict thermophysical properties of various systems quickly. The most important advantages of the model may be that it does not rely on the experimental data. This method was extended to predict solubility,600 activity coefficient,601 and Henry’s law constant.602 Zhang et al.603 employed the COSMO-RS method to predict the Henry’s law constant of CO2 in 408 ILs at 298.2 K and found that the ILs with the anion tris(pentafluoroethyl)trifluorophosphate ([FEP]) have high CO2 absorption capacities, which were agreeable with the next experimental measurement.74 GC is simple and capable of giving a reasonable accuracy when all the necessary group increments are obtained from the experimental data.590 The properties of ILs, such as density,604 electrical conductivity, and thermal conductivity, viscosity, electrical conductivity, and thermal conductivity,605 can be predicted by GC. Gardas et al.605 have successfully regressed the parameters in different equations for the calculation of viscosity, conductivity, and compressibility of ILs. Mousazadeh et al.589 related the surface tensions with the reduced temperatures using a linear regression based on the corresponding states principle. Ge and Hardacre606 predicted gas and liquid Cp using the Joback GC. However, the ideal gas Cp of ILs cannot be tested using experimental data; thus, the prediction places a level of uncertainty on the prediction of the liquid Cp. Zhang et al.99 proposed a fragment contribution-corresponding states (FC− CS) method to predict properties of ILs. 46 fragments were specially classified from more than 500 ILs, and the fragment increments were determined using the experimental density data. The density and surface tension were predicted by six and two different equations based on the calculated critical properties. The results showed that the FC−CS was reasonable with an average absolute relative deviation less than 4%. The FC−CS were also available for density prediction of the mixtures of ILs. The ANN method is a powerful tool for predicting the properties such as Tm,592 ρ,607 and thermal conductivity.608

ΔGD = ΔGDe , c =

∑ GDke ,out − ∑ GDke ,in k

k

(52)

where ΔGDe,c refers to the GD value of total energy consumption. 7.2.2. LCA Analysis. Life cycle assessment (LCA) is a method to quantitatively analyze the potential impact of products and processes over the whole life cycle, from the resource, manufacture, use, disposal, to recycle.612−614 Such an assessment considers several indicators; the global warming potential, the photochemical ozone creation potential, and the complicated results may confuse the decider. However, the difference between the LCA and GD methods is that the LCA method does not consider the selected indicators in one according to the GD formula; thus, it is simple and easy to understand. At present the LCA has been successfully applied to several chemical processes.609−611,615 Besides their negligible vapor pressures of ILs, many aspects, such as abiotic depletion and ozone layer depletion, should be considered. Zhang et al.616 used the LCA method to compare the environmental impact of [C4C1im][BF4] and the molecular solvents. The results showed that the process with the IL [C4C1im][BF4] is likely to have a larger environmental impact than the traditional methods. Kralisch and Stark et al.617 gave an energetic, environmental, and economic evaluation for [C4 C1 im][BF 4 ] synthesis and application. The toxicity, persistency in the environment, and bioaccumulation were selected to be the environmental impact evaluation criteria because of the limited data. The results showed that dimethylsulfate and frequently used fluorinated anions are highly toxic, and the fluorinated anions can persistently exist in the environment, making it important to develop alternatives. Zhang et al.616 estimated such data, such as energy and organic solvent consumption, and some of the data, such as the amount of solvent emitted. Aspen process simulation614 was performed to obtain some data for process integration. 6668

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

7.3. Objective Optimization Analysis

distributed. Both of the distribution indices for crossover and mutation operators are set as 20. The results of the optimization are shown in Figure 34A and B, which describes the solutions at generations 1 and 100, respectively. It is found that the solutions at generation 1 are randomly scattered, which seeks design and selection of optimal ILs for satisfying the demands of different aspects of performance.

The objective optimization method is a powerful tool to design and construct a cost-effective, environmentally friendly, and socially beneficial energy system and is widely applied to the design, operation, and supply chain optimization of chemical processes.618−620 Specifically, the nonlinear programming (NLP) model621 was proposed for the rigorous optimization of a hydrocarbon biorefinery that integrated the techno-economic analysis and LCA assessment by a multiobjective optimization framework.622 The method has been used to optimize the bioprocesses based on economic and energy criteria.623−625 Liu et al.626 proposed a LCA multiobjective biofuel supply chain model with multiconversion paths to provide a good guide for the design of biofuel supply chains in China. Gerber et al.627 developed a methodology for the sustainable design and synthesis of energy systems coupled with process integration, life cycle assessment, and industrial ecology. More recently, the studies focus on the design of optimal ILs for a variety of industrial applications. Karunanithi et al.628 presented a computer-aided tailor-made IL design approach on the basis of semiempirical QSPR models, which considered several promising IL candidates. A logical and systematic methodology for the simultaneous design of IL entrainers and azeotropic separation processes was introduced.629 It integrated a GC solubility parameter model, a GC property models, and a UNIFAC model to predict the minimum IL concentration needed to break the desired azeotrope. The results showed that the process using a designed IL enhanced the material and energy efficiency compared with the experimentally selected IL. A MILP model computational molecular design for using ILs within environmentally benign absorption refrigeration systems also was proposed by McLeese et al.630 The method was used to screen the ILs for enhancing CO2 capture.631,632 However, such a model is only a single objective optimization. Zhang et al. carried out a multiobjective optimization of the 1,3-butadiene separation process using the IL-containing acetonitrile (ACN).611 The optimization involves maximizing the GD, the purity, and the recovery of 1,3-butadiene, and minimization of the energy consumption. After giving importance to the effect of each variable on the objective functions, the industrial requirements, and the ease, the variables were chosen as the decision variables in the process: ST1, RRT1, ST2, RRT2, RRT8, and DFT10, as shown in the specification of the variables in Table 3. These decision variables can usually be varied in a narrow operability range, which is determined based on the design limitations and sensitivity analysis. The nondominated sorting genetic algorithm (NSGA-II), which uses elitism. And a crowded comparison operator is adopted to keep diversity without specifying any additional parameters, and generates solution that are both nondominated and uniformly

7.4. Process Simulations for IL Key Applications

7.4.1. CO2 Capture and Conversion. IL-based carbon capture and conversion have been investigated by process simulations. The aim is to tailor IL reactivity, self-assembly, and selectivity for CO2 capture and conversion; finally, the technologies can be scaled up and become economically feasible. It was reported33 that CO2 exhibits good solubility in the IL [C6C1im][PF6] with increasing 0.0881% mass after the absorption. Process simulations were performed for CO2 absorption in ILs633−637 and in traditional monoethanolamine (MEA) aqueous solvent.638−640 Figure 35 shows the flow sheets of the two process simulations. To achieve the desired CO2 purity, a cryogenic column was deployed in the IL flow sheet and another flash before the cryogenic column unit was added to the sheet to remove the water in the cryogenic column (Figure 35B). A compressor was also included in the sheet, helping to establish suitable operating conditions for the cryogenic column. The simulation results showed that the IL process has a higher removal of 99.9% and achieved a higher CO2 purity of 95.1%, and the annual capacity was 47000 t. In addition, the IL-based process has a reduced 16% energy consumption so that the investment of such a process will be 11% lower than that of the MEA-based process.637 Thus, the IL solvents brought potential promises in CO2 capturing technology with high absorption capacity and low energy consumption when compared to the traditional MEA technology. Some new ILs, such as tetraethylphosphonium benzimidazolide ([P2222][Bnim]), were synthesized to capture CO2 from the flue gas.641 This process can make good use of the heat generated by this phase change to reduce parasitic power consumption, and the simulation shows 55% lower than that of the MEA process in energy consumption. By using the ILs [C2C1im][Ac] and [C 4C 1im][Ac], Krupiczka et al.642 found the ILs have comparable CO2 absorption capacities with 15% MEA solution, but the cost of electricity is 28% lower than that of the MEA solution process in a packed bed column. The mathematical model636 based on Peng−Robinson (PR) equation of state analysis found that the energy requirement of the stripper used to flash the CO2 out of the IL is much less. The IL-based CO2 capture process is also environmentally friendly. The GD analysis finds that the GD value of ILS is the highest among the PWS, MAS, and ILS technologies;609 that is, the consumptions of ILS and PWS are about 50% lower than that of MAS, indicating ILS is an environmentally friendly process. Based on the COSMO-RS model, Liu et al.643 compared the IL process with the MEA process for CO2 capture, and found that the single-stage and multistage process alternatives would reduce the total energy consumption by 42.8% and 66.04%. The IL− amine hybrid solvents showed more effective aspects such as high solubility and selectivity in the CO2 capture process. Three ILs, 1-butyl-3-methylimidazolium dicyanamide ([C4C1im][DCA]), [C4C1im][BF4], and [C4pyr][BF4] mixing solutions with MEA,644 were used, and the [C4pyr][BF4]-MEA process can reduce 15% in regeneration heat duty, 12% in the whole energy

Table 3. Variables of Multiobjective Optimization and the Corresponding Bound Ranges for the IL-Containing 1,3Butadiene Separation Process Variables

Variable Specification

Lower Bound

Upper Bound

ST1 RRT1 ST2 RRT2 RRT8 DFT10

Solvent ratio of T1 Mass reflux ratio of T1 Solvent ratio of T2 Mass reflux ratio of T2 Mass reflux ratio of T8 Distillate flow of T10

8.3 1.5 0.12 1.2 2 2700 kg/h

9.9 2.7 0.18 2.8 6 3545.8 kg/h 6669

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

Figure 34. Evolution of the optimum solution for the IL process at generation 1 (A) and at generation 100 (B).

Figure 35. CO2 capture flow sheets using MEA (A) and ILs (B).634 (A) CO2 removal = 0.9349, flow (kg CO2 captured/h) = 101.55, CO2 purity = 0.9740, and heat consumption (kW) = 155. (B) CO2 removal = 0.9999, flow (kgCO2 captured/h) = 106.47, CO2 purity = 0.9514, and heat consumption (kW) = 45. Reproduced with permission from ref 634. Copyright 2015 American Chemical Society.

Figure 36. Schematic two-step reaction and technologies for cycloaddition of CO2 and EO to produce EC, EG, and DMC.

hydroxy, and carbonyl, or some hydrophilic group anions, such as halide, [BF4]−, [CF3COO]−, [C3F7COO]−, and [CF3SO3]−, were also designed. Many processes using these traditional catalysts were optimized in history. Union Carbide and Scientific Design (UCSD) used KI or tetraethylammonium bromide as catalyst.684 Asahi Chemical developed a process by a pretreating strongly basic anionic exchanger resin Dowex MSA1. Jefferson Chemical (now part of Texaco Chemical)685 optimized a process by adding 10−15% triethylamine and ethyl bromide into EC. OMEGA developed a process catalyzed by tetrabutylphosphoniumiodide and K2CO3 dissolved in EG in the presence of H2O.686,687 For a new IL catalytic system, Tian et al. performed a process simulation to optimize the flow sheet.688 Figure 37 shows the flow sheets of the simulations. EO and CO2 are added in the reactor R101 to produce EC; then the EC and remaining mixture of EO and CO2 entered into the tube reactor R102 for further usage. The mixture from R102 was then sent to the flash tank V101. Here, the light component, CO2, from the top of the tank returned to R101 for recycling through the compressor K101. The residue of EC and the additional catalyst was delivered into the top of the reaction distillation tower (C201) for the countercurrent contacting with methanol from the bottom (MeOH:EC ≈ 8:1 (in molar ratio)). The light components, DMC and DME, standing at the top of the tower were delivered

penalty, and 13.5% in capture cost when compared with the MEA process. In recent years, conversion of CO2 into EG has been paid much attention in developing greener processes. As shown in the reaction scheme in Figure 36, CO2 and EO are converted into EC in step 1 by IL catalysts; then the EC produces ethylene glycol (EG) in step 2 by hydration. In step 1, the ILs as catalyst exhibit good catalytic activity,645,646 that especially correlates with Hbonding interactions;647,648 thus, some new IL catalysts649,650 have been exploited, such as imidazolium [C4C1im], pyridinial n[C4Pyr] salts, and amino-functionalized and super base-derived PILs, which can promote cycloaddition reaction of EO and CO2 to EC.651,652 EC resulting from the coupling of EO with CO2 is generally regarded as the “vector” for the ethylglycol (EG) in industry. EG is an important organic compound,653,654 and produced by hydration of EC.655−658 Numerous catalysts have also been exploited for the step 2 reaction,, such as molecular sieves,659 metal-based systems,660−670 anion-exchange resins,671 functionalized chitosan,672,673 polymer systems,674 silica,675 alkali metals, 676 guaternary ammonium salts,677 and IL catalysts.678−682 Especially, the 1-butyl-4-azo-1-azoniabicyclo[2.2.2]octane hydroxide IL exhibits high activity and has as much as 95% conversion at 70 °C and 4 h−1 with 10 times excess methanol to EC.679 Functional ILs,680,683 such as carboxyl, amidogen, 6670

DOI: 10.1021/acs.chemrev.6b00776 Chem. Rev. 2017, 117, 6636−6695

Chemical Reviews

Review

Figure 37. Process simulation flow sheets for IL catalytic coproduction EG/DMC.78 (A) Modeling flow sheet of EG/DMC production (R101, carbonation reactor; V101, flash drum; C201, reactive distillation column; C202, azeotrope separation column; C203, dimethyl carbonate column; C204, ethylene glycol column). (B) The flow sheet has added an EO absorption column C101. The raw materials enter R101 and mix with CO2, which offered the main product EC and a small amount of byproducts, EG, and diethylglycol (DEG). After one-step flash operation, the gas enters C101 to remove the unreacted EO, and the resulting reaction mixture was then directly transferred into C201 for the alcoholysis with MeOH. (C) The flow sheet has an additional EC purification column C102 between two main reactors. The real product from step 1 contains not only the expected ester EC but also its hydrolysis products, EG and DEG. In consideration of an equilibrium reaction characteristic of the second alcoholysis, the presence of EG and DEG in the feed material will decrease the conversion. It is expected reasonably to improve the efficiency of alcoholysis by supplying the purified EC from its mixture with EG and DEG. (D) The flow sheet for removal of water from EO was added, and its purification was initiated going through the washing column C101. In column C101, the impurity gas was removed, while the aqueous material (90% EO) at the bottom was subsequently sent by the pump P101 to the upper section of the dehydration tower C102. After a multistage distillation, EO of 99% purity could be obtained at the top of C102, which was then sent to R101. The remaining sections were the same as (C). Reproduced with permission from ref 78. Copyright 2015 Royal Society of Chemistry.

Table 4. Effect of Temperature on Biomass Dissolution Process by ILs699 Biomass Solubility at given temperature in (wt/%) or (g/mol) Noa

ILs

40 °C

type

50 °C

60 °C

70 °C

e

80 °C

90 °C

MCC 11.5 12.5 13 15 NG 1 [C4C1im][OAc] 2 [C4C1im][HSCH2COO] MCC