Crossover from tumbling to tank-treading-like motion in dense simulated suspensions of red blood cells

Abstract

Via computer simulations, we provide evidence that the shear rate induced red blood cell tumbling-to-tank-treading transition also occurs at quite high volume fractions, where collective effects are important. The transition takes place as the ratio of effective suspension stress to the characteristic cell membrane stress exceeds a certain value and does not explicitly depend on volume fraction or cell deformability. This value coincides with that for a transition from an orientationally less ordered to a highly ordered phase. The average cell deformation does not show any signature of transition, but rather follows a simple scaling law independent of volume fraction.

---

1 Introduction

Despite the large interest in a better understanding of the circulatory system and related diseases, there are still many open issues regarding the microscopic mechanisms which ultimately determine the rheology of suspensions of red blood cells (RBCs), the major particulate constituent of blood. Focusing on purely mechanical aspects, an RBC can be considered as a thin and flexible incompressible biconcave membrane which encloses an internal fluid (haemoglobin solution) and resists shearing and bending.

Depending on the ambient shear rate [small gamma, Greek, dot above], viscosity contrast (ratio between internal and external fluid viscosities), membrane viscoelasticity and other parameters, one observes tumbling, swinging, or tank-treading motion of isolated RBCs.^1–4 While in the case of single vesicles the dynamical phase space has been investigated thoroughly,^5–11 there are only a few studies on dense suspensions of deformable particles or RBCs.^12–15 These studies investigate the influence of concentration, deformability and aggregation tendency on suspension rheology.

In the present work, we focus on a dynamic phenomenon in dense RBC suspensions. Via computer simulations, we provide evidence for the transition from a tumbling to a tank-treading-like state (henceforth called TB–TT transition) upon increasing the capillary number Ca. The capillary number is defined as the ratio of fluid stress to the characteristic membrane stress, where the latter is related to the shear elasticity of the cell. Transcending previous studies, we find that the TB–TT transition occurs not only in the well-known case of a dilute suspension, but up to haematocrit values as high as Ht = 65% (haematocrit, or volume fraction, is the ratio of the volume occupied by RBCs to the total suspension volume). It is shown that—for all studied shear rates, cell deformabilities, and volume fractions—this transition is characterised by an effective capillary number Ca* (ratio between effective suspension stress and the characteristic membrane stress) rather than by the bare capillary number Ca. A detailed analysis of the average RBC inclination angle [small theta, Greek, macron] (a measure of average cell orientation) and the corresponding order parameter Q_> is also provided. When plotted versus Ca*, all values for Q_> and [small theta, Greek, macron] tend towards a master curve as Ca* exceeds a certain threshold value Ca^*_cr which, remarkably, coincides with that for the onset of the TB–TT transition for a single cell.

Our results provide evidence that the change from tumbling to tank-treading-like dynamics is accompanied by a phase transition from weak to strong orientational ordering of the cells. Interestingly, the average deformation of a cell, quantified by the Taylor deformation parameter D_a, changes continuously in the investigated parameter range, without any signature of the observed transition.

The article is organised as follows. The setup of the simulations is discussed in Section 2. Section 3 contains the results and discussions of the TB–TT transition (Section 3.1), the orientational order (Section 3.2), and the RBC deformation (Section 3.3). We conclude our work and provide an outlook in Section 4.

2 Simulation setup

2.1 Numerical model

A hybrid immersed-boundary-lattice-Boltzmann-finite-element model has been used.^16,17 The liquid phase is fully considered and solved by the lattice-Boltzmann method,^18,19 while the fluid–particle interaction is realised by the immersed-boundary method.²⁰ The red blood cells are represented by a finite-element mesh consisting of 1620 triangular facets. Physically, the cell membrane is made of a lipid bilayer (being essentially an incompressible 2D fluid), giving rise to a finite resistance against bending and changes of surface area, and a cytoskeleton, leading to a finite shear resistance.^21–23 The ensuing effects on the collective and rheological properties of the RBC suspension can be captured by assuming a total membrane energy of the form

E = E_S + E_B + E_A + E_V, (1)

