A theoretical expression is derived for the mechanical contribution to the mean swimming speed of a vibrating two-sphere, consisting of two spheres connected by an elastic spring, and immersed in a viscous incompressible fluid. The spring provides a harmonic potential for oscillations about a mean distance between centers. The system is made to oscillate at a chosen frequency by activating forces which sum to zero. The mechanical contribution to the resulting mean swimming velocity is calculated from the mechanical equations of motion and the corresponding impedance matrix of linear response. The frequency-dependent pair friction coefficients are found from approximate expressions derived earlier. The mechanical contribution is calculated to second order in the amplitude of stroke as a function of the scaling number, a dimensionless combination of size, frequency, and kinematic viscosity. Retarded friction and added mass determine the functional behavior.

## I. INTRODUCTION

The swimming velocity of a deformable body in a fluid can be decomposed as a sum of two contributions. The first contribution follows from the equations of classical mechanics' for the individual body parts. These equations involve friction and added mass effects. The second contribution follows from the Navier–Stokes equations for fluid velocity and pressure. Both contributions are a result of rectification due to nonlinear interactions. We call the first contribution mechanical and the second one hydrodynamic. In the present article, we study the first mechanical contribution for a vibrating two-sphere, consisting of two spheres, connected by an elastic spring, and immersed in a viscous incompressible fluid. A picture of a two-sphere is shown in Fig. 1.

The work of Klotsa *et al.*^{1} established the vibrating two-sphere as an important model in the theory of swimming. The model is surprisingly complex but sufficiently simple to allow some mathematical analysis. Its partial solution can improve understanding and elucidate principles. This will be helpful in the analysis of more sophisticated swimmers.

In the case of the two-sphere, the body parts are the two spheres. The pair mobilities can be evaluated with various degrees of accuracy.^{2} Hubert *et al.*^{3} calculated the mean swimming speed of a two-sphere with Oseen type interactions for large distance between centers. For small amplitude of stroke, we generalized this result to shorter distance^{4} and included added mass effects in dipole approximation. We used a kinematic approach in which the relative distance between centers is prescribed to oscillate harmonically. The mean speed was obtained as a function of the oscillation frequency.

In the following, we study the two-sphere chain with harmonic elastic interaction, set into motion by periodic activating forces, which sum to zero. First, we consider instantaneous hydrodynamic interactions without specifying the distance dependence and derive a novel expression for the mean speed of small amplitude swimming in terms of the impedance matrix of linear response theory. For Oseen type hydrodynamic interactions, with or without dipolar added mass effects, the expression for the mean speed reduces to that derived earlier by expansion in harmonics.^{4}

Next, we generalize the theory to frequency-dependent friction coefficients, as first developed by Stokes for a single sphere on the basis of the linearized Navier–Stokes equations.^{5} It was realized by Basset^{6} that this corresponds to retarded friction, as described by a memory term in the equation of motion for the sphere.^{6,7} We use also retarded hydrodynamic pair interactions^{8} in the impedance matrix.

The second contribution to the swimming velocity results from flow generated by the swimmer via the movement of the spheres and via Reynolds stress. Since we need only the time-averaged swimming velocity, it suffices to solve the steady-state Stokes equations.^{9} The generated net flow is called the effect of steady streaming.^{10} We postpone the study of its contribution to the mean swimming velocity to a separate article.

Here, we consider specifically two different systems, as described by four models of increasing sophistication. In the first system, the two spheres are neutrally buoyant. This is indicated by *N*. In the second system, the larger sphere is again neutrally buoyant, but the smaller sphere is loaded such that it has twice the mass of the larger sphere. This system is indicated by *L*. In the first model, the hydrodynamic interactions between spheres are treated in Oseen-approximation. This is indicated by *O* and is the model studied by Hubert *et al.*^{3} In the second model, we take approximate account of added mass effects by including the single sphere added mass in the equations of motion. The model is called the Oseen^{*}-model and is indicated by *S*. The third model takes account of pair added mass effects in dipole approximation.^{11} This model is indicated by *D*. In the fourth model, the friction coefficient of a single sphere is replaced by the frequency-dependent coefficient introduced by Stokes,^{5} and the pair hydrodynamic interactions are also made frequency-dependent, corresponding to retarded hydrodynamic interactions.^{8} In addition, the single sphere friction coefficient is corrected for a pair effect such that the added masses are identical with those of the dipole model. This is called the retarded friction model and is indicated by *R*. Thus, *NR* refers to the neutrally buoyant system described by the retarded friction model.

