Swimming of Spermatozoa in a Maxwell Fluid

It has been suggested that the swimming mechanism used by spermatozoa could be adopted for self-propelled micro-robots in small environments and potentially applied to biomedical engineering. Mammalian sperm cells must swim through a viscoelastic mucus layer to find the egg cell. Thus, understanding how sperm cells swim through viscoelastic liquids is significant not only for physiology, but also for the design of micro-robots. In this paper, we developed a numerical model of a sperm cell in a linear Maxwell fluid based on the boundary element slender-body theory coupling method. The viscoelastic properties were characterized by the Deborah number (De), and we found that, under the prescribed waveform, the swimming speed decayed with the Deborah number in the small-De regime (De < 1.0). The swimming efficiency was independent of the Deborah number, and the decrease in the swimming speed was not significantly affected by the wave pattern.


Introduction
The swimming of spermatozoa is responsible for the reproduction process, especially for finding egg cells. The human spermatozoon is about 50-100 µm [1], and has one long flagellum. By beating the flagellum, the sperm can propel into surrounding liquid, with an average swimming speed on the order of O(10 −5 ) to O(10 −4 ) m/s [1,2]. Recently, sperm-inspired micro-robots have been developed by researchers in the field of biomedical engineering and biomimetics [3,4]. Thus, understanding how sperm swim is important not only for physiology, but also for the fabrication of micro-robots.
Over small scales, the viscous effect can dominate when the inertial effects become sufficiently small. Thus, the locomotion of sperm cells is governed by micro-fluid dynamics. Several researchers have studied how sperm swim from the perspective of fluid mechanics. Smith et al. [5] numerically investigated sperm cells swimming near a planar wall. They showed that the cells accumulated at the wall and, due to the conical envelope of the flagellar trajectory, became aligned with it. The equilibrium height and stability were analyzed by Ishimoto and Gaffney [6]. Gadêlha et al. [7] investigated the fluid-structure interaction of flagellar motions and reported that nonlinear flagellar dynamics are essential for breaking the symmetry of the waveform. Kantsler et al. [8] experimentally observed the rheotaxis of sperm cells, which is self-organized upwards, swimming in the fluid flow. This mechanism can be explained in terms of fluid mechanics [9,10]. More detail on the fluid mechanics of sperm cells is provided in [1,11].
Mammalian sperm cells must swim through mucus layers in the genitals, such as the oviduct and the uterus. Mucus has viscoelastic properties [1]. Ishimoto et al. [12] investigated the locomotion of sperm cells through a viscoelastic fluid. They observed small yawing when the sperm were swimming through a viscoelastic medium, related to the greater density of the regularized point-force singularity along the flagellum. Wróbel et al. [13] developed a model of a flagellar swimmer moving through a viscoelastic network, on a length scale comparable to that of the swimmer. They concluded that stiffer networks enhance the velocity of the swimmer. Although a number of researchers have investigated sperm cells swimming through a viscoelastic liquid, many questions remain unanswered. In particular, the efficiency and locomotive power of sperm in a viscoelastic fluid have not been studied. As viscoelastic fluid can easily be found in the human body, knowledge about the efficiency of locomotion through such media should be improved to further our understanding of reproduction and our ability to fabricate micro-robots capable of operating in biological environments.
In this study, we carried out a numerical investigation of sperm cells swimming through viscoelastic fluid. The viscoelastic fluid was modeled by a linear Maxwell fluid. We evaluated the swimming speed and efficiency whilst varying the Deborah number. In Section 2, we explain the equations governing the motion and our numerical algorithm. Our results are presented in Section 3. Finally, we conclude the paper in Section 4.

Governing Equations and Numerical Methods
In this section, we explain the equations governing the motion of the sperm as well as our numerical methods. Details of the fluid mechanics of a Maxwell fluid in the low-Reynolds number regime can be found in [14,15], and we explain it briefly here.

