Swirling Instability of the Microtubule Cytoskeleton

David B. Stein
These authors contributed equally
Center for Computational Biology, Flatiron Institute, 162 5th Ave., New York, NY 10010
Gabriele De Canio
These authors contributed equally
Department of Applied Mathematics and Theoretical Physics, Centre for Mathematical
Sciences, University of Cambridge, Cambridge CB3 0WA, United Kingdom
Eric Lauga
Department of Applied Mathematics and Theoretical Physics, Centre for Mathematical
Sciences, University of Cambridge, Cambridge CB3 0WA, United Kingdom
Michael J. Shelley
Center for Computational Biology, Flatiron Institute, 162 5th Ave., New York, NY 10010
Courant Institute, New York University, 251 Mercer St., New York, NY 10012
Raymond E. Goldstein
Department of Applied Mathematics and Theoretical Physics, Centre for Mathematical
Sciences, University of Cambridge, Cambridge CB3 0WA, United Kingdom

Swirling Instability of the Microtubule Cytoskeleton

In the cellular phenomena of cytoplasmic streaming, molecular motors carrying cargo
along a network of microtubules entrain the surrounding fluid.
The piconewton forces produced by individual motors are sufficient to deform
long microtubules, as are the collective fluid flows generated by many moving motors.
Studies of streaming during oocyte development in the fruit
fly D. melanogaster have shown a transition from a spatially-disordered
cytoskeleton, supporting flows with only
short-ranged correlations, to an ordered state with a cell-spanning vortical flow.
To test the hypothesis that this transition is driven by fluid-structure interactions
we study a discrete-filament model and a coarse-grained continuum theory
for motors moving on a deformable cytoskeleton, both of which
are shown to exhibit a swirling instability to spontaneous large-scale rotational
motion, as observed.

A striking example of fluid-structure interactions within cells
[1] occurs in oocytes of the fruit fly
Drosophila melanogaster [2].  These develop over
a week from a single cell through repeated rounds of cell division,
differentiation and growth, ultimately reaching hundreds of microns across.
This pathway has historically
been divided into $14$ stages, and it is in stages $9-11$, at days $6.5-7$ [3],
that fluid motion is most noticeable.  In stage $9$ (Fig. <ref>), microtubules (MTs) reach
inward from the oocyte periphery, forming a dense assembly along which
molecular motors (kinesins) move at tens of nm/sec, carrying messenger
RNAs and other nanometric particles.  This motion
entrains the surrounding fluid, producing cytoplasmic streaming [4, 5]
that can be visualized several ways: in brightfield by the motion of endogenous
particles [6, 7, 8], via their
autofluorescence [9, 10], and through a combination of particle image velocimetry
and fluorescently labelled microtubules
[11, 12, 13].
Previous work [7, 11]
revealed that these flows initially take the form of transient, recurring
vortices and jets whose correlation length is a fraction of the cell scale, with
no long-range order.  But by stage $11$, a dramatic reconfiguration of the cytoskeleton occurs, coincident with the appearance of a vortex
spanning the entire cell [6, 7, 10, 14].

Kinesin motors move from the minus ends of microtubules (attached to the
oocyte periphery) to the plus ends (free in the interior).
Transport of cargo through the network
depends on motor-microtubule binding details [15, 16] and
on the mesh architecture [17, 18].  As a motor pulls cargo toward the plus
end the filament experiences a localized minus-end-directed compressive force,
as in Euler buckling.  For a filament of length $L$ and
bending modulus $A$ [19],
the buckling force is $\sim \pi A/L^2\sim 60\,{\rm pN}/L^2$,
where $L$ is measured in microns.  Thus, a
kinesin's force of several pN [20] can buckle
MTs $10-40\,\mu$m long.

Cytoplasmic streaming flows in the Drosophila oocyte.  (a) Experimental
flow field [28]
and schematic of the disordered swirling flows and microtubule organization in early stages
of development.  (b) Later flows organize into a single vortex as MTs lie parallel to the cell

The coupled filament-motor problem is richer than
Euler buckling because a motor exerts a “follower force" [21]
that is aligned with the filament.  This feature breaks the
variational structure of the problem and allows a filament pinned at its minus end to
oscillate even at zero Reynolds number [22, 23, 24].
By exerting a force on the fluid a motor induces long-range flows
which, if compressive, can further deform
filaments [25, 26].

Discrete filament computations.  (a) $N$ equally spaced elastic filaments clamped at their
attachment points, reach inward from a no-slip spherical shell. Each has a continuous distribution of
tangential point forces (red)
that (b) exert a force ${\bf \Gamma}$ on the fluid and an equal and opposite
compressive force on the filament. Synchronous oscillations ($N=7$,
$\sigma=1700$), (d) steady, bent configuration ($N=9$, $\sigma=500$) and swirling flow field.

It has been hypothesized [10, 14] that the transition from disordered
streaming flows to a single vortex in stage $11$ is a consequence of the kinds of fluid-structure
interactions described above, facilitated by a decrease in cytoplasmic viscosity that accompanies
the disappearance of a coexisting network of the biopolymer f-actin.
Here, through a combination of direct computations on the coupled filament-flow problem
[23]
and studies of a recent continuum theory for dense filament suspensions
[27], we confirm this hypothesis
by showing the existence of a novel swirling instability of the cytoskeleton.

