Description of Ion Properties Using Molecular Orbital Energy Levels

1. Introduction

The vast majority of activity parameters (descriptors) characterize the catalytic activity of heterogeneous catalysts. [1,2] Researchers have proposed a lot of electronic parameters that quite accurately describe the activity of various materials: centers of transition metal atomic bands, [3–5] orbital occupancy for transition metal oxides, [6] Fermi level, [7,8] effective coordination numbers [9] and associated electron density for metals, [10–12] density of states at the Fermi level for graphene, [13,14] etc. Such studies normally consider industrial reactions, for example, oxygen reduction, [14,15] hydrogen oxidation, [13,16] nitrogen fixation, [17,18] CO₂ conversion, [19,20] etc.

The activity parameters for molecules, [21–25] which are mainly used in molecular catalysis, design of biologically active substances, etc., are somewhat different. First and foremost in this list of parameters are the values of ionization potential (IP) and electron affinity (EA) and the associated energies of the highest occupied molecular orbital (HOMO) IP=-ε_HOMO and the lowest unoccupied molecular orbital (LUMO) EA=-ε_LUMO. [23] There are also a number of derived parameters based on HOMO and LUMO: [21–23] orbital energy difference (ε_HOMO-LUMO), chemical potential (µ), hardness (η), softness (S), electronegativity (χ), and global electrophility index (ω). These molecular descriptors find applications in catalysis, [21,24] corrosion studies, [25] evaluation of biological activity of compounds, [24] etc. In addition, complex descriptors for neutral molecules [26] and ions [26,27] were proposed for characterizing some parameters of complexes (dissociation constants, reaction rate constants, etc.). These parameters describe specific properties of complexes of molecules or ions based on a whole range of parameters of the system under study (dipole moment, molecular volume, hydrogen bonding strength, etc.).

The vast majority of studies devoted to searching for activity descriptors usually consider a specific reaction proceeding on a variety of materials or active centers. [6,13–15,17–20] The main goal of such studies is to select materials on which the considered reaction has the highest rate. Despite the current trends, it is also appropriate to consider the reverse problem, where the study concerns the interaction of a variety of absorbents on the surface of one catalyst or molecule. Such studies can be extremely important for controlling the surface composition and surface activity of various applications in various industrial processes. In addition to catalysis, an example of an applied problem that can be solved is selection of the optimal electrolyte [28] for a particular electrode material in batteries or supercapacitors. Another practical application is binding of toxic impurities from solutions, especially binding of heavy metal ions. [29–31] Such descriptors will make it possible to evaluate the efficiency of a particular sorbent for a range of ions.

The existence of such descriptors is confirmed by the data of refs, [32–34] which show correlations of the adsorption energy of neutral molecules on the oxide surface with the HOMO and LUMO energy level and other parameters. Refs. [35,36] reveals correlations of the solvation energy with hardness. These facts suggest the possibility of a simplified description of the dependence of adsorption and some physicochemical and electrochemical parameters of molecules or ions on the energy level of the orbitals.

For better clarity and potential practical application, the present work analyzes the adsorption coupling of various ions on the surface of a dielectric (aluminum oxide (0001) surface), a semiconductor (graphene) and a metal (gold (111) surface). The intermolecular interaction is analyzed by studying the binding of water and dimethyl carbonate (solvents frequently used in chemical and electrochemical practice) to ions of molecules. Correlations between some physicochemical and electrochemical parameters of ions and their electronic structure are also analyzed. The systems under study are a number of surface inactive and active ions, which are used in electrochemical investigations of surface phenomena [37,38] and electron transfer, [37,39] organic electrochemistry [40–43] and potentially applicable in chemical current sources, [28,44–51] such as batteries, supercapacitors, etc.

2. Methods

Quantum-chemical calculations of ion adsorption on slabs were performed in the Quantum Espresso software package [52] within the generalized gradient approximation by the Perdew-Burke-Ernzerhof (PBE) functional [53] using standard PAW pseudopotentials. [54,55] Spin polarization was taken into account in the calculations and the maximum energy limit was 600 eV. The Monkhorst-Pack special point method [56] with an electron smearing of 0.1 eV and a k-grid density of 5x5x1 by the Methfessel-Paxton method [57] was selected for integration over the Brillouin zone. The energy convergence criterion of the self-consistent field method was set at 10^- [5] eV. The van der Waals interactions were considered within the Grimme method (DFT-D3). [58]

The relaxation of atoms in the slabs was carried out by the Broyden-Fletcher-Goldfarb-Shanno (BFGS) algorithm, the convergence criterion for the energy optimization of the atom geometry was set at 10^– [4] eV and the maximum force was 0.01 eV/Å.

A 72-atom graphene, a 60-atom Au(111), and a 60-atom α-Al₂O₃ (0001) slabs with cross-sectional dimensions of 6x6, 4x4, and 4x4 unit cells, respectively, were used as the model electrodes. The adsorption energy was calculated as the energy of the optimized slab with the ion on the surface minus the energy of the slab with the ion at a distance of 9−11 Å from the surface.

The cations in the study were the following systems contained in ionic liquids, [42,43] electrolytes for power sources, [48] including aqueous batteries that have been popular in the last decade, [28,47] as well as supporting electrolytes used in organic electrochemistry: [40,41] alkali metal cations (Li⁺, Na⁺, K⁺, Cs⁺, Rb⁺), tetraalkylammonium cations R₄N⁺ (R = -Me, -Et, -Pr, -Bu), derivatives of the pyridinium cation RPy⁺ (R = -H, -Me, -Et), alkylmethylpyrrodilidone cation RMPyr⁺ (R = -Et, -Bu) and alkylmethylimidazolium cation RMI⁺ (R = -Et, -Bu, -Hexcyl (H), -Octyl (O), -Dodecyl (D), -Ph), cetyltrimethylammonium Me₃NCet⁺ (CTA⁺) and cetylpyridinium CetPy⁺ (CP⁺) cations.

The anions were selected in a similar way: [47–49] OH^-, ClO₄^-, NO₃^-, and dycyanamide (DCA^-) grouped into a separate category for convenience; halides (F^-, Cl^-, Br^-); RSO₃^- derivatives (R = -NH₂, -Me, -Et, -Cet, -Ph) and HSO₄^- (R = -OH) added here; RCOO^- carboxylates (R= -Me, -Et, -Pr, i-Pr and -Bu) and fluorine-containing anions (CF₃COO^-, tetrafluoroborate BF₄^-, hexaflurophosphate PF₆^-, trifluromethanesulfonate TfO^- and bis(trifluromethylsulfonyl) imide TFSI^-).

In addition to the above systems, we considered the products of water autoprotholysis: hydroxide (OH^-) and hydroxonium (H₃O⁺) ions.

The geometry optimization of ion complexes with water and dimethyl carbonate (DMC) molecules was carried out in the ORCA software package [59] with the B3LYP functional [60,61] and the 6-31++G(d,p) basis set. The free energy of ion-molecule binding was calculated as the difference between the free energy of the optimized aqueous complex and the original ions and water or DMC molecule in vacuum.

The dipole moment, polarizability, solvation energy and redox potential were calculated in the 6-311++G(3df,3pd) basis set. The solvation and redox potentials were determined within the PCM model [62–64] taking into account the dielectric parameters of the respective solvent; H₂O or DMC was used as the environment. The dielectric and other physical parameters of DMC for the PCM model construction were taken from refs. [65–67] For simplicity, the solvation energy (G_solv) was calculated as the difference between the free energy of the ion in the PCM model (G_PCM) and that of the ion in vacuum (G_vac): G_solv = G_PCM – G_vac+G^°→*, where G^°→*=-k_BT∙ln(24.4509) at T = 298 K is the term describing the process transition of the solute from gas 1 atm to solution 1 mole/L. The redox potential was calculated for the cations from the reduction reaction Ion⁺ + e^- = Ion [0] and for the anions from the oxidation reaction Ion^- – e^- = Ion [0] as the difference between the energy of the initial ion (G_PCM(Ion^±)) and that of the reduced or oxidized ion (GPCM(Ion0)): [28,68] E_pot = – [G_PCM(Ion⁺) – G_PCM(Ion [0] )]/F or [G_PCM(Ion [0] ) – G_PCM(Ion^-)]/F, where F was the Faraday constant. The absolute value of the redox potential (E_pot) was converted to a value relative to the standard hydrogen electrode (SHE): E_pot (vs. SHE) =E_pot - 4.44 V. [69] The values of the electronic levels of the ions were calculated by the B3LYP functional in the 6-311++G(3df,3pd) basis set in vacuum.

The distance between the center of mass of the ion and the surface (z₀) or the center of mass of the solvent molecule (r₀) was calculated using the Python library of the Atomic Simulation Environment (ASE). [70] A number of adsorbent and aqueous complex parameters, such as adsorption energies (E_ads), binding energies (G_bind) and distances (z₀ and r₀), were statistically averaged taking into account different conformations of the ions and solvent molecules.

3. Results and Discussion

3.1. Relation of Electronic Levels of Ions to Adsorption and Binding Energy

Let us evaluate how the energy level of an ion's orbitals can affect binding into an adsorption complex with solvent molecules or an electrode surface.

