Julien Sigüenza
email: j.siguenza@sim-and-cure.com

Simon Mendez

Franck Nicoud

How should the optical tweezers experiment be used to characterize the red blood cell membrane mechanics?

Keywords: Red blood cells, Optical tweezers, Membrane modeling, Cytoskeleton, Lipid bilayer, Fluid-structure interactions, Immersed boundary method

Stretching red blood cells using optical tweezers is a way to characterize the mechanical properties of their membrane by measuring the size of the cell in the direction of the stretching (axial diameter) and perpendicularly (transverse diameter). Recently, such data have been used in numerous publications to validate solvers dedicated to the computation of red blood cell dynamics under flow. In the present study, different mechanical models are used to simulate the stretching of red blood cells by optical tweezers. Results first show that the mechanical moduli of the membranes have to be adjusted as a function of the model used. In addition, by assessing the area dilation of the cells, the axial and transverse diameters measured in optical tweezers experiments are found to be insufficient to discriminate between models relevant to red blood cells or not. At last, it is shown that other quantities such as the height or the profile of the cell should be preferred for validation purposes since they are more sensitive to the membrane model.

Introduction

Blood is a complex substance consisting in a suspension of platelets, white blood cells and red blood cells (RBCs) in a Newtonian fluid, the plasma. The RBCs, which typically represent 40-45% of the whole blood volume, are composed of a membrane enclosing an internal fluid, the cytoplasm. The RBC membrane is a composite structure composed of a lipid bilayer and a two-dimensional elastic cytoskeleton, both linked through temporary tethering sites thanks to transmembrane proteins embedded in the lipid bilayer. This complex structure confers to the RBC membrane very specific mechanical properties: The cytoskeleton provides a resistance to shear solicitations and slightly resists to area dilatation, while the lipid bilayer provides to the membrane its bending stiffness and quasi-incompressibility. The RBCs have a biconcave discocyte shape at rest with a remarkable deformability, because of the excess of surface area enclosing the inner volume. RBCs are thus able to undergo very large deformation preserving their area, squeezing through capillaries with inner diameter less than 3 µm, although the average large diameter of a RBC is about 8 µm. As mentioned by [START_REF] Mohandas | Red cell membrane: past, present, and future[END_REF], the normal RBC can deform with linear extensions of up to 250%, but a 3-4% increase in surface area results in cell lysis.

So far, there is no universal model to describe the mechanical behavior of the RBC membrane. The local elasticity of the RBC membrane is generally described using either continuum models [START_REF] Le | An implicit immersed boundary method for three-dimensional fluidmembrane interactions[END_REF][START_REF] Klöppel | A novel two-layer, coupled finite element approach for modeling the nonlinear elastic and viscoelastic behavior of human erythrocytes[END_REF][START_REF] Farutin | 3D numerical simulations of vesicle and inextensible capsule dynamics[END_REF][START_REF] Sinha | Dynamics of a single red blood cell in simple shear flow[END_REF] or network models [START_REF] Li | Spectrin-level modeling of the cytoskeleton and optical tweezers stretching of the erythrocyte[END_REF][START_REF] Dao | Molecularly based analysis of deformation of spectrin network and human erythrocyte[END_REF][START_REF] Pivkin | Accurate coarse-grained modeling of red blood cells[END_REF]Fedosov et al. 2010aFedosov et al. , b, 2014;;[START_REF] Chen | Investigation of membrane mechanics using spring networks: application to red-blood-cell modelling[END_REF], which can be complemented with other global models to treat the quasi-incompressibility of the lipid bilayer [START_REF] Pivkin | Accurate coarse-grained modeling of red blood cells[END_REF]Fedosov et al. 2010a, b). Detailed experi-123 mental investigations of the RBC mechanics are nonetheless needed in order to: (1) characterize and validate a numerical model of the RBC membrane and (2) once validated, determine the mechanical parameters of the model.

To gain insight into the mechanical behavior of RBCs, experimental techniques were developed for measurements of the RBC membrane properties [START_REF] Abkarian | On the importance of the deformability of red blood cells in blood flow[END_REF]. Micropipette aspiration [START_REF] Evans | New membrane concept applied to the analysis of fluid shear-and micropipette-deformed red blood cells[END_REF]) and optical tweezers [START_REF] Hénon | A new determination of the shear modulus of the human erythrocyte membrane using optical tweezers[END_REF][START_REF] Mills | Nonlinear elastic and viscoelastic deformation of the human red blood cell with optical tweezers[END_REF]) are the most popular ones and were notably used to determine the shear modulus of the RBC membrane. The optical tweezers experiment by [START_REF] Mills | Nonlinear elastic and viscoelastic deformation of the human red blood cell with optical tweezers[END_REF] provides a useful means for the analysis of the single cell mechanics under a variety of well-controlled stress states, where stretching of an isolated RBC is generated by means of attached silica microbeads and optical trap. Using a continuum model of the RBC membrane to solve the deformation of the RBC subjected to optical stretching, [START_REF] Yeoh | Some forms of the strain energy function for rubber[END_REF] successfully matched the force-extension data obtained from the experiment, thus enabling the extraction of the shear modulus of the RBC membrane.

A recent work of [START_REF] Dimitrakopoulos | Analysis of the variation in the determination of the shear modulus of the erythrocyte membrane: effects of the constitutive law and membrane modeling[END_REF] showed that large differences of shear modulus reported in various studies may be explained based on the different membrane models used to fit the experimental data. Theoretically investigating continuum models under uniaxial extension and local area incompressibility, he showed that the only constitutive law able to properly match the wide variety of experimental data available in the literature is the Skalak law, specifically developed by [START_REF] Skalak | Strain energy function of red blood cell membranes[END_REF] to represent the in-plane elasticity of the RBC membrane. Based on this finding, Dimitrakopoulos stated that [START_REF] Mills | Nonlinear elastic and viscoelastic deformation of the human red blood cell with optical tweezers[END_REF] found the shear modulus that represents the Yeoh law, but not the true shear modulus of the RBC membrane. This purely theoretical work was specifically dedicated to the response of different membrane laws under small, moderate and large shear strains. A more realistic configuration where the whole RBC is stretched as in the [START_REF] Mills | Nonlinear elastic and viscoelastic deformation of the human red blood cell with optical tweezers[END_REF] experiment was, however, not considered by [START_REF] Dimitrakopoulos | Analysis of the variation in the determination of the shear modulus of the erythrocyte membrane: effects of the constitutive law and membrane modeling[END_REF]; this is done computationally in the present paper.

