Taming active turbulence with patterned soft interfaces

Active matter embraces systems that self-organize at different length and time scales, often exhibiting turbulent flows apparently deprived of spatiotemporal coherence. Here, we use a layer of a tubulin-based active gel to demonstrate that the geometry of active flows is determined by a single length scale, which we reveal in the exponential distribution of vortex sizes of active turbulence. Our experiments demonstrate that the same length scale reemerges as a cutoff for a scale-free power law distribution of swirling laminar flows when the material evolves in contact with a lattice of circular domains. The observed prevalence of this active length scale can be understood by considering the role of the topological defects that form during the spontaneous folding of microtubule bundles. These results demonstrate an unexpected strategy for active systems to adapt to external stimuli, and provide with a handle to probe the existence of intrinsic length and time scales.

M any examples of either living or activated entities display intriguing modes of non-equilibrium behavior that cover a wide range of length scales and obey vastly different operational mechanisms. This general framework defines what is referred to as the realm of active matter 1,2 . Flows in active systems are intrinsically autonomous and their lack of spatial and temporal coherence seem to preclude any sort of external controlling intervention. This vision is supported by the emergence of the flow regime known as active turbulence [3][4][5][6][7][8] , which determines the collective behavior in bacterial colonies [9][10][11] or cytoskeletal networks [12][13][14] , and has a major role in chemical transport and nutrient distribution 9 .
In spite of recent theoretical [15][16][17] and experimental [18][19][20] attempts, less attention has been dedicated to analyzing the role of confining conditions on active fluids. This is an important issue in a cell biology context, for instance in relation to cytoplasmic streaming in eukaryotes 21 . In fact, recent experiments with active biofluids have revealed that strong coupling with substrates can lead to the steering of otherwise chaotic flows 11,22 . Boundary conditions can determine the configuration of soft materials, such as colloidal suspensions and liquid crystals. Yet, a fundamental question is how active materials, endowed with nonequilibrium length and time scales, adapt to confinement, and whether external constraints can be used to modify the behavior of active systems in a controlled manner, for instance in a microfluidic lab-on-a-chip context.
In vitro reconstitutions of the cytoskeletal network 12 provide an excellent model active material to address such fundamental issues without the intricacies of a living system. A recently developed tubulin-kinesin gel has been shown to organize as an ordered active material when condensed at a soft interface. This preparation has been given the name active nematic (AN) 12 . Under continuous supply of adenosine triphosphate (ATP), a layer of ordered filamentous protein bundles self assembles featuring a turbulent-like behavior. The currently accepted paradigm for active turbulence [3][4][5][6][7] in ANs assumes that vortical flows originate from local stress sources that are associated to orientational gradients within the active material [23][24][25][26] . Recent numerical studies 7 suggest that, in spite of the apparently chaotic flow behavior, this regime features an intrinsic length scale, l α $ ffiffiffiffiffiffiffiffi ffi K=α p . This scale arises from a balance between active and passive stresses [23][24][25] , and may be experimentally varied by tuning the activity (α) or the rigidity (K) of the AN.
Within this framework, here we study the mechanism with which an active nematic constrained to flow in a quasi-twodimensional geometry adapts to lateral confinement exerted by a structured liquid interface. Our analysis is based on providing a quantitative assessment of the geometry of the active flows, comparing the AN in contact with an isotropic oil with the AN flowing in contact with a reconfigurable interface featuring a well-defined mesostructure. In this work, we will call the former the unconstrained configuration while the latter will be the constrained configuration. We begin by presenting clear experimental evidence for an exponential distribution of vortex sizes in the unconstrained active turbulence regime. Next, we show that this geometry is reversibly altered when the active flows are locally captured by mesoscopic circular domains, which impose their scale-free size distribution on the AN. By modifying the activity and the rigidity of the AN, we demonstrate that the active length scale reemerges in the constrained regime as a cutoff domain size for the effective capture of the active flow. Finally, we demonstrate that the prevalence of l α is due to its intimate relation with the dynamics of topological defects that underlie the active flows.