To consider the simplest scenario of the interaction of an ion with a solvent molecule through a one-electron operator, let us represent their wave functions in terms of Gaussian orbitals for a solvent molecule (assuming that the solvent consists of sufficiently light atoms) ψ_solv=A_solv∙rⁿ∙exp(-α_solv∙r [2]), where n=0 for the s-orbital and n=1 for the p-orbital and ψ_ion=A_ion∙rⁿ∙exp(-α_ion∙r [2]) for the s-, p- or d-orbital for the ion (n=2 for the d-orbital). A_solv and A_ion are normalizing coefficients, α_solv=0.1 a.u.^- [2] and α_ion are exponential coefficients, and all the wave functions are normalized $〈{\mathsf{\psi }}_{\mathrm{solv}}|{\mathsf{\psi }}_{\mathrm{solv}}〉=〈{\mathsf{\psi }}_{\mathrm{ion}}|{\mathsf{\psi }}_{\mathrm{ion}}〉=1$. The overlap of the wave functions is given by the expression $\mathrm{S}=〈{\mathsf{\psi }}_{\mathrm{solv}}|{\mathsf{\psi }}_{\mathrm{ion}}〉$.

The eigenvalues of the orbitals of the solvent molecule and ion [71,72] will look as follows:

$${\mathsf{\epsilon }}_{\mathrm{solv}}=〈{\mathsf{\psi }}_{\mathrm{solv}}\left|\stackrel{\wedge }{\mathrm{H}}\right|{\mathsf{\psi }}_{\mathrm{solv}}〉$$

(1)

$$\mathrm{and}{\mathsf{\epsilon }}_{\mathrm{ion}}=〈{\mathsf{\psi }}_{\mathrm{ion}}\left|\stackrel{\wedge }{\mathrm{H}}\right|{\mathsf{\psi }}_{\mathrm{ion}}〉$$

(2)

where

$\stackrel{\wedge }{\mathrm{H}}={\mathrm{H}}_{\mathrm{kin}}+{\mathrm{H}}_{\mathrm{en}}=-{\displaystyle \frac{{\mathrm{h}}^2}{{2\mathrm{m}}_{\mathrm{e}}}}{\nabla }^2+{\displaystyle \sum _{\mathrm{i}}\frac{{\mathrm{e}}_0^2{\mathrm{q}}_{\mathrm{i}}}{{\mathsf{\epsilon }}_0\left|\mathrm{r}-{\mathrm{r}}_{\mathrm{i}}\right|}}$ is the one-electron perturbation operator, which is the sum of the electron kinetic energy (H_kin) and electron-nuclear interaction (H_en); ${\nabla }^2$ is the Laplace operator; $\left|\mathrm{r}-{\mathrm{r}}_{\mathrm{i}}\right|=\sqrt{{\left(\mathrm{x}-{\mathrm{x}}_{\mathrm{i}}\right)}^2+{\left(\mathrm{y}-{\mathrm{y}}_{\mathrm{i}}\right)}^2+{\left(\mathrm{z}-{\mathrm{z}}_{\mathrm{i}}\right)}^2}$ is the distance between the electron density (product of the wave functions) at an arbitrary point x, y, z and a point charge q_i with x_i, y_i, and z_i coordinates of the i-atom; h is the Planck constant; m_e and e₀ are the mass and charge of the electron; ε₀ is the vacuum dielectric constant; i is the serial number of an atom in the solvent molecule.

The matrix element (V) for the orbitals is expressed through operator $\stackrel{\wedge }{\mathrm{H}}$:

$$\mathrm{V}=〈{\mathsf{\psi }}_{\mathrm{solv}}\left|\stackrel{\wedge }{\mathrm{H}}\right|{\mathsf{\psi }}_{\mathrm{ion}}〉$$

(3)

The binding energy between the ψ_solv and ψ_ion orbitals is described by the expression: [4]

$${\mathrm{E}}_{\mathrm{bind}}{=2\mathsf{\epsilon }}_{+}-{\mathsf{\epsilon }}_{\mathrm{solv}}-{\mathsf{\epsilon }}_{\mathrm{ion}}$$

(4)

where ${\mathsf{\epsilon }}_{+}$ is the solution to the Schrödinger equation for the ground state,

$${\mathsf{\epsilon }}_{+}={\displaystyle \frac{\left({\mathsf{\epsilon }}_{\mathrm{solv}}{+\mathsf{\epsilon }}_{\mathrm{ion}}\right)-\sqrt{4{\mathrm{V}}^2+\left(1-{\mathrm{S}}^2\right){\left({\mathsf{\epsilon }}_{\mathrm{solv}}-{\mathsf{\epsilon }}_{\mathrm{ion}}\right)}^2}}{2\left(1-{\mathrm{S}}^2\right)}}$$

(5)

The energy level of the ψ_ion orbital can be changed by varying the α_ion parameter in the range of 0.001 – 0.5 a.u.^- [2].

Now let us consider the adsorption of ions on the electrode using the Newnes-Andersen (NA) model. [3,73] According to the NA model, E_ads is determined by the sum of the contributions: [4,5,74,75]

$${\mathrm{E}}_{\mathrm{ads}}{=\mathrm{E}}_{\mathrm{hyb}}+{\mathrm{E}}_{\mathrm{ortho}}$$

(6)

The first term ${\mathrm{E}}_{\mathrm{hyb}}={\displaystyle \underset{-\infty }{\overset{+\infty }{\int }}\mathsf{\epsilon }\mathrm{f}\left(\mathsf{\epsilon }\right)}{\mathsf{\rho }}_{\mathrm{a}}\left(\mathsf{\epsilon }\right)\mathrm{d}\mathsf{\epsilon }-{2\mathsf{\epsilon }}_{\mathrm{a}}$represents the hybridization energy due to the overlap between the adsorbent orbitals and the electronic states of the electrode, the second term ${\mathrm{E}}_{\mathrm{ortho}}=-2\left({\mathrm{n}}_{\mathrm{a}}+\mathrm{f}\right)\mathrm{S}\left|{\mathrm{V}}_{\mathrm{ak}}\right|$ is the orthogonalization energy (Pauli repulsion),

where f(ε) is the Fermi function, ε_a is the energy level of the adsorbent valence orbital, n_a is the filling of the adsorbent valence orbital determined by the expression ${\mathrm{n}}_{\mathrm{a}}={\displaystyle \underset{-\infty }{\overset{+\infty }{\int }}\mathrm{f}\left(\mathsf{\epsilon }\right)}{\mathsf{\rho }}_{\mathrm{a}}\mathrm{d}\mathsf{\epsilon }$, and $\mathrm{f}={\displaystyle \underset{-\infty }{\overset{+\infty }{\int }}\mathrm{f}\left(\mathsf{\epsilon }\right)}{\mathsf{\rho }}_{\mathrm{el}}\mathrm{d}\mathsf{\epsilon }/{\displaystyle \underset{-\infty }{\overset{+\infty }{\int }}{\mathsf{\rho }}_{\mathrm{el}}\mathrm{d}\mathsf{\epsilon }}$ is the degree of filling of the electrode electronic level.

For simplicity, we assume that V is proportional to S=-α∙|V|, where α is a coefficient, and thus arrive at ${\mathrm{E}}_{\mathrm{ortho}}=2\left({\mathrm{n}}_{\mathrm{a}}+\mathrm{f}\right)\mathsf{\alpha }{\mathrm{V}}_{\mathrm{a}}^2$.

The density of states on the adsorbent is given by the expression

$${\mathsf{\rho }}_{\mathrm{a}}\left(\mathsf{\epsilon }\right)={\displaystyle \frac{1}{\mathsf{\pi }}}{\displaystyle \frac{{\mathrm{\Delta }}^2(\mathsf{\epsilon })}{{\left(\mathsf{\epsilon }-{\mathsf{\epsilon }}_{\mathrm{a}}-\mathrm{\Lambda }(\mathsf{\epsilon })\right)}^2+{\mathrm{\Delta }}^2(\mathsf{\epsilon })}}$$

(7)

where Δ and Λ are the chemisorption functions described by the expressions [76] $\mathsf{\Delta }\left(\mathsf{\epsilon }\right)=\mathsf{\pi }{\displaystyle \sum _{\mathrm{k}}{\left|{\mathrm{V}}_{\mathrm{k}}\right|}^2\mathsf{\delta }\left(\mathsf{\epsilon }-{\mathsf{\epsilon }}_{\mathrm{k}}\right)}$, and Λ is the Hilbert transform of Δ: $\mathrm{\Lambda }\left(\mathsf{\epsilon }\right)={\displaystyle \frac{1}{\mathsf{\pi }}}\mathrm{P}{\displaystyle \int \frac{\left|\mathsf{\Delta }(\mathsf{\epsilon }\prime )\right|}{\mathsf{\epsilon }-\mathsf{\epsilon }\prime }}$, where P denotes the principal part. For simplicity, instead of the sum of values for all the k-vectors, we represent Δ(ε) in the form $\mathsf{\Delta }\left(\mathsf{\epsilon }\right)\approx \mathsf{\pi }{\mathsf{\rho }}_{\mathrm{el}}\left(\mathsf{\epsilon }\right){\left|\mathrm{V}\right|}^2$. [5,72,77]

The adsorbent nature in the NA model, eq. (6−7), is given directly only by parameters ε_a and V [2,5] S is indirectly related to the matrix element. In case of the ions, as the data show (Supplementary Material of ref. [72]) for O₂ and intermediates of oxygen reduction OH and OOH, the matrix element increases significantly due to the contribution of the electron-nuclear interaction. That is why we will evaluate the influence of, first of all, V [2] and then ε_a on E_ads.

