Introduction

The quest on how to navigate or steer to optimally reach a target is important, e.g., for airplanes to save fuel while facing complex wind patterns on their way to a remote destination, or for the coordination of the motion of the parts of a space-agent to safely land on the moon. These classical problems are well-explored and are usually solved using optimal control theory1. Likewise, navigation and search strategies are frequently encountered in a plethora of biological systems, including the foraging of animals for food2, or of T cells searching for targets to mount an immune response3. Very recently there is a growing interest also in optimal navigation problems and search strategies4,5,6,7,8,9 of microswimmers10,11,12,13 and “dry” active Brownian particles14,15,16,17,18. These active agents can self-propel in a low-Reynolds-number solvent, and might play a key role in tomorrow’s nanomedicine as recently popularized, e.g. in ref. 19. In particular, they might become useful for the targeted delivery of genes20 or drugs21,22 and other cargo23,24 to a certain target (e.g. a cancer cell) through our blood vessels, requiring them to find a good, or ideally optimal, path toward the target avoiding, e.g., obstacles and unfortunate flow field regions. In the following, we refer to the general problem regarding the optimal trajectory of a microswimmer which can freely steer but cannot control its speed toward a predefined target (point-to-point navigation) as “the optimal microswimmer navigation problem”.

The characteristic differences between the optimal microswimmer navigation problem and conventional optimal control problems for macroagents like airplanes, cruise-ships, or moon-landers root in the presence of a low-Reynolds-number solvent in the former problem only. They comprise (i) overdamped dynamics, (ii) thermal fluctuations, and (iii) long-ranged fluid-mediated hydrodynamic interactions with interfaces, walls, and obstacles, all of which are characteristic for microswimmers15. Notice in particular that the non-conservative hydrodynamic forces which microswimmers experience call for a distinct navigation strategy than the conservative gravitational forces acting, e.g. on space vehicles. Recent work has explored optimal navigation problems of dry active particles (and particles in external flow fields) accounting for (i) and partly also for (ii): Specifically, the very recent works4,5,8,9,25,26,27,28,29,30,31,32,33 have pioneered the usage of reinforcement learning34,35,36, e.g. to determine optimal steering strategies of active particles to optimally navigate toward a target position4,5,8,9 or to exploit external flow fields to avoid getting trapped in certain flow structures by learning smart gravitaxis25. Meanwhile, refs. 5,6,37 have used (deep) reinforcement learning to explore microswimmer navigation problems in mazes and obstacle arrays assuming global5 or only local6 knowledge of the environment. Very recent analytical approaches7,8 to optimal active particle navigation complement these works and allow testing machine-learned results8,9. (In addition, note that a significant knowledge exists on the complementary problem of optimizing body-shape deformation of deformable swimmers with optimal control theory; see e.g. refs. 38,39,40.) Despite this remarkable progress in recent years, (iii), and its interplay with (ii), remains an important open problem to understand the optimal microswimmer navigation strategy.

To fill this gap, in the present work, we systematically explore the optimal microswimmer navigation problem in the presence of walls or obstacles, where hydrodynamic microswimmer–wall interactions are well known to occur41,42,43,44,45,46,47,48,49, but whose impact on optimal microswimmer navigation is essentially unknown. Combining an analytical approach with numerical simulations, we find that in the presence of remote obstacles or walls, the shortest path is not fastest for microswimmers, even in the complete absence of external force or flow fields. Thus, unlike dry active particles (or light rays following Fermat’s princple), in the presence of remote obstacles microswimmers generically have to take excursions to reach their target fastest. In the presence of fluctuations, the “optimal” navigation strategy also protects microswimmers from fluctutions and can drastically decrease the traveling time as compared to cases where microswimmers head straight toward the target. This offers a promising perspective on the motion of microorganisms near surfaces or interfaces: it suggests that microorganisms might have a survival advantage when actively regulating their distance to remote walls in order to approach a food source via a strategic detour, rather than directly heading toward it. Besides their possible biological implications, our findings might provide a benchmark for future research on optimal navigation strategies of active particles.

Results and discussion

Before introducing our detailed model, let us illustrate the consequences of the finding that the shortest path is not fastest for microswimmers: Consider a microswimmer which can freely control its swimming direction (but not its speed) and aims to reach a predefined target in the presence of two obstacles (Fig. 1): While in the absence of hydrodynamics (dry active particle), the shortest path is optimal (blue), an actual microswimmer takes a qualitatively different path to reach its target fastest (red and green curves) because it produces a flow field which is reflected by the obstacles and changes its speed. In particular, “source dipole” microswimmers (specified below) which produce a flow field which slows them down near walls (Fig. 1b), take an excursion (red curve in Fig. 1a) to avoid the obstacles. Other source dipole microswimmers which produce an analogous but sign-reversed flow field (Fig. 1c), which speeds them up near walls, take a qualitatively different excursion to benefit from the proximity of both obstacles (green curve in Fig. 1a). More generally, we find that the role of hydrodynamics on optimal microswimmer routes can be subtle and lead to counterintuitive trajectory-shapes: while for source dipole swimmers the direction of the flow field plays a decisive role for the resulting trajectory-shapes, as have just seen (Fig. 1), for force-dipole swimmers, we will see that pushers and pullers lead to identical optimal trajectories.

Fig. 1: Optimal microswimmer navigation near spherical obstacles.
figure 1

a Curves represent optimal trajectories for microswimmers (red and green) and for a dry active particle (blue), i.e., in the absence of hydrodynamic interactions (HIs) with obstacles (shown in gray). The swimmers are micron sized and the swimming trajectories are assumed to take place in the plane passing through the centers of the spherical obstacles. Here, σ is the source dipole coefficient and x, z are spatial coordinates. Panels b, c show the flow field streamlines induced by a source dipole at position r = (−3, 3) μm in the presence of a spherical obstacle with radius R = 3 μm. Black arrows indicate the orientation of the swimmers. The fluid-mediated hydrodynamic interactions with the spherical boundary induce a deceleration of the swimming agent for σ > 0 leading to larger speeds as the swimmer gets away from the obstacle. In contrast to that, hydrodynamic interactions cause an acceleration for σ < 0 yielding an increased speed near the obstacle.

Full size image

Model

