# Self-organized Beating and Swimming of Internally Driven Filaments

## Abstract

We study a simple two-dimensional model for motion of an elastic filament subject to internally generated stresses and show that wave-like propagating shapes which can propel the filament can be induced by a self-organized mechanism via a dynamic instability. The resulting patterns of motion do not depend on the microscopic mechanism of the instability but only of the filament rigidity and hydrodynamic friction. Our results suggest that simplified systems, consisting only of molecular motors and filaments could be able to show beating motion and self-propulsion.

Cilia and flagella are hair-like appendages of many cells which generate motion and are used for self-propulsion and to stir the surrounding fluid. They all share the characteristic architecture of their core structure, the axoneme, a common structural motive that was developed early in evolution. It is characterized by nine parallel pairs of microtubules, which are long and rigid protein filaments, that are arranged in a circular fashion together with a large number of dynein molecular motors . In the presence of a fuel which is ATP, the dynein motors attached to the microtubules generate relative forces while acting on neighboring microtubules. The resulting internal stresses induce the bending and a wave-like motion of the axoneme.

These biological systems are complex, they consist of a large number of different components and various patterns of motion have been observed. Attempts to model their behavior are either based on the assumption that some unknown control system generates oscillatory motor activity  or that a self-organized mechanism is at work . Generically, the latter involves a dynamical instability. Theoretical studies of simple models for collective action of molecular motors have demonstrated the possibility of such instabilities . Several examples of oscillatory motion of biological many-motor systems are known. Recently, it was suggested that spontaneous oscillations observed in muscles could be a property of the motor-filament system alone . This idea is supported by the fact that the oscillations continue to exist after all regulatory systems are removed  but also by the observation that an in vitro motor-filament system shows the signature of a dynamic transition . Furthermore, the observations that flagellar dyneins are able to generate oscillatory motion on microtubules , and that isolated and de-membranated flagella in solution containing ATP above a threshold concentration swim with a simple wave-like motion  support the idea that basic types of flagellar beating could result from a dynamic instability.

In this article, we introduce a simple two-dimensional model which reveals many physical aspects of the motion of an elastic filament driven by internal forces, that should be relevant for flagellar beating. Our approach is inspired by studies of semiflexible filaments subject to external forces , however, in our case all motion is induced by internal stresses. Our model consists of two incompressible but elastic filaments of length $`L`$ arranged at constant distance $`aL`$ and rigidly attached only at one end which we call the head. A large number of molecular motors and passive elements holding the filament pair together are assumed to generate a coarse-grained force per unit length $`f`$ which is an internal stress acting in opposite directions on the two filaments and induces the bending of the filament pair. The dynamic equations of this model define patterns of beating motion resulting from the internal forces which are assumed to oscillate. More interestingly, we show that characteristic wave-like patterns which propagate along the filament are generated most naturally by a dynamic instability of the motor-filament system, see Fig. 1 for examples. As we show below, the qualitative shapes of these patterns do not depend on the microscopic mechanism of force generation but only on the elastic properties of the filaments and on hydrodynamic friction. We demonstrate that these patterns lead to self-propulsion of the system and calculate the velocity of motion.

In order to define our model and to derive the dynamic equations, we start from the enthalpy functional

$$G=_0^L\left[\frac{\kappa }{2}(C(s))^2+f(s)\mathrm{\Delta }(s)+\mathrm{\Lambda }(s)(_s\stackrel{}{r})^2\right]𝑑s$$
(1)   
where $`\stackrel{}{r}(s)`$ is a parameterization of the shape of the filament pair by the arclength $`s`$, $`\kappa `$ is the bending rigidity, and $`C=\stackrel{}{n}_s^2\stackrel{}{r}`$ is the local curvature with the filament normal $`\stackrel{}{n}`$. The internal force density $`f`$ couples to the relative local displacement of the filament pair $`\mathrm{\Delta }(s)=a_0^sC(s^{})𝑑s^{}`$ . In order to impose the constraint of local incompressibility $`(_s\stackrel{}{r})^2=1`$, we have introduced the Lagrange multiplier $`\mathrm{\Lambda }(s)`$. The equation of motion can be written as

$$_t\stackrel{}{r}=\left(\frac{1}{\xi _{}}\stackrel{}{n}\stackrel{}{n}+\frac{1}{\xi _{}}\stackrel{}{t}\stackrel{}{t}\right)\frac{\delta G}{\delta \stackrel{}{r}}$$
(2)   
where $`\stackrel{}{t}\stackrel{}{t}`$ and $`\stackrel{}{n}\stackrel{}{n}`$ are projectors on the filament tangent and normal, and we assume local anisotropic friction with tangent and normal coefficients $`\xi _{}`$ and $`\xi _{}`$, respectively. The Lagrange multiplier is determined by the condition $`(_s\stackrel{}{r})_t_s\stackrel{}{r}=0`$ .

