# On the key role of electrolyte-electrode van der Waals interactions in the simulation of ionic liquids-based supercapacitors

Camille Bacon, Alessandra Serva, Céline Merlet, Patrice Simon, Mathieu Salanne

<figure><img src="image_1.png" /></figure>

## To cite this version:

Camille Bacon, Alessandra Serva, Céline Merlet, Patrice Simon, Mathieu Salanne. On the key role of electrolyte-electrode van der Waals interactions in the simulation of ionic liquids-based supercapacitors. Electrochimica Acta, 2023, 455, pp.142380. 10.1016/j.electacta.2023.142380. hal-04131515

HAL Id: hal-04131515

https://hal.science/hal-04131515

Submitted on 16 Jun 2023

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# On the key role of electrolyte-electrode van der Waals interactions in the simulation of ionic liquids-based supercapacitors

Camille Bacon,<sup>1,2</sup> Alessandra Serva,<sup>1,2</sup> Céline Merlet,<sup>3,2</sup> Patrice Simon,<sup>3,2</sup> and Mathieu Salanne<sup>1,2,4, a)</sup>

<sup>1)</sup>Sorbonne Université, CNRS, Physicochimie des Électrolytes et Nanosystèmes Interfaciaux, F-75005 Paris, France

<sup>2)</sup>Réseau sur le Stockage Electrochimique de l'Energie (RS2E), FR CNRS 3459, 80039 Amiens Cedex, France

<sup>3)</sup>CIRIMAT, Université de Toulouse, CNRS, 31062 Toulouse Cedex 9, France

<sup>4)</sup>Institut Universitaire de France (IUF), 75231 Paris, France

(Dated: 7 March 2023)

The performance of supercapacitors is governed by the structure and dynamics of ions at the solid/liquid interface. At the molecular scale, these properties result from a subtle combination of electrolyte-electrolyte and electrolyte-electrode interactions. Although the former are well captured by conventional force fields, validated against experiments, the latter are much more difficult to parameterize accurately. In this work, by using constant potential classical molecular dynamics, we investigate the effect of the strength of the electrode-electrolyte van der Waals interactions on the interfacial properties for a system composed of the 1-ethyl-3-methylimidazolium bis(trifluoromethylsulfonyl)imide ionic liquid and a graphite electrode. We show that stronger van der Waals interactions lead to a decrease in the exchange of co-ions by counter-ions with the increase of potential difference and, thus, to a lower capacitance of the devices. The ion exchange dynamics is strongly affected, but the charging rate remains constant over the whole range of studied parameters. This study emphasizes the need for a careful parameterization of force fields for electrode materials in future classical molecular dynamics studies.

## I. INTRODUCTION

Supercapacitors (SCs) are electrochemical energy storage devices that operate by adsorption and reorganization of ions from an electrolyte on the electrode surface.¹ SCs exhibit high power density and long cycle life, but lower energy density than batteries.² The two devices can complement each other since they are used in different applications. Several experimental and theoretical studies have been conducted over the years to increase the energy density stored in SCs.

In order to obtain devices with both high energy and power density, a deep understanding of the electrode-electrolyte interface at a molecular scale, complementing electrochemical experiments, is required. In this respect, molecular dynamics (MD) simulations appear as a promising tool. In particular, classical MD simulations are employed in the case of SCs because they allow to simulate large systems (> 10,000 atoms) for long times (several nanoseconds).³⁻⁵ One of the crucial aspects of classical molecular dynamics is the force field used to describe the electrolyte-electrolyte and electrode-electrolyte interactions. This force field has to be carefully chosen in order to accurately account for the complex interactions in the system and produce realistic structural and dynamical properties of the interface.

For the bulk electrolytes, generic force fields such as the ones developed by Canongia Lopes, Pádua and coworkers for Ionic Liquids (ILs)⁶⁻⁸ were shown to reproduce accurately the structure and thermodynamics of the

system. The dynamical properties are generally more difficult to capture, a problem that could be solved by introducing rescaled charges⁹ or polarization effects.⁸,¹⁰ The interactions between the electrolyte and the electrode are less well known and studied, despite their importance for accurately determining the physical properties of the interface, but the same force fields are generally used without further reparameterization.¹¹⁻¹⁴

However, previous simulations have shown that modifying the short-range water-carbon interaction energy changes the contact angle of a water droplet on a graphite surface,¹⁵ and similar results have been found for ILs on graphene.¹⁶ Other studies have shown that the van der Waals (vdW) strength modulates the structure of the electrochemical double layer (EDL), as well as the capacitance of aqueous electrolytes or pure ILs in contact with planar electrodes.¹⁷⁻¹⁹ In particular, Yang *et al.* have found a non-monotonic evolution of the capacitance when increasing the strength of the ion-electrode attraction for an aqueous electrolyte.¹⁷ Jin *et al.* have observed a decrease of the capacitance when the strength of the ion-electrode short-range interaction is reduced for the BMIM-PF₆ IL at 450 K. However, both studies used constant charges on the electrodes, instead of a more realistic constant potential approach,²⁰ which could lead to different interfacial properties. Kondrat and co-workers have also widely studied the case of so-called “ionophobic” pores, showing that such systems could improve both the energy and the power density of supercapacitors, but their simulations were not explicitly considering carbon materials.²¹,²²