Let us consider a self-propelling active particle interacting with a 3D fluctuating environment. The particle’s center of mass evolves as \(\dot{{\bf{r}}}(t)=({v}_{x}({\bf{r}},\psi ,t),{v}_{y}({\bf{r}},\psi ,t),{v}_{z}({\bf{r}},\psi ,t))+\sqrt{2D}\ {\boldsymbol{\eta }}(t)\), wherein vx, vy, and vz are the components of the deterministic swimming velocity in Cartesian coordinates, which depend on the hydrodynamic swimmer–wall interactions as well as on the propulsion direction of the swimmer \((\cos \psi \cos \phi ,\cos \psi \sin \phi ,\sin \psi )\). Here, D is the diffusion coefficient which determines the strength of thermal fluctuations; for now, we choose D = 0 and consider a two-dimensional motion in the xz-plane and thus set ϕ = 0. We will discuss the effect of fluctuations later. Given a predefined initial r(t = 0) = rA and terminal point r(t = T) = rB in the xz-plane, we ask for the optimal connecting trajectory minimizing the travel time T, when the swimmer is allowed to steer freely. (This may represent, e.g., biological swimmers which steer through body shape deformations or synthetic swimmers controlled by external feedback.) That is, ψ(t) can be freely chosen so as to minimize the travel time of the swimmer. This is a well-defined optimal control problem determining the optimal trajectory, the navigation protocol {ψ(t)∣t ∈ [0, T]}, and T. It resembles classical navigation problems, e.g., of an airplane, which can steer freely and move at a speed which is determined by the wind (assuming some favorable constant engine power). Interestingly, however, while for such macroagents or dry active particles in constant external fields the shortest path is optimal7,50, for microswimmers excursions can pay off, as we will see in the following. Note that the considered model can be tested with programmed active colloids51,52,53 but should also be relevant for biological microswimmers which often swim at an almost constant speed and are able to control their self-propulsion direction on demand. These swimmers are often in contact with (remote) walls or interfaces which drastically influence their overall swimming speed and direction of motion.

As a model microswimmer, in the following, we employ the multipole description of swimming microorganisms. Accordingly, the self-generated flow field is decomposed in the far-field limit, to a good approximation, into a superposition of higher-order singularities of the Stokes equations. This model has favorably been verified experimentally for swimming Escherichia coli54 and Chlamydomonas algae55. Even though the theory is based on a far-field description of the hydrodynamic flow, it has been demonstrated using boundary integral simulations that, in some cases, the far-field approximation is surprisingly accurate, all the way down to a microswimmer–wall distance of one-tenth of a swimmer length56.

Source dipole swimmers

To develop an elementary understanding of optimal microswimmer navigation, let us first consider a source dipole microswimmer (e.g. Paramecium or active colloids with uniform surface mobility) aiming to reach a target in the presence of a distant hard wall infinitely extended in the xy-plane, for an initial and target position in the xz-plane, yielding56 \({v}_{x}=\left[{v}_{0}-\sigma /\left(4{z}^{3}\right)\right]\cos \psi\), vy = 0, and \({v}_{z}=\left({v}_{0}-\sigma /{z}^{3}\right)\sin \psi\) wherein the deviation from v0 is due to hydrodynamic swimmer–wall interactions.

To reduce the parameter space to its essential dimensions, we choose the length unit as \(\ell ={(\left|\sigma \right|/{v}_{0})}^{1/3}\), which represents the swimmer–wall distance at which the swimmer displacement per time due to hydrodynamic interactions and due to self-propulsion become comparable. For the time unit, we chose the associated time scale \(\tau ={\left|\sigma \right|}^{1/3}/{v}_{0}^{4/3}\). In reduced units, the noise-free equations of motion for a source dipole swimmer then read,

$${\dot{x}}^{* }=\left(1-\frac{s}{4{{z}^{* }}^{3}}\right)\cos \psi ,$$
(1)
$${\dot{z}}^{* }=\left(1-\frac{s}{{{z}^{* }}^{3}}\right)\sin \psi ,$$
(2)

where x* = x/ℓ, z* = z/ℓ, and s = sgn(σ). Accordingly, by expressing the equations of motion together with the underlying boundary conditions in reduced units, the optimal microswimmer trajectories only depend on the sign of the source dipole parameter, but not on its strength or on the swimmer speed. For microswimmers achieving self propulsion through surface activity (ciliated microorganisms like Paramecium, active colloidal particles with uniform surface mobility), one expects σ > 0, i.e. s = 1 whereas σ < 0 (s = −1) applies to some non-ciliated microswimmers with flagella57.

To solve the optimal navigation problem, we first eliminate ψ from the equations of motion (“Methods”). We then obtain the travel time as \({T}^{* }=\mathop{\int}\nolimits_{{x}_{A}^{* }}^{{x}_{B}^{* }}| {\dot{x}}^{* }{| }^{-1}\ {\rm{d}}{x}^{* }=\mathop{\int}\nolimits_{{x}_{A}^{* }}^{{x}_{B}^{* }}{{\mathcal{L}}}_{{\rm{SD}}}\left({x}^{* },{z}^{* }({x}^{* }),{{z}^{* }}^{\prime}({x}^{* })\right) {\rm{d}}{x}^{* }\), wherein z*(x*) represents the path which we optimize and \({\dot{z}}^{* }={z}^{* ^{\prime} }{\dot{x}}^{* }\) with \({z}^{* ^{\prime} }=\partial {z}^{* }/\partial {x}^{* }\). To find the optimal path which minimizes T*, we now determine the Lagrangian as (see “Methods” for details)

$${{\mathcal{L}}}_{{\rm{SD}}}^{* }:= {\left|{\dot{x}}^{* }\right|}^{-1}={\left(\frac{1}{{\left(1-\frac{s}{4{{z}^{* }}^{3}}\right)}^{2}}+\frac{{{z}^{* }}^{\prime}}{{\left(1-\frac{s}{{{z}^{* }}^{3}}\right)}^{2}}\right)}^{1/2}$$
(3)

and then numerically solve the Euler–Lagrange equation for \({{\mathcal{L}}}_{{\rm{SD}}}^{* }\) as a boundary value problem, with the boundary conditions \({z}^{* }({x}_{A}^{* })={z}_{A}^{* }\) and \({z}^{* }({x}_{B}^{* })={z}_{B}^{* }\), using shooting methods.

Force dipole and force quadrupole swimmers

Similarly, the translational swimming velocities due to force dipolar hydrodynamic interactions (E. coli, Chlamydomonas) with a planar hard wall reads54,56 \(\dot{x}={v}_{0}\cos \psi +3\alpha \sin (2\psi )/(8{z}^{2})\) and \(\dot{z}={v}_{0}\sin \psi +3\alpha \left[1-3\cos (2\psi )\right]/(16{z}^{2})\) where α is the force dipole coefficient. In units of \(\ell ={(| \alpha | /{v}_{0})}^{1/2}\) and \(\tau ={\left|\alpha \right|}^{1/2}/{v}_{0}^{3/2}\), the equations of motion read

