Long-range hydrodynamic effect due to a single vesicle in linear flow

A long-range hydrodynamic interaction between objects immersed in a fluid is recognized as one of the main ingredients contributing to the effective viscosity of suspensions of meso-scale particles. The nonlocal nature of the hydrodynamic interaction is a key physical factor that makes this problem so rich, but at the same time so difficult. The velocity distribution in the fluid near an immersed object is affected by the presence of another object at a finite distance. It results in modification of the fluid stresses at the surface of the first object and correspondingly of its translational and rotational motion. The effect critically depends on the nature and elastic properties of the immersed objects.

In this respect, central questions are: what is the characteristic spatial scale of a long-range hydrodynamic interaction and whether and how does it depend on the nature and elastic properties of objects immersed in a fluid? In the case of solid spherical particles in a fluid subjected to a linear flow field, this problem was first addressed theoretically in ref. [1] and then studied in many papers both theoretically and numerically (see, e.g., [2,3] and references therein). In the case of vesicles, the long-range hydrodynamic vesicle pair interaction and the corresponding spatial scale was investigated recently both experimentally [4,5] as well as numerically [6]. In this paper we address the problem directly related to the long-range hydrodynamic interaction of a pair of vesicles in a planar linear flow, namely: How strong are the perturbations of a linear velocity field caused by a vesicle and how far away from it are they observed? What is the scaling of the perturbations of the linear velocity field with the distance from the vesicle? We explore these questions both experimentally and numerically and compare the results.

A unilamellar liposome vesicle is a drop of fluid encapsulated in a lipid bi-layer membrane, suspended in either the same or different fluid. The impermeability and the inextensibility of the lipid membrane dictate the conservation of both volume and surface area of the vesicle. Vesicles undergo a tank-treading motion in a planar linear flow at internal-to-external viscosity contrast $\lambda=\eta_{\textit{in}}/\eta_{\textit{out}}= 1$, vorticity-to-strain ratio $\omega/s<2$, and excess area $\Delta = S/R_0^2 -4\pi <1$ [7]. Here $R_0=\left(V\bigl/\frac{4\pi}{3}\right)^{^1\!/_3}$ is the vesicle effective radius, V is the vesicle volume, and S is the vesicle surface area; ω and s are the vorticity and strain rates, defined following ref. [8]. At a fixed viscosity contrast, the inclination angle θ decreases with increasing Δ and $\omega/s$ [8–11].

The details of the experimental methods were presented in ref. [5]. Below we provide a brief overview of the velocity field measurements in the presence of a vesicle, which are in the focus in this work. The experiments were conducted in a microfluidic four-roll mill apparatus [7,12]. The key component of this apparatus is the dynamical trap, which allows long observation times compared to the vesicles orbit period in the flow. The flow was driven by gravity and the control parameter $\omega/s$ was set via setting the pressure drop. Time-lapse velocity fields were inferred by means of 2D micro Particle Image Velocimetry (μPIV) at the mid-plane of the dynamical trap. The aqueous vesicles suspension was seeded with fluorescent particles, allowing imaging of the flow in the trap. The images were then filtered by a Laplace filter using Gaussian second derivatives [13], and processed using Gpiv [14] at interrogation windows of 32 × 32 pixels (corresponding to about $8.6\times8.6\ {\mu}\text{m}^2$) with 50% overlap [5].

The 2D numerical model used here was introduced in ref. [6]. A single vesicle with membrane bending rigidity κ, initial length L₀, initial area A₀, and $R_0 = L_0/2\pi$ being the vesicle effective radius, is studied in a wall-bounded shear flow with a shear rate $\dot{\gamma} = \omega/2 = s/2$, which is a specific case of a planar linear flow. Periodic boundary conditions are used in the x-direction; the upper and the lower walls slide along this direction at velocities $\pm u_{\textit{wall}}$, respectively, resulting in a linear flow profile $(V_x, V_y) = (\dot{\gamma}y,0)$, with shear rate $\dot{\gamma} = 2u_{\textit{wall}}/H_y$, H_y being the flow cell width. The system size is $H_x\times H_y =18.85R_0\times 8.12R_0$. The steady velocity field without a vesicle is taken in a frame of size $L_x\times L_y = 10.47R_0\times H_y$ centred at $(H_x/2,H_y/2)$. When the vesicle is embedded in the solvent, the steady velocity pattern is recorded in a frame of size $L_x \times L_y$ centred at $(x_{cm},H_y/2)$, where $(x_{cm}, y_{cm})$ denotes the position of the vesicle centre of mass. The bending rigidity is set to $\kappa=6.58k_B TR_0$, where k_BT is the thermal energy. This value of κ gives rise to a similar amplitude of undulation modes as for lipid bilayer membranes in 3D (where $\kappa\simeq 10k_BT$). The reduced shear rate is $\dot{\gamma}^* = \dot{\gamma}\tau_c=3$, where $\tau_c = \eta_{\textit{out}} R_0^3/\kappa$ is the relaxation time of the vesicle, so that the Reynolds number is $Re = 0.08$ allowing to neglect inertial effects. The viscosity contrast is set to $\lambda = 1$.

