Active micromachines: Microfluidics powered by mesoscale turbulence

See allHide authors and affiliations

Science Advances  08 Jul 2016:
Vol. 2, no. 7, e1501854
DOI: 10.1126/sciadv.1501854


Dense active matter, from bacterial suspensions and microtubule bundles driven by motor proteins to cellular monolayers and synthetic Janus particles, is characterized by mesoscale turbulence, which is the emergence of chaotic flow structures. By immersing an ordered array of symmetric rotors in an active fluid, we introduce a microfluidic system that exploits spontaneous symmetry breaking in mesoscale turbulence to generate work. The lattice of rotors self-organizes into a spin state where neighboring discs continuously rotate in permanent alternating directions due to combined hydrodynamic and elastic effects. Our virtual prototype demonstrates a new research direction for the design of micromachines powered by the nematohydrodynamic properties of active turbulence.

  • Mesoscale turbulence
  • active matter
  • microrotor array
  • self-organised spin-state
  • activity-powered micromachines
  • biological motors


A hallmark of active fluids, such as bacterial suspensions and filament-motor protein mixtures, is the emergence of large-scale chaotic, turbulent-like structures known as mesoscale (or active) turbulence (14). The intrinsic disorder of these flow fields acts as a barrier to the extraction of useful power from active suspensions. Predesigned microfluidic systems have been built to direct cell motion (58), pump fluid (9), or drive microrotors (10, 11). These impose symmetry breaking through geometrical constraints or external biases. In particular, it is possible to obtain a persistent rotation for specially designed microscopic gears submersed in bacterial baths (1214). However, it has been argued that an asymmetric gear shape is required to achieve a persistent rotation and that bacterial “traps” in the form of symmetric teeth are essential for generating directed motion of gears (13). Here, we argue that asymmetry of the gears is not a prerequisite and propose a means of exploiting spontaneous symmetry breaking of mesoscale turbulence to turn many symmetric rotors.

We use numerical simulations to demonstrate that mesoscale turbulence can be exploited to generate work by considering a two-dimensional (2D) square array of cylindrical rotors that are free to turn in response to mesoscale flows in the solution. The array of discs achieves a state of permanent rotation below a given interdisc spacing (Fig. 1A). The spin state of the rotors is fixed, with neighboring discs rotating in opposite directions. Hence, the disc array rectifies the active turbulence into a steady power source.

Fig. 1 Emergence of antiferromagnetic spin state in an array of rotating discs.