$${\dot{x}}^{* }=\cos \psi +\frac{3s\sin (2\psi )}{8{{z}^{* }}^{2}},$$
(4)
$${\dot{z}}^{* }=\sin \psi +\frac{3s\left[1-3\cos (2\psi )\right]}{16{{z}^{* }}^{2}}.$$
(5)

Here, s = sgn(α) is the sign of the singularity coefficient. After some algebra, as detailed in the “Methods” section, the resulting Lagrangian follows as

$${{\mathcal{L}}}_{{\rm{FD}}}=\left|\frac{72{{z}^{* }}^{2}{{z}^{* }}^{\prime}}{{r}_{\pm }^{2}-16{{z}^{* }}^{4}-27}\right|,$$
(6)

where r± are the roots of a lengthy quartic polynomial, the coefficients of which are explicitly known functions of z* and \({{z}^{* }}^{\prime}\) (see “Methods”). The optimal swimming trajectories then result again from solving the Euler–Lagrange equation as a boundary value problem using shooting methods. Interestingly, \({{\mathcal{L}}}_{{\rm{FD}}}\) is sign invariant in the force dipole coefficient α, which means that pushers and pullers show identical optimal trajectories, albeit they require different navigation strategies for this. We will discuss this further in the “Results” section.

Finally, we also calculate the Lagrangian for force quadrupole microswimmers as detailed in the “Methods” section.

Optimal microswimmer trajectories

Flat walls

As shown in Fig. 2, in the presence of an infinitely extended and distant flat wall, we find that all considered swimmers (source and force dipole swimmers as well as force-quadrupole swimmers) follow significant excursions to reach the target fastest. That is, the shortest path is not fastest. For instance, Fig. 2c shows that source dipoles with s = 1 (which slow down when approaching the wall) follow a parabola bended away from the wall, whereas those with s = −1 prefer reducing their distance to the wall which speeds them up. The corresponding steering angles required for the optimal navigation strategy are shown in Fig. 2d. In contrast to source dipole swimmers, perhaps surprisingly, for force dipole swimmers, the shape of the resulting parabola depends only on the force dipole strength but not on the sign of the flow field (Fig. 2a, b). This pusher–puller-identity is generic not only for planar walls but also applies for spherical obstacles and obstacle landscapes, as can be directly seen from the independence of the Lagrangian of s (Eq. (6)). Interestingly, however, the required steering protocol, i.e. the temporal evolution of the optimal value of the control variable ψ(t), which the swimmer has to choose to realize the optimal path is different for pushers and for pullers (Fig. 2b). Force quadrupolar microswimmers describing small microswimmers with elongated flagella57,58, can also be solved using the Lagrangian approach; their swimming trajectories and steering angles are presented in Fig. 2e and f, respectively. The resulting parabolic curves are bent toward or away from the wall depending on the sign of the force quadrupole coefficient.

Fig. 2: Optimal 2D microswimmer trajectories near a hard and infinitely extended wall.
figure 2

Exemplary trajectories minimizing traveling times between \({{\bf{r}}}_{A}^{* }=(0,0,3)\) and \({{\bf{r}}}_{B}^{* }=(10,0,3)\) are presented together with the corresponding steering angles for force dipole swimmers (a and b), source dipole swimmers (c and d), and force quadrupole swimmers (e and f). Lengths are made dimensionless by ℓ representing the distance at which the effect of hydrodynamic interactions becomes important. The swimming trajectories take place in the xz-plane. Red solid and green dashed lines correspond to positive (s = 1) and negative (s = −1) singularity coefficient, respectively.

Full size image

Spherical obstacles and complex landscapes

Based on these results, we can now understand why hydrodynamic interactions with obstacles can have a drastic impact on the required navigation strategy to cross an obstacle field fastest. As shown in Fig. 1 without hydrodynamic interactions, the agent takes the shortest path (blue curve), whereas a source dipole microswimmer takes qualitatively different path, which depends on the sign of the singularity coefficient. This is because source dipole swimmers with s = 1 (σ > 0) tend to avoid flat walls (red solid lines in Fig. 2c) as well as spherical obstacles (red solid lines in Fig. 3) and are faster when staying at a certain distance to the obstacles. This explains that the red path, which is longer than the blue one, is faster than the blue one for s = 1-source dipole swimmers in Fig. 1. Conversely, swimmers with s = −1 (σ < 0) speed up near flat or spherical walls (green dashed lines in Figs. 2c and 3), which explains why they manage to cross the obstacle field in Fig. 1 faster when following the green trajectory than the shorter blue trajectory. These observations demonstrate that the optimal microswimmer navigation strategy qualitatively differs from the optimal strategy of dry active particles or macroagents.

Fig. 3: Zoom of optimal swimming trajectories close to a spherical boundary.
figure 3

The trajectories are directed from (0,0,13) to (5,0,11) for source dipole microswimmers near a spherical boundary of scaled radius R* = 10 positioned at the center of the system of coordinates. Here, lengths are made dimensionless by ℓ and s = ±1 is the sign of the source dipole coefficient. Physically, the singularity coefficients characterize the hydrodynamic flow induced by microswimmers in the far-field description.

Full size image

Fluctuating environments

In the world of microswimmers, fluctuations often play an important role. Besides Brownian noise which significantly displaces small biological microorganisms or active colloids on their way to a target, steering errors (or delay effects59) can effectively lead to fluctuations even in larger microswimmers. We now exemplarically consider source dipole microswimmers and set \({D}^{* }:= D/\left({\ell }^{2}/\tau \right)=D/({v}_{0}^{2/3}| \sigma {| }^{1/3}) \, \ne \, 0\) assuming that D* does not depend on space for simplicity. (Note that accounting for rotational diffusion, e.g. to represent imperfect steering, does not qualitatively change the following results.) Let us now compare the following two different navigation strategies: The first one, which we call the “straight swimming strategy” is to steer straight toward the target at each instant of time. An alternative strategy is to re-calculate the optimal path of the underlying noise-free problem at each point in time, using the present position as a starting point, and to steer in the correspondingly determined direction. We refer to this as the “optimal swimming strategy”. While the latter strategy is of course expected to be better at weak noise, for strong noise, one might expect the opposite.

