Discrete and Continuum Analyses for Confinement Effects of an Ionic

Jan 3, 2019 - Discrete and Continuum Analyses for Confinement Effects of an Ionic Liquid on the EDL Structure and the Pressure Acting on the Wall...
0 downloads 0 Views 4MB Size
Subscriber access provided by UNIV OF LOUISIANA

C: Physical Processes in Nanomaterials and Nanostructures

Discrete and Continuum Analyses for Confinement Effects of an Ionic Liquid on the EDL Structure and the Pressure Acting on the Wall Yu Dong Yang, Gi Jong Moon, Jung Min Oh, and In Seok Kang J. Phys. Chem. C, Just Accepted Manuscript • DOI: 10.1021/acs.jpcc.8b09630 • Publication Date (Web): 03 Jan 2019 Downloaded from http://pubs.acs.org on January 10, 2019

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

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

Page 1 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(a) MD simulation snapshot without wall charges (red-anion, blue-cation). The simulation domain (gray) shrinked from the initial domain (dashed), and yellow dotted line shows the centerline of each nanoslit. (b) Electric potential from the electrode under the nanoslit & bulk condition. Black solid line shows the potential at the bulk condition while the red line shows the nanoslit condition (H = 50 Å). Transparent boxes show the electrode (blue), nanoslit region (red), and bulk region (green). 157x57mm (300 x 300 DPI)

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Comparison of the MD and continuum results in the bulk case (H = 20 d) (a) Ion concentration distribution near the electrode. (b) Concentrations of ions (bar graph) and charge (point graph) in each layer (c) Effective ion concentration (up) and electrostatic potential (down). 171x40mm (300 x 300 DPI)

ACS Paragon Plus Environment

Page 2 of 39

Page 3 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Interference of the EDL structure predicted by bulk EDL structure. (upper) EDL layer near the left and right wall (lower) Superposition of ‘Left’ and ‘Right’ EDLs and the MD simulation result at (a) H = 2.0 d (b) H = 2.5 d (c) H = 3.0 d. 160x92mm (300 x 300 DPI)

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(a) Pressure exerted on the wall varying the nanoslit width. The inset figures are for the enlarged pictures when 0.5 d ≤ H < 5.5 d. (b) Average ion density in a nanoslit. Inset graph shows the average charge density in a nanoslit. 156x54mm (300 x 300 DPI)

ACS Paragon Plus Environment

Page 4 of 39

Page 5 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Ion distribution in the when H = 2.0 ~ 5.0 d. (upper) Ion distribution snapshots at the last time-step. (lower) Ion distribution graphs. 160x146mm (300 x 300 DPI)

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Ion distributions inside the nanoslit region when H = 1.0, 1.5 and 2.0 d. Ion distribution snapshots in the xy-plane (upper) and in the xz-plane (middle) and ion distribution graphs (lower). Inset figures are for the enlarged ion distribution graphs. 157x79mm (300 x 300 DPI)

ACS Paragon Plus Environment

Page 6 of 39

Page 7 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

A snapshot of the MD simulation model for the H = 200 Å case. Orange boxes show the wall atoms for calculating the stress. (red-anion, blue-cation) 58x84mm (300 x 300 DPI)

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Ion distribution graph varying the analysis bin width when H = 5.0 d. 160x81mm (300 x 300 DPI)

ACS Paragon Plus Environment

Page 8 of 39

Page 9 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Ion distribution graph inside a negatively charged nanoslit. 160x87mm (300 x 300 DPI)

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Effects of the parameters used in the continuum model. (a) Steric effect. (b) Correlation effect. (c), (d) Individual results with concept figures.

ACS Paragon Plus Environment

Page 10 of 39

Page 11 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Interaction between the ions inside a nanoslit and the wall atoms (a) Interaction force between the wall atoms and ions inside a nanoslit and counter-ions in the first ion layer. (b) Interaction energy between the wall atoms and ions inside a nanoslit and ions in the first & second layer. 157x58mm (300 x 300 DPI)

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 39

Discrete and Continuum Analyses for Confinement Effects of an Ionic Liquid on the EDL Structure and the Pressure Acting on the Wall Yu Dong Yang†1, Gi Jong Moon†1, Jung Min Oh1,2*, In Seok Kang1*

1

Department of Chemical Engineering, Pohang University of Science and Technology, 77 Cheongam-Ro, Nam-Gu, Pohang, Gyeongbuk, 37673, South Korea

2 Center

for Soft and Living Matter, Institute for Basic Science (IBS), 50 UNIST-gil, Ulju-gun, Ulsan, 44919, South Korea

†These

authors contributed equally

*Corresponding Authors Dr. In Seok Kang E-mail: [email protected] Phone: +82-54-279-2273 Dr. Jung Min Oh E-mail: [email protected] Phone: +82-10-9459-1910

1

ACS Paragon Plus Environment

Page 13 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

 Abstract We use both MD and continuum models to analyze the EDL structure and pressure acting on the wall of a symmetrically valenced ionic liquid in a nanoconfinement. The nanoslit width varies from 1 layer of ion size to 20 layers thick, in which the BSK model may not be compatible. However, the comparison with the MD model shows that the continuum model can successfully predict the averaged charge density and decay of the pressure. The continuum model also shows the exact location of the 1st and 2nd layer of ions from the confinement wall when considering both steric and correlation effects. However, the continuum model does not account for the oscillatory behaviors of the ion density and the pressure due to the typical characteristics of ionic discreteness, and it slightly overestimates the counter ion concentrations compared with the MD analysis. In conclusion, the continuum model can be used to estimate the overall tendency of the pressure acting on the wall when ionic liquid is confined in nanometer thickness though the oscillatory behavior of the pressure, which is originated from the discreteness of ions, cannot be predicted exactly. It can provide useful information in designing nano-porous structures for various electrochemical applications.

2

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 39

