Dipolar Hard Spheres: Comprehensive Data from Monte Carlo

Feb 1, 2019 - The simulation results further confirm noteworthy chain and ring structures of the dipolar hard-sphere particles at lower densities and ...
1 downloads 0 Views 1MB Size
Download PDF
Article Cite This: J. Chem. Eng. Data XXXX, XXX, XXX−XXX

pubs.acs.org/jced

Dipolar Hard Spheres: Comprehensive Data from Monte Carlo Simulations Marc Theiss* and Joachim Gross*

J. Chem. Eng. Data Downloaded from pubs.acs.org by WEBSTER UNIV on 02/05/19. For personal use only.

Institute of Thermodynamics and Thermal Process Engineering, University of Stuttgart, Pfaffenwaldring 9, 70569 Stuttgart, Germany ABSTRACT: In this work, we report comprehensive molecular simulation data for the Helmholtz energy, the internal energy, the constant-volume heat capacity, and the relative permittivity of pure nonpolarizable dipolar hard spheres. The dimensionless density (ρ*) thereby covers a range from 0.05 to 1.0, and dimensionless squared dipole moments (μ*2) range from 0.04 to 7.0. An empirical polynomial correlation of the Helmholtz energy is parametrized, which can be used to substitute molecular simulations within this range of density and dipole moment. The correlation and its derivatives may facilitate the development of new theories and equation-of-state models of dipolar fluids. The simulation results further confirm noteworthy chain and ring structures of the dipolar hard-sphere particles at lower densities and at high dipole moments, which are known to have an impact on the thermodynamic properties of a system.

■

INTRODUCTION Results from molecular simulation of model fluids are important for developing and evaluating fluid theories, because the underlying molecular models of theory and simulation can be identical, allowing the unambiguous identification of shortcomings of theoretical approaches. Dipolar and ionic fluids are particularly relevant for such comparison, because new theoretical development is needed for such fluids. The relative permittivity (ε) is of central interest in the thermodynamics of purely dipolar fluids and of electrolyte solutions. This property can be associated with a solvent’s capability of shielding ions, present in electrolyte solutions, by aligning the dipole molecules around the ions to effectively suppress the long-range electrostatic interactions.1 The relative permittivity is important especially in the McMillan−Mayer theory of solutions. Here, one applies a simplifying approach based on the assumption the solvent molecules are smeared out to form a uniform background with a relative permittivity ε (also termed the primitive model). Famous theoretical models are based on the McMillan−Mayer approach, e.g., the Debye− Hü ckel theory2 and the mean spherical approximation theories.3,4 A more detailed representation of electrolyte solutions is possible with nonprimitive models, where the solvent molecules are explicitly considered. In this model class, both ions and solvent molecules are treated explicitly with a molecular model. There are some simulation studies of nonprimitive model electrolyte solutions, which are most commonly described by a mixture of charged and dipolar hard spheres, focused on the Helmholtz energy, the compressibility, and the internal energy.5−11 Drunsel and Gross12 calculated the chemical potential of electrolyte solutions for several liquid-like densities © XXXX American Chemical Society

with varying ionic charges, dipole moments, and ion compositions. Some scientists computed vapor−liquid equilibria13−15 and static dielectric properties16−20 of the Stockmayer fluid, whereas some previous simulation studies focused on dipolar hard spheres.21−28 Depending on the density, temperature, or dipole moment, a dipolar fluid can form different kinds of structures such as rings, chains, and branches. For example, Levesque and Weis29 observed that dipolar hard-sphere particles can assemble into chainlike structures, especially at low densities (0.05 ≤ ρ* ≤ 0.2) and high dipole moments. Additionally, at lower temperatures and densities, Rovigatti et al.30 found a stabilization of rings compared to chains. These structures have a significant impact on the thermodynamic properties of a system,31−33 for example, on the phase behavior of dipolar soft spheres. Wei and Patey34 considered higher densities (ρ* ≥ 0.6) of dipolar soft spheres and observed that an orientationally ordered fluid phase can be formed for sufficiently low temperatures. For dipolar hard spheres, no liquid−vapor transitions have been observed so far.30,33 In this work, we present Monte Carlo (MC) simulation data for the Helmholtz energy, the internal energy, the constantvolume heat capacity, and the relative permittivity of pure nonpolarizable dipolar hard-sphere fluids. Only the configurational part of the energies is considered here. On the basis of these data, we correlate an analytic three-dimensional Helmholtz energy surface that has some particularly noteworthy properties. First, the Helmholtz energy surface could be used to Received: December 7, 2018 Accepted: January 17, 2019

A

DOI: 10.1021/acs.jced.8b01169 J. Chem. Eng. Data XXXX, XXX, XXX−XXX