Another way to estimate adsorption is to represent the effective orbital of an electrode as an s-function distributed along the z-axis perpendicular to the electrode, as in the jellium model: [78,79]

ψ_el=A_el∙exp(-α_el∙z[2]) (8)

Such approach has been used in numerous works [71,72,80,81,82,83,84] to estimate the matrix element of the overlap between the electrode and redox system orbitals. In the following calculations, α_el was set at 0.01 a.u.^- [2] and the normalization factor A_el in (8) was calculated at the electrode surface (z=0) and at some distance from it using the expression $1={\displaystyle \underset{\mathrm{z}=0}{\overset{\mathrm{z}=+\infty }{\int }}{\mathsf{\psi }}_{\mathrm{el}}^2}$.

In this method, the ψ_el function substitutes for ψ_solv in expressions (1−3). After that, all further calculations are made using equations (4, 5), and the final expression for the adsorption energy can be rewritten as

$${\mathrm{E}}_{\mathrm{ads}}{=2\mathsf{\epsilon }}_{+}-{\mathsf{\epsilon }}_{\mathrm{el}}-{\mathsf{\epsilon }}_{\mathrm{ion}}$$

(9)

where $\mathrm{S}=〈{\mathsf{\psi }}_{\mathrm{solv}}|{\mathsf{\psi }}_{\mathrm{ion}}〉$,${\mathsf{\epsilon }}_{\mathrm{el}}=〈{\mathsf{\psi }}_{\mathrm{el}}\left|\stackrel{\wedge }{\mathrm{H}}\right|{\mathsf{\psi }}_{\mathrm{el}}〉$,${\mathsf{\epsilon }}_{\mathrm{ion}}=〈{\mathsf{\psi }}_{\mathrm{ion}}\left|\stackrel{\wedge }{\mathrm{H}}\right|{\mathsf{\psi }}_{\mathrm{ion}}〉$,$\mathrm{V}=〈{\mathsf{\psi }}_{\mathrm{el}}\left|\stackrel{\wedge }{\mathrm{H}}\right|{\mathsf{\psi }}_{\mathrm{ion}}〉$ and the ground state energy ${\mathsf{\epsilon }}_{+}={\displaystyle \frac{\left({\mathsf{\epsilon }}_{\mathrm{el}}{+\mathsf{\epsilon }}_{\mathrm{ion}}\right)-\sqrt{4{\mathrm{V}}^2+\left(1-{\mathrm{S}}^2\right){\left({\mathsf{\epsilon }}_{\mathrm{el}}-{\mathsf{\epsilon }}_{\mathrm{ion}}\right)}^2}}{2\left(1-{\mathrm{S}}^2\right)}}$.

The results of the calculations of binding of the orbitals of cations and anions with the solvent orbital using the model of eq. (1−5) give a set of graphs that fully or partially form two volcano-shaped dependencies (Figure 1 A,B). The graphs show that the orbital energies for the anions are mostly lower than those for the cations. It means that the centers of the volcano-shaped dependencies on the energy scale are shifted to the left for the anions and to the right for the cations.

In the regions of decreasing or increasing E_bind values, the graphs are quite close to linear or exponential dependencies. In case of s/s (schematic notation: solvent molecule orbital/ion orbital), s/p_z and ${\mathrm{s}/\mathrm{d}}_{{\mathrm{z}}^2}$cation orbital interactions, the maximum region becomes sloped, which is caused by a large overlap of the orbitals due to their shape (σ-bond formation), lowering the binding energy value. When the solvent s-orbitals are combined with the ${\mathrm{d}}_{\mathrm{xy}}$and ${\mathrm{d}}_{\mathrm{xz}}$ orbitals of the ions, rather weak π-bonds are formed. In case of the solvent p-orbitals, similar trends are observed for p_x/p_z and p_z/${\mathrm{d}}_{{\mathrm{z}}^2}$ and, to a lesser extent, for p_z/d_xz and p_x/p_x. The steeper slope of p_z/d_xz and p_x/p_x compared to s/p_x, etc. is explained by the larger extent of the p-orbital compared to the s-orbital, which results in a strong overlap and a decrease in the E_bind value. The specific shape of p_z/${\mathrm{d}}_{{\mathrm{z}}^2}$ is caused by the orientation of p_z-orbital directly towards the point charge of the ion. As a result, the contribution of the H_ne operator makes the bond stronger.

The most interesting trend is observed when the p_z-orbital of the anion participates in the binding, which causes a split in the left side of the branch (Figure 1 A,B). A special case are the results for the combination of s/s, and s/${\mathrm{d}}_{{\mathrm{z}}^2}$ orbitals at q=-1; p_x/p_x and p_x/p_y at q= -1 and p_x/p_x, p_x/p_y, and p_z/d_xy at q= +1, for which stable solutions are obtained only for one of the branches, the shape of which is close to a linear function, and in some cases (s/${\mathrm{d}}_{{\mathrm{z}}^2}$) to an exponential one.

The dependence of E_ads(ε_a) according to the NA model, eq. (6, 7), is almost linear and does not depend significantly on coefficient α (Figure 1C), while in case of E_ads(V [2]) an exponential dependence is observed (Figure 1D). In the region of moderately high values V [2]> 0.5 eV [2], the dependence of E_ads(V [2]) is almost linear. Notably, the approximation of specific results of DFT calculations of DOS adsorbents, for example, [74,75] shows that each bond type or band involved in the adsorption can have its own set of parameters ε_a and ε_c (band center). Therefore, for sufficiently high adsorption energies and, respectively, high V [2], the shape of the E_ads(V [2]) or E_ads(ε_a) dependencies may become more complicated due to the contribution of a larger number of electrode bands and adsorbent orbitals involved in the adsorption process.

The calculations for the ions carried out by eq. (9) for the ${\mathrm{p}}_{\mathrm{x}}$, ${\mathrm{d}}_{\mathrm{xy}}$, ${\mathrm{d}}_{{\mathrm{x}}^2{-\mathrm{y}}^2}$orbitals lead to two volcano-shaped dependencies (Figure 1E). Their side branches are close to linear dependencies and the interaction of the orbitals with the electrode is rather weak in this case. In case of the ${\mathrm{p}}_{\mathrm{z}}$ and${\mathrm{d}}_{{\mathrm{z}}^2}$ orbitals, a quickly decreasing exponential dependence is observed. For the cations involving these orbitals, the E_ads dependence represents a function that decreases as the ε_ion value becomes higher. In general, this result is very close to the data obtained by the first model by eq. (1−5).

Since the HOMO or LUMO orbitals of real ions represent a superposition of a number of MOs, the total contribution of such composite orbitals to the final adsorption energy or energy of binding can become more complicated. If we consider the interaction between the solvent molecules and the ions, the dependence of the binding energy on the ion energy level becomes volcano-shaped. When the p-orbitals of the ion are oriented along the z-axis (perpendicular to the MO of the solvent molecule), the left side of the volcano-shaped E_bind dependence splits into two branches. In addition to the splitting, a scenario of parallel dependencies caused by the simultaneous contribution of multiple MOs of the ion is also possible. Within the NA model, the adsorption dependencies can be linear or exponential. However, if we consider the adsorption on the electrode through the effective wave function by eq. (8, 9), we can expect the same trends as those obtained in the first model by eq. (1−5), i.e. formation of volcano-like dependencies or an exponential or linear decrease/increase in the adsorption energy.

3.2. Trends in Ion Adsorption on α-Al₂O₃ (0001), Graphene and Au (111)

Summarizing the data on the adsorption of ions (Figure S1.1−S1.12, Table S2.3 in Supplementary Information) and their binding into complexes (Figure S2.1−S2.4, Table S3.3), we can conclude that it is most convenient to evaluate the energy levels of the HOMO (ε_HOMO) and LUMO (ε_LUMO) orbitals (Table S4.1, S4.2). In some cases, second highest (HOMO-1) and second lowest (LUMO+1) molecular orbitals can also be used. In addition, there are descriptors that combine HOMO and LUMO with HOMO-1 and LUMO+1 orbitals, e.g. hyperhardness or hypersoftness. [21] In case of strong correlations, combining parameters into descriptors may well make practical sense in some problems. Working directly with ε_HOMO and ε_LUMO levels allows us to describe the adsorption and distance from the electrode to the ion, avoiding additional calculations of descriptors made up of number of parameters. Therefore, we will mainly focus on the analysis of trends in adsorption and binding depending on one orbital only.

In case of Al₂O₃, quantum chemical calculations show that the value of the cation adsorption energy decreases exponentially (Figure 2A) as the ε_HOMO value becomes lower and almost linearly when the ε_LUMO value goes down (Figure 2B). The strongest binding to the surface is observed for alkali metal cations and H₃O⁺ and the weakest binding for EMI⁺. The distance from the electrode surface to the cation also follows the E_ads trends, where z₀ becomes lower as the ε_HOMO (Figure 2A) and ε_LUMO (Figure 2B) values decrease. Only H₃O⁺ is slightly out of the general E_ads(ε_HOMO) trend, which may be caused by the formation of strong adsorption bonds between the hydrogen of H₃O⁺ and the oxygen of the Al₂O₃ surface.

In case of z₀(ε_LUMO), only PyMe⁺ does not follow the linear trend. Since the z₀ distance is determined at the ion center of mass, the elongated shape of the PyMe molecule and, consequently, the orbitals may cause overestimation of the z₀ value.