In this work, we employ constant potential classical MD to investigate the effect of the electrode-

<footnote>a) Electronic mail: mathieu.salanne@sorbonne-universite.fr</footnote>
<figure><img src="image_3.png" /><figcaption>FIG. 1. Snapshot of the simulated system: neat EMIM-TFSI ionic liquid in contact with planar carbon electrodes (electrodes in light green, cations in blue, anions in red).</figcaption></figure>

electrolyte vdW interactions on the capacitance, interfacial structure and dynamics of a neat IL (1-ethyl-3-methylimidazolium bis(trifluoromethylsulfonyl)imide, EMIM-TFSI) in contact with planar carbon electrodes. To focus on ions-electrode interactions, planar electrode have been chosen. Our results show that the effect of the electrode-electrolyte vdW interactions is far from being negligible since all interfacial properties are modified: ion packing, capacitance and ion adsorption dynamics vary noticeably with the vdW interaction strength.

## II. METHODS

Constant potential MD simulations of neat EMIM-TFSI between two graphite electrodes have been performed using the MetalWalls software.²³ A representative snapshot of the system is shown in Figure 1. The electrolyte consists in 322 ion pairs and two electrodes made of three rigid parallel graphene sheets separated by 0.34 nm. The box dimensions are $L_x = 3.4$ nm, $L_y = 3.7$ nm and $L_z = 11.2$ nm, with 2D-periodic boundary conditions applied in the $x$ and $y$ directions. In the constant potential method, the carbon atoms of the electrode carry a Gaussian charge distribution of width of 0.40 Å. The value of the charge is allowed to fluctuate in order to meet the constant potential condition.

The all-atom CL&P⁶,²⁴ force field has been adopted to represent EMIM-TFSI. Atomic charges are scaled by a factor of 0.8 in order to mimic polarizable effects, as suggested in previous studies of bulk liquids.²⁵⁻²⁷ VdW interactions are modelled using a Lennard-Jones (LJ) potential, with two parameters for each interaction, $\epsilon$ and $\sigma$, according to the following equation:

$$U_{ij} = 4\epsilon_{ij} \left[ \left( \frac{\sigma_{ij}}{r_{ij}} \right)^{12} - \left( \frac{\sigma_{ij}}{r_{ij}} \right)^{6} \right] \quad (1)$$

The Lorentz-Berthelot mixing rules are used to obtain cross LJ parameters between different atom types.

For the carbon atoms of the electrodes, three sets of LJ parameters are used in this work. The first set of parameters ($\epsilon_{cc} = 0.23$ kJ mol⁻¹, $\sigma_{cc} = 3.37$ Å) is taken from Ref. 28, and we will refer to it as “weak” (interactions)

set. These parameters were adjusted for graphite in contact with a noble gas, and have already been used for previous simulations of supercapacitors.²⁹,³⁰ For the two other sets of LJ parameters, $\sigma_{cc}$ is kept unchanged while $\epsilon_{cc}$ is increased to 0.51 kJ mol⁻¹ and 1.00 kJ mol⁻¹, respectively. We will refer to these sets of LJ parameters as “moderate” (interactions) and “strong” (interactions), respectively. Since Lorentz-Berthelot mixing rules are used to calculate the ‘cross’ pair interactions ($\sigma_{ij}$ and $\epsilon_{ij}$), modifying $\epsilon_{cc}$ from one set to another means that all the interactions between the graphene carbon atoms and the electrolyte are modified. Figure 2 shows the resulting LJ potentials for the interaction between the carbon atom of the electrode and the nitrogen atom of the EMIM⁺ cation, as an example.

<figure><img src="image_4.png" /><figcaption>FIG. 2. Lennard-Jones potential curves describing the van der Waals interactions between a carbon atom of the electrode and a nitrogen atom of the EMIM⁺ cation for the three values of $\epsilon_{cc}$ studied in this work: 0.23 kJ mol⁻¹ (blue), 0.51 kJ mol⁻¹ (orange) and 1.00 kJ mol⁻¹ (green).</figcaption></figure>

The starting configuration of the system was generated using the PACKMOL software.³¹ In order to obtain the correct liquid density, the electrodes are allowed to move in the $z$ direction for few nanoseconds using a piston, i.e. by applying a pressure of 0 bar on both sides. The temperature applied is 298 K. When the electrode positions are stabilized, a bulk density of 1.576 g cm⁻³ is reached, very close to the experimental value at ambient temperature.³²

