Flow-induced clustering and alignment of vesicles and red blood cells in microcapillaries

Abstract

The recent development of microfluidic devices allows the investigation and manipulation of individual liquid microdroplets, capsules, and cells. The collective behavior of several red blood cells (RBCs) or microcapsules in narrow capillaries determines their flow-induced morphology, arrangement, and effective viscosity. Of fundamental interest here is the relation between the flow behavior and the elasticity and deformability of these objects, their long-range hydrodynamic interactions in microchannels, and thermal membrane undulations. We study these mechanisms in an in silico model, which combines a particle-based mesoscale simulation technique for the fluid hydrodynamics with a triangulated-membrane model. The 2 essential control parameters are the volume fraction of RBCs (the tube hematocrit, H(T)), and the flow velocity. Our simulations show that already at very low H(T), the deformability of RBCs implies a flow-induced cluster formation above a threshold flow velocity. At higher H(T) values, we predict 3 distinct phases: one consisting of disordered biconcave-disk-shaped RBCs, another with parachute-shaped RBCs aligned in a single file, and a third with slipper-shaped RBCs arranged as 2 parallel interdigitated rows. The deformation-mediated clustering and the arrangements of RBCs and microcapsules are relevant for many potential applications in physics, biology, and medicine, such as blood diagnosis and cell sorting in microfluidic devices.

Fig. 1.

Conformations of elastic vesicles for various densities and mean flow velocities v*₀, for a single vesicle in the simulation cylinder (n_ves = 1) with various lengths L*_ves. The hematocrit is given by H_T = n_vesV_ves/πL_zR_cap² = 0.28/L*_ves, where V_ves is the volume of single vesicle. (A) Average asphericity α, which measures the deviation from a spherical shape, as a function of the mean flow velocity v*₀. Here, the asphericity is given by α = [ (λ₁ − λ₂)² + (λ₂ − λ₃)² + (λ₃ − λ₁)²]/2R_g⁴, with the eigenvalues λ₁ ≤ λ₂ ≤ λ₃ of the gyration tensor and the squared radius of gyration R_g² = λ₁ + λ₂ + λ₃. RBCs transits from the discocyte (with α ≃ 0.15) to the parachute shape (with α ≲ 0.05). Simulation snapshots show the parachute (Lower Left) and bowl (Upper Right) shapes for v*₀ = 7.7 with L*_ves = 2.25 and L*_ves = 1, respectively. (B) Average inclination angle θ for a discocyte vesicle for low flow velocities v*₀ < v*_c, as a function of the vesicle distance L_ves. The inclination angle θ measures the deviation of the vesicle symmetry axis (determined by the eigenvector associated with the minimum eigenvalue of the gyration tensor) with the flow direction (z axis).

Fig. 2.

Streamlines (blue) and velocity field (red arrows) of the flow between 2 vesicles, in the comoving frame with the vesicle velocity at L*_ves = 2 (H_T = 0.14) and v*₀ = 10. A sliced snapshot (black line) of the vesicle is also shown. A flow vortex (bolus) is seen between vesicles.

Fig. 3.

Snapshots of n_ves = 6 elastic vesicles in the simulation channel. (A) Disordered-discocyte phase for L*_ves = 0.875 (H_T = 0.32) and v*₀ = 2.5, where vesicles appear usually as discocytes; the degree of shape deformation increases with increasing v*₀. (B) Aligned-parachute phase for L*_ves = 0.875 (H_T = 0.32) and v*₀ = 10. (C) Zigzag-slipper phase for L*_ves = 0.75 (H_T = 0.37) and v*₀ = 10.

Fig. 4.

Phase diagram and pressure drop of dense RBC suspensions in flow. (A) Phase behavior as a function of average vesicle distance L*_ves and mean flow velocity v*₀ for n_ves = 6. The hematocrit varies between H_T = 0.22 and H_T = 0.45, because H_T = 0.28/ L*_ves. Symbols represent the disordered-discocyte (*), aligned-parachute (■), and zigzag-slipper (●) phases, respectively. The phase boundaries are drawn to guide the eye. (B) Pressure drop ΔP*_drp per vesicle for the aligned-parachute phase (simulations with n_ves = 1) and the zigzag-slipper phase (simulations with n_ves = 6) at the same volume fraction (L*_ves = 0.75, corresponding to H_T = 0.37). The pressure drop is given by ΔP*_drp = ΔP_drpR_cap/η₀n_vesv_m = 8 (v₀ − v_m)L_z/n_vev_mR_cap, where v_m is mean fluid velocity in the presence of vesicles.

Fig. 5.

Sequential snapshots of 6 elastic vesicles in dilute suspension (H_T = 0.084) at v*₀ = 7.7.

Fig. 6.

Pair distribution function G(z*_nb) (A) and cluster-size probabilities P(n_cl) (B) in dilute suspension with H_T = 0.084. The black and red solid lines represent data v*₀ = 7.7 and 10, respectively. Results are obtained from simulations with n_ves = 6 vesicles. The dotted line in A represents the pair distribution function of only the nearest-neighbor vesicles at v*₀ = 10. In B, vesicles closer than 2.5 R_cap are defined to belong to the same cluster.