The swirling instability can be understood in a
simplified model of the oocyte: a rigid sphere of radius $R$ containing a fluid of viscosity
$\mu$, with $N$ elastic filaments reaching inwards from clamped attachment points
equally spaced around the equator.  A slice in the filament plane
(Fig. <ref>(a)) appears like a confocal slice of the oocyte (Fig. <ref>).
The filaments have a radius $r$, a constant length $L$, bending modulus $A$ and a uniform line density $f$ of
follower forces (Fig. <ref>(b)).  Although free microtubules have a complex
dynamics of growth and decay, recent evidence [29] for
`superstable' cortically-bound microtubules in stages displaying unidirectional
streaming justifies the constant-length approximation.

Microtubules are the quintessential slender bodies [30] of biophysics, with aspect
ratios $\varepsilon=r/L$ of ${\cal O}(10^{-3})$.  As their self-interactions are weak, we use local slender-body theory [31, 32] to obtain the dynamics.
In an arclength parameterization $s$,
the $j^{\rm th}$ filament ${\bf r}^j(s,t)$ evolves as

\begin{equation}
\eta\left(\r_t^j - \U^j\right) = (\mathbb{I} + \r_s^j\r_s^j)\left(-A\r_{4s}^j + (\Lambda^j\r_s^j)_s + f\r^j_s\right),
\label{eom1}
\end{equation}

where $\r^j_s$ is the unit tangent, $\eta=8\pi\mu/c$, with $c=|\ln (e\varepsilon^2)|$, and the
Lagrange multiplier
$\Lambda^j$ enforcing inextensibility obeys a second-order PDE [33].
In the background flow ${\bf U}^j={\bf u}^j+{\bf u}^{i\to j} +
{\bf v}^{i\to j}$, ${\bf u}^j$ is that produced by the motors on
$j$, ${\bf u}^{i\to j}$ is due to the motors on
$i\ne j$, and ${\bf v}^{i\to j}$ is due to motion of filaments $i\ne j$.
For example, the induced flow due to the $j$th fiber is ${\bf u}(\x)=\int_0^L ds
f \r_s^j(s)\cdot {\bf G}(\x-\r^j(s))$ (see Supplemental Material [34, 35]), with ${\bf G}$ the Greens
function appropriate to the interior of a no-slip sphere [36].
Filament clamping at the sphere implies that $\r^j(0,t)$ remains fixed and that
$\r_s^j(0,t)$ is the inward sphere normal at the attachment point. The free end is torque- and force-free: $\r^j_{ss}(L,t)=\r^j_{sss}(L,t)=\Lambda(L,t)=0$.

A single fiber clamped at a flat wall displays
a supercritical
Hopf bifurcation which, expressed in terms of the dimensionless motor force $\sigma\equiv fL^3/A$, occurs
at $\sigma^*\simeq 124.2$, beyond which the
filament exhibits steady oscillations
whose amplitude grows as $\sqrt{\sigma-\sigma^*}$ [23].
When several filaments interact within the sphere (<ref>c) they
also oscillate, but with their motions synchronized
in phase, very much like eukaryotic
flagella [37].  The dynamical model (<ref>) contains
two ingredients often found necessary for such synchronization [38]: hydrodynamic interactions and
the ability of a filament to change shape and thereby
adjust its phase in response to those flows.

As the filament density and motor strength are increased we find the swirling instability: a
transition to a steady configuration of bent filaments whose distal parts are oriented almost parallel
to the wall (Fig. <ref>(d)).  The bent configuration is maintained
by azimuthal flows,
induced by the motors, that generate drag along the
distal part of the filaments, and thus a torque opposing the bending torques near the filament base.
As with any such spontaneous breaking of symmetry,
both left- and right-hand configurations are possible; the choice between the two is dictated by
initial conditions.  This transition is reminiscent of the self-organized
rotation of cytoplasmic droplets extracted from plants [39] and
the spiral vortex state of confined bacterial suspensions [40], both
modeled as suspensions of stresslets [41, 42, 43].

While direct computations on denser arrays of discrete filaments are possible [44], cortically
bound oocyte microtubules are so tightly packed, with an inter-fiber spacing
$\delta\ll L$ [10, 11, 12, 13], that a continuum approach is
justified.  The description we use [27], in which microtubules form an
anisotropic porous medium, is based on the map $\X=\r(\balpha)$,
where the Lagrangian coordinate $\balpha=(\alpha,s)$ encodes the location $\alpha$
of the minus ends of the microtubules and arclength $s$.
In a system of units made
dimensionless by the length $L$ and elastic relaxation time $\eta L^4/A$,
we obtain a continuum version of (<ref>),

\begin{equation}
\r_t - \u|_{\r(\balpha)}  = (\mathbb{I}+\hatt\hatt)\cdot\left(-\r_{ssss} + (\Lambda\hatt)_s
- \sigma\hatt\right)\,.
\label{eq:continuum:fiber}
\end{equation}

The fluid velocity
$\u$ arises from the force distribution along the filaments and is evaluated at the
Eulerian position $$ according to an inhomogeneous Stokes equation,
\begin{equation}
-\nabla^2\u + \grad p = \chi_\textnormal{mt}\rho\left[\mathcal{J}^{-1}(-\r_{ssss} + (\Lambda\r_s)_s)\right]|_{\r^{-1}(\x)}\,,
\label{eq:continuum}
\end{equation}
subject to the incompressibility constraint $·=̆ 0$. The indicator function $χ_mt$ is
supported where the MT array is present (Fig.~\ref{fig3}a).  Here, $ρ=8πρ_0 L^2/c$ is the
rescaled areal number density of microtubules, expressible as $ρ=ϕ(L/δ)^2$, where the constant
$ϕ$ depends only on the MT slenderness and packing geometry at the wall;
$ϕ≈4$ when $c≈10$ and the MTs are hexagonally packed.
The quantity $𝒥=[∂/̊∂]$ measures the change in
microtubule density due to deformations of the array; $𝒥^-1$ increases as fibers move closer together.

\begin{figure}[t]
\includegraphics[clip=true,width=\columnwidth]{Swirling_figure3.pdf}
\caption{Continuum model in planar geometry. (a) Bi-infinite array of
MTs, whose minus ends are clamped vertically at a no-slip boundary. Results of full
computations at (b) $\rho=4.65$, $\sigma=70$ and (c) $\sigma=39$. Colors
denote time, from cyan (early) to pink (late).}
\label{fig3}
\end{figure}

The simplest geometry is an infinite planar array of MTs, as shown in
Fig.~\ref{fig3}(a).  As in the discrete model, the MTs are normally
clamped to a no-slip wall and are
force- and torque-free at their plus ends. At a distance $H$ above the wall,
no-penetration and zero-tangential stress conditions are imposed on the fluid. For dynamics homogeneous
along $x$, the fluid flow is unidirectional and constant above the MTs, so $H$ plays no role.
Nonlinear computations [45] reveal both oscillatory dynamics and the emergence of steady streaming.
Fig.~\ref{fig3}(b) shows the dynamics when $ρ=4.65$ and $σ=70$: self-sustaining oscillations of the
MT array are observed, similar to those in Fig.~\ref{fig2}(c). Note that while Fig.~\ref{fig3}(b) shows only a single filament, it represents the common dynamics of
{\it all} of the collectively beating filaments in the array. When $σ$ is decreased to $≈39$, the
MT array deforms and stabilizes into a \textit{steady}
bent state (Fig.~\ref{fig3}(c)). This
represents the continuum description of the swirling
transition, with similar dynamics to those observed in the discrete results.

An equilibrium of the system occurs when filaments are aligned straight along $z$, with $=̆0$ and $Λ=-σ(1-z)$. For $σ>0$, the motor-force is
compressive and buckling may occur.
A small transverse perturbation in fiber shape of the form  $=ẑ+ϵg(z)x̂$ ($ϵ≪1$) evolves as
\begin{equation}
g_t = -g_{zzzz} -\sigma\left[(1-z)g_{z}\right]_z + \rho
\left[\sigma(1-z)g+g_{zz}\right].
\label{linstab}
\end{equation}
The first two terms are like those of an elastic filament under an aligned gravitational load, with an internal tension varying linearly
from one end to the other [46, 47]. The third is the fiber forces filtered through the
non-local Stokes operator, capturing hydrodynamic interactions within the fiber array (and hence the $ρ$ prefactor). That this term is local is both fortunate and surprising, and follows from the simplicity of the
Stokes flow in this case. The term $ρg_zz$ captures the additional resistance to bending from flow: if a MT is to bend, it must move the
fluid around it, bending other MTs; the term $ρσ(1-z)g$ is destabilizing: if a
MT is to remain straight, it must resist the fluid motions generated by MTs around it.

While the planar geometry reproduces all qualitative features of the streaming transition [34], to capture
the key feature of confined hydrodynamic interactions in the oocyte we extend the analysis to a cylindrical domain,
where the no-flow steady state is an array of straight MTs
pointing inwards. Fig.~\ref{fig4} shows the results of a linear stability analysis for an experimentally
relevant ratio of cylinder diameter to MT length of $10:1$.

\begin{figure*}[t]
\includegraphics[clip=true,width=2.0\columnwidth]{Swirling_figure4.pdf}
\caption{Continuum model in cylindrical geometry. (a) Results of linear stability analysis about the radially aligned
state, with $R=5L$.
(b) Steady-state fiber deformations and velocity field for $\sigma=150$ and $\rho=80$.
Density of visualized fibers corresponds to the physical density. Top inset shows deformed MTs and the dynamics of
a representative one (see also supplemental video [34]). Bottom inset shows the azimuthal velocity field as a function of $r$.
(c) Dimensional streaming velocities in parameter space; hatched region is consistent
with \emph{in vivo} estimates of $100-400\,$nm/s. Yellow dot denotes simulation shown in (b).}
\label{fig4}
\end{figure*}

For $ρ≪1$, the continuum model behaves like isolated fibers with negligible collective fluid entrainment.
For small $σ$, straight fiber arrays are stable (regions I $&$ II, with region II having oscillatory decay to
equilibrium), but with increasing $σ$ there is a Hopf bifurcation to a state that nonlinear simulations
show has oscillations (cf. Fig.~\ref{fig2}(c)). For $ρ≳2.8$ ($δ≲1.2L$), a new region of
instability (IV) appears, with real and positive eigenvalues; nonlinear simulations show this
leads to collective MT bending and swirling flows.

Figure \ref{fig4}(b) shows a nonlinear simulation of the transition to streaming
in region IV.  The upper inset shows the development of the instability,
with successive MTs bending over to form a dense canopy above their
highly curved bases. At steady state, the concentrated motor forces within
the canopy are azimuthally aligned, almost a $δ$-function a distance $∼L/4$ above the
wall, and drive the large-scale
streaming flow. The ooplasmic flow beneath the MT canopy is nearly a linear shear flow,
transitioning above to solid body rotation, the solution to Stokes flow forced at a
cylindrical boundary.

We now estimate ranges of density and force that are consistent with
observed streaming speeds $u≈100-400 $nm/s (Fig.~\ref{fig1} and [29, 14]). Taking $L=20 μ$m,
$μ=1 $ Pa s [11] and $A=20 $pN$μ$m$^2$, we obtain a velocity scale $A/ηL^3 ≈1 $nm/s and a force-density scale $A/L^3≈2.5 $fN/$μ$m.
Figure \ref{fig4}c shows the streaming speeds calculated from nonlinear
simulations in region IV. Those with maximum speeds falling in the experimental range
lie in the hatched area. Increasing $ρ$ only marginally increases streaming speeds, and so to increase flow speed while remaining in region IV requires increasing both $ρ$ and $σ$. The minimum value of $ρ≈20$ that is consistent with observed streaming velocities corresponds to $δ≲0.4 L$, a more stringent constraint than that required for the streaming transition. The force densities consistent with streaming speeds are $f∼ 0.1-0.6 $pN/$μ$m. Speeds on the higher end of the physical range approach the $≈ 700 $nm$/$s of kinesin-1 under negligible load [20], while cargo speeds on oocyte MTs are $200-500 $nm$/$s [14, 29, 48]. Assuming a linear force-velocity relation and a stall force of $6 $pN [20] gives a single motor force of $≈ 2 $pN; approximately $1-6$ kinesins are required per $20  μ$m MT to generate these force densities.

It may be surprising that the streaming speed only weakly depends on $ρ$. A heuristic argument for the flow speeds views the cytoskeleton as a porous medium of permeability
$k∼δ^2$, in which speed $u∼(k/μ)∇p$, where the pressure gradient (force/volume) from motors is $f/δ^2$, yielding $u∼f/μ∼(A/ηL^3)(8π/c)σ$, independent of $ρ$.  This relationship is surprisingly accurate [34].

% We now estimate ranges of density and force that are consistent with
% observed streaming speeds $u\approx 100-400\,$nm/s (Fig.~\ref{fig1} and [29, 14]). Taking $L=20\,\mu$m,
% $\mu=1\,$ Pa s [11] and $A=20\,$pN$\mu$m$^2$, we obtain a velocity scale $A/\eta L^3 \approx 1\,$nm/s and a force-density scale $A/L^3\approx 2.5\,$fN/$\mu$m.
% Figure \ref{fig4}c shows the streaming speeds calculated from nonlinear
% simulations in region IV. Those with maximum speeds falling in the experimental range
% lie in the hatched area.
% A heuristic argument for the flow speeds views the cytoskeleton as a porous medium of permeability
% $k\sim \delta^2$, in which speed $u\sim (k/\mu)\nabla p$, where the pressure gradient (force/volume)
% from motors is $f/\delta^2$, yielding $u\sim f/\mu \sim (A/\eta L^3)(8\pi/c)\sigma$,
% independent of $\rho$.  This relationship is surprisingly accurate [34];
% increasing $\rho$ only marginally increases streaming speeds, and so to
% increase flow speed while remaining in region IV requires increasing both $\rho$ and $\sigma$. The minimum
% value of $\rho\approx20$ that is consistent with observed streaming velocities corresponds to
% $\delta\lesssim 0.4 L$, a more stringent constraint than that required for the streaming transition.
% The force densities consistent with streaming speeds are $f\!\sim\! 0.1-0.6\,$pN/$\mu$m.
% Speeds on the higher end of the physical range approach the $\approx\! 700\,$nm$/$s
% of kinesin-1 under negligible load [20], while cargo speeds on oocyte MTs
% are $200-500\,$nm$/$s [14, 29, 48]. Assuming a linear force-velocity
% relation and a stall force of $6\,$pN
% [20] gives a single motor force of $\approx\! 2\,$pN; approximately $1-6$ kinesins
% are required per $20\, \mu$m MT to generate these force densities.

When the density $ρ$ is sufficiently high, the swirling instability first appears for force densities $σ$
substantially smaller than those that induce buckling instabilities in a single filament. Thus this transition
must be driven by the additional hydrodynamic destabilization that neighboring fibers impart (in the simplest geometry,
given by the term $ρσ(1-z)g$ in Eq.~\ref{linstab}). This observation motivates a simple heuristic argument
for the instability, in which a filament is bent by the flow produced primarily
by its upstream neighbor, whose distal half is nearly parallel to the wall.
Seen from a distance, that bent portion acts on the fluid like a point force [49]
$F∼(fL/2)_̊s(L)$ oriented along its distal tangent vector
(Fig.~\ref{fig5}), displaced a distance $h∼L/2$
from the surface.  Near a no-slip wall, the far-field flow
along $x$ due to a force $F∥x̂$ a distance $δ$ upstream is simple shear
[50, 51],
\begin{equation}
{\bf U}(x,z)=\dot\gamma z {\bf \hat e}_x\,,
\label{shearflow}
\end{equation}
where $γ̇=3h F/2πμδ^3$. Self-consistency requires the magnitude of the force driving the shear be
given by the projection of $F$ along $x$, so $γ̇→γ̇sin(θ(L))$.

\begin{figure}[b]
\includegraphics[clip=true,width=0.95\columnwidth]{Swirling_figure5.pdf}
\caption{Self-consistent model.  An upstream point force ${\bf F}$
parallel to the distal end of a filament produces shear flow that
deflects the filament. Two variants of the model: (i) rigid rod
with a torsional spring at its base, (ii) a clamped elastic filament.}
\label{fig5}
\end{figure}

The very simplest model to illustrate the self-consistency condition
is a rigid MT with a torsional spring at its base  that provides a restoring
torque $-kθ$ (Fig.~\ref{fig5}(i)).
With $z(s)=scosθ$ and $ηn̂n̂·U$
the local normal force on a segment, the local torque about the point $s=0$ is
$ηγ̇s^2 cos^2θ$ which, when integrated along the filament and balanced
against the spring torque, yields the
self-consistency condition
\begin{equation}
\theta=B\sin\theta\cos^2\!\theta\,,
\end{equation}
where $B=ηγ̇L^3/k$.
For $B<1$ (slow flow or a stiff spring) $θ=0$ is the only fixed point, while for
$B≳1$ two mirror-image swirling solutions appear through a pitchfork bifurcation,

To study the interplay between filament oscillations and swirling we use \eqref{shearflow}
in the filament dynamics \eqref{eom1}, where the control parameter for the shear flow is [25, 26]
\begin{equation}
M=\frac{\eta\dot\gamma L^4}{A}\sim \frac{3\sigma}{c}\left(\frac{\rho}{\phi}\right)^{3/2},
\label{Mparameter}
\end{equation}
where the second relation uses the above estimates for $F$ and $h$.
Since a clamped elastic filament
behaves like a torsional spring with a spring constant $k=A/L$, we see consistency with the
parameter $B$ defined above.
A numerical self-consistent calculation confirms the existence of a swirling instability [34].

Through simplified discrete and continuum models we have demonstrated the existence of
a novel swirling instability of arrays of elastic filaments, thus lending support to the hypothesis
[14]
that cytoplasmic streaming flows in \textit{Drosophila}
oocytes are tied to self-organization of the microtubule cytoskeleton.  Future studies could shed light on the detailed mechanism involved in the untangling of the
\textit{Drosophila} oocyte cytoskeleton when it transitions to the vortical state,
and the possibility of reproducing this transition \textit{in vitro}. Lastly, this study highlights the role of active force dipoles in the self-organization of fluid-biopolymer systems [41, 42, 43].

\begin{acknowledgments}
We are indebted to Maik Drechsler and Isabel Palacios for sharing the data in Fig.~\ref{fig1} and to them, Daniel St Johnston and Vladimir Gelfand for discussions on \textit{Drosophila} streaming.
This work was supported in part by ERC Consolidator grant 682754 (EL),
Wellcome Trust Investigator Award 207510/Z/17/Z,
Established Career Fellowship EP/M017982/1 from the
Engineering and Physical Sciences Research Council, and the Schlumberger Chair Fund (REG). MJS acknowledges the support of NSF Grant DMS-1620331.

\end{acknowledgments}
\begin{thebibliography}{99}

[1]
D. Needleman and M.J. Shelley, The stormy fluid dynamics of the living cell,
\href{https://doi.org/10.1063/PT.3.4292}{Physics Today {\bf 72}, 32 (2019)}.

[2] A short primer on the subject is:
R. Bastock and D. St. Johnston, \textit{Drosophila} oogenesis,
\href{https://doi.org/10.1016/j.cub.2008.09.011}{Curr. Biol. {\bf 18}, R1082 (2008)}.

[3] L. He, X. Wang, and D.J. Montell, Shining light on {\it Drosophila} oogenesis:
live imaging of egg development,
\href{https://doi.org/10.1016/j.gde.2011.08.011}{Curr. Op. Gen. \& Dev. {\bf 21}, 612 (2011)}.

[4] The phenomenon of streaming was first discovered in plants, as reported by
B. Corti, {\it Osservazione Microscopische sulla Tremella e sulla Circulazione del Fluido
in Una Planto Acquaguola} (Appresso Giuseppe Rocchi, Lucca, Italy, 1774).

[5]
R.E. Goldstein and J.-W. van de Meent,
A Physical Perspective on Cytoplasmic Streaming,
\href{https://doi.org/10.1098/rsfs.2015.0030}{Interface Focus {\bf 5}, 20150030 (2015)}.

[6] H. Gutzeit and R. Koppa, Time-lapse film analysis of cytoplasmic streaming during
late oogenesis of {\it Drosophila},
\href{https://doi.org/}{J. Embryol. Exp. Morphol. {\bf 67}, 101 (1982)}.

[7] W. Theurkauf, S. Smiley, M. Wong, and B. Alberts, Reorganization of the cytoskeleton
during {\it Drosophila} oogenesis: implications for axis specification and intercellular
\href{https://doi/org/}{Development {\bf 115}, 923 (1992)}.

[8] W.E. Theurkauf, Premature microtubule-dependent cytoplasmic streaming in
cappuccino and spire mutant oocytes,
\href{https://doi/org/10.1126/science.8091233}{Science {\bf 265} 2093 (1994)}.

[9] I.M. Palacios and D. St. Johnston, {\it Kinesin light chain}-independent
function of the {\it Kinesin Heavy Chain} in cytoplasmic streaming and posterior localisation in
the {\it Drosophila} oocyte,
\href{https://doi/org/10.1242/dev.00119}{Development 129:5473–5485 (2002)}.

[10] L.R. Serbus, B.J. Cha, W.E. Theurkauf, W.M. Saxton, Dynein and the actin
cytoskeleton control kinesin-driven cytoplasmic streaming in {\it Drosophila} oocytes,
\href{https://doi/org/10.1242/dev.01956}{ Development {\bf 132}, 3743 (2005)}.

[11] S. Ganguly, L.S. Williams, I.M. Palacios, and R.E. Goldstein,
Cytoplasmic Streaming in \textit{Drosophila} Oocytes Varies with Kinesin Activity
and Correlates With the Microtubule Cytoskeleton Architecture,
\href{https://doi.org/10.1073/pnas.1203575109}{Proc. Natl. Acad. Sci. USA {\bf 109}, 15109 (2012)}.

[12] M. Drechsler, F. Giavazzi, R. Cerbino, and I.M. Palacios,
Active diffusion and advection in \textit{Drosophila} oocytes results from the interplay
of actin and microtubules,
\href{https://doi.org/10.1038/s41467-017-01414-6}{Nat. Comm. {\bf 8}, 1520 (2017)}.

[13] M. Drechsler, L.F. Lang, L. Al-Khatib, H. Dirks, M. Burger, C.-B. Sch{\"o}nlieb and
I.M. Palacios, Optical flow analysis reveals that Kinesin-mediated advection impacts the orientation of
microtubules in the {\it Drosophila} oocyte,
\href{https://doi.org/10.1091/mbc.E19-08-0440}{Mol. Biol. Cell {\bf 31}, 1246 (2020)}.

[14] C.E. Monteith, M.E. Brunner, I. Djagaeva, A.M. Bielecki, J.M. Deutsch,
and W.M. Saxton, A Mechanism for Cytoplasmic Streaming: Kinesin-Driven Alignment of
Microtubules and Fast Fluid Flows,
\href{https://doi.org/10.1016/j.bpj.2016.03.036}{Biophys. J. {\bf 110}, 2053 (2016)}.

[15]
P. Khuc Trong, J. Guck, and R.E. Goldstein,
Coupling of Active Motion and Advection Shapes Intracellular Cargo Transport,
\href{https://doi.org/10.1103/PhysRevLett.109.028104}{Phys. Rev. Lett. {\bf 109}, 028104 (2012)}.

[16] L.S. Williams, S. Ganguly, P. Loiseau, B.F. Ng and I.M. Palacios,
The auto-inhibitory domain and the ATP-independent microtubule-binding
region of Kinesin Heavy Chain are major functional domains for
transport in the \textit{Drosophila} germline,
\href{https://doi.org/10.1242/dev.097592}{Development {\bf 141}, 176 (2014)}.

[17]
P. Khuc Trong, H. Doerflinger, J. Dunkel, D. St. Johnston, and R.E. Goldstein,
Cortical Microtubule Nucleation Can Organise the Cytoskeleton of
\textit{Drosophila} Oocytes to Define the Anteroposterior Axis,
\href{https://doi.org/10.7554/eLife.06088}{eLife {\bf 4}, e06088 (2015)}.

[18] W. Lu, M. Lakonishok, A.S. Serpenskaya, D. Kirchenb{\"u}echler,
S.-C. Ling
and V.I. Gelfand, Ooplasmic flow cooperates with transport and anchorage in {\it Drosophila} oocyte
posterior determination,
\href{https://doi.org/10.1083/jcb.201709174}{J. Cell. Biol. {\bf 217}, 3497 (2018)}.

[19] F. Gittes, B. Mickey, J. Nettleton and J. Howard, Flexural Rigidity of
Microtubules and Actin Filaments Measured from Thermal Fluctuations in Shape,
\href{https://doi.org/10.1083/jcb.120.4.923}{J. Cell Biol. {\bf 120}, 923 (1993)}.

[20] K. Visscher, M.J. Schnitzer, and S.M. Block, Single Kinesin Molecules
Studied With a Molecular Force Clamp,
\href{https://doi.org/10.1038/22146}{Nature {\bf 400}, 184 (1999)}.

[21] G. Herrmann and R.W. Bungay RW, On the stability of
elastic systems subjected to nonconservative forces,
\href{https://doi.org/10.1115/1.3629660}{J. Appl. Mech. {\bf 31}, 435 (1964)}.

[22] P.V. Bayly and S.K. Dutcher, Steady dynein forces induce flutter
instability and propagating waves in mathematical models of flagella,
\href{https://doi.org/10.1098/rsif.2016.0523}{J. R. Soc. Interface {\bf 13}, 20160523 (2016)}.

[23] G. De Canio, E. Lauga, and R.E. Goldstein, Spontaneous oscillations of elastic
filaments induced by molecular motors,
\href{https://doi.org/10.1098/rsif.2017.0491}{J. R. Soc. Interface {\bf 14}, 20170491 (2017)}.

[24] F. Ling, H. Guo, and E. Kanso, Instability-driven oscillations of elastic microfilaments,
\href{https://doi.org/10.1098/rsif.2018.0594}{J. R. Soc. Interface {\bf 15}, 149 (2018).}

[25] Y.-N. Young and M.J. Shelley, Stretch-Coil Transition and
Transport of Fibers in Cellular Flows,
\href{https://doi.org/10.1103/PhysRevLett.99.058303}{Phys. Rev. Lett. {\bf 99}, 058303

[26] V. Kantsler and R.E. Goldstein, Flucutations, Dynamics, and the
Stretch-Coil Transition of Single Actin Filaments in Extensional Flows,
\href{https://doi.org/10.1103/PhysRevLett.108.038103}{Phys. Rev. Lett.
{\bf 108}, 038103 (2012)}.

[27]
D.B. Stein and M.J. Shelley,
Coarse graining the dynamics of immersed and driven fiber assemblies,
\href{https://doi.org/10.1103/PhysRevFluids.4.073302}{Phys. Rev. Fluids {\bf 4}, 073302 (2019)}.

[28] I.M. Palacios and M. Drechsler, private communication (2020), based on methods detailed
earlier [11, 16].

[29] W. Lu, M. Winding, M. Lakonishok, J. Wildonger and V.I. Gelfand, Microtuble-microtubule
sliding by kinesin-1 is essential for normal cytoplasmic streaming in {\it Drosophila} oocytes,
\href{https://doi.org/10.1073/pnas.1522424113}{Proc. Natl. Acad. Sci. USA {\bf 113}, E4995 (2016)}.

[30] J.B. Keller and S.I. Rubinow, Slender-body theory for slow viscous flow,
\href{https://doi.org/10.1017/S0022112076000475}{J. Fluid Mech. {\bf 75}, 705 (1976)}.

[31] J. Gray and G.J. Hancock, The propulsion of sea-urchin spermatozoa,
J. Exp. Biol.{\bf 32}, 802 (1955).

[32] A.K. Tornberg, and M.J. Shelley, Simulating the dynamics and interactions of flexible fibers in Stokes flows, \href{https://doi.org/10.1016/j.jcp.2003.10.017}{J. Comput. Phys. {\bf 196}, 1 (2004)}.

[33] R.E. Goldstein and S.A. Langer, Nonlinear dynamics of stiff polymers,
\href{https://doi.org/10.1103/PhysRevLett.75.1094}{Phys. Rev. Lett. {\bf 75}, 1094 (1995)}.

[34]  See Supplemental Material at \url{http://link.aps.org/supplemental/xxx} for
further details and results.

[35] Further examples and details are in:
G. De Canio, Motion of filaments induced by molecular motors:
from individual to collective dynamics, PhD thesis, University of Cambridge (2018).

[36] C. Maul and S. Kim, Image systems for a Stokeslet inside a rigid
spherical container,
\href{https://doi.org/10.1063/1.868223}{Phys. Fluids {\bf 6}, 2221 (1994)}.

[37] R.E. Goldstein, Green Algae as Model Organisms for Biological Fluid Dynamics,
\href{https://doi.org/10.1146/annurev-fluid-010313-141426}{Annu. Rev. Fluid Mech. {\bf 47}, 343 (2015)}.

[38] T. Niedermayer, B. Eckhardt and P. Lenz, Synchronization, phase locking, and metachronal
wave formation in ciliary chains,
\href{https://doi.org/10.1063/1.2956984}{Chaos {\bf 18}, 037128 (2008)}.

[39] Y. Yotsuyanagi, Recherches sur les ph{\'e}énomen{\`e}s moteurs dans les fragments de
protoplasme isol{\'e}s. I. Mouvement rotatoire et le processus de son apparition,
\href{https://doi.org/10.1508/cytologia.18.146}{Cytologia {\bf 18}, 146 (1953)};
Recherches sur les ph{\'e}énomen{\`e}s moteurs dans les fragments de
protoplasme isol{\'e}s. II. Mouvements divers d{'e}termin{'e}s par la condition de milieu,
\href{https://doi.org/10.1508/cytologia.18.202}{Cytologia {\bf 18}, 202 (1953)}.

[40] H. Wioland, F.G. Woodhouse, J. Dunkel, J.O. Kessler, and R.E. Goldstein,
Confinement Stabilizes a Bacterial Suspensions into a Spiral Vortex,
\href{https://doi.org/10.1103/PhysRevLett.110.268102}{Phys. Rev. Lett. {\bf 110}, 268102 (2013)}.

[41] D. Saintillan and M.J. Shelley, Instabilities and Pattern Formation in Active Particle Suspensions:
Kinetic Theory and Continuum Simulations,
\href{https://doi.org/10.1103/PhysRevLett.100.178103}{Phys. Rev. Lett. 100, 178103 (2008).}

[42] D. Saintillan, M. J. Shelley, A. Zidovska, Extensile motor activity drives coherent motions in a model of interphase chromatin, \href{https://doi.org/10.1073/pnas.1807073115}{Proc. Natl. Acad. Sci. USA {\bf 115}, 11442 (2018).}

[43] F.G. Woodhouse and R.E. Goldstein, Spontaneous Circulation of
Confined Active Suspensions,
\href{https://doi.org/10.1103/PhysRevLett.109.168105}{Phys. Rev. Lett. {\bf 109}, 168105 (2012)}.

[44] E. Nazockdast, A. Rahimian, D. Zorin, and M.J. Shelley, A fast platform for simulating semi-flexible fiber suspensions applied to cell mechanics, \href{https://doi.org/10.1016/j.jcp.2016.10.026}{J. Comput. Phys. {\bf 329} (2017)}.

[45] Differentiation of \eqref{eq:continuum:fiber} yields an equivalent equation
for the tangent-vector field $\hatt$ [27], which we have found to be
numerically more stable.

[46] L.D. Landau and E.M. Lifshitz, \textit{Theory of Elasticity},
2nd ed. (Pergamon Press, Oxford, 1970), p. 99, Problem 7.

[47] R.E. Goldstein, P.B. Warren and R.C. Ball, Shape of a Ponytail and the
Statistical Physics of Hair Fiber Bundles,
\href{https://doi.org/10.1103/PhysRevLett.108.078101}{Phys. Rev. Lett. {\bf 108}, 078101 (2012)}.

[48] P. Loiseau,, R. Davies, L.S. Williams, M. Mishima and I.M. Palacios,
\textit{Drosophila} PAT1 is required for Kinesin-1 to transport cargo and to maximize its
motility, \href{https://doi.org/10.1242/dev.048108}{Development {\bf 137}, 2763 (2010)}.

[49] D.R. Brumley, K.Y. Wan, M. Polin and R.E. Goldstein, Flagellar synchronization through
direct hydrodynamic interactions,
\href{https://doi.org/10.7554/eLife.02750}{eLife {\bf 3}, e02750 (2014)}.

[50] J.R. Blake, A note on the image system for a stokeslet in a no-slip boundary,
\href{https://doi.org/10.1017/S0305004100049902}{Math. Proc. Camb. Phil. Soc. {\bf 70}, 303 (1971)}.

[51] J.-B. Thomazo, E. Lauga, B. Le R{\'e}v{\'e}rend, E. Wandersman,
and A.M. Prevost, Collective stiffening of soft hair assembles,
\href{https://www.arxiv.org/abs/2002.02834}{arXiv:2002.02834}.

\end{thebibliography}

\vfil
\eject

\section{Supplemental Material}

This file contains calculational details for computations of discrete filaments
and the continuum model, further results from linear stability analyses in the
latter, and quantitative analysis of streaming speeds in the swirling state.

\setcounter{equation}{0}
\setcounter{figure}{0}
\setcounter{table}{0}
\setcounter{page}{1}
\makeatletter
\renewcommand{\theequation}{S\arabic{equation}}
\renewcommand{\thefigure}{S\arabic{figure}}
\renewcommand{\bibnumfmt}[1]{[S#1]}
\renewcommand{\citenumfont}[1]{S#1}
%%%%%%%%%% Prefix a "S" to all equations, figures, tables and reset the counter %%%%%%%%%%

\section{Calculational Details for Discrete Filaments}

\subsection{Image system inside a sphere}

By modelling the motor-cargo ensemble with fluid-entraining
follower forces and
approximating the system geometry by a spherical container,
we take advantage of the analytical solution [52] for the velocity field
of a point force inside a sphere with a no-slip wall.
These simplifications allow us to compute the
background flow directly,
instead of numerically solving the
Stokes equations for the fluid flow.
For a point force of magnitude $F$ located at $_0$ inside a sphere of unit radius,
the m$^th$ component of the velocity is
\begin{equation}
u_{\text m} = \frac{F_{\text j}}{8 \pi \mu} \left[ G_{\text{jm}}
\left(\x - \x_0 \right)  +
\bar {G}_{\text{jm}} \left(\x \right) \right] \, ,
\label{vel_sphere}
\end{equation}
where $μ$ is the dynamic viscosity
of the fluid, $G_jm$ is the Green's function
\begin{equation}
G_{\text{jm}} \left( \x-\x_0 \right)  =
\frac{\delta_{\text{jm}}}{r} + \frac{( x_{\text j} - x_{0_{\text j}} )
\left(  x_{\text m} - x_{0_{\text m}} \right)}{r^3} \, ,
\end{equation}
\begin{widetext}
\begin{align}
\bar {G}_{\text{jm}} \left(\x\right)   =  & \ \hat e_{\text m} \hat e_{\text k}
\Bigg\{ \frac{1-3R^2_0}{2R^3_0} \  G_{\text{jk}} \left(\x - \bar\x_0 \right)
-\frac{1-R^2_0}{R^4_0} \ \hat e_{\text l} \   G_{\text{jk},l} \left(\x - \x_0 \right)
-\frac{(1-R^2_0)^2}{4R^5_0} \ \nabla^2  G_{\text{jk}} \left(\x - \bar\x_0 \right) \Bigg\} \notag \\
& + \left( \delta_{\text{km}} -\hat e_{\text k} \hat e_{\text m} \right)
\Bigg\{  \frac{3R^2_0-5}{2R^3_0} \  G_{\text{jk}} \left(\x - \bar\x_0 \right)
+\frac{(1-R^2_0)^2}{4R^5_0} \ \nabla^2  G_{\text{jk}} \left( \x - \bar\x_0 \right) \Bigg\} \notag \\
& + \frac{1-R^2_0}{R^4_0} \ \hat e_{\text k}  \left( \delta_{\text{lm}} -\hat e_{\text l} \hat e_{\text m} \right)
G_{\text{jk},l} \left(\x - \bar\x_0 \right)
-\frac{3(R^2_0-1)}{R^3_0}  \frac{\left( \delta_{\text{jm}} -\hat e_{\text j} \hat e_{\text m} \right)}{\bar r}
+\left(R^2 - 1 \right) \left( \delta_{\text{km}} -\hat e_{\text k} \hat e_{\text m} \right)
\frac{\partial \varphi_{\text k}}{\partial x_{\text j}} \, ,
\end{align}
\end{widetext}
\begin{equation}
\varphi_{\text k} = - \frac{3(R^2_0 - 1)}{2R^3_0}\ \frac{x_{\text k}}{\bar r} \ \frac{R-\bar R_0 \cos \alpha +
\bar r \cos \alpha}{R \bar R^2_0 \sin^2\alpha} \, .
\end{equation}
Here, $R_0$ and $R̅_0$ are the norms of $_0$ and its image
$_0$, respectively, $ = _0/R_0$ is the unit vector in the axial direction,
$r = |-_0|$, $r̅ = |-_0|$,
$R$ is the norm of $$, and $α$ is the
angle between $$ and $$.

\subsection{Geometrical formulation}

For both the discrete-filament computations and the self-consistent model, we use the tangent angle representation of the
filament.  For a filament of length $L$, parameterized by arclength $s$,
anchored at wall in the
$x-y$ plane, with $z$ orthogonal to the plane, into the fluid, the position vector is $r(s)=(x(s),z(s))$,
and the unit tangent and normal vectors $t̂(s)=r_s$ and $n̂(s)$
\begin{eqnarray}
{\bf \hat t}(s)&=&\cos\theta(s){\bf \hat e}_x + \sin\theta(s){\bf \hat e}_z\\
{\bf \hat n}(s)&=&\sin\theta(s){\bf \hat e}_x -\cos\theta(s){\bf \hat e}_z~,
\end{eqnarray}
where $θ(s)$ is the angle the tangent vector makes with respect to the $x$-axis.
The curve is obtained by integration of $t̂$, and assuming $(x(0),z(0))=(0,0)$, we have
\begin{equation}
x(s)=\int_0^s \! ds' \cos\theta(s')~, \ \ \ z(s)=\int_0^s\! ds' \sin\theta(s')~.
\end{equation}
It follows from the above that
\begin{equation}
\frac{\partial{\bf \hat t}}{\partial t}=-\theta_t{\bf \hat n}~,
\end{equation}
from which we can obtain the equation of motion for $θ$ given that for $r$.

The equation of motion is that of elastohydrodynamics within resistive force theory,
\begin{equation}
\left(\zeta_{\parallel} {\bf \hat t}{\bf \hat t}+\zeta_{\perp}{\bf \hat n}{\bf \hat n}\right)
\cdot\left({\bf r}_t-{\bf U}\right)=-A{\bf r}_{ssss}-\left(\Lambda{\bf r}_s\right)_s~,
\label{ehd}
\end{equation}
where $A$ is the filament bending modulus, $Λ$ is the Lagrange multiplier, and
$U=Uê_x+Vê_z$ is the background flow experienced by the filament.
We express the derivatives in \eqref{ehd} in terms of the curvature $κ=θ_s$, using the
Frenet-Serret equations $t̂_s=-κn̂$ and $n̂_s=κt̂$,
and rescale the equation of motion via
\begin{align}
{\bf r}&=L{\bf r}'~, \ \ \ z=Lz'~, \ \ \ s=Ls'~, \ \ \ \kappa = \frac{1}{L} \kappa' ~, \nonumber \\
\Lambda &=\frac{A}{L^2} \Lambda' ~,  \ \ \ t= \frac{\zeta_\perp L^4}{A} t' ~, \ \ \ U= \frac{\zeta_\perp L^3}{A} U'.
\end{align}
If we let $ζ_⊥/ζ_∥≃2$, and use
$ê_x=cosθt̂ + sinθn̂$ and
$ê_z=sinθt̂ - cosθn̂$, and then drop the primes, we obtain
\begin{widetext}
\begin{equation}
{\bf r}_t=\left(U\cos\theta +V\sin\theta\right){\bf \hat t}+\left(U\sin\theta-V\cos\theta\right) {\bf \hat n}
+\left(\kappa_{ss}-\kappa^3+\kappa\Lambda\right){\bf \hat n}
+\eta\left(3\kappa\kappa_s-\Lambda_s\right){\bf \hat t}~.
\label{eom2}
\end{equation}
Differentiating \eqref{eom2} with respect to $s$ to obtain the equation of motion of
the tangent angle $\theta$ and that for $\Lambda$:
\begin{equation}
\theta_t=-\theta_{ssss}+3\theta_s^2\theta_{ss}-\left(\theta_s\Lambda\right)_s
\label{final1}
\end{equation}
while enforcing inextensibility leads to the equation for $\Lambda$,
\begin{equation}
\left(\partial_{ss}-\frac{1}{2}\theta_s^2\right)\Lambda=
\label{final2}
\end{equation}
\end{widetext}
The boundary conditions at the attachment point are those of a clamped filament,
\begin{align}
&\theta(0)=\frac{\pi}{2}~, \ \ \  \theta_{sss}(0)-\theta_s(0)^3+
\theta_s(0)\Lambda(0)=0~, \\ \nonumber
\end{align}
whereas at the free end we have
\begin{equation}
\theta_s(1)=0~, \ \ \ \theta_{ss}(1)=0~,\ \ \ \Lambda(1)=0~. \ \ \
\end{equation}

\subsection{Self-consistent model}

If we take the typical bent filament shape in the streaming regime (Fig. 2(d) of main text) and place
it near a flat no-slip wall, we can understand in the simplest situation the flow it produces.  As shown
in the streamlines and velocity colormap of Figure \ref{flow_near_filament} there is a boundary-layer
phenomenon near the wall that arises from the combination of tangential forcing of the flow by the
motors and the no-slip condition at the wall.  The most prominent part of the flow is a lobe of high speeds
emanating from the bent portion of the filament, directed downstream.  This feature forms the basis
of the self-consistent model for the swirling transition.

\begin{figure}[b]
\includegraphics[width=0.98\columnwidth]{Swirling_figureS1.pdf}
\caption{Flow field near a filament in the swirling regime.  A filament in the bent shape found in the
spherical geometry of the main text was placed near a flat, no-slip wall and the resulting flow field
computed numerically.  Note the shear flow downstream.}
\label{flow_near_filament}
\end{figure}

As Fig. \ref{flow_near_filament} shows, the flow downstream from a bent filament is approximately
simple shear.  This can be seen directly in Blake's analysis
[53] of the flow due to a point force near a no-slip wall.  For the purposes of a
simple self-consistent model, we consider only the asymptotic form of that flow evaluated at
$x=(x,y,z)$ due to a point force $F$ at $(0,0,h)$:
\begin{align}
u_i\simeq&\frac{F_k}{8\pi\mu}\biggl[\frac{12hx_ix_{\alpha}x_3\delta_{k\alpha}}{\vert{\bf x}\vert^5}\\ \nonumber
& +h^2\delta_{k3}\left(-\frac{(12+6\delta_{i3})x_ix_3}{\vert{\bf x}\vert^5}+
\frac{30x_ix_3^3}{\vert{\bf x}\vert^7}\right)\biggr]~.
\end{align}
Taking the leading-order term, and $F$ along $x$, we obtain
\begin{equation}
{\bf u}(x,z)=\dot\gamma(x) z {\bf \hat e}_x
\end{equation}
where the shear rate is
\begin{equation}
\dot\gamma(x)=\frac{3Fh}{2\pi\mu x^3}~.
\end{equation}
In the notation used in the tangent angle dynamics \eqref{final1}, we have $U=Mz(s)$ and $V=0$, where
\begin{equation}
M=\frac{\eta\dot\gamma L^4}{A}\sim \frac{3\sigma}{c}\left(\frac{\rho}{\phi}\right)^{3/2}~.
\end{equation}
a control parameter that appears also in the dynamics elastic filaments in extensional flows
[54, 55].  Here, the second relation uses the definitions of the dimensionless
control parameters $σ$ and $ρ$ from the main text.
To complete the self-consistency calculation, we set the $x$-component of the force $F$ to be
proportional to the projection of the tangent vector at the distal tip, so
$M →Msin(θ(L))$

\begin{figure}[t]
\includegraphics[width=0.98\columnwidth]{Swirling_figureS2.pdf}
\caption{Results of the self-consistent calculation. Relaxation of a filament to steady bent configuration,
at $\sigma=45$ and $\rho=5$, with colors interpolating between highly bent initial condition (cyan) to
final state (magenta) over a (dimensionless) time of $0.05$.}
\label{filament_relaxation_sc}
\end{figure}

The tangent angle dynamics \eqref{final1} and \eqref{final2}
were solved by a Crank-Nicholson method
via a pentadiagonal matrix for the fourth derivative,
treating all the nonlinearities
explicitly, and discretizing the $Λ$ problem as a tridiagonal matrix in which
the diagonal term involving $θ_s^2$ was updated after each time step.  The initial
condition was $θ(0,s)=π/2-asin(πs/2)^2$, with $a=2.5$ and the time step
$dt=0.1 ds^4$, with $ds=1/(N-1)$, where $N$ is the total number of
grid points ($N=41$ is sufficient).

Figure \ref{filament_relaxation_sc} shows a typical result of the self-consistent
calculation in the regime of stable streaming, illustrating how a highly bent initial condition relaxes
to a conformation like that seen in the calculations described in the main text (compare Figs. 3(c) from
discrete filament calculations and
4(b) from continuum model).

% A coarse parameter sweep results in the
% numerical stability diagram shown in Fig. \ref{stability_diagram_sc}, which displays the same
% features as that from the continuum model, including a steep vertical asymptote at
% $\sigma\sim 25$, a tongue of stable streaming at intermediate values, and collective oscillations
% at higher $\sigma$.  Given the simplifying assumptions made in the model \textemdash replacement of
% a fraction of the filament by a single stokeslet, approximation of the flow field as pure shear, etc.
% \textemdash this level of agreement lends support to the heuristic picture of the swirling instability
% embodied in the self-consistent model.

%	\begin{center}
%	\includegraphics[width=0.98\columnwidth]{Swirling_figureS3.pdf}
%	\caption{Stability diagram (placeholder for final version).}
%	\end{center}
%	\label{stability_diagram_sc}
% \end{figure}

\section{Further Results from Continuum Model}

\subsection{Dynamics of the transition to swirling}

In Supplemental Video 1, we show the dynamics of the transition from the unstable equilibrium of radially oriented filaments to stable streaming flow, corresponding to Fig. 4b (main text). In this parameter regime, the filament array initially buckles, and large azimuthal flows are generated by rapid deformation of the array. As the filament motion slows, the filaments become oriented parallel to the cylinder walls. Azimuthal flows are continuously generated by the streaming mechanism, and the drag from these flows moving past the stationary filaments is sufficient to pin them into the observed conformation.

\subsection{Effect of cylindrical confinement}

\begin{figure}[b]
\includegraphics[width=0.9\columnwidth]{Swirling_figureS4.pdf}
\caption{Results of the stability analysis in the planar geometry. Regions are colored and coded as in the main text.}
\label{flat_phase}
\end{figure}

\begin{figure*}[t]
\includegraphics[width=1.8\columnwidth]{Swirling_figureS5.pdf}
\caption{Results of the stability analysis in the confined geometry. The confinement ratios are 2.2:1 (Panel a), 4:1 (Panel b), and 20:1 (Panel c). See Fig.~\ref{flat_phase} for a description of the different regions.}
\label{cylindrical_phase}
\end{figure*}

In the main text, we present the evolution equation for a small perturbation to the equilibrium solution in a planar geometry. This equation:
\begin{equation}
g_t = -g_{zzzz} -\sigma\left[(1-z)g_{z}\right]_z + \rho
\left[\sigma(1-z)g+g_{zz}\right],
\label{linstabSM}
\end{equation}
remains local due to the simplicity of the Stokes flow, allowing the individual terms to be easily identified and interpreted. In order to better capture the real geometry of the oocyte, in the main text we present the results of the linear stability analysis in a cylindrical geometry (Fig 4a, main text). The same analysis, over the same range of effective densities ($ρ$) and effective motor force-densities ($σ$) is shown in Fig. \ref{flat_phase}, computed for the planar geometry. These results are qualitatively the same as those shown in the main text in Fig 4a, indicating that the essential structure of the swirling transition is due to the interactions between nearby filaments in the fiber bed, and does not depend in a critical way on the specific geometry of the oocyte.

In, Fig. \ref{cylindrical_phase} we present the same analysis for different confinement ratios --- 2.2:1 (Fig \ref{cylindrical_phase}a), 4:1 (Fig \ref{cylindrical_phase}b), and 20:1, (Fig \ref{cylindrical_phase}c). Note that the ratio 2.2:1 indicates a cylinder diameter of $2.2L$, which is extremely confined, with filaments reaching $≈91%$ of the way to the center of the cylinder. Even under such strong confinement, the qualitative structure remains, although the quantitative value at which the bifurcations occur changes significantly. The effect of strong confinement is to suppress the appearance of streaming at larger forcing values. This makes sense: when forcing is large, the tendency of lone fibers is to oscillate. In the planar system, when fibers are sufficiently dense, the strong flows caused by motors walking on neighboring fibers is sufficient to pin the fiber into a stable deformed state, suppressing oscillations (in fact, the drag from these flows produces extensile forces on the fibers, reducing the net contractile forcing generated by the motors). Under high confinement, the strength of such flows near the fiber tips is reduced, and the flow can no longer suppress the tendency of the fiber to oscillate. When the confinement is not so extreme, even the quantitative values change little, from 10:1 (Fig 4a, main text), to 20:1 (Fig \ref{cylindrical_phase}c), to the flat geometry (Fig \ref{flat_phase}).

It may surprise the reader that when $ρ→0$, the behavior (as a function of $σ$) at all confinement ratios converges to \emph{exactly} the same thing. This is because the continuum model imposes boundary conditions on the coarse-grained flow only; as $ρ→0$ this flow goes to $0$ with it and the dynamics are those of a fiber moving in a quiescent background.

\begin{figure}[b]
\includegraphics[width=0.9\columnwidth]{Swirling_figureS6.pdf}
\caption{Dependence of streaming speed on the force density ($\sigma$) generated by kinesin-1 motors. When $\rho$ is large, the streaming speed approaches the estimate $v_\textnormal{streaming}=2\sigma$.}
\label{streaming_speeds}
\end{figure}

\subsection{Streaming Speeds in Dense Arrays}

By solving for $-_̊ssss + (Λ_̊s)_s$ in Eq. 2 (main text), substituting into Eq. 3 (main text), defining $ξ≡χ_MTρ𝒥^-1$, and ignoring frame transformations for the sake of simplicity, we find that:
\begin{equation}
-\grad^2\u + \grad p + \xi(\mathbb{I}-\r_s\r_s/2)(\u-\r_t) = \xi\sigma\r_s.
\label{eq:simple_balance}
\end{equation}
This is a forced Brinkman equation with an anisotropic permeability. When $ξ$ is large, the skeletal drag term dominates the left-hand side, and at steady state $_̊t=0$. We may thus approximate the steady streaming velocity as simply $≈̆σ(𝕀+_̊s_̊s)_̊s = 2σ_̊s$. For homogeneous flows in the planar geometry, and axisymmetric flows in the cylindrical geometry, the flow is purely in the $x̂$ or $θ̂$ (azimuthal) direction, respectively. Since $_̊s$ is not purely coincident with $x̂$ or $θ̂$, some of the forcing generated by the molecular motors must be absorbed into the pressure gradient, and so we take this estimate to be an upper bound: $v_streaming≤2σ$. The maximum steady-state streaming speed as a function of $σ$, as computed by the continuum model in 10:1 confinement, is shown in Fig.~\ref{streaming_speeds}, for $ρ=10, 20, …, 50$. For smaller values of $ρ$, the fiber is less deformed and drag contributes less to the balance in Eq.~\ref{eq:simple_balance}, leading to a relationship closer to $v_streaming=σ$. When $ρ$ is larger, drag dominates Eq.~\ref{eq:simple_balance}, the fiber is nearly azimuthally aligned, and the streaming speed approaches the bound $v_streaming=2σ$.

In dimensional units, this bound is $v_streaming≤2σA/ηL^3$. In the main text, we used a simple argument based on Darcy's law to estimate $v_streaming=(σA/ηL^3)(8π/c)$. The leading factor of two in the more refined estimate presented in this SI comes from taking into account fiber anisotropy and the nearly azimuthal alignment of the fibers, while the disappearance of the factor of $8π/c$ arises from properly accounting for the fibers aspect ratio. In either case, the disappearance of the MT density from the estimate of streaming velocity arises due to the fact that both the forcing ($ξσ_̊s$) and the skeletal drag ($ξ(𝕀-_̊s_̊s/2))$ scale with the density $ρ$.

\begin{thebibliography}{99}

[52] C. Maul and S. Kim, Image systems for a Stokeslet inside a rigid
spherical container,
\href{https://doi.org/10.1063/1.868223}{Phys. Fluids {\bf 6}, 2221 (1994)}.

[53] J.R. Blake, A note on the image system for a stokeslet in a no-slip boundary,
\href{https://doi.org/10.1017/S0305004100049902}{Math. Proc. Camb. Phil. Soc. {\bf 70}, 303 (1971)}.

[54] Y.-N. Young and M.J. Shelley, Stretch-Coil Transition and
Transport of Fibers in Cellular Flows,
\href{https://doi.org/10.1103/PhysRevLett.99.058303}{Phys. Rev. Lett. {\bf 99}, 058303

[55] V. Kantsler and R.E. Goldstein, Flucutations, Dynamics, and the
Stretch-Coil Transition of Single Actin Filaments in Extensional Flows,
\href{https://doi.org/10.1103/PhysRevLett.108.038103}{Phys. Rev. Lett.
{\bf 108}, 038103 (2012)}.

\end{thebibliography}
\end{document}