Pairwise hydrodynamic interactions and diffusion in a vesicle suspension
Abstract
The hydrodynamic interaction of two deformable vesicles in shear flow induces a net displacement, in most cases an increase of their distance in the transverse direction. The statistical average of these interactions leads to shearinduced diffusion in the suspension, both at the level of individual particles which experience a random walk made of successive interactions, and at the level of suspension where a nonlinear downgradient diffusion takes place, an important ingredient in the structuring of suspension flows. We make an experimental and computational study of the interaction of a pair of lipid vesicles in shear flow by varying physical parameters, and investigate the decay of the net lateral displacement with the distance between the streamlines on which the vesicles are initially located. This decay and its dependency upon vesicle properties can be accounted for by a simple model based on the well established law for the lateral drift of a vesicle in the vicinity of a wall. In the semidilute regime, a determination of selfdiffusion coefficients is presented.
pacs:
82.70.Uv, 83.80.Lz, 83.50.HaNow at ]Laboratoire de Biomécanique et de Bioingénierie, UMR CNRS 7338,
Université de Technologie de Compiègne, France.
Now at ]Zentralinstitut für Medizintechnik, Technische Universität München, Germany
I Introduction
Liquid suspensions of deformable particles are the focus of permanent interest due to their ubiquity in life science and applications, from emulsions to blood, a dense suspension of red blood cells. It is well known since Batchelor Batchelor and Green (1972); Zinchenko (1984) that the viscosity of a semidilute suspension of rigid spheres departs from the classic linear Einstein law of viscosity for volume fractions of particles above a few percents, due to the additional dissipation induced by hydrodynamic interactions between particles.
In addition to their influence on the effective viscosity at significant volume fractions, these hydrodynamic interactions (which are sometimes called binary collisions in the literature, although still mediated by hydrodynamics) can lead to irreversible perturbations of the particle trajectories which result in an effective random walk of individual particles in the suspension. This shearinduced diffusion has two main consequences: enhanced mixing and transport even at low Reynolds number Goldsmith (1971); Goldsmith and Marlow (1979); Cha and Beissinger (2001); Higgins et al. (2009); Zhao and Shaqfeh (2011); Tan, Le, and Chiam (2012), and a modification of the structure of suspensions via diffusion along gradients of concentration of the particles Hudson (2003); Rusconi and Stone (2008); Podgorski et al. (2011); Grandchamp et al. (2013).
In shear flow, two identical particles located on different streamlines and moving towards each other will generally experience irreversible drift in the shear and vorticity directions after they have interacted. However, for smooth rigid spherical particles in a dilute regime dominated by pairwise interactions, the crossstream lateral displacement is expected to be negligible at low Reynolds number since trajectories must be symmetric due to the flowreversal symmetry of the Stokes equation Bretherton (1962) and the symmetry of the geometrical configuration. Symmetry breaking can be obtained by considering rough particles Da Cunha and Hinch (1996); Blanc, Peters, and Lemaire (2011) or deformable twofluid systems such as bubbles Van Wijngaarden and Jeffrey (1976) or drops Loewenberg and Hinch (1997); Guido and Simeone (1998); Singh and Sarkar (2009); Le and Chiam (2011). More recently, systems made of closed membranes have been investigated numerically or experimentally. Hydrodynamic interaction between elastic capsules were studied numerically in several papers Lac, Morel, and BarthèsBiesel (2007); Lac and BarthèsBiesel (2008); Pranay et al. (2010); Omori et al. (2013). During the interaction, net displacement of the capsules can be coupled with wrinkling or buckling of the membranes, whose tension strongly increases during interaction.
Set  Solution  (mPa.s)  

1 


1.0  
2 


3.8  
3 


