Maxime Joos

Chengjie Ding

Vivien Loo

Guillaume Blanquer

Elisabeth Giacobino

Alberto Bramati

Valentina Krachmalnicoff

Quentin Glorieux
email: quentin.glorieux@lkb.upmc.fr

Polarization Control of Linear Dipole Radiation Using an Optical Nanofiber

teaching and research institutions in France or abroad, or from public or private research centers.

I. INTRODUCTION

The wave nature of light implies that radiation phenomena are the combined manifestation of a source and its environment. Thus, properties of light sources such as power, frequency, or directivity can be tuned by confining the source in an appropriate environment. Nanophotonics aims at controlling these properties using the interaction with wavelength-scaled structures [START_REF] Benson | Assembly of hybrid photonic architectures from nanophotonic constituents[END_REF] such as cavities, antennas [START_REF] Novotny | Antennas for light[END_REF] or waveguides [START_REF] Guo | High-bit rate ultra-compact light routing with mode-selective on-chip nanoantennas[END_REF]. These structures rely on the discretization of the optical modes when the spatial confinement is of the order of the wavelength but until recently, they did not explicitly make use of the vectorial nature of light. In the last decade, exploiting the strong confinement of vectorial fields has proven to be a very rich aspect of nanophotonics. Indeed, the sharp transverse confinement of vectorial fields breaks the symmetry of the intensity profile and modifies the polarization structure of the light [START_REF] Le Kien | Field intensity distributions and polarization orientations in a vacuum-clad subwavelength-diameter optical fiber[END_REF]. When the transverse electric field varies significantly over a wavelength -as it is the case in evanescent waves -a longitudinal field component emerges and oscillates in phase quadrature with the transverse field. The local polarization becomes dependent on the propagation direction of light, an effect referred to as spin-momentum or polarization-direction locking of light [START_REF] Van Mechelen | Universal spin-momentum locking of evanescent waves[END_REF] [START_REF] Lodahl | Chiral quantum optics[END_REF]. In this context, directional emission [START_REF] Francisco | Near-field interference for the unidirectional excitation of electromagnetic guided modes[END_REF][START_REF] Petersen | Chiral nanophotonic waveguide interface based on spin-orbit interaction of light[END_REF][START_REF] Söllner | Deterministic photonemitter coupling in chiral photonic circuits[END_REF], non-reciprocal light-matter interaction [START_REF] Sayrin | Nanophotonic optical isolator controlled by the internal state of cold atoms[END_REF] and developments in nano-polarization devices have been achieved [START_REF] Francisco | Universal method for the synthesis of arbitrary polarization states radiated by a nanoantenna[END_REF][START_REF] Espinosa-Soria | On-chip optimal stokes nanopolarimetry based on spinorbit interaction of light[END_REF].

Optical nanofibers are privileged interfaces to investigate the interaction of spin-momentum locked light with matter. In such cylindrical waveguides, the light extends outside the fiber in the form of an evanescent field and enables easy interfacing with emitters such as atoms, quantum dots [START_REF] Yalla | Efficient channeling of fluorescence photons from single quantum dots into guided modes of optical nanofiber[END_REF] or plasmonic scatterers [START_REF] Sadgrove | Polarization response and scaling law of chirality for a nanofibre optical interface[END_REF].

In this work, we exploit the strong confinement of light in optical nanofibers to demonstrate that a linear (1D) dipole can emit elliptically polarized light in the transverse plane when coupled to spin-momentum locked modes. We study how the guided light polarization depends on the local orientation and position of the dipole and further show that the whole Poincaré sphere is accessible with this nanophotonics system. For example, we observed emission of quasi-circularly polarized light by a linear dipole oriented at 46°with respect to the propagation direction of light. We remark that this effect does not involve any birefringent component. This system constitutes a step toward a complete control of light emission at the nanoscale.

II. EXPERIMENTAL SYSTEM

A. Setup Our experimental physical system is composed of a single gold nanorod deposited on the surface of an air-clad optical nanofiber as shown in figure 1a. A linearly polarized, focused laser beam illuminates the nanorod particle that behaves as a linear dipole scatterer. The polarization of the light collected in the guided mode of the nanofiber is analyzed as a function of the position of the nanorod on the fiber and of its orientation, which are measured by means of a Scanning Electron Microscope (SEM).

We use optical nanofibers produced by a standard heat-and-pull process (see [START_REF] Hoffman | Ultrahigh transmission optical nanofibers[END_REF]) in order to locally stretch a commercial single-mode fiber to a subwavelength-diameter of around 300 nm while maintaining over 95% optical transmission. At our working wavelength of 637 nm, the nanofiber only guides the fundamental mode HE 11 studied in [START_REF] Le Kien | Field intensity distributions and polarization orientations in a vacuum-clad subwavelength-diameter optical fiber[END_REF].

B. Model

We now consider in more detail how an electric dipole moment couples to the fundamental mode of the nanofiber. The dipole is assumed to lie flat on the nanofiber surface at an azimuthal position denoted by the angle α and forms an angle θ with respect to the propagation axis z. We use a coordinate system (x y z) such that the dipole is located on the y -axis, i.e at x = 0. The x and y axes are obtained by rotation of axes x and y by an angle α around the z-axis (see figure 1). We will express the guided electric field E g emitted by the dipole as a linear combination of the quasi-linearly polarized hybrid modes envelopes [START_REF] Le Kien | Field intensity distributions and polarization orientations in a vacuum-clad subwavelength-diameter optical fiber[END_REF][START_REF] Reitz | Nanofiber-based double-helix dipole trap for cold neutral atoms[END_REF]:

HE 11,x = e x (r, ϕ) + ie z (r) cos ϕ, HE 11,y = e y (r, ϕ) + ie z (r) sin ϕ (1)
where r and ϕ are the polar coordinates in the x y -plane, e x (r, ϕ), e y (r, ϕ) are the real-valued transverse components and e z (r) cos ϕ, e z (r) sin ϕ are the longitudinal components of the guided modes with main polarization along x and y respectively (see Appendix A for details). The guided field can then be expressed as:

E g = [A HE 11,x + B HE 11,y ] e i(βz-ωt) (2)
FIG. 1. Experimental setup. A single gold nanorod is deposited on the surface of an optical nanofiber. A focused laser beam contained in the yz-plane illuminates the particle that scatters light in the fundamental mode of the nanofiber. The nanorod lies in a plane tangent to the nanofiber surface. The normal to this plane (y ) makes an angle α with the vertical (y) axis. In this plane the nanorod forms an angle θ with the nanofiber (z) direction. The guided light is then analyzed with a polarimeter allowing to measure Stokes parameters.

where β is the propagation constant and ω the angular frequency of the light field.

E g = [C sin θ HE 11,x + i D cos θ HE 11,y ] e i(βz-ωt) (3)
where C and D are real (see Appendix B). One first consequence of the coupling of the dipole to the nanofiber is that, whatever the orientation θ of the dipole, the guided field (3) never vanishes. When aligned along the nanofiber (θ = 0°), the dipole can even radiate light in its axis direction -thanks to the nanofiber guided modes-in strong contrast to the dipole radiation pattern in free-space. Considering now the emitted polarization in the nanofiber: the terms of the right-hand side of (3) oscillate in phase quadrature and give rise in general, to elliptical polarization in the transverse plane. This can be understood as follow: the dipole has non-zero overlap only with the longitudinal z -component of HE 11,y and with the transverse x '-component of HE 11,x . The dipole excites HE 11,y through its longitudinal component and HE 11,x through its (transverse) x '-component. For spin-momentum locked light however, transverse and longitudinal components are in phase quadrature so that the hybrid modes HE 11,x and HE 11,y will also oscillate in phase quadrature.

Let us now illustrate this mechanism with two cases: when the dipole is aligned (θ = 0°) or perpendicular (θ = 90°) to the nanofiber. For θ = 0°, the electric field (3) of the guided light reduces to E g ∝ HE 11,y which corresponds to quasi-linearly polarized light along y . For θ = 90°, light is quasi-linearly polarized along x . Between this two limiting cases, the continuous range of polarization ellipses with axes along x and y is accessible. A remarkable intermediate dipole orientation, for which |A| = |B| gives rise to quasi-circularly polarized light. This can be numerically solved and leads to θ circ = ±43°for typical experimental conditions and a nanofiber diameter of 305 nm. The dipole orientation hence defines the ellipticity of the guided light.

The azimuth α of the dipole on the other hand determines the major axis direction of the polarization ellipse. This can be easily seen considering the situation θ = 0°. As already mentioned, the guided light polarization is quasi-linear along the y direction which is, by definition, the azimuth of the dipole.

Hence, all polarization states are deterministically accessible through the following mapping: the nanorod azimuthal position α ∈ [-90°, 90°] determines the polarization ellipse orientation ψ or equivalently the longitude 2α on the Poincaré sphere as shown on figure 2. The nanorod orientation θ ∈ [-θ circ , θ circ ] defines the polarization ellipticity or equivalently the latitude on the Poincaré sphere according to a non trivial mapping f (θ) ≈ 90°θ circ θ.

III. EXPERIMENT A. Polarization measurement