Journal of Chemical & Engineering Data

Article

Table 1. Numerical Values of the Internal Energy, the Residual Helmholtz Energy, the Constant-Volume Heat Capacity, and the Relative Permittivity for 0.05 ≤ ρ* ≤ 0.5 and 0.04 ≤ μ*2 ≤ 7.0a ρ*

μ*2

0.05 0.05 0.05 0.05 0.05 0.05 0.05 0.05 0.15 0.15 0.15 0.15 0.15 0.15 0.15 0.15 0.25 0.25 0.25 0.25 0.25 0.25 0.25 0.25 0.35 0.35 0.35 0.35 0.35 0.35 0.35 0.35 0.45 0.45 0.45 0.45 0.45 0.45 0.45 0.45

0.04 0.5 1.0 2.0 3.0 4.0 6.0 7.0 0.04 0.5 1.0 2.0 3.0 4.0 6.0 7.0 0.04 0.5 1.0 2.0 3.0 4.0 6.0 7.0 0.04 0.5 1.0 2.0 3.0 4.0 6.0 7.0 0.04 0.5 1.0 2.0 3.0 4.0 6.0 7.0

⟨βU ⟩ N

−0.0001(0) −0.0182(1) −0.0745(2) −0.3558(6) −1.121(2) −3.074(10) −10.16(2) −13.12(2) −0.0004(0) −0.0556(2) −0.2171(5) −0.8873(10) −2.2(3) −4.379(10) −10.33(2) −13.21(1) −0.0007(0) −0.0949(3) −0.3519(5) −1.28(1) −2.785(3) −4.912(6) −10.37(1) −13.2(1) −0.001(0) −0.1354(4) −0.479(6) −1.603(1) −3.22(3) −5.305(4) −10.45(1) −13.21(1) −0.0014(0) −0.1772(4) −0.6032(5) −1.887(1) −3.583(3) −5.639(5) −10.58(1) −13.31(2)

β(A − A 0) N

−0.0001(0) −0.009(0) −0.0365(1) −0.1586(1) −0.4226(3) −0.9745(12) −3.517(6) −5.304(6) −0.0002(0) −0.0282(1) −0.1104(2) −0.4376(3) −1.02(1) −1.929(1) −4.799(4) −6.614(5) −0.0003(0) −0.0489(1) −0.1848(2) −0.6829(4) −1.467(1) −2.545(1) −5.539(3) −7.36(4) −0.0005(0) −0.0708(2) −0.2604(2) −0.9085(5) −1.848(1) −3.048(1) −6.143(2) −7.959(4) −0.0007(0) −0.0943(1) −0.337(3) −1.122(0) −2.194(1) −3.499(1) −6.695(3) −8.53(3)

CV kN

ε

ρ*

μ*2

0.0(0) 0.019(1) 0.081(3) 0.54(2) 2.5(1) 8.1(14) 6.9(19) 3.1(12) 0.0(0) 0.052(3) 0.22(1) 0.98(4) 3.0(2) 5.9(11) 8.1(23) 3.9(10) 0.001(0) 0.08(4) 0.31(2) 1.1(1) 2.7(2) 4.7(6) 5.4(13) 4.9(15) 0.001(0) 0.12(1) 0.36(2) 1.1(1) 2.6(2) 3.7(3) 4.7(10) 4.5(12) 0.001(0) 0.16(1) 0.43(2) 1.1(1) 2.0(1) 3.2(4) 4.8(9) 4.4(13)

1.008(0) 1.11(3) 1.228(6) 1.492(13) 1.866(25) 2.631(40) − − 1.025(1) 1.338(9) 1.746(19) 2.887(49) 4.626(97) − − − 1.042(1) 1.621(14) 2.416(40) 4.705(93) 8.43(16) − − − 1.06(1) 1.916(26) 3.194(57) 6.86(14) 12.87(31) 23.24(46) − − 1.077(2) 2.257(31) 4.063(86) 9.66(21) 19.07(46) − − −

0.1 0.1 0.1 0.1 0.1 0.1 0.1 0.1 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.4 0.4 0.4 0.4 0.4 0.4 0.4 0.4 0.5 0.5 0.5 0.5 0.5 0.5 0.5 0.5

0.04 0.5 1.0 2.0 3.0 4.0 6.0 7.0 0.04 0.5 1.0 2.0 3.0 4.0 6.0 7.0 0.04 0.5 1.0 2.0 3.0 4.0 6.0 7.0 0.04 0.5 1.0 2.0 3.0 4.0 6.0 7.0 0.04 0.5 1.0 2.0 3.0 4.0 6.0 7.0

β(A − A 0) N

⟨βU ⟩ N

