1. Introduction
Swimming is one of the oldest survival strategies in the microscopic theatre of life. From
Escherichia coli
(Berg & Brown Reference Berg and Brown1972; Sourjik & Berg Reference Sourjik and Berg2002; Saragosti et al. Reference Saragosti, Calvez, Bournaveas, Perthame, Buguin and Silberzan2011) foraging to spermatozoa navigating (Suarez & Pacey Reference Suarez and Pacey2006; Smith et al. Reference Smith, Gaffney, Gadêlha, Kapur and Kirkman-Brown2009; Denissenko et al. Reference Denissenko, Kantsler, Smith and Kirkman-Brown2012; Miki & Clapham Reference Miki and Clapham2013), microorganisms tens of micrometres in size operate in a world where inertia is negligible and viscous forces reign – Reynolds numbers as low as
$O(10^{-5}){-}O(10^{-2})$
(Brennen & Winet Reference Brennen and Winet1977; Purcell Reference Purcell1977). Many aquatic microorganisms, including zooplankton and certain bacteria, evade predators by achieving a swimming Reynolds number of approximately
$O(1)$
(Wickramarathna et al. Reference Wickramarathna, Noss and Lorke2014). This is influenced by fluid inertia on their movement and predator avoidance, which are both crucial for survival (Kiorboe et al. Reference Kiorboe, Jiang and Colin2010; Hamel et al. Reference Hamel, Fisch, Combettes, Dupuis-Williams and Baroud2011; Wang & Ardekani Reference Wang and Ardekani2012b
; Scholz et al. Reference Scholz, Jahanshahi, Ldov and Lowen2018).
In the past seventy years, a solid framework for studying microorganisms swimming in Newtonian fluids has been developed (Taylor Reference Taylor1951; Lighthill Reference Lighthill1952; Tuck Reference Tuck1968; Blake Reference Blake1971a , Reference Blakeb ; Wang & Ardekani Reference Wang and Ardekani2012a ; Khair & Chisholm Reference Khair and Chisholm2014; Ouyang et al. Reference Ouyang, Lin and Phan-Thien2022). In fact, real biological environments are rarely ideal Newtonian fluids. The mammalian reproductive tract, biofilm, DNA solution, plant vascular bundles, and cell cultures in industrial fermenters are all commonly rich in polymer or biomacromolecular networks that exhibit remarkable properties of viscoelasticity and anisotropy (Olsvik & Kristiansen Reference Olsvik and Kristiansen1994; Suarez & Pacey Reference Suarez and Pacey2006; Wang et al. Reference Wang, Liu, Huang and Ma2014; Duchesne et al. Reference Duchesne, Rainville and Galstian2015; Tung et al. Reference Tung, Lin, Harvey, Fiore, Ardon, Wu and Suarez2017; Chi et al. Reference Chi, Gavrikov, Berlyand and Aranson2022; Gmach & Opdenbosch Reference Gmach and Opdenbosch2022). These findings take the study of microorganisms’ propulsion into the realm of non-Newtonian complex fluids. Lauga (Reference Lauga2007) extended the work of Taylor (Reference Taylor1951) into Oldroyd-B viscoelastic fluids, which suggests that Oldroyd-B swimmers swim slower than Newtonian swimmers. For flexible swimmers, fluid elasticity can lead to a faster swimming speed (Espinosa-Garcia et al. Reference Espinosa-Garcia, Lauga and Zenit2013; Quan et al. Reference Quan, Yu, Lin, Ouyang, Xia and Lin2025a ). Furthermore, Quan et al. (Reference Quan, Zhang, Yu, Lin, Xia and Lin2025b ) found that the swimming speed of microorganisms in viscoelastic fluids exhibits a non-monotonic trend with Reynolds number, markedly different from the monotonic decrease observed in Newtonian fluids (Tuck Reference Tuck1968), which demonstrates the interplay of fluid inertia and viscoelasticity.
In anisotropic environments, swimmers exhibit more complex behaviours. The undulatory swimmer can swim either faster or slower in nematic liquid crystal than in a Newtonian fluid in different parameter regimes, and can even move in the opposite direction relative to that in an isotropic fluid if the rotational viscosity is large enough (Krieger et al. Reference Krieger, Spagnolie and Powers2014, Reference Krieger, Spagnolie and Powers2015, Reference Krieger, Spagnolie and Powers2019; Lintuvuori et al. Reference Lintuvuori, Würger and Stratford2017; Shi & Powers Reference Shi and Powers2017). Recent studies of lyotropic liquid crystal polymers (LCPs) by Lin et al. (Reference Lin, Chen and Gao2021, Reference Lin, Yu, Li and Gao2022, Reference Lin, Wang, Quan, Yu, Gao and Chen2025) examined alignment-dependent swimming and near-wall motions by combining asymptotic analysis with direct simulations. Here, lyotropic LCPs refer to a solution or suspension of rod-like LCPs described by a Q-tensor model based on Doi’s theory (Doi Reference Doi1981; Feng & Leal Reference Feng and Leal1997). Within this model, the equilibrium state can be either isotropic or nematic depending on the strength of molecular interactions (Lin et al. Reference Lin, Chen and Gao2021, Reference Lin, Yu, Li and Gao2022). In this bottom-up description, anchoring conditions are not required to impose molecular alignment, which differs from most studies based on classical Ericksen–Leslie or Landau–de Gennes descriptions of small-molecule liquid crystals (De Gennes & Prost Reference De Gennes and Prost1993; Larson Reference Larson1999). The nematic liquid crystal considered in Quan et al. (Reference Quan, Xia, Ouyang, Yu, Lin and Lin2026) refers to the nematic phase of the same Q-tensor model, and that work uses immersed-boundary simulations to study finite-length, large-amplitude swimmers with two asymmetric strokes, kicker and burrower, in nematic LCPs at finite Reynolds numbers, focusing on inertia-induced swimming dynamics and reorientation. The present work addresses a complementary question by developing a finite-Reynolds-number Taylor swimming sheet analysis for the same Q-tensor model. This analysis provides a unified theoretical framework for undulatory swimming in LCPs at finite Reynolds numbers, allowing us to isolate the coupling between finite-Reynolds-number effects and the phase-dependent microstructure of LCPs, and to clarify the different responses of the isotropic and nematic regimes to fluid inertia.
These considerations motivate the present study of undulatory swimming in lyotropic LCPs at finite Reynolds numbers by combining asymptotic analysis and direct numerical simulations. The structure of the paper is as follows. Section 2 introduces the mathematical model. Section 3 presents the analytical and numerical results for free swimming in lyotropic LCPs at finite Reynolds number. The findings are summarised in § 4, with Appendices A, B and C providing supplementary details regarding the numerical algorithm and asymptotic analysis.
2. Mathematical model
To study undulatory swimming in LCPs, we examine an infinite sheet, i.e. Taylor’s swimming sheet, described by
$y(x,t)=A\sin (kx-\omega t)$
, where
$A$
,
$k$
and
$\omega$
represent the amplitude, wavenumber and angular frequency of the wave, respectively. We rescale it as
$y(x,t)=\varepsilon \sin (x-t)$
with
$\varepsilon =Ak$
, by choosing
$1/k$
as the length scale, and
$1/\omega$
as the time scale. For the fluid phase, we build upon previous works (Lin et al. Reference Lin, Chen and Gao2021, Reference Lin, Yu, Li and Gao2022, Reference Lin, Wang, Quan, Yu, Gao and Chen2025) and use a Q-tensor model based on Doi’s theory (Doi Reference Doi1981; Doi & Edwards Reference Doi and Edwards1988) to characterise lyotropic LCPs. For non-dimensionalisation, we choose
$2\nu k_B T$
as the LCP’s stress scale, which leads to the non-dimensional system
Here,
$\boldsymbol{u}$
represents the fluid velocity,
$p$
is the fluid pressure, and
${\tau }_p$
is the polymer stress. The parameter
$\zeta$
indicates the strength coefficient of interaction between molecules of liquid crystal (Maier & Saupe Reference Maier and Saupe1958). The LCP’s field is described by a probability distribution function
$\varPsi (\boldsymbol x, \boldsymbol p, t)$
for concentrated micro-rods, defining their centre-of-mass position
$\boldsymbol{x}$
and orientation
$\boldsymbol{p}$
. The second-moment configuration tensor is defined as
${\unicode{x1D63F}}(\boldsymbol{x},t) = \int _{\mathbb{S}^{1}}\boldsymbol{p}\boldsymbol{p}\,\varPsi (\boldsymbol{x},\boldsymbol{p},t)\,\mathrm{d}\boldsymbol{p}$
, and the fourth-moment tensor as
${\unicode{x1D64E}}(\boldsymbol{x},t) = \int _{\mathbb{S}^{1}}\boldsymbol{p}\boldsymbol{p}\boldsymbol{p}\boldsymbol{p}\,\varPsi (\boldsymbol{x},\boldsymbol{p},t)\,\mathrm{d}\boldsymbol{p}$
, where the integrals are taken over all orientations
$\boldsymbol{p}$
on the unit circle
$\mathbb{S}^{1}$
; the strain-rate tensor is
${{\unicode{x1D640}}}=({\boldsymbol{\nabla }}\boldsymbol{u}+{\boldsymbol{\nabla }}\boldsymbol{u}^{\rm T})/2$
. The Reynolds number is defined as
$\textit{Re}=\rho _{\!f}\omega /(\mu _sk^2)$
, with
$\rho _{\!f}$
denoting the fluid density, and
$\mu _s$
the fluid solvent viscosity. The Ericksen number
$\textit{Er}=2\nu k_B T/(\mu _s\omega )$
characterises the coupling between the LCPs and the viscous fluid, where
$\nu \sim l^{-3}$
(rod molecule length) represents the effective volume fraction,
$k_B$
is the Boltzmann constant, and
$T$
is the absolute temperature. The two Péclet numbers,
$\textit{Pe}=\omega /{4d_r}$
and
$\textit{Pe}_t= \omega /(d_tk^2)$
, represent the ratios of the time scales for the rod’s rotation and transport, arising from the translational and rotational diffusivity coefficients
$d_t$
and
$d_r$
, respectively. Notably,
$ \textit{Pe}_t \sim O(10){-}O(100)$
is chosen to be at least one order of magnitude greater than
$\textit{Pe} \sim O(10^{-1}){-}O(10)$
, so that the translational diffusion effect is negligible. In addition, a small effective crowdedness coefficient
$ \beta = \beta _0 \omega /({4 d_r (\nu l^3)^{2}}) \sim O(10^{-2}){-}O(10^{-3})$
is chosen in the extra stress. The detailed numerical algorithm and full derivation of the asymptotic analysis in the present work are given in Appendices A and B.
3. Results and discussion
Schematic of the three configurations considered in this work: (a) isotropic LCPs, (b) nematic LCPs with
$\parallel$
alignment, and (c) nematic LCPs with
$\perp$
alignment. The pink solid line represents Taylor’s swimming sheet.