For each set of LJ parameters, the system is charged by applying a potential difference at $t = 0$ of 0 V, 1 V, 2 V and 4 V between the two electrodes. The simulations are performed at such constant applied potential in the NVT ensemble using a timestep of 2 fs. A first equilibration of 8 ns is carried out until the total charge on the electrode reaches a plateau (see Figure 3). This equilibration is followed by a production run of 60 ns. The temperature is set to 298 K, and maintained constant using a Nosé-Hoover thermostat chain of length 5 with a relaxation
<figure><img src="image_5.png" /><figcaption>FIG. 3. Time evolution of the total charge. The vertical dashed line shows the transition between the equilibration and the production periods.</figcaption></figure>

time of 500 fs. The cut-off for the vdW interactions is set to 1.48 nm. The Ewald summation method is used to compute long-range interactions, the cut-off distance for real space interactions being of 1.70 nm.

### III. RESULTS AND DISCUSSION

To investigate the effect of the carbon-electrolyte vdW interactions on the interfacial properties, we first want to assess the variation of the electrode charge as a function of the applied potential difference. We start by determining the potential profile across the simulation cell for all sets of interactions. This is done along the *z* direction, perpendicular to the electrode surface (see eq.1 and Figure S1 in SI). This allows us to determine the potential of the bulk liquid, which is used as the reference potential. In the following, $\Delta\Psi = \Psi_{\text{electrode}} - \Psi_{\text{bulk}}$ is defined as the potential drop at each of the interfaces.

Figure 4 displays the total charge at equilibrium on the electrode with respect to $\Delta\Psi$, extracted from the production phase as defined on Figure 3. The strength of the ion-electrode short-range interactions affects the amount of accumulated charge. The average decrease of the total electrode charge when going from $\epsilon_{cc} = 0.23 \text{ kJ mol}^{-1}$ to $\epsilon_{cc} = 0.51 \text{ kJ mol}^{-1}$ is around 11%, and of up to 20% for the larger value of $\epsilon_{cc}$, i.e. 1.00 kJ mol$^{-1}$. Although the difference looks smaller for the lowest applied potential, this may be attributed to the relatively large error on the potential drop value at such potentials.

The differential capacitance of a single electrode is defined as the derivative of this electrode charge with respect to the corresponding potential drop,

<figure><img src="image_6.png" /><figcaption>FIG. 4. Total electrode charge as a function of the difference of potential between the electrode and the bulk, for the three sets of ion-electrode vdW interactions. The standard deviations are smaller than the symbols size.</figcaption></figure>

$$C_{\text{diff}}(\Delta\Psi) = \frac{\partial Q_{\text{electrode}}}{\partial \Delta\Psi} \quad (2)$$

Although it would be possible to do polynomial fits of the $Q_{\text{electrode}} = f(\Delta\Psi)$ plot, in practice, due to the insufficient number of $\Delta\Psi$ values sampled, it is difficult to extract a value for each point. However, by approximating the curve with two different linear functions for the negative and positive $\Delta\Psi$ values, we can extract average $C_{\text{diff}}^{+}$ and $C_{\text{diff}}^{-}$ values for the positive and negative electrodes, respectively. These are reported in Table I.

<table><caption>TABLE I. Average differential capacitances on the positive and negative electrode depending on the $\epsilon_{CC}$ LJ parameter.</caption><tr><td>ε<sub>CC</sub> (kJ mol<sup>-1</sup>)</td><td>0.23</td><td>0.51</td><td>1.00</td></tr><tr><td>C<sup>-</sup><sub>diff</sub> (μF cm<sup>-2</sup>)</td><td>6.9</td><td>7.1</td><td>5.6</td></tr><tr><td>C<sup>+</sup><sub>diff</sub> (μF cm<sup>-2</sup>)</td><td>2.4</td><td>2.2</td><td>1.8</td></tr></table>

We observe a decrease of the differential capacitance when the LJ parameter becomes stronger. Albeit, the values for weak interactions are consistent with previous simulations$^{29}$, it is difficult to compare properly with experiments since the latter results show important discrepancies: as discussed by Tu and McDaniel$^{33}$, the experimental results depend on the measurement method and on the model used to extract the differential capacitance values. However, we note that numbers ranging between 6 and 20 μF cm$^{-2}$ were obtained for planar electrodes and the EMIM-TFSI ionic liquid$^{34-37}$. The agreement is therefore better for the weaker interaction potential, albeit the numbers obtained for the positive electrode are certainly too low. This points towards the existence of additional interactions between the electrode and the
TFSI anion, and will be the subject of future studies. This result shows that one must use properly parameterized LJ parameters to describe the interactions between the ions and the electrode in order to improve the accuracy of capacitance estimations. In addition, Tu et al. have recently shown that accounting for the surface roughness can lead to significant increase of the capacitance due to anion-anion pairing.³⁸ Here we used perfectly flat surfaces, which is probably not the case in experimental studies.