As a consequence, the numerical results of [START_REF] Mills | Nonlinear elastic and viscoelastic deformation of the human red blood cell with optical tweezers[END_REF] were successfully matched to the force-extension data obtained from optical tweezers using the Yeoh law, whereas a proper modeling of the RBC membrane should rather rely on the Skalak law. This reveals the simplistic nature of these experimental data, which was also pointed out by [START_REF] Dimitrakopoulos | Analysis of the variation in the determination of the shear modulus of the erythrocyte membrane: effects of the constitutive law and membrane modeling[END_REF]. Despite this observation, optical tweezers data continue to be used as a way to validate numerical models of the RBC membrane [START_REF] Li | Spectrin-level modeling of the cytoskeleton and optical tweezers stretching of the erythrocyte[END_REF][START_REF] Dao | Molecularly based analysis of deformation of spectrin network and human erythrocyte[END_REF][START_REF] Pivkin | Accurate coarse-grained modeling of red blood cells[END_REF][START_REF] Le | An implicit immersed boundary method for three-dimensional fluidmembrane interactions[END_REF]Fedosov et al. 2010aFedosov et al. , b, 2014;;[START_REF] Klöppel | A novel two-layer, coupled finite element approach for modeling the nonlinear elastic and viscoelastic behavior of human erythrocytes[END_REF][START_REF] Chen | Investigation of membrane mechanics using spring networks: application to red-blood-cell modelling[END_REF][START_REF] Farutin | 3D numerical simulations of vesicle and inextensible capsule dynamics[END_REF][START_REF] Sinha | Dynamics of a single red blood cell in simple shear flow[END_REF], notably to probe the accuracy of solvers dedicated to the study of the RBC dynamics under flow. However, a proper validation test case needs to be selective to discriminate between appropriate and inappropriate models. There is a suspicion that computing optical tweezers experiment does not constitute a true validation test case.

The present paper constitutes a numerical study which first aims at emphasizing previous findings of Dimitrakopoulos (2012), highlighting the limitations of the optical tweezers experiment for characterizing the mechanics of the RBC membrane. Theoretical investigations of [START_REF] Dimitrakopoulos | Analysis of the variation in the determination of the shear modulus of the erythrocyte membrane: effects of the constitutive law and membrane modeling[END_REF] are here complemented with detailed simulations of the optical tweezers experiment by [START_REF] Mills | Nonlinear elastic and viscoelastic deformation of the human red blood cell with optical tweezers[END_REF], using a numerical method dedicated to the simulation of the dynamics of RBCs under flow. After a brief description of this numerical method, an easy-to-implement computational setup is presented and validated against the numerical results of [START_REF] Mills | Nonlinear elastic and viscoelastic deformation of the human red blood cell with optical tweezers[END_REF]. Then, different continuum membrane models are investigated, based on various combinations of strain, area conservation and bending energies. If the membrane incompressibility can be easily imposed theoretically [START_REF] Dimitrakopoulos | Analysis of the variation in the determination of the shear modulus of the erythrocyte membrane: effects of the constitutive law and membrane modeling[END_REF], it is rarely done in models of red blood cells. The membrane is generally modeled using a mechanical resistance to area dilatation, which enables some small area variation of the membrane [START_REF] Li | Spectrin-level modeling of the cytoskeleton and optical tweezers stretching of the erythrocyte[END_REF][START_REF] Dao | Molecularly based analysis of deformation of spectrin network and human erythrocyte[END_REF][START_REF] Pivkin | Accurate coarse-grained modeling of red blood cells[END_REF]Fedosov et al. 2010aFedosov et al. , b, 2014;;[START_REF] Chen | Investigation of membrane mechanics using spring networks: application to red-blood-cell modelling[END_REF][START_REF] Sinha | Dynamics of a single red blood cell in simple shear flow[END_REF]. In the present study, the impact of this area dilatation resistance is carefully investigated, restraining the area variation of the membrane either locally or globally. Detailed analysis of the shape of the stretched RBC are also carried out in order to identify which kind of additional experimental data could be helpful to better characterize the mechanics of the RBC membrane.

Numerical method

The present numerical method is very similar to the one developed by [START_REF] Mendez | An unstructured solver for simulations of deformable particles in flows at arbitrary Reynolds numbers[END_REF] and [START_REF] Sigüenza | Validation of an immersed thick boundary method for simulating fluid-structure interactions of deformable membranes[END_REF] for fluid-structure interactions (FSI) of deformable membranes, and is based on the immersed boundary method (IBM) introduced by [START_REF] Peskin | The immersed boundary method[END_REF]. Two independent meshes are considered to discretize the RBC membrane and the fluid. The RBC membrane is discretized by a moving Lagrangian mesh, and the fluid is discretized by a fixed Eulerian unstructured mesh. The different steps of the present method are the following:

(1) The membrane force -→ F is calculated on the Lagrangian mesh, which depends on the membrane deformation and on the models used to represent the membrane rheology.

(2) The forces exerted by the membrane on the fluid are represented by the fluid volumetric force -→ f , calculated on the Eulerian mesh by regularizing the membrane force -→ F such as

-→ f -→ x , t = s -→ F -→ X , t δ -→ x - -→ X dX,
where -→ x and -→ X , respectively, denote the coordinates vectors of the Eulerian fluid nodes and Lagrangian nodes, s denotes the solid domain defining the RBC membrane and δ is the well-known Dirac function.

(3) The fluid velocity -→ v is calculated on the Eulerian mesh by solving the Navier-Stokes equations (forced by the source term -→ f ). (4) The membrane velocity -→ V is calculated on the Lagrangian mesh by interpolating the fluid velocity -→ v such as

-→ V -→ X , t = f -→ v -→ x , t δ -→ x - -→ X dX,
where f denotes the fluid domain.