−0.0003(0) −0.0365(2) −0.1466(5) −0.6449(11) −1.768(3) −3.92(8) −10.28(5) −13.28(2) −0.0005(0) −0.0752(3) −0.2847(5) −1.095(1) −2.527(2) −4.683(6) −10.37(2) −13.21(2) −0.0008(0) −0.1147(3) −0.4168(7) −1.449(1) −3.017(3) −5.11(6) −10.41(1) −13.18(3) −0.0012(0) −0.156(4) −0.5427(7) −1.749(1) −3.409(2) −5.484(6) −10.52(1) −13.24(1) −0.0016(0) −0.1988(4) −0.6628(6) −2.017(1) −3.753(3) −5.814(4) −10.68(1) −13.35(1)

−0.0001(0) −0.0184(1) −0.0734(2) −0.3032(3) −0.7488(6) −1.526(1) −4.297(5) −6.092(5) −0.0003(0) −0.0384(1) −0.1476(2) −0.5633(5) −1.254(1) −2.26(1) −5.203(3) −7.013(5) −0.0004(0) −0.0598(1) −0.2226(2) −0.7975(4) −1.662(1) −2.806(1) −5.853(3) −7.665(4) −0.0006(0) −0.0824(2) −0.2987(3) −1.016(0) −2.024(1) −3.277(1) −6.424(2) −8.247(3) −0.0008(0) −0.1066(2) −0.3755(3) −1.224(1) −2.359(1) −3.713(1) −6.959(2) −8.811(3)

CV kN

ε

0.0(0) 0.036(2) 0.16(1) 0.82(4) 2.7(2) 7.0(10) 9.9(64) 4.3(13) 0.001(0) 0.069(4) 0.26(1) 1.1(1) 2.7(2) 5.4(8) 7.3(18) 6.0(21) 0.001(0) 0.1(0) 0.34(2) 1.2(1) 2.5(2) 3.9(5) 5.4(10) 6.0(25) 0.001(0) 0.14(1) 0.43(2) 1.2(1) 2.0(1) 3.8(4) 5.0(13) 5.0(13) 0.001(0) 0.17(1) 0.46(2) 1.2(1) 2.0(1) 2.8(3) 5.1(12) 4.9(10)

1.017(0) 1.221(6) 1.486(12) 2.101(27) 3.094(54) 5.07(11) − − 1.034(1) 1.481(13) 2.074(30) 3.723(71) 6.59(13) 11.16(25) − − 1.051(1) 1.761(19) 2.807(43) 5.92(13) 10.93(28) − − − 1.068(2) 2.078(31) 3.673(74) 8.43(21) 16.54(41) 28.44(73) − − 1.086(2) 2.431(39) 4.543(98) 11.56(31) 23.65(62) − − −

a

The uncertainty in the last digit is shown in parentheses.

validate equations of state or other theoretical models with respect to the Helmholtz energy or its derivatives. Second, it can be used to estimate results from molecular simulations under various conditions (0.05 ≤ ρ* = ρσ3 ≤ 1.0 and 0.04 ≤ μ*2 = βμ2/ σ3 ≤ 7.0, respectively, where ρ is the particle density, β = 1/kT, k is Boltzmann’s constant, T is the temperature, σ is the hardsphere diameter, and μ denotes the dipole moment). Through the analytic Helmholtz energy correlation, costly computer simulations involving free energy methods, such as thermodynamic integration, can be avoided.

■

A − A0 = =

∫0 ∫0

1

1

∂A ∂λ

∂Uλ ∂λ

l

≈

∑ Δλi i=1

dλ λ

dλ NVT , λ

∂Uλ ∂λ

NVT , λi

(1)

with Uλ = U0 + λUdd

MOLECULAR SIMULATION

(2)

where λ is a coupling parameter defined between 0 and 1 and U0 represents the total configurational energy of the reference system, i.e., the hard-sphere fluid, while Udd = ∑i,j>i udd(ij) represents the total dipole potential, which is computed using

Thermodynamic integration (TI) is one of the simplest methods for calculating Helmholtz energy differences from molecular simulation.35 Here, we use the so-called λ-TI method B

DOI: 10.1021/acs.jced.8b01169 J. Chem. Eng. Data XXXX, XXX, XXX−XXX

Journal of Chemical & Engineering Data

Article

Table 2. Numerical Values of the Internal Energy, the Residual Helmholtz Energy, the Constant-Volume Heat Capacity, and the Relative Permittivity for 0.55 ≤ ρ* ≤ 1.0 and 0.04 ≤ μ*2 ≤ 7.0a ρ*

μ*2