Results
Geometry of unconstrained active turbulence. The material we consider here is constituted by fluorescently labeled, micron-sized microtubules (MTs), bundled together by the depleting action of  Fig. 1 Structure of the active turbulent flow in contact with an isotropic oil. a Confocal fluorescence micrograph (375 × 375 μm 2 ) of the AN in the turbulent regime when in contact with an isotropic oil. An example of the proliferating +1/2 (blue) and −1/2 (red) defects is sketched on top of the image (Supplementary Movie 1). b Instantaneous flow field (vector plot) and computed Okubo-Weiss parameter field (density plot, arbitrary units), which is employed to determine the location and size of vortices. c Binary image corresponding to the Okubo-Weiss field. d Statistical analysis of the distribution of vortex sizes accumulated for a series of one thousand snapshots during the active turbulent flow. The solid line is an exponential fit to the data. The range of vortex sizes is limited by the field of view polyethylene glycol (PEG), and cross-linked by clusters of kinesin motor proteins 12 . Under continuous supply of ATP, the tens-ofmicrons-long filamentous bundles are subjected to permanent internal stresses due to the action of the motor complexes (see Methods). This aqueous preparation spontaneously selfassembles into a flow-permeated active gel 27-29 that condenses as an AN layer at surfactant-decorated oily interfaces. As a reference experiment, a layer of AN is interfaced with an isotropic oil in a custom-built cell (see Methods). The AN textures that form upon continuous stretching and folding of MT bundles are punctured by disclinations that proliferate under high activity conditions, giving rise to a two-dimensional active turbulent steady state 3-8 ( Fig. 1a; also see Supplementary Movie 1). Defects can be easily identified in the micrographs as moving regions devoid of MTs around which active filaments are organized either with a parabolic (+1/2 defects) or a hyperbolic (−1/2 defects) arrangement (Fig. 1a). Because of their symmetry, moving +1/2 defects actively organize the flow, while −1/2 defects are simply advected. This configuration constitutes the active turbulent regime, whose geometry we have characterized by analyzing the statistical distribution of vortex sizes. To locate and measure the area of each vortex, the local instantaneous velocity of the active flow is evaluated from a sequence of micrographs. Velocimetry data are used to obtain values of the Okubo-Weiss parameter, which provides with a standard criterion for vortex location by considering the extension of each vortex to be bound by the condition OW < 0 7, 30 (Fig. 1b, c). A statistical analysis reveals an exponential distribution of vortex sizes (Fig. 1d), from which we can extract a characteristic vortex area, A * , that determines de decay of population as a function of vortex size. For the experimental conditions of the experiment in Fig. 1 ([ATP] = 700 μM, [PEG] = 0.8% (wt/wt), oil viscosity = 12.5 Pa s) we estimate A * =585±17 μm 2 . We can, therefore, extract a characteristic length scale for the unconstrained flow geometry l α $ ffiffiffiffiffi ffi A Ã p that is known to increase when the AN rigidity increases or when its activity decreases 7 . We can also extract a characteristic time scale for this system by analyzing the lifetime of individual vortices. To this purpose, we have analyzed the time oscillations of the OW parameter in regions of size A * (Supplementary Fig. 1), and obtained a well-defined frequency of 25 mHz, which yields a typical lifetime of~40 s for vortices in this regime. The experiments reported here have been performed with the AN flowing in contact with a silicon oil with a dynamic viscosity of the same order as the low-shear viscosity of the anisotropic oil described below to better appraise the effect that imposing lateral order has in the active flow geometry.
Active flow constrained by a patterned interface. In order to impose lateral order on the AN flow, we take advantage of the hydrodynamic coupling that has been demonstrated to exist between the active and passive fluids across the water/oil interface 31 . Moreover, we take advantage of the giant viscous anisotropy of Smectic-A (SmA) liquid crystal (LC) oils 32 to establish a patterning in the boundary conditions for the active flow. Note that a patterning imposed by means of a solid surface is unpractical in the present case, since the AN has been observed to only form in contact with a liquid interface 12 . The chosen mesogen, octyl-cyanobiphenyl (8CB), features mesophases at temperatures compatible with protein activity (see Methods). At room temperature, free energy minimization results in the SmA phase spontaneously self-assembling into polydisperse domains, known as toroidal focal conic domains (TFCDs) 32 , that organize into a fractal tiling known as Apollonian gasket 32,33 (Fig. 2a). At the water/LC interface, TFCDs have a circular footprint, and are formed by concentric SmA planes perpendicular to the interface. Thus, 8CB molecules, which are both parallel to the water layer due to the interaction with the used surfactant, and perpendicular to the SmA planes, orient radially in concentric rings ( Fig. 2a; Supplementary Fig. 2). Crucially, molecules in the SmA phase can diffuse freely within a given smectic plane but their transport is severely hindered in the direction perpendicular to the planes. Although 8CB is, from a macroscopic point of view, a liquid, the local interfacial shear stress probed by the active material is This results in the size distribution of swirls to become commensurate with the Apollonian gasket. The original exponential vortex size distribution in the unconstrained turbulent AN is now transformed into a scale-free power law, matching the distribution of underlying domains, for which we find a better fit with an exponent −0.7 (Fig. 2d). Our experiments reveal the emergence of a cutoff size above which TFCDs can locally imprint a swirl on the AN, while small domains can, at most, scatter moving defects. To provide with a quantitative assessment of the flow confinement due to a given TFCD, we define, for each +1/2 defect that hovers within the domain's limits, a winding number, Q = Δα/2π, where Δα is the accumulated angle of rotation around the TFCD center until the defect is either annihilated or it escapes from the domain. For the reference experimental conditions [ATP] = 700 μM, [PEG] = 0.8% (wt/wt), we have plotted the average value of Q in Fig. 2d (Supplementary Fig. 4a). Our data reveal that domains smaller than 25 μm in diameter do not significantly alter the defect trajectories, exerting an average confinement Q < 0.5, while domains larger than 40 μm in diameter exert a confinement Q > 1, effectively trapping each passing defect. When comparing the unconstrained and the constrained active flow, we observe that the characteristic vortex size measured from the exponential distribution in the former (indicated by a line in Fig. 2d) matches with the cutoff domain size observed in the latter.
The effectiveness of the active flow confinement depends both on the size of each TFCD and on the state of the AN. We have quantified the influence on the average winding number of the concentration of ATP, which will directly modify the activity parameter α~log[ATP], and the bending rigidity (K) of the active bundles, which is expected to increase, albeit with an unknown dependency, with the concentration of the depleting agent PEG 34  . For a TFCD with diameter 60 μm (Fig. 3a), for which Q > 1 under the reference experimental conditions, passing defects are forced to circulate around the domain center, following the circular easy flow directions defined by the underlying SmA phase (Fig. 3a). Accumulated defect trajectories reveal an inner circular region of diameter~20 μm forbidden for circulating defects. A decrease in α (Fig. 3b, d) or an increase in K (Fig. 3c, e) lead to a shift to smaller Q values for a given domain size, i.e., to an increase in the critical size for domains to capture the active flow.
Intrinsic length and time scales in the constrained flow.
To understand the origin of the threshold TCFD size for the effective capture of moving filaments, we focus on how the AN self-organizes into swirls. As a characteristic feature of active turbulence, defects are permanently renovated through instabilities of the elongated filaments 23-25 , even their number may change to a small extend, continuously assembling and disassembling the core structure of each rotating swirl. From a topological perspective, the confined rotating filaments constitute a singularity of total charge S = + 1, which corresponds to a full rotation around the domain center. However, the spontaneous folding of the extensile AN filaments can only create +1/2 defects (at the tip of the fold) and −1/2 defects (at the tail of the fold). This topological mismatch poses an additional constraint on the ensemble of AN defects that may evolve within a TFCD. At all times, a dynamic balance must be established such that the arithmetic sum of defect charges inside a single domain adds up to one, as dictated by topology ( Fig. 4; also see Supplementary Movie 4). Clearly, the minimum number of semi-integer defects required to organize a rotating swirl is two S = + 1/2 defects (Fig. 4a). Since defect spacing is regulated by l α 7, 8, 25, 35, 36 , the above topological arguments establish a minimum size for a rotating swirl. Considering that the threshold domain size for effective flow trapping scales with material parameters consistently with l α , we conclude that it is indeed this active length scale that also determines the crossover from turbulent to laminar flow under constrained conditions. The ability of constraining the active flow allows to easily probe intrinsic time scales in the system. It is known that the extensile nature of the studied AN material makes a configuration of parallel elongated filaments intrinsically unstable, leading to a defect-forming bending instability 23,24 . As observed in Figs. 2c and 3a, the distribution of the trapped active filaments around large TFCDs is heterogeneous, with a central region mostly devoid of MTs surrounded by a denser region occupied by the outer arms of +1/2 defects. Since filaments are locally parallel within this region, the inner structure of organized S = + 1 swirls is unstable, being periodically disrupted and reconstituted. An example of this process is depicted in Fig. 5 (Supplementary  Movie 5), where we show a swirl, which has been assembled by a TFCD of approximate diameter 200 μm. The vortices in the original turbulent regime have been transformed into swirls of circulating laminar flow. In the first micrograph (Fig. 5a), the AN is structured as an annular band of circularly-aligned MT bundles. This configuration is unstable (Fig. 5b), and the filaments develop a radial buckling instability that leads to the formation of semi-integer defect pairs, moving inwards and outwards of the domain (Fig. 5c, d). Some of these new defects will eventually annihilate with existing inner ones (Fig. 5e, f), finally reorganizing the original circularly-aligned state (Fig. 5g). It is remarkable that the handedness of flow rotation inside each swirl, which is randomly selected upon self-assembly, is preserved through these episodes of structural instability. This reproducible sequence of dynamic events is also exhibited as a time periodic modulation of the velocity of the active flows constrained inside domains. In Fig. 6a, we show the trajectory of a micro-particle trapped at the SmA/water interface, moving along the circular easy-flow directions of the TFCD, and propelled by coupling with the AN. An analysis of the time evolution of the particle's trajectory, which acts as an active flow tracer, displays a remarkable periodicity (Fig. 6b, c; also see Supplementary Movie 6), thus revealing an intrinsic time scale of the AN 22,25 . The average period of the oscillations, T~25 s, is of the same order as the lifetime of vortices in the unconstrained AN flow.