We first conduct an asymptotic analysis, and expand all variables under the small-amplitude condition
$\varepsilon \ll 1$
, deriving the leading-order solutions at
$O(\varepsilon ^2)$
. In the present two-dimensional Q-tensor LCP model, the interaction parameter
$\zeta$
determines the equilibrium phase, with the isotropic–nematic transition occurring at the critical value
$\zeta _c=4$
: the isotropic regime corresponds to
$\zeta \lt \zeta _c$
, whereas the nematic regime corresponds to
$\zeta \gt \zeta _c$
. For the nematic regime, we further consider the deep-nematic case
$\zeta \gg 1$
in deriving the corresponding asymptotic results. It is important to mention that the fluid velocity close to the sheet is of the order of
$\sqrt {\textit{Re}}\,\varepsilon$
when
$\textit{Re}$
is large, and there exists a boundary-layer structure with thickness
$O(\textit{Re}^{-1/2})$
in relation to the wavelength (Childress Reference Childress2008). Our analytical findings should be interpreted as asymptotic in
$\varepsilon$
for fixed
$\textit{Re}$
, and are not applicable uniformly for large
$\textit{Re}$
. Therefore, in the present work, we focus on the cases
$\textit{Re}\sim O(1)$
to ensure that the velocity near the sheet is uniform to one order in
$\varepsilon$
. The full derivation can be found in Appendix B. Meanwhile, we perform direct simulations of an infinite-length swimmer undergoing a small amplitude (
$\varepsilon \sim O(10^{-2})$
) in LCPs at finite Reynolds numbers to further validate our analytical predictions. These simulations solve the equations using the immersed-boundary method described in Appendix A. The infinite sheet is represented in a periodic computational domain containing one wavelength of the prescribed travelling wave, and the mean swimming speed is obtained from the time-averaged translational velocity of the sheet after the swimming motion reaches a steady state. As shown in figure 1, the three swimming configurations considered in this work include the isotropic regime and two equilibrium cases in the nematic regime; in figures 1(b,c), the symbols
$\parallel$
and
$\perp$
denote nematic directors parallel and perpendicular to the swimming direction, respectively.
(a) Asymptotic results and (b) numerical results of mean-speed ratio
$U_{\textit{LC}}/U_{\textit{LC}, \textit{Re}=0}$
as a function of Reynolds number
$\textit{Re}$
, with the inset in (a) showing the initial increase in swimming speed in the small-
$\textit{Re}$
region. (c) The Reynolds number that corresponds to the maximum speed (
$\textit{Re}_{\textit{max}}$
) for varied Péclet number in the isotropic regime. These parameters are fixed:
$\varepsilon =0.02\pi$
,
$\zeta =2$
in the isotropic regime.

In the isotropic regime, we obtain the time-averaged swimming speed as
\begin{equation} U_{\textit{LC}}=\left\{\frac {F}{2}-\dfrac {\dfrac {\textit{Re}G(2F+1)}{(F+1)^2+G^2}-\textit{Er}\textit{Pe}\,\dfrac {2G(F+1)(4-\zeta )\,\textit{Pe}-8F^2\,\textit{Pe}^2}{16\,\textit{Pe}^2+(4-\zeta )^2}}{(4+\textit{Er}\textit{Pe})F}\right\}\varepsilon ^2+o(\varepsilon ^2)\quad \mathrm{(isotropic)}, \end{equation}
where
$F$
and
$G$
are constants involving
$\textit{Er}$
,
$\textit{Pe}$
and
$\zeta$
, with their explicit forms given in Appendix B. To examine the effects of fluid inertia on swimming in LCPs, we compare the swimming speed with that obtained in the Stokesian regime, i.e.
$U_{\textit{LC}, \textit{Re}=0}$
, introducing the speed ratio
$U_{\textit{LC}}/U_{\textit{LC}, \textit{Re}=0}$
as a function of
$\textit{Re}$
. Figures 2(a,b) demonstrate good agreement between analytical and numerical predictions in the isotropic regime, confirming the reliability of our asymptotic solutions and numerical algorithm. The results in figures 2(a,b), with the inset in figure 2(a), show a viscoelastic-like inertial response similar to that observed in Oldroyd-B fluids (Quan et al. Reference Quan, Zhang, Yu, Lin, Xia and Lin2025b
): the swimming speed first increases and then decreases with increasing
$\textit{Re}$
, indicating a non-monotonic dependence of the swimming speed on the Reynolds number. For small
$\textit{Er}$
or
$\textit{Pe}$
, the initial speed increase occurs in a small-
$\textit{Re}$
range shown in the inset, whereas for larger
$\textit{Er}$
, the subsequent decrease is shifted to larger Reynolds numbers, typically at
$\textit{Re}\gt 1$
. More specifically, for larger values of
$\textit{Er}$
and
$\textit{Pe}$
, the swimming speed can be consistently enhanced as
$\textit{Re}$
increases, as shown by the case with
$\textit{Er}=5.0$
and
$\textit{Pe}=1.0$
, where the enhancement reaches approximately 10 %. This same non-monotonic dependence can be understood from (3.1) by incorporating an effective polymer viscosity
$\mu _p$
(Doi & Edwards Reference Doi and Edwards1988; Feng & Leal Reference Feng and Leal1997) as
$\mu _p/\mu _s=\alpha\, \textit{Er}\textit{Pe}$
. Thus
$\textit{Er}\textit{Pe}$
defines the effective polymer contribution to the zero-shear-rate viscosity, and
$\alpha$
is an order-dependent coefficient. Particularly for the special case
$\zeta =0$
(the limit without molecular orientational interactions), we estimate
$\alpha =1$
(Lin et al. Reference Lin, Chen and Gao2021). In this case, (3.1) has the same viscosity ratio structure as the finite-Reynolds-number Oldroyd-B swimming-speed result, given in (2.22) of Quan et al. (Reference Quan, Zhang, Yu, Lin, Xia and Lin2025b
). This correspondence can be further seen by rewriting
$4/(4+\textit{Er}\textit{Pe})$
as
$4\mu _s/(4\mu _s+\mu _p)$
, which identifies the effective viscosity ratio in the isotropic LCP. Together, these observations indicate that the behaviour of the swimmer is still dominated by the fluid viscoelasticity in the isotropic regime.
To further understand the above results, we expand (3.1) in terms of
$\textit{Re}$
as
The coefficients
$H_i$
(
$i=0,1,2$
) in (3.2) depend on
$\textit{Er}$
,
$\textit{Pe}$
and
$\zeta$
, with their explicit forms given in Appendix C. The zeroth-order term recovers the Stokesian-regime result of our previous work (Lin et al. Reference Lin, Chen and Gao2021). The leading-order inertial correction is the
$O(\textit{Re})$
term, whose coefficient is always positive (i.e.
$H_1\gt 0$
); therefore, the linear response is always positive, and the swimming speed initially increases with increasing
$\textit{Re}$
. We associate this initial speed enhancement with the coupled effect of fluid inertia and the LCP properties in the isotropic regime, which corresponds to a viscoelastic-like response (Quan et al. Reference Quan, Zhang, Yu, Lin, Xia and Lin2025b
). The
$O(\textit{Re}^2)$
contribution is negative for moderate polymer effects, and can offset the positive linear contribution, leading to a subsequent decrease in swimming speed at finite
$\textit{Re}$
. For larger
$\textit{Er}$
and
$\textit{Pe}$
, however, this term can become positive, so that the speed enhancement persists over a wider range of Reynolds numbers. This explains why fluid inertia can still consistently increase the swimming speed at
$\textit{Re}\leqslant 1$
for the case with
$\textit{Er}=5.0$
and
$\textit{Pe}=1.0$
in figure 2(a). Thus the non-monotonic dependence observed in the finite-Reynolds-number result reflects the competition between the initial inertial enhancement and higher-order finite-
$\textit{Re}$
corrections.
The above discussion shows that in the isotropic regime, fluid inertia can contribute to an enhancement of the swimming speed over an appropriate range of
$\textit{Re}$
, similar to the behaviour observed in viscoelastic fluids. It is intriguing to explore the Reynolds numbers
$\textit{Re}_{\textit{max}}$
that correspond to the maximum speed at different Péclet number, as demonstrated in figure 2(c). It can be observed that as the Péclet number increases, the swimming speed reaches its maximum at larger Reynolds numbers. This variation of
$\textit{Re}_{\textit{max}}$
reflects how
$\textit{Er}$
and
$\textit{Pe}$
modify the inertial response of the swimmer through the polymer contribution. When
$\textit{Er}\rightarrow 0$
or
$\textit{Pe}\rightarrow 0$
(polymeric effect vanishes), the Q-tensor model reduces to the Newtonian model, and the maximum speed is recovered in the Stokesian regime (Tuck Reference Tuck1968; Childress Reference Childress2008). As
$\textit{Er}$
or
$\textit{Pe}$
increases, the positive inertial correction in (3.2) persists over a wider range of
$\textit{Re}$
, so that the maximum speed is shifted to finite Reynolds numbers. This does not mean that the absolute swimming speed increases monotonically with
$\textit{Er}$
or
$\textit{Pe}$
; rather, the normalised inertial enhancement becomes more pronounced when the polymeric contribution is stronger. To quantify this effect further, we examine the maximum speed ratio
$U_{\textit{LC}, \textit{max}}/U_{\textit{LC},\textit{Re}=0}$
corresponding to
$\textit{Re}_{\textit{max}}$
at different
$\textit{Er}$
and
$\textit{Pe}$
. For larger
$\textit{Er}$
, fluid inertia contributes more strongly to the relative speed enhancement, leading to a larger value of
$U_{\textit{LC}, \textit{max}}/U_{\textit{LC},\textit{Re}=0}$
. The maximum speed ratio also increases monotonically with
$\textit{Pe}$
(not shown here). In this sense, the polymer elastic contribution of the isotropic LCPs does not directly imply a larger absolute swimming speed; Lin et al. (Reference Lin, Chen and Gao2021) showed that the swimming speed in isotropic LCPs decreases with increasing
$\textit{Pe}$
in the Stokesian regime. Related low-Reynolds-number studies also found that polymer stresses can reduce Taylor-sheet swimming in Oldroyd-B fluids, or enhance flexible-swimmer locomotion in viscoelastic fluids (Lauga Reference Lauga2007; Teran et al. Reference Teran, Fauci and Shelley2010; Espinosa-Garcia et al. Reference Espinosa-Garcia, Lauga and Zenit2013). The present work therefore isolates how fluid inertia modifies the swimming speed relative to the corresponding Stokesian-regime case, showing that the polymer elastic contribution amplifies the inertial enhancement at finite Reynolds numbers.
In nematic LCPs, the swimming speeds for the two equilibrium states are obtained as
\begin{equation} U_{\textit{LC}}=\begin{cases} \left \{\dfrac {1+f}{4f}+{\dfrac {\textit{Er}\textit{Pe}}{\zeta }\dfrac {2f^3+2f^2+f-1}{4f^3}}\right \}\varepsilon ^2+o(\varepsilon ^2)&{(\parallel )},\\[15pt] \left \{\dfrac {1+f}{4f}{-\dfrac {\textit{Er}\textit{Pe}}{\zeta }\dfrac {(f-1)(2f^2-1)}{4f^3}}\right \}\varepsilon ^2+o(\varepsilon ^2)&{(\perp )} \end{cases}\quad \mathrm{(nematic)}. \end{equation}
In the above expression, the coefficient
$f$
depends on
$\textit{Re}$
, as given in Appendix B. Likewise, figures 3(a,b) show good agreement between analytical and numerical results in the nematic regime, demonstrating the reliability of our prediction. The speed ratio
$U_{\textit{LC}}/U_{\textit{LC}, \textit{Re}=0}$
of perpendicular gait is slightly lower than that for the parallel case. Besides, the swimming speed of parallel and perpendicular gaits decreases monotonically with Reynolds number. Surprisingly, it is similar to behaviours observed in Newtonian fluids (Tuck Reference Tuck1968; Childress Reference Childress2008). We then expand swimming speed in the nematic regime (3.3) in terms of
$\textit{Re}$
for further understanding:
\begin{align} U_{\textit{LC}}=\begin{cases} \left \{\dfrac {1}{2}+\dfrac {\textit{Er}\textit{Pe}}{\zeta }-\dfrac {1}{32}\left (1+\dfrac {\textit{Er}\textit{Pe}}{\zeta }\right )\textit{Re}^2\right \}\varepsilon ^2+o(\varepsilon ^2\,\textit{Re}^2)& {(\parallel )},\\[15pt] \left \{\dfrac {1}{2}-\dfrac {1}{32}\left (1+\dfrac {\textit{Er}\boldsymbol{Pe}}{\zeta }\right )\textit{Re}^2\right \}\varepsilon ^2+o(\varepsilon ^2\,\textit{Re}^2)&{(\perp )} \end{cases}\quad \mathrm{(nematic)}. \end{align}
The speed expansion in (3.4) contains only even powers of
$\textit{Re}$
, with no odd-order contributions. The first order in
$\textit{Re}$
therefore vanishes, and the leading-order term begins at
$O(\textit{Re}^2)$
and is always negative, leading to the swimming speed decreasing with Reynolds number for both parallel and perpendicular cases. It resembles the asymptotic solution derived for Newtonian fluids,
$U_{N}=\{({1}/{2})-({1}/{32})\,\textit{Re}^2\}\varepsilon ^2+o(\varepsilon ^2\,\textit{Re}^2)$
(Tuck Reference Tuck1968; Childress Reference Childress2008). This alignment-dependent difference is consistent with previous studies of swimming in nematic liquid crystals, where swimmer behaviour can be strongly affected by the nematic director configuration and liquid crystal anisotropy (Krieger et al. Reference Krieger, Spagnolie and Powers2015; Lintuvuori et al. Reference Lintuvuori, Würger and Stratford2017; Shi & Powers Reference Shi and Powers2017; Lin et al. Reference Lin, Yu, Li and Gao2022). The present results further show that in the finite-Reynolds-number cases considered here, both parallel and perpendicular nematic swimming exhibit a Newtonian-like inertial response.
(a) Asymptotic results and (b) numerical results of mean-speed ratio
$U_{\textit{LC}}/U_{\textit{LC}, \textit{Re}=0}$
as a function of Reynolds numbers
$\textit{Re}$
. (c) The rescaled
$\Delta U_{\textit{LC}}$
as a function of Reynolds numbers
$\textit{Re}$
. These parameters are fixed:
$\varepsilon =0.02\pi$
,
$\zeta =50$
in the nematic regime.