However, in our simulations we find that the optimal swimming strategy notably outperforms the straight swimming strategy over the entire considered noise regime (Fig. 4a), i.e. from D* = 0 up to D* ≈ 0.15. Interestingly, the difference between the two strategies increases with the noise strength, such that the choice of the swimming strategy gets more and more important for a microswimmer as fluctuations become important. This finding might be relevant, e.g., for microswimmers when trying to reach a food source: they do much better when seeking the proximity of nearby walls first (or getting into some reasonable distance), rather than greedily heading straight toward the target.

Fig. 4: Competition of navigation strategies in a fluctuating environment.
figure 4

a Relative travel-time difference between the straight swimming strategy (travel time T1) and the optimal swimming strategy (T2) for source dipole swimmers as a function of D*. Here, s is the sign of source dipole singularity. Lengths and times are made dimensionless by ℓ and τ, respectively. In addition, D and D* denote the dimensional and dimensionless diffusion coefficients, respectively. Panels b–e: Probability distribution (averaged over 5000 trajectories) of the coarse-grained position of a microswimmer navigating from rA = (0, 0, 1.2) to rB = (1, 0, 1.2) for D* = 0.025 and s = +1 (b, c) or s = −1 (d, e). The target is counted as reached when the microswimmer enters a spherical domain centered at the target point rB of dimensionless radius of r* = 0.025. Fluctuations are considered here in 3D.

Full size image

To understand these observations, let us first consider the case s = −1 where optimal swimming tends to reduce the microswimmer–wall distance and guides the swimmer to locations where hydrodynamic interactions are comparatively important and speed up the swimmer (Fig. 4d, e). Thus, for s = −1, the swimmer can steadily approach the target for comparatively large D*-values. In contrast, when following the straight swimming strategy, nothing stops fluctuations from transferring the swimmer to regions where it is very slow (Fig. 4e). The swimmer is then dominated by noise at comparatively low D*-values and might reach the target only after following a long and winding path.

Let us now discuss the case s = 1, where optimal swimming reduces traveling times over the whole range of explored D*-values, although the above mechanism does not apply, because swimmers slow down when they are close to the wall. To see the strategic advantage of optimal swimming also here, note that when following the straight swimming strategy, fluctuations may accidentally displace the swimmer to locations close to the wall, where it is slow. In contrast, the optimal swimming strategy makes the swimmer stay away from the wall (Fig. 4b, c) and prevents it from getting trapped in regions where it is slow and dominated by noise.

Time-dependent microswimmers

We finally complement our discussion of the optimal microswimmer navigation problem by an exploration of time-dependent cases. This is inspired by microswimmers moving by body-shape deformations such as, e.g., the algae Chlamydomona reinhardtii, which alternatively moves forward (stroke) and backward (recovery stroke) and creates an oscillatory flow field60. We exemplary consider a time-dependent source dipole swimmer with \(\dot{x}=[{v}_{t}-{\sigma }_{t}/(4{z}^{3})]\cos \psi , \dot{z}= [{v}_{t}-{\sigma }_{t}/{z}^{3}]\sin \psi\), where vt ≔ v0g1(t) and σt ≔ σg2(t) are explicitly time-dependent functions. Using again length- and time-units \(\ell ={(\left|\sigma \right|/{v}_{0})}^{1/3}\), \(\tau ={\left|\sigma \right|}^{1/3}/{v}_{0}^{4/3}\), this translates to \({\dot{x}}^{* }= [{g}_{1}-s{g}_{2}/(4{{z}^{* }}^{3})]\cos \psi ,\ {\dot{z}}^{* }=[{g}_{1}-s{g}_{2}/{{z}^{* }}^{3}]\sin \psi\). While the case g1 = g2 yields the same trajectories as the time-independent case (not shown), for g1 ≠ g2 nontrivial trajectories occur (Fig. 5a, b). Necessary conditions for these trajectories can be determined based on Pontryagin’s maximum principle from optimal control theory1,61,62,63 as detailed in the “Methods” section. Choosing for illustrative purposes e.g. \({g}_{1}=1+\lambda \sin (\omega t)=1+\lambda \sin ({\omega }^{* }{t}^{* })\), where ω* = ωτ, and \({g}_{2}={g}_{1}^{2}/2\), leads to optimal trajectories (Fig. 5a, b) which feature a characteristic step-plateau-like structure. Following such a trajectory, the microswimmer mainly changes its distance from the wall in phases where it is slow, essentially to “improve” its distance from the wall for subsequent phases. When ω* increases, the plateau length decreases and for ω* → ∞ the optimal trajectory approaches a parabola (purple curve in Fig. 5a) which differs from the optimal trajectory for λ = 0, because we have 〈g2(t)〉 ≠ 〈g1(t)〉 for the time-averages of g1, g2.

Fig. 5: Time-dependent microswimmer navigation.
figure 5

Optimal trajectories and travel times of time-dependent source dipole microswimmers for fixed amplitude λ = 1 and different dimensionless frequencies ω* = ωτ shown in the key (a, c) and for fixed ω* = 6π and different λ (b, d). The shown travel times, T* = T/τ, for the optimal trajectory are always shorter than when choosing the shortest path to the target. Here, lengths and times are made dimensionless by ℓ and τ, respectively.

Full size image

The resulting travel time, monitored as a function of frequency (Fig. 5c), features a sequence of extrema occurring at frequencies where the swimmer reaches the target before completing a full driving cycle. For example, the global minimum corresponds to ω* ≈ π where the swimmer reaches the target at maximum speed without experiencing a phase where the swimmer is slower than its average speed. The travel time also depends non-monotonously on λ (Fig. 5d); it features a local minimum around λ = 0.5, where the time-average 〈g2(t)〉 is smallest, and a local maximum at λ = 1. To understand the decrease of the travel time for λ > 1, note that for λ > 1 the velocity temporarily changes sign. Since the swimmer can freely steer, it immediately turns and swims forward again with an effective speed of \({v}_{0}| 1+\lambda \sin (\omega t)|\). This leads to an average swimmer speed which increases with λ, yielding the observed decrease of T. These exact results exemplify the complexity of finding the optimal strategy in time-dependent cases and might serve as useful reference calculations to challenge corresponding machine-learning-based approaches.

Parameter regimes

