Short-time dynamics of colloidal suspensions in confined geometries

We analyze the short-time dynamical behavior of a colloidal suspension in a confined geometry. We analyze the relevant dynamical response of the solvent, and derive the temporal behavior of the velocity autocorrelation function, which exhibits an asymptotic negative algebraic decay. We are able to compare quantitatively with theoretical expressions, and analyze the effects of confinement on the diffusive behavior of the suspension. @S1063-651X~99!09204-1#


I. INTRODUCTION
The dynamics of fluids in confined geometries is a subject of interest, both for its relevance in chemical engineering and environmental science, and because of the new phenomena that arise as a consequence of the competition between the dynamical intrinsic features of the fluid and the geometrical restrictions introduced by the walls.Such an interaction may give rise to a qualitatively different behavior from those predicted in unbounded systems.For example, it is well known that in equilibrium new thermodynamic phases may be induced by the geometrical constraints.Here we will analyze the effects on the dynamics of suspended particles suspended as a consequence of the modification of the solvent hydrodynamics.
In a previous paper ͓1͔, we showed how the presence of constraining walls modifies the dynamics of colloidal particles.We indicated there that the effect of the surfaces on the dynamics of the fluid depends on the specific momentum exchange mechanisms at the interfaces, and not only on the geometrical constraints.For example, for slip boundary conditions we saw that the fluid behaved as an effective medium of dimensionality d* equal to the dimensions of the system that were unconstrained by the presence of the boundaries.On the other hand, for stick surfaces, the propagation of sound modes in the fluid was modified.Depending on the geometry of the constraining medium, the effective sound velocity may become imaginary, developing into a diffusive sound wave.This effect, in turn, qualitatively changes the short-time dynamical behavior of a suspended particle.The asymptotic decay is controlled by the diffusive sound mode, rather than by vorticity diffusion, as happens in unbounded fluids ͓2͔.
The effect that the coupling of a fluid to an elastic matrix has on the dynamical properties of the fluid has been considered previously following the model introduced by Biot ͓3͔.
According to this model, the dynamics of the solid and fluid phases are treated on the same footing, and all the dissipation arises due to the momentum exchange between the fluid and the solid when they are in relative motion.This dissipation is included as a frictional force acting on the fluid, with a friction coefficient that is either left undetermined or is computed assuming that the solid frame is a porous medium and the fluid is filling it.One of the predictions of Biot's theory is that this coupling induces two sound modes in the fluid, one of which has a velocity smaller than the corresponding sound velocities in the two media.Subsequently, it was pointed out that the velocity of the slow sound could give rise to a diffusive sound mode ͓4͔.Analogous ideas were used to analyze the dynamics of a thin fluid layer adsorbed on a solid substrate ͓5͔, as well as the dynamics of gels ͓6͔, where it was argued that, due to the characteristic small sizes of the pores, such a diffusive mode could play an important role in the low frequency response of the gel.
In this paper we will restrict ourselves to a simple geometry, namely, a fluid between plane walls, where an exact hydrodynamic calculation can be carried out without making any assumption about the force that the solid exerts on the fluid.It then becomes apparent that it is only the role of the boundary conditions that determines the hydrodynamic modes of the fluid, and from which the diffusive sound mode can be completely characterized.Although a careful hydrodynamic analysis of the modes of a fluid layer confined between two parallel plates has previously been carried out ͓7͔, it was restricted to high frequencies.We will concentrate on the low frequency regime, and in particular we will consider the behavior of the velocity autocorrelation function ͑VACF͒, defined as for N suspended particles with velocities ͕v ជ i (t)͖, where the average is carried out on an equilibrium fluid over initial conditions.The VACF contains many of the features that characterize the dynamics of the suspension, and it is the simplest example where we can understand the macroscopic properties of the suspension in terms of its microscopic dynamics, since the time integral of the VACF is related to the self-diffusion coefficient of the suspension.We will analyze *Present address: Department of Physics and Astronomy, Univer- the short-time regime for the colloidal dynamics, a time scale in which colloid displacements can be neglected, but in which the solvent hydrodynamics has time to develop.
In Sec.II we will theoretically consider the VACF of a particle suspended in a fluid.To this end, we will determine the relevant fluid modes characteristic of the confined solvent.In Sec.III we present the simulation results, and compare quantitatively with the theoretical predictions.We look in detail at the behavior of the fluid, and we study both the translational and rotational motion of a suspended particle.In Sec.IV we use the results of the preceding sections to study the effect of confinement on the diffusion of colloidal suspensions.We end with a discussion of our results.In Appendix A we give some of the quantities used in Sec.II, and in Appendix B we present an intuitive derivation of the results of Sec.II.