Figure 3. Long description
Panel A: A line graph titled ‘Analytical’ shows the mean-speed ratio as a function of Reynolds numbers. The x-axis is labeled ‘Re’ and the y-axis is labeled ‘U_LC/U_LC,Re = 0’. Multiple lines represent different values of Er Pe (0.1, 1.0, 5.0) and orientations (parallel and perpendicular), along with a line for Newtonian fluid. Panel B: A line graph titled ‘Numerical’ also shows the mean-speed ratio as a function of Reynolds numbers. The axes and legends are the same as in Panel A. Panel C: A line graph titled ‘Analytical’ shows the rescaled parameter as a function of Reynolds numbers. The x-axis is labeled ‘Re’ and the y-axis is labeled ‘ΔU_LC/(Er Pe (ε^2 δ))’. Lines represent parallel case, perpendicular case, and Newtonian fluid.
Furthermore, a comparison is made between behaviours in nematic LCPs and Newtonian fluids. We find that the swimming speed of perpendicular swimming is slightly lower than that of parallel swimming, with Newtonian swimming speed lying between the two, which is similar to those observed in the Stokesian regime (Lin et al. Reference Lin, Chen and Gao2021). To quantify the effect of fluid inertia on LCPs, we propose the inertia-induced anisotropic speed
$\Delta U_{\textit{LC}}={U_{\textit{LC}}-U_{N}}$
to distinguish them from Newtonian fluids. Surprisingly, by rescaling this quantity with
$\textit{Er}\textit{Pe}$
, the swimming speed of parallel motions (or perpendicular motions) is unified, merging into a single line, which highlights the additional impact arising specifically from the coupling between the inertial effect and the anisotropic property. In figure 3(c), the
$\Delta U_{\textit{LC}}$
of parallel swimming decreases with Reynolds numbers, but is positive. This corresponds to the enhanced swimming speed of parallel swimming (
$U_{\textit{LC}}\gt U_{N}$
), and indicates that fluid inertia dampens this increase in speed. This is consistent with the trend of speed decreasing monotonically with Reynolds number. As for the perpendicular case, it is negative (except for case
$\textit{Re}=0$
) and also decreases with increasing Reynolds numbers, which denotes the retarded swimming motion at finite Reynolds numbers, i.e.
$U_{\textit{LC}}\lt U_{N}$
. This also reflects that perpendicular swimming speed decays more rapidly with Reynolds numbers under the coupled effect of fluid inertia and the LCP.
The pressure force acting on the sheet as a function of Reynolds number
$\textit{Re}$
: (a) isotropic regime; (b) nematic regime. (c) The rescaled
$\Delta F^{pres}_{\textit{LC}}$
as a function of Reynolds numbers
$\textit{Re}$
.