0.28 
The dynamics and rheology of suspensions of lipid vesicles have recently been the focus of several studies, due to their relevance to the understanding of blood flows, considering giant vesicles as models of red blood cells, and the challenging theoretical questions they pose as a consequence of their rich microscopic dynamics. Vesicles are closed lipid bilayers with mechanical properties similar to those of living cells. A key property is the membrane inextensibility, which leads to local area conservation, while volume conservation is generally obtained once osmotic equilibrium is reached. The vesicles mechanical response, as well as their shapes Helfrich (1973), are governed by a bending energy of order a few , where is Boltzmann’s constant and the temperature. These particular properties, especially the nonlinearities due to the constraint of local area conservation, are responsible for the various dynamics of single vesicles in shear flow Deschamps, Kantsler, and Steinberg (2009); Farutin, Biben, and Misbah (2010); Biben, Farutin, and Misbah (2011). The phase diagram of microscopic dynamics has a signature on the rheology of vesicle suspensions Danker and Misbah (2007), but there is still disagreement, especially in the semidilute regime where two experimental studies show contradictory results Kantsler, Segre, and Steinberg (2008); Vitkova et al. (2008). In an effort to resolve this contradiction, Kantsler et al. Kantsler, Segre, and Steinberg (2008) and Levant et al. Levant et al. (2012) have investigated the influence of interactions on fluctuations and correlations of the inclination angles of interacting vesicles and suggest that they may be responsible for discrepancies between theories in the dilute regime and experimental measurements of the effective viscosity, which are often made in a semidilute regime for sensitivity reasons. On the analytical side, the trajectories of interacting vesicles have been recently studied in the limit were they are initially very distant from each other Gires, Danker, and Misbah (2012). This study has been later on refined in the case of vesicles located in the same shear plane Farutin and Misbah (2013). Very recently, such trajectories have been calculated numerically by Zhao and Shaqfeh Zhao and Shaqfeh (2013). They also calculated the rheology of a semidilute suspension and found good agreement with the experiments by Vitkova et al. Vitkova et al. (2008), a strong indication that, in that concentration regime, interactions between vesicles cannot explain the contradiction between the latter experiments and the one by Kantsler et al. Kantsler, Segre, and Steinberg (2008)
Along with their influence on rheology, hydrodynamic interactions significantly affect the structure of suspensions, especially in confined flows where a balance between migration away from walls and shearinduced diffusion due to repulsive interactions leads to the formation of a nonhomogeneous distribution of vesicles Podgorski et al. (2011). During heterogeneous interactions of vesicles or capsules with different mechanical or geometrical characteristics, asymmetric displacements take place, which leads to segregation or margination Podgorski et al. (2011); Zhao and Shaqfeh (2011); Narsimhan, Zhao, and Shaqfeh (2012); Kumar and Graham (2012); Fedosov, Fornleitner, and Gompper (2012), a phenomenon also observed in blood flows Goldsmith and Spain (1984); Tilles and Eckstein (1987); Uijttewaal et al. (1993); Yeh and Eckstein (1994).
In this paper, we report on our experimental and numerical investigation of the interaction of two identical vesicles in shear flow, with a focus on the net lateral displacement as a function of initial configuration and vesicle properties.
With a good agreement between experiments and simulations, the amplitude of the lateral displacement is found to be weakly dependent on vesicle deflation and viscosity ratio, at least in the tanktreading regime to which we restrict our study. Thanks to the simulations, we also discuss to which extent the discrepancies between the ideal case of two identical and neutrally buoyant vesicles, placed in the shear plane of an infinite simple shear flow, and the realistic case of channel flow, influence the final result.
Finally, from the numerical results, an evaluation of the selfdiffusion coefficient, obtained by averaging displacements over all initial configurations, is proposed and compared to the recent results of Zhao and Shaqfeh Zhao and Shaqfeh (2013).
Ii Experimental setup
Fluid vesicles are prepared by following the electroformation method Angelova et al. (1992), which produces vesicles of various size and deflation (that is, the surface to volume ratio). They are made of a dioleoylphosphatidylcholine (DOPC) lipid bilayer. We consider three sets of outer and inner solutions in order to vary the viscosity ratio between the inner and the outer fluids (see table 1). The different additives (sugars and dextran) used for inner and outer solutions provide an optical index contrast which is convenient for phase contrast microscopy.
We wish to observe interactions in simple shear flow between vesicles located in the same plane, where is the flow direction and the shear direction. To that end, the vesicle suspension is injected in a standard polydimethylsiloxane (PDMS) microfluidic device. The observation channel is 184 µm wide ( direction) and 100 µm deep ( direction). The imposed flow is along the axis (see Fig. 1). Before the observation section, vesicles flow in a channel of several centimeters long, so that centering in the direction is generally rather well achieved Coupier et al. (2008), as confirmed by the location of all vesicles within a focal plane of thickness of order 5 µm. The interacting vesicles were followed manually translating the stage. The observation window is 477x358 , with a resolution of µm/pixel. The use of a channel flow, rather than a fourroll mill deviceKantsler, Segre, and Steinberg (2008); Levant et al. (2012), allows to measure the final lateral displacement due to the interaction, a key parameter in the discussion of diffusion phenomena.
As measurable interactions only occur when vesicles have initial separation not larger than 2 radii, it appeared necessary to favor such an initial condition by adding a flow focusing device at the entrance of the observation channel. Two lateral inlets were then added, where vesiclefree fluid was injected in order to focus the suspension in a narrow area. This area is located at around one fourth of the total width of the channel, that is, far from the wall and far from the center, where the flow can be considered as a simple shear flow, in a first approximation to be discussed later. Vesicles stay at this favorable position thanks to the balance between lift forces Abkarian, Lartigue, and Viallat (2002); Callens et al. (2008); Coupier et al. (2008) and gravity. In addition, dilution by the lateral inlets decreases the probability of perturbation of the interaction trajectories by other vesicles.
In the observation window, at most 3 or 4 vesicles (including the two studied vesicles) are present at the same time. In the selected interaction sequences, the additional vesicles of non negligible size are always at a distance from the pair larger than 5 radii and are located almost on the same streamlines, so that they will not come close to the pair within the duration of the studied interaction process. As in Fig. 2, very small vesicles may come closer, but the induced perturbation is expected to be negligible: from an asymptotic approachGires, Danker, and Misbah (2012), we can expect the velocity perturbations induced by a vesicle 4 times smaller than the studied ones to be smaller than the one coming from the interacting vesicles by a factor .
Once an appropriate pair of vesicles is chosen, the vesicles are followed along their trajectories and the coordinates of the vector linking their geometrical centers are determined, as well as their shapes. We denote by the initial position and by the final one. By convention, and . An example of selected snapshots taken along a trajectory is shown in Fig. 2. As we only have access to their twodimensional crosssection in the plane, we characterize the 2D shapes by the effective radius , , defined by , where is the crosssection perimeter, and by a reduced area , where is the crosssectional area. These two parameters are evaluated before the vesicles strongly interact, at which point outofplane deformations occur.
In this study, we focus on pairs of vesicles of similar size and deflation (within maximal variations of 10 percent for the radii and 5 percent for the reduced area). We denote by and the arithmetic averages of the radii and reduced areas of the two interacting vesicles. lies between 5 and 19 µm, and between 0.73 and 1. The flow velocity is set so that the capillary number lies typically between 10 and 100. This capillary number qualitatively represents the ratio between the magnitude of the liquid viscous stresses exerted on a membrane, and its resisting bending stresses, controlled by bending rigidity . is the shear rate and the suspending fluid viscosity. is the effective radius of the vesicle, defined from its volume by . Note that, due to vesicle volume conservation, this 3D effective radius is constant and characteristic of the considered vesicle, while the observed 2D radius depends on the applied flow. As can only be roughly estimated in the experiments, we only have access to estimated values for .
From the obtained trajectories, we extract the main information, that is the lateral displacement as a function of initial lateral separation . Both distances are rescaled by . Several initial positions are scanned either by considering different pairs of vesicles, or by making a given pair going back and forth thanks to flow reversal.
Iii Experimental results
We first focus on vesicles with no viscosity contrast. Results for as a function of are shown in Fig. 3. Initial and final positions are measured by averaging over several positions long before and after interaction. Error bars are associated to the fluctuations in these positions due to the presence of other small vesicles or flow variations due to channel roughness. Such events are likely to occur because of the large ratio between the relative velocity along the axis between the considered vesicles and the other vesicles or the wall, and the velocity along the axis. The studied pairs are split into two subpopulations according to their reduced areas. Vesicles with undergo negligible deviations which are not measurable within experimental errors. This result is expected for spherical particles and allows to check that no uncontrolled drift alters the experimental results. All other pairs of vesicles yield comparable deviations whatever the reduced area in the range . Data scattering can be due to variations in reduced areas (including within a pair), sizes, capillary numbers, but also to non complete colocation in the same plane. In addition, displacements might be affected by the flow perturbation induced by the walls, which depends on the lateral position of the vesicles, a parameter that varies from one pair to another.
Lateral deviation is a decreasing function of initial lateral separation, and becomes negligible for initial lateral separations greater than one diameter. An empirical estimate for this deviation can be obtained by considering that, in the reference frame of the bottom vesicle, the displacement of the top vesicle is due to the interaction with a wall of finite extent in the direction, whose role is played by the bottom vesicle. The lift velocity of a vesicle in a simple shear flow and at a distance from a wall was experimentally shown Callens et al. (2008) to agree with the scaling law suggested or confirmed by several theoretical works Olla (1997); Sukumaran and Seifert (2001); Vlahovska and Serral Gracia (2007); Farutin and Misbah (2013), where is a dimensionless parameter that depends on viscosity ratio and reduced volume. The reduced volume is defined as the ratio between the vesicle volume and the volume of the sphere having the same area; due to volume conservation and membrane incompressibility, this is a constant parameter that characterizes the deflation of the vesicle. The top vesicle flows with relative velocity , so that . Interaction takes place on a finite distance of order . Integrating the latter equation on this distance for and between and for , one finally finds
(1) 
where contains the interaction details. is the maximal displacement, obtained for . Taking an estimate of from Olla’s work Olla (1997), we have for a prolate ellipsoid with , , which is close to the maximal displacement seen in Fig. 3, where a full fit of the whole data set with Eq. 1 yields . Note that this is a single parameter fit, so that the distance at which interaction becomes negligible is fully determined by the maximal deviation obtained for quasialigned vesicles. In particular, vesicles initially distant by one diameter in the direction will deviate only by from their initial trajectories.
From this law, we can estimate how the final displacement should vary with reduced volume. Following for instance the recent study by Farutin and Misbah Farutin and Misbah (2013), scales as , so that the maximal displacement scales as . This sharp increase around explains the strong difference between the quasispherical vesicles () and the more deflated ones () seen on Fig. 3. On the other hand, from Olla’s results Olla (1997), is multiplied by a factor 2.7 between prolate vesicles of reduced areas 0.98 and 0.82, respectively. The maximum displacement for vanishing should then be multiplied by . Such a tiny variation is within the scattering and error of experimental data.
Similarly, when the viscosity ratio is varied, no significant displacement variation is observed, as shown in Fig. 4. Once again, this is consistent with our empirical law, since from Olla’s results again, for , decreases only by 2% between vesicles of viscosity ratio 0.28 and 1, and by 12% between vesicles of viscosity ratio 1 and 3.8. According to Zhao and ShaqfehZhao and Shaqfeh (2013), the maximum displacement drops by about 30% between viscosity ratio 1 and 7.
In the next section, we address the same questions with full 3D numerical simulations restricted to the case , following our discussion on the weak influence of the viscosity ratio in the previous section. We then confront the numerical results with the experimental ones.
Iv Model and numerical method
iv.1 Liquid and membrane
In this section we outline the model and numerical method. The internal and external liquids are modeled as incompressible, homogeneous, Newtonian fluids. We restrict our study to the case where their densities, as well as their viscosities, are equal. Both liquids flow in the creeping regime.
The membranes are modeled by two dimensional surfaces. As for the liquids, their inertia is negligible. Their areas stay locally constant. They resist bending with an energy , given by Zhongcan and Helfrich (1989)
(2) 
where is the membrane surface, the bending rigidity, and the mean curvature. The sign convention for the curvatures is taken so that the mean curvature of a sphere is negative.
The resulting surface force density that the membrane exerts on the fluids is
where is the unit normal vector pointing outward, the Gaussian curvature, and a Lagrange multiplier that enters the total energy, obtained by adding to (2) . It ensures local membrane incompressibility and satisfies:
(3) 
where is the surface gradient operator and is the membrane velocity.
Note that we don’t include in our model any other small range interaction than the hydrodynamic forces within the lubricating film, described as squeezed between athermal membranes. As a first approximation, we considered that theses stresses grow fast enough so that the minimal distance between the membranes, that we denote , remains higher than a typical distance under which other type of interactions become significant. The first one that would appear is linked to the inhibition of thermal fluctuationsde Haas et al. (1997), which leads to an entropic repulsion pressure. It is of order . It would balance the imposed pressure, estimated as , that tends to push the two vesicles towards each other, if
(4) 
Using the typical value , for the smallest shear rate in our experiments , one finds that reaches values in the range of nm. We checked that, in the trajectories we investigated, remains higher than the previous estimate. The facts that the entropic force is repulsive, and that, on the contrary to some rigid particles, there are no heterogeneities on the phospholipid membranes that can facilitate the drainage of the lubricating film, support even more our approximation.
iv.2 Boundary conditions
The membranes are supposed to be at osmotic equilibrium and are modeled as impermeable. Together with the no slip boundary condition, this leads to an advection of the membranes with the local velocity of the flow.
A force balance on the membrane yields
(5) 
where is the liquid stress tensor with a or superscript respectively for the external and internal fluids, defined as . Far from the vesicles, the imposed simple shear flow is recovered.
We denote by the vector linking the centers of mass of the two vesicles. We shall study the evolution of as a function of , that is, the trajectory of vesicle 2 in the frame centered on vesicle 1. Different initial positions will be scanned, with initial longitudinal distance much larger than the vesicles radii. A sketch of the initial state of the system is presented in Fig. 5.
iv.3 Numerical method
The full set of equations in the Stokes regime can be converted into a boundary integral formulation Pozrikidis (1992). The integral equation (recalled below) is solved numerically in three dimensions following the work by Biben et al.Biben, Farutin, and Misbah (2011). The new elements of the present study are the extensions to two vesicles and, in a second time, to the presence of a wall that turns out to be a relevant ingredient when confronting the numerical results with the experimental one. We shall first study the situation without wall, which is the main goal of the paper.
The integral equation provides the expression of the membrane velocities as a function of boundary integrals and reads
(6) 
where is the position vector of a membrane point, the boundaries present in the system under consideration, which are in the present case the two vesicle membranes, and the Green’s function of an incompressible fluid following Stokes equation. As we consider an unbounded domain, an appropriate choice is the Green’s function associated to a point force in an infinite liquid, such that , where Pozrikidis (1992)
(7) 
For most simulations, the vesicles are meshed by vertices, and the time step is . We checked that for a typical trajectory (), results were relatively independent from a reduction of the mesh size and time step: increasing the number of vertices to and reducing the time step by a factor 2 led to relative changes in the transverse migration of . A challenge is to achieve an evolution of the membrane shapes ensuring a local conservation of the area. We present in appendix A details showing that our study conserves the area with a good approximation. The simulations start with both vesicles having the steady inclination angle of an isolated vesicle in shear flow, obtained from a preliminary simulation.
V Numerical results
v.1 Identical vesicles in the same shear plane
We start with the case of identical vesicles, with the typical parameters , in the same shear plane of an infinite simple shear flow. The capillary number is taken in the set , so that the whole range of possible experimental values is covered. We plot in Fig. 6 the interaction curve (that is the difference between the final and initial positions), with initial and final distances corresponding to and .
For all values of , we recover the decrease of from around to zero. All deviations become smaller than for . There is a good agreement with the simple model based on the law for the lift of a vesicle near a wall, that was presented in the preceding section. It thus validates this model as a convenient tool to anticipate the dependency of the lift with the mechanical properties of the vesicles. Note however that, since the shape in Olla’s model is prescribed, no dependency on can arise from it.
Overall, considering that there are no fitting parameters (but some experimental uncertainty on ), the agreement is rather satisfactory. However, the experiments lead to smaller displacements, as the numerical curve passes through the error bars of only about of the experimental points. This discrepancy may be explained by differences between the experimental configuration and the ideal unbounded simple shear flow on three aspects. First, the suspension is slightly polydisperse, both in shape and size. Second, the centering in the direction might not be perfect. Indeed, the depth of focus of the microscope is about µm, which allows to differ from by amounts up to . Third, the balance between wallinduced lift forces and sedimentation in the direction is perturbed during the interaction, and may also not be fully reached before interaction starts, because of preceding interactions. All those effects could be non negligible. We use the numerical model to study their relative importance.
v.2 Departure from interaction of two identical vesicles in a shear plane of a simple shear flow
Influence of polydispersity
Regarding the influence of polydispersity, we computed several sets of interaction curves, with , first changing the radius ratio so that , and then both reduced volumes, so that . We find relatively small effects, not sufficient to explain all the data scattering : for instance, for , the maximal variation in is . Such a small effect was expected from the qualitative discussion presented in Sec. III.
Influence of
We plot the interaction curve for . The result is reported in Fig. 7.
As expected, the deviation decreases with . Considering that the vesicles can be initially shifted in the vorticity direction by the maximal distance allowed by the focal depth of the microscope, a better agreement between experimental data and simulations is found (about 70% of experimental points, for vesicles of radii 10 µm).
Influence of the bottom channel wall
We consider now the influence of an imbalance between walloutward migration and sedimentation. For simplicity and since gravity acts similarly on both vesicles, we only consider the wall migration effect. As lift is a decreasing function of the distance to the wall Coupier et al. (2008), we expect the upper vesicle to migrate less relatively to the wall, so that the distance between the two vesicles is indirectly reduced due to that wallinduced lift forces.
We compare the whole trajectory obtained by our code with the one corresponding to the experiment shown on Fig. 2. The geometrical input parameters of the simulation are the reduced volume and the 3D effective radius , in contrast with the experimentally measured reduced area and the 2D effective radius . From the study of isolated vesicles in simple shear flow, we find that, for , vesicles having same 2D crosssections as the vesicles of Fig. 2 are characterized respectively by and .
In order to quantify the bottom wall effect, we adopt the Green’s function corresponding to a semiinfinite fluid Pozrikidis (1992); Kaoui et al. (2009), and include the quadratic part of the flow in the plane of shear, so that the imposed flow is . The initial distance of vesicle 1 from the wall is µm. The comparison between the experiment and the numerical study is presented in Fig. 8, without and with wall, for vesicles in the same shear plane (). A possible shift is also considered together with the presence of the wall.
As expected, lift by the wall leads to a slight initial attraction (a decrease of for ), which results in a slightly smaller final lateral displacement when . It appears however that this correction is too small to account for the remaining discrepancy between the simulations and the experiments, for which the initial attraction of around 1 µm, that is seen on Fig. 2, appears on most trajectories.
Anyhow, this secondorder effect is most probably linked with the presence of the wall, as suggested by recent simulations by Narsimhan et al., where interacting red blood cells in the vicinity of a wall are studied Narsimhan, Zhao, and Shaqfeh (2013). They show that, for particles close enough to a wall, the relative lateral distance might decrease before interaction (see trajectories on their Fig. 15(b)), sometimes even leading to a completely different scenario of interaction involving swapping trajectories where particles do not cross. As shown by the authors, the presence of the bottom wall in the y direction induces the formation of a recirculation vortex behind the first particle, which is mainly responsible for the initial attraction. It is likely that when walls are also present in the z direction, as is the case in the experiment, the strength and extension of this recirculation are larger, leading to the stronger attraction observed in the first stage of experimental trajectories.
To sum it up, starting from comparable results for experiments and simulations, we have shown that a shift in the vorticity direction and the contribution of walls, both being inherent to the experiment, lead to a decrease of the repulsion, thus to some scattering in the experimental data, that all lie right below the ideal curve given by the simulations.
v.3 Deflection in the vorticity direction
As an extension to experimental results, the model also allows to investigate the effect of the interaction on the deflection .
In Fig. 9, we present the interaction curves and , for , the other parameters remaining the same as previously.
We find that there is a range of initial transverse positions for which the interaction leads to a transverse attraction between the vesicles, mostly in the vorticity direction. A similar phenomenon has been predicted for the interaction of capsules Lac and BarthèsBiesel (2008), but not for drops Loewenberg and Hinch (1997). An asymptotic study, for vesicles in the far field interacting regime, also predicts such an attraction Gires, Danker, and Misbah (2012). However, here the vesicles become close during the interaction, so a qualitative interpretation of the predicted attraction may involve the description of the flow of the thin liquid film between the two tanktreading membranes, as done for drops Loewenberg and Hinch (1997).
Vi Hydrodynamic diffusion
From the numerical study, one can expect to deduce results about the hydrodynamic diffusion properties of vesicle suspensions, in a regime where the solution is concentrated enough so that interaction effects are not negligible, but dilute enough so that pairwise interactions dominate over threebody interactions. We mostly study the case of selfdiffusion, a phenomenon related to the average transverse motion of a single vesicle. We also find that an estimation of the collective diffusion coefficients is not possible only considering twovesicle interactions, due to the long range of hydrodynamic interactions.
vi.1 Selfdiffusion
Theoretical background
We consider a homogeneous suspension of vesicles, described at a mesoscopic level by a volume fraction . For a given initial state, if this suspension is sheared by an imposed flow , a given vesicle will interact with the others and, as a result, will undergo a net displacement from its original streamline. In an unstructured semidilute suspension, the transverse motion of the vesicle is expected to be a random walk due to successive interactions with different vesicles. At long times, its meansquared displacement is described by the selfdiffusion coefficients , defined by
with , being an ensemble average over all possible initial states of the suspensionLoewenberg and Hinch (1997).
As detailed by Da Cunha and Hinch Da Cunha and Hinch (1996), assuming that only twovesicle interactions occur leads to the following expression for :
(8) 
with
(9) 
where is the deviation of a test vesicle in the laboratory frame after one interaction. For identical vesicles, . In the latter integral and from now on, all lengths are expressed in units.
Analysis of the formal convergence of the diffusion coefficient
As the hydrodynamic interaction between two vesicles slowly decreases, the question of the convergence of the previous integral arises. As all displacements are bounded, the convergence of the expression 9 is linked to the contribution of the integration domain . We analyze this contribution by using an asymptotic study of two interacting quasispherical vesicles remaining very distant from each other, that was recently proposed by Gires et al. Gires, Danker, and Misbah (2012). We first need to determine the domain in the space for which this asymptotic study is valid, a discussion that was not provided in the original paper. For the asymptotic study to be valid, vesicles must remain far enough along the whole trajectory, so that, at all times, . As , one could expect this criterion to be always satisfied. However, let us consider . If there was no interaction, both vesicles would flow with the same velocity. But, since the velocity field induced by one vesicle is radial, vesicle collision may occur, which is inconsistent with the asymptotic approach. These considerations hint to the fact that the asymptotic study may not be valid for .
In order to get a more accurate validity criterion, we assume the asymptotic study to be valid for all times, and check that the intervesicle distance remains large. We expect that this approach can be used as each vesicle is not in the vicinity of a bifurcation phenomenon, such as the transition between the tanktreading and vacillatingbreathing modes.
As detailed in Gires et al.Gires, Danker, and Misbah (2012), within this asymptotic approach with respect to the intervesicle distance, the trajectory of with respect to is of the form:
(10)  
and
(11) 
where , and are constants linked to the perturbation of the velocity field induced by the vesicles.
As the symmetry of the system does not depend on the reduced volume, we expect these scalings to be valid for vesicles of arbitrary deflation in the tanktreading regime, the dependency on the reduced volume being accounted for by the tensor .
As , the distance between the vesicles will remain large if initially large. However, as , problems may arise at small , as discussed earlier. In this case, the prevalent term in Eq. 10 is the term proportional to . If , this could lead to a minimal distance in the vorticity direction of the form , with . In order that the asymptotic approach remains valid starting with , we impose the condition that , where is a constant. This criterion can be achieved if , with .
(12)  
(13) 
It is clear from these expressions that the integral of Eq. 9, restricted to the region where the asymptotic expression is valid, is convergent.
As for the region of large and small with , where the asymptotic expression is not valid, since is bounded by its maximal value and the integral of on this region is finite, its contribution to the integral in Eq. 9 is bounded, and finally the whole integral is convergent.
Numerical determination of the selfdiffusion coefficient
We now estimate the value of (Eq. 9) that enters the expression of the diffusion coefficient (Eq. 8). For that purpose we need to run several simulations by starting with different initial position in the plane (which is the plane orthogonal to the flow direction) and determine by how much the initial relative positions and have varied (by amounts and ) after the two vesicles have interacted. We discretize the domain of initial values ) by considering the following domain size (in units of vesicle radius ). The discretized lattice of initial positions is shown in Fig. 10(a) with dark gray disks (blue online). Since the interaction is important only when the two vesicles are separated by about 2 or 3 radii, the lattice has a wide enough periodicity far away from , whereas in the vicinity further refinements are chosen in order to gain numerical precision. More precisely, the domain is decomposed into three regions , and , consisting in , and . The lighter gray disks (red online) in Fig. 10(a) show the final relative positions and . We note that in region C the effect is weak, while it becomes more and more pronounced in region and . The contributions of the integral involved in Eq. 9 on the different subdomains , and are then evaluated using a trapezoidal rule. The results are given in Table 2.