No correlation is observed between the anion adsorption and ε_HOMO (Figure 2E), while in case of E_ads(ε_LUMO), there is a volcano-like dependence with a maximum at ε_LUMO = 0.245 a.u. (Figure 2F). The strongest binding to the surface is observed for the RCOO^-, OH^-, and F^- systems, while that for the PF₆^- anion is the weakest. OH^- and F^- are out of the E_ads(ε_LUMO) trend, with the adsorption energies for the above ions getting as high as -4.0 and -4.3 eV, respectively; the adsorption energy in these cases is close to that of a strong chemical bond with the Al atoms.

The shape of z₀(ε_HOMO) and z₀(ε_LUMO) becomes more complicated (Figure 2G,H), with the dependencies splitting into two curves that are reasonably well approximated by linear and exponential functions. In case of z₀(ε_HOMO), the linear dependence is formed by PF₆^-, BF₄^- and F^-, while the exponential one is represented by the other anions (TFSI^-, TfO^-; HSO₄^-, MeSO₃^-; ClO₄^-, NO₃^-, OH^-; HCOO^-, AcO^-, EtCOO^- and Cl^-, Br^-). z₀(ε_LUMO) also splits into two dependencies, with the first linear one consisting of OH^-; HSO₄^-, MeSO₃^-; HCOO^-, AcO^-, and EtCOO^- and the second exponential one represented by the TFSI^-, TfO^-, BF₄^-, PF₆^-; ClO₄^-, and NO₃^- anions and halides. In all the cases, z₀ decreases smoothly as the HOMO and LUMO energy values increase.

As the analysis (Figure S1.13) shows, the distance from the ions to the electrode surface for a number of systems correlates with the adsorption energy. For the anions, there is no E_ads(z₀) correlation, while for the cations a linear dependence is observed (Figure S1.13 A, B).

It can be concluded that the adsorption of cations on aluminum oxide is the result of contributions of both the HOMO and the LUMO orbitals. The adsorption of the anions is determined by the LUMO, HOMO-1, and LUMO+1 levels (B, D in Figure S1.1, S1.2) and is independent of the HOMO one.

Let us continue our consideration of adsorption by studying a semiconductor as an example. As in the previous system, the alkali metal cations are strongly adsorbed on the graphene surface, the PyMe⁺ cation is weakly bound to the surface. The E_ads value for the cations on graphene decreases linearly as the ε_HOMO value becomes lower (Figure 3A) and in case of E_ads(ε_LUMO), the dependence (Figure 3B) is a volcano-shaped plot with an energy maximum in the region where ε_LUMO ~ -0.10 a.u. At the minimum distance from the electrode are the compact alkali metal cations, whereas the farthest from the surface are the voluminous TPA⁺ and TBA⁺ ions. As the ε_HOMO value becomes higher, the z₀ distance from the electrode increases exponentially (Figure 3C), and at an increase in the ε_LUMO value, z₀ grows linearly (Figure 3D).

The E_ads(ε_HOMO) dependence for the anions is almost linear; in case of E_ads(ε_LUMO) two branches are observed, the first linear one is represented by the ClO₄^-, NO₃^-, OH^-; RSO₃^-; TFSI^-, TfO^-, BF₄^-, PF₆^-; and RCOO^- anions and the second one consisting of halides is close to an exponential dependence. Only i-PrCOO^- is out of the trend. The presence of an i-Pr- voluminous group probably leads to a rather strong repulsion due to the van der Waals interactions, making the adsorption energy value for this ion positive. The lowest adsorption energy values are observed for BF₄^-, PF₆^- and NO₃^-, the highest energies are registered for some of the halides (Cl^-, Br^-), EtSO₃^- and i-PrCOO^-.

z₀(ε_HOMO) for the anions splits into two linear dependencies, the first one represented by the BF₄^-, PF₆^-, F^- and OH^- ions, the second one – by the other ions; i-PrCOO^- and p-TsO^- are noticeably out of the trend. In case of p-TsO^-, such spread of values may also be caused by the total volume and the associated asymmetry of the system due to the Ph- group, leading to an overestimation of the z₀ distances relative to the trend. The i-PrCOO^- distance from the electrode is also associated with the voluminous i-Pr group. The z₀(ε_LUMO) plot, as the one for the ions on Al₂O₃ (Figure 2H), is formed by a linear dependence consisting of the OH^-; HSO₄^-, MeSO₃^-, EtSO₃^-, NH₂SO₃^-, p-TsO^-; HCOO^-, AcO^-, EtCOO^-, PrCOO^-, and i-PrCOO^- ions and an exponential dependence consisting of TFSI^-, TfO^-, BF₄^-, PF₆^-; ClO₄^-, NO₃^- and halides.

As in case of Al₂O₃, the E_ads(z₀) dependence on graphene is observed only for the cations (Figure S1.13 C, D) and has a volcano shape with a maximum at z₀ ~ 4 Å and E_ads ~ -0.5 eV. According to the above dependencies (Figure 3) and the data of (Figure S1.5 B, D), the adsorption of the cations is determined by the HOMO orbital (except HOMO-1), as well as LUMO and LUMO+1. In case of the anions (Figure 3, Figure S1.7 B, D), the HOMO, HOMO-1 and LUMO orbitals are involved in the adsorption.

The strongest binding to the gold surface is observed for hydroxonium, which is followed by alkali metal cations and EMPyr⁺. PyMe⁺, as in case of graphene, is weakly bound to the surface. For the Au (111) surface, the voluminous organic PyMe⁺ and EMPyr⁺ cations drop out of the exponential dependence of E_ads(ε_HOMO) and E_ads(ε_LUMO) (Figure 4A,B). In addition, a significant shift towards negative energies is observed for H₃O⁺. In case of the alkali metal, TMA⁺ and EMI⁺ cations, the dependencies are almost exponential, the E_ads value in both cases decreases with an increase in ε_HOMO and ε_LUMO.

z₀ grows exponentially with an increase in ε_HOMO and ε_LUMO (Figure 4С,D). Quite far from the z₀(ε_HOMO) dependence are only TMA⁺, EMPyr⁺ and PyMe⁺. The alkali metal cations, H₃O⁺, and EMI⁺ fit quite well into the exponential dependence (Figure 4 С), and all the ions considered fit into the linear z₀(ε_LUMO) dependence (Figure 4 D).

Among the anions, the strongest interaction with the surface is observed for halides, which are followed by OH^- and HCOO^-. The weakest interaction with the (111) Au surface is observed for the nitrate ion, followed by TfO^- and MeSO₃^-. The anions form a characteristic volcano-shaped E_ads(ε_HOMO) dependence (Figure 4E) with a maximum at ε_HOMO ~ -0.25 a.u., while no correlation is found for E_ads(ε_LUMO) (Figure 4F). The z₀(ε_HOMO) dependence is represented by a linear relationship consisting of BF₄^-, PF₆^- and OH^- ions and an exponential relationship with the other systems (Figure 4G). Two dependencies are observed for z₀(ε_LUMO) (Figure 4H), with the first linear one formed by OH^-; HSO₄^-, MeSO₃^-, EtSO₃^-, NH₂SO₃^-, p-TsO^-; HCOO^-, AcO^-, EtCOO^-, PrCOO^-, and i-PrCOO^- ions and the second exponential one made up by TFSI^-, TfO^-, BF₄^-, PF₆^-; ClO₄^-, NO₃^- and halides.

Some degradation of the correlations on the metal and deviation of a number of cations (on average, cations have stronger binding to Au (111) than anions) can be caused by the increased complexity of the adsorption process of bulk systems on the metal. In addition to a bigger population of electronic states in the metal and, consequently, a bigger number of electrode wave functions involved in the adsorption, there are also local perturbations of the electron density in the vicinity of the ions. In particular, we can talk about image charges, [79,85,86] which can make a significant contribution to the adsorption energy [87] when the adsorbent approaches the metal electrode at a sufficiently close distance. [85]

Unlike the previous model electrodes, on Au (111) the E_ads(z₀) dependence is observed for both the cations (an exponential one, Figure S1.13 E, except for PyMe⁺, EMPyr⁺ and H₃O⁺) and the anions (a linear one, Figure S1.13 F).

To sum it up (Figure 4, Figure S1.9 B, D and Figure S1.11 B, D), the adsorption of the cations, except for some systems, and the anions on Au (111) is formed by both HOMO, HOMO-1 and LUMO, LUMO-1. Such conclusion is expectable considering the high density of electronic levels in the metal involved in the adsorption process.

The proposed descriptors characterize the adsorption parameters of cations on Al₂O₃ and graphene quite well (R [2]=0.64−0.86 for E_ads and R [2]=0.87−0.96 for z₀). On Au (111), the correlations for the cations look slightly worse because some organic cations fall out of the dependencies. The correlation coefficients of adsorption and distances for the cations increase as follows: Au (111) < graphene < Al₂O₃ (0001). In case of the anions, the correlations for E_ads are somewhat worse R [2]=0.35−0.55, but still acceptable for qualitative and semi-quantitative evaluations. And the correlations for the z₀ parameter are approximately at the same level R [2]=0.64−0.95 as those for the cations.