II. THEORETICAL ANALYSIS
In order to study the VACF, we start by analyzing the velocity generated by a localized force perturbation in a pure solvent in the presence of walls, which will show clearly the connection between the dynamics of the suspended particle and the relevant hydrodynamic modes of the solvent.We will consider a fluid confined between two parallel plates, located at yϭϮL/2.For simplicity, we restrict ourselves to the two-dimensional ͑2D͒ situation, in which the confining walls are slits, but it is straightforward to generalize the calculation to higher dimensions.
Since the force pertubation is weak, the fluid density and velocity fields v ជ , obey the linearized continuity and Navier-Stokes equations, respectively, where c s is the speed of sound, and and ⌫ are the kinematic shear and bulk viscosities, respectively.The fluid is initially at rest, perturbed only by the external weak force F ជ applied in the center of the slit and pointing along the direction of the slit, x ˆ, F ជ ϭ f 0 x ˆ␦(t)␦(yϪy 0 ).In Eqs.͑2͒ and ͑3͒ we have disregarded thermal effects, assuming the process to be isothermal, but its generalization to adiabatic processes will not modify the features presented in this paper.Since the fluid is confined between two walls, in order to study the dynamics of the flow generated by a given force it is necessary to specify the boundary conditions at the walls.We will assume that stick boundary conditions are satisfied at the boundaries, which is the usual situation for solid interfaces.In this case, the velocity of the fluid is equal to that of the wall, namely, v ជ (yϭϮL/2)ϭ0.
In order to solve Eqs.͑2͒ and ͑3͒, we Fourier transform them.Defining we can rewrite Eqs.͑2͒ and ͑3͒ as a set of ordinary differential equations where the vectors V and f ជ are defined as The hydrodynamic matrix D ˜, which is written down in Appendix A, has, as eigenvalues, From Eq. ͑5͒ we can express the velocity field at any point in the fluid as

͑9͒
with D being the diagonalized hydrodynamic matrix and M the corresponding eigenvector matrix.Explicit expressions for both of them are given in Appendix A. Finally, C ជ is a vector of constants that should be specified by imposing the boundary conditions in Eq. ͑9͒.Although the derivation so far has been focused on a localized force, we view it as a ͑small͒ colloidal particle.In the limit of a vanishing radius, the identification will be exact.
It is well established that the dynamics of suspended particles can be understood in terms of couplings between the particle variables and the hydrodynamic modes ͓8͔.Since we consider times at which the displacement of the colloid can be neglected, we can take advantage of Onsager's regression hypothesis ͓9͔, and can relate the decay of the velocity of the fluid at point y 0 ͑where the force has been applied͒ with the VACF of a particle placed at that position.The velocity of the fluid at the point where we have initially applied the force is given by and the constants c 1 , . . .,c 4 are the components of the vector C ជ .Substituting their appropriate expressions for stick boundary conditions, one obtains where we have introduced the function H(k x , 1 , 2 ,L,y 0 ), which is written down explicitly in Appendix A. If we now transform back to real space and time, we find