Direct measurement of the polarization, i.e. determination of the Stokes parameters in the nanofiber region is hard to implement and we rather choose to measure the polarization state after light exits from one end of the fiber. While light propagates in the transition region, the adiabatic tapering enables almost perfect conversion from the fundamental mode of the nanofiber to the LP 01 mode of the standard fiber. The polarization however is not maintained in practice, because of the significant birefringence in fibers. In order to still map the hybrid mode basis onto an accessible paraxial mode basis, we use a uniaxial birefringent plate mounted in the form of a Bereck compensator to compensate the birefringence effect of the fiber (see Appendix C).

B. Nanorods as linear dipoles

In our experiment, we use single gold nanorods that, under suitable illumination, behave as dipole scatterers with dipole orientations along the rods longitudinal axes. The deposition of a single nanorod on the nanofiber surface is performed by touching it with a small droplet of a commercial colloidal particle dispersion (Nanocomposix); this technique has an accuracy of a few hundred microns in the position of the particles along the nanofiber. The nanorods were chosen so that their longitudinal surface plasmon resonance (LSPR) matches our laser wavelength of 637 nm when deposited on the nanofiber. We hence work with nanorods of average aspect ratio of 2.6 yielding a LSPR centered around our working wavelength and with a size of 17 ± 1 nm in diameter and 45 ± 6 nm in length, considered small compared to the nanofiber geometry.

We model a nanorod as an anisotropic scatterer on which an exciting electric field E exc = (E L , E T ) induces a dipole moment p = (α L E L , α T E T ) where α L and α T are the longitudinal and transverse polarizabilities of the rod, respectively. Close to the LSPR, Im(α L )

Im(α T ) and the transverse component of the dipole moment is strongly suppressed compared to the longitudinal component. We can hence suppose that, as long as the excitation field does not oscillate too perpendicularly with respect to the rod, in which case α L E L might be comparable to α T E T , the dipole is induced along the nanorod.

To verify this crucial assumption experimentally, we rotated the excitation polarization while monitoring the power scattered by the nanorod into the nanofiber. We observed a sinusoidal dependence corresponding to a Malus law as shown in figure 3. The beam polarization angle for which the scattered signal is maximum is referred to as χ max . In contrast to the collected power into the nanofiber, the guided polarization is almost unchanged when turning the excitation polarization. This validates the assumption that the dipole tends to align along the rod, independently of the excitation polarization. In the following experiments, we will always maximise the scattered signal by adjusting the beam polarization orientation to the value χ max prior to acquiring the guided polarization state. The intrinsic diffusion of the excitation beam by the bare nanofiber is lower than 0.5 % of the overall scattered signal and can hence be neglected.

FIG. 3. Nanorod optical response to rotating linear polarization excitation. The lower graph represents the normalised power in the nanofiber as the excitation polarization angle χ rotates over 80°, starting at χ = χmax for which the scattered signal is maximum. The measured collected power is fitted to a Malus law. Inset: SEM image of the measured nanorod on the nanofiber. The upper plots represent the measured polarization ellipses for different excitation polarization angles. The guided polarization is almost independent of χ which constitutes a consistency check of our assumption that the dipole tends to form along the rod.

C. Experimental procedure

The experimental demonstration of the dependence between the nanorod geometry and the polarization of the guided light requires an ensemble of nanorods with various different positions and orientations. The experimental procedure is the following: (1) We deposit a first particle, (2) we compensate the birefringence of the fiber portion between the particle and the detection system, (3) we align the excitation polarization along the nanorod and (4) we record the guided polarization state. Next to the first particle, we then deposit a second particle and iterate the experimental procedure for the second particle. After that, we deposit a third particle and so on.

Since they are several particles on the same nanofiber, the light scattered by the nanorod under investigation and guided through the nanofiber towards the detection system could be scattered by the other particles, which could alter the polarization outcome. To prevent this issue, each new particle is deposited on the nanofiber between the previously deposited one and the detection system. This ensures that the light emitted by the last deposited nanorod encounters no obstacles when guided towards the detection system. A special care is also dedicated to have sufficient space between the different particles in order for the excitation beam to address nanorods individually. We also record the exact relative positions of the particles with respect to each other. This is meant to facilitate the finding of the nanorods when we ultimately observe the nanofiber with a SEM.

On the SEM, we identify the various depositions, disregard the clusters or the rods with odd shape to record only the azimuth α and the orientation θ of proper single nanorods. The orientation is determined from analysing SEM images as the one shown in the inset of figure 3. From this specific image, we extracted θ = 5 ± 3°. The uncertainty associated to the orientation is typically around a few degrees and is mainly due to the quality of the SEM images. We note that experimentally, we never observed rods with orientation |θ| 60°which we understand as the tendency for the nanorod to maximise its contact surface with the nanofiber.

D. Results & discussion

In our model of a linear dipole coupled to nanofiber modes, we predicted that the dipole orientation determines the polarization ellipticity. In order to express this dependence and compare our measurement to the model, we consider the normalised Stokes parameter S 3 as a function of the dipole/nanorod orientation θ as shown in figure 4a. S 3 ranges from -1 to 1 and provides us immediate information on the degree of circular polarization -defined as |S 3 | -and on the sense of circulation of the electric field in the transverse plane, related to the sign of S 3 (see Appendix D for a discussion of the uncertainties associated to the measured Stokes parameters). The figure also shows the measured corresponding polarization ellipses for five characteristic nanorods, the orientation of which varies from -40°to +46°. As input parameters for the model, we took the average nanofiber diameter extracted from SEM images (2a = 305 ± 7 nm), the working wavelength of 637 nm, 1.457 and 1.000 for the refractive index of silica and air respectively and we assumed 9 nm, which corresponds to the radius of a nanorod for the distance of the dipole to the nanofiber surface.

Our measurements of S 3 fit well our dipole model and present the characteristic behaviour presented above namely that, in the range |θ| θ circ , the degree of circular polarization grows with |θ| as the two components in (3) reach comparable amplitude. The case of circular polarization is expected at θ circ = ±43°and experimentally, we observed almost purely circularly polarized light for a nanorod oriented at θ = 46 ± 3°. In contrast, nanorods nearly aligned along the nanofiber (small θ) induce a low degree of circular polarization and hence a narrow polarization ellipse. We also notice that the measured handedness of the polarization fits to our model as positive (negative) θ gives rise to counter-clockwise (clockwise) circulation of the electric field in the transverse plane.

We now consider the relation between the azimuthal position α of the nanorod and the polarization ellipse orientation ψ. This is easily observable for rods giving rise to narrow polarization ellipses, i.e. for which the orientation θ forms a small angle with the nanofiber. We hence consider five rods with orientations θ lying between -20°and 20°. Figure 4b shows the measured ellipse orientations as a function of the nanorod azimuth. The measured orientations of the polarization ellipses tend to follow the azimuthal positions of the nanorods.

Thus our experiments show that the measured polarization states as a function of the orientation and position of the nanorod are in good agreement with our linear dipole model. This constitutes an experimental demonstration that this system can generate all possible polarization states.

This realisation opens the possibilities for polarization control of single photon emitters coupled to nanofibers. Single photon sources with fully controllable polarization are of crucial importance for quantum cryptography protocols and a device based on our prototype can fulfil the requirements for this major application.

Indeed, the coupling of single photon sources with nanofibers has been reported [START_REF] Yalla | Efficient channeling of fluorescence photons from single quantum dots into guided modes of optical nanofiber[END_REF][START_REF] Liebermeister | Tapered fiber coupling of single photons emitted by a deterministically positioned single nitrogen vacancy center[END_REF], as well as the polarized single photon emission of colloidal quantum dots [START_REF] Vezzoli | Exciton fine structure of cdse/cds nanocrystals determined by polarization microscopy at room temperature[END_REF][START_REF] Manceau | Effect of charging on cdse/cds dot-in-rods single-photon emission[END_REF]. Additionally, orientation of nanoparticules has been achieved using optical tweezers [START_REF] Onofrio | Optical trapping and manipulation of nanostructures[END_REF]. Our results bring the last step for the realisation of a potential device with tunable polarization by showing the direct mapping between the spatial degrees of freedom of the source and the photon polarization.

Combining this mapping with an optical trapping and alignment technique [START_REF] Tong | Alignment, rotation, and spinning of single plasmonic nanoparticles and nanowires using polarization dependent optical forces[END_REF], it is therefore possible to align the emitter to a desired orientation by sending inside the nanofiber a tweezer beam with a controlled polarization during the deposition procedure [START_REF] Tong | Alignment, rotation, and spinning of single plasmonic nanoparticles and nanowires using polarization dependent optical forces[END_REF][START_REF] Pelton | Optical trapping and alignment of single gold nanorods by using plasmon resonances[END_REF]. Finally, once the emitter is correctly aligned, the same mapping (but in the reciprocal way) guarantees the desired polarization at the output of the fiber, as shown previously.

IV. CONCLUSION

In conclusion, we used gold nanorods deposited on a nanofiber to demonstrate that a linear dipole can radiate elliptically polarized light when coupled to spin-momentum locked modes. We showed that this system can, in principle, generate all possible polarization states. The azimuth position of the rod controls the inclination of the polarization ellipse while its orientation with respect to the nanofiber defines the ellipticity. This constitutes a demonstration of the mapping between purely geometrical degrees of freedom of a light source and polarization states. This system opens a new way for controlling the polarization of light sources at the nanoscale without involving birefringent components or magnetic fields.

Random errors