1. Introduction In the past two decades, there has been increasing interest in room temperature ionic liquids (RTILs) due to their unique characteristics such as tunable ion design, low vapor pressure, and no necessity of a solvent to dissociate them into cations and anions. Especially, with their high charge density, RTILs offer a remarkable advantage as an electrolyte in various devices such as batteries and super-capacitors1-6. When an electrolyte is confined to a nanostructure, however, several physical phenomena far from the existing bulk can occur such as electric double layer (EDL) overlapping or electro-capillary phenomenon 7-11. A number of researchers have been investigating the physical phenomena inside the EDL12-16. Especially, the theoretical interpretation of a nano-confined electrolyte has received considerable attention in various fields which are extremely difficult or time-consuming to analyze in experiments. The atomistic simulations (e.g., molecular dynamics or Monte Carlo simulations) are the preferable tools to investigate an electrolyte in a nanoconfinement17-26. Several researchers have investigated discontinuous effects on a small length scale, such as the topology effect26-27, edge effect28-30, and confinement effect13, 31-35, which are far from the physical phenomena in the bulk. Since Chmiola et al. reported the anomalous increase in carbon capacitance at pore sizes less than 1 nm36, the confinement effect, i.e., size effect of the confinement, has come into the spotlight. Many studies have shown that the confinement effect on capacitance leads to a decaying oscillatory behavior of the capacitance or internal ion density with respect to the confinement width in slit-like pores. This tendency becomes prominent as the nanoslit width approaches the ion size12, 31-40, and researchers have attributed this oscillatory behavior to interference from the EDLs31,33. However, molecular dynamics (MD) simulations are too expensive in terms of the computation time, especially when a porous electrode with networked confinements has a multi-scale domain (from nm to sub-μm). Additionally, the analysis of dynamic processes such as ion diffusion or electrochemical interactions require a rather long time simulation (from s to ms), which is almost infeasible in MD simulation at least until now. The continuum simulations can be a good alternative, which could greatly reduce the computation time. Continuum models have been constantly developed starting from the basic concept of the EDL proposed by Helmholtz in the early 20th century41-46. The existence of a diffuse double layer 3

ACS Paragon Plus Environment

Page 15 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

due to the rapid thermal motion was considered in the Gouy-Chapman model, which is also known as the Poisson-Boltzmann (PB) equation47-48. Recently, the electric double layer model of ionic liquids under a solvent-free condition was further advanced by the Bazant, Storey and Kornyshev (BSK model)41, 49, which accounts for over-screening and crowding effects. Recently, our group has been focused on the osmotic pressure of an electrolyte in a nanoconfinement, which is important in physicochemical processes such as capillarity in pores50-51. Additionally, it can be directly related to the total stress acting on a wall, which could induce deformation of porous electrodes52. Several studies have reported that the pressure also changes non-monotonically as the confinement size changes50-51, 53-54. As the nanoslit width decreases from bulk, the osmotic pressure in the slit centerline increases when the slit size becomes comparable to the ion size, and it sharply decreases resulting in a negative (joining) pressure when the nanoslit width is very small. In this domain scale, the continuum hypothesis should be checked by comparing with a molecular-level analysis55-59. Because the osmotic pressure at the nanoslit centerline is exactly the same as the total pressure exerted on the electrode surface50-51, 60, the osmotic pressure tendency can be directly compared with the MD result which can directly calculate the pressure exerted on the wall. In this study, the continuum model and MD model were directly compared for the same problem domain. The pressure acting on the wall and the ion distribution in a nanoslit were mainly observed. A primitive ion model is the simplest model to consider the finite ion size. By simplifying the RTIL under solvent free conditions, we can focus on the ionic behavior and distribution, not on the solvent and ion shape effect. The model has been widely used to study physical properties such as capacitance or EDL structures24, 38-41, 59. While varying the nanoslit width, the tendency of the pressure acting on the wall (MD results) and the osmotic pressure at the nanoslit centerline (continuum results) were compared to check the applicable range of the continuum hypothesis. Different from other MD studies, which fixed the internal ion density, ion reservoirs were included in this simulation which can take into account the ion density change inside a nanoslit31-33, 38. The results part consists of three themed sections. In the first section, the continuum and MD models are compared in the bulk scale to verify both models. The BSK model well predicted the internal ion and charge distribution up to the second ionic layer with the steric 4

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 39

factor and correlation length calculated from the MD conditions49. In the second section, the interference from the EDLs is predicted to anticipate the confinement effect using the bulk EDL structure. By comparing the predicted EDL structure and simulation results, we inferred that there might be other effects besides the EDL overlapping as the nanoslit width changed. In the third section, the pressure acting on the wall and the ion distribution in a nanoslit is analyzed by varying the nanoslit width. The continuum and MD results are compared, and the difference is explained in the scope of the EDL overlapping effect and ionic discreteness due to the packing structure. By clarifying the validity of the MD and continuum model, several future theoretical studies can use the continuum model instead of the MD model, which can largely decrease the calculation cost. Also, the prediction of the force exerted onto the confinement wall can suggest how strong the confinement should be or the durability of the system.

2. Theoretical Methods 2.1. MD Simulation Methodology. The cations and anions were both treated as singly charged (z = ± 1) Lennard-Jones (LJ) spheres to simulate monovalent ions with the same size. The diameter of both ions (d) is taken to be 10 Å with a short range repulsive LJ potential17-18,

 d 12  d 6  uLJ  k BT        r   r  

(1)

(kB: Boltzmann constant; T: Temperature) with a cutoff length of 21/6 d. The polarizability of the ionic liquid is modeled by scaling the electrostatic interactions with the dielectric constant that equals 2.0. The method of considering the IL polarizability by screening the electrostatic force has been widely used in other studies17-19, 61-62. The wall atoms were also modeled as LJ spheres with a diameter of 0.1 d, which were much smaller than the ion size, to consider the nanoslits as flat surfaces reducing the surface roughness. The simulation domain consisted of two ion reservoirs and two nanoslits shown in Figure 1a. Two thousand ion pairs were randomly packed in each 200 × 350 × 50 Å cuboid reservoir. Each electrode is spaced apart by a distance of H between the two ion reservoirs, and the area of each electrode is 250 Å × 50 Å. In a nanoslit, each electrode is spaced apart by a 5

ACS Paragon Plus Environment

Page 17 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

distance of H through which ions could freely move. The upper nanoslit is centered at 250 Å while the other slit is centered at 100 Å which is enough spacing, considering the cutoff length. The nanoslit spacing H is changed from 10 Å to 200 Å (1 d to 20 d), and each atom of the electrode had an electron charge of 0.006 e; thus, the charge density across each electrode (σw) is ± 9.613 μC/cm2. To satisfy the overall electroneutrality of the system, the upper slit is negatively charged while the lower slit is charged with a positive charge. The case for H > 10 d is calculated with a different geometry, which is explained in the supporting information. For non-dimensional analysis, the length scale in the subsequent analysis is expressed in terms of the ion diameter (d = 10 Å).

Figure 1. (a) MD simulation snapshot without wall charges (red-anion, blue-cation). The simulation domain (gray) shrinked from the initial domain (dashed), and yellow dotted line shows the centerline of each nanoslit. (b) Electric potential from the electrode under the nanoslit & bulk condition. Black solid line shows the potential at the bulk condition while the red line shows the nanoslit condition (H = 50 Å). Transparent boxes show the electrode (blue), nanoslit region (red), and bulk region (green).

For conducting the MD simulations, we used the Packmol63 program to obtain the initial configuration of the ionic liquids, and the open source software LAMMPS (Large-scale Atomic/Molecular Massively Parallel Simulator)64 for the computation. After the MD computation, the VMD (Visual Molecular Dynamics)65 tool was used for the visualization. The long-range Coulombic interactions were computed by the particle-particle particle-mesh 6

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 39