Discussion
It is well accepted in the literature that the spatial arrangement of ANs is determined by an intrinsic length scale l a , which establishes spatial features in the material such as the steady-state defect separation and the bending radius of the extensile MT bundles. Earlier simulations revealed that this length scale also determines the geometry of the active turbulent regime 7 , where the vortex size distribution follows an exponential whose decay length is related to l α . Our experiments provide quantitative evidence to confirm these predictions when the AN evolves in contact with an isotropic oil (Fig. 1). On the other hand, when flowing in contact with the highly anisotropic SmA phase of 8CB, the active turbulence is reorganized in a scale-free power law distribution of swirl sizes that adapt to the underlying SmA tiling. Even in this regime, an intrinsic length scale resurges in the form of a minimum swirl size inside of which AN flow is laminar. We have provided topological arguments that justify the existence of this minimum size, since at least two +1/2 defects (parabolic folds) of the active filaments are required to organize a circulating swirl, and the defect separation is regulated by l α 7,8,25,35,36 . Moreover, we have found that the scaling of this characteristic size with the material parameters α and K is consistent with the one expected for l α . This evidence, along with the topological arguments presented above, let us summarize our findings as a manifestation of the different roles played by the same active length scale l α , depending on the confinement constraints.
In addition, the ability to entrap the flows of active materials provides with a useful tool to constrain their dynamics and study new organizational features arising in bounded conditions. In this work, by confining the AN in circulating swirls, we have put into evidence the metastable nature of the extensile active material and revealed the existence of periodic rearrangements in the aligned AN. Active stress builds up, and it is released in periodic bursts of defect creation, followed by annihilation and the recovery of the aligned state. Our experiments reveal the existence of a welldefined time scale for these oscillations, which is compatible with the lifetime of vortices in the unconstrained active turbulence regime. One could define a time scale, τ, by combining the average +1/2 defect speed, v, and the intrinsic active length scale, so that τ $ l / =v. Although we have shown that the same intrinsic length scale emerges regardless on the confinement constraints, it is unclear whether the defect speed should be affected by confinement and, therefore, what should be the expected scaling behavior of τ.
The ability to harness active turbulence has allowed in this work to organize well-defined rotating flows, which might have the potential to act as microscale machines. As a proof-ofconcept, we have seeded the active material with passive polystyrene micro-beads, which are dragged by the active flow (Supplementary Movie 6). As they are not bound to a circular trajectory, tracers inside rotating swirls drift outwards while performing spiral trajectories. However, particles can be adsorbed at the LC phase, thus describing perfect circular trajectories that follow the geometry of the SmA planes.
In summary, our experiments have demonstrated that constrainment with an anisotropic soft interface can lead to the modification of the nature and geometry of active flows, whose structure is nevertheless always determined by the same single intrinsic length scale. It is the role of this length scale, and not its value, that is affected by the confinement conditions. This observation overturns our idea of active flows as it confers them with a minimal regularity, and permits to envisage strategies of control vastly different from those applying to their passive counterparts. In addition, the adaptability to soft confinement provides with a handle to probe the intrinsic length and time scales of active materials, and paves the way for further studies comparing classical and active two-dimensional turbulence. Moreover, our findings may unveil an alternative mechanism by which cells self-regulate and control cytoplasmic flows through reconfigurable membranes. In fact, the role of topology Kinesin expression. Drosophila melanogaster heavy chain kinesin-1 K401-BCCP-6His (truncated at residue 401, fused to biotin carboxyl carrier protein (BCCP) and labeled with six histidine tags) has been expressed in Escherichia coli by using the plasmid WC2 from The Gelles Laboratory (Brandeis University) and purified with a nickel column. After dialysis against 500 mM Imidazole buffer, kinesin concentration is estimated by means of absorption spectroscopy and stored at a specific concentration in a 40% (wt/vol) sucrose solution at −80°C 40 .
Assembly of the active gel. MTs are mixed with the motor clusters that act as cross-linkers, and with ATP (Sigma; A2383) that drives the activity of the gel. The aqueous dispersion contains a non-adsorbing polymer (PEG, 20 kDa) (Sigma; 95172) that promotes the formation of filament bundles through depletion interaction ( Supplementary Fig. 2). To maintain a constant concentration of ATP during the experiments, an enzymatic ATP-regenerator system is used, consisting on phosphoenolpyruvate (Sigma; P7127) that fuels pyruvate kinase/lactate dehydrogenase (PK/LDH) (Invitrogen; 434301) to convert ADP back into ATP. Several antioxidant components are also included in the solution to avoid protein denaturation, and to minimize photobleaching during characterization by means of fluorescence microscopy. The PEG-based triblock copolymer surfactant Pluronic F-127 (Sigma; P-2443) is added at 1% (wt/wt) (final concentration) to procure a biocompatible water/oil interface in subsequent steps.
Experimental Setup. The studied AN is formed at the interface with either an isotropic silicone oil (BlueStar Silicones; BlueSil ® V47:12500) or 8CB. The interface is prepared in a cylindrical well of diameter 5 mm and depth 2 mm, manufactured with a block of poly-dimethylsiloxane (Sylgard ® , Dow Corning) using a custom mold. The block is glued onto a bioinert and superhydrophilic polyacrylamide (PAA)-coated glass, which is prepared following reported protocols 41 (Supplementary Fig. 5). In brief, clean and activated glass is first silanized with an acidified ethanolic solution of 3-(Trimethoxysilyl)propylmethacrylate (Sigma; 440159), which will act as polymerization seed. The silanized substrates are subsequently immersed in a aqueous solution of acrylamide monomers (Sigma; 01697) (2% (wt/vol), for at least 2 h) in the presence of the initiator ammonium persulphate (Sigma; A3678) and N,N,N′,N′-Tetramethylethylenediamine (Sigma; T7024), which catalyzes both initiation and polymerization of acrylamide. The cavity is first filled with~30 μl of oil and, subsequently, 1 μl of the water-based active gel is injected between the hydrophobic oil and the superhydrophilic glass plate (Supplementary Fig. 5). Pluronic F-127 stabilizes the interface and avoids direct contact between the oil and the protein-based active mixture. Moreover, this surfactant ensures a planar alignment of 8CB molecules at the water/LC interface. When 8CB is used, samples are heated up to 35°C to promote transition to the less viscous nematic phase of the LC, which facilitates the spreading of the active gel onto the PAA-coated substrate. Temperature is controlled by using a thermostatic oven built with Thorlabs SM1 tube components and foil heater, which is regulated with a Thorlabs TC200 controller. After several minutes at room temperature, the active material in the gel spontaneously condenses onto the flat water/oil interface, leading to the formation of the AN layer. Unlike conventional flow cells, in which a layer of the active gel is confined in a thin gap between two glass plates 12 , our open setup enables us to prepare the interface using oils with a high viscosity.
Imaging. Routine observations of the AN were performed by means of conventional epifluorescence microscopy. We used a custom-built inverted microscope with a white led light source (Thorlabs MWWHLP1) and a Cy5 filter set (Edmund Optics). Image acquisition was performed with a QImaging ExiBlue cooled CCD camera operated with ImageJ μ-Manager open-source software. For sharper imaging of the interfacial region, we used a Leica TCS SP2 laser-scanning confocal microscope equipped with a photomultiplier as detector and a HeNe-633-nm laser as light source. A ×20 oil immersion objective was employed. We performed confocal acquisition in fluorescence and reflection modes, simultaneously.
Image analysis. Tracer-free velocimetry analysis of the AN was performed with a public domain particle image velocimetry (PIV) program implemented as an ImageJ plugin 42 . Flows were also traced by dispersing PEGylated spherical polystyrene microbeads of diameter 1.7 μm (Micromod; 08-56-173). Manual Tracking ImageJ plugin was used to manually track trajectories of particles or defects in the AN. Further analysis of velocimetry data was performed with custom-written MatLab codes. In order to quantify the vortex distribution in the turbulent AN, the Okubo-Weiss (OW) parameter was mapped from the velocity fields obtained by image velocimetry (Fig. 1b). The vortex areas were quantified from thresholded OW field binary images in ImageJ (Fig. 1c). Vortices on the edges of the field of view were excluded.