## II. MECHANICS WITH SIMPLE FRICTION

We consider two spheres of radii *a* and *b* and mass densities *ρ _{a}* and

*ρ*immersed in a viscous incompressible fluid of shear viscosity

_{b}*η*and mass density

*ρ*. The fluid is of infinite extent in all directions. We assume that at all times, the centers are located on the

*z*axis of a fixed Cartesian system of coordinates. The dynamics of the system is governed by a harmonic interaction potential $ V int = 1 2 k ( z 2 \u2212 z 1 \u2212 d ) 2$, depending on the instantaneous configuration of centers $ z 1 ( t ) , z 2 ( t )$, and by actuating forces $ E 1 ( t )$ and $ E 2 ( t )$ in the

*z*direction and varying harmonically in time with frequency

*ω*. We assume in this section that friction with the fluid can be described by a 2 × 2 friction matrix with elements which depend on distance $ z ( t ) = z 2 ( t ) \u2212 z 1 ( t )$ but are independent of frequency

*ω*. Added mass effects are described by a 2 × 2 mass matrix $ m ( z )$, calculated from potential theory.

In the course of time, the motion of the system becomes periodic with period $ T = 2 \pi / \omega $ and with a steady mean swimming velocity $ U \xaf s w$, and the average being over a period. In the following, we derive a simple expression for $ U \xaf s w$ for the case of small amplitude harmonic variation of the activating forces $ E 1 ( t ) = \phi 1 \u2009 cos ( \omega t )$ and $ E 2 ( t ) = \u2212 E 1 ( t )$. The mean velocity will be calculated to the second order in the amplitude $ \phi 1$.

*z*. The dynamics of the system is assumed to be governed by the approximate equations of motion,

^{12}

*V*are invariant under translations.

_{int}We abbreviated $ E = ( E 1 , E 2 )$. We assume $ E 2 ( t ) = \u2212 E 1 ( t )$, thus ensuring that the total applied force vanishes at any time. Another choice would imply that the center of mass of the swimmer is swung back and forth by an applied force during its motion. The activating forces cause internal motions of the body without acceleration of its center of mass.

We evaluate these quantities from the linear response.

## III. CORRELATIONS

*z*=

*d*. For oscillating activating forces at frequency

*ω*, it is convenient to use complex notation. The equation for the amplitudes $ \psi \omega = ( q 1 \omega , q 2 \omega , U 1 \omega , U 2 \omega )$ then takes the form

*d*,

*z*=

*d*. In complex notation, we use that for two real quantities

*A*(

*t*) and

*B*(

*t*) oscillating at frequency

*ω*,

For the Oseen^{*} and Oseen-Dipole models, Eq. (3.11) agrees with the expression for the mean swimming velocity derived by a different method.^{4} The mean swimming velocity must be supplemented with the additional contribution from the motion of the fluid. First, we show how the above derivation must be modified when the friction coefficients depend on frequency and find a general expression for the moment in Eq. (3.11).

## IV. FREQUENCY-DEPENDENT FRICTION

It was first shown by Stokes^{5} on the basis of the linearized Navier–Stokes equations that the flow pattern about an oscillating sphere differs from that about a steadily moving sphere. Correspondingly, the friction coefficient of a single sphere depends on frequency. For a sphere of radius *a* with no-slip boundary condition, it takes the value $ \zeta 0 = 6 \pi \eta a$ at zero frequency, but at higher frequency, the complex coefficient $ \zeta ( \omega )$ differs. It was shown by Basset^{6} that the frequency-dependent friction corresponds to a retardation effect in the hydrodynamic force exerted by the fluid. If the fluid is viscoelastic, then that is an additional reason for retardation. In the following, we discuss how the retardation effect can be taken into account in the theory of swimming.

^{7}

*m*. The added masses are defined by

_{ij}*m*

_{1}and

*m*

_{2}are the bare masses. These definitions imply that the friction coefficients do not have an imaginary part proportional to

*ω*at high frequency. Since only the impedance matrix $ z ( \omega )$ is relevant, we could alternatively use the bare masses only and include the added masses in the friction coefficients.

^{12}

^{,}