After an initial transient, vesicles reach a tank-treading steady state while diffusing, due to thermal noise, across the centre line of the channel. Values for the reduced area $A^* = A_0/\pi R_0^2$, related to the excess length $\Delta_L$ by $A^*=(1+\Delta _L/2\pi)^{-2}$, were chosen to be comparable with the experimental ones: At the initial simulation time t₀ it is $A^*(t_0) = 1$, corresponding to a nearly spherical vesicle, which in the flow settles at $A^*=0.984$, and $A^*(t_0) = 0.89$ for the deflated one, reaching $A^*=0.885$ in the steady flow. The time averages of the measured vesicle inclination angles are $\langle \theta \rangle = (0.226\pi, 0.162\pi)\ \text{rad}$ with standard deviations $\sigma_\theta=(0.048\pi, 0.025\pi)\ \text{rad}$, for $A^* = 0.984, 0.885$, respectively.

A few notes related to the numerical conservation of vesicle length and area are worth adding. In the numerical model used here, the membrane is modelled using beads successively connected by elastic bonds to form a closed ring (see ref. [15] for details). The bond stiffness is the control parameter which maintains the membrane length L in proximity to its initial value L₀, both locally and globally. To respect the physical constraints, the vesicle area is kept close to the initial value A₀ by introducing a quadratic potential depending on a compression modulus [15]. Here, the deviations from the values A₀ and L₀ are below 0.6% and 0.5%, respectively. However, the values of bond stiffness and compression modulus cannot be arbitrarily increased due to numerical instability. This limitation is well known and appears in other membrane dynamics modelling, erythrocytes for example [16], in particular at larger shear stresses.

As a quantitative measure for the strength of the vesicle effect on the velocity field, we chose to study the decay of the disturbance field magnitude $\| \delta \vec{V} \|$ with the distance from the vesicle. To this end, the following analysis procedure was applied to both the experimental and the numerical disturbance fields: A ridge detection algorithm, explained in detail as part of ref. [17] (the code is available online), based on a differential geometric descriptor [18], was applied to the $\|\delta \vec{V}\|$ landscape. Position coordinates residing on a ridge were sorted based on their location in the frame of reference whose origin is at the vesicle centre of mass. Then, linear regression was applied to find the best fitting lines approximating the principal ridge directions within a distance of up to 2.5 vesicle effective radii. To analyse $\| \delta \vec{V} \|$ along these directions, for every sampled x coordinate within the data, a corresponding y coordinate from the data grid was chosen by minimal distance from the best fitting line.

We characterise and quantify the effect of a single vesicle on the surrounding velocity field which drives it, via the analysis of the disturbance field $\delta \vec{V}$ it induces; here $\delta \vec{V}$ is the difference between the velocity field in the absence and presence of a vesicle. Figure 1 shows the resulting disturbance fields, where the panels on the left show the experimental results and the right ones show the numerical simulations; quasi-spherical vesicles are presented in the top part, while the deflated ones can be found at the bottom. This ordering is respected throughout figs. 1–3.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