To understand the origin of the variation in the amount of charge stored on the electrode, the number of ions adsorbed at the interface has been determined. An ion is considered adsorbed if at least one of its atoms is closer than 5 Å from the electrode surface. This value is chosen because it corresponds to the distance for which all the ions in contact with the electrode are included but not the ions in the second layer (defined as ions separated from the electrode by at least another ion, see Figure 5). Using such an atom-based criterion is more reliable than selecting the ions only depending on the position of their center of mass, as is usually done,³⁹,⁴⁰ because the latter definition leaves out some ions which only have a small portion in contact with the electrode (Figure S2 in SI).

The top panels of Figure 6 show the evolution of the number of cations and anions adsorbed on the electrode as a function of the electrode potential for the three sets of LJ interaction parameters. In all cases, the charging mechanism observed corresponds mainly to an exchange of co-ions with counter-ions. This means that when the electrode becomes negatively charged, some anions desorb and are replaced by cations, in an approximately equivalent number. The first layer of ions thus

becomes positively charged, with the total charge on the electrode almost exactly counterbalancing the charge of the electrode. An opposite phenomenon takes place when the electrode becomes positively charged. For the more attractive electrodes, fewer exchanges are detected, i.e. fewer co-ions are replaced by counter-ions leading to the observed reduction of the total charge of the first layer of ions.

In addition to the difference in the number of co-ions and counter-ions exchanged, the strength of the ion-electrode vdW interactions also influences the orientation of the adsorbed species with respect to the surface of the electrode. To characterize this, a representative angle $\theta$ is defined for both ions. For the cation, it is the angle between a vector perpendicular to the imidazolium ring plane and the z direction (perpendicular to the electrode surface). For the anion, the angle is defined between the vector connecting the two carbon atoms and the z direction. Figure S3 in SI illustrates the definition of these angles. Following the probability distribution of the angle of the adsorbed ions (Figure S4 in the SI), the cation is considered to be oriented perpendicular when $\theta$ is comprised between 45° and 135°, and parallel otherwise. The anion is considered as parallel for angle values in the range 70°-110° and perpendicular otherwise. This interval allows to take into account angle variations of the anion due to the rotation of the S-N-S-C dihedral, i.e cis and trans configurations.

Percentages of ions in the perpendicular orientation for all the ion-electrode interactions and all potentials studied are shown in the bottom panels of Figure 6. It is clear from this figure that large re-orientation effects take place when the potential of the electrode is changed. The percentage of cations perpendicular to the electrode decreases for negative potentials and increases for positive potentials. On the contrary, the percentage of anions perpendicular to the electrode increases for negative potentials and decreases for positive ones. This behavior is in agreement with previous theoretical and experimental studies on the same system.³⁹,⁴¹⁻⁴³

Remarkably, the ions orientation is also strongly affected by the strength of the attraction between the ions and the electrode. For all potential differences applied, less ions are oriented perpendicular to the electrode for stronger attractions. This variation in the orientation is more pronounced for the cations. Similarly, for the BMIM–PF₆ IL, Jin et al. see more cations oriented parallel to the surface when the force of the interactions between the ions and the electrode is increased at both the positive and negative surfaces.¹⁸ The structure of such systems consists in successive layers of dense fluid, up to a few nanometers, as was evidenced in many studies (see Figure S5 in the SI).³⁹,⁴¹,⁴²,⁴⁴,⁴⁵ For EMIM-TFSI here, the change in orientation leads to an increase in the intensity observed for the first peak in the ion density profiles. This is expected as increasing $\epsilon_{cc}$ enhances the attraction between the atoms of the ions and the atoms of the electrode, and a larger number of atoms can be found in

<figure><img src="image_7.png" /><figcaption>FIG. 5. Snapshot of one configuration of the studied system with a zoom at one interface. Carbon atoms are represented in light blue, adsorbed ions in purple and ions in the other layers in orange.</figcaption></figure>
<figure><img src="image_8.png" /><figcaption>FIG. 6. Top panels: number of cations (left) and anions (right) adsorbed on the electrode per nm² as a function of the potential difference between the electrode and the bulk, for the three values of ε<sub>cc</sub>. Bottom panels: percentage of adsorbed cations (left) and anions (right) oriented perpendicular to the electrode with respect to the potential difference between the electrode and the bulk. Error bars shown correspond to the standard deviations.</figcaption></figure>

contact with the electrode if the ions are oriented paral-
lel. It should be noted however that the number of ions
parallel to the electrode does not depend only on vdW
interactions. It results from a balance between vdW in-
teractions, electrostatic interactions and steric effects.