(A) An array of counter-rotating discs in mesoscale turbulence generated by active matter. Rotors have a radius R = 5 and are fixed on a square lattice with cell edge of 100. The director field is shown by solid green lines and is overlaid by streamlines (black lines with arrows). Discs are colored according to their angular velocity as shown in the colorbar (red, anticlockwise; blue, clockwise). A 3 × 3 portion of an 8 × 8 lattice of rotors is shown (movies S1 and S2). (B) Cross-correlations, 〈ωi(tj(t)〉, between neighboring rotors (i, j) as a function of the minimum gap size d between rotors for (i) homeotropic anchoring, (ii) no anchoring, and (iii) planar anchoring of the director at the rotor surface. The root mean square (rms) angular velocity, Embedded Image, is shown in the inset. (C) Time correlations of rotor angular velocity 〈ω(t)ω(0)〉, normalized by its mean value 〈ω(0)2〉 for different lattice sizes as a function of scaled time T = t/(2π/ωrms) for the case with homeotropic boundary conditions on the rotor surface. Inset: The decorrelation time τ is shown as a function of gap size for d > 10.

To determine how the antiferromagnetic spin state of the rotors is established, we consider a solitary disc immersed in an active bath. Protracted, but nonpermanent, rotation can be achieved by synergistic hydrodynamic and elastic effects, which both depend on the properties of the active flow and on the orientational anchoring of the active fluid at the disc surface. We find an exponential dependence of the rotational persistence time on rotor size, which suggests that the decorrelation of rotor angular velocity is controlled by an activated barrier-crossing process. In addition, by placing a single rotor within an immobile circular cavity, we explain how rotation stabilization results from fluid confinement between neighboring discs in the lattice configuration and demonstrate that the characteristic length of activity-induced vortices sets the optimal length scales for the microfluidic design.


A square lattice of 64 rotors is deployed in simulations with periodic boundary conditions on domain boundaries and with either homeotropic (perpendicular to the boundary), planar (parallel to the boundary), or no anchoring of the director (characterizing the orientation field) at the rotor surface. Drag on the rotors due to spontaneous active flows causes them to spin (Fig. 1A). Each disc preferentially rotates in the opposite direction to its nearest neighbors. Although Fig. 1A shows only a subset of rotors, the alternating pattern repeats throughout the entire system (movies S1 and S2). This strong anticorrelation between angular velocities occurs when the gap size d between of the nearest neighbors is small. However, if the spacing is increased, the spins of neighboring rotors decorrelate, and the rotors spin randomly (Fig. 1B).

Moreover, the rotation of a single disc in the array is extremely enduring (Fig. 1C and movie S3). When the separation between rotors is small, it becomes constant on the time scale of our simulations (107 simulation time steps, which covers more than 1000 rotations of the rotors), but when the spacing increases, the angular velocity autocorrelation function of a single rotor, 〈ω(t)ω(0)〉, approaches zero over time (Fig. 1C). The autocorrelation function describes the propensity of the rotor to continue rotating in its current direction and is characterized by the time scale τ. Even for the largest spacings, the rotors maintain a nonzero persistence time of rotation.

Streamlines in Fig. 1A demonstrate that local flows move in the same direction as each rotor spins for small gap sizes. The array of spinning rotors stabilizes the active fluid flow, which would otherwise exhibit mesoscale turbulence. Simulations that consider the extreme case of infinite friction (by disallowing disc rotation) show the same streamlines and director fields as the arrays of free rotors, although with reduced flow speeds. This indicates that the self-organized spin state should be expected in experimental realizations of microfluidic systems of rotors that include finite rotational friction or have been designed to act as active turbines.

For a free rotor to change its direction of rotation in the small separation regime, it must overcome the local flow field, which is not only established through the boundary conditions of the disc in question but also entrenched by the flow field of neighboring rotors. Although the decorrelation time of the angular velocity increases with decreasing gap size, the rotation rate is nonmonotonic with a sharp maximum (Fig. 1B, inset). The maximum is set by the competition between gap size and active nematic fields. To better understand the mechanisms leading to steady anticorrelated rotation and the impact of rotor spacing, we systematically investigate the dynamics of a solitary rotor and a single rotor within a concentric cavity. These single-rotor results demonstrate how the permanent spin state arises in the array and explain how the optimal rotation rate relates to the intrinsic characteristic vorticity length scale of the active turbulence.

A solitary disc in an active fluid begins to rotate once activity-driven flow fields are established (Fig. 2). The surrounding active fluid exhibits self-sustained mesoscale turbulence (Fig. 2A) even at zero Reynolds number (15). The chaotic vorticity field is driven by a “reverse energy cascade” produced at microscopic length scales and escalating to larger scales (16, 17). In active nematic fluids, the continuous fluid motion is coupled to the formation of lines of strong bend distortions in the orientation field, referred to as walls (18, 19), and their unzipping through the creation and annihilation of topological defect pairs (20, 21). Notably, the resulting stochastic velocity and vorticity fields have been noted to exhibit many features of turbulent flow (1, 15, 22) and can be characterized by the vorticity-vorticity correlation length scale.

Fig. 2 Rotation of a single rotor immersed in a turbulent flow of active matter with homeotropic anchoring.

(A) Color maps show the vorticity contours, and solid lines represent streamlines. Vorticity is normalized to its maximum value, and blue and red represent clockwise and anticlockwise rotation, respectively. An animated version is available (movie S3). (B) Director field around the rotating disc. Red dots and yellow triangles denote +1/2 and −1/2 topological defects, respectively. Inset: Schematic of the wall in the director field that encircles each rotor and the resulting active flow of fluid along the wall.

The rotational dynamics are controlled by both the nature of the mesoscale turbulence and the rotor characteristics. Spontaneous symmetry breaking randomly establishes clockwise or counterclockwise spin, and, thereafter, the rotor sporadically switches direction. Hence, as shown in Fig. 3B, the rotation of the disc in response to the flow is propulsive, at short timesEmbedded Image(1)where ωrms is the root mean square angular velocity. The short-time behavior is referred to as propulsive rather than ballistic because rotor inertia is viscously damped and rotation is due to the nematohydrodynamic propulsion from the active flow field. At longer times, there is a crossover to a diffusive regimeEmbedded Image(2)where Dr is the effective rotational diffusion constant. The diffusive regime arises as active vortices impinge on the rotor and is therefore controlled by the characteristic decay time of the vorticity-vorticity correlation function Dr ~ τ−1. The decay time τ is substantially prolonged by increasing the size of the rotor, as is apparent from both Fig. 3B and the decay of the angular velocity autocorrelation function.

Fig. 3 Rotational dynamics of a single rotor in mesoscale turbulence.

(A) Mean-squared angular displacement of rotors versus time for different rotor radii R. Director boundary conditions on the rotor surface are homeotropic. (B) Persistence time as a function of R comparing homeotropic and no anchoring director boundary conditions on the rotor surface.

The rotation is driven by two synergistic mechanisms: (i) entrainment of an active turbulence vortex around the rotor and (ii) wall-driven flow at the surface of the rotor (Fig. 2B, inset). The second mechanism arises due to the orientational boundary condition on the surface of a disc. For homeotropic anchoring, bend deformations in the director field appear naturally at the surface (23, 24). The bend deformations result in the formation of a wall around the perimeter of each rotor (Fig. 2B, inset) (18, 19). Simulations of colloids pulled through active nematic fluids have previously shown that orientational distortion can lead to local flows that cause non-Stokesian dynamics and even negative drag (25). Likewise, in this system, the rotor-encircling wall establishes a unidirectional flow (18) that generates localized vorticity (movie S3). Supporting this argument, if, instead of homeotropic anchoring, no anchoring boundary conditions are implemented, the rotational persistence time of the rotor is significantly diminished because no wall is formed at the disc surface (Fig. 3C and movie S4). The wall-generated flows spontaneously establish a direction of rotation, and any corotational vortices enhance the rotational velocity when the characteristic vortex size of the mesoscale turbulence is commensurate to the rotor size. Thus, the unidirectional flow generated within the wall and the swirling flows within active vortices collude to give a persistent rotation of the disc.

Extracting the decay times from the autocorrelation functions and plotting them against R shows that τ grows exponentially with the radius of the rotor (Fig. 3C). The exponential dependence of rotational persistence time on rotor radius suggests an activated barrier crossing, in which flipping a rotor’s spin is a stochastic process over an energy barrier that is linearly dependent on rotor size. The barrier that must be overcome by local fluctuations in the vorticity field to decorrelate the propulsive rotation is the reorientation of the wall around each rotor, with a total energy proportional to the perimeter ~ R. In this way, spin flipping is analogous to a nonthermal activated Kramers’ process with τ ~ eβR, where β is a nonthermal constant. Once the wall flips direction, it drives a reversal of the unidirectional flow in the immediate vicinity of the surface. In this way, the intrinsic active elastic properties of the nematic establish a barrier to active turbulence that flips the rotor spin. As the activity is increased, the rotation rate increases linearly (Fig. 4A). However, the increasingly energetic fluctuations in the mesoscale turbulence are more able to propel the system over the energy barrier between different directions. Hence, the correlation time decreases with activity (Fig. 4B).

Fig. 4 Characterization of a single rotor in mesoscale turbulence.

(A) Increasing activity results in an enhanced angular velocity. (B) Increasing activity shortens the rotational persistence time.

A solitary symmetric rotor ensnares a vortex and establishes a barrier to decorrelation, leading to large persistence times but not as persistent a spin state as observed for arrays of discs. We now ask why the array of rotors can maintain such a stable spin state. We consider a single disc within a fixed circular cavity (Fig. 5). Several authors have shown that confinement can stabilize active flows and leads to spontaneous circulations (23, 26, 27). Even in the absence of a prescribed orientation at the surfaces, the active fluid self-organizes within the confinement of the annulus to induce a unidirectional flow, which leads to an enduring rotation of the inner cylindrical rotor (Fig. 5A). As the gap size D between the rotor and confinement boundary increases from small values, the angular velocity of the rotor increases (Fig. 5). The rotation rate continues to increase with gap size until defects and distortions appear within the orientation field (Fig. 5B), at which point the rotation rate crests and begins to decrease. The maximum corresponds to the gap size becoming equal to the characteristic length scale of vorticity in bulk active turbulence, measured from the decay of the vorticity autocorrelation function (Fig. 5 and movie S5).

Fig. 5 Confinement-induced rotation of a disc immersed in a bath of active matter showing the transition from unidirectional to turbulent flow.

(A) The main panel shows the variation of the angular velocity of the rotor with increasing gap size. The vertical dashed line corresponds to the characteristic length scale of bulk mesoscale turbulence measured from the vorticity-vorticity correlation function in bulk. The insets show the transition from unidirectional flow to active turbulence with increasing gap size D: (i) D = 8; (ii) D = 30; (iii) D = 70. Shown is the director field (dashed green lines), superposed by streamlines (solid black lines with arrows) and topological defects (red and blue dots denoting +1/2, and −1/2 defects, respectively). (B) Mean square angular displacement of the rotor as a function of simulation time, where the transition from propulsive to diffusive behavior is seen when gap size is > 70.

The maximum in the rotation rate is also marked by the appearance of topological defects in the director field. These leave the streamlines relatively unperturbed, with the defects advected in a laminar manner (Fig. 5B). However, for still larger gaps, obvious vortices appear in the flow field (Fig. 5C), and the rotor rotation rate slowly decreases with increasing D (Fig. 5).

The rotational persistence time of the rotor is finite at large gap sizes with multiple vortices but is substantially enhanced when only a single stabilized vortex fits within the confinement. On the other hand, dissipation is increased when the gap size is smaller than the natural size of the vortices. From this, we understand the persistent counter-rotating spin state of a disc array (Fig. 1A): each rotor establishes a wall-induced active shear flow, while neighboring rotors do the same, simultaneously creating a confinement that stabilizes the local flow. This process is most efficient if the spacing between rotors is comparable to the characteristic vorticity length scale of the mesoscale turbulence (Fig. 5 and movie S5). In typical bacterial suspensions in 2D geometries, the characteristic length scale is ~ 40 μm (15).


We have demonstrated the persistent, correlated rotation of an array of symmetric rotors immersed in a dense suspension of active matter, such as bacterial suspensions (1), filament-motor protein solutions (28), or cellular monolayers (29). The results show that harnessing useful work from the disordered and chaotic flows of active materials does not require specifically designed asymmetric microscopic gears or external fields.

The persistent spin state of the rotor array is antiferromagnetic, similar to the recently observed spin states of bacteria suspensions confined within interconnected cavities (30). However, the physical mechanisms leading to these particular spin states are markedly different. Whereas the rotation of Bacillus subtilis in interconnected circular cavities is dictated by a counter-rotating single-cell layer that is controlled by the size of the intersection compared to the cell size (30), the antiferromagnetic spin state of rotors studied here arises through a combination of an activated barrier-crossing process and of mooring active nematic vortices to each rotor. The array of spinning rotors stabilizes the active fluid flow, which would otherwise be in a state of mesoscale turbulence.

The way the rotor lattice constrains active excitations of the vorticity field, which in turn causes nonzero torques to arise on the rotors themselves, bears a strong conceptual similarity to Casimir forces (as defined in the general sense and not just for the electromagnetic field fluctuations) (31). The notion of using this concept to design micromachines has been explored in the context of the traditional Casimir forces, revealing interesting possibilities afforded because of a lack of physical contact (32, 33). A key aspect of the Casimir effect is that it depends primarily on geometric characteristics of the system and, as such, offers a simple route to tuning the behavior of the micromachine. We note that Casimir forces have been postulated to exist in active systems (3436).

By investigating the dynamics of solitary rotors, we determine how different length scales, namely, the rotor size, the lattice spacing, and the gap between the rotors, affect the behavior of the array of rotors. We find that in larger gaps, decorrelating vortices can be excited between rotors, whereas in smaller spaces, counter-rotating vortices around neighboring discs viscously dissipate energy. Additionally, the spin state of each rotor is locked by an energy barrier due to the thin active nematic wall that forms on its surface, and larger rotors have longer perimeters and higher barriers that further stabilize the spin state. Thus, the choices of configuration and size of rotors are specific to a microfluidic design, and these choices can be optimized using vorticity correlation length as the characteristic length scale of the fluid.

Our findings demarcate a promising direction toward harnessing power from active matter. The ability to exploit the continuous injection of energy from constituent elements of active matter to obtain an ordered lattice of counter-rotating rotors has implications in many applications, from providing driving power for microelectromechanical devices to acting as novel microfluidic mixers.


We use the equations of active nematohydrodynamics based on the theory of liquid crystals, which has proven successful in describing spatiotemporal dynamics of active matter systems, including bacterial suspensions (37), microtubule bundles (2, 38, 39), and cellular monolayers (29, 40). The orientational order of microscopic active particles is represented by the nematic tensor Embedded Image, where q is the magnitude of the orientational order, n is the director, and I is the identity tensor, which evolves due to a corotation term and relaxation through rotational diffusivity Γ of the molecular field, which accounts for both the Landau–de Gennes bulk free energy and the nematic distortion free energy (41), assuming a single elastic constant K. The total density and the velocity field u of the active matter obey the incompressible Navier-Stokes equations, with a stress tensor Π that must account for contributions from the viscosity η, the elastic stress (which includes the pressure P) (42), and the active contribution to the stress Πact = − ζQ (43). Any gradient in Q generates a flow field through the active stress, with the strength determined by the activity coefficient ζ.

The equations of active nematohydrodynamics are solved using a hybrid lattice Boltzmann technique (21) that does not include thermal fluctuations. The momentum equation is solved using the lattice Boltzmann method to resolve the hydrodynamics. The method of lines (21) is implemented to determine the order parameter, in which a finite-difference approach is used for spatial discretization, and the temporal evolution is obtained through a Euler integration scheme. For details, see the studies of Thampi et al. (21), Marenduzzo et al. (24), and Fielding et al. (44).

As usual in lattice Boltzmann schemes, discrete space and time steps are chosen as unity, and all quantities can be converted to physical units in a material-dependent manner (37, 45, 46). The parameters used in the simulations are Γ = 0.34, K = 0.01, ζ = 0.01, and η = 2/3, in lattice units. These parameters result in active turbulence in the bulk. Unless otherwise stated, simulations were performed in a 2D domain of size 200 × 200. Microrotors are modeled as discs discretized on the simulation lattice. Each rotor is fixed in space but allowed to spin freely with a moment of inertia I = 103 for all radii. For the parameters used here, viscous damping dominates over inertia in approximately ~ IR2 ≈ 100 time steps. No-slip boundary conditions allow flows to apply torques to the rotors.


Supplementary material for this article is available at

movie S1. Antiferromagnetic spin state in an array of rotating discs.

movie S2. Antiferromagnetic spin state in an array of rotating discs.

movie S3. Rotation of a single rotor immersed in active turbulence.

movie S4. Rotation of a single rotor immersed in active turbulence.

movie S5. Confinement-induced rotation of a disc immersed in an active nematic bath.

This is an open-access article distributed under the terms of the Creative Commons Attribution-NonCommercial license, which permits use, distribution, and reproduction in any medium, so long as the resultant use is not for commercial advantage and provided the original work is properly cited.


Acknowledgments: We are grateful to M. Pushkin for insightful discussions. Funding: This work was funded by ERC Advanced Grant 291234 MiCE and was supported by EMBO funding to T.N.S. (ALTF181-2013). Author contributions: J.M.Y. conceived and supervised the project. S.P.T. designed and performed the simulations. All authors contributed to interpreting the results and writing the manuscript. Competing interests: The authors declare that they have no competing interests. Data and materials availability: All data needed to evaluate the conclusions are present in the paper and the Supplementary Materials. Additional data related to this paper may be requested from the authors.
View Abstract

Navigate This Article