The Dirac function δ used in the procedures of regularization and interpolation of steps (2) and ( 4) is numerically represented by a smooth discrete Dirac function, which is adapted to unstructured meshes using the Reproducing Kernel particle method [START_REF] Pinelli | Immersed-boundary methods for general finite-difference and finite-volume Navier-Stokes solvers[END_REF][START_REF] Mendez | An unstructured solver for simulations of deformable particles in flows at arbitrary Reynolds numbers[END_REF][START_REF] Sigüenza | Validation of an immersed thick boundary method for simulating fluid-structure interactions of deformable membranes[END_REF]. Interpolation of the fluid velocity on the membrane Lagrangian mesh leads to small mass conservation errors. A specific algorithm has been developed to perfectly conserve the volume of the RBC during the calculations [START_REF] Mendez | An unstructured solver for simulations of deformable particles in flows at arbitrary Reynolds numbers[END_REF][START_REF] Sigüenza | Validation of an immersed thick boundary method for simulating fluid-structure interactions of deformable membranes[END_REF].

Membrane forces computation

In the present method, the RBC membrane is considered to be infinitely thin and is represented by a triangulated surface. The membrane force is derived from a combination of strain, area conservation and bending energies. Resistances to shear and area dilatation are modeled thanks to a hyperelastic strain energy function W , which is written as a function of the local in-plane principal values of strain λ 1 and λ 2 , following the method of [START_REF] Charrier | Free and constrained inflation of elastic membranes in relation to thermoforming nonaxisymmetric problems[END_REF], [START_REF] Eggleton | Large deformation of red blood cell ghosts in a simple shear flow[END_REF], [START_REF] Sui | Dynamic motion of red blood cells in simple shear flow[END_REF] and [START_REF] Doddi | Lateral migration of a capsule in a plane Poiseuille flow in a channel[END_REF]. Several hyperelastic models are investigated in the present study:

• The neo-Hookean law,

W NH = E s 2 λ 2 1 + λ 2 2 + λ -2 1 λ -2 2 -3 , (1)
where E s stands for the membrane in-plane shear modulus.

• The Yeoh law,

W YE = E s 2 λ 2 1 + λ 2 2 + λ -2 1 λ -2 2 -3 + C 3 λ 2 1 + λ 2 2 + λ -2 1 λ -2 2 -3 3 , (2)
which is an extension of the previous neo-Hookean law, with the addition of a nonlinear term driven by the nonlinear modulus C 3 . • The law introduced by [START_REF] Skalak | Strain energy function of red blood cell membranes[END_REF] for red blood cells,

W SK = E s 4 λ 2 1 + λ 2 2 -2 2 + 2 λ 2 1 + λ 2 2 -λ 2 1 λ 2 2 -1 + E a 4 λ 2 1 λ 2 2 -1 2 , (3)
where shear resistance and area dilatation resistance are separately taken into account through the shear modulus E s and the area dilatation modulus E a , respectively. It can also be written with the ratio of the area dilatation modulus to the shear modulus, C = E a /E s ,

W SK = E s 4 λ 2 1 +λ 2 2 -2 2 + 2 λ 2 1 +λ 2 2 -λ 2 1 λ 2 2 -1 + C λ 2 1 λ 2 2 -1 2 . ( 4
)
Although the Skalak law can be used to control area variations of the RBC membrane, another approach consists in using a global area conservation energy:

E S = κ S 2 (S -S 0 ) 2 S 0 , (5)
with κ S the area modulus, S the area of the membrane and S 0 its target area. This energy is actually already used in other formulations based on discrete approaches [START_REF] Pivkin | Accurate coarse-grained modeling of red blood cells[END_REF]Fedosov et al. 2010b) or in shape predictions by energy minimization [START_REF] Lim | Stomatocyte-discocyteechinocyte sequence of the human red blood cell: evidence for the bilayer-couple hypothesis from membrane mechanics[END_REF][START_REF] Lim | Red blood cell shapes and shape transformations: Newtonian mechanics of a composite membrane, soft matter, vol lipid bilayers and red blood cells[END_REF].

Conveniently, the force applied by the membrane on the fluid associated with the energy term E S can be expressed explicitly:

-→ F S = -2κ S (S -S 0 ) S 0 H -→ n , (6)
with H the mean curvature and -→ n the outward normal vector to the surface.

In addition, the bending resistance of the membrane can be represented using the bending energy E b , proposed by [START_REF] Helfrich | Elastic properties of lipid bilayers: theory and possible experiments[END_REF]:

E b = κ b 2 S (2H -c 0 ) 2 dS, (7)
123

Table 1 Different energies available to model the RBC membrane and associated notations of mechanical moduli with κ b = 2.0 × 10 -19 N m [START_REF] Lim | Stomatocyte-discocyteechinocyte sequence of the human red blood cell: evidence for the bilayer-couple hypothesis from membrane mechanics[END_REF][START_REF] Lim | Red blood cell shapes and shape transformations: Newtonian mechanics of a composite membrane, soft matter, vol lipid bilayers and red blood cells[END_REF] the bending modulus, and c 0 a possible spontaneous curvature (which is set to zero in the present study). The bending force applied by the membrane on the fluid reads:

-→ F b = κ b (2H -c 0 ) 2H 2 -2K + c 0 H + 2 LB H ] -→ n , (8)
where LB denotes the surface Laplacian operator (Zhongcan and Helfrich 1989) (also called the Laplace-Beltrami operator) and K is the local Gaussian curvature of the surface.

The terms of the bending force are calculated by local fitting of a quadratic approximation of the surface. The method is similar to the one used by [START_REF] Farutin | 3D numerical simulations of vesicle and inextensible capsule dynamics[END_REF]. Table 1 summarizes the three energies introduced, with the associated parameters. Every combination of these energies (W , E S , E b ) can be used to model the RBC membrane.

Navier-Stokes equations solver

The fluid inside and outside the RBC is supposed to be incompressible and Newtonian. The YALES2BIO flow solver is used [START_REF] Mendez | An unstructured solver for simulations of deformable particles in flows at arbitrary Reynolds numbers[END_REF][START_REF] Chnafa | Image-based large-eddy simulation in a realistic left heart[END_REF][START_REF] Sigüenza | Validation of an immersed thick boundary method for simulating fluid-structure interactions of deformable membranes[END_REF][START_REF] Zmijanovic | About the numerical robustness of biomedical benchmark cases: interlaboratory FDA's idealized medical device[END_REF] to solve the forced Navier-Stokes equations over the Eulerian unstructured mesh by using a projection method [START_REF] Chorin | Numerical solution of the Navier-Stokes equations[END_REF]). The momentum conservation equations reads:

∂ -→ v ∂t + -→ ∇ . -→ v ⊗ -→ v = - -→ ∇ p ρ + ν -→ v + -→ f ρ , (9)
where -→ v and p are the velocity vector and the pressure, ρ the density and ν the kinematic viscosity. For an incompressible fluid, the mass conservation equation becomes:

-→ ∇ . -→ v = 0 (10)
The fluid velocity is advanced using a fourth-order centered scheme in space and a fourth-order Runge-Kutta scheme in time. A divergence-free velocity field is obtained at the end of the time step by solving a Poisson equation for pressure and correcting the predicted velocity. A deflated preconditioned conjugate gradient (DPCG) algorithm is used to solve this Poisson equation [START_REF] Moureau | Design of a massively parallel CFD code for complex geometries[END_REF][START_REF] Malandain | Optimization of the deflated conjugate gradient algorithm for the solving of elliptic equations on massively parallel machines[END_REF].

The YALES2BIO solver was validated in several test cases where reference data (either experimental, analytical or numerical) are available. This is described in previous publications, where the reader can also find additional implementation details [START_REF] Martins Afonso | On the damped oscillations of an elastic quasi-circular membrane in a two-dimensional incompressible fluid[END_REF][START_REF] Mendez | An unstructured solver for simulations of deformable particles in flows at arbitrary Reynolds numbers[END_REF][START_REF] Sigüenza | Characterisation of a dedicated mechanical model for red blood cells: numerical simulations of optical tweezers experiment[END_REF][START_REF] Sigüenza | Validation of an immersed thick boundary method for simulating fluid-structure interactions of deformable membranes[END_REF][START_REF] Zmijanovic | About the numerical robustness of biomedical benchmark cases: interlaboratory FDA's idealized medical device[END_REF].

Optical tweezers modeling

The purpose of this section is to establish a computational setup allowing the computation of the optical tweezers experiment by [START_REF] Mills | Nonlinear elastic and viscoelastic deformation of the human red blood cell with optical tweezers[END_REF]. The computational setup presented in this section is built heavily on the one developed by [START_REF] Dao | Mechanics of the human red blood cell deformed by optical tweezers[END_REF], which has also been used by [START_REF] Mills | Nonlinear elastic and viscoelastic deformation of the human red blood cell with optical tweezers[END_REF] to simulate the optical tweezers experiment.

Figure 1a illustrates the experimental setup used in [START_REF] Mills | Nonlinear elastic and viscoelastic deformation of the human red blood cell with optical tweezers[END_REF] to perform the stretching of the RBC. Two silica microbeads, of diameter 4.12 µm, are attached to the cell at diametrically opposite points. The left bead is anchored to the surface of a glass slide while the right bead is trapped by a laser beam. The trapped bead remaining at rest, moving the slide and attached left bead stretches the cell. Then, the axial diameter D A (in the direction of the stretching), and the transverse diameter D T (orthogonal to the stretching direction) are measured on the stretched RBC.

Computational setup

The analytical model of the RBC biconcave shape proposed by [START_REF] Evans | Improved measurements of the erythrocyte geometry[END_REF] is used to define the RBC geometry:

z = ±0.5R 0 1 - x 2 + y 2 R 2 0 ×   A 1 + A 2 x 2 + y 2 R 2 0 + A 3 x 2 + y 2 R 2 0 2   ( 11
)
where R 0 = 3.91 µm is the average RBC radius, A 1 = 0.207161, A 2 = 2.002558 and A 3 = -1.122762.

Rather than explicitly solving the contact between the beads and the RBC (as [START_REF] Dao | Mechanics of the human red blood cell deformed by optical tweezers[END_REF][START_REF] Mills | Nonlinear elastic and viscoelastic deformation of the human red blood cell with optical tweezers[END_REF]), most of the works simulating the optical tweezers experiment consider pure Neumann loading conditions to simulate the RBC stretching, applying a constant stretching force F over a certain percentage of nodes at the extremities of the RBC mean positions of A stretching force F is applied over the two circular edges delimitating the contact areas between the RBC and the beads, with a contact size d c = 2 µm [START_REF] Le | An implicit immersed boundary method for three-dimensional fluidmembrane interactions[END_REF][START_REF] Farutin | 3D numerical simulations of vesicle and inextensible capsule dynamics[END_REF][START_REF] Chen | Investigation of membrane mechanics using spring networks: application to red-blood-cell modelling[END_REF][START_REF] Fedosov | Multiscale modeling of blood flow: from single cells to blood rheology[END_REF][START_REF] Sinha | Dynamics of a single red blood cell in simple shear flow[END_REF]. The drawback of this approach was nonetheless pointed out by [START_REF] Klöppel | A novel two-layer, coupled finite element approach for modeling the nonlinear elastic and viscoelastic behavior of human erythrocytes[END_REF]: The rigidity of the beads is not properly taken into account, leading to a larger axial diameter (D A ), and thus a higher estimation of the in-plane shear modulus. An alternative methodology which mimics the beads rigidity is introduced in what follows, within a three-step strategy:

• The contact areas between the beads and the RBC are properly defined following the procedure of [START_REF] Dao | Mechanics of the human red blood cell deformed by optical tweezers[END_REF]. As shown in Fig. 1b, these contact areas are defined by intersecting the surface of the RBC with two opposite planes perpendicular to the stretching direction. The position of these planes is chosen such that the contact size between the beads and the RBC is d c = 2 µm [START_REF] Dao | Mechanics of the human red blood cell deformed by optical tweezers[END_REF]. • Rather than applying the stretching force F over all the nodes of the contact areas, the force is applied only to the nodes located on the edges delimiting the contact areas (see Fig. 1b). • Instead of evaluating the axial diameter (D A ) as the distance between the extremities of the stretched RBC, the axial diameter is determined by calculating the mean position of each loaded edge, which are deformed during the RBC stretching (as sketched in Fig. 1b).