In order to keep the description simple, we consider small deformations of a filament parallel to the $`x`$-axis, $`\stackrel{}{r}(s)=(s+u(s),h(s))`$, which we describe by an expansion in the transverse and longitudinal displacements $`h`$ and $`u`$. To quadratic order in $`_xh(x)`$ we can write

$$G_0^L\left[\frac{\kappa }{2}(_x^2h)^2+af(x)(_xh(x)_xh(0))\right]𝑑x,$$
(3)   
where we use the $`x`$-coordinate as parameter. We first discuss transverse motion which for small deformations is independent of longitudinal forces  and satisfies the equation $`\xi _{}_th=\kappa _x^4h+a_xf`$ together with two boundary conditions at the head with $`x=0`$ and two conditions at the tail for $`x=L`$. We assume a free tail which implies $`_x^2h(L)=0`$ and $`\kappa _x^3h(L)=af(L)`$. At the head, we distinguish three different cases as shown in Fig. 1: (A) clamped head with $`h(0)=0`$ and $`_xh(0)=0`$; (B) fixed head with $`h(0)=0`$ and $`\kappa _x^2h(0)=a_0^Lf(x)𝑑x`$; and (C) a viscous load at $`x=0`$ with friction coefficient $`\eta `$ for which the condition on $`h(0)`$ in (B) is replaced by $`\eta _th(0)=af(0)\kappa _x^3h(0)`$.

We demonstrate the basic properties of this model, by first assuming that an oscillating force density with constant amplitude is generated by some unspecified mechanism: $`f_m(t)=\mathrm{Re}(\stackrel{~}{f}_0e^{i\omega t})`$. The total force density $`f`$ acting on the filament pair is the sum of the force $`f_m`$, internal dissipative forces and in general the forces of elastic elements which locally connect the filaments. Introducing the complex Fourier amplitude $`\stackrel{~}{h}`$, where $`h(x,t)=\mathrm{Re}(\stackrel{~}{h}(x)e^{i\omega t})`$, we can express the total force density as $`\stackrel{~}{f}=\chi \stackrel{~}{v}+\stackrel{~}{f}_0`$, where $`\stackrel{~}{v}i\omega a(_x\stackrel{~}{h}(x)_x\stackrel{~}{h}(0))`$ is the complex amplitude of the local sliding velocity $`v=_t\mathrm{\Delta }`$. The coefficient $`\chi =(\lambda +K/i\omega )`$ describes a viscoelastic response of the material between the filaments with dissipation coefficient $`\lambda `$ and elastic modulus $`K`$. The oscillating state is characterized by

$$\kappa _x^4\stackrel{~}{h}a^2i\omega \chi _x^2\stackrel{~}{h}+\xi _{}i\omega \stackrel{~}{h}=0.$$
(4)   
The homogeneous active force $`\stackrel{~}{f}_0`$ only enters via boundary conditions. Eq. (4) and boundary conditions represent an inhomogeneous linear system which is solved by $`\stackrel{~}{h}=Ae^{kx/L}`$ leading to four complex values of $`k`$. The corresponding coefficients $`A`$ are adjusted to satisfy the boundary conditions which leads to a solution with an amplitude proportional to the internal force $`\stackrel{~}{h}(x)\stackrel{~}{f}_0`$. We can distinguish two different regimes: (i) hydrodynamic friction dominates $`|\chi |^2\kappa \xi _{}/\omega a^4`$; (ii) internal viscoelasticity dominates $`|\chi |^2\kappa \xi _{}/\omega a^4`$. We can neglect in Eq. (4), $`\chi `$ in case (i) and $`\xi _{}`$ in case (ii). Fig. 2 shows examples of the amplitude $`H`$ and the gradient of the phase $`\varphi `$ of $`\stackrel{~}{h}(x)=H(x)e^{i\varphi (x)}`$ for $`\chi =0`$ and different boundary conditions as dashed lines. The corresponding time dependent solutions

$$h(x,t)=H(x)\mathrm{cos}(\omega t\varphi (x))$$
(5)   
are propagating wave-like shapes qualitatively similar to those shown in Fig. 1. The sign of the local propagation velocity $`v_p=\omega /_x\varphi `$ of the phase allows us to determine the direction of apparent wave propagation.