wheredescribes the energy penalty against shear and area dilation,is the energy change due to bending, whileandpenalise changes of membrane surface area A and cell volume V over their undeformed counterparts, A⁽⁰⁾ and V⁽⁰⁾, respectively. Eqn (2) represents Skalak's energy density²⁴ with κ_S and κ_α being the shear and area resistance, respectively. I₁ and I₂ denote the in-plane strain invariants, which are related to the local membrane deformation tensor (see Krüger et al.¹⁶ for more details). Note that both the parameters κ_α and κ_A penalise changes in the membrane surface area (physically, κ_α is related to the cytoskeleton, while κ_A arises from the incompressibility of the liquid bilayer). The bending energy in eqn (3) has the classical Helfrich form,²⁵ of which a simplified versionis employed in the simulations.²⁶ Here, κ_B is the bending modulus and the sum runs over all pairs of neighbouring facets of the tessellated RBC surface with mutual equilibrium normal-to-normal angles θ^eq_ij. We remark that, in the present case, the biconcave shape of the RBC is not a result of the minimization of the bending energy subject to the constraint of fixed volume,²³ but rather specified as input to the simulation via the construction of the membrane mesh (and ensured by the presence of the θ^eq_ij in eqn (6)).

In our simulations, a fixed number (494, 635, 776, and 917) of neutrally buoyant RBCs with equilibrium shape of a biconcave disk of radius r = 9 (all quantities in lattice units) and equal internal and external viscosities have been distributed throughout a fixed fluid volume (L_x × L_y × L_z = 100 × 100 × 160) resulting in volume fractions Ht = 35, 45, 55, and 65%, respectively. To improve numerical stability and avoiding particle overlap, a repulsive force inspired by Feng and Michaelides²⁷ is introduced. It is zero for distances larger than one lattice constant and behaves like 1/r² for smaller distances. The effect of the membrane energy, eqn (1), is realised in terms of local forces derived from the principle of virtual work.¹⁶

The shear flow (periodic in x- and y-directions) has been achieved by moving the bottom and top walls at z = 0 and z = 160 in opposite directions along the x-axis with velocities ±u_w, giving rise to an average shear rate of [small gamma, Greek, dot above] = 2u_w/L_z. Inertia is neglected. One layer of RBCs has been stuck to each wall to avoid slip. A typical simulation snapshot illustrating the geometry is shown in Fig. 1(a).

Fig. 1 (a) Snapshot of a simulation with average haematocrit (volume fraction) Ht = 55% and capillary number Ca = 0.011. The bottom (top) wall is located at z = 0 (z = L_z = 160) and moves to the left (right) with velocity ±u_w (black arrows). (b) Layer-resolved haematocrit, Ht, shear stress, σ_xz, and flow velocity, u_x, for the same control parameters as in panel (a). Each curve is obtained as an average over steady state and independent runs. In the central region (grey, from L_z/4 to 3L_z/4), the volume fraction is nearly constant, the velocity profile is linear (red solid line), and the stress is constant.

The elasticity parameter values κ_S and κ_B in simulation units can be found by matching the relevant dimensionless parameters of the problem (in this case the capillary number η₀[small gamma, Greek, dot above]r/κ_S and the reduced bending modulus κ_B/(κ_Sr²)). Two different RBC rigidities have been considered. The softer RBCs, (κ_S, κ_B, κ_α, κ_A, κ_V) = (0.2, 0.004, 1, 1, 1), correspond to healthy RBCs (κ_S = 5 × 10⁻⁶ Nm⁻¹ and κ_B = 2 × 10⁻¹⁹ Nm).²³ The more rigid RBCs obey (κ_S, κ_B, κ_α, κ_A, κ_V) = (0.06, 0.012, 1, 1, 1). The former are denoted as ‘s’, the latter as ‘r’. The choice of κ_α = κ_A = κ_V = 1 ensures that local area and volume deviations are restricted to a few percent.

In order to improve the statistics, all simulations have been repeated with various random initial RBC configurations. Ten and five independent runs have been performed for the softer and more rigid RBCs, respectively. All reported observables are averaged over independent runs and time in the steady state. Due to the absence of thermal fluctuations in the present model, all observed effects are shear-induced.