0.55 0.55 0.55 0.55 0.55 0.55 0.55 0.55 0.65 0.65 0.65 0.65 0.65 0.65 0.65 0.65 0.75 0.75 0.75 0.75 0.75 0.75 0.75 0.75 0.85 0.85 0.85 0.85 0.85 0.85 0.85 0.85 0.95 0.95 0.95 0.95 0.95 0.95 0.95 0.95

0.04 0.5 1.0 2.0 3.0 4.0 6.0 7.0 0.04 0.5 1.0 2.0 3.0 4.0 6.0 7.0 0.04 0.5 1.0 2.0 3.0 4.0 6.0 7.0 0.04 0.5 1.0 2.0 3.0 4.0 6.0 7.0 0.04 0.5 1.0 2.0 3.0 4.0 6.0 7.0

⟨βU ⟩ N

−0.0018(0) −0.2206(4) −0.7221(8) −2.143(1) −3.92(2) −5.983(3) −10.78(1) −13.55(2) −0.0022(0) −0.2642(5) −0.8372(7) −2.387(1) −4.248(2) −6.331(4) −11.11(2) −13.83(1) −0.0027(1) −0.3093(4) −0.9507(9) −2.622(1) −4.566(2) −6.71(3) −11.6(2) −14.3(2) −0.0032(1) −0.3551(5) −1.063(1) −2.852(1) −4.9(3) −7.17(10) −12.2(1) −14.79(1) −0.0038(1) −0.4019(5) −1.173(1) −3.083(1) −5.296(9) −7.775(7) −13.24(1) −15.81(11)

β(A − A 0) N

−0.0009(0) −0.1189(2) −0.4143(3) −1.326(1) −2.52(1) −3.924(1) −7.232(2) −9.096(3) −0.0011(0) −0.1451(2) −0.4925(3) −1.526(1) −2.837(1) −4.337(1) −7.78(2) −9.7(3) −0.0014(0) −0.1723(2) −0.572(3) −1.722(1) −3.146(1) −4.749(1) −8.373(3) −10.36(0) −0.0017(0) −0.2009(2) −0.6527(4) −1.918(1) −3.455(1) −5.17(2) −9.021(2) −11.09(0) −0.002(0) −0.2307(2) −0.7342(4) −2.115(1) −3.771(1) −5.635(2) −9.786(3) −12.0(0)

CV kN

ε

ρ*

μ*2

0.002(0) 0.18(1) 0.45(2) 1.1(1) 1.8(1) 2.8(3) 4.6(13) 5.0(11) 0.002(0) 0.2(1) 0.48(2) 1.1(1) 1.7(1) 2.2(2) 4.9(17) 4.4(10) 0.003(0) 0.22(1) 0.51(2) 0.97(6) 1.6(1) 2.3(2) 3.5(18) 5.2(16) 0.003(0) 0.24(1) 0.54(3) 1.0(7) 1.5(2) 2.7(8) 3.0(7) 3.3(6) 0.004(0) 0.26(1) 0.55(3) 1.1(1) 2.6(7) 2.2(4) 2.3(5) 23.0(171)

1.093(2) 2.585(44) 5.2(11) 13.5(32) 26.55(69) − − − 1.11(3) 2.928(54) 6.36(15) 17.98(46) 43.5(12) − − − 1.131(3) 3.443(65) 7.86(17) 25.4(63) 61.5(15) − − − 1.15(4) 3.969(64) 9.48(22) 37.94(84) − − − − 1.167(4) 4.409(88) 11.76(26) 64.2(15) − − − −

0.6 0.6 0.6 0.6 0.6 0.6 0.6 0.6 0.7 0.7 0.7 0.7 0.7 0.7 0.7 0.7 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0

0.04 0.5 1.0 2.0 3.0 4.0 6.0 7.0 0.04 0.5 1.0 2.0 3.0 4.0 6.0 7.0 0.04 0.5 1.0 2.0 3.0 4.0 6.0 7.0 0.04 0.5 1.0 2.0 3.0 4.0 6.0 7.0 0.04 0.5 1.0 2.0 3.0 4.0 6.0 7.0

⟨βU ⟩ N

−0.002(0) −0.2422(5) −0.7801(7) −2.265(1) −4.081(2) −6.153(4) −10.94(1) −13.64(2) −0.0025(1) −0.2864(4) −0.8949(8) −2.504(1) −4.405(2) −6.519(3) −11.37(2) −14.07(2) −0.003(1) −0.3318(6) −1.006(1) −2.736(1) −4.731(2) −6.911(4) −11.89(1) −14.54(1) −0.0035(1) −0.3788(6) −1.118(1) −2.968(1) −5.082(3) −7.478(11) −12.5(1) −15.12(1) −0.0042(0) −0.4257(3) −1.228(0) −3.205(1) −5.569(8) −8.079(5) −13.69(1) −16.93(0)