Consistent with the numerical framework described in Sect. 2, the computation of the RBC stretching consists in solving a transient fluid-structure interaction problem until stabilization of the shape. The RBC is immersed in a fluid box extended from -4R 0 to 4R 0 in the x direction (direction of the stretching), from -2R 0 to 2R 0 in the y direction (direction orthogonal to the stretching), and from -R 0 to R 0 in the z direction (direction perpendicular to the plane of the RBC). The fluid mesh is composed of 881 992 tetrahedral elements, with a constant mesh resolution of R 0 /12.5. The RBC membrane is composed of 6 434 nodes, with a constant mesh resolution of R 0 /25.

The stretching force is applied on the RBC membrane as an external force, with a time-dependent ramp ranging from 0 to the desired value of F. This external force is seen by the fluid which starts moving, and deforms the RBC. After a transient phase, the mechanical forces inside the membrane and the applied external force balance, and a steady deformation is obtained. The choice of the fluid properties and the size of the computational domain may affect the transient phase, but have no influence on the steady deformation of the RBC and calculated axial and transverse diameters. Only the final stabilized shapes are postprocessed.

Validation

With the aim of validating the present computational setup, the optical tweezers experiment by [START_REF] Mills | Nonlinear elastic and viscoelastic deformation of the human red blood cell with optical tweezers[END_REF] is simulated, and the present simulations are compared with the numerical simulations performed by [START_REF] Mills | Nonlinear elastic and viscoelastic deformation of the human red blood cell with optical tweezers[END_REF]. Two cases are simulated, corresponding to different modeling of the RBC membrane. These two cases are summarized in Table 2. For both cases, only the local in-plane elasticity is considered. The membrane is assumed to follow the neo-

123

Table 2 Cases simulated with the present computational setup and compared with the results of Mills et al. ( 2004) Hookean law (Eq. 1) in case 1, and the Yeoh law (Eq. 2) in case 2.

W E S E b Case 1 NH: E s = 7.3 µN/m X X Case 2 YE: E s = 7.3 µN/m X X C 3 = E s /30
Figure 2 shows both axial (D A ) and transverse (D T ) diameters of the RBC stretched by optical tweezers, as a function of the applied force, for cases 1 and 2. As the cell is more and more elongated when increasing the stretching force, it is seen that the axial diameter (D A ) increases. The elongation of the cell leads to its contraction in the orthogonal direction, resulting in a decrease of the transverse diameter (D T ).

When using pure Neumann loading conditions to simulate the RBC stretching [START_REF] Le | An implicit immersed boundary method for three-dimensional fluidmembrane interactions[END_REF][START_REF] Farutin | 3D numerical simulations of vesicle and inextensible capsule dynamics[END_REF][START_REF] Chen | Investigation of membrane mechanics using spring networks: application to red-blood-cell modelling[END_REF][START_REF] Fedosov | Multiscale modeling of blood flow: from single cells to blood rheology[END_REF][START_REF] Sinha | Dynamics of a single red blood cell in simple shear flow[END_REF], the rigidity of the beads used in the optical tweezers experiment is not taken into account, which is known to strongly influence the deformation of the stretched RBC, especially the estimation of the axial diameter (D A ) (Klöppel and Wall 2011). The present results, however, show that it is possible to mimic the beads rigidity using a customized computational setup based on pure Neumann loading conditions, which is seen to faithfully reproduce the numerical results obtained by [START_REF] Mills | Nonlinear elastic and viscoelastic deformation of the human red blood cell with optical tweezers[END_REF], who explicitly solved the contact between the beads and the RBC.

As pointed out by [START_REF] Mills | Nonlinear elastic and viscoelastic deformation of the human red blood cell with optical tweezers[END_REF], comparison of the numerical results of case 1 with the experimental data shows that the neo-Hookean law is not adapted to describe the behavior of the RBC membrane. Indeed, experimental trends are well captured over the range of 0-88 pN. However, the model deviates gradually for loadings higher than 88 pN, showing a strain-softening behavior under large deformation [START_REF] Barthès-Biesel | Effect of constitutive laws for two-dimensional membranes on flow-induced capsule deformation[END_REF]. Conversely, the Yeoh law provides accurate predictions of diameters over the entire range of experimental data. The strain-hardening behavior of RBCs under large deformation is thus well transcribed by the model. Regarding the mechanical response of the stretched RBC in terms of axial (D A ) and transverse (D T ) diameters, the membrane modeling corresponding to case 2, using the Yeoh law, provides a good description of the membrane mechanical behavior.

In order to investigate the influence of the mesh resolution, two meshes were constructed from the mesh used in Fig. 2: A coarse mesh whose resolution is twice coarser than the reference mesh resolution and a fine mesh whose resolution is twice finer than the reference mesh resolution. Axial (D A ) and transverse (D T ) diameters obtained from these three meshes are compared in Table 3 with the diameters obtained from numerical simulations of [START_REF] Mills | Nonlinear elastic and viscoelastic deformation of the human red blood cell with optical tweezers[END_REF] for the largest loading F = 193 pN. The mesh resolution has almost no influence on the prediction of the axial diameter (D A ), and only small influence on the prediction of the transverse diameter (D T ). This indicates that the reference mesh is sufficiently refined, and can thus be used in the remainder of this study.

Figure 3 shows the deformation of the RBC for different values of the stretching force F, which ranges from 0 to 193 pN. A detailed analysis of the shape of the RBC shows that as the cell is elongated when increasing the force, a large fold is appearing, as also observed in the numerical simulations of [START_REF] Mills | Nonlinear elastic and viscoelastic deformation of the human red blood cell with optical tweezers[END_REF]. Occurence of such a folding is, however, not investigated in the experiment. Fig. 3 Visualization of the red blood cell deformation over the entire range of stretching force, for both cases 1 and 2. Only half of the cell is displayed

Influence of the membrane modeling