2.2 Characterisation of the flow

Since we are interested in bulk properties rather than wall-induced effects, we first examine RBC concentration, suspension velocity, and shear stress profiles to see whether a bulk-like behaviour may be expected.

Fig. 1(b) reveals that, within the inner 50% of the volume between the walls (henceforth called the central region), the volume fraction is nearly constant and the velocity is linear, defining a constant shear rate. The shear stress is confirmed to be constant throughout the entire volume as expected for simple shear flow.¹⁷ Thus, an effective bulk shear viscosity η([small gamma, Greek, dot above]) ≡ σ_xz([small gamma, Greek, dot above])/[small gamma, Greek, dot above] can be defined as the ratio of shear stress σ_xz and shear rate [small gamma, Greek, dot above] using the data in the central region (Fig. 2). For the present purposes, the viscosity is only needed to define an effective capillary number (see eqn (10) below), which will turn out to be the central quantity governing the cell dynamics. A detailed investigation of the rheology will be presented in another publication. The results analysed in the central region are expected to be characteristic of the bulk properties of the system. In the following, system properties will thus be determined in the central region where all studied observables have been found to be rather independent of the transverse position z.

Fig. 2 Relative suspension viscosity η/η₀ as a function of bare capillary number Ca = η₀[small gamma, Greek, dot above]r/κ_S. η₀ is the reference viscosity of the suspending fluid without red blood cells. For all volume fractions, the suspension is shear thinning. Additionally, the viscosity increases with increasing volume fraction. One observes that the data for soft and hard RBCs overlap, when plotted versus Ca. Error bars are of the order of the symbol size. It is also noteworthy that the viscosities can be approximated by a Herschel–Bulkley law (dashed lines) of the form η/η₀ = a + b × Ca^c with the parameters (a,b,c) = (0,2.0, 0.78) for 25%, (0.004, 2.4, 073) for 45%, (0.015, 3.4, 0.72) for 55%, and (0.043, 5.0, 0.71) for 65%.

3 Results and discussion

3.1 Transition from tumbling to tank-treading-like dynamics

We investigate the average RBC tumbling (rigid body-like) frequency [small omega, Greek, macron] by tracking the orientation of the cells' inertia tensor in time. We recall the average tumbling frequency of a single rigid ellipsoid in a simple viscous shear flow,where p = a/c (a ≥ b ≥ c are the semi-axes of the ellipsoid, a and c are aligned with the shearing plane, see Fig. 3).²⁸ For a rigid sphere (p = 1), the tumbling frequency is . One finds for the inertia ellipsoid of an undeformed RBC a = b = 1.1r and c = 0.36r, and Jeffery's solution²⁸ predicts [small omega, Greek, macron]/[small gamma, Greek, dot above] = 0.30. In contrast, for purely tank-treading cells, the average is [small omega, Greek, macron] = 0, as they do not show any tumbling activity.

Fig. 3 Cross-section through an undeformed red blood cell (red) together with its equivalent inertia ellipsoid (black). r is the large RBC radius (dotted), a,b are the two large semiaxes (a = b = 1.1r) and c is the small semiaxis (c = 0.36r) of the inertia ellipsoid. The vector ô is perpendicular to the a–b-plane and is given by the normalised eigenvector of the inertia ellipsoid corresponding to the smallest eigenvector (i.e., the shortest semiaxis, c). The intrinsic membrane orientation vector î, on the other hand, is computed from the positions of the vertices of the RBC mesh and connects the centre of masses (green dots) of the top and bottom halves of the undeformed membrane. The identity of a mesh vertex of being in the “top” or “bottom” half remains unaltered during the course of the simulation, making î an intrinsic measure of the RBC orientation that is not available from the inertia tensor.

Additionally, we define the instantaneously reduced angular velocity of an RBC as

where ac is approximately the RBC cross-sectional area in the shearing plane andis the surface-averaged rotational component of the membrane velocity (r and v are location and velocity of an RBC surface element relative to the current RBC centre, respectively, and A is the RBC surface area). Note that [small omega, Greek, macron]^* is sensitive to any form of membrane rotation, even if the cell's inertia tensor is not rotating in space (i.e., if there is no tumbling).