͑12͒
This equation contains the complete information on the dynamics of the particle.We will now focus on its asymptotic behavior, in the situation of strong confinement, i.e., 1 LӶ1 and 2 LӶ1.The opposite limit has been examined in the literature ͓7͔, in the context of interfacial hydrodynamic modes.In that case, the observed modes differ from the usual bulk ones in the sense that, although qualitatively they are not modified, the damping increases due to the presence of the solid walls.
From the structure of Eq. ͑11͒, one can see that the dynamics of the particle is determined by the relevant hydrodynamic modes.The expression of the denominator shows that different dynamical behaviors will exist at different time scales.For example, the confinement will induce resonances of the propagating modes at the time scale in which sound propagates the width of the system.At low frequencies, however, a diffusive mode develops and controls the asymptotic relaxation of fluid perturbations.In this time regime, Eq. ͑11͒ reduces to , which is a diffusive mode.Therefore, density perturbations in the tube will decay diffusively, and we can define an effective diffusion coefficient for sound, D*, as As already pointed out, sound diffusion corresponds to the relevant low frequency dynamical response of the fluid in the tube.In this case we have derived an expression for the diffusion coefficient, which comes from the exact hydrodynamic analysis, and which presents a simple dependence on the geometry considered.It clearly shows that sound diffusion arises purely from the dissipation at the walls, without making any assumption about the momentum exchange.Nevertheless, the expression we obtain coincides with the predictions made on the basis of permeability ͓4͔, which in this case is L 2 /12.The previous analysis is also easily ex-tended to explore the effect of different kinds of boundary conditions on the hydrodynamic response of the fluid.For example, for slip boundary conditions, the hydrodynamic modes correspond to those of a fluid with an effective dimensionality d*ϭDϪ1 ͓1͔.In this case sound is always a propagative mode.
From Eqs. ͑12͒ and ͑13͒, we can derive the asymptotic decay of the velocity of a suspended particle,

͑15͒
The asymptotic behavior is algebraic, and the velocity reverses its direction during its decay.Also for unbounded fluids an algebraic long-time tail is observed ͓2͔.In this case, the slow decay is due to the coupling of the particle velocity to vorticity diffusion.The velocity does not change sign, and the predicted exponent is different.In fact, on the basis of vorticity, one would have expected an exponential decay of the velocity of the suspended particle in the present situation ͓10͔.As we will show in Sec.III, vorticity cannot develop at long times due to the interaction with the walls.It was this fact which led previous authors to predict an exponential decay for the VACF of a particle between walls ͓10͔.Our result shows that the particle velocity couples to the slowest decaying mode in the fluid.The solid walls qualitatively modify the low frequency modes with respect to those characteristics of an unbounded fluid, and the coupling to the diffusing compressible mode leads to a qualitatively different decay for the VACF.
We have restricted our analysis to a two-dimensional fluid, but it can be easily generalized to higher dimensions.In all cases a diffusive sound mode is present, and, accordingly, a power law characterizes the long-time tail of the VACF of a suspended particle.The value of the exponent can be easily related to the number of dimensions of the fluid that are not constrained, d* ͓1͔.In Appendix B we present a more intuitive derivation of the negative algebraic long-time tail.It gives a simple physical picture of the role of the walls on the development of sound modes, and predicts the corresponding exponent of the long-time tail for any dimension.
So far we have presented a hydrodynamic analysis.If the colloidal particle is diffusing, the amplitude of the long-time tail will be modified.It is easy to account for it by making use of mode-coupling theory ͓11͔.In this case the algebraic decay is not modified, and the amplitude of the long-time tail will depend both on and the diffusion coefficient of the particle in the usual way.Finally, we have always assumed that the solid surfaces are completely rigid, which ensures a complete time scale separation between the propagation of the fluid and wall modes.If the surfaces are deformable, sound modes can partially propagate through the walls, leading also to a decrease in the amplitude of the long-time tail.

III. COMPUTER SIMULATIONS
We have performed computer simulations to gain a better understanding of the dynamical effects induced by the confining walls, using the lattice Boltzmann model to simulate a fluid ͓12͔.This method is a preaveraged version of a latticegas cellular automaton model of a fluid.The basic dynamic quantity is the fraction of particles moving in a given direction at a certain lattice node.With this technique it is easy to simulate the dynamics of colloidal particles.These are introduced as surfaces where the collision rules of the populations in the neighboring nodes are modified in order to ensure the appropriate boundary conditions.Therefore, it is also simple to introduce bounding walls for the fluid.We use an implicit updating scheme for the moving boundaries ͓13͔ to avoid instabilities, and to allow us to simulate buoyant particles, while the original scheme could only deal with heavy particles ͓12͔.We have taken the lattice spacing as the unit of length, and the time step as the unit of time.In these units, the speed of sound of the fluid is c s ϭ1/ͱ2, and the kinematic viscosity used is ϭ 1 2 unless otherwise stated.In order to compute the VACF we should, in principle, perform equilibrium averages for the velocities of the particles, taking into account that they move in a fluctuating fluid.Since we are interested in the short-time dynamics, as far as their displacements can be neglected, we can make use of an Onsager regression hypothesis ͓9͔ and study the decay of their velocities from an initial perturbation.In this case, we can then disregard the fluctuations of the fluid, with the corresponding improvement in the simulation performance ͓14͔.In all the simulations we will consider times such that sound has not had time to travel the length of the walls.Therefore, the system can be regarded as infinite, and we do not have to consider finite-size effects.
Throughout this section we will constrain ourselves to the simplest geometry.In this sense, the confining walls will be either planes ͑in 2D slits͒ or a cylindrical tube.When we look at particle dynamics, they will always be spheres ͑disks in two-dimensional͒.