(PPPM) technique66. A Nose-Hoover thermostat and barostat with a relaxation time constant of 0.1 ps was used to control the temperature and pressure67-68. Every initial configuration was equilibrated for 10 ns with neutral electrodes under an NPT ensemble where the external pressure equaled 1 atm. The dynamic simulations ran with charged wall atoms for 9 ns, equilibration for 7 ns and statistical analysis for 2 ns with a NVT ensemble. A time step of 1 fs was used, and the velocity was integrated using the Verlet algorithm. While varying the nanoslit width, the total ion number and electrode surface area were kept constant. The charge and number density distribution of the ionic liquid were calculated every 0.1 Å from the electrode wall with uniform bins along the perpendicular direction to the electrode. A bin thickness of 0.1 Å was the marginal value without losing accuracy, and the details are explained in the supporting information. The ion concentration was analyzed in the lower nanoslit (between the positively charged electrodes) for every step and averaged for 2 ns. The ion distribution between the two negatively charged electrodes (upper region) showed a symmetric result, and it is included in the supporting information. The pressure acting on the wall was calculated using the virial stress, which was the symmetric per-atom stress tensor. The virial stress acting on each wall atom in the middle region (150 × 50 Å) only was calculated to neglect the edge effect. By integrating the normal direction (y-direction) stresses and divide them by the wall area, the pressure could be obtained. The details for calculating the virial stress are explained in the supporting information.

2.2. Continuum Simulation Methodology. To compare the result of the MD simulation, a study based on the continuum approach should be done. In this paper, we adopted the model proposed by Bazant, Storey, and Kornyshev as a model equation for ionic liquids (BSK model)49 because ionic liquids have not only a steric effect but also a correlation effect due to their finite ion size. The meaning of the steric effect and the correlation effect are discussed in the supporting information.

 1  lc2 2   2  2nb ze

sinh  ze / k BT 

1   cosh  ze / k BT   1

(2)

This governing equation was solved with two different boundary conditions (bulk 7

ACS Paragon Plus Environment

Page 19 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

region and nanoslit region) shown in Figure 1b. In our previous research51, a fixed electric potential was applied. However, in this study, a fixed charge density condition is applied, which is the same value as the wall charge density in the MD simulation. From the charge conservation condition69, 70, nˆ  D  nˆ   1  lc2 2     w ,

(3)

and the correlation effect is neglected near the surface ( lc  0 at x  0 ). Therefore, the boundary condition (3) is divided into nˆ  

w 

    0 

w , 

nˆ    2    0     0   0 .

(4) (5)

This governing equation was solved with two different boundary conditions (bulk region and nanoslit region) shown in Figure 1b. In the case of the bulk boundary condition, the semi-infinite condition is used. Therefore, the electric potential, 1st derivative and 3rd derivative of the electric potential should be zero at the bulk end,

 ,  ,    0 as X   .

(6)

For practicality, X = 500 Å, which is large enough to be considered as bulk. In the case of the nanoslit, both electrodes are symmetrically charged, and the electric field at the center is 0. Additionally, the nanoslit boundary conditions can be represented as

   H / 2   0,    0      H   0 .

(7)

At first, we solved the governing equations using the bulk boundary conditions, and then, we got the electric potential curve and ion concentration distribution result. After that, we solved the same equations again but used the nanoslit boundary conditions. In the case of the nanoslit boundary conditions, the result depends on the nanoslit width H. Other simulation parameters are listed in Table 1. After equilibrating for 10 ns, the bulk concentration can be calculated from the MD result in the bulk case (H = 200 Å). The steric factor (γ) is the ratio of the bulk ion density to the maximum possible density, and it can 8

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 39

be calculated with the ion size and bulk concentration. We set γ = 2nb(π/6)(d/2)3/Φmax, for random close packing of the spheres at a volume fraction Φmax of 0.6349. The correlation length is set to be the same as the ion radius. The dielectric constant of the background medium equals 10.0 to represent the solvent-free ionic liquid condition, while the MD simulation considers the dielectric constant as 2 to screen the electrostatic interactions. This method was widely used for comparing the MD and continuum simulations under the bulk ionic liquid system and showed a well matched result49. Electric permittivity

  10 0

Bulk concentration

nb  0.65 mol/L

Steric factor

  0.65

Bulk region length

X  500Å

Correlation length

lc  5Å

Nanoslit width

H  10 ~ 200Å

Charge valence

z1

Charge density

 w   9.613  C/cm 2

< Table 1 > Simulation parameters used in continuum model The bulk surface pressure and the nanoslit surface pressure should be calculated prior to the osmotic pressure at the centerline by the following equations51. From the BSK model, the electric potential distribution is calculated, and then, we can get the first and second derivatives of the electric potential. Additionally, the osmotic pressure gradient can be represented as

 lc2  2 2  1      f         E  E        2 2 4    1   2   3



(8)



2 2 2 and we already know that  f   1  lc    from the Gauss law. By substituting this into

(8) we can get

 lc2 d       (4)       dx 2

(9)

The osmotic pressure at the bulk region can be expressed by integrating the osmotic pressure gradient. The bulk surface osmotic pressure 

A

is given as 9

ACS Paragon Plus Environment

Page 21 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

  x  Bulk   0 

 lc2  2         ,   0  Bulk   A 2 2

(10)

The bulk osmotic pressure  0 can be treated as the integration constant, and it is calculated from50

0 

 2n0 k BT  1  1    ln    1 ln      ln   2      

(11)

Without considering the steric factor, this result is the same as that calculated by van’t Hoff equation, but if the steric factor is considered, then, the bulk osmotic pressure slightly increases50, 71. Moreover, we can represent the nanoslit surface osmotic pressure 

B

by

adapting the same method and take the reference point as the centerline point (x = h).

  x  Nanoslit    h  

l2  2    c    ,   0  Nanoslit   B 2 2

(12)

To calculate the osmotic pressure at the nanoslit region, however, the electric correlation effect should be considered. In Eqn. (8),  1 and  2 are a function of the electric potential; thus, it will be constant,  1A   1B and  2A   2B , along the same electrode surface. The only difference between the bulk surface and the nanoslit surface    B

A

is  3B   3A , due to the

EDL overlapping. Therefore, the osmotic pressure at the nanoslit surface can be represented by

   B

A

2  lc2   2 

  4  n 2  

Nanoslit

  2   2   n 

2

   Bulk 

(13)

Now, the osmotic pressure at the nanoslit centerline can be obtained by using the osmotic pressure value at point B. From Eqn. (11), we have

 lc2  2   h               2 2 C

B

(14)

As discussed above, the osmotic pressure at the nanoslit centerline is the same as the normal stress exerted on the wall surface50-51, 60. After solving the governing equation, we 10

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 22 of 39