Random errors are dominated by fluctuations of the power scattered by the nanorod, principally due to vibrations of the nanofiber holding the nanorod in the focus of the excitation laser beam and, to a smaller extent, to intrinsic laser intensity fluctuations. These power fluctuations introduce a dispersion on the measured normalised Stokes parameters as low as σ = 0.01.

Systematic error

The uncertainty on the measured state of polarization of the light emitted by the nanorod is in fact largely dominated by a systematic error. It comes from the non-perfect birefringence compensation described above. Ideally, once the Berek compensator has been adjusted correctly, quasi-linearly polarized light in the nanofiber region should be mapped onto linearly polarized light at the polarimeter location. In practice, we observe a residual ellipticity, which implies a systematic error on the measured normalized Stokes parameters. Typical errors on the normalized Stokes parameters range from 0.01 to 0.24 with an average of 0.1.

  A and B are complex coefficients determined by the projection of the dipole moment amplitude d onto the complex envelope of the hybrid modes (1) evaluated at the dipole position (r = a, ϕ = π/2) [8, 17], where a is the nanofiber radius: A = d * • HE 11,x (a, π/2) and B = d * • HE 11,y (a, π/2). For a linear dipole expressed in the (x , y , z) coordinate system, d is proportional to (sin θ, 0, cos θ). It follows that A = d • e x (a, π/2) reduces to a single real term and B = d • [e y (a, π/2) + ie z (a)] reduces to a purely imaginary term since e y (a, π/2) has no x -component and d no y -component, leading to the electric field:

FIG. 2 .

 2 FIG. 2. Mapping between dipole geometry and emitted polarization. (a) Cut of the nanofiber. The azimuth of the dipole α defines the orientation of the polarization ellipse ψ = α. (b) Poincaré sphere: a dipole with geometrical parameters (α, θ) gives rise to a guided polarization represented by a point with coordinate (2α, f (θ)) on the Poincaré sphere. V, 45°and LCP represent vertically linearly, linearly at ψ = 45°and left-circularly polarized light respectively.

FIG. 4 .

 4 FIG. 4. Measured state of polarization as a function of the orientation θ and azimuth α of the nanorod. (a) Ellipticity as a function of θ. Graph : measured Stokes parameter S3 for 9 different nanorods compared to the expected dipole model. We plot the measured polarization ellipse for characteristic nanorods. (b) Measured ellipse orientation ψ for 5 nanorods with different azimuth. The inset represents a cut of the nanofiber showing the azimuth of a specific nanorod and the associated polarization ellipse. (c) Poincaré sphere representation of the measured polarization for 4 nanorods from (b).

V. ACKNOWLEDGMENTS

The authors thanks Clément Sayrin for his helpful suggestions throughout this work, Imène Estève for her decisive support on the SEM platform and Sébastien Bidault for fruitful discussions. This work is supported by the Emergence program from Ville de Paris, by PSL Research University in the framework of project COSINE and by the Cai Yuanpei program.

Appendix A: Hybrid modes An optical nanofiber in the single mode regime only guides the fundamental mode associated to a unique propagation constant β where n 1 k > β > n 2 k, k is the wavenumber of light in vacuum, n 1 and n 2 are the refractive indices of the fiber and surrounding medium respectively. General expressions for the electric field of the fundamental mode can be found in [START_REF] Le Kien | Field intensity distributions and polarization orientations in a vacuum-clad subwavelength-diameter optical fiber[END_REF]. We consider here the orthogonal basis of quasi-linearly polarized modes along x and y where we use the cartesian coordinate system (x y z) introduced in the main text. For the region outside of the nanofiber (r > a) and light propagating in the positive z-direction, the complex enveloped of the modes are:

where r, ϕ are the polar coordinates in the x y -plane. N is a normalization constant, a is the nanofiber radius,

the Bessel functions of the first kind and the modified Bessel functions of the second kind respectively, and the prime stands for the derivative.

Appendix B: Modal amplitudes for a linear dipole coupled to HE11

At the dipole (nanorod) position (r = ρ, ϕ = π/2), the hybrid modes take the value:

where

K1(qa) K 1 (qρ) and ρ is the radial position of the dipole. The dot products between the dipole moment d * = d [sin θ, 0, cos θ] and the hybrid modes read:

where d is real. Emission of quasi-circularly polarized light occurs when |A| = |B| which leads to:

This can be numerically solved and leads to θ circ = ±43°for our experimental parameters: ρ = a+9 nm, 2a = 305 nm, λ = 637 nm, n 1 = 1.457 and n 2 = 1.000. At the output of the fiber, i.e. at the detection system, the normalised Stokes parameters are given by:

where

, it is equal to 1 for purely polarized light. In our experiment, we measured an average degree of polarization of 0.98 ± 0.01 for the guided light emitted by the nanorod. S 1 and S 2 are dependent of the coordinate system whereas S 3 is rotation invariant. In the coordinate system (x y ) introduced in the main text, S 2 = 0.

In the main text, we showed that a dipole with geometrical parameters (α, θ) gives rise to a guided polarization represented by a point with coordinate (2α, f (θ)) on the Poincaré sphere. For purely polarized light and θ ∈ [-θ circ , θ circ ], we have:

where the error made with the approximation does not exceed 3% in typical experimental conditions.

Appendix C: Birefringence compensation

In order to find the adjustment of the Berek compensator that compensates the birefringence introduced by the fiber setup, we observe the Rayleigh scattering [START_REF] Vetsch | Nanofiber-based optical trapping of cold neutral atoms[END_REF] in the nanofiber region when light is sent from the polarimeter toward the nanofiber, through the compensator and the fiber setup. The light scattering is assumed to originate from imperfections in the nanofiber. Because Rayleigh scattering preserves the polarization, its analysis can give information on the polarization of the propagating fundamental mode in the nanofiber. Especially, Rayleigh scattering enables to align the polarization in the nanofiber quasi-linearly along the axis of observation or orthogonal to it. Since we can recognize an orthogonal hybrid basis in the nanofiber region, the compensation procedure consists in finding the Berek compensator adjustment that map the orthogonal polarization basis at the polarimeter location with the hybrid basis in the nanofiber.

Appendix D: Uncertainty on the measured Stokes parameters

The measurement of the state of polarization of the guided light emitted by the nanorod is subject to random errors and systematic errors.


Joao Lebel

A Trovao

R German
email: ronan.german@univ-lille1.fr

N Solis

L Reyes

L Silva

F.-A Lebel

João P Trovão

A Bouscayrol

Coupled Electric and Thermal Batteries Models using Energetic Macroscopic Representation (EMR) for Range Estimation in Electric Vehicles

HAL is

Introduction

Electric Vehicles (EVs) are vehicles of the future. But they suffer of a lack of autonomy, a long charging time, and a high cost [START_REF] Chan | Overview of electric, hybrid and fuel vell vehicles[END_REF]. The range anxiety is key issue for EV's drivers [START_REF] Neubauer | The impact of range anxiety and home, workplace and public charging infrastructure on simulated battery EV lifetime utility[END_REF]. Therefore, estimation of the correct driving range is an important matter.

In some cold winter countries, EV manufacturers noticed reduced driving range (autonomy) when temperature goes under 0 °C. For example, a Canadian study dealt with the influence of low temperature on battery with sub-zero temperatures. The result was a capacity decrease up to 35 % at -20 °C compared to 25 °C [START_REF] Jaguemont | Low Temperature Discharge Cycle Tests for a Lithium Ion Cell[END_REF]. Other electric parameters (such as series resistance for instance…) also depend on the temperature [START_REF] Lin | Online Parameterization of Lumped Thermal Dynamics in Cylindrical Lithium Ion Batteries for Core Temperature Estimation and Health Monitoring[END_REF]. These variations of battery behavior are non-permanent and the original driving range is reached when the temperature goes higher.

Moreover, cold and hot weathers both have an important effect on battery ageing rate. The result of ageing is a permanent decrease of the discharge capacitance (reduced driving range) and an increase of the battery resistance (increase of Joule losses). For subzero temperature, the ageing rate is drastically increased [START_REF] Jaguemont | Lithium-Ion Battery Aging Experiments at Subzero Temperatures and Model Development for Capacity Fade Estimation[END_REF]. For temperatures from 10 °C to 60 °C it is commonly considered that the ageing rate doubles when the temperature increases by 10 °C [START_REF] Baghdadi | Lithium battery aging model based on Dakin's degradation approach[END_REF]. Thus, estimating the temperature of the battery during one cycle is also very important for long term driving range.

The goal of this paper is to build an Electro-Thermal (ET) model of the battery with electric parameters thermal dependence to better estimate the driving range of the EVs under various temperature conditions. This model is intended to be used for vehicle simulation and hardware in the loop for future works. Therefore, it must be as simple as possible to limit the calculation time. Some papers have already proposed ET models focused on the battery impedance parameters evolution. For example, papers [START_REF] Lin | Online Parameterization of Lumped Thermal Dynamics in Cylindrical Lithium Ion Batteries for Core Temperature Estimation and Health Monitoring[END_REF] and [START_REF] Lin | A lumped-parameter electro-thermal model for cylindrical batteries[END_REF] include a temperature dependence on the series resistance of the battery. In this paper, the capacity thermal dependence (based on the experimental results of [START_REF] Jaguemont | Low Temperature Discharge Cycle Tests for a Lithium Ion Cell[END_REF]) is also considered in order to better estimate the EV range.