A. Diffusive sound modes
The simplest way to show the diffusive character of sound modes in a confined geometry is by analyzing the temporal decay of a localized density perturbation in the absence of any solid particle.Initially, the density of the fluid is homogeneous except for a node located in the center of a tube.We compute then the moments of the density distribution along the direction of the tube, averaging across the transverse direction, i.e., ͗x(t) n ͘ϭ͐ Ϫϱ ϱ dx͐ ϪL/2 L/2 ((x,y,t)Ϫ 0 ) n dy, with nϭ1,2, . . ., where 0 is the equilibrium density and L is the width of the tube.We have analyzed the behavior of the diffusion coefficient obtained from the second moment of the density distribution, i.e.D*ϭ 1 2 (d/dt)͗x(t) 2 ͘ at long times, as a function of the viscosity and the radius of the tube.In Fig. 1 we display the diffusion coefficient measured as a function of the inverse of the shear viscosity.The dependence of the diffusion coefficient of sound on the viscosity predicted in Eq. ͑14͒ is clearly recovered.We have analogously checked its dependence on the width of the slit.Both scalings also allow us to deduce that the diffusion coefficient measured in the simulations behaves as which agrees quantitatively with Eq. ͑14͒, except for the factors 1 2 and 1, which are due to lattice artifacts.In order to understand how the diffusive regime settles down, in Fig. 2 we show the second and fourth moments of the distribution as a function of time, normalized by the corresponding values for Gaussian diffusion, i.e., ͗x 2 (t)͘ ϭ2D*t and ͗x 4 (t)͘ϭ12(D*t) 2 .We use as a value for D* the one given in Eq. ͑16͒.At long times, both moments approach 1, which means that the diffusive regime settles in.The diffusive regime is reached on a time scale in which momentum can diffuse the width of the tube, ϭL 2 /.One can also see that the second moment approaches the Gaussian behavior faster than the fourth moment.This effect is readily seen by looking at the second cumulant ␣ 2 (t) ϵ͗x 4 (t)͘/"3͗x 2 (t)͘ 2 …Ϫ1, which characterizes the non-Gaussianity of the diffusion process, which is also displayed in Fig. 2. One can see that it is a slowly decaying function.In fact, a study of its asymptotic behavior shows that it vanishes as ␣ 2 (t)ϳ1/t, which implies than, even when the diffusive regime is achieved, still diffusion will be non-Gaussian until longer times.
Finally, in Fig. 3 we display the time evolution of the density profile generated by the initially localized density perturbation.At short times two peaks start to displace in opposite directions at the speed of sound.As time proceeds, these peaks are progressively damped.On the other hand, the region between the two peaks does not relax to the equilibrium density.The inhomogeneity generated in this region as a result of the successive collisions of the initial density perturbation with the solid wall, develops subsequently into a Gaussian and decays diffusively.It is on this time scale that the dynamics is controlled by the diffusion coefficient D* derived in Sec.II.