First we address the case of a quasi-spherical vesicle having a very small excess surface area Δ (excess length $\Delta_L$ in the 2D numerical simulations). In the experiments, the vesicle was located in the vicinity of the stagnation point of the dynamical trap and moving on a small orbit in a planar linear flow of ω = 0.46 s⁻¹ and $\omega/s=1.8$ [5]. The second figure in ref. [5] shows the reference field (top panel), from which the measured velocity field in the presence of a vesicle (middle panel) was subtracted, resulting in the disturbance field (bottom panel); an adapted plot of the disturbance field is provided here in fig. 1(a). It allows to study the back-reaction of the vesicle on the flow and to estimate the effective distance of this back-reaction. Note that the characteristic disturbance field magnitude is smaller than the velocities in the reference field by more than an order of magnitude (see fig. 2 in [5]). Two pronounced qualitative features of the disturbance velocity field are evident already at this stage: a) the flipping of strain directions in the vicinity of the vesicle, which is significant up to a few vesicle effective radii along the strain main axes, and b) the presence of four vortices generated in the vicinity of the vesicle. The former can be viewed as a result of the tank-treading motion of the membrane, which moves at a uniform tangential velocity throughout the membrane, thus transferring momentum from the shear direction to the perpendicular one. The uniformity of the membrane tangential velocity is due to the vesicle's volume and area conservation, which opposes extension and compression. These main features of the vesicle back-reaction on a shear flow are captured in the 2D numerical simulations, as can be seen from the time-average of the disturbance field around a nearly spherical vesicle provided in fig. 1(b). They were also found and presented in a snapshot of the velocity field around a nearly spherical vesicle in fig. 1 of ref. [15], as well as in snapshots of the disturbance field generated by a solid spherical particle in a shear flow obtained in 3D numerical simulations [3]. It seems that dimensionality does not play a key role in determining the flow features of interest to us in this work. A similar disturbance field and streamlines caused by a torque-free solid sphere in a simple shear flow were presented in fig. 9 of ref. [19], where the results of numerical simulations for $Re=0,1,10$ were shown. However, one notices that in the Re = 0 case the streamlines are straight and the vortices are absent [19]. At Re = 1 the streamlines are slightly curved and two vortices in the compressed direction become visible in sharp contrast to the perturbed velocity field caused by a solid spherical particle in ref. [3] and by a spherical vesicle at $Re \ll 1$ in both the experimental results and numerical simulations reported here (see fig. 1(a), (b)), where four vortices are evident.

Next, we proceed with the experimental results of a deflated vesicle moving on a small orbit in a plane linear flow with ω = 0.46 s⁻¹ and $\omega/s=1.4$; due to the comparable size of the vesicle and its orbit, we performed partial time averaging including only the sections highlighted by (orange) thicker lines in fig. 1(c), adapted from fig. 3 of ref. [5]. The corresponding numerical results are presented in fig. 1(d). These are qualitatively consistent with the experimental ones, fig. 1(c). A quantitative comparison to the nearly spherical case, as well as between the experimental and the numerical results, can be found in what follows.

In order to quantitatively analyse the disturbance field, we focus on its magnitude landscape $\| \delta \vec{V} \|$, plotted in fig. 2 as heat maps. The main features of these landscapes were extracted using the application of a ridge detection algorithm [17], as explained in the previous section. The arrow heads in fig. 2 designate the detected ridge loci, pointing perpendicular to the ridge (they are directed along the eigenvectors of the local Hessian matrix, corresponding to the least principal curvature); colours correspond to the grouping by stretching/compression directions.

We next wanted to study the functional dependence of the disturbance magnitude decay on the distance from a vesicle. The principal directions of the ridges were extracted using a linear approximation at the vicinity of the vesicle within the dimensionless distance from the vesicle center of mass of $d / R_0 \lesssim 2.5$ (corresponding ridge loci are circled in fig. 2, the linear regression results are plotted as dashed lines), and extrapolated based on these best fitting lines farther on up to $d/R_0\approx 6$. Based on the experimental data for the vesicles at $\lambda=1$, the angles formed by the approximating lines with respect to the flow direction are $\varphi_\pm \approx (0.237\pi, - 0.254\pi)\ \text{rad}$ in the nearly spherical case, and $\varphi_\pm \approx (0.166\pi, -0.325\pi)\ \text{rad}$ in the deflated case. The pronounced difference between the deflated vesicle and the nearly spherical one is in the deviation of the angles $\varphi_\pm$ from the strain principal directions: In the case of the deflated vesicle we study here, the stretching direction deviates from the strain direction by about $\pi/12\ \text{rad}$. The angles resulting from the analysis of the numerical data are $\varphi_\pm \approx (+0.239\pi, - 0.270\pi)\ \text{rad}$ and $\varphi_\pm \approx (+0.186\pi, -0.311\pi)\ \text{rad}$ for reduced areas $A^*=0.984$ and 0.885, respectively. Note that the principal directions are nearly orthogonal to each other, that is $\varphi_+ -\varphi_- \approx \pi/2\ \text{rad}$, and that the positive angle $\varphi_+$ is consistent with the average inclination angle $\langle \theta \rangle$ for both values of the reduced area.

