The Density Functional Theory for Electrolyte Solutions

For an electrolyte solution close to a charged surface with temperature , total volume , and chemical potential of each species specified, the grand potential, , is written as

$$\Omega[{\rho_i(\boldsymbol{r})},\psi(\boldsymbol{r})] = F^\text{id}[{\rho_i(\boldsymbol{r})}] + F^\text{exc}[{\rho_i(\boldsymbol{r})}]+ F^\text{coul}[{\rho_i(\boldsymbol{r})},\psi(\boldsymbol{r})]+ \sum_i \int_{V} d \boldsymbol{r} [V_i^{(\text{ext})}(\boldsymbol{r})-\mu_i] \rho_i(\boldsymbol{r})+ \int_{\partial V}d \boldsymbol{r} \sigma(\boldsymbol{r}) \psi(\boldsymbol{r})$$

where $\rho_i(\boldsymbol{r})$ is the local density of the component i, $\psi(\boldsymbol{r})$ is the local electrostatic potential, and $V^\text{ext}_{i}$ is the external potential. The volume of the system is V and $\partial V$ is the boundary of the system.

The ideal-gas contribution $F^\text{id}$ is given by the exact expression

$$F^\text{id}[{\rho_i (\boldsymbol{r})}] = k_B T\sum_i \int_{V} d\boldsymbol{r}\ \rho_i(\boldsymbol{r})[\ln(\rho_i (\boldsymbol{r})\Lambda_i^3)-1]$$

where is the Boltzmann constant, is the absolute temperature, and is the well-known thermal de Broglie wavelength of each ion.

The Coulomb's free-energy is obtained by the addition of the electric field energy density and the minimal-coupling of the interaction between the electrostatic potential and the charge density , and it can be written as

$$F^\text{coul}[{\rho_{i}(\boldsymbol{r})},\psi(\boldsymbol{r})] = -\int_V d\boldsymbol{r}\ \frac{\epsilon_0 \epsilon_r}{2} |\nabla{\psi(\boldsymbol{r})}|^2 + \int_{V} d\boldsymbol{r}\ \sum_i Z_i e \rho_{i}(\boldsymbol{r}) \psi(\boldsymbol{r})$$

where is the valence of the ion i, is the elementary charge, is the vacuum permittivity, and is the relative permittivity.

The excess Helmholtz free-energy, , is the free-energy functional due to particle-particle interactions splitted in the form

$$F^\text{exc}[{\rho_i(\boldsymbol{r})}] = F^\text{hs}[{\rho_i(\boldsymbol{r})}] + F^\text{ec}[{\rho_i(\boldsymbol{r})}]$$

where $F^{\textrm{hs}}$ is the hard-sphere excess contribution and $F^{\textrm{ec}}$ is the electrostatic correlation excess contribution.

The hard-sphere contribution, $F^{\textrm{hs}}$, represents the hard-sphere exclusion volume correlation and it can be described using different formulations of the fundamental measure theory (FMT) as

• Rosenfeld Functional (RF) - Rosenfeld, Y., Phys. Rev. Lett. 63, 980–983 (1989)
• White Bear version I (WBI) - Yu, Y.-X. & Wu, J., J. Chem. Phys. 117, 10156–10164 (2002); Roth, R., Evans, R., Lang, A. & Kahl, G., J. Phys. Condens. Matter 14, 12063–12078 (2002)
• White Bear version II (WBII) - Hansen-Goos, H. & Roth, R. J., Phys. Condens. Matter 18, 8413–8425 (2006)

The electrostatic correlation can be described using different approximations as

• Mean-Field Theory (MFT) -
• Bulk Fluid Density (BFD) - Kierlik and Rosinberg, Phys.Rev.A 44, 5025 (1991); Y. Rosenfeld, J. Chem. Phys. 98, 8126 (1993)
• functionalized Mean Spherical Approximation (fMSA) - Roth and Gillespie, J. Phys.: Condens. Matter 28, 244006 (2016)
• Reference Fluid Density (RFD) - Gillespie, D., Nonner, W. & Eisenberg, R. S., J. Phys. Condens. Matter 14, 12129–12145 (2002); Gillespie, D., Valiskó, M. & Boda, D., J. Phys. Condens. Matter 17, 6609–6626 (2005)

Finally, The chemical potential for each ionic species is defined as $\mu_i = \mu_i^\text{id} + \mu_i^\text{exc}$, where superscripts id and exc refer to ideal and excess contributions, respectively.

The thermodynamic equilibrium is obtained by the minimum of the grand-potential, $\Omega$, which can be obtained by the functional derivatives, such that, the equilibrium condition for each charged component is given by

$$\left. \frac{\delta \Omega}{\delta \rho_i(\boldsymbol{r})} \right|_{{\mu_k},V,T} = k_B T \ln[\Lambda_i^3 \rho_i(\boldsymbol{r})] + \frac{\delta F^\text{exc}}{\delta \rho_i(\boldsymbol{r})} + Z_i e \psi(\boldsymbol{r}) + V^{\text{ext}}_i(\boldsymbol{r}) - \mu_i =0$$

and for the electrostatic potential it is

$$\left. \frac{\delta \Omega}{\delta \psi(\boldsymbol{r})} \right|_{{\mu_k},V,T} = \epsilon_0 \epsilon_r\nabla^2{\psi(\boldsymbol{r})} + \sum_i Z_i e \rho_i(\boldsymbol{r}) =0$$

valid in the whole volume V, this is the well-known Poisson's equation of the electrostatic potential with the boundary conditions

$$\left. \frac{\delta \Omega}{\delta \psi(\boldsymbol{r})} \right|{{\mu_k},V,T} = \left. \epsilon_0 \epsilon_r \boldsymbol{\hat{n}}(\boldsymbol{r}) \cdot \boldsymbol{\nabla}{\psi(\boldsymbol{r})} \right|{\partial V} + \sigma(\boldsymbol{r}) = 0$$

valid on the boundary surface $\partial V$, where $\boldsymbol{\hat{n}}(\boldsymbol{r})$ is denoting the vector normal to the surface pointing inward to the system.

Fig.1 - The ionic density profiles of an 1:1 electrolyte solution with c_+= 0.01 M and σ = -0.5C/m².	Fig.2 - The electrostatic potential profile of an 1:1 electrolyte solution with c_+= 0.01 M and σ = -0.5C/m².

Fig.3 - The ionic density profiles of an 1:1 electrolyte solution with c_+= 1.0 M and σ = -0.5C/m².	Fig.4 - The electrostatic potential profile of an 1:1 electrolyte solution with c_+= 1.0 M and σ = -0.5C/m².