Effect of the shape of microswimmers and slip boundary conditions on the dynamic characteristics of near-wall microswimmers

Abstract

Although the interaction between microswimmers and walls during near-wall swimming has been extensively studied, the effect of microswimmer shapes and slip boundary conditions on the dynamic characteristics of near-wall microswimmers has received less attention. In this study, elliptical microswimmer models have been developed with various aspect ratios based on circular microswimmers. The lattice Boltzmann method has been used for the numerical simulation of the dynamic behaviour of microswimmers near walls. Under slip boundary conditions, the escape or capture of microswimmers by the walls is influenced by the swimming Reynolds number (Re_s), wall slip length (l_s) and the aspect ratio (C_ab) of a microswimmer. Changes in the C_ab value of a microswimmer considerably affect its swimming state, especially for puller-type microswimmers. The tendency of pullers to be captured by the wall increases with increasing C_ab. Moreover, changes in l_s within the slip boundary condition of a puller can induce a transition in its movement state from a wall oscillation state to a stable sliding state and eventually to a wall lock-up state, a process influenced by the C_ab value of the puller. Pusher-type microswimmers show a considerably increased tendency to escape from walls with increasing C_ab and no wall lock-up state is observed, which is opposite to the case of pullers. Pushers and pullers show an increased tendency to be captured by the wall with increasing initial swimming angle of the microswimmer. The findings of this study enhance our understanding of the swimming patterns of natural microswimmers near walls and are of substantial importance for the design of artificial microswimmers and microfluidic devices.

1. Introduction

The swimming of microorganisms in natural and biological environments is considerably influenced by confined boundaries within enclosed spaces (Li 2023). These boundaries guide microswimmer behaviour through surface–microorganism interaction. Bacteria are stably captured, which accumulate on walls, facilitating the formation of biofilms (Li et al. 2011; Dey, Saha & Chakraborty 2020). Walls also have a remarkable impact on the movement of microswimmers near boundaries. The vaginal walls of mammals guide sperms to the egg (Denissenko et al. 2012), and freely swimming male and female brown algae have a probability of encountering each other near boundaries (Kinoshita et al. 2016). Microswimmers exhibit a variety of behaviours near solid boundaries, such as various wall responses of Chlamydomonas (Li & Ardekani 2014), the reverse movement of sperm along walls (Miki & Clapham 2013), distinct motion patterns of Volvox (Drescher et al. 2011) and circular swimming in different directions by flagellated bacteria (Lauga et al. 2006; Di Leonardo et al. 2011). The advancement of microfluidic technology has aided the application of artificial microswimmers in areas such as laboratory devices (Bechinger et al. 2016), cell manipulation (Bunea & Taboryski 2020) and micromotor propulsion (Poddar, Bandopadhyay & Chakraborty 2019; Chen et al. 2021). Furthermore, the design of complex boundaries in microfluidic devices (protrusions and obstacles) has developed various artificial microswimmers (Simmchen et al. 2016).

The squirmer is a typical propulsion model for microswimmers (Lighthill 1952; Blake 1971a,b), where self-propulsion is achieved by applying a specific tangential slip velocity on its surface. This propulsion is categorized into front-propelled pullers and rear-propelled pushers and is used to describe the movement of various representative microorganisms such as Bacillus subtilis, Chlamydomonas reinhardtii, Escherichia coli and Volvox. It is extensively used to explore the hydrodynamic properties of microswimmers, including the behaviour of squirmers in vertical pipes (Guan, Lin & Nie 2022; Ying et al. 2022; Nie et al. 2023), their movement in different background flow fields (More & Ardekani 2021; Liu, Ouyang & Lin 2022; Guan & Lin 2023), swimming in non-Newtonian fluids (De Corato, Greco & Maffettone 2015; Nganguia & Pak 2018; Ouyang, Lin & Ku 2018) and their behaviour near walls (Kuron et al. 2019; Chaithanya & Thampi 2021; Ishimoto, Gaffney & Smith 2023). The diffusion (Ishikawa & Pedley 2007), nutrient absorption (Magar & Pedley 2005) and rheological properties (Uspal et al. 2015) of squirmers have also been thoroughly studied.

Previous studies have primarily focused on spherical or circular squirmers. However, ~80 % of bacteria and microorganisms are elongated, with an average aspect ratio of ~3. The motility of squirmers increases with their aspect ratio (Dusenbery 2009). The body shape of phytoplankton elongates as their volume increases (Gibson, Atkinson & Gordon 2007) in organisms such as C. reinhardtii, Paramecium and E. coli. Zöttl & Stark (2013) extended the model of spherical squirmers to elongated squirmers in Poiseuille flow. They reported that the frequency of oscillation and tumbling depended on the aspect ratio although their trajectories were similar. Kyoya et al. (2015) studied Stokes flow around ellipsoidal squirmers, noting that the collective movement of squirmers is primarily influenced by their aspect ratio. This highlights the role of the squirmer shape in collective motion. Zantop & Stark (2020) introduced the squirmer rod model and compared it with spherical squirmers to find the generation of complex fluid dynamic multipole moments by squirmer rods. Ouyang & Lin (2021) established a squirmer rod model comprising multiple spherical squirmers and studied their movement in fluid. They found that the swimming speed and the energy consumption of pullers and pushers vary with the number of spherical squirmers and their proximity. Zantop (2023) investigated the collective dynamical behaviour of neutral-, puller- and pusher-type squirmer rods and reported various states ranging from disordered to clogged states. This study indicated the role of flow field types and squirmer rod density in the transition of squirmer rod movement states. Ashtari et al. (2022) used a no-slip condition and studied the movement of elliptical squirmers in a planar channel made of flexible membranes. They found that the aspect ratio of an elliptical squirmer has a significant impact on the movement and transport of squirmers, and demonstrated how elliptical squirmers differ from the spherical squirmers in a fluid dynamics environment. Compared with spherical squirmers, elliptical squirmers move faster at specific Reynolds numbers, providing important insights for the design of artificial microswimmers.

The above studies are based on no-slip boundary conditions on the wall. The slip phenomenon at the solid–fluid interface is a complex behaviour resulting from the interaction of various physicochemical parameters such as wettability, surface roughness, impurities and dissolved gases, which is commonly observed on hydrophobic surfaces. The slip length, defined as the hypothetical distance from the solid–fluid interface to the point where the fluid velocity extrapolates to zero, is a key metric for measuring the slip phenomenon. There have been some studies on slip boundary conditions. Choi, Westin & Breuer (2003) discovered that the slip length (~30 nm) of water on hydrophobic surfaces is approximately linearly related to the shear rate. Cottin-Bizonne et al. (2005) confirmed the existence of boundary slip and proposed that a slip length of 20 nm is consistent with theories and simulations for non-wetting smooth surfaces. Joseph & Tabeling (2005) increased experimental precision and found that slip lengths are generally <100 nm. Huang et al. (2008) used molecular dynamics simulations to study the hydrodynamic slip of water on hydrophobic surfaces, finding that slip lengths ranged from nanometres to tens of nanometres, thereby resolving previous disputes regarding these measurements. In addition to nanometre-scale slip lengths, micrometre-scale slip lengths can exist under certain physical and chemical conditions. Zhu & Granick (2001) pointed out that shear-induced bubble nucleation could form a gas film, which significantly increased the slip length. Tretheway & Meinhart (2002) discovered that applying a hydrophobic monolayer coating on microchannel surfaces produced a slip length of ~1 mm. Thus, the slip at the solid–fluid interface can be modelled as a uniform partial slip at the microscale level. Currently, wall slip conditions are widely applied in research on autonomous underwater vehicles (Guo & Hou 2023) and in exploring the stability of liquid flow (Samanta 2017).

