1. Introduction
The dynamics of slender fibres in viscous flow is key to our understanding of many complex phenomena encountered in diverse areas ranging from biology to physics and engineering (Du Roure et al. Reference Du Roure, Lindner, Nazockdast and Shelley2019). Advances in microfluidics have made it possible to probe various fluid–structure interaction problems. Recent progress has centred around the dynamics of freely transported filaments in steady flows. Examples include morphological dynamics in steady shear (Schroeder et al. Reference Schroeder, Teixeira, Shaqfeh and Chu2005; Harasim et al. Reference Harasim, Wunderlich, Peleg, Kröger and Bausch2013; Kuei et al. Reference Kuei, Słowicka, Ekiel-Jeżewska, Wajnryb and Stone2015; Liu et al. Reference Liu, Chakrabarti, Saintillan, Lindner and du Roure2018; Żuk et al. Reference Żuk, Słowicka, Ekiel-Jeżewska and Stone2021) and hyperbolic extensional flows in the vicinity of stagnation points (Schroeder et al. Reference Schroeder, Babcock, Shaqfeh and Chu2003; Kantsler & Goldstein Reference Kantsler and Goldstein2012; Manikantan & Saintillan Reference Manikantan and Saintillan2015; Chakrabarti et al. Reference Chakrabarti, Liu, LaGrone, Cortez, Fauci, du Roure, Saintillan and Lindner2020). However, many biophysical processes and industrial lab-on-a-chip applications, such as cardiovascular transport (Tarbell et al. Reference Tarbell, Shi, Dunn and Jo2014), filtration of cells and clog mitigation (Cheng et al. Reference Cheng, Ye, Ma, Xie and Wang2016; Lee et al. Reference Lee, Kim, Li, Kim and Kim2018), or sorting and mixing of particles (Dincau, Dressaire & Sauret Reference Dincau, Dressaire and Sauret2020), involve unsteady or time-periodic flows. Elucidating the physics behind the transport and deformation of flexible fibres in such conditions is thus of paramount importance yet remains nascent.
In the present work, we focus on studying the dynamics and morphologies of elastic Brownian filaments in oscillatory shear flow. We use individual F-actin filaments as an experimental model system. Actin is a semi-flexible biopolymer usually found in the cell cytoskeleton and plays a variety of important roles in sub-cellular processes (Pollard & Borisy Reference Pollard and Borisy2003). In physiological conditions, actin polymerizes into filaments whose length 
$L$ is of the order of their persistence length 
$\ell _p$. As a result, they possess a degree of flexibility intermediate between the limits of entropy-dominated long-chain polymers such as DNA (
$\ell _p\ll L$) and stiff, rigid rods (
$\ell _p\gg L$) such as microtubules. For filaments with 
$L \sim \ell _p$, the bending and thermal energies have similar magnitudes, and together they determine the emergent dynamics under viscous loading. The combination of filament elasticity coupled with their rotation (Jeffery Reference Jeffery1922) in flow is central to a plethora of rich dynamics. In particular, actin filaments can undergo buckling instabilities when viscous forces overcome their bending rigidity, much like in the Euler buckling of a macroscopic elastic beam loaded mechanically at its ends. This results in a series of morphological transitions in steady shear flow (Harasim et al. Reference Harasim, Wunderlich, Peleg, Kröger and Bausch2013; Liu et al. Reference Liu, Chakrabarti, Saintillan, Lindner and du Roure2018; Słowicka, Stone & Ekiel-Jeżewska Reference Słowicka, Stone and Ekiel-Jeżewska2020; Żuk et al. Reference Żuk, Słowicka, Ekiel-Jeżewska and Stone2021), which are responsible for the so-called stretch-coil transition (Young & Shelley Reference Young and Shelley2007; Wandersman et al. Reference Wandersman, Quennouz, Fermigier, Lindner and Du Roure2010; Kantsler & Goldstein Reference Kantsler and Goldstein2012), and can lead to the formation of helicoidal structures (Chakrabarti et al. Reference Chakrabarti, Liu, LaGrone, Cortez, Fauci, du Roure, Saintillan and Lindner2020) in compressional flow. In parallel, numerical explorations suggest that these structural instabilities strongly impact the statistical properties of tumbling dynamics under shear (Munk et al. Reference Munk, Hallatschek, Wiggins and Frey2006; Lang, Obermayer & Frey Reference Lang, Obermayer and Frey2014), leading to deviations in classical scaling laws of flexible and stiff polymers (Schroeder et al. Reference Schroeder, Babcock, Shaqfeh and Chu2003). These microscopic instabilities also play a role in the macroscopic rheological properties of polymeric suspensions, where they give rise to shear-thinning and positive normal stress differences in shear flow (Becker & Shelley Reference Becker and Shelley2001; Tornberg & Shelley Reference Tornberg and Shelley2004; Chakrabarti et al. Reference Chakrabarti, Liu, du Roure, Lindner and Saintillan2021).
In such steady flows, the dynamics of a filament is primarily characterized by the dimensionless elastoviscous number 
$\bar {\mu }_{m}$ that compares the time scale of bending relaxation to the shear rate and provides an effective measure of the strength of viscous forces compared to elastic resistance. In steady simple shear flow, compressive viscous forces can lead to an Euler-like buckling instability at a critical value 
$\bar {\mu }_{m}^c \approx 306.8$ (Becker & Shelley Reference Becker and Shelley2001). This gives rise to deformed C-shaped configurations. In even stronger flows, higher-order buckling modes are triggered along with a series of conformational transitions that have been characterized previously (Harasim et al. Reference Harasim, Wunderlich, Peleg, Kröger and Bausch2013; Liu et al. Reference Liu, Chakrabarti, Saintillan, Lindner and du Roure2018).
Time-periodic oscillatory flows introduce another time scale in the problem, namely the period 
$T$ of oscillation, making way for a variety of unexplored morphological dynamics. Here, we characterize systematically filament dynamics in an oscillatory shear flow through a combination of fluorescence microscopy experiments in microfluidic channels, theoretical modelling, and Brownian dynamics simulations. Our experiments use vertical Hele-Shaw cells to drive a horizontal oscillatory flow generated by a pair of pressure controllers connected to the channel inlets. For all experimental conditions, we checked that the Womersley number 
$W\!o$ (Womersley Reference Womersley1955) satisfied 
$W\!o\ll 1$, meaning that oscillation-induced inertial effects were negligible compared to viscous forces. As a consequence, the flow is essentially proportional to the instantaneous pressure gradient and well approximated at the length scale of a filament as a simple shear flow in the horizontal plane with 
$\dot {\gamma }(t) = \dot {\gamma }_m \sin {(2{\rm \pi} t/T)}$. The additional time scale is described through a dimensionless number 
$\rho = \dot {\gamma }_m T$ that compares the maximum shear rate 
$\dot {\gamma }_m$ to the time period 
$T$ of the imposed flow.
In this work, we focus on three central questions. First, we characterize filament morphologies and possible dynamics over a half time period 
$T/2$ of the imposed shear flow. Our findings reveal that, in contrast to steady shear, the emergent dynamics is sensitive to initial filament orientations and are strongly affected by the time period 
$T$ of the flow. Second, we show that the coupling of orientational dynamics of the filaments with the flow results in novel higher-order buckling modes absent in steady shear. Finally, we find that unique to the periodic forcing of the flow is a possibility of suppression of buckling instabilities, even at maximum shear rates that exceed the buckling threshold for steady flow (
$\bar {\mu }_{m} \gg \bar {\mu }_{m}^c$). We corroborate this remarkable dynamical behaviour through a weakly nonlinear Landau theory of buckling. The paper is organized as follows. In § 2, we discuss the experimental set-up, the modelling framework and governing equations of the problem. We present our key findings from experiments, simulations and theory in § 3, and we conclude in § 4.
2. Experimental methods, modelling and governing equations
2.1. Materials and methods
We follow a well-controlled and reproducible protocol for the synthesis of actin filaments, the details of which are described extensively in Liu et al. (Reference Liu, Chakrabarti, Saintillan, Lindner and du Roure2018). The filaments in our experiments are labelled fluorescently and stabilized with phalloidin. They have typical lengths 
$L \sim 5\unicode{x2013}25\,\mathrm {\mu }\mathrm {m}$, diameter 
$a \sim 8\,\mathrm {nm}$ and persistence length 
$\ell _p=17\,\mathrm {\mu }$m, as measured by analysis of thermal-fluctuation-induced conformational changes (Brangwynne et al. Reference Brangwynne, Koenderink, Barry, Dogic, MacKintosh and Weitz2007; Liu Reference Liu2018).
A sketch of the experimental set-up is shown in figure 1. The experiments are conducted in an inverted microscope equipped with an LED light source (Ziess Colibri 7) for dye excitation and a high numerical aperture (NA) water-immersion objective (Zeiss 63x C-Apochromat/1.2 NA), which provides large working distances of up to 
$280\,\mathrm {\mu }$m and a depth of focus around 
$1\,\mathrm {\mu }$m.
Figure 1. Sketch of the set-up used to carry out oscillatory shear flow experiments on actin filaments. (a) Lateral three-dimensional view of the microchannel and imaging system. (b) Top view of the microchannel along with the system for generating pressure-driven oscillatory flows. (c) Imposed pressure (blue symbols) and local shear rate experienced by the filament (red symbols) in a typical oscillatory experiment. Solid black lines are sinusoidal fits. (d) Raw and processed image of a filament. Overlaid to the reconstructed centreline are the compression and extension quadrants of the flow. According to our notation, the mean orientation angle is measured clockwise from the negative 
$x$-axis. This results in compressive viscous stresses for 
$0^{\circ }< \theta <90^{\circ }$, and extensile stresses for 
$90^{\circ }<\theta <180^{\circ }$ when the shear rate is positive, as in the case shown.
The flow is driven in vertical Hele-Shaw poly(dimethylsiloxane) (PDMS) channels of length 
$L_{ch} =30$ mm, height 
$H =500\,\mathrm {\mu }$m and width 
$W = 150\,\mathrm {\mu }$m (figure 1a). Thanks to the high NA of the objective, we visualize the filaments at distance 
$z_{obs}\approx W = 150\,\mathrm {\mu }$m from the coverslip, where the flow profile is Poiseuille-like in the horizontal 
$x\unicode{x2013}y$ plane and has a velocity plateau along the channel height (Liu et al. Reference Liu, Chakrabarti, Saintillan, Lindner and du Roure2018). Besides, in order to avoid wall interaction effects and the change of sign of 
$\dot \gamma$ near the channel centre, we consider only filaments that are located inside the dark blue regions in figure 1(b). This also ensures that we can neglect the out-of-plane shear component of the flow, as in these regions, the in-plane shear rate 
$\dot \gamma _y$ is much larger than the shear rate 
$\dot \gamma _z$ in the perpendicular direction (Liu Reference Liu2018). Therefore, the filaments are deformed mainly by a nearly two-dimensional shear flow, provided that they remain in the observation plane during the experiment. This can be checked systematically by monitoring the filament contour length over time. As shown in figure 8(c), our analysis considers only filaments whose length remains nearly constant with time, further confirming that three-dimensional deformations can be neglected. Finally, the channel width 
$W$ being much larger than the typical dimension of the deformed filament (
${\sim }10\,\mathrm {\mu }$m), the flow gradient can be approximated as constant over that length scale, with the 
$y$-position determining the actual shear rate experienced by the filament (Liu et al. Reference Liu, Chakrabarti, Saintillan, Lindner and du Roure2018).
Since we observe filaments far away from the coverslip, we face two major problems. First, the intensity of the fluorescent light emitted by the filament is low. Second, a slight refractive index mismatch between the solution and the PDMS can result in optical aberrations. To circumvent these problems, we add 45.5 % (wt/vol) sucrose to the actin solution, which allows us to match the refractive index of the PDMS channel (
$n=1.41$). This increases the suspension viscosity to 5.6 mPa s (Liu et al. Reference Liu, Chakrabarti, Saintillan, Lindner and du Roure2018).
Images are acquired using a CMOS camera (HAMAMATSU ORCA flash 4.0LT) at relatively large exposure times in the range 
$t_{exp} = 40\unicode{x2013}60$ ms. This, together with the high NA of the objective, allows the collection of as much light as possible from the filaments, therefore enhancing image quality, but limits the maximum frequency of acquisition to 16–25 fps and sets an upper boundary to the maximum flow speed where filaments can be observed without image blur. Images are processed using a custom-made MATLAB routine, which involves a Hessian-based multiscale filtering, noise reduction using Gaussian blurring, binarization, skeletonization and B-spline reconstruction; see figure 1(d) for an example.
The periodic pressure-driven flow is imposed through a pair of commercially available pressure controllers (Fuigent Lineup Flow EZ) connected to the channel inlets (figure 1b). Each unit has a typical response time 50 ms and a pressure drop in the range 1–25 mbar with 0.1 mbar precision, and is controlled by dedicated software. In order to apply a sinusoidal flow with zero offset, one channel inlet is provided with constant pressure of magnitude 
$p_0$, while a time-dependent pressure 
$p(t) = p_0 + \delta p\sin {(2{\rm \pi} t/T)}$ is imposed in the opposite inlet (
$\delta p=1.3\unicode{x2013}12$ mbar). This results in a total pressure drop 
${\rm \Delta} p(t) = \delta p\sin {(2{\rm \pi} t/T)}$. Typical data are reported in figure 1(c) in blue.
We apply maximum flow rates in the range 
$Q = 0.8\unicode{x2013}11$ nl s
$^{-1}$, corresponding to typical maximum velocities 
$u_x = 15\unicode{x2013}175\,\mathrm {\mu }$m s
$^{-1}$ in the observation plane, filament Reynolds numbers around 
$10^{-5}$–
$10^{-4}$, and maximum shear rates 
$\dot \gamma _m=0.4\unicode{x2013}3.7\,\mathrm {s}^{-1}$. In order to obtain these shear rate values, we have to reduce the pressure drop in the channel. This is achieved as shown in figure 1(b) by inserting before each channel inlet a very thin tubing (length 
$l\approx 10$ cm, inner diameter 
$d= 75\,\mathrm {\mu }$m) with flow resistance 
$R_{tub}$ much larger than the channel resistance 
$R_{ch}\sim 1.5 \times 10^{12}$ Pa s m
$^{-3}$, which is estimated as 
$R_{ch} = 12\mu L_{ch}/[W^3H (1-0.63W/H)]$ (Bruus Reference Bruus2008). Typically, this yields a total flow resistance 
$R_{tot} = 2R_{tub}+R_{ch}$, with 
$R_{tot}/R_{ch}\sim 100$.
Unlike in low-
${{Re}}$ steady flows where only viscous effects dominate, ensuring laminar conditions, inertia can be important in a time-dependent flow and possibly lead to deviations from the steady velocity profile as well as phase shifts between the flow and the pressure gradient, even for 
${Re}\ll 1$. We compute the Womersley number 
$W\!o$ (Womersley Reference Womersley1955), a dimensionless parameter that compares the relative importance of transient inertial effects to viscous forces:
\begin{equation} W\!o = D\left(\frac{2{\rm \pi}\varrho}{T\mu}\right)^{1/2}, \end{equation}
where 
$D$ is a characteristic length scale of the flow (here, taken as the largest transverse size of the channel, 
$D\sim H$), and 
$\varrho$ is the fluid density. In our experiments, the Womersley number is set by the oscillation frequency, and decreases from 
$W\!o\sim 0.5$ to 
${\sim }0.15$ as we vary the period 
$T$ in the range 1–10 s. Since in all cases 
$W\!o<1$, inertia is not relevant, whereas viscous resistance dominates. As shown in Appendix A, the use of small Womersley numbers guarantees that a steady Poiseuille-like velocity profile has time to develop during each cycle, and that the flow is nearly in phase with the instantaneous external pressure (Dincau et al. Reference Dincau, Dressaire and Sauret2020).
In this frequency regime, therefore, the time-dependent shear rate experienced by the filament essentially follows the imposed pressure. This is demonstrated clearly in figure 1(c), where the reported shear rate (red symbols), reconstructed from the filament centre-of-mass velocity along 
$x$ (refer to Appendix A), is well fitted by a sinusoidal function 
$\dot {\gamma }(t) = \dot {\gamma }_m \sin {(2{\rm \pi} t/T)}$.
2.2. Slender-body theory for a Brownian filament
We choose to model the slender filaments as one-dimensional space curves described by their centreline, which is identified by a Lagrangian marker 
${\boldsymbol x}(s,t)$ parametrized by arc length 
$s \in [-L/2,L/2]$. Hydrodynamics is captured using local slender-body theory (SBT) (Tornberg & Shelley Reference Tornberg and Shelley2004), in which the centreline position evolves as
\begin{equation} 8 {\rm \pi}\mu\left[\partial_{t} {\boldsymbol x}(s, t)-{\boldsymbol u}_{\infty}(t)\right]={-}\boldsymbol{\varLambda} \boldsymbol{\cdot} {\boldsymbol f}(s, t). \end{equation}
Here, 
$\mu$ is the fluid viscosity and 
${\boldsymbol u}_\infty (t) = (\dot {\gamma }(t) y,0,0)$ is the background shear flow with oscillatory shear rate 
$\dot {\gamma }(t) = \dot {\gamma }_m \sin (2 {\rm \pi}t/T)$. The force per unit length exerted by the filament on the fluid is modelled as 
$\boldsymbol {f}=B \boldsymbol {x}_{s s s s}-(T \boldsymbol {x}_{s})_{s}+\boldsymbol {f}^{b}$, where 
$B$ is the bending rigidity, 
$T(s)$ is a Lagrange multiplier that enforces inextensibility of the filament and can be interpreted as a line tension, and 
${\boldsymbol f}^b$ is the Brownian force density obeying the fluctuation–dissipation theorem. The local mobility operator 
$\boldsymbol {\varLambda }$ accounts for drag anisotropy and is given by
\begin{equation} \boldsymbol{\varLambda} \boldsymbol{\cdot} \boldsymbol{f}=\left[(2-c) \boldsymbol{\mathsf{I}}-(c+2) \boldsymbol{x}_{s} \boldsymbol{x}_{s}\right] \boldsymbol{\cdot} \boldsymbol{f}, \end{equation}
where 
$c = \ln (\epsilon ^2 e) < 0$ is an asymptotic geometric parameter depending on the aspect ratio 
$\epsilon = a/L \ll 1$. Note that the formulation of (2.2)–(2.3) neglects long-ranged hydrodynamic interactions between distant parts of the filament: these interactions could be accounted for using non-local slender-body theory as in our past work on steady shear flow (Liu et al. Reference Liu, Chakrabarti, Saintillan, Lindner and du Roure2018), where we found that they have a negligible effect on the dynamics. We scale lengths by 
$L$, time by the characteristic relaxation time of bending deformations 
$\tau _r = 8 {\rm \pi}\mu L^4/B$, elastic forces by the bending force scale 
$B/L^2$, and Brownian forces by 
$\sqrt {L/\ell _p}\,B/L^2$. The dimensionless equation of motion then reads
\begin{equation} \partial_{t} \boldsymbol{x}(s, t)=\bar{\mu}_{m}\,\boldsymbol{u}_{\infty}(t)-\boldsymbol{\varLambda} \boldsymbol{\cdot}\left[\boldsymbol{x}_{s s s s}-\left(T \boldsymbol{x}_{s}\right)_{s}+\sqrt{L / \ell_{p}}\, \boldsymbol{\zeta}\right], \end{equation}
where 
$\boldsymbol {\zeta }$ is a Gaussian random vector with zero mean and unit variance. There are three dimensionless numbers that govern the evolution of the filament. The first is the elastoviscous number. As mentioned in § 1, this serves as the measure of the effective hydrodynamic forcing and is defined as
\begin{equation} \bar{\mu}_{m} = \frac{8 {\rm \pi}\mu \dot{\gamma}_m L^{4}}{B\,|c|}. \end{equation}
The second number is the dimensionless time period 
$\rho = \dot {\gamma }_m T$, which enters the dimensionless time-periodic external flow given as 
${\boldsymbol u}_\infty (t) = [\sin (2 {\rm \pi}\bar {\mu }_{m}\,|c|\,t/\rho ) y,0,0]$. Finally, the third number is 
$L/\ell _p$, which captures the strength of thermal shape fluctuations. Consistent with findings from our previous work (Liu et al. Reference Liu, Chakrabarti, Saintillan, Lindner and du Roure2018; Chakrabarti et al. Reference Chakrabarti, Liu, LaGrone, Cortez, Fauci, du Roure, Saintillan and Lindner2020), we will see that thermal shape fluctuations have little direct effect on the filament shape dynamics, and that their primarily role is to trigger instabilities and smooth out sharp deterministic bifurcations.
We solve (2.4) following the numerical algorithm outlined in Tornberg & Shelley (Reference Tornberg and Shelley2004) and Liu et al. (Reference Liu, Chakrabarti, Saintillan, Lindner and du Roure2018). For a freely suspended filament, we use the boundary conditions 
${\boldsymbol x}_{ss} = {\boldsymbol x}_{sss} = \boldsymbol {0}$ and 
$T = 0$ at 
$s=\pm 1/2$. The unknown Lagrange multiplier 
$T(s)$ is obtained by using the constraint of inextensibility 
${\boldsymbol x}_s \boldsymbol {\cdot } {\boldsymbol x}_s = 1$. In our experiments, we have access to filament conformations in the plane of the microscope. To facilitate comparison, we thus restrict ourselves to two-dimensional simulations where the filament is confined to the plane of the flow.
2.3. Conformation characterization
To quantify the morphologies and orientations of the filaments in flow, we introduce the two-dimensional gyration tensor (Liu et al. Reference Liu, Chakrabarti, Saintillan, Lindner and du Roure2018) defined as
\begin{equation} G_{ij}(t) = \frac{1}{L}\int_{{-}L/2}^{L/2} \left [x_i(s,t) - \bar{x}_i(t)\right] \left[x_{j}(s,t) - \bar{x}_{j}(t) \right] {\mathrm{d}} s, \end{equation}
where 
$\bar {{\boldsymbol x}}$ is the instantaneous position of the centre of mass of the filament. The eigenvalues of the gyration tensor (
$\lambda _1,\lambda _2$) are combined to compute a sphericity parameter 
$\omega =1-4\lambda _1\lambda _2/(\lambda _1+\lambda _2)^2$, which quantifies filament deformations: it varies between 
$\omega \approx 1$ (
$\lambda _1\gg \lambda _2\approx 0$) for a straight undeformed configuration, to 
$\omega \approx 0$ (
$\lambda _1 \approx \lambda _2$) for a nearly isotropic, hence strongly deformed, shape. The eigenvector corresponding to the largest eigenvalue is used to obtain the mean orientation angle 
$\theta$ with respect to the negative flow direction, as illustrated in figure 1(d).
3. Results and discussion
3.1. Classification of dynamics over a half-period
In an attempt to illustrate the distinct features of an oscillatory shear flow compared to its steady counterpart, let us first consider the more straightforward case of a rigid non-Brownian rod in the limit of vanishing aspect ratio 
$\epsilon \rightarrow 0$. A rigid rod will rotate and perform a periodic tumbling motion about its centre in a steady, planar shear flow, while translating with the fluid. The tumbling motion is described by Jeffery's model (Jeffery Reference Jeffery1922), which predicts that the tumbling frequency is proportional to the shear rate, 
$\nu _t \sim \dot {\gamma }$ (Jeffery Reference Jeffery1922; Tornberg & Shelley Reference Tornberg and Shelley2004). When the rotational diffusion of the rod is taken into account, the tumbling dynamics is also affected by the additional time scale 
$D_r^{-1}$. As reported in previous works (Puliafito & Turitsyn Reference Puliafito and Turitsyn2005; Kobayashi & Yamamoto Reference Kobayashi and Yamamoto2010; Harasim et al. Reference Harasim, Wunderlich, Peleg, Kröger and Bausch2013), as long as 
$\dot \gamma \gg D_r$, thermal fluctuations dominate over the effects of shear only in a restricted angular region around the flow axis, where the effect of the noise is to drive the rod across the line 
$\theta =0$ in a finite time and thus initiate a new tumbling cycle. In a steady shear flow, therefore, a Brownian filament alternates between fast deterministic phases, in which it undergoes both compressional and extensional viscous stresses, and relatively long diffusive phases dominated by thermal fluctuations around the flow-aligned state. During each tumble, any deformation in the compressive quadrant of the flow is followed by relaxation in the extensional part before a new rotational cycle begins. Therefore, each tumbling phase starts with a thermally fluctuating, flow-aligned straight conformation with 
$\theta _0 \approx 0$.
The situation is altered fundamentally in an oscillatory flow. The mean filament orientation now oscillates with a frequency set by the flow's characteristic time period 
$T$. For a rigid non-Brownian rod in the plane of the flow, Jeffery's equation reads
\begin{equation} \frac{{\rm d} \theta}{{\rm d} t} ={-} \sin^2 \theta(t) \sin (2 {\rm \pi}t/ \rho). \end{equation}
The orientational amplitude and frequency are thus governed by the single dimensionless number 
$\rho =\dot {\gamma }_m T$, which, unlike in steady shear, can be varied independently by changing either the maximum shear rate 
$\dot {\gamma }_m$ or the time period 
$T$. In addition, the periodic flow reversal allows the filaments to explore different orientations at the beginning of each oscillation. For a fixed time period and shear rate, this initial orientation determines the relative fraction of time that the filament spends in the compressive or extensile quadrants before the flow reverses at 
$t = T/2$. As long as the initial filament orientation, as well as its orientation at reversal, remains sufficiently far from the flow direction, the influence of rotational diffusion should remain negligible. We thus expect the deformation dynamics to be sensitive to both 
$\rho$ and the initial orientation 
$\theta _0 \equiv \theta (t=0)$ at the start of the cycle and to be well described by Jeffery dynamics.
Also, the orientational dynamics for the flexible filaments is affected by 
$\bar {\mu }_{m}$ with the possibility of large deformations. These deformations lead to a non-reciprocal conformational evolution. As a result, the flow reversal at 
$t=T/2$ does not generally lead to retracing of the morphological dynamics of the first half-period, which implies that 
$\theta (T) \neq \theta _0$ in general. Thus it is helpful to first restrict our analysis to the half-period 
$0< t< T/2$, where we have precise control over all the dynamical parameters of the problem. Note that we consider only filaments that are not deformed at 
$t=0$. We reserve the discussion of the filament evolution over an entire time period for § 3.4.
Figure 2 illustrates the role of the initial orientation 
$\theta _0$ on the emergent morphological dynamics over the half-period. The mean orientation angle 
$\theta$ in the 
$x$–
$y$ shear plane is measured clockwise from the negative 
$x$-axis; see figure 1(d). For 
$0< t< T/2$, the shear rate is positive, which results in compressive viscous stresses for 
$0^{\circ }< \theta <90^{\circ }$, and extensile stresses for 
$90^{\circ }<\theta <180^{\circ }$. In all these examples, we fixed 
$\ell _p/L = 1.3$ and 
$\rho = \dot {\gamma }_m T = 13$. The elastoviscous number was fixed at 
$\bar {\mu }_{m} = 10^4$, which is significantly higher than the buckling threshold of 
$\bar {\mu }_{m}^c = 306.8$ for steady shear. In figure 2(a), we show snapshots of filament conformations over the half time period from both experiments and simulations, which show very good qualitative agreement. Quantitative agreement is observed in figures 2(b,c), where for each horizontal pair, we respectively plot the time evolution of 
$\omega$ and 
$\theta$ from experiments (coloured symbols) and simulations (solid lines in corresponding colours). We also report in figure 2(c) the evolution of the angle 
$\theta (t)$ expected for a non-Brownian rigid rod obeying Jeffery's equation (3.1) (black dashed line). The good agreement between experiments, simulations and the rigid-rod model indicates that Jeffery's model for orientational dynamics works well even in cases where the polymers are deformed. This also confirms that orientational fluctuations play only a secondary role, in agreement with the fact that we focus on relatively small values of 
$\rho$ (in the range 1–20), for which the filaments do not have enough time to get aligned with the flow axis where Brownian effects are expected to dominate.
Figure 2. The role of initial orientation 
$\theta _0$ on the emerging dynamics of filaments over a half time period. We classify four distinct dynamical behaviours (a) by the characteristic evolution of the anisotropy parameter 
$\omega$ (b) and orientation 
$\theta (t)$ defined in terms of the end-to-end vector (c). The evolution of 
$\theta (t)$ is further compared to a rigid non-Brownian rod model (RR (3.1)) as shown by the black dashed line. The symbols correspond to experiments (E) and compare well with the direct simulations (S) indicated by solid lines. Yellow/grey regions in (c) correspond to compression and extension, respectively. In (b), 
$\omega _c=0.95$ is the value of 
$\omega$ below which a deformation is classified as a buckling event. Parameter values: 
$\bar {\mu }_{m} = 10^4$, 
$\rho = 13$ and 
$\ell _p/L = 1.3$.
Depending on the initial orientation 
$\theta _0$ of the filament, we identify four types of dynamical behaviours, each characterized by distinct morphological dynamics over the course of the half-period. Their characteristics can be summarized as follows.
(i) Continuous buckling (CB). For a small initial angle close to the flow axis,
$\theta _0 \sim 14^{\circ }$ (blue parts of figure 2), the filament buckles in the presence of compressive viscous loading. Here, the filament is contained entirely in the compressional quadrant over the half-period, and never enters the extensional quadrant. After initiation of buckling, it keeps getting compressed, and $\omega$ reaches a minimum when $t=T/2$.(ii) Buckling then partial stretching (BTPS). Upon slightly increasing
$\theta _0 \sim 16^{\circ }$ (magenta parts of figure 2), we transition to a regime where the filament rotates into the extensional quadrant of the flow, as illustrated by the angle $\theta$ in figure 2(c). At this point, the filament stops deforming and gets partially stretched out. This can be seen from the final conformation as well as from the parameter $\omega$.(iii) Buckling then stretching (BTS). Upon further increase of
$\theta _0 \sim 25^{\circ }$ (orange parts of figure 2), we observe a dynamics similar to BTPS with a key difference: since the filament spends less time in the compressional quadrant, its deformation is smaller when entering the extensional quadrant (at around $t=T/4$), so that it has enough time to get stretched entirely before the end of the half-period. This is indicated by a distinct minimum in the evolution of $\omega$ followed by an evolution towards a straight conformation with $\omega \approx 1$.(iv) Tumbling (T). For a large
$\theta _0 \sim 46^{\circ }$ (yellow parts of figure 2), the filament spends little time in the compressive quadrant of the flow. As a result, it is found to tumble without any observable deformation, in a way akin to a rigid Brownian rod. This is quantified further by the evolution of the anisotropy parameter $\omega$ that remains close to unity.
We emphasize that we need to distinguish buckling-induced large deformations from fluctuation-induced conformational changes in both our experiments and simulations. This becomes a challenging task for 
$L \sim \ell _p$. In all of the examples shown here, we have chosen a threshold of 
$\omega _c = 0.95$ (dash-dotted lines in figure 2b) below which we classify a filament shape as buckled. As we will see in the following subsections, this choice does not significantly alter the major conclusions of our analysis.
Figure 3 indicates the same dynamical transitions as a function of the time period 
$\rho$, for a fixed 
$\theta _0$ and 
$\bar {\mu }_{m}$. As 
$\rho$ is increased, the orientation amplitude increases (figure 3c). This results in a transition from tumbling for the smallest 
$\rho$ (yellow) to the occurrence of CB (blue) and eventually BTS (orange) dynamics for the intermediate and largest values of 
$\rho$, respectively.
Figure 3. The role of the dimensionless time period 
$\rho$ on the emerging dynamics over a half time period. The initial orientation in all these examples is fixed at 
$\theta _0 \approx 25^{\circ }$. Parameters: 
$\ell _p/L = 1.3$ and 
$\bar {\mu }_{m} = 5 \times 10^3$. (a) Snapshots from experiments (E) and simulations (S) for increasing values of 
$\rho$. (b) Evolution of the anisotropy parameter 
$\omega (t)$ and (c) of the orientation 
$\theta (t)$. Colour panels and legends are the same as in figure 2.
As highlighted in figures 2 and 3, the two critical aspects of an oscillatory flow are the possibilities of (i) exploring many different initial orientations of the filament, and (ii) controlling independently the angular amplitude of the filament motion. In very strong flows, we can show that these two aspects beget higher-order buckling modes absent in steady shear. Figure 4 illustrates morphological dynamics as a function of the elastoviscous number 
$\bar {\mu }_{m}$ for a fixed orientation and time period. Identical to steady shear, we find that filaments can have a characteristic C-shaped Euler buckling as 
$\bar {\mu }_{m}$ is increased. However, unique to the oscillatory forcing is the emergence of S- and W-shaped configurations at large 
$\bar {\mu }_{m}$ (figures 4b,c). These higher-order buckling modes along with global tumbling of the filament backbone are not observed in steady shear where the filaments are always aligned with the flow direction at the beginning of any time period, and as a result undergo a tank-treading motion with hairpin configurations for sufficiently large 
$\bar {\mu }_{m}$ (Liu et al. Reference Liu, Chakrabarti, Saintillan, Lindner and du Roure2018).
Figure 4. Various buckling modes and conformational dynamics observed with increasing values of 
$\bar {\mu }_{m}$, for a fixed initial orientation 
$\theta _0 \approx 20^{\circ }$ and dimensionless time period 
$\rho \approx 13$.
3.2. Suppression of buckling instabilities in strong flows
A remarkable feature in figure 3 is the absence of buckling instabilities at small 
$\rho$ (yellow), even though the filament is contained entirely in the compressional quadrant of the flow and the elastoviscous number exceeds the theoretical buckling threshold, 
$\bar {\mu }_{m} \gg \bar {\mu }_{m}^c$. A linear stability analysis in the spirit of Becker & Shelley (Reference Becker and Shelley2001) would suggest the emergence of deformed conformations. However, in the presence of a high-frequency periodic forcing, we observe suppression of finite-size deformations. Interestingly, as illustrated in the yellow parts of figure 2, filament buckling is absent even at large 
$\rho$ for initial orientations close to 
$45^{\circ }$, for which the compression rate is maximum.
These two examples highlight the subtle interplay of two nonlinear effects beyond the predictions of a linear stability theory. First, filament orientational dynamics dictate the duration and strength of the effective compressive viscous loading acting along the filament backbone, depending on both filament orientation and instantaneous flow strength. This is sensitive to the initial orientation 
$\theta _0$ and the time period 
$\rho$. Second is the coupling between the duration of filament exposure to compressive forces and the characteristic time scale over which deformations grow to be detectable over thermal fluctuations. In § 3.3, we will develop a weakly nonlinear Landau theory that accounts for all of these effects, but first, we focus on their quantitative characterization from experiments and simulations.
To this end, we performed simulations by systematically varying 
$\rho \in [2,18]$ and 
$\theta _0 \in [0,{\rm \pi} /2)$ while keeping 
$\bar {\mu }_{m} = 2\times 10^4 \gg \bar {\mu }_{m}^c$ and 
$\ell _p/L = 1.3$ fixed. For each combination of 
$(\rho,\theta _0)$, we performed 50 numerical simulations with different random seeds over one half-period for statistical averaging. This then allows us to define and estimate a probability of buckling as
\begin{equation} \mathcal{P}_B(\rho,\theta_0) = \left.\frac{\mathcal{N}(\omega(t) \le \omega_c)}{\mathcal{N}_{{tot}}}\right|_{(\bar{\mu}_{m},\ell_p/L)}, \end{equation}
where 
$\mathcal {N}(\omega (t) \le \omega _c)$ is the number of cases for which the filament buckled at some point during the half-period, and 
$\mathcal {N}_{{tot}} = 50$. It is important to point out that the buckled filament conformations will arise due to the various morphological dynamics characterized in the previous subsection. In contrast to the previous characterization, here we focus on understanding the role of flow frequency and filament orientation in the initiation of the buckling instability, independent of the precise nature of the underlying morphological dynamics. The probability of observing specific types of morphological dynamics for a particular set of parameters and initial condition is discussed in Appendix B.
Figure 5(a) displays the buckling probability as a function of the initial filament orientation 
$\theta _0$ for three different values of the dimensionless time period 
$\rho$. For 
$\theta _0 \gtrsim 70^{\circ }$, instability growth is not sufficient to lead to detectable deformations, and irrespective of 
$\rho$, the buckling probability 
$\mathcal {P}_B$ is small. For intermediate initial angles, buckling probabilities are significant, and we note that for a given 
$\theta _0$, the probability of buckling 
$\mathcal {P}_B$ increases with increasing 
$\rho$. When initial angles get close to 
$\theta _0 \sim 0^{\circ }$, a decrease in the buckling probability is again observed.
Figure 5. (a) Buckling probability 
$\mathcal {P}_B$ as a function of 
$\theta _0$ as obtained from simulations (50 independent runs) for three different dimensionless time periods 
$\rho$. (b) The ensemble-averaged minimum value 
$\omega _m$ of the deformation parameter, corresponding to the maximum deformation, from the same simulations. (c,d) Phase charts of 
$\mathcal {P}_B$ and 
$\omega _{m}$ versus 
$(\rho,\theta _0)$ for a larger set of dimensionless periods 
$\rho \in [2,18]$. In (d), the maximum deformation from simulations (isolines) is also compared to experimental data (symbols) under identical conditions. The vertical coloured lines in (c,d) indicate the constant 
$\rho$ slices shown in (a,b). Parameter values: 
$\bar {\mu }_{m} = 2\times 10^4$, 
$\ell _p/L = 1.3$.
Figure 5(b) shows the ensemble-averaged minimum of 
$\omega (t)$ reached during the evolution over the half-period and serves as a measure of the deformation. The existence of strongly deformed filament conformations with a small value of 
$\langle \omega _m \rangle$ correlates with a probability of buckling close to unity. These behaviours are further corroborated in the phase plots shown in figures 5(c,d), where we varied systematically the dimensionless period in the range 
$\rho =2\unicode{x2013}18$. In figure 5(d), superimposed on the simulations, we also indicate experimental results under identical flow conditions and find good quantitative agreement. Note that here, due to intrinsic difficulties in gathering experimental data with identical values of the elastoviscous number, dimensionless period and initial orientation, the ensemble-averaged simulations (isolines) are compared with only single experiments (symbols).
The observations in figure 5 are in qualitative agreement with the nonlinear effects discussed above. For initial orientations larger than 
$\theta _0 \gtrsim 70^{\circ }$, the time spent by the filament in the compressional quadrant is not sufficient for the instability to grow, regardless of the value of 
$\rho$. For small initial angles, the growth of the instability is limited by the combined effect of the reduced rotational excursion for these initial angles, which impedes the filament from reaching maximum compression, and the time-varying flow.
3.3. Weakly nonlinear Landau theory
In the previous subsection, we highlighted that even in strong flows, filaments may not undergo detectable buckling instabilities, in contrast to linear stability theory. To explain this finding quantitatively, we propose here a weakly nonlinear Landau theory for filament deformations.
Our theory neglects Brownian fluctuations, and we scale time with the inverse shear rate 
$\dot {\gamma }_m^{-1}$. In this case, the non-Brownian SBT equation becomes
\begin{equation} \bar{\mu}_{m}\,\partial_{t} \boldsymbol{x}(s, t)=\bar{\mu}_{m}\, \boldsymbol{u}_{\infty}(t)-\boldsymbol{\varLambda} \boldsymbol{\cdot}\left[\boldsymbol{x}_{s s s s}-\left(T \boldsymbol{x}_{s}\right)_{s}\right], \end{equation}
where 
${\boldsymbol u}_\infty (t) = [\sin (2 {\rm \pi}t/\rho ) y, 0, 0 ]$. We start our analysis by recapitulating the linear stability analysis (Becker & Shelley Reference Becker and Shelley2001). The base state of the filament is an undeformed straight rod. Suppose that at any instant, the rod makes an angle 
$\theta (t)$ with the negative 
$x$-axis. We describe the instantaneous orientation of the rod by the director 
${\boldsymbol p}(t) = -\cos \theta \,\hat {\boldsymbol {x}} + \sin \theta \,\hat {\boldsymbol {y}}$. We also define the vector orthogonal to the director 
${\boldsymbol p}^\perp (t) = \sin \theta \,\hat {\boldsymbol {x}} + \cos \theta \,\hat {\boldsymbol {y}}$. The orientation of the rod is assumed to evolve according to Jeffery's equation (3.1). Associated with this is a time-varying parabolic tension profile given by
\begin{equation} T_0(s,t) = \frac{\bar{\mu}_{{inst}}(t)}{8}\left(s^2 - \frac{1}{4}\right), \end{equation}where
\begin{equation} \bar{\mu}_{{inst}}(t) = \bar{\mu}_{m} \sin (2 {\rm \pi}t/ \rho) \sin (2 \theta) \end{equation}
is the instantaneous elastoviscous number experienced by the filament as it rotates in the unsteady flow. In order to first understand the linear stability, we perturb the filament in the transverse direction. The perturbed state is given by 
${\boldsymbol x} = s {\boldsymbol p} + h_\perp (s,t)\,{\boldsymbol p}^\perp$, where 
$h_\perp (s,t) \sim {O}(\varepsilon )$ with 
$\varepsilon \ll 1$. Linearizing around the base state results in a linear equation for 
$h_\perp (s,t)$ given by 
$\bar {\mu }_{m}\,\partial _t h_\perp = \bar {\mu }_{{inst}}(t)\,\mathcal {L}[h_\perp ]/2 - \partial _s^4 h_\perp$, where 
$\mathcal {L}$ is the differential operator
\begin{equation} \mathcal{L}[f] = \left[1 + s \partial_s + \tfrac{1}{4}\left(s^2-\tfrac{1}{4}\right) \partial_s^2 \right]f(s). \end{equation}
Following Becker & Shelley (Reference Becker and Shelley2001), we use the ansatz of normal modes and write 
$h(s,t) = \phi (s) \exp (\sigma t)$, where 
$\sigma$ is the growth rate. This leads to an eigenvalue problem for the growth rate 
$\sigma$ and the eigenfunction 
$\phi (s)$, from which we find that there is an Euler buckling instability when the elastoviscous number exceeds 
$\bar {\mu }_{m}^c \approx 306.8$. At the onset of this instability, we have 
$\sigma = 0$ and 
$\phi = \phi ^c(s)$. This means that the eigenfunction at the critical threshold satisfies
\begin{equation} \mathcal{L}[ \phi^c] = \frac{2}{\bar{\mu}_{m}^c}\,\partial_s^4 \phi^c. \end{equation} We now proceed to describe the nonlinear dynamics of the filament away from the instability threshold. For that purpose, we develop a weakly nonlinear theory in the vicinity of the first bifurcation. We expand the transverse perturbation on the basis of the first mode of deformation 
$\phi ^c(s)$ as 
$h_\perp (s,t) = A(t)\,\phi ^c(s)$, where 
$A(t) \sim {O}(\varepsilon )$. Transverse perturbations induce changes in length along the axis of the filament. As we will see, these length changes are higher order in 
$A(t)$. We account for this by introducing axial perturbations 
$h_{\parallel }(s,t)$. The perturbed filament conformation can then be written as
\begin{equation} {\boldsymbol x}(s,t) = \left(s + h_\parallel\right) {\boldsymbol p} + h_\perp {\boldsymbol p}^\perp. \end{equation}
Similarly, we perturb the tension around its base state 
$T_0(s,t)$ as
\begin{equation} T(s,t) = T_0(s,t) + \tilde{T}(s,t). \end{equation}
We choose to describe the filament backbone and its mean orientation with the rigid-rod model, which is supported by the observations of figure 2. So at any instant, we assume that the orientation 
$\theta (t)$ follows (3.1). We next outline the set of steps to obtain the evolution equation for the amplitude 
$A(t)$.
(i) Inextensibility. Filament inextensibility implies that
${\boldsymbol x}_s \boldsymbol {\cdot } {\boldsymbol x}_s = 1$, which yields
(3.10)We first note that in the absence of any axial perturbations, the inextensibility condition was satisfied up to\begin{equation} {\boldsymbol x}_s \boldsymbol{\cdot} {\boldsymbol x}_s = 1 + \left[2 \partial_s h_{{\parallel}} + (\partial_s h_\perp)^2 \right] + (\partial_s h_{{\parallel}})^2. \end{equation}${O}(\varepsilon ^2)$. This points to the fact that $h_{\parallel } \sim {O}(\varepsilon ^2)$. We can now satisfy inextensibility up to ${O}(\varepsilon ^4)$ and obtain $h_{\parallel }(s,t)$ by setting the second term in (3.10) to zero. Thus we obtain
(3.11)where we have defined\begin{equation} h_\parallel(s,t) ={-}\frac{1}{2}\int_{{-}1/2}^s \left[\partial_{s'} h_\perp(s',t) \right]^2 \mathrm{d} s'= {A(t)^2}\,v(s), \end{equation}$v(s) = -\frac {1}{2}\int _{-1/2}^{s} \phi ^c_{s'}(s')^2\,{\mathrm {d}} s'$.(ii) Tension perturbation. In order to solve for the tension, we first derive an auxiliary equation from the SBT using the fact that
$\partial _t {\boldsymbol x}_s \boldsymbol {\cdot } {\boldsymbol x}_s = 0$. This provides a second-order ordinary differential equation for the tension $T(s,t)$ as
(3.12)Retaining terms up to\begin{equation} 2 T_{s s}- ({\boldsymbol x}_{ss}\boldsymbol{\cdot} {\boldsymbol x}_{ss}) T={-}\bar{\mu}_{m} {\boldsymbol x}_s \boldsymbol{\cdot} \partial_s {\boldsymbol u} -6 {\boldsymbol x}_{s s s} \boldsymbol{\cdot} {\boldsymbol x}_{s s s}-7 {\boldsymbol x}_{s s} \boldsymbol{\cdot} {\boldsymbol x}_{s s s s} \equiv \mathcal{R}. \end{equation}$A(t)^3$, the equation for the tension perturbation becomes
(3.13)The above linear system can be solved numerically along with the boundary condition\begin{align} 2 \tilde{T}_{ss} &= \left(\bar{\mu}_{m}\,A(t) \cos(2 \theta)\,\phi^c_s + A(t)^2 \left[-\bar{\mu}_{m} \sin(2 \theta)\,(\phi^c_s)^2 - 6 (\phi^c_{sss})^2 - 7 \phi^c_{ss} \phi^c_{ssss} \right] \right) \nonumber\\ &\quad \times\sin(2 {\rm \pi}t/ \rho)+ A(t)^2\,(\phi^c_{ss})^2\,T_0(s,t). \end{align}$\tilde {T}(s=\pm 1/2) = 0$.(iii) Amplitude equation at the linear order. Before deriving the nonlinear amplitude equation, we first consider the nature of the linearized amplitude equation. At the linear order, we have
(3.14)On using (3.7), we can simplify the above equation and obtain\begin{equation} \bar{\mu}_{m}\,\frac{{\mathrm{d}} A}{{\mathrm{d}} t}\,\phi^c(s) = \left(\bar{\mu}_{{inst}}(t)\,\mathcal{L}[\phi^c]/2 - \partial^4_s \phi^c \right)A(t). \end{equation}(3.15)where\begin{equation} \frac{{\mathrm{d}} A}{{\mathrm{d}} t} = \left\langle \phi^c,\partial^4_s \phi^c \right\rangle \left[\frac{\bar{\mu}_{{inst}}(t) - \bar{\mu}_{m}^c}{\bar{\mu}_{m}\bar{\mu}_{m}^c} \right] A(t), \end{equation}$\langle f,g \rangle$ is the inner product of two functions. The above relation amounts to recasting the results of the linear stability analysis, and predicts exponential growth of the amplitude $A(t)$ whenever $\bar {\mu }_{{inst}} > \bar {\mu }_{m}^c$.(iv) Landau equation for the amplitude. In order to derive a nonlinear amplitude equation, we make use once again of the SBT equation. Using our ansatz for the perturbations, we project the equations onto the transverse direction. Retaining all the terms up to order
$A(t)^3$ and using the solution from the perturbed tension equation, we can obtain a nonlinear amplitude equation of the form
(3.16)In this equation,\begin{equation} \frac{{\mathrm{d}} A}{{\mathrm{d}} t} = \left\langle \phi^c,\partial^4_s \phi^c \right\rangle \left[\frac{\bar{\mu}_{{inst}}(t) - \bar{\mu}_{m}^c}{\bar{\mu}_{m}\bar{\mu}_{m}^c} \right] A(t) + \frac{1}{\bar{\mu}_{m}}\left\langle \phi^c, {\varOmega} \right\rangle. \end{equation}$\varOmega$ captures the effect of nonlinearities and is given by
(3.17)It is worth pointing out that the operator in (3.7) used for the stability analysis is not self-adjoint, and as a result the eigenfunctions of (3.7) are not orthonormal. However, this subtlety does not show up in the present analysis since we are interested only in the first mode of deformation. It would be interesting to formulate a theory that expands upon this analysis and incorporates the dynamics of higher-order buckling modes. This would require the use of the adjoint eigenmodes as done in our past work (Chakrabarti et al. Reference Chakrabarti, Liu, LaGrone, Cortez, Fauci, du Roure, Saintillan and Lindner2020).\begin{align} {\varOmega} = A(t)^2\,v(s) \sin^2 \theta + A(t) \left[2 \phi^c_s \tilde{T}_s + \tilde{T} \phi^c_{ss}\right] - A(t)^3 \left[\phi^c_{s} v_{ssss} + (\phi^c_s)^2 \phi^c_{ssss}\right]. \end{align}(v) Numerical solution. We integrate (3.16) numerically using a fourth-order Runge–Kutta scheme. The amplitude equation is complemented with (3.1), which provides the instantaneous orientation of the filament. At every time step, we solve the linear system for the tension perturbation
$\tilde {T}(s)$ to compute the nonlinear terms ${\varOmega }$.
In our numerical explorations, we varied 
$\rho$ and 
$\theta _0$, and studied the evolution of filament deformation predicted from the nonlinear theory. We initiated the simulations by perturbing the filament backbone as 
${\boldsymbol x}(s,t=0) = s {\boldsymbol p} + A_0\,\phi ^c(s)\,{\boldsymbol p}^\perp$, where 
$A_0$ is the initial perturbation amplitude. For all our simulations, we chose 
$A_0 = 0.02$ for the initial amplitude, which is the typical magnitude of transverse fluctuations of a freely suspended actin filament with 
$\ell _p/L \sim {O}(1)$ as considered in our experiments. We characterize an event as buckling in our theory whenever the amplitude exceeds 
$A \approx 0.04$, which corresponds to an anisotropy parameter 
$\omega _c \approx 0.95$. We show the phase chart of minimum 
$\omega _m$ computed from the theory in figure 6(a). It predicts that buckling is favoured only for certain ranges of 
$\theta _0$, and that its probability increases with increasing 
$\rho$, an observation consistent with our experiments and Brownian simulations shown in § 3.2.
Figure 6. (a) Phase plot of the maximum deformation for 
$\bar {\mu }_{m} = 3.5 \times 10^3$ as computed from the nonlinear Landau theory outlined in § 3.3. The two vertical coloured lines indicate slices of 
$\rho = 6$ (b) and 
$\rho = 16$ (c), respectively, for which we show the evolution of the growth rate 
$\sigma (t)$ and the perturbation amplitude 
$A(t)$ as functions of the instantaneous orientation of the filament, for several values of the initial angle 
$\theta _0$.
In order to explain this behaviour, we now plot the evolution of the growth rate 
$\sigma (t)$ and the amplitude of the perturbation 
$A(t)$ as functions of the orientation angle of the filament 
$\theta (t)$ from the nonlinear model. Typical data for two values of 
$\rho$ and several initial orientations 
$\theta _0$ are reported respectively in figures 6(b) (
$\rho =6$) and 6(c) (
$\rho =16$).
A number of observations can be made from these figures. First, the range of orientation angles traversed by the filament varies strongly as a function of the initial orientation, even for a fixed value of 
$\rho$ (compare, for example, the blue and red curves in figure 6b). This is linked to the fact that the Jeffery dynamics is slow close to alignment with the negative flow direction, but faster for increasing initial angles. Filaments starting aligned close to the negative flow direction thus do not rotate very far compared to filaments with initial angles closer to 45
$^\circ$. Second, filaments that start close to 90
$^\circ$ will spend little time in the compressive quadrant, where positive growth rates can exist. In a time-independent shear flow, all filaments will experience maximum compression and thus a maximum growth rate at 45
$^\circ$. In the case of an oscillatory flow, the flow forcing itself is time-dependent and the maximum growth rate results from a combination of the orientation angle and the time-dependent flow intensity, as is clear from the figures where maximum growth rates are not necessarily observed at 45
$^\circ$, even for filaments sweeping by this orientation. Comparing figures 6(b) (
$\rho =6$) and 6(c) (
$\rho =16$), it is obvious that larger values of 
$\rho$ lead to larger growth rates over more extended ranges of angles. It is also interesting to note that negative growth rates are observed even in the compressive quadrant, for situations where 
$\bar {\mu }_{inst}$ remains below the buckling threshold.
As pointed out above, to observe a buckling instability, it is not sufficient for the filament to experience a positive growth rate; the perturbation also needs sufficient time to grow so as to lead to a measurable signature in the anisotropy parameter 
$\omega$. This can be seen from the bottom plots of figures 6(b,c), which show the instantaneous amplitude resulting from the cumulated growth of the perturbation during rotation of the filament through the compressive quadrant. These examples thus underscore a key aspect of our experimental measurements: linear instability does not always correspond to detectable deformations of the filament.
This critical distinction between linear instability and measurable deformations also applies to steady shear flows, even if to a lower extent, but seems to have been unreported. We elaborate further on this point in Appendix C.
3.4. Dynamics over one full period
So far, we have focused on the deformation dynamics in the first half of the forcing period, considering only the particular case of a straight filament at 
$t=0$. This has helped us to understand the role of the time period 
$\rho$ and the initial orientation 
$\theta _0$ that can be both controlled independently and varied systematically. We now discuss briefly the filament dynamics over one entire period of the oscillatory flow. The key aspect here is that the orientation and the deformation of the filament at the beginning of the second half of the period is governed by its deformation and orientational dynamics during the first half, giving rise to new, interesting dynamics. In particular, one also has to consider cases in which the filament is already deformed at the beginning of the second half-cycle. The interplay between oscillation period, flow strength, initial orientation and deformation of the filament leads to a very rich spectrum of dynamical modes. Typical examples are reported in figure 7. In figure 7(a), the filament undergoes a continuous buckling (CB) event during the first half-cycle. The filament remains oriented in the extensional quadrant when the flow is reversed. As a result, it is stretched by the flow during the second half-cycle, and the anisotropy parameter increases continuously. As 
$\rho$ is increased (figure 7b), the dynamics in the first half of the period transitions to a BTPS event. Consequently, when the shear is reversed, the filament first experiences compression, deforming from an already buckled state. As the filament enters the extensional quadrant, it is then stretched out by the flow.
Figure 7. Dynamical features over a complete period of tumbling. During the second half-period, the filament can be either continuously stretched by the flow (a) or compressed once again (b). Parameter values: 
$\ell _p/L = 1.3$ and 
$\bar {\mu } = 1 \times 10^4$.
These examples give a glimpse of the complex evolution of the filament deformation over multiple cycles. Since the deformation in general brings about irreversible orientational dynamics (
$\theta (T) \neq \theta _0$), the asymptotic filament dynamics is expected to be characterized by the emergence of a chaotic behaviour, even in the athermal limit and at small Reynolds number.
4. Conclusion
Using stabilized actin filaments as a model system, we have explored systematically the dynamics of semi-flexible polymers in oscillatory shear flows. Our microfluidic experiments were shown to be in excellent agreement with Brownian simulations based on Euler–Bernoulli elasticity and local slender-body theory hydrodynamics, and were also well described by a weakly nonlinear Landau model. Through a systematic variation of the filament length and of the shear rate and frequency of the imposed flow, we characterized the orientational and morphological dynamics of the filaments over a half oscillation period of the flow, and identified four distinct modes of deformation that are unique to periodic forcing.
For a given flow strength (elastoviscous number 
$\bar {\mu }_{m}$), these morphological dynamics were shown to be sensitive to the dimensionless time period 
$\rho$ of the flow and initial orientation 
$\theta _0$ of the filament, which can take on any values at the start of an oscillation cycle. These two variables govern the sign, magnitude and duration of the viscous stresses experienced by the filament as it rotates through the applied flow, giving rise to a rich variety of behaviours from rigid tumbles to elastoviscous buckling followed by varying levels of stretching. In the case of buckling, a wide range of deformation modes (C, S and W modes) was also found to emerge. These findings are in contrast with the case of steady shear, where every filament tumble initiates from the same flow-aligned stretched conformation, resulting in a more constrained set of deformation modes (C modes and tank-treading dynamics) that are governed entirely by the constant applied flow strength (Liu et al. Reference Liu, Chakrabarti, Saintillan, Lindner and du Roure2018).
Another striking feature of oscillatory shear flow is the absence of observable buckling in strong flows at high oscillation frequencies (low values of 
$\rho$), even when the maximum elastoviscous number greatly surpasses the buckling threshold for steady shear (
${\bar {\mu }_{m}\gg \bar {\mu }_{m}^c}$). This is a consequence of the short duration of time during which compressive stresses are applied before the flow changes direction: the half-period is too short for flow-induced deformations to grow sufficiently to be detected over thermal shape fluctuations, even when buckling is predicted to occur. Our simulations were used to calculate a probability of buckling, which peaks at a value close to 1 for 
$\theta _0\approx 25^\circ$ in strong flows, but drops sharply when 
$\rho \lesssim O(1)$.
Our study has focused on characterizing the dynamics over a half-period of the flow, and only briefly addressed the case of a full period in § 3.4. There, the dynamics is even richer as the conformation of the filament at the end of the first half-period may involve significant deformation. The evolution of filament conformations over multiple time periods remains an open question of great interest: for instance, it is not known whether the filament dynamics would cycle, perhaps chaotically, through a variety of morphological transitions, or if the system may end up selecting specific deformation modes corresponding to preferred initial orientations after multiple oscillations (Agrawal & Mitra Reference Agrawal and Mitra2022). Some of the findings from our study also hint at interesting rheological behaviours for suspensions of semi-flexible polymers in time-periodic flows. In steady shear, our previous work demonstrated that the onset of buckling instabilities is accompanied by enhanced shear-thinning and an increase in normal stress differences (Chakrabarti et al. Reference Chakrabarti, Liu, du Roure, Lindner and Saintillan2021). In an oscillatory flow, we expect the same effect to occur but also anticipate a non-trivial dependence on the flow period, with a possible increase in viscosity at high frequencies due to the suppression of buckling in that regime. Some of these outstanding questions will be addressed in future work.
Acknowledgements
We are grateful to the team of G. Romet-Lemonne and A. Jégou for providing purified actin, and to Michael Shelley for fruitful discussions. This work has benefited from the technical contribution of the joint service unit CNRS UAR 3750.
Funding
A.L. and F.B. acknowledge support from the ERC Consolidator Grant PaDyFlow under grant agreement 682367. We thank Institut Pierre-Gilles de Gennes (Investissements d'avenir ANR-10-EQPX-34).
Declaration of interests
The authors report no conflict of interest.
Author contributions
F.B. and B.C. contributed equally to this work.
Appendix A. Velocity profile and shear rate in the observation plane
In this appendix, we show that the flow profile in the channel horizontal plane is always parabolic over the course of a typical oscillatory experiment, and highlight the procedure used to obtain the time-dependent, local shear rate experienced by the filament.
Figure 8(a) shows the measured velocities perpendicular (
$v_{cm}$) and parallel (
$u_{cm}$) to the flow of 
$n=5$ filaments that have been tracked during a typical oscillatory experiment with 
$T=5$ s. The intermittent lack of experimental points in the data for some filaments is due to the fact that these filaments are located close to the borders of the image, and so were lost temporarily during the oscillation. First, we note that all the measured transverse velocities fluctuate only around zero. This indicates the absence of any significant drift perpendicular to the flow. As expected, periodic oscillations are observed for the velocities in the streamwise direction. Fitting to a sinusoidal function (solid black lines) yields, for each filament, an independent estimate of the period 
$T$, phase 
$\phi$ and amplitude 
$\delta u_{cm}$ of the oscillation. We get 
$T=5.094\pm 0.003$ s and 
$\phi =(-0.09\pm 0.02)$ rad (
${\rm mean}\pm {\rm std}$) for the period and phase, respectively. Since the camera and the pressure controllers are synchronized, 
$\phi$ corresponds to the actual phase shift between the imposed pressure and the flow. It indicates that the flow lags behind the input pressure only slightly. Figure 8(b) shows the velocity amplitude 
$\delta u_{cm}$ inferred from the fit as a function of the 
$y$-position of the filament 
$y_{cm}$. As seen, the profile is well described by a Poiseuille flow (black solid line). The fit also yields channel width 
$W=147.5\pm 16\,\mathrm {\mu }$m, in good agreement with the expected one. Finally, we report in figure 8(c) typical data showing the time evolution of the filament length 
$L$, as inferred from image post-processing, during half-periods of oscillation. Data and colours correspond to the first three examples from the top of figure 2 in the main text, and are normalized by the length at the beginning of the cycle 
$L_0 = L(t=0)$. The small fluctuations in filament length are due mainly to experimental noise in the reconstructed images and demonstrate clearly that the filaments remain in the horizontal shear plane during the half-cycle, allowing us to neglect three-dimensional effects. Consequently, as the flow is parabolic at any time, and considering the flow gradient as constant over the length scale of the filament (§ 2.1), the instantaneous shear rate is computed as
\begin{equation} \dot\gamma(t) = \frac{u_{cm}(t)\left[W-2\,{\rm \Delta} y_{cm}\right]}{{\rm \Delta} y_{cm} \left[W-{\rm \Delta} y_{cm}\right]}, \end{equation}
where 
${\rm \Delta} y_{cm}=\vert \langle y_{cm}(t)\rangle -y_{wall}\vert$ is the absolute distance between the average centre-of-mass position 
$\langle y_{cm}(t)\rangle$ of the filament and the channel wall.
Figure 8. (a) Parallel (
$u_{cm}$) and perpendicular (
$v_{cm}$) centre-of-mass velocities of filaments located in different regions of the observation plane. (b) The velocity amplitude 
$\delta u_{cm}$ versus 
$y$-position of the filament is parabolic (solid line), which attests to a Poiseuille-like flow profile. (c) Evolution of the filament contour length over a half-cycle of oscillation for the first three examples (from the top) in figure 2. The data are normalized by the filament length at the beginning of the oscillation.
Appendix B. Statistical nature of the morphological dynamics
In § 3.1, the different dynamics over the half-period are classified based on the orientational dynamics and the conformation parameter (
$\omega$) of the trajectories. As stated in the main text, due to Brownian fluctuations, a trajectory with same set of parameters and initial conditions may exhibit different morphological dynamics in different realizations. In figure 9, we show an example of how an ensemble of initial conditions and parameters spans the space of various morphological transitions. We report the probability of the different deformation modes versus 
$(\theta _0,\rho )$, as obtained from a set of 50 simulations with 
$\bar {\mu }_{m} = 10^4$ and 
$\ell _p/L=50$. As illustrated, the different morphological dynamics peak at different locations in the phase space, although one can clearly observe some overlap between the distributions; the CB probability is higher for high oscillation frequencies and flow-aligned states, while the most probable deformation mode is BTS, upon increasing 
$\rho$ and the initial orientation. As expected, the region in which one observes the largest probability of simple tumbling, T, coincides with the smallest probability for the filament to buckle (see figure 5). Note that these distributions are in good agreement with the data in figures 2 and 3, despite the larger value of 
$\ell _p/L$ used here.
Figure 9. The distribution of the various deformation modes as a function of the dimensionless period and initial filament orientation, estimated from numerical simulations. Parameter values: 
$\bar {\mu }_{m} = 10^4$ and 
$\ell _p/L = 50$.
Appendix C. Revisiting the linear stability analysis in steady simple shear flow
In the main text, we have explained how the buckling dynamics depends on the orientation and the time period of the flow. One of the highlights of our present work is the absence of buckling at large flow rates. Here, we revisit the classical linear stability analysis in steady simple shear flow (Becker & Shelley Reference Becker and Shelley2001), and delineate the dependence of buckling instability and filament orientations. For steady shear flow, the eigenvalue problem for the growth rate 
$\sigma$ is given by
\begin{equation} \bar{\mu}_{m} \sigma \phi =\tfrac{1}{2} \bar{\mu}_{{inst}}(t)\,\mathcal{L}[\phi] - \partial_s^4 \phi, \end{equation}
where 
$\bar {\mu }_{{inst}}(t) = \bar {\mu }_{m} \sin 2 \theta (t)$. The orientation obeys Jeffery's equation as before: 
$\dot {\theta } = -\sin ^2 \theta$. As mentioned previously, the buckling instability sets in at 
$\bar {\mu }_{m}^c \approx 306.8$. Since the effective flow strength varies with orientation, so does the growth rate 
$\sigma$. At the onset of the instability, the growth rate 
$\sigma$ is positive only at 
$\theta = {\rm \pi}/4$. Figure 10(a) shows the variation of the growth rate 
$\sigma (t)$ and the associated orientation angle 
$\theta (t)$ as functions of time for 
$\bar {\mu }_{m} = 700$. We are showing the compressional quadrant of the problem where the growth rate is positive for only a small range of angle 
$\theta (t)$ for which 
$\bar {\mu }_{{inst}}(t) > \bar {\mu }_{m}^c$. The angle 
$\theta (t)$ at which the growth rate first becomes positive is denoted by 
$\theta _0^c$. We now use our weakly nonlinear theory to find the evolution of the amplitude 
$A(t)$ of a perturbation for this case. Figure 10(b) shows 
$A(t)$ for three different initial orientations 
$\theta (t=0) = \theta _0$. For 
$\theta _0 < \theta _0^c$, the growth rate 
$\sigma (t)$ for the initial period is less than zero and the perturbation amplitude decays. Thus even though the growth rate becomes positive when 
$\theta > \theta _0^c$, buckling cannot be detected from the filament conformation. If 
$\theta _0 > \theta _0^c$, then the perturbation can grow for a while before 
$\sigma (t)$ becomes negative as 
$\theta (t)$ enters the compressional quadrant. As expected, maximum growth is observed for 
$\theta _0 = \theta _0^c$ and provides us with the possibility of detecting buckling from filament conformations. This example highlights two important points. First, it re-emphasizes the role of initial filament orientation in detecting buckling from filament conformations even in steady shear. Second, it highlights the crucial role of Brownian fluctuations. Unlike in the non-Brownian Landau theory, in Brownian simulations, the filament backbone is perturbed continuously as it rotates, which makes it more susceptible to buckling instabilities.
Figure 10. (a) Variation of the growth rate 
$\sigma (t)$ and orientation angle 
$\theta (t)$ of a filament in steady flow as a function of time. Here, 
$\theta _0^c$ indicates the critical angle at which the eigenvalue becomes positive for the first time. The markers indicate initial angles 
$\theta _0$ for the amplitude evolution shown in (b). (b) Evolution of the deformation amplitude 
$A(t)$ for three different initial orientation angles 
$\theta _0$ as computed from the nonlinear Landau equation.
Open access