To further clarify the mechanism underlying the inertial swimming dynamics in LCPs, we next examine the hydrodynamic force acting on the sheet. In Appendix B, the time-averaged
$x$
-direction contributions from the pressure, viscous stress and polymer stress are calculated separately, as given in (B33) and (B73), and are denoted by
$F_{x}^{pres/visc/poly}$
. We then project these components along the swimming direction of the swimmer, i.e.
$F^{pres/visc/poly}_{\textit{LC}}=-F_{x}^{pres/visc/poly}$
, and refer to the three projected components below as the pressure force, viscous force and elastic force associated with the polymer stress, respectively. In the isotropic regime, the magnitude of the pressure force increases with Reynolds number, as shown in figure 4(a). Besides, it is always positive and thereby functions as a thrust force and propels the sheet forwards. The viscous and elastic forces are always negative and therefore act as resistance, counterbalancing the thrust. It should be noted that their magnitudes still depend on respective effective viscosity ratios, i.e.
$F^{visc}_{\textit{LC}}=-({4}/({\textit{Er}\textit{Pe}+4}))F^{pres}_{\textit{LC}},\ F^{poly}_{\textit{LC}}=-( {\textit{Er}\textit{Pe}}/({\textit{Er}\textit{Pe}+4}))F^{pres}_{\textit{LC}}$
, indicating that the isotropic regime exhibits a mechanical response similar to that in viscoelastic fluids (Quan et al. Reference Quan, Zhang, Yu, Lin, Xia and Lin2025b
). For small values of
$\textit{Er}\textit{Pe}$
, the elastic force is close to 0. With increasing
$\textit{Er}$
or
$\textit{Pe}$
, the magnitude of elastic force increases, accompanied by an increase in the equilibrium pressure.
As shown in figure 4(b), in the nematic regime, the pressure force acting on the parallel motions equals that in the perpendicular case when higher-order terms are neglected. In both cases, the magnitude of the pressure force increases with
$\textit{Re}$
, and at fixed
$\textit{Re}$
, it increases with
$\textit{Er}\textit{Pe}$
. Its value is always positive, representing a thrust. Conversely, the viscous force is negative, acting as a resistive force to counterbalance it. Especially, the magnitude of the elastic force exactly equals
$0$
when parallel (or perpendicular) swimming, which means
$F^{visc}_{\textit{LC}}=-F^{pres}_{\textit{LC}}$
ignoring higher-order terms. This indicates that the swimmer is primarily subject to the influence of pressure and viscous forces, similar to those in Newtonian fluids. Therefore, we believe that the phase transition of the LCP fundamentally alters its mechanical response. Corresponding to the swimming speed in figure 3(c), we introduce the inertia-induced anisotropic force
$\Delta F^{pres}_{\textit{LC}}={F^{pres}_{\textit{LC}}-F^{pres}_{N}}$
in figure 4(c) to distinguish between the pressure force experienced in nematic LCP swimming and that in Newtonian fluids, where
$F^{pres}_{N}$
denotes the pressure force acting on the swimmer in Newtonian fluids. We rescale this quantity with
$\textit{Er}\textit{Pe}$
, and focus on the coupled effect of fluid inertia and LCP on the hydrodynamic force. The results in figure 4(c) demonstrate that
$\Delta F^{pres}_{\textit{LC}}$
values acting on the parallel and perpendicular swimming are unified. Besides, it is always positive and increases with Reynolds number, which indicates that swimmers (both parallel and perpendicular motions) are subjected to stronger thrust in nematic LCPs.
The flow patterns around the swimming sheet at
$t=0.5T_{c}$
, for (a–c)
$\zeta =2$
, (d–f)
$\zeta =50$
. The black and purple solid lines represent the infinite-length swimming sheet and streamlines. The background colour corresponds to the magnitude of polymer stress.

Figure 5. Long description
Panel A: A heat map showing flow patterns around a swimming sheet with Er = 2.0, Pe = 0.1, Re = 0. The background color indicates the magnitude of polymer stress, ranging from 0 to 0.2. Black and purple solid lines represent the infinite-length swimming sheet and streamlines, respectively. The heat map is isotropic. Panel B: A heat map showing flow patterns around a swimming sheet with Er = 2.0, Pe = 0.1, Re = 1.0. The background color indicates the magnitude of polymer stress, ranging from 0 to 0.2. Black and purple solid lines represent the infinite-length swimming sheet and streamlines, respectively. The heat map is isotropic. Panel C: A heat map showing flow patterns around a swimming sheet with Er = 2.0, Pe = 1.0, Re = 1.0. The background color indicates the magnitude of polymer stress, ranging from 0 to 0.6. Black and purple solid lines represent the infinite-length swimming sheet and streamlines, respectively. The heat map is isotropic. Panel D: A heat map showing flow patterns around a swimming sheet with Er = 2.0, Pe = 0.1, Re = 0. The background color indicates the magnitude of polymer stress, ranging from 0 to 0.0150. Black and purple solid lines represent the infinite-length swimming sheet and streamlines, respectively. The heat map is nematic. Panel E: A heat map showing flow patterns around a swimming sheet with Er = 2.0, Pe = 0.1, Re = 1.0. The background color indicates the magnitude of polymer stress, ranging from 0 to 0.0150. Black and purple solid lines represent the infinite-length swimming sheet and streamlines, respectively. The heat map is nematic. Panel F: A heat map showing flow patterns around a swimming sheet with Er = 2.0, Pe = 1.0, Re = 1.0. The background color indicates the magnitude of polymer stress, ranging from 0 to 0.10. Black and purple solid lines represent the infinite-length swimming sheet and streamlines, respectively. The heat map is nematic.
The swimming efficiency ratio
$\eta _{\textit{LC}}/\eta _{\textit{LC},\textit{Re}=0}$
of an infinite-length sheet as a function of Reynolds number
$\textit{Re}$
. (a,c) Asymptotic solutions of Taylor’s swimming sheet, with the inset in (a) magnifying the initial efficiency enhancement in the small-
$\textit{Re}$
region. (b,d) Results of numerical simulations. These parameters are fixed:
$\varepsilon =0.02\pi$
,
$\zeta =2(50)$
in the isotropic (nematic) regime.

Subsequently, based on our analytical results, we examine the flow patterns around the swimming sheet for various
$\textit{Pe}$
and
$\textit{Re}$
, which are described in figure 5. As
$\textit{Re}=0$
, i.e. Stokes flow, the flow pattern appears simple and symmetric in both the isotropic and nematic regimes. By comparing figure 5(a) with figure 5(d), it can be seen that the polymer stress distributions differ markedly between the isotropic and nematic regimes. In the isotropic regime, high-stress regions are uniformly distributed around the sheet, whereas in the nematic regime, they concentrate predominantly at locations of maximum
$y$
-direction velocity of the sheet (i.e. at
$y = 0$
). A comparison between figures 5(a,b) and 5(d,e) reveals that fluid inertia mainly changes the flow patterns of the flow field distribution, causing the streamlines to be deflected. Meanwhile, the high-polymer-stress regions also undergo a spatial shift, clustering in the wake behind the swimmer. In addition, the magnitude of polymer stress increases with increasing
$\textit{Pe}$
, as shown by comparing figures 5(b,c) and 5(e, f), corresponding to the results in figure 4. Taking fluid inertia into account, we conclude that the flow field lags behind, to some extent, the direction in which the swimmer is swimming, which exhibits nonlinear characteristics. These results indicate that the nonlinear regime at finite Reynolds numbers disrupts the symmetry in the original Stokes regimes.
Finally, we analyse the swimming efficiency based on its time-averaged speed and power, which are detailed in Appendices A and B. The swimming efficiency ratio
$\eta _{\textit{LC}}/\eta _{\textit{LC},\textit{Re}=0}$
is introduced in figure 6. Here,
$\eta _{\textit{LC}}$
denotes the mechanical swimming efficiency in LCPs at finite Reynolds numbers, while
$\eta _{\textit{LC},\textit{Re}=0}$
denotes the corresponding efficiency in the Stokesian regime. Therefore,
$\eta _{\textit{LC}}/\eta _{\textit{LC},\textit{Re}=0}$
measures the inertial contribution to swimming efficiency. Figures 6(a,b) and 6(c,d) demonstrate that the analytical solutions are in excellent agreement with the numerical predictions. In the isotropic regime, the efficiency shows the same inertial response as the swimming speed, and first increases before decreasing with increasing
$\textit{Re}$
, indicating a non-monotonic dependence on the Reynolds number similar to swimming in viscoelastic fluids. For small
$\textit{Er}$
or
$\textit{Pe}$
, the efficiency enhancement occurs in the small-
$\textit{Re}$
region highlighted by the inset in figure 6(a), whereas for larger
$\textit{Er}$
, the enhancement can persist to
$\textit{Re}\gt 1$
. Thus fluid inertia can contribute to an enhancement of the efficiency over an appropriate range of
$\textit{Re}$
in isotropic LCPs, rather than producing a monotonic increase. However, in the nematic regime, the swimming efficiency decreases monotonically with increasing Reynolds number, like the behaviour in Newtonian fluids, with the efficiency of perpendicular motions slightly lower than that of parallel swimming, and Newtonian swimming efficiency lying between the two. These findings suggest that different phase states of lyotropic LCPs give rise to markedly distinct responses of the undulatory swimmer to fluid inertia.
4. Conclusion
In the present work, we extend previous studies (Lin et al. Reference Lin, Chen and Gao2021, Reference Lin, Yu, Li and Gao2022) of undulatory swimming in lyotropic liquid crystal polymers (LCPs) to the inertial regime. Combining asymptotic analysis and direct numerical simulations, we show that fluid inertia can enhance the swimming speed and efficiency of an undulatory swimmer in isotropic LCPs at finite Reynolds numbers, similar to behaviour in viscoelastic fluids. By contrast, in nematic LCPs, the swimming speed and efficiency decrease monotonically with Reynolds number for both parallel and perpendicular swimming, exhibiting Newtonian-like behaviour. These results connect the inertial swimming behaviour in LCPs with existing studies of complex-fluid swimming. The isotropic phase exhibits viscoelastic-like swimming behaviour through effects of polymer, whereas the nematic phase retains alignment-dependent differences but shows a Newtonian-like inertial response. Analysis of the hydrodynamic force acting on the sheet suggests that the phase state of an LCP fundamentally alters its mechanical response. More complex phenomena during anisotropic swimming in the inertial regime will be a key focus of our future work, motivated by the observation that fluid inertia can disrupt the symmetric flow patterns observed at low Reynolds numbers.
Funding
The work was supported by the National Natural Science Foundation of China (grant nos 12132015, 12332015, 12302340).
Declaration of interests
The authors report no conflict of interest.
Appendix A. Numerical method
The numerical simulations used to validate the asymptotic results solve (2.1)–(2.4) with an immersed-boundary representation of the swimming sheet. In the numerical implementation, the momentum equation in (2.2) is augmented by an immersed-boundary body force
$\boldsymbol{f}_{\!e}$
, which transfers the elastic force generated by the prescribed sheet motion to the surrounding fluid. The overall algorithm follows the immersed-boundary formulation and time-stepping strategy used by Quan et al. (Reference Quan, Xia, Ouyang, Yu, Lin and Lin2026); it is summarised here to specify the ingredients needed for the present validation.
The swimmer considered in this work is a classic Taylor’s swimming sheet
$y=A{}\sin (kx-\omega t)$
, initially positioned along the
$x$
axis, which in the Lagrangian frame creates a rightward wave-like curvature on the target body as
where
$A$
is amplitude,
$k$
is wavenumber and
$\omega$
is angular frequency. Its kinematics in the parametric form is given as
$\boldsymbol{X}(s,t)$
, with local arc length
$s\in [0, L_s]$
. To simulate its motion, we set it in a rectangular box with periodic boundary conditions in the
$x$
-direction. According to the immersed-boundary approach (Peskin Reference Peskin2002), the Lagrangian body force
$\boldsymbol{F}_{\!e}$
resulting from the imposed actuation (A1) can be derived by taking the variational derivative of the elastic energy
$E$
(Fauci & Peskin Reference Fauci and Peskin1988):
The elastic body force
$\boldsymbol{F}_{\!e}$
is calculated explicitly at each time step. The elastic energy is defined as
\begin{equation} E[\boldsymbol{X}(s,t)]=\frac {\sigma _{s}}{2}\int _{\varOmega _{L}}\left (\left |\frac {\partial \boldsymbol{X}}{\partial s}\right |-1\right )^{2}\mathrm{d}s+\frac {\sigma _{b}}{2}\int _{\varOmega _{L}}\left (\frac {\partial ^{2}\boldsymbol{X}}{\partial s^{2}}\boldsymbol{\cdot }\boldsymbol{n}-\kappa _{0}\right )^{2}\mathrm{d}s, \end{equation}
where
$\boldsymbol{n}$
denotes the local normal direction. The total elastic energy has two parts: one from stretching with modulus
$\sigma _s$
, and the other from bending with modulus
$\sigma _b$
. The Lagrangian force is spread to the Eulerian grid, and the sheet velocity is interpolated back from the surrounding fluid velocity, using the same compact bilinear kernel as in Quan et al. (Reference Quan, Zhang, Yu, Lin, Xia and Lin2025b
, Reference Quan, Xia, Ouyang, Yu, Lin and Lin2026). The numerical swimming efficiency is calculated from the numerical time-averaged speed
$U_{\textit{LC}}$
and the mechanical power:
\begin{equation} \eta _{\textit{LC}}=\frac {U_{\textit{LC}}^{2}}{\big \langle ({1}/{L_{s}})\int _{{L}_s}{\boldsymbol{U}}\boldsymbol{\cdot }{\boldsymbol{F}}_{\!e}\,{\mathrm{d}}s\big \rangle }. \end{equation}
For the fluid phase, a fractional-step projection method (Yu Reference Yu2005) with a second-order half-staggered finite-difference discretisation is used in the fixed Eulerian domain; the resulting Helmholtz and Poisson problems are solved by the ADI technique (Namiki Reference Namiki1999). For the configuration tensor constitutive equation, we use Cholesky decomposition (Vaithianathan & Collins Reference Vaithianathan and Collins2003) to preserve the symmetric positive definiteness of the configuration tensor, followed by a fourth-order Runge–Kutta update. An arctangent change of variables is further applied to the Cholesky factors so that the diagonal components of the configuration tensor remain bounded, i.e.
$0\leqslant D_{11},D_{22}\leqslant 1$
. The detailed discrete equations are the same as those in Quan et al. (Reference Quan, Xia, Ouyang, Yu, Lin and Lin2026) and are therefore not repeated here.
Convergence studies of an infinite swimming sheet free swimming in the isotropic regime with the time-dependent velocity
$\boldsymbol{U}(U_x, U_y)$
by varying: (a) domain size, (b) Eulerian grid width. These parameters are fixed:
$\textit{Re}\approx 1.0$
,
$\textit{Er}=2.0$
,
$\textit{Pe}=1.0$
,
$\zeta=2$
.