β(A − A 0) N

−0.001(0) −0.1318(2) −0.4534(3) −1.427(1) −2.679(1) −4.132(1) −7.503(2) −9.387(4) −0.0013(0) −0.1587(2) −0.5322(4) −1.624(1) −2.992(1) −4.543(1) −8.068(2) −10.02(0) −0.0015(0) −0.1866(2) −0.6119(4) −1.82(1) −3.3(1) −4.957(1) −8.687(2) −10.71(0) −0.0018(0) −0.2157(2) −0.6932(4) −2.016(1) −3.611(1) −5.394(2) −9.368(4) −11.49(0) −0.0021(0) −0.2461(2) −0.7757(2) −2.213(0) −3.941(1) −5.888(2) −10.27(1) −12.64(1)

CV kN

ε

0.002(0) 0.18(1) 0.47(2) 1.1(1) 1.9(1) 2.6(3) 3.7(5) 4.6(18) 0.002(0) 0.21(1) 0.52(2) 1.1(1) 1.5(1) 2.1(2) 3.8(17) 4.0(14) 0.003(0) 0.23(1) 0.53(2) 0.98(6) 1.4(1) 2.4(3) 4.1(16) 3.1(8) 0.003(0) 0.25(1) 0.52(2) 0.97(6) 1.5(2) 2.5(7) 2.9(6) 2.8(6) 0.004(0) 0.25(1) 0.54(2) 1.1(1) 2.4(5) 1.9(3) 1.9(8) 1.6(3)

1.102(3) 2.801(51) 5.76(11) 15.53(45) 32.35(74) − − − 1.124(3) 3.212(50) 7.07(16) 20.24(51) 51.7(14) − − − 1.139(4) 3.573(62) 8.56(20) 31.39(73) − − − − 1.158(4) 4.106(88) 10.72(27) 46.5(12) 278.5(78) − − − 1.178(5) 4.71(88) 13.18(33) 89.7(22) − − − −

a

The uncertainty in the last digit is shown in parentheses.

standard Ewald summation.36 By introducing artificial intermediate states λi, we can compute the desired energy differences between the fluid of interest A = Aλ=1 and the reference fluid A0 = Aλ=0. We define 20 equally spaced λ values (l = 20) and apply Simpson’s rule for carrying out the integration. Note that we use the so-called ion extended dipole model in our computer simulations, which means that the dipole moments are represented as point charges, separated by a small distance d of 0.1σ. Drunsel and Gross37 showed that this molecular model mimics hard-sphere fluids with point dipoles sufficiently well. We used an Ewald parameter of 5.6 together with a sum over seven reciprocal vectors. Computations with nine reciprocal vectors led to indistinguishable results, which we validated for two high-dipole density states (ρ* = 0.8; μ*2 = 2.0; μ*2 = 4.0). We consider N = 1000 particles in Monte Carlo simulations and collect samples, e.g., for the internal energy from 105 × N

trial moves in the production period. Note that we use at least the same number for the equilibration phase. In our Monte Carlo simulations, we solely conduct displacement and rotation moves of individual particles. No special algorithms (e.g., the aggregation-volume-bias Monte Carlo algorithm38,39 or the unbonding and bonding algorithm40) were implemented to increase the sampling efficiency for those thermodynamic states in which chain and ring structures occur. In a postprocessing step, we apply the method of Shirts and Chodera41 to remove correlated data together with a resampled data size of 200 to compute the uncertainty of the mean with a 95% confidence interval using a bootstrap method.42 In the canonical ensemble, the (configurational) isochoric heat capacity can be obtained by calculating the fluctuations of the configurational electrostatic energy as CV /kN = [⟨(βU )2 ⟩ − ⟨βU ⟩2 ] /N C

(3)

DOI: 10.1021/acs.jced.8b01169 J. Chem. Eng. Data XXXX, XXX, XXX−XXX

Journal of Chemical & Engineering Data

Article

where β = 1/kT, k is Boltzmann’s constant, T is the temperature, and N is the particle number. A reliable way of calculating relative permittivity ε is to run a computer simulation using the standard Ewald summation without the surface term to consider a good conductor as a surrounding medium and then to calculate the fluctuation of total dipole moment M.43 We calculate ε as 4 ε = 1 + πρβ(⟨M2⟩ − ⟨M⟩2 )/N (4) 3 where M = ∑Ni=1 μi and μi = μiμ̂i, with μ̂ i being the orientational unit vector of the dipole of particle i.

■