part A 
0.028  0.002 
part B 
0.003  0.004 
part C 
0.001  0.005 

We find that the contributions for are decreasing. As we proved the convergence of the expression, we expect the contributions of the remaining part of the plane to be at most of the order of the contribution of the subdomain C, and thus get the following estimation for :
(14) 
The uncertainty of is a rough estimate coming from a study of the sensitivity of the code to some numerical parameters, like a tension parameter used to preserve locally the area of the membrane.
For , we do not get decreasing contributions, due to the slow decrease of with when . A similar study has been presented by Zhao and Shaqfeh Zhao and Shaqfeh (2013), who calculated for and . Using the effective radius based on the surface as a length scale (, where is the vesicle membrane area) they find . With the same convention instead of our choice of radius based on the volume, we find , for and , which is a consistent result since lateral displacement increases with (Fig. 4). Zhao and Shaqfeh also estimated the value of , restricting to the integration domain : their value matches ours on the same region. However, the present study shows that restricting the integration to this domain is not sufficient to get a quantitative value of , due to the slow decrease of the attraction with for vesicles characterized by .
We are not aware of experimental measures of to which we could compare our estimation. On the basis of studies on suspensions of spheres, the assumption of considering only twovesicle interactions could be a good approximation up to volume fractions of about Larson (1999).
Discussion
From simulated trajectories, Loewenberg and Hinch Loewenberg and Hinch (1997) calculated and for pairs of drops as a function of viscosity ratio and capillary number. They evoke the scaling at long distance , which is similar to ours, to prove the convergence of the integral of Eq. 9. It appears that in the case of drops, restricting the integration domain to A+B is sufficient, for as for . For , was found to be around , depending on the capillary number, a result close to ours. A more quantitative comparison is precluded by the dependency with capillary number and the difference in nature between the elastic restoring forces involved in drops and vesicles. Interestingly, Loewenberg and Hinch Loewenberg and Hinch (1997) find , while we already find a result 3 times larger by integration over A+B+C. We can conclude that anisotropy in selfdiffusion is lower for vesicles than for drops. This weaker anisotropy is also stated by Lac and BarthèsBiesel in their study of capsules collisions, though and are not calculated Lac and BarthèsBiesel (2008).
vi.2 Downgradient diffusion
The collective diffusion property of a vesicle suspension can be modeled in the following way: we consider a suspension of vesicles with a linearly increasing concentration given by , sheared by an imposed flow . As a result of the hydrodynamic interactions between the vesicles, we expect a collective diffusion of the vesicles to appear, consisting of a transverse flux of vesicles. As for molecular diffusion due to thermal motion, is expected to be in the opposite direction of , of order . Thus, for , we expect that , with . As done by Da Cunha and Hinch in the case of rough spheres Da Cunha and Hinch (1996), we tried to estimate assuming only twovesicle interactions. This leads to an expression involving the integral of over the plane. However, as and the integral of over , with , is divergent, it turns out that the estimated expression is not convergent. A renormalization procedure, analogous to the one used by Batchelor Batchelor (1972); Batchelor and Green (1972), and followed by Wang et al.Wang, Mauri, and Acrivos (1998) in the case of the study of the hydrodynamic diffusion properties of a suspension of spheres, may lead to a convergent expression. It is hoped to investigate this matter further in a future work.
Vii Conclusion
We performed an experimental and numerical study of the trajectory deviations of identical vesicles interacting in shear flow. In experiments, restricted to pairs of vesicles in the same shear plane, the amplitude of the net displacement decreases quickly when the initial lateral distance increases and becomes negligible when this distance is larger than approximately two vesicle radii.
A simplified model based on the well established law for the lift of a vesicle near a wall was proposed, which allows to estimate quantitatively how the displacement should vary with the mechanical properties of the vesicles.
With no fitting parameter, the deviations are found to be in rather good agreement with our 3D simulations, even if smaller deviations are found experimentally. We found than the main part of this discrepancy can be due to differences between the experimental configuration and the ideal case of unbounded shear flow where the two vesicles would be perfectly coplanar. The effect of walls, recently highlighted by Narsimhan et al. Narsimhan, Zhao, and Shaqfeh (2013), would need to be quantified thoroughly in complementary experiments where our requirements of similar deflation within the pair of vesicles could be loosened for simplicity, since the effect of deflation has been characterized and shown to be weak. We also indicate that, according to partial results not shown here, the requirement of identical size within a pair may be released, as rescaling of the displacements by the average radius of two vesicles of different size lead to a similar curve for as a function of .
In addition, displacements in the vorticity direction were explored through the simulations and found to be about an order of magnitude lower than in the shear direction, with a range of initial distances leading to a weak attraction of vesicles.
Shearinduced diffusion coefficients can be obtained by a proper averaging of the net displacement over all initial configurations. The selfdiffusion, related to the random walk of vesicles in a suspension, can be quantified using a discrete integration over a relatively small domain for the diffusivity in the shear direction, and could be determined in the vorticity direction if a larger integration area was considered, due to the slower decrease of the amplitude of displacement in that direction. Note that this integration would not have been possible in 2D Li and Pozrikidis (2000) where the displacements would scale like instead of .
An estimation of the downgradient diffusivities as defined by Da Cunha and Hinch Da Cunha and Hinch (1996) was not possible due to the long range of hydrodynamic interactions, leading to a divergence of the integrals. In this case, the dilute limit assumption breaks down and one can no longer consider only pair interactions as is the case for rough spheres with short range interactions Da Cunha and Hinch (1996).
Appendix A Local conservation of the area
We present in this appendix results about the local conservation of the area of a vesicle, during a typical trajectory. The parameters chosen were . First, we plot in Fig. 11 the maximal relative variation of the area of the mesh faces, between two time steps, if they were advected by the full velocity field (in the simulation, the vertices are only advected by the normal component of the velocity field).
We find that this maximum can reach values higher than one. However, the proportion of faces corresponding to such values stays lower than , as shown in Fig. 12.
Acknowledgements.
The authors thank T. Biben (LPMCN, Lyon, France) for providing the initial simulation code and CNES (Centre National d’Etudes Spatiales) and ESA (European Space Agency) for financial support. P.Y.G. also acknowledges funding from the French Ministry of Higher Education and Research, and the CNRS (Centre National de la Recherche Scientifique).Footnotes
 thanks: both first authors equally contributed to the production of the data presented in this paper. P.Y. G. produced the theoretical and computational results and A. S. performed the experiments.