B. Velocity correlation function for a particle
We will start by analyzing the VACF of a sphere of radius rϭ2.5 in the center of a cylinder of radius Rϭ4.5, as shown in Fig. 4. The initial velocity of the particle is parallel to the axis of the cylinder.The decay coincides with the one obtained for an unbounded fluid until the initial velocity perturbation has reflected back from the wall on the particle.This deviation will scale linearly with the tube width, because it is controlled by sound propagation.At longer times, superposition of sound reflections modifies the decay of the velocity which becomes negative, exhibiting a minimum.Such a minimum appears, roughly, at a time c ϭL/c s , when momentum has propagated the tube width.It is, subsequently, at times when momentum has diffused the tube width that the decay becomes algebraic.As shown in the inset of Fig. 4, the exponent of the algebraic decay is Ϫ 3 2 in this case, which coincides with the theoretical prediction for a 2D fluid of the preceding section.
In order to test the prediction for the amplitude of the long-time tail, we have studied the asymptotic decay of a localized velocity perturbation, and of particles of radius r ϭ0.5 and 2.5, all of them located at the center of a twodimensional slit of width 16.In Fig. 5 we compare the amplitude of the algebraic tail predicted by Eq. ͑15͒ with the simulation results in all three cases.Due to the discreteness of the lattice, there exists an uncertainty about the actual location of the solid boundaries ͓12͔.Both for a point force and for the particle of radius 0.5 this indeterminacy is negligible.For a particle of radius 2.5 we have fitted the curve multiplying the distance from the center of the particle to the wall by a factor 1.1.This assumes a 10% error in the location of the interfaces, leading to errors in the distance that are smaller than one lattice spacing, and which are in agreement with previous estimates of the uncertainty in the particle-wall separations ͓15͔.
We have focused on motions of the particles parallel to the directions that are not confined.If the particle initially moves perpendicularly to the walls, an exponential decay is obtained as could be expected, since in this case the density perturbation generated by the motion of the particle cannot propagate through the system.Therefore, even if in a suspension particles move at random, the long-time behavior will be controlled by the components parallel to the walls.
One can also analyze the behavior of the angular velocity autocorrelation function ͑AVACF͒.The decay is always exponential as long as the particle is in equidistant from the walls.This is due to the fact that the angular velocity does not couple to any compressible mode, and its decay is controlled by vorticity.However, when it is not equidistant, then it will induce a translational velocity, and therefore it will decay algebraically at long times, with the same exponents as those corresponding to the VACF.Again, if slip boundary conditions are considered the situation is qualitatively different.In this case the fluid has an effective dimension d* ϭDϪ1, vorticity can diffuse in this effective medium, and the AVACF decays algebraically with the corresponding algebraic power, without changing its direction of motion during its decay ͓16͔.
If we would have considered the VACF of a particle in a closed container, one expects that eventually the decay will be exponential.For a roughly isotropic container, in which all its dimensions are of the same magnitude, either the decay of the VACF is purely exponential or vorticity develops, when the system is large enough, showing a later crossover to the final exponential decay ͓17͔.If we consider an elongated axisymmetric container, even if vorticity cannot diffuse, if the long length L ͉͉ and the short one L Ќ are related in such a way that momentum can diffuse the short distance before travelling the long one, i.e., L ͉͉ /c s ӷL Ќ 2 /, then diffusion of sound along L ͉͉ will have time to develop.In this case, for times L ͉͉ /c s ϾtϾL Ќ 2 / the decay of the VACF of a particle will be algebraic, and only at later times it will become exponential.In this case, the algebraic decay is not the asymptotic behavior of the short-time dynamics, as seen in Fig. 6.