RESULTS We provide tabulated values of thermodynamic properties of the dipolar hard-sphere fluid obtained from MC simulations for densities (ρ*) in the range from 0.05 to 1.0 and dipole moments (μ*2) from 0.04 to 7.0. Tables 1 and 2 summarize values for the residual Helmholtz energy, the internal energy, the isochoric heat capacity, and the relative permittivity. On the basis of the data listed in Tables 1 and 2, we provide a polynomial correlating the Helmholtz energy of a pure fluid for 0.05 ≤ ρ* ≤ 1.0 and 0.04 ≤ μ*2 ≤ 7.0. We use the polynomial of the form β

Figure 1. Residual Helmholtz energy of a dipolar hard-sphere fluid: reference data obtained from MC simulations using λ-thermodynamic integration (black symbols) and results obtained from the correlation (eq 5) (colored surface).

A − A0 = c1μ*2 + c 2(μ*2 )2 + c3(μ*2 )3 + c4μ*2 ρ* + c5(μ*2 )2 N ρ* + c6μ*2 (ρ*)2 + c 7(μ*2 )2 (ρ*)2 + c8(μ*2 )3 ρ* + c 9μ*2 (ρ*)3 + c10(μ*2 )3 (ρ*)2 + c11(μ*2 )2 (ρ*)3 + c12(μ*2 )3 (ρ*)3 + c13(μ*2 )3 (ρ*)4 + c14(μ*2 )4 (ρ*)3 + c15(μ*2 )4 (ρ*)4

(5)

The correlated coefficients ci of eq 5 are listed in Table 3. The coefficient of determination R2 = 0.9999. Figure 1 compares data

Figure 2. Snapshot of a NVT ensemble MC simulation with N = 1000 particles at ρ* = 0.05 and μ*2 = 6.0.

Table 3. Numerical Values of the Coefficients in Equation 5a c1 c2 c3 c4 c5 c6 c7 c8

0.06601 0.00126 −0.01499 −0.13228 −0.87347 −0.81392 0.88407 0.07715

c9 c10 c11 c12 c13 c14 c15

functions, change considerably when slowly relaxing structures appear. The quality of the relative permittivity can be determined by calculating ε in two different ways. In one approach, eq 4 is used directly; in a second approach, one can use the fact that the total averaged dipole moment is zero for isotropic fluids. We calculate ε in both ways, first with ⟨M⟩2 as obtained from molecular simulation and second with ⟨M⟩2 = 0. Larger deviations between both calculation methods indicate that slowly relaxing nonisotropic fluid structures have evolved in a simulation, indicating that the computed value for ε is less reliable. Figure 3 gives the relative deviations between the two calculation methods. The green and yellow labeled symbols indicate relative deviations of ≤1 and ≤10%, respectively. The red labeled symbols highlight relative deviations of >10%. The numerical values of the relative permittivity corresponding to the green labeled data points are highlighted in Figure 4 and are represented in Tables 1 and 2.

0.56577 −0.01446 −0.52534 0.03364 −0.01048 −0.01072 0.0068

a

The coefficient of determination R2 = 0.9999.

for the residual Helmholtz energy from molecular simulations (black symbols; see Tables 1 and 2) to the proposed correlation. It is important to note that our simulations also show chain and ring formations, most pronounced at low densities and at high dipole moments. Figure 2 shows a snapshot configuration of dipolar hard spheres at ρ* = 0.05 and μ*2 = 6.0. Chain structures are clearly observed in this configuration, and an eight-membered ring is seen at the bottom left of of Figure 2. These ring and chain structures are important for conducting molecular simulations, because these structures are slowly relaxing so that true ensemble averages require (prohibitively) long simulations. The relative permittivity is particularly demanding for systems with slowly relaxing chain structures. Also, for the development of fluid theories, such chain structures are important, because both macroscopic properties, such as heat capacities, are affected by the appearance of chain structures and microscopic quantities, such as n-body distribution

■

CONCLUSION This study reports comprehensive data for dipolar hard-sphere fluid obtained from molecular simulations. Data for the Helmholtz and internal energy, for the constant-volume heat capacity, and for the relative permittivity are reported for reduced densities (ρ*) ranging from 0.05 to 1.0 and for reduced dipole moments (μ*2) ranging from 0.04 to 7.0. The Helmholtz energy was correlated using a polynomial in density and dipole moment. The correlation can be used to substitute raw data from molecular simulations in the given range of density and D

DOI: 10.1021/acs.jced.8b01169 J. Chem. Eng. Data XXXX, XXX, XXX−XXX

Journal of Chemical & Engineering Data

■

REFERENCES