Surprisingly, the total number of ions adsorbed on the
electrode (see Figure S6 in the SI) remains constant, on
average, for any value of ε<sub>cc</sub>. This is due to the pres-
ence of more ions parallel to the surface, and means that
the enhancement of the short-range attraction strength
between the ions and the electrode changes the relative
balance between electrostatic and steric interactions at
the interface.

Apart from the changes in the first layer of ions ad-
sorbed on the electrode surface, it is interesting to see
whether modifying the strength of the ions–electrode
short-range interactions also affects the layering adopted

by the liquid at longer range. Upon application of a non-
zero potential on the electrode, the first layer becomes
enriched in counter-ions, while the second one becomes
enriched in co-ions and so on, up to typically 3 or 4 layers
of liquid depending on the system. Ionic density profiles
for anions and cations for applied potential differences of
1 V and 4 V are shown in Figure 7. In the case of the
weak interaction potential, the ionic density oscillations
at both interfaces are enhanced when the applied poten-
tial difference increases. A similar behavior is observed
at the negative interface for the moderate and strong in-
teraction sets. Interestingly, it is much less pronounced
at the positive electrode. In fact, for the moderate and
strong interactions the oscillations remain almost iden-
tical independently of the potential at the positive in-
terface. Such results are surprising as short-range vdW
interactions vanish after less than a nanometer (see Fig-
<figure><img src="image_9.png" /></figure>

FIG. 7. Ionic density profiles in the direction perpendicular to the electrode (z) for the cations (top panels) and for the anions (bottom panels), for applied potential differences of 1 V (left panels) and 4 V (right panels). The negative electrode is on the left, with the graphite sheet in contact with the electrolyte at z = 0 nm. The positive electrode is on the right. See Figure S7 for the same quantities at applied potentials of 0 and 2 V.

ure 2), so they are not expected to perturb ionic densities at a longer range.

The layered structure of the electrolyte at the interface is mostly due to packing constraints, which induce a preferential organization of the ions in terms of orientation. Figure 8 shows the profile of the mean value of $cos θ$ along the z axis at 4V. The oscillations in the orientation profile, which were already observed in previous studies,⁴¹ resemble the oscillations of the ionic density profiles reported in Figure 7, which means that the layers observed correspond to a perturbation in the orientation of the ions due to the interface. When comparing the orientation profiles obtained for the three sets of LJ parameters, one can see that the oscillations in the case of the cations are reduced at the positive interface when the ion-electrode interactions are stronger, while they remain similar at the negative interface. This behavior is similar to the one observed for the density profiles. On the contrary, the orientation profiles for the anions seem to be independent on the strength of the ion-electrode interactions.

Finally, the dynamic properties have also been investi-

<figure><img src="image_10.png" /></figure>

FIG. 8. Mean orientation profiles along the z direction for the cations (top panels) and for the anions (bottom panels) for an applied potential difference of 4 V. $cos(θ) = 1$, respectively 0, corresponds to the ions oriented parallel, respectively perpendicular, to the electrode.

gated. Figure 9 shows the time evolution of the number of cations and anions adsorbed along the simulation for a potential difference of 0 V. The curve is the running av-

<figure><img src="image_11.png" /></figure>

FIG. 9. Time evolution of the number of cations (left panel) and anions (right panel) adsorbed at the interface for an applied potential of 0 V. The running average over 1 ns is displayed in order to reduce the variations due to ions that shortly leave or enter the cut-off distance of 5 Å.

erage over 1 ns in order to remove some of the noise due to ions that shortly leave or enter the cut-off distance
of 5 Å used to define the adsorbed layer. For the same simulation time, fewer adsorption or desorption events are detected when the strength of the attraction between the ions on the electrode is raised. Note that in this case the number of adsorbed ions does not oscillate around a mean value. However, this has no significant impact on the total charge on the electrodes which fluctuates around an equilibrated mean value (see Figure 3). As charging storing is due to adsorption and desorption of ions, the reduction of the dynamics of adsorption is expected to slow down the charging dynamics. However, we do not provide charging times, as in order to obtain reliable values it would be necessary to run several MD simulations in the charging region.

## IV. CONCLUSIONS

In this work we have studied the effect of the strength of electrolyte-electrode vdW interactions on the structural, dynamical and capacitive properties of the interface between pure EMIM-TFSI and a planar carbon electrode. The decrease of capacitance resulting from the increase of the ion-electrode attraction observed in this study is related to the limitation of exchanges between co-ions and counter-ions adsorbed at the interface when a difference of potential between the electrodes is applied. In addition, modifying the ion-electrode short-range interactions impacts the structure of the ions at a distance where these interactions are expected to have vanished. In fact, the orientation of the ions adsorbed is modified due to the increase of the adsorption strength, and this modification perturbs the orientation in the subsequent layers of ions.