^{13,14}We call the triangle as the impedance triangle. With 3-vectors $ P 0 = ( P 0 , 0 ) , \u2009 P 1 = ( P 1 , 0 ) , and P 2 = ( P 2 , 0 )$, and difference vectors $ S 1 = P 1 \u2212 P 0$ and $ S 2 = P 2 \u2212 P 0$, the value $ N ( \omega )$ is the third component of the cross product $ N ( \omega ) = S 1 \xd7 S 2$. The vector $ N ( \omega )$ can be regarded as a generalized torque corresponding to the impedance matrix $ z ( \omega )$.

*s*. It follows from Eq. (3.9) that $ d 2 ( 1 ) = \u2212 d 1 ( 1 )$, so that it suffices to consider $ ( d 1 ( 1 ) ) 2 \xaf$. For $ E 1 ( t ) = \phi k a \u2009 cos ( \omega t )$, the mean square displacement is calculated as

*k*. The dependence on

*k*has been fully absorbed in the expression (4.18) for the amplitude factor $\phi $.

## V. EFFECTIVENESS

*et al.*

^{16}

^{,}$ R = M 2$. We define the dimensionless effectiveness

*V*(

*R*) as

*ε*at a scale number

*R*. Clearly, one wishes to optimize this value. We recall that the efficiency, as defined by Shapere and Wilczek,

^{16}involves both the mean swimming velocity and the power dissipated during a period. The effectiveness of the stroke is defined regardless of the power.

This shows an interesting symmetry between mass and friction.

^{4}

*R*and

_{xM}*V*. The function $ V M 2 ( R )$ has maximum value

_{xM}*V*at

_{xM}*R*, and the curve has half-width $ \Delta R = 2 3 R x M$. We call a curve in the

_{xM}*VR*-plane with behavior described by Eq. (5.7) a

*VR*2-curve.

*R*

^{4}

^{4}we derived explicit expressions for the parameters

*R*and

_{xM}*V*for three specific models. The Oseen-model is based on the Oseen-approximation for the hydrodynamic pair mobility, and added mass effects are neglected. The Oseen

_{xM}^{*}-model uses the same approximation for the pair mobility and includes the single sphere added mass effect. The Oseen-Dipole model describes in addition the pair added mass in dipole approximation. We note that for the Oseen

^{*}-model $ M = m 1 * + m 2 *$, where $ m j * = m j + 1 2 m f j$ with $ m f j = ( 4 \pi / 3 ) \rho a j 3$. For the Oseen-Dipole model, $ M = m 11 + 2 m 12 + m 22$ with

*m*given below. The factors $ Z 0 d$ and

_{ij}*Z*

_{0}are the same for the three models. The Oseen mobility matrix corresponds to the friction matrix,

^{11}

^{,}

*z*to terms of order $ z \u2212 6$, this agrees with the expression given by Lamb.

^{21}The diagonal elements include the single sphere added mass effect.

^{4}except for a copying error in the coefficient

*r*in Eq. (6.16) of Ref. 4. This should read instead

_{b}*et al.*

^{3}for large

*d*/

*a*. Their Eq. (18) reads in our notation

^{*}-model, the same equations hold with

*m*replaced by $ m i * = m i + ( 2 \pi / 3 ) \rho a i 3$ or with

_{i}*σ*replaced by $ \sigma i * = \sigma i + 1 2$.

_{i}The numerator of the last factor in Eq. (5.17) for fixed $ \sigma a , \xi , \delta $ can switch sign as a function of *σ _{b}*. For example, for $ \sigma a = 1 , \u2009 \xi = 1 / 2 , and\u2009\u2009 \delta = 3$, the expression switches sign at $ \sigma b = 8 / 3$. This corresponds to a change of swimming direction.

*V*depends on the mass densities $ \rho a , \rho b$ only via the ratio $ r = \rho b / \rho a = \sigma b / \sigma a$. The expression is, from Eq. (5.11),

_{xMO}^{*}-model, replace

*r*by $ r * = \sigma b * / \sigma a *$. For the position of the maximum, one finds

*R*. These simple properties do not hold for the Oseen-Dipole model. For typical values of $ \xi , \delta $, the value of

_{xMS}*V*ranges from positive to negative values. From Eq. (5.19), one sees that in the Oseen-model, the sign switch occurs at

_{xMO}*η*and mass density

*ρ*. The mass densities of the spheres can be varied and take values $ \rho a = \sigma a \rho $ and $ \rho b = \sigma b \rho $. We consider first neutrally buoyant spheres with $ \sigma a = 1 and \sigma b = 1$. In Fig. 2, we show the effectiveness