(1) Robinson, R. A.; Stokes, R. H. Interpretation of Conductance, Chemical Potential and Diffusion in Solutions of Simple Electrolytes, 2nd ed.; Butterworths Scientific Publication: London, 1959. (2) Hückel, E.; Debye, P. The theory of electrolytes: I. lowering of freezing point and related phenomena. Phys. Z. 1923, 24, 185−206. (3) Blum, L. Mean spherical model for asymmetric electrolytes: I. Method of solution. Mol. Phys. 1975, 30, 1529−1535. (4) Blum, L.; Høye, J. Mean spherical model for asymmetric electrolytes. 2. Thermodynamic properties and the pair correlation function. J. Phys. Chem. 1977, 81, 1311−1316. (5) Drunsel, F.; Gross, J. Theory of model electrolyte solutions: Assessing the short- and long-ranged contributions by molecular simulations. Fluid Phase Equilib. 2016, 430, 195−206. (6) Chan, K.-Y.; Gubbins, K. E.; Henderson, D.; Blum, L. Monte Carlo and simple theoretical calculations for ion-dipole mixtures. Mol. Phys. 1989, 66, 299−316. (7) Eggebrecht, J.; Ozler, P. Multipolar electrolyte solution models. I. Computer simulation of the charged and dipolar hard sphere mixture. J. Chem. Phys. 1990, 93, 2004−2015. (8) Eggebrecht, J.; Ozler, P. Multipolar electrolyte solution models. IV. Thermodynamic perturbation theory. J. Chem. Phys. 1993, 98, 1552−1565. (9) Eggebrecht, J.; Peters, G. H. Multipolar electrolyte solution models. II. Monte Carlo convergence and size dependence. J. Chem. Phys. 1993, 98, 1539−1545. (10) Lo, W. Y.; Chan, K.-Y.; Henderson, D. Improved Monte Carlo simulations of the structure of ion-dipole mixtures. Mol. Phys. 1993, 80, 1021−1029. (11) Caillol, J.; Levesque, D.; Weis, J. Monte Carlo simulation of an ion-dipole mixture. Mol. Phys. 1990, 69, 199−208. (12) Drunsel, F.; Gross, J. Chemical potential of model electrolyte solutions consisting of hard sphere ions and hard dipoles from molecular simulations. Fluid Phase Equilib. 2016, 429, 205−213. (13) Van Leeuwen, M. E.; Smit, B.; Hendriks, E. Vapour-liquid equilibria of stockmayer fluids. Mol. Phys. 1993, 78, 271−283. (14) Smit, B.; Williams, C.; Hendriks, E.; Leeuw, S. D. Vapour-liquid equilibria for Stockmayer fluids. Mol. Phys. 1989, 68, 765−769. (15) Bartke, J.; Hentschke, R. Phase behavior of the Stockmayer fluid via molecular dynamics simulation. Phys. Rev. E 2007, 75, 061503. (16) Adams, D.; Adams, E. Static dielectric properties of the Stockmayer fluid from computer simulation. Mol. Phys. 1981, 42, 907−926. (17) Neumann, M.; Steinhauser, O.; Pawley, G. S. Consistent calculation of the static and frequency-dependent dielectric constant in computer simulations. Mol. Phys. 1984, 52, 97−113. (18) Pollock, E.; Alder, B. Static dielectric properties of stockmayer fluids. Phys. A 1980, 102, 1−21. (19) Gray, C.; Sainger, Y.; Joslin, C.; Cummings, P.; Goldman, S. Computer simulation of dipolar fluids. Dependence of the dielectric constant on system size: A comparative study of Ewald sum and reaction field approaches. J. Chem. Phys. 1986, 85, 1502−1504. (20) Levesque, D.; Weis, J. On the calculation of dielectric properties from computer simulations. Phys. A 1984, 125, 270−274. (21) Patey, G. N.; Valleau, J. P. Dipolar hard spheres: A Monte Carlo study. J. Chem. Phys. 1974, 61, 534−540. (22) Levesque, D.; Patey, G.; Weis, J. A Monte Carlo study of dipolar hard spheres The pair distribution function and the dielectric constant. Mol. Phys. 1977, 34, 1077−1091. (23) Patey, G.; Levesque, D.; Weis, J. On the theory and computer simulation of dipolar fluids. Mol. Phys. 1982, 45, 733−746. (24) Leeuw, S. W. D.; Perram, J. W.; Smith, E. R. Simulation of Electrostatic Systems in Periodic Boundary Conditions. III. Further Theory and Applications. Proc. R. Soc. London, Ser. A 1983, 388, 177− 193. (25) de Leeuw, S. W.; Perram, J. W.; Smith, E. R. Simulation of electrostatic systems in periodic boundary conditions. I. Lattice sums and dielectric constants. Proc. R. Soc. London, Ser. A 1980, 373, 27−56.

Figure 3. Relative deviations of ε from eq 4 compared to ε from eq 4 with ⟨M⟩2 = 0. Relative deviations of ≤1 and ≤10% are colored green and yellow, respectively. The color red highlights relative deviations of >10%.