Slip boundaries do influence the swimming behaviour of microswimmers near walls. The swimming pattern of E. coli near walls forms circular trajectories, which is related to the counter-rotation of its body and flagellar bundle. Lemelle et al. (2013) explained the relation between the circular motion of E. coli and slip at the solid–fluid interface, suggesting that biopolymers and surfactants controlling bacterial motility could be a novel therapeutic strategy. Hu et al. (2015) reported that the motion of E. coli shifted from clockwise circular trajectories to linear or counterclockwise circular trajectories with an increase in the slip length. Lopez & Lauga (2014) used far-field hydrodynamic methods to simulate the swimming of helical flagellated bacteria near walls and found that the slip boundary reoriented the cells parallel to the surface and attracted them to the surface. Their research has potential applications in the separation of individual cells. Ketzetzi et al. (2020) noted that the propulsion speed of microswimmers near walls was affected by the slip length of the wall. Poddar, Bandopadhyay & Chakraborty (2020) further explored the mechanism affecting the swimming speed of microswimmers near walls by the slip length, demonstrating changes in the trajectories of puller- and pusher-type swimmers. Ghosh & Poddar (2023) studied the coupling between wall slip and external flow velocity gradients, indicating that the slip length significantly affects the rheological properties of microswimmers.

However, the aforementioned studies on the motion of microswimmers near walls under slip boundary conditions have primarily focused on circular or spherical microswimmers. The geometric shape of microswimmers affects their stable swimming states in general flow fields. Hence, it can be inferred that the geometric shape of microswimmers also influences their swimming states near walls. This paper aimed to study the motion characteristics of elliptical microswimmers near walls under slip boundary conditions. (i) Elliptical microswimmer models with different aspect ratios were constructed based on circular microswimmers. (ii) The surface wettability of walls was simulated by combining modified reflection and mirror reflection conditions, thereby reflecting slip boundary conditions. (iii) The lattice Boltzmann method (LBM) was used to numerically simulate the dynamic behaviour of microswimmers near walls. Based on this foundation, the study analysed the swimming characteristics of puller- and pusher-type microswimmers near walls with different geometric configurations. It explored the impact of key factors such as the slip length, swimming Reynolds number, aspect ratio and swimming type parameters on swimming characteristics. Additionally, the study presents the dependency of microswimmer swimming modes on the initial swimming angle.

2. Numerical method

2.1. The lattice Boltzmann method (LBM)