*V*as a function of

_{M}*R*for the Oseen-model, the Oseen

^{*}-model, and the Oseen-Dipole model. In the notation introduced at the end of Sec. I, the three models are indicated, respectively, as the

*NO*-model, the

*NS*-model, and the

*ND*-model. In Table I, we list the values of the parameters

*R*and

_{xM}*V*for the three models.

_{xM}$ ( \rho a , \rho b )$ . | Model . | $ R 2 x M$ . | $ 10 6 V 2 x M$ . | $ Z 0 / ( \eta a )$ . | $ M / ( \rho a 3 )$ . |
---|---|---|---|---|---|

$ ( \rho , \rho )$ | Oseen | 4.571 | 2480 | 21.5423 | 4.712 |

$ ( \rho , \rho )$ | $ Osee n *$ | 3.048 | 2480 | 21.5423 | 7.069 |

$ ( \rho , \rho )$ | Dipole | 3.072 | 2538 | 21.5423 | 7.012 |

$ ( \rho , \rho )$ | RetFr | 5.255 | 1298 | 21.7598 | 7.014 |

$ ( \rho , 16 \rho )$ | Oseen | 1.714 | –7440 | 21.5423 | 12.566 |

$ ( \rho , 16 \rho )$ | $ Osee n *$ | 1.444 | –5874 | 21.5423 | 14.923 |

$ ( \rho , 16 \rho )$ | Dipole | 1.449 | –5879 | 21.5423 | 14.866 |

$ ( \rho , 16 \rho )$ | RetFr | 4.409 | –2827 | 21.7598 | 14.868 |

$ ( \rho a , \rho b )$ . | Model . | $ R 2 x M$ . | $ 10 6 V 2 x M$ . | $ Z 0 / ( \eta a )$ . | $ M / ( \rho a 3 )$ . |
---|---|---|---|---|---|

$ ( \rho , \rho )$ | Oseen | 4.571 | 2480 | 21.5423 | 4.712 |

$ ( \rho , \rho )$ | $ Osee n *$ | 3.048 | 2480 | 21.5423 | 7.069 |

$ ( \rho , \rho )$ | Dipole | 3.072 | 2538 | 21.5423 | 7.012 |

$ ( \rho , \rho )$ | RetFr | 5.255 | 1298 | 21.7598 | 7.014 |

$ ( \rho , 16 \rho )$ | Oseen | 1.714 | –7440 | 21.5423 | 12.566 |

$ ( \rho , 16 \rho )$ | $ Osee n *$ | 1.444 | –5874 | 21.5423 | 14.923 |

$ ( \rho , 16 \rho )$ | Dipole | 1.449 | –5879 | 21.5423 | 14.866 |

$ ( \rho , 16 \rho )$ | RetFr | 4.409 | –2827 | 21.7598 | 14.868 |

For our choice of parameters $ \xi = 1 / 2 , \u2009 \delta = 3$, we have $ R xMNO = 32 / 7 = 4.571$ and $ V xMNO = 5 / 2016 = 0.002 \u2009 480$. As a function of $ log 10 R$, the curve resembles a Gaussian. For the *NS*-model, the corresponding numbers read: $ R xMNS = 64 / 21 = 3.048$ and $ V xMNS = 5 / 2016 = 0.002 \u2009 480$. We compare with the curve $ V MND ( R )$ for the *ND*-model, characterized by the same friction matrix Eq. (5.14), but with mass matrix given by Eq. (5.15). The function $ V MND ( R )$ has the form Eq. (5.7) with parameter values *R _{xMND}* = 3.063 and

*V*= 0.00 2538. Thus, the two $ V R 2 \u2212$curves are nearly identical.

_{xMND}For the *NS*-model, the points *P _{j}* of the impedance triangle, defined below Eq. (4.16), have coordinates $ P 0 = ( 5.386 \eta a , 0 ) , \u2009 P 1 = ( 21.542 \eta a , \u2212 19.149 \omega \rho a 3 ) , \u2009 and P 2 = ( 10.771 \eta a , \u2212 2.394 \omega \rho a 3 )$ at $ R = R xMO$. The value of the coefficient $ N \u2032$ is $ N \u2032 = 21.15 \eta \rho a 4$. For the

