(574e) Binding Debye-Hückel Equation of State for Charged Hard Sphere Fluids

The thermodynamic modeling of electrolyte solutions is of great importance in many fields of science and engineering, including chemical engineering, material science, and biophysics. Electrolyte solutions are ubiquitous in many natural and industrial processes, such as electrochemical energy storage devices, chemical synthesis, and biological systems. Therefore, the development of accurate thermodynamic models for electrolyte solutions is essential for advancing our understanding of many natural and industrial processes and for designing new materials and processes.

Accurately predicting the properties of electrolyte solutions requires an understanding of the long-range ion-ion interactions. For the ion-ion interactions in electrolyte solutions, two theories are mostly used: Debye-Hückel¹ and Mean Spherical Approximation^2,3 theories. In these theories, it has been assumed that charged hard spheres are in continuum medium of solvent with a finite static permittivity. Furthermore, it has been assumed that electrolytes fully dissociate in the solution. However, it has been demonstrated that the second assumption may not hold true under all conditions^4–7. Therefore, it has been shown that accurate property predictions for certain systems or under certain conditions of the system necessitate the consideration of ion-ion association.

Moreover, the development of a more accurate equation of state (EOS) for charged hard sphere fluids has significant implications for understanding the behavior of many complex systems, including ionic liquids and electrolyte solutions. The new EOS provides a more accurate representation of ion-ion association, which is a critical factor in the thermodynamic properties of electrolyte solutions.

In this study, we formulated a novel equation of state (EOS) for a charged hard sphere fluid, accounting for ion-ion association. To develop this model, we have used the Wertheim two-density formulation⁸ (total (ρ_i⁰) and unbounded number densities (ρ_i = α_i ρ_i⁰)) for association which has proven to be effective in modeling hydrogen bonding (eq. ( 1 )).

βA^Assoc/N_tot = (1 / ρ_tot ) Σ_i [ ρ_i⁰ (Ln (α_i) + 0.5(1 - α_i) )] (1)

In this equation, β, N_tot, ρ_tot, and α_i are the are the Boltzmann factor (1/k_BT), total number of charged hard spheres, total number density, and fraction unbounded ions. A^Assoc is also the contribution of ion-ion association to the Helmholtz free energy. Moreover, we have employed the reference cavity approximation suggested by Stell et al.⁹ to calculate the fraction of unbounded ions. In this approach, the number density of ion pairs is calculated from eq. ( 2 ).

ρ_mn = K°_mn y_mn [ ρ_m⁰ - ∑_l ρ_ml][ρ_n⁰ - ∑_l ρ_ln] (2)

K°_mn = ∫ r² exp (2q/r) dr (3)

In this equation, K°_mn is the association constant at infinite dilution between cation m and anion n (eq. ( 3 )), N is number of anions, and M is number of cations in the solution¹⁰. The lower limit in the integral in eq. (3) is the contact distance between ion pairs (σ_mn) and the upper limit is the distance from the center of ions in which counter-ions considered as ion pairs (l_mn). y_mn is the reference cavity function between cation m and anion n calculated from eq. (17) in ref.¹¹.

The reference cavity function includes an electrostatic contribution (Ln(y^E) = β(μ^r_m + μ^r_n - μ^r_mn(ρ_mn=0))) and an excluded volume contribution (g^HS) which is the hard-sphere radial distribution function between ion pairs at contact distance^12,13. μ^r_m, μ^r_n, and μ^r_mn are the residual chemical potential of cation, anion, and ion pair due to electrostatic interactions in the no association limit, respectively. For the ions, the residual chemical potential has been calculated from the Debye-Hückel theory. For the ion pairs, the residual chemical potential has been calculated from the Debye-Hückel theory (if |Z_m|≠|Z_n|) and ion-dipole interactions from the Kirkwood equation. g^HS is also calculated from the hard-sphere models of Boublik¹² and Mansoori et al.¹¹ Finally, after calculation of number density of unbounded ions from eq. ( 2 ), the residual Helmholtz energy of the charged hard sphere fluid can be from eq. ( 4 ):

βA^r/N_tot = β(A^HS + A^DH + A^Assoc + A^Born) / N_tot (4)