First evidence for the onset of the TB–TT transition is provided in Fig. 4 where the reduced angular membrane velocity [small omega, Greek, macron]^*/[small gamma, Greek, dot above] and the reduced tumbling frequency [small omega, Greek, macron]/[small gamma, Greek, dot above] are plotted versus the effective capillary number

which gives the relative strength of the suspension stress σ_xz = η([small gamma, Greek, dot above])[small gamma, Greek, dot above], acting on the cells and the characteristic membrane stress, κ_S/r. The approximately constant behaviour of [small omega, Greek, macron]^* in panel (a) indicates that rotational motion is always present in the studied Ca*-range. Despite this fact, we observe a rather sharp decrease of the tumbling frequency at Ca^*_cr ≃ 0.2 in panel (b). This clearly shows that, for Ca^* > Ca^*_cr, a type of rotational motion is dominant which allows for a non-rotating tensor of inertia. Tank-treading-like dynamics provides such an alternative. Interestingly, for an isolated RBC we see a TB–TT transition in the same region (Fig. 7), which is in qualitative agreement with previous works.^29–31 Furthermore, at σ ≈ 0.4Pa, which corresponds to Ca* ≈ 0.2 for healthy RBCs, the transition has been observed experimentally for single RBCs.¹⁰ Note that all data for Ca^* > Ca^*_cr collapse onto a single master curve, which is not the case when plotted as a function of the bare capillary number Ca (not shown). The decay of the tumbling frequency in the region Ca* ≈ 0.2–0.3 can be captured by a simple power law, [small omega, Greek, macron]/[small gamma, Greek, dot above] ∝ Ca^*−3. We have to emphasise, however, that this fit is entirely ad-hoc. The establishment of such a power-law requires a more careful analysis which goes beyond the scope of the present work.

Fig. 4 Reduced (a) angular velocity [small omega, Greek, macron]^*/[small gamma, Greek, dot above] and (b) tumbling frequency [small omega, Greek, macron]/[small gamma, Greek, dot above] versus effective capillary number Ca*. The legend applies to both panels. The grey area denotes the TB–TT transition. As seen in panel (b), for Ca* > 0.1, the average tumbling frequency of a single isolated cell is zero, whereas [small omega, Greek, macron]/[small gamma, Greek, dot above] approximately decays like [small omega, Greek, macron]/[small gamma, Greek, dot above] ∝ Ca^*−3 (inclined solid line) in the dense suspensions. In contrast, the angular membrane velocity in (a) is rather independent of Ca*. Note that [small omega, Greek, macron]^*/[small gamma, Greek, dot above] and [small omega, Greek, macron]/[small gamma, Greek, dot above] are nearly identical for Ca* ≲ 0.1. The analytic values of [small omega, Greek, macron]/[small gamma, Greek, dot above] for a rigid sphere and ellipsoid with aspect ratio p = a/c = 1.1/0.36 are also shown in (b).

The picture of a TB–TT transition characterised by a critical capillary number without explicit dependence on haematocrit is further supported by a survey of the dynamics at the level of the individual cells. Fig. 5 exemplifies the typical temporal behaviour of an individual RBC in a moderately dense suspension, from which one can clearly identify tumbling motion at low effective capillary numbers, whereas tank-treading is performed at large Ca*. Interestingly, even at this comparatively large value of haematocrit (45%), the typical rotational motion of an individual cell does not appear to be significantly different from that of a single cell, apart from occasional hydrodynamical collisions with its neighbours. While a more detailed analysis of these aspects will be presented in a separate work, we point out here that a comparison of the rate of flips³² of the equivalent inertia ellipsoid on the one hand and of a suitable intrinsic membrane orientation vector† on the other hand provides an independent means to distinguish tumbling from tank-treading behaviour. This is demonstrated in Fig. 6, where we plot the ratio of the total number of flips (obtained by integrating the flip rate over several decades in inverse shear rate) of the two orientation vectors. Consistent with the expectations from the analysis of the angular velocities, we find that for low effective capillary numbers, both the inertia ellipsoid and intrinsic membrane orientation vector flip at the same rate, indicative of tumbling-like motion. In contrast, at higher Ca*, only the intrinsic membrane orientation vector performs flips while the orientation of the inertia ellipsoid remains virtually fixed, as would be expected for tank-treading-like motion.