Let us now briefly discuss the generic relevance of our results for typical microswimmers. The force dipole coefficient of pushers and pullers is expected to scale as α ~ a2v064, where a is the body size of the microswimmer. Thus, for v0 ~ 10 μm/s we have ℓ ≈ a. For instance, consider microswimmers which use beating elements (limbs) with an overall extension of ℓ to create a flow field which advects the body of size a with a speed v0. If the extension of the limbs is of the same order of magnitude as the body size, i.e. ℓ ~ a, we obtain \(\ell =\sqrt{| \alpha | /{v}_{0}} \sim a\) and hence α ~ a2v0. Following Figs. 2b and 3, this means that even at a wall distance of 2–3 body length, which commonly occurs in the life of many microswimmers, the deviation of the optimal path from the shortest one is significant. These estimates can be specified for E. coli bacteria where the force dipole coefficient has been measured in various experiments and amounts to α = 8 − 75 μm3/s64 yielding \(\ell =\sqrt{| \alpha | /{v}_{0}} \sim 0.5-2\)μm for v0 ~ 20 μm/s, which again is comparable to the length scale of the swimmer and means that the influence of hydrodynamic wall interactions on the required navigation strategy to reach a target fastest is highly significant. A similar discussion also applies to source dipole microswimmers, where the source dipole coefficient is expected to scale as σ ~ a3v0.

Regarding the generic relevance of our findings in the presence of noise, let us now estimate the typical value of \({D}^{* }=D/(| \sigma {| }^{1/3}{v}_{0}^{2/3})\) for source dipole swimmers to compare with Fig. 4a. Using σ ~ a3v0 as well as the Stokes–Einstein relation D = kBT/(6πηa), where η ≈ 10−3 Pa s is the viscosity of water and kBT is the thermal energy, at room temperature, we obtain D* ~ 10−2 − 10−1, depending on the size of the considered microswimmer. This roughly coincides with the parameter regime shown in Fig. 4a, which means that typical microswimmers can save a significant fraction of their traveling time to reach a food source or another target lying one or a few body lengths away from a wall, when strategically regulating their distance to the wall rather than greedily heading straight toward the target.

Conclusions

The message of this work is that, to reach their target fastest, microswimmers require navigation strategies which qualitatively differ from those used to optimize the motion of dry active particles or motile macroagents like airplanes. This finding hinges on hydrodynamic interactions between microswimmers and remote boundaries, which oblige the swimmers to take significant detours to reach their target fastest, even in the absence of external fields. Such strategic detours are particularly useful in the presence of (strong) fluctuations: they effectively protect microswimmers against fluctuations and allow them to reach a food source or another target up to twice faster than when greedily heading straight toward it. This suggests that strategically controlling their distance to remote walls might benefit the survival of motile microorganisms—which serves as an alternative to the common viewpoint, that the microswimmer–wall distance is a direct (i.e. non-actively-regulated) consequence of hydrodynamic interactions.

Our results might be relevant for future studies on microswimmers in various complex environments involving hard walls or obstacle landscapes65,66,67, penetrable boundaries68,69, or external (viscosity) gradients70,71,72. For such scenarios, our results (or generalizations based on the same framework) can be used as reference calculations, e.g., to test machine-learning-based approaches to optimal microswimmer navigation5,6 and perhaps also to help programming navigation systems for future microswimmer generations. They should also serve as a useful ingredient for future works on microswimmer navigation problems in environments which are not globally known but subsequently discovered by the microswimmers. Finally, for future work, it would also be interesting to explicitly solve the Hamilton-Jacobi-Bellman equation for the present problem with noise to compare the discussed navigation strategies which are optimal in the absence of noise and highly useful in the presence of noise, with the optimal navigation strategy following from this equation.

Methods

Here we discuss details regarding the two approaches used to solve the optimal microswimmer navigation problem based on a Lagrangian approach and on Pontryagin’s maximum principle respectively. Both approaches lead to identical results but have been found to be advantageous in different situations: the Lagrangian approach leads to a boundary value problem which is more immediate to implement, numerically simpler, and more robust than the corresponding higher-dimensional problem resulting from Pontryagin’s principle. The latter in turn allows for solving more general problems applying e.g. also to explicitly time-dependent microswimmers.

Lagrangian approach for source dipole microswimmers

To find the optimal path, we write the path connecting the starting point (\({x}_{A}^{* },{z}_{A}^{* }\)) and the terminal point (\({x}_{B}^{* },{z}_{B}^{* }\)) as a function z*(x*) and write the traveling time as \({T}^{* }=\mathop{\int}\nolimits_{{x}_{A}^{* }}^{{x}_{B}^{* }}| {\dot{x}}^{* }{| }^{-1}\ {\rm{d}}{x}^{* }=\mathop{\int}\nolimits_{{x}_{A}^{* }}^{{x}_{B}^{* }}{{\mathcal{L}}}^{* }\left({x}^{* },{z}^{* }({x}^{* }),{{z}^{* }}^{\prime}({x}^{* })\right)\ {\rm{d}}{x}^{* }\). Following the Lagrangian optimization approach, a necessary condition for minimizing T* follows from the Euler–Lagrange equation

$$\frac{{\rm{d}}}{{\rm{d}}{x}^{* }}\frac{\partial {{\mathcal{L}}}^{* }}{\partial {{z}^{* }}^{\prime}}-\frac{\partial {{\mathcal{L}}}^{* }}{\partial {z}^{* }}=0,$$
(7)

where the Lagrangian, \({{\mathcal{L}}}^{* }\), depends on the microswimmer under consideration. First, considering source-dipolar hydrodynamic interactions with a planar interface, an explicit expression for the Lagrangian can readily be obtained. It follows from Eqs. (1) and (2) that \(\cos \psi ={x}^{* }/\left(1-s/\left(4{{z}^{* }}^{3}\right)\right)\) and \(\sin \psi ={z}^{* }/\left(1-s/{{z}^{* }}^{3}\right)\). By enforcing the relation \({\cos }^{2}\psi +{\sin }^{2}\psi =1\) and using the fact that \({\dot{z}}^{* }={z}^{* ^{\prime} }{\dot{x}}^{* }\) with \({z}^{* ^{\prime} }=\partial {z}^{* }/\partial {x}^{* }\), the Lagrangian can explicitly be obtained as

$${{\mathcal{L}}}_{{\rm{SD}}}^{* }:= {\left|{\dot{x}}^{* }\right|}^{-1}={\left({\left(1-\frac{s}{4{{z}^{* }}^{3}}\right)}^{-2}+{{{z}^{* }}^{\prime}}^{2}{\left(1-\frac{s}{{{z}^{* }}^{3}}\right)}^{-2}\right)}^{\frac{1}{2}}.$$
(8)

Inserting this Langrangian into the Euler–Lagrange Eq. (7) shows that the optimal swimming trajectory is governed by the following second-order differential equation