In this equation, A^HS, A^DH, A^Assoc, and A^Born are the hard-sphere, Debye-Hückel, ion-ion association (eq. ( 1 )), and Born contributions to the residual Helmholtz energy. Then, the thermodynamic properties of charged hard sphere fluids can be calculated from the mole number, volume, and temperature derivatives of the residual Helmholtz energy.

To validate the new EOS, we have compared the predicted mean ionic activity coefficient and osmotic coefficient by the BiDH EOS for 1:1, 1:2, 2:2, and 3:1 electrolytes with the Monte Carlo (MC) simulations reported in the literature. These MC simulations are from Lamperski¹⁴ (L1-L5), Gutiérrez-Valladares et al.¹⁵ (G1-G3), and Abbas et al.¹⁶ (S1-S104). The results show that the new EOS accurately predicts the thermodynamic properties of the system. The model is also compared with the Mean Spherical Approximation, Binding Mean Spherical Approximation¹⁷, and Debye-Hückel models, demonstrating its superiority in accurately predicting the behavior of charged hard sphere fluids.

In addition, the performance of the EOS has been studied by comparing mean and individual ionic activity coefficient calculated by implicit-solvent molecular dynamic simulations reported by Saravi et al.¹⁸. In this comparison, the Leonard-Jones contribution to the residual Helmholtz free energy has been calculated from eq. (10) in ref.¹⁹.

For this comparison, we have used exactly the same parameters (ionic diameter and Leonard-Jones potential). We have shown that the mean ionic activity coefficient predicted by the BiDH EOS is in good agreement with MD simulations. Nonetheless, we also noted a discrepancy between the MD simulations and the BiDH model's prediction of individual activity coefficients.

In conclusion, the results of this study provide a significant contribution to the understanding of the thermodynamic properties of charged hard sphere fluids and demonstrate the effectiveness of the developed EOS in accurately predicting the behavior of these systems. This investigation holds substantial implications for advancing current electrolyte solution equations of state, which prove highly beneficial in the design and optimization of processes within the chemical industries.

Acknowledgment

We thank the ERC Advanced Grant Project No. 832460.

References

1 P. Debye and E. Hückel, Physikalische Zeitschrift, 1923, 24, 185–206.

2 L. Blum, Mol Phys, 1975, 30, 1529–1535.

3 L. Blum and J. S. Høye, Journal of Physical Chemistry, 1977, 81, 1311–1316.

4 Y. Marcus and G. Hefter, Chem Rev, 2006, 106, 4585–4621.

5 S. Naseri Boroujeni, B. Maribo-Mogensen, X. Liang and G. M. Kontogeorgis, Journal of Molecular Liquids (submitted).

6 S. Naseri Boroujeni, B. Maribo-Mogensen, X. Liang and G. M. Kontogeorgis, Fluid Phase Equilib, 2023, 567, 113698.

7 S. Naseri Boroujeni, X. Liang, B. Maribo-Mogensen and G. M. Kontogeorgis, Ind Eng Chem Res, 2022, 61, 3168–3185.

8 M. S. Wertheim, J Chem Phys, 1986, 85, 2929–2936.

9 G. Stell and Y. Zhou, Fluid Phase Equilib, 1992, 79, 1–20.

10 N. Bjerrum, Kgl. Danske Videnskab. Selskab, 1926, 7, 1–48.

11 S. Naseri Boroujeni, B. Maribo-Mogensen, X. Liang and G. M. Kontogeorgis, unpublished work.

12 G. A. Mansoori, N. F. Carnahan, K. E. Starling and T. W. Leland, J Chem Phys, 1971, 54, 1523–1525.

13 T. Boublík, J Chem Phys, 1970, 53, 471–472.

14 S. Lamperski, Mol Simul, 2007, 33, 1193–1198.

15 E. Gutiérrez-Valladares, M. Lukšič, B. Millán-Malo, B. Hribar-Lee and V. Vlachy, Condens Matter Phys, 2011, 14, 1–15.

16 Z. Abbas, E. Ahlberg and S. Nordholm, Journal of Physical Chemistry B, 2009, 113, 5905–5916.

17 O. Bernard and L. Blum, Journal of Chemical Physics, 1996, 104, 4746–4754.

18 S. H. Saravi and A. Z. Panagiotopoulos, Journal of Chemical Physics, 2021, 155, 184501.

19 J. Gross and G. Sadowski, Ind Eng Chem Res, 2001, 40, 1244–1260.