We perform convergence tests for an infinite-length sheet free swimming in LCPs at finite Reynolds numbers as shown in figure 7, where we examine the time-dependent velocity components, and demonstrate good agreement by varying the grid width, time step and domain size separately, confirming the reliability of the numerical approach. The detailed computational parameters are given in table 1.
Computational parameters.

Table 1. Long description
A table with two column and sixteen rows listing computational parameters. The table includes the following parameters: Domain size, Eulerian grid width, Time step, Swimmer grid width, Reynolds number, Ericksen number, Péclet number (rotation), Péclet number (transport), Strength coefficient, Effective crowdness coefficient, Swimmer length, Amplitude, Wavenumber, Angular frequency, Bending stiffness, and Stretching stiffness. Each parameter is listed with its corresponding value or range of values.
Appendix B. Asymptotic analysis
To determine the effect of fluid inertia on microorganism swimming in LCPs, we employ the classical model of Taylor’s swimming sheet, which is undulating with a propagating wave of function
$y(x,t)=A\sin (kx-\omega t)$
. We rescale it as
$y(x,t)=\varepsilon \sin (x-t)$
, by choosing
$1/k$
as the length scale,
$1/\omega$
as the time scale, and
$\omega /k$
as the velocity scale. Following previous work (Lin et al. Reference Lin, Chen and Gao2021, Reference Lin, Yu, Li and Gao2022, Reference Lin, Wang, Quan, Yu, Gao and Chen2025), we assume a small amplitude
$\varepsilon = Ak\ll 1$
.
Recalling the forced Navier–Stokes equation
the polymer stress, when neglecting the crowdedness effect (i.e.
$\beta =0$
), is written as
By ignoring the translational diffusion (i.e.
$\textit{Pe}_t^{-1}\to 0$
), the configuration tensor evolution equation is written as
For the two-dimensional problem, we adopt a stream function
$\varphi (x,y,t)$
to compute the velocity as
Correspondingly, the boundary conditions for
$\varphi (x, y, t)$
in the moving frame of the swimmer arise from conditions at infinity and on the undulatory sheet with a steady speed
$-U_{\textit{LC}}\,\hat {\boldsymbol{e}}_x$
:
After applying curl on both sides of (B1) and substituting the stream function, we have
B.1. Isotropic cases
The fluid phase exhibits isotropy where
$0\leqslant \zeta \lt \zeta _{c}$
(
$\zeta _{c}=4$
). In this regime, we close the fourth-moment tensor
${\unicode{x1D64E}}$
as
which leads to
Next, we expand all the variables in terms of
$\varepsilon$
up to second order as
Expanding the boundary conditions of (B5), the first-order problem satisfies
while the boundary conditions for the stream function at order
$\varepsilon ^{2}$
are
\begin{equation} \boldsymbol {\nabla }\varphi ^{(2)}\big |_{(x,+\infty )}=U_{\textit{LC}}^{(2)}\,\hat {\boldsymbol{e}}_{y},\quad \boldsymbol {\nabla }\varphi ^{(2)}\big |_{(x,0)}=-\sin (x-t)\,\boldsymbol{\nabla }\left .\left (\frac {\partial \varphi ^{(1)}}{\partial y}\right )\right |_{(x,0)}. \end{equation}
B.1.1. Swimming speed
Since the motion of the sheet is of order
$\varepsilon$
, the solutions are homogeneous at the zeroth order, i.e.
$\varphi ^{(0)} = 0$
,
$U_{\textit{LC}}^{(0)} = 0$
,
${\tau }_p^{(0)} = \boldsymbol{0}$
and
${{\unicode{x1D63F}}}^{(0)} = {{\unicode{x1D644}}}/2$
.
At
$O(\varepsilon )$
order, (B6) is derived as
\begin{equation} \boldsymbol{\nabla }\times \big ( \boldsymbol{\nabla }\boldsymbol{\cdot }{\tau }_p^{(1)}\big)=\dfrac {1}{\textit{Er}}\left ({{\nabla} ^{4}}\varphi ^{(1)}-\textit{Re}\,\dfrac {\partial {{\nabla} ^2}\varphi ^{(1)}}{\partial t}\right )\!. \end{equation}
Substituting (B2) into (B3), we obtain
where
$\alpha ={4}/[{\textit{Er}\,(4-\zeta )}]$
and
$\beta =1/(\textit{Er}\textit{Pe})$
. After taking the divergence and applying curl to both sides of (B16), we get
\begin{equation} \left (\alpha \frac {\partial }{\partial t}+\beta \right )\left ({{\nabla} ^{4}}\varphi ^{(1)}-\textit{Re}\,\dfrac {\partial {{\nabla} ^2}\varphi ^{(1)}}{\partial t}\right )+\dfrac {1}{4}{{\nabla} ^{4}}\varphi ^{(1)}=0. \end{equation}
By applying separation of variables (Childress Reference Childress2008), we can obtain the first-order solution of the stream function:
where
\begin{equation} F=\mathrm{Re} \left(\sqrt {\dfrac {4\beta +1-4\alpha\, \textit{Re}+(4\alpha +4\beta\, \textit{Re})\mathrm{i}}{4\beta +1+4\mathrm{i}\alpha }} \right)\!,\quad G=\mathrm{Im} \left(\sqrt {\dfrac {4\beta +1-4\alpha\, \textit{Re}+(4\alpha +4\beta\, \textit{Re})\mathrm{i}}{4\beta +1+4\mathrm{i}\alpha }} \right)\!. \end{equation}
Here, Re and Im are the operators taking the real and imaginary parts of a complex number, respectively. Then substituting
$\varphi ^{(1)}(x,y,t)$
into (B1) and (B16), we can obtain the first-order pressure
$p^{(1)}$
and polymer stress
${\tau }_p^{(1)}$
. Using stream function
$ \varphi ^{(1)}$
and boundary conditions of (B13), we derive that the first-order flow field does not exert net swimming speed, i.e.
$U_{\textit{LC}}^{(1)}=0$
.
At
$O(\varepsilon ^2)$
order, (B6) is expressed as
\begin{align} \boldsymbol{\nabla }\times \big (\boldsymbol{\nabla }\boldsymbol{\cdot }{\tau }_p^{(2)}\big)=\dfrac {1}{\textit{Er}}\left ({{\nabla} ^{4}}\varphi ^{(2)}-\textit{Re}\left (\dfrac {\partial {{\nabla} ^{2}}\varphi ^{(2)}}{\partial t}+\dfrac {\partial \varphi ^{(1)}}{\partial y}\dfrac {\partial {{\nabla} ^{2}}\varphi ^{(1)}}{\partial x}-\dfrac {\partial \varphi ^{(1)}}{\partial x}\dfrac {\partial {{\nabla} ^{2}}\varphi ^{(1)}}{\partial y}\right )\right )\!. \end{align}
After taking the divergence and applying curl to both sides of the equation, the governing equation for
$\varphi ^{(2)}$
can be obtained as
\begin{align} \begin{aligned} &\left (\alpha \dfrac {\partial }{\partial t}+\beta \right )\left [{{\nabla} ^{4}}\varphi ^{(2)}-\textit{Re}\left (\dfrac {\partial {{\nabla} ^{2}}\varphi ^{(2)}}{\partial t}+\dfrac {\partial \varphi ^{(1)}}{\partial y}\dfrac {\partial {{\nabla} ^{2}}\varphi ^{(1)}}{\partial x}-\dfrac {\partial \varphi ^{(1)}}{\partial x}\dfrac {\partial {{\nabla} ^{2}}\varphi ^{(1)}}{\partial y}\right )\right ]+\dfrac {1}{4}{{\nabla} ^{4}}\varphi ^{(2)}\\[5pt] &= \boldsymbol{\nabla }\times \left \{ \boldsymbol{\nabla }\boldsymbol{\cdot }\left [-2{{\unicode{x1D640}}}^{\kern1pt(1)}:{{\unicode{x1D64E}}}^{(1)}-\boldsymbol{u}^{(1)}\boldsymbol{\cdot }\boldsymbol {\nabla }{{\unicode{x1D63F}}}^{(1)}+\big({{\unicode{x1D63F}}}^{(1)}\boldsymbol{\cdot } \boldsymbol{\nabla }\boldsymbol{u}^{(1)}+ \boldsymbol{\nabla }\boldsymbol{u}^{(1){\rm T}}\boldsymbol{\cdot }{{\unicode{x1D63F}}}^{(1)}\big )-\dfrac {\zeta }{1-\zeta /4}\dfrac {\partial {{\unicode{x1D643}}}}{\partial t}\right ]\right \}\!, \end{aligned} \end{align}
where
${{\unicode{x1D643}}}={{\unicode{x1D63F}}}^{(1)}\boldsymbol{\cdot }{{\unicode{x1D63F}}}^{(1)}-{{\unicode{x1D63F}}}^{(1)}:{{\unicode{x1D64E}}}^{(1)}$
. To get the steady swimming speed, it is convenient to consider its time-averaged solution
$\langle \varphi ^{(2)}\rangle$
. After performing the time-averaging on (B21), we obtain
\begin{equation} \begin{aligned} &\left (\beta +1/4\right )\frac {{{\rm d}}^{4}}{{\mathrm{d}}y^{4}}\langle \boldsymbol \varphi ^{(2)}\rangle -\beta\, \textit{Re}\left \langle \frac {\partial \boldsymbol \varphi ^{(1)}}{\partial y}\frac {\partial {\nabla} ^{2}\boldsymbol \varphi ^{(1)}}{\partial x}-\frac {\partial \boldsymbol \varphi ^{(1)}}{\partial x}\frac {\partial {\nabla} ^{2}\boldsymbol \varphi ^{(1)}}{\partial y}\right \rangle \\[3pt] &=\langle \boldsymbol{\nabla }\times \{\boldsymbol{\nabla }\boldsymbol{\cdot }[-2{{\unicode{x1D640}}}^{\kern1pt(1)}:{{\unicode{x1D64E}}}^{(1)}-\boldsymbol{u}^{(1)}\boldsymbol{\cdot }\boldsymbol{\nabla }{{\unicode{x1D63F}}}^{(1)}+({{\unicode{x1D63F}}}^{(1)}\boldsymbol{\cdot }\boldsymbol{\nabla }\boldsymbol{u}^{(1)}+\boldsymbol{\nabla }\boldsymbol{u}^{(1){\rm T}}\boldsymbol{\cdot }{{\unicode{x1D63F}}}^{(1)})]\}\rangle . \end{aligned} \end{equation}
We integrate both sides of (B22) to get
With boundary conditions given by the time-average of (B14),
we obtain the second-order swimming speed
\begin{equation} U_{\textit{LC}}^{(2)}=C_{0}=\dfrac {F}{2}-\dfrac {\dfrac {\textit{Re}\,G(2F+1)}{(F+1)^2+G^2}-\textit{Er}\textit{Pe}\,\dfrac {2G(F+1)(4-\zeta )\,\textit{Pe}-8F^2\,\textit{Pe}^2}{16\,\textit{Pe}^2+(4-\zeta )^2}}{(4+\textit{Er}\textit{Pe})F}. \end{equation}
Finally, we obtain the swimming speed up to second order:
B.1.2. Hydrodynamic forces
We analyse the hydrodynamic force acting on the sheet in the
$x$
-direction (Katz Reference Katz1974; Ives & Morozov Reference Ives and Morozov2017), averaged over time:
where the unit outer normal to the sheet is given by
$\boldsymbol{n}=\boldsymbol{e}_y-\varepsilon \cos (x-t)\,\boldsymbol{e}_x$
.
After some manipulations of (B27) with Taylor expansion, we can obtain
$O(\varepsilon )\mathrm{\,order}$
\begin{equation} \left .\left \langle \frac {\partial u_{x}^{(1)}}{\partial y}+\textit{Er}\,\tau _{p,12}^{(1)}\right \rangle \right |_{y=0}=0, \end{equation}
and
$O(\varepsilon ^2)\mathrm{\,order}$
\begin{equation} \begin{gathered} \left \langle \big [\cos (x-t)\,p^{(1)}\big ]+\left [\frac {\partial u_{x}^{(2)}}{\partial y}+\sin (x-t)\frac {\partial ^{2}u_{x}^{(1)}}{\partial y^{2}}\right .\right .-2\cos (x-t)\left .\frac {\partial u_{x}^{(1)}}{\partial x}\right ] \\{}+\textit{Er}\left [\tau _{p,12}^{(2)}+ \left .\left .\sin (x-t)\frac {\partial \tau _{p,12}^{(1)}}{\partial y}-\cos (x-t)\,\tau _{p,11}^{(1)}\right ]\right \rangle \right |_{y=0}=0. \end{gathered} \end{equation}
The first-order flow field does not exert a net force on the sheet, the same as the Newtonian flow field (Katz Reference Katz1974; Ives & Morozov Reference Ives and Morozov2017). And the second-order flow field contributes a net force
\begin{align} F_{visc,x}^{(2)}&=2\left \langle \dfrac {\partial u_{x}^{(2)}}{\partial y}+\sin (x-t)\dfrac {\partial ^{2}u_{x}^{(1)}}{\partial y^{2}}-2\cos (x-t)\left .\dfrac {\partial u_{x}^{(1)}}{\partial x}\right \rangle \right |_{y=0}\nonumber \\[3pt]&=-\frac {\beta }{\beta +1/4}\frac {(-F^2+G^2-F)(4\alpha ^2+4\beta ^2+\beta )-\alpha (2\textit{FG}+G)}{4(\alpha ^2+\beta ^2)}=-\frac {4}{\textit{Er}\textit{Pe}+4}F_{\textit{pres},x}^{(2)}, \end{align}
\begin{align} F_{\textit{poly},x}^{(2)} &=2\,\textit{Er}\left .\left \langle \tau _{p,12}^{(2)}+\sin (x-t)\dfrac {\partial \tau _{p,12}^{(1)}}{\partial y}-\cos (x-t)\,\tau _{p,11}^{(1)}\right \rangle \right |_{y=0}\nonumber \\[3pt]&=-\frac {1/4}{\beta +1/4}\frac {(-F^{2}+G^{2}-F)(4\alpha ^{2}+4\beta ^{2}+\beta )-\alpha (2\textit{FG}+G)}{4(\alpha ^{2}+\beta ^{2})}=-\frac {\textit{Er}\textit{Pe}}{\textit{Er}\textit{Pe}+4}F_{_{\textit{pres},x}}^{(2)}. \end{align}
Thus we can gain three components of hydrodynamic forces in the isotropic regime as
B.1.3. Swimming efficiency
Next, we analyse the swimming efficiency based on the average rate of work done by the sheet. The time-averaged rate of work (Lauga Reference Lauga2007; Teran et al. Reference Teran, Fauci and Shelley2010) can be calculated as
After some manipulations, we can obtain
where the leading-order term is
\begin{align} \begin{aligned} P_{\textit{LC}}^{(2)}&=2\left .\left \langle u_{x}^{(1)}\left (\left (\frac {\partial u_{x}^{(1)}}{\partial y}+\frac {\partial u_{y}^{(1)}}{\partial x}\right )+\textit{Er}\,\tau _{p,12}^{(1)}\right )+u_{y}^{(1)}\left (2\frac {\partial u_{y}^{(1)}}{\partial y}-p^{(1)}-\textit{Er}\,\tau _{p,11}^{(1)}\right )\right \rangle \right |_{y=0}\\[3pt]&=\frac {(-F^2+G^2-F)(4\alpha ^2+4\beta ^2+\beta )-\alpha (2\textit{FG}+G)}{4(\alpha ^2+\beta ^2)}. \end{aligned} \end{align}
Then the mechanical efficiency (Lauga Reference Lauga2007; Teran et al. Reference Teran, Fauci and Shelley2010) is further obtained as
B.2. Nematic cases
In the nematic regime where
$\zeta \gt \zeta _{c}$
, we focus on the solutions of parallel and perpendicular gaits. In this regime, we adopt the quadratic closure to approximate the fourth-moment tensor
${\unicode{x1D64E}}$
(Doi & Edwards Reference Doi and Edwards1988) as
Here, we consider the cases
$\zeta \gg 1$
following previous work (Lin et al. Reference Lin, Yu, Li and Gao2022, Reference Lin, Wang, Quan, Yu, Gao and Chen2025), which tends to be more accurate in deep nematic. Then we expand all the variables with
$\varepsilon$
to second order, and
$\delta =\zeta ^{-1}$
to first order, i.e.
After some manipulations of boundary conditions, the
$O(\varepsilon )$
problem satisfies
Correspondingly, the boundary conditions for the stream function at order
$\varepsilon ^{2}$
are
\begin{equation} \boldsymbol{\nabla }\varphi ^{(20)}\big |_{(x,+\infty )}=U_{\textit{LC}}^{(20)}\,\hat {\boldsymbol{e}}_{y},\quad \boldsymbol {\nabla }\varphi ^{(20)}\big |_{(x,0)}=-\sin (x-t)\,\boldsymbol{\nabla }\left .\left (\frac {\partial \varphi ^{(10)}}{\partial y}\right )\right |_{(x,0)}, \end{equation}
\begin{equation} \boldsymbol{\nabla }\varphi ^{(21)}\big |_{(x,+\infty )}=U_{\textit{LC}}^{(21)}\,\hat {\boldsymbol{e}}_{y},\quad \boldsymbol {\nabla }\varphi ^{(21)}\big |_{(x,0)}=-\sin (x-t)\, \boldsymbol{\nabla }\left .\left (\frac {\partial \varphi ^{(11)}}{\partial y}\right )\right |_{(x,0)}. \end{equation}
B.2.1. Swimming speed
For convenience, we denote in the following text
At zeroth order, the solutions are homogeneous, i.e.
$\varphi ^{(0i)} = 0$
,
$U_{\textit{LC}}^{(0i)} = 0$
,
${{\tau }_p}^{(0i)} = \boldsymbol{0}{} (i=0,1)$
. For the configuration tensor, we solve for the equilibrium solution through (B2):
\begin{equation} {{\unicode{x1D63F}}}^{(00)}=\begin{pmatrix}1&0\\0&0\end{pmatrix}\!,\quad {{\unicode{x1D63F}}}^{(01)}=\begin{pmatrix}-\dfrac {1}{2}&0\\0&\dfrac {1}{2}\end{pmatrix}\!,\quad {{\unicode{x1D63F}}}^{(02)}=\begin{pmatrix}-\dfrac {1}{4}&0\\0&\dfrac {1}{4}\end{pmatrix}\quad {(\parallel )}, \end{equation}
\begin{equation} {{\unicode{x1D63F}}}^{(00)}=\begin{pmatrix}0&0\\0&1\end{pmatrix}\!,\quad {{\unicode{x1D63F}}}^{(01)}=\begin{pmatrix}\dfrac {1}{2}&0\\0&-\dfrac {1}{2}\end{pmatrix}\!,\quad {{\unicode{x1D63F}}}^{(02)}=\begin{pmatrix}\dfrac {1}{4}&0\\0&-\dfrac {1}{4}\end{pmatrix}\quad{(\perp )}. \end{equation}
Correspondingly, (B6) can be expressed as
At
$O( \varepsilon \delta ^0)$
, we substitute
${{\unicode{x1D63F}}}^{(00)}$
into (B2) and derive configuration tensor
${{\unicode{x1D63F}}}^{(10)}$
as
which yields
${\tau }_p^{(10)}=\mathbf{0}$
. Then the governing equation for the
$O(\varepsilon ^{1}\delta ^{0})$
order stream function can be obtained from (B50) as
Using the boundary conditions in (B43) on the sheet, we obtain the solution
where
$f=\mathrm{Re} (\sqrt {1+\mathrm{i}\,\textit{Re}} ),\ g=\mathrm{Im} (\sqrt {1+\mathrm{i}\,\textit{Re}} )$
, which further leads to the solution of
$p^{(10)}$
and
${\tau }^{(10)}$
. Given by the boundary condition at infinity (B43), we can get the steady speed at
$O( \varepsilon \delta ^0)$
order, i.e.
$U_{\textit{LC}}^{(10)}=0.$
At
$O(\varepsilon \delta )$
order, substituting (B51) into (B2), we can obtain
Similarly, using (B3), we can derive the solution of
${\tau }_p^{(11)}$
, which yields the following equation for
$\varphi ^{(11)}$
:
\begin{equation} \begin{aligned}{\nabla} ^{4}\varphi ^{(11)}-\textit{Re}\left(\frac {\partial {\nabla} ^{2}\varphi ^{(11)}}{\partial t}\right)={}&\textit{Er}\textit{Pe}\bigl\{\mathrm e^{-y}\left [a_{1}\cos (x-t)+a_{2}\sin (x-t)\right ]\\&{}+\mathrm e^{-fy}\left [a_{3}\cos (x-t+gy)+a_{4}\sin (x-t+gy)\right ]\bigr\}. \end{aligned} \end{equation}
Next, using the boundary conditions in (B44), we get the mean swimming speed
$U_{\textit{LC}}^{(11) }=0$
and the solution
\begin{align} \varphi ^{(11)}=\textit{Er}\textit{Pe} \left\{\begin{matrix} y\,\mathrm e^{-y}\left [b_{1}\cos (x-t)+b_{2}\sin (x-t)\right ]+y\,\mathrm e^{-fy}\left [b_{3}\cos (x-t+gy)+b_{4}\sin (x-t+gy)\right ] \\[5pt]{}+\mathrm e^{-y}\left [b_{5}\cos (x-t)+b_{6}\sin (x-t)\right ]+\mathrm e^{-fy}\left [b_{7}\cos (x-t+gy)+b_{8}\sin (x-t+gy)\right ] \end{matrix} \right\}\!, \end{align}
where
$a_m\ (m=1,\dots ,4)$
and
$b_n\ (n=1,\dots ,8)$
are constants, detailed in Appendix C.
It is evident from the calculations presented above that the sheet does not exhibit a net velocity in the first order of
$\varepsilon$
. Hence we calculate the steady velocity to the second order, i.e.
$O(\varepsilon ^2)$
.
At
$O(\varepsilon ^2\delta ^0)$
, we substitute
${{\unicode{x1D63F}}}^{(10)}$
and
${{\unicode{x1D63F}}}^{(11)}$
into (B2)–(B3), then obtain
${{{\tau }_p}^{(20)}}=0.$
This leads to the governing equation for
$\varphi ^{(20)}$
:
\begin{equation} {{\nabla} ^{4}}\varphi ^{(20)}-\textit{Re}\left (\frac {\partial {{\nabla} ^{2}}\varphi ^{(20)}}{\partial t}+\frac {\partial \varphi ^{(10)}}{\partial y}\frac {\partial {{\nabla} ^{2}}\varphi ^{(10)}}{\partial x}-\frac {\partial \varphi ^{(10)}}{\partial x}\frac {\partial {{\nabla} ^{2}}\varphi ^{(10)}}{\partial y}\right )=0. \end{equation}
We perform time-averaging on both sides of (B57):
\begin{equation} \frac {{{\rm d}}^{4}}{{\mathrm{d}}y^{4}}\langle \varphi ^{(20)}\rangle -\textit{Re}\left \langle \frac {\partial \varphi ^{(10)}}{\partial y}\frac {\partial {{\nabla} ^{2}}\varphi ^{(10)}}{\partial x}-\frac {\partial \varphi ^{(10)}}{\partial x}\frac {\partial {{\nabla} ^{2}}\varphi ^{(10)}}{\partial y}\right \rangle =0. \end{equation}
Given the boundary conditions in (B45), we get the non-zero swimming velocity at
$O(\varepsilon ^2\delta ^0)$
:
Here, the result
$U_{\textit{LC}}^{(20)}$
covers the results in Newtonian fluid with inertia (Tuck Reference Tuck1968; Childress Reference Childress2008).
At
$O(\varepsilon ^2\delta )$
, we derive the equations via (B2)–(B3) and get
\begin{equation} \tau _{p,12}^{(21)}=\begin{cases} 2D_{12}^{(10)}D_{11}^{(12)}& {(\parallel )},\\[11pt] -2D_{12}^{(10)}D_{11}^{(12)}&(\perp ). \end{cases} \end{equation}
Again, we use the same manipulation in (B50) and take the time average to get the equation of
$\varphi ^{(21)}$
:
\begin{align} \begin{aligned} \frac {{{\rm d}}^{4}}{{\mathrm{d}}y^{4}}\langle \varphi ^{(21)}\rangle &-\textit{Re}\left \langle \frac {\partial \varphi ^{(10)}}{\partial y}\frac {\partial {\nabla} ^{2}\varphi ^{(11)}}{\partial x}+\frac {\partial \varphi ^{(11)}}{\partial y}\frac {\partial {\nabla} ^{2}\varphi ^{(10)}}{\partial x}-\frac {\partial \varphi ^{(10)}}{\partial x}\frac {\partial {\nabla} ^{2}\varphi ^{(11)}}{\partial y}-\frac {\partial \varphi ^{(11)}}{\partial x}\frac {\partial {\nabla} ^{2}\varphi ^{(10)}}{\partial y}\right \rangle \\&=-\textit{Er}\frac {\partial ^{2}}{\partial y^{2}} \big\langle \tau ^{(21)}_{p,12} \big\rangle . \end{aligned} \end{align}
We integrate on both sides of (B61) and use the boundary conditions in (B46) to get
\begin{equation}U_{\textit{LC}}^{(21)}=\begin{cases} \dfrac {\textit{Er}\textit{Pe}\,(2f^3+2f^2+f-1)}{4f^3}&(\parallel ),\\[11pt] -\dfrac {\textit{Er}\textit{Pe}\,(f-1)(2f^2-1)}{4f^3}&(\perp ). \end{cases} \end{equation}
Finally, we get the swimming speed in the nematic regime as
B.2.2. Hydrodynamic forces
In the next step, we analyse the hydrodynamic force acting on the sheet in the
$x$
-direction, averaged over time:
where the unit outer normal to the sheet is given by
$\boldsymbol{n}=\boldsymbol{e}_y-\varepsilon \cos (x-t)\,\boldsymbol{e}_x$
.
After some manipulations of (B64) with Taylor expansion, we can obtain
$O(\varepsilon \delta ^i)\mathrm{\,order}$
\begin{equation} \left .\left \langle \frac {\partial u_{x}^{(1i)}}{\partial y}+\textit{Er}\,\tau _{p,12}^{(1i)}\right \rangle \right |_{y=0}=0, \end{equation}
and
$O(\varepsilon ^2\delta ^i)\mathrm{\,order}$
\begin{equation} \begin{gathered} \left \langle \big [\cos (x-t)\,p^{(1i)}\big ]+\left [\frac {\partial u_{x}^{(2i)}}{\partial y}+\sin (x-t)\frac {\partial ^{2}u_{x}^{(1i)}}{\partial y^{2}}\right .\right .-2\cos (x-t)\left .\frac {\partial u_{x}^{(1i)}}{\partial x}\right ] \\[5pt] {}+\textit{Er}\left [\tau _{p,12}^{(2i)}+ \left .\left .\sin (x-t)\frac {\partial \tau _{p,12}^{(1i)}}{\partial y}-\cos (x-t)\,\tau _{p,11}^{(1i)}\right ]\right \rangle \right |_{y=0}=0. \end{gathered} \end{equation}
Similarly, the net force occurs at the second-order flow field, i.e.
$O(\varepsilon ^2)$
:
$O(\varepsilon ^2\delta ^0)\mathrm{\,order}$
\begin{align} F_{visc,x}^{(20)}&=2\left \langle \dfrac {\partial u_{x}^{(20)}}{\partial y}+\sin (x-t)\dfrac {\partial ^{2}u_{x}^{(10)}}{\partial y^{2}}-2\cos (x-t)\left .\dfrac {\partial u_{x}^{(10)}}{\partial x}\right \rangle \right |_{y=0}\nonumber \\&=f^2+f-g^2, \\[-12pt]\nonumber \end{align}
\begin{align} F_{\textit{poly},x}^{(20)} &=2\,\textit{Er}\left .\left \langle \tau _{p,12}^{(20)}+\sin (x-t)\dfrac {\partial \tau _{p,12}^{(10)}}{\partial y}-\cos (x-t)\,\tau _{p,11}^{(10)}\right \rangle \right |_{y=0}=0, \end{align}
and
$O(\varepsilon ^2\delta ^1)\mathrm{\,order}$
\begin{align} F_{\textit{pres},x}^{(21)}&=2 \big \langle \cos (x-t)\,p^{(11)}\big \rangle \big |_{y=0}\nonumber \\[5pt]&=\textit{Er}\textit{Pe}\,\left\{\begin{matrix} -b_7g^3+(-3b_8f+3b_4)g^2+(3b_7f^2+\textit{Re}\,b_8-6b_3f-2b_7)g+b_8f^3\\[5pt] -3b_4f^2+(-\textit{Re}\,b_7-2b_8)f+(-b_5+b_1+b_3)\,\textit{Re}-b_6-b_2+2b_4 \end{matrix}\right\}\!, \end{align}
\begin{align} F_{visc,x}^{(21)}={}&-2 \big\langle \cos (x-t)\,p^{(11)}\big \rangle \big |_{y=0}\nonumber \\[5pt]={}&-\textit{Er}\textit{Pe}\,\left\{\begin{matrix} -b_7g^3+(-3b_8f+3b_4)g^2+( 3b_7f^2+\textit{Re}\,b_8-6b_3f-2b_7)g\\[5pt] +b_8f^3-3b_4f^2+(-\textit{Re}\,b_7-2b_8)f+( -b_5+b_1+b_3)\,\textit{Re}-b_6-b_2+2b_4 \end{matrix}\right\}\!, \end{align}
\begin{align} F_{\textit{poly},x}^{(21)} &=2\,\textit{Er}\left .\left \langle \tau _{p,12}^{(21)}+\sin (x-t)\dfrac {\partial \tau _{p,12}^{(11)}}{\partial y}-\cos (x-t)\,\tau _{p,11}^{(11)}\right \rangle \right |_{y=0}=0. \end{align}
Then we can gain three components of hydrodynamic forces as
B.2.3. Swimming efficiency
Finally, we analyse the swimming efficiency of the sheet in the nematic regime. Here, its rate of work can be expressed as
where each component is
\begin{align} \begin{aligned} P_{\textit{LC}}^{(20)}=\,&2\left .\left \langle u_{x}^{(10)}\left (\left (\frac {\partial u_{x}^{(10)}}{\partial y}+\frac {\partial u_{y}^{(10)}}{\partial x}\right )+\textit{Er}\,\tau _{p,12}^{(10)}\right )+u_{y}^{(10)}\left (2\frac {\partial u_{y}^{(10)}}{\partial y}-p^{(10)}-\textit{Er}\,\tau _{p,11}^{(10)}\right )\right \rangle \right |_{y=0}\\=\,&-(f^2+f-g^2),\\[3pt] P_{\textit{LC}}^{(21)}=\,&2\left .\left \langle u_{x}^{(11)}\left (\left (\frac {\partial u_{x}^{(10)}}{\partial y}+\frac {\partial u_{y}^{(10)}}{\partial x}\right )+\textit{Er}\,\tau _{p,12}^{(10)}\right )+u_{y}^{(11)}\left (2\frac {\partial u_{y}^{(10)}}{\partial y}-p^{(10)}-\textit{Er}\,\tau _{p,11}^{(10)}\right )\right \rangle \right |_{y=0}\\&+2\left .\left \langle u_{x}^{(10)}\left (\left (\frac {\partial u_{x}^{(11)}}{\partial y}+\frac {\partial u_{y}^{(11)}}{\partial x}\right )+\textit{Er}\,\tau _{p,12}^{(11)}\right )+u_{y}^{(10)}\left (2\frac {\partial u_{y}^{(11)}}{\partial y}-p^{(11)}-\textit{Er}\,\tau _{p,11}^{(11)}\right )\right \rangle \right |_{y=0}\\=\,&\textit{Er}\textit{Pe}\,\Bigl\{-b_7g^3+(-3b_8f+3b_4)g^2+(3b_7f^2+\textit{Re}\,b_8-6b_3f-2b_7)g\\&+b_8f^3-3b_4f^2+(-\textit{Re}\,b_7-2b_8)f+(-b_5+b_1+b_3)\,\textit{Re}-b_6-b_2+2b_4\Bigr\}. \end{aligned}\end{align}
According to swimming speed (B63) and work power (B74), we can finally get the mechanical efficiency.
Appendix C. Variables and coefficients
In this appendix, the detailed variables and coefficients in the asymptotic analysis are offered as follows.
C.1. Isotropic cases
\begin{equation} H_0=\frac {1}{2}\left (1-\frac {\textit{Er}\textit{Pe}}{\textit{Er}\textit{Pe}+4}\frac {16\,\textit{Pe}^2}{16\,\textit{Pe}^2+\left (4-\zeta \right )^2}\right )\!, \end{equation}
\begin{equation} M_{1}=\frac {\left \{\begin{aligned}&F^{5}\alpha ^{2}+(-G\alpha \beta +4\alpha ^{2})F^{4}+2((G^{2}+3\alpha ^{2}-G\alpha \beta )F^{3}\\&{}+2((2G\,\textit{Re}\,\beta +3G^{2}+2)\alpha ^{2}-G^{3}\alpha \beta +2\beta ^{3}G\,\textit{Re})F^{2}\\&{}+((G^{4}-8G\,\textit{Re}\,\beta +6G^{2}+1)\alpha ^{2}+2\beta (-G^{3}+G)\alpha - 8\beta ^{3}G\,\textit{Re})F\\&-G(2(-2G^{2}\,\textit{Re}\,\beta -G^{3}+2\,\textit{Re}\,\beta -G)\alpha ^{2}+\beta (G^{4}-1)\alpha -4G^{2}\,\textit{Re}\,\beta ^{3}+4\,\textit{Re}\,\beta ^{3})&\end{aligned}\right \}}{((F^{2}+G^{2}+1)^{2}-4F^{2})(8\beta +2)(\alpha ^{2}+\beta ^{2})}, \end{equation}
\begin{equation} M_{2}=\frac {\left \{\begin{aligned} &F^{5}\alpha \beta +F^{4}G\alpha ^{2}+2((2\,\textit{Re}\,\beta +G)\alpha ^{2}+\beta (G^{2}-1)\alpha +2\,\textit{Re}\,\beta ^{3})F^{3}\\&{}+2(G^{3}\alpha ^{2}-G^{2}\alpha \beta )F^{2}+(G-1)(2(2\,\textit{Re}\,\beta +G)\alpha ^{2} \end{aligned}\right \}}{((F^{2}+G^{2}+1)^{2}-4F^{2})(8\beta +2)(\alpha ^{2}+\beta ^{2})}, \end{equation}
C.2. Nematic cases


∥
⊥
ULC/ULC,Re=0
Re
Re
Remax
ε=0.02π
ζ=2
ULC/ULC,Re=0
Re
ΔULC
Re
ε=0.02π
ζ=50
Re
ΔFLCpres
Re
t=0.5Tc
ζ=2
ζ=50
ηLC/ηLC,Re=0
Re
Re
ε=0.02π
ζ=2(50)
U(Ux,Uy)
Re≈1.0
Er=2.0
Pe=1.0
ζ=2