Green's Function for a Linear Maxwell Fluid
We assumed that a sperm cell was immersed in an unbounded linear Maxwell fluid. We estimated the Reynolds number of the sperm swimming through the fluid to be much smaller than unity [1]. Hence, we assumed that the fluid motion was a Stokes flow. The deviatoric stress tensor of the Maxwell fluid is defined as: where τ is the deviatoric stress of the fluid, λ is the relaxation time, and µ is the viscosity. D is the rate of deformation tensor, which is defined as: where v is the velocity of the fluid. By introducing the invertible linear operator L = µ(1 + λ∂/∂t) −1 , the momentum equation of the Maxwell fluid with a point force F can be written as [15]: where p is the pressure and δ(x,t) is the delta distribution function in space and time, so that the point force F is concentrated at a certain point x s and time t s . By using Green's functions, the singularity solution of Equation (3) is given by [15]: where P(x,t; x s ,t s ) and J(x,t; x s ,t s ) are the pressure kernel and the Stokeslet for the Maxwell fluid, respectively, which are singular at x s and t s : and where P 0 (x; x s ) and J 0 (x; x s ) are the Green's functions of a Newtonian fluid: and where r = |r| and r = xx s .

Boundary Integral Equation
To describe the fluid dynamics of a sperm cell in viscoelastic fluid, we now introduce the boundary integral equation. We were required to define an initial condition due to the presence of the time derivative in Equation (4). We assumed, as a starting condition, that the dynamics were trivial for t < 0 and that we had Stokes dynamics for t = 0 [15]. For t = t s (> 0), the flow field was calculated by integrating the singularity solution with respect to time and space: where q is the traction force per unit area, and S represents the surface of the cell. In this study, the slender-body theory [16] was applied to the motion of the flagellum. The flagellum was modeled as a thin, curved rod with a length-radius ratio a/L = 3.57 × 10 −3 [1], where a is the flagellar radius and L is the length. Then, the flow field can be described by the following boundary integral-the slender-body coupling formulation: where f is the force per unit length of the flagellum, and K 0 is the slender-body kernel [16]. The delta distribution gives the limiting value of the time integral, and the boundary integral of the Maxwell fluid is reduced to the Newtonian kernel with the leading time derivative:

Flagellar Waveform and Head Geometry
To accurately describe the flagellar waveform, we introduce an orthonormal body frame, ξ i , with the origin at the head-tail junction. We assumed a rigid connection between the head and the flagellum, and the minus ξ 1 -axis was set to coincide with the orientation vector. The waveform was specified a priori for a given transformation in terms of the following sinusoidal function [17]: where T is the period of the flagellar waveform, η is the wavelength (which can vary with respect to ξ 1 ), and α is the chirality parameter. These parameters were selected such that they coincided with experimental observations. An experimental study of the waveform of the human sperm flagellum was reported in [18,19]. In low-viscosity liquid, high apparent curvatures appeared along with whole flagellum, and the head made repeated yawing motions. In high-viscosity liquid, on the other hand, a steeping wave appeared near the tip of the flagellum and the yawing motion of the head was less pronounced. We defined two different beating patterns to express these different wave modes, as shown in Figure 1, by varying the wavelength η [17]: η = 1 for Mode 1, and η = 0.9613 − 0.038 tan −1 (5ξ 1 − 4) for Mode 2. The amplitudes of the two waves were adjusted to match the swimming speed in a Newtonian fluid. If the chirality parameter, α, is also set to α = 0.25 [18], then the sperm makes rolling motions as it swims.
Micromachines 2019, 10, x FOR PEER REVIEW 4 of 9 The amplitudes of the two waves were adjusted to match the swimming speed in a Newtonian fluid. If the chirality parameter, α, is also set to α = 0.25 [18], then the sperm makes rolling motions as it swims. Both human and bull sperm are similar to asymmetric ellipsoids [1]. We mimicked the elliptical sperm head using the following mapping function [10]: where X sp is a point on the head of the sperm, X is a point on the sphere with radius , and c1, c2, and c3 are non-dimensional shape parameters. The parameters were set to L = 4.17 × 10 −2 , c1 = 3.0, c2 = 2.0, and c3 = 4.0 to ensure that the morphology accurately reflected a human sperm cell.