compared the result of the MD simulation and continuum approach by changing the nanoslit width H. By comparing the normal stress exerted on the wall, we can find the length scale where the continuum model is still valid.

3. Results and Discussion 3.1. Bulk Case Verification. First, the ion and charge distribution near the electrode surface were compared for the bulk case. Here, the bulk case represents a system in which a 1D flat electrode is placed in a semi-infinite medium. When the nanoslit width is wide enough (H = 20 d), the net charge density and electric potential at the slit centerline become 0; thus, this case can be considered as bulk.

Figure 2. Comparison of the MD and continuum results in the bulk case (H = 20 d) (a) Ion concentration distribution near the electrode. (b) Concentrations of ions (bar graph) and charge (point graph) in each layer (c) Effective ion concentration (up) and electrostatic potential (down).

When the steric factor and correlation length are fitted from the MD simulation conditions, Figure 2a shows that both cases could similarly predict the location of the first counter-ion layer (the first peak of the red lines from the wall) and the adjacent co-ion layer (the first peak of the blue lines). However, the ion distribution near the wall appears to have different shapes when the ion distribution is analyzed by each data point (ρion). The difference arises from the difference in the calculation method. In the case of the MD simulation, the ions are assumed to have charges and masses concentrated at the center point, while the ion size is defined by the van der Waals interaction. Therefore, the discrete ion concentration is calculated 11

ACS Paragon Plus Environment

Page 23 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

through the center of mass (CoM) distribution of ions. On the other hand, because the probability of finding an ion in each lattice is used to calculate the ion concentration in the continuum model, it shows the continuous ion distribution. Because the method of obtaining the ion concentration is different between the two models, it is not appropriate to compare the CoM distribution (MD) and probability of ion distribution (continuum) in the discrete perspective as in Figure 2a. One of the methods to compare both models in the continuous perspective is through the average of the ions constituting each EDL21, 49, and similar EDL structures are observed up to the second ionic layer in Figure 2b. The region inside a nanoslit is divided into layers with a thickness of 0.7 d to separate the EDL layer-by-layer. By integrating the pointwise concentration (ρion) in each layer, the ion concentration (ρint) and the difference between the cation and anion concentration, that is, the relative charge density (σint) in each layer, are obtained. It is meaningful that the continuum model can predict the first and second ion layer similar to the MD model because the scale of the total interaction between the wall and the internal ion is largely determined by the interaction with the first and second layers, as explained in the supporting information. As a result, the pressure acting on the wall (Pbulk) and the bulk ion density (ρbulk) are also similarly predicted in both models. The other method is to compare the effective ion concentration (ρeff), which is calculated from the CoM distribution. Assuming that the mass and charge of an ion are uniformly distributed in the spherical shape, the effective ion concentration can be calculated by taking into account the partial volume of the spherical ion inside a unit slab at each calculation point. Additionally, the electrostatic potential can be derived by solving the Poisson’s equation, and the both results are well matched. The detailed analysis methodologies for Figure 2b and 2c are explained in the supporting information. In summary, it is appropriate to compare the continuum and MD model through the steric factor and correlation length directly obtained by the MD simulation condition.

3.2. Prediction of the EDL Overlapping in Nanoslit Cases. As the nanoslit width decreases, EDLs from each electrode side may interfere with each other. The confinement effect in the nano-scale has a large effect on the capacitance oscillation, and it is directly related to the EDL 12

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 24 of 39

overlapping31, 33. Similar to the superposition of a wave, the ion density may increase or decrease due to the overlapping. When the ion layer coincides with the same charged ion layers, the internal ion density increases (constructive interference), while it decreases (destructive interference) when the ion layer coincides with the oppositely charged ion layers. The interference from the EDLs in a nanoslit are predicted using the EDL structure in the bulk case. The EDL layer near the left and right wall are obtained using the ion distribution in the bulk case, respectively (‘Left’ and ‘Right’ in Figure 3). By shifting the positions of the peaks, the EDL overlapping phenomenon is predicted in Figure 3. When H = 2.0 and 2.5 d, both EDLs seem to constructively interfere with each other. The maximum constructive superposition is observed in Figure 3b (H = 2.5 d), while the destructive superposition is shown in Figure 3c (H = 3.0 d). However, in the MD results, the increase in the co-ion density due to constructive superposition does not appear to vary significantly between H = 2.0 d and 2.5 d. Instead, the cation and anion densities were found to be maximum at H = 2.0 d in the MD results. The differences are clearly shown in later analysis with a detailed explanation. It suggests that the confinement effect cannot be predicted only from the EDL overlapping.

Figure 3. Interference of the EDL structure predicted by bulk EDL structure. (a-c) EDL layer near the left and right wall (d-f) Superposition of ‘Left’ and ‘Right’ EDLs and the MD 13

ACS Paragon Plus Environment

Page 25 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

simulation result at (a, d) H = 2.0 d (b, e) H = 2.5 d (c, f) H = 3.0 d.

The ion density change due to the EDL overlapping may also significantly affect the pressure, or normal stress acting on the wall. The idea came from the numerical analysis done by Aguilar-Pindea et al., decomposing the normal stress on the wall into a contact and an electrical component53. The normal stress was calculated from the sum of both components. The contact component means the force generated by the collision between the ions and the wall, which acts in the outward direction on the wall surface. On the other hand, the electrical component is the force due to the electrostatic interaction between the ions and the charged wall atoms, which is similar to the Maxwell stress affected by the electric field from the wall. Because the counter-ions mainly exist near the wall surface, the electrical component mainly acts to pull the wall inward. Because both components acting on the opposite directions tend to increase as the internal ion concentration increases, predicting which the force will increase more requires complicated calculations. Therefore, in the following section, the relationship between the pressure acting on the wall and the internal ion density change are analyzed using the MD and continuum simulation.

3.3. Confinement Effect on the Pressure Acting on the Wall / Internal Ion Density. The change of pressure acting on the wall with respect to the nanoslit width is shown in Figure 4a. Both the MD and continuum results show that the pressure acting on the wall converges to the bulk pressure as H increases as expected. Interestingly, the plot shows that a thickness of 5 ~ 6 d layer, which is too thin to be valid for the continuum postulate, is enough to get the bulk pressure. It implies the possibility that the continuum model can be used where the continuum postulate may be violated.

14

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 26 of 39

Figure 4. (a) Pressure exerted on the wall varying the nanoslit width. The inset figures are for the enlarged pictures when 0.5 d ≤ H < 5.5 d. (b) Average ion density in a nanoslit. Inset graph shows the average charge density in a nanoslit.