$$A({z}^{* }){z}^{* ^{\prime\prime} }({x}^{* })+B({z}^{* }){z}^{* ^{\prime} }{({x}^{* })}^{2}+C({z}^{* })=0,$$
(9)

where we have defined the coefficients

$$A({z}^{* })={z}^{* }({x}^{* })\left({z}^{* }{({x}^{* })}^{3}-s\right){\left(4{z}^{* }{({x}^{* })}^{3}-s\right)}^{3},$$
(10)
$$B({z}^{* })=-3s\left(2{z}^{* }{({x}^{* })}^{3}+s\right){\left(4{z}^{* }{({x}^{* })}^{3}-s\right)}^{2},$$
(11)
$$C({z}^{* })=48s{\left({z}^{* }{({x}^{* })}^{3}-s\right)}^{3}.$$
(12)

Equations (9) through (12) subject to Dirichlet boundary condition of imposed vertical distance at the start and end points can readily be solved numerically using a computer algebra software by means of a standard shooting method.

Lagrangian approach for force dipole microswimmers

Next, for force-dipolar hydrodynamic interactions, we first solve Eq. (5) for the orientation angle ψ, which yields four distinct solutions. They are given by

$${\psi }_{1,2}=\arctan \left({A}_{+},\pm {B}_{+}\right),$$
(13)
$${\psi }_{3,4}=\arctan \left({A}_{-},\pm {B}_{-}\right),$$
(14)

where we have defined the arguments

$${A}_{\pm }=\frac{s}{9}\left(-4{{z}^{* }}^{2}\pm {\phi }_{1}\right),$$
(15)
$${B}_{\pm }=\frac{s}{9}{\left({\phi }_{2}\pm 8{{z}^{* }}^{2}{\phi }_{1}\right)}^{1/2},$$
(16)

with

$${\phi }_{1}={\left(72s{{z}^{* }}^{2}{\dot{z}}^{* }+27+16{{z}^{* }}^{4}\right)}^{1/2},$$
(17)
$${\phi }_{2}=2\left(-36s{{z}^{* }}^{2}{\dot{z}}^{* }+27-16{{z}^{* }}^{4}\right).$$
(18)

Note that, for \(a,b\in {\mathbb{R}}\), the function \(\arctan (b,a)\) returns the principal value of the argument of the complex number c = a + ib, i.e.73,

$$\arctan (b,a)=-{\rm{i}}\mathrm{ln}\,\left(\frac{c}{| c| }\right)\in (-\pi ,\pi ],$$
(19)

where \(| c| ={\left({a}^{2}+{b}^{2}\right)}^{1/2}\). We note that only real values of the steering angle should be considered. Now inserting Eqs. (13) and (14) into Eq. (4), setting \({\dot{z}}^{* }={{z}^{* }}^{\prime}{\dot{x}}^{* }\), and solving the resulting equations for \({\dot{x}}^{* }\), the Lagrangian is obtained as

$${{\mathcal{L}}}_{{\rm{FD}}}:= {\left|{\dot{x}}^{* }\right|}^{-1}=\left|\frac{72{{z}^{* }}^{2}{{z}^{* }}^{\prime}}{{r}_{\pm }^{2}-16{{z}^{* }}^{4}-27}\right|,$$
(20)

where r± are the roots of the quartic polynomial

$${P}_{\pm }(Z)={a}_{0}\pm {a}_{1}Z+{a}_{2}{Z}^{2}+\pm {a}_{3}{Z}^{3}+{a}_{4}{Z}^{4},$$
(21)

the coefficients of which are explicitly given by

$${a}_{0}= \, 9{\left(16{{z}^{* }}^{4}+27\right)}^{2}+256{{{z}^{* }}^{\prime}}^{2}{{z}^{* }}^{4}\left(16{{z}^{* }}^{4}-81\right),\\ {a}_{1}= \, -64{{{z}^{* }}^{\prime}}^{2}{{z}^{* }}^{2}\left(16{{z}^{* }}^{4}+81\right),\\ {a}_{2} = \, -6\left(3+2{{{z}^{* }}^{\prime}}^{2}\right)\left(16{{z}^{* }}^{4}+27\right),\\ {a}_{3}= \, 32{{{z}^{* }}^{\prime}}^{2}{{z}^{* }}^{2},\\ {a}_{4} = \, 9+4{{{z}^{* }}^{\prime}}^{2}.$$

The nature of the roots of the quartic polynomial is primarily determined by the sign of the discriminant74. Assuming that z* is a weakly-varying function about the value h > 0, such that z*(x*) = h + ϵf(x*), where ∣ϵ∣ ≪ h, the discriminants Δ± of the polynomial function given by Eq. (21) can be expanded to leading order about ϵ = 0 as

$${{{\Delta }}}_{\pm } \sim K\left(27+16{h}^{4}\right)\left(3s+4{h}^{2}\right)\left(3s-4{h}^{2}\right){{{z}^{* }}^{\prime}}^{4},$$

where K is a positive real number. In the far-field limit, we expect that h2 ≫ 3/4, and thus Δ± < 0. Accordingly, the polynomial functions has two distinct real roots and two complex conjugate non-real roots75.

If we denote by r1 and r2 the real roots of P+, then it can readily be noticed that −r1 and −r2 are the real roots of P− since P−(−Z) = P+(Z). Consequently, the system admits two possible Lagrangians, as can be inferred from Eq. (20).

The roots r1 and r2 can be obtained via computer algebra systems. They are not listed here due to their complexity and lengthiness.

Physically, the Lagrangian yielding the shortest traveling time is the one that needs to be considered7.

Lagrangian approach for force-quadrupolar microswimmers

Finally, we investigate the optimal swimming due to force-quadrupolar hydrodynamic interactions with the interface. In this case, the translational swimming velocities read56,58

$$\dot{x}={v}_{0}\cos \psi +\frac{\nu \cos \psi }{32{z}^{3}}\left[27\cos \left(2\psi \right)-13\right],$$
(22)
$$\dot{z}={v}_{0}\sin \psi +\frac{\nu \sin \psi }{8{z}^{3}}\left[9\cos \left(2\psi \right)+5\right],$$
(23)

where ν is the force quadrupolar coefficient. In units of \(\ell ={(\left|\nu \right|/{v}_{0})}^{1/3}\) and \(\tau ={\left|\nu \right|}^{1/3}/{v}_{0}^{4/3}\), Eqs. (22) and (23) can be expressed in a dimensionless form as