The present computational setup is now used to investigate different continuum models of the RBC membrane. With the present numerical method, the different mechanical properties of the RBC membrane can be modeled by a combination of strain, area conservation and bending energies. Four new cases are summarized in Table 4. Note that the bending stiffness of the lipid bilayer was neglected in cases 1 and 2, but is accounted for in the others. Using the Yeoh law (Eq. 2) to describe the local in-plane elasticity of the RBC membrane was seen to provide a good agreement with the optical tweezers experiment (see Fig. 2b). Case 3 thus appears to be a first obvious candidate to model the mechanics of the RBC membrane. As stated by Dimitrakopoulos (2012), the RBC membrane should rather be modeled by the Skalak law instead of the Yeoh law. Cases 4 and 5 are thus introduced, with two different values of the ratio C (low value in case 4, and high value in case 5). Note, however, that when using the Skalak law to model the local in-plane elasticity of the RBC membrane, a high value of C should be considered to restrain the area variations of the RBC membrane, thus modeling the quasi-incompressibility of the lipid bilayer. Consequently, case 4 does not constitute a potential candidate to model the mechanics of the RBC membrane, but is only introduced to investigate the influence of the ratio C on the mechanical response of the RBC subjected to optical stretching. Finally, case 6 proposes a hybrid modeling of the RBC membrane, dissociating the cytoskeleton and the lipid bilayer: The Skalak law with low ratio C is used to model the local in-plane elasticity of the cytoskeleton, allowing local area changes of the cytoskeleton; on top of this, the global area conservation energy is used to model the reorganisation of the quasi-incompressible lipid bilayer, sliding along the cytoskeleton. It is noticed that a twice smaller shear modulus E s is considered when using the Skalak law in cases 4, 5 and 6, as compared to case 3. This factor of 2 is explained in the work of [START_REF] Dimitrakopoulos | Analysis of the variation in the determination of the shear modulus of the erythrocyte membrane: effects of the constitutive law and membrane modeling[END_REF] by the fact that the Yeoh and Skalak laws behave differently at moderate and high deformation. It is thus required to multiply the shear modulus E s by 2 when considering the Yeoh law as compared to the Skalak law, in order to have a good comparison with the optical tweezers experiment in the large deformation range. Note that this results in an underestimation of the cell deformation for low stretching forces with the Yeoh law, as illustrated in the next section.

Comparison of axial and transverse diameters

Figure 4 shows the numerical predictions of the axial (D A ) and transverse (D T ) diameters for the different modeling cases introduced in Table 4. All cases provide a good comparison with the experimental results of [START_REF] Mills | Nonlinear elastic and viscoelastic deformation of the human red blood cell with optical tweezers[END_REF]. Cases 5 and 6 are in a slightly better agreement with the experiment, especially regarding the transverse diameter (D T ) in the higher range of imposed stretching force. However, differences between all the modeling cases are contained within the experimental error bars.

It is interesting to note that increasing the resistance to area dilatation of the RBC membrane between case 4 and case 5 (by increasing the ratio C) has only a marginal influence on the predictions of the axial (D A ) and transverse (D T ) diameters, which was also observed in previous works [START_REF] Sigüenza | Characterisation of a dedicated mechanical model for red blood cells: numerical simulations of optical tweezers experiment[END_REF][START_REF] Sinha | Dynamics of a single red blood cell in simple shear flow[END_REF]. In addition, restraining the area variation of the RBC membrane either locally (in case 5) or globally (in case 6) leads to almost identical predictions of the axial (D A ) and transverse (D T ) diameters.

Characterization of the RBC shape

The deformation of the stretched RBC at different stretching forces is displayed in Fig. 5. First, it is seen that the shapes obtained in case 3 differ from the ones obtained in case 2 (see Fig. 3), which also uses the Yeoh law to model the local in-plane elasticity of the RBC membrane. The large fold which appears during the RBC stretching in case 2 is 123 Table 4 Summary of different continuum models of the RBC membrane investigated by means of optical tweezers simulations (see Table 2 for cases 1 and 2) 4 restrained in case 3 by the bending stiffness of the lipid bilayer, modeled by the bending energy (neglected in case 2). The fold is, however, still visible during the stretching, but much smoother. In case 4, when switching the hyperelastic model to the Skalak law, the RBC tends to lose its biconcave shape with increasing stretching. This phenomenon is even more pronounced and faster in case 5, when the area dilatation resistance is increased, leading to a more rounded shape at maximum stretching. Finally, case 6 exhibits a very similar behavior of case 5, with a faster transition from the biconcave to the rounded shape (see shapes at F = 67 pN in Fig. 5), and a more circular shape at maximum stretching. Note that simulations have also been performed combining the Yeoh law with the global area conservation energy, showing the same transition from the biconcave to the rounded shape (not shown). This indicates that this mechanical behavior does not come from the use of the Skalak law itself, but from the area variation restriction of the RBC membrane, achieved either using the Skalak law or the global area conservation energy. 4 In the light of these observations, it appears relevant to introduce two additional lengths measured on the deformed RBC: The in-plane length L P , defined as being the height in the direction perpendicular to the plane of the RBC (see Fig. 5); the folding length L F , also aligned with the direction perpendicular to the plane of the RBC, but evaluated at the fold location (see Fig. 5). As shown in Fig. 6, the discrimination between the different modeling cases is more obvious when analyzing the evolution of the in-plane (L P ) and folding (L F ) lengths than the classical analysis made on the axial (D A ) and transverse (D T ) diameters (in Fig. 4). Previous observations of Fig. 5 can be highlighted: In case 3, the in-plane (L P ) and folding (L F ) lengths show parallel evolutions, meaning that the RBC keeps its biconcave shape for the whole range of stretching force; in case 4, lengths get closer with increasing stretching force, showing that the RBC progressively loses its biconcave shape when subjected to stretching; in cases 5 and 6, a transition from a biconcave folded shape to a rounded shape occurs when the two lengths become identical (for F = 109 pN in case 5, and F = 88 pN in case 6), and the shape of the RBC becomes more and more circular as the lengths increase with the stretching force.

W E S E b Case 3 YE: E s = 7.3 µN/m X κ b = 2.0 × 10 -19 N.m C 3 = E s /30

Area variation

The ability of the quasi-incompressible lipid bilayer to restrain area variations during the RBC deformation is an important mechanical feature of the RBC membrane [START_REF] Mohandas | Red cell membrane: past, present, and future[END_REF]. Figure 7 shows the evolution of the global area variation of the RBC membrane during stretching for the different modeling cases introduced in Table 4. In case 3, the area increase reaches 28%, since the Yeoh law is not designed to restrain area variations of the RBC membrane. Using the Skalak law in case 4 enables to restrain the area variation to a maximum value of 12%. Area variations are even more restrained when increasing the resistance to area dilatation in cases 5 (0.3%) and 6 (0.4%).