The pressure has the maximum value when H = 2.0 d, but when the nanoslit width becomes narrower, the pressure rapidly decreases resulting in a negative (joining) pressure. The pressure increase at H = 2.0 d and the negative pressure at the narrower confinement are also predicted by the continuum model. It may seem unrealistic for the negative pressure to appear, but it has been shown in several studies, and it was verified thermodynamically that liquids can exist in absolute negative pressure regimes72,73. However, only the MD model predicts the oscillatory behavior of the pressure acting on the wall at 2.0 d ≤ H < 5.0 d (transparent red region). The oscillation shows a gradual decay and eventually converges to the bulk pressure, similar to the capacitance or pressure tendency observed in other studies31, 33, 53. Figure 4b shows the average ion density and average charge density with respect to the nanoslit width. The average densities are obtained by integrating the pointwise ion density in a nanoslit region. The average charge density in the nanoconfinement can be fairly accurately predicted by the continuum model, but the ion density shows a slight difference. In the previous section, the ion density was expected to be maximum at H = 2.5 d because the EDLs from the opposite electrodes are supposed to be overlapped constructively. However, a higher density at H = 2.0 d was observed than at H = 2.5 d, which means there may be other effects besides EDL overlapping. Additionally, the continuum model seems to slightly overpredict the counter-ion density. These differences come from the discreteness due to the packing structure 15

ACS Paragon Plus Environment

Page 27 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

of the ion layers. Let us imagine that in a tightly packed slit width with n-layers of ions and the width slightly increases, but not enough to fill (n+1)-layers. Then, the ion layers become more loosely packed, which means the ion density and the pressure may be decreased. Moreover, if the slit width increases a little bit more to fill (n+1)-layers tightly, the charge density and the pressure will increase again due to the tight packing. This will repeat as the slit width increases. It leads to the oscillatory behaviors of ion density and the pressure according to the nanoslit width, and the amplitude of the oscillation decays as the number of ion layers increases. It is also notable that the decay is so fast that the thickness of 20 d is large enough to be neutral at the given condition. The oscillatory decay of the pressure is a typical characteristic of the discreteness of ions, which cannot be explained by the continuum model. It is because the lattice size is assumed not to change with respect to the slit width in the continuum model. However, the effective lattice size that the ion can occupy may change with respect to the slit width in the discrete (MD) model. Actually, this difference causes the difference in the prediction of the average ion density in the nanoslit, and we may say the continuum model makes an upper bound of the ion density. When 2.0 d ≤ H ≤ 5.0 d, the discrete effect greatly influences the ion distribution and the pressure acting on the wall because only a few ion layers are formed in a nanoslit. For a detailed analysis, an ion distribution snapshot at the last time-step in the MD results and the ion distribution graphs in a nanoslit are plotted in Figure 5. It directly shows the loosely packed ionic layers (H = 3.0, and 4.0 d) and tightly packed ionic layers (H = 2.0, 3.5 d). For tightly packed structures, the counter-ion peak next to the wall is higher than that of the loosely packed cases. Because the ions are closely packed to each other, there is not enough space where the ions can freely move, and it results in a higher ion density. From this result, we can confirm that the size of the available occupancy sites for ions changes as the packing structure changes in the MD simulations.

16

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 28 of 39

Figure 5. Ion distribution in the when H = 2.0 ~ 5.0 d. (upper) Ion distribution snapshots at the last time-step. (lower) Ion distribution graphs.

Because both electrodes are charged with the same sign, the total number of ion layers between the electrodes is observed to be odd. The case of H = 3.0 d is expected to have 4 ion layers, and the slit width is not large enough for having 5 ion layers. Therefore, it results in 3 layers with a loosely packing structure, with a large free space in which the ions can freely move. In the case of H = 2.0 and 3.5 d, however, the nanoslit widths are suitable (well-fitted) for 3 layers and 5 layers, respectively. The pressure and ion concentration increases when the number of ion layers increases so that the layers are tightly packed again. This discrete effect 17

ACS Paragon Plus Environment

Page 29 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

gradually becomes negligible as the number of ion layers increases, and thus, the pressure acting on the wall seems to converge to the bulk value after H = 5.0 d. In the cases of H < 2.0 d, however, the co-ion density abruptly decreases in a nanoslit. As discussed above, the normal stress can be decomposed into the electrical and contact components. As co-ions are pushed out of the confinement, the inward electrical component becomes more prominent. When the inward electrical component becomes larger than the outward contact component, the sum of two components results in an inward force (negative pressure) applied on the wall. The negative pressure applied to the wall, which is a result of a competition between the electrical and contact components, is similarly predicted in both the MD and continuum results. Figure 6 shows the ion distribution snapshot at the last time-step in the MD simulation and the ion distribution graphs for the H = 1.0, 1.5 and 2.0 d cases. For the snapshots, the spheres in the solid red color mean anions above the centerline while the transparent red ones depict the ions below the centerline. When H = 1.5 d, ions from one side hinders the ion layer formation on the other side due to the ion size. Thus, the average counter-ion density is low, and a weak pressure is applied on the wall. At this point, the MD simulation shows that no coion exists in a nanoslit while the continuum model does. In the MD simulation, the co-ion density in the confinement rapidly converges to zero as the confinement size decreases, whereas it gradually decreases in the continuum simulation. When H = 1.0 d, the co-ion is completely pushed out of the confinement in the continuum simulation. Despite both models showing a difference in predicting the decreasing tendency of the co-ion concentration, they similarly predict the scale of the negative pressure acting on the wall.

18

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 30 of 39

Figure 6. Ion distributions inside the nanoslit region when H = 1.0, 1.5 and 2.0 d. Ion distribution snapshots in the xy-plane (upper) and in the xz-plane (middle) and ion distribution graphs (lower). Inset figures are for the enlarged ion distribution graphs.

4. Conclusion We investigated the confinement effect of a solvent-free RTIL inside a nanoslit using the MD and continuum (BSK) model. The predicted internal ion density, charge density and pressure acting on the wall were analyzed to find the tendency with respect to the confinement size and the validity of the continuum model in such a small scale. When the confinement size is sufficiently large (H > 5.0 d), the pressure, average ion and charge density in a nanoslit are well matched in the MD and continuum simulations. Both simulations predicted the same pressure acting on the wall, average ion and charge density inside a nanoslit. Additionally, the EDL structures are well matched up to the second ionic layer in both the discrete and the continuous perspective. Thus, it was confirmed to be a great choice to use the continuum model using the steric factor and correlation factor obtained by the MD simulation condition. When the size of the nanoslit becomes comparable to the ion size, it is hard to predict the oscillatory behavior with the continuum model due to the ionic discreteness (whether the 19

ACS Paragon Plus Environment

Page 31 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

ion layer structure is packed tightly or loosely), and the continuum model slightly overpredicts the counter-ion density. However, to some extent, it can predict the overall tendency of the internal charge or ion density with respect to the slit width. In conclusion, the continuum model can still be useful to estimate the order of the pressure of an ionic liquid in a nanoconfinement even where the continuum postulate is not valid. In the cases of H < 2.0 d, however, the co-ion exclusion phenomenon greatly influences the pressure and EDL structure. The continuum model can predict the same scale of negative (joining) pressure as the MD model. However, it shows a slower decrease of the co-ion density as the confinement becomes narrower compared to the MD simulation. In other words, the continuum model can work better on this scale to predict the pressure acting on the wall if the co-ion density is neglected in a small confinement (using the counter-ion only assumption). This study can be a guideline for using the continuum model to analyze a nanoconfined liquid system, which is directly related to the calculation cost. Additionally, this study will give inspiration to researchers in various fields who are interested in the pressure acting on the wall for designing stable nanoporous electrode systems. However, we found no signatures for the emergence of long decay ranges of force as seen in the surface force experiments74. Further studies on the edge effect or solvent effect will help to understand more precisely about the validity of the model. In particular, the solvent effect has been found to affect the dynamics or capacitance of ionic liquids75-78. However, the scale of solvent effect may be controversial and further research is needed.