In the study of electrochemical systems, and in particular when the electrolyte is a pure ionic liquid, it is commonly admitted that Coulomb interactions play a major role in the structural and dynamic properties, with vdW interactions having a minor role. This common assumption has led to put a lot of effort to improve the description of the electrostatic interactions in such system. For example, modeling the electrodes using constant potential²⁰ method is nowadays widely used, as well as, the development of models including explicit polarization for the ions of the electrolyte.³³,⁴⁶

Yet, this work shows that the intensity of the ion-electrode interactions also has a significant impact on the electrode-electrolyte interface, in particular in the determination of the capacitance. As a consequence, the development of a dedicated force field not only for the electrolyte but also for the electrode is of a particular importance in the simulation of electrolyte-electrode interfaces if one wants to determine the capacitance at a quantitative level. Challenges are associated with this endeavour. Modern force fields are generally parameterized using electronic structure calculations at the Density Functional theory (DFT) level, yet in the case of carbon surfaces the vdW interactions are notoriously difficult

to catch.⁴⁷,⁴⁸ The use of higher level methods, such as Quantum Monte-Carlo, can help to benchmark the DFT calculations and ultimately lead to improved classical force fields between the carbon and ionic species, as was shown recently for lithium.⁴⁹ However the case of organic ions such as the ones constituting typical ionic liquids are even more challenging and thus may require significant amounts of work in the future years.

## SUPPLEMENTARY MATERIAL

See the supplementary material for Figures S1-S7 as described in the manuscript.

## ACKNOWLEDGMENTS

This project has received funding from the European Research Council under the European Union's Horizon 2020 research and innovation programme (grant agreements No. 771294 and 714581). This work was supported by the French National Research Agency (Labex STORE-EX, Grant ANR-10-LABX-0076). The authors acknowledge HPC resources granted by GENCI (resources of IDRIS, Grant No. A0120910463 and A0110911061) and by the SACADO service unit of Sorbonne University.

¹J. M. Griffin, A. C. Forse, W.-Y. Tsai, P.-L. Taberna, P. Simon, and C. P. Grey, "In situ NMR and electrochemical quartz crystal microbalance techniques reveal the structure of the electrical double layer in supercapacitors," Nat. Mater. **14**, 812-819 (2015).

²J. R. Miller and P. Simon, "MATERIALS SCIENCE: Electrochemical capacitors for energy management," Science **321**, 651-652 (2008).

³K. Xu, H. Shao, Z. Lin, C. Merlet, G. Feng, J. Zhu, and P. Simon, "Computational insights into charge storage mechanisms of supercapacitors," Energy Environ. Mat. **3**, 235-246 (2020).

⁴S. Park and J. G. McDaniel, "Interference of electrical double layers: Confinement effects on structure, dynamics, and screening of ionic liquids," J. Chem. Phys. **152**, 074709 (2020).

⁵G. Jeanmairet, B. Rotenberg, and M. Salanne, "Microscopic simulations of electrochemical double-layer capacitors," Chem. Rev. **122**, 10860-10898 (2022).

⁶J. N. C. Lopes and A. A. H. Pádua, "Molecular force field for ionic liquids composed of triflate or bistriflimide anions," J. Phys. Chem. B **108**, 16893-16898 (2004).

⁷S. W. Coles, M. Mishin, S. Perkin, M. V. Fedorov, and V. B. Ivanistsev, "The nanostructure of a lithium glyme solvate ionic liquid at electrified interfaces," Phys. Chem. Chem. Phys. **19**, 11004-11010 (2017).

⁸K. Goloviznina, J. N. Canongia Lopes, M. Costa Gomes, and A. A. H. Padua, "Transferable, polarizable force field for ionic liquids," J. Chem. Theory Comput. **15**, 5858-5871 (2019).

⁹C. Schröder and O. Steinhauser, "Simulating polarizable molecular ionic liquids with drude oscillators," J. Chem. Phys. **133**, 154511 (2010).

¹⁰J. G. McDaniel, E. Choi, C. Y. Son, J. R. Schmidt, and A. Yethiraj, "Ab initio force fields for imidazolium-based ionic liquids," J. Phys. Chem. B **120**, 7024-7036 (2016).

¹¹T. Méndez-Morales, J. Carrete, M. Pérez-Rodríguez, O. Cabeza, L. J. Gallego, R. M. Lynden-Bell, and L. M. Varela, "Molecular dynamics simulations of the structure of the graphene-ionic
liquid/alkali salt mixtures interface,” Phys. Chem. Chem. Phys. **16**, 13271–13278 (2014).

¹²V. Ivanistsev, T. Méndez-Morales, R. M. Lynden-Bell, O. Cabeza, L. J. Gallego, L. M. Varela, and M. V. Fedorov, “Molecular origin of high free energy barriers for alkali metal ion transfer through ionic liquid–graphene electrode interfaces,” Phys. Chem. Chem. Phys. **18**, 1302–1310 (2016).