The LBM with a single relaxation time mode (Qian, d'Humières & Lallemand 1992) is used to simulate the swimming characteristics of microswimmers near a slip wall. The discrete lattice Boltzmann relation is shown in (2.1)

(2.1)\begin{equation}{f_i}(\boldsymbol{x} + {\boldsymbol{e}_i}\Delta t,t + \Delta t) - {f_i}(\boldsymbol{x},t) ={-} \frac{1}{\tau }[\,{f_i}(\boldsymbol{x},t) - f_i^{(0)}(\boldsymbol{x},t)],\end{equation}

where f_i(x, t) is the distribution function for the microscopic velocity, e_i is the ith direction, Δt is the time step of the simulation, τ is the relaxation time related to the fluid viscosity, $f_i^{(0)}(\boldsymbol{x},t)$ is the equilibrium distribution function expressed as (2.2)

(2.2)\begin{equation}f_i^{(0)}(\boldsymbol{x},t) = {w_i}\rho \left[ {1 + \frac{{3{\boldsymbol{e}_i}\boldsymbol{\cdot }\boldsymbol{u}}}{{{c^2}}} + \frac{{9{{({\boldsymbol{e}_i}\boldsymbol{\cdot }\boldsymbol{u})}^2}}}{{2{c^4}}} - \frac{{3{\boldsymbol{u}^2}}}{{2{c^2}}}} \right],\end{equation}

where c = Δx/Δt represents the lattice speed, Δx is the lattice spacing, the speed of sound has the relation ${c_s} = c/{(3)^{1/2}}$, lattice spacing Δx and the time step Δt are fixed at 1, weights are w ₀ = 4/9, w _1–4 = 1/9 and w _5–8 = 1/36, the density and velocity of the fluid are expressed as (2.3)

(2.3a,b)\begin{equation}\rho = \sum\limits_i {{f_i}} ,\quad \rho \boldsymbol{u} = \sum\limits_i {{f_{\boldsymbol{i}}}{\boldsymbol{e}_i}} .\end{equation}

The popular D2Q9 (i.e. nine discrete velocities in two dimensions) lattice model with nine velocities is used. The discrete velocity vectors are shown in (2.4)

(2.4)\begin{equation}{\boldsymbol{e}_i} = \begin{cases} {(0,0),}&{\textrm{for}\;i = 0,}\\ {({\pm} 1,0)c,\;(0, \pm 1)c,}&{\textrm{for}\;i = 1\;\textrm{to}\;4,}\\ {({\pm} 1, \pm 1)c,}&{\textrm{for}\;i = 5\;\textrm{to}\;8.} \end{cases}\end{equation}

The macroscopic mass and momentum, expressed as (2.5) and (2.6), respectively, in the low Mach number limit, can be recovered by applying the Chapman–Enskog expansion

(2.5)\begin{gather}\frac{{\partial \rho }}{{\partial t}} + \boldsymbol{\nabla }\boldsymbol{\cdot }(\rho \boldsymbol{u}) = 0,\end{gather}

(2.6)\begin{gather}\frac{{\partial (\rho \boldsymbol{u})}}{{\partial t}} + \boldsymbol{\nabla }\boldsymbol{\cdot }(\rho \boldsymbol{uu}) ={-} \boldsymbol{\nabla }p + \rho \nu {\nabla ^2}\boldsymbol{u}.\end{gather}

2.2. Boundary conditions

The proper treatment of a moving boundary in the LBM is a critical factor in ensuring accurate simulations. Lallemand & Luo (2003) proposed an improved version of the bounce-back scheme including interpolation, which is stable, robust and easy to implement. Figure 1 shows a square and circle representing a solid and fluid node within and outside the squirmer boundary, respectively, and a triangle representing a boundary node that rebounds on the squirmer boundary. All fluid and solid nodes might change after the collision. Hence, the equilibrium distribution function from the solid node to the fluid node needs to be recalculated based on the accurate position of the boundary node (triangle presented in figure 1). The distribution function ( f ₅) of the fluid node (A) was originally calculated using the solid node (B), which was inside the curved boundary. Thus, f ₅ was not determined. A rebound scheme can be used to satisfy the no-slip boundary condition by setting the f ₅ = f ₇ condition. However, this method assumes that the boundary node (F) is exactly located at the midpoint of the line connecting two lattice nodes and is not universal. Therefore, f ₅ needs to be recalculated, and a parameter q = |AF|/|AB| is introduced to determine the position of the boundary node. Finally, f ₅ is recalculated based on nearby nodes (B, C and D) using an interpolation method using (2.7)

(2.7) \begin{align}{f_5}(A) = \begin{cases} {q(1 + 2q){f_7}(B) + (1 - 4{q^2}){f_7}(A) - q(1 - 2q){f_7}(C) - 2{w_5}\rho \dfrac{{{\boldsymbol{e}_5}\boldsymbol{\cdot }{\boldsymbol{u}_F}}}{{c_s^2}}}&{q < \dfrac{1}{2},}\\ {\dfrac{1}{{q(2q + 1)}}{f_7}(B) + \dfrac{{2q - 1}}{q}{f_5}(C) - \dfrac{{2q - 1}}{{2q + 1}}{f_5}(D) - \dfrac{{2{w_5}}}{{q(2q + 1)}}\dfrac{{{\boldsymbol{e}_5}\boldsymbol{\cdot }{\boldsymbol{u}_F}}}{{c_s^2}}}&{q \ge \dfrac{1}{2}.} \end{cases}\end{align}

The forces and torques applied to the body of microswimmers are calculated using momentum exchange methods via scanning every fluid–solid boundary link (Lallemand & Luo 2003). We used the method proposed by Aidun, Lu & Ding (1998) during the motion of the microswimmer involving the coverage and creation of fluid nodes to calculate these additional forces and torques (including the repulsive force in (2.11)). The total force and torque are obtained by integrating the obtained force and torque on the squirmer surface. Therefore, the motion trajectory of the microswimmer is determined by solving Newton's equations of motion based on the resultant forces and torques.

Figure 1. Schematic of the bounce-back scheme in the LBM.

2.3. Squirmer model

2.3.1. Circle squirmer

The squirmer model can be achieved by adding oscillatory fluctuations in the radial and tangential directions to the boundaries of particles using (2.8) (Blake 1971a,b)

(2.8)\begin{equation}{\boldsymbol{u}^s}(\boldsymbol{\theta }) = \sum\limits_{n = 0}^\infty {{A_n}\cos (n\theta )\hat{\boldsymbol{r}}} + \sum\limits_{n = 1}^\infty {{B_n}\sin (n\theta )\hat{\boldsymbol{\theta }}} .\end{equation}

Here, $\hat{\boldsymbol{r}}$ and $\hat{\boldsymbol{\theta }}$ are radial and tangential unit vectors at one point on the surface of the squirmer, as shown in figure 2(a) and A_n and B_n are the corresponding time-dependent amplitudes, respectively. The radial component is usually ignored in problems related to the squirmer motion. Hence, a simplified squirmer model represented by a second-order truncated tangential velocity (n ≤ 2) is expressed using (2.9)

(2.9)\begin{equation}{u^s}(\theta ) = {B_1}\sin \theta + {B_2}\sin \theta \cos \theta .\end{equation}

Figure 2. Schematic diagram of squirmers moving to the right. (a) Circle squirmer and (b) ellipse squirmer.

The coefficient B ₁ determines the self-propulsion strength of the microswimmer, thus determining its steady-state swimming velocity (U_s) under Stokes flow conditions of B ₁/2. The coefficient (B ₂) influences the vorticity strength around the microswimmer and β = B ₂/B ₁ is typically introduced as a key parameter to characterize different types of microswimmers: (i) puller-type microswimmers (β > 0), such as E. coli and B. subtilis, which obtain thrust from their rear end, (ii) pusher-type microswimmers (β < 0), such as C. reinhardtii, which obtain thrust from their front end (Koch & Subramanian 2011) and (iii) neutral-type microswimmers (β = 0), which are associated with vortex-free symmetric flows (Fadda, Molina & Yamamoto 2020), such as Volvox.

2.3.2. Ellipse squirmer

Inspired by the work of Zantop & Stark (2020) who took a squirmer rod as a combination of several circular squirmers with equal diameters, in our model of elliptical squirmer, the velocity at each point (e.g. E in figure 2b) on the surface of elliptical squirmer is equal to the velocity at the point on the surface of circular squirmer that is tangent to that point, as shown in figure 2(b).

In the computation, we choose enough computed nodes (e.g. E in figure 2b) on the surface of elliptical squirmer and sufficiently small spacing between adjacent squirmers to ensure high accuracy and continuity of tangential velocity at the computed nodes on the surface of the elliptical squirmer.

To demonstrate the effectiveness of the above model, we compare the swimming speed computed using the model of elliptical squirmer with that given by analytical solution, $2Us/\{ {B_1}\ast [1 + ({C_{ab}} - 1)/({C_{ab}} + 1)]\}$, as shown in figure 6(b), obtained from the conformal theorem (mapping a circle to a ellipse) (Liu et al. 2022). The two results are consistent.

For guaranteeing the continuity of surface velocities on the elliptical squirmer, before officially starting the computation, we conducted a trial run and choose enough points on the surface of the elliptical squirmer for computation. Besides, the specific implementation method for the elliptical squirmer is as follows: according to the method outlined in § 2.3, it is possible to determine the boundary nodes and the normal vectors of the elliptical squirmer, as shown at point E and vector n in figure 1. Then, an inscribed circle, circle 1, which shares the same node E and normal vector n, can be determined, and the tangential velocity at point E on circle 1 is mapped onto the ellipse. This process is repeated until the tangential velocities at all boundary points on the ellipse are derived from the circle.

2.4. Repulsive force model

Collisions might occur when swimmers approach a wall. Thus, we adopted the repulsive force model proposed by Glowinski et al. (2001) to resolve collision issues. This model has been widely applied in many studies (Ouyang et al. 2018; Ying et al. 2022; Nie et al. 2023) and has been proven to be effective and reliable for simulating interactions between swimmers and the wall. In this model, the collision between swimmer and the wall is assumed to be smooth, i.e. precaution will be taken to avoid overlapping of the regions occupied by the swimmer and wall when a swimmer collides with the wall surface. To achieve this goal, a short-range repulsive force is given

(2.10) \begin{equation}F_{ij}^w = \begin{cases} 0&{\textrm{if}\;{d^{\prime}_{ij}} > {R_i} + {R_j} + \xi }\\ {\dfrac{{{c_{ij}}}}{{{\varepsilon_w}}}{{\left( {\dfrac{{{R_i} + {R_j} + \xi - {d^{\prime}_{ij}}}}{\xi }} \right)}^2}\dfrac{{({\boldsymbol{G}_i} - {{\boldsymbol{G^{\prime}}}_j})}}{{{d^{\prime}_{ij}}}}}&{\textrm{if}\;{d^{\prime}_{ij}} \le {R_i} + {R_j} + \xi } \end{cases}.\end{equation}

Here, R_i represents the radius of a squirmer with the centroid located at G_i, and ${\boldsymbol{G^{\prime}}_j}$ is the centroid position of a hypothetical squirmer located inside the wall, as shown in figure 3. Also, ${d^{\prime}_{ij}} = |{\boldsymbol{G}_i} - {\boldsymbol{G}^{\prime}_j}|$ is the distance between two centroids, ξ is a cutoff distance below which the repulsive force would come into effect, which is set to 1.5 Δx. The reason for setting 1.5Δx as the cutoff distance instead of other values is that it can achieve the required calculation accuracy (referring to 3.4) while avoiding the need for too much calculation time. The parameter ε_w is a given positive stiffness coefficient set to 10⁻³ and c_ij is a force scaling factor and is defined as the difference between the gravitational and buoyancy forces of the squirmer in the simulations.

Figure 3. Schematic diagram of the repulsive force model. (a) Circle and (b) ellipse.

The original repulsive force in (2.10) was designed for circular particles. Therefore, we modified this formula to accommodate collisions involving elliptical particles. The modified expression of the repulsive force is given as (2.11). The label N_s presented in figure 3(b) represents the point on the elliptical particle closest to the wall, while N_w (figure 3b) is the position on the wall corresponding to that point. Thus, |N_s−N_w| represents the minimum distance from the elliptical particle to the wall. Unlike circular particles, the repulsive force from the wall to the elliptical particle might not necessarily point toward its centroid. This repulsive force generates a torque on the elliptical particle. Hence, this torque must be considered in calculations as expressed in (2.11)

(2.11)\begin{equation}F_{ij}^{we} = \begin{cases} 0&{\textrm{if}\;|{\boldsymbol{N}_s} - {\boldsymbol{N}_w}|> \xi }\\ {\dfrac{{{c_{ij}}}}{{{\varepsilon_w}}}{{\left( {\dfrac{{|{\boldsymbol{N}_s} - {\boldsymbol{N}_w}|- \xi }}{\xi }} \right)}^2}\dfrac{{({\boldsymbol{N}_s} - {\boldsymbol{N}_w})}}{{|{\boldsymbol{N}_s} - {\boldsymbol{N}_w}|}}}&{\textrm{if}\;|{\boldsymbol{N}_s} - {\boldsymbol{N}_w}|\le \xi } \end{cases}.\end{equation}

2.5. Slip boundary conditions

For the slip boundary conditions of the wall we combined the methods of reflection and mirror reflection conditions (Wang et al. 2018). Figure 4 illustrates the application of this method using the D2Q9 model; j presented in figure 4 represents the jth layer of grids and s (figure 4) is the force driving the fluid flow along the x-axis with the y-axis in the vertical direction. The distribution functions f ₀, f ₁, f ₃, f ₄, f ₇ and f ₈ for nodes can be determined during the flow process, but f ₂, f ₅ and f ₆ remain unknown. These unknown distribution functions are defined as (2.12)–(2.14) (Succi 2002)

(2.12)\begin{gather}{f_2} = {f_4},\end{gather}

(2.13)\begin{gather}{f_5} = {r_1}\,{f_7} + (1 - {r_1}){f_8} + 2{r_1}\rho {\omega _i} \cdot {u_w}/2c_s^2,\end{gather}

(2.14)\begin{gather}{f_6} = {r_1}\,{f_8} + (1 - {r_1}){f_7} + 2{r_1}\rho {\omega _i} \cdot {u_w}/2c_s^2,\end{gather}

where u_w represents the velocity of the wall and r ₁ is the adjustment coefficient . The value of r ₁ is calculated using (Wang et al. 2018)

(2.15)\begin{equation}{r_1} = \frac{1}{{1 + \dfrac{{{l_s}}}{\tau }}},\end{equation}

where l_s is the slip length, τ is the is the relaxation time related to the fluid viscosity and r ₁ is a crucial parameter linking the slip length l_s with the boundary slip effect. Using (2.15), l_s in the code can be adjusted to match the experimental result.

Figure 4. Boundary grid for the slip boundary conditions using the D2Q9 model.

3. Flow parameters and validation of method

3.1. Flow field parameters

Figure 5 presents a schematic of a microswimmer approaching a wall with slip boundary conditions. The numerical simulations have the following fluid domains: length (L), width (H), density (ρ) and viscosity (ν). The density of the microswimmer is denoted as ρ_s. The diameter of the circular microswimmer is d, while the major and minor axes of the elliptical microswimmer are 2a and 2b (2b = d), respectively, with an aspect ratio C_ab = a/b. The angle between the swimming direction of the microswimmer and the x-axis is represented by θ. In order to simulate an infinitely long area in the x-direction, a method of moving computational domain is used, i.e. both the flow and the position of the microswimmer shift to the left by one lattice spacing when the microswimmer moves beyond a distance of X ₀ + Δx from its initial position. Using this method can ensure that microswimmers do not escape the computational domain, while significantly reducing the required number of computational grids and improving resolution. A boundary condition of fully developed non-equilibrium extrapolation is applied to the upper boundary. Using this boundary condition can ensure that the upper boundary does not affect the flow near the wall surface.

Figure 5. Schematic diagram of squirmer towards a slip boundary.

The swimming Reynolds number (Re_s) is introduced to describe the self-propulsion strength of the microswimmer, as defined in (2.15), where U_s = B ₁/2. Variations in the self-propulsion strength of the microswimmer can be achieved by changing the value of Re_s. Unless specifically stated, the values of some parameters in numerical simulations are as mentioned: L = 12d, H = 8d, ρ = 1, ρ_s = 1, Re_s ranges from 0.01 to 1 and d = 40. The initial position of the microswimmer is set to X ₀ = 0.5L and Y ₀ = h = 0.4H with an initial swimming angle of θ ₀ = −0.1${\rm \pi}$

(3.1)\begin{equation}R{e_s} = \frac{{{U_s}d}}{\nu },\end{equation}

3.2. Squirmer model validation

The squirmer model presented in § 2.3 (2.9) was validated via numerical simulations of a single squirmer freely swimming in a 25d × 25d square region with all the boundaries of the region set as periodic boundary conditions. The simulated squirmers were both circular and elliptical with an aspect ratio of C_ab = 2.0. Their swimming speeds were normalized using U_s/(B ₁/2) for the circular squirmer and U_s/[(B ₁/2) × (1 + (C_ab−1)/(C_ab + 1))] for the elliptical squirmer. Figure 6 shows the results of the swimming speeds of neutral (β = 0), puller (β = 5) and pusher (β = −5) squirmers tending toward their theoretical speeds. It is evident that both circular and elliptical swimming speeds eventually converge to the theoretical speeds. This indicates that the circular and elliptical squirmer models used in this paper are credible. Terminal swimming velocity in figure 6 is based on Re_s = 0.01 rather than a larger value, the reason is that we want to confirm whether the terminal swimming velocity of microswimmers reaches the theoretical value B ₁/2 under the Stokes flow condition (extremely low Reynolds numbers) when using the squirmer model. Due to the fact that validation is carried out under specific conditions such as Stokes flow, choosing Re_s = 0.01 is typical.

Figure 6. Terminal swimming velocity for a freely swimming squirmer (Re_s = 0.01 and d = 40); (a) C_ab = 1.0 and (b) C_ab = 2.0.

We also performed a grid independence verification. The results of the swimming speeds at different mesh densities are shown in figure 7, where the diameter d of the squirmer is 5, 10, 20, 40 and 80. When d ≥ 20, the stable swimming speed of the squirmer reaches the theoretical speed, and the speed curve is smooth and stable. Considering factors such as accuracy and computational time, we chose d = 40 in the simulation.

Figure 7. Swimming speed for a freely swimming squirmer under different mesh divisions.

3.3. Slip boundary validation

This section validates the slip boundary conditions. Figure 8(a) shows the two-dimensional shear flow field features of a top plate moving at a uniform speed of U ₀ = 0.04, while the bottom plate remains stationary under slip boundary conditions. Periodic boundary conditions are applied at the inlet and outlet; N_x and N_y represent the number of grids along the X and Y directions, respectively. Figure 8(b) shows the computational results for different grid numbers (N_x × N_y) under slip boundary conditions. The results from different grid numbers are almost identical, providing a basis for the selection of grid numbers in the simulations. Additionally, this demonstrates that the calculated l_s closely matches the theoretical result (l_{s 0}).

Figure 8. Validation of slip boundary conditions. (a) Calculation area and (b) grid independence verification.

3.4. Cutoff distance validation

In § 2.4, we mentioned the cutoff distance, ξ. To determine ξ, we conducted more detailed numerical simulations for ξ = 1.2Δx, 1.5Δx, 2.0Δx and 2.5Δx, and the trajectories of the squirmer at different ξ are shown in figure 9, where the impact of ξ on the trajectories of the squirmer is very small, so we chose ξ = 1.5Δx because it can achieve the required calculation accuracy while avoiding the need for too much calculation time.

Figure 9. Trajectories of the squirmer at different cutoff distances ξ(Re_s = 0.08, β = 3, l_s = 3).

4. Results and discussion

4.1. Movement of squirmers with different C_ab values

Figure 10 displays the instantaneous streamlines around squirmers with different aspect ratios. The puller (β = 1) propels by pulling the fluid from the front and experiences restricted movement near the wall, thereby affecting the structure of streamlines. The pusher (β = −1) propels by pushing the fluid from the rear and has a different streamline distribution near the wall than that of the puller. Thus, different propulsion methods result in distinct streamlined structures around pullers and pushers. The distribution of streamlines around neutral (β = 0) squirmers resembles passive particles. A comparison of figures 10(a)–10(c) shows that the surrounding streamlines and velocity distributions vary for squirmers with different C_ab values. A recirculation area at the head of the puller increases with increasing aspect ratio C_ab when β = 1. Conversely, the recirculation area at the tail of the pusher gradually decreases when β = −1. This reflects the differences in the self-propulsion methods of squirmers and the degree of interaction between squirmers with different C_ab values and the wall.

Figure 10. Instantaneous streamlines and the velocity field around a squirmer (Re_s = 0.08, θ ₀ = −0.1${\rm \pi}$, h = 0.2H); (a) C_ab = 0.8, (b) C_ab = 1.0 and (c) C_ab = 1.2.

4.2. The interaction between a puller and the wall

Figure 11 presents the phase diagram of the final states of a puller, revealing four distinct motion states: escape, oscillation, sliding along the wall and wall lock-up. The escape state occurs when the self-propulsion strength (β) is relatively low. The puller cannot escape from the wall once β exceeds a certain critical value and is found in one of the other three states. The puller exhibits a transitional change near the wall as the β value changes, oscillating to sliding along the wall and finally locking up. Figure 12 shows the trajectories of the puller in these four states, illustrating that, apart from escape, the other three states occur near the wall.

Figure 11. Phase diagram of the final states of a puller with different C_ab values, solid circle: wall escape, solid square: oscillation near the wall, hollow square: sliding along the wall and solid triangle: wall lock-up (the region from ${\rm \pi}$/2 to 7${\rm \pi}$/6: Re_s = 0.08, the region from 7${\rm \pi}$/6 to 11${\rm \pi}$/6: Re_s = 0.8, the region from 11${\rm \pi}$/6 to ${\rm \pi}$/2: Re_s = 0.4, the radii of the circle, from innermost to outermost, represent l_s values of 0, 0.1, 0.5, 1, 5, 10, respectively); (a) C_ab = 0.8, (b) C_ab = 1.0, (c) C_ab = 1.2 and (d) C_ab = 2.0.

Figure 12. Movement trajectories of a puller in the four different states.

4.2.1. Wall-escape state of the puller

Figure 13(a) shows the variation of the horizontal velocity (V_x) of a puller during swimming as a function of the relative distance h/H from the wall under different l_s. The puller initially moves toward the wall from a distance of 0.4H, moves away from the wall under its influence and returns to a position of 0.4H away from the wall. The value of V_x gradually increases during the approach toward the wall; V_x first decreases and then increases on moving toward and away from the wall, respectively. The trend in the change of V_x remains the same as l_s increases, but the puller moves close to the wall (small h/H values). There is a significant velocity gradient in the narrow gap between the puller and the wall when the puller approaches the wall. A large l_s weakens this velocity gradient. This conclusion is confirmed by the instantaneous streamline and vorticity diagrams, as shown in figure 14. Vorticity between the puller and wall with a l_s of 10 is smaller than that without a slip length (l_s = 0). In addition, horizontal velocity of the puller increases with increasing aspect ratio C_ab. The reason is that, on the one hand, the shape of the microswimmer transitions from circular to elliptical as C_ab increases, becoming more streamlined and thereby reducing the flow resistance and allowing for an increase in V_x. On the other hand, the tangential velocity and its related V_x on the surface of an elliptical squirmer increases as C_ab increases. At l_s = 10, there's a noticeably distinct behaviour around h/H = 0.1 compared with other l_s values. On the one hand, the closer to the wall, the greater the effect of l_s, h/H = 0.1 is closest to the wall compared with other positions, so there is a noticeably distinct behaviour around h/H = 0.1. On the other hand, the larger the value of l_s, the greater the effect of reducing the velocity gradient between the wall and the microswimmer, which allows the microswimmer to swim closer to the wall at l_s = 10, thereby causing changes in speed and swimming angle and triggering different motion behaviours.

Figure 13. Variation in the horizontal velocity, torque experienced and swimming angle of a puller relative to its distance from the wall under different l_s values (θ ₀ = −0.1${\rm \pi}$, Re_s = 0.08 and β = 1). (a) Horizontal velocity of pullers, (b) torque experienced of pullers and (c) swimming angle of pullers.

Figure 14. Instantaneous streamline and vorticity diagrams around a puller at different l_s values (θ ₀ = −0.1${\rm \pi}$, Re_s = 0.08, β = 1, C_ab = 0.8, h/H = 0.2); (a) l_s = 0 and (b) l_s = 10.

Figure 13(b) shows that the torque experienced by the puller decreases with an increase in l_s at the same relative distance of h = 0.2H from the wall, inducing the puller to swim close to the wall when l_s values are high. This explains the variation in the swimming angle of the puller at different l_s values, as shown in figure 13(c). The larger the aspect ratio C_ab, the greater the torque variation during the interaction between the puller and the wall.

The escape state of a puller occurs when the values of β and Re_s are relatively low, as previously mentioned and illustrated in figure 11. Ishimoto & Gaffney (2013) reported that the puller would not generate a hydrodynamic rotation toward the wall when the value of β was not very high. However, the puller would produce additional rotation and form a movement pattern near the wall when β exceeds a certain critical value. The phase diagram presented in figure 15 elaborates on the influence of the l_s and C_ab values of the puller on its movement patterns to investigate the mechanism behind this phenomenon (Re_s = 0.08 and β ≤ 2). The grey area (figure 15) represents the escape state of the puller, while the yellow area (figure 15) indicates the state where the puller is captured by the wall. The corresponding initial angle is θ ₀ = −0.1${\rm \pi}$, as the attraction of the wall on the puller can be neglected for θ ₀ > 0 (Li & Ardekani 2014). Figure 15 shows that the puller is captured (yellow area) only when both β and l_s are relatively large for the three values of C_ab because an increased l_s changes the characteristics of the interaction between the puller and wall. The puller is in the escape state (grey area, figure 15) below the critical value (β_c) (β < β_c). The value of β_c increases with the increase of C_ab, suggesting that the pullers with a high C_ab have a high probability of escaping from the wall at Re_s = 0.08 and β ≤ 2. The reason is that the swimming speed of the microswimmer increases with increasing C_ab, while the higher speed gives the microswimmer a greater repulsive force after colliding with the wall, which makes it easier for the microswimmer to escape from the wall.

Figure 15. Phase diagram of the final states of a puller. Grey area: escape state and yellow area: wall capture state (Re_s = 0.08); (a) C_ab = 0.8, (b) C_ab = 1.0 and (c) C_ab = 1.2.

4.2.2. Periodic oscillation state of the wall capture of pullers

An increase of l_s value would prompt the puller to change from an escape state to a periodic oscillating wall capture state. Figure 16(a) shows the puller transition from an escape state to an oscillating state when l_s exceeds a certain critical value. Figure 16(b) shows the swimming angle (θ) of the puller deflecting upward as it approaches and collides with the wall, generating an upward velocity component that drives the puller away from the wall. The deflection of θ is the result of the combined effect of factors such as self-propulsion, wall effect, repulsive force and the pressure distribution around the puller. As time progresses and the distance between the puller and wall increases, the influence of the wall gradually diminishes. The propulsive force driving the puller toward the wall is dominant, reversing θ and the direction of the puller to propel toward the wall again. The limit cycle in figure 16(b) is a special pattern of system behaviour, in which the state of the system (such as the position and angle of the microswimmer) periodically returns to its initial state over time. In the pattern diagram, it appears as a closed loop.

Figure 16. Transition of puller's swimming state at different l_s (C_ab = 1.0, Re_s = 0.08 and β = 1.6). (a) Trajectory and (b) angle height.

Figure 16(a) shows that the periodic oscillatory trajectory of the puller changes with an increase of l_s. The trajectory of the puller undergoes small-scale and large-amplitude oscillations before being captured by the wall, reaching its maximum amplitude during the first departure from the wall. This maximum amplitude decreases as l_s increases. Figure 17 demonstrates that the amplitude (h_a) and frequency ( f) of the periodic oscillation of the puller change with variations in l_s. We can see that, as l_s increases, h_a decreases and f increases. The reason is that there is a velocity gradient in the narrow gap between the puller and the wall when the puller approaches the wall, the larger the value of l_s, the greater the effect of reducing the velocity gradient, which allows the puller to swim closer to the wall, thereby causing the oscillation of puller to be suppressed by the wall, i.e. decreased in amplitude h_a. While h_a decreases, the oscillation frequency f increases. The curve fitting presented in figure 17 indicates that, when l_s increases logarithmically, the changes in h_a and f are not in a straight line. From a mathematical perspective, this fitting relationship can be used to relatively easily calculate h_a and f when l_s is known.

Figure 17. Periodic oscillation state of wall capture for a puller at different l_s (C_ab = 1.0, Re_s = 0.4, β = 1).

Figure 18(a) shows the change in the periodic oscillatory trajectory of the puller with its C_ab. The trajectory of the puller undergoes small-scale large-amplitude oscillations with the maximum amplitude increasing as C_ab increases. The maximum amplitude exceeds the initial distance between the puller and wall at C_ab = 1.2. The range of large-amplitude oscillations also significantly increases as C_ab increases, which is opposite to the case of increasing l_s, as shown in figure 16(a). In figure 18(b), we display the response of the h_a value and frequency f of the puller as C_ab changes from 0.8 to 1.2. The results show that, as C_ab increases, the growth of h_a is quite noticeable, particularly when C_ab increases from 0.8 to 1.0, where h_a accelerates quickly, then the rate of increase slows down towards C_ab = 1.2. This indicates that the h_a of the puller is sensitive to increases in C_ab. In contrast, the changes in f are relatively stable, overall demonstrating insensitivity to changes in C_ab.

Figure 18. Periodic oscillation state of wall capture for a puller at different C_ab (Re_s = 0.4, β = 1, l_s = 0.1). (a) Trajectory and (b) h_a, f vs. C_ab.

4.2.3. Stable sliding state of pullers

The swimming state of the puller changes from a periodic oscillatory wall capture state to a stable wall-sliding state as l_s increases, where the puller slides parallel to the wall in a fixed configuration. Figure 19(a) shows a periodic oscillation state of the puller when l_s ≤ 0.3, forming a limit cycle in the phase diagram. Oscillation dampens and eventually becomes a stationary point in the phase diagram when l_s surpasses a critical value. This transition can be attributed to the unique pressure distribution in the flow around the puller. Figure 19(a) illustrates a large θ value of the puller as it approaches the wall with an increase of l_s, enhancing the interaction between the low-pressure area near the head of the puller and wall. The ability of the puller to reorient diminishes with the strengthening of this interaction, causing the swimming mode to change from periodic oscillation to wall sliding.

Figure 19. Swimming angle, streamline and pressure distribution of the flow field around a puller when C_ab = 0.8(Re_s = 0.08 and β = 5). (a) Swimming angle and (b) streamline motion and pressure distribution (l_s = 0.5).

Significant differences in the swimming behaviour of pullers with different C_ab values in the steady sliding state are observed. Figure 20(a) displays the steady sliding state of a puller when C_ab = 1.0. A series of magnified repeated line segments show clear periodicity. The differences in the swimming behaviour of pullers with varying C_ab values are presented in figures 19(b) and 20(b), where the streamline and pressure distributions in the steady sliding state are shown for C_ab = 0.8 and C_ab = 1.0, respectively. Although the movement modes of the pullers are similar in both cases (sliding is parallel to the wall at a fixed height), their configurations are different. The puller is inclined toward the wall with asymmetric pressure and streamline distributions when C_ab = 0.8, allowing a stable slide near the wall. Alternately, the pressure and streamline distributions are almost symmetric with slight asymmetry causing the puller to slide slowly near the wall when C_ab = 1.0.

Figure 20. Swimming angle, streamline flow and pressure distribution of the flow field around a puller when C_ab = 1.0 (Re_s = 0.08 and β = 5). (a) Swimming angle and (b) streamline and pressure distribution (l_s = 0.8).

Figure 21(a) shows the trajectory of a puller in the steady sliding state when C_ab = 0.8. Thus, the puller gradually moves closer to the wall as l_s increases. The inset in figure 21(a) indicates that the tilt angle of the puller enlarges with increasing l_s with a decreasing speed. Figure 21(b) shows that the oscillation amplitude of the tilt angle of the puller decreases as l_s increases when C_ab = 1.0. This is reflected in both the maximum and minimum values of the swimming angle of the puller decreasing and eventually stabilizing within a range of ~−90° (perpendicular to the lower wall) ± 1.5°, with an amplitude of ~3°. Thus, l_s significantly influences the swimming angle of the puller near the wall, affecting its swimming speed and subsequent movement patterns.

Figure 21. State of wall steady sliding of pullers with different C_ab (Re_s = 0.4 and β = 5.0); (a) C_ab = 0.8 and (b) C_ab = 1.0.

4.2.4. Wall lock-up state of pullers

The parameter C_ab plays a crucial role in the formation of the wall lock-up state for pullers. Figure 11 shows that the wall lock-up state occupies only a small area when C_ab = 0.8 and 1.0, this state covers more than two thirds of the area when C_ab = 1.2. The mode almost spans the entire area when C_ab = 2.0. Ishimoto (2017) indicated that the asymmetry of the geometric shape of the puller causes an additional torque owing to spatial obstruction interactions as C_ab increases. In this study, a large C_ab indicates a high torque of the puller upon collision with the wall, as shown in figure 13(b). This increased torque rapidly changes the swimming angle of the puller. The low-pressure area between the puller and the wall intensifies when |θ| exceeds a certain critical value, leading to the capture of the puller by the wall. Consequently, the puller loses its swimming speed in the X-direction and becomes firmly locked at a specific location on the wall (figure 22a). Moreover, increasing l_s brings the puller close to the wall. The closer the puller is to the wall, the greater the torque generated by the interaction with the wall, causing effects similar to increasing C_ab. Figure 22(b) shows that the maximum torque induced by the collision of the puller with the wall increases with the increase of l_s. This further confirms that l_s affects the swimming mode of the puller, especially when the torque is generated due to the interaction with the wall.

Figure 22. Characteristics of the wall lock-up state for a puller (C_ab = 1.2, Re_s = 0.4 and β = 5). (a) Streamline and pressure distribution (l_s = 5.0) and (b) maximum torque.

4.2.5. Impact of the initial swimming angle of the pullers

This section examines the final swimming state of the puller at different initial swimming angles (θ ₀). Figure 23 presents the phase diagrams of the final state of the puller for four initial swimming angle (θ ₀) values, illustrating how Re_s, β and l_s jointly affect the final state. The puller exhibits relatively uniform swimming states at different parameter settings. When β is large (β ≥ 3), the impact of θ ₀ on the final state of the puller is the minimum. The puller primarily exhibits wall oscillation and wall-sliding states when C_ab = 0.8 and 1.0, whereas the puller is predominantly in the wall lock-up state when C_ab = 1.2, 1.5 and 2.0. The final state of the puller can shift from a wall-escape state to the wall lock-up state as θ ₀ increases for small β values (β < 3). Particularly, pullers with a high C_ab and small β tend to escape from the wall at small θ ₀ values, indicating that pullers are more likely to detach from the wall at small θ ₀.

Figure 23. Phase diagrams of the final state of the puller at different θ ₀ values at Re_s = 0.08. Solid circles: wall escape, solid squares: oscillation near the wall, hollow squares: sliding along the wall and solid triangles: wall lock-up (the region from ${\rm \pi}$/10 to ${\rm \pi}$/2: C_ab = 1.0, the region from ${\rm \pi}$/2 to 9${\rm \pi}$/10: C_ab = 0.8, the region from 9${\rm \pi}$/10 to 13${\rm \pi}$/10: C_ab = 2.0, the region from 13${\rm \pi}$/10 to 17${\rm \pi}$/10: C_ab = 1.5, the region from 17${\rm \pi}$/10 to ${\rm \pi}$/10: C_ab = 1.2, the radii of the circle, from innermost to outermost, represent l_s values of 0, 0.1, 0.5, 1, 5, 10, respectively); (a) θ ₀ = −0.05${\rm \pi}$, (b) θ ₀ = −0.1${\rm \pi}$, (c) θ ₀ = −0.2${\rm \pi}$ and (d) θ ₀ = −0.3${\rm \pi}$.

4.3. Interaction between pushers and the walls

The swimming behaviour of pushers also exhibits characteristics of wall capture and wall-escape states. Unlike pullers, pushers show a tendency to escape from the wall even at a high swimming Re_s, and the likelihood of wall escape increases with the increase in Re_s. However, increasing l_s induces pushers to transform from a state of wall escape to that of wall capture when Re_s and β are the same. Figure 24 illustrates that pushers primarily exhibit three types of movement states: wall escape, oscillation and steady sliding states. Pushers do not exhibit the wall lock-up state observed in pullers, which is related to the differences in pressure and streamline distribution around the pusher and puller in the flow field. Figures 27 and 28 will further illustrate this.

Figure 24. Phase diagram of the final states of a pusher with different C_ab (solid circle: wall escape; solid square: oscillation near the wall; hollow square: sliding along the wall); (a) C_ab = 0.8, (b) C_ab = 1.0, (c) C_ab = 1.2 and (d) C_ab = 2.0.

4.3.1. Wall-escape state of pushers

Figure 25(a) shows that the V_x value of the pusher gradually increases as it approaches the wall; V_x increases with C_ab, similar to the situation of the puller depicted in figure 13(a). However, V_x first increases and then decreases when the pusher approaches the wall and begins to escape, which is opposite to the puller. This difference is attributed to the distinct torque experienced by the microswimmer near the wall and the swimming angle that determines whether it can escape, and this can be explained through the surrounding pressure distribution and streamline structure. Pullers first experience a negative torque when they approach the wall because the head of a puller is typically in a low-pressure state, thus the head is attracted by the wall, causing the swimming direction to gradually turn towards the wall. The phenomenon can be found in figure 13(c), where the swimming angle of the puller decreases initially. In contrast to pullers, pushers first experience a positive torque as they approach the wall. As shown in figure 25(c), the swimming angle of the pusher initially increases because the high-pressure region is located at the head and tail of the pusher, which leads to a tendency for the head to move away from the wall. In addition, the horizontal velocity of the pusher increases with increasing aspect ratio C_ab. Figure 25(b) shows that the negative torque the pusher experiences while swimming is significantly higher than the positive torque, a contrast to the situation of the puller shown in figure 13(b). The larger the aspect ratio C_ab, the greater the torque variation during the interaction between the pusher and the wall. Furthermore, the trend in the change of the swimming angle of the pusher shown in figure 25(c) also differs from that of the puller (figure 13c).

Figure 25. Variation in horizontal velocity, torque experienced and swimming angle of a pusher relative to its distance from the wall under different l_s (θ ₀ = −0.1${\rm \pi}$, Re_s = 0.08, and β = −1). (a) Horizontal velocity of pushers, (b) torque experienced by pushers and (c) swimming angle of pushers.

The different physical mechanisms during the swimming of pullers and pushers are shown in figure 26, which identifies four representative points in the process of approaching and escaping the wall for pullers and pushers. Figures 27 and 28 display the streamline and pressure distribution around the puller and pusher, respectively, at these four points. It can be observed that two recirculation zones appear around the head of the puller during the initial stage of approach to the wall with high pressure on the sides and low pressure around the head and tail of the puller (figure 27a). This unique pressure distribution causes the low-pressure area at the puller's head to be easily captured by the wall, resulting in what is called a ‘wall-locking state’. In contrast, the pressure and streamline distribution around pushers is the exact opposite of pullers (figure 28a). The pressure features around pushers are high pressure at the head and tail and low pressure at the sides, a distribution that does not favour the formation of a stable captured state near the wall. Therefore, pushers do not exhibit the wall-locking state observed in pullers. The presence of a recirculation zone on the side close to the wall disrupts the force balance on either side of pullers and pushers. Thus, the θ value of the front-driven pullers tends to tilt toward the wall, while the rear-driven θ of the pusher deviates away from the wall. The influence of the wall on the swimming intensifies as the puller and pusher approach the wall, enhancing the tendency of the θ value of the pusher to deviate away from the wall and altering the direction of the θ value of the puller toward the wall. When the puller and pusher begin to move away from the wall, the influence of the wall on their swimming diminishes, allowing them to leave the wall at a fixed θ (figures 27d and 28d). Additionally, the pusher requires a large swimming angle |θ| to escape from the wall owing to the presence of a low-pressure area between the side of the pusher and wall (figure 28c).

Figure 26. Four key points for pullers and pushers in the wall-escape state (points 1 and 4: h/H = 0.3, point 2: maximum T, point 3: minimum T) (C_ab = 1.0, Re_s = 0.08, β = 1.0, and l_s = 5.0). (a) Puller and (b) pusher.

Figure 27. Distribution of streamlines and pressure around the puller at four key points. (a) Point 1, (b) point 2, (c) point 3 and (d) point 4.

Figure 28. Distribution of streamlines and pressure around the pusher at four key points. (a) Point 1, (b) point 2, (c) point 3 and (d) point 4.

The θ value of the pusher continuously increases as it approaches the wall owing to the difference in propulsion methods between pullers and pushers, ensuring an angle of the pusher to move away from the wall near it. Figure 28(c) marks the forces on the pusher, revealing a balance between the high and low pressure on the tail and side, respectively. However, the angle at which the pusher moves away from the wall facilitates its escape from the influence of the wall. The phase diagram presented in figure 24 suggests that pushers have the probability of being in a wall-escape mode; l_s significantly influences the θ of the pusher as it approaches the wall.

Figure 29 illustrates the variation of the swimming angle θ for pushers with different C_ab as a function of l_s. The value of θ decreases with a decrease of l_s, and the rate of decrease of θ diminishes with increasing C_ab. This behaviour occurs because an increase in l_s weakens the velocity gradient in the narrow gap between the pusher and wall by bringing the pusher closer to the wall and altering interactions between the pusher and wall.

Figure 29. The variation of the θ of the pusher with l_s for different C_ab (h/H = 0.12); (a) Re_s = 0.08 and β = −1, (b) Re_s = 0.4 and β = −3 and (c) critical l_s.

The changing l_s changes the swimming states of the pusher when |β| is large (β = |3|). Figure 29(b) shows that a pusher with a C_ab of 0.8 has a stable sliding state along the wall and θ decreases as l_s increases. Alternately, the θ value of pushers decreases to a certain threshold value with C_ab values of 1.0 and 1.2, owing to the increase of l_s; the pusher is captured by the wall and begins to oscillate near the wall. With further decrease in θ, the oscillatory motion of the pusher eventually transforms into a stable sliding state. The critical l_s value at which the pusher is captured by the wall is dependent on C_ab during this process; as shown in figure 29(c), the critical value of l_s gradually decreases as C_ab increases.

4.3.2. Periodic oscillation state during the wall capture of the pusher

Figure 30 shows the distance from the wall and θ of pushers and pullers in the state of periodic oscillation. There are differences between the pusher and puller. The initial direction of the movement of the pusher is away from the wall, while that of the puller is toward the wall. The distance of the pusher from the wall and the amplitude of its directional oscillation are more than those of the puller with the same parameters. This indicates that the dynamic behaviour of pushers and pullers in their interaction with the wall is different, especially in terms of adjusting their swimming direction to cope with the constraints of the wall.

Figure 30. Changes in the distance from the wall and the swimming angle of microswimmers over time in the state of periodic oscillation (C_ab = 1.0, Re_s = 0.08, and l_s = 1.0). (a) Pusher (β = −3) and (b) puller (β = 3).

The h_a and f values of the trajectories of pushers and pullers with changing l_s in the periodic oscillation state are shown in figure 31. We can see that, as l_s increases, h_a decreases and f increases. Under the same parameters, the h_a value of pushers is higher than that of pullers. A large h_a allows a wide range of periodic oscillation states. The slip boundary conditions on the wall significantly impact the swimming states of microswimmers, and the extent of this impact varies between different propulsion types.

Figure 31. Variation of trajectory amplitude (h_a) and frequency ( f) with slip length (l_s) in the oscillatory state of microswimmers (C_ab = 1.0 and Re_s = 0.08). (a) Pusher (β = −3) and (b) puller (β = 3).

Figure 32 presents the swimming behaviour of pushers with three different C_ab values near the wall, revealing that all pushers exhibit an oscillatory state. However, significant differences are observed in the oscillation characteristics between the pushers with varying C_ab values. The oscillation amplitude of a pusher with C_ab = 0.9 (figure 32a) remains constant with an initial maximum amplitude (h_ima) of 31.09. Alternately, the oscillation amplitudes of the pushers with C_ab = 1.0 (figure 32b) and C_ab = 1.2 (figure 32c) gradually decay with h_ima values of 30.14 and 22.24, respectively. A comparison the f of pushers with three different C_ab values shows that f is inversely proportional to C_ab. The trajectory of the pusher with C_ab = 1.2 reaches a stable oscillation state in the shortest time, indicating that the C_ab value of microswimmers significantly affects their oscillation amplitude and frequency near walls, with a more pronounced effect in pushers.

Figure 32. Periodic oscillations of pushers with different C_ab (Re_s = 0.08, β = −3, and l_s = 0.1); (a) C_ab = 0.9, (b) C_ab = 1.0 and (c) C_ab = 1.2.

4.3.3. Stable sliding state of pushers

Figure 33 presents the streamlines and pressure distribution of pushers sliding along the wall with three C_ab values. Pushers are captured by the wall owing to the low-pressure areas between them and the wall, inducing their sliding movement along it. The hydrodynamic impact on the flow field varies with the C_ab value of pushers, resulting in distinct flow line and pressure distributions. The arrows inside the pushers (figure 33) indicate their direction of motion, showing an increased upward angle as C_ab increases. Figure 34 further illustrates the relation between the θ value of pushers and their C_ab during the wall-sliding state. The value of θ sharply increases from 4° to 30° when C_ab increases from 0.8 to 1.2. As C_ab further increases from 1.2 to 1.5, θ slightly decreases (from 30° to 27°). This increase in angle θ means that the ‘head’ or ‘front end’ of the pusher is oriented more towards directions away from the wall. Therefore, as C_ab increases, the trajectory of the pusher tends to form a larger angle with the wall, thereby increasing the likelihood of escaping from the wall. In other words, the increased swimming angle helps the pusher to be less constrained by the wall when near it, making it easier to deviate from its original path and ultimately escape from the wall.

Figure 33. Streamlines and pressure distribution of pushers with three different C_ab in the wall-sliding state (Re_s = 0.08, β = −5, and l_s = 5.0); (a) C_ab = 0.8, (b) C_ab = 1.0 and (c) C_ab = 1.2.

Figure 34. Relation between θ and C_ab of pushers in the wall-sliding state (Re_s = 0.08, β = −5, l_s = 5.0).

Figure 35 shows the relationship between θ and l_s for pushers with different C_ab during wall-sliding state. It is evident that, as l_s increases, pushers at different C_ab exhibit distinct trends in θ. When C_ab = 0.8, the θ value of the pusher is smaller and less sensitive to l_s. When C_ab = 1.2 and 1.5, the θ values are larger and decrease slightly with l_s, When C_ab = 1.0, the θ value of the pusher significantly decreases as l_s increases, indicating that the θ value of spherical pushers is more significantly affected by l_s. Since the θ of a pusher influences its behaviour and state near the wall, controlling the shape of the microswimmer can govern its behaviour and state near the wall. This provides guidance for the design of synthetic active particles in practical applications.

Figure 35. Relation between θ and l_s of pushers in the wall-sliding state (Re_s = 0.4, β = −5).

4.3.4. Impact of the initial swimming angle of the pushers

Figure 36 shows the phase diagrams of the final states of the pusher with four initial swimming angles (θ ₀), which are similar to the results of pullers presented in figure 23. Pushers with a large C_ab show noticeable differences in their final states with an increase in θ ₀. Pushers are more inclined to escape from the wall when θ ₀ is small. However, they gradually move toward a state of sliding along the wall instead of detaching from it as θ ₀ increases. This suggests that θ ₀ significantly impacts the final swimming states of pushers, especially those with a large C_ab.

Figure 36. Phase diagrams of the final state of the pusher at different θ ₀ values. Solid circles: wall escape, solid squares: oscillation near the wall, hollow squares: sliding along the wall (Re_s = 0.08); (a)θ ₀ = −0.05${\rm \pi}$, (b) θ ₀ = −0.1${\rm \pi}$, (c) θ ₀ = −0.2${\rm \pi}$ and (d) θ ₀ = −0.3${\rm \pi}$.

5. Conclusion

We used the squirmer model to describe the locomotion of microswimmers in this study. We developed elliptical squirmer models with varying aspect ratios based on the spherical squirmer and used the LBM for numerical simulations of the squirmer dynamics near walls with hydrodynamic slip. We discussed the effects of the swimming Reynolds number (Re_s), squirmer-type factor (β), wall slip length (l_s) and the aspect ratio (C_ab) of the microswimmers on their near-wall locomotion behaviour and states with the following conclusions.

1. (i) Microswimmers exhibited three fundamental locomotion states: escaping from the wall, periodic oscillation near the wall and stable sliding near the wall. The stable sliding state of pullers near the wall might change into a fourth state, wall locking.

2. (ii) The velocity gradient in the gap between the microswimmer and wall decreased when the boundary condition of the wall included slip, increasing the probability of the microswimmer being captured by the wall. This induced a possible transition of the microswimmer from an escaping state to a wall-captured state. Pullers captured by the wall exhibited a transition from a periodic oscillation state to a stable sliding state and even to a wall-locking state depending on the values of l_s and β. The frequency ( f) of the microswimmer increased with increasing l_s, but both amplitude (h_a) and swimming direction (θ) in the stable sliding state decreased with increasing l_s. Hence, the consideration of slip in the wall boundary condition and the extent of slip directly affect the locomotion behaviour and the states of microswimmers near the wall.

3. (iii) The C_ab value of the microswimmer significantly impacted the swimming speed, trajectory and final state of the microswimmer. The swimming speed of the microswimmer increased with an increase in C_ab. A large C_ab easily disrupted the symmetry of the pressure distribution in the flow field surrounding the microswimmer, generating a high torque and influencing the dynamic characteristics of the motion of the swimmer. Pullers and pushers exhibited distinct responses to changes in their C_ab. Pullers were more inclined to move away from the wall with an increase in C_ab. Conversely, an increase in the C_ab of pushers reduced the critical l_s required for them to be captured by the wall. The value of C_ab significantly affected the f of pushers but with less impact on the f of pullers. The swimming angle of pushers decreased with an increase in C_ab as microswimmers entered a wall-sliding state. Both pushers and pullers with an increased initial swimming angle were more likely to be captured by the wall.

4. (iv) The research findings revealed the diversity of microswimmer behaviours under various swimming parameter conditions. This increased our understanding of the motion behaviours and the states of natural microswimmers near confined boundaries. They provide design guidance for the fabrication of artificial active particles in practical applications.

5. (v) In the present study, we mainly focus on the effect of Re_s, l_s and C_ab on the dynamic characteristics of squirmers. The geometric scales involved in Re_s, l_s and C_ab are one-dimensional scales, and the impact of two- and three-dimensional properties on Re_s, l_s and C_ab is not so significant, so the results can provide valuable physical insights into the dynamic characteristics of squirmers. Real biological squirmers exist and move in a three-dimensional environment after all, therefore, the study of the dynamic characteristics of squirmers in a three-dimensional environment will be conducted in the future.