Supporting Information. Explanations of bulk case MD model, methods to calculate the pressure, verification of specific bin width, ion distribution near the negatively charged electrodes, the physical meaning of the parameters used in the BSK model and the wall-ion interaction.

ACKNOWLEDGMENT 20

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 32 of 39

This work was supported by Basic Science Research Program through the National Research Foundation of Korea (Grant Number: 2017R1D1A1B05035211) and by the BK 21 Plus program of Korea. The research at IBS Center for Soft and Living Matter was supported by IBS-R020-D1. REFERENCES

1. Armand, M.; Endres, F.; MacFarlane, D. R.; Ohno, H.; Scrosati, B. Ionic-Liquid Materials for the Electrochemical Challenges of the Future. Nat. Mater. 2009, 8, 621-629. 2. Haerens, K.; Van Deuren, S.; Matthijs, E.; Van Der Bruggen, B. Challenges for Recycling Ionic Liquids by Using Pressure Driven Membrane Processes. Green Chem. 2010, 12, 21822188. 3. Singh, M. P.; Singh, R. K.; Chandra, S. Ionic Liquids Confined in Porous Matrices: Physicochemical Properties and Applications. Prog. Mater. Sci. 2014, 64, 73-120. 4. Fedorov, M. V.; Kornyshev, A. A. Ionic Liquids at Electrified Interfaces. Chem. Rev. 2014, 114, 2978-3036. 5. Perkin, S. Ionic Liquids in Confined Geometries. Phys. Chem. Chem. Phys. 2012, 14, 50525062. 6. Hayes, R.; Warr, G. G.; Atkin, R. Structure and Nanostructure in Ionic Liquids. Chem. Rev. 2015, 115, 6357-6426. 7. Baldessari, F. Electrokinetics in Nanochannels. Part I. Electric Double Layer Overlap and Channel-to-Well Equilibrium. J. Colloid Interface Sci. 2008, 325, 526-538. 8. Bratko, D.; Daub, C. D.; Leung, K.; Luzar, A. Effect of Field Direction on Electrowetting in a Nanopore. J. Am. Chem. Soc. 2007, 129, 2504-2510. 9. Lu, W.; Han, A.; Kim, T.; Punyamurtula, V. K.; Chen, X.; Qiao, Y. Field-Responsive Ion Transport in Nanopores. Appl. Phys. Lett. 2009, 94, 023106. 10. Kong, X.; Qiao, Y. An Electrically Controllable Nanoporous Smart System. J. Appl. Phys. 2006, 99, 064313. 11. Joseph, S.; Aluru, N. R. Why Are Carbon Nanotubes Fast Transporters of Water? Nano Lett. 2008, 8, 452-458. 21

ACS Paragon Plus Environment

Page 33 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

12. Wu, P.; Huang, J.; Meunier, V.; Sumpter, B. G.; Qiao, R. Complex Capacitance Scaling in Ionic Liquids-Filled Nanopores. ACS Nano 2011, 5, 9044-9051. 13. Xing, L.; Vatamanu, J.; Borodin, O.; Bedrov, D. On the Atomistic Nature of Capacitance Enhancement Generated by Ionic Liquid Electrolyte Confined in Subnanometer Pores. J. Phys. Chem. Lett. 2013, 4, 132-140. 14. Osti, N. C.; Gallegos, A.; Dyatkin, B.; Wu, J.; Gogotsi, Y.; Mamontov, E. Mixed Ionic Liquid Improves Electrolyte Dynamics in Supercapacitors. J. Phys. Chem. C 2018, 122, 1047610481. 15. Bo, Z.; Li, C.; Yang, H.; Ostrikov, K.; Yan, J.; Cen, K. Design of Supercapacitor Electrodes Using Molecular Dynamics Simulations. Nano-Micro Lett. 2018, 10, No. 33. 16. Ghosal, S. Electrokinetic Flow and Ion Transport in Nanochannels. In Encyclopedia of MicrofluIdics and Nanofludics; Li, D., Ed.; Springer: New York, 2008; pp 496-506. 17. Fedorov, M. V.; Kornyshev, A. A. Ionic Liquid near a Charged Wall: Structure and Capacitance of Electrical Double Layer. J.Phys. Chem. B 2008, 112, 11868-11872. 18. Fedorov, M. V.; Kornyshev, A. A. Towards Understanding the Structure and Capacitance of Electrical Double Layer in Ionic Liquids. Electrochim. Acta 2008, 53, 6835-6840. 19. Kirchner, K.; Kirchner, T.; Ivaništšev, V.; Fedorov, M. V. Electrical Double Layer in Ionic Liquids: Structural Transitions from Multilayer to Monolayer Structure at the Interface. Electrochim. Acta 2013, 110, 762-771. 20. Docampo-Álvarez, B.; Gómez-González, V.; Montes-Campos, H.; Otero-Mato, J. M.; Méndez-Morales, T.; Cabeza, O.; Gallego, L. J.; Lynden-Bell, R. M.; Ivaništšev, V. B.; Fedorov, M. V.; Varela, L. M. Molecular Dynamics Simulation of the Behaviour of Water in Nano-Confined Ionic Liquid-Water Mixtures. J. Phys.: Condens. Matter 2016, 28, 464001. 21. Jiang, X.; Huang, J.; Zhao, H.; Sumpter, B. G.; Qiao, R. Dynamics of Electrical Double Layer Formation in Room-Temperature Ionic Liquids under Constant-Current Charging Conditions. J. Phys.: Condens. Matter 2014, 26, 284109. 22. García, G.; Atilhan, M.; Aparicio, S. Interfacial Properties of Double Salt Ionic Liquids: A Molecular Dynamics Study. J. Phys. Chem. C 2015, 119, 28405-28416. 23. Arneson, C.; Wawrzyniakowski, Z. D.; Postlewaite, J. T.; Ma, Y. Lithiation and Delithiation Processes in Lithium–Sulfur Batteries from Ab Initio Molecular Dynamics Simulations. J. Phys. Chem. C 2018, 122, 8769-8779. 22

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 34 of 39