Fig. 5 Exemplary cross-sections of a rotating RBC in a 45%-suspension at (a) Ca* = 0.05 and (b) Ca* = 0.49. The cross-sections correspond to the shearing plane and run through the cell centre. The shear flow is in the clockwise direction. In the left panel, a flipping event (rigid body-like rotation) is shown (temporal order: red, blue, green; time difference between snapshots is 0.14/[small gamma, Greek, dot above]). The right panel illustrates the irregular tank-treading-like rotation (temporal order: red, blue, green; time difference between snapshots is 0.35/[small gamma, Greek, dot above]). Shape fluctuations due to collisions with neighbouring particles are clearly visible. Arrows indicate the local surface velocity.

Fig. 6 Ratio of the number of flips (i.e., half-turns) of the inertia ellipsoid, N_flips,IE, and intrinsic cell orientation vector (see Fig. 3), N_flips,orient, of RBCs in a suspension at different haematocrit values. At low effective capillary numbers, both inertia ellipsoid and orientation vector flip at the same rate (the flip number ratio being close to unity), characteristic for tumbling-like motion. In contrast, at high Ca* only the intrinsic orientation vector performs flips, while the flip rate of the inertia ellipsoid is strongly reduced, indicating tank-treading-like motion.

Fig. 7 This sequence shows the time evolution of the rotational behaviour of a single RBC in an external shear flow with rate [small gamma, Greek, dot above]. The cross-section is parallel to the shearing plane, and the vorticity of the shear flow is clockwise. The evolution as a function of dimensionless time [small gamma, Greek, dot above]t is shown for three different capillary numbers, Ca = 0.1 (loosely dashed line), 0.2 (densely dashed line), and 0.5 (solid line). For the smallest value (Ca = 0.1), the RBC performs a tumbling motion since it is not sufficiently deformed to undergo a tank-treading motion. However, already for Ca = 0.2, the RBC can rotate without tumbling, but shape oscillations are still observed. The tank-treading motion is fully developed for the highest capillary number (Ca = 0.5). Note that in the dilute limit Ca* ≈ Ca holds.

It deserves special notice that, while the TB–TT transition is fairly well understood in terms of flow gradient effects and cell deformability (it is expected to occur when the fluid stress is strong enough to push the cells through the maximum of the elastic energy),² its occurrence in a dense suspension, where collective effects play a major role, is a more complex phenomenon. Our results suggest that, even at a haematocrit of Ht = 65%, RBCs can perform a tank-treading-like motion in an effective medium whose viscosity is no longer the viscosity of the ambient fluid but the significantly higher effective suspension viscosity. Also it has to be noted that, due to RBC collisions and related complicated dynamics, an unperturbed swinging or tank-treading motion for Ca ≳ Ca^*_cr is not expected. Instead, the cells' inclination and deformation fluctuate about their time averages, which can also be inferred from Fig. 5. Therefore, the term tank-treading-like is used to emphasise that the cells are not tumbling, but neither show a perfect, unperturbed tank-treading motion.

3.2 Cell alignment and orientational ordering

The onset of the TB–TT transition seems to be accompanied by a transition in the orientational order parameter Q_> of the system which is defined as the largest eigenvalue of the order tensor³⁴

The RBC orientation vector ô is defined as the inertia tensor eigenvector corresponding to the shortest semi-axis c (see Fig. 3). The related eigenvector of Q is called the director n indicating the average orientation of the cells.

As illustrated in Fig. 8, both Q_> and the average inclination angle [small theta, Greek, macron] (average angle between the director and the flow axis) change rapidly at Ca^* ∼ Ca^*_cr. We have also investigated the time evolution and spatial dependence of n and Q_>. Both quantities show only slight fluctuations about their averages. Note that RBC ordering at high shear rates has been observed experimentally.¹ Similarly to the behaviour of the tumbling frequency, all simulated data follow a master curve for Ca^* ≥ Ca^*_cr and can be described by one single rather than by two independent parameters (Ca* instead of Ca and Ht).