Most of the dependencies on the orbital energy are exponential, in some cases they are close to the linear type, as the theoretical estimates suggest (Figure 1). Linear dependencies are also reported for the adsorption of neutral molecules on metal oxides in refs. [32,33,34] Volcano-shaped dependencies are observed for the anions on Al₂O₃ (Figure 2F) and cations on graphene (Figure 3B), which is also in agreement with the theory. The splitting of the dependencies for E_ads into two branches is observed only for the anions on graphene (Figure 3F). There is a split in all the z₀ anion dependencies on ε_HOMO and ε_LUMO. One cause of the division of the anions into two groups may be the presence of several MOs with the same energy (Table S6.1). This effect will be discussed in more detail in the next paragraph.

Some limitations on the applicability of ε_HOMO and ε_LUMO as adsorption descriptors may arise in case of voluminous organic ions whose interaction with the electrode is poorly described by the MO energy level due to steric perturbations. Another limiting factor may be strong chemisorption involving the formation of a chemical bond with the surface, comparable in energy to the interatomic bond in a solid. In general, the descriptors provided in the manuscript are best suited for characterizing semiconductors and dielectrics.

3.3. Ion Binding in Complexes with Water and Dimethyl Carbonate

The free energy of cation binding into an aqueous complex is described well by an exponential function of ε_HOMO (Figure 5A), the TBA⁺ cation falls out of the trend, and a polynomial dependence is observed for G_bind(ε_LUMO) (Figure 5B). The distance from the center of mass of the cation and water is reasonably well described by the exponential r₀(ε_HOMO) and linear r₀(ε_LUMO) functions (Figure 5C,D). Only the pyridinium derivatives (RPy⁺) fall out of the r₀(ε_LUMO) trend; the distance between the centers of mass for RPy⁺ according to quantum chemical calculations is almost 1 Å longer than the obtained trend. The cause of this discrepancy may be greater asymmetry of RPy⁺ compared to the other systems. The lowest binding energies and shortest distances are registered for the alkali metal cations, while the highest values are observed for the organic cations.

For the anions, the dependence of G_bind on ε_HOMO is divided into linear and exponential curves (Figure 5E), while there is no correlation between G_bind and ε_LUMO (Figure 5F). The lower linear G_bind(ε_HOMO) correlation (Figure 5E) is formed by BF₄^-, PF₆^-, F^- and OH^- ions, the upper one is made up by ClO₄^-, NO₃^-, DCA^-; TFSI^-, TfO^-; Cl^-, Br^-; RSO₃^-; and RCOO^- ions. The ion-water distances (Figure 5G, H) depend only on ε_LUMO, and no correlations are observed for ε_HOMO. The linear dependence of r_o(ε_LUMO) is represented by OH^-; HSO₄^-, MeSO₃^-, EtSO₃^-, NH₂SO₃^-, p-TsO^-; HCOO^-, AcO^-, EtCOO^-, PrCOO^-, and i-PrCOO^- ions, the exponential one – by TFSI^-, TfO^-, BF₄^-, PF₆^-; ClO₄^-, NO₃^-; F^-, Cl^-, and Br^- ions. Only DCA^- deviates slightly from the exponential trend. The lowest values of binding energies and distances are observed for F^- and OH^-, while the highest values are observed for the organic cations. The highest values of binding energies are observed for fluorine-containing anions TFSI^-, TfO^-, BF₄^-, and PF₆^-, as well as DCA^- and RSO₃^-. The longest r₀ distances are found for BF₄^-, PF₆^- and DCA^-.

For the cations, G_bind is exponentially related to r₀ (Figure S2.9 A), while for the anions there is no such correlation (Figure S2.9 B).

It can be concluded that the HOMO, LUMO (Figure 5A,B) and HOMO-1, LUMO+1 orbitals (Figure S2.1. B, D) are actively involved in the cation binding to water. The anionic aqueous complex is only made up of HOMO and HOMO-1 orbitals (Figure 5E, Figure S2.3 C), but not the LUMO orbitals.

The dependence of the free energy of cation binding to DMC on ε_HOMO is approximated by an exponential function (Figure 6A) and on ε_LUMO by a polynomial function with a maximum at -0.09 a.u. (Figure 6B). The distances from the center of mass of the cation and DMC grow exponentially depending on ε_HOMO and inversely exponentially depending on ε_LUMO (Figure 6C,D). Only PyMe⁺ and EMI⁺ drop out of the r₀(ε_LUMO) dependence, with the predicted r₀ values tending to be almost 1 Å lower than the calculated r₀ ones. The strongest binding to DMC is observed for the alkali metal cations, while the weakest interaction is registered for TEA⁺ and EMPyr⁺.

The plot of the anion binding energy vs ε_HOMO represents two linear dependencies (Figure 6E), the first one consisting of PF₆^-, BF₄^- and F^-, the second one — of ClO₄^-, NO₃^-; MeSO₃^-; TfO^-, TFSI^-, halides and RCOO^- ions. G_bind(ε_LUMO) is formed (Figure 6F) by a linear relationship composed of MeSO₃^- anions, halides and RCOO^- and an exponential one consisting of ClO₄^-, NO₃^-; PF₆^-, BF₄^-, TfO^-, and TFSI^-. Like G_bind(ε_HOMO), r₀(ε_HOMO) for the anions also represents two dependencies (Figure 6G), while r₀(ε_LUMO) is a linear relationship (Figure 6H). The strongest interaction with the solvent is observed for F^- and the weakest one — for the TFSI^- and PF₆^- anions.

For the cations, the G_bind(r₀) data are approximated by an exponential function (Figure S2.9 C) and for the anions — by a linear function (Figure S2.9 D).

The HOMO, LUMO (Figure 6A,B), HOMO-1, and LUMO+1 orbitals (Figure S2.6B,D) are actively involved in the binding of cations to DMC. However, the formation of the anionic complex with DMC involves only the HOMO, LUMO (Figure 6E,B), and HOMO-1 (Figure S2.7 B) orbitals.

In addition to describing the adsorption on the electrode, the energy level of the HOMO and LUMO is an excellent descriptor of the energy of binding into a complex (R [2]=0.70−0.95) and the distance between the solvent molecule and the ions (R [2]=0.36−0.98). In almost all the cases, the correlation coefficients are much better for the complexes than for the electrodes. Reducing the number of orbitals involved in the binding makes the complex formation easier and closer to theoretical models (paragraph 1).

The splitting of the dependencies is mainly due to the orbital type of the ion involved in the complex formation. For example, the dependence of the energy of solvent binding to anions on the ε_HOMO parameter (Figure 5E, Figure 6E) is split into two curves: one consisting of OH^-, F^-, BF₄^-, and PF₆^- and the second one comprised of the other ions. The MO analysis shows (Table S6.1) that the ions from the first dependence possess multiple HOMOs of the same energy. A similar division into two dependencies can be seen in z₀(ε_HOMO), for example, for the complexes with DMC (Figure 6G) and for the adsorption of anions on Al₂O₃ (0001) (Figure 2G), graphene (Figure 3G), and Au (111) (Figure 4G). The z₀(ε_LUMO) plots are also represented by two dependencies: a linear one consisting of OH^-; RSO₃^-; and RCOO^- and an exponential one composed of ClO₄^-, NO₃^-; BF₄^-, PF₆^-, TfO^-, TFSI^-, and halides. There is a similar trend for the adsorption on electrodes (Figures 2−4H) and for complexation with water (Figure 5H), while it is not observed for DMC (Figure 6H). Theoretical evaluations can shed light on the nature of such dependencies by establishing the distance between the ion and the electrode or solvent molecule at which the total energy of the system would be minimized. Taking into account the available data, we can assume that the division is caused by the specific composition of the LUMO orbital and does not depend on the number of LUMO orbitals with the same energy. Taking into account the data on DMC, it is quite natural to conclude that the molecular orbitals of the solvent molecule influence the complex formation.

3.4. Correlations with Physicochemical and Electrochemical Parameters of Ions

Judging by the dependence for dipole moment [88] and the more complicated one for polarization [89] (since the perturbation of the system by the external field is calculated), the wave functions are explicitly included in all the expressions of these parameters. In several models, the wave functions are also included in some terms of the solvation energy expression. [63] Notably, in addition to a number of summands responsible for the environment polarization, [63,64] cavity formation and other contributions, [90] the solvation energy is determined by the first solvation shell, [91,92] which is formed by the solvent directly bound to the ion. As paragraph 3 shows, the energy of binding to water and the intermolecular distances are related to the HOMO and LUMO. The redox potential, in turn, contains the solvation energies of the reactants [93,94] in the original entry and, therefore, depends on the internal rearrangement of the bonds and the environment response to the charge. Such facts and some literature data [35] suggest a possible correlation of the mentioned parameters with the orbital energy level. In general, the mentioned physicochemical parameters, as in works, [26,27] are likely to depend on a whole range of descriptors, but for simplicity of the analysis we will focus on the correlations with MO levels only..

According to the calculations (Figure 7A), the dependence of the dipole moment for the cations on ε_HOMO is nearly sigmoidal, for compact one-atom ions (alkali metals) and symmetric ions (TMA⁺) in the range ε_HOMO = -3 – -0.7 a.u. the dipole moment is equal to zero. Starting from ε_HOMO > -0.62 a.u. a sharp increase in the dipole moment value is observed, with the maximum values corresponding to the voluminous and asymmetric CP⁺, CTA⁺ and DMI⁺ cations. The dependence of the dipole moment of anions on ε_LUMO (Figure 7B) decreases exponentially in the range ε_LUMO = 0.08 − 0.25 a.u., the maximum dipole moment is observed for voluminous DS^-, followed by BuCOO^-, p-TsO^-. At ε_LUMO > 0.25 a.u. the dipole moment of the anions becomes close to zero (for halides and symmetric NO₃^-, ClO₄^-, etc.).