Numerical Procedure
To express the locomotion of the sperm cell, we assumed that the cell moves rigidly. One motion of the flagellum is defined by Equation (12). The velocity at a point xs on the cell can be decomposed as: , (14) where V is the translational velocity, Ω is the angular velocity, r = x s − x g , and xg is the center-of-mass of the head. The velocity of the flagellum wave is v fla (xs,t), a function of the point of interest and time. When xs is located on the cell head, then v fla is equal to zero. When xs is located on the flagellum, then v fla is determined by the time derivative of Equation (12). We solved the resistance problem defined in Equations (11) and (14) with respect to the unknowns V, Ω and the tractions q and f. The head and flagellum were discretized into 1280 triangular meshes and 200 curved segments by spline interpolation. The surface integral of Equation (11) was solved using Gaussian numerical integration, while the finite time difference scheme was applied to the time derivative. We assumed the following force-free and torque-free conditions so that we could express the free-swimming sperm cell as: Both human and bull sperm are similar to asymmetric ellipsoids [1]. We mimicked the elliptical sperm head using the following mapping function [10]: where X sp is a point on the head of the sperm, X is a point on the sphere with radius , and c 1 , c 2 , and c 3 are non-dimensional shape parameters. The parameters were set to /L = 4.17 × 10 −2 , c 1 = 3.0, c 2 = 2.0, and c 3 = 4.0 to ensure that the morphology accurately reflected a human sperm cell.