Fig. 8 (a) Order parameter and (b) average inclination angle versus effective capillary number Ca*. The curves obtained for one single cell are also shown. The grey area denotes the TB–TT transition. The legend applies to both panels. (c) Inclination angle probability distribution for the softer cells with Ht = 55% and selected values of Ca*. The analytic curve for a single rigid ellipsoid (aspect ratio p = a/c = 1.1/0.36) is also shown. The single cell in panel (b) approaches 90° for decreasing Ca*, as expected from the analytical solution.

In qualitative agreement with analytical predictions^2,3,35 and previous simulations,^36,37 the average inclination angle approaches 90° in the limit of vanishing capillary number²⁸ and grows with increasing capillary number toward the onset of tank-treading (it is bounded from above by 135°). In the tank-treading regime, the average inclination angle is expected to decrease again,^35,36,38 as the capsule becomes more elongated and aligns further with the flow. For the limited range of capillary numbers accessible to present simulations, however, we do not observe such a behaviour. Rather, the inclination angle is found to remain roughly constant or even slightly increases (which is, however, in line with the predictions of Skotheim and Secomb²).

A typical distribution of inclination angles is shown in Fig. 8(c). For comparison, the probability of finding the rigid inertia ellipsoid of an RBC (aspect ratio p = a/c) with a given orientation angle in shear flow is

This relation can be extracted from the known expression for the angular velocity of the ellipsoid.²⁸

As [small gamma, Greek, dot above] increases, the centre of P(θ) is shifted towards larger angles and its shape becomes narrower. While the former results in an increase of [small theta, Greek, macron], the latter leads to a higher orientational order Q_> since deviations from a given cell orientation become restricted to a narrower range. It is striking that, for Ca^* ≳ Ca^*_cr, the average inclination angle equals that of the isolated cell (Fig. 8(b)).

3.3 Deformation behaviour

In Fig. 9(a), some typical deformation probability distributions P(D_a) are shown. Here,is the Taylor deformation parameter as a measure for the RBC asymmetry in the a–b-plane, where the index zero refers to the undeformed state and a₀ = b₀ holds for an RBC. Using these data, we determine D^max_a, the deformation parameter with maximum probability. It is found (Fig. 9(b)) that the data on D^max_a collapse onto a single master curve when plotted as a function of Ca* for all studied values of control parameters.

Fig. 9 (a) Deformation probability distribution for the softer cells at Ht = 65% for different effective capillary numbers Ca*. For each curve, an ellipsoid is shown whose deformation parameter equals D^max_a, the deformation parameter with maximum probability. (b) D^max_aversus Ca*. The power-law fit, a/b|_max − 1 = 0.89 × Ca*^0.9, is shown as a dashed line. There is a one-to-one correspondence between a/b and D_a. The grey area denotes the TB–TT transition. Note the excellent agreement with the experimental data obtained for young RBCs reported in Fig. 2 in Pfafferott et al.³³ (we converted the shear stress reported therein to Ca* by assuming a shear elasticity of κ_S = 5 μN m⁻¹).²³

In marked contrast to the behaviour of tumbling frequency and orientational order parameter, D^max_a shows no signature of the TB–TT transition. A similar observation for dilute RBC suspensions has been reported before.¹⁰ The entire observed range of capillary numbers can be parameterised by a simple power law, a/b|_max − 1 ∝ Ca*^0.9, as shown in Fig. 9(b).

4 Conclusions and outlook

The focus of the present study is on the tumbling-to-tank-treading transition (TB–TT transition) in suspensions of aggregation-free red blood cells in the limit of high volume fractions where collective effects become dominant. We provide evidence that this transition occurs when the ratio of the effective suspension stress to the characteristic membrane stress (effective capillary number) exceeds a threshold value Ca^*_cr. Note that all dependence on volume fraction is implicit in Ca* through its dependence on the effective stress. For a single cell, Ca^*_cr corresponds to the stress where the cell is driven through the maximum of its elastic energy, thus allowing tank-treading-like dynamics.