C. Hydrodynamic fields
We will now look in detail at the flow fields that are generated in a fluid confined between two parallel plates as a response of an initially moving sphere in an otherwise quiescent fluid.We will consider a particle of radius rϭ2.5 lattice units, placed in the center of a tube of radius Rϭ8.5, in a fluid of shear viscosity ϭ 1 6 .We focus on the features of the flow, vorticity, and density fields for a particle diffusing along the axis of the tube, during the initial decay of the VACF, i.e., for times at which the asymptotic decay has not been reached yet.
The total simulation time corresponds to t/ Ϸ1.The flow, vorticity, and density fields depend on the three spatial coordinates.The fields are examined in the zϭ0 plane.The velocity field and density are directly calculated in the lattice-Boltzmann scheme.The vorticity, , at a lattice point (x,y), which is a local quantity is given by ͑x,y ͒ϭ This is a two-dimensional discretization of ϭ 1 2 ٌϫv, the three-dimensional definition of vorticity.The time evolution of the various fields is shown in Figs.7 and 8.
Initially, the velocity field contains the vortex pair that is responsible for the positive long-time tail in the VACF for an unbounded fluid, although as soon as the perturbation induced by the particle reaches the walls vortices creeping at the wall of the tube are also seen.On the other hand, the motion of the particle produces an increase of the density in the front, and a decrease in the back, as expected.In Fig. 7 it is already seen that the two vortices cannot grow diffusively, and that meanwhile the vortices that were creeping along the wall are already outside the picture.The density plot indicates that in front of the particle there is a sizable density increase.At later times, when the particle velocity reverses direction, a clear change in the qualitative behavior of the FIG. 5. Amplitude of the long-time tail measured in the simulation A sim relative to the theoretical prediction A th given in Eq. ͑15͒ as a function of the distance from the center of the disturbance or from the particle center y to the center of the two-dimensional tube, y c , expressed in lattice units.The tube width is Lϭ16, and the viscosity ϭ 1 2 .
FIG. 6. Decay of the VACF for a disk of radius rϭ2.5 in the center of a box of sides L y ϭ9 and L x ϭ199.The coefficient ␣ used to show the long-time exponential decay has been fitted.We have used ␣ϭ0.006.Lengths are expressed in lattice units.flow field is observed.The two vortices have been damped out, which shows that the perturbation induced by vorticity decays exponentially due to the presence of the walls ͓10͔, and over a significant area the velocity in the x direction is opposed to the direction at which the particle was initially moving.This feature becomes more salient at later times, as can be seen in Fig. 8, which corresponds to the minimum of the VACF.Due to the back flow around the particle, the vorticity near the colloidal particle has almost disappeared.At later times the negative flow field keeps growing in size, at the end spanning the total plotted region.The velocity field resembles the Poiseuille flow field for these times, which is the parabolic steady state flow profile for flow in a tube, with a pressure drop, and the density field becomes more one dimensional.