Figure 8 shows the local area variation of the RBC membrane for modeling cases 5 and 6. In case 5, the use of the Skalak law with high ratio C allows very small local area variations of the RBC membrane. In case 6, the quasi-incompressibility of the lipid bilayer is independently modeled using the global area conservation energy, whereas the Skalak law with lower ratio C is used to model the own area dilatation resistance of the cytoskeleton. This results in higher local area variations, which correspond to the defor- Biconcave stress-free shape Quasi-spherical stress-free shape Mills experiment Fig. 9 Influence of the stress-free shape of the RBCs on the evolution of the axial (D A ) and transverse (D T ) diameters, for the modeling case 6 mation of the cytoskeleton. In both cases, the maximum local area variations are obtained at the extremities of the cell, near to the bead/RBC contact areas. These regions of high stretching may thus be the locations where the RBC is the most prone to lysis. Note that variation of cytoskeleton area was measured by [START_REF] Discher | Molecular maps of red cell deformation: hidden elastic and in situ connectivity[END_REF] in a micropipette aspiration experiment, but the authors are not aware of similar measurements in optical tweezers experiment.

RBC stress-free shape

Recent studies suggest that RBCs have a quasi-spherical stress-free shape [START_REF] Lim | Stomatocyte-discocyteechinocyte sequence of the human red blood cell: evidence for the bilayer-couple hypothesis from membrane mechanics[END_REF][START_REF] Khairy | Minimum-energy vesicle and cell shapes calculated using spherical harmonics parameterization[END_REF][START_REF] Cordasco | Comparison of erythrocyte dynamics in shear flow under different stress-free configurations[END_REF][START_REF] Peng | Erythrocyte responses in lowshear-rate flows: effects of non-biconcave stress-free state in the cytoskeleton[END_REF][START_REF] Peng | Stability of the tank treading modes of erythrocytes and its dependence on cytoskeleton reference states[END_REF][START_REF] Dupire | A simple model to understand the effect of membrane shear elasticity and stress-free shape on the motion of red blood cells in shear flow[END_REF], meaning that the well-known biconcave shape of the RBCs (Eq. 11) is pre-stressed. This initial pre-stress has not been taken into account so far in the present study, but could eventually play a significant role. [START_REF] Klöppel | A novel two-layer, coupled finite element approach for modeling the nonlinear elastic and viscoelastic behavior of human erythrocytes[END_REF] recently investigated the influence of an initial prestressed biconcave shape of a RBC subjected to stretching deformation, and almost no influence of this initial pre-stress was observed. They concluded that when investigating static deformation of RBCs, the biconcave initial shape of the RBCs can be assumed as being stress-free.

In this section, the influence of the quasi-spherical stressfree shape of the RBCs is investigated. Figure 9 compares previous simulations of the modeling case 6, assuming a biconcave stress-free shape of the RBC, with simulations where the stress-free shape of the RBC is a quasi-spherical shape having a reduced volume V /V 0 = 0.98 (with V 0 the volume of a sphere having the same surface area). As usually done when modifying the stress-free shape [START_REF] Cordasco | Comparison of erythrocyte dynamics in shear flow under different stress-free configurations[END_REF][START_REF] Peng | Erythrocyte responses in lowshear-rate flows: effects of non-biconcave stress-free state in the cytoskeleton[END_REF], the spontaneous curvature c 0 (Eq. 7) is adjusted so that the equilibrium shape is similar to the para- Biconcave stress-free shape Quasi-spherical stress-free shape Fig. 10 Influence of the stress-free shape of the RBCs on the evolution of the in-plane (L P ) and folding (L F ) lengths, for the modeling case 6 metric biconcave shape described by Eq. ( 11). In the present case, the spontaneous curvature is set to c 0 = 4.6 × 10 6 m -1 . Consistent with the conclusions made by [START_REF] Klöppel | A novel two-layer, coupled finite element approach for modeling the nonlinear elastic and viscoelastic behavior of human erythrocytes[END_REF], it is seen that there is no significant effect of the stress-free shape regarding the evolution of the axial (D A ) and transverse (D T ) diameters. Regarding the latter quantity, a slightly better agreement with the experiment is nonetheless observed when considering a quasi-spherical stress-free shape.

The evolution of the in-plane (L P ) and folding (L F ) lengths is displayed in Fig. 10, showing a more significant effect of stress-free shape. Indeed, it is seen that the in-plane length (L P ) is higher for the quasi-spherical stress-free shape before the transition from the biconcave to the rounded shape occurs.

Discussion

In the present paper, the optical tweezers experiment by Mills et al. ( 2004) is simulated using a numerical method dedicated to the simulation of the dynamics of RBCs under flow. A computational setup for simulating the RBC stretching is presented, which is seen to perfectly reproduce the numerical results obtained by [START_REF] Mills | Nonlinear elastic and viscoelastic deformation of the human red blood cell with optical tweezers[END_REF]. Influence of the RBC membrane modeling is then investigated, introducing different continuum models to describe the membrane mechanics.

Comparison of the numerical results with the forceextension data provided by the experiment (i.e., the axial (D A ) and transverse (D T ) diameters of the stretched RBC) shows that all the considered modeling approaches are able to reproduce the mechanical response of the RBC subjected to optical stretching (see Fig. 4). An adjustment of the shear modulus E s is, however, required depending if the RBC membrane is described using the Yeoh law or the Skalak law (E s is twice smaller when using the Skalak law). It is also seen that some of these models allow non-physiological area variations of the RBC membrane during stretching (see Fig. 7), especially the Yeoh law which was considered in previous works as a suitable model of the RBC membrane [START_REF] Mills | Nonlinear elastic and viscoelastic deformation of the human red blood cell with optical tweezers[END_REF][START_REF] Suresh | Connections between single-cell biomechanics and human disease states: gastrointestinal cancer and malaria[END_REF]. Consistent with the findings of [START_REF] Dimitrakopoulos | Analysis of the variation in the determination of the shear modulus of the erythrocyte membrane: effects of the constitutive law and membrane modeling[END_REF], this indicates that the Yeoh law should not be used to describe the mechanical behavior of the RBC membrane. This also indicates that the single analysis of the axial (D A ) and transverse (D T ) diameters of the stretched RBC is not sufficient for characterizing the mechanics of the RBC membrane, and cannot be used alone to validate numerical models of the RBC membrane.