The average tumbling frequency and cell inclination angle signal the onset of the TB–TT transition and a scaling collapse above Ca^*_cr when expressed in dependence of the effective capillary number. Remarkably, the most probable cell deformation changes continuously in the investigated parameter range following a scaling law, ∝ Ca*^0.9, without any signature of the observed transition. The value of Ca^*_cr coincides with that for a transition from a less ordered cell orientation distribution to a highly ordered phase.

The above findings support the following conclusions: (i) the cell dynamics is dominated by Ca* and therefore the suspension stress. The cells in the tank-treading-like state behave more like isolated particles experiencing their environment only via the suspension stress. (ii) The large degree of orientational ordering at Ca^* ≥ Ca^*_cr is related to the onset of tank-treading-like dynamics. (iii) Cell deformation alone does not seem to be a relevant factor for the TB–TT transition.

Interestingly, scaling of static and dynamic quantities in terms of an effective, concentration-corrected dimensionless number is not uncommon in soft matter systems and has been observed previously, for instance, in polymer solutions.^39–41 An important difference between these and the present system is, however, the volume and area conservation of the RBCs, which implies that excluded volume effects should play a much more dominant role for RBC suspensions. This might be one of the reasons why many quantities lack scaling with Ca* in the tumbling regime, where, in contrast to the tank-treading case, mutual disturbance of the cells is expected to be important.

During tank-treading, particles explore a smaller volume as compared to tumbling, thus providing a less strong hindrance to the motion of other particles. This is expected to affect other quantities, such as stresses and diffusivity, which are sensitive to the interactions between particles. Indeed, it is commonly expected that the suspension viscosity is reduced as the number of tank-treading particles increases.^1,10,14,42

Furthermore, our work opens a route to equipping coarse-grained blood models^43,44 with proper constitutive laws for the tumbling/tank-treading probability or the average particle deformation given an ambient stress value.

Finally, while we focused on the limit of high Ca*, many interesting open questions remain regarding the opposite limit of small shear rates such as the possible existence of a yield stress or flow heterogeneity.^45,46

Acknowledgements

This work is financially supported by the DFG-project VA205/5-1. We thank B. Kaoui, S. Frijters, J. Harting, and M. Ripoll for valuable discussions. We are also grateful for the computational time granted by the Juelich Supercomputing Centre (Project ESMI17). ICAMS acknowledges funding from its industrial sponsors, the state of North-Rhine Westphalia and the European Commission in the framework of the European Regional Development Fund (ERDF).

References