We have thus developed the framework to calculate and analyze wave-propagating solutions of our model and can now study motion generated by the properties of the motor-filament system via a Hopf bifurcation. We assume that the material between the two filaments which contains both molecular motors and passive elements has properties which can be characterized on a coarse-grained level by a nonlinear history-dependent response function. We will study the instability of a non-moving solution $`h(x)=0`$ towards wave-like patterns. For this case it is sufficient to consider only small amplitudes, $`|_xh|1`$ as described above. Furthermore, in this regime the local sliding velocity $`v`$ is small and we can ignore nonlinearities in $`v`$ and restrict ourselves to the frequency dependent linear response $`\stackrel{~}{f}=\chi \stackrel{~}{v}`$. Here, we have set the artificially introduced force $`\stackrel{~}{f}_0=0`$ and characterize both passive and active material properties by the complex response function $`\chi (\omega ,\mathrm{\Omega })`$ which can e.g. be calculated explicitly for a simple model  or measured experimentally . The out-of-equilibrium nature of the system is characterized by the parameter $`\mathrm{\Omega }`$ which can for example be identified with the ATP concentration. Note, that for an active system $`\chi `$ can have unusual behaviors which formally correspond to a negative friction ($`\mathrm{Re}(\chi )<0`$) or a negative elastic response ($`\mathrm{Im}(\chi )>0`$).

In the case $`\stackrel{~}{f}_0=0`$, Eq. (4) and boundary conditions becomes a homogeneous linear system which always has the solution $`\stackrel{~}{h}(x)=0`$ and which can now be reinterpreted as an eigenvalue problem for $`\chi `$. Spontaneous motion corresponds to nontrivial solutions to this problem. A discrete set of such solutions $`\stackrel{~}{h}_i`$ exists, each $`\stackrel{~}{h}_i`$ corresponds to a complex eigenvalue $`\chi =\overline{\chi }_i(\omega )`$, $`i=1,2\mathrm{}\mathrm{}`$ which we order according to $`|\overline{\chi }_i(\omega )||\overline{\chi }_{i+1}(\omega )|`$.

Consider now a system initially at equilibrium with $`\mathrm{\Omega }=0`$. If $`\mathrm{\Omega }`$ is increased, an instability occurs as soon as a critical value $`\mathrm{\Omega }_c`$ is reached for which $`\chi (\omega _c,\mathrm{\Omega }_c)=\overline{\chi }_i(\omega _c)`$ for a frequency $`\omega _c`$. In the vicinity of this point, the system develops for $`\mathrm{\Omega }>\mathrm{\Omega }_c`$ motion with this frequency and a shape characterized by the nontrivial solution $`\stackrel{~}{h}_i(x)`$. This scenario applies to a supercritical bifurcation. Nonlinear terms of the response function and nonlinear corrections to the simple Monge representation can become important for larger $`\mathrm{\Omega }`$, or they could change the nature of the bifurcation to subcritical. Typically, the instability occurs for the smallest value $`\chi =\overline{\chi }_1(\omega )`$ since larger $`|\chi |`$ require larger values of $`\mathrm{\Omega }`$ which correspond to more system activity. Note, that the resulting pattern of motion is independent of the microscopic mechanism which leads to the instability. It is sufficient that the active material is capable to generate the response $`\chi =\overline{\chi }_i`$ .

Fig. 2 displays examples of the amplitude and the gradient of the phase of $`\stackrel{~}{h}_1(x)`$ for boundary conditions (A) and (B), snapshots of the corresponding motion are shown in Fig. 1. The boundary conditions play an essential role in selecting different types of motion. Observing the sign of $`_x\varphi `$ which determines the direction of the phase velocity we find that for clamped head (A) the wave propagates from the tail towards the head while for case (B) it propagates in the opposite direction. The amplitude $`H(x)`$ also differs significantly between cases (A) and (B), see Fig. 2. The case (C) of a free head with viscous load $`\eta `$ is similar to case (B) and therefore not shown in Fig. 2. For this example, the qualitative properties of motion induced by the dynamic instability are the same as those of the system driven by a homogeneous force $`\stackrel{~}{f}_0`$, see Fig. 2. In fact, for the parameters chosen, $`|\chi _1|^2\kappa \xi _{}/(\omega a^4)`$ and the corresponding solution is not far from the solution for $`\chi =0`$. The case of homogeneous force $`\stackrel{~}{f}_0`$ is simple and allows us to explain the effect of boundary conditions. A homogeneous internal force $`\stackrel{~}{f}_0`$ can be rewritten as boundary terms in the expression of the energy: $`Gaf_0h(L)af_0h(0)aLf_0_xh(0)+_0^L𝑑x(_x^2h)^2\kappa /2`$. Its action is equivalent to two opposite transverse forces $`af_0`$ acting at both ends together with a torque $`aLf_0`$ applied at the head. In the case of a clamped head this apparent force and torque are suppressed and the system is driven by a virtual force at the tail, propagating the wave towards the head . If the head is not clamped, the virtual oscillating torque at the head can propagate a wave in the opposite direction.