IV. DIFFUSION COEFFICIENTS
The self-diffusion of the suspended particle can be related to the integral of the VACF by means of a Green-Kubo expression.We can then obtain macroscopic information about the behavior of the suspension, from the knowledge gained in the analysis of its macroscopic behavior.
In Fig. 9 we show the values for the diffusion coefficient for a sphere of radius 2.5 at the center of a cylindrical tube of variable radius, normalized by the value of the diffusion coefficient in the absence of the tube.We also plot the prediction of the center-line approximation ͓18͔, calculated assuming a purely incompressible fluid, i.e., neglecting sound effects.We obtain a perfect agreement with the theoretical predictions.This can be understood, because the algebraic decay we have described in the preceding sections arises from the coupling to the compressible modes of the fluid, and these vanish in the stationary limit.These modes do not affect the integral of the VACF since their amplitude is in fact proportional to the frequency itself.The diffusion coefficient is determined by the decay of transverse velocity perturbations ͓10͔ which, in a confined system, is exponential.So, while compressibility effects dominate the long-time dynamics, they still do not contribute to the diffusion coefficient.However, it is important to resolve all the time scales in order to obtain the proper diffusion coefficient from the integral of the VACF.In Fig. 10, we show the time dependent diffusion coefficient, which is the integral of the VACF from 0 to time t.One can see that it exhibits a maximum.If we would have assumed an exponential decay for the VACF we would have obtained an overestimate for the value of the diffusion coefficient.
Due to the flexibility of the lattice Boltzmann method to deal with complex geometries, we can use the technique as a means to derive profiles for the diffusion coefficient of suspended particles in general geometries.In Fig. 11, we show the values of the diffusion coefficient for a spherical particle of radius 2.5 in a cylinder of radius 9.As we move offcenter, the three components of the diffusion coefficient are not equal.In particular, the component of the diffusion parallel to the tube (D x in the figure͒, exhibits a maximum as FIG. 7. Flow field, vorticity, and density, from top to bottom, respectively (t/ ϭ0.47).The velocity field, which is at the top, is scaled with sϭMax v x in the x direction, whereas for the y direction 4s is used to scale the arrows.For the sake of clarity one half of the velocity vectors in both directions are omitted from the picture.The vorticity field is obtained by applying Eq. ͑17͒ to the velocity field, and is shown in the middle.It is scaled such that ten isovorticity lines are shown, at heights which vary linearly between the maximum and minimum vorticities.The largest vorticity is depicted by the darkest color.The same procedure was performed for the density, the bottom picture.we move toward the plate.A perturbative calculation ͓19͔ implies the existence of a maximum in the parallel component of the diffusion for the motion inside a cylinder, although we do not know of any explicit calculation of the value of the maximum.The same behavior persists for a two-dimensional fluid, as shown in Fig. 12.However, the diffusion coefficient parallel to the the walls is a monotonous function for a three-dimensional fluid between two parallel plates.In Fig. 13 we show the corresponding profile, and we compare it with the value for the density predicted theoretically if we neglect the hydrodynamic interactions between the two walls.In this case, the diffusion coefficient of the sphere is the sum of the diffusion coefficients induced by each wall independently, corrected by the value at the center of the layer.For the diffusion coefficient of a sphere in the presence of a single wall a theoretical expression exists that is sufficiently accurate except close to contact ͓20͔.The plot shows that the agreement is very good, even for a relatively narrow fluid layer as the one considered, in which the aspect ratio is ϵ2r/L y ϭ1/3, with r being the radius of the colloid and L y the width of the tube.The deviations that can be expected close to contact, where lubrication will be important, cannot be resolved in a lattice Boltzmann calculation with the small sphere considered in these simulations.

Diffusion coefficient for a confined suspension
We have also considered the VACF for a suspension in a slit.The time dependence observed is analogous to the obtained in the preceding sections for a single particle.By in-FIG.10.Dependence of the diffusion coefficient with time for a particle of radius rϭ2.5 in a slit of width L y ϭ7.t w ϭ(L/2) 2 / is the time needed the vorticity to diffuse the width of the tube.FIG.11.Profile of the diffusion coefficient of a particle of radius 2.5 in a tube of width 9 lattice spacings, as a function of its position off-center.x is the direction of the axis of the cylinder, and we move off-center in the z direction.D refers to the average local diffusion coefficient defined as the mean of its three components.tegrating, we obtain the short-time diffusion coefficient for a colloidal suspension, and we analyze its behavior as a function of the aspect ratio .It is clear that the wider the tube the larger the diffusion coefficient, since the mobility is hindered due to the hydrodynamic interaction of the particles with the walls.This is clearly depicted in Fig. 14, where we show the decay of the value of the diffusion coefficient for a particle equidistributed in the tube, as a function of .However, if we normalize the diffusion coefficients of the suspension by the average mobility at zero volume fraction, then the curves for the diffusion coefficients of suspensions corresponding to different look much more similar to each other, as shown in Fig. 15.In fact, a significant deviation from the general behavior is only observed for strong confinement.As soon as roughly two layers of suspended particles fit between the walls, their relative motion becomes dominant, and the walls can be seen as if they provide only the mean friction felt by the colloidal particles.

V. CONCLUSIONS
In this paper we have analyzed the short-time dynamics of confined colloidal suspensions.We have restricted ourselves to the simplest geometry, namely, a tube, where exact analytical results can be derived for point particles and one can focus on the basic dynamical features of the short-time colloidal dynamics.
We have shown that the presence of stick walls qualitatively modifies the hydrodynamic modes that characterize the response of the fluid.The coupling of the fluid to a solid elastic medium ͑which is a typical example where stick boundary conditions are fulfilled͒ induces the development of a diffusive sound mode at low frequencies.We have also shown that this diffusive mode is initially non-Gaussian, and it is only on longer time scales that it becomes Gaussian.This diffusive mode decays slower than vorticity, which cannot propagate diffusively due to the presence of the walls, and it becomes the leading hydrodynamic response at low frequencies.The dynamics of the suspension will couple to the longest-lived solvent mode.Therefore, the long-time decay of the velocity autocorrelation function is controlled by sound, contrary to what happens in an unbounded fluid, where it is controlled by vorticity.The diffusive decay of the initial density perturbation induces a change in the direction of the suspended objects, as was clearly seen from the pictures of the hydrodynamic fields.Also, the angular velocity autocorrelation function differs from the usual behavior in unbounded fluids.Since its decay is due to vorticity, it always exhibit an exponential decay, insofar as it is not coupled with the velocity.Due to the presence of the walls, this is an unrealistic requirement for colloidal suspensions.As soon as the particles are not equidistant from the walls, an initial angular velocity will induce a translational velocity of the particle ͓21͔.This coupling term will decay slower and will control the subsequent decay of the AVACF.
From the analysis of the VACF we have extracted the values for the diffusion coefficient of a confined suspension.In order to obtain proper values, it is necessary to take into account the algebraic decay of the VACF, although the final values for the diffusion coefficients agree with the predictions for incompressible fluids.We have looked at the local values of the diffusion coefficient for different geometries.Due to the presence of the walls the diffusion is no longer isotropic, and we have focused on the component of the diffusion parallel to the walls.We have seen that it may exhibit a maximum, although it is not a generic feature.For example, in the case of two parallel plates it has a monotonous behavior, and we have been able to compare our values with the predictions obtained from the additivity assumption, which neglects the hydrodynamic interactions between the two walls.This approximation is seen to work quite well even for quite narrow tubes.The flexibility of the lattice Boltzmann method in dealing with general boundaries converts it into a useful technique to obtain maps of the diffusion coefficient in general geometries, that can be used as inputs for other simulations techniques that work in the particle diffusion time scale, such as Brownian dynamics ͓22͔.
The development of the diffusive sound mode can be understood on a simple physical basis.For time and length scales in which the fluid is interacting with the solid interfaces, momentum is no longer a conserved variable.Under these conditions, only mass is conserved, and has therefore a diffusive character.This interpretation provides also an analogy with a Lorentz gas.This model is known to exhibit an algebraic velocity decay at long times with an exponent ϪD/2ϩ1 ͓23͔ for a D-dimensional system, and which is precisely the one we have obtained from solving the Navier-Stokes equations, taking D as the effective dimension in which the density is the only conserved variable.This simple connection explains the robustness of the results presented in this paper.We have studied the decay of the VACF for particles moving in a porous medium built up with fixed random particles, and also the decay of a particle moving in a regular array of fixed rectangles, obtaining in all cases a negative algebraic decay, with an exponent that can be understood as the number of unconstrained dimensions.In these two systems, the solid walls do not completely restrict any of the spatial dimensions, and we find, consistently, d*ϭD, D being the dimensionality of the system.This also implies that the same kind of dynamical behavior we have described here should apply to the dynamics of particles embedded in gels or in suspended membranes, where in the latter case the possibility of developing diffusive sound modes has also been pointed out ͓24͔.
Finally, if we would have considered that the fluid is contained between slip walls, then the fluid would develop after a transient into an effective fluid of dimension d*.This clearly shows that the dynamical behavior of the suspensions in confined geometries will be sensitive to the specific momentum fluxes that occur at the boundaries.
which after diagonalizing gives The eigenvector matrix M corresponding to D, and which is needed to derive the appropriate expression for the velocity field from Eq. ͑9͒, has the form Finally, the function H(k x , 1 , 2 ,L,y 0 ) introduced in Eq. ͑11͒ and which specifies the velocity at the point where an initial perturbation is applied in the fluid is given by PRE 59 4467 SHORT-TIME DYNAMICS OF COLLOIDAL . . .

FIG. 1 .
FIG. 1.The effective diffusion coefficient of density perturbations D* as a function of 1/.Results were obtained in a twodimensional slit of width Lϭ9.The points denote the simulation results, and the line is a guide to the eye.D* and are expressed in lattice units.

FIG. 2 .
FIG. 2. Second and fourth moments, and second cumulant ␣ 2 of of the density distribution as a function of time for an initially localized density perturbation in a 2D slit of width Lϭ10 and kinematic viscosity ϭ0.5.The moments are normalized by the corresponding moment for Gaussian diffusion.Dimensional quantities are expressed in lattice units.

FIG. 12 .
FIG.12.Profile for the component of the diffusion coefficient parallel to the walls in a two-dimensional slit for two different tube widths, L y .The particle has always radius rϭ2.5 lattice units.
FIG.15.Normalized diffusion coefficient as a function of the volume fraction for a confined suspension, as a function of the aspect ratio.In all cases the particles have a radius rϭ2.5 lattice units.