The polarization dependence on ε_HOMO for the cations (Figure 7C) is split into two curves: an exponential one represented by H₃O⁺, alkali metal cations, RMPyr⁺, and R₄N⁺, and a linear one represented by CTA⁺, CP⁺, RMI⁺ (except PhMI⁺, which is out of the trend), and RPy⁺. The dependence of P on ε_LUMO, as the one for the dipole moment, decreases exponentially (Figure 7D), with TFSI^- falling out of the trend. The maximum polarization values are obtained for the voluminous and heavy DS^- and p-TsO^- cations, the lowest ones for the light F^-, OH^- and BF₄^- ions.

The value of the solvation energy for the cations decreases exponentially with an increase in the ε_HOMO or ε_LUMO values (Figure 7E,F). The lowest solvation values are observed for Li⁺, which is followed by the other alkali metals, with the maximum values obtained for the organic cations. The G_solv(ε_HOMO-1) dependence for the anions has a characteristic volcano shape with a maximum at ε_HOMO-1 = -0.42 a.u. (Figure 7G). The left part of the dependence is formed by the halides, the right part – by the other ions. A significant deviation from the trend is observed for OH^-, the energy value of which is significantly underestimated. TFSI^-, the solvation energy of which is slightly overestimated, is out of the trend. The lowest values of G_solv are observed for F^-, OH^- and the highest values for TFSI^-.

A dependence of the oxidation potentials on ε_HOMO is only observed for the anions (Figure 7H). E_SHE decreases exponentially as ε_HOMO goes up.

The present paragraph clearly shows that in addition to describing the adsorption of ions or their binding to the solvent, the HOMO and LUMO parameters are also well suited for characterizing a number of physicochemical and electrochemical parameters. The latter, first of all, include dipole moment and polarization, whose correlation coefficients of the approximations reach R [2]=0.8−0.94 and the cation solvation energy R [2]=0.95. Slightly worse are R [2]=0.74−0.78, but there are clear dependencies on the solvation energy and oxidation potential of anions.

3.5. Practical Evaluations

In this paragraph, we will show how orbital values can be applied to solving a practical problem. For this purpose, in continuation of the study of aqueous electrolytes, [28,44,45,46,47,50,51] let us consider the ability of cations and anions to displace water from the inner layer of the electric double layer (EDL). The point of this task is to displace the electroactive component of the electrolyte water as far away from the electrode surface as possible [28,44,45,46] and, as a result, to increase the electrochemical potential window and energy capacity of the batteries. Since such problem considers a potential electrolyte for a battery, we will use a Li⁺ ion as the cation. The parameters of the anions will be calculated using the approximations obtained previously (Table S2.3, S3.3, S5.3). We will use graphene as the electrode because it can serve as a simplified model system representing different carbon materials for supercapacitors [95,96] and batteries. [97,98] By varying the ε_HOMO and ε_LUMO parameters in the range characteristic of anions, we will find systems with the minimal water concentration in the inner layer of EDL.

Under equilibrium conditions, the reversible process of displacement of ions (i) by water molecules (w) in the inner layer from the surface into the bulk ${\mathrm{i}}^{\mathrm{b}\mathrm{u}\mathrm{l}\mathrm{k}}+{\mathrm{w}}^{\mathrm{s}\mathrm{u}\mathrm{r}\mathrm{f}}\rightleftarrows {\mathrm{i}}^{\mathrm{s}\mathrm{u}\mathrm{r}\mathrm{f}}+{\mathrm{n}}_{\mathrm{i}}{\mathrm{w}}^{\mathrm{b}\mathrm{u}\mathrm{l}\mathrm{k}}$is described by a multicomponent adsorption isotherm: [99,100,101]

$$\ln {\displaystyle \frac{{\mathrm{\theta }}_{\mathrm{i}}}{{\mathrm{\theta }}_{\mathrm{W}}^{{\mathrm{n}}_{\mathrm{i}}}}}+\ln {\displaystyle \frac{{\mathrm{f}}_{\mathrm{i}}^{\mathrm{s}\mathrm{u}\mathrm{r}\mathrm{f}}}{{\left({\mathrm{f}}_{\mathrm{w}}^{\mathrm{s}\mathrm{u}\mathrm{r}\mathrm{f}}\right)}^{{\mathrm{n}}_{\mathrm{i}}}}}=\ln \left\{{\mathrm{\beta }}_{\mathrm{i}}{\displaystyle \frac{{\mathrm{f}}_{\mathrm{i}}^{\mathrm{b}\mathrm{u}\mathrm{l}\mathrm{k}}{\mathrm{с}}_{\mathrm{i}}^{\mathrm{b}\mathrm{u}\mathrm{l}\mathrm{k}}}{{\left({\mathrm{f}}_{\mathrm{w}}^{\mathrm{b}\mathrm{u}\mathrm{l}\mathrm{k}}{\mathrm{с}}_{\mathrm{w}}^{\mathrm{b}\mathrm{u}\mathrm{l}\mathrm{k}}\right)}^{{\mathrm{n}}_{\mathrm{i}}}}}\right\}-{\displaystyle \frac{{\mathrm{q}}_{\mathrm{i}}{\mathrm{e}}_0{\mathrm{\phi }}_{\mathrm{d}}}{{\mathrm{k}}_{\mathrm{B}}\mathrm{T}}}-\left\{{\displaystyle \frac{{\mathrm{q}}_{\mathrm{i}}{\mathrm{e}}_0\left({\mathrm{l}}_{\mathrm{e}\mathrm{f}\mathrm{f}}-{\mathrm{r}}_{\mathrm{i}}\right)}{\mathrm{l}{\mathrm{k}}_{\mathrm{B}}\mathrm{T}}}+{\mathrm{b}}_{\mathrm{i}}-{\mathrm{n}}_{\mathrm{i}}{\mathrm{b}}_{\mathrm{w}}\right\}\mathrm{\Delta }\mathrm{\phi }+\left({\mathrm{c}}_{\mathrm{i}}-{\mathrm{n}}_{\mathrm{i}}{\mathrm{c}}_{\mathrm{w}}\right)\mathrm{\Delta }{\mathrm{\phi }}^2,$$

(10)

and parameters

$${\mathrm{b}}_{\mathrm{j}}={\mathrm{\mu }}_{\mathrm{j}}/{\mathrm{l}}_{\mathrm{e}\mathrm{f}\mathrm{f}}{\mathrm{k}}_{\mathrm{B}}\mathrm{T}(\mathrm{j}=\mathrm{i},\mathrm{w})$$

(11)

$${\mathrm{с}}_{\mathrm{j}}={\mathrm{\epsilon }}_0{\mathrm{\epsilon }}_{\mathrm{j}}{\mathrm{A}}_{\mathrm{j}}/2{\mathrm{l}}_{\mathrm{e}\mathrm{f}\mathrm{f}}{\mathrm{k}}_{\mathrm{B}}\mathrm{T},$$

(12)

where θ_i and θ_w are the degrees of surface coverage by ions and water molecules, $\mathrm{\Delta }\mathrm{\phi }={\displaystyle \frac{\left({\mathrm{\sigma }}_0+{\mathrm{\sigma }}_{\mathrm{i}\mathrm{o}\mathrm{n}\mathrm{s}}\right){\mathrm{l}}_{\mathrm{e}\mathrm{f}\mathrm{f}}}{{\mathrm{\epsilon }}_{\mathrm{e}\mathrm{f}\mathrm{f}}{\mathrm{\epsilon }}_0}}$ is the potential drop in the inner layer of the EDL, σ₀ is the surface charge density on the electrode, ${\mathrm{\sigma }}_{\mathrm{i}\mathrm{o}\mathrm{n}\mathrm{s}}={\displaystyle \sum _{\mathrm{i}}\frac{{\mathrm{\theta }}_{\mathrm{i}}{\mathrm{q}}_{\mathrm{i}}{\mathrm{e}}_0}{\mathrm{\pi }{\mathrm{r}}_{\mathrm{i}}^2}}$ is the surface charge density of adsorbed ions in the inner layer, ${\mathrm{\phi }}_{\mathrm{d}}={\displaystyle \frac{2\mathrm{R}\mathrm{T}}{\mathrm{F}}}\arcsin \mathrm{h}\left({\displaystyle \frac{-({\mathrm{\sigma }}_0+{\mathrm{\sigma }}_{\mathrm{i}\mathrm{o}\mathrm{n}\mathrm{s}})}{\sqrt{8\mathrm{R}\mathrm{T}{\mathrm{\epsilon }}_{\mathrm{b}}{\mathrm{\epsilon }}_0{\mathrm{c}}_0\cdot 1000}}}\right)$is the potential drop in the diffuse layer according to the Gui-Chapman model, [102] ε_b is the dielectric permittivity of the solution in the bulk (since dilute aqueous solutions are considered, the value is approximated as 78 [103]), c₀ is the concentration of the supporting electrolyte in mol/L, q_i is the i ion charge, ${\mathrm{n}}_{\mathrm{i}}={\mathrm{r}}_{\mathrm{i}}^2/{\mathrm{r}}_{\mathrm{w}}^2$is the ratio of areas of the i ion and water molecule, r_i and r_w are the effective radius of the ion or water molecule. [28]