Can these beating patterns propel the filament? Time-reversal symmetry has to be broken, $`h(x,t)h(x,t)`$, for propulsion to be possible . According to Eq. (5), this requirement is fulfilled since $`_x\varphi 0`$. Because of the symmetry $`h(x,t)=h(x,t+\pi /\omega )`$, there can be no net motion in the transverse direction. In order to estimate longitudinal motion, we have to study the displacement $`u(x)`$. To second order in $`_xh`$, we can write

$$u(x)u(0)\frac{1}{2}_0^x(_xh)^2𝑑x^{},$$
(6)   
indicating that the dynamics of $`u(x)`$ is governed by the motion $`h(x,t)`$. Note that $`u(x)u(0)`$ is small, but the filament displacement $`u(0)`$ can become large. The longitudinal component $`f_l`$ of the hydrodynamic force density $`(\xi _{}\stackrel{}{n}\stackrel{}{n}+\xi _{}\stackrel{}{t}\stackrel{}{t})_t\stackrel{}{r}`$ acting locally on the filament is given by $`f_l(\xi _{}\xi _{})_xh_th\xi _{}_tu(x,t)`$ in our approximation. The velocity of motion $`V`$ is the time-average of $`_tu(0)`$ and follows from the condition that the total longitudinal force vanishes. If an isotropic viscous load is attached to the head, this condition is $`_0^Lf_l𝑑x+\eta _tu(0)=0`$ and we find $`V=V_0/(1+\eta /\xi _{}L)`$ where

$$V_0=\left(\frac{\xi _{}}{\xi _{}}1\right)\frac{\omega }{2L}_0^LH(x)^2_x\varphi dx$$
(7)   
is the no-load velocity. If the head is not permitted to move, the filament generates a force $`F=V_0\xi _{}L`$ at the head. Note, that for isotropic friction both $`V`$ and $`F`$ vanish. For a semiflexible rod-like filament $`\xi _{}/\xi _{}2`$  and the direction of motion is opposite to the direction of phase propagation.

The parameters chosen in Figs. 1 and 2 correspond e.g. to $`L40\mu `$m, $`a20`$nm, $`\xi _{}210^3`$Ns/m<sup>2</sup>, $`\kappa 410^{22}`$Nm<sup>2</sup>, which is the elastic modulus of about $`20`$ Microtubules , and a frequency $`\omega /(2\pi )30`$s<sup>-1</sup>. For this choice, we find a critical value $`\chi _1(10+20i)`$Ns/m<sup>2</sup>. We choose an amplitude of $`\stackrel{~}{h}_1`$ with maximal value $`H/L0.1`$. In this case, the maximal local sliding velocity is $`v8\mu `$m/s. In axoneme, dynein motors are spaced every $`24`$nm along the microtubules. Assuming that only one microtubule pair is active, we estimate that for this choice a force per motor of $`4`$pN corresponds to the critical value $`f|\chi _1|v`$. This is a typical force created by molecular motors. Larger forces could be necessary to generate beating with larger amplitudes. Our result suggests that in this case several microtubule pairs could be active at the same time thus allowing for smaller forces per motor. Using the motion $`\stackrel{~}{h}_1`$ obtained for boundary conditions (C) with a viscous load $`\eta =510^8Ns/m`$, we find a spontaneous velocity of lateral motion $`V40\mu `$m/s, which is significantly larger than local sliding velocities and not far from experimentally observed values.

We have demonstrated that a pair of elastic filaments held together by an active, force-generating material, can induce wave-like patterns by a dynamic instability of the system. This study is motivated by biological flagella such as those of sperms which use such motion for self-propulsion. Our model suggests that the boundary conditions imposed at the ends select the type of beating pattern. This could be tested by micro-manipulation experiments which apply external forces and torques at the ends of beating flagella. We have restricted our study to a two-dimensional system, small deformations and the linear regime of the instability. Nonlinearities could play an important role in the regime of rapid self-propulsion with large amplitudes. Furthermore, the possibility of torsional motion in three-dimensional systems could allow for new types of behavior.

The observation that wave-like motion can be generated in a self-organized way raises the idea that sophisticated control mechanisms may have evolved after the development of the basic axonemal structure in order to fine-tune the system and to create more complex types of motion. This concept suggests that artificially constructed systems consisting only of motors and filaments could already undergo beating motion and self-propulsion.

We thank R. Everaers for stimulating collaborations, M. Bornens and H. Delacroix for introducing us to cilia and flagella, A. Parmeggiani and C. Wiggins for useful discussions and A. Maggs for a critical reading of the manuscript.