A dynamic density functional theory for particles in a flowing solvent

J. Chem. Phys. 127, 244906 (2007); https://doi.org/10.1063/1.2806094 127, 244906 © 2007 American Institute of Physics. A dynamic density functional theory for particles in a flowing solvent Cite as: J. Chem. Phys. 127, 244906 (2007); https://doi.org/10.1063/1.2806094 Submitted: 10 September 2007 . Accepted: 12 October 2007 . Published Online: 28 December 2007 Markus Rauscher, Alvaro Domínguez, Matthias Krüger, and Florencia Penna ARTICLES YOU MAY BE INTERESTED IN Dynamical density functional theory and its application to spinodal decomposition The Journal of Chemical Physics 121, 4246 (2004); https://doi.org/10.1063/1.1778374 Driven colloidal suspensions in confinement and density functional theory: Microstructure and wall-slip The Journal of Chemical Physics 140, 094701 (2014); https://doi.org/10.1063/1.4866450 Dynamic density functional theory of fluids The Journal of Chemical Physics 110, 8032 (1999); https://doi.org/10.1063/1.478705 A dynamic density functional theory for particles in a flowing solvent Markus Rauschera兲 Max-Planck-Institut für Metallforschung, Heisenbergstr. 3, D-70569 Stuttgart, Germany and Institut für Theoretische und Angewandte Physik, Universität Stuttgart, Pfaffenwaldring 57, D-70569 Stuttgart, Germany Alvaro Domínguez Física Teórica, Universidad de Sevilla, Apdo. 1065, E-41080 Sevilla, Spain Matthias Krüger Max-Planck-Institut für Metallforschung, Heisenbergstr. 3, D-70569 Stuttgart, Germany and Institut für Theoretische und Angewandte Physik, Universität Stuttgart, Pfaffenwaldring 57, D-70569 Stuttgart, Germany Florencia Penna Universidad Autónoma de Madrid, E-28049 Madrid, Spain 共Received 10 September 2007; accepted 12 October 2007; published online 28 December 2007兲 We present a dynamic density functional theory 共dDFT兲which takes into account the advection of the particles by a flowing solvent. For potential flows, we can use the same closure as in the absence of solvent flow. The structure of the resulting advected dDFT suggests that it could be used for nonpotential flows as well. We apply this dDFT to Brownian particles 共e.g., polymer coils兲in a solvent flowing around a spherical obstacle 共e.g., a colloid兲and compare the results with direct simulations of the underlying Brownian dynamics. Although numerical limitations do not allow for an accurate quantitative check of the advected dDFT both show the same qualitative features. In contrast to previous works which neglected the deformation of the flow by the obstacle, we find that the bow wave in the density distribution of particles in front of the obstacle as well as the wake behind it are reduced dramatically. As a consequence, the friction force exerted by the 共polymer兲 particles on the colloid can be reduced drastically. © 2007 American Institute of Physics. 关DOI: 10.1063/1.2806094兴 I. INTRODUCTION The generalization of classical density functional theory 共DFT兲to nonequilibrium states has become a valuable tool to study the dynamics of directly interacting Brownian particles. This dynamic DFT 共dDFT兲for the ensemble averaged density was proposed in Refs. 1and 2and studied in a meanfield model in Ref. 3. More recently, dDFT has been extended to mixtures4,5and anisotropic particles.6Hydrodynamics is known to play a crucial role in the dynamics of suspensions, but up to now, hydrodynamic interactions have been treated only in a mean-field manner7and the solvent has been assumed to be at rest. dDFT has also been used to investigate the distribution of solute particles around a strongly repulsive potential moving through the suspension.8 The intention was to model a colloidal particle moving through a polymer solution. A similar model has been used to study the depletion interaction between two colloidal particles moving through a polymer solution.9,10 In all these studies, the hydrodynamic flow of the solvent around the colloid was neglected and the solvent effectively passed through the colloid. A real colloid would displace the solvent as it moves, as shown for the case of a small and a large colloid in Fig. 1. For a spherical colloidal particle of radius R dragged through an unbounded incompressible viscous Newtonian solvent with velocity cat low Reynolds number, the flow field u共r兲共in a frame of reference comoving with the colloid兲is given by the solution of the Stokes equation,11 u共r兲=3R 4r 冉 1+ R2 3r2 冊 c+3R 4r3r共r·c兲 冉 1−R2 r2 冊 −c.共1兲 For large distances from the colloid, rⰇR, the flow field is well approximated by u共r兲=−c. For solute particles 共e.g., polymers or other colloids兲which only feel this far field, the model presented in Refs. 8–10 is a reasonable approximation. This is the case for large solute particles with a radius dⰇR. Their centers can approach the dragged colloidal particle only up to a distance D=R+d⬇d, see Fig. 1共a兲, and thus will feel a flow field u共r兲⬇−c, if one neglects the additional effect of the solute particles on the solvent flow; in other words, the hydrodynamic interaction between the solute particle and the colloid is neglected. This is a reasonable approximation for polymer coils but certainly a bad one for solid solute particles. However, small solute particles of radius dⰆRcan get much closer to the colloid and feel the distortion of the solvent velocity field, as illustrated in Fig. 1共b兲. The solute particles will be deviated from the colloid by the flow field and this will reduce significantly the bow-wave effect in front of the colloid presented in Ref. 8and the strength of the nonequilibrium depletion force discussed in Refs. 9and 10. Extremely small solute particles would not a兲Electronic mail: [email protected]. THE JOURNAL OF CHEMICAL PHYSICS 127, 244906 共2007兲 0021-9606/2007/127共24兲/244906/8/$23.00 © 2007 American Institute of Physics127, 244906-1 show this effect at all since they would behave like solvent molecules. In this limit, however, the basis of the theory discussed here, i.e, the description of the solute particles as overdamped Brownian particles, is no longer valid because it is based on a separation of the length and time scales associated with the solvent molecules and the solute particles. In this paper, we present a generalization of the dDFT derived in Refs. 1and 2to the case of Brownian solvent particles advected by a flow, thereby incorporating some aspects of the hydrodynamics of the solvent into the theory. However, we do not model hydrodynamic interactions between the solute particles as well as the back reaction of the solute particles on the flow field, e.g., by a concentration dependent viscosity, or by a reduced mobility of the solute particles in the vicinity of the colloid. However, the latter can be included in a straightforward manner as we discuss in the conclusions in Sec. IV. In Sec. II, we derive the advected dynamic density functional theory using the method described in Refs. 12 and 13. In Sec. III, we discuss two sample cases, namely, ideal solute particles and Gaussian solute particles which stress the importance of taking into account the solvent flow. II. ADVECTED dDFT We start with the Langevin equation of an ensemble of N advected interacting Brownian particles confined to a finite volume Vin the overdamped limit, dri dt =u共ri兲−⌫ⵜi 冋 U共ri兲+兺 j=1 N V共兩ri−rj兩兲 册 + ␩ i共t兲,共2兲 with the pair interaction potential between the particles V共r兲 and an external potential U共r兲. Both the flow field u共r兲and the external potential U共r兲can depend on time. For example, a time dependent potential has been used in Ref. 14 to model the oscillating cavities containing soft particles. For clarity of notation we will not make this dependence explicit in the equations. The flow field is not necessarily divergence-free 共i.e., the solvent can be compressible兲.ⵜidenotes the gradient with respect to ri. We approximate the noise generated by the thermal motion of the solvent particles by a Wiener process, 具 ␩ i共t兲典 =0,共3兲 具 ␩ i ␣ 共s兲 ␩ j ␤ 共t兲典 =2T⌫ ␦ ij ␦ ␣␤ ␦ 共t−s兲,共4兲 with the temperature Tmeasured in units of energy 共setting kB=1兲and the mobility coefficient ⌫⬎0. The mobility coefficient has to appear in the correlation in Eq. 共4兲in order to fulfill the fluctuation-dissipation theorem and to get the correct equilibrium distribution for u共r兲=0. The boundaries of Vare impermeable for the particles or periodic 共or a mixture of both兲and, therefore, the number of particles is conserved. The Fokker-Planck equation corresponding to the Langevin equation 共2兲gives the time evolution of the probability density W共r1,...,rN,t兲for finding the particles at time tat the positions r1,...,rN,15,16 ⳵ W ⳵ t=−兺 i=1 N ⵜi· 再 ⌫ 冋 u共ri兲 ⌫−ⵜiU共ri兲 −兺 j=1 N ⵜiV共兩ri−rj兩兲 −Tⵜi 册 W 冎 .共5兲 For a potential flow, the velocity field can be written as the gradient of a scalar field, u共r兲=−⌫ⵜ⌽共r兲,17 such that the external potential and the effect of the flow field can be combined into a modified external potential U*共r兲=U共r兲+⌽共r兲. If U*is time independent, one can find a stationary probability density Weq *共r1,...,rN兲which fulfills the detailed balance condition for Eq. 共5兲, i.e., the term in curly brackets is zero for each i=1,...,N, 冋 ⵜiU*共ri兲+兺 j=1 N ⵜiV共兩ri−rj兩兲 +Tⵜi 册 Weq *=0. 共6兲 The solution is Weq *共r1, ...,rN兲=1 Z*e−共1/T兲兺i=1 N 关 U*共ri兲+兺j=1 NV共兩ri−rj兩兲 兴 ,共7兲 normalized with the sum of states Z*such that 冕冕 VN d3r1, ...,d3rNWeq *共r1, ...,rN兲=1. 共8兲 For such a situation, the whole apparatus of equilibrium statistical mechanics can be used in order to calculate expectation values and correlations in a stationary nonequilibrium situation. However, this is restricted to cases when the detailed balance condition holds, which implies necessarily a potential flow. Even in the Stokes flow 共1兲or in a simple shear flow 共e.g., in Couette or Poiseuille flow兲this is not true because ⵜ⫻u⫽0. For flows with a finite vorticity, there is no detailed balance in a strict sense, see 关Ref. 15, Eq. 共5.3.4c兲兴. From Eq. 共5兲, we can calculate the time evolution of the noise averaged particle density ␳ 共r,t兲, namely, the expectation value of the density operator ␳ ˆ共r,t兲=兺i=1 N ␦ 共r−ri共t兲兲, FIG. 1. 共Color兲Cross section of the flow field given in Eq. 共1兲in a plane parallel to the direction of motion around 共a兲a small spherical colloid and 共b兲a big one 共full circles, radius R兲. The dashed circles of radius Dmark the points of closest approach of solute’s centers 共point in the center of open circles兲to the colloid. The solute’s diameter in its mutual interaction is ␴ and, for nonadditive mixtures, not necessarily equal to its diameter 2d =2共D−R兲共indicated by the dotted circle兲in the interaction with the colloid. The component of the flow field normal to the dashed circle is larger for the small colloid 共a兲than for the large colloid 共b兲. ␴ and dare the same in both figures. 244906-2 Rauscher et al. J. Chem. Phys. 127, 244906 共2007兲 ⳵␳ ⳵ t+ⵜ·共 ␳ u兲 =ⵜ·⌫ 冋 ␳ ⵜU+Tⵜ ␳ +ⵜ 冕 V d3r⬘V共兩r−r⬘兩兲 ␳ 共2兲 ⫻共r,r⬘,t兲 册 ,共9兲 with the mean density ␳ 共r,t兲=N 冕冕 VN−1 d3r2, ...,d3rNW共r,r2, ...,rN,t兲,共10兲 and the nonequilibrium density-density correlation function ␳ 共2兲共r,r⬘,t兲=N共N−1兲 ⫻ 冕冕 VN−2 d3r3, ...,d3rNW共r,r⬘,r3, ...,rN,t兲. 共11兲 Equation 共9兲is the starting point of a hierarchy of Nevolution equations which connect the time derivative of the n-point density correlation function to the 共n+1兲-point density correlation function, similar to the BBGKY hierarchy for deterministic systems with inertia or the BGY hierarchy for equilibrium correlation functions. In order to find a closed equation for the time evolution of ␳ 共r,t兲, we approximate the interaction term in Eq. 共9兲by its value in an equilibrium system with the same interaction potential V共r兲. Let us first restrict our considerations to the case that detailed balance holds 共so that, in particular, u =−⌫ⵜ⌽兲. We modify our system of Brownian particles by applying an external potential ⌿共r兲so as to create a system whose equilibrium density distribution is ␳ eq ⌿共r兲= ␳ 共r,t兲. This new potential ⌿共r兲depends on ␳ 共r,t兲and will be different for each tas long as ␳ 共r,t兲is not stationary. The FokkerPlanck equation for the modified system will be Eq. 共5兲but with U共r兲being replaced by U共r兲+⌿共r兲. The equilibrium probability density Weq ⌿of the modified system is given by Eq. 共7兲but with U*共r兲being replaced by U*共r兲+⌿共r兲.Ifwe integrate the detailed balance condition 共6兲for Weq ⌿over N −1 positions, we get u ␳ eq ⌿=⌫ 冋 ␳ eq ⌿ⵜ共U+⌿兲+Tⵜ ␳ eq ⌿ +ⵜ 冕 V d3r⬘V共兩r−r⬘兩兲 ␳ ⌿ 共2兲共r,r⬘兲 册 ,共12兲 with the equilibrium pair correlation function ␳ ⌿ 共2兲共r,r⬘兲for the modified system in the external potential ⌿. From equilibrium density functional theory, one knows that the equilibrium density distribution in the grand canonical ensemble is the minimum of the grand canonical functional ⍀关 ␳ 兴=Fex关 ␳ 兴+ 冕 V d3r兵T ␳ 关ln共 ␳ ⌳3兲−1兴+共Uext − ␮ 兲 ␳ 其, 共13兲 with the thermal wavelength ⌳, the chemical potential ␮ , and the sum of all external potentials Uext=U+⌽+⌿. The excess free energy Fex关 ␳ 兴summarizes the effect of the particle interactions and it is not known exactly in general. We take the gradient of the Euler-Lagrange equation following from the functional 共13兲: since in thermal equilibrium the chemical potential is constant across the whole system, we get u ⌫=ⵜ共U+⌿兲+T ␳ eq ⌿ⵜ ␳ eq ⌿+ⵜ冏 ␦ Fex关 ␳ 兴 ␦ ␳ 冏 ␳ eq ⌿ .共14兲 If we compare Eq. 共12兲with Eq. 共14兲we can see that ⵜ 冕 V d3r⬘V共兩r−r⬘兩兲 ␳ ⌿ 共2兲共r,r⬘兲= ␳ eq ⌿ⵜ冏 ␦ Fex关 ␳ 兴 ␦ ␳ 冏 ␳ eq ⌿ .共15兲 Note that the right hand side does not depend on the velocity potential ⌽while the dependence on ⌿enters only through ␳ eq ⌿. We will use Eq. 共15兲as a closure to the hierarchy of equations starting with Eq. 共9兲: Hereby we assume that the density correlations at time tin the nonequilibrium system with mean density ␳ 共r,t兲are the same as those in an equilibrium system with the additional potential ⌿and with equilibrium mean density ␳ eq ⌿共r兲= ␳ 共r,t兲. We then get ⳵␳ ⳵ t+ⵜ·共 ␳ u兲=ⵜ· 冉 ⌫ ␳ ⵜ ␦ F关 ␳ 兴 ␦ ␳ 冊 ,共16兲 with the free energy functional F关 ␳ 兴=Fex关 ␳ 兴+ 冕 V d3r兵T ␳ 关ln共 ␳ ⌳3兲−1兴+ ␳ U其.共17兲 In thermodynamic equilibrium for u=0 and time independent U, the equilibrium density distribution given by ␮ =兩 ␦ F关 ␳ 兴/ ␦ ␳ 兩 ␳ eq is a stationary solution of Eq. 共16兲.An H-theorem, ⳵ ⳵ t 冕 V d3rF关 ␳ 兴=− 冕 V d3r⌫ ␳ 冉 ⵜ ␦ F关 ␳ 兴 ␦ ␳ 冊 2 艋0, 共18兲 guarantees that the time evolution actually converges to the equilibrium distribution. 关The dynamics in Eq. 共2兲together with the boundary conditions taken for Vimply that the particle current through the system boundaries is zero and, therefore, the surface terms from partial integration vanish.兴 The final chemical potential is then determined by the conserved number of particles in the system. As discussed above, a system in a potential flow corresponds to an equilibrium system with a modified external potential U*=U +⌽. Equation 共16兲can then be written in the form ⳵␳ ⳵ t=ⵜ 冉 ⌫ ␳ ⵜ ␦ F*关 ␳ 兴 ␦ ␳ 冊 ,共19兲 with the modified free energy functional F*关 ␳ 兴=F关 ␳ 兴 +兰Vd3r⌽ ␳ . Thus, we have an H-theorem for F*关 ␳ 兴instead of 244906-3 Dynamic density functional theory J. Chem. Phys. 127, 244906 共2007兲 F关 ␳ 兴and the “equilibrium state” is determined by ␦ F*关 ␳ 兴/ ␦ ␳ = ␮ *. The right hand side of Eq. 共15兲is completely independent of the flow field and one could be tempted to use it as a closure to Eq. 共9兲for the most general case 共i.e., nonpotential flows兲. However, this would mean approximating density correlations in a driven nonequilibrium system where detailed balance cannot be achieved by thermal equilibrium correlations. While we could argue that close to equilibrium, Eq. 共15兲may be a reasonable approximation if detailed balance still holds, there is no such argument for the most general case that detailed balance is violated. The study in Ref. 8 addresses a situation where the approximation is found to be good although there is no detailed balance since a net particle current is driven through the system. In the next section we consider some examples with the purpose of assessing 共i兲 the effect of a more realistic Stokes flow as discussed in the Introduction, and 共ii兲the validity of the approximate Eq. 共16兲 for a nonpotential flow. III. EXAMPLES As an example, we study a solution of polymers 共radius d, density ␳ 0兲in an incompressible Newtonian solvent flowing around a spherical colloidal particle 共radius R兲in a stationary situation. We model the polymer coils as pointlike particles from the point of view of the solvent, but with a finite interaction range ␴ concerning other polymer coils and D=d+Rconcerning the colloidal particle, the interaction with the latter being of hard-wall type 共see Fig. 1兲. The velocity field of the solvent is given by the Stokes flow 关Eq. 共1兲兴and we choose c=ce ˆz. Measuring lengths in terms of D, the dimensionless parameters determining the system are the Péclet number c*=cD/共⌫T兲, the colloid radius R* =R/D, and the polymer’s mutual interaction range ␴ /D. A. Ideal polymers For ideal solute particles in an incompressible solvent with U=0, the stationary condition ⳵␳ / ⳵ t=0 from Eq. 共16兲 reads u·ⵜ ␳ =⌫T⌬ ␳ ,共20兲 where ⌫Tis the diffusion constant of the solute particles. The hard interaction with the colloid is written as a boundary condition for the current density of the solute particles j =u ␳ −⌫Tⵜ ␳ at r/D=1, 兩共e ˆr·j兲兩r/D=1 =0. 共21兲 We expand the density field ␳ 共r兲in spherical harmonics up to an order Nand obtain a system of N+1 ordinary differential equations for the 兩r兩-dependent coefficients which we solve numerically with AUTO 2000.18 AUTO 2000 is a software which solves autonomous boundary value problems for systems of ordinary differential equations by continuation, i.e., by starting from a known solution for a specific set of problem parameters 共for c=0, we have ␳ = ␳ 0兲and changing paFIG. 2. Contour plots of the density of ideal polymers for a flow velocity c*=10. The white circle at the origin is the colloidal particle with radius R, the black circle is the annulus of thickness d, and outer diameter Dwhich is unaccessible to the polymer centers due to the hard-wall interaction, see Fig. 1.共a兲corresponds to a uniform flow u共r兲=−ce ˆz共i.e., to R*=0兲.The maximum density in front of the colloid is ␳ 共r兲/ ␳ 0=6.31. 共b兲corresponds to R*=0.9. The bow-wave effect is reduced drastically. The maximum density in front of the colloid is ␳ 共r兲/ ␳ 0=1.05. FIG. 3. Density of ideal polymers at the point x=y=0, z/D=1, i.e., right in front of the forbidden zone around the colloid. 共a兲shows ␳ 共0,0,D兲as a function of the colloid size for c*=1 and c*=10. 共b兲shows ␳ 共0,0,D兲as a function of c*for different values of R*. 244906-4 Rauscher et al. J. Chem. Phys. 127, 244906 共2007兲 rameters 共in our case c兲continuously until the desired value is reached. As demonstrated for R*=0 and R*=0.9 in Figs. 2共a兲and 2共b兲, respectively, the bow-wave effect is large when the colloid is small compared to the polymers and does not distort the flow too much, reaching its maximum as R*→0. This is the case investigated in Ref. 8. When the colloid is large compared to the polymers, the bow-wave effect is small. In the limit d/R→0, the effect vanishes completely since the polymers behave like solvent molecules. Figures 3共a兲and 3共b兲show the density right in front of the colloid as a function of R*and of c*, respectively. The density of ideal solute particles scales almost linearly with the velocity c*. B. Gaussian polymers Here, we address the case of interacting polymers. We consider the same polymer-polymer interaction potential studied in Ref. 8, namely, V共r兲=Texp关−共r/ ␴ 兲2兴.共22兲 The interaction of the polymers with the colloidal particle is modeled as an external potential of the form U共r兲=10Texp关−共r/a兲6兴.共23兲 This potential rises steeply up to 10T, thus resembling a hard wall. The length amust be related to the radius of the forbidden zone Daround the colloidal particle. We conventionally set the value of aby the condition U共D兲=2T, giving a ⬇0.924D. We take ␴ =2d共i.e., additive mixture of polymers and colloidal particle兲and R=1.7 ␴ , leading to R*⬇0.77, ␴ /D⬇0.46. Finally, we also considered the choice R*=0, ␴ /D⬇0.46, which represents a hard particle that does not distort the uniform flow 关R=0 in Eq. 共1兲兴in a nonadditive mixture 共 ␴ ⫽2d兲: this was the model addressed in Ref. 8. We ran Brownian dynamics 共BD兲simulations of this system for two values of the flow velocity corresponding to the polymer Péclet numbers 共 ␴ /D兲c*=1 and 10 studied in Ref. 8共i.e., c*⬇2.2 and 22兲. We considered a colloidal particle at the center of a box of dimensions Lx=Ly=12 ␴ and Lz=24 ␴ with periodic boundary conditions. The box contained N=3456 polymers, corresponding to a mean polymer number density ␳ 0 ␴ 3=1. We took a time step of 0.003 ␴ 2⌫/T for the discretized Langevin dynamics. The system was allowed to relax for 105time steps, after which collection of data was carried out during 106time steps. Even though the simulated system is finite, we used the analytically known flow field around a sphere in an infinite medium 关Eq. 共1兲兴. The error due to the truncation of this flow by the boundary of the simulation box is the largest 共about 20%兲at the midplane of the colloid 共z=0兲. This introduces effectively a discontinuity in the flow velocity field at the boundary which we discuss later. We also solved numerically the dDFT in the random phase approximation 共a mean-field model兲, i.e., with3,8,9 Fex关 ␳ 兴=1 2 冕冕 V2 d3rd3r⬘V共兩r−r⬘兩兲 ␳ 共r兲 ␳ 共r⬘兲.共24兲 The time evolution given in Eq. 共16兲of an initially homogeneous density was solved in cylindrical coordinates on a grid spanning the domain −120⬍z/ ␴ ⬍24, 0艋r⬜/ ␴ ⬍60, where r⬜=冑x2+y2. The grid constant was 0.0125 ␴ near the colloid, i.e., for 兩z兩,r⬜⬍6 ␴ , and 0.1 ␴ in the rest of the domain. For details on the numerical procedure see Ref. 8. The boundary condition at the domain border was ␳ = ␳ 0, also in this case we used the flow field given in Eq. 共1兲. The error introduced here is smaller than that in the BD simulations since the integration domain for the dDFT is larger than the BD simulation box. Figure 4presents the density field ␳ ¯ 共z兲, spatially averaged over thin disks of radius ␴ and thickness 2⌬z=0.05 ␴ centered at the zaxis, i.e., ␳ ¯ 共z兲ª 1 ␴ 2⌬z 冕 z−⌬z z+⌬z dz⬘ 冕 0 ␴ dr⬜r⬜ ␳ 共r⬜,z⬘兲.共25兲 The results of both the BD simulations and the dDFT illustrate the dramatic effect of advection by the Stokes flow 共1兲. In particular, at the higher velocity c*=22 and R*=0 共uniform flow兲there is a marked accumulation of polymers in front of the particle and a strong depletion behind it which are hardly observable for R*=0.77. In general, the effect of the Stokes flow is to weaken the influence of the colloid on FIG. 4. Plots of ␳ ¯ 共z兲defined in Eq. 共25兲as provided by the numerical solution of dDFT 共lines兲and as measured in BD simulations 共symbols兲.共a兲 corresponds to a velocity c*=2.2 and 共b兲to c*=22. In each plot the results for both uniform flow 共R*=0兲and Stokes flow 共R*=0.77兲are presented. 244906-5 Dynamic density functional theory J. Chem. Phys. 127, 244906 共2007兲 the density profile, as the polymers tend to be advected by the stream and to travel around the particle. Actually, for c*=2.2 共and smaller兲the deformation by the Stokes flow is so tiny that the dDFT profile ␳ ¯ 共z兲in Fig. 4共a兲is indistinguishable from the equilibrium profile 共i.e., c*=0兲. We notice a discrepancy between the density profiles measured in the BD simulation and those calculated numerically in the dDFT. We attribute this to two finite-size effects in the simulation which have been confirmed by performing BD simulations in a smaller box at the same polymer number density 共Lx=Ly=8 ␴ ,Lz=16 ␴ ,N=1024兲. First, the numerical solution of the dDFT for a Stokes flow with c*=22 exhibits a long 共⬇15 ␴ 兲tail of slight polymer depletion 共 ␳ ⬇0.98 ␳ 0兲 behind the colloidal particle. The tail is longer than the length of the BD simulation box in the zdirection and, as a consequence of the periodic boundary conditions, the inflowing density far ahead of the particle is smaller than ␳ 0. This screening effect is very noticeable when the flow is approximated as uniform because the depletion of polymers behind the colloidal particle is very large 关see Fig. 4共b兲兴. However, in the case of the Stokes flow, this effect seems to be less important compared to the second effect: the discontinuity of the normal component of the flow field at the lateral boundaries of the simulation box leads to a nonvanishing divergence of the flow 关we remind that the flow described by Eq. 共1兲has ⵜ·u=0兴. Figure 5represents the density profile averaged over thin cylindrical shells of height ␴ and radial thickness ⌬r=0.05 ␴ coaxial with the zaxis, i.e., ␳ ˆ共r⬜,zc兲=1 ␴ 冉 r⬜+⌬r 2 冊 ⌬r ⫻ 冕 zc−共 ␴ /2兲 zc+ ␴ /2 dz⬘ 冕 r⬜ r⬜+⌬r dr⬜ ⬘r⬜ ⬘ ␳ 共r⬜ ⬘,z⬘兲.共26兲 The periodic boundary conditions imply that ⵜ·u⬍0 effectively at the side boundaries located upstream, where therefore the density is enhanced: Even though we expect the density to decrease toward ␳ 0as the radial distance r⬜to the colloid increases, we find instead an increase of the density at the boundary of the simulation box. At the side boundaries located downstream, on the other hand, ⵜ·u⬎0 and the region near the boundary of the simulation box becomes depleted of polymers. As expected this effect is enhanced for reduced box size 共see the inset in Fig. 5兲. The comparison of the BD results with the dDFT calculation in Fig. 4indicates that the overall consequence of these effects is a density enhancement near the colloidal particle. In view of these important finite-size effects, we cannot quantify the validity of the approximation given by Eq. 共16兲for a realistic flow. When there is only uniform flow 共R*=0兲, however, the finite-size effects are much less pronounced and we find a good agreement between BD simulations and dDFT calculations, in concordance with Ref. 8. C. Drag force on colloids We have also measured the force in the zdirection exerted by the polymers on the particle. This force is additional to the Stokes drag force FStokes exerted by the flowing solvent. If ⌫in Eq. 共2兲is assumed to be given by the StokesEinstein relation for spherical polymers of diameter ␴ , the Stokes drag for a colloid of radius Rin the same solvent is given by FStokes=2c*共R/ ␴ 兲共T/D兲. For example, for the colloid radius R*=0.77 considered in the BD simulations this gives FStokes=3.4c*T/D. In the BD simulations, the force excerted on the colloid by the polymers can be measured directly. In the case of ideal particles discussed in Sec. III A, we use the ideal gas law p=T ␳ in order to calculate the local pressure on the colloid surface. Integrating the local pressure over the surface yields the force on the colloid. Table Icollects the mean force for different types of flow 共R*=0 for uniform flow and R*=0.77 for Stokes flow兲and velocities. The results confirm the necessity to take into account the solvent flow: the mean force in the case of Stokes flow is markedly smaller 共at most of the order of FStokes兲than in the case of a homogeneous flow and the dependence on c* is milder. This can be understood in terms of the reduction of the bow-wave effect in the density profile around the particle by the Stokes flow in the solvent which advects the polymers. The forces in the BD simulations are of the same order of magnitude as in the ideal case, but in the case of the Stokes flow they have a weaker dependence on c*. Because the sizes involved are of the order of the microFIG. 5. The plot represents ␳ ˆ共r⬜,zc兲, see Eq. 共26兲, in a Stokes flow 共R* =0.77兲with c*=22 upstream 共at zc=1.5 ␴ , triangles up兲and downstream 共at zc=−1.5 ␴ triangles down兲of the colloid, as measured in BD simulations. The inset shows the results for a smaller simulation box. TABLE I. Mean force exerted by the polymers measured in the BD simulations for different types of flow 共R*=0 for the uniform flow, and R* =0.77 for the Stokes flow兲, compared to the force exerted by ideal polymers and to the Stokes friction of the colloid. The forces are given in units of T/D. Type of flow c*兩具Fz典兩 Ideal FzFStokes Uniform 2.2 42.5 41.9 Stokes 2.2 6.38 3.14 7.48 Uniform 22 215 290 Stokes 22 10.3 18.0 74.8 244906-6 Rauscher et al. J. Chem. Phys. 127, 244906 共2007兲 scopic length D, the variance of the force measured in the Brownian dynamics is relatively large. However, we find that it is not affected by the flow type and velocity and it coincides with the variance of Fzwe have measured in the equilibrium state 共c*=0兲. For comparison, Table II collects the force measured in the BD simulation in the smaller box. The increased force is consistent with the density enhancement near the particle caused by the finite-size effects. We also attribute the increase of the fluctuations to these effects. IV. CONCLUSIONS We have proposed a dDFT 关Eq. 共16兲兴for interacting Brownian particles in a flowing solvent under the assumption that detailed balance holds 共which requires, in particular, a curl-free flow兲. We get the same equation as already derived in Ref. 1, but with the partial time derivative replaced by the total 共material兲time derivative. The whole effect of the flow field can be summarized into a modified external potential, allowing application of the whole machinery of equilibrium statistical mechanics. Thus, we are able to find an H-theorem for a modified free energy. In this paper, we include the displacement of the solvent by the colloid, but the hydrodynamic interactions, both, among the solute particles and between these and the colloid, were not taken into account. The treatment of the first kind of interactions is a highly nontrivial and still open problem. The equivalent of Eq. 共9兲including hydrodynamic interactions in a pairwise approximation has been derived in Refs. 19 and 20 but it contains additional terms, one involving pair correlations in a form such that the closure relation in Eq. 共15兲 cannot be used directly, and one containing the three-point density-density correlation function. In the context of sedimentation, this kind of hydrodynamic interactions has been included in terms of a density-dependent mobility ⌫共 ␳ 兲.7The hydrodynamic interactions of the solute particles with the colloid, however, can be included in a straightforward manner by replacing the mobility ⌫by a space dependent and symmetric mobility tensor ⌫共ri兲in Eqs. 共2兲and 共4兲. Thereby the noise becomes multiplicative and the appropriate calculus has to be considered such that it leads to the FokkerPlanck equation 共5兲with ⌫replaced by ⌫共ri兲. Then the equilibrium distribution is not changed. For spherical particles in the vicinity of planar walls, the mobility tensor can be calculated in the limit of large distances.21 This result has been extended to surfaces with a partial slip boundary condition in Ref. 22. Due to the translational symmetry of the system, ⌫ is diagonal. While the mobility perpendicular to the wall increases with distance, the distance dependence of the mobility parallel to the wall depends on the slip condition. For no-slip it increases while for total slip it decreases with the distance to the wall. The hydrodynamic interaction between two spheres has been calculated, e.g., in the Rotne-Prager approximation.23,24 The derivation of Eq. 共16兲essentially carries through with the only exception that the quotient u共ri兲/⌫in Eq. 共5兲has to be replaced by a field u ˜ 共ri兲with ⌫共ri兲·u ˜ 共ri兲=u共ri兲. In order to absorb the flow field into a modified external potential, u ˜ 关and not only u共ri兲兴 has to be curl-free with u ˜ =−ⵜ⌽. Instead of Eq. 共16兲we then get ⳵␳ ⳵ t+ⵜ·共 ␳ u兲=ⵜ· 冉 ␳ ⌫共r兲·ⵜ ␦ F关 ␳ 兴 ␦ ␳ 冊 .共27兲 Recently, dDFT has been used to describe the dynamics of mixtures4,5as well as anisotropic particles with orientational degrees of freedom.6The derivation presented in this paper can be generalized to mixtures in a straightforward manner. For anisotropic particles, the coupling of the vorticity of u共r;t兲to the orientational degrees of freedom has to be taken into account. Nevertheless, we see no obvious reason why it should not be possible to obtain a dDFT also in this case. In Ref. 8, the polymer distribution was studied in a polymer solution flowing uniformly through a spherical particle which is hard only for the polymers. In spite of the violation of detailed balance by the boundary conditions, the comparison between simulations and the numerical solution of the proposed dDFT was good. In this paper, we have considered the more realistic case of a Stokes flow 共1兲around the particle. The exact solution of the ideal case 共no polymerpolymer interaction兲, the numerical solution of the interacting case as well as the corresponding Brownian dynamics simulations evidence all the dramatic effect by advection on the properties of the stationary solution. We conclude that the approximation of uniform flow, as employed in Refs. 8–10, is quantitatively bad. We have found discrepancies in the density distribution of polymers as measured in the simulations and as computed numerically in the framework of the dDFT. However, the discrepancies could be rationalized in terms of finite-size effects in the simulations due to the slow decay of the Stokes flow far from the obstacle. Thus, although a quantitative check of the validity of the approximations leading to the dDFT in Eq. 共16兲and of its validity for nonpotential flows was not possible the results are encouraging. ACKNOWLEDGMENTS The authors thank S. Dietrich for financial support and fruitful discussions. A.D. acknowledges financial support from the Junta de Andalucía 共Spain兲through the program “Retorno de Investigadores.” M.R. acknowledges funding by the Deutsche Forschungsgemeinschaft within the priority program SPP 1164 “Microand Nanofluidics” under Grant No. RA 1061/2-1. 1U. M. B. Marconi and P. Tarazona, J. Chem. Phys. 110,8032共1999兲. 2U. M. B. Marconi and P. Tarazona, J. Phys.: Condens. Matter 12, A413 共2000兲. TABLE II. Mean and variance of the force exerted by the polymers in a Stokes flow measured in BD simulations of different boxsizes 共see Sec. III B兲. The forces are given in units of T/D. Simulation box size c*兩具Fz典兩 冑具Fz 2典−具Fz典2 Small 2.2 11 63.4 Large 2.2 6.4 103 Small 22 19 65.6 Large 22 10 103 244906-7 Dynamic density functional theory J. Chem. Phys. 127, 244906 共2007兲 3J. Dzubiella and C. N. Likos, J. Phys.: Condens. Matter 15, L147 共2003兲. 4A. J. Archer, J. Phys.: Condens. Matter 17, 1405 共2005兲. 5A. J. Archer, P. Hopkins, and M. Schmidt, Phys. Rev. E 75, 040501共R兲 共2007兲. 6M. Rex, H. H. Wensink, and H. Löwen, Phys. Rev. E 76, 021403 共2007兲. 7C. P. Royall, J. Dzubiella, M. Schmidt, and A. van Blaaderen, Phys. Rev. Lett. 98, 188304 共2007兲. 8F. Penna, J. Dzubiella, and P. Tarazona, Phys. Rev. E 68, 061407 共2003兲. 9J. Dzubiella, H. Löwen, and C. N. Likos, Phys. Rev. Lett. 91, 248301 共2003兲. 10 M. Krüger and M. Rauscher, J. Chem. Phys. 127, 034905 共2007兲. 11 L. D. Landau and E. M. Lifshitz, Fluid Mechanics, Course of Theoretical Physics Vol. 6, 2nd. ed. 共Butterworth-Heinemann, Elsevier, 2005兲. 12 A. J. Archer and R. Evans, J. Chem. Phys. 121, 4246 共2004兲. 13 A. J. Archer and M. Rauscher, J. Phys. A 37, 9325 共2004兲. 14 M. Rex, C. N. Likos, H. Löwen, and J. Dzubiella, Mol. Phys. 104, 527 共2006兲. 15 C. W. Gardiner, Handbook of Stochastic Methods for Physics, Chemistry and the Natural Sciences, Springer Series in Synergetics Vol. 13, 1st ed. 共Springer, Berlin, 1983兲. 16 H. Risken, The Fokker-Planck Equation, Springer Series in Synergetics Vol. 18 共Springer, Berlin, 1984兲. 17 Note that ⌽has to be uniquely defined up to a constant, also for periodic boundaries of V. For example, for the uniform flow field u=共ux,0,0兲 with periodic boundaries in the xdirection, this is not the case. 18 See http://sourceforge.net/projects/auto2000/ for the source code. 19 S. Harris, J. Phys. A 9, 1895 共1976兲. 20 B. U. Felderhof, J. Phys. A 11,929共1978兲. 21 J. Happel and H. Brenner, Low Reynolds Number Hydrodynamics 共Prentice-Hall, Englewood Cliffs, NJ, 1965兲. 22 E. Lauga and T. M. Squires, Phys. Fluids 17, 103102 共2005兲. 23 D. J. Jeffrey and Y. Onishi, J. Fluid Mech. 139, 261 共1984兲. 24 J. K. G. Dhont, An Introduction to Dynamics of Colloids, Studies in Interface Science Vol. II 共Elsevier, Amsterdam, 1997兲. 244906-8 Rauscher et al. J. Chem. Phys. 127, 244906 共2007兲