Effective forces in square well and square shoulder fluids
EEffective forces in square welland square shoulder fluids
D. Fiocco, † G. Pastore, ‡ and G. Foffi ∗ , † Institute of Theoretical Physics (ITP), Ecole Polytechnique Fédérale de Lausanne (EPFL),CH-1015 Lausanne, Switzerland, and Dipartimento di Fisica dell’Università di Trieste andCNR-IOM, Strada Costiera 11, 34151 Trieste, Italy
E-mail: giuseppe.foffi@epfl.ch
Abstract
We derive an analytical expression for the effective force between a pair of macrospheresimmersed in a sea of microspheres, in the case where the interaction between the two unlikespecies is assumed to be a square well or a square shoulder of given range and depth (or height).This formula extends a similar one developed in the case of hard core interactions only. Qual-itative features of such effective force and the resulting phase diagram are then analyzed in thelimit of no interaction between the small particles. Approximate force profiles are then ob-tained by means of integral equation theories (PY and HNC) combined with the superpositionapproximation and compared with exact ones from direct Monte Carlo simulations.
Colloidal systems are ubiquitous and have attracted a growing interest in the last decades. Theimpact on everyday life is vast and ranges from food to materials, from biology to photonic ∗ To whom correspondence should be addressed † Ecole Polytechnique Fédérale de Lausanne (EPFL) ‡ Dipartimento di Fisica dell’Università di Trieste and CNR-IOM a r X i v : . [ c ond - m a t . s t a t - m ec h ] J u l rystals. On a more fundamental level, experiments on colloids have been the testing ground forseveral physical phenomena, from thermodynamics to the glass transition. In general, the in-terest comes from the possibility of tuning macroscopic properties by encoding them at the levelof the interactions between the constituents of colloidal solutions. For these reasons these systemshave dominated the scientific discussions in these fields during the last years. In particular, thepossibility of engineering systems with tunable microscopic properties allowed the observation ofexotic phenomena like the Alder transition for hard spheres, the metastable liquid-liquid phaseseparation, the existence of two distinct glassy phases for short ranged attractive potentials andthe emergence of a thermodynamically stable cluster phase. As most soft-matter systems, colloids are characterized by a large number of degrees of freedomspanning several length and time scales. In several circumstances, however, it is possible to inte-grate out some of these degrees of freedom to gain more insight. In fact, colloids are one of theprototypical systems where a coarse-grained approach can be extremely fruitful. Typically, theextra degrees of freedom are mapped into some effective interaction among colloidal particles, enabling the use of the full arsenal of statistical mechanics. In particular, most of the results forsimple liquids (integral equations, perturbation theory, mode coupling theory, etc.) can be success-fully used to address questions of great relevance for colloidal systems.A classical example of such success is the case of depletion interactions. In the case of a binarymixture of colloidal hard spheres, when the diameters of the two species are very different, anentropic force starts to set in. This results in a net attraction between the colloids belonging tothe largest species. Asakura and Oosawa and, independently, Vrij derived an effective interac-tion potential casting the problem of a binary mixture into an effective single-component system.This very simple result was the beginning of a series of investigations (see for instance (16,17)).It was possible, for the first time, to investigate a system of particles interacting with an attractionwhose range could be tuned, something that is not feasible for atomic or molecular systems. Theprofound effects both on thermodynamics and the dynamics of a short range attraction are nowwell established. 2he use of effective interactions, far from being restricted to colloidal systems, has important im-plications also for proteins, especially for what concerns crystal nucleation and phase behavior. Bymodeling proteins as short range attractive colloids it has been possible, for example, to rationalizethe enhancement of crystal nucleation in the proximity of the liquid-liquid critical point. Thisextremely simplified modeling has been also useful in the case of binary mixtures of eye-lens pro-teins, whose experimental behavior was successfully modeled using a simple hard sphere potentialwith specific square well attractions.
Here, a fundamental role was played by the intensity ofthe interspecies attraction: by varying it, it was possible to reduce the instability with respect todemixing due to depletion. Apart from the medical and biological implications, this work showedhow a modification of the mutual interaction among the components can alter the phenomenologyof the whole mixture.The aim of this paper is to rationalize this fact in terms of effective interactions. We investigatedthe effective interactions between the larger component of an asymmetric binary mixture of hardspheres where the potential between the two components has the form of a square well or a squareshoulder. To derive our results we followed a reasoning similar to Attard’s derivation of the effec-tive force between two hard spheres immersed in a sea of smaller ones. We generalized Attard’sargument to the case of square well (or shoulder) interspecies interaction. The calculation of theactual force, as in (21), requires the evaluation of the density profile of the smaller spheres arounda pair of larger ones. We performed this step following different methods. We shall first proposea simple approximation that treats the smaller component as an ideal gas. In this way, we shallfind a generalization of the Asakura-Oosawa (AO) analytical expression of the force to the case inwhich the mutual interaction potential is an attractive square well or a repulsive shoulder. Withinsuch approximation we’ll explore, by means of first order thermodynamic perturbation theory, thequalitative changes in the phase diagram determined by the mutual interaction. The analogy withthe case of binary mixtures of eye-lens proteins are evident. As for the purely repulsive case, theapproximation above holds only when the interactions between the smaller particles are negligi-ble. Thus we improved the results of the force calculation using an integral equation approach3ombined with the superposition approximation for the case of hard sphere interaction betweenthe smaller species. All the aforementioned theoretical results have been tested with respect to the“exact values” of the force as calculated directly from new Monte Carlo simulations.The paper is organized as follows: in Section 2 we present our generalization of Attard’s deriva-tion of the effective force in the case of square well/shoulder interactions between the species ofa binary mixture; in Section 3 the methods used to estimate the density profile needed to com-pute the effective force are described. Details about the parameters used are also provided and adescription of the qualitative effect of the effective interactions on the phase diagram is given; inSection 4 quantitative results obtained with the different methods are shown and discussed. Finallyin Section 5 conclusions of this work are drawn.
Consider an asymmetric binary mixture (we label with 1 the larger component and with 2 the otherone) and let the interaction potential between two particles of the larger species φ ( r ) and theinterspecies potential φ ( r ) have a square-well or square shoulder form, that is φ i j ( r ) = ∞ if r < σ i j − ε i j if σ i j ≤ r < λ i j σ i j λ i j σ i j ≤ r , (1)where σ i j denotes the range of the hard core interaction, ε i j is a positive (negative) energy mea-suring the depth (height) of the well (shoulder) and λ i j determines the width of the well (shoulder)with respect to the range of the hard core. We fix ε = φ ( r ) unspecified. We now wish to calculate the effective force between twoparticles of type 1 immersed in a sea of particles of type 2.The problem in the particular case of ε = in the past. If two spheres of type 1 are fixed in 000 (the origin) and RRR respectively, the4alue f of the radial component of the effective force between each other can be written as (21): f ( R ) = − πβ σ (cid:90) π ρ ( SSS σ ; R ) cos θ sin θ d θ , (2)where β is the inverse temperature 1 / k B T , SSS σ is a vector of length σ whose tail is in 000 and θ is such that SSS σ · RRR = R σ cos θ (the meaning of the parameters that appear in Eq. (2) is furtherelucidated in Appendix A and Figure 1). Thus the force depends on the density ρ ( rrr ) of smallparticles around the pair of macrospheres.In Appendix A a derivation similar to Attard’s is used to obtain the expression of the effectiveforce between two macrospheres when the interaction between the species 1 and 2 is of the kinddescribed by Eq. (1): f ( R ) = − πβ (cid:18) σ (cid:90) π ρ ( SSS σ ; R ) cos θ sin θ d θ +( − e β ε ) λ σ (cid:90) π ρ ( SSS λ σ ; R ) cos θ sin θ d θ (cid:19) , (3)where the shorthand notation λ σ = λ σ has been used, SSS λ σ is defined in a fashion similarto SSS σ and a representation of all the parameters is given in Figure 1.Our expression differs from the hard sphere case ( Eq. (2)) because of the presence of a secondterm. This is due to the introduction of the square well (shoulder) of finite depth (height) | ε | andhas the same form of the previous term but a weight 1 − e β ε , which is negative for square wellpotentials. This dependence of the force on the temperature is a by-product of the non-zero valueof ε , which introduces a natural energy scale that is absent in the athermal hard sphere system.Thus, the force is not entirely entropic as in Eq. (2) but has also an energetic origin. The expressionof the force in Eq. (3) correctly reduces (as it should) to Attard’s formula in the ε = σ = λ σ in the limit ε → − ∞ (infinitelyhigh shoulder). We emphasize the fact that Eq. (3) yields the exact force if the correct form ofthe density is known regardless the interaction potential between the small particles. In this paper5e shall restrict ourselves to the two cases in which the small particles are either non-interacting( φ ( r ) =
0) or hard spheres (so that φ ( r ) will take the form of Eq. (1) with ε = In order to determine the effective force the density ρ ( rrr ) was obtained using three different meth-ods: the dilute gas approximation, integral equations together with the superposition approximationand Metropolis Monte Carlo simulations. If species 2 is dilute we can approximate the density around a pair of macrospheres with that of anideal gas: ρ DL ( rrr ; R ) = ρ e − β V ( rrr ; R ) , (4)where V ( rrr ; R ) is the potential energy of a particle of type 2 centered in rrr due to the presence of thefixed pair of type 1 (see also Eq. (21)).This is equivalent to the Asakura-Oosawa approximation in the case ε =
0. We expect suchapproximation to be less accurate in more coupled regimes but still able to give qualitative infor-mation about the force profile. In App. Appendix B we make use of Eq. (4) and obtain an analytic expression of the effective force. For a wide choice of the parameters ε , σ , λ σ it has featureslike those appearing in Figure 2. The profile obtained when ε = The effect of negative valuesof ε is that of increasing both the range and the strength of the effective force, approaching the ε = − ∞ solution mentioned above. In this case the qualitative features are essentially the same asin the ε = Positive values of ε determine a completely different behavior: increasinglylarge values determine the onset of a repulsion at short distances (maximum repulsion occurringat R = σ ) and of a strong attractive force at intermediate distances (with maximum attractionfound at R = σ + λ σ ). 6he physical origin of such form of the effective interaction can be understood (at least quali-tatively) by means of the argument below.For sake of simplicity, let’s consider the ε = repel each other when they come into contact, i.e. when they are at a distance r = σ . Itfollows that if a macrosphere is inserted in a sea of small ones, it will experience on average a forcedepending on the density of small spheres at distance r = σ . If a macrosphere is isolated fromothers the density distribution of small ones will be spherically symmetric so it won’t experienceany net force. If two macrospheres come close together, however, the hemispheres at r = σ fac-ing each other become depleted of microspheres: the “push” from the two external hemispheres isnot balanced anymore and the result is an effective attraction. The situation is depicted in Figure 3.Now consider the case when a square well or shoulder is also present. We have to take intoaccount that a macroparticle and a microparticle experience an attraction or a repulsion dependingon the sign of ε when their separation is λ σ . Thus, the effect of a well is to make a macrospherebeing attracted by the density of small spheres localized at a distance r = λ σ , while a shoulderwill make the macrosphere be pushed away from it. Taking this into account we can understandthe effect of ε on f DL ( R ) .It’s easy to realize (see also Figure 3) that the case ε < ε =
0, with the novel contribution due to the density localized at the outer rim of the shoulderand a weaker contribution from the hard wall due to lower density of small particles inside theshoulder.The situation is more complicated if ε > σ + λ σ < R < λ σ part of the outer rim of the well is inside the well of the neighboring bigparticle and there is a strong “pull” due to the high density of small particles here (green arrows inFigure 4(b)). The result is an effective attraction. At lower R such attraction is counterbalanced bya strong “push” due to the fact that the density at the hard wall is higher where the two cores faceeach other (see blue arrows in Figure 4(c)). At the same time the attractive “pull” at the outer rimof the well (green arrows on the right in Figure 4(c)) is weaker because part of it is in the region not7ccessible to the small spheres ( ρ = R and big enough ε repulsiondominates over attraction ( Figure 4(d)). Qualitative dependence of the phase diagram from ε We studied the effect of the presence of the smaller species on the phase diagram of the largercomponent, to see if the stability of the system can be tuned by means of varying the parameter ε .Effects of the mutual attraction on the stability of this class of systems has previously been reportedin the literature. In the semi-grand canonical ensemble we can write for the thermodynamicpotential F ( N , V , z ) e − β F = N ! Λ Tr e − β H eff , (5)where N is the number of large particles, V is the volume of the system, z is the fugacity of thesmaller species, Tr denotes the integration operator over the coordinates of the particles of thelarger species (cid:82) dddRRR . . . dddRRR NNN , Λ ≡ h / (cid:112) π m / β its thermal wavelength and H eff is the effectivepotential. H eff in turn can be written as H eff = H + Ω , (6)where H is the sum of the pair interactions between the larger particles ∑ Ni (cid:54) = j φ ( R ) and Ω is thegrand potential of the smaller species in a fixed configuration of the large particles. It can be shownthat Ω can be expressed as a sum of n -body terms Ω = ∑ n Ω n (7)and it is straightforward (see Appendix C) to prove that the 0- and 1-body terms are linear in thedensity of the effective component ρ and therefore do not alter the phase behavior in the caseexamined here. In what follows the 2-body term is assumed to be given by8 = N ∑ i (cid:54) = j V DL ( R i j ) , (8)where V DL ( R ) = − (cid:90) f DL ( R ) dR + c , (9)where c is chosen to guarantee that lim r → ∞ V DL ( r ) = can thus be employed to explore the phasebehavior of the effective component, in a way similar to that applied by Gast and coworkers in theAO case ( ε = We thus write the free energy per particle of the effective component F eff / N with β F eff / N = β F HS / N + β ρ (cid:90) V DL ( r ) g HS ( r ) π r dr , (10)where F HS / N and g HS are respectively the free energy per particle and the pair distributionfunction associated to the hard sphere reference system. The former is approximated here bymeans of the Carnahan-Starling expression and the Verlet-Weis description is used for thelatter. A check of the validity of the first-order approach can be performed by verifying that theBarker-Henderson second-order correction to F eff / N is reasonably smaller than the secondterm on the RHS of Eq. (10). The chemical potential µ and the pressure p can be expressed with β µ = ∂∂ ρ ( ρ β F eff / N ) , (11) β p = ρ β µ − ρ β F eff / N . (12)Coexistence lines are found by plotting the value of the volume fraction η of the bigger spheresat which the parametric curve in the ( p , µ ) plane self-intersects for various values of the volumefraction η of the smaller particles. 9 .2 Integral equations and superposition approximation (PY + S and HNC+ S) methods Another approximate way to determine the density is using integral equation theory for a fluidmixture in order to get the unlike (big-small) pair distribution g ( r ) . This in turn can be pluggedinto the superposition approximation for the density profile around the pair of macrospheres ρ S ( rrr ; R ) = g ( r ) g ( | rrr − RRR | ) ρ . (13)Note how the superposition approximation amounts to saying that the density around the pair isequal to the product of the densities that the particles would have around themselves if they wereisolated.The pair distribution function g in Eq. (13) can be determined solving the two-componentOrnstein-Zernike equation h i j ( r ) = c i j ( r ) + ∑ k = (cid:90) ρ k c ik ( (cid:12)(cid:12) rrr − rrr (cid:48)(cid:48)(cid:48) (cid:12)(cid:12) ) h k j ( r (cid:48) ) dddrrr (cid:48)(cid:48)(cid:48) , (14)and the closure equation g i j ( r ) = e − β φ i j ( r )+ h i j ( r ) − c i j ( r ) − E i j ( r ) , (15)where i , j , k ∈ { , } , h i j , c i j , E i j ( r ) are respectively the indirect correlations, direct correlationsand bridge functions, ρ k is the density of the species k and ∑ k ρ k = ρ . In order to solve theintegral equations we implemented a version of Gillan algorithm. Special care was taken to treatdiscontinuous potentials of the form of Eq. (1).The well-known Percus-Yevick closure (PY) E i j ( r ) = ln [ + h i j ( r ) − c i j ( r )] − h i j ( r ) + c i j ( r ) (16)10nd the hypernetted-chain closure (HNC) E i j ( r ) = The exact density can be sampled using Metropolis Monte Carlo simulations. The force can then beobtained with the method used in the two-sphere studies in (22) (a possibly more efficient alterna-tive method is that described in (29)). We briefly summarize such method here. Two macrospheresof diameter σ are inserted at fixed positions ( , , ) and ( R , , ) in a cell whose dimensions are H along the x direction and L in the y and z directions. The same cell contains N smaller particleswhich interact with each other with a potential of the form of Eq. (1) such that ε = H and L arechosen so that the density profile is flat away from the pair of macrospheres (i.e. in what we canconsider to be the bulk ) and is equal to ρ . To keep ρ constant for each value of R , the number N istuned to compensate the variation in the volume accessible to the small spheres. Periodic boundaryconditions and the NVT ensemble are used.At the beginning of the simulation microparticles are placed randomly. At each MC step a mi-crosphere is selected at random and a random displacement is attempted. The displacement isaccepted or rejected according to the Metropolis scheme, the displacement is tuned to reach an ac-ceptance ratio of 0.25 and cell lists are used in order to improve the efficiency of the algorithm. Configuration samples are taken after the mean square displacement of the microspheres equals σ or after a number of moves sufficient to decorrelate the total energy per particle (in the SW andSS case). The value of the integrals that appear in Eq. (3) can be obtained using the relation2 π v (cid:90) π ρ ( SSS vvv ; R ) sin θ cos θ d θ ≈ (cid:10) dr ∑ v < v i < v + dv cos θ i (cid:11) , (18)11here we sum all the cos θ i whose distance v i from the center of the macrosphere at position ( R , , ) is in the interval [ v , v + dv ] and the angle brackets denote an average over all the samplescollected during a simulation. To obtain a better estimate of the integrals in Eq. (3) a third-orderpolynomial was fitted with the RHS of Eq. (18) corresponding to different v ’s and extrapolated thecurves to v = σ + and v = λ σ + . The errors of the estimates are evaluated on the basis of the errorsof the fit parameters. Data points in the fits were weighted according to the statistical uncertaintiesof the measured values of the RHS of Eq. (18). Plots of the density
It’s useful to obtain plots of the density obtained with the different methods in order to highlightany differences between their predictions. Such plots can be obtained from MC simulations bydividing the simulation box in cells, counting the number of particles contained in each of them,and averaging on multiple configurations. This data can be projected in 2D afterwards exploitingthe symmetries of the system (for example ρ ( rrr ; R ) ≡ ρ ( r , θ , φ ; R ) = ρ ( r , θ ; R ) if φ defines a rotationaround the direction of RRR ).Plots of the density associated to the integral equation and superposition method can be ob-tained using Eq. (13) by means of sampling it on a very fine mesh in real space. Down-samplingto the same mesh used with MC data produces results that can be compared to those obtained withthe MC method.
We used the same geometries analyzed by Dickman et al. in 22, i.e. two size ratios ξ ≡ σ / σ = 5, 10, taking φ to be of the form of Eq. (1) with σ = ε = η = πρσ / = . , . , . λ σ = . , . ξ = ,
10) and whose depth (height) | ε | = / β = H = L =
16 for ξ = H = L = ξ = reduced force f ∗ MC ( R ) = β f MC ( R ) / ( πρσ ) for values of R ranging from σ to σ + . σ at regular intervals of 0 . σ were obtained. Each data point took about 10-20 hoursof CPU time on an Intel Xeon 3.2 GHz processor to be determined. They are shown in Figure 6and Figure 7. Data taken from (22) (relative to the case ε =
0) are also plotted for comparison.Table 1: Values of the parameters used to obtain the f ∗ MC ( R ) profiles. By a MC step here we meana single particle displacement attempt. The averages have been evaluated by using the reportedsampling frequency over the total number of MC steps shown in the last column. In all cases atleast 10 equilibration MC steps were performed. N Sampling frequency (MC steps) Production (MC steps) ξ = , η = . ≈ ξ = , η = . ≈ · · ξ = , η = . ≈ · · ξ = , η = . ≈ ξ = , η = . ≈ · · ξ = , η = . ≈ · · In the PY+S method values of g ( r ) were sampled on an equispaced ( dr = .
01) mesh of 4096or 8192 points in r -space respectively for ξ = ,
10. Values for any r were obtained from the dis-crete sample through linear interpolation (linear extrapolation was used to obtain the values neardiscontinuous points). These in turn allowed us to obtain ρ PY + S ( rrr ; R ) for any rrr through Eq. (13).The f ∗ PY + S ( R ) could finally be obtained performing a numerical integration of the RHS side ofEq. (3) for various values of R . The force profiles are shown in Figure 6, Figure 7. Each requireda few minutes of CPU time on a desktop computer. The same procedure was carried out using theHNC+S method, and results are shown in Figure 6 for the case ξ = , η = . η , η ) was obtained in the particular case of ξ = ε and is shown in Figure 5(b). The associated effective potentials (shown inFigure 5(a)) were obtained by numerical integration of samples of the force in Eq. (34). The firstorder perturbation in in Eq. (10) was also found by numerical integration, while the derivative re-quired by in Eq. (11) was obtained via numerical differentiation. The discretization of real spaceneeded to carry out such operations was performed on a grid sufficiently fine so that further refine-ments had no appreciable effect on the scale of the plots.In addition we obtained further information about the density of microspheres for the case ξ = η = . [ ρ MC ( r , θ ; R ) − ρ PY + S ( r , θ ; R )] / ρ PY + S ( r , θ ; R ) at R = . ρ MC ( r , θ ; R ) ≡ ρ MC ( r , − θ ; R ) when θ <
0, and noise at θ ≈ = π r sin θ drd θ ). The MC method allows to measure the exact effective force in the various cases. Our data arenot always in good agreement with those found in (22) for the case ε = ε = ε = ± ε = ρ MC and ρ PY + S . In this case PY+S onlyslightly underestimates the particle density in the zone close to both spheres (see Figure 8). Thefact that the density match is very good everywhere but in this region suggests that this discrepancyis due to the superposition approximation. In the ε = ± ξ = , η = . η = .
116 of Figure 6, where the PY+S and HNC+S methods bothfail to describe accurately the force profile at short ranges.The phase diagram (see Figure 5(b)) obtained from first order thermodynamic perturbation theorywithin the DL approximation shows that for small positive values of ε one needs to move athigher densities of the smaller species as ε increases in order to observe phase separation. Thisis due to the fact that in this regime an increase in ε corresponds to an additional repulsive termin the effective potential (see Figure 5(a)).Such “stabilizing effect” due to mutual attraction doesn’t hold at higher values of ε , whereincrements in ε correspond to progressively stronger effective attraction and lowering of thecoexistence line in the ( η , η ) plane. The critical density in the case under examination is higherin the low ε regime than in the high ε one. 15 Conclusions
In this paper we have studied the effective forces between two big hard spheres dispersed into afluid of smaller particles. In particular we studied the effect of interspecies interactions. To thisaim we have extended the force derivation done by Attard to the case in which these take the shapeof a square well or shoulder, that is in the case where both entropic and energetic effects play arole. As in the case of the original formula, it is sufficient to know the density profile of the smallparticles around the big ones to obtain the exact value of the effective force. The determination ofthe density however is a difficult task and one has to rely on some approximations. Here we haveproposed an equivalent of the Asakura-Oosawa approximation, i.e. the assumption that the smallparticles behave like an ideal gas. This approximation leads to an analytical expression that cap-tures, at least qualitatively, the correct behavior as predicted by our MC simulations. An approachbased on integral equations and a superposition approximation improves the results but stills failsin estimating the density in the neighborhood of the two large spheres regardless the choice of theclosure used (at least in the cases examined here). We have shown how the addition of a shoulderor of a well in the interspecies interaction can change qualitatively the well known depletion forceprofile observed in hard sphere binary mixtures. The presence of a well, for example, causes theonset of a repulsion at short ranges and, in the case of deep wells, of an attraction at intermediateranges. This happens because when a well is present the larger particles are surrounded by a layerof smaller ones. This layer makes a close contact between the large particles unfavorable, but atthe same time if two bigger spheres share part of their surrounding layers the small particles sit-ting between them act as “glue”, stabilizing a configuration of intermediate distance between thepair. Such description qualitatively agrees with what has been observed in molecular simulationsof binary mixtures of eye-lens proteins of hard-core potentials with a Yukawa tail and couldbe confirmed by experimental determinations of the force on colloids interacting with a squareshoulder/potential (polymer grafted colloids might be a valid candidate member to this class) us-ing optical tweezers. The phase diagram as studied with simple thermodynamic perturbation theory within the DL ap-16roximation shows that the introduction of shallow well (low ε ) has the effect of pushing tohigher densities of the smaller component the phase separation. Increasing further ε , however,one reaches the point where deepening the well has the opposite effect.The effect of steric hindrance of the small particles in the purely hard sphere case has alreadybeen claimed to be used to stabilize solutions and foods. It’s clear that taking into account thepossibility of tuning the interspecies interactions could broaden even more the possible routes forstabilization. Summarizing, we confirm, in agreement with previous work, that the mutualattraction could be an extra parameter to play with when tuning the stability of a binary mixture.The present work provides a qualitative and quantitative analysis of the resulting changes.
Acknowledgement
We would like to thank Phil Attard for useful discussions and Simone Belli who wrote with one ofus (D.F.) the Fortran implementation of the Gillan algorithm of solution of the IE. We thank alsoFrancesco Varrato and Nicolas Dorsaz for comments on the manuscript. D.F. and G.F. acknowl-edge support by the Swiss National Science Foundation (grant no. PP0022_119006).
A Derivation of the effective force
It is shown in (21) that the force can be expressed as f ( R ) = − ∂ φ ∂ R − ∂ F ∂ R = − ∂ φ ∂ R + β Z ∂ Z ∂ R , (19)where φ is the interaction potential between particles of type 1, F is the free energy of theparticles of type 2 that move in the potential generated by the fixed pair of type 1, Z the partitionfunction associated to it and β = / k B T .Following (21) we also have that ∂ Z ∂ R = − (cid:90) ∂∂ R (cid:16) − e − β V ( rrr ) (cid:17) e β V ( rrr ) ρ ( rrr ) Z dddrrr , (20)17here V ( rrr ) is the potential energy of a single particle of type 2 due to the presence of the fixedcouple of type 1. V ( rrr ) can be written as V ( rrr ) = ∞ if r < σ or | rrr − RRR | < σ − ε if σ ≤ r < λ σ xor σ ≤ | rrr − RRR | < λ σ − ε if σ ≤ r < λ σ and σ ≤ | rrr − RRR | < λ σ . (21)We can write the resulting Mayer function as:1 − e − β V ( rrr ) = ( − e β ε ) (cid:2) H λ σ ( rrr ) + H λ σ ( rrr −−− RRR ) (cid:3) ++ ( e β ε − e β ε − ) H λ σ ( rrr ) H λ σ ( rrr −−− RRR )++ ( e β ε − e β ε ) (cid:2) H σ ( rrr ) H λ σ ( rrr −−− RRR ) + H λ σ ( rrr ) H σ ( rrr −−− RRR ) (cid:3) + − e β ε H σ ( rrr ) H σ ( rrr −−− RRR )++ e β ε [ H σ ( rrr ) + H σ ( rrr −−− RRR )] , (22)where the geometric definition of the support of V ( rrr ) has been encoded using the characteristicfunctions H D ( rrr ) defined as: H D ( rrr ) = r < D r ≥ D . (23)Plugging Eq. (22) inside Eq. (20) and rearranging yields:18 Z ∂ R = − ( − e β ε ) (cid:90) RRRR · rrr −−− RRR | rrr −−− RRR | δ ( | rrr −−− RRR | − λ σ ) e β V ( rrr ) ρ ( rrr ) Z dddrrr + − ( e β ε − e β ε − ) (cid:90) H λ σ ( rrr ) RRRR · rrr −−− RRR | rrr −−− RRR | δ ( | rrr −−− RRR | − λ σ ) e β V ( rrr ) ρ ( rrr ) Z dddrrr + − ( e β ε − e β ε ) (cid:90) H σ ( rrr ) RRRR · rrr −−− RRR | rrr −−− RRR | δ ( | rrr −−− RRR | − λ σ ) e β V ( rrr ) ρ ( rrr ) Z dddrrr + − ( e β ε − e β ε ) (cid:90) H λ σ ( rrr ) RRRR · rrr −−− RRR | rrr −−− RRR | δ ( | rrr −−− RRR | − σ ) e β V ( rrr ) ρ ( rrr ) Z dddrrr ++ e β ε (cid:90) H σ ( rrr ) RRRR · rrr −−− RRR | rrr −−− RRR | δ ( | rrr −−− RRR | − σ ) e β V ( rrr ) ρ ( rrr ) Z dddrrr + − e β ε (cid:90) RRRR · rrr −−− RRR | rrr −−− RRR | δ ( | rrr −−− RRR | − σ ) e β V ( rrr ) ρ ( rrr ) Z dddrrr , (24)where we have used the shorthand notation λ σ to refer to λ σ .The integrals in Eq. (24) can be evaluated in the three dimensional space performing the changeof variables sss ≡ rrr −−− RRR , using spherical coordinates centered in
RRR with θ defined by sss · RRR = sR cos θ and exploiting the azimuthal symmetry of the problem.The first one thus reads (cid:90) RRRR · ssss δ ( s − λ σ ) e β V ( sss +++ RRR ) ρ ( sss +++ RRR ) Z dddsss == (cid:90) ∞ (cid:90) π (cid:90) π RRRR · ssss δ ( s − λ σ ) e β V ( sss +++ RRR ) ρ ( sss +++ RRR ) Z s sin θ ds d θ d φ == πλ σ (cid:90) π cos θ e β V ( SSS λσ ) ρ ( SSS λ σ ) Z sin θ d θ , (25)where in the last line SSS λ σ is defined with the notation SSS v ≡ vvv +++ RRR with | vvv | = v . (26)The sixth integral in Eq. (24) has the same form and its value is2 πσ (cid:90) π cos θ e β V ( SSS σ ) ρ ( SSS σ ) Z sin θ d θ . (27)19s for the second (cid:90) H λ σ ( rrr ) RRRR · rrr −−− RRR | rrr −−− RRR | δ ( | rrr −−− RRR | − λ σ ) e β V ( rrr ) ρ ( rrr ) Z drrr == πλ σ (cid:90) π H λ σ ( SSS λ σ ) cos θ e β V ( SSS λσ ) ρ ( SSS λ σ ) Z sin θ d θ . (28)whereas the fourth (cid:90) H λ σ ( rrr ) RRRR · rrr −−− RRR | rrr −−− RRR | δ ( | rrr −−− RRR | − σ ) e β V ( rrr ) ρ ( rrr ) Z drrr == πσ (cid:90) π H λ σ ( SSS σ ) cos θ e β V ( SSS σ ) ρ ( SSS σ ) Z sin θ d θ . (29)while the third and the fifth ones vanish because wherever ρ is nonzero H σ is zero and viceversa.The H ’s can be eliminated introducing the auxiliary functions ρ λ σ ( SSS λ σ ) = ρ ( SSS λ σ ) if 0 ≤ θ < θ λ σ , ρ λ σ ( SSS λ σ ) = ρ ( SSS λ σ ) if θ λ σ ≤ θ < θ λ σ ≤ π , ρ σ ( SSS σ ) = ρ ( SSS σ ) if 0 ≤ θ < θ σ , ρ σ ( SSS σ ) = ρ ( SSS σ ) if θ σ ≤ θ < θ σ ≤ π , (30)20here θ λ σ = arccos (cid:18) R λ σ R (cid:19) θ λ σ = arccos (cid:18) λ σ − R − σ − σ R (cid:19) θ σ = arccos (cid:18) σ − R − λ σ − λ σ R (cid:19) θ σ = arccos (cid:18) R σ R (cid:19) (31)are the angles at which the spheres of radii σ and λ σ centered in 000 and RRR intersect. Defining ρ ( SSS λ σ ) = ρ λ σ ( SSS λ σ ) + ρ λ σ ( SSS λ σ ) ρ ( SSS σ ) = ρ σ ( SSS σ ) + ρ σ ( SSS σ ) (32)and substituting these inside the integrals and some algebra eventually yields Eq. (3). B Expression of the force in the dilute limit
Use of Eq. (4) together with Eq. (21) allows to rewrite the densities that appear in Eq. (3) as ρ ( SSS λ σ ) = ρ if 0 ≤ θ < θ λ σ ρ e β ε if θ λ σ ≤ θ < θ λ σ , ρ ( SSS σ ) = ρ e β ε if 0 ≤ θ < θ σ ρ e β ε if θ σ ≤ θ < θ σ . (33)21here the angles are the same defined in Eq. (31). Substitution of such density inside Eq. (3) yieldsa piecewise expression of the force valid for R ∈ [ σ , + ∞ ] :if σ < R ≤ σ : f DL ( R ) = − πρβ (cid:40) e β ε (cid:32) − (cid:18) R + σ − λ σ R (cid:19) + σ (cid:33) ++ e β ε (cid:32) − (cid:18) R (cid:19) + (cid:18) R + σ − λ σ R (cid:19) (cid:33) ++( − e β ε ) (cid:34)(cid:32) − (cid:18) R (cid:19) + λ σ (cid:33) ++ e β ε (cid:32) − (cid:18) R + λ σ − σ R (cid:19) + (cid:18) R (cid:19) (cid:33)(cid:35)(cid:41) ;if 2 σ < R ≤ σ + λ σ : f DL ( R ) = − πρβ (cid:40) e β ε (cid:32) − (cid:18) R + σ − λ σ R (cid:19) + σ (cid:33) ++ e β ε (cid:32) − σ + (cid:18) R + σ − λ σ R (cid:19) (cid:33) ++( − e β ε ) (cid:34)(cid:32) − (cid:18) R (cid:19) + λ σ (cid:33) ++ e β ε (cid:32) − (cid:18) R + λ σ − σ R (cid:19) + (cid:18) R (cid:19) (cid:33)(cid:35)(cid:41) ;if σ + λ σ < R ≤ λ σ : f DL ( R ) = − πρβ (cid:40) ( − e β ε ) (cid:34)(cid:32) − (cid:18) R (cid:19) + λ σ (cid:33) ++ e β ε (cid:32) − λ σ + (cid:18) R (cid:19) (cid:33)(cid:35)(cid:41) ;22f 2 λ σ < R : f DL ( R ) = . (34) C Linear dependence of volume terms with ρ we calculate the volume terms (0- and 1-body) of Ω in the effectivepotential of Eq. (6).The first term can be interpreted as the grand potential of a pure system of small particles at fugacity z enclosed in a volume V Ω = z β (cid:90) V dddrrr = z V β , (35)while the second is Ω = ∑ N z β (cid:90) V f i dddrrr (36) = ∑ N z β (cid:18) π σ (( e − β ε − )( λ − ) − ) (cid:19) (37) = ρ V z β × const , (38)where we have used the definition f i = e − β φ ( RRR iii −−− rrr ) −
1. From Eq. (35) and Eq. (38) we see that ( Ω + Ω ) / V is linear with respect to ρ and thus does not alter the results of the phase diagramconstruction that leads to Figure 5(a). References (1) Mezzenga, R.; Schurtenberger, P.; Burbidge, A.; Michel, M.
Nature Materials , , 729–740.(2) Evans, D.; Wennerström, H. The colloidal domain: where physics, chemistry, biology, andtechnology meet ; Wiley-VCH, 1994. 233) Pusey, P.; Van Megen, W.
Nature , , 340–342.(4) Gast, A.; Hall, C.; Russel, W. Journal of Colloid and Interface Science , , 251–267.(5) Foffi, G.; McCullagh, G.; Lawlor, A.; Zaccarelli, E.; Dawson, K.; Sciortino, F.; Tartaglia, P.;Pini, D.; Stell, G. Physical Review E , , 31407.(6) Alder, B.; Wainwright, T. Journal of Chemical Physics , , 1208.(7) Gast, A.; Russel, W. Physics Today , , 24–31.(8) Pham, K.; Puertas, A.; Bergenholtz, J.; Egelhaaf, S.; Moussaid, A.; Pusey, P.; Schofield, A.;Cates, M.; Fuchs, M.; Poon, W. Science , , 104.(9) Dawson, K.; Foffi, G.; Fuchs, M.; Götze, W.; Sciortino, F.; Sperl, M.; Tartaglia, P.; Voigt-mann, T.; Zaccarelli, E. Physical Review E , , 11401.(10) Stradner, A.; Sedgwick, H.; Cardinaux, F.; Poon, W.; Egelhaaf, S.; Schurtenberger, P. Nature , , 492–495.(11) Sciortino, F.; Mossa, S.; Zaccarelli, E.; Tartaglia, P. Physical Review Letters , , 55701.(12) Likos, C. Physics Reports , , 267.(13) Frenkel, D. Physica A , , 26–38.(14) Asakura, S.; Oosawa, F. Journal of Chemical Physics , , 1255–1256.(15) Vrij, A. Pure and Applied Chemistry , , 471.(16) Cinacchi, G.; Martínez-Ratón, Y.; Mederos, L.; Navascués, G.; Tani, A.; Velasco, E. TheJournal of Chemical Physics , , 214501.(17) Lajovic, A.; Tomšiˇc, M.; Jamnik, A. The Journal of Chemical Physics , , 104101.(18) Wolde, P.; Frenkel, D. Science , , 1975.2419) Stradner, A.; Foffi, G.; Dorsaz, N.; Thurston, G.; Schurtenberger, P. Physical Review Letters , , 198103.(20) Dorsaz, N.; Thurston, G. M.; Stradner, A.; Schurtenberger, P.; Foffi, G. Journal of PhysicalChemistry B , , 1693–1709.(21) Attard, P. Journal of Chemical Physics , , 3083.(22) Dickman, R.; Attard, P.; Simonian, V. Journal of Chemical Physics , , 205.(23) Dijkstra, M.; Brader, J.; Evans, R. Journal of Physics: Condensed Matter , , 10079–10106.(24) Hansen, J.; McDonald, I. Theory of simple liquids ; Academic Press, 2006.(25) Verlet, L.; Weis, J.
Physical Review A , , 939–952.(26) Barker, J.; Henderson, D. Journal of Chemical Physics , , 2856.(27) Caccamo, C. Physics Reports , , 1–105.(28) Gillan, M. J. Molecular Physics , , 1781–1794.(29) Wu, J.; Bratko, D.; Blanch, H.; Prausnitz, J. The Journal of Chemical Physics , ,7084.(30) Frenkel, D.; Smit, B. Understanding molecular simulation: from algorithms to applications ;Academic Pr, 2002.(31) Louis, A.; Allahyarov, E.; Löwen, H.; Roth, R.
Physical Review E , , 61407.(32) Crocker, J.; Matteo, J.; Dinsmore, A.; Yodh, A. Physical Review Letters , , 4352–4355.(33) Wasan, D.; Nikolov, A.; Henderson, D. AIChE Journal , , 550–556.2534) Xu, W.; Nikolov, A.; Wasan, D.; Gonsalves, A.; Borwankar, R. Journal of Food Science , , 183–188. 26 ρ ( r ) S σ σ S λσ λσ θ Figure 1: Scheme of the geometrical parameters that appear in Eq. (2) and Eq. (3). The solidwhite area indicates the region inaccessible to the centers of the microparticles due to the hard partof the potential φ ( r ) . The rippled area represents the structuration of the density ρ ( rrr ) due to thepresence of the macrospheres. 27 σ + λσ A tt r a c t i on ← − f D L − → R epu l s i on σ λσ R (cid:15) (cid:28) (cid:15) < (cid:15) = 0 (cid:15) > (cid:15) (cid:29) Figure 2: Qualitative force profiles for various values of ε in the dilute limit. This qualitativebehavior is observed in a quite broad region in the space of the geometrical parameters. Theconfigurations corresponding to R = σ , σ , σ + λ σ , λ σ are shown. The smaller disksrepresent the hard spheres, the white coronas the volume forbidden to the centers of the micro-spheres, and the bigger coronas are the square well/shoulder. The size of a microsphere is that ofthe red disks. The geometrical parameters used in the drawing are σ = σ = λ σ = . a) (b)(c) (d) Figure 3: Figure 3(a), Figure 3(b): Plot of the density of small particles around a big sphere (theblue disk) at large R (isolated particle case) and when another sphere (not shown) is present atsmall R in the case ε = R and small R in the case ε <
0. The anisotropy of the density determines the onset of an attractive force as inFigure 3(b), this time due to the contribution of the density both at the surfaces of radius σ (hardcore) and λ σ (outer rim of the shoulder). 29 a) (b)(c) (d) Figure 4: Density of small particles around a macrosphere for decreasing R in the case ε > σ (much stronger due to the higher density inside the well) and the “pull” at λ σ (outer rim of the well). 30 . . . . . R − − − − − V D L / ρ ξ = 5 (cid:15) = 0 . (cid:15) = 0 . (cid:15) = 0 . (cid:15) = 1 . (cid:15) = 1 . (cid:15) = 1 . (a) . . . . η . . . . . . . η ξ = 5 (cid:15) = 0 . (cid:15) = 0 . (cid:15) = 0 . (cid:15) = 1 . (cid:15) = 1 . (cid:15) = 1 . (b) Figure 5: Figure 5(a): Effective potentials obtained by numerical integration of the force in Eq. (34)using the parameters reported in the text relative to the case ξ = ε . Fig-ure 5(b): Coexistence curves obtained by first order perturbation theory using the potentials inFigure 5(a). 31 . . . . . . . R − − − − − − − f ∗ ( R ) η = 0 . ξ = 5 . MC: H = 22 , L = 16 , (cid:15) = 1 . AO: (cid:15) = 1 . PY + S: (cid:15) = 1 . HNC + S: (cid:15) = 1 . MC: H = 22 , L = 16 , (cid:15) = 0 . AO: (cid:15) = 0 . MC: [5]
PY + S: (cid:15) = 0 . HNC + S: (cid:15) = 0 . MC: H = 22 , L = 16 , (cid:15) = − . AO: (cid:15) = − . PY + S: (cid:15) = − . HNC + S: (cid:15) = − . . . . . . . . R − − − − − f ∗ ( R ) η = 0 . ξ = 5 . MC: H = 22 , L = 16 , (cid:15) = 1 . AO: (cid:15) = 1 . PY + S: (cid:15) = 1 . MC: H = 22 , L = 16 , (cid:15) = 0 . AO: (cid:15) = 0 . MC: [5]
PY + S: (cid:15) = 0 . MC: H = 22 , L = 16 , (cid:15) = − . AO: (cid:15) = − . PY + S: (cid:15) = − . . . . . . . . R − − − − − f ∗ ( R ) η = 0 . ξ = 5 . Figure 6: Force profiles obtained with the various methods in the case ξ = η = . , . , .
341 for different values of ε . In the case η = .
116 the force profiles obtainedwith the HNC closure and the superposition approximation are also presented. The legend thatholds for all the other force profiles in the article is that shown for the case η = . . . . . . . . R − − − − − f ∗ ( R ) η = 0 . ξ = 10 . . . . . . . . R − − − − − f ∗ ( R ) η = 0 . ξ = 10 . . . . . . . . R − − − − − f ∗ ( R ) η = 0 . ξ = 10 . Figure 7: Same as Figure 6, but using ξ =
10. MC data are obtained using a simulation box whose H =
30 and L =
24. 33igure 8: Plot of the difference of the densities of smaller particles around a pair of large particlesobtained via MC and PY+S for different values of ε12