Detailed analysis of the shape of the stretched RBC reveal different behaviors among the investigated models (see Fig. 5). A transition of the RBC shape from a biconcave folded shape to a rounded shape is observed when restraining the area variations of the RBC membrane, either locally or globally. This observation may be due to the fact that the RBC tends to lose its biconcave shape when subjected to optical stretching, to prevent area variations of the RBC membrane. Note that such ellipsoidal shapes were also reported in previous numerical studies [START_REF] Li | Spectrin-level modeling of the cytoskeleton and optical tweezers stretching of the erythrocyte[END_REF][START_REF] Klöppel | A novel two-layer, coupled finite element approach for modeling the nonlinear elastic and viscoelastic behavior of human erythrocytes[END_REF][START_REF] Farutin | 3D numerical simulations of vesicle and inextensible capsule dynamics[END_REF][START_REF] Sigüenza | Characterisation of a dedicated mechanical model for red blood cells: numerical simulations of optical tweezers experiment[END_REF].

This transition from a biconcave folded shape to a rounded shape can be characterized by introducing two additional measurements in the direction perpendicular to the plane of the RBC: The in-plane length L P and the folding length L F (see Fig. 6). Experimental measurements of such lengths could thus be of prime interest to make the optical tweezers experimental setup more helpful to characterize the mechanics of the RBC membrane. Indeed, these quantities reveal to be more sensitive to the area variation restriction of the RBC membrane than the usual force-extension data and could thus enable to better investigate the mechanical behavior of the membrane. Such data are also expected to be sensitive to the bending stiffness of the RBC membrane. However, the latter is seen to mainly influence the shape of the fold of the stretched RBC, but is not at the origin of the transition from the biconcave to the rounded shape. Indeed, this transition may occur even when the bending stiffness of the RBC membrane is not considered [START_REF] Sigüenza | Characterisation of a dedicated mechanical model for red blood cells: numerical simulations of optical tweezers experiment[END_REF]. The lengths L P and L F are thus expected to be of interest to qualitatively challenge RBC modeling, while quantitative comparisons should account for the possible influence of the bending stiffness of the RBC membrane.

More sophisticated measurements of the shape of the stretched RBC must, however, be performed with reasonable experimental uncertainties. One of the main sources of uncertainty is expected to come from the contact areas between the beads and the RBC, which may vary from one experiment to the other. In the present computational setup, these contact areas are defined by the contact size d c which is initially chosen to be d c = 2 µm, as in the computations of [START_REF] Mills | Nonlinear elastic and viscoelastic deformation of the human red blood cell with optical tweezers[END_REF]. Figure 11 shows the influence of this contact size on the numerical predictions of the axial (D A ) and transverse (D T ) diameters of the stretched RBC (using the modeling case 6), when the contact size is successively set to d c = 1 µm, d c = 2 µm and d c = 3 µm. It is seen that the contact size strongly influence the prediction of the axial diameter (D A ), showing a more rigid behavior with increasing d c , but has no influence on the prediction of the transverse diameter (D T ). This may explain the large and increasing error bars obtained by [START_REF] Mills | Nonlinear elastic and viscoelastic deformation of the human red blood cell with optical tweezers[END_REF] in the experimental measurements of the axial diameter (D A ), as compared to the smaller and monotonous error bars obtained for the transverse diameter (D T ). This finding sheds doubt on the meaningfulness of the use of the axial diameter (D A ) for the determination of the shear modulus E s . Indeed, the choice of the contact size d c may strongly influence the determined value of the shear modulus E s . Conversely, simulations should rather be fitted to the transverse diameter (D T ) which is less sensitive to the choice of this contact size.

Figure 12 shows that the contact size d c has only a little influence on the predictions of the in-plane (L P ) and folding (L F ) lengths, which means that comparison between computed and measured values of these quantities would be robust to the uncertainties related to the bead/RBC contact areas. The authors hope that these findings will arouse an interest for updated optical tweezers experiments.

  Fig. 1 a Illustration of the experimental setup of Mills et al. (2004). The axial (D A ) and transverse (D T ) diameters of the stretched RBC are measured. b Computational setup used to simulate the optical tweezers

Fig. 2

 2 Fig. 2 Axial (D A ) and transverse (D T ) diameters of the RBC stretched by optical tweezers. Comparison with the experimental and numerical data from Mills et al. (2004). a The RBC membrane is assumed to follow the neo-Hookean law, corresponding to case 1. b The RBC membrane is assumed to follow the Yeoh law, corresponding to case 2

Fig. 4

 4 Fig. 4 Comparison of the axial (D A ) and transverse (D T ) diameters of the RBC stretched by optical tweezers for the different modeling cases introduced in Table 4

Fig. 5 Fig. 6

 56 Fig. 5 Visualization of the red blood cell deformation over the entire range of stretching force, for the different modeling cases introduced in Table 4. Only half of the cell is displayed

Fig. 7

 7 Fig. 7 Global area variation of the RBC membrane for the different modeling cases introduced in Table 4

Fig. 8

 8 Fig. 8 Comparison of the local area variations of the red blood cell membrane at the maximum stretching force F = 193 pN, for the modeling approaches of cases 5 and 6

Fig. 11 Fig. 12

 1112 Fig. 11 Influence of the bead/RBC contact areas on the axial (D A ) and transverse (D T ) diameters, for the modeling case 6

Table 3

 3 Influence of the mesh resolution for case 2, at the maximum imposed force of 193 pN

		D A (µm)	D T (µm)
	Mills simulation	16.14	4.90
	Coarse mesh	15.92	4.94
	Reference mesh	15.93	4.81
	Fine mesh	15.93	4.72
	Case 1	Case 2
	0 pN
	67 pN
	130 pN
	193 pN

Acknowledgements V. Moureau and G. Lartigue from the CORIA lab, and the SUCCESS scientific group are acknowledged for providing the YALES2 solver which constitutes the basis of the YALES2BIO tool.

Funding This study was performed with supports from ANR (FORCE project ANR-11-JS09-0011), from BPIfrance (DAT@DIAG Project No. I1112018W) and from the NUMEV Labex (ANR-10-LABX-20).

Compliance with ethical standards

Conflict of interest

The authors declare that they have no conflict of interest.