Figure 4. Molecular simulation data of the relative permittivity (see eq 4) for 0.05 ≤ ρ* ≤ 1.0 and 0.04 ≤ μ*2 ≤ 4.0. The ε values correspond to the green labeled data points in Figure 3.

dipole moment. The simulation data and the Helmholtz−energy correlation are intended to facilitate the development of fluid theories for dipolar (and ionic) fluids. In agreement with previous studies on the subject, our simulations indicate the occurrence of chain and ring structures at small densities and high dipole moments.

■

Article

AUTHOR INFORMATION

Corresponding Authors

*E-mail: [email protected]. *E-mail: [email protected]. ORCID

Joachim Gross: 0000-0001-8632-357X Funding

The authors acknowledge support by the state of BadenWürttemberg through bwHPC and the German Research Foundation (DFG) through Grant INST 40/467-1 FUGG. Notes

The authors declare no competing financial interest. E

DOI: 10.1021/acs.jced.8b01169 J. Chem. Eng. Data XXXX, XXX, XXX−XXX

Journal of Chemical & Engineering Data

Article

(26) De Leeuw, S.; Perram, J.; Smith, E. Simulation of electrostatic systems in periodic boundary conditions. II. Equivalence of boundary conditions. Proc. R. Soc. London, Ser. A 1980, 373, 57−66. (27) De Leeuw, S.; Perram, J. W.; Smith, E. Computer simulation of the static dielectric constant of systems with permanent electric dipoles. Annu. Rev. Phys. Chem. 1986, 37, 245−270. (28) Adams, D. Computer simulation of highly polar liquids: The hard sphere plus point dipole potential. Mol. Phys. 1980, 40, 1261−1271. (29) Levesque, D.; Weis, J. J. Orientational and structural order in strongly interacting dipolar hard spheres. Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top. 1994, 49, 5131−5140. (30) Rovigatti, L.; Russo, J.; Sciortino, F. Structural properties of the dipolar hard-sphere fluid at low temperatures and densities. Soft Matter 2012, 8, 6310−6319. (31) Sindt, J. O.; Camp, P. J. Simulations of dipolar fluids using effective many-body isotropic interactions. J. Chem. Phys. 2015, 143, 024501. (32) Rovigatti, L.; Tavares, J. M.; Sciortino, F. Self-Assembly in Chains, Rings, and Branches: A Single Component System with Two Critical Points. Phys. Rev. Lett. 2013, 111, 168302. (33) Ganzenmüller, G.; Patey, G.; Camp, P. J. Vapour-liquid phase transition of dipolar particles. Mol. Phys. 2009, 107, 403−413. (34) Wei, D.; Patey, G. N. Ferroelectric liquid-crystal and solid phases formed by strongly interacting dipolar soft spheres. Phys. Rev. A: At., Mol., Opt. Phys. 1992, 46, 7783−7792. (35) Bruckner, S.; Boresch, S. Efficiency of alchemical free energy simulations. I. A practical comparison of the exponential formula, thermodynamic integration, and Bennett’s acceptance ratio method. J. Comput. Chem. 2011, 32, 1303−1319. (36) Ewald, P. P. Die Berechnung optischer und elektrostatischer Gitterpotentiale. Ann. Phys. 1921, 369, 253−287. (37) Drunsel, F.; Zmpitas, W.; Gross, J. A new perturbation theory for electrolyte solutions. J. Chem. Phys. 2014, 141, 054103. (38) Chen, B.; Siepmann, J. I. A novel Monte Carlo algorithm for simulating strongly associating fluids: Applications to water, hydrogen fluoride, and acetic acid. J. Phys. Chem. B 2000, 104, 8725−8734. (39) Chen, B.; Siepmann, J. I. Improving the Efficiency of the Aggregation- Volume- Bias Monte Carlo Algorithm. J. Phys. Chem. B 2001, 105, 11275−11282. (40) Wierzchowski, S.; Kofke, D. A. A general-purpose biasing scheme for Monte Carlo simulation of associating fluids. J. Chem. Phys. 2001, 114, 8752−8762. (41) Shirts, M. R.; Chodera, J. D. Statistically optimal analysis of samples from multiple equilibrium states. J. Chem. Phys. 2008, 129, 124105. (42) Efron, B. Computers and the Theory of Statistics: Thinking the Unthinkable. SIAM Rev. 1979, 21, 460−480. (43) Allen, M. P.; Tildesley, D. J. Computer simulation of liquids; Oxford University Press, 1989.

F

DOI: 10.1021/acs.jced.8b01169 J. Chem. Eng. Data XXXX, XXX, XXX−XXX