*ND*-model, the points are $ P 0 = ( 5.386 \eta a , \u2212 0.089 \omega \rho a 3 ) , \u2009 P 1 = ( 21.542 \eta a , \u2212 19.308 \omega \rho a 3 ) , and P 2 = ( 10.771 \eta a , \u2212 2.413 \omega \rho a 3 )$, and $ N \u2032 = 21.47 \eta \rho a 4$.

As a second example, we consider the *L*-system with the same geometry, but with mass densities $ \rho a = \rho $ and $ \rho b = 16 \rho $, so that $ m 2 = 2 m 1$. For this geometry $ r 0 = 8 / 3$, so that the *r*-values of the two systems, *r* = 1 and *r* = 16, lie on opposite sides of *r*_{0}. Hence, the two systems swim in opposite directions. We list the corresponding values in Table I.

## VI. RETARDED FRICTION AND ADDED MASS

In this section, we study the effect of frequency-dependent friction. We do not aim at a complete calculation based on exact frequency-dependent friction coefficients but limit ourselves to the dominant effects for not too small distance between centers. At very large distance, the interaction tends to zero, and the friction matrix is diagonal with single sphere terms given by the frequency-dependent friction coefficient first calculated by Stokes.^{5}

The mutual mobility will be approximated by a result,^{8} based on a single Green function between the spheres and on Faxén theorems. We have used the approximation earlier in a study of velocity relaxation of a pair of spheres.^{17} In the limit of zero frequency, the approximation reduces to that of Rotne and Prager,^{18} as generalized to spheres of different sizes.^{19} The Rotne–Prager approximation is popular in soft matter modeling.^{20} The self-mobility is given by the single sphere friction coefficient, corrected by a distance-dependent term, as detailed below, chosen such that the mass matrix is identical with that of the Oseen-Dipole model, given by Eq. (5.15).

*A*and

*B*is, respectively,

^{7}

^{,}

*z*between centers reads

^{8}

*d*. With the above expressions, we find

*d*,

To the order shown, this is identical with the exact result.^{19}

In Fig. 3, we show the mechanical effectiveness $ V M ( R )$ for the *RF*-model as specified above. For the *N*-system, the maximum is at $ R mMNR = 5.573$ and has value *V _{mMNR}* = 0.001 300. We compare with the

*VR*2-curve for these parameter values. It is evident that the retardation of friction has a significant effect. In a plot of

*V*vs

*R*, the half-width of the actual curve is $ \Delta R = 73.28$, whereas for the $ V R 2 \u2212$ curve, $ \Delta R 2 = 2 3 R mMNR = 19.30$. The values for the minimum for the

*L*-system with $ m 2 = 2 m 1$ are

*R*= 4.389 and $ V mMLR = \u2212 0.002 \u2009 828$. The half-width of the spectral curve is $ \Delta R = 64.02$, whereas for the $ V R 2 \u2212$ curve, $ \Delta R 2 = 2 3 R mMLR = 15.20$. Hence, in both cases, the retarded friction causes a shift to higher frequency of the extremum and a strong broadening of the spectral curve. In Fig. 4, we show the analogue of Fig. 2 for the

_{mMLR}*L*-system, and in Fig. 5, we show the analogue of Fig. 3.

## VII. DISCUSSION

In the earlier, we analyzed the mechanical contribution to the mean swimming velocity of a two-sphere swimmer. We showed that this contribution can be expressed in terms of the impedance matrix of the two-sphere. The latter embodies the linear response of the two-sphere to applied oscillatory forces acting on either sphere. The response depends strongly on frequency. We studied this dependence for various, more or less realistic choices of the hydrodynamic interactions incorporating the effects of friction and added mass.

Part of the desired information on the mean swimming velocity is missing. We must yet add the hydrodynamic contribution due to fluid flow generated by the oscillatory motion of the two spheres and affected by the Reynolds stress in the fluid. Many aspects of the required calculations were discussed by Derr *et al.*^{15} In separate work, we shall compare the hydrodynamic and mechanical contributions. Once a reliable calculation of the mean swimming velocity of the vibrating two-sphere is established, we can apply the method to other systems of interest, for example to a vibrating chain of three spheres.

## AUTHOR DECLARATIONS

### Conflict of Interest

The authors have no conflicts to disclose.

### Author Contributions

**B. U. Felderhof:** Conceptualization (lead); Writing – review & editing (lead).

## DATA AVAILABILITY

The data that support the findings of this study are available within the article.

## REFERENCES

*Microhydrodynamics: Principles and Selected Applications*

*The Shoelace Book*