First the classic electric and thermal models of Li-ion battery are presented. Then, they are associated using Energetic Macroscopic Representation (EMR) [START_REF] Bouscayrol | Graphic Formalisms for the Control of Multi-Physical Energetic Systems: COG and EMR[END_REF]. As EMR has already been used to model traction systems of EVs and HEVs (Hybrid Electric Vehicles) [START_REF] Castaings | Real-time energy management strategies of a battery/supercapacitors system for Electric Vehicle under limitations[END_REF]- [START_REF] Chen | Energetic Macroscopic Representation and Inversion-based control: Application to an Electric Vehicle with an electrical differential[END_REF], such a formalism will enable an easy integration of this ET battery model for vehicle simulation. The experimental thermal dependence of the capacity and the series resistance of the battery is showed through literature results [START_REF] Jaguemont | Low Temperature Discharge Cycle Tests for a Lithium Ion Cell[END_REF], [START_REF] Lin | Online Parameterization of Lumped Thermal Dynamics in Cylindrical Lithium Ion Batteries for Core Temperature Estimation and Health Monitoring[END_REF]. The electric and thermal models are coupled using the EMR formalism. The electric parameters of the model are calculated with the thermal dependence extracted from the literature results. In a second time the model is validated by experimental results [START_REF] Lin | Online Parameterization of Lumped Thermal Dynamics in Cylindrical Lithium Ion Batteries for Core Temperature Estimation and Health Monitoring[END_REF]. Finally, the range of an EV is tested at different temperatures.

Li-ion Batteries classical modeling

Classic Li-ion Electrical Model and EMR Organization

Different models can be found for batteries. The complexity of the models depends on the application. Some models are called electrochemical models. They are based on the kinetic of the energy storage reactions inside the battery. This type of model is the nearest from the chemistry of the system [START_REF] Ahmed | Reduced-Order Electrochemical Model Parameters Identification and State of Charge Estimation for Healthy and Aged Li-Ion Batteries-Part II: Aged Battery Model and State of Charge Estimation[END_REF]. Therefore they can be used to monitor battery (ageing for example).

Equivalent circuit models are another kind of model. The battery electric behavior is represented by equivalent impedances and sources. Each electric parameter is linked with a phenomenon happening in the battery. For example, ionic transfer delays can be modeled by RC circuits [START_REF] Baghdadi | Electro-thermal model of lithium-ion batteries for electrified vehicles applications[END_REF] and ionic diffusion can be modelled by a constant phase element (CPE) [START_REF] Eddahech | Thermal characterization of a high-power lithium-ion battery: Potentiometric and calorimetric measurement of entropy changes[END_REF].

From the EV point of view, the battery is considered as a power source subsystem. Therefore, the equivalent circuit model is adapted for EV application. For a matter of simplicity we neglect ionic transfer and diffusion. The cell dispersion (i.e. the slight differences for electric parameters between cells grouped in a battery [START_REF] Lebel | A Lithium-Ion Battery Electro-Thermal Model of Parallellized Cells[END_REF]) is not considered as we just model a battery composed of one cell for this study. Therefore a simple Thevenin model is chosen. It consists in an open circuit voltage (OCV) dependent on the state-of-charge (SoC) and a series resistance Rbat (see Fig. 2).

Where SoCinit is the initial SOC (at t=0 s), ibat is the discharge current and t is the time.

The EMR is a causal method to organize various models [START_REF] Bouscayrol | Graphic Formalisms for the Control of Multi-Physical Energetic Systems: COG and EMR[END_REF]. Four kinds of pictograms (see appendix) are used to describe basic energetic functions (energy source, accumulation, distribution and conversion). Each pictogram highlights the causal inputs and outputs of a component. The product of the input and output is the power transiting through the function described by the related pictograms. The battery is considered as a voltage source and a conversion element. The voltage source (output OCV) takes as input the battery current ibat. The power flowing out the voltage source is obtained by the product of OCV and ibat. The resistance is described by an energy conversion element.

