15982
J. Phys. Chem. C 2007, 111, 15982-15988
Integral Equations for the Pair Structure: An Efficient Method for Studying the Potential of Mean Force in Strongly Confined Colloids† S. Amokrane,* A. Ayadim, and J. G. Malherbe Physique des Liquides et Milieux complexes, Faculte´ des Sciences et Technologie, UniVersite´ Paris Est (Cre´ teil), 61 AV. du Ge´ ne´ ral de Gaulle, 94010 Cre´ teil Cedex, France ReceiVed: May 18, 2007; In Final Form: September 4, 2007
The domain of validity of a recently proposed method for studying inhomogeneous colloidal mixtures is discussed. In this method, bridge functions derived from Rosenfeld’s hard-sphere bridge functional are used in the reference hypernetted chain closure of the Ornstein-Zernike equations, the inhomogeneity being treated by following the technique of Henderson, Abraham, and Barker. To better assess the accuracy of this method, we present new theoretical and simulation results for the wall-macroparticle potential of mean force in demanding conditions that combine the effect of strong confinement and of very short-range attractions between the colloidal particles. We show that this method satisfactorily predicts an efficient role of very short-range solvation forces, in a pseudobinary mixture confined in a slit pore that can accommodate only one big macroparticle in the direction normal to the walls.
1. Introduction The study of colloids in confined geometry is a question of both practical and fundamental importance. Indeed, these systems find applications in domains ranging from biology to nanotechnology1,2 besides providing stringent tests for the theory of complex fluids. The theoretical study of such supramolecular assemblies of nanoscopic objects is, however, a difficult task. One aspect of this difficulty is the presence of widely different length and time scales. When one is interested in those equilibrium properties that require structural information at the scale of the size of the objects forming the suspension, such as that embodied in the pair distribution function of the center of mass of the colloidal particles, the appropriate theoretical framework is the statistical mechanics of an assembly of classical objects interacting through some effective (coarsegrained) pair interaction potential φij. Once these potentials are given (for recent studies of this question, see, for example, ref 3), the thermodynamics of the system can be investigated by standard methods of liquid state theory such as Monte Carlo (MC) or molecular dynamics simulations and the “analytical” methods, such as the density functional theory (DFT) or the Ornstein-Zernike (OZ) integral equations4 for the pair structure. When studying, for example, the phase behavior of colloidal dispersions, one is tempted in a first stage to reduce the complicated actual mixture to an assembly of interacting structureless spheres whose sole characteristic is the large difference in size, the so-called binary mixture of hard spheres (HS). One may expect such a view to be appropriate, for instance, to sterically stabilized dispersions or to charge stabilized ones in the regime of high salt concentration. In both cases, the colloidal particles are deemed to behave essentially as hard-sphere-like objects. However, the study of even such a largely coarse-grained model still involves considerable technical difficulties. This is even more true when short-range tails are 5 grafted to the hard-sphere potentials φHS ij . The necessity to †
Part of the “Keith E. Gubbins Festschrift”. * Corresponding author. E-mail:
[email protected].
introduce such “residual” interactions arises for example from measurements of the potential of mean force (PMF) (see for example refs 6 and 7) or the phase diagrams8,9). This potential of mean force10 is the main quantity that we will discuss in the present paper. In recent years some progress has been made in the description of pure HS mixtures due to the development of reliable functionals in the DFT, especially that using Rosenfeld’s fundamental measure functional (FMF)11 and its variants.12,13 The question is still open when the interaction potentials φij involve non-HS contributions, due to the lack of a generic free energy functional. The alternative route using the OrnsteinZernike equations is, on the contrary, generic and should thus be used instead. However, it requires an additional relation between the interaction potential and the correlation functions, the so-called closure relation. This relation is not known exactly, even for simple fluids. In this case, the reference hypernetted chain14 (RHNC) closure with a simple parametrized hard-sphere bridge function15 proved to be very successful. Kahl et al. have shown16 that when bridge functions deduced from the HS bridge functional11 are used, this closure is quite accurate for mixtures of particles with equal (or slightly different) sizes and various interaction potentials (it will be referred to here as the RHNCFMF closure). The situation of colloids we are interested in (that is, the highly asymmetric regime) has, however, received less attention. In a previous work,17 we have shown that this method is quantitative or at worse semiquantitative for homogeneous dispersions with quite general interactions. More recently,18 we have extended its analysis to the case of inhomogeneous mixtures such as colloids in a slit pore. Its efficiency was demonstrated in situations that are beyond the reach of standard DFT. A similar treatment of a Lennard-Jones fluid in a slit pore has been done by Oettel.19 In this work, we illustrate it further by computing the PMF for a big macroparticle in a colloidal suspension under strong confinement. These conditions are close to those prevailing in the study of confined colloids in an external field1 or in biological nanopores.2 But rather than the intrinsic interest of these situations, our motiva-
10.1021/jp073834d CCC: $37.00 © 2007 American Chemical Society Published on Web 10/03/2007
Integral Equations for the Pair Structure
J. Phys. Chem. C, Vol. 111, No. 43, 2007 15983
tion here is to evaluate the efficiency of the RHNC-FMF closure in the strong confinement regime. This paper is organized as follows: In section II, the theoretical methods we used to compute the potential of mean force for inhomogeneous mixtures are briefly recalled. In section III, we assess these methods by comparing the theoretical prediction for the PMF with MC simulations for pure hardsphere potentials and in a typical situation involving non-HS interactions: a solvated macroparticle in a slit pore. We end this paper by a brief conclusion. 2. Theory Our implementation of the OZ equations with the FMF bridge function in the RHNC closure is detailed in our previous papers (ref 17 for the bulk phase, and ref 18 for inhomogeneous systems). We thus briefly recall here the treatment of inhomogeneous mixtures in a slit-like pore because the use of the RHNC-FMF closure in this context is less familiar. 2.1. Potential of Mean Force from the OZ Integral Equations. We consider here the PMF between two objects a and b that are present at infinite dilution in a mixture of n components. The objects can be two macroparticles (possibly different) or one macrosphere and a wall or a pore. This special case is obtained by letting the radius (or the radii for a pore) of object a go to infinity following the method of Henderson, Abraham, and Barker20 (HAB). For a homogeneous mixture, it is convenient to use the packing fractions ηi ) (π/6)FiDi3 with Fi the bulk number density and Di a suitably defined hard core diameter for the ith component. Equivalently, the mixture may be thought of as being in equilibrium with a reservoir that fixes the chemical potentials µi. The corresponding set of bulk densities {Fj} in the reservoir are used as the independent variables (for example, in an open pore). The (excess) potential of mean force Wab between these objects10 is most directly obtained from the pair distribution function (pdf) at infinite dilution of both objects:
gab(r;Fa,bf0,µk*a,b) ) exp[-β(φab(r) + Wab(r;Fa,bf0,µk*a,b))] (1) φab being the direct interaction between the objects (β ) 1/kBT, with T the temperature). The pdf gab is computed by taking the limit Fa,b f 0 of the OZ equation:4
γij )
∑k Fkcik X hkj
(3)
which require in practice some approximation for the bridge function bij. From eqs 1 and 3, the PMF is
Wab ) kBT(-γab + bab)
HS HS bti[{Fi(r);r}] ) β(µi,ex [{Figti(r);r}] - µi,ex ({ Fi})) +
X htj ∑j Fjc(2),HS ij
(4)
The use of these equations to study the PMF in colloids is well-known (for earlier work see, for example, refs 21 and 22, and for more recent work see refs 17 and 23-29). The difference between the various studies lies in the approximation for the bridge functions bij. As in ref 17, we used the bij obtained from Rosenfeld’s test particle consistent DFT,11 as detailed in ref 18
(5)
where µi,ex[{Fi(r);r}] ) δFex[{Fi(r)}]/δFi(r) is the excess (with respect to the ideal gas) chemical potential functional and c(2),HS is the corresponding HS direct correlation function (Fi is ij the density of particles of type i far from the inhomogeneity). Notice that in this approach, one obtains naturally the bridge functions for inhomogeneous systems (as generated by the HAB method) by inserting the inhomogeneous density profiles in the bridge functional (eq 5). Under the assumption of the universality of the HS bridge functional11,16 the bridge functions computed by inserting in eq 5 the actual correlation functions can be used in the RHNC closure for general interactions φij. One needs in practice an explicit expression of the excess HS free energy functional FHS ex . In several variants of the weighted density approximation,30 it is taken as βFHS ex [{Fi(r)}] ) ∫dx Φ[{nR(x)}], where {nR(x)} is a set of weighted densities constructed from the actual densities {Fi(r)} as nR(x) ) ΣiFi X ω(R) with weight functions ω(R) i i . As in ref 18 we used the modified version proposed in refs 12 and 13 that is more accurate than Rosenfeld’s original one.11 One has also to enforce symmetrization for unlike particles: both Rosenfeld’s ansatz31
bti ) (xibti + xtbit)/(xi + xt)
(6)
and the alternative32
bti ) (xibit + xtbti)/(xi + xt)
(7)
were used (xi is the concentration of component i). In the situation of infinite dilution of the bigger objects considered here, these are used as the natural test particle; that is, we take bit when xt ) 0 in eq 7. Finally, to apply the method to soft potentials, one needs a criterion for determining the hard-sphere radii σij to be used in eq 5. For bulk systems, one usually uses Lado’s criterion33 that minimizes the free energy with respect to the hard-sphere diameter. For mixtures, this criterion is16
(2)
where γij ) hij - cij and X designates a convolution product: [fXg](r) ) ∫dr′ f(r′) g(r′-r). hij ) gij - 1 is the total correlation function, and cij is the direct correlation function (dcf). The OZ equations are supplemented by the closure relations
gij ) exp{-βφij + hij - cij - bij}
for the inhomogeneous case. The expression of the FMF bridge functions follows ref 11 when the density profile equation obtained by minimizing the grand potential is recast in the form of a closure of the OZ eq 3:
∂bij[{gij};r]
∑ij FiFj∫dr [gij(r) - grefij (r)]
∂σkk
)0
k ) 1, 2 (8)
When taking the infinite dilution limit, only the terms with Fi * 0 survive in eq 8. For the PMF of macroparticles at infinite dilution in a one-component “solvent”, for example, one is left with a single equation corresponding to the pure fluid of small spheres. One may also use either bit or bti, depending also on the type of test particle that is used.17 The fact that in this procedure one ends up with the optimum diameter for the pure fluid of small spheres, even for strongly inhomogeneous systems, suggests that this criterion might be inappropriate in the HAB limit. An alternative one has been proposed for inhomogeneous situations by Oettel19 but this point has not yet received a clear answer. We are presently investigating this point, but it will not be detailed further here because we will consider, in what follows, only situations where there is no need, a priori, to optimize the diameters σij.
15984 J. Phys. Chem. C, Vol. 111, No. 43, 2007
Amokrane et al.
2.2. Case of Confined Colloids. We consider here the case of a pseudobinary mixture of two differently sized colloids. The suspending medium is assumed to play no explicit role because of the much smaller size of the solvent molecules. The system is confined between two parallel structureless planar hard walls. For simplicity, we take here the infinite dilution limit for the bigger colloidal particles. The goal is to investigate the way attractive forces between the big and the small colloids modify the wall-macroparticle depletion potential. Until now, the latter has indeed been studied mostly in the context of pure HS interactions (see for example refs 34 and 35). A system corresponding to a fluid near a wall or in a pore was generated from the HAB20 limit of the homogeneous OZ equations. In eqs 2, one lets the radius (or the radii for a pore) of the object a go to infinity after having taken the limit Fa f 0 in a binary mixture (indexes C, s, and w are for the big colloids, the small ones, and the wall, respectively). To distinguish between one wall and a pore, the latter will be labeled as w2 : Rw/Rsf ∞ and R ≡ RC/Rs > 1. The relevant equations are (see also ref 36) as follows: (i) the pure fluid (small colloids) OZ equation and the associated closure
γss ) Fscss X hss
(9)
gss(r) ) exp(-βφss(r) + γss - bss)
(10)
(ii) the big-small colloids OZ equation and the associated closure
γsC ) Fscss X hsC
(11)
gCs(r) ) exp(-βφCs(r) + γCs - bCs)
(12)
(iii) the fluid-wall (pore) OZE and its closure
γsw ) Fscss X hsw
(13)
gsw(z) ) exp(-βφsw(z) + γsw - bsw)
(14)
This is next inserted into the big colloid-wall OZ equation
γCw(z) ) FscCs X hsw
(15)
to obtain without iteration the wall-colloid PMF as
WCw(z) ) kBT(-γCw(z) + bCw(z))
(16)
To obtain the total colloid-wall effective interaction, one simply adds to WCw the direct interaction φCw.The colloid wall distribution function is then gCw(z) ) exp[-β(φCw + WCw)]. Notice that the PMF obtained from this method is completely free of adjustable parameters. In the case of a slit pore (two walls) one has to use in eq 14
φsw2(z) ) φsw(z) + φsw(H - z)
(17)
instead of φsw(z), where H is the pore width. In these equations, the bridge functions bss(r), bCs(r), bsw(z), and bCw(z) are obtained by inserting the appropriate density profiles Figti(r) and dcf c(2),HS with t, i ) s, w, C in eq 5. The fact that all these bridge ij functions follow from the same theory is one significant feature of this approach. The method used to solve these equations has been detailed in our previous work18 (see also ref 37). To compute WCw for a real system, one should specify the solvent-solvent, solvent-colloid, and solvent-wall potentials
Figure 1. Schematic view of the confined system. z ) 0 corresponds to the distance of closest approach to the wall of the small colloids. The zero for zCW corresponds to z ) RC - Rs.
φij (solvent possibly meaning the small colloids in a pseudobinary mixture). Notice that until specification of the interaction potentials, the equations recalled above do not distinguish between the situation we are interested in and that of a confined simple fluid (see, however, below). The different confinement situations that have been considered in our previous work18 to test this method will not be reconsidered here. Though additionnal tests would be useful, our main concern in this paper is to assess its limits of validity in a situation close to certain experimental conditions, that is for a colloidal mixture confined in a slit pore having a width comparable to the colloidal diameter. To better see how it performs in this specific situation, it is useful to recall how the various parameters of the model enter the theory. One should consider to this end the expression of the colloid-wall PMF given in eq 16. The most relevant ones are the pore width H and the size ratio DC/Ds. The former is the truly specific parameter in the present study. It determines the density profile of the confined pure fluid (without the big colloids). The second one determines the big-small colloids correlation function cCs for the bulk mixture and the colloidwall bridge function. As discussed in ref 18, the situation H >> Ds and DC/Ds >> 1 in which the pore width is much larger than the size of the smaller colloids adds little to what is known about the behavior of the PMF for a macroparticle at a single wall. As we are interested in the strong confinement regime, we need to take H/Ds ∼ 1 and hence also DC/Ds ∼ 1. In these conditions cCs in the bulk from the RHNC equation should be very accurate. This does not also constitute per se a difficult condition for bCw(z). The accuracy of the method in the present conditions depends thus largely on that of hsw(z) for the pure fluid. It would of course be interesting to perform a systematic evaluation of the RHNC and alternative integral equations for such strongly confined fluids (see for example, ref 36) but this point is disconnected with the question of the influence on the accuracy of the RHNC of the size ratio in the colloidal mixture. One should also notice that the accuracy of the integral equations generally improves when the fluid is less confined. The example shown below illustrating an excellent behavior of hsw(z) from the RHNC equation in the strong confinement regime, H/Ds ∼ 1, is therefore sufficient for the scope of this paper (note also that for bulk fluids, the quality of the RHNC method has been tested previously32 at various densities and size ratios, up to DC/Ds ) 20 with hard spheres). 3. PMF for a Strongly Confined Macroparticle 3.1. Model Systems. As we are interested in evaluating the method outlined above in the condition of strong inhomogeneity, we take a macroparticle of diameter DC ) 2Ds in a pore of width H ) 2.8Ds, for the small colloids (see Figure 1). The
Integral Equations for the Pair Structure
J. Phys. Chem. C, Vol. 111, No. 43, 2007 15985
coordinate z of the center of mass of the small colloids has access to the range [0, H]. If we assume the same hard core repulsion for the big colloids, their center of mass has access to the range [RC - Rs, H + Rs - RC] where Ri refers to the radii. With H ) 2.8Ds, the pore can accommodate only one big colloidal particle in the direction normal to the walls. We will first consider the case of pure HS potentials and the situation of a “solvated” colloid: the interaction between the big and the small colloids has an attractive Yukawa tail:
φCs(r) )
{
r < DCs ∞ (18) -Cs exp{-λCs(r - DCs)}/r r gDCs
The other interactions are modeled by the hard-sphere potential: φsw, φss ) φHS ij with
φHS ij (reDij) ) +∞
φHS ij (r>Dij) ) 0
where Dij is the HS diameter. Here DCC ≡ DC and Dss ≡ Ds (two indexes are required for the cross diameter DCs). Specifying this cross diameter is an intrinsic aspect of eq 18: we take here the simplest model DCs ) (DC + Ds)/2. A different question is the choice of the optimum HS diameter in the HS bridge function (eqs 5 and 8). In the context of Rosenfeld’s FMF, there is no need to optimize the diameters for the small-small, smallbig colloids, and the big colloid-wall bridge functions because (i) the relevant interactions are pure HS ones and (ii) the bridge functional in eq 5 is relative to additive hard spheres. We thus have DOpt Cs ) (DC + Ds)/2. It is conceivable that a more flexible functional with an optimized cross diameter might be preferable when the cross interaction involves a non-HS contribution (Cs * 0) but it is not easy to devise a simple modification of the FMF to this end (see, for example, ref 38). One important remark now concerns the assumption of a flat (structureless) wall. This is certainly more justified in the present context than for a molecular fluid confined in a slit pore whose width is comparable to the molecular diameter. In the latter case, the molecular structure of the wall and the resulting threedimensional structuring of the confined fluid should be taken into account. This is usually done by averaging in the direction parallel to the wall the total pair interactions between the wall particles and one fluid particle. A typical example is the Steele potential39 used in adsorption studies. The resulting wall-fluid particle potential Wwp(z) that depends only on the distance z to the wall is strictly speaking consistent only with a mean field treatment of the wall fluid coupling. In our case, this remark would also hold, but only for the suspending medium. Because the colloids/molecular solvent size ratio exceeds the big/small colloid one by at least 2 orders of magnitude, the suspending medium plays no role here, as specified above. Neglecting the actual structuring in the direction parallel to the walls is hence certainly more justified here. This is the reason why the model we consider is indeed appropriate to a pseudobinary mixture of colloids, even though we take DC/Ds ∼ DC/H ∼ 1 to emphasize the confinement effect. Another distinction with respect to molecular fluids and relative to the interaction range will be made below. 3.2. Hard-Sphere Potentials. We show in Figure 2 the macroparticle-wall PMF for a series of small colloids densities Fs. Comparison with simulations shows that the method is accurate even at relatively high values of Fs. To compare the theoretical PMF with the simulation data, one must establish the correspondence between the control variables: for the theory one uses Fs and in the simulation the number of small colloids and the volume of the simulation box. For the latter, one must
Figure 2. Colloid-wall potential of mean force versus reduced distance from contact for pure hard-sphere potentials. Solid lines: RHNC-FMF integral equations. Symbols: MC simulation data for small colloids reduced density FsDs3 ) 0.3 (circles), 0.57 (squares), 0.675 (triangles). Size ratio DC/Ds ) 2, pore width H ) 2.8Ds. The zero for zCW corresponds to z ) RC - Rs in Figure 1.
check that for the given pore width, the surface of the wall must be large enough to ensure in practice the infinite dilution of the macroparticles (see below). During the simulation, one computes the small colloids chemical potential in the pore µs. The corresponding bulk density Fs(µs) is determined here (for HS colloids) from the Carnahan-Starling equation of state4 (in a more general situation, one has to perform a simulation for the pure bulk fluid to obtain µs). A second problem is the absolute magnitude of the PMF which is difficult to obtain by simulation. For the present purpose we simply shifted the MC data so as to match the theoretical value of the PMF at the pore center. In Figure 2, we observe the expected deepening of the depletion well, and the usual structuring due to the small-small and big-small colloids correlations. There is also a progressive buildup of a secondary minimum at the middle of the pore. Besides this interesting feature that will be discussed below, the main observation is the expected effect of depletion forces already discussed for bulk colloids and colloids at a wall (see references, for example, in refs 34 and 35). We may now ask a question that is suggested by the shape of the PMF in the central region: is it possible to deepen further this secondary minimum (with respect to the height of the maximum) without changing the density of the small colloids? 3.3. Solvated Macroparticle. A simple way18 to act on the PMF is to “solvate” the macroparticle by grafting an attractive Yukawa tail to φCs(r) (Cs * 0 in eq 18). To mimic the situation of hard-sphere-like colloids in a narrow pore, we must take a very short interaction range, thus λCs-1