The parameters included in coefficients b_j and c_j (11,12) have the following physical meaning: μ_j is the dipole moment perpendicular to the electrode surface; ${\mathrm{\epsilon }}_{\mathrm{j}}={1+{\mathrm{P}}_{\mathrm{j}}/{\mathrm{\epsilon }}_0{\mathrm{A}}_{\mathrm{j}}\mathrm{l}}_{\mathrm{e}\mathrm{f}\mathrm{f}}$is the dielectric constant of component j, calculated according to the model described in ref. [104] (P_j is the polarizability); ${\mathrm{A}}_{\mathrm{j}}=2\mathrm{\pi }{\mathrm{r}}_{\mathrm{j}}^2$is the area occupied by component j; $\mathrm{\epsilon }={\displaystyle \sum _{\mathrm{j}=\mathrm{i},\mathrm{w}}{\mathrm{\epsilon }}_{\mathrm{j}}{\mathrm{\theta }}_{\mathrm{j}}}$and ${\mathrm{l}}_{\mathrm{e}\mathrm{f}\mathrm{f}}={\displaystyle \sum _{\mathrm{j}=\mathrm{i},\mathrm{w}}{\mathrm{z}}_{\mathrm{j}}{\mathrm{\theta }}_{\mathrm{j}}}$are the effective values of the dielectric constant and thickness of the inner layer.

The ${\mathrm{f}}_{\mathrm{i}}^{\mathrm{s}\mathrm{u}\mathrm{r}\mathrm{f}}/{\left({\mathrm{f}}_{\mathrm{w}}^{\mathrm{s}\mathrm{u}\mathrm{r}\mathrm{f}}\right)}^{{\mathrm{n}}_{\mathrm{i}}}$ratio of the activity coefficients of the electrolyte components in the inner layer according to the mean field approximation [99,100] is generally determined by the expression

$$\begin{array}{l}\ln {\displaystyle \frac{{\mathrm{\theta }}_{\mathrm{m}}}{{{\mathrm{f}}_{\mathrm{W}}}^{{\mathrm{r}}_{\mathrm{m}}}}}=\left\{{\mathrm{\theta }}_{\mathrm{w}}\left(1-{\mathrm{\theta }}_{\mathrm{m}}\right)-{\mathrm{r}}_{\mathrm{m}}{\mathrm{\theta }}_{\mathrm{m}}\left(1-{\mathrm{\theta }}_{\mathrm{w}}\right)\right\}{\mathrm{A}}_{\mathrm{m}\mathrm{k}+1}+{\displaystyle \sum _{\mathrm{g}\ne \mathrm{m},\mathrm{k}+1}{\mathrm{\theta }}_{\mathrm{g}}\left(1-{\mathrm{\theta }}_{\mathrm{m}}+{\mathrm{r}}_{\mathrm{m}}{\mathrm{\theta }}_{\mathrm{m}}\right){\mathrm{A}}_{\mathrm{g}\mathrm{m}}}\\ -{\displaystyle \sum _{\mathrm{g}\ne \mathrm{m},\mathrm{k}+1}\left\{{\mathrm{\theta }}_{\mathrm{w}}{\mathrm{\theta }}_{\mathrm{g}}+{\mathrm{r}}_{\mathrm{m}}{\mathrm{\theta }}_{\mathrm{g}}\left(1-{\mathrm{\theta }}_{\mathrm{w}}\right)\right\}{\mathrm{A}}_{\mathrm{g}\mathrm{k}+1}}+{\displaystyle \sum _{\mathrm{g}<\mathrm{h}}{\mathrm{\theta }}_{\mathrm{g}}{\mathrm{\theta }}_{\mathrm{h}}\left({\mathrm{r}}_{\mathrm{m}}-1\right){\mathrm{A}}_{\mathrm{g}\mathrm{h}}}\end{array}$$

(13)

For a 1:1 electrolyte (A stands for cations, B designates anions) and a solvent, equation (13) looks as follows:

$$\ln {\displaystyle \frac{{\mathrm{f}}_{\mathrm{A}}}{{{\mathrm{f}}_{\mathrm{W}}}^{{\mathrm{n}}_{\mathrm{A}}}}}={\mathrm{\theta }}_{\mathrm{B}}\left\{1+\left({\mathrm{n}}_{\mathrm{A}}-1\right){\mathrm{\theta }}_{\mathrm{A}}\right\}{\mathrm{A}}^{\mathrm{A}\mathrm{B}}+\left\{{\mathrm{\theta }}_{\mathrm{W}}\left(1-{\mathrm{\theta }}_{\mathrm{A}}\right)-{\mathrm{\theta }}_{\mathrm{A}}\left(1-{\mathrm{\theta }}_{\mathrm{W}}\right){\mathrm{n}}_{\mathrm{A}}\right\}{\mathrm{A}}^{\mathrm{A}\mathrm{W}}-{\mathrm{\theta }}_{\mathrm{B}}\left\{{\mathrm{n}}_{\mathrm{A}}-\left({\mathrm{n}}_{\mathrm{A}}-1\right){\mathrm{\theta }}_{\mathrm{W}}\right\}{\mathrm{A}}^{\mathrm{B}\mathrm{W}}$$

(14)

$$\ln {\displaystyle \frac{{\mathrm{f}}_{\mathrm{B}}}{{{\mathrm{f}}_{\mathrm{W}}}^{{\mathrm{n}}_{\mathrm{B}}}}}={\mathrm{\theta }}_{\mathrm{A}}\left\{1+\left({\mathrm{n}}_{\mathrm{B}}-1\right){\mathrm{\theta }}_{\mathrm{B}}\right\}{\mathrm{A}}^{\mathrm{A}\mathrm{B}}+\left\{{\mathrm{\theta }}_{\mathrm{W}}\left(1-{\mathrm{\theta }}_{\mathrm{B}}\right)-{\mathrm{\theta }}_{\mathrm{B}}\left(1-{\mathrm{\theta }}_{\mathrm{W}}\right){\mathrm{n}}_{\mathrm{B}}\right\}{\mathrm{A}}^{\mathrm{B}\mathrm{W}}-{\mathrm{\theta }}_{\mathrm{A}}\left\{{\mathrm{n}}_{\mathrm{B}}-\left({\mathrm{n}}_{\mathrm{B}}-1\right){\mathrm{\theta }}_{\mathrm{W}}\right\}{\mathrm{A}}^{\mathrm{A}\mathrm{W}}$$

(15)

If k≠m, the A_km coefficients between components k and m are written as

$${\mathrm{A}}_{\mathrm{k}\mathrm{m}}={\mathrm{z}}_{\mathrm{C}\mathrm{N}}\left\{{\mathrm{E}}_{\mathrm{k}\mathrm{m}}-\left({\mathrm{E}}_{\mathrm{k}\mathrm{k}}+{\mathrm{E}}_{\mathrm{m}\mathrm{m}}\right)/2\right\}/{\mathrm{k}}_{\mathrm{B}}\mathrm{T}$$

(16)

where z_CN is the effective coordination number in the two-dimensional lattice, E is the energy of the interaction between the electrolyte components, when k=i and m=w: E_iw is the energy of binding of the ion and the water molecule (taken equal to the free binding energy of the components in the following estimations), E_ww is the energy of binding of the water molecules into a dimer, E_ii is the energy of binding of the ions with each other.

For simplicity, E_ii is represented as a Coulomb interaction of charges q_i separated from each other by the distance of their own effective ionic radii: ${\mathrm{E}}_{\mathrm{A}\mathrm{A}}={\displaystyle \frac{1}{4\mathrm{\pi }{\mathrm{\epsilon }}_0{\mathrm{\epsilon }}_{\mathrm{e}\mathrm{f}\mathrm{f}}}}{\displaystyle \frac{{\mathrm{q}}_{+}^2{\mathrm{e}}_0^2}{(2{\mathrm{r}}_{+})}}$ for two cations, ${\mathrm{E}}_{\mathrm{B}\mathrm{B}}={\displaystyle \frac{1}{4\mathrm{\pi }{\mathrm{\epsilon }}_0{\mathrm{\epsilon }}_{\mathrm{e}\mathrm{f}\mathrm{f}}}}{\displaystyle \frac{{\mathrm{q}}_{-}^2{\mathrm{e}}_0^2}{(2{\mathrm{r}}_{-})}}$ for two anions, and ${\mathrm{E}}_{\mathrm{A}\mathrm{B}}={\displaystyle \frac{1}{4\mathrm{\pi }{\mathrm{\epsilon }}_0{\mathrm{\epsilon }}_{\mathrm{e}\mathrm{f}\mathrm{f}}}}{\displaystyle \frac{{\mathrm{q}}_{+}{\mathrm{q}}_{-}{\mathrm{e}}_0^2}{({\mathrm{r}}_{+}+{\mathrm{r}}_{-})}}$ for one cation and one anion.

β_i is the parameter which determines the difference between the chemical potentials of the reaction components on the surface and in the bulk.

$$\ln \left({\mathrm{\beta }}_{\mathrm{i}}\right)=\left\{\left({\mathrm{\mu }}_{\mathrm{i}}^{0,\mathrm{b}\mathrm{u}\mathrm{l}\mathrm{k}}-{\mathrm{\mu }}_{\mathrm{i}}^{0,\mathrm{s}\mathrm{u}\mathrm{r}\mathrm{f}}\right)-{\mathrm{n}}_{\mathrm{i}}\left({\mathrm{\mu }}_{\mathrm{w}}^{0,\mathrm{b}\mathrm{u}\mathrm{l}\mathrm{k}}-{\mathrm{\mu }}_{\mathrm{w}}^{0,\mathrm{s}\mathrm{u}\mathrm{r}\mathrm{f}}\right)\right\}/{\mathrm{k}}_{\mathrm{B}}\mathrm{T}$$

