Phase Behavior and Structure of Binary Hard-Sphere Mixtures

VOLUME 81, NUMBER 11 PHYSICAL REVIEW LETTERS 14 SEPTEMBER 1998 Phase Behavior and Structure of Binary Hard-Sphere Mixtures Marjolein Dijkstra, René van Roij, and Robert Evans H. H. Wills Physics Laboratory, University of Bristol, Bristol BS8 1TL, United Kingdom (Received 17 April 1998) By integrating out the degrees of freedom of the small spheres in a binary mixture of large and small hard spheres, we derive an explicit effective Hamiltonian for the large spheres. Using the two-body (depletion potential) contribution to this effective Hamiltonian in simulations, we find stable fluid-solid and both metastable fluid-fluid and solid-solid coexistence in a mixture with size ratio q ­ 0.1. For q ­ 0.05 the solid-solid coexistence becomes stable. [S0031-9007(98)07074-4] PACS numbers: 61.20.Gy, 64.70. – p, 82.70.Dd Understanding the structure and phase equilibria of binary hard-sphere mixtures is a long-standing problem in liquid state physics. These idealized systems provide a natural reference system for determining the properties of more realistic models of mixtures of simple (atomic) fluids, of colloids and polymers, and of other colloidal systems. A contentious issue, which attracts much attention, is whether fluid-fluid phase separation occurs in this model system. A classic study [1], based on the PercusYevick approximation, showed that hard spheres mix at all state points, for any ratio of diameters q ; ss ysl . In 1991, improved integral equation studies provided evidence for a spinodal instability when q # 0.2 [2]. The main reason for the subsequent interest resides in the fact that any mechanism for a demixing transition in hardsphere systems must be purely entropic. In Ref. [2] the depletion effect was identified as the mechanism behind the possible instability. This effect, which was first invoked to explain phase separation in colloid-polymer mixtures [3], is based on the idea that clustering of the large spheres allows more free volume for the small ones which may lead to an increase of the entropy. A scaled particle theory for the free volume yielded a fluid-fluid spinodal [4]. The weakness of the integral equation and the free volume theories lies in the sensitivity of the existence and location of the spinodal instability to fine details of the theory [5]. Moreover, experimental work on colloidal systems indicates that any demixing is strongly coupled to the freezing transition [6,7], whereas both types of theories are not designed to deal with solid phases—see Ref. [8], however. One might suppose that the issue could be settled by computer simulations; however, direct simulations of a highly asymmetric binary mixture are prohibited by slow equilibration. It therefore remains an open question as to whether or not (meta)stable fluid-fluid phase separation does occur and, indeed, just what the phase diagram of the binary hardsphere mixture is when the size ratio q is small. In this Letter, we take advantage of the large size asymmetry, and integrate out the degrees of freedom of the small spheres to obtain an effective Hamiltonian for the large ones. The depletion effect is now described in terms of effective potentials between the large spheres, and for q # 0.1 we expect the pairwise contribution to dominate. The pairwise (depletion) potential is essentially attractive, and arises from an unbalanced osmotic pressure of the “sea” of small spheres when the surface-surface separation of two large spheres is # ss [3]. Provided the attraction is sufficiently strong and of sufficient range, it might drive fluid-fluid phase separation—in keeping with the classical van der Waals picture of vapor-liquid separation in a simple fluid. On the other hand, when the range of the attraction is much smaller than that of the repulsion, it is known that the vapor-liquid transition becomes metastable with respect to (w.r.t.) the fluidsolid [9] and that for very short-ranged attraction an isostructural solid-solid transition can appear in the phase diagram of a simple model fluid [10]. Since the range of the attraction in the depletion potential is # ss , and that of the repulsion is sl , one might hope to find solidsolid, in addition to fluid-fluid and fluid-solid coexistence, in hard-sphere mixtures with q # 0.1. We investigate all of these possibilities using Monte Carlo simulations for the effective one-component fluid. We consider Nl large and Ns small hard spheres with diameter ratio q in a macroscopic volume V at temperature T . The total Hamiltonian consists of (trivial) kinetic energy contributions and interaction terms H ­ Hll 1 Hls 1 Hss . It is convenient to consider the system in the sNl , zs , V , T d ensemble, in which the fugacity zs of the small spheres is fixed. The Helmholtz free energy F of this system can be written as expf2bFg ­ Trl expf2bHeff g, where Heff ­ Hll 1 V is the effective Hamiltonian of the large spheres and b ­ 1ykB T. Here, V is the grand potential of the fluid of small spheres in the external field of a fixed configuration of Nl large spheres with coordinates P hRi j; i ­ 1, 2, . . . , Nl , and is given by expf2bVg ­ `Ns ­0 zsNs Trs expf2bsHls 1 n times the Hss dg. The trace Trn is short for 1yNn ! L3N n volume integral over the coordinates of the particles of species n, where Ln is the thermal wavelength. Once V, and thus Heff , are known for all values of zs , the thermodynamics and the phase behavior of the mixture can be determined. To this end, we expand 2268 © 1998 The American Physical Society 0031-9007y98y81(11)y2268(4)$15.00 VOLUME 81, NUMBER 11 PHYSICAL REVIEW LETTERS expf2bVg in terms of the Mayer functions associated with the pair-potentials fls and fss . After taking the logarithm, and using standard diagrammatic techniques [11,12], the resulting terms of the diagrammatic expansion of 2bV can be classified according to the number n ­ 0, 1, 2, . . . , Nl of large hard spheres that interact simultaneously with the sea of small spheres, so that PNl bVn . We give expressions for bVn for bV ­ n­0 n ­ 0, 1, and 2, and argue that higher order (3-body, 4-body, etc.) terms should not be important for highly asymmetric binary hard-sphere mixtures. V0 ­ 2pszs dV is the grand potential of a pure system of small spheres at fugacity zs in a volume V , where pszs d is the pressure of that system. We can also show that V1 yNl is the grand-potential difference between a system in a volume V at fugacity zs with and without a large sphere at the origin. An accurate approximation for V1 is given in [13]. This consists of a volume, a surface, and a term K which is independent of sl : V1 yNl ­ pszs dpsl3 y6 1 gszs dpsl2 1 Kszs d, where gszs d is the surface tension of the small hard-sphere fluid at a hard spherical wall of diameter sl . Explicit scaled particle results are given for p, g, and K in Ref. [13]. Within the same formalism, PNl V2 can be written as a sum of pairpotentials V2 ­ i,j fdep sRij ; zs d, where we can show that fdep is the grand potential difference between a sea of small spheres at fugacity zs containing two large spheres separated by a distance Rij ; jRi 2 Rj j and by infinite distance. It follows that fdep can be identified with the standard definition of the depletion potential [14,15]. Several theories and approximations exist for fdep and these are summarized in Ref. [16]. Recently, Mao et al. [17] used a virial expansion to calculate fdep , within the Derjaguin approximation. For q ­ 0.1 and for values of hsr as high as 0.37 their results at third order in hsr were in remarkable agreement with those of a simulation [15] for separations in the important depletion range: sl , Rij , sl 1 ss . Here hsr is the packing fraction of a reservoir of small hard spheres at fugacity zs . We use a simpler third order expression, derived in [16], which provides an equally accurate account of the simulation data: 11q f3l2 hsr 1 s9l 1 12l2 d shsr d2 bfdep sRij d ­ 2 2q 2 1 s36l 1 30l d shsr d3 g for 21 , l , 0 , (1) where l ­ Rij yss 2 1yq 2 1. This form describes a deep and very narrow potential well close to the surface of the large sphere, whose depth increases with increasing hsr , followed by a small repulsive barrier. We set fdep ­ 0 for Rij . sl 1 ss , neglecting any weak oscillations in this range [15,17]; we expect these to be unimportant for the phase behavior of the mixture. Note that the first term in (1) corresponds to the Asakura-Oosawa approximation to fdep . 14 SEPTEMBER 1998 In what follows, we set Vn ­ 0 for n $ 3. The neglect of 3-body and higher potentials can be supported by geometric arguments for q , 0.154, since then three or more nonoverlapping large spheres cannot simultaneously overlap with a small one [15,18]. Moreover, analysis of simulation data [15] implies that pairwise additivity should be an excellent approximation for q # 0.1, even for high packing fractions hl of the large spheres. We thus arrive at P the effective one-component Hamiltonian Nl feff sRij d, where H0 ­ 2pszs d s1 2 Heff ­ H0 1 i,j 2 hl dV 1 gszs dpsl Nl 1 KNl is irrelevant for the phase equilibria, although it does contribute to the pressure of the mixture. The effective pair-potential is feff ­ fll 1 fdep . Note that the formalism of mapping the two-component system onto an effective one-component system is not restricted to hard spheres. At first sight, one might think that the phase behavior of this effective one-component system can be determined by standard perturbation theory based on the pure hard-sphere reference system. Indeed this was the approach adopted in earlier studies [3,18] of colloid-polymer mixtures based on the Asakura-Oosawa results for fdep . Using first order theory for q ­ 0.1, we do not find any indication for a fluid-fluid spinodal (see also [16]). However, when simulations are performed using fdep , we find that the radial distribution function gsrd differs enormously from that of the reference hard-sphere fluid. This is illustrated in Fig. 1, where we plot gsrd for hl ­ 0.35 and hsr ­ 0.25 [19]. We find that gssl d , 42, which should be compared with the much lower contact value, ,3, for the hard-sphere reference system at the same hl . Similar large contact values, or strong tendencies for clustering, have been observed in previous simulation and integral equation studies [15,20], as well as in experiments on colloidal hard-sphere mixtures [21]. The vast 3 50 40 30 2 g(r/σl) 20 10 0 0 1 2 3 4 5 6 1 with repulsive barrier without 0 1.0 1.5 2.0 2.5 r/σl FIG. 1. The radial distribution function gsrysl d for the effective one-component system with packing fractions hl ­ 0.35, hsr ­ 0.25, and size ratio q ­ 0.1 using fdep (1), with and without the small repulsive barrier. 2269 VOLUME 81, NUMBER 11 PHYSICAL REVIEW LETTERS difference between gsrd of the reference hard-sphere system and that of the effective system signals the breakdown of perturbation theory, and we thus resort to full numerical simulations for the free energy F of the effective system. Before describing the results, we compare gsrd for a depletion potential with and without the repulsive barrier. Figure 1 shows that the contact value and most other features are not sensitive to the barrier. The small well near r ­ 1.07sl does reflect the presence of the barrier, but the free energies calculated from the two potentials differ only slightly. We conclude that small differences in the choice of depletion potential should not have a drastic effect on the resulting phase equilibria. We calculate F from Monte Carlo simulations using the thermodynamic integration technique [9]. For a given zs , i.e., a given hsr , the integration path starts at a hard-sphere fluid or solid (fcc) at the required hl , proceeds by gradually switching on the depletion potential (1), and finishes at the full effective one-component system. For the free energy of the hard-sphere reference system, we use the CarnahanStarling expression for the fluid, and the equation of state proposed by Hall [9] for the solid phase. In the latter case, an integration constant was determined such that the known fluid-solid coexistence of the pure hard-sphere system is recovered [22]. Simulations were performed in a similar fashion to those in Ref. [9]. In Fig. 2, we plot F as a function of hl at several hsr for q ­ 0.1. For hsr . 0.06 we find that the solid branch of F becomes nonconvex, indicative of a spinodal instability. For hsr . 0.29 another spinodal instability is found on the fluid branch. This instability can be seen clearly in the inset of Fig. 2, where hsr ­ 0.31. Note that subtracting a linear function of hl does not affect the common tangent construction. This simultaneous existence of a fluid-fluid and a solid-solid spinodal instability has not been observed before for binary hard-sphere mixtures. We fitted polynomials to F and computed the pressure and chemical potential at each hl . The densities of two coexisting phases can then be determined by equating the pressures and chemical potentials in the two phases. In Fig. 3, we show the resulting phase diagram. At hsr ­ 0, we find the usual freezing transition of the pure hard-sphere system. As hsr increases, a widening of the fluid-solid transition occurs, implying that a fluid with a low packing fraction of large spheres coexists with a dense crystal. This observation is consistent with results of a perturbation theory [18] and with experiments on colloids, where adding small amounts of small spheres induces a rapid decrease in the lattice constant of the crystal [7]. The results of free energy calculations reveal that both the isostructural solid-solid and the fluid-fluid transitions are metastable with respect to the wide fluid-solid coexistence, although the critical point of the solid-solid spinodal near hl ­ 0.63, hsr ­ 0.06 is close to the stable fluid-solid phase boundary. In Fig. 3, we also plot experimental state points for a colloidal hard-sphere system with q , 0.1075 [7]. In order to convert the experimental packing fraction 2270 6 4 * βF 14 SEPTEMBER 1998 r ηs = 0.09 r ηs = 0.10 r ηs = 0.12 r ηs = 0.13 r ηs = 0.31 2 0 0 -5 -10 -15 -20 0.0 0.2 0.4 ηl 0.6 0.8 FIG. 2. Reduced free energy bF p ­ bfF 2 sV0 1 3 V1 dgsl yV versus hl for q ­ 0.1 at several hsr . The curves for hl $ 0.54 are the solid branches, while the curves for lower hl are the fluid branches. Note the difference in scale for hsr ­ 0.31. For clarity, we subtracted a linear fit in hl to the data for hsr ­ 0.31 (see inset). of small spheres to hsr , we used Ns ­ 2≠bFy≠ log zs , 2≠bsV0 1 V1 dy≠ log zs , and scaled particle expressions for pszs d and gszs d [13]. (The details and accuracy of this procedure will be discussed elsewhere.) The following correlations are striking: (i) The crosses [denoting a (meta)stable fluid state] and the dashed line (the experimental estimate of the binodal) are, for hsr , 0.12, close to our stable fluid phase boundary. At higher hsr , there is a substantial deviation. Further details and a discussion of the origin of this deviation will be given elsewhere [12]; (ii) the open squares, denoting state points that correspond to fluid-solid coexistence, are well inside our fluid-solid coexistence region, and these appear to extend to large values of hl provided there is no intervening metastable fluidfluid or solid-solid binodal; (iii) the triangles, denoting the glassy states, are all within or close to the metastable fluidfluid or solid-solid binodal. A similar correlation between the formation of nonequilibrium phases and the presence of metastable phase coexistence has been observed in experiments on colloid-polymer mixtures [23]. We also note that the state point (hl ­ 0.1, hsr , 0.367) investigated in the effective one-component simulations of Ref. [15] lies well inside our metastable fluid-fluid binodal, which might explain the observed two-stage demixing dynamics. It is tempting to argue that the rapid clustering at the first stage reflects the fluid-fluid binodal, while the subsequent slow relaxation of clusters reflects the ultimate crystallization process. VOLUME 81, NUMBER 11 PHYSICAL REVIEW LETTERS In conclusion, we have shown that our effective onecomponent treatment predicts that a binary hard-sphere mixture with size ratio q ­ 0.10 exhibits a stable fluidsolid, a metastable fluid-fluid, and a metastable solid-solid coexistence. For q ­ 0.05, the solid-solid coexistence becomes stable w.r.t. fluid-solid near hsr ­ 0.05, while the fluid-fluid coexistence remains metastable. This work was supported by Grants No. ERBFM-BICT972446, No. EPSRC GRyL89013, and No. ERBFM-BICT971869. R. E. has benefited from many discussions with B. Götzelmann. 0.35 F+F 0.30 0.25 r ηs 0.20 0.15 F+S 0.2 0.10 0.1 F+F S+S F+S F S+S 0.05 0.0 0.00 0.00 0.00 F 0.25 S S 0.50 0.25 14 SEPTEMBER 1998 0.75 ηl 0.50 0.75 FIG. 3. Phase diagram for q ­ 0.1 and q ­ 0.05 (inset) in the hsr 2 hl plane. F and S denote the stable fluid and solid (fcc) phase. F 1 S, F 1 F, and S 1 S denote, respectively, the stable fluid-solid, the metastable fluid-fluid, and (meta)stable solid-solid coexistence regions. The triangles, open squares, and crosses are experimental state points taken from Ref. [7], representing glassy states, fluid-solid demixing, and (meta)stable fluid, respectively. The thin dashed line denotes the fluid branch of the experimental binodal [7]. We used the same procedures to compute the phase diagram for q ­ 0.05, which is also shown in Fig. 3. The most striking feature is the downward shift of the solid-solid binodal compared to that of the fluid-solid. This gives rise to a stable solid-solid coexistence in a small regime, where hsr , 0.05, whose critical point is shifted towards close packing. Such a trend is consistent with studies [10] for a square-well fluid in which the width of the well is reduced. Note that decreasing q is equivalent to decreasing the range of the attraction. Because the fluid-solid coexistence also shifts down as q is decreased, fluid-fluid coexistence remains metastable at q ­ 0.05, and we find that for hsr $ 0.11 coexistence occurs between an extremely dilute fluid and a solid cp whose density approaches that of close packing (hl ­ 0.7405). Note that for q ! 0 the phase diagram should approach that of the sticky-sphere model [15]. This model shows for finite T , equivalent to hsr . 0, coexistence of an infinitely dilute gas and a close-packed solid, and a metastable solid-solid transition at close packing [10]. Finally, we note that recent simulations for the actual two-component system have been carried out [20] for several q for the single state point hl ­ hs ­ 0.1215, which we convert to hsr , 0.136. The results revealed no, weak, and strong tendencies to demix for q ­ 0.1, 0.05, and 0.033, respectively. According to Fig. 3 for q ­ 0.05, the clustering observed in [20] should be associated with crystallization rather than with fluid-fluid demixing, since the state point falls outside the fluid-fluid but inside the fluid-solid binodal. [1] J. L. Lebowitz and J. S. Rowlinson, J. Chem. Phys. 41, 133 (1964). [2] T. Biben and J. P. Hansen, Phys. Rev. Lett. 66, 2215 (1991). [3] S. Asakura and F. Oosawa, J. Chem. Phys. 22, 1255 (1954); A. Vrij, Pure Appl. Chem. 48, 471 (1976). [4] H. N. W. Lekkerkerker and A. Stroobants, Physica (Amsterdam) 195A, 387 (1993). [5] For recent results and a summary of earlier work, see T. Coussaert and M. Baus, Phys. Rev. Lett. 79, 1881 (1997). [6] A. D. Dinsmore, A. G. Yodh, and D. J. Pine, Phys. Rev. E 52, 4045 (1995). [7] A. Imhof and J. K. G. Dhont, Phys. Rev. Lett. 75, 1662 (1995); Phys. Rev. E 52, 6344 (1995). [8] W. C. K. Poon and P. B. Warren, Europhys. Lett. 28, 513 (1994). [9] M. H. J. Hagen and D. Frenkel, J. Chem. Phys. 101, 4093 (1994). [10] P. Bolhuis and D. Frenkel, Phys. Rev. Lett. 72, 2211 (1994); C. F. Tejero et al., Phys. Rev. Lett. 73, 752 (1994). [11] J. P. Hansen and I. R. McDonald, Theory of Simple Liquids (Academic, London, 1986). [12] M. Dijkstra, R. van Roij, and R. Evans (to be published). [13] J. R. Henderson, Mol. Phys. 50, 741 (1983). [14] P. Attard, J. Chem. Phys. 91, 3083 (1989). [15] T. Biben, P. Bladon, and D. Frenkel, J. Phys. Condens. Matter 8, 10 799 (1996). [16] B. Götzelmann, R. Evans, and S. Dietrich, Phys. Rev. E 57, 6785 (1998). These authors discuss limitations of the virial expansion for fdep . [17] Y. Mao, M. E. Cates, and H. N. W. Lekkerkerker, Physica (Amsterdam) 222A, 10 (1995). [18] A. P. Gast, C. K. Hall, and W. B. Russel, J. Colloid Interface Sci. 96, 251 (1983). [19] Note that this state point falls inside the fluid-solid coexistence region [see Fig. 3 (below)]. The peaks in gsrd near rysl ­ 1.74 and 2 are similar to those found at state point “C” of Ref. [15]. [20] A. Buhot and W. Krauth, Phys. Rev. Lett. 80, 3787 (1998). Phase equilibria were not determined here. [21] E. K. Hobbie, Phys. Rev. E 55, 6281 (1997). [22] W. G. Hoover and F. H. Ree, J. Chem. Phys. 49, 3609 (1968). [23] W. C. K. Poon, A. D. Pirie, M. D. Haw, and P. N. Pusey, Physica (Amsterdam) 235A, 110 (1997). 2271