References
 G. Batchelor and J. Green, “The determination of the bulk stress in a suspension of spherical particles to order ,” J. Fluid Mech. 56, 401–427 (1972).
 A. Z. Zinchenko, “Effect of hydrodynamic interactions between the particles on the rheological properties of dilute emulsions,” J. Applied Math. Mech. 48, 198–206 (1984).
 H. L. Goldsmith, “Red cell motions and wall interactions in tube flow,” Fed. Proc. 30, 1578 (1971).
 H. L. Goldsmith and J. C. Marlow, “Flow behaviour of erythrocytes,” J. Colloid Interface Sci. 71, 383–407 (1979).
 W. Cha and R. L. Beissinger, “Evaluation of shearinduced particle diffusivity in red cell ghosts suspensions,” Korean J. Chem. Eng. 18, 479 (2001).
 J. M. Higgins, D. T. Eddington, S. N. Bhatia, and L. Mahadevan, “Statistical dynamics of flowing red blood cells by morphological image processing,” PLoS Computational Biology 5, e1000288 (2009).
 H. Zhao and E. S. G. Shaqfeh, “Shearinduced platelet margination in a microchannel,” Phys. Rev. E 83, 061924 (2011).
 M. H.Y. Tan, D.V. Le, and K.H. Chiam, “Hydrodynamic diffusion of a suspension of elastic capsules in bounded simple shear flow,” Soft Matter 8, 2243–2251 (2012).
 S. D. Hudson, “Wall migration and shearinduced diffusion of fluid droplets in emulsions,” Phys. Fluids 15, 1106–1113 (2003).
 R. Rusconi and H. A. Stone, “Shearinduced diffusion of platelike particles in microchannels,” Phys. Rev.Lett. 101, 254502 (2008).
 T. Podgorski, N. Callens, C. Minetti, G. Coupier, F. Dubois, and C. Misbah, “Dynamics of vesicle suspensions in shear flow between walls,” Microgravity Sci. Technol. 23, 263–270 (2011).
 X. Grandchamp, G. Coupier, A. Srivastav, C. Minetti, and T. Podgorski, “Lift and downgradient shearinduced diffusion in red blood cell suspensions,” Phys. Rev. Lett. 110, 108101 (2013).
 F. P. Bretherton, “The motion of rigid particles in shear flow at low reynolds number,” J. Fluid Mech. 14, 284–304 (1962).
 F. Da Cunha and E. Hinch, “Shearinduced dispersion in a dilute suspension of rough spheres,” J. Fluid Mech. 309, 211–223 (1996).
 F. Blanc, F. Peters, and E. Lemaire, “Experimental signature of the pair trajectories of rough spheres in the shearinduced microstructure in noncolloidal suspensions,” Phys. Rev. Lett. 107, 208302 (2011).
 L. Van Wijngaarden and D. J. Jeffrey, “Hydrodynamic interaction between gas bubbles in liquid,” J. Fluid Mech. 77, 27– 44 (1976).
 M. Loewenberg and E. Hinch, “Collision of two deformable drops in shear flow,” J. Fluid Mech. 338, 299 (1997).
 S. Guido and M. Simeone, “Binary collision of drops in simple shear flow by computerassisted video optical microscopy,” J. Fluid Mech. 357, 1–20 (1998).
 R. K. Singh and K. Sarkar, “Effects of viscosity ratio and three dimensional positioning on hydrodynamic interactions between two viscous drops in a shear flow at finite inertia,” Phys. Fluids 21, 103303 (2009).
 D.V. Le and K.H. Chiam, “Hydrodynamic interaction between two nonspherical capsules in shear flow,” Phys. Rev. E 84, 056322 (2011).
 E. Lac, A. Morel, and D. BarthèsBiesel, “Hydrodynamic interaction between two identical capsules in simple shear flow,” J. Fluid. Mech. 573, 149 (2007).
 E. Lac and D. BarthèsBiesel, “Pairwise interaction of capsules in simple shear flow: threedimensional effects,” Phys. fluids 20, 040801 (2008).
 P. Pranay, S. G. Anekal, J. P. HernandezOrtiz, and M. D. Graham, “Pair collisions of fluidfilled elastic capsules in shear flow: Effects of membrane properties and polymer additives,” Phys. Fluids 22 (2010).
 T. Omori, T. Ishikawa, Y. Imai, and T. Yamaguchi, “Membrane tension of red blood cells pairwisely interacting in simple shear flow,” J. Biomech. 46, 548 (2013).
 W. Helfrich, “Elastic Properties of Lipid Bilayers: Theory and Possible Experiments,” Z. Naturforsch. 28c, 1–12 (1973).
 J. Deschamps, V. Kantsler, and V. Steinberg, “Phase Diagram of Single Vesicle Dynamical States in Shear Flow,” Phys. Rev. Lett. 102 (2009).
 A. Farutin, T. Biben, and C. Misbah, “Analytical progress in the theory of vesicles under linear flow,” Phys. Rev. E 81, 061904 (2010).
 T. Biben, A. Farutin, and C. Misbah, “Threedimensional vesicles under shear flow: Numerical study of dynamics and phase diagram,” Phys. Rev. E 83, 031921 (2011).
 G. Danker and C. Misbah, “Rheology of a dilute suspension of vesicles,” Phys. Rev. Lett. 98, 088104 (2007).
 V. Kantsler, E. Segre, and V. Steinberg, “Dynamics of interacting vesicles and rheology of vesicle suspension in shear flow,” Europhys. Lett. 82, 58005 (2008).
 V. Vitkova, M.A. Mader, B. Polack, C. Misbah, and T. Podgorski, “Micromacro link in rheology of erythrocyte and vesicle suspensions,” Biophys. J. 95, 33 (2008).
 M. Levant, J. Deschamps, E. Afik, and V. Steinberg, “Characteristic spatial scale of vesicle pair interactions in a plane linear flow,” Phys. Rev. E 85, 056306 (2012).
 P.Y. Gires, G. Danker, and C. Misbah, “Hydrodynamic interactions between two vesicles in a linear shear flow: asymptotic study,” Phys. Rev. E 86, 011408 (2012).
 A. Farutin and C. Misbah, “Analytical and numerical study of three main migration laws for vesicles under flow,” Phys. Rev. Lett. 110, 108104 (2013).
 H. Zhao and E. S. G. Shaqfeh, “The dynamics of a nondilute vesicle suspension in a simple shear flow,” J. Fluid Mech. 725, 709–731 (2013).
 V. Narsimhan, H. Zhao, and E. S. G. Shaqfeh, “Shearinduced particle migration and margination in a cellular suspension,” Phys. Fluids 24, 011902 (2012).
 A. Kumar and M. D. Graham, “Mechanism of margination in confined flows of blood and other multicomponent suspensions,” Phys. Rev. Lett. 109, 108102 (2012).
 D. A. Fedosov, J. Fornleitner, and G. Gompper, “Margination of white blood cells in microcapillary flow,” Phys. Rev. Lett. 108, 028104 (2012).
 H. L. Goldsmith and S. Spain, “Margination of leukocytes in blood flow through small tubes,” Microvas. Res. 27, 204 – 222 (1984).
 A. W. Tilles and E. C. Eckstein, “The nearwall excess of plateletsized particles in blood flow: Its dependence on hematocrit and wall shear rate,” Microvasc. Res. 33, 211 – 223 (1987).
 W. S. Uijttewaal, E. J. Nijhof, P. J. Bronkhorst, E. Den Hartog, and R. M. Heethaar, “Nearwall excess of platelets induced by lateral migration of erythrocytes in flowing blood,” Am. J. Physiol. Heart Circ. Physiol. 264, H1239–H1244 (1993).
 C. Yeh and E. Eckstein, “Transient lateral transport of plateletsized particles in flowing blood suspensions,” Biophys. J 66, 1706 – 1716 (1994).
 M. I. Angelova, S. Soleau, P. Meleard, J. F. Faucon, and P. Bothorel, “Preparation of giant vesicles by external ac electric fields. kinetics and applications,” Progr. Colloid. Polym. Sci 89, 127–131 (1992).
 G. Coupier, B. Kaoui, T. Podgorski, and C. Misbah, “Noninertial lateral migration of vesicles in bounded Poiseuille flow,” Phys. Fluids 20, 111702 (2008).
 M. Abkarian, C. Lartigue, and A. Viallat, “Tank treading and unbinding of deformable vesicles in shear flow: Determination of the lift force,” Phys. Rev. Lett. 88, 068103 (2002).
 N. Callens, C. Minetti, G. Coupier, M. Mader, F. Dubois, C. Misbah, and T. Podgorski, “Hydrodynamic lift of vesicles under shear flow in microgravity,” Europhys. Lett. 83, 24002 (2008).
 P. Olla, “The lift on a tanktreading ellipsoidal cell in a shear flow,” J. Phys. II France 7, 1533–1540 (1997).
 S. Sukumaran and U. Seifert, “Influence of shear flow on vesicles near a wall: a numerical study,” Phys. Rev. E 64, 011916 (2001).
 P. M. Vlahovska and R. Serral Gracia, “Dynamics of a viscous vesicle in linear flows,” Phys. Rev. E 75, 016313 (2007).
 O. Zhongcan and W. Helfrich, “Bending energy of vesicle membranes : general expressions for the first, second and third varaition of the shape energy and applications to spheres and cylinders,” Phys. Rev. A 39, 1280 (1989).
 K. H. de Haas, C. Blom, D. van den Ende, M. H. G. Duits, B. Haveman, and J. Mellema, “Rheological behavior of a dispersion of small lipid bilayer vesicles,” Langmuir 13, 6658 (1997).
 C. Pozrikidis, Boundary integral and singularity methods for linearized viscous flow (Cambridge University Press, Cambridge, 1992).
 B. Kaoui, G. Coupier, C. Misbah, and T. Podgorski, ‘‘Lateral migration of vesicles in microchannels: effects of walls and shear gradient,” Houille Blanche 5, 112–119 (2009).
 V. Narsimhan, H. Zhao, and E. S. G. Shaqfeh, “Coarsegrained theory to predict the concentration distribution of red blood cells in wallbounded couette flow at zero reynolds number,” Phys. Fluids 25, 061901 (2013).
 R. Larson, The structure and rheology of complex fluids (Oxford University Press, Oxford, 1999).
 G. Batchelor, “Sedimentation in a dilute suspension of spheres,” J. Fluid Mech. 52, 245–268 (1972).
 Y. Wang, R. Mauri, and A. Acrivos, “Transverse shearinduced gradient diffusion in a dilute suspension of spheres,” J. Fluid Mech. 357, 279–287 (1998).
 X. Li and C. Pozrikidis, “Wallbounded shear flow and channel flow of suspensions of liquid drops,” Int. J. Mult. Flow 26, 1247 (2000).