¹³A. M. Sampaio, G. F. L. Pereira, M. Salanne, and L. J. A. Siqueira, “Comparing the performance of sulfonium and phosphonium ionic liquids as electrolytes for supercapacitors by molecular dynamics simulations,” Electrochim. Acta **364**, 1371–1380 (2020).

¹⁴P. H. L. Ferreira, A. M. Sampaio, and L. J. A. Siqueira, “Energy and power performances of binary mixtures of ionic liquids in planar and porous electrodes by molecular dynamics simulations,” Electrochim. Acta **410**, 13998–1402 (2022).

¹⁵T. Werder, J. H. Walther, R. L. Jaffe, T. Halicioglu, and P. Koumoutsakos, “On the water-carbon interaction for use in molecular dynamics simulations of graphite and carbon nanotubes,” J. Phys. Chem. B **107**, 1345–1352 (2003).

¹⁶R. Burt, G. Birkett, M. Salanne, and X. S. Zhao, “Molecular dynamics simulations of the influence of drop size and surface potential on the contact angle of ionic-liquid droplets,” J. Phys. Chem. C **120**, 15244–15250 (2016).

¹⁷H. Yang, Z. Bo, J. Yang, J. Yan, and K. Cen, “Towards understanding the effects of van der Waals strengths on the electric double-layer structures and capacitive behaviors,” J. Power Sources **366**, 218–225 (2017).

¹⁸W. Jin, X. Liu, Y. Han, S. Li, and T. Yan, “Effects of repulsive interaction on the electric double layer of an imidazolium-based ionic liquid by molecular dynamics simulation,” Phys. Chem. Chem. Phys. **17**, 2628–2633 (2015).

¹⁹R. Sundararaman, D. Vigil-Fowler, and K. Schwarz, “Improving the accuracy of atomistic simulations of the electrochemical interface,” Chem. Rev. **122**, 10651–10674 (2022).

²⁰C. Merlet, C. Péan, B. Rotenberg, P. A. Madden, P. Simon, and M. Salanne, “Simulating supercapacitors: Can we model electrodes as constant charge surfaces?” J. Phys. Chem. Lett. **4**, 264–268 (2012).

²¹S. Kondrat, P. Wu, R. Qiao, and A. A. Kornyshev, “Accelerating Charging Dynamics in Subnanometre Pores,” Nat. Mater. **13**, 387–393 (2014).

²²S. Kondrat and A. Kornyshev, “Pressing a spring: what does it take to maximize the energy storage in nanoporous supercapacitors?” Nanoscale Horiz. **1**, 45–52 (2016).

²³A. Marin-Laflèche, M. Haefele, L. Scalfi, A. Coretti, T. Dufils, G. Jeanmairet, S. Reed, A. Serva, R. Berthin, C. Bacon, S. Bonella, B. Rotenberg, P. Madden, and M. Salanne, “Metal-Walls: A classical molecular dynamics software dedicated to the simulation of electrochemical systems,” J. Open Source Softw. **5**, 2373 (2020).

²⁴J. N. C. Lopes, J. Deschamps, and A. A. H. Pádua, “Modeling ionic liquids using a systematic all-atom force field,” J. Phys. Chem. B **108**, 11250–11250 (2004).

²⁵V. Chaban, “Polarizability versus mobility: atomistic force field for ionic liquids,” Phys. Chem. Chem. Phys. **13**, 16055 (2011).

²⁶B. L. Bhargava and S. Balasubramanian, “Refined potential model for atomistic simulations of ionic liquid [bmim][PF6],” J. Chem. Phys. **127**, 114510 (2007).

²⁷Y. Zhang and E. J. Maginn, “A simple AIMD approach to derive atomic charges for condensed phase simulation of ionic liquids,” J. Phys. Chem. B **116**, 10036–10048 (2012).

²⁸M. W. Cole and J. R. Klein, “The interaction between noble gases and the basal plane surface of graphite,” Surf. Sci. **124**, 547–554 (1983).

²⁹C. Merlet, M. Salanne, B. Rotenberg, and P. A. Madden, “Imidazolium ionic liquid interfaces with vapor and graphite: Interfacial tension and capacitance from coarse-grained molecular simulations,” J. Phys. Chem. C **115**, 16613–16618 (2011).

³⁰C. Merlet, B. Rotenberg, P. A. Madden, P.-L. Taberna, P. Simon, Y. Gogotsi, and M. Salanne, “On the molecular origin of supercapacitance in nanoporous carbon electrodes,” Nat. Mater. **11**, 306–310 (2012).

³¹L. Martínez, R. Andrade, E. G. Birgin, and J. M. Martínez, “PACKMOL: A package for building initial configurations for molecular dynamics simulations,” J. Comput. Chem. **30**, 2157–2164 (2009).