24. Lu, H.; Nordholm, S.; Woodward, C. E.; Forsman, J. Ionic Liquid Interface at an Electrode: Simulations of Electrochemical Properties Using an Asymmetric Restricted Primitive Model. J. Phys.: Condens. Matter 2018, 30, 074004. 25. Pal, T.; Beck, C.; Lessnich, D.; Vogel, M. Effects of Silica Surfaces on the Structure and Dynamics of Room-Temperature Ionic Liquids: A Molecular Dynamics Simulation Study. J. Phys. Chem. C 2018, 122, 624-634. 26. Vatamanu, J.; Cao, L.; Borodin, O.; Bedrov, D.; Smith, G. D. On the Influence of Surface Topography on the Electric Double Layer Structure and Differential Capacitance of Graphite/Ionic Liquid Interfaces. J. Phys. Chem. Lett. 2011, 2, 2267-2272. 27. Xing, L.; Vatamanu, J.; Smith, G. D.; Bedrov, D. Nanopatterning of Electrode Surfaces as a Potential Route to Improve the Energy Density of Electric Double-Layer Capacitors: Insight from Molecular Simulations. J. Phys. Chem. Lett. 2012, 3, 1124-1129. 28. Yuan, W.; Zhou, Y.; Li, Y.; Li, C.; Peng, H.; Zhang, J.; Liu, Z.; Dai, L.; Shi, G. The Edgeand Basal-Plane-Specific Electrochemistry of a Single-Layer Graphene Sheet. Sci. Rep. 2013, 3, 2248. 29. Pak, A. J.; Paek, E.; Hwang, G. S. Impact of Graphene Edges on Enhancing the Performance of Electrochemical Double Layer Capacitors. J. Phys. Chem. C 2014, 118, 2177021777. 30. Yang, H.; Yang, J.; Bo, Z.; Zhang, S.; Yan, J.; Cen, K. Edge Effects in Vertically-Oriented Graphene Based Electric Double-Layer Capacitors. J. Power Sources 2016, 324, 309-316. 31. Jiang, D.-e.; Jin, Z.; Wu, J. Oscillation of Capacitance inside Nanopores. Nano Lett. 2011, 11, 5373-5377. 32. Feng, G.; Qiao, R.; Huang, J.; Sumpter, B. G.; Meunier, V. Ion Distribution in Electrified Micropores and Its Role in the Anomalous Enhancement of Capacitance. ACS Nano 2010, 4, 2382-2390. 33. Feng, G.; Cummings, P. T. Supercapacitor Capacitance Exhibits Oscillatory Behavior as a Function of Nanopore Size. J. Phys. Chem. Lett. 2011, 2, 2859-2864. 34. Kondrat, S.; Wu, P.; Qiao, R.; Kornyshev, A. A. Accelerating Charging Dynamics in Subnanometre Pores. Nat. Mater. 2014, 13, 387-393. 35. Rajput, N. N.; Monk, J.; Singh, R.; Hung, F. R. On the Influence of Pore Size and Pore Loading on Structural and Dynamical Heterogeneities of an Ionic Liquid Confined in a Slit 23

ACS Paragon Plus Environment

Page 35 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Nanopore. J. Phys. Chem. C 2012, 116, 5169-5181. 36. Chmiola, J.; Yushin, G.; Gogotsi, Y.; Portet, C.; Simon, P.; Taberna, P. L. Anomalous Increase in Carbon Capacitance at Pore Sizes Less Than 1 Nanometer. Science 2006, 313, 1760-1763. 37. Centeno, T. A.; Sereda, O.; Stoeckli, F. Capacitance in Carbon Pores of 0.7 to 15 Nm: A Regular Pattern. Phys. Chem. Chem. Phys. 2011, 13, 12403-12406. 38. Jiang, D.-e.; Jin, Z.; Henderson, D.; Wu, J. Solvent Effect on the Pore-Size Dependence of an Organic Electrolyte Supercapacitor. J. Phys. Chem. Lett. 2012, 3, 1727-1731. 39. Kalluri, R. K.; Konatham, D.; Striolo, A. Aqueous NaCl Solutions within Charged CarbonSlit Pores: Partition Coefficients and Density Distributions from Molecular Dynamics Simulations. J. Phys. Chem. C 2011, 115, 13786-13795. 40. Kalluri, R. K.; Ho, T. A.; Biener, J.; Biener, M. M.; Striolo, A. Partition and Structure of Aqueous NaCl and CaCl2 Electrolytes in Carbon-Slit Electrodes. J. Phys. Chem. C 2013, 117, 13609-13619. 41. Kornyshev, A. A. Double-Layer in Ionic Liquids: Paradigm Change? J. Phys.Chem. B 2007, 111, 5545-5557. 42. Maggs, A. C.; Podgornik, R. General Theory of Asymmetric Steric Interactions in Electrostatic Double Layers. Soft Matter 2016, 12, 1219-1229. 43. Han, Y.; Huang, S.; Yan, T. A Mean-Field Theory on the Differential Capacitance of Asymmetric Ionic Liquid Electrolytes. J. Phys.: Condens. Matter 2014, 26, 284103. 44. Kilic, M. S.; Bazant, M. Z.; Ajdari, A. Steric Effects in the Dynamics of Electrolytes at Large Applied Voltages. I. Double-Layer Charging. Phys. Rev. E 2007, 75, 021502. 45. Gavish, N.; Elad, D.; Yochelis, A. From Solvent-Free to Dilute Electrolytes: Essential Components for a Continuum Theory. J. Phys. Chem. Lett. 2018, 9, 36-42. 46. Gavish, N.; Yochelis, A. Theory of Phase Separation and Polarization for Pure Ionic Liquids. J. Phys. Chem. Lett. 2016, 7, 1121-1126. 47. Israelachvili, J. Intermolecular and Surface Forces; Academic Press: New York, 2011. 48. Masliyah, J. H.; Bhattacharjee, S. Electrokinetic and Colloid Transport Phenomena; John Wiley & Sons: Hoboken, NJ, 2006. 49. Bazant, M. Z.; Storey, B. D.; Kornyshev, A. A. Double Layer in Ionic Liquids: Overscreening versus Crowding. Phys. Rev. Lett. 2011, 106, 046102. 24

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 36 of 39