bat bat Ri - ) SoC ( OCV U = ( 2 )
To simplify, the traction part of the vehicle is considered as a current source in this EMR (equivalent traction source).

Generally, the parameters of these electrical models are not considered temperature dependent even though the energy storage principle of the battery is dependent on the temperature. Thus, the accuracy of the electric model is reduced if the impact of temperature is neglected.

Li-ion Battery Thermal Modeling

For a matter of simplicity, strong hypotheses have been taken for thermal modeling [START_REF] Forgez | Thermal modeling of a cylindrical LiFePO4/graphite lithium-ion battery[END_REF] (see Fig. 3):

• The battery cell is considered cylindrical with a length far bigger than the diameter. It is placed in open air (at Tamb temperature). • The battery is composed of two parts (core and surface) perfectly homogeneous with associated thermal capacitances Ccore and Csurf (in J/K). The Biot number [START_REF] Mahamud | Spatial-resolution, lumped-capacitance thermal model for cylindrical Li-ion batteries under high Biot number conditions[END_REF] quantifying the thermal gradient in the cell is neglected. Therefore, the temperature in the core is considered homogeneous. The thermal resistances (Rcond and Rconv (in K/W)) are located at the interfaces between two parts (core and surface for instance). • The heat source is located at the center of the core. The series resistance Rbat is the only source of heat. Reversible heat sources linked to the entropic effect [START_REF] Eddahech | Thermal characterization of a high-power lithium-ion battery: Potentiometric and calorimetric measurement of entropy changes[END_REF] will be neglected.

• The only heat transfer mechanisms considered are conduction for solids and convection for solid to gas. Radiations are neglected as the self-heating of batteries is moderate in EV. These hypotheses will be validated in part 4.1 with experimental thermal results coming from ref. [START_REF] Lin | Online Parameterization of Lumped Thermal Dynamics in Cylindrical Lithium Ion Batteries for Core Temperature Estimation and Health Monitoring[END_REF].The equivalent thermal circuit is given in Fig. 4. In equivalent thermal circuit power is equivalent to electric current and temperature is equivalent to voltage. The joule losses are a power generator (Pheat) and the ambient air is a source of temperature (Tamb). From a thermal point of view, the resistance Rbat is a component generating heat when a current goes through (Fig. 4). Thus, Rbat is described by a source of entropy flow qS1 in EMR (Fig. 5).

core 2 bat bat 1 s T i . R q = ( 3 )
In the thermal domain, the flow variable is the entropy flow qS (expressed in W/°K [START_REF] Horrein | Impact of heating system on the range of an electric vehicle[END_REF]) and the effort variable is the temperature (°K). Like in the electric domain, the product of the flow variable and the effort variable leads to the power. The battery thermal model can be considered as a source with the entropy flow as input and the inside battery temperature Tcore as output. The EMR of the thermal model is given in Fig. 5. The thermal capacitances are energy accumulation elements and the thermal resistances are energy conversion elements. The air is considered as a temperature source (Tamb).

Air

R conv R cond T core T amb T core q s1 q s3 q s4 q s7 T surf q s6 C surf T surf Heat Source C core R bat

Fig. 5. EMR of the thermal model

The following equations are describing the behavior of the thermal model. As EMR is a causal representation, all the equations are expressed in a causal way (i.e. without derivative).

)) T ln( dt ) q - q ( C 1 exp( = T init core t 0 3 s 1 s core core + ∫ ( 4 ) core cond surf core T . R T - T = 3 qs ( 5 ) surf cond surf core T . R T - T = 4 qs ( 6
)

)) T ln( dt ) q - q ( C 1 exp( = T init surf t 0 6 s 4 s surf surf + ∫ ( 7 ) surf conv amb surf T . R T - T = 6 qs ( 8
)
amb conv amb surf T . R T - T = 7 qs ( 9
)
Where Tcore init and Tsurf init are the initial temperatures of the core and the surface of the battery. For the next parts, we consider that the battery is at rest for a long time at t=0s. Consequently, the initial core and surface temperatures are equal to Tamb.

Li-ion Batteries Electro-thermal modeling organized using EMR

Influence of Temperature on Li-ion Battery Electric Parameters

The experimental results presented in Fig. 6 and Fig. 7 are directly taken from [START_REF] Jaguemont | Low Temperature Discharge Cycle Tests for a Lithium Ion Cell[END_REF] and [START_REF] Lin | Online Parameterization of Lumped Thermal Dynamics in Cylindrical Lithium Ion Batteries for Core Temperature Estimation and Health Monitoring[END_REF]. They show that the electric parameters of the batteries are influenced by the temperature in the operating temperature range of the batteries: We notice that :

• Cbat increases linearly with temperature (impact on the EV range) [START_REF] Jaguemont | Low Temperature Discharge Cycle Tests for a Lithium Ion Cell[END_REF].

• Rbat decreases exponentially with temperature (impact on the EV efficiency and range) [START_REF] Lin | Online Parameterization of Lumped Thermal Dynamics in Cylindrical Lithium Ion Batteries for Core Temperature Estimation and Health Monitoring[END_REF].

EMR for Electro-Thermal (ET) Coupling

The common element between electric and thermal model is the resistance of the battery Rbat (see Fig. 8). The electric and thermal models are coupled by Rbat which is now a multi-domain (electro/thermal ET) coupling element. The temperature of the battery core Tcore is used to update the value of Cbat and Rbat.

Air

C core R conv R cond T core T amb T core q S1 q s3 q s4 q s7 T surf q s6 C surf T surf

Tract

. Validation Using a Highly Dynamic Urban Assault Driving Cycle

The validation of the model is achieved with experimental results from the literature [START_REF] Lin | Online Parameterization of Lumped Thermal Dynamics in Cylindrical Lithium Ion Batteries for Core Temperature Estimation and Health Monitoring[END_REF]. The tested cells are Lithium Iron Phosphate cells (Fig. 9) which is a technology used in EVs. The electric and thermal parameters are given in Table I. The thermal dependence is plotted in Fig. 6 and Fig. 7. The convection resistance is identified by [START_REF] Lin | Online Parameterization of Lumped Thermal Dynamics in Cylindrical Lithium Ion Batteries for Core Temperature Estimation and Health Monitoring[END_REF] with a forced airflow (fan blowing) in a thermal chamber.

The value extracted will vary as a function of the airflow.

Parameters Values

Ubat nom 3.3 V Cnom @ 25°C 2.5 A.h Rnom @ 25°C

Table I: Parameters of the studied battery [4]

The model has been tested with the highly dynamic urban assault cycle (UAC) [START_REF] Lin | Online Parameterization of Lumped Thermal Dynamics in Cylindrical Lithium Ion Batteries for Core Temperature Estimation and Health Monitoring[END_REF]. This cycle is expressed in Crate in Fig 10 (1 C= 2.5 A). Then, the temperature given by the model has been compared with experimental results from [START_REF] Lin | Online Parameterization of Lumped Thermal Dynamics in Cylindrical Lithium Ion Batteries for Core Temperature Estimation and Health Monitoring[END_REF]. The maximal error is lower than 2°C (i.e. 4 %) validating the accuracy of the E/T model (see Fig. 11).

Test of the Cells with the NEDC

The ET model is suited for any cycle with lower dynamics than the UAC. For a first EV range estimation, the New European Driving Cycle (NEDC) is considered. This cycle is composed of an urban part and an extra urban part (Fig. 12). The studied vehicle is the commercial Tazzari Zezo (2seat vehicle of 600 kg) [START_REF] Zero | [END_REF]. The parameters of the pure electric model of the battery are classically identified at 25 °C [START_REF]A123 Systems[END_REF].

Fig. 13 presents the evolution of the estimated SoC during the NEDC for an ambient temperature of -20°C. Significant differences can be observed on the SOC estimation. For -20°C and one NEDC cycle the SOC variation of the ET model is + 22 % bigger than the SOC variation of the pure electric model. The process has been repeated for other temperatures (25°C and 40°C). Results are summarized in Table II. With a pure electric battery model the SoC variation is underestimated by 22 %. The driving range will thus be overestimated by a significant factor. Thus, the ET model with capacity temperature dependence is very useful for cold weather conditions. For 25°C the electric model gives satisfactory SoC estimation for NEDC cycle. Thus, at 25 °C temperature conditions the pure electric model is satisfactory for EV driving range estimation. For higher temperature (+40 °C) the driving range is under estimated by the electric model by a significant factor. Therefore, a part of the energy inside the battery is unused if the pure electric model is used for SoC estimation. Overall the ET model seems to be necessary for EV driving estimation as -20°C and + 40°C temperatures are encountered in many countries.

Conclusion

An electric model of a battery is coupled with a thermal model. This coupling enables an update of the electric parameters from the temperature evolution. This ET model has been validated by experimental results using a high-dynamics driving cycle. This model has then been used for a standard driving cycle for a small commercial EV. At low temperature, the driving range using a classical electric model is overestimated by about 20% and for high ambient temperature (40 °C) it is underestimated by about 10 %. The use of the ET model will thus improve the estimation of the driving range. This study is focused on one cell temperature. As perspective, the temperature interaction between cells in a battery pack can be considered [START_REF] Hijazi | Electro-thermal sizing of supercapacitor stack for an electrical bus: Bond graph approach[END_REF].

Fig. 1 .

 1 Fig. 1. Electric model for battery Fig. 2. EMR organization of electric model The SoC (in %) is frequently used as to quantify the energy stored in the batteries. It is calculated with the capacity of the battery Cbat (in Ah) which represents the current necessary to fully discharge the battery in one hour. It represents the remaining capacity stored in the battery over its maximum capacity Cbat.

Fig. 3 .Fig. 4 .

 34 Fig. 3. Heat transfer phenomena in the battery

5 Fig. 6 . 3 ]Fig. 7 .

 5637 Fig.[START_REF] Baghdadi | Lithium battery aging model based on Dakin's degradation approach[END_REF]. Battery capacity evolution with temperature[START_REF] Jaguemont | Low Temperature Discharge Cycle Tests for a Lithium Ion Cell[END_REF] Fig.7. Battery resistance evolution with temperature[START_REF] Lin | Online Parameterization of Lumped Thermal Dynamics in Cylindrical Lithium Ion Batteries for Core Temperature Estimation and Health Monitoring[END_REF]

Fig. 8 .

 8 Fig. 8. EMR of the E/T model of the battery 4. Driving range estimation 4.1. Validation Using a Highly Dynamic Urban Assault Driving Cycle

Fig. 9 .

 9 Fig.9. Studied Li-ion battery[START_REF]A123 Systems[END_REF]

Fig. 10 .

 10 Fig. 10. C rate of the urban assault cycle (UAC) Fig. 11. Temperature evolutions for the UAC

Fig. 12 .

 12 Fig. 12. Crate during the NEDC Fig. 13. SoC evolution of the NEDC

Table II : Difference between electric and ET model for different ambient temperature after one NEDC (SoC init = 100%)

 II

	Ambient temperature	Electric model final SoC	ET model final SoC	Error on Δ SoC
	-20 °C		88 %	+ 22 %
	25 °C	90.6 %	90.65 %	-0.5%
	40 °C		91.3 %	-8 %

Acknowledgements

This work was supported by the french general consulate in Quebec (réalisé avec le soutien du Service de à Québec).


Simon Mendez

Manouk Abkarian

In-plane elasticity controls the full dynamics of red blood cells in shear flow

The problem of understanding the movement of red blood cells (RBCs) is at the heart of hemorheology. It has thus motivated an extensive body of experimental and numerical works, which showed that RBCs display a rich dynamical behavior in pure shear flow. However, a clear physical understanding of the coupling between cells orientation, membrane deformations and circulation is still not emerging, notably due to the lack of a comprehensive and tractable model to serve as the theoretical foundations for data analysis. Here, we propose a new low-order model which, combined with detailed simulations and existing experimental data, demonstrates how membrane in-plane deformations and elasticity are the essential ingredients responsible for RBCs dynamics at low shear stresses. Our approach demonstrates that out-of-plane deformations and fluid inertia are not necessary to explain the RBC dynamics and underlines the importance of membrane stress-free shape. By reproducing all the details of known RBCs dynamics in shear flow, the new low-order provides a single framework to understand the full dynamics of RBCs at low shear stresses.

Red blood cells (RBCs) are the essential oxygen carriers of our body and compose almost half of blood's entire volume. Understanding their behavior in flow is therefore a key step to characterize both physiological and pathological blood circulatory conditions, and to propose non-invasive diagnostic tools by testing their impaired dynamics. The apparent structural simplicity of RBCs contrasts with the complexity of predicting their dynamics [START_REF] Abkarian | Fluid-structure interactions in low-Reynolds-number flows[END_REF]. Indeed, a RBC is a micron-scale deformable bag, of the shape of a biconcave disk at rest, containing a Newtonian solution of hemoglobin of viscosity η i , typically five times the viscosity of water and higher than the outer viscosity η o in physiological conditions. This solution is enclosed by a thin composite membrane made of an incompressible fluid lipid bilayer, scaffolded inside by a two-dimensional network of elastic filaments of spectrin [START_REF] Mohandas | Red cell membrane: past, present, and future[END_REF] with a typical in-plane shear modulus G s . Understanding the dynamics of a RBC comes down to two main issues: how its shape adapts to viscous stresses and how the resulting internal flow impacts cell behavior. This fluid-structure problem is highly non-linear and non-local in nature and despite decades of research, the full phase diagram of the different dynamical states of a RBC in pure shear flow is still not available, as demonstrated by the recent discovery of new dynamical states observed in physiological conditions of strong flows [START_REF] Lanotte | Red cells' dynamic morphologies govern blood shear thinning under microcirculatory flow conditions[END_REF].

However, experiments realized during the last decades on single RBCs for shear stresses η o γ lower than 0.1 Pa [START_REF] Abkarian | Swinging of red blood cells under shear flow[END_REF][START_REF] Dupire | Full dynamics of e red blood cell in shear flow[END_REF] identified numerous dynamics of orientation where membrane out-of-plane deformations are small and cells maintain their discocyte shape. For instance, RBCs present a dual solid-liquid behavior for increasing η o γ, depending on the viscosity ratio λ = η i /η o [START_REF] Goldsmith | Flow behaviour of erythrocytes. I. rotation and deformation in dilute suspensions[END_REF][START_REF] Fischer | The red cell as a fluid droplet: Tank treadlike motion of the human erythrocyte membrane in shear flow[END_REF]. For λ > ∼ 3, the cell flips like a solid coin. Its axis of symmetry precesses around the vorticity direction, describing a closed trajectory. However, unlike a solid ellipsoid [START_REF] Jeffery | The motion of ellipsoidal particles immersed in a viscous fluid[END_REF], this orbit progressively drifts for increasing shear rate γ until the symmetry axis aligns with the vorticity direction: the cell rolls [START_REF] Dupire | Full dynamics of e red blood cell in shear flow[END_REF][START_REF] Goldsmith | Flow behaviour of erythrocytes. I. rotation and deformation in dilute suspensions[END_REF][START_REF] Bitbol | Red blood cell orientation in orbit C=0[END_REF][START_REF] Yao | Low viscosity Ektacytometry and its validation tested by flow chamber[END_REF]. For λ < ∼ 3, the orbital drift to rolling is maintained [START_REF] Dupire | Full dynamics of e red blood cell in shear flow[END_REF][START_REF] Bitbol | Red blood cell orientation in orbit C=0[END_REF], but for higher γ, rolling is not stable and is replaced by a modified droplet-like behavior called swinging [START_REF] Abkarian | Swinging of red blood cells under shear flow[END_REF][START_REF] Dupire | Full dynamics of e red blood cell in shear flow[END_REF]: the membrane circulates around the center of mass like a tank-tread, while the orientation of one of the cell's largest axis of inertia oscillates around an average angle whose value decreases with γ. A rich corpus of numerical simulations reproduced part of these motions [START_REF] Lanotte | Red cells' dynamic morphologies govern blood shear thinning under microcirculatory flow conditions[END_REF][START_REF] Cordasco | Orbital drift of capsules and red blood cells in shear flow[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 low-shear-rate flows: effects of non-biconcave stress-free state in the cytoskeleton[END_REF][START_REF] Tsubota | Elastic behavior of a red blood cell with the membrane's nonuniform natural state: equilibrium shape, motion transition under shear flow, and elongation during tank-treading motion[END_REF][START_REF] Sinha | Dynamics of a single red blood cell in simple shear flow[END_REF][START_REF] Mauer | Complex shapes and dynamics of red blood cells in shear flow under physiological conditions[END_REF], but even a simple question as the essential ingredients yielding this complex dynamics remains unanswered. Dupire et al. [START_REF] Dupire | Full dynamics of e red blood cell in shear flow[END_REF] hypothesized for instance that in-plane elasticity may be responsible for orbital drift, but the only test provided was to rigidify the RBCs, indeed preventing drift but simultaneously precluding in-plane and out-of-plane deformations and inner flow. Whether out-of-plane deformations are necessary to explain RBC dynamics at low shear stress, or solely the viscosity ratio is unknown. Moreover, inertia, which is known to promote rolling motion for solid ellipsoids [START_REF] Rosén | Effect of fluid and particle inertia on the rotation of an oblate spheroidal particle suspended in linear shear flow[END_REF], may play a role. Inertia and out-of-plane deformations are small in this flow regime, but the orbital drift is a very slow process, occurring over hundreds of time units γ-1 . This suggests that even small effects may contribute to the long-term dynamics of RBCs.

Important breakthroughs in the physical understanding of RBC dynamics sprung from the development of theoretical models where the axis of symmetry of the cell is bound to the shear plane. Keller and Skalak [START_REF] Keller | Motion of a tank-treading ellipsoidal particle in a shear flow[END_REF] understood the fundamental role played by λ in the transition between flipping and tank-treading. Later, the existence of the membrane stress-free shape different from the sphere was added to this approach [START_REF] Abkarian | Swinging of red blood cells under shear flow[END_REF][START_REF] Skotheim | Red blood cells and other nonspherical capsules in shear flow: Oscillatory dynamics and the tank-treading-to-tumbling tran-sition[END_REF], inspired by the observations of in-plane elasticity and shape memory for the RBC membrane [START_REF] Fischer | The stress-free shape of the red blood cell membrane[END_REF][START_REF] Fischer | Shape memory of human red blood cells[END_REF]. These models predicted the transition between flipping and swinging with γ [START_REF] Abkarian | Swinging of red blood cells under shear flow[END_REF]. The membrane elements have to overcome an energy barrier to circulate around the RBC and make the cell tanktread/swing. In-plane elasticity was thus identified as a key element controlling dynamics for this planar configuration. However, such conclusions have not been generalized to 3-D dynamics: the absence of comprehensive approaches coupling arbitrary orientation, membrane circulation and stress, isolating the role played by the geometry and by each structural element, has been a true limitation for the understanding of RBC dynamics.

In this letter, we propose a three-dimensional low-order model for RBCs dynamics that we compare to detailed simulations and existing experimental data. We demonstrate the essential role played by both membrane inplane deformations and viscosity ratio. Our approach underlines the importance of membrane stress-free shape and reproduces all the dynamics observed both experimentally and numerically, simultaneously ruling out outof-plane deformations and fluid inertia as necessary elements to explain the behavior of RBCs under shear. The problem geometry is shown in Fig. 1. xi (unit vectors êi ) denote the coordinates of a fixed frame R, and x i (unit vectors e i ) the coordinates of a body frame R of same origin, whose axes correspond to the principal axes of the cell. To define the orientation of the cell with respect to the fixed frame, the Euler angles θ, ϕ and ψ are used (Fig. 1(a) and supplementary material [44]). An external shear flow is imposed: û0 1 = γ x2 . The relative strength of the shear flow is defined by the capillary number Ca = η o γa/G s , with a = 2.82 µm a characteristic size of the RBC [START_REF] Cordasco | Comparison of erythrocyte dynamics in shear flow under different stress-free configurations[END_REF].

The small axis (symmetry axis) of the RBC is along x 3 . The Euler angle ϕ is the angle between x 3 and the vorticity direction x3 , and thus characterizes the orbits of the cell. Membrane circulation is tracked in a Lagrangian way following the location of the small axis of the stressfree shape of the cytoskeleton in the body frame, denoted by vector X (Fig. 1b). When X is not along x 3 , the membrane has circulated and stored in-plane elastic energy.

The present analysis is based on models of decreasing complexity. We first present simulations of the full fluidstructure interaction problem, including all modes of deformations and small yet non-zero fluid inertia. Simulations provide reference data for the subsequent low-order modeling. Simulations are performed using YALES2BIO [START_REF] Lanotte | Red cells' dynamic morphologies govern blood shear thinning under microcirculatory flow conditions[END_REF][START_REF] Mauer | Complex shapes and dynamics of red blood cells in shear flow under physiological conditions[END_REF][START_REF] Gibaud | Numerical simulation of red blood cells flowing in a blood analyzer[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][START_REF] Sigüenza | How should the optical tweezers experiment be used to characterize the red blood cell membrane mechanics?[END_REF][START_REF] Sigüenza | Fluidstructure interaction of a pulsatile flow with an aortic valve model: A combined experimental and numerical study[END_REF] (see supplementary material [44]), a solver for the incompressible Navier-Stokes equation based on the immersed boundary method. The infinitely thin membrane resists shear and area changes following the Skalak model [START_REF] Skalak | Strain energy function of red blood cell membranes[END_REF]. Bending resistance is modeled by the Helfrich curvature energy [START_REF] Helfrich | Elastic properties of lipid bilayers: Theory and possible experiments[END_REF][START_REF] Zhong-Can | Bending energy of vesicle membranes: General expressions for the first, second, and third variation of the shape energy and applications to spheres and cylinders[END_REF]. The membrane is prestressed at rest, its stress-free shape being a quasispherical oblate ellipsoid [START_REF] Cordasco | Comparison of erythrocyte dynamics in shear flow under different stress-free configurations[END_REF][START_REF] Peng | Erythrocyte responses in low-shear-rate flows: effects of non-biconcave stress-free state in the cytoskeleton[END_REF][START_REF] Lim | Stomatocyte-discocyte-echinocyte 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[END_REF], with the same surface area as the RBC. Different stress-free ellipsoids are used, defined by their reduced volume v ref , the volume of the ellipsoid divided by the volume of the sphere having the same surface area. This reference shape is first deflated until an equilibrium biconcave shape is reached, which is then subjected to shear. The higher v ref , the lower the elastic energy barrier for membrane circulation [START_REF] Peng | Erythrocyte responses in low-shear-rate flows: effects of non-biconcave stress-free state in the cytoskeleton[END_REF][START_REF] Fischer | The stress-free shape of the red blood cell membrane[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]. In all simulations, G s = 2.5 µN.m -1 [START_REF] Lim | Red blood cell shapes and shape transformations: Newtonian mechanics of a composite membrane[END_REF][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], and the bending modulus is κ b = 3.0 × 10 -19 J [START_REF] Mauer | Complex shapes and dynamics of red blood cells in shear flow under physiological conditions[END_REF].

As a first result, we show an example of orbit selection in Fig. 1(c: after a transient flipping, the cell initially in the shear plane (ϕ = π/2) reorients and flips over a specific limit orbit, characterized by ϕ orbit , the minimum value of ϕ after reorientation. At higher external stress, the cell may reorient to rolling: the RBC small axis aligns with the vorticity direction (ϕ → 0, Fig. 1d). During the reorientation, the reference shape circulates along the RBC without major cell deformation (Fig. 1e).

The orbital drift of RBCs as a function of Ca is first examined for λ = 1.0 (Fig. 2): RBCs flip over orbits of decreasing value of ϕ orbit when Ca increases, until rolling is reached (ϕ orbit = 0), as seen in experiments [START_REF] Dupire | Full dynamics of e red blood cell in shear flow[END_REF]. For higher Ca, RBCs return to the shear plane (ϕ orbit = π/2) and tank-tread/swing [START_REF] Abkarian | Swinging of red blood cells under shear flow[END_REF][START_REF] Dupire | Full dynamics of e red blood cell in shear flow[END_REF]. This reorientation from rolling to tank-treading is performed without flipping: the RBC takes a fixed inclination and spins without angle oscillations, a motion called 'frisbeeing' or 'hovering', reported in simulations [START_REF] Cordasco | Comparison of erythrocyte dynamics in shear flow under different stress-free configurations[END_REF][START_REF] Sinha | Dynamics of a single red blood cell in simple shear flow[END_REF] and transiently in experiments with a step increase of shear rate [START_REF] Dupire | Full dynamics of e red blood cell in shear flow[END_REF].

In Fig. 2, results at λ = 1 show that the stress-free shape of the cell determines the stress at which orbital drift occurs. The energy barrier for membrane circulation being higher at v ref = 0.96 than at 0.997, the reorientation occurs for higher Ca. As a consequence, the frisbee motion disappears. Instead, the RBC buckles [START_REF] Lanotte | Red cells' dynamic morphologies govern blood shear thinning under microcirculatory flow conditions[END_REF][START_REF] Mauer | Complex shapes and dynamics of red blood cells in shear flow under physiological conditions[END_REF] and becomes a stomatocyte [START_REF] Mauer | Complex shapes and dynamics of red blood cells in shear flow under physiological conditions[END_REF]: at such higher stresses, the whole RBC deforms. Note that results align well with the experiments [START_REF] Dupire | Full dynamics of e red blood cell in shear flow[END_REF] for v ref = 0.96 (Fig. 2, inset). The simulations highlight the importance of the membrane stress-free shape 13, 32], and thus its in-plane elasticity, as a key parameter controlling the 3-D dynamics of RBCs. However, the roles of inertia or shape deformations, inherently present both in simulations and in experiments, remain unclear. In order to separate the different mechanisms possibly leading to orbit selection and rolling, we introduce a theoretical model in which inertia and shape deformations are absent. This is inspired by the works of Jeffery [START_REF] Jeffery | The motion of ellipsoidal particles immersed in a viscous fluid[END_REF], Keller and Skalak [START_REF] Keller | Motion of a tank-treading ellipsoidal particle in a shear flow[END_REF] and their extensions for elastic membranes [START_REF] Abkarian | Swinging of red blood cells under shear flow[END_REF][START_REF] Skotheim | Red blood cells and other nonspherical capsules in shear flow: Oscillatory dynamics and the tank-treading-to-tumbling tran-sition[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].

The model considers an axisymmetric ellipsoid suspended in an unbounded shear flow û0 1 = γ x2 , at zero Reynolds number. The internal fluid can flow inside the ellipsoid of prescribed fixed shape. The membrane can circulate along this shape. In the existing models for fluid ellipsoids [START_REF] Abkarian | Swinging of red blood cells under shear flow[END_REF][START_REF] Keller | Motion of a tank-treading ellipsoidal particle in a shear flow[END_REF][START_REF] Skotheim | Red blood cells and other nonspherical capsules in shear flow: Oscillatory dynamics and the tank-treading-to-tumbling tran-sition[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], one principal axis was lying in the shear plane. Here, the orientation is left free to adapt to the flow and the membrane circulation is solved in three dimensions. The tank-treading velocity is defined by an unknown vector Ω around which the membrane circulates (Fig. 1b): the membrane surface velocity v m i in the body frame reads v m i = a i ijk Ωj x k /a k [START_REF] Keller | Motion of a tank-treading ellipsoidal particle in a shear flow[END_REF], where ijk is the Levi-Civita symbol in three dimensions and Ωj is the tank-treading rate around axis x j .

Since the total moment acting on a freely suspended neutrally buoyant particle vanishes, it may be shown that the spins of the ellipsoid around its axes ( ω1 = θ sin ϕ sin ψ + φ cos ψ, ω2 = θ sin ϕ cos ψ -φ sin ψ, ω3 = θ cos ϕ + ψ) read:

                     ω1 = - γ 2 sin ϕ sin ψ + a 2 2 -a 2 3 a 2 2 + a 2 3 e 0 23 - 2a 2 a 3 a 2 2 + a 2 3 Ω1 , ω2 = - γ 2 sin ϕ cos ψ + a 2 3 -a 2 1 a 2 3 + a 2 1 e 0 31 - 2a 3 a 1 a 2 3 + a 2 1 Ω2 , ω3 + Ω3 = - γ 2 cos ϕ, (1)
with a 1 , a 2 and a 3 the semi-axes of the ellipsoid (a 1 = a 2 for the axisymmetric ellipsoid considered). Ω3 cannot be distinguished from ω3 because x 3 is an axis of symmetry. The e 0 ij are the components of the strain rate tensor of the undisturbed flow in the body frame: (2) The tank-treading rates are obtained by writing the balance between the stress jump across the membrane and the membrane elastic load (calculated assuming that the membrane is a prestressed Kelvin-Voigt viscoelastic material [START_REF] Abkarian | Swinging of red blood cells under shear flow[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]):

           Ω1 = -2f 3 f 2 -λ ef f f 1 e 0 23 + γβC Ca ξ 2 ξ 3 , Ω2 = 2f 3 f 2 -λ ef f f 1 e 0 31 + γβC Ca ξ 1 ξ 3 . (3)
The f i and β are geometric factors [START_REF] Abkarian | Fluid-structure interactions in low-Reynolds-number flows[END_REF][START_REF] Keller | Motion of a tank-treading ellipsoidal particle in a shear flow[END_REF]. λ ef f is an effective viscosity ratio accounting for membrane viscosity (here we use λ ef f = λ). C is a non-dimensional constant accounting for membrane prestress [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] (the larger the reduced volume of the reference shape, the lower C). ξ i are the coordinates of the small axis of the reference shape mapped on a sphere: ξ i = X i /a i (no summation). X i and thus ξ i , are advanced from the tank-treading velocity v m i . The model is consistent with existing models for solid and fluid ellipsoids [START_REF] Abkarian | Swinging of red blood cells under shear flow[END_REF][START_REF] Jeffery | The motion of ellipsoidal particles immersed in a viscous fluid[END_REF][START_REF] Keller | Motion of a tank-treading ellipsoidal particle in a shear flow[END_REF]. The system of non-linear ordinary differential equations is solved numerically. The input data are the geometry of the ellipsoid (here, a 1 = a 2 = 4.2375 µm and a 3 = 1.2511 µm to match the volume and surface of the simulated RBC), γ, Ca, C, λ and the initial orientation. Results of the model are displayed in the fixed frame R by projecting on a unit sphere the evolution of the small axes of both the ellipsoid with the vector e 3 R and of the stress-free shape with the vector ξ R .

Figure 3 shows that the model reproduces all the typical motions of RBCs reported by experiments and simulations at low shear stress [START_REF] Abkarian | Fluid-structure interactions in low-Reynolds-number flows[END_REF]. For λ > ∼ 3, the model yields a flipping dynamics over orbits (Fig. 3a) determined by Ca, the orbit shrinking around the vorticity axis when Ca increases, until rolling is reached (Fig. 3b). For higher values of Ca and λ > ∼ 3, rolling is the only stable motion predicted. For λ < ∼ 3 and increasing Ca, the model first predicts the flipping and the orbital drift to rolling as for λ > ∼ 3, then the frisbee [START_REF] Cordasco | Comparison of erythrocyte dynamics in shear flow under different stress-free configurations[END_REF][START_REF] Sinha | Dynamics of a single red blood cell in simple shear flow[END_REF] dynamics (Fig. 3c). For higher Ca, the RBC goes to the shear plane to perform swinging [START_REF] Abkarian | Swinging of red blood cells under shear flow[END_REF], possibly with a transient 'kayaking' oscillation around the shear plane (Fig. 3d) [START_REF] Cordasco | Orbital drift of capsules and red blood cells in shear flow[END_REF][START_REF] Dupont | Stable equilibrium configurations of an oblate capsule in simple shear flow[END_REF]. In addition to qualitative agreement with the documented behaviors of RBCs, the model aligns well with full simulations (Fig. 2). For v ref = 0.997, the prestress coefficient C is adjusted in the model (Eq. 3) using data at λ = 5 and kept constant to predict results at λ = 1 and λ = 0.2. The model is able to predict both the orbital drift and the reorientation from rolling to tank-treading in the frisbee state. The comparison with simulations is excellent, despite some discrepancies, in particular when Ca increases. This comparison supports that the in-plane elasticity, the only mechanical property included in the model, is the main effect driving the reorientation. Out-of-plane deformation and inertia only have a second-order effect on the RBC dynamics.

Figures 3(a-d) also show the interplay between the dynamics of the ellipsoid and of its membrane: the dynamics depends on the ability of the membrane to circulate over the ellipsoid, which results in the separation between the small axes of the ellipsoid and of the reference shape. To further understand the results, we remark that Eq. ( 1) expresses an equilibrium between moments of forces [START_REF] Keller | Motion of a tank-treading ellipsoidal particle in a shear flow[END_REF]. The spins of the ellipsoid result from the sum of the moments of forces exerted by the external flow and those associated to membrane circulation (terms in Ωi ). The latter can be divided into two contributions (Eq. 3), a purely viscous one and an elastic one, in 1/Ca. Equation (1) could thus be formally expressed as ω = ωExt.F low + ωV isc.Circ. + ωElast. , where the terms of the right-hand side are associated with the external flow, the viscous circulation and the elastic contribution to the circulation, respectively.

Figure 4(a) displays a comparison of the trajectory of the small axis of the ellipsoid in the solid case (Ca = 0, Jeffery orbits) and in the purely viscous case (Ca = +∞). Vectors ωV isc.Circ. are displayed to analyze how circulation of the membrane exerts an additional moment which modifies Jeffery's trajectory. It is worth noting that torque vectors pointing inwards slow down the trajectory without changing the orbit. On the other hand, torques in the same (resp. opposite) direction as the trajectory deviate the trajectory towards the inside (resp. outside). Figure 4(a) thus shows how the additional moment associated with the circulation of the membrane makes the orbit thinner for Ca = +∞. Orbits get actually thinner as λ decreases (not shown). Note also that ωV isc.Circ. has the same symmetries as ωExt.F low .

Figure 4(b) displays the trajectory of e 3 R and ξ R in a case with membrane elasticity. Contrary to ωV isc.Circ. , the elastic torque ωElast. is not symmetric with respect to x1 and x2 : the orbit thus becomes tilted. It is easily shown that the elastic torque is actually proportional to (e 3 •ξ)e 3 ×ξ: its direction is thus normal to e 3 and ξ, while its norm is determined by the angle between e 3 and ξ, as evidenced by Fig. 4(b). The elastic torque is of course generated by the circulation of the membrane, which is responsible for the storage of in-plane elastic energy. However, a second mechanism determines the direction of the elastic torque: the ellipsoid spins around its symmetry axis (see Eq. 1). This spinning is responsible for a rotation of ξ around e 3 , and consequently of a reorientation of the restoring elastic torque. This is seen in Fig. 4 by the difference in the orientation of ωV isc.Circ. and ωElast. . The spinning of the ellipsoid thus makes the membrane restore its elastic energy with a torque oriented differently from the viscous torque which made the membrane circulate. This process is summarized in a schematic way in Fig. 4(c). This element is essential for the understanding of how elasticity controls the dynamics of red blood cells, in particular by modifying their orbits, which are not Jeffery's orbits.

Finally, a full study of the dynamics of RBCs at low shear stresses is made possible by our model: it gathers the different dynamics of RBCs and capsules that were studied before in experiments [START_REF] Dupire | Full dynamics of e red blood cell in shear flow[END_REF][START_REF] Bitbol | Red blood cell orientation in orbit C=0[END_REF][START_REF] Yao | Low viscosity Ektacytometry and its validation tested by flow chamber[END_REF] and simulations [START_REF] Lanotte | Red cells' dynamic morphologies govern blood shear thinning under microcirculatory flow conditions[END_REF][START_REF] Cordasco | Orbital drift of capsules and red blood cells in shear flow[END_REF][START_REF] Cordasco | Comparison of erythrocyte dynamics in shear flow under different stress-free configurations[END_REF][START_REF] Sinha | Dynamics of a single red blood cell in simple shear flow[END_REF][START_REF] Dupont | Stable equilibrium configurations of an oblate capsule in simple shear flow[END_REF][START_REF] Omori | Reorientation of a nonspherical capsule in creeping shear flow[END_REF], but only partially explained until now. By proposing a model able to reproduce the flipping dynamics over orbits, the orbital drift to rolling, frisbeeing, kayaking, swinging in shear flow depending on λ and Ca, we demonstrate that the whole dynamics at low shear rates is mainly controlled by the storage and release of in-plane elastic energy in the membrane.

The behaviors of fluid and rigid ellipsoids are different: Jeffery [START_REF] Jeffery | The motion of ellipsoidal particles immersed in a viscous fluid[END_REF] demonstrated that a rigid ellipsoid rotates indefinitely in shear flow, its symmetry axis precessing around the vorticity direction, describing closed orbits. These orbits are degenerate since an infinite number of them are allowed, their selection depending only on the initial orientation of the particle. This degeneracy, proper to Stokes flows, is removed in real flows. Ellipsoids have been shown to display a rich variety of behaviors depending on both particle and fluid Reynolds numbers [START_REF] Rosén | Effect of fluid and particle inertia on the rotation of an oblate spheroidal particle suspended in linear shear flow[END_REF][START_REF] Rosén | The dynamical states of a prolate spheroidal particle suspended in shear flow as a consequence of particle and fluid inertia[END_REF]. For instance, neutrally buoyant oblate ellipsoids have been shown to generally roll on their edge for moderate Reynolds numbers [START_REF] Rosén | Effect of fluid and particle inertia on the rotation of an oblate spheroidal particle suspended in linear shear flow[END_REF], but present a steady inclination for higher inertia. Another cause of orbit selection is the visco-elasticity of the outer fluid [START_REF] Gauthier | Particle motions in non-newtonian media. I: Couette flow[END_REF][START_REF] Maffettone | Particle dynamics in viscoelastic liquids[END_REF], though its origin remains unclear. Here, we have identified a new cause for removing the degeneracy of Jeffery's orbits: the energy barrier associated with the circulation of the membrane over the RBC. When the RBC is not in the shear plane, the sliding of its membrane around its ellipsoidal shape generates a lateral torque that makes the particle deviate from Jeffery's orbit.

The role of in-plane elasticity is thus identified and can be now separated from the other effects thanks to our model. We show how the circulation of the stress-free shape coupled to spinning, together with the viscosity contrast, are the two main ingredients to capture the full RBCs phase diagram of movement at low shear stresses. Our model is easy to implement and runs orders of magnitude faster than full simulations. Most importantly, it will allow the stability analysis of the different regimes of movement and a better characterization of their transitions. In addition, the present model paves the way not only to the development of theoretical predictions for blood rheology at low volume fractions, but enables predictions of the dynamics of RBCs in general linear microflows as will be detailed in a full-length publication.

This work was supported by the Labex Numev Convention grants ANR-10-LABX-20. YALES2BIO (http: //www.math.univ-montp2.fr/ ~yales2bio) simulations have used HPC resources from GENCI-CINES (grants 2016-c2016037194 and 2017-A0020307194). Profs. M. Hillairet, T. Mignon and F. Nicoud (IMAG Montpellier) are gratefully acknowledged for fruitful discussions.

FIG. 1 .

 1 FIG. 1. Geometry and variables of interest (a,b) and examples of dynamics in shear flow (c-e). (a) The orientation of the body frame (xi) with respect to the fixed frame (xi) is defined by the Euler angles [44]. (b) The membrane tank-treading is defined in the body frame by the tank-treading rate vector Ω and the small axis of the reference shape X. The RBC is colored by the reference shape radius (radial distance to the cell center), from blue (low values, small axis) to red (large values, long axes). (c) Simulation of an orbit selection (λ = 1, Ca = 0.016, v ref = 0.997): ϕ along time. The minimum angle in the stationary state, ϕ orbit , characterizes the limit orbit. (d-e) Simulation from flipping to rolling (λ = 1, Ca = 0.025, v ref = 0.997): ϕ along time (d) and corresponding sequence (e) of the simulated RBC colored as in (b). The + in (d) mark the instants when the shapes are displayed in (e).

FIG. 2 .

 2 FIG. 2. Orbital drift of a RBC with Ca in shear flow varying λ and membrane energy barrier. Comparison of ϕ orbit from simulations (filled symbols) and low-order model (solid lines). The energy barrier is modified in the simulations by changing v ref (0.997 and 0.96) and adjusted in the model (Eq. 3) by changing C (0.00625 and 0.0175, respectively) based on the data at λ = 5 and λ = 1, respectively. Frisbee angles are also displayed for v ref = 0.997 and λ = 1 (open symbols for simulations and dashed line for the model). Inset: angle α between the flow axis and the small axis when the RBC is on the edge while flipping: experiments for different RBCs in PBS-dextran solutions with dextran molecular weight of 10 5 g/mol (solid lines) [5] and simulations at λ = 1 (symbols).

22

  cos ψ sin 2θ sin 2ϕ + cos 2θ sin ψ sin ϕ , sin ψ sin 2θ sin 2ϕ -cos 2θ cos ψ sin ϕ .

FIG. 3 .

 3 FIG. 3. Dynamics of a fluid ellipsoid (a1 = a2 = 4.2375 µm, a3 = 1.2511 µm, C = 0.01), in shear flow shown by the trajectories projected in the shear plane of e3 R (solid red line) and ξ R (blue dotted line). The circle of radius 1 is also shown; black dots mark the initial orientation. Dynamics (see supplementary movies [44]) for different values of λ and Ca: flipping with orbit selection (a), rolling (b) frisbee (c), transient kayaking and swinging (d).

FIG. 4 .

 4 FIG. 4. Example of closed orbits (followed in the clockwise direction) for a fluid ellipsoid in shear flow (a1 = a2 = 4.2375 µm, a3 = 1.2511 µm, C = 0.01, λ = 5), projected in the shear plane. The circle of radius 1 is shown in dotted lines. (a) Cases Ca = 0 (i.e. solid ellipsoid, black dashed line) and Ca = +∞ (fluid ellipsoid without membrane elasticity, solid magenta line). Vectors ωV isc.Circ. are displayed at 16 evenly spaced instants over the period. (b) Case Ca = 0.013: trajectory of e3 R (solid red line) and ξ R (blue dashed line). Vectors ωElast. are displayed. (c) Schematic of the drift due to the elastic torque.