³²J. L. Trenzado, Y. Rodríguez, A. Gutiérrez, A. Cincotti, and S. Aparicio, “Experimental and molecular modeling study on the binary mixtures of [emim][bf4] and [emim][tfsi] ionic liquids,” J. Mol. Liq. **334**, 116049 (2021).

³³Y.-J. Tu and J. G. McDaniel, “Structure capacitance relationships of graphene/ionic liquid electrolyte double layers,” J. Phys. Chem. C **125**, 20204–20218 (2021).

³⁴D. J. Bozym, B. Uralcan, D. Limmer, M. A. Pope, N. J. Szamreta, P. G. Debenedetti, and I. A. Aksay, “Anomalous capacitance maximum of the glassy carbon- ionic liquid interface through dilution with organic solvents,” J. Phys. Chem. Lett. **6**, 2644–2648 (2015).

³⁵M. Jitvisate and J. R. T. Seddon, “Direct measurement of the differential capacitance of solvent-free and dilute ionic liquids,” J. Phys. Chem. Lett. **9**, 126–131 (2018).

³⁶J. M. Klein, E. Panichi, and B. Gurkan, “Potential dependent capacitance of [emim][tfsi], [n1114][tfsi] and [pyr 13][tfsi] ionic liquids on glassy carbon,” Phys. Chem. Chem. Phys. **21**, 3712–3720 (2019).

³⁷C. Gu, L. Yin, S. Li, B. Zhang, X. Liu, and T. Yan, “Differential capacitance of ionic liquid and mixture with organic solvent,” Electrochim. Acta **367**, 137517 (2021).

³⁸Y.-J. Tu, S.-C. Wu, and J. G. McDaniel, “Importance of anion–anion pairing for capacitance of carbon/ionic liquid interfaces,” J. Phys. Chem. C **126**, 20213–20225 (2022).

³⁹J. Vatamanu, O. Borodin, D. Bedrov, and G. D. Smith, “Molecular dynamics simulation study of the interfacial structure and differential capacitance of alkylimidazolium bis(trifluoromethanesulfonyl)imide [cnmim][TFSI] ionic liquids at graphite electrodes,” J. Phys. Chem. C **116**, 7940–7951 (2012).

⁴⁰C. Merlet, D. T. Limmer, M. Salanne, R. van Roij, P. A. Madden, D. Chandler, and B. Rotenberg, “The electric double layer has a life of its own,” J. Phys. Chem. C **118**, 18291–18298 (2014).

⁴¹N. N. Rajput, J. Monk, and F. R. Hung, “Structure and dynamics of an ionic liquid confined inside a charged slit graphitic nanopore,” J. Phys. Chem. C **116**, 14504–14513 (2012).

⁴²A. Fang and A. Smolyanitsky, “Simulation study of the capacitance and charging mechanisms of ionic liquid mixtures near carbon electrodes,” J. Phys. Chem. C **123**, 1610–1618 (2019).

⁴³S. Baldelli, “Surface structure at the ionic liquid-electrified metal interface,” Acc. Chem. Res. **41**, 421–431 (2008).

⁴⁴G. Feng, S. Li, W. Zhao, and P. T. Cummings, “Microstructure of room temperature ionic liquids at stepped graphite electrodes,” AIChE J **61**, 3022–3028 (2015).

⁴⁵R. Atkin, S. Z. El Abedin, R. Hayes, L. H. S. Gasparotto, N. Borisenko, and F. Endres, “AFM and STM studies on the surface interaction of [BMP]TFSA and [EMIm]TFSA ionic liquids with Au(111),” J. Phys. Chem. C **113**, 13266–13272 (2009).

⁴⁶A. Coretti, C. Bacon, R. Berthin, A. Serva, L. Scalfi, I. Chubak, K. Goloviznina, M. Haefele, A. Marin-Laflèche, B. Rotenberg, S. Bonella, and M. Salanne, “MetalWalls: Simulating electrochemical interfaces between polarizable electrolytes and metallic electrodes,” J. Chem. Phys. **157**, 184801 (2022).

⁴⁷Y. Al-Hamdani, D. Alfe, and A. Michaelides, “How strongly do hydrogen and water molecules stick to carbon nanomaterials?” J. Chem. Phys. **146**, 094701 (2017).

⁴⁸J. G. Brandenburg, A. Zen, D. Alfe, and A. Michaelides, “Interaction between water and carbon nanostructures: How good are current density functional approximations?” J. Chem. Phys. **151**, 164707 (2019).

⁴⁹M. Ruggeri, K. Reeves, T.-Y. Hsu, G. Jeanmairet, M. Salanne, and C. Pierleoni, “Multi-scale simulation of the adsorption of
lithium ion on graphite surface: From quantum monte carlo to
molecular density functional theory," J. Chem. Phys. 156, 094709

(2022).