The values of the reduced disturbance field magnitude $\|\delta\vec{V}\| /\dot{\gamma} R_0$ in the vicinity of the regression lines are plotted in fig. 3 as a function of the reduced distance $d/R_0$ from the vesicle centre of mass. The insets in fig. 3 present these data sets on a log-log scale, indicating that there is an intermediate regime where the decay can be approximated by a power law scaling $\propto (d/R_0)^{-3/2}$ to a reasonable degree. The experimental data for the deflated vesicle exhibit a larger scatter, yet show no contradiction to the power law decay with the same exponent. This larger scatter is attributed mainly to the averaging over a smaller data set, with respect to the one of the quasi-spherical vesicle. In the deflated case, the vesicle orbit was comparable to its size, hence averaging was taken over a sub-sample where the vesicle spent the larger fraction of time.

Finally, we wanted to assert that the numerical approach employed in this work captures the expected spherical solid-body behaviour at high enough values of viscosity contrast. Analysing the numerical results of a quasi-spherical vesicle at $\lambda=20$ resulted in the following angles for the principal directions: $\varphi_\pm \approx (+0.250\pi, - 0.256\pi)\ \text{rad}$. Indeed, this is close to what is expected for a solid sphere. The disturbance magnitude decay with distance, presented in fig. 4 for comparison, exhibits a scaling behaviour similar to the one previously described.

Zoom In Zoom Out Reset image size

In this work we have studied the back-reaction of a single vesicle, undergoing tank-treading motion, on the driving planar linear velocity field which surrounds it. Here we focused on the disturbance field induced by the vesicle, which results from the subtraction of the velocity field when the vesicle is present, from the reference field without it. The approach presented was applied to experimental results as well as to numerical simulations at comparable vesicle and flow configurations. The comparison between the wet experiment and the computerised one reveals an impressive agreement.

The landscape of the disturbance field magnitude was analysed to achieve quantitative results. We tracked ridges in this landscape and found two principal directions of vesicle effect, opposing the strain compression and stretching directions. The two principal directions are nearly perpendicular to each other and approximately follow the average tilt angle of the vesicle.

The decay of the disturbance field magnitude with the distance d from the vesicle along these principal directions reveals an intermediate regime which can be characterised as a power law decay $\|\delta \vec{V} \| /\dot{\gamma}R_0 \propto (d/R_0)^{-3/2}$. This result seems to be supported by our data for a nearly spherical vesicle, as well as for a moderately deflated one. The reported scaling for the decay is different from the theoretical prediction for long-range hydrodynamic effect at large distances $\sim 1/d +1/d^2 +O(1/d^3)$ at $d/R_0 \gg 1$, resulting from a membrane in an external shear flow and obtained as a special solution of the inhomogeneous Stokes equation for the velocity field [20].

The effect of a single vesicle undergoing tank-treading motion in a planar linear flow remains significant up to a distance of about 4 effective vesicle radii. This value is consistent with the separation distance between two tank-treading vesicles in the same flow, at which their dynamics is altered [5]. Unfortunately, an attempt to reproduce numerically the experimental result on the dependence of the root-mean-squared fluctuations of the deflated vesicle inclination angle as a function of distance from a spherical one, was not successful so far and resulted in very scattered data. Nevertheless, the results we found in numerical simulations showed no contradiction to the scaling found experimentally [5].

To summarise, we found an indication for the long-range hydrodynamic effects of a single vesicle on the velocity field, which are indeed present in a planar linear flow; this conclusion is supported both by experimental results and numerical simulations. Therefore, this calls for the incorporation of such long-range hydrodynamic effects into future models for improved predictions of the rheology of a vesicle suspension.

This work was partially supported by grants from Israel Science Foundation and Volkswagen Foundation through the Lower Saxony Ministry of Science and Culture Cooperation (Germany). Part of this research was performed while the first-listed author (EA) was a postdoc in Prof. David Harel's group at the Weizmann Institute of Science.