$${\dot{x}}^{* }=\cos \psi +\frac{s\cos \psi }{32{{z}^{* }}^{3}}\left[27\cos \left(2\psi \right)-13\right],$$
(24)
$${\dot{z}}^{* }=\sin \psi +\frac{s\sin \psi }{8{{z}^{* }}^{3}}\left[9\cos \left(2\psi \right)+5\right],$$
(25)

where we have used the abbreviation \(s=\mathrm{sgn}\,(\nu )\).

Solving Eq. (24) for ψ yields three possible distinct values:

$${\psi }_{1}=\arccos \left({\phi }_{+}\right),$$
(26)
$${\psi }_{2}=\pi -\arccos \left(\frac{1}{2}\ {\phi }_{+}-{\rm{i}}\ \frac{\sqrt{3}}{2}\ {\phi }_{-}\right),$$
(27)
$${\psi }_{3}=\pi -\arccos \left(\frac{1}{2}\ {\phi }_{+}+{\rm{i}}\ \frac{\sqrt{3}}{2}\ {\phi }_{-}\right),$$
(28)

where we have defined for convenience the abbreviations \(Z={{z}^{* }}^{3}/s\) and \(E={\left(27{\dot{x}}^{* }Z+{f}^{\frac{1}{2}}\right)}^{\frac{1}{3}}\). Moreover,

$$f=64{Z}^{3}+3\left(243{\dot{x}}^{2}-80\right){Z}^{2}+300Z-125,$$
(29)

and

$${\phi }_{\pm }=\frac{2}{9}\left(E\pm \frac{5-4Z}{E}\right).$$
(30)

Substituting the expression of ψ1 given by Eq. (26) into Eq. (25), using the fact that \({\dot{z}}^{* }={{z}^{* }}^{\prime}{\dot{x}}^{* }\), and solving the resulting equation for \({\dot{x}}^{* }\) yields the expression of the Lagrangian

$${{\mathcal{L}}}_{{\rm{FQ}}}:= {\left|{\dot{x}}^{* }\right|}^{-1}=\frac{54{{z}^{* }}^{3}}{\left|{\rho }^{3}-s{\rho }^{-3}{\left(4{{z}^{* }}^{3}-5s\right)}^{3}\right|},$$
(31)

where ρ denotes the roots of explicitly known polynomial of degree 12. In the physical range of parameters, this polynomial admits either one real radical having the order of multiplicity four or three distinct radicals also having the order of multiplicity four.

It turned out that the same set of radicals are obtained when making substitution with ψ2 of ψ3. Again, only real values of the steering angle should physically be considered.

In order to proceed further, we evaluate numerically the Lagrangians and accurately fit the results using a standard nonlinear bivariate hypothesis of the form

$${\mathcal{L}}({z}^{* },{{z}^{* }}^{\prime})=\mathop{\sum }\limits_{m=0}^{N}\mathop{\sum }\limits_{n=0}^{N}{a}_{ij}{{z}^{* }}^{m}{{{z}^{* }}^{\prime}}^{n},$$
(32)

where aij are fitting parameters. Here, we have taken N = 5 but checked that taking larger values does not alter our results.

Hydrodynamic interactions near spherical boundaries

The translational swimming velocities resulting from source dipolar hydrodynamic interactions with a solid sphere of radius R positioned at the origin of coordinates (see Fig. 6) can be decomposed into two terms

$${\bf{V}}=\hat{{\bf{e}}}+{{\bf{v}}}^{{\rm{HI}}},$$
(33)

with \(\hat{{\bf{e}}}=\cos \theta \ \hat{{\bf{n}}}+\sin \theta \ \hat{{\bf{t}}}\) denoting the instantaneous orientation angle of the swimmer, and

$${{\bf{v}}}^{{\rm{HI}}}=P\sin \theta \ \hat{{\bf{t}}}+Q\cos \theta \ \hat{{\bf{n}}},$$
(34)

quantifies the effect of hydrodynamic interactions with the spherical boundary. This contribution can readily be determined from the Green’s function near a rigid sphere76. Here, we have defined for convenience

$$P=-\frac{\sigma R\left(3{h}^{2}+6hR+8{R}^{2}\right)}{{h}^{3}{\left(h+2R\right)}^{3}},$$
(35)
$$Q=\frac{\sigma R\left(3{h}^{2}+6hR+4{R}^{2}\right)\left({h}^{2}+2hR-2{R}^{2}\right)}{4{h}^{3}{\left(h+2R\right)}^{3}{\left(h+R\right)}^{2}},$$
(36)

where, again, we have scaled lengths by a characteristic length scale of the swimmer L, and velocities by the bulk swimming speed v0.

Fig. 6: Graphical illustration of a microswimmer moving near a spherical boundary.
figure 6

Here, \({\hat{{\bf{x}}}}_{0}\) and \({\hat{{\bf{z}}}}_{0}\) are the unit vectors of the system of Cartesian coordinates, R denotes the radius of the spherical obstacle, \(\hat{{\bf{e}}}\) is the instantaneous orientation vector of the swimmer, \(\hat{{\bf{n}}}\) and \(\hat{{\bf{t}}}\) are unit vectors in the microswimmer’s frame of reference, θ is an inclination angle, and h defines the distance between the microswimmer and the obstacle.

Full size image

Without loss of generality, we consider motion in the plane y = 0.

To obtain the translational swimming velocities near two obstacles, as illustrated in Fig. 1a, we use the commonly employed superposition approximation77,78,79. The latter conjectures that the solution for the Green’s function near two widely separated obstacles can conveniently be approximated by superimposing the contributions due to each obstacle independently. Accordingly,

$${\bf{V}}={\hat{{\bf{e}}}}_{1}+{{\bf{v}}}_{1}^{{\rm{HI}}}+{{\bf{v}}}_{2}^{{\rm{HI}}},$$
(37)

where

$${{\bf{v}}}_{i}^{{\rm{HI}}}={P}_{i}\sin {\theta }_{i}\ {\hat{{\bf{t}}}}_{i}+{Q}_{i}\cos {\theta }_{i}\ {\hat{{\bf{n}}}}_{i},\quad i\in \{1,2\}.$$
(38)

In addition,

$$\cos {\theta }_{1}\ {\hat{{\bf{n}}}}_{1}+\sin {\theta }_{1}\ {\hat{{\bf{t}}}}_{1}=:{\hat{{\bf{e}}}}_{1}\equiv {\hat{{\bf{e}}}}_{2}:= \cos {\theta }_{2}\ {\hat{{\bf{n}}}}_{2}+\sin {\theta }_{2}\ {\hat{{\bf{t}}}}_{2}.$$
(39)