(17)

To simplify the problem, let us write this parameter as a function of the adsorption energies of the components on the electrode

$$\ln \left({\mathrm{\beta }}_{\mathrm{i}}\right)=\left\{(-{\mathrm{E}}_{\mathrm{i},\mathrm{a}\mathrm{d}\mathrm{s}})-{\mathrm{n}}_{\mathrm{i}}(-{\mathrm{E}}_{\mathrm{w},\mathrm{a}\mathrm{d}\mathrm{s}})\right\}/{\mathrm{k}}_{\mathrm{B}}\mathrm{T}=\left({\mathrm{n}}_{\mathrm{i}}{\mathrm{E}}_{\mathrm{w},\mathrm{a}\mathrm{d}\mathrm{s}}-{\mathrm{E}}_{\mathrm{i},\mathrm{a}\mathrm{d}\mathrm{s}}\right)/{\mathrm{k}}_{\mathrm{B}}\mathrm{T}$$

(18)

As work [28] suggests, the effective radius of cations should be estimated taking into account the solvation water. Keeping that in mind, we used the following parameters as the initial data for the lithium cation with a correction for solvation water: E_ads = -1.27 eV, E_bind,w = -0.33 eV, ${\mathrm{r}}_{\mathrm{L}{\mathrm{i}}^{+}}$= 3.00 Å (including the solvation water (Table S7.1)), z₀ = 3.45 Å, μ and P are equal to zero; the parameters for the water molecule are E_ads = -0.11 eV, r_w = 1.41 Å, z₀ = 3.29 Å, μ = 2.19 D, P =7.21 Å [3], E_ww = -0.13 eV. We described the dipole moment and polarizability of all the anions by the expressions µ_D=0.55598+353.79894∙exp(-31.72148∙ε_LUMO) and P=194.4018∙exp(-ε_LUMO/0.03634) 12.32756∙ε_LUMO+7.63034. The effective coordination number was taken to be z_CN = 6.

In accordance with the specific features of the dependencies of adsorption of the anions and their interaction with water on ε_HOMO and ε_LUMO, we divided all the systems into three sets (Figure 8A–C).

The first set of systems (A) consists of F^-, BF₄^- and PF₆^-; the adsorption of these ions is described by the expressions: E_ads=1.93766∙ε_HOMO+0.27664 and z₀=1.54334-5.459∙ε_HOMO, the aqueous complex parameters – by the expressions G_bind=-3.43293∙ε_HOMO-1.47558 and r₀=6.13597-16.63311∙ε_LUMO. According to the obtained solutions (Figure 8A) to equations (10−18), the lowest water coverage degree (θ_w < 0.08) is observed for the BF₄^- and PF₆^- ions. The only practical limitation on using BF₄^- and PF₆^- as aqueous electrolytes is their slow hydrolysis, [105,106,107,108,109] which is why they have not been mentioned as salts for aqueous systems. Figure 8А does not show F^-, which lies much higher at ε_LUMO= 0.54 a.u., but this system also extends to the region of low values θ_w < 0.08. In contrast to BF₄^- and PF₆^-, the fluoride appears in studies [97,110] as a component of an aqueous electrolyte with a wide potential window.

The second set of systems (B) is represented by ions whose adsorption is described in a way similar to the systems in (Figure 8A), except for the parameter z₀=1.92102-7.83216∙ε_HOMO, and the aqueous complex parameters are expressed as G_bind=0.03092∙exp(-9.00895∙ε_HOMO)-0.46743 and r₀=6.13597-16.63311∙ε_LUMO. All the systems (Figure 8B) lie in the region θ_w < 0.08. Almost all of the above systems are reported in the literature as aqueous electrolytes with a wide potential window: TFSI^-, [44,45,46] TfO^-, [111] ClO₄^-, [112,113] and NO₃^-. [114,115] The exceptions are Cl^- and Br^- halides, whose electroactivity is quite high, especially that of Br^-, due to the relatively low oxidation potential, which limits their use in high voltage electrolytes. Only KCl is employed in supercapacitors but its use is limited too. [116,117,118]

The third set of systems (C) is represented by OH^-, HSO₄^-, CF₃COO^-, RSO₃^-, and RCOO^- anions whose adsorption is described by expressions similar to those in (Figure 8B), the parameters of interaction with water are described by the equations G_bind=0.03092∙exp(-9.00895∙ε_HOMO)-0.46743 and r₀=6.37246∙exp(-6.99633∙ε_LUMO)+2.3211. The lowest water coverage θ_w < 0.08 (Figure 8C) is observed for the CF₃COO^- system as well as for HSO₄^- and NH₂SO₃^-. CF₃COO^- is mentioned in refs. [119,120] In general, the application of such system is limited by their cost and toxicity. HSO₄^- in an aqueous medium is characterized by dissociation with a marked pH decrease, [121] which significantly narrows the potential window by shifting the cathodic region in the acidic medium. NH₂SO₃^- is electrochemically unstable as the experiments and calculations show. [28] The MeSO₃^-, EtSO₃^-, and HCOO^- systems lie in the moderate water coverage degree region 0.08 < θ_w < 0.16, AcO^- lies at the boundary of the region θ_w ~ 0.16. OH^- itself is quite electroactive and shifts the potential window; however, in concentrated alkalis the potential window slightly increases to ~1.6 V, [116,117,122,123] which is higher than in pure water (1.23 V [115]). Such systems are used in supercapacitors. [116,117,123] In practice, sulfonate is more often used in its acid form as an additive to electrolytes in electroplating and hydrometallurgy, [124] and as an anion in redox flow batteries. [125] EtSO₃^- has hardly ever been mentioned as an aqueous electrolyte in the literature, except for some references where it was mentioned as an electrolyte for zinc-ion batteries. [126] HCOO^- as an electrolyte was mentioned in ref. [127]; in general, alkaline salts of carboxylic acids are often reported as potential electrolytes. [28,127,128,129,130]

The low values of the electrode water coverage in electrolytes are in most cases consistent with the electrochemical stability of the system, confirming that the coverage is one of the main parameters [28] determining the potential window. In addition to the known systems, electrolytes with anions represented by sulfonic acid derivatives are potentially interesting as objects of further investigation. Of course, a comprehensive analysis of the ion suitability for electrolytes requires, in addition to the water coverage, evaluation of the electrochemical stability of the anion itself and its potential intermediates during electrochemical decomposition, as it was done in case of NH₂SO₃^-. [28]

The use of the ε_HOMO and ε_LUMO parameters makes it possible to quickly evaluate the water coverage of the electrode and select the optimal anions, as it has been shown by a particular example. Due to the dependencies of the other electrochemical parameters (solvation energy, redox potential, etc.) on MO energy levels, such analysis can be performed in a comprehensive way, taking into account all the necessary parameters of cations and anions at the same time.

4. Conclusions

In contrast to a lot of other works that study the effect of the electrode material or solvent type on the interactions of ions with the electrode surface and solvent molecules, the present work shows the dependence of these interactions on the ion nature. The calculations carried out using theoretical models of adsorption and complex formation for ions showed that the energy of adsorption on the electrode or the energy of binding into a complex with the solvent depends on the energy of the ion valence orbital. Most of the obtained dependencies are close to linear and exponential types, in some cases the solution leads to volcano-shaped dependencies.

A comparison of a number of ion parameters obtained by quantum chemical calculations revealed clear correlations with the HOMO and LUMO energies of the ion orbitals. The first parameter to be considered was the adsorption energies on model electrodes represented by aluminum oxide (0001), graphene and Au (111). The best correlations were observed for graphene and aluminum oxide; in case of gold, some deviations were found for the organic cations. Volcano-like correlations were revealed between the adsorption energy of anions on Al₂O₃ (0001) and graphene surfaces with ε_LUMO, and on Au (111) with ε_HOMO. Strong correlations were found for the energies of ion binding to solvent molecules represented by water and dimethyl carbonate. Volcano-like dependencies were also observed for the energy of binding of cations to H₂O and DMC on ε_LUMO. In addition, a clear relationship was observed between the dipole moment, polarization and solvation energy of ions, on the one hand, and the HOMO and LUMO energy level. For the anions, the oxidation potential was found to be dependent on the HOMO energy level.

The described dependencies of the ion parameters on the energy levels of molecular orbitals were approximated by expressions that allow us to use the energy level of orbitals as a descriptor in applied problems. The use of such parameters makes the search for ions with desired properties a lot easier. First of all, such descriptors can be relevant in the field of catalysis, especially for catalyst surface composition control and in surface chemistry, for example, when searching for ions with the highest adsorption for a given material to solve the problem of heavy metal binding by different sorbents. An example of possible practical application of the descriptors was clearly demonstrated by aqueous electrolytes for batteries. Li⁺ was used as the cation, the anion parameters were estimated by the previously obtained expressions, the task was to select an anion the use of which would minimize the degree of water coverage of the electrode surface. Reduction in the water coverage increases the working window potential of the electrolyte. The obtained estimates for the anions widely applied in electrochemical practice agree with the experimental data, new electrolytes can be potentially found among alkaline salts of sulfonic acid derivatives.