1. H. Schmid-Schönbein and R. Wells, Science, 1969, 165, 288–291.
2. J. Skotheim and T. Secomb, Phys. Rev. Lett., 2007, 98, 78301.
3. M. Abkarian, M. Faivre and A. Viallat, Phys. Rev. Lett., 2007, 98, 188302.
4. Y. Sui, Y. Chew, P. Roy, Y. Cheng and H. Low, Phys. Fluids, 2008, 20, 112106.
5. S. Keller and R. Skalak, J. Fluid Mech., 1982, 120, 27–47.
6. J. Beaucourt, F. Rioual, T. Séon, T. Biben and C. Misbah, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2004, 69, 011906.
7. H. Noguchi and G. Gompper, Phys. Rev. Lett., 2004, 93, 258102.
8. V. Kantsler and V. Steinberg, Phys. Rev. Lett., 2006, 96, 036001.
9. C. Misbah, Phys. Rev. Lett., 2006, 96, 028104.
10. A. Forsyth, J. Wan, P. Owrutsky, M. Abkarian and H. Stone, Proc. Natl. Acad. Sci. U. S. A., 2011, 108, 10986.
11. A. Yazdani, R. Kalluri and P. Bagchi, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2011, 83, 046305.
12. S. Doddi and P. Bagchi, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2009, 79, 046318.
13. J. R. Clausen, D. A. Reasor and C. K. Aidun, J. Fluid Mech., 2011, 685, 1–33.
14. D. A. Fedosov, W. Pan, B. Caswell, G. Gompper and G. E. Karniadakis, Proc. Natl. Acad. Sci. U. S. A., 2011, 108, 11772–11777.
15. D. A. Reasor, J. R. Clausen and C. K. Aidun, J. Fluid Mech., 2013, 726, 497.
16. T. Krüger, F. Varnik and D. Raabe, Comput. Math. Appl., 2011, 61, 3485–3505.
17. T. Krüger, F. Varnik and D. Raabe, Philos. Trans. R. Soc., A, 2011, 369, 2414–2421.
18. S. Succi, The Lattice Boltzmann Equation for Fluid Dynamics and Beyond, Oxford University Press, 2001, p. 368.
19. D. Aidun and J. Clausen, Annu. Rev. Fluid Mech., 2010, 42, 439–472.
20. C. Peskin, Acta Numerica, 2002, 11, 479–517.
21. N. Mohandas and E. Evans, Annu. Rev. Biophys. Biomol. Struct., 1994, 23, 787–818.
22. S. Svetina, D. Kuzman, R. Waugh, P. Ziherl and B. Zeks, Bioelectroch., 2004, 62, 107–113.
23. G. Gompper and M. Schick, Soft Matter: Lipid Bilayers and Red Blood Cells, Wiley-VCH, 2008.
24. R. Skalak, A. Tozeren, R. Zarda and S. Chien, Biophys. J., 1973, 13, 245–264.
25. W. Helfrich, Z. Naturforsch., C: Biochem., Biophys., Biol., Virol., 1973, 28, 693–703.
26. G. Gompper and D. Kroll, J. Phys. I, 1996, 6, 1305–1320.
27. Z.-G. Feng and E. Michaelides, J. Comput. Phys., 2004, 195, 602–628.
28. G. Jeffery, Proc. R. Soc. London, Ser. A, 1922, 102, 161–179.
29. C. Pozrikidis, J. Fluid Mech., 2001, 440, 269–291.
30. C. Pozrikidis, Ann. Biomed. Eng., 2003, 31, 1194–1205.
31. J. Zhang, P. Johnson and A. Popel, Microvasc. Res., 2009, 77, 265–272.
32. F. Janoschek, F. Mancini, J. Harting and F. Toschi, Philos. Trans. R. Soc., A, 2011, 369, 2337–2344.
33. C. Pfafferott, G. Nash and H. Meiselman, Biophys. J., 1985, 47, 695–704.
34. A. Saupe, Angew. Chem., Int. Ed., 1968, 7, 97.
35. P. M. Vlahovska, Y.-N. Young, G. Danker and C. Misbah, J. Fluid Mech., 2011, 678, 221–247.
36. P. Bagchi and R. Kalluri, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2009, 80, 016307.
37. J. R. Clausen and C. K. Aidun, Phys. Fluids, 2010, 22, 123302.
38. D. Barths-Biesel, J. Fluid Mech., 1980, 100, 831853.
39. J. S. Hur, E. S. G. Shaqfeh, H. P. Babcock, D. E. Smith and S. Chu, J. Rheol., 2001, 45, 421.
40. C.-C. Huang, R. G. Winkler, G. Sutmann and G. Gomppper, Macromolecules, 2010, 43, 10107.
41. S. P. Singh, D. A. Fedosov, A. Chatterji, R. G. Winkler and G. Gompper, J. Phys.: Condens. Matter, 2012, 24, 464103.
42. S. Chien, Science, 1970, 168, 977–979.
43. F. Janoschek, F. Toschi and J. Harting, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2010, 82, 056710.
44. S. Melchionna, Macromol. Theory Simul., 2011, 20, 548561.
45. F. Varnik, L. Bocquet, J. Barrat and L. Berthier, Phys. Rev. Lett., 2003, 90, 95702.
46. S. Mandal, M. Gross, D. Raabe and F. Varnik, Phys. Rev. Lett., 2012, 108, 098301.

---

Footnote

† The intrinsic membrane orientation vector is computed at the beginning of the simulation as the vector connecting the centre of masses of the top and bottom halves of the undeformed RBC, see Fig. 3.

---