We then project Eq. (37) along the unit vectors \({\hat{{\bf{t}}}}_{1}\) and \({\hat{{\bf{n}}}}_{1}\) to obtain two scalar equations. Next, we project Eq. (39) along the unit vectors \({\hat{{\bf{t}}}}_{2}\) and \({\hat{{\bf{n}}}}_{2}\) to obtain two additional equations. Subsequently, the unknown quantities \(\cos {\theta }_{1}\), \(\sin {\theta }_{1}\), \(\cos {\theta }_{2}\), and \(\sin {\theta }_{2}\) can be expressed in terms of x, z, and \(z^{\prime}\) by solving the resulting linear system of four equations, upon using the fact that \(\dot{z}={z}^{\prime}\dot{x}\). By solving the identity \({\cos }^{2}{\theta }_{1}+{\sin }^{2}{\theta }_{1}=1\) for \(\dot{x}\), the Lagrangian can then be obtained explicitly as \({{\mathcal{L}}}_{{\rm{SD}}}={\left|\dot{x}\right|}^{-1}\). The expression of the resulting Lagrangian is rather lengthy and complicated and thus omitted here. Finally, the Euler–Lagrange equation can be solved numerically in Matlab using the \(\tt ode45\) routine.

Optimal control theory for microswimmer navigation

We now use Pontryagin’s principle from optimal control theory1,61,62,63 as an alterntive approach to solve the optimal navigation problem for microswimmers. While the Lagrangian approach is numerically rather convenient for the present set of problems, Pontryagin’s approach provides a more general framework to determine optimal microswimmer trajectories. In particular, it allows us to discuss optimal trajectories in higher spatial dimensions or to calculate the optimal path for microswimmers with a time-dependent propulsion speed as we discuss now.

Let us consider the cost J[r(t), p(t), ψ(t), t] = T which we want to minimize subject to the equations of motion for microswimmers interacting with remote walls (Eqs. (1) and (2)) and the boundary conditions r(t = 0) = rA and r(t = T) = rB where r ≔ (x, z). Here ψ is the control variable, p is the costate, and T is the (unknown) traveling time corresponding to the fastest route. As the cost can be determined from the endpoint (or the “endtime”) alone, it can be considered as a pure endpoint cost, whereas the running cost is zero. Allowing for explicitly time-dependent microswimmers with a time-dependent source dipole strength σt and/or a time-dependent self-propulsion speed vt, this is a non-autonomous, fixed endpoint, free endtime optimal control problem. To solve it, we first construct the optimal control Hamiltonian as H(r, p, ψ, t) = F ⋅ p, where F is given by \(\dot{{\bf{r}}}=\left([{v}_{t}-{\sigma }_{t}/\left(4{z}^{3}\right)]\cos \psi ,[{v}_{t}-{\sigma }_{t}/{z}^{3}]\sin \psi \right)=:{\bf{F}}\). We now minimize H with respect to the (unconstrained) control variable ψ by solving ∂ψH = 0 yielding \(\tan {\psi }^{* }={u}_{z}/{u}_{x}\) where ψ* is the optimal control and where we have defined ux ≔ px(vt − σt/(4z3)) and uz := pz(vt − σt/z3). Using \({\mathcal{H}}({\bf{r}},{\bf{p}},t)\equiv H({\bf{r}},{\bf{p}},{\psi }^{* },t)\), after a few lines of straightforward algebra, we obtain the minimized (lower) Hamiltonian as \({\mathcal{H}}({\bf{r}},{\bf{p}},t)=\pm \sqrt{{u}_{x}^{2}+{u}_{z}^{2}}\). In the following, we use only the plus–branch, as the minus–branch does not seem to result in sensible trajectories. The optimal trajectory now follows from the Hamilton equations of motion \(\dot{{\bf{r}}}={\partial }_{{\bf{p}}}{\mathcal{H}}\) and \(\dot{{\bf{p}}}=-{\partial }_{{\bf{r}}}{\mathcal{H}}\) which have to be solved as a boundary value problem, such that r(t = 0) = rA and r(t = T) = rB where T is the unknown traveling time. Thus, up to now, we have determined four first-order differential equations which are complemented by four boundary conditions. Since this problem depends on the unknown T, we need one additional boundary condition. This is given by the Hamiltonian endpoint condition (or transversality condition) which determines the value of \({\mathcal{H}}({\bf{r}}(T),{\bf{p}}(T),T)\). This condition can be determined from the endpoint Lagrangian \(\bar{L}=T+{\boldsymbol{\nu }}\cdot {\bf{e}}\) where T is again the endpoint cost; ν is the (constant) endpoint costate and e = (r(T) − rB) (such that e = 0 represents the boundary conditions in standardized form). The Hamiltonian endpoint condition simply yields \({\mathcal{H}}({\bf{r}}(T),{\bf{p}}(T),T)=-{\partial }_{T}\bar{L}=-1\) for the considerd minimal time problem. (The remaining two transversality conditions for the costate, \({\bf{p}}(T)=-{\partial }_{{\bf{r}}(T)}\bar{L}\), just relate the two unknowns p(T) to two other unknowns v and hence do not provide any additional information.)

To numerically solve the Hamilton equations together with r(0) = rA, r(T) = rB and \({\mathcal{H}}({\bf{r}}(T),{\bf{p}}(T),T)=-1\), it is convenient to rescale time via \(t=Tt^{\prime}\) such that, in primed units, the terminal time is T = 1 and the endpoint boundary conditions are explicitly known as \({\bf{r}}(t^{\prime} =1)={{\bf{r}}}_{B}\); i.e. they no longer depend on the unknown variable T, which instead shows up in the equations of motion. In addition, to ease the usage of standard shooting methods to solve this boundary value problem, it is convenient to treat T as a dynamical variable, i.e. we introduce T =: a(t), leading to the additional equation of motion \(\dot{a}=0\), which we solve together with the Hamilton equations and the five mentioned boundary conditions. We have compared the results from this approach for time-independent microswimmers in the presence of flat walls with those from the Lagrangian approach and find the same trajectories. For time-dependent microswimmers, corresponding results are presented.

Let us finally mention that an alternative approach to treat the present non-autonomous optimal control problem would be to first transform it to an autonomous problem, by introducing an extra variable \(\dot{\tau }=1\) with τ(t = 0) = 0 such that τ(t) = t. The corresponding approach leads to the same numerical results to those presented in Fig. 5.