50. Rajni; Oh, J. M.; Kang, I. S. Ion Size Effects on the Osmotic Pressure and Electrocapillarity in a Nanoslit: Symmetric and Asymmetric Ion Sizes. Phys. Rev. E 2016, 93, 063112. 51. Moon, G. J.; Ahn, M. M.; Kang, I. S. Osmotic Pressure of Ionic Liquids in an Electric Double Layer: Prediction Based on a Continuum Model. Phys. Rev. E 2015, 92, 063020. 52. Kowalczyk, P.; Furmaniak, S.; Gauden, P. A.; Terzyk, A. P. Methane-Induced Deformation of Porous Carbons: From Normal to High-Pressure Operating Conditions. J. Phys. Chem. C. 2012, 116, 1740-1747. 53. Aguilar-Pineda, G. E.; Jiménez-Ángeles, F.; Yu, J.; Lozada-Cassou, M. Van der WaalsLike Isotherms in a Confined Electrolyte by Spherical and Cylindrical Nanopores. J. Phys. Chem. B 2007, 111, 2033-2044. 54. Hsiao, Y.-W.; Hedström, M. Swelling Pressure in Systems with Na-Montmorillonite and Neutral Surfaces: A Molecular Dynamics Study. J. Phys. Chem. C 2017, 121, 26414-26423. 55. Lee, J. W.; Nilson, R. H.; Templeton, J. A.; Griffiths, S. K.; Kung, A.; Wong, B. M. Comparison of Molecular Dynamics with Classical Density Functional and PoissonBoltzmann Theories of the Electric Double Layer in Nanochannels. J. Chem. Theory Comput. 2012, 8, 2012-2022. 56. Kim, D.; Darve, E. Molecular Dynamics Simulation of Electro-Osmotic Flows in Rough Wall Nanochannels. Phys. Rev. E 2006, 73, 051203. 57. Joly, L.; Ybert, C.; Trizac, E.; Bocquet, L. Hydrodynamics within the Electric Double Layer on Slipping Surfaces. Phys. Rev. Lett. 2004, 93, 257805. 58. Celebi, A. T.; Beskok, A. Molecular and Continuum Transport Perspectives on Electroosmotic Slip Flows. J. Phys. Chem. C 2018, 122, 9699-9709. 59. Giera, B.; Henson, N.; Kober, E. M.; Shell, M. S.; Squires, T. M. Electric Double-Layer Structure in Primitive Model Electrolytes: Comparing Molecular Dynamics with LocalDensity Approximations. Langmuir 2015, 31, 3553-3562. 60. Lee, J. A.; Kang, I. S. Electrocapillarity of an Electrolyte Solution in a Nanoslit with Overlapped Electric Double Layer: Continuum Approach. Phys. Rev. E 2014, 90, 032401. 61. Dommert, F.; Wendler, K.; Berger, R.; Delle Site, L.; Holm, C. Force Fields for Studying the Structure and Dynamics of Ionic Liquids: A Critical Review of Recent Developments. ChemPhysChem 2012, 13, 1625-1637. 62. Schmidt, J.; Krekeler, C.; Dommert, F.; Zhao, Y.; Berger, R.; Delle Site, L.; Holm, C. Ionic 25

ACS Paragon Plus Environment

Page 37 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Charge Reduction and Atomic Partial Charges from First-Principles Calculations of 1,3Dimethylimidazolium Chloride. J. Phys. Chem. B 2010, 114, 6150-6155. 63. Martínez, L.; Andrade, R.; Birgin, E. G.; Martínez, J. M. Packmol: A Package for Building Initial Configurations for Molecular Dynamics Simulations. J. Comput. Chem. 2009, 30, 21572164. 64. Plimpton, S. Fast Parallel Algorithms for Short-Range Molecular Dynamics. J. Comput. Phys. 1995, 117, 1-19. 65. Humphrey, W.; Dalke, A.; Schulten, K. VMD: Visual Molecular Dynamics. J. Mol. Graphics 1996, 14, 33-38. 66. Hockney, R. W.; Goel, S. P.; Eastwood, J. W. Quiet High-Resolution Computer Models of a Plasma. J. Comput. Phys. 1974, 14, 148-158. 67. Nosé, S. A Unified Formulation of the Constant Temperature Molecular Dynamics Methods. J. Chem. Phys. 1984, 81, 511-519. 68. Hoover, W. G. Canonical Dynamics: Equilibrium Phase-Space Distributions. Phys. Rev. A 1985, 31, 1695-1697. 69. Bazant, M. Z.; Kilic, M. S.; Storey, B. D.; Ajdari, A. Towards an Understanding of InducedCharge Electrokinetics at Large Applied Voltages in Concentrated Solutions. Adv. Colloid Interface Sci., 2009, 152, 48-88. 70. Storey, B. D.; Bazant, M. Z. Effects of Electrostatic Correlations on Electrokinetic Phenomena. Phys. Rev. E 2012, 86, 056303. 71. Hua, C. K.; Kang, I. S.; Kang, K. H.; Stone, H. A. Steric Effects of Ions in the ChargeRelated Wetting Phenomena. Phys. Rev. E 2010, 81, 036314. 72. Imre, A. R.; Martinás, K.; Rebelo, L. P. N. Thermodynamics of Negative Pressures in Liquids. J. Non-Equilib. Thermodyn. 1998, 23, 351-375. 73. Tas, N. R.; Mela, P.; Kramer, T.; Berenschot, J. W.; van den Berg, A. Capillarity Induced Negative Pressure of Water Plugs in Nanochannels. Nano Lett., 2003, 3, 1537-1540. 74. Gebbie, M. A.; Smith, A. M.; Dobbs, H. A.; Lee, A. A.; Warr, G. G.; Banquy, X.; Valtiner, M.; Rutland, M. W.; Israelachvili, J. N.; Perkin, S.; Atkin, R. Long Range Electrostatic Forces in Ionic Liquids. Chem. Commun., 2017, 53, 1214-1224. 75. Bo, Z.; Yang, H.; Zhang, S; Yang, J.; Yan, J.; Cen, K. Molecular Insights into Aqueous NaCl Electrolytes Confined within Vertically-oriented Graphenes. Sci. Rep. 2015, 5, 14652. 26

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 38 of 39

76. Burt, R.; Breitsprecher, K.; Daffos, B; Taberna, P.-L.; Simon, P.; Birkett, G.; Zhao, X.

S.; Holm, C.; Salanne, M. Capacitance of Nanoporous Carbon-Based Supercapacitors Is a Trade-Off between the Concentration and the Separability of the Ions. J. Phys. Chem. Lett. 2016, 7, 4015-4021. 77. Yang, H.; Zhang, X.; Yang, J.; Bo, Z.; Hu, M.; Yan, J.; Cen, K. Molecular Origin

of Electric Double-Layer Capacitance at Multilayer Graphene Edges. J. Phys. Chem. Lett., 2017, 8, 153-160. 78. Ma, K.; Wang, X.; Forsman, J.; Woodward, C. E. Molecular Dynamic Simulations

of Ionic Liquid’s Structural Variations from Three to One Layers inside a Series of Slit and Cylindrical Nanopores. J. Phys. Chem. C 2017, 121, 13539-13548.

27

ACS Paragon Plus Environment

Page 39 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

TOC Graphic

28

ACS Paragon Plus Environment