Numerical Procedure
To express the locomotion of the sperm cell, we assumed that the cell moves rigidly. One motion of the flagellum is defined by Equation (12). The velocity at a point x s on the cell can be decomposed as: v( where V is the translational velocity, Ω is the angular velocity,r = x s − x g , and x g is the center-of-mass of the head. The velocity of the flagellum wave is v fla (x s ,t), a function of the point of interest and time. When x s is located on the cell head, then v fla is equal to zero. When x s is located on the flagellum, then v fla is determined by the time derivative of Equation (12). We solved the resistance problem defined in Equations (11) and (14) with respect to the unknowns V, Ω and the tractions q and f . The head and flagellum were discretized into 1280 triangular meshes and 200 curved segments by spline interpolation. The surface integral of Equation (11) was solved using Gaussian numerical integration, while the finite time difference scheme was applied to the time derivative. We assumed the following force-free and torque-free conditions so that we could express the free-swimming sperm cell as: The above equations were also discretized so that Gaussian numerical integration could be performed. Once V and Ω were determined, the point of interest was updated using the second-order Runge-Kutta method. For more detail about the numerical method, please refer to our previous study [10]. We then introduced an important parameter, the Deborah number, which represents the elasticity and viscosity ratio of the viscoelastic fluid, and is characterized by the relaxation time λ. In this study, the Deborah number is defined as De = ωλ, where ω is the angular frequency of the flagellar beat. In the small-Deborah number regime, De << 1, the mechanical properties of a Maxwell fluid asymptotically approach those of a Newtonian fluid. Due to the numerical instability, the maximum Deborah number is limited to be smaller than unity in this study. Although maximum De is limited, a weak viscoelastic liquid can be easily found in vitro situations [1,2,12]. For example, in an in vitro experiment of sperm swimming in 1% methylcellulose [12], the Deborah number was estimated to be De ≈ 0.2. We thus decided to vary De in the range 0 ≤ De ≤ 1.

Swimming Motion of a Sperm Cell in a Newtonian Fluid
Firstly, we investigated a sperm cell swimming through a Newtonian fluid. In the Stokes flow regime, the viscosity of the fluid simply functions as a multiplier for the force and traction. We can then take µ to be unity without loss of generality. Cell profiles from one period, with different flagellar beats, are shown in Figure 2a,b. In the case of Mode 1, the flagellar motion was similar to a traveling sinusoidal wave, whereas in Mode 2, the flagellar beat was more forceful at the tip. The amplitude of the beat in the middle part of the flagellum was comparable between the two modes. In both cases, the sperm cell exhibited a rolling motion before the trajectories began to follow a left-handed helicoid. The trajectories of the center-of-mass of the head in the (x,y)-plane are shown for 20 beats in Figure 2c. The initial swimming orientation was set to (1,0,0) for both beat modes. The head oscillated less in the direction perpendicular to the swimming motion. Although the direction in which the sperm swam was almost straight, small deviations in the y-direction were opposite in the two modes. This difference was caused by the initial condition. To facilitate quantitative discussion of the swimming speed, we define two types of average swimming speed. The swimming speed U was defined as the time average of the magnitude of the instantaneous translational velocity: where n represents the number of periods T. The latter swimming speed V was defined in terms of the displacement of the sperm cell over a sufficiently long time, shown with dotted arrows in Figure 2, that is,

Swimming Motion of a Sperm Cell in a Maxwell Fluid
Next, we considered a sperm cell swimming through a linear Maxwell liquid. Since the impact of the initial conditions decays over the timescale of the relaxation time λ, we limit our discussion to the physical values observed after the first period T, t ≥ T, which was sufficiently large for the initial conditions to decay.
The temporal flow fields and pressure distributions at time t = T with De = 0.0 and De = 1.0 are shown in Figure 3. Note that the pressure field was determined by the kernel shown in Equation (5). As shown in Figure 3, high pressure was observed in the forepart of the flagellum/head movement, while the negative pressure appeared at the rear in all cases. In the case of wave mode 1 with De = 0.0 (Figure 3a), the pressure was enhanced within the region of r/L < 0.1 along with the whole flagellum, that is, the red region in Figure 3, where r is the distance from the nearest flagellum/head surface and L is the flagellar length. By increasing the Deborah number, the high-pressure region shrank to near the surface (r/L < 0.05, see Figure 3b). Accordingly, the flow generated by the sperm's swimming diminished as the Deborah number increased. The same tendency was observed in wave Mode 2 (Figure 3c,d), although both flow and pressure fields were different between the two modes. In this study, the flagellum waveform was prescribed, that is, the beating speed was defined, and the traction was determined based on the time derivatives according to Equation (11). The traction force decreased with the flagellar velocity. As small traction forces cannot produce strong fluid flow, the flow was weak in the high-De regime. The swimming speeds U and V of the two beat modes are shown in Figure 2d. As the trajectories were almost straight, the values of U and V were similar between both beat modes, as shown in Figure 2d.

Swimming Motion of a Sperm Cell in a Maxwell Fluid
Next, we considered a sperm cell swimming through a linear Maxwell liquid. Since the impact of the initial conditions decays over the timescale of the relaxation time λ, we limit our discussion to the physical values observed after the first period T, t ≥ T, which was sufficiently large for the initial conditions to decay.
The temporal flow fields and pressure distributions at time t = T with De = 0.0 and De = 1.0 are shown in Figure 3. Note that the pressure field was determined by the kernel shown in Equation (5). As shown in Figure 3, high pressure was observed in the forepart of the flagellum/head movement, while the negative pressure appeared at the rear in all cases. In the case of wave mode 1 with De = 0.0 (Figure 3a), the pressure was enhanced within the region of r/L < 0.1 along with the whole flagellum, that is, the red region in Figure 3, where r is the distance from the nearest flagellum/head surface and L is the flagellar length. By increasing the Deborah number, the high-pressure region shrank to near the surface (r/L < 0.05, see Figure 3b). Accordingly, the flow generated by the sperm's swimming diminished as the Deborah number increased. The same tendency was observed in wave Mode 2 (Figure 3c,d), although both flow and pressure fields were different between the two modes. In this study, the flagellum waveform was prescribed, that is, the beating speed was defined, and the traction was determined based on the time derivatives according to Equation (11). The traction force decreased with the flagellar velocity. As small traction forces cannot produce strong fluid flow, the flow was weak in the high-De regime. Given that the swimming modes differ, why do the swimming speeds vary in the same way? In the next section, we considered the power and efficiency so that we could identify the mechanism that caused the swimming speed to decay.

Power and Swimming Efficiency
To see the effect of the Deborah number in more detail, we investigated the power generated by the locomotion. The power P due to the cellular locomotion is:  Figure 4. The value was averaged in the range 10T ≤ t ≤ 20T and normalized by the swimming speed through a Newtonian fluid U 0 . In the small-Deborah number region De << 1, the swimming speed converged to U 0 . Interestingly, we obtained similar results for the two different swimming modes, with both values decaying with De. Given that the swimming modes differ, why do the swimming speeds vary in the same way? In the next section, we considered the power and efficiency so that we could identify the mechanism that caused the swimming speed to decay. Given that the swimming modes differ, why do the swimming speeds vary in the same way? In the next section, we considered the power and efficiency so that we could identify the mechanism that caused the swimming speed to decay.

Power and Swimming Efficiency
To see the effect of the Deborah number in more detail, we investigated the power generated by the locomotion. The power P due to the cellular locomotion is:

Power and Swimming Efficiency
To see the effect of the Deborah number in more detail, we investigated the power generated by the locomotion. The power P due to the cellular locomotion is: The time-averaged powers for the different wave modes are shown in Figure 5a. For all values of De tested, the power P of Mode 1 was always higher than that of Mode 2. However, the two curves matched when we normalized the value by the power when De → 0 (i.e., at P 0 ), as shown in Figure 5b. We then investigated the efficiency of the locomotion. The efficiency was defined by the swimming speed per unit power, µfL 2 U/P, where µ is the viscosity, f is the beat frequency, and L is the flagellar length. The results are shown in Figure 6, where we can clearly see that the efficiency was almost independent of De, and the value became constant. When the efficiency was fixed, the tendency of the swimming speed to decay with De was less dependent on the wave pattern. Thus, the swimming speed varied in the same way in both wave modes.
Micromachines 2019, 10, x FOR PEER REVIEW 8 of 9 The time-averaged powers for the different wave modes are shown in Figure 5a. For all values of De tested, the power P of Mode 1 was always higher than that of Mode 2. However, the two curves matched when we normalized the value by the power when De → 0 (i.e., at P0), as shown in Figure 5b. We then investigated the efficiency of the locomotion. The efficiency was defined by the swimming speed per unit power, μfL 2 U/P, where μ is the viscosity, f is the beat frequency, and L is the flagellar length. The results are shown in Figure 6, where we can clearly see that the efficiency was almost independent of De, and the value became constant. When the efficiency was fixed, the tendency of the swimming speed to decay with De was less dependent on the wave pattern. Thus, the swimming speed varied in the same way in both wave modes.

Conclusion
In this study, we numerically investigated a sperm cell swimming in a linear Maxwell fluid in the small-De regime (De < 1.0). We found that, for the given waveform, the efficiency of the motion remained constant as the Deborah number varied. The fixed efficiency diminished the effect of the wave pattern on the decrease in swimming speed with De. These findings could be relevant to one-way fluid-structure interaction and may be helpful to researchers designing micro-swimmers in viscoelastic fluids with prescribed velocity conditions. For further study, to understand the physiology of sperm swimming, we must consider full fluid-structure interactions by developing a mechanical model of the inner structure of the flagellum.
The time-averaged powers for the different wave modes are shown in Figure 5a. For all values of De tested, the power P of Mode 1 was always higher than that of Mode 2. However, the two curves matched when we normalized the value by the power when De → 0 (i.e., at P0), as shown in Figure 5b. We then investigated the efficiency of the locomotion. The efficiency was defined by the swimming speed per unit power, μfL 2 U/P, where μ is the viscosity, f is the beat frequency, and L is the flagellar length. The results are shown in Figure 6, where we can clearly see that the efficiency was almost independent of De, and the value became constant. When the efficiency was fixed, the tendency of the swimming speed to decay with De was less dependent on the wave pattern. Thus, the swimming speed varied in the same way in both wave modes.

Conclusion
In this study, we numerically investigated a sperm cell swimming in a linear Maxwell fluid in the small-De regime (De < 1.0). We found that, for the given waveform, the efficiency of the motion remained constant as the Deborah number varied. The fixed efficiency diminished the effect of the wave pattern on the decrease in swimming speed with De. These findings could be relevant to one-way fluid-structure interaction and may be helpful to researchers designing micro-swimmers in viscoelastic fluids with prescribed velocity conditions. For further study, to understand the physiology of sperm swimming, we must consider full fluid-structure interactions by developing a mechanical model of the inner structure of the flagellum.

Conclusions
In this study, we numerically investigated a sperm cell swimming in a linear Maxwell fluid in the small-De regime (De < 1.0). We found that, for the given waveform, the efficiency of the motion remained constant as the Deborah number varied. The fixed efficiency diminished the effect of the wave pattern on the decrease in swimming speed with De. These findings could be relevant to one-way fluid-structure interaction and may be helpful to researchers designing micro-swimmers in viscoelastic fluids with prescribed velocity conditions. For further study, to understand the physiology of sperm swimming, we must consider full fluid-structure interactions by developing a mechanical model of the inner structure of the flagellum.

Funding:
The authors acknowledge the support of JSPS KAKENHI (Grants No. 18K18354 and 17H00853). We also acknowledge Mr. Eisaku Shimomura, who took the data in early stages of this study.