Rosenfeld functional for non-additive hard spheres

Matthias Schmidt¹

Soft Condensed Matter Group, Debye Institute, Utrecht University, Princetonplein 5, 3584 CC Utrecht, The Netherlands

¹On leave from: Institut für Theoretische Physik II, Heinrich-Heine-Universität Düsseldorf, Universitätsstraße 1, D-40225 Düsseldorf, Germany

Received 18 June 2004
Published 16 July 2004

Density-functional theory (DFT) is a powerful approach to study equilibrium properties of inhomogeneous systems, including dense liquids and solids of single- and multi-component substances [1]. Its practical applicability depends on the quality of the approximation to the central object, the (Helmholtz) excess free energy functional arising from the interparticle interactions. The specific model of additive hard sphere mixtures constitutes the reference system par excellence to describe mixtures governed by steric repulsion, and Rosenfeld's fundamental-measure theory (FMT) [2-5] is arguably the best available approximation to tackle inhomogeneous situations. A rapidly increasing number of applications to interesting physical problems can be witnessed [6].

The more general non-additive hard sphere mixture is defined through pair potentials between particles of species i and j, given as for r<σ_ij and 0 otherwise, where r is the centre-centre distance between the two particles, and σ_ij is the distance of minimal approach between particles of species i and j. In a binary mixture non-additivity is measured conventionally through the parameter Δ  =  2σ₁₂/(σ₁₁ + σ₂₂) - 1. The physics of non-additive hard sphere mixtures is considerably richer than that of the additive case. In particular the case of Δ>0 is striking, as small values of Δ are known to be already sufficient to induce stable fluid-fluid demixing into phases with different chemical compositions (for recent studies see e.g. [7-10]).

The treatment of general non-additivity is elusive in the FMT framework. The author is aware of successful studies only in four special cases. First, for the Asakura-Oosawa-Vrij (AOV) model [11, 12], where species 1 represents colloidal hard spheres and species 2 (with σ₂₂  =  0) represents non-interacting polymer coils with radius of gyration equal to σ₁₂ - (σ₁₁/2), an excess free energy functional was given [13]. Second, a free energy functional for the Widom-Rowlinson (WR) model, where σ₁₁  =  σ₂₂  =  0, was obtained [14]. Third, the depletion potential between two big spheres immersed in a sea of smaller spheres was obtained through `Roth's trick' of working on the level of the one-body direct correlation functional [15-17]. In this case the functional for the additive case is sufficient to obtain results, but the approach is limited to small concentration of big spheres. Fourth, in the FMT of Lafuente and Cuesta for lattice hard core models, due to an odd-even effect of the particle sizes (measured in units of lattice constants), non-additivity of the size of one lattice spacing arises [18]. This effect, however, is specific to lattice models and vanishes in the continuum limit.

The aim of the present letter is to generalize FMT for hard spheres to the case of general positive and moderate negative non-additivity and arbitrary size asymmetry. The proposed extended framework accommodates, in the respective limits, the Rosenfeld functional for additive hard sphere mixtures [2], the DFT for the extreme non-additive AOV model [13], and the exact virial expansion up to second order in densities. The structure of the theory, however, goes qualitatively beyond that of either limit.

The excess (over ideal) Helmholtz free energy functional is expressed as

where ρ_i(r) is the one-body density distributions of species i  =  1,2 dependent on position r, k_BT is the thermal energy, Φ_αβ for α,β  =  0,1,2,3 is the free energy density depending on the sets of weighted densities {n_ν⁽ⁱ⁾} for i  =  1,2, and the kernels K_αβ⁽¹²⁾(r) are a means to control the range of non-locality between unlike components and depend solely on distance r. The weighted densities are built in the usual way [2] through convolution of the respective bare density profile with appropriate weight functions:

where ν  =  0,1,2,3 labels the type of weight function, and R_i  =  σ_ii/2 is the particle radius of species i  =  1,2. The (fully scalar) Kierlik-Rosinberg form [19, 20] of the w_ν(r,R) is used in the following, as this renders the determination of the K_αβ⁽¹²⁾(r) more straightforward. The w_ν(r,R) are

where R  =  R_i, the prime denotes the derivative w.r.t. the argument, δ(·) is the Dirac distribution, and Θ(·) is the step function. Alternatively, in Fourier space the weight functions are and are given as

with the abbreviations s  =  sin(kR) and c  =  cos(kR). The kernels K_αβ⁽¹²⁾(r) in (1) can be viewed as αβ-components of a second-rank tensor

where indexing is such that the top row contains K₀₀⁽¹²⁾,...,K₀₃⁽¹²⁾, etc, and † distinguishes different elements. All K_αβ⁽¹²⁾(r) possess a range of R₁₂  =  σ₁₂ - (σ₁₁ + σ₂₂)/2, i.e. vanish for values of r beyond that distance (see figure 1 for an illustration of the length scales). The dimension of K_αβ⁽¹²⁾ is (length) ^- α - β, andhence the dimension of w_γ is (length)^γ - 3. The elements of are defined, with R  =  R₁₂>0, through (3), and furthermore

with the derivatives δ^(γ)(x)  =  d^γδ(x)/d x^γ for γ  =  3,4,5. Again, we also give the Fourier space representation (being together with (4) also valid for R  =  R₁₂<0), which reads

In order to express the dependence of the free energy density, Φ_αβ in equation (1), on the weighted densities (2) we introduce ansatz functions A_αγ⁽ⁱ⁾ for species i  =  1,2 that possess the dimension of (length)^α - 3 and the order γ in density (i.e. contain γ factors n_τ⁽ⁱ⁾). Explicit expressions for the non-vanishing terms are

The excess free energy density is then constructed as

where is the γth derivative of the zero-dimensional excess free energy as a function of the average occupation number η [3], _0d(η)  =  (1 - η)ln(1 - η) + η, and for γ  =  0. The specific form (10) ensures both that the terms in the sum in (1) possesses the correct dimension of (length) ^- 6 and that the prefactor of _0d^(γ) in (10) is of the total order γ in densities, as is common in FMT. This completes the prescription for the functional; a full account of all details, also for multi-component mixtures and for lower spatial dimensionality, will be given elsewhere.

Here we discuss some of the properties of the theory. For small densities it is straightforward to show that the correct virial expansion up to second order in densities is obtained, , where the Mayer functions, f_ij(r)  =  exp( - V_ij(r)/(k_BT)) - 1, are recovered through

where * denotes the convolution, . In the limit of an additive mixture, and hence , one finds that if and otherwise. This leads to a cancellation of one spatial integration in (1) and yields the Rosenfeld functional for a binary additive hard sphere mixture [2] with radii R₁ and R₂. In the AOV limit, , one finds that if , and otherwise. The integration over in (1) together with the kernel and the fact that the density appears linearly in Φ_αβ, see A₀₁⁽²⁾ in (8), plays the same role that building weighted densities for the polymer species in the AOV case does. The resulting functional is equal to that for the AOV model [13]. However, in the WR limit, in contrast to [14], terms higher than on the second virial level vanish. For Δ  =   - 1 the two species decouple, and which is not obeyed by the present approximation, limiting its applicability to small values of Δ if Δ<0.

We next turn to an investigation of bulk properties of the theory. To assess structure, direct correlation functions can be obtained via

which can be shown to feature Percus-Yevick- (PY-) like behaviour: c_ij⁽²⁾(r>σ_ij)  =  0. Inverting the Ornstein-Zernike (OZ) relations permits us to calculate partial structure factors, S_ij(k), and partial pair correlation functions, g_ij(r). We have carried out Monte Carlo (MC) computer simulations in the canonical ensemble with 1024 particles and 10⁵ MC moves per particle; histograms of all distances between particles yield benchmark results for g_ij(r). We have chosen an intermediate size ratio of σ₂₂/σ₁₁  =  0.5 and have considered various values of Δ from  - 0.3 to 0.5 and a range of statepoints characterized by packing fractions, η_i  =  πρ_iσ_ii³/6 for i  =  1,2. For Δ  =  0, the current DFT reproduces the solution of the PY integral equation, as the functional reduces to the Rosenfeld case (which is known to yield the same c_ij⁽²⁾(r) as the PY approximation). Results for the representative case Δ  =  0.2 at two different statepoints are shown in figure 2. The core condition, g_ij(r<σ_ij)  =  0, is only approximately fulfilled, but the overall agreement between results from theory and simulation is quite remarkable.

In principle, one could envisage that this approach permits us to study the depletion potential, V_depl⁽¹¹⁾(r), between particles of species 1 being generated by the immersion into a `sea' of particles 2 through V_depl⁽¹¹⁾(r)  =   - k_BTlng₁₁(r) for , and ρ₂  =  const. However, for the (relevant) case of small size ratios (e.g.  , see [15, 16]) already in both limits of additive hard spheres and the AOV model the results are only of rather moderate accuracy, underestimating the strength of the depletion attraction [13], similar to results from the PY approximation. However, results from the present theory obtained through the OZ route (not shown) cross over smoothly between the additive hard sphere case and the AOV case, similar to the correct behaviour [15, 16]. Hence one can conclude that the pair structure predicted by the current DFT is similar to that of the PY approximation. This is a remarkable property, and one can anticipate test-particle calculations to yield superior results.

Evaluating (1) at constant density fields yields an analytic expression for the bulk excess free energy for fluid states, . The total Helmholtz free energy is then , where Λ_i is the (irrelevant) de Broglie wavelength of species i, and V is the system volume. Via Taylor expanding F_exc in both densities one can show that it features the exact second virial coefficients (consistent with the correct incorporation of f_ij(r) on the second virial level) and also the exact third virial coefficients (see e.g. [7]) provided 2σ₁₂>max(σ₁₁,σ₂₂).

The fluid-fluid demixing spinodal can be obtained from (numerical) solution of , and the location of the critical point can be determined from minimizing one of the chemical potentials, μ₁ or μ₂, along the spinodal. Such results are compared in figure 3 to those from simulations for σ₁₁  =  σ₂₂, performed in the semi-grand ensemble by Jagannathan and Yethiraj [10] and by Góźdź [9], the latter study including a finite size analysis, for a variety of non-additivities ranging from Δ  =  0.1-1. For the highly asymmetric case of σ₂₂  =  0.1σ₁₁ results from Gibbs ensemble simulations were obtained by Dijkstra [7]. For both size ratios the strong decrease of the total critical packing fraction with increasing values of Δ, as well as the overall functional dependence, are very well described by the theory. However, the precise value at given Δ is underestimated. This behaviour is not uncommon for mean-field-like theories and is also present in the AOV case. A benefit of working on the level of the density functional is that the structure is consistent with the free energy. In figure 4 partial structure factors are shown for a range of values of Δ evaluated at the fluid-fluid critical point obtained from the free energy, and indeed .

In conclusion, having demonstrated the good accuracy of the predictions of the current theory for bulk fluid properties of the non-additive hard sphere mixture, we are confident that it is well suited to study interesting and relevant interfacial situations, such as the structure and tension of interfaces between demixed phases, wetting at substrates [21], and more. Note that any colloidal mixture interacting with soft repulsive forces, as e.g. present in charge-stabilized dispersions, can be mapped (e.g. by the Barker-Henderson procedure) onto an effective non-additive hard sphere system. Hence one can anticipate experimental consequences of the structure and phase separation predicted by the present theory. The treatment of freezing [8] requires additional contributions to the free energy functional [3, 4].

H Löwen, R Evans, R Blaak, K Jagannathan and J A Cuesta are thanked for useful comments. Support by the SFB TR6 of the DFG is acknowledged. This work is part of the research programme of FOM, that is financially supported by the NWO.

References

[1]

Evans R 1992 Fundamentals of Inhomogeneous Fluids ed D Henderson (New York: Dekker) chapter 3 p 85 

[2]

Rosenfeld Y 1989 Phys. Rev. Lett. 63 980 
CrossRefPubMed

[3]

Rosenfeld Y, Schmidt M, Löwen H and Tarazona P 1997 Phys. Rev. E 55 4245 
CrossRef

[4]

Tarazona P 2000 Phys. Rev. Lett. 84 694 
CrossRefPubMed

[5]

Cuesta J A, Martinez-Raton Y and Tarazona P 2002 J. Phys.: Condens. Matter 14 11965 
IOPscience

[6]

See e.g. the special issue on DFT of liquids, 2002 J. Phys.: Condens. Matter 14 

[7]

Dijkstra M 1998 Phys. Rev. E 58 7523 
CrossRef

[8]

Louis A A, Finken R and Hansen J P 2000 Phys. Rev. E 61 R1028 
CrossRef

[9]

Góźdź W T 2003 J. Chem. Phys. 119 3309 
CrossRef

[10]

Jagannathan K and Yethiraj A 2003 J. Chem. Phys. 118 7907 
CrossRef

[11]

Asakura S and Oosawa F 1954 J. Chem. Phys. 22 1255 

[12]

Vrij A 1976 Pure Appl. Chem. 48 471 
CrossRef

[13]

Schmidt M, Löwen H, Brader J M and Evans R 2000 Phys. Rev. Lett. 85 1934 
CrossRefPubMed

[14]

Schmidt M 2001 Phys. Rev. E 63 010101(R) 
CrossRef

[15]

Roth R and Evans R 2001 Europhys. Lett. 53 271 
IOPscience

[16]

Roth R, Evans R and Louis A A 2001 Phys. Rev. E 64 051202 
CrossRef

[17]

Louis A A and Roth R 2001 J. Phys.: Condens. Matter 13 L777 
IOPscience

[18]

Lafuente L and Cuesta J A 2002 J. Phys.: Condens. Matter 14 12079 
IOPscience

[19]

Kierlik E and Rosinberg M L 1990 Phys. Rev. A 42 3382 
CrossRef

[20]

Phan S, Kierlik E, Rosinberg M L, Bildstein B and Kahl G 1993 Phys. Rev. E 48 618 
CrossRef

[21]

Brader J M, Evans R, Schmidt M and Löwen H 2002 J. Phys.: Condens. Matter 14 L1 
Brader J M, Evans R and Schmidt M 2003 Mol. Phys. 101 3349 
Wessels P P F, Schmidt M and Löwen H